@@ -41,7 +41,7 @@ TEST(NormalEstimation, OCC) {
4141 Bnd_Box B;
4242 BRepBndLib::AddClose (face, B);
4343
44- double d = 0.05 ;
44+ double d = 0.01 ;
4545 double x0, y0, z0, x1, y1, z1;
4646 B.Get (x0, y0, z0, x1, y1, z1);
4747
@@ -63,61 +63,60 @@ TEST(NormalEstimation, OCC) {
6363
6464 std::cout << " center " << (*it).format () << std::endl;
6565
66- for ( int i = 1 ; i <= 10 ; i += 2 ) {
66+ int i = 2 ;
6767
68- visitor<26 > vis;
69- vis.max_depth = 0.1 * i / d;
68+ visitor<26 > vis;
69+ vis.max_depth = 0.1 * i / d;
7070
71- std::cout << " md " << (0.1 * i) << std::endl;
71+ std::cout << " md " << (0.1 * i) << std::endl;
7272
73- std::vector<float > coords;
73+ std::vector<float > coords;
7474
75- auto selection = storage->empty_copy ();
75+ auto selection = storage->empty_copy ();
7676
77- vis ([&coords, &selection](const tagged_index& pos) {
78- if (pos.which == tagged_index::VOXEL) {
79- coords.push_back (pos.pos .get (0 ));
80- coords.push_back (pos.pos .get (1 ));
81- coords.push_back (pos.pos .get (2 ));
82- selection->Set (pos.pos );
83- }
84- else {
85- throw std::runtime_error (" Unexpected" );
86- }
87- }, storage, *it);
88-
89- /* {
90- std::ofstream ofs("m2.obj");
91- ((chunked_voxel_storage<bit_t>*)selection)->obj_export(ofs);
77+ vis ([&coords, &selection](const tagged_index& pos) {
78+ if (pos.which == tagged_index::VOXEL) {
79+ coords.push_back (pos.pos .get (0 ));
80+ coords.push_back (pos.pos .get (1 ));
81+ coords.push_back (pos.pos .get (2 ));
82+ selection->Set (pos.pos );
9283 }
84+ else {
85+ throw std::runtime_error (" Unexpected" );
86+ }
87+ }, storage, *it);
9388
94- storage->boolean_subtraction_inplace(selection);
95-
96- {
97- std::ofstream ofs("m1.obj");
98- storage->obj_export(ofs);
99- }*/
89+ /* {
90+ std::ofstream ofs("m2.obj");
91+ ((chunked_voxel_storage<bit_t>*)selection)->obj_export(ofs);
92+ }
10093
101- std::cout << " neighbours " << (coords. size () / 3 ) << std::endl ;
94+ storage->boolean_subtraction_inplace(selection) ;
10295
103- Eigen::MatrixXf points = Eigen::Map<Eigen::MatrixXf>(coords.data (), 3 , coords.size () / 3 ).transpose ();
96+ {
97+ std::ofstream ofs("m1.obj");
98+ storage->obj_export(ofs);
99+ }*/
104100
105- Eigen::MatrixXf centered = points.rowwise () - points.colwise ().mean ();
106- Eigen::MatrixXf cov = centered.adjoint () * centered;
107- Eigen::SelfAdjointEigenSolver<Eigen::MatrixXf> eig (cov);
101+ std::cout << " neighbours " << (coords.size () / 3 ) << std::endl;
108102
109- auto estimated = eig. eigenvectors (). col ( 0 );
103+ Eigen::MatrixXf points = Eigen::Map<Eigen::MatrixXf>(coords. data (), 3 , coords. size () / 3 ). transpose ( );
110104
111- std::cout << estimated << std::endl;
105+ Eigen::MatrixXf centered = points.rowwise () - points.colwise ().mean ();
106+ Eigen::MatrixXf cov = centered.adjoint () * centered;
107+ Eigen::SelfAdjointEigenSolver<Eigen::MatrixXf> eig (cov);
112108
113- auto angle = std::acos (estimated.dot (norm));
114- if (angle > M_PI / 2 ) {
115- angle = M_PI - angle;
116- }
117- auto angle_degrees = angle / M_PI * 180 .;
109+ auto estimated = eig.eigenvectors ().col (0 );
118110
119- std::cout << " angle " << angle_degrees << " d " << std::endl;
111+ std::cout << estimated << std::endl;
120112
121- ASSERT_LT (angle_degrees, 90 .);
113+ auto angle = std::acos (estimated.dot (norm));
114+ if (angle > M_PI / 2 ) {
115+ angle = M_PI - angle;
122116 }
117+ auto angle_degrees = angle / M_PI * 180 .;
118+
119+ std::cout << " angle " << angle_degrees << " d" << std::endl;
120+
121+ ASSERT_LT (angle_degrees, 1 .);
123122}
0 commit comments