OBB在第二次通过"Assertion Failed"时崩溃的特征向量计算

Question

抱歉，如果标题不够清晰，我真的不确定如何总结这个问题，因为它非常具体和奇怪。基本上，我有一个我绘制的网格，我正在尝试计算一个定向边界框。为了做到这一点，我使用了主成分分析。计算特征向量后，我找到了旋转矩阵，该矩阵将分别将特征向量转换为 x、y 和 z 轴。然后，我通过此旋转矩阵转换网格，计算轴对齐的边界框，然后通过旋转矩阵的反转来旋转网格和框。当我运行该程序时，它成功地运行了一次边界框函数，虽然没有框可见，但在第二次运行显示函数时，它在尝试计算边界框时崩溃了。我已经确定了导致崩溃的线，它是特征向量计算的一部分。但是，因为我正在使用一个库来计算特征向量，所以我不完全确定函数中发生了什么，我无法确定为什么它不起作用。当它崩溃时，它会给出错误"断言失败"，后跟我的程序的文件路径和行。

我计算协方差的函数在这里：

std::array<std::array<double, 3>, 3> covarianceCalc2()
{
   std::array<std::array<double, 3>, 3> sum = {{{0, 0, 0}, {0, 0, 0}, {0, 0, 0,}}};
   std::array<double, 3> tempVec;
   double mean;
   for(int i = 0; i < meshVertices.size(); i++)
   {
     mean = (meshVertices[i].x + meshVertices[i].y + meshVertices[i].z)/3;
     tempVec[0] = meshVertices[i].x - mean;
     tempVec[1] = meshVertices[i].y - mean;
     tempVec[2] = meshVertices[i].z - mean;
     sum = matrixAdd(sum, vectorTranposeMult(tempVec));
   }
 sum = matrixMultNum(sum, 1/(meshVertices.size()));
 return sum;
}

以下是计算特征向量并旋转它们的代码：

 Compute_EigenV(covarianceCalc2(), eigenValues, eigenVectors_1, eigenVectors_2, eigenVectors_3);
 std::array<double, 3> x = {1, 0, 0};
 std::array<double, 3> y = {0, 1, 0};
 std::array<double, 3> z = {0, 0, 1};
 std::array<std::array<double, 3>, 3> transformX = findRotation(eigenVectors_1, x);
 std::array<std::array<double, 3>, 3> transformY = findRotation(eigenVectors_2, y);
 std::array<std::array<double, 3>, 3> transformZ = findRotation(eigenVectors_3, z);

 memcpy(transformX, axisRotateX, sizeof(double) * 3 * 3);
 memcpy(transformY, axisRotateY, sizeof(double) * 3 * 3);
 memcpy(transformZ, axisRotateZ, sizeof(double) * 3 * 3);

本节的第一行是导致崩溃的行。Compute_Eigen函数如下：

void Compute_EigenV(std::array<std::array<double, 3>, 3> covariance, std::array<double, 3> eigenValues, std::array<double, 3> eigenVectors_1, std::array<double, 3> eigenVectors_2, std::array<double, 3> eigenVectors_3)
{
    printf("Matrix Stuffn");
    MatrixXd m(3, 3);
    m << covariance[0][0], covariance[0][1], covariance[0][2],
         covariance[1][0], covariance[1][1], covariance[1][2],
         covariance[2][0], covariance[2][1], covariance[2][2];
    // volving SVD
    printf("EigenSolvern");
    EigenSolver<MatrixXd> solver(m);
    MatrixXd all_eigenVectors = solver.eigenvectors().real();
    MatrixXd all_eigenValues = solver.eigenvalues().real();
    // find the max index
    printf("Find Max Indexn");
    int INDEX[3];
    double max;
    max=all_eigenValues(0,0);
    int index=0;
    for (int i=1;i<3;i++){
        if (max<all_eigenValues(i,0)){
            max=all_eigenValues(i,0);
            index=i;
        }
    }
    INDEX[0]=index;
    // find the min index
    printf("Find Min Indexn");
    double min;
    min=all_eigenValues(0,0);
    index=0;
    for (int i=1;i<3;i++){
        if (min>all_eigenValues(i,0)){
            min=all_eigenValues(i,0);
            index=i;
        }
    }
    INDEX[1]=3-index-INDEX[0];
    INDEX[2]=index;
    // giave eigenvalues and eien vectors to matrix
    printf("Give values and vector to matrixn");
    eigenValues[0]=all_eigenValues(INDEX[0],0);
    printf("1");
    eigenValues[1]=all_eigenValues(INDEX[1],0);
    printf("1n");
    eigenValues[2]=all_eigenValues(INDEX[2],0);
    printf("Vector 1n");
    VectorXd featureVector_1 = all_eigenVectors.col(INDEX[0]);
    eigenVectors_1[0]=featureVector_1(0);
    eigenVectors_1[1]=featureVector_1(1);
    eigenVectors_1[2]=featureVector_1(2);
    printf("Vector 2n");
    VectorXd featureVector_2 = all_eigenVectors.col(INDEX[1]);
    eigenVectors_2[0]=featureVector_2(0);
    eigenVectors_2[1]=featureVector_2(1);
    eigenVectors_2[2]=featureVector_2(2);
    printf("Vector 3n");
    VectorXd featureVector_3 = all_eigenVectors.col(INDEX[2]);
    eigenVectors_3[0]=featureVector_3(0);
    eigenVectors_3[1]=featureVector_3(1);
    eigenVectors_3[2]=featureVector_3(2);
}

