(57)【要約】【課題】より正確に対象物をモデリングすることが可能なモデリング技術を提供する。【解決手段】モデリング装置は、所定の評価関数を最適化するように対象物に関する標準モデルを変形して、計測データに基づき対象物のモデルを生成する。所定の評価関数は、標準モデルの制御点Ｃｊと計測データにおいて当該制御点Ｃｊに対応する対応点との位置関係を反映したものである。対応点決定に際しては、計測データに含まれる複数の部分領域Ｂｒ１，Ｖｒ１，...のそれぞれについて対象物らしさを示す指標値が算出され、複数の部分領域Ｂｒ１，Ｖｒ１，...のうち、指標値に関する所定の基準を満たす部分領域Ｂｒ１に属する点の中から対応点Ｑｂ１が選択される。【選択図】図２１
2007312838.tif 000002
【特許請求の範囲】
【請求項１】
対象物のモデルを生成するモデリング装置であって、
前記対象物に関する計測データを取得するデータ取得手段と、
所定の評価関数を最適化するように前記対象物に関する標準モデルを変形して、前記計測データに基づき前記対象物のモデルを生成するモデル生成手段と、
を備え、
前記所定の評価関数は、前記標準モデルの制御点と前記計測データにおいて前記制御点に対応する対応点との位置関係を反映したものであり、
前記モデル生成手段は、前記計測データに含まれる複数の部分領域のそれぞれについて前記対象物らしさを示す指標値を算出し、前記複数の部分領域のうち前記指標値に関する所定の基準を満たす部分領域である基準充足領域を求め、前記基準充足領域に属する点の中から前記対応点を選択することを特徴とするモデリング装置。
【請求項２】
請求項１に記載のモデリング装置において、
前記計測データは、複数のスライス画像の集合として構成される３次元画像であり、
前記モデル生成手段は、前記計測データの各スライス画像に含まれる複数の部分領域のそれぞれについて前記対象物らしさを示す指標値を算出し、当該複数の部分領域のうち前記指標値に関する所定の基準を満たす部分領域である基準充足領域を各スライス画像ごとに求め、前記複数のスライス画像の前記基準充足領域に属する点の中から前記対応点を選択することを特徴とするモデリング装置。
【請求項３】
請求項２に記載のモデリング装置において、
各部分領域の前記指標値は、その部分領域に属する画素の階調値のばらつき度合いを反映した値であることを特徴とするモデリング装置。
【請求項４】
請求項２または請求項３に記載のモデリング装置において、
各部分領域の前記指標値は、その部分領域に属する画素の最大階調値を反映した値であることを特徴とするモデリング装置。
【請求項５】
請求項２から請求項４のいずれかに記載のモデリング装置において、
各部分領域の前記指標値は、その部分領域の面積を反映した値であることを特徴とするモデリング装置。
【請求項６】
請求項２から請求項４のいずれかに記載のモデリング装置において、
各部分領域の前記指標値は、その部分領域の円形度を反映した値であることを特徴とするモデリング装置。
【請求項７】
請求項１から請求項６のいずれかに記載のモデリング装置において、
前記モデル生成手段は、前記基準充足領域に属する点であって且つ前記制御点に最も近い点を、前記対応点として決定することを特徴とするモデリング装置。
【請求項８】
対象物のモデルを生成するモデリング方法であって、
（ａ）前記対象物に関する計測データを取得するステップと、
（ｂ）所定の評価関数を最適化するように前記対象物に関する標準モデルを変形して、前記計測データに基づき前記対象物のモデルを生成するステップと、
を備え、
前記所定の評価関数は、前記標準モデルの制御点と前記計測データにおいて前記制御点に対応する対応点との位置関係を反映したものであり、
前記ステップ（ｂ）においては、前記計測データに含まれる複数の部分領域のそれぞれについて前記対象物らしさを示す指標値が算出され、前記複数の部分領域のうち前記指標値に関する所定の基準を満たす部分領域である基準充足領域が求められ、前記基準充足領域に属する点の中から前記対応点が選択されることを特徴とするモデリング方法。
【請求項９】
コンピュータに、
（ａ）対象物に関する計測データを取得するステップと、
（ｂ）前記対象物に関する標準モデルを所定の評価関数を最適化するように変形して、前記計測データに基づき前記対象物のモデルを生成するステップと、
を実行させるためのプログラムであって、
前記所定の評価関数は、前記標準モデルの制御点と前記計測データにおいて前記制御点に対応する対応点との位置関係を反映したものであり、
前記ステップ（ｂ）においては、前記計測データに含まれる複数の部分領域のそれぞれについて前記対象物らしさを示す指標値が算出され、前記複数の部分領域のうち前記指標値に関する所定の基準を満たす部分領域である基準充足領域が求められ、前記基準充足領域に属する点の中から前記対応点が選択されることを特徴とするプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、モデルフィッティング手法によって計測データと標準モデルとに基づいて対象物のモデルを生成する技術に関する。
【背景技術】
【０００２】
近年、医療分野において、Ｘ線ＣＴ(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座標の相関値が閾値を越えない。
【００７０】
条件ＣＤ１は、末梢血管を排除し、大動脈内の１つの点を開始点として選択することを意図するものである。
【００７１】
また、図１０に示すように、血管領域の最大ＣＴ値（各スライス画像における血管領域内に含まれる全画素のＣＴ値のうちの最大値）は、（一般的に骨領域の最大ＣＴ値よりも小さく）所定の閾値（例えばＴＨ９）以下であることが多い。条件ＣＤ２はこのような特性を利用して、骨領域を排除して血管領域を抽出することを意図している。なお、図１０は、或る３次元画像について、全領域のＣＴ値の度数分布と、血管領域のＣＴ値の度数分布とを示す図である。
【００７２】
また、血管領域のＣＴ値の標準偏差は（一般的に骨領域のＣＴ値の標準偏差よりも小さく）所定の閾値以下であることが多い。条件ＣＤ３は、このような特性を利用して、骨領域を排除して血管領域を抽出することを意図している。
【００７３】
また、背骨に沿って走る大動脈の血管領域Ｖｒ１は、各スライス画像において円形に近い。円形度に関する条件ＣＤ４，ＣＤ５は、このような特性を利用して血管領域Ｖｒ１を選択することを意図するものである。
【００７４】
そして、部分領域がこれらの条件ＣＤ１〜ＣＤ５の全てを満たす場合には当該部分領域の中心座標を求め、当該中心座標を開始点として決定して開始点決定処理を終了する。
【００７５】
一方、部分領域がこれらの条件ＣＤ１〜ＣＤ５のいずれか１つでも満たさないときには、次の部分領域を検索して同様の判定動作を再び行う。また、１つのスライス画像に上記の条件ＣＤ１〜ＣＤ５の全てを満たす部分領域が１つも存在しない場合には、別のスライス画像について同様の処理を行う。
【００７６】
このような処理によって、上記の条件を満たす開始点ＳＰが決定される。なお、図６等に示すように、背骨領域Ｂｒと血管領域Ｖｒとは全てのスライス画像で接しているのではなく、両領域Ｂｒ，Ｖｒが互いに接していないスライス画像も存在する。上記の条件ＣＤ１〜ＣＤ５を考慮することによれば、両領域Ｂｒ，Ｖｒが互いに接していないスライス画像における血管領域Ｖｒ内の点が開始点として決定される。
【００７７】
＜３−２．第１の領域拡張処理（ステップＳ３）＞
上述のように、ステップＳ３の第１の領域拡張処理では、領域拡張法を採用して、血管領域Ｖｒと背骨領域Ｂｒとの双方を含む混在領域を暫定的に抽出する。
【００７８】
図１１は、領域拡張法を用いた暫定領域抽出処理（ステップＳ３）の詳細を示すフローチャートである。図１２は、スライス画像ＰＭ上の画素Ｍ１を中心として隣接する画素付近を拡大表示した図である。図１２においては、スライス画像ＰＭの画素Ｍ１の隣接画素として、同一階層のスライス画像ＰＭにおける隣接８画素Ｍ２〜Ｍ９に加えて、スライス画像ＰＵにおける隣接９画素Ｕ１〜Ｕ９と、スライス画像ＰＤにおける隣接９画素Ｄ１〜Ｄ９とが示されている。なお、スライス画像ＰＵは、スライス画像ＰＭに対して＋ｚ方向に隣接するスライス画像（換言すればスライス画像ＰＭの直上層のスライス画像）であり、スライス画像ＰＤは、スライス画像ＰＭに対して−ｚ方向に隣接するスライス画像（換言すればスライス画像ＰＭの直下層のスライス画像）である。
【００７９】
まず、ステップＳ３１において、ステップＳ２で決定された開始点ＳＰ１が拡張開始画素（注目画素とも称する）に設定される。
【００８０】
次に、ステップＳ３２において、注目画素に隣接する画素が拡張候補画素として検出される。詳細には、ステップＳ３１において注目画素として画素Ｍ１が選択されたとすると、画素Ｍ１に隣接する画素全てが拡張候補画素となる。つまり、スライス画像ＰＭ上の画素Ｍ２〜Ｍ９、スライス画像ＰＵ上の画素Ｕ１〜Ｕ９、及びスライス画像ＰＤ上の画素Ｄ１〜Ｄ９の計２６画素が拡張候補画素となる。
【００８１】
次に、ステップＳ３３において拡張候補画素がステップＳ３２で検出されたか否かが判定され、ステップＳ３２で拡張候補画素が検出されている場合には、ステップＳ３４へ移行する。
【００８２】
ステップＳ３４では、検出された拡張候補画素の階調値が所定の閾値ＴＨ１以上であるという条件を満たすか否かが判定される。拡張候補画素の階調値が所定の閾値ＴＨ１以上であると判定されると、ステップＳ３５に移行する。
【００８３】
ステップＳ３５では、当該拡張候補画素を、拡張後の領域（拡張領域）内の画素とする領域拡張を行う。
【００８４】
一方、ステップＳ３４において、拡張候補画素の階調値が所定の閾値ＴＨ１以上でないと判定されると、当該拡張候補画素の領域拡張は行われない。
【００８５】
次に、ステップＳ３６では、拡張候補画素の階調値判定工程（ステップＳ３４）をまだ経ていない拡張候補画素が存在するか否かが判断される。
【００８６】
階調値判定工程（ステップＳ３４）を経ていない拡張候補画素が存在する場合には、当該未終了の拡張候補画素に対してステップＳ３４〜Ｓ３６の処理が実行される。
【００８７】
一方、全ての拡張候補画素に対して階調値判定工程（ステップＳ３４）が終了している場合には、ステップＳ３２に戻る。そして、ステップＳ３２においては、ステップＳ３５において拡張された画素（換言すれば、領域内の画素であるとして特定された画素）が新たな注目画素として設定されるとともに、当該注目画素にさらに隣接する画素であって、領域内に存在するか否かが未だ判定されていない画素が、新たな拡張候補画素として検出される。その後、新たに検出された拡張候補画素に対して、上述のステップＳ３３〜Ｓ３６の工程が実行される。ステップＳ３２〜Ｓ３６の工程は、新たな拡張候補画素を検出することが可能な限り繰り返され、新たな拡張候補画素が検出できなくなると、暫定領域抽出処理は終了する（ステップＳ３３）。
【００８８】
このようにステップＳ３２〜Ｓ３６の工程が繰り返し実行されることによって開始点ＳＰから徐々に領域拡張が行われ、背骨領域Ｂｒと血管領域Ｖｒとが混在した領域（混在領域）が暫定的に抽出される。
【００８９】
図１３は、混在領域抽出後のスライス画像ＳＣを示す図である。図１３においては、抽出された領域に斜線が付されて示されている。なお、ステップＳ３終了時点における抽出状態は、図７においても（立体的に）示されている。
【００９０】
＜３−３．標準モデルの初期位置決定処理（ステップＳ４）＞
次に、モデルフィッティング手法における標準モデルの初期位置を決定する初期位置決定処理（ステップＳ４）について図１３および図１４を参照しながら説明する。
【００９１】
この実施形態においては、次のステップＳ５のモデルフィッティング処理用の背骨ＢＮの標準モデルＳＢとして、図１４に示すような比較的簡略化された標準モデルＳＢ１を採用する。この標準モデルＳＢ１は、背骨の各椎骨の概略形状を円柱とみなしてモデル化されたものである。標準モデルＳＢ１は、複数の椎骨部分を表現する各円柱ＣＬと、当該各円柱ＣＬ間を接続する曲面ＡＰとによって構成されている。また、当該円柱部分（背骨領域Ｂｒ１に相当）以外の部分、具体的には、当該円柱部分から外側に突出した突起部分（背骨領域Ｂｒ２に相当）はモデル化の対象から除外されている。
【００９２】
なお、血管領域Ｖｒと接するのは、椎骨の略円柱状部分に対応する背骨領域Ｂｒ１であることが多く、他の背骨領域Ｂｒ２〜Ｂｒ４と血管領域Ｖｒとが接することはほとんど無いという事情がある。このような事情のため、簡易な標準モデルＳＢ１を用いてモデルフィッティングを行って背骨領域Ｂｒ１を抽出し、当該背骨領域を排除するような第２の領域拡張処理（ステップＳ７）を行うことによっても、血管領域Ｖｒのみを抽出することが可能である。
【００９３】
さて、標準モデルＳＢ１は、背骨領域Ｂｒ１に対応する椎骨部分を中心にモデル化されたものであるため、背骨領域Ｂｒ１に近い位置に標準モデルＳＢ１を初期配置することが好ましい。また、背骨領域Ｂｒ２は突起部分を有しているため背骨領域Ｂｒ２の位置（例えば背骨領域Ｂｒ２の下端点ＰＴ２）を特定することは比較的容易である。そこで、当該背骨領域Ｂｒ２の位置を利用することに加えて、背骨領域Ｂｒ１と背骨領域Ｂｒ２との標準的な位置関係をもさらに利用して、標準モデルＳＢ１の初期位置を決定する。
【００９４】
図１３に示すように、領域Ｂｒ２は、混在領域のうち、スライス画像内において最も下方に存在している。したがって、混在領域のうち、各スライス画像内において最も下方に存在する点が領域Ｂｒ２の下端ＰＴ２の候補となる。
【００９５】
より詳細には、各スライス画像について、ステップＳ３で抽出された混在領域のＹ方向における最下端位置を検出する。そして、図１５に示すように、当該最下端位置が複数のスライス画像にわたって連続的に変化するとき、連続的変化区間Ｋｉの中央のスライス画像における最下端位置（×印）を、当該区間Ｋｉに対応する突出部分の基準点の位置として決定する。また、このような突出部分の基準点（例えばＰＴ２（図１３））と当該突出部分に対応する椎骨の基準点（例えばＰＴ１（図１３））との標準的な位置関係に基づいて、各椎骨の基準点（例えば位置ＰＴ１）の推定座標値が特定される。このような推定座標値は複数の椎骨のそれぞれについて求められる。その後、各椎骨の基準点の推定座標値に対して、標準モデルＳＢ１の各対応制御点（標準モデルＳＢ１において各椎骨の基準点に対応する制御点（後述））の位置が一致するように、標準モデルＳＢ１を移動する。
【００９６】
このようにして標準モデルＳＢ１の初期配置が行われる。
【００９７】
＜３−４．モデルフィッティング処理（ステップＳ５）＞
ステップＳ５のモデルフィッテイング処理は、初期配置された「標準モデル」を、計測データにおけるモデリング対象領域（除外領域）に関する情報（形状等）を用いて変形する処理である。なお、本出願においては、モデルフィッティング処理による変形後の標準モデル（換言すれば、計測データの情報が反映された標準モデル）を「個別モデル」とも称するものとする。
【００９８】
図１６は、モデルフィッティング処理（ステップＳ５）の詳細を示すフローチャートである。図１７は、標準モデルＳＢおよびその各制御点Ｃｊを示す概念図であり、図１８は、計測データ中に配置された標準モデルＳＢを示す概念図である。図１９は、図１８の領域ＲＢ付近の拡大図であり、モデルフィッティングにおける制御点の移動の様子を示している。また、図２０は、制御点間を仮想バネで繋いだ様子を示す概念図である。
【００９９】
なお、モデリング対象物体（ここでは背骨）の標準モデルは、例えば、微小な多角形（例えば、三角形）（「ポリゴン」とも称する）で構成され、記憶部２２等に保存されている。ポリゴンで構成された標準モデルは、各ポリゴンの頂点の３次元座標によってそのモデルの表面形状を表現することができる。また、標準モデルにおけるポリゴンの頂点のうち、全頂点あるいは代表的な幾つかの頂点を「制御点」とも称する。そして、この制御点を代表点としてモデルフィッティング処理を行う。
【０１００】
以下では、図１６〜図２０を用いてモデルフィッティング処理について詳述する。このモデルフィッティング処理においては、標準モデルＳＢを１制御点ごとに微小量Ｌずつ徐々に変形させることによって、個別モデルが作成される。
【０１０１】
まず、ステップＳ５１（図１６）では、計測データにおけるモデリング対象物体ＯＢａ（図１８参照）の輪郭（境界）を示す画素（境界点）のうち、標準モデルＳＢの各制御点Ｃｊに対応する点（以下、「対応点」とも称する）Ｑｊを決定する処理が行われる（図１８参照）。すなわち、標準モデルＳＢにおける各制御点（代表点）Ｃｊと計測データにおける対応点Ｑｊとの対応関係が決定される。例えば、図１８および図１９に示すように、標準モデルの初期配置後のデータ空間（３次元空間）において、制御点Ｃｊの最も近傍に存在する境界点（画素）Ｑｊが、制御点Ｃｊの対応点として求められる。なお、このステップＳ５１の処理については後に詳述する。
【０１０２】
次に、ステップＳ５２では、標準モデルＳＢの制御点Ｃｊのうち任意の一点（以下、「移動対象点」とも称する）（例えば制御点Ｃ１）が一方向（例えばＡ１方向）に微小量Ｌ移動される（図１９参照）。
【０１０３】
さらに、ステップＳ５３では、ステップＳ５２において移動対象点を移動させ一時的に変形させた状態のモデル（以下、「一時変形モデル」とも称する）の総合エネルギーＵｅが算出される。
【０１０４】
総合エネルギーＵｅは、式（１）に示されるように、制御点Ｃｊと対応点Ｑｊとの距離に関する外部エネルギー項Ｆｅと、過剰な変形を避けるための内部エネルギー項Ｇｅとの和で表される。なお、外部エネルギー項（単に外部エネルギーとも称する）Ｆｅ、および内部エネルギー項（単に内部エネルギーとも称する）Ｇｅについては後述する。
【０１０５】
【数１】
2007312838.tif 000003
【０１０６】
次に、ステップＳ５４において、移動対象点が全ての方向に移動されたか否かを判定する。例えば、３次元空間の全方位における移動対象点の移動方向を２６方向（当該移動対象点（ボクセル）に隣接する２６画素へ向かう方向）とすると、２６方向全ての方向に制御点Ｃ１が移動されたか否かを判定する。
【０１０７】
移動対象点（制御点Ｃ１）を２６方向全てに移動させた一時変形モデルの作成が終了していなければ、移動対象点の移動方向を変更して再び微小量Ｌ移動させ、異なる移動方向パターンの一時変形モデルを作成し、各一時変形モデルの総合エネルギーＵｅを算出する（ステップＳ５２、Ｓ５３）。
【０１０８】
そして、ステップＳ５４において全方向の移動が終了したと判定されると、ステップＳ５５へ移行する。
【０１０９】
ステップＳ５５では、作成された全パターンの一時変形モデルの中から、総合エネルギーＵｅを最小とする一時変形モデルが選択される。換言すれば、或る制御点Ｃｊを２６方向に移動させて生成された各一時変形モデルのうち、総合エネルギーＵｅを最小化する一時変形モデルが選択される。例えば、後述する２つの項Ｆｅ，Ｇｅのうち外部エネルギーＦｅのみを考慮する場合には、その制御点が対応点に近づく方向へと移動する一時変形モデルが選択されることになる。また、内部エネルギーＧｅを考慮する場合には、その制御点が対応点に近づく方向とは異なる方向へと移動する一時変形モデルが選択されることもある。
【０１１０】
次に、ステップＳ５６において全制御点Ｃｊの移動が終了したか否かを判定する。具体的には、標準モデルＳＢの全ての制御点Ｃｊについて微小量Ｌの移動が終了したか否かを判定し、終了していない制御点（未了点とも称する）が存在すれば、移動対象点を当該制御点（未了点）に変更してステップＳ５２〜Ｓ５５の動作を繰り返し、全ての制御点Ｃｊの移動を完了した一時変形モデルを作成する。一方、全ての制御点の移動が終了していれば、ステップＳ５７へ移行する。
【０１１１】
ステップＳ５７では、モデルフィッティング処理を終了するか否かを判定する。具体的には、全ての制御点Ｃｊの移動を完了させた一時変形モデルにおける複数の制御点とその対応点との距離の平均値が所定値以下であることを条件とし、当該条件を満たす場合に、モデルフィッティング処理を終了するようにすればよい。これによれば、各制御点が対応点に所定程度近づいた場合に同処理を終了することができる。または、これに加えてあるいはこれに代えて、（全ての制御点Ｃｊの移動を完了させた）前回の一時変形モデルと今回の一時変形モデルとの総合エネルギーＵｅの変化量が所定量以下であるか否かを終了判定の基準として用いてもよい。これによれば、制御点を移動しても総合エネルギーがあまり変化しなくなった場合に処理を終了することができる。
【０１１２】
ステップＳ５７においてモデルフィッティング処理を終了しないと判定される場合は、ステップＳ５８へと移行する。
【０１１３】
ステップＳ５８では、上述のステップＳ５２〜Ｓ５７で実行される処理を単位処理ループとして当該単位処理ループが所定回数Ｗ（例えば１０回）実行されたか否かを判定する。所定回数Ｗ実行されていなければ、所定回数Ｗ実行されるまで再びステップＳ５２〜Ｓ５７の処理ループを繰り返し、所定回数Ｗ実行されていれば、ステップＳ５９へと移行する。すなわち、全ての制御点の移動を完了させた一時変形モデルが所定回数Ｗ作成されるまで、単位処理ループが繰り返される。
【０１１４】
ステップＳ５９では、ステップＳ５１で決定した対応点の更新が行われる。具体的には、上述の処理ループによって所定回数Ｗ移動した後の各制御点Ｃｊの周辺領域ＲＢ内を再設定し、当該周辺領域ＲＢ内に存在する点（ボクセル）の中から新しい対応点を選択する対応点の更新動作が行われ、再びステップＳ５２〜Ｓ５９の処理が繰り返し行われる。このような「対応点の更新」によれば、制御点の移動に伴って制御点の位置が変わる場合にも、対応点の適正化を図ることができる。
【０１１５】
一方、ステップＳ５７においてモデルフィッティング処理を終了すると判定される場合は、最終的に得られた一時変形モデル（変形後の標準モデル）が、計測データの情報を反映したモデル、すなわち個別モデルとして決定される。そして、当該個別モデルの計測データ内での存在位置が背骨領域Ｂｒとして抽出される。
【０１１６】
＜３−５．  総合エネルギー（評価関数）＞
ここで、総合エネルギーＵｅを構成する外部エネルギーＦｅ及び内部エネルギーＧｅについて説明する。
【０１１７】
外部エネルギーＦｅは、式（２）に示されるように、標準モデルにおける制御点（代表点）Ｃｊと計測データにおける対応点Ｑｊとの距離に関する要素を有している。なお、αは定数、Ｎｔは制御点の数を表すものとする。
【０１１８】
【数２】
2007312838.tif 000004
【０１１９】
このような外部エネルギーＦｅが大きくなるような一時変形モデル、すなわち、制御点Ｃｊと対応点Ｑｊとの距離が移動前より大きくなった一時変形モデルは、その総合エネルギーＵｅが大きくなるため、上述のステップＳ５５（総合エネルギーＵｅを最小とする制御点の移動を採用する工程）において採用されにくくなる。逆に言えば、外部エネルギーＦｅを考慮することによって、各制御点Ｃｊと対応点Ｑｊとの距離が移動前よりも小さくなる（すなわち各制御点Ｃｊが対応点Ｑｊに近づく）一時変形モデルが選択されやすくなる。
【０１２０】
外部エネルギーＦｅは、標準モデルＳＢの各制御点Ｃｊと計測データにおける各対応点Ｑｊとの位置関係を反映したエネルギー項であるとも表現される。あるいは、外部エネルギーＦｅは、標準モデルの各制御点Ｃｊと計測データにおける各対応点Ｑｊとの距離が近づくときに最適化されるエネルギー項であるとも表現される。端的に言えば、この外部エネルギーＦｅは、標準モデルＳＢの形状を計測データにおける輪郭に近づけようとする役割を果たす。
【０１２１】
また、内部エネルギーＧｅは、例えば、図２０に示されるように、制御点Ｃｊ間が仮想バネＳＰＲ（ＳＰＲ１，ＳＰＲ２，ＳＰＲ３，...）によって繋がれていると想定すると、式（３）のように表される。
【０１２２】
【数３】
2007312838.tif 000005
【０１２３】
但し、βは定数、Ｋは仮想バネのバネ係数、Ｎｈは仮想バネの本数、ｗｉは第ｉ番目の仮想バネの自然長からの変位量を表すものとする。
【０１２４】
式（３）によると、各制御点Ｃｊの過剰な移動は、仮想バネＳＰＲに蓄えられるエネルギーの増大として表現される。例えば、１つの制御点Ｃｚが、或る点Ｖｚへと移動し他の制御点との相対変位が増大したとすると、仮想バネＳＰＲ１，ＳＰＲ２及びＳＰＲ３には、各仮想バネの伸びによるエネルギーが蓄えられ内部エネルギーＧｅひいては総合エネルギーＵｅが大きくなる。
【０１２５】
このような過剰変形を伴う一時変形モデルは、その内部エネルギーＧｅが大きくなり、その総合エネルギーＵｅも大きくなるため、上述のステップＳ５５（総合エネルギーＵｅを最小とする制御点の移動を採用する工程）において採用されにくくなる。
【０１２６】
換言すれば、内部エネルギーＧｅを減少させて総合エネルギーＵｅを減少させるような一時変形モデルがステップＳ５５で選択されることによって、各制御点Ｃｊの移動による過剰な変形を防止する作用を得ることができる。
【０１２７】
つまり、このような内部エネルギーＧｅを導入することによって、標準モデルＳＢの形状すなわち標準モデルＳＢを構成する各ポリゴンの形状を損なわない制御点Ｃｊの移動が可能となる。なお、内部エネルギーＧｅを考慮する場合には、各制御点Ｃｊが対応点Ｑｊに近づく方向とは異なる方向へと移動する一時変形モデルが選択されることもある。
【０１２８】
また、ステップＳ５１において対応点Ｑｊを求める際の境界抽出処理が正確でないなどの理由によって、対応点Ｑｊが正しいものでない場合が存在する。このような場合において、外部エネルギーＦｅのみを考慮するときには、或る制御点が、不正確な対応点に向けて移動するため、不正確なモデルフィッティングが行われることになる。これに対して、内部エネルギーＧｅを考慮することによれば、このような場合においても、標準モデルＳＢの形状を損なわないような変形動作が行われることになるため、より正確なモデリング処理が可能になる。
【０１２９】
この内部エネルギーＧｅは、標準モデルＳＢの各制御点Ｃｊ相互の関係が一定の関係（すなわち、全仮想バネが自然長を有する状態）に近づくときに最適化されるエネルギー項であるとも表現される。端的に言えば、この内部エネルギーＧｅは、標準モデルＳＢを一定の形状に維持しようとする役割を果たす。
【０１３０】
以上のようなモデルフィッティングを実行することによって、領域抽出装置１は背骨領域Ｂｒを高精度に抽出することが可能である。
【０１３１】
＜３−６．  対応点決定処理（ステップＳ５１）＞
ここで、ステップＳ５１の対応点決定処理について詳細に説明する。
【０１３２】
上述のように、或る制御点Ｃｊに対応する対応点Ｑｊを決定する際には、上記のステップＳ３で抽出された混在領域の境界上の点のうち、制御点Ｃｊに最も近い点を、制御点Ｃｊの対応点として求めることが可能である。
【０１３３】
ただし、上述のように混在領域には、モデリング対象物体の領域（背骨領域）とモデリング対象物体以外の領域（血管領域）との双方が包含されている。そのため、背骨に関する標準モデルＳＢの初期配置の状況次第では、モデリング対象領域ではない別の物体領域（例えば血管領域Ｖｒ１）の境界の方が、本来のモデリング対象領域（背骨領域Ｂｒ１）の境界よりも制御点Ｃｊの近傍に存在する場合がある。図２１は、血管領域Ｖｒ１の境界が、背骨領域Ｂｒ１の境界よりも、標準モデルＳＢの制御点Ｃｊの近くに存在する状況を示す概念図である。この場合、背骨領域Ｂｒ１と血管領域Ｖｒ１との双方を含む混在領域の境界上の点のうち制御点に最も近い境界点であることのみを基準（条件）にして、制御点Ｃｊの対応点Ｑｊを決定すると、モデリング対象領域ではない血管領域Ｖｒ１の境界上の点（例えば点Ｑｖ１）が誤って対応点として決定されてしまう。そして、このような誤った対応関係に基づいてモデルフィッティング処理が行われると、当該モデルフィッティング処理の結果に悪影響を与えることがある。
【０１３４】
そこで、この実施形態においては、計測データにおける複数の部分領域のうち最もモデリング対象物体らしい部分領域を選択し、選択された部分領域に属する点の中から対応点Ｑｊを選択するものとする。またここでは、上記と同様に、制御点Ｃｊに最も近い境界点であることをも、対応点選択の条件として採用する。
【０１３５】
具体的には、まず、計測データにおけるモデリング対象物体の領域と別物体の領域とを含む複数の領域（部分領域とも称する）のそれぞれについて、モデリング対象物体らしさを示す指標値ＤＸを算出する。詳細には、各スライス画像ＳＤにおける混在領域を連続領域ごとに区分し、それぞれの連続領域を部分領域として検出し、各部分領域の「対象物らしさ」（ここでは「背骨らしさ」）を示す指標値ＤＸ（後述）を算出する。
【０１３６】
次に、各スライス画像ＳＤにおける複数の部分領域のうち、指標値ＤＸに関する所定の基準（条件）を満たす部分領域（「対象物らしさ」が高い部分領域）を「基準充足領域」として求め、当該基準充足領域をモデリング対象物体の領域として特定する。これによって、例えば、対応点の探索対象から血管領域Ｖｒ１が除外され、背骨領域Ｂｒ１が対応点の探索対象として選択されることになる。
【０１３７】
また、このような各スライス画像でのモデリング対象領域を決定する動作を、複数のスライス画像について繰り返し行うことによって、複数のスライス画像のそれぞれにおける基準充足領域を特定する。これによって、各スライス画像においてモデリング対象領域（背骨領域Ｂｒ１）が特定される。
【０１３８】
その後、標準モデルＳＢ１の複数の制御点Ｃｊについて、それぞれ、対応点Ｑｊが選択されて決定される。詳細には、複数のスライス画像において背骨領域Ｂｒ１であると判定された領域の境界上の点のうち、さらに別の基準（条件）（ここでは、制御点Ｃｊに（３次元的に）最も近いこと）に基づいて選択された点（例えば点Ｑｂ１（図２１参照））が、或る制御点Ｃｊに対応する対応点Ｑｊとして決定される。なお、「制御点Ｃｊに３次元的に最も近い」という基準を用いる場合には、制御点Ｃｊとはその階層位置（Ｚ方向位置）が異なるスライス画像内の点が対応点Ｑｊとして決定されることもある。
【０１３９】
以上の処理によれば、モデリング対象物体の境界の近くに別物体の境界が存在する場合においても、当該別物体の境界上の点を対応点として誤検出することを回避し、本来のモデリング対象物体の境界上の点を対応点として求めることができる。したがって、対応点をより正確に求めて、モデルフィッティング処理をより正確に行うことが可能になる。
【０１４０】
ここで、指標値ＤＸおよび当該指標値ＤＸに関する条件について説明する。指標値および当該指標値ＤＸに関する条件としては、それぞれ、例えば次のような内容（（１）〜（４）参照）のものを採用することができる。
【０１４１】
（１）分散値ＤＸ１が閾値以上：
例えば、指標値ＤＸとして、各部分領域に属する画素の階調値のばらつき度合いを反映した値、具体的には各部分領域に属する画素のＣＴ値（階調値）の分散値ＤＸ１、を採用することができる。また、当該指標値ＤＸ１に関する条件としては、
条件ＣＮ１：分散値ＤＸ１が閾値以上であること、
を採用することができる。
【０１４２】
背骨領域Ｂｒ１のＣＴ値の標準偏差は（一般的に血管領域のＣＴ値の標準偏差よりも大きく）所定の閾値以上であることが多い。すなわち、背骨領域Ｂｒ１のＣＴ値は比較的ばらつき度合いが大きい。条件ＣＮ１は、このような特性を利用して、血管領域を排除して背骨領域Ｂｒ１を抽出することを意図している。
【０１４３】
このように、分散値ＤＸ１に関する条件ＣＮ１を利用することによって、複数の部分領域の中から本来のモデリング対象領域（背骨領域Ｂｒ１）を特定することができる。
【０１４４】
（２）最大値ＤＸ２が閾値以上：
或いは、指標値ＤＸとして、各部分領域に属する画素の最大階調値を反映した値、具体的には、各部分領域に属する画素のＣＴ値（階調値）の最大値ＤＸ２、を採用することができる。また、当該指標値ＤＸ２が満たすべき所定の条件としては、
条件ＣＮ２：最大値ＤＸ２が閾値以上であること、
を採用することができる。
【０１４５】
背骨領域Ｂｒ１の最大ＣＴ値（各スライス画像における背骨領域Ｂｒ１内に含まれる全画素のＣＴ値のうちの最大値）は、（一般的に血管領域の最大ＣＴ値よりも大きく）所定の閾値（例えばＴＨ９（図１０参照））以上であることが多い。条件ＣＮ２はこのような特性を利用して、血管領域を排除して背骨領域Ｂｒ１を抽出することを意図している。
【０１４６】
このように、最大値ＤＸ２に関する条件ＣＮ２を利用することによって、複数の部分領域の中から本来のモデリング対象領域（背骨領域Ｂｒ１）を特定することができる。
【０１４７】
（３）面積値ＤＸ３が閾値以上：
或いは、指標値ＤＸとして、各部分領域の面積を反映した値、例えば、スライス画像における各部分領域の面積値ＤＸ３自体、を採用することができる。また、当該指標値ＤＸ３が満たすべき条件としては、
条件ＣＮ３：面積値ＤＸ３が閾値以上であること；
を採用することができる。
【０１４８】
条件ＣＮ３は、各スライス画像における背骨領域Ｂｒ１の面積は血管領域の面積よりも大きく且つ所定の閾値よりも大きいこと、に基づく条件である。条件ＣＮ３はこのような特性を利用して、血管領域を排除して背骨領域Ｂｒ１を抽出することを意図している。
【０１４９】
このように、面積値ＤＸ３に関する条件ＣＮ３を利用することによって、複数の部分領域の中から本来のモデリング対象領域（背骨領域Ｂｒ１）を特定することができる。
【０１５０】
（４）円形度ＤＸ４が閾値以下：
或いは、指標値ＤＸとして、各スライス画像における各部分領域の円形度を反映した値、具体的には、円形度を表す式（４）の値ＤＸ４、を採用することができる。ただし、値Ｓａはその部分領域の面積を表し、値ｌはその部分領域の周囲長を表す。
【０１５１】
【数４】
2007312838.tif 000006
【０１５２】
また、当該指標値ＤＸ４が満たすべき条件としては、
条件ＣＮ４：円形度ＤＸ４が閾値以下であること、
を採用することができる。
【０１５３】
背骨領域Ｂｒ１は、各スライス画像において血管領域Ｖｒ１ほど円形に近くはない。円形度に関する条件ＣＮ４は、このような特性を利用して背骨領域Ｂｒ１を選択することを意図するものである。
【０１５４】
このように、円形度ＤＸ４に関する条件ＣＮ４を利用することによって、複数の部分領域の中から本来のモデリング対象領域（背骨領域Ｂｒ１）を特定することができる。
【０１５５】
以上（（１）〜（４））のようにして、各スライス画像における複数の部分領域のうち、その指標値ＤＸが所定の基準（条件）を満たす部分領域をモデリング対象物体の領域（モデリング対象領域）として特定することができる。
【０１５６】
なお、上記において、２つ以上の部分領域が上記条件（例えばＣＮ１）を満たす場合には、当該２つ以上の部分領域のうち、その指標値ＤＸが最も高い部分領域（最高指標値を有する部分領域）をモデリング対象領域として特定するようにしてもよい。
【０１５７】
また、上記においては、条件ＣＮ１〜ＣＮ４をそれぞれ独立に用いてモデリング対象領域を特定する場合を例示しているが、これに限定されず、各条件ＣＮ１〜ＣＮ４の全てあるいは一部を組み合わせて用いるようにしてもよい。例えば、条件ＣＮ１と条件ＣＮ２との双方を満たす部分領域を、モデリング対象領域として特定するようにしてもよい。あるいは、３つの条件ＣＮ１，ＣＮ３，ＣＮ４を満たす部分領域を、モデリング対象領域として特定するようにしてもよい。
【０１５８】
さらに、上記の条件ＣＮ１〜ＣＮ４に限定されず、その他の条件を用いるようにしてもよい。例えば、指標値ＤＸ１と指標値ＤＸ２との線形和を新たな指標値として算出し、当該新たな指標値が或る閾値以上であること、という基準を用いて、モデリング対象領域を特定するようにしてもよい。また、上記（１）では分散値ＤＸ１を例示したが、標準偏差を用いるようにしてもよい。さらには、全く異なる観点による指標値を用いるようにしてもよい。
【０１５９】
ところで、上記のステップＳ５１においては、必ずしも全てのスライス画像において背骨領域Ｂｒ１が適切に特定されるのではなく、一部のスライス画像においては背骨領域Ｂｒ１が適切に特定できないことがある。例えば、ステップＳ３の完了時点では、背骨領域Ｂｒ１と血管領域Ｖｒ１とが接しているスライス画像ＳＣにおいて、背骨と血管との両物体が混在する領域としての部分領域が抽出されているに過ぎず、背骨領域Ｂｒ１が血管領域Ｖｒ１から未だ分離されていない。そのため、このようなスライス画像ＳＣ等においては「モデリング対象領域」（背骨領域Ｂｒ１）が正しく特定されないことがある。
【０１６０】
しかしながら、図２２に示すように、背骨領域Ｂｒ１と血管領域Ｖｒ１とが接しているスライス画像（例えば画像ＳＣ）は、３次元画像ＴＤＰ１に含まれる複数のスライス画像のうちの少数であり、他の大部分のスライス画像ＳＥにおいては、背骨領域Ｂｒ１と血管領域Ｖｒ１とが分離された状態で撮影されている。そして、この実施形態では、大部分のスライス画像ＳＥにおいては、上述のような対応点検索によって、「モデリング対象領域」が適切に特定される。したがって、多くの制御点に関して、大部分のスライス画像における本来の「モデリング対象領域」の境界上の点が、各制御点に対応する対応点として適切に特定されるため、モデルフィッティング処理は大きな悪影響を受けずに済む。
【０１６１】
すなわち、仮に少数の制御点に対して不適切な点が対応点として選択されたとしても、他の多数の制御点に関しては適切な点が対応点として選択され適切な対応関係が求められるため、多数の適切な対応関係に基づく外部エネルギーの作用と、内部エネルギーの作用とによって、適切なモデリング動作が行われる。
【０１６２】
＜３−７．境界設定処理（ステップＳ６）＞
次に、ステップＳ６の境界設定処理について説明する。
【０１６３】
ステップＳ６においては、上記のステップＳ５で抽出された背骨領域Ｂｒ１の境界を、「拡張限界」として設定する。
【０１６４】
図２３は、ステップＳ５の終了時点におけるスライス画像ＳＣを示す図である。図２３においては、モデルフィッティング処理によって抽出された領域Ｂｒ１にクロスハッチングが付されて示されている。なお、ステップＳ５終了時点における抽出状態は、図８においても（立体的に）示されている。
【０１６５】
また、図２４は、図２３の一部を拡大して示す図である。図２４に示されるように、抽出された背骨領域Ｂｒの境界（境界面ないし境界線）ＢＤが、第２の領域拡張法（ステップＳ７）における「拡張限界」として設定される。具体的には、変形後の標準モデル（個別モデル）の計測データにおける存在位置が背骨領域Ｂｒであることから、変形後の標準モデル（個別モデル）ＳＢの境界上に存在する画素（ボクセル）ＶＸが、「拡張限界」を表すものとして設定される。
【０１６６】
＜３−８．第２の領域拡張処理（ステップＳ７）＞
ステップＳ７においては、再び領域拡張法を用いた領域抽出処理が実行される。ＳＰＳ７の領域拡張処理（第２の領域拡張処理）にも基本的にはステップＳ３の領域拡張処理（第１の領域拡張処理）と同様の手法が用いられる。ただし、この第２の領域拡張処理は、背骨領域Ｂｒに基づく「拡張限界」の制約の下で行われる点で、第１の領域拡張処理と相違する。以下では、ステップＳ３との相違点を中心に図１１を参照しながら説明する。
【０１６７】
具体的には、「拡張限界」が設定されているため、ステップＳ３２における処理が異なる。
【０１６８】
ステップＳ３２においては、或る画素が新たな拡張候補画素として検出されるための条件としては、注目画素にさらに隣接する画素であること、および領域内に存在するか否かが未だ判定されていない画素であることだけが求められるのではなく、「背骨領域Ｂｒ１の境界画素でないこと」も加重的に求められる。このように、第２の領域拡張法は、「拡張限界」の制約の下で行われる。
【０１６９】
この結果、ステップＳ７においては、血管領域Ｖｒ内の開始点ＳＰから拡張された血管領域Ｖｒのみが抽出される。ステップＳ７の領域拡張に際して、背骨領域Ｂｒの境界まで到達するとそれより先には進めない（すなわち背骨領域Ｂｒ内には進めない）。換言すれば、背骨領域Ｂｒの境界が「拡張限界」（あるいは拡張障壁）として機能して、当該境界を超えて領域を拡張することができない。したがって、背骨領域Ｂｒを抽出することなく、血管領域Ｖｒのみを抽出することが可能である。
【０１７０】
また、ステップＳ３４の判定処理においては、検出された拡張候補画素の階調値が所定の閾値ＴＨ２以上であるという条件を満たすか否かが判定される。
【０１７１】
閾値ＴＨ２が閾値ＴＨ１と同一の値であるときには、図２５に示すように、ステップＳ３で抽出された混合領域から背骨領域Ｂｒを除外した領域がステップＳ７において血管領域Ｖｒとして抽出されることになる。図２５は、ステップＳ７終了時点の（すなわち最終的な血管領域抽出後における）スライス画像ＳＣを示す図である。図２５においては、最終的に抽出された血管領域Ｖｒに斜線が付されて示されている。なお、ステップＳ７終了時点における抽出状態は、図９においても（立体的に）示されている。
【０１７２】
なお、ステップＳ３での領域拡張法に用いられる閾値ＴＨ１とステップＳ７での領域拡張法に用いられる閾値ＴＨ２とは、互いに異なる値に設定されてもよい。例えば、閾値ＴＨ１を閾値ＴＨ２よりも大きな値に設定して（すなわち比較的厳しい基準で）ステップＳ３の領域抽出を行うようにしてもよい。これによれば、ノイズをできるだけ排除した状態で背骨領域Ｂｒ１等を抽出することができ、ひいてはステップＳ５のモデルフィッティング処理において良好に背骨領域を抽出することが可能になる。また、これは、逆に言えば、閾値ＴＨ２を閾値ＴＨ１よりも小さな値（低い値）に設定して（すなわち比較的緩やかな基準で）ステップＳ７の領域抽出を行うことを意味する。したがって、図２６に示すように、「ステップＳ３で抽出された混合領域から背骨領域Ｂｒを除外した領域」よりも多くの領域（血管領域Ｖｒ１１等をも含む領域）、すなわち多くの血管領域Ｖｒを抽出することが可能になる。図２６においては、本来存在している血管領域Ｖｒ１１がさらに抽出されている。このように、細部の血管をも検出することが可能になる。
【０１７３】
＜４．変形例等＞
以上、この発明の実施の形態について説明したが、この発明は上記説明した内容のものに限定されるものではない。
【０１７４】
たとえば、上記実施形態においては、第１の領域拡張法を用いて、背骨領域Ｂｒを含む混在領域を抽出する場合を例示したが、これに限定されず、第１の領域拡張法とは別の手法によって、背骨領域Ｂｒを含む混在領域を抽出するようにしてもよい。例えば、各スライス画像ＳＤにおいて閾値ＴＨ１を用いた二値化処理を行って、背骨領域Ｂｒを含む混在領域を抽出するようにしてもよい。この場合においても、除外領域（背骨領域Ｂｒ）に関する標準モデルＳＢを用いて、計測データから除外領域Ｂｒをモデルフィッティング手法によって抽出することが可能である。したがって、除外領域Ｂｒの境界を拡張限界とする領域拡張処理（ステップＳ７）をさらに行うことによって、計測データから血管領域Ｖｒ（抽出対象領域）を抽出することが可能である。
【０１７５】
なお、上述の第１の領域拡張法を用いることによれば、領域の連続性を用いて混在領域を抽出するので、二値化処理を用いて混在領域を抽出する場合に比べて、ノイズ領域（誤検出領域）の量を低減することが可能である。
【０１７６】
また、上記実施形態においては、開始点ＳＰ１としては１つの点を指定する場合を例示したが、これに限定されず、複数の点を開始点として設定するようにしてもよい。ただし、情報処理量減少の観点からは、開始点の数を減少させることが好ましい。
【０１７７】
また、上記実施形態においては、背骨領域Ｂｒ２〜Ｂｒ４と血管領域Ｖｒとが接することはほとんど無いという事情に基づいて、簡易な標準モデルＳＢ１を用いたモデルフィッティング処理を行うことによって背骨領域Ｂｒ１を抽出する場合を例示した。しかしながら、これに限定されず、複雑な形状のより正確な標準モデルを用いたモデルフィッティング処理を行うことによって、背骨領域Ｂｒをより厳密に（例えば、背骨領域Ｂｒ１〜Ｂｒ４の全てを）抽出するようにしてもよい。
【０１７８】
また、上記実施形態の拡張候補画素の検出工程（ステップＳ３２）においては、拡張開始画素に隣接する全画素（２６画素）を拡張候補画素とし、領域拡張処理を実行していたが、これに限定されない。例えば、拡張開始画素Ｍ１に対して±ｘ方向に隣接する画素（Ｍ２、Ｍ３）、±ｙ方向に隣接する画素（Ｍ４、Ｍ７）及び±ｚ方向に隣接する画素（Ｕ１、Ｄ１）の計６画素を拡張候補画素として検出してもよい。
【０１７９】
また、上記実施形態においては、モデルフィッティング処理（ステップＳ５）の期間中において式（２）の定数αおよび式（３）の定数βをそれぞれ常に一定の値とする場合を例示したが、これに限定されない。例えば、ステップＳ５のモデルフィッティング処理を前半段階および後半段階の２つの段階に区分し、前半段階と後半段階とで定数α，βの値を変更するようにしてもよい。具体的には、前半段階においては、内部エネルギーＧｅを重視し、比較的大きな値β（あるいは比較的小さな値α）を用いて総合エネルギーＵｅを収束させるとともに、後半段階においては、外部エネルギーＦｅを重視し、比較的大きな値α（あるいは比較的小さな値β）を用いて総合エネルギーＵｅを収束させるようにしてもよい。これによれば、前半段階では標準モデルの形状をできるだけ維持した状態で標準モデルと計測データとの位置合わせ等を行うとともに、後半段階では計測データにあわせて標準モデルの形状を適切に変化させることが可能である。
【０１８０】
また、上記実施形態においては、図１６のフローチャートに示すようにして総合エネルギーＵｅを最小化（エネルギー関数を最適化）する場合を例示したが、エネルギー関数最適化手法は、上記手法に限定されない。例えば、準ニュートン法或いはニュートン法等の数学的手法を用いて総合エネルギーＵｅを最適化するようにしてもよい。
【０１８１】
また、上記実施形態においては、一部で接する複数の領域Ｂｒ１，Ｂｒ２，...，Ｖｒ１，Ｖｒ２，...から一の領域Ｂｒ１を抽出する場合について例示したが、本発明はこれに限定されない。例えば、比較的容易に互いに分離可能な複数の領域のうちの一の領域をモデルフィッティング手法を用いたモデリングによって抽出する場合に、本発明の思想を適用するようにしてもよい。その場合、上述のような領域拡張法を用いた処理（ステップＳ３，Ｓ７）は必要に応じて用いれば良い。
【図面の簡単な説明】
【０１８２】
【図１】Ｘ線ＣＴ装置によって取得される３次元画像を示す図である。
【図２】スライス画像を示す図である。
【図３】領域抽出装置の概要を示す図である。
【図４】領域抽出装置の各種機能を示すブロック図である。
【図５】領域抽出装置の全体動作を示すフローチャートである。
【図６】３次元画像内の状態（初期状態）を示す概念図である。
【図７】第１の領域拡張処理終了時点での状態を示す概念図である。
【図８】モデルフィッティング処理終了時点での状態を示す概念図である。
【図９】第２の領域拡張処理終了時点での状態を示す概念図である。
【図１０】全領域および血管領域の各ＣＴ値の度数分布を示す図である。
【図１１】第１の領域拡張処理の詳細を示すフローチャートである。
【図１２】隣接画素を示す図である。
【図１３】混在領域抽出後のスライス画像を示す図である。
【図１４】比較的簡略化された標準モデルを示す図である。
【図１５】突起部の最下端位置がスライス画像ごとに変化する様子を示す図である。
【図１６】モデルフィッティング処理を示すフローチャートである。
【図１７】標準モデルおよびその各制御点を示す図である。
【図１８】計測データ中に配置された標準モデルを示す図である。
【図１９】図１８の領域ＲＢ付近の拡大図である。
【図２０】制御点間を仮想バネで繋いだ模式図である。
【図２１】血管領域の境界が背骨領域の境界よりも制御点の近くに存在する状況を示す図である。
【図２２】多数のスライス画像においては、背骨領域と血管領域とが分離されている状況を示す図である。
【図２３】モデルフィッティング処理終了時点におけるスライス画像を示す図である。
【図２４】図２３の一部拡大図である。
【図２５】第２の領域拡張処理終了時点におけるスライス画像を示す図である。
【図２６】変形例に係るスライス画像を示す図である。
【符号の説明】
【０１８３】
１  領域抽出装置
ＢＮ  背骨
Ｂｒ，Ｂｒ１〜Ｂｒ４  背骨領域
Ｃｊ  制御点
Ｆｅ  外部エネルギー
Ｇｅ  内部エネルギー
ＰＤ，ＰＭ，ＰＵ，ＳＣ，ＳＤ  スライス画像
Ｑｊ  対応点
ＳＢ，ＳＢ１  標準モデル
ＴＤＰ１  ３次元画像
Ｕｅ  総合エネルギー
ＶＬ，Ｖｒ１，Ｖｒ２  血管領域
【図１】
2007312838.tif 000007
【図２】
2007312838.tif 000008
【図３】
2007312838.tif 000009
【図４】
2007312838.tif 000010
【図５】
2007312838.tif 000011
【図６】
2007312838.tif 000012
【図７】
2007312838.tif 000013
【図８】
2007312838.tif 000014
【図９】
2007312838.tif 000015
【図１０】
2007312838.tif 000016
【図１１】
2007312838.tif 000017
【図１２】
2007312838.tif 000018
【図１３】
2007312838.tif 000019
【図１４】
2007312838.tif 000020
【図１５】
2007312838.tif 000021
【図１６】
2007312838.tif 000022
【図１７】
2007312838.tif 000023
【図１８】
2007312838.tif 000024
【図１９】
2007312838.tif 000025
【図２０】
2007312838.tif 000026
【図２１】
2007312838.tif 000027
【図２２】
2007312838.tif 000028
【図２３】
2007312838.tif 000029
【図２４】
2007312838.tif 000030
【図２５】
2007312838.tif 000031
【図２６】
2007312838.tif 000032