メインコンテンツに移動

基本概念

CGALの基本概念

カーネル

幾何形状は, 距離空間の中の点と線から構成されます. 次元，座標系や測度によって様々なものがありますが，CGALは，それをカーネルサンダースと呼びます．代表例は

  • CGAL::Simple_cartesian<double>  普通の直交座標系ですねー
    • 欠点は, 無限遠を識別できないところです．遠方はええ加減な表記しかできません．
  • CGAL::Simple_homogeneous<long> 同次座標系を定義
  • CGAL::Exact_predicates_exact_constructions_kernel
    • 表記はみんな大好きCartesian座標系です
    • 座標数値は，みんな大好きdoubleです
    • exact geometric predicates と exact geometric constructions をサポート
  • CGAL::Exact_predicates_inexact_constructions_kernel
    • 表記はみんな大好きCartesian座標系です
    • 座標数値は，みんな大好きdoubleです
    • exact geometric predicates と inexact geometric constructions をサポート

・・・で, predicate と construction とは．

predicate: 多分, 座標点の内部表現みたいなやつだな. 

construction: なんか要素を生み出す関数のこと. 例えば distance() は距離の2乗を生む. sqrt(distance()) もあった方が良い? いやいや, それは丸め誤差を生み出す inexact な construction だ, と考える. らしい.

カーネルを構成する要素は

Point_2 Point_3 (x,y)座標を表す
Vector_2 Vector_2 ベクトルなんじゃないでしょうか
Direction_2 Direction_3 方向なんじゃあねえかな，普通に考えると方向以外考えられない
Line_2 Line_3 無限直線
Ray_2 Ray_3 半無限直線
Segment_2 Segment_3 線分

ポリゴン

ポリゴンとは多角形です．多面体ではありません．多角形とは, 一連の線分(edge)が閉ループを構成するものです．

ポリヘドラ

ポリヘドラとは多面体です．多面体の構成要素は

  • vertices 頂点ですよね．
  • halfedge 辺(edge)ですが, 行きと帰りの2通りを考えることが多いので, セットにした片割れを edge の半分, halfedge と呼びます．
  • facet 面です．faceではありません(なんのこっちゃ)

三角錐 tetraphedron

三角錐を作りましょう. 一発で作れるのは三角錐と三角形だけなんですけどね:

Kernel::Point_3 p( 1.0, 0.0, 0.0);
Kernel::Point_3 q( 0.0, 1.0, 0.0);
Kernel::Point_3 r( 0.0, 0.0, 1.0);
Kernel::Point_3 s( 0.0, 0.0, 0.0);
CGAL::Polyhedron_3<Kernel> P;
P.make_tetrahedron( p, q, r, s);
//登録した順にイテレータが使える
for (auto v = P.vertices_begin(); v != P.vertices_end(); ++v)
   std::cout << v->point() << std::endl;

三角錐の辺は6本ですので, halfedgeは12本あります. halfedgeの逆方向のやつは, opposite()で取得できます．vertex()->point()で行き先がわかりますので, 結局, どの点からどの点に向かうhalfedgeであるかをゲットできます：

int n=0;
for (auto h = P.halfedges_begin(); h != P.halfedges_end(); ++h) {
   std::cout << n << "th halfedge " << h->opposite()->vertex()->point()
                                     << " --> " << h->vertex()->point() <<  std::endl;
   n++;
}

上記のhalfedgeのイテレータは, 全てのhalfedgeが一覧されますが, 幾何学的形状を反映しているわけではありません．H0の反対向きがH1, H0の次がH2, H2の次がH4,...ですが, H8からH10までジャンプしてしまいます．

tetra1

このために, halfedgeのnext()があります．next()は, ある面を巡回するようにhalfedgeを進めます．opposite()すると, 別の面を巡回し始めます．

tatra2

面を一覧することもできます:

auto triangle_area=[](const CGAL::Polyhedron_3<Kernel>::Facet& f){
   return Kernel::Compute_area_3()(f.halfedge()->vertex()->point(),
       f.halfedge()->next()->vertex()->point(),f.halfedge()->opposite()->vertex()->point() );
};
auto triangle_normal=[](const CGAL::Polyhedron_3<Kernel>::Facet& f){
   return Kernel::Construct_normal_3()(f.halfedge()->vertex()->point(),
   f.halfedge()->next()->vertex()->point(),f.halfedge()->opposite()->vertex()->point() );
};
n=0;
for (auto f = P.facets_begin(); f != P.facets_end(); ++f) {
   std::cout << "Facets " << ++n;
   if (f->is_triangle())
     std::cout << " is triangle area normal to " << triangle_normal(*f)
     << " with area= " << triangle_area(*f);
   std::cout << std::endl;
}

この例では, 面要素が三角形であれば, 法線ベクトルと面積を出力します.

三角錐に限定した CGAL::Tetraphedrom<Kernel> には, 体積を求めるvolume()があります. これはP1-P2-P3平面に対しP4が正側にあるかどうかで正負がついた体積を戻しますので注意（上の図では, P1-P2-P3の循環に対して右ねじ方向の逆側にP4があるので, マイナス三角錐ですね）：

std::vector<Kernel::Point_3> ARRY;
for (auto v = P.vertices_begin(); v != P.vertices_end(); ++v) ARRY.push_back(v->point());
CGAL::Tetrahedron_3<Kernel> T(ARRY[0],ARRY[1],ARRY[2],ARRY[3]);
std::cout << "Volume= " << std::abs(T.volume()) << std::endl;

多面体の追加工

作った多面体はEuler操作とかいうので加工できます．

  • split_facet
  • join_facet
  • split_vertex
  • join_vertex
  • split_loop
  • join_loop

では先ほどの三角錐を変形してみましょう：

auto h=halfedge;//Polyphedron_3を作った時, 辺を一つ返却してくるのでとりあえずそれな
auto g=h->next()->opposite()->next();
P.split_edge(h->next());P.split_edge(g->next());P.split_edge(g);
h->next()->vertex()->point()=Kernel::Point_3(1,0,1);
g->next()->vertex()->point()=Kernel::Point_3(0,1,1);
g->opposite()->vertex()->point()=Kernel::Point_3(1,1,0);
print_vertex(P);
print_facets(P);

追加工は別料金

辺を分割すると頂点は自動的に増えますが, 面の数は不変です．折れ曲がった異常な面が出来上がります：

これがfacet

これは, halfedgeを2個指定してsplit_facetすれば良い. これは, 第一引数の終点と第二引数の終点を結ぶように面を分割し, 第一引数のnextとなるhalfedgeを返却する. というわけで, 新しいhalfedgeの根元に頂点を作成して, 残りの角の値に設定すれば良い．

auto f=P.split_facet(g->next(),g->next()->next()->next());
auto e=P.split_edge(f);
e->vertex()->point()=Kernel::Point_3(1,1,1);;

a

もう一度 split_facet()すれば, 立方体が構成できる．

P.split_facet(e,f->next()->next());

chromeのクソが手間とらせやがって