JP4910478B2 - モデリング装置、モデリング方法およびプログラム - Google Patents

Description

本発明は、計測データに基づき対象物のモデルを生成するモデリング技術に関する。

近年、医療分野において、Ｘ線ＣＴ(Computed Tomography)装置或いはＭＲＩ（Magnetic Resonance Imaging)装置等によって撮影された生体内の３次元画像が広く用いられている。これにより、体内の臓器等の情報が視覚的に把握可能となり診断精度の向上が期待できる。しかし、その一方で、診断に用いられる３次元画像は数十枚〜数百枚のスライス画像から構成されているため、このような膨大な情報量の中から診断に必要な情報のみを得ることは読影医師にとって大きな負担となっている。

そこで、計算機を援用した定量的或いは自動的診断の要望が強まり、計算機による診断支援（ＣＡＤ：Computer−aided diagnosis）システムの研究が盛んに行われている。

計算機による診断支援を行うには、診断に必要な情報、つまり臓器領域等を正確に抽出することが重要な課題となる。

臓器の領域抽出手法としては、予め用意した標準モデルをエネルギー最小化原理に基づいて変形させ、目的の輪郭を見つけるモデルフィッティング（Model Fitting）手法が提案されている（非特許文献１）。

ツァガーン・バイガルマ、清水昭伸、小畑秀文、「３次元可変形状モデルによる腹部ＣＴ像からの腎臓領域抽出法の開発」、電子情報通信学会論文誌、D-II、2002年１月、Vol.J85-D-II、No.１、pp.140-148

しかしながら、抽出対象とする領域（抽出対象領域）と抽出対象以外の他の領域とが、互いに接している場合には、両者の接触部分において境界が区別できないため、モデルフィッティング手法を用いても両者の接触部分を正確にモデリングできないという問題がある。

そこで、本発明は上記課題に鑑みてなされたものであり、抽出対象領域と他の領域とが接触する場合においても、当該接触部分における抽出対象領域のモデリングを精度良く行うことが可能な技術を提供することを課題とする。

上記の課題を解決するため、請求項１の発明は、３次元計測データに基づいて対象物に関するモデリングを行うモデリング装置であって、前記対象物に関する標準モデルを取得する手段と、前記３次元計測データから、前記対象物に対応する領域である対象領域と、当該対象領域以外の他の領域とを含む暫定領域を抽出する暫定領域抽出手段と、モデルフィッティング手法により前記標準モデルを変形することによって、前記対象物に関するモデルを作成するモデル作成手段とを備え、前記モデル作成手段は、前記暫定領域に対して縮退処理を施す縮退手段と、前記暫定領域が前記縮退処理により分断される際の分断点を特定する分断点特定手段と、前記分断点を含み、前記対象領域と前記他の領域とのつながり方向に基づき規定される仮想面を設定する設定手段と、前記仮想面と前記暫定領域とが交差する画素を特定することによって、前記暫定領域において前記対象領域と前記他の領域とが接触する部分の仮想的な接触点を取得する接触点取得手段と、前記標準モデル表面上の複数の代表点の中から、前記接触点に対応する対応代表点を決定する対応点決定手段と、前記接触点と前記対応代表点との距離に関する要素を含む評価関数を最適化するように前記標準モデルを変形させる変形手段とを有することを特徴とする。

また、請求項２の発明は、請求項１の発明に係るモデリング装置において、前記モデル作成手段は、前記分断点近傍の周辺情報を取得する手段をさらに有し、前記設定手段は、前記分断点と前記周辺情報とに基づいて、前記分断点における前記対象領域と前記他の領域とのつながり方向を表すベクトルを算出し、前記ベクトルを法線ベクトルとし且つ前記分断点を含む仮想平面を前記仮想面として設定することを特徴とする。

また、請求項３の発明は、請求項１または請求項２の発明に係るモデリング装置において、前記暫定領域抽出手段は、領域拡張手法により前記暫定領域を抽出することを特徴とする。

また、請求項４の発明は、計測データに基づいて対象物に関するモデリングを行うモデリング装置であって、前記対象物に関する標準モデルを取得する手段と、前記計測データから、前記対象物に対応する領域である対象領域を含む暫定領域を抽出する暫定領域抽出手段と、モデルフィッティング手法により前記標準モデルを変形することによって、前記対象物に関するモデルを作成するモデル作成手段とを備え、前記モデル作成手段は、前記暫定領域に対して縮退処理を施す縮退手段と、前記暫定領域が前記縮退処理により分断される際の分断点を特定する分断点特定手段と、前記標準モデル表面上の複数の代表点の中から、前記分断点に対応する対応代表点を決定する対応点決定手段と、前記分断点と前記対応代表点との距離に関する要素を含む評価関数を最適化するように前記標準モデルを変形させる変形手段とを有することを特徴とする。

また、請求項５の発明は、３次元計測データに基づいて対象物に関するモデリングを行うモデリング方法であって、前記対象物に関する標準モデルを取得する工程と、前記３次元計測データから、前記対象物に対応する領域である対象領域と、当該対象領域以外の他の領域とを含む暫定領域を抽出する暫定領域抽出工程と、モデルフィッティング手法により前記標準モデルを変形することによって、前記対象物に関するモデルを作成するモデル作成工程とを備え、前記モデル作成工程は、前記暫定領域に対して縮退処理を施す縮退処理工程と、前記暫定領域が前記縮退処理により分断される際の分断点を特定する分断点特定工程と、前記分断点を含み、前記対象領域と前記他の領域とのつながり方向に基づき規定される仮想面を設定する工程と、前記仮想面と前記暫定領域とが交差する画素を特定することによって、前記暫定領域において前記対象領域と前記他の領域とが接触する部分の仮想的な接触点を取得する接触点取得工程と、前記標準モデル表面上の複数の代表点の中から、前記接触点に対応する対応代表点を決定する対応点決定工程と、前記接触点と前記対応代表点との距離に関する要素を含む評価関数を最適化するように前記標準モデルを変形させる変形処理工程とを有することを特徴とする。

また、請求項６の発明は、コンピュータに、対象物に関する標準モデルを取得する工程と、３次元計測データから、前記対象物に対応する領域である対象領域と、当該対象領域以外の他の領域とを含む暫定領域を抽出する暫定領域抽出工程と、モデルフィッティング手法により前記標準モデルを変形することによって、前記対象物に関するモデルを作成するモデル作成工程とを実行させるためのプログラムであって、前記モデル作成工程は、前記暫定領域に対して縮退処理を施す縮退処理工程と、前記暫定領域が前記縮退処理により分断される際の分断点を特定する分断点特定工程と、前記分断点を含み、前記対象領域と前記他の領域とのつながり方向に基づき規定される仮想面を設定する工程と、前記仮想面と前記暫定領域とが交差する画素を特定することによって、前記暫定領域において前記対象領域と前記他の領域とが接触する部分の仮想的な接触点を取得する接触点取得工程と、前記標準モデル表面上の複数の代表点の中から、前記接触点に対応する対応代表点を決定する対応点決定工程と、前記接触点と前記対応代表点との距離に関する要素を含む評価関数を最適化するように前記標準モデルを変形させる変形処理工程とを有することを特徴とする。

請求項１から請求項６に記載の発明によれば、暫定的に抽出された暫定領域において対象領域と対象領域以外の他の領域とが接触する部分の仮想的な接触点を取得するとともに、標準モデル表面上の代表点の中から、接触点に対応する対応代表点を決定し、接触点と対応代表点との距離に関する要素を含む評価関数を最適化するように標準モデルを変形するので、対象領域と対象領域以外の他の領域とが接触する部分におけるモデリングを精度良く行うことが可能になる。

以下、本発明の実施形態を図面に基づいて説明する。

＜１．計測データ＞
この実施形態では、Ｘ線ＣＴ装置によって取得される３次元画像（ボリュームデータないし３次元計測データとも称する）を計測データとして取得し、当該計測データから所望の領域（「抽出対象領域」とも称する）を抽出する処理について説明する。なお、本発明は、これに限定されず、他の計測データ（例えば、ＭＲＩ装置によって取得される３次元画像等）に対して適用することも可能である。

図１は、Ｘ線ＣＴ装置によって取得される３次元画像ＴＤＰ１を示す図である。

図１に示されるように、３次元画像ＴＤＰ１は、物体（人体）を輪切りにした断面を示す複数（例えば、数十枚から数百枚）のスライス画像（断層画像）ＳＤで構成されている。各スライス画像ＳＤは、各点（各画素）でのＸ線の吸収量（ＣＴ値）を濃淡表示して可視化した画像である。図１においては、３次元画像ＴＤＰ１を構成する複数のスライス画像ＳＤのうちの所定位置におけるスライス画像ＳＣも併せて示されている。

この実施形態に係る３次元画像ＴＤＰ１は、血管に造影剤を注入した状態で撮影されており、各スライス画像ＳＤにおいては、血管に対応する領域（血管領域とも称する）と背骨に対応する領域（背骨領域とも称する）とがその他の臓器に比べて高いＣＴ値を有している。

