等長リメッシング (Isotropic Remeshing)†
三角形サーフェスメッシュを、辺の長さがどこでもほぼ等しく、正三角形に近い三角形だけで構成されたメッシュに作り直す操作を等長リメッシング(isotropic remeshing)と呼びます。
CADから出力したSTL、医用画像やボクセルデータにマーチングキューブスを適用して得たメッシュ、スキャンデータなどは、三角形の大きさと形が極端に不均一になりがちです。このようなメッシュは、そのままでは数値解析にも後段のメッシュ操作にも使いにくいため、形状を保ったまま三角形を作り直します。
なぜ等長メッシュが必要か†
- 細長く潰れた三角形は補間誤差を大きくし、剛性行列の条件数を悪化させます
- 陽解法では最小辺長が時間刻みの上限を決めるため、1本の極端に短い辺が計算全体の速度を支配します
- スムージング、細分割、デローニー分割法による体積メッシュ生成など、後段の処理は入力の三角形品質にそのまま影響されます
- 可視化やレイキャストでも、潰れた三角形は法線が不安定になり陰影に破綻が出ます
アルゴリズム†
Botsch と Kobbelt が提案した incremental remeshing を用います。目標辺長 L をひとつ決め、次の5段階を1反復として繰り返します。
1. 分割 (Split)†
L の 4/3 より長い辺を中点で分割します。粗すぎる領域に点を増やす操作です。
2. 統合 (Collapse)†
L の 4/5 より短い辺を中点に潰します。ただし、潰した結果 4/3 L を超える辺ができる場合と、メッシュの位相が壊れる場合は実行しません。密すぎる領域から点を減らす操作です。
4/5 と 4/3 という2つの閾値は目標辺長を挟むように置かれており、分割と統合が互いの結果を打ち消し合って振動することを防いでいます。
3. 反転 (Flip)†
隣り合う2つの三角形が作る四辺形の対角線を張り替え、頂点の価数(その頂点に集まる辺の数)が6に近づく方向に辺を反転します。境界上の頂点では4を目標にします。正三角形が敷き詰められた状態の価数が6であることに対応します。
4. 接平面方向の平滑化 (Tangential relaxation)†
各頂点を、周囲の頂点の重心(一様ラプラシアン)へ向けて動かします。ただし移動量を頂点の接平面へ射影してから適用するため、形状を膨らませたり縮めたりせずに、面内での点の配置だけを均します。
5. 元形状への投影 (Projection)†
平滑化で面から浮いた頂点を、元の入力サーフェス上の最近点へ投影して形状を戻します。最近点探索には vtkStaticCellLocator を用います。
弊社の実装 vtkIsotropicRemeshingFilter†
上記のアルゴリズムを VTK 向けのフィルタとして実装しました。vtkPolyDataAlgorithm を継承しているため、通常のVTKフィルタと同じようにパイプラインに接続できます。
論文のみに基づく独自実装で、外部のジオメトリ処理ライブラリには依存していません。内部ではハーフエッジ構造を構築し、分割・統合・反転をすべて局所操作として処理します。平滑化と投影は vtkSMPTools で並列化しています。
ソースコードはGitHubで公開しています。ライセンスはBSD-3-Clauseです。
パラメータ†
| パラメータ | 既定値 | 説明 |
| TargetEdgeLength | 0 (自動) | 目標とする辺の長さ。0以下を指定すると入力メッシュの平均辺長を使います |
| NumberOfIterations | 3 | リメッシュの反復回数。増やすほど均一になります |
| NumberOfRelaxationSteps | 1 | 1反復あたりの接平面平滑化の回数 |
| DoProject | true | 元のサーフェスへ投影して形状を保つかどうか |
| InterpolateAttributes | true | 入力の点データ・セルデータを出力へ補間するかどうか |
| ProtectFeatures | false | 鋭い稜線を保護するかどうか |
| FeatureAngle | 60度 | 稜線と判定する二面角の閾値。ProtectFeatures が有効なときに使います |
使い方†
#include "vtkIsotropicRemeshingFilter.h" vtkNew<vtkIsotropicRemeshingFilter> remesh; remesh->SetInputData(surface); // vtkPolyData remesh->SetTargetEdgeLength(0.5); // モデルの大きさに合わせて決める remesh->SetNumberOfIterations(5); remesh->ProtectFeaturesOn(); // CAD形状の稜線を残したいとき remesh->SetFeatureAngle(60.0); remesh->Update(); vtkPolyData* result = remesh->GetOutput();
コマンドラインツール†
ファイルを直接リメッシュする RemeshSurface も用意しています。
RemeshSurface input.{vtp|stl|ply|obj} output.{vtp|stl|ply} [目標辺長] [反復回数] [特徴角度]
目標辺長を省略するか0にすると入力の平均辺長を使います。特徴角度を与えると稜線の保護が有効になります。
効果の確認†
形状の比較†
vtkSphereSource で作った半径1の球です。緯度経度に沿ってテセレーションされるため、極の周辺に細長い三角形が集中し、赤道付近との大きさの差も大きくなっています。これを目標辺長0.08、5反復でリメッシュしました。
極の特異点が解消され、球面全体がほぼ同じ大きさの三角形で覆われています。
メッシュ品質の定量評価†
上の例の数値です。
| 指標 | 入力 | リメッシュ後 |
| 点数 | 1,058 | 2,125 |
| 三角形数 | 2,112 | 4,246 |
| 辺長 最小 | 0.0178 | 0.0443 |
| 辺長 平均 | 0.1294 | 0.0827 |
| 辺長 最大 | 0.1888 | 0.1153 |
| 辺長 標準偏差 | 0.0401 | 0.0078 |
| 辺長 最大/最小 | 10.60 | 2.60 |
| 三角形の最小角 | 7.4度 | 36.1度 |
| 三角形の最小角の平均 | 31.2度 | 54.2度 |
辺長のばらつきは標準偏差で約5分の1になり、最大と最小の比は10.6倍から2.6倍まで縮まりました。三角形の最小角は最悪値が7.4度から36.1度へ改善し、平均は正三角形の60度に近い54.2度になっています。
入力の分布が離散的な棒状になっているのは、緯度経度に沿った規則的なテセレーションだからです。リメッシュ後は目標辺長を中心とした単峰の分布に変わっています。
鋭い稜線の保護†
平滑化と投影は形状を滑らかにする方向に働くため、そのままでは立方体のような角張った形状の稜線が削れてしまいます。ProtectFeatures を有効にすると、二面角が FeatureAngle を超える辺を境界と同じ拘束付き扱いにします。この辺は反転されず、辺上の頂点は稜線に沿ってしか動かず、稜線の交点は固定されます。
一辺2の立方体を3回細分割したメッシュを、目標辺長0.15、5反復でリメッシュした結果です。
保護なし(中央)では稜線と頂点が丸まって形状が崩れています。保護あり(右)では稜線と8つの角が保たれたまま、各面が等長メッシュになっています。
境界と属性の扱い†
境界の保存†
開いたサーフェスの境界辺は、設定によらず常に拘束されます。境界上の頂点は境界の折れ線に沿ってのみ移動し、幾何学的な角にあたる頂点は固定されます。板や薄膜のような形状をリメッシュしても外形が変わりません。
属性の補間†
InterpolateAttributes が有効(既定)のとき、入力が持つ点データとセルデータを出力に引き継ぎます。点データは、元のサーフェス上の最近点における重心座標で補間します。セルデータは、出力三角形の重心に最も近い元セルからコピーします。
最終的な投影によって出力頂点は元のサーフェス上に載るため、点データの補間は事実上、厳密な線形補間になります。実際、z座標を点スカラーとして与えた球で確認したところ、補間誤差の最大値は 2.22e-16 でした。
計算時間†
球を入力として、目標辺長を変えながら5反復で計測しました。
| 入力 | 目標辺長 | 出力 | 時間 |
| 1,682点 / 3,360三角形 | 0.05 | 5,628点 / 11,252三角形 | 35.7 ms |
| 19,602点 / 39,200三角形 | 0.02 | 34,441点 / 68,878三角形 | 405 ms |
| 79,202点 / 158,400三角形 | 0.01 | 141,593点 / 283,182三角形 | 1.83 s |
| 318,402点 / 636,800三角形 | 0.005 | 568,323点 / 1,136,642三角形 | 7.55 s |
最大ケースでの段階別の内訳です。
| 段階 | 時間 |
| ハーフエッジ構築 | 798 ms |
| ロケータ構築 | 34 ms |
| 分割 | 830 ms |
| 統合 | 2,111 ms |
| 反転 | 1,546 ms |
| 平滑化 | 191 ms |
| 投影 | 2,044 ms |
| 合計 | 7,554 ms |
計測環境は AMD Ryzen 7 2700X (8コア16スレッド)、メモリ32GB、VTK 9系のReleaseビルドです。統合と投影が支配的で、出力100万三角形規模でも10秒以内に収まっています。
回帰テスト†
実装には次の4つのテストを用意しています。
- 球: 異方的にテセレーションされた球をリメッシュし、閉じた多様体の三角形メッシュであること、辺長がほぼ均一であること、面積が保たれることを確認します
- 平面: 矩形板をリメッシュし、境界の折れ線と4隅が保存されること、全頂点が平面内に留まることを確認します
- 立方体: 特徴保護を有効にしてリメッシュし、全頂点が立方体表面上にあること、8つの角が保存されることを確認します
- 属性補間: 点スカラーとセルスカラーを与え、補間結果が元の値を再現することを確認します
制限事項†
- 入力は多様体で、かつ向きが一貫している必要があります。向きが揃っていない場合は事前に vtkPolyDataNormals を ConsistencyOn で通してください。三角形以外のセルは内部で三角形化し、重複点は統合します
- 非多様体な辺を含む入力は扱えません
- 目標辺長は全体で1つの値です。曲率に応じて粗密を変える適応的なリメッシュは行いません
- 投影によって頂点は元の面上に載りますが、区分的に平面である以上、面積はわずかに変化します。半径1の球の例では 12.3641 が 12.3004 になり、0.5%程度の減少でした
参考文献†
- M. Botsch, L. Kobbelt, "A Remeshing Approach to Multiresolution Modeling", Symposium on Geometry Processing, 2004.
- Kenichiro-Yoshimi/IsotropicRemeshing (ソースコード)
- VTK
- ParaView


