http://doc.cgal.org/latest/Point_set_processing_3/index.html#Point_set_processing_3NormalOrientation を少し変更して実行した。
法線推定
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h> #include <CGAL/pca_estimate_normals.h> #include <CGAL/mst_orient_normals.h> #include <CGAL/property_map.h> #include <CGAL/IO/read_off_points.h> #include <CGAL/IO/write_xyz_points.h> #include <cstring> // to use strcmp #include <string> #include <utility> #include <list> #include <fstream> #include <iostream> #include <cstdlib> // to use EXIT_XXX macros typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; typedef Kernel::Point_3 Point; typedef Kernel::Vector_3 Vector; // Point with normal vector stored in a std::pair. typedef std::pair<Point, Vector> PointVectorPair; // Concurrency #ifdef CGAL_LINKED_WITH_TBB typedef CGAL::Parallel_tag Concurrency_tag; #else typedef CGAL::Sequential_tag Concurrency_tag; #endif int main(int argc, char* argv[]) { std::list<PointVectorPair> points; if ( argc != 3 || !std::strcmp(argv[1], argv[2]) ) { std::cerr << "usage: ./a.out input.off output.off" << std::endl << " input and output must have different names" << std::endl; return EXIT_FAILURE; } const std::string input_filename = argv[1]; const std::string output_filename = argv[2]; std::ifstream input(input_filename.c_str()); if (!input || !CGAL::read_off_points(input , std::back_inserter(points) , CGAL::First_of_pair_property_map<PointVectorPair>())) { std::cerr << "Error: cannot read " << input_filename << std::endl; return EXIT_FAILURE; } // estimate normals const int nb_neighbors = 18; CGAL::pca_estimate_normals<Concurrency_tag>(points.begin(), points.end(), CGAL::First_of_pair_property_map<PointVectorPair>(), CGAL::Second_of_pair_property_map<PointVectorPair>(), nb_neighbors); // orient normals std::list<PointVectorPair>::iterator unoriented_points_begin = CGAL::mst_orient_normals(points.begin(), points.end(), CGAL::First_of_pair_property_map<PointVectorPair>(), CGAL::Second_of_pair_property_map<PointVectorPair>(), nb_neighbors); std::ofstream output(output_filename.c_str()); if (!output || !CGAL::write_xyz_points_and_normals(output , points.begin(), points.end() , CGAL::First_of_pair_property_map<PointVectorPair>() , CGAL::Second_of_pair_property_map<PointVectorPair>())) { std::cerr << "Error: cannot write " << output_filename << std::endl; return EXIT_FAILURE; }<img src="http://f.st-hatena.com/images/fotolife/k/kuyata/20160422/20160422014516.png" width=350> return EXIT_SUCCESS; }
コンパイルはclang++ -std=c++11 -I/home/vagrant/cgal/cgal/usr/local/include/ -I/usr/include/eigen3/ -I/usr/include/ -lCGAL -lboost_thread -lboost_system -lgmp -lmpfr
- rglで可視化する
camel <- read.table("outputcamel.xyz") body <- 50 * camel[,1:3] #(適当に拡大する) normal <- camel[,4:6] plot3d(rgl) # 9770頂点の内900頂点に法線をプロットする for (i in 1:900) { ori <- body[10*i,] iro <- body[10*i,] + normal[10*i,] x <- rbind(ori,iro) segments3d(as.matrix(x), add=TRUE, c='blue') }
この段階では形は保存されているように見えるが、足の部分の法線は方向がバラバラになっている。
- 法線だけプロットしたもの
poisson reconstruction
CGALのサンプルを実行する - tak0kadaの何でもノートと同様にした
結果
またしてもうまく行かず、後ろ足が切れてしまった。
分析
元データcamel.offには法線データがないのでkitten.xyzとoni.xyzで分析してみる。
# データ読み込み、ctrlが元データ、testが推定データ ctrl <- read.table("ctrl.xyz")<img src="http://f.st-hatena.com/images/fotolife/k/kuyata/20160422/20160422014516.png" width=350> test <- read.table("test.xyz") # 座標を比較する # TRUE any(ctrl[, 1:3] == test[, 1:3]) # 法線のずれを調べる diff_norm <- apply(test[,4:6] - ctrl[,4:6] ,1 ,function(x){x[1]^2 + x[2]^2 + x[3]^2}) # tagcloudを使ってグラデーションをつけることができる # colorは差が小さと1→青、大きいと0→赤に対応(分かりにくいが) library(tagcloud) color <- (max(diff_norm) - diff_norm) / (max(diff_norm) - min(diff_norm)) plot3d(ctrl[,1:3], col = smoothPalette(color, pal="RdBu"))
比較的形が保存されているように見えるネコでもメッシュが異なる。鬼は口と角でメッシュの違いが大きい。青と赤の点で表示したのは法線の差の絶対値の2乗で、ネコはほぼ均一、鬼はメッシュの違いが表れた口と角でズレが大きい。元のデータで相当数の頂点が含まれているがメッシュ化後では頂点数がかなり減少している。
法線の差をプロットすると、凹凸になっている部分で外側への成分が過小評価されているのが分かる。