図２は、スライス画像ＳＣを示す図である。図２に示されるように、スライス画像ＳＣにおいては、血管領域Ｖｒ（詳細には、Ｖｒ１，Ｖｒ２，...）と背骨領域Ｂｒ（詳細には、Ｂｒ１，Ｂｒ２，Ｂｒ３，Ｂｒ４）とが存在する。

スライス画像ＳＣにおいては、背骨領域Ｂｒ１と血管領域Ｖｒ１とが部分的に接している。そのため、通常の手法では、背骨領域Ｂｒと血管領域Ｖｒとをその他の領域から抽出することは可能であるが、背骨領域Ｂｒと血管領域Ｖｒとを自動的に区別することは難しい。例えば、血管領域Ｖｒを抽出対象領域として設定し、領域拡張法で血管領域Ｖｒを抽出することを試みると、血管領域Ｖｒ（Ｖｒ１，Ｖｒ２，...）だけでなく背骨領域Ｂｒ１をも含む領域（混在領域とも称する）が抽出されてしまうことになる。また、血管をモデル化することが困難なため、血管の標準モデルを用いたモデルフィッティング手法によって血管を抽出することも困難である。

そこで、この実施形態では、モデル化が困難な対象物（血管）とモデル化が可能な対象物（背骨）とが一部において接する状態を有する計測データから、後述のような手法を用いて、モデル化が困難な対象物（血管）を、モデル化が可能な対象物（背骨）と分離して抽出する場合について例示する。これによれば、背骨付近の血管を正確に抽出することが可能である。

また特に、この実施形態においては、血管領域を最終的に抽出するために、背骨領域の抽出処理が「背骨」のモデリング処理を伴って行われる（後述）。この「背骨」のモデリング処理は、本発明の思想に基づいて実行されるものである。換言すれば、この実施形態においては、本発明の思想が「背骨」のモデリング処理に適用される場合について例示されている。

＜２．装置構成＞
図３は、本発明の実施形態に係る領域抽出装置１の概要を示す図である。領域抽出装置１は、計測データ（詳細には３次元計測データ（立体計測データとも称する））から所望の対象物の領域（詳細には立体領域）を抽出する装置である。この実施形態では「血管」が立体的な抽出対象領域として抽出される。なお、この領域抽出装置１は、最終的な血管領域抽出処理に先立って行われる背骨領域抽出処理に際して、「背骨」のモデルを作成するモデリング処理をも行うことから、モデリング装置とも称せられる。

図３に示すように領域抽出装置１は、パーソナルコンピュータ（以下、単に「パソコン」と称する）２と、モニター３と、装着部５と、操作部４とを備えている。

パソコン２は、制御部２０、入出力Ｉ／Ｆ２１、及び記憶部２２を備えている。

入出力Ｉ／Ｆ２１は、モニター３、操作部４および装着部５とパソコン２との相互間でデータの送受信を行うためのインターフェイス（Ｉ／Ｆ）であり、制御部２０との間でデータの送受信を行う。

記憶部２２は、例えばハードディスクなどで構成されており、後述の各種処理を実行するためのソフトウェアプログラム（以下、単に「プログラム」と称する）ＰＧ等を格納している。

制御部２０は、主にＣＰＵ、ＲＯＭ２０ａ及びＲＡＭ２０ｂ等を有し、パソコン２の各部を統括制御する部位である。

モニター３は、例えば、ＣＲＴで構成され、制御部２０で生成される表示用画像を可視的に出力する。

操作部４は、キーボード及びマウス等から構成され、使用者（ユーザ）に各種操作にしたがって各種信号を入出力Ｉ／Ｆ２１に送信する。

また、装着部５は、メモリカード５１等の記憶媒体を着脱自在に装着することができる。そして、装着部５に装着されたメモリカード５１に格納される各種データ又はプログラム等を入出力Ｉ／Ｆ２１を介して制御部２０或いは記憶部２２に取り込むことができる。

図４は、領域抽出装置１の各種機能を示すブロック図である。

図４に示されるように、領域抽出装置１は、計測データ入力部１１と開始点決定部１２と暫定領域抽出部１３と初期位置決定部１４とモデルフィッティング部１５と境界設定部１６と領域抽出部１７と抽出領域出力部１８とモデル格納部１９とを備えている。

これらの各種機能は、制御部２０内のＣＰＵ等の各種ハードウェアを用いて所定のプログラムＰＧを実行することによって実現される。

計測データ入力部１１は、ＣＴ装置あるいはＭＲＩ装置などによって取得された３次元計測データ（ボリュームデータ）を入力する機能を有している。

開始点決定部１２は、領域拡張法における開始点を決定する機能を有している。開始点決定部１２は、第１の領域拡張法（後述）における開始点と第２の領域拡張法（後述）における開始点との双方を決定する。なお、この実施形態においては第１および第２の領域拡張法において、同一の開始点が用いられる。

暫定領域抽出部１３は、領域拡張法（第１の領域拡張法とも称する）を用いて、抽出対象物体（領域）を暫定的に抽出する機能を有している。暫定領域抽出部１３によって抽出された暫定領域は、上述の血管と背骨との双方を含む混在領域として求められる。

モデル格納部１９は、モデルフィッティング処理に用いられる「標準モデル」（図１４，図１７等参照）ＳＢを格納している。領域抽出装置１は、モデル格納部１９から読み出すことによって「標準モデル」を取得することができる。ここで、「標準モデル」は、モデルフィッティング処理（後述）用に予め準備されたモデルであり、或る物体についての一般的（標準的）な３次元形状データ等を有するものである。後述するように、この実施形態では、本来の抽出対象物体以外の物体、具体的には、「背骨」の標準モデルが用いられる。換言すれば、抽出対象物体の領域（抽出対象領域）として抽出されるべき領域から除外すべき領域（「除外領域」ないし「非対象領域」とも称する）に関する標準モデルが用いられる。

初期位置決定部１４は、暫定領域抽出部１３によって暫定的に抽出された領域を用いて、標準モデルの初期位置を決定する機能を有している。

モデルフィッティング部１５は、モデルフィッティング手法を用いて所定の評価関数（後述する総合エネルギーＵｅ）を最適化するように標準モデルを計測データ内において変形させ、「除外領域」のモデルを生成する機能を有している。換言すれば、モデルフィッティング部１５は、「除外領域」に関する標準モデルを用いて、計測データから「除外領域」をモデルフィッティング手法によって抽出する機能を有している。

詳細には、このモデルフィッティング部１５は、上述の暫定領域に対して縮退処理を施す機能、当該縮退処理によって暫定領域が分断される際の分断点（後述）を特定する機能、当該分断点近傍の周辺情報を取得する機能、暫定領域において対象領域と当該対象領域以外の他の領域とが接触する部分の仮想的な接触点（後述）を分断点および周辺情報に基づいて取得する機能、標準モデル表面上の複数の代表点の中から、接触点に対応する対応代表点（後述）を決定する機能、ならびに接触点と対応代表点との距離に関する要素を含む評価関数（総合エネルギーＵｅ）を最適化するように標準モデルを変形させる機能等を有している。

境界設定部１６は、計測データにおいて、モデルフィッティング後の「除外領域（背骨）」のモデルが存在する領域の境界を、領域抽出部１７による領域拡張法（第２の領域拡張法とも称する）での「拡張限界」として設定する。

領域抽出部１７は、除外領域の境界を拡張限界とする領域拡張法（第２の領域拡張法）によって、計測データから抽出対象領域（血管）を抽出する機能を有している。

抽出領域出力部１８は、抽出した対象物体（血管）をモニター３に表示出力する機能を有している。また、抽出領域出力部１８は、抽出した対象物体（血管）に関するデータをファイル出力する機能をも有している。

＜３．動作＞
次に、領域抽出装置１における動作の概要について図５〜図９等を参照しながら説明する。

図５は、領域抽出装置１の全体動作を示すフローチャートである。また、図６〜図９は、それぞれ、３次元画像ＴＤＰ１を概念的に示す図である。図６〜図９においては、３次元画像ＴＤＰ１内に含まれる各部位、具体的には、背骨領域Ｂｒ１を含む背骨ＢＮと、血管領域Ｖｒ１を含む血管ＶＬとが模式的に示されている。また、図示の簡略化のため、背骨領域Ｂｒ２，Ｂｒ３，...および血管領域Ｖｒ２，Ｖｒ３，...は適宜省略されている。

図５に示すように、まず、ステップＳ１においては、処理対象の計測データ（３次元画像ＴＤＰ１）が領域抽出装置１に入力される。図６には、この初期状態における未処理データが概念的に示されている。なお、この入力直後の時点においては、領域抽出装置１は、計測データを、スライス画像ごとに区分された画素（ボクセル）の集合体として認識しているに過ぎない。

