Level Set 法 3次元への適用

in English

はじめに

以前こちらでLevel Set法を紹介しました。今回これを3次元に拡張したのであらためて紹介します。
（プログラムにバグがありましたので以下の文章を書き直しました。この修正により、計算速度が大幅に改善されました。2013年1月20日）

デモ

こんな感じです。
デモ1:

2つの球体をその内部に持つように立方体（初期波面）を配置します。Level Set法を実行すると、立方体が徐々に縮まり、2つの球面に分裂する様子を見ることが出来ます（画面キャプチャなし視点変更なしなら、70秒ほどで計算を終えます）。
デモ2:

立方体（初期波面）を2つの球体の外に置いた場合です。立方体は球体を呑み込みながら成長していき、ワイヤフレームで示した大きな球面に衝突して止まります。言い換えると、小さな立方体は2つの球面を内部に含む大きな球面に成長します（画面キャプチャなし視点変更なしなら、120秒ほどで計算を終えます）。

どちらのデモも空間の大きさは200x200x200pixlesです。描画にはOpenGLを使いました。

開発環境

1. Mac OS X 10.8.2
2. プロセッサ：3.06 GHz Intel Core 2 Duo
3. メモリ：4GB
4. Xcode4.5.2 with Apple LLVM 4.1（C++ Language Dialect → C++11, C++ Standard Library → libc++）
5. boost-1.51.0（上記コンパイラでコンパイルしたもの。こちら）
6. opencv-2.4.3（上記コンパイラでコンパイルしたもの。こちら）

注意：lambda式やthreadなど、C++11の新しい機能を使用しています。

ソース

こちらです。上記の開発環境で動作確認しました。

実行方法

実行ファイ名は、LevelSetMethod2です。引数なしで実行すると以下を出力します。 各引数の意味は以下の通りです。

1. dim:次元です。2か3を指定します。
2. verbose:途中経過を出力したいときに指定します。
3. input: 2次元のとき入力ファイルへのパスを、3次元のとき"pattern0"か"pattern1"を渡します。3次元の場合はプログラム内で既定の立体を作成します。任意の立体を与えられるようにする予定です（こちらへ）。
4. wband: Narrow Band Level Set法を採用しています。波面（境界面）の計算領域の幅を設定します（後述）。
5. wreset: Level Set法の精度を高めるため、ある計算ステップごとに波面の初期化を行います。その引き金となる領域の幅を設定します（後述）。
6. time_step: Level Set法では波面を時間発展させます。その時間幅を設定します。
7. gain: 波面の速度は定数項と波面の曲率に依存する項から成ります。後者の割合です（後述）。
8. constant_speed: 波面速度の定数項です（後述）。
9. speed_threshold: この速度以下の場合、波面速度を0とみなします。
10. left, right, top, bottom, front, back: 最初に与える立方体のサイズを指定します。frontとbackは３次元の場合に使います。

以下の図はwbandとwresetを示したものです。

wbandのうち緑で示した領域の幅がwresetです。波面が緑の領域に到達すると、波面の再初期化が実行されます。 本アプリで使用している速度関数 は以下で定義されます。

ここで、 は曲率です。constant_speedは に相当します。閉曲面を内側に向かって収縮させるなら負の値を、その逆なら正の値を設定します。また、gainは に相当します。この項は曲面を滑らかにする効果を持ちます。

①３次元の場合:
以下のコマンドを実行します。最初に示したデモ1が始まります。 ここでキーボードから以下の操作が可能です。

1. ESC: 終了します。
2. p: 一時停止します。
3. f: 波面(front)の描画を有効・無効にします。
4. o: 物体(object)の描画を有効・無効にします。
5. <: 視点を物体に近づけます。
6. >: 視点を物体から遠ざけます。
7. u: x軸周りに回転します（右ねじ）。
8. n: x軸周りに逆回転します。
9. h: y軸周りに回転します（右ねじ）。
10. j: y軸周りに逆回転します。

引数に--verboseを指定したので、コマンドライン上に波面（front）を構成するピクセルに関する量が表示されます。 その意味は以下の通りです。

1. zero speed rate: 成長速度が0となったピクセルの割合。1のとき波面の成長は完全に停止します。
2. front size: 全ピクセル数。
3. total speed: ピクセルの成長速度の絶対値の和。

本アプリは、total speedが３回連続して同じ値になると終了します。

②２次元の場合:
以下のコマンドを実行します。 引数に--verboseを指定しない場合、波面に関する上記3つの量は表示されません。 なお、入力画像は１チャンネルのグレー画像でなければなりません。こんな感じになります( 画面キャプチャなしなら0.8秒程で計算を終えます)。画像サイズは256x256pixelsです。（追記：こちらに２次元のデモを追加しました。）

ここで、キーボードから以下の操作が可能です。

1. ESC: 終了します。
2. p: 一時停止します。
3. f: 波面(front)の描画を有効・無効にします。
4. o: 物体(object)の描画を有効・無効にします。

ソース概説

main関数は以下の通りです。
//main.cpp クラスCommandLineInterfaceには以下2つのstatic関数が定義されています。

1. execute_level_set_method_in_2d
2. execute_level_set_method_in_3d

例えばexecute_level_set_method_in_3dではクラスLevelSetMethodViewer3dを呼び出しています。
//CommandLineInterface.cpp 9行目でクラスLevelSetMethodViewer3dのインスタンスを作成し、そのあと、viewerの初期化（19行目のinitialize_viewer）とLevel Set法の初期化（20行目のinitialize_level_set_method）を行っています。後者の関数の中身は以下の通りです。
// LevelSetMethodViewer.h 13行目にあるクラスLevelSetMethod<Dimension>がLevel Set法を実行します。これはテンプレートクラスです。
// LevelSetMethod.h テンプレート引数Dには次元を表現するクラスを与えます。
// DimensionTypes.h

参考文献

1. コンピュータビジョン 最先端ガイド1, 第1章
2. Fast Level Set Methodの提案とビデオ画像の移動物体のリアルタイム追跡, 情報処理学会論文誌, Vol.44, No.8, Aug, 2003
3. Fast Level Set Methodを用いた複数移動物体の三次元追跡, 日本ロボット学会誌, Vol.23, No.7, pp.813-820, 2005
4. Level Set Methods and Fast Marching Methods, J.A.Sethian, Cambridge University Press, ISBN 0-521-64557-3

追記

こちらにも2Dへの適用例を掲載しました。