5 . 1 序論
これまでの変形可能曲面モデルは，球面と同相な物体への利用を前提としており，物体の 統合や分際によって生じる位相変化には対応できない.したがって，位相に関する事前知識 をもつことなく任意の位相に適応して柔軟に全周型モデノレを生成する手法が望まれる.
ところで，数値解析の分野では， OsherとSethian[58，59]が，平均曲率流方程式に従って 動く曲線(面)の運動を解析するため，等高面法(LevelSet Methods)を提案した.等高面法 の基本的な考え方は，時間の経過とともに動く曲面を一つ高い次元で定義されたある補助関 数のゼロ等高面とみなすことである.この補助関数の運動方程式を用いると，元の曲面の位 相を陽に扱う必要はなくなるため，曲面の位相に依存しない時間大域的な追跡が可能となる.
Malladiら[35]は，この等高面法を輪郭線抽出問題に応用し，複数の領域や穴を含む形状の 動的輪郭線抽出を行なった.そこで，位相適応型の三次元モデリングを実現するため，この 動的輪郭線抽出手法を三次元に拡張することが考えられる.
このときに問題となるのは，運動方程式に従って動く曲面が実曲面付近を通過しつつある とき，曲面の移動をいかにして停止させるかということである.上述の動的輪郭線抽出では，
曲面点が対象輪郭を越えるのを避けるため，濃度勾配が急激に変化するところに隙間のない 連結した停止領域を設けている.ところが，多視点距離画像のような離散データに対し，こ
のような停止領域を確保することは困難である.
本章では，複数レベノレの手続きで曲面を追跡することにより任意の位相に適応しながら全 周型物体モデリングを行なう手法を提案する.まず，各データ点から比較的遠い距離内にあ る点をすべて停止点とすることで，隙間のない停止領域を設定する.そして，曲面を等速成 長させることにより，停止領域への曲面当てはめを行なう.次に，停止領域を狭くしながら，
この曲面当てはめの過程をモデル曲面が実物体と同相になるまで繰り返す.最後に，上で得 られた曲面を改めて初期曲面とし，データとの近接性に依存する速度で曲面を移動させるこ とにより，原形状に忠実な物体モデリングを達成する.
本章では，まず，等高面法の概要について述べる.次に，三次元曲面再構成問題に等高面 法を応用する際に生じるデータ隙間の問題を指摘し，この問題を解決した位相適応型物体モ
5.2. 等高面法の基礎理論と位相適応性 70
デリングについて述べる.最後に，合成データを用いた実験の結果を示し，提案する手法の 有効性を明らかにする.
5 . 2 等高面法の基礎理論と位相適応性
等高面法で三次元曲面を扱うため，与えられた初期閉曲面f(t
= 0 )
に対し，速度F
で法 ベクトル方向に動く曲面f(t)の運動方程式を導出する.ここで，Fは曲率や法ベクトル方向 などを変数とする速度関数である.等高面法の基本的な考え方は，時間の経過とともに動く 曲面f(t)をより高次の関数 u =ゆ(x，
y，
z，
t)のゼロ等高面とみなして，それを埋め込むことで ある.いま，三次元位置x=(x
，
y，
z)に対し，補助関数ゆ(x，t = 0)をゆ(x
，
t = 0) =士d (5.1)と定義する.ここで， d はxから初期曲面f(t= 0)までの最短距離である.また，右辺の符号 は， xが関曲面f(t
=
0)の外{ftlJにあるときに正，内側にあるときに負の符号がとられる.こ れにより，f(t
=
0)=
{x 14>(x，
t=
0)=
O} (5.2) を満足する補助関数ゆ( x
，t=
0) :冗3→
7之が得られる.この簡単な例として，三次元閉曲面f(t)を二次元閉曲線 γ(t)におきかえてみる.初期閉 曲線γ(0)としてX Y平面上の円が与えられたとき (Fig.5.1(a))，これを式(5.2)のように，
三次元空間7之3における連続な初期関数 z=ゆ(x
，
y，
t=
0)のゼロ等高線{ゆ=O}とみること ができる (Fig.5.1(b)).そして，この閉曲線が時間の経過とともに動くと (Fig.5.1(c))，こ れをゼロ等高線{ゆ=O}にもつ関数z=ゆ(x，
y， t )
も変形していく (Fig.5.1(d)).本節の目的は， f(りをゼロ等高面にもつ補助関数ゆ(x，t)の方程式を導出することである.
いま，x(t)を曲面 f(t) 上のある点が曲面変形に伴って移動していく経路とすると，初期位 置x(t
=
0)は初期曲面f(t=
0)上の点であり，関数ゅのゼロ等高面が曲面f(t)と常に一致す るための必要条件はゆ(x(t)
，
t) = 0 (5.3) と表される.式(5.3)の両辺をtで微分してゆ
t+
マゆ( x(
t)，
t) .x ' (
t) = 0 (5.4)5.2. 等高面法の基礎理論と位相適応性
y
x
y
z=φ(x， y，七=0)
z=φ(x， y，七)
y
γ(七)=Level Setφ=0
X
Fig.5.1:等高面法による運動方程式の定式化
71
5.2. 等高面法の基礎理論と位相適応性 72
Z z=φ(x，y，七=0) Z z=φ(x，y，七)
x
，
y= = >
¥ γ ! 七 ) ヌ
X(c)
( d )
Fig.5.2:等高面法による位相変化への適応
5.3. 数値解法 73
となる.ここで，外向き法ベクトル方向の移動速度Fは，三次元位置
x ( t )
の移動速度x ' ( t )
と曲面
r
の外向き単位法ベクトル n =マゆ1 /
マゆ!とを用いて F=x ' ( t ) .
nと表せる.よって，ゅに関する運動方程式(5.4)は
ゆt
+
FIY'利 =0 (5.5)と表せる.この方程式(5.5)を等高面方程式と呼ぶ[58，59].
等高面法の好ましい特徴の一つは，速度関数Fが滑らかのとき，曲面f(t)が変形の途中 で位相変化を生じた場合でも，曲面の追跡が可能なことである.簡単な例として，二次元曲線 の場合をとりあげる.Fig.5.2(c)， (d)に示すように，単一の滑らかな初期閉曲線γ(0)が縮小し ていく過程でちぎれて2個の閉曲線に分離した場合を考える.このような動きをする曲線の 追跡に通常の媒介変数表示を用いたのでは，分裂する際に生じる特異点のため，分裂後の曲 線γの追跡が不可能となる.一方，曲線γに対する補助関数z=ゆ(x
，
仏t)は， Fig.5.2( a)， (b) に示すように，その問も位相が不変な単一の関数のままであるから， γをゆ(x，
y，
t)のゼロ等 高面とみなすことにより，その追跡が可能となる.このことは三次元曲面を追跡する場合も 全く同じことがいえる.5 . 3 数値解法
本節では，まず，動く曲面を追跡するための等高面方程式の数値計算法を説明する.次に，
等高面法を適用した動的輪郭線抽出手法を紹介する.最後に，これを三次元モデリングに適 用しようとする際に生じる問題点について述べる.
5.3.1 等 高 面 方 程 式 の 離 散 化
ゆ
( x ( t ) ， t )
は，曲面の位相が変化しても，単一の関数のままであるから， xの定義域を離散化 し，空間と時間に関する微分を有限微分近似に置き換えることが可能である.例えば，三次元空 間全体を辺の長さ九の立法格子であるセノレC(iふた)の集合に分割し，何j，kを解ゆ(ih，
jh，
kh，
nムt )
の近似値とすると，等高面方程式(5.5)から差分方程式 ゆ
f ; ; ‑ P
t，j，~
+
(F)Iマ
i，j，k何
j，kl=
0ムt (5.6)
が得られる.ここで，ムtは時間きざみ幅，nは時間ステップ数，マi，j，kは空間微分を適当な差 分演算子で置き換えたものである.量子化誤差や雑音の影響を受けにくくするため，本章で は差分演算子マi，j，kに対し前向き差分スキーム
[ 6 0 ]
を採用している.また，何j，kのゼロ等高面
η r
を求めるための具体的な手順を以下に示す.5.3. 数値解法 74
(1)不等式
max(
佼
j，k'・・1ゆじ1，j，k+l，ゆじ1，j+l，k+l)くO (5.7) を満たすセルC ( i
，j，k )
からなる連結領域Lm
の集合{ L l
ぃ・・，LM}
を抽出する.ただし，Lm nLm
，=
日(m"#m')とする.(2)各領域
Lm
を閉曲面r
nの内部領域とみなして，Lm
に外部から隣接するセルを成分と する連結領域r
弘を求める.(3) 和集合 U~=lr弘をゼロ等高面?とする.
ここでは，全空間領域を構成するセルの三次元配列を曲面に属しているか否かで分類する ことにより三次元曲面を表現していることに注意されたい.このような占有表現法は，解像 度を上げると多くの記憶容量が必要になるが，切断面の生成や符号化などの計算が容易にな
るという特徴がある.
5.3.2 輪 郭 線 抽 出 の 場 合 の 運 動 方 程 式
等高面法を用いた輪郭線抽出手法において， Malladiら[35]は，速度
F = ‑1 ‑eκ
( 5 . 8 )
で等速成長する滑らかな曲線を利用しようと考えた.ここで， κは曲線ゆ
( x ，
y， t )
= 0の外向 き法ベクトル方向の平均曲率を表す.ここで，第一項は，曲線を一定速度 ‑1で一様に縮小さ せようとする移流項である.第二項一日は，曲率の大きなところを平滑化する拡散項であり，薄板モデルにおける内部変形エネルギー項などと同じ正則化の効果を与える.すなわち，式
( 5 . 8 )
は等速成長方程式と平均曲率流方程式の二つの効果をあわせもつ曲線の運動方程式[ 6 1]
といえる.また， eはこれら二項の相対的重みを制御する正定数である.
しかし，このままでは曲線の移動が止まらず，一点に絡んでしまう.したがって，速度関 数
( 5 . 8 )
に対し，曲線の移動を輪郭付近で停止させる速度関数を画像データからどのようにし て定義するかが問題となる.これに対し， Malladiらはフィルタk I ( Z J ) = 1 (5.9)
を式
( 5 . 8 )
に掛けて，速度をF=k](‑l‑e
κ) (5.10)5.3.数値解法 75
と定義した.式(5.9)における Gσ*1は，標準偏差σのガウス関数による画像Iのたたみ込み積 分を表し，画像をガウス関数でぽかす操作に相当する.そのため， Gσ*1を微分したマGσ*1は 光強度が急激に変化するところ以外ではほとんど0となる.これより，フィルタ
k I ( X ，
y)は境 界から離れたところでは1となり，光強度が鋭く変化するところ，すなわち停止領域ではOとなるため，式(5.10)から曲線が対象輪郭の付近で移動を停止することになる.
ここで，上の速度関数(5.10)は，ゼロ等高線，すなわちゅ(X
，
y，
t) = 0上の点でしか意味 をもたないことに注意されたい.他方，等高面方程式は，全体領域を定義域とする関数ゅに 対して記述されなければならない.この問題に対処するため， Malladiら[35]は，ゼロ以外の 等高面ゆ(X，
y，
t)=
C， Ci ‑
0上の点の速度Fに，その点から最も近い曲線ゆ(x，
y，
t)=
0上の 点のFの値を設定した (Fig.5.3).この方法には，等高面同士の衝突を回避できるという利 点がある.この輸事s線抽出手法を実濃淡画像に適用した結果をFig.5.4に示す.対象物は二台 の車である.まず，対象物体群を包含するように初期輪郭を設定する (Fig.5.4(a)).この閉 曲線は，速度(5.10)に従って縮小する (Fig.5.4(b)，(c)).その後，物体問で分裂を生じ，最 終的に二つの対象物の輪郭を抽出している (Fig.5.4(d)) . Fig.5.4(b)とFig.5.4(c)に，この 閉曲線が次第に縮小し，対象輪郭付近に到達した点が移動を停止する様子を示す.Fig..5.4(d) に，この閉曲線が二つの閉曲線に分離し，対象輪郭付近で移動を停止した最終結果を示す.5.3.3 三 次 元 曲 面 再 構 成 の 場 合 の 問 題 点
本節では，上述の二次元形状モデリング手法の三次元物体への拡張を試みる.そのため，
式(5.10)の類推により，開曲面ゆ(x
，
y，
z，
t) = 0の外向き移動速度Fを F = k R ( ‑1 ‑eK)と定義する.ここで，Kは三次元平均曲率を表し，
K= (
必+ゅ;+必)一3ベ(仰+ゆzz)必+(伽+ゆzz) ゅ~ +
ゆ(xx+
ゆ川必‑2 [ 州 内 + ゆzゅ ん + ゆ 仰 の
z ] }
のように与えられる.また， eは式(5.10)と同じ正定数である.
(5.11)
(5.12)
ここで，曲面の移動が物体曲面のところで止まるようなフィルタkRを定義しなければな らない.二次元の輪郭線抽出の場合，隙間のない停止領域の設定が可能なため，式(5.10)に 従って曲線上のどの点も境界に到達した時点で移動を停止することができた.しかし，三次 元の場合，計測データは隊散的な三次元点からなる集合であり，使用センサの精度や走査性 能の限界のため，データのない領域が常に存在する.したがって，ここで隙間のない停止領
5.4. 最適化に基づく曲面変形 76
域を設定することは，二次元の場合よりもいっそう困難である.例えば，データ点が存在す る位置でkRがゼロになるように定めても， Fig.5.5(a)に示すように，曲面が欠損部を通過し たりデータ点をとび越えるなどして実曲面付近で移動が止まらない.したがって，移動する 曲面を停止させるため，さらに洗練された方法を開発しなければならない.
5 . 4 最適化に基づく曲面変形
本節では，三次元データに存在する隙間の問題を克服するため，新たに運動方程式を定式 化する.
式(5.11)からフィルタkRを取り除き，最初の等速移動の項‑1をあるデータ近接項に置 き換える.いま， 曲面点
X=(X ，
y， Z ) T
から三次元データ集合までの自乗距隣をd(x)
=
(x ‑Xc)2+
(y ‑Yc)2+
(z ‑zc)2 (5.13) と定義する.ただし， Xc=
(xc，
Yc，
zc)T はxから最も近いデータ点である.ここで，モデル 表面点xに対し， d(x)の値が小さくなる方向，すなわちベクトル‑¥l d = (‑dx
，
‑dy，
‑dz)T= (‑2(x ‑xc)
，
‑2(y ‑yふ
‑2(z‑zc))Tの方向に引力を作用させる.したがって，移動速度Fは F = 9 ‑eK によって与えられ， 9は
g =一α
マ
d.n(5.14)
(5.15)
(5.16) となる.ここで， eとKは式(5.11)と同じ記号， αは正定数である.また， nは単位法ベク
トノレ
(5.17)
であるから，式(5.16)は
9 =
‑ s {
(x ‑xc)CTx+
(y ‑Yc)cty+
(z ‑zc)CTz } (5.18)と表せる.ここで，
s
は式(5.16)のαを￨マゆ￨で割った正定数である.式(5.15)はデータ近 接項gと平滑化項‑eKの両方を含んでいるため，正則化アプローチとみなすことができる.この運動方程式により，モデ、ノレ曲面は実曲面付近に到達すると，ゆるやかに自動停止するよ うになり，停止領域を陽に設定することも不要となる.(Fig.5.6) .