次のステップＳ２，Ｓ３においては、第１の領域拡張法による領域抽出処理が行われる。具体的には、ステップＳ２において開始点ＳＰ１（図７参照）が設定され、ステップＳ３において３次元の領域拡張手法によって領域が抽出される。ステップＳ３においては、開始点ＳＰ１からその周辺へと向かって検索対象範囲を徐々に広げていき、所定の閾値ＴＨ１以上の階調値を有する画素（ボクセル）を順次抽出領域に追加して当該抽出領域を拡張していく処理が行われる。

図７は、ステップＳ３終了時点における抽出状態を示す図である。図７は、開始点ＳＰ１から拡張された領域（斜線）が抽出されている様子を示している。背骨領域Ｂｒと血管領域Ｖｒとが一部において接しているため、図７に示すように、開始点ＳＰ１が血管領域Ｖｒ１内に設定されたとしても、ステップＳ３の領域拡張法によって、血管領域Ｖｒ（Ｖｒ１，Ｖｒ２，...）だけでなく背骨領域Ｂｒ１も抽出される。したがって、ステップＳ３終了時点では、背骨領域Ｂｒと血管領域Ｖｒとが混在した領域（混在領域）が抽出されることになる。

この実施形態では、血管ＶＬに隣接する背骨ＢＮのモデル化が可能であることに着目して、背骨領域Ｂｒ（除外領域）を除外して血管領域Ｖｒ（本来の抽出対象領域）のみを抽出する処理を行う。

具体的には、まず、ステップＳ３で抽出された混在領域と、背骨ＢＮに関する標準モデルＳＢとを用いたモデルフィッティング手法によって、計測データから背骨領域Ｂｒを抽出する（ステップＳ４，Ｓ５）。詳細には、モデルフィッティング手法における背骨ＢＮの標準モデルＳＢ（図８参照）の初期位置を決定し（ステップＳ４）、実際にモデルフィッティング手法を適用して背骨領域Ｂｒを抽出する（ステップＳ５）。ステップＳ５では、初期位置に配置された標準モデルがモデルフィッティング手法によって変形され、変形後の標準モデルの位置に存在する領域が、対象物の領域として計測データから抽出される。図８は、ステップＳ５終了時点における抽出状態を示す図である。図８は、背骨ＢＮの標準モデルＳＢを用いたモデルフィッティング手法によって、背骨領域Ｂｒが抽出された状態を示している。

次に、抽出された背骨領域Ｂｒの境界を第２の領域拡張法における「拡張限界」として設定し（ステップＳ６）、「拡張限界」の制約の下で第２の領域拡張法によって血管領域Ｖｒが抽出される（ステップＳ７）。ステップＳ７においては、開始点ＳＰ２（ここでは開始点ＳＰ１と同一点）からその周辺へと向かって検索対象範囲を徐々に広げていき、所定の閾値ＴＨ２以上の階調値を有する画素（ボクセル）を順次抽出領域に追加して当該抽出領域を拡張していく処理が行われる。ステップＳ７の領域拡張処理においては、拡張限界が設定されているため、背骨領域Ｂｒを含まない状態で拡張領域が生成されることになる。

図９は、ステップＳ７終了時点における抽出状態を示す図である。図９においては、血管領域Ｖｒ（Ｖｒ１）が、背骨領域Ｂｒから分離された状態で抽出されている様子が示されている。なお、図９においては図示を省略しているが他の血管領域Ｖｒ２等をも含む血管領域Ｖｒが抽出される。

その後、ステップＳ８では、抽出された血管領域Ｖｒが抽出対象物体としてモニター３等に出力される。また、抽出された血管領域Ｖｒのデータがファイル出力される。

上記のようにして、モデル化が困難な血管領域Ｖｒを、背骨領域Ｂｒから分離して、計測データから正確に抽出することができる。

以下、ステップＳ２〜Ｓ７の各処理について順次に詳述する。

＜３−１．開始点決定処理（ステップＳ２）＞
上述したように、ステップＳ２においては、領域拡張法における開始点ＳＰ１および開始点ＳＰ２（総称する場合は開始点ＳＰとも称する）が決定される。開始点ＳＰ１（および開始点ＳＰ２）は、操作者が、モニター３に表示された３次元画像ＴＤＰ１を見ながら、マウス等を操作して手動で決定するようにしてもよいが、領域抽出装置１が自動的に決定するようにしてもよい。ここでは自動決定処理の一例について説明する。

混在領域抽出用の開始点ＳＰ１は、血管領域Ｖｒ内の点であってもよく背骨領域Ｂｒ内の点であってもよい。また、血管領域抽出用の開始点ＳＰ２は血管領域Ｖｒ内の点であることを要する。ここでは、簡単化のため、血管領域Ｖｒ内の１つの点を、開始点ＳＰ２として設定するとともに開始点ＳＰ１としても設定するものとする。すなわち、血管領域Ｖｒ内の１つの点を、両開始点ＳＰ１，ＳＰ２として兼用するものとする。

また、血管の連続性を考慮すれば、開始点ＳＰとして１つの点を指定すれば十分である。ここでは、下記のような厳しい条件を満たすような１つの点を選択することによって、血管領域Ｖｒ１内の点であることが確実な点を開始点として得るものとする。

具体的には、まず、１つのスライス画像内において、所定の閾値ＴＨ１以上のＣＴ値を有する領域を、連続した部分ごとに区分した部分領域を抽出し、当該部分領域が次の条件ＣＤ１〜ＣＤ５を満たすか否かを判定する。

・条件ＣＤ１：部分領域のサイズが閾値以上；
・条件ＣＤ２：部分領域内の最大ＣＴ値が閾値を越えない；
・条件ＣＤ３：部分領域内のＣＴ値の標準偏差が閾値を越えない；
・条件ＣＤ４：部分領域に外接する矩形領域（Bounding Box）に占める当該部分領域の面積割合が閾値以上；
・条件ＣＤ５：部分領域のx座標, y座標の相関値が閾値を越えない。

条件ＣＤ１は、末梢血管を排除し、大動脈内の１つの点を開始点として選択することを意図するものである。

また、図１０に示すように、血管領域の最大ＣＴ値（各スライス画像における血管領域内に含まれる全画素のＣＴ値のうちの最大値）は、（一般的に骨領域の最大ＣＴ値よりも小さく）所定の閾値（例えばＴＨ９）以下であることが多い。条件ＣＤ２はこのような特性を利用して、骨領域を排除して血管領域を抽出することを意図している。なお、図１０は、或る３次元画像について、全領域のＣＴ値の度数分布と、血管領域のＣＴ値の度数分布とを示す図である。

また、血管領域のＣＴ値の標準偏差は（一般的に骨領域のＣＴ値の標準偏差よりも小さく）所定の閾値以下であることが多い。条件ＣＤ３は、このような特性を利用して、骨領域を排除して血管領域を抽出することを意図している。

また、背骨に沿って走る大動脈の血管領域Ｖｒ１は、各スライス画像において円形に近い。円形度に関する条件ＣＤ４，ＣＤ５は、このような特性を利用して血管領域Ｖｒ１を選択することを意図するものである。

そして、部分領域がこれらの条件ＣＤ１〜ＣＤ５の全てを満たす場合には当該部分領域の中心座標を求め、当該中心座標を開始点として決定して開始点決定処理を終了する。

一方、部分領域がこれらの条件ＣＤ１〜ＣＤ５のいずれか１つでも満たさないときには、次の部分領域を検索して同様の判定動作を再び行う。また、１つのスライス画像に上記の条件ＣＤ１〜ＣＤ５の全てを満たす部分領域が１つも存在しない場合には、別のスライス画像について同様の処理を行う。

このような処理によって、上記の条件を満たす開始点ＳＰが決定される。なお、図６等に示すように、背骨領域Ｂｒと血管領域Ｖｒとは全てのスライス画像で接しているのではなく、両領域Ｂｒ，Ｖｒが互いに接していないスライス画像も存在する。上記の条件ＣＤ１〜ＣＤ５を考慮することによれば、両領域Ｂｒ，Ｖｒが互いに接していないスライス画像における血管領域Ｖｒ内の点が開始点として決定される。

＜３−２．第１の領域拡張処理（ステップＳ３）＞
上述のように、ステップＳ３の第１の領域拡張処理では、領域拡張法を採用して、血管領域Ｖｒと背骨領域Ｂｒとの双方を含む混在領域を暫定的に抽出する。

図１１は、領域拡張法を用いた暫定領域抽出処理（ステップＳ３）の詳細を示すフローチャートである。図１２は、スライス画像ＰＭ上の画素Ｍ１を中心として隣接する画素付近を拡大表示した図である。図１２においては、スライス画像ＰＭの画素Ｍ１の隣接画素として、同一階層のスライス画像ＰＭにおける隣接８画素Ｍ２〜Ｍ９に加えて、スライス画像ＰＵにおける隣接９画素Ｕ１〜Ｕ９と、スライス画像ＰＤにおける隣接９画素Ｄ１〜Ｄ９とが示されている。なお、スライス画像ＰＵは、スライス画像ＰＭに対して＋ｚ方向に隣接するスライス画像（換言すればスライス画像ＰＭの直上層のスライス画像）であり、スライス画像ＰＤは、スライス画像ＰＭに対して−ｚ方向に隣接するスライス画像（換言すればスライス画像ＰＭの直下層のスライス画像）である。

