IGLライブラリー
2018年発進のIGLを利用しましょう. 作者はコレ ですがIGLはCGALかCORKと繋げないと動かないのですが，CGALの方はLinuxでは(コンパイラが古すぎて？)動かない．CORKの方は計算がおかしい．という問題がありますね．
うさちゃん
君は，うさぎを見たか．
見たことがなければ，ここの例を実行．面倒な人はダウンロードしてください．Macですので, libglad.a, libglfw.a の他に IOKit/CoreVideo/Cocoa/GLUT/OpenGL のframeworkを追加して再生する：
うさちゃん
まあこれで，いちいちBlenderを使わなくても, OFFやOBJなどの3Dモデルデータをプレビューできるわけですね．
で, うさぎの構成要素は血と肉ではなく, 座標と面であるわけです．座標は行列 Eigen::MatrixXd V です.
Eigen::MatrixXd(行数, 列数)
列数を3にすると3次元座標で, 点の数だけ行を設定します．うさちゃんの面は三角形なので, 点数3と, 三角形を構成する座標点の点番号[0 行番号-1]が羅列した整数配列 igen::MatrixXi F ですね．
まあこれですね．バカじゃねえんだから，Fで座標点の並べ順は正規なように並べろよ，ということでした．だからFの1行目が0 1 2 なら2行目は1 2 3 ではダメで，2行目は, 1 3 2 ではないといけねえのよ．でもこれで，わかりますよね，立体だって何だって描けますわ．
三角錐
まずは三角錐ですよね.
void make_tetrahedron(Eigen::MatrixXd& V,Eigen::MatrixXi& F){
int number_of_vertices=4;
int number_of_faces=4;
V.resize(number_of_vertices,3);
V(0,0)=0.; V(0,1)=0.; V(0,2)=0.; //0
V(1,0)=1.; V(1,1)=0.; V(1,2)=0.; //1
V(2,0)=0.; V(2,1)=1.; V(2,2)=0.; //2
V(3,0)=0.; V(3,1)=0.; V(3,2)=1.; //3
F.resize(number_of_faces,3);
F(0,0)=0; F(0,1)=1; F(0,2)=3;
F(1,0)=1; F(1,1)=0; F(1,2)=2;
F(2,0)=2; F(2,1)=0; F(2,2)=3;
F(3,0)=1; F(3,1)=2; F(3,2)=3;
}
できまぴた：
で, ファイルへの保存は
#include <igl/writeOFF.h>
#include <igl/writeOBJ.h>
....
void save_obj(std::string fname,Eigen::MatrixXd& V,Eigen::MatrixXi& F){
igl::writeOBJ(fname,V,F);
};
void save_off(std::string fname,Eigen::MatrixXd& V,Eigen::MatrixXi& F){
igl::writeOFF(fname,V,F);
}
こんな感じでまあ十分でしょう．
基本的には3角形の集合を記述するので, 立方体でもこんな感じです：
void make_cube(Eigen::MatrixXd& V,Eigen::MatrixXi& F){
int number_of_vertices=8;
int number_of_faces=12;
V.resize(number_of_vertices,3);
V(0,0)=0.; V(0,1)=0.; V(0,2)=0.; //0
V(1,0)=0.; V(1,1)=0.; V(1,2)=1.; //1
V(2,0)=0.; V(2,1)=1.; V(2,2)=0.; //2
V(3,0)=0.; V(3,1)=1.; V(3,2)=1.; //3
V(4,0)=1.; V(4,1)=0.; V(4,2)=0.; //4
V(5,0)=1.; V(5,1)=0.; V(5,2)=1.; //5
V(6,0)=1.; V(6,1)=1.; V(6,2)=0.; //6
V(7,0)=1.; V(7,1)=1.; V(7,2)=1.; //7
F.resize(number_of_faces,3);
F(0,0)=0; F(0,1)=6; F(0,2)=4;
F(1,0)=0; F(1,1)=2; F(1,2)=6;
F(2,0)=0; F(2,1)=3; F(2,2)=2;
F(3,0)=0; F(3,1)=1; F(3,2)=3;
F(4,0)=2; F(4,1)=7; F(4,2)=6;
F(5,0)=2; F(5,1)=3; F(5,2)=7;
F(6,0)=4; F(6,1)=6; F(6,2)=7;
F(7,0)=4; F(7,1)=7; F(7,2)=5;
F(8,0)=0; F(8,1)=4; F(8,2)=5;
F(9,0)=0; F(9,1)=5; F(9,2)=1;
F(10,0)=1; F(10,1)=5; F(10,2)=7;
F(11,0)=1; F(11,1)=7; F(11,2)=3;
}
で，出来上がりは
Eigenでは, 行列Vの行数は V.rows() 列数は V.cols() で取得可能です．列ベクトルや行ベクトルは V.row(i)とかV.col(j)ですよ.
for(int i=0;i<V.rows();i++) {すると次の結果が出ます．流石に, ベクトルを立てベクトルで出力するようにはなってねえ?
auto x=V.row(i); std::cout << x << std::endl;
}
Column 8 x 3
0 0 0
0 0 1
...
いやいや, それは貴様がボケているのだ．列ベクトルを出したのだから，これは横ベクトルですよね．transpose()しねえとダメだボケ：
for(int i=0;i<V.rows();i++) {auto x=V.row(i).transpose(); std::cout << x << std::endl << "--" << std::endl;}
これで
0
0
2
--
0
0.5
0
--
こんな感じに出力される．
座標点が線形代数ライブラリEigenに入っているので, まあ, Affine変換は普通に計算すれば良いですね.
Eigen::Affine3d SC(Eigen::Scaling(w[0],w[1],w[2])); //スケール変換部分
SC.translation()=Eigen::Vector3d(fr[0],fr[1],fr[2]); //平行移動部分
SC.linear().col(1)=vv; //変換行列は linear()で取得でき
for(int i=0;i<V_u.rows();i++) //V_uの頂点を変換
V.row(i)=SC*Eigen::Vector3d(V_u.row(i));る
で，こうなる：
Boolean Operations (Cork Version)
CGALのNef図形の論理演算をIGLで行うことができる．
- 本格的なのが嫌いな人に！
- 計算が爆速なだけで，おかしいじゃねーか，というプロ人が多い
- おかしい理由は無限集合の扱いがダメだ，ってことらしい．
- じゃ，十分なのでは．
#include <igl/copyleft/cgal/mesh_boolean.h>
...
igl::MeshBooleanType boolean_type=static_cast<igl::MeshBooleanType>(0);
igl::copyleft::cork::mesh_boolean(VA,FA,VB,FB,boolean_type,V,F);
0が演算を表す．以下のものが可能である
| 0 | 1 | 2 |
| Union | Intersect | Minus |
|
|
 |
