2016-12-29 29 views
0

ポイントのセットに最適なフィッティング平面を見つけることを試みており、SVDを使用してax+by+cz+d=0の平面方程式を計算しています。SVDを使用して最良のフィッティング平面方程式を計算する際にエラーが発生する

私はSVDを実装して正常に飛行機に乗ることができましたが、dを計算することができません。

掘り出した後、計算した重心を式に戻してdと計算しましたが、間違った値になっています。私はこれをRANSACメソッドと比較しているので、これは間違った値であると確信しています。

私は取得しています結果は

pcl::ModelCoefficients normal_extractor::plane_est_svd(pcl::PointCloud<pcl::PointXYZ>::ConstPtr point_cloud) 
{ 
    Eigen::MatrixXd points_3D(3,point_cloud->width); 
    //assigning the points from point cloud to matrix 
    for (int i=0;i<point_cloud->width;i++) 
    { 
     points_3D(0,i) = point_cloud->at(i).x; 
     points_3D(1,i) = point_cloud->at(i).y; 
     points_3D(2,i) = point_cloud->at(i).z; 
    } 
    // calcaulating the centroid of the pointcloud 
    Eigen::MatrixXd centroid = points_3D.rowwise().mean(); 
    //std::cout<<"The centroid of the pointclouds is given by:\t"<<centroid<<std::endl; 
    //subtract the centroid from points 
    points_3D.row(0).array() -= centroid(0); 
    points_3D.row(1).array() -= centroid(1); 
    points_3D.row(2).array() -= centroid(2); 
    //calculate the SVD of points_3D matrix 
    Eigen::JacobiSVD<Eigen::MatrixXd> svd(points_3D,Eigen::ComputeFullU); 
    Eigen::MatrixXd U_MAT = svd.matrixU(); 
    //std::cout<<"U matrix transpose is:"<<U_MAT<<std::endl<<std::endl<<"U matrix is:"<<svd.matrixU()<<std::endl; 
    /********************************************************************************************* 
    * caculating d by sybstituting the centroid back in the quation 
    *  aCx+bCy+cCz = -d 
    ********************************************************************************************/ 
    //double d = -((U_MAT(0,2)*points_3D(0,1))+ (U_MAT(1,2)*points_3D(1,1)) + (U_MAT(1,2)*points_3D(1,2))); 
    double d = -((U_MAT(0,2)*centroid(0))+ (U_MAT(1,2)*centroid(1)) + (U_MAT(1,2)*centroid(2))); 

    pcl::ModelCoefficients normals; 
    normals.values.push_back(U_MAT(0,2)); 
    normals.values.push_back(U_MAT(1,2)); 
    normals.values.push_back(U_MAT(2,2)); 
    normals.values.push_back(d); 
    return(normals); 

} 

を次のように私のコードの実装がある

RANSAC方法:

a = -0.0584306 b = 0.0358117 c = 0.997649 d = -0.161604 

SVD方法:

a = 0.0584302 b = -0.0357721 c = -0.99765 d = 0.00466139 

から結果、私法線方向が逆転しているにもかかわらず、法線はちょうど良いと考えられますが、dの値は正しくありません。私はどこが間違っているのか分からない。どんな助けでも本当に感謝しています。

事前に感謝..

+1

場合 ' d'が正しく計算されていない場合は、使用している式を確認してください。最後の言葉( 'U_MAT(1,2)* centroid(2)')は正しくありません( 'U_MAT(2、2)')? – 1201ProgramAlarm

+0

はい、もちろん、ありがとうございます。何らかの理由で、私はそれに選択的に盲目的でした。 @ 1201ProgramAlarmもう一度ありがとうございます。 – spacemanspiff

答えて

1

1201ProgramAlarmは右である、U_MAT(1,2)*centroid(2)にタイプミスがあります。

なタイプミスより良い書き込みを回避するには、次の

d = -centroid.dot(U_MAT).col(2); 

をあなたも簡素化することができます。今後の参考のために

points_3D.colwise() -= centroid; 

を、ここでは自己完結型の例である:

#include <iostream> 
#include <Eigen/Dense> 
using namespace Eigen; 
using namespace std; 

int main() 
{ 
    int n = 10; 
    // generate n points in the plane centered in p and spanned bu the u,v vectors. 
    MatrixXd points_3D(3,n); 
    Vector3d u = Vector3d::Random().normalized(); 
    Vector3d v = Vector3d::Random().normalized(); 
    Vector3d p = Vector3d::Random(); 
    points_3D = p.rowwise().replicate(n) + u*VectorXd::Random(n).transpose() + v*VectorXd::Random(n).transpose(); 
    MatrixXd initial_points = points_3D; 

    Vector3d centroid = points_3D.rowwise().mean(); 
    points_3D.colwise()-=centroid; 
    JacobiSVD<MatrixXd> svd(points_3D,ComputeFullU); 
    Vector3d normal = svd.matrixU().col(2); 
    double d = -normal.dot(centroid); 

    cout << "Plane equation: " << normal.transpose() << " " << d << endl; 
    cout << "Distances: " << (normal.transpose() * initial_points).array() + d << endl; 
} 
+0

あなたの答えをありがとう。はい、それは単にタイプミスでした。私はあなたが示唆した方法で重心の 'd 'と減算の計算を簡素化しようとしましたが、エラー' YOU_TRIED_CALLING_A_VECTOR_METHOD_ON_A_MATRIX''を投げます。私はなぜそれがわかりません。また、 'd = -centroid.dot(U_MAT).col(2);'をd = -centroid.dot(U_MAT.col(2)); ..に変更しました - ありがとう – spacemanspiff

+0

私は自己完結型の例を追加しました。私は問題は、あなたがベクトル型の代わりにジェネリック 'MatrixXd'を使ったことだと思います。 – ggael

関連する問題