1

我在编程生活中遇到了最奇怪的问题。我正在对大约 800x800x600 个元素的 3D 体积进行图像分析，从每个像素中提取一个基于灰度的粗麻布矩阵，并在检测体积中描述的纤维的对齐之后计算一些代数（使用Eigen ）。

为了构建卷，我使用了一个 double[x][y][z] 数组，该数组是从 nhdr+raw 文件中读取的。到目前为止，对于域问题的介绍非常有用。

我执行将大体积分成子体积的分析（总执行时间约为 14 小时，但这是意料之中的，就像科学分析软件一样）。在某个时刻，当我处于数组的 x=(max-2) （我在边界前 2 个像素处停止分析）索引处时，对于几个元素（[792][247][76] 和 [ 792][84][56]，但在那些之后可能还有其他人，可能仍然是 x=max），我的特征函数突然失败。

Debuggin 我发现要构建 Hessian 矩阵的偏导数，我访问元素 [x+2][y][z] 并从中获得 NaN，因此显然 Eigen 在对其进行操作时发疯了。

最奇怪的是，如果在软件开始时（加载卷之后），如果我打印那个确切元素的值，它就存在并且它也有一个有意义的值。这怎么可能？？我多次运行该软件，对于相同的两个像素，相同位置的相同错误，所以即使猜测它可能是 RAM 错误，它不应该由于我的 PC 中发生的其他事情而以某种方式波动和改变位置吗？

我进一步进行了测试。我在子卷上循环，一切都很好（我继续跟踪一个固定的 volume[][][] 元素的值，故障本身就是在这个元素上表现出来的）像素值在以下函数之外保持不变，那个像素一个像素地分析子体积。据我所知，我感兴趣的像素值是 51941（并且在进入给定子卷的以下函数之前），我在值发生变化时设置了警戒线。

这里发生了什么 > 子卷中的像素值开始：51941 元素：5982；636,260,62; 循环中的像素值：1.65031e-22

在 5982 次循环（整个子卷需要大约 300 万次）之后，值发生了变化，并且在下面的代码中我没有触及它！什么可能导致这样的事情？

Matrix3d HesseAnalysis(int subX, int subY, int subZ, int maxX, int maxY, int maxZ){

//int maxX = subX+deltaX;
//int maxY = subY+deltaY;
//int maxZ = subZ+deltaZ;

int counter=0;

for (int x=subX;x<(maxX-2);x++){
    for (int y=subY;y<(maxY-2);y++){
        for (int z=subZ;z<(maxZ-2);z++){


            if(volume[792][247][76]!=51941){
                cout << "Element :" << counter << "; " << x << "," << y << "," << z << ";" << endl;
                cout << "Pixel value in loop: " << volume[792][247][76]<< endl;
                exit(0);
            }

            fxx=((volume[x+2][y][z]-volume[x][y][z])/2-(volume[x][y][z]-volume[x-2][y][z])/2)/2;
            fyy=((volume[x][y+2][z]-volume[x][y][z])/2-(volume[x][y][z]-volume[x][y-2][z])/2)/2;
            fzz=((volume[x][y][z+2]-volume[x][y][z])/2-(volume[x][y][z]-volume[x][y][z-2])/2)/2;
            fxy=((volume[x+1][y+1][z]-volume[x+1][y-1][z])-(volume[x-1][y+1][z]-volume[x-1][y-1][z]));
            fxz=((volume[x+1][y][z+1]-volume[x+1][y][z-1])-(volume[x-1][y][z+1]-volume[x-1][y][z-1]));
            fyz=((volume[x][y+1][z+1]-volume[x][y+1][z-1])-(volume[x][y-1][z+1]-volume[x][y-1][z-1]));

            //compose hessian matrix for the pixel, remember that  
 hessian <<     fxx, fxy, fxz,
    fxy, fyy, fyz,
    fxz, fyz, fzz;

//extract eigenvalues and choose the eigenvector related to the smallest eigenvalue,
//and do the outer product of it with itself

EigenSolver<Matrix3d> solver(hessian);
int minorEigen = minorEigenvalue(solver.eigenvalues().real());
Vector3d v3 =solver.eigenvectors().col(minorEigen).real();
V3outerV3 = v3*v3.transpose();
OuterProducts[(x-subX)-2][(y-subY)-2][(z-subZ)-2]= V3outerV3;

counter++;

            }
        }
    }
4

2 回答 2

1

我找到了。我仍然没有时间探索真正的原因，但这取决于 Eigen 的一个可能的错误。

这一行是导致混乱的行：

Vector3d v3 =solver.eigenvectors().col(minorEigen).real();

出于某种原因，在某个时刻，它决定侵入已经为我的 volume[ ][ ][ ] 数组分配的内存，这是最奇怪的部分，我用于测试的像素，volume[792] [247][76] 值根据变量 minorEigen 的值（即简单地返回函数中最小特征值的索引）而改变。

因此，例如，minorEigen = 0 像素变为 1.65031e-22，minorEigen = 1 变为 2.1402e+5... 我将调查 Eigen 错误报告，检查更新的版本，或者为矩阵实现我自己的 eigensolver。

于 2013-10-16T15:02:06.750 回答
1

通常 C/C++ 使用从 0 开始的索引。如果你不希望你x的超出范围，你应该停在x=max-3，所以索引可以保持在范围 from 0to max-1？

于 2013-10-15T11:05:57.740 回答