程序崩溃的特定部分是：

    eigenValues[0]=all_eigenValues(INDEX[0],0);
    eigenValues[1]=all_eigenValues(INDEX[1],0);
    eigenValues[2]=all_eigenValues(INDEX[2],0);

在 上使用 转到定义 all_eigenValues 将显示以下代码部分：

/** returns a reference to the coefficient at given the given row and column.
  *
  * sa operator[](Index)
  */
EIGEN_STRONG_INLINE Scalar&
operator()(Index row, Index col)
{
  eigen_assert(row >= 0 && row < rows()
      && col >= 0 && col < cols());
  return derived().coeffRef(row, col);
}

不幸的是，我对 C 没有足够的经验来真正弄清楚导致这里错误的原因。因为它可能是相关的，所以这是我计算旋转矩阵的函数：

std::array<std::array<double, 3>, 3> findRotation(std::array<double, 3> vec1, std::array<double, 3> vec2)
{
std::array<std::array<double, 3>, 3> tempMatrix, tempMatrix2;
std::array<double, 3> crossProd = crossProduct(vec1, vec2);
double dotProd = dotProduct(vec1, vec2);
double sinAngle = modVector(crossProd)/(modVector(vec1) * modVector(vec2));
std::array<std::array<double, 3>, 3> xProdMatrix = {{{0, crossProd[2] * -1, crossProd[1]}, {crossProd[2], 0, crossProd[0] * -1}, {crossProd[1] * -1, crossProd[0], 0}}}; 
tempMatrix = matrixAdd(iMatrix, xProdMatrix);
tempMatrix2 = matrixMult(xProdMatrix, xProdMatrix);
tempMatrix2 = matrixMultNum(tempMatrix2, (1 - dotProd)/pow(sinAngle, 2));
tempMatrix = matrixAdd(tempMatrix, tempMatrix2);
return tempMatrix;
}

如果您需要更多代码段，我可以为他们提供它们，我已经遇到了这个问题，并且在几个小时内没有取得任何进展，因此任何帮助将不胜感激。

我假设由于特征计算来自一个库，因此错误是我的而不是函数的错误，如果我不得不猜测我的错误在哪里，它会在旋转矩阵计算中，但我无法自己发现错误，如果那是它所在的地方。

编辑：好的，我已经解决了我希望的协方差计算的问题，但不幸的是，同样的崩溃仍然发生。我还根据另一位评论者的建议将所有数组更改为 std：：arrays。我已经在这里待了 7 个小时了。

Answer 1

我相信

我发现了你的问题。

你需要再次检查理论！我可能还记得，理论中的协方差定义为：

C=1/M Sum ( (p-pmean)*(p-pmean)^t )

好吧，你可能会注意到C是一个3x3矩阵，而不是一个值。因此，当您调用Compute_EigenV并尝试制作3x3矩阵的 SVD 时，它会崩溃。

编辑看起来你解决了这部分，并解释了一点，非常好。该断言是eigen库的一部分，如您所见：

  eigen_assert(row >= 0 && row < rows()   &&  col >= 0 && col < cols());

它只是确保输入是 3x3 矩阵（在您的情况下）。不知何故，这并没有发生。我很确定Compute_EigenV代码无需修改即可工作，因为好吧，我的同事写了它并且我已经尝试过了，因此您声明归咎于函数的变量的方式一定有问题。

---

同样，该代码看起来很熟悉;)。