まず、ステップＳ３１において、ステップＳ２で決定された開始点ＳＰ１が拡張開始画素（注目画素とも称する）に設定される。

次に、ステップＳ３２において、注目画素に隣接する画素が拡張候補画素として検出される。詳細には、ステップＳ３１において注目画素として画素Ｍ１が選択されたとすると、画素Ｍ１に隣接する画素全てが拡張候補画素となる。つまり、スライス画像ＰＭ上の画素Ｍ２〜Ｍ９、スライス画像ＰＵ上の画素Ｕ１〜Ｕ９、及びスライス画像ＰＤ上の画素Ｄ１〜Ｄ９の計２６画素が拡張候補画素となる。

次に、ステップＳ３３において拡張候補画素がステップＳ３２で検出されたか否かが判定され、ステップＳ３２で拡張候補画素が検出されている場合には、ステップＳ３４へ移行する。

ステップＳ３４では、検出された拡張候補画素の階調値が所定の閾値ＴＨ１以上であるという条件を満たすか否かが判定される。拡張候補画素の階調値が所定の閾値ＴＨ１以上であると判定されると、ステップＳ３５に移行する。

ステップＳ３５では、当該拡張候補画素を、拡張後の領域（拡張領域）内の画素とする領域拡張を行う。

一方、ステップＳ３４において、拡張候補画素の階調値が所定の閾値ＴＨ１以上でないと判定されると、当該拡張候補画素の領域拡張は行われない。

次に、ステップＳ３６では、拡張候補画素の階調値判定工程（ステップＳ３４）をまだ経ていない拡張候補画素が存在するか否かが判断される。

階調値判定工程（ステップＳ３４）を経ていない拡張候補画素が存在する場合には、当該未終了の拡張候補画素に対してステップＳ３４〜Ｓ３６の処理が実行される。

一方、全ての拡張候補画素に対して階調値判定工程（ステップＳ３４）が終了している場合には、ステップＳ３２に戻る。そして、ステップＳ３２においては、ステップＳ３５において拡張された画素（換言すれば、領域内の画素であるとして特定された画素）が新たな注目画素として設定されるとともに、当該注目画素にさらに隣接する画素であって、領域内に存在するか否かが未だ判定されていない画素が、新たな拡張候補画素として検出される。その後、新たに検出された拡張候補画素に対して、上述のステップＳ３３〜Ｓ３６の工程が実行される。ステップＳ３２〜Ｓ３６の工程は、新たな拡張候補画素を検出することが可能な限り繰り返され、新たな拡張候補画素が検出できなくなると、暫定領域抽出処理は終了する（ステップＳ３３）。

このようにステップＳ３２〜Ｓ３６の工程が繰り返し実行されることによって開始点ＳＰから徐々に領域拡張が行われ、背骨領域Ｂｒと血管領域Ｖｒとが混在した領域（混在領域）が暫定的に抽出される。

図１３は、混在領域抽出後のスライス画像ＳＣを示す図である。図１３においては、抽出された領域に斜線が付されて示されている。なお、ステップＳ３終了時点における抽出状態は、図７においても（立体的に）示されている。

＜３−３．標準モデルの初期位置決定処理（ステップＳ４）＞
次に、モデルフィッティング手法における標準モデルの初期位置を決定する初期位置決定処理（ステップＳ４）について図１３および図１４を参照しながら説明する。

この実施形態においては、次のステップＳ５のモデルフィッティング処理用の背骨ＢＮの標準モデルＳＢとして、図１４に示すような比較的簡略化された標準モデルＳＢ１を採用する。この標準モデルＳＢ１は、背骨の各椎骨の概略形状を円柱とみなしてモデル化されたものである。標準モデルＳＢ１は、複数の椎骨部分を表現する各円柱ＣＬと、当該各円柱ＣＬ間を接続する曲面ＡＰとによって構成されている。また、当該円柱部分（背骨領域Ｂｒ１に相当）以外の部分、具体的には、当該円柱部分から外側に突出した突起部分（背骨領域Ｂｒ２に相当）はモデル化の対象から除外されている。

なお、血管領域Ｖｒと接するのは、椎骨の略円柱状部分に対応する背骨領域Ｂｒ１であることが多く、他の背骨領域Ｂｒ２〜Ｂｒ４と血管領域Ｖｒとが接することはほとんど無いという事情がある。このような事情のため、簡易な標準モデルＳＢ１を用いてモデルフィッティングを行って背骨領域Ｂｒ１を抽出し、当該背骨領域を排除するような第２の領域拡張処理（ステップＳ７）を行うことによっても、血管領域Ｖｒのみを抽出することが可能である。

さて、標準モデルＳＢ１は、背骨領域Ｂｒ１に対応する椎骨部分を中心にモデル化されたものであるため、背骨領域Ｂｒ１に近い位置に標準モデルＳＢ１を初期配置することが好ましい。また、背骨領域Ｂｒ２は突起部分を有しているため背骨領域Ｂｒ２の位置（例えば背骨領域Ｂｒ２の下端点ＰＴ２）を特定することは比較的容易である。そこで、当該背骨領域Ｂｒ２の位置を利用することに加えて、背骨領域Ｂｒ１と背骨領域Ｂｒ２との標準的な位置関係をもさらに利用して、標準モデルＳＢ１の初期位置を決定する。

図１３に示すように、領域Ｂｒ２は、混在領域のうち、スライス画像内において最も下方に存在している。したがって、混在領域のうち、各スライス画像内において最も下方に存在する点が領域Ｂｒ２の下端ＰＴ２の候補となる。

より詳細には、各スライス画像について、ステップＳ３で抽出された混在領域のＹ方向における最下端位置を検出する。そして、図１５に示すように、当該最下端位置が複数のスライス画像にわたって連続的に変化するとき、連続的変化区間Ｋｉの中央のスライス画像における最下端位置（×印）を、当該区間Ｋｉに対応する突出部分の基準点の位置として決定する。また、このような突出部分の基準点（例えばＰＴ２（図１３））と当該突出部分に対応する椎骨の基準点（例えばＰＴ１（図１３））との標準的な位置関係に基づいて、各椎骨の基準点（例えば位置ＰＴ１）の推定座標値が特定される。このような推定座標値は複数の椎骨のそれぞれについて求められる。その後、各椎骨の基準点の推定座標値に対して、標準モデルＳＢ１において各椎骨の基準点に対応する点の位置が一致するように、標準モデルＳＢ１を移動する。

このようにして標準モデルＳＢ１の初期配置が行われる。

＜３−４．モデルフィッティング処理（ステップＳ５）＞
ステップＳ５のモデルフィッティング処理は、初期配置された「標準モデル」を、計測データにおけるモデリング対象領域に関する情報（形状等）を用いて変形する処理である。なお、本出願においては、モデルフィッティング処理による変形後の標準モデル（換言すれば、計測データの情報が反映された標準モデル）を「個別モデル」とも称するものとする。

なお、モデリング対象物体（ここでは背骨）の標準モデルＳＢは、例えば、微小な多角形（例えば、三角形）（「ポリゴン」とも称する）で構成され、記憶部２２等に保存されている。ポリゴンで構成された標準モデルは、各ポリゴンの頂点の３次元座標によってそのモデルの表面形状を表現することができる。また、標準モデルにおけるポリゴンの頂点のうち、全頂点あるいは代表的な幾つかの頂点を「制御点」とも称する。そして、この制御点を代表点としてモデルフィッティング処理を行う。

本実施形態におけるモデルフィッティング処理では、上述のように簡易な標準モデルＳＢ１を用いて、混在領域ＫＲから背骨領域Ｂｒ１が抽出される。具体的には、背骨領域Ｂｒ１と血管領域Ｖｒ１との接触部分において接触点Ｔｋを仮想的に設定し、標準モデルＳＢ１表面上の制御点Ｃｊの中から接触点Ｔｋに対応する制御点（「対応制御点」或いは「対応代表点」とも称する）Ｃｋを決定する。そして、接触点Ｔｋと対応制御点Ｃｋとの距離が近づけば良化する要素を含んだ評価関数を最適化するように標準モデルを変形させ、背骨領域Ｂｒ１（抽出対象物体）に関する個別モデルを最終的に作成する。

以下では、このような背骨領域Ｂｒ１に対して行われるモデルフィッティング処理について説明する。なお、当該説明においては簡単化のため背骨領域Ｂｒ１と背骨領域Ｂｒ１に接する血管領域Ｖｒ１とに注目して説明する。