|
| 3 | 4 | |
| XOR | RESOLVE |
で，コンパイルするときには libcork.a をりんくするのぢゃよ
Corkは計算は早い. CGALの数十倍である. しかし誤差は・・・有効数字2-3桁がいいところである. 計算できないはずがない場合でも, なんかエラーする.
Boolean Operations (CGAL Version)
CGALのNef図形の論理演算をIGLで行うことができる
- Macではできたが, Linuxでは, libigl-cgal の結合部分がおかしくてビルドできねえ
- だから，まあ普通に使っていこう，と言う人向きではない
- 開集合（無限集合）にも強い（いらねえ〜〜〜〜機能だ）
#define IGL_STATIC_LIBRARY
#include <igl/copyleft/cgal/mesh_boolean.h>
...
igl::MeshBooleanType boolean_type=static_cast<igl::MeshBooleanType>(0);
igl::copyleft::cgal::mesh_boolean(VA,FA,VB,FB,boolean_type,V,F);で, コンパイルするときにはこんな感じ：
どういうわけか CGALは libigl_cgal.a に取り込まれている．
にもかかわらず, mpfr とか gmp は取り込まれていない．ので， dylib が必要である．
にゃお，libboost_thread も内部で使っているみたいだよ．
三角分解
IGLは, 論理演算はCGALの借り物で高速に計算を行った. 三角分解ももしや借り物か？というと, そのトーリ. Tetgenを利用する．そうだ(マニュアル)．こいつは
- PLC(直線の辺を持つセルから構成される3Dモデル)を三角化する.
- 主に頂点を使うが, Steiner点を追加しなければならないこともある. Steiner点を追加しても三角分解できないPLC領域は存在しないことが証明されている
- 境界面や曲率に合わせてメッシュを高密度に分布させる（いらねーよ等分割しろよ）
で，こんな感じで
#define IGL_STATIC_LIBRARY
#include <igl/copyleft/tetgen/tetrahedralize.h>
...
void make_tetra(Eigen::MatrixXd V,Eigen::MatrixXi F,
Eigen::MatrixXd& TV,Eigen::MatrixXi& TT,Eigen::MatrixXi& TF){
std::string SW;
SW="pcQ";//PLC input, keep convex hull, quiet
igl::copyleft::tetgen::tetrahedralize(V,F,SW.c_str(), TV,TT,TF);
}
あーそれから，どういうわけか
- libiglでは内部で #include "tetgen.h" とかボケた書き方をしているので, XCodeのHeader Search Path に /usr/local/include/tetgen を追加してくだせえ
- libigl_tetgen.a(IGL→tetgenバインド), libtet.a (tetgen本体)のリンクじゃボケ
で，こうなる
凸物体よし!
ではこれは？ 
アウト！！！！！
やっぱ，事前にConvex Decomposition が必要なのだ. そらなー
でも, tetgenの使用例では非凸物体も三角分解している．使い方の問題か．と思ったら，違う考え方があるのである．
サーフェス縛り
これは，「前にサーフェスであったものは，後でもサーフェスでなければならぬ」という拘束プレイである．なるほど，元々閉領域であったのだから，この縛りの下では穴が埋まることができない．
SW="pYQ";//PLC input, keep surface, quiet
igl::copyleft::tetgen::tetrahedralize(V,F,SW.c_str(), TV,TT,TF);
これで
おう！無事だ．
体積
体積ももちろん計算できる：
#include <igl/volume.h>
#include <igl/volume.cpp> // インライン関数が定義されている
double make_volume(Eigen::MatrixXd& TV,Eigen::MatrixXi& TT){
Eigen::MatrixXd VOL;
igl::volume(TV, TT, VOL);
double sum=0.;
for(int i=0;i<TT.rows();i++) sum+=VOL(i,0);
return sum;
};
これで上の例を計算すると:
|A| = 1, |B|=1.25 |A-B|=0.776574, |B-A|=1.02657, |A⋂B|=0.223426, |A⋃B|=2.02657
合ってるのか？
|A⋂B| + |A-B| = |A| はOK, |A⋂B| + |B-A| = |B| はOK, |A⋂B|+ |A-B| + |B-A| = |A⋃B| もOK.
できたわ．
| 添付 | サイズ |
|---|---|
| 06 - libigl.pdf | 14.15 MB |
| TetGen マニュアル | 2.14 MB |