図１６および図１７は、モデルフィッティング処理（ステップＳ５）の詳細を示すフローチャートである。図１８は、暫定的に抽出された混在領域ＫＲに含まれる背骨領域Ｂｒ１と血管領域Ｖｒ１とを示す図であり、図１９は、背骨領域Ｂｒ１と血管領域Ｖｒ１とに縮退処理を施した図である。図２０および図２１は、縮退処理による分断直前の分断周辺情報を示す図である。図２０および図２１では、分断点を中心とする或る平面内の５×５の画素が示されている。図２２は、混在領域ＫＲに設定された仮想平面ＩＦを示す図である。図２３は、計測データ中に配置された標準モデルＳＢ１を示す図であり、図２４は、図２３の領域ＲＢ１付近の拡大図である。図２５は、計測データ中に配置された標準モデルＳＢ１を示す図であり、図２６は、図２５の領域ＲＢ２付近の拡大図であり、モデルフィッティングにおける制御点の移動の様子を示している。また、図２７は、制御点間を仮想バネで繋いだ様子を示す概念図である。なお、暫定的に抽出された混在領域ＫＲを含む抽出（計測）データは２値化されており、暫定的に抽出された混在領域ＫＲの各画素の画素値は「１」として表され、当該領域ＫＲに含まれない各画素の画素値は「０」として表される。

まず、ステップＳ５１（図１６）では、ステップＳ３において暫定的に抽出された混在領域ＫＲ（図１８参照）に対して、縮退処理が施される。縮退処理は、暫定的に抽出された混在領域の輪郭を示す境界画素を１画素ずつ削除する処理であり、図１９に示されるように混在領域が背骨領域Ｂｒ１と血管領域Ｖｒ１とに分断されるまで実行される。

混在領域が分断されたか否かは、例えば、混在領域を示す画素のうち、削除された画素に隣接する画素を拡張開始画素として領域拡張処理を実行して混在領域を再度抽出し、削除前後の混在領域の画素数を比較することにより判断することができる。具体的には、削除後の混在領域の画素数が削除前の混在領域の画素数より（約）１画素少ない場合には、混在領域は未だ分断されていないと判断し、削除前後の画素数が相当数減少している場合には、混在領域は分断されたと判断する。

また、縮退処理が実行されると、縮退処理実行後の混在領域を表すデータ（混在領域データ）が保存される。

なお、本実施形態では、注目画素に隣接する画素全てが当該注目画素に連続し得る画素である。例えば、図１２に示される画素Ｍ１を注目画素とした場合には、注目画素Ｍ１に隣接する計２６画素が、注目画素と連続し得る画素となる。そして、隣接画素と注目画素とがともに混在領域の画素である場合（具体的には、二値化された画像データにおいて、注目画素および隣接画素が同じ画素値「１」を有する場合）には、２つの画素がつながっている（連続する）と判断される。

ステップＳ５２において、混在領域が分断されたと判断されると、ステップＳ５３に移行する。

ステップＳ５３では、計測データにおける分断された箇所（「分断箇所」、「分断点」とも称する）ＢＰの３次元座標が特定される。また、分断される直前の混在領域データから分断箇所周辺の情報（「分断周辺情報」とも称する）が保存される。分断周辺情報としては、例えば、分断箇所ＢＰを中心とした５×５×５の画素情報を保存すればよい。

次のステップＳ５４では、分断周辺情報に基づいて、背骨領域Ｂｒ１と血管領域Ｖｒ１との接触部分における接触点が取得される。

ここで、接触点の取得手法について詳述する。背骨領域Ｂｒ１と血管領域Ｖｒ１との接触点は、分断点（分断画素）ＢＰを含む仮想的な平面（仮想平面）ＩＦとステップＳ３において暫定的に抽出された混在領域とが交差（重複）する画素を算出することにより求められる。換言すれば、背骨領域Ｂｒ１と血管領域Ｖｒ１との接触点は、混在領域内において仮想平面ＩＦ上に存在する画素を取得することにより求められる。

例えば、ステップＳ５３において図２０に示されるような分断箇所（分断画素）ＢＰを中心とする分断周辺情報が取得されたとすると、まず、分断箇所ＢＰにおける背骨領域Ｂｒ１と血管領域Ｖｒ１とのつながり方向（「接続方向」とも称する）を表すベクトルＶＴが算出される。当該ベクトルＶＴは、分断箇所ＢＰから混在領域を表す各画素に向かう各単位ベクトルの成分に対して、主成分分析 (Principal Component Analysis, ＰＣＡ)を行うことによって算出される。そして、当該ベクトルＶＴを法線ベクトルとし、かつ、分断箇所ＢＰを含む仮想平面ＩＦが設定される（図２１参照）。次に、このように設定された仮想平面ＩＦと混在領域とが交差する画素が接触点Ｔｋ（図２２中×印）として取得（算出）される。なお、ここでは分断点ＢＰは接触点Ｔｋの１つでもある。

次のステップＳ５５（図１６）では、初期配置された標準モデルＳＢ１の各制御点Ｃｊの中から、接触点Ｔｋに対応する制御点（対応制御点）Ｃｋを決定する処理が行われる。具体的には、図２３および図２４に示すように、標準モデルＳＢ１の初期配置後のデータ空間（３次元空間）において、接触点Ｔｋの最も近傍に存在する制御点Ｃｊが、接触点Ｔｋの対応点として求められる。例えば、接触点Ｔ１０の対応制御点Ｃｋとして、接触点Ｔ１０に最も近い制御点Ｃ２０を、接触点Ｔ１０の対応点として求めることができる（図２４参照）。

なお、このような接触点Ｔｋに対応する対応点の決定は、接触点Ｔｋごとに行われ、異なる接触点Ｔｋ同士が共通の制御点Ｃｋを対応点とすることもある。例えば、接触点Ｔ１５の対応点は、接触点Ｔ１０の対応制御点と同じ制御点Ｃ２０となり、接触点Ｔ２０の対応点は制御点Ｃ２１となる（図２４）。

次のステップＳ５６（図１７）では、暫定的に抽出された混在領域ＫＲの輪郭（境界）を示す境界画素（境界点）のうち、標準モデルＳＢ１の各制御点Ｃｊに対応する点（以下、「対応点」とも称する）Ｑｊを決定する処理が行われる（図２５参照）。すなわち、標準モデルＳＢ１における各制御点（代表点）Ｃｊと計測データにおける対応点Ｑｊとの対応関係が決定される。ここでは、混在領域ＫＲの輪郭（境界）を示す画素（境界点）のうち、制御点Ｃｊの最も近傍に存在する画素を当該制御点Ｃｊの対応点とする手法を用いる。例えば、図２６においては、制御点Ｃ３１の対応点は、制御点Ｃ３１の最も近傍に存在する境界点（画素）Ｑ３１となる。

次に、ステップＳ５７では、標準モデルＳＢ１の制御点Ｃｊのうち任意の一点（以下、「移動対象点」とも称する）（例えば制御点Ｃ３１）が一方向（例えばＡ１方向）に微小量Ｌ移動される（図２６参照）。

さらに、ステップＳ５８では、ステップＳ５７において移動対象点を移動させ一時的に変形させた状態のモデル（以下、「一時変形モデル」とも称する）の総合エネルギーＵｅが算出される。

総合エネルギーＵｅは、式（１）に示されるように、制御点Ｃｊと対応点Ｑｊとの距離に関する外部エネルギー項Ｆｅと、過剰な変形を避けるための内部エネルギー項Ｇｅと接触点Ｔｋと当該接触点Ｔｋに対応する制御点Ｃｋとの距離に関する拘束エネルギー項Ｈｅとの和で表される。なお、外部エネルギー項（単に外部エネルギーとも称する）Ｆｅ、内部エネルギー項（単に内部エネルギーとも称する）Ｇｅ、および拘束エネルギー項（単に拘束エネルギーとも称する）Ｈｅについては後述する。

次に、ステップＳ５９において、移動対象点が全ての方向に移動されたか否かを判定する。例えば、３次元空間の全方位における移動対象点の移動方向を２６方向（当該移動対象点（ボクセル）に隣接する２６画素へ向かう方向）とすると、２６方向全ての方向に制御点Ｃ１が移動されたか否かを判定する。

移動対象点（制御点Ｃ１）を２６方向全てに移動させた一時変形モデルの作成が終了していなければ、移動対象点の移動方向を変更して再び微小量Ｌ移動させ、異なる移動方向パターンの一時変形モデルを作成し、各一時変形モデルの総合エネルギーＵｅを算出する（ステップＳ５７、Ｓ５８）。

そして、ステップＳ５９において全方向の移動が終了したと判定されると、ステップＳ６０へ移行する。

ステップＳ６０では、作成された全パターンの一時変形モデルの中から、総合エネルギーＵｅを最小とする一時変形モデルが選択される。換言すれば、或る制御点Ｃｊを２６方向に移動させて生成された各一時変形モデルのうち、総合エネルギーＵｅを最小化する一時変形モデルが選択される。例えば、後述する３つの項Ｆｅ，Ｇｅ，Ｈｅのうち外部エネルギーＦｅのみを考慮する場合には、その制御点が対応点に近づく方向へと移動する一時変形モデルが選択されることになる。また、内部エネルギーＧｅを考慮する場合には、その制御点が対応点に近づく方向とは異なる方向へと移動する一時変形モデルが選択されることもある。

次に、ステップＳ６１において全制御点Ｃｊの移動が終了したか否かを判定する。具体的には、標準モデルＳＢ１の全ての制御点Ｃｊについて微小量Ｌの移動が終了したか否かを判定し、終了していない制御点（未了点とも称する）が存在すれば、移動対象点を当該制御点（未了点）に変更してステップＳ５７〜Ｓ６０の動作を繰り返し、全ての制御点Ｃｊの移動を完了した一時変形モデルを作成する。一方、全ての制御点の移動が終了していれば、ステップＳ６２へ移行する。

ステップＳ６２では、モデルフィッティング処理を終了するか否かを判定する。具体的には、全ての制御点Ｃｊの移動を完了させた一時変形モデルにおける複数の制御点とその対応点との距離の平均値が所定値以下であることを条件とし、当該条件を満たす場合に、モデルフィッティング処理を終了するようにすればよい。これによれば、各制御点が対応点に所定程度近づいた場合に同処理を終了することができる。または、これに加えてあるいはこれに代えて、（全ての制御点Ｃｊの移動を完了させた）前回の一時変形モデルと今回の一時変形モデルとの総合エネルギーＵｅの変化量が所定量以下であるか否かを終了判定の基準として用いてもよい。これによれば、制御点を移動しても総合エネルギーがあまり変化しなくなった場合に処理を終了することができる。

ステップＳ６２においてモデルフィッティング処理を終了しないと判定される場合は、ステップＳ６３へと移行する。

ステップＳ６３では、上述のステップＳ５７〜Ｓ６２で実行される処理を単位処理ループとして当該単位処理ループが所定回数Ｗ（例えば１０回）実行されたか否かを判定する。所定回数Ｗ実行されていなければ、所定回数Ｗ実行されるまで再びステップＳ５７〜Ｓ６２の処理ループを繰り返し、所定回数Ｗ実行されていれば、ステップＳ６４へと移行する。すなわち、全ての制御点の移動を完了させた一時変形モデルが所定回数Ｗ作成されるまで、単位処理ループが繰り返される。

ステップＳ６４では、ステップＳ５６で決定した対応点の更新が行われる。具体的には、上述の処理ループによって所定回数Ｗ移動した後の各制御点Ｃｊの最も近傍に存在する境界画素（境界点）を各制御点それぞれの新しい対応点とする対応点の更新が行われ、再びステップＳ５７〜Ｓ６４の処理が繰り返し行われる。このような「対応点の更新」によれば、制御点の移動に伴って制御点の位置が変わる場合にも、対応点の適正化を図ることができる。

一方、ステップＳ６２においてモデルフィッティング処理を終了すると判定される場合は、最終的に得られた一時変形モデル（変形後の標準モデル）が、計測データの情報を反映したモデル、すなわち個別モデルとして決定される。そして、当該個別モデルの計測データ内での存在位置が背骨領域Ｂｒとして抽出される。

＜３−５． 総合エネルギー（評価関数）＞
ここで、総合エネルギーＵｅを構成する外部エネルギーＦｅ、内部エネルギーＧｅおよび拘束エネルギーＨｅについて説明する。

外部エネルギーＦｅは、式（２）に示されるように、標準モデルにおける制御点（代表点）Ｃｊと計測データにおける対応点Ｑｊとの距離に関する要素を有している。なお、αは定数、Ｎｔは制御点の数を表すものとする。

このような外部エネルギーＦｅが大きくなるような一時変形モデル、すなわち、制御点Ｃｊと対応点Ｑｊとの距離が移動前より大きくなった一時変形モデルは、その総合エネルギーＵｅが大きくなるため、上述のステップＳ６０（総合エネルギーＵｅを最小とする制御点の移動を採用する工程）において採用されにくくなる。逆に言えば、外部エネルギーＦｅを考慮することによって、各制御点Ｃｊと対応点Ｑｊとの距離が移動前よりも小さくなる（すなわち各制御点Ｃｊが対応点Ｑｊに近づく）一時変形モデルが選択されやすくなる。

外部エネルギーＦｅは、標準モデルＳＢの各制御点Ｃｊと計測データにおける各対応点Ｑｊとの位置関係を反映したエネルギー項であるとも表現される。あるいは、外部エネルギーＦｅは、標準モデルの各制御点Ｃｊと計測データにおける各対応点Ｑｊとの距離が近づくときに最適化されるエネルギー項であるとも表現される。端的に言えば、この外部エネルギーＦｅは、標準モデルＳＢの形状を計測データにおける輪郭に近づけようとする役割を果たす。

また、内部エネルギーＧｅは、例えば、図２７に示されるように、制御点Ｃｊ間が仮想バネＳＰＲ（ＳＰＲ１，ＳＰＲ２，ＳＰＲ３，...）によって繋がれていると想定すると、式（３）のように表される。

但し、βは定数、Ｋは仮想バネのバネ係数、Ｎｈは仮想バネの本数、ｗｉは第ｉ番目の仮想バネの自然長からの変位量を表すものとする。

式（３）によると、各制御点Ｃｊの過剰な移動は、仮想バネＳＰＲに蓄えられるエネルギーの増大として表現される。例えば、１つの制御点Ｃｚが、或る点Ｖｚへと移動し他の制御点との相対変位が増大したとすると、仮想バネＳＰＲ１，ＳＰＲ２及びＳＰＲ３には、各仮想バネの伸びによるエネルギーが蓄えられ内部エネルギーＧｅひいては総合エネルギーＵｅが大きくなる。

このような過剰変形を伴う一時変形モデルは、その内部エネルギーＧｅが大きくなり、その総合エネルギーＵｅも大きくなるため、上述のステップＳ６０（総合エネルギーＵｅを最小とする制御点の移動を採用する工程）において採用されにくくなる。

換言すれば、内部エネルギーＧｅを減少させて総合エネルギーＵｅを減少させるような一時変形モデルがステップＳ６０で選択されることによって、各制御点Ｃｊの移動による過剰な変形を防止する作用を得ることができる。

つまり、このような内部エネルギーＧｅを導入することによって、標準モデルＳＢの形状すなわち標準モデルＳＢを構成する各ポリゴンの形状を損なわない制御点Ｃｊの移動が可能となる。なお、内部エネルギーＧｅを考慮する場合には、各制御点Ｃｊが対応点Ｑｊに近づく方向とは異なる方向へと移動する一時変形モデルが選択されることもある。

また、ステップＳ５６において対応点Ｑｊを求める際の境界抽出処理が正確でないなどの理由によって、対応点Ｑｊが正しいものでない場合が存在する。このような場合において、外部エネルギーＦｅのみを考慮するときには、或る制御点が、不正確な対応点に向けて移動するため、不正確なモデルフィッティングが行われることになる。これに対して、内部エネルギーＧｅを考慮することによれば、このような場合においても、標準モデルＳＢの形状を損なわないような変形動作が行われることになるため、より正確なモデリング処理が可能になる。

この内部エネルギーＧｅは、標準モデルＳＢの各制御点Ｃｊ相互の関係が一定の関係（すなわち、全仮想バネが自然長を有する状態）に近づくときに最適化されるエネルギー項であるとも表現される。端的に言えば、この内部エネルギーＧｅは、標準モデルＳＢを一定の形状に維持しようとする役割を果たす。

また、拘束エネルギーＨｅは、式（４）に示されるように、接触点Ｔｋと当該接触点Ｔｋに対応する制御点Ｃｋとの距離に関する要素を有している。なお、γは定数、Ｎｓは接触点の数を表すものとする。

拘束エネルギーＨｅは、計測データにおける２つの領域の接触部分において仮想的に設定された接触点Ｔｋと標準モデルＳＢの対応制御点Ｃｋとの位置関係を反映したエネルギー項であり、接触点Ｔｋと標準モデルＳＢの対応制御点Ｃｋとの距離が近づくときに最適化されるエネルギー項である。端的に言えば、この拘束エネルギーＨｅは、標準モデルＳＢにおいて接触部分に対応する部分の形状を計測データにおける接触点Ｔｋに近づけようとする役割を果たす。つまり、このような拘束エネルギーＨｅによれば、接触部分におけるモデル変形を精度よく行うことが可能となる。

なお、標準モデルＳＢ１が初期配置された際における接触点Ｔｋと対応制御点Ｃｋとの距離が所定値以上ある場合には、当該接触点Ｔｋは接触点として適正でないとして排除しても良い。これによれば、例えば接触部分以外で仮想平面ＩＦと混在領域ＫＲとが交差する場合に、接触点として不適切な点（本来の接触部分以外に存在する点）を拘束エネルギーＨｅの算出から有効に排除することが可能になる。

＜３−６．境界設定処理（ステップＳ６）＞
次に、ステップＳ６の境界設定処理について説明する。

ステップＳ６においては、上記のステップＳ５で抽出された背骨領域Ｂｒ１の境界を、「拡張限界」として設定する。

図２８は、ステップＳ５の終了時点におけるスライス画像ＳＣを示す図である。図２８においては、モデルフィッティング処理によって抽出された領域Ｂｒ１にクロスハッチングが付されて示されている。なお、ステップＳ５終了時点における抽出状態は、図８においても（立体的に）示されている。

また、図２９は、図２８の一部を拡大して示す図である。図２９に示されるように、抽出された背骨領域Ｂｒの境界（境界面ないし境界線）ＢＤが、第２の領域拡張法（ステップＳ７）における「拡張限界」として設定される。具体的には、変形後の標準モデル（個別モデル）の計測データにおける存在位置が背骨領域Ｂｒであることから、変形後の標準モデル（個別モデル）ＳＢの境界上に存在する画素（ボクセル）ＶＸが、「拡張限界」を表すものとして設定される。

＜３−７．第２の領域拡張処理（ステップＳ７）＞
ステップＳ７においては、再び領域拡張法を用いた領域抽出処理が実行される。ＳＰＳ７の領域拡張処理（第２の領域拡張処理）にも基本的にはステップＳ３の領域拡張処理（第１の領域拡張処理）と同様の手法が用いられる。ただし、この第２の領域拡張処理は、背骨領域Ｂｒに基づく「拡張限界」の制約の下で行われる点で、第１の領域拡張処理と相違する。以下では、ステップＳ３との相違点を中心に図１１を参照しながら説明する。

具体的には、「拡張限界」が設定されているため、ステップＳ３２における処理が異なる。

ステップＳ３２においては、或る画素が新たな拡張候補画素として検出されるための条件としては、注目画素にさらに隣接する画素であること、および領域内に存在するか否かが未だ判定されていない画素であることだけが求められるのではなく、「背骨領域Ｂｒ１の境界画素でないこと」も加重的に求められる。このように、第２の領域拡張法は、「拡張限界」の制約の下で行われる。

この結果、ステップＳ７においては、血管領域Ｖｒ内の開始点ＳＰから拡張された血管領域Ｖｒのみが抽出される。ステップＳ７の領域拡張に際して、背骨領域Ｂｒの境界まで到達するとそれより先には進めない（すなわち背骨領域Ｂｒ内には進めない）。換言すれば、背骨領域Ｂｒの境界が「拡張限界」（あるいは拡張障壁）として機能して、当該境界を超えて領域を拡張することができない。したがって、背骨領域Ｂｒを抽出することなく、血管領域Ｖｒのみを抽出することが可能である。

また、ステップＳ３４の判定処理においては、検出された拡張候補画素の階調値が所定の閾値ＴＨ２以上であるという条件を満たすか否かが判定される。

閾値ＴＨ２が閾値ＴＨ１と同一の値であるときには、図３０に示すように、ステップＳ３で抽出された混合領域から背骨領域Ｂｒを除外した領域がステップＳ７において血管領域Ｖｒとして抽出されることになる。図３０は、ステップＳ７終了時点の（すなわち最終的な血管領域抽出後における）スライス画像ＳＣを示す図である。図３０においては、最終的に抽出された血管領域Ｖｒに斜線が付されて示されている。なお、ステップＳ７終了時点における抽出状態は、図９においても（立体的に）示されている。

なお、ステップＳ３での領域拡張法に用いられる閾値ＴＨ１とステップＳ７での領域拡張法に用いられる閾値ＴＨ２とは、互いに異なる値に設定されてもよい。例えば、閾値ＴＨ１を閾値ＴＨ２よりも大きな値に設定して（すなわち比較的厳しい基準で）ステップＳ３の領域抽出を行うようにしてもよい。これによれば、ノイズをできるだけ排除した状態で背骨領域Ｂｒ１等を抽出することができ、ひいてはステップＳ５のモデルフィッティング処理において良好に背骨領域を抽出することが可能になる。また、これは、逆に言えば、閾値ＴＨ２を閾値ＴＨ１よりも小さな値（低い値）に設定して（すなわち比較的緩やかな基準で）ステップＳ７の領域抽出を行うことを意味する。したがって、図３１に示すように、「ステップＳ３で抽出された混合領域から背骨領域Ｂｒを除外した領域」よりも多くの領域（血管領域Ｖｒ１１等をも含む領域）、すなわち多くの血管領域Ｖｒを抽出することが可能になる。図３１においては、本来存在している血管領域Ｖｒ１１がさらに抽出されている。このように、細部の血管をも検出することが可能になる。

以上のように、領域抽出装置１は、ステップＳ５のモデルフィッティング処理において、暫定的に抽出された混在領域ＫＲに縮退処理を施し、背骨領域Ｂｒ１と血管領域Ｖｒ１とが接触する部分の接触点Ｔｋを取得するとともに、標準モデルＳＢ表面上の制御点Ｃｊの中から、接触点Ｔｋに対応する対応制御点Ｃｋを決定する（図２２〜図２４等参照）。そして、接触点Ｔｋと対応制御点Ｃｋとの距離に関する要素を含む評価関数を最適化するように標準モデルＳＢを変形するので、背骨領域Ｂｒ１と血管領域Ｖｒ１とが接触する部分におけるモデリングを精度良く行うことが可能になる。この結果、背骨領域Ｂｒ１を高精度に抽出することが可能になり、血管領域Ｖｒ１も高精度に抽出することが可能になる。

＜４．変形例等＞
以上、この発明の実施の形態について説明したが、この発明は上記説明した内容のものに限定されるものではない。

例えば、上記実施形態においては、分断箇所（分断点）ＢＰを含む複数の接触点Ｔｋを取得し、複数の接触点Ｔｋに基づいて拘束エネルギーＨｅを算出していたが、これに限定されない。具体的には、分断点ＢＰ以外の接触点Ｔｋを算出せずに、分断点ＢＰのみに基づいて拘束エネルギーＨｅを算出してもよい。

なお、複数の接触点Ｔｋに基づいて拘束エネルギーＨｅを算出する方が、分断点ＢＰのみに基づいて拘束エネルギーＨｅを算出する場合よりも、２つの領域が接触する部分におけるモデル変形をさらに精度良く行うことが可能となる。

具体的には、図３２を用いて説明する。図３２は、複数の接触点Ｔｋに基づいて拘束エネルギーＨｅが算出される場合において変形された標準モデルＳＢ１（ＳＢ２とも表記する）、および分断点ＢＰのみに基づいて拘束エネルギーＨｅが算出される場合において変形された標準モデルＳＢ１（ＳＢ３とも表記する）を示す図である。

複数の接触点Ｔｋに基づいて拘束エネルギーＨｅが算出される場合には、複数の接触点Ｔｋのうちのいずれかの対応制御点Ｃｋが互いに異なる可能性が高くなり、対応制御点Ｃｋが複数存在する可能性が高くなる。そして、これらの複数の対応制御点Ｃｋは、モデル変形の際に、拘束エネルギーＨｅを小さくするようにそれぞれ移動される。このため、複数の接触点Ｔｋに基づいて拘束エネルギーＨｅが算出される場合には、図３２に示されるように、標準モデルＳＢ２の所定部分が接触部分ＴＬに添うように変形される可能性が高くなる。

一方、分断点ＢＰのみに基づいて拘束エネルギーＨｅが算出される場合には、分断点ＢＰに対応付けられる対応制御点Ｃｋは一点となり、当該対応制御点Ｃｋのみが、拘束エネルギーＨｅを小さくするように移動される。このため、分断点ＢＰのみに基づいて拘束エネルギーＨｅが算出される場合には、図３２に示されるように、標準モデルＳＢ３の一部（対応制御点）が接触部分ＴＬに近づくように変形される。

このように、複数の接触点Ｔｋに基づいて拘束エネルギーＨｅを算出する方が、分断点ＢＰのみに基づいて拘束エネルギーＨｅを算出する場合よりも、２つの領域が接触する部分におけるモデル変形をさらに精度良く行うことが可能となる。

また、上記実施形態におけるステップＳ６４では、ステップＳ５６で決定した対応点の更新を行っているが、これに加えてさらに接触点Ｔｋの対応制御点Ｃｋの更新を行ってもよい。具体的には、所定回数Ｗ移動した後の各制御点Ｃｊの中から、接触点Ｔｋの最も近傍に存在する制御点を接触点Ｔｋの対応制御点Ｃｋとする更新動作が行われ、再びステップＳ５７〜Ｓ６４の処理が繰り返し行われる。このような「対応点制御点の更新」によれば、制御点の移動に伴って制御点の位置が変わる場合にも、対応制御点の適正化を図ることができる。

また、上記実施形態においては、モデル化が困難な抽出対象領域（血管領域）を抽出するため、抽出対象領域以外の除外領域（背骨領域）についてモデルフィッティングを行う場合を例示したが、これに限定されない。具体的には、抽出対象領域のモデル化が可能であれば、抽出対象領域そのものを上述のモデルフィッティング手法を用いて抽出してもよい。

また、上記実施形態においては、開始点ＳＰ１としては１つの点を指定する場合を例示したが、これに限定されず、複数の点を開始点として設定するようにしてもよい。ただし、情報処理量減少の観点からは、開始点の数を減少させることが好ましい。

また、上記実施形態の拡張候補画素の検出工程（ステップＳ３２）においては、拡張開始画素に隣接する全画素（２６画素）を拡張候補画素とし、領域拡張処理を実行していたが、これに限定されない。例えば、拡張開始画素Ｍ１に対して±ｘ方向に隣接する画素（Ｍ２、Ｍ３）、±ｙ方向に隣接する画素（Ｍ４、Ｍ７）及び±ｚ方向に隣接する画素（Ｕ１、Ｄ１）の計６画素を拡張候補画素として検出してもよい。

また、上記実施形態においては、図１６のフローチャートに示すようにして総合エネルギーＵｅを最小化（エネルギー関数を最適化）する場合を例示したが、エネルギー関数最適化手法は、上記手法に限定されない。例えば、準ニュートン法或いはニュートン法等の数学的手法を用いて総合エネルギーＵｅを最適化するようにしてもよい。

Ｘ線ＣＴ装置によって取得される３次元画像を示す図である。 スライス画像を示す図である。 領域抽出装置の概要を示す図である。 領域抽出装置の各種機能を示すブロック図である。 領域抽出装置の全体動作を示すフローチャートである。 ３次元画像内の状態（初期状態）を示す概念図である。 第１の領域拡張処理終了時点での状態を示す概念図である。 モデルフィッティング処理終了時点での状態を示す概念図である。 第２の領域拡張処理終了時点での状態を示す概念図である。 全領域および血管領域の各ＣＴ値の度数分布を示す図である。 第１の領域拡張処理の詳細を示すフローチャートである。 隣接画素を示す図である。 混在領域抽出後のスライス画像を示す図である。 比較的簡略化された標準モデルを示す図である。 突起部の最下端位置がスライス画像ごとに変化する様子を示す図である。 モデルフィッティング処理を示すフローチャートである。 モデルフィッティング処理を示すフローチャートである。 暫定的に抽出された混在領域に含まれる背骨領域と血管領域とを示す図である。 背骨領域と血管領域とに縮退処理を施した図である。 分断直前の分断周辺情報を示す図である。 分断直前の分断周辺情報を示す図である。 混在領域に設定された仮想平面を示す図である。 計測データ中に配置された標準モデルを示す図である。 図２３の領域ＲＢ１付近の拡大図である。 計測データ中に配置された標準モデルを示す図である。 図２５の領域ＲＢ２付近の拡大図である。 制御点間を仮想バネで繋いだ様子を示す概念図である。 モデルフィッティング処理終了時点におけるスライス画像を示す図である。 図２８の一部拡大図である。 第２の領域拡張処理終了時点におけるスライス画像を示す図である。 変形例に係るスライス画像を示す図である。 標準モデルの変形の様子を示す図である。

符号の説明

１ 領域抽出装置
ＢＮ 背骨
Ｂｒ，Ｂｒ１〜Ｂｒ４ 背骨領域
ＶＬ，Ｖｒ１，Ｖｒ２ 血管領域
ＢＮ 背骨
ＢＰ 分断点（分断箇所）
Ｑｊ 対応点
Ｃｊ 制御点
Ｃｋ 対応制御点
Ｔｋ 接触点
ＩＦ 仮想平面
ＫＲ 混在領域
ＴＬ 接触部分

Claims (6)

1. ３次元計測データに基づいて対象物に関するモデリングを行うモデリング装置であって、
   前記対象物に関する標準モデルを取得する手段と、
   前記３次元計測データから、前記対象物に対応する領域である対象領域と、当該対象領域以外の他の領域とを含む暫定領域を抽出する暫定領域抽出手段と、
   モデルフィッティング手法により前記標準モデルを変形することによって、前記対象物に関するモデルを作成するモデル作成手段と、
   を備え、
   前記モデル作成手段は、
   前記暫定領域に対して縮退処理を施す縮退手段と、
   前記暫定領域が前記縮退処理により分断される際の分断点を特定する分断点特定手段と、
   前記分断点を含み、前記対象領域と前記他の領域とのつながり方向に基づき規定される仮想面を設定する設定手段と、
   前記仮想面と前記暫定領域とが交差する画素を特定することによって、前記暫定領域において前記対象領域と前記他の領域とが接触する部分の仮想的な接触点を取得する接触点取得手段と、
   前記標準モデル表面上の複数の代表点の中から、前記接触点に対応する対応代表点を決定する対応点決定手段と、
   前記接触点と前記対応代表点との距離に関する要素を含む評価関数を最適化するように前記標準モデルを変形させる変形手段と、
   を有することを特徴とするモデリング装置。
2. 請求項１に記載のモデリング装置において、
   前記モデル作成手段は、
   前記分断点近傍の周辺情報を取得する手段、
   をさらに有し、
   前記設定手段は、前記分断点と前記周辺情報とに基づいて、前記分断点における前記対象領域と前記他の領域とのつながり方向を表すベクトルを算出し、前記ベクトルを法線ベクトルとし且つ前記分断点を含む仮想平面を前記仮想面として設定することを特徴とするモデリング装置。
3. 請求項１または請求項２に記載のモデリング装置において、
   前記暫定領域抽出手段は、領域拡張手法により前記暫定領域を抽出することを特徴とするモデリング装置。
4. 計測データに基づいて対象物に関するモデリングを行うモデリング装置であって、
   前記対象物に関する標準モデルを取得する手段と、
   前記計測データから、前記対象物に対応する領域である対象領域を含む暫定領域を抽出する暫定領域抽出手段と、
   モデルフィッティング手法により前記標準モデルを変形することによって、前記対象物に関するモデルを作成するモデル作成手段と、
   を備え、
   前記モデル作成手段は、
   前記暫定領域に対して縮退処理を施す縮退手段と、
   前記暫定領域が前記縮退処理により分断される際の分断点を特定する分断点特定手段と、
   前記標準モデル表面上の複数の代表点の中から、前記分断点に対応する対応代表点を決定する対応点決定手段と、
   前記分断点と前記対応代表点との距離に関する要素を含む評価関数を最適化するように前記標準モデルを変形させる変形手段と、
   を有することを特徴とするモデリング装置。
5. ３次元計測データに基づいて対象物に関するモデリングを行うモデリング方法であって、
   前記対象物に関する標準モデルを取得する工程と、
   前記３次元計測データから、前記対象物に対応する領域である対象領域と、当該対象領域以外の他の領域とを含む暫定領域を抽出する暫定領域抽出工程と、
   モデルフィッティング手法により前記標準モデルを変形することによって、前記対象物に関するモデルを作成するモデル作成工程と、
   を備え、
   前記モデル作成工程は、
   前記暫定領域に対して縮退処理を施す縮退処理工程と、
   前記暫定領域が前記縮退処理により分断される際の分断点を特定する分断点特定工程と、
   前記分断点を含み、前記対象領域と前記他の領域とのつながり方向に基づき規定される仮想面を設定する工程と、
   前記仮想面と前記暫定領域とが交差する画素を特定することによって、前記暫定領域において前記対象領域と前記他の領域とが接触する部分の仮想的な接触点を取得する接触点取得工程と、
   前記標準モデル表面上の複数の代表点の中から、前記接触点に対応する対応代表点を決定する対応点決定工程と、
   前記接触点と前記対応代表点との距離に関する要素を含む評価関数を最適化するように前記標準モデルを変形させる変形処理工程と、
   を有することを特徴とするモデリング方法。
6. コンピュータに、
   対象物に関する標準モデルを取得する工程と、
   ３次元計測データから、前記対象物に対応する領域である対象領域と、当該対象領域以外の他の領域とを含む暫定領域を抽出する暫定領域抽出工程と、
   モデルフィッティング手法により前記標準モデルを変形することによって、前記対象物に関するモデルを作成するモデル作成工程と、
   を実行させるためのプログラムであって、
   前記モデル作成工程は、
   前記暫定領域に対して縮退処理を施す縮退処理工程と、
   前記暫定領域が前記縮退処理により分断される際の分断点を特定する分断点特定工程と、
   前記分断点を含み、前記対象領域と前記他の領域とのつながり方向に基づき規定される仮想面を設定する工程と、
   前記仮想面と前記暫定領域とが交差する画素を特定することによって、前記暫定領域において前記対象領域と前記他の領域とが接触する部分の仮想的な接触点を取得する接触点取得工程と、
   前記標準モデル表面上の複数の代表点の中から、前記接触点に対応する対応代表点を決定する対応点決定工程と、
   前記接触点と前記対応代表点との距離に関する要素を含む評価関数を最適化するように前記標準モデルを変形させる変形処理工程と、
   を有することを特徴とするプログラム。

Images