JP5954846B2 - 形状データ生成プログラム、形状データ生成方法及び形状データ生成装置 - Google Patents

Description

形状データ生成プログラム、形状データ生成方法及び形状データ生成装置に関する。

医療分野においては、手術シミュレータ及び臓器シミュレータ等のシミュレータが、治療方針の決定、診断、術後予測、医薬品及び医療機器の開発等に利用されている。このようなシミュレータによるシミュレーションにおいては、臓器の三次元形状データを使用するが、臓器の三次元形状データを生成するのは容易ではないことが多い。臓器は生体内部に存在しているため目視や直接計測を行うことができない上に、臓器の形状は元来複雑であるためである。

従来、ＣＴ（Computed Tomography）画像及びＭＲＩ（Magnetic Resonance Imaging）画像等の断層画像の輝度値に基づき対象の領域を抽出することで、シミュレーションに利用する三次元形状データを生成する技術（例えばリジョングローイング法）が存在する。例えば、図１に示した心臓の断層画像のデータから、リジョングローイング法により、図２において実線２１で囲まれた領域（右室流体の領域）及び実線２２で囲まれた領域（左室流体の領域）を抽出し、右室流体及び左室流体の三次元形状データを生成したいとする。ところが、図１に示した心臓の断層画像データに対して実際にリジョングローイング法を実行すると、図３に示すような結果が得られる。図３に示す抽出結果においては、右室流体の一部の領域を、造影剤の影響により抽出できていない。また、左室流体において乳頭筋に該当する領域も抽出できていない。このように、輝度値のみに基づいた領域の抽出には限界がある。

また、輪郭を示す画素であるか否かをその特徴量に基づき評価する評価関数を求め、その評価関数を用いて、領域内の画素が輪郭を示す画素であるかを示す評価値を算出し、評価値に基づき輪郭を決定する技術が存在する。

また、誤差拡散法によって２値化された画像において、細かいエッジのつぶれを最小限にとどめ且つ階調性が高い鮮明な多値画像へ変換する技術が存在する。具体的には、入力された２値画像における注目画素の少なくともエッジ強度を検出し、検出結果に基づき画像の明るさの閾値を設定する。複数の多値化フィルタ内の黒画素及び白画素の割合と閾値とを比較することにより多値化フィルタの大きさを決定する。そして、決定された大きさの多値化フィルタ内の白画素及び黒画素の割合に基づき２値化画像データを多階調画像データに変換する。

しかし、上記のような技術を用いても、対象の領域を正確に抽出できないことがあるため、精度が高い三次元形状データを生成することはできない。

特開２００７−３０７３５８号公報 特開平８ー１９１３８７号公報

A. Tremeau, N. Borel, A Region Growing and Merging Algorithm to Color Segmentation, Pattern recognition, Volume 30, Issue 7, July 1997, Pages 1191-1203

従って、１つの側面では、本発明の目的は、高精度の三次元形状データを生成するための技術を提供することである。

本発明に係る形状データ生成方法は、複数の断層画像のうち対象物体が占める領域が指定されている断層画像を用いて、対象物体を表す三次元のボクセルデータを生成し、対象物体が占める領域が指定されている断層画像から、指定されている領域に含まれるボクセルの輝度値を抽出し、抽出された当該輝度値を用いて、対象物体が占める領域に含まれるボクセルである確率を算出するための関数を生成し、複数の断層画像を包含するボクセル空間における各ボクセルについて、当該ボクセルの輝度値と関数とを用いて確率を算出し、対象物体の形状を表す三次元のボクセルデータとボクセル空間における各ボクセルについて算出された確率とを用いて、対象物体の形状を表す三次元のボクセルデータを再生成する処理を含む。

高精度の三次元形状データを生成できるようになる。

図１は、心臓の断層画像を示す図である。 図２は、心臓の断層画像にシミュレーションに利用する領域を重ねて表示した図である。 図３は、心臓の断層画像に対してリジョングローイング法を実行した結果を示す図である。 図４は、形状データ生成装置の機能ブロック図である。 図５は、形状データ生成部の機能ブロック図である。 図６は、確率マップ生成部の機能ブロック図である。 図７は、ＤＩＣＯＭデータに含まれるヘッダ情報の一例を示す図である。 図８は、ＤＩＣＯＭデータに含まれる輝度値情報の一例を示す図である。 図９は、ピクセルの値の記述方法について説明するための図である。 図１０は、セグメント画像データ格納部に格納されているデータの一例を示す図である。 図１１は、メインの処理フローを示す図である。 図１２Ａは、初期値生成処理の概要について説明するための図である。 図１２Ｂは、初期値生成処理の処理フローを示す図である。 図１３は、形状データ生成処理の処理フローを示す図である。 図１４は、ＴＰＳＷａｒｐについて説明するための図である。 図１５Ａは、ボクセル化処理の処理フローを示す図である。 図１５Ｂは、ボクセル化処理の途中で生成されるデータの一例を示す図である。 図１５Ｃは、ボクセル化処理の途中で生成されるデータの一例を示す図である。 図１６Ａは、ボクセル化処理の処理フローを示す図である。 図１６Ｂは、ボクセル化処理の途中で生成されるデータの一例を示す図である。 図１７は、ボクセル化処理により生成されるデータを示す図である。 図１８は、境界特定処理の処理フローを示す図である。 図１９は、ぼかし処理により生成されるデータを示す図である。 図２０は、確率マップ生成処理の処理フローを示す図である。 図２１は、エロージョン処理について説明するための図である。 図２２は、輝度値格納部に格納されているデータの一例を示す図である。 図２３は、確率マップ生成処理の処理フローを示す図である。 図２４は、内側領域についてのヒストグラム及び外側領域についてのヒストグラムを示す図である。 図２５は、確率マップ生成処理の結果を特定のｚ座標について示した図である。 図２６は、初期値生成処理の結果を特定のｚ座標について示した図である。 図２７は、領域抽出処理の結果を特定のｚ座標について示した図である。 図２８は、表示装置に表示される画像の一例を示す図である。 図２９は、表示装置に表示される画像の一例を示す図である。 図３０は、表示装置に表示される画像の一例を示す図である。 図３１は、コンピュータの機能ブロック図である。

図４に、本実施の形態における形状データ生成装置の機能ブロック図を示す。形状データ生成装置１は、断層画像データ格納部１０と、入力処理部１１と、セグメント画像データ格納部１２と、形状データ生成部１３と、確率マップ生成部１４と、マージ部１５と、初期値格納部１６と、領域抽出部１７と、抽出結果格納部１８と、出力部１９と、標準形状データ格納部２０とを含む。

入力処理部１１は、断層画像データ格納部１０に格納されている断層画像データに対する、シミュレーションに用いる領域の指定をユーザから受け付ける。そして、入力処理部１１は、領域指定された断層画像データ（以下、セグメント画像データと呼ぶ）をセグメント画像データ格納部１２に格納する。形状データ生成部１３は、セグメント画像データ格納部１２に格納されているデータ及び標準形状データ格納部２０に格納されているデータを用いて、後述する形状データ生成処理を実行し、形状データを生成する。確率マップ生成部１４は、セグメント画像データ格納部１２に格納されているデータ及び断層画像データ格納部１０に格納されているデータを用いて、後述する確率マップ生成処理を実行し、確率マップを生成する。マージ部１５は、形状データ生成部１３により生成された形状データ及び確率マップ生成部１４により生成された確率マップを用いて初期値を生成し、初期値格納部１６に格納する。領域抽出部１７は、初期値格納部１６に格納されている初期値を用いて領域抽出処理を実行し、処理結果を抽出結果格納部１８に格納する。出力部１９は、抽出結果格納部１８に格納されているデータを表示装置（例えばディスプレイ）に表示させる。

図５に、形状データ生成部１３の機能ブロック図を示す。形状データ生成部１３は、変形処理部１３１と、第１形状データ格納部１３２と、ボクセル化処理部１３３と、ボクセル化データ格納部１３４と、境界特定部１３５と、特定結果格納部１３６と、ぼかし処理部１３７と、第２形状データ格納部１３８とを含む。

変形処理部１３１は、セグメント画像データ格納部１２に格納されているデータを用いて、標準形状データ格納部２０に格納されている標準形状のデータ（本実施の形態においては、左室流体（乳頭筋部分を含む）の標準形状のデータ）を変形し、変形後の形状データを第１形状データ格納部１３２に格納する。ボクセル化処理部１３３は、第１形状データ格納部１３２に格納されている形状データに対して、後述するボクセル化処理を実行し、処理結果をボクセル化データ格納部１３４に格納する。境界特定部１３５は、ボクセル化データ格納部１３４に格納されているデータを用いて左室流体と左室流体ではない部分との境界に該当するボクセルを特定する処理を実行し、処理結果を特定結果格納部１３６に格納する。ぼかし処理部１３７は、ボクセル化データ格納部１３４に格納されているデータ及び特定結果格納部１３６に格納されているデータを用いてぼかし処理を実行し、処理結果を第２形状データ格納部１３８に格納する。

図６に、確率マップ生成部１４の機能ブロック図を示す。確率マップ生成部１４は、正規化部１４１と、第１データ格納部１４２と、モルフォロジ処理部１４３、第２データ格納部１４４、輝度値抽出部１４５及び輝度値格納部１４９を含むヒストグラム生成部１５０と、ヒストグラムデータ格納部１４６と、マッピング部１４７と、確率マップ格納部１４８とを含む。

正規化部１４１は、断層画像データ格納部１０に格納されている断層画像データの輝度値を正規化する処理を実行し、処理結果を第１データ格納部１４２に格納する。モルフォロジ処理部１４３は、セグメント画像データ格納部１２に格納されているセグメント画像データのうち、領域指定された部分に対してエロージョン処理を実行し、処理結果を第２データ格納部１４４に格納する。輝度値抽出部１４５は、第１データ格納部１４２に格納されているデータ、第２データ格納部１４４に格納されているデータを用いて処理を行い、処理結果を輝度値格納部１４９に格納する。また、輝度値抽出部１４５は、輝度値格納部１４９に格納されているデータを用いてヒストグラムを生成し、ヒストグラムデータ格納部１４６に格納する。マッピング部１４７は、ヒストグラムデータ格納部１４６に格納されているデータを用いて、断層画像データ格納部１０に格納されている断層画像データを包含するボクセル空間における各ボクセルに確率データをマッピングする処理を実行し、処理結果を確率マップ格納部１４８に格納する。

断層画像データ格納部１０には、例えば図１に示すような断層画像データが複数枚分格納されている。断層画像データのフォーマットは、例えばＤＩＣＯＭ（Digital Imaging and COmmunication in Medicine）である。ＤＩＣＯＭのデータは、ヘッダ情報と輝度値の情報とを含む。ヘッダ情報には、例えば図７に示すように、番号毎に患者情報、医療用測定装置の種類及び断層画像の解像度等の情報が含まれる。図７において、「０．４５６」はピクセル１つ分の距離（ミリ）を表しており、「０．５」は断層画像間の距離（ミリ）を表している。また、輝度値情報には、例えば図８に示すように、Ｒ（赤）、Ｇ（緑）及び青（Ｂ）の各々についての値が含まれる。ＣＴの場合には、ＣＴ値（単位はHounsfield Unit）と呼ばれる値が含まれ、最小値は−１０００であり、最大値は３０９５である。図８における各ピクセルの値は、図９に示すように、断層画像の左上隅から右方向に記述され、右端に到達すると、次の行の左端から右方向へ同様に記述される。なお、ＣＴ及びＭＲＩ等によって取得された複数枚の断層画像を積層すると、三次元の臓器の形状を包含するデータを生成することができる。このデータは、ボリュームデータと呼ばれる。

図１０に、セグメント画像データ格納部１２に格納されているデータの一例を示す。図１０の例は、図１に示した断層画像データに、シミュレーションに用いる部分であるとしてユーザによって指定された領域（本実施の形態においては、左室流体（乳頭筋部分を含む）の部分）に所定のラベル値が付与された画像データである。セグメント画像データ格納部１２には、このようなセグメント画像データが複数枚分格納されている。ユーザは、断層画像データ格納部１０に格納されている断層画像データの一部（例えば、３０枚に１枚）に、手動で領域指定を行っていく。なお、領域指定はセグメンテーションと呼ばれることもある。

標準形状データ格納部２０には、左室流体の標準的な形状（以下、標準形状と呼ぶ）のデータが格納されている。具体的には、形状の頂点の情報及び頂点の接続情報を含むＳＴＬ（Standard Triangulated Language）データが格納されており、形状は三角形メッシュからなる。但し、データのフォーマットはＳＴＬデータのフォーマットに限られるわけではない。

次に、図１１乃至図３０を用いて、図４に示した形状データ生成装置１の動作について説明する。

まず、入力処理部１１は、断層画像データ格納部１０に格納されている断層画像データの一部（例えば、３０枚に１枚）を表示装置に順に出力する（図１１：ステップＳ１）。これにより、医師等であるユーザに領域指定を促す。ユーザは、各断層画像データにおける目的の領域（本実施の形態においては、左室流体である領域）に対し、マウス等を操作してラベル値を付与する。これに応じ、入力処理部１１は、領域指定を受け付け（ステップＳ３）、ユーザによる領域指定により生成されたセグメント画像データをセグメント画像データ格納部１２に格納する。

そして、形状データ生成装置１は、初期値生成処理を実行する（ステップＳ５）。初期値生成処理については、図１２Ａ乃至図２６を用いて説明する。

まず、図１２Ａを用いて、初期値生成処理の概要について説明する。初期値生成処理により生成される初期値は、ボクセル空間上のデータであり、後述する反応拡散フィルタの初期値になる。初期値は、形状データと確率マップとをマージすることにより生成される。形状データは、標準形状データをセグメント画像データにより特定される形状に合うように変形することにより得られる。確率マップは、セグメント画像データを用いて作成された輝度値のヒストグラムを用いて、確率データを断層画像データにおける各ボクセルにマッピングすることにより得られる。

以下では、初期値生成処理の処理フローについて説明する。

まず、形状データ生成部１３は、形状データ生成処理を実行する（図１２Ｂ：ステップＳ１１）。形状データ生成処理については、図１３乃至図１９を用いて説明する。

変形処理部１３１は、標準形状データ格納部２０から標準形状のデータを読み出すと共に、セグメント画像データ格納部１２からセグメント画像データを読み出す。そして、変形処理部１３１は、標準形状のデータ及びセグメント画像データを含むランドマーク設定画面を表示装置に表示させる（図１３：ステップＳ２１）。

ユーザは、表示装置に表示されたランドマーク設定画面を見て、標準形状とセグメント画像データにより特定される形状（以下、目標形状と呼ぶ）との大まかな位置合わせを行う。具体的には、（１）標準形状における所定の部分にソースランドマークを設定する。（２）目標形状における、ソースランドマークを配置した位置に対応する部分にターゲットランドマークを設定する。なお、所定の部分とは、例えば左室流体における特徴的な部分である。

変形処理部１３１は、ソースランドマーク及びターゲットランドマークの設定を受け付け、ソースランドマーク及びターゲットランドマークのデータ（例えば三次元座標）をメインメモリ等の記憶装置に格納する（ステップＳ２３）。

変形処理部１３１は、後述するＴＰＳＷａｒｐ（Thin Plate Spline Warp）等の手法により標準形状を変形する処理を実行する（ステップＳ２５）。そして、変形処理部１３１は、処理結果の形状のデータを第１形状データ格納部１３２に格納する。なお、処理結果の形状データのフォーマットは、ＳＴＬデータのフォーマットである。

ここで、ＴＰＳＷａｒｐについて説明する。図１４に示すように、ＴＰＳＷａｒｐにおいては、ソースランドマーク及び当該ソースランドマークに対応するターゲットランドマークが与えられると、ソースランドマークがターゲットランドマークに重なるように変形を行う。なお、ＴＰＳＷａｒｐの詳細については、例えば「Fred L. Bookstein, "Principal Warps: Thin-Plate Splines and the Decomposition of Deformations",IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. 11, NO. 6, PP. 567-585, JUNE 1989」を参照のこと。

このようにすれば、標準形状を、ユーザによって指定された領域に合うように変形できるようになる。

そして、ボクセル化処理部１３３は、第１形状データ格納部１３２に格納されている形状データを用いてボクセル化処理を実行する（ステップＳ２７）。ボクセル化処理については、図１５Ａ乃至図１７を用いて説明する。

はじめに、ボクセル化処理を説明するにあたって前提となる事項を説明する。

をボクセル空間における座標（ｘ，ｙ，ｚ）におけるボクセル値とする。初期値は−１とする。

をボクセル空間におけるボクセルとする。ｖＮは、ボクセル空間に含まれるボクセルの数である。ボクセル空間において、ｘ軸方向におけるボクセルの数をｘＮ、ｙ軸方向におけるボクセルの数をｙＮ、ｚ軸方向におけるボクセルの数をｚＮとする。従って、ｘＮ×ｙＮ×ｚＮ＝ｖＮである。

をボクセル空間の下限（左下）の座標とする。

をボクセル空間の上限（右上）の座標とする。

をボクセルのサイズ（ｘｓ×ｙｓ×ｚｓ）とする。

第１形状データ格納部１３２に格納されている形状データをＳとし、形状データＳを以下のように定める。

ここで、｛ｐ_ｉ｝_{ｉ＝１，．．．，Ｎ}は形状データＳにおける頂点を表し、ベクトル｛ｍ_ｊ｝_{ｊ＝１，．．．，Ｍ}は形状データＳにおける三角形メッシュを表す。Ｎは形状データＳにおける頂点の数である。Ｍは形状データＳにおける三角形メッシュの数である。

メッシュｍ_ｊに含まれる頂点を以下のように定める。

なお、ボクセル化処理によって生成されるボクセル空間は、断層画像データ格納部１０に格納されている断層画像データによって特定されるボクセル空間と一致するボクセル空間である。すなわち、下限及び上限の座標、各軸方向のボクセル数並びにボクセルサイズが同じであるとする。

以下では、ボクセル化処理の処理フローについて説明する。

まず、ボクセル化処理部１３３は、三角形メッシュを識別するための変数であるｊをｊ＝１と設定する（図１５Ａ：ステップＳ３１）。

ボクセル化処理部１３３は、ボクセル空間におけるｘ座標を表す変数であるｒをｒ＝ｒ_ｍｉｎと設定する（ステップＳ３２）。ｒ_ｍｉｎは、以下のような値である。

ここで、ガウス記号中の分子の第１項は、ｘ_{ｍ（ｊ），１}（（ｊ）はｊがｍの下付きの文字であることを表す）の最小値を求める関数である。

ボクセル化処理部１３３は、ボクセル空間におけるｙ座標を表す変数であるｓをｓ＝ｓ_ｍｉｎと設定する（ステップＳ３３）。ｓ_ｍｉｎは、以下のような値である。

ボクセル化処理部１３３は、ボクセル空間におけるｚ座標を表す変数であるｔをｔ＝ｔ_ｍｉｎと設定する（ステップＳ３４）。ｔ_ｍｉｎは、以下のような値である。

ボクセル化処理部１３３は、（ｒ，ｓ，ｔ）が形状データＳの内部に有るか判断する（ステップＳ３５）。ステップＳ３５においては、以下の式を満たすか判断する。

ここで、
とする。

数式１１において、（（ｖｓ・ｖ＋Ｌ）−ｐ_{ｍ（ｊ），１}）（（ｊ）はｊがｍの下付きの文字であることを表す）は、始点をｐ_{ｍ（ｊ），１}とし終点をボクセル（ｒ，ｓ，ｔ）とするベクトルを表す。（（ｐ_{ｍ（ｊ），２}−ｐ_{ｍ（ｊ），１}）×（ｐ_{ｍ（ｊ），３}−ｐ_{ｍ（ｊ），１}））は、三角形メッシュを含む平面を形状の内側から外側の方向へ貫く法線ベクトルを表す。前者のベクトルと法線ベクトルとの内積が０より小さいか否かを判断することにより、ボクセル（ｒ，ｓ，ｔ）が形状データＳの内部に有るか否かを判断することができる。但し、ステップＳ３５の処理における判断の対象となるボクセルは、左室流体と左室流体でない部分との境界付近のボクセルに限られている。

（ｒ，ｓ，ｔ）が形状データＳの内部に有る場合（ステップＳ３５：Ｙｅｓルート）、ボクセル化処理部１３３は、Ｖ（ｒ，ｓ，ｔ）＝１と設定する（ステップＳ３６）。一方、（ｒ，ｓ，ｔ）が形状データＳの内部に無い場合（ステップＳ３５：Ｎｏルート）、ボクセル化処理部１３３は、Ｖ（ｒ，ｓ，ｔ）＝０と設定する（ステップＳ３７）。

ボクセル化処理部１３３は、ｔ≧ｔ_ｍａｘであるか判断する（ステップＳ３８）。ｔ_ｍａｘは、以下のような値である。

ここで、ガウス記号中の分子の第１項は、ｚ_{ｍ（ｊ），１}（（ｊ）はｊがｍの下付きの文字であることを表す）の最大値を求める関数である。

ｔ≧ｔ_ｍａｘではない場合（ステップＳ３８：Ｎｏルート）、ボクセル化処理部１３３は、ｔを１インクリメントし（ステップＳ３９）、ステップＳ３５の処理に戻る。一方、ｔ≧ｔ_ｍａｘである場合（ステップＳ３８：Ｙｅｓルート）、ボクセル化処理部１３３は、ｓ≧ｓ_ｍａｘであるか判断する（ステップＳ４０）。ｓ_ｍａｘは、以下のような値である。

ｓ≧ｓ_ｍａｘではない場合（ステップＳ４０：Ｎｏルート）、ボクセル化処理部１３３は、ｓを１インクリメントし（ステップＳ４１）、ステップＳ３４の処理に戻る。一方、ｓ≧ｓ_ｍａｘである場合（ステップＳ４０：Ｙｅｓルート）、ボクセル化処理部１３３は、ｒ≧ｒ_ｍａｘであるか判断する（ステップＳ４２）。ｒ_ｍａｘは、以下のような値である。

ｒ≧ｒ_ｍａｘではない場合（ステップＳ４２：Ｎｏルート）、ボクセル化処理部１３３は、ｒを１インクリメントし（ステップＳ４３）、ステップＳ３３の処理に戻る。一方、ｒ≧ｒ_ｍａｘである場合（ステップＳ４２：Ｙｅｓルート）、ボクセル化処理部１３３は、ｊ≧Ｍであるか判断する（ステップＳ４４）。ｊ≧Ｍではない場合（ステップＳ４４：Ｎｏルート）、ボクセル化処理部１３３は、ｊを１インクリメントし（ステップＳ４５）、ステップＳ３２の処理に戻る。一方、ｊ≧Ｍである場合（ステップＳ４４：Ｙｅｓルート）、処理は端子Ａを介して図１６ＡのステップＳ５１に移行する。

ここまでの処理によって生成されるボクセルデータの一例を図１５Ｂに示す。図１５Ｂは、特定のｚ座標におけるボクセル値を示している。図１５Ｂにおいては、ボクセル値が１であるボクセルを白で示し、ボクセル値が０であるボクセルを灰色で示し、ボクセル値が−１であるボクセルを黒で示している。ここまでの処理では、左室流体と左室流体ではない部分との境界付近のボクセル以外のボクセルについては、初期値である−１が設定されたままである。

図１５Ｃは、図１５Ｂにおいて白線１５５で囲まれた領域を拡大した図である。図１５Ｃに示すように、境界付近のボクセルのうち形状データＳの内部にあるボクセルのボクセル値は１であり、境界付近のボクセルのうち形状データＳの外部にあるボクセルのボクセル値は０である。

図１６Ａの説明に移行し、ボクセル化処理部１３３は、ボクセル空間に対して、シード点を（０，０，０）としてリジョングローイング法を実行し（図１６Ａ：ステップＳ５１）、処理結果をメインメモリ等の記憶装置に格納する。リジョングローイング法は、シード点に連結する点を取り込むことで領域を拡張し、対象領域を抽出する方法である。リジョングローイング法については、例えば「A. Tremeau, N. Borel, A Region Growing and Merging Algorithm to Color Segmentation, Pattern recognition, Volume 30, Issue 7, July 1997, Pages 1191-1203」を参照のこと。なお、ステップＳ５１においてはシード点を（０，０，０）としているが、左室流体ではない部分（すなわち、形状データＳの外側の部分）のボクセルであればシード点として利用することができる。

ボクセル化処理部１３３は、ｘ座標を識別するための変数ｘをｘ＝１と設定する（ステップＳ５２）。

ボクセル化処理部１３３は、ｙ座標を識別するための変数ｙをｙ＝１と設定する（ステップＳ５３）。

ボクセル化処理部１３３は、ｚ座標を識別するための変数ｚをｚ＝１と設定する（ステップＳ５４）。

ボクセル化処理部１３３は、ボクセル（ｘ，ｙ，ｚ）がリジョングローイング法によって抽出されたか判断する（ステップＳ５５）。抽出された場合（ステップＳ５５：Ｙｅｓルート）、ボクセル（ｘ，ｙ，ｚ）は形状データＳの外側に存在するので、ボクセル化処理部１３３は、Ｖ（ｘ，ｙ，ｚ）＝０と設定する（ステップＳ５６）。抽出されていない場合（ステップＳ５５：Ｎｏルート）、ボクセル化処理部１３３は、Ｖ（ｘ，ｙ，ｚ）が−１であるか判断する（ステップＳ５７）。

Ｖ（ｘ，ｙ，ｚ）が−１である場合（ステップＳ５７：Ｙｅｓルート）、ボクセル（ｘ，ｙ，ｚ）は形状データＳの内側にあり且つ値が初期値のままなので、ボクセル化処理部１３３は、Ｖ（ｘ，ｙ，ｚ）＝１と設定する（ステップＳ５８）。一方、Ｖ（ｘ，ｙ，ｚ）が−１ではない場合（ステップＳ５７：Ｎｏルート）、Ｖ（ｘ，ｙ，ｚ）＝０又は１と既に設定されているので、ステップＳ５９の処理に移行する。

ボクセル化処理部１３３は、ｚ≧ｚＮであるか判断する（ステップＳ５９）。ｚ≧ｚＮではない場合（ステップＳ５９：Ｎｏルート）、ボクセル化処理部１３３は、ｚを１インクリメントし（ステップＳ６０）、ステップＳ５５の処理に戻る。

ｚ≧ｚＮである場合（ステップＳ５９：Ｙｅｓルート）、ボクセル化処理部１３３は、ｙ≧ｙＮであるか判断する（ステップＳ６１）。ｙ≧ｙＮではない場合（ステップＳ６１：Ｎｏルート）、ボクセル化処理部１３３は、ｙを１インクリメントし（ステップＳ６２）、ステップＳ５４の処理に戻る。

ｙ≧ｙＮである場合（ステップＳ６１：Ｙｅｓルート）、ボクセル化処理部１３３は、ｘ≧ｘＮであるか判断する（ステップＳ６３）。ｘ≧ｘＮではない場合（ステップＳ６３：Ｎｏルート）、ボクセル化処理部１３３は、ｘを１インクリメントし（ステップＳ６４）、ステップＳ５３の処理に戻る。一方、ｘ≧ｘＮである場合（ステップＳ６３：Ｙｅｓルート）、元の処理に戻る。

図１６Ｂを用いて、ステップＳ５６の処理について説明する。本実施の形態においては、シード点を（０，０，０）としている（すなわち、形状データＳの外側の部分としている）ため、リジョングローイング法においては、形状データＳの外側の部分が抽出される。よって、ステップＳ５６において、リジョングローイング法によって抽出されたボクセルの値を０とすると、図１６Ｂのようになる。図１６Ｂは、図１５Ｂと同様、特定のｚ座標におけるボクセル値を示しており、ボクセル値が１であるボクセルを白で示し、ボクセル値が０であるボクセルを灰色で示し、ボクセル値が−１であるボクセルを黒で示している。但し、ステップＳ３７の処理においてＶ（ｘ，ｙ，ｚ）＝０と設定されたボクセルについては、灰色の上に黒の縞を付している。そして、図１６Ｂにおいて形状データＳの内側の黒色部分に該当するボクセルの値は、ステップＳ５８の処理により１に置換される。従って、形状データＳの内側のボクセル値が１となり、形状データＳの外側のボクセル値が０となる。

以上のようなボクセル化処理を実行することにより、ＳＴＬデータであった形状データを、ボクセル空間上のデータとする。ボクセル化処理を実行すると、例えば図１７に示すようなデータを生成することができる。図１７においては、ボクセル値が１である部分を明るい色（白）で示しており、ボクセル値が０である部分を黒で示している。この段階では、左室流体と左室流体ではない部分との境界のボクセル値は０又は１であるため、境界は滑らかではなく凸凹がある。

図１３の説明に戻り、ボクセル化処理（ステップＳ２７）が終了すると、境界特定部１３５は、境界特定処理を実行する（ステップＳ２８）。境界特定処理については、図１８を用いて説明する。

まず、境界特定部１３５は、ボクセルを識別するための変数ｉをｉ＝１と設定する（図１８：ステップＳ７１）。

境界特定部１３５は、ｖ_ｉは左室流体と左室流体ではない部分との境界に該当するか判断する（ステップＳ７３）。ステップＳ７３においては、以下の式を満たすか判断する。

ｖ_ｉが境界に該当する場合（ステップＳ７３：Ｙｅｓルート）、境界特定部１３５は、ｖ_ｉ＝（ｘ，ｙ，ｚ）を境界ｂ＝｛ｖ_ｂ｝に追加し（ステップＳ７５）、ｖ_ｉ＝（ｘ，ｙ，ｚ）の座標を特定結果格納部１３６に格納する。ｖ_ｉが境界に該当しない場合（ステップＳ７３：Ｎｏルート）、ステップＳ７７の処理に移行する。

境界特定部１３５は、ｉ≧ｖＮであるか判断する（ステップＳ７７）。ｉ≧ｖＮではない場合（ステップＳ７７：Ｎｏルート）、境界特定部１３５は、ｉを１インクリメントし（ステップＳ７９）、ステップＳ７３の処理に戻る。一方、ｉ≧ｖＮである場合（ステップＳ７７：Ｙｅｓルート）、元の処理に戻る。

このように、境界に該当するボクセルを特定しておけば、ぼかし処理に利用できるようになる。

図１３の説明に戻り、ぼかし処理部１３７は、ボクセル化データ格納部１３４に格納されているデータにより特定される形状の境界に対してぼかし処理を実行し、処理結果を第２形状データ格納部１３８に格納する（ステップＳ２９）。そして元の処理に戻る。

ステップＳ２９においては、特定結果格納部１３６にデータが格納されているボクセルのボクセル値の平均が０．５（すなわち、０と１との平均）になるように、ガウシアンフィルタを実行する。ガウシアンフィルタを実行すると、ボクセル値の分布は図１９に示すようになる。図１９において、ボクセル値は０から１までの間の値を取るようになっており、ボクセル値が大きいほど明るい色で示している。境界においては、ボクセル値は１から０までに徐々に変化するようになっている。

図１２Ｂの説明に戻り、形状データ生成処理が終了すると、確率マップ生成部１４は、確率マップ生成処理を実行する（ステップＳ１３）。確率マップ生成処理については、図２０乃至図２５を用いて説明する。

まず、正規化部１４１は、断層画像データ格納部１０に格納されている断層画像データの輝度値を正規化する処理を実行し、処理結果を第１データ格納部１４２に格納する（図２０：ステップＳ８１）。ステップＳ８１においては、例えば、輝度値を−１０００から１０００の間の値に正規化する。但し、このような値の範囲に限られるわけではない。

ヒストグラム生成部１５０におけるモルフォロジ処理部１４３は、セグメント画像データ格納部１２から未処理のセグメント画像データを１枚分特定する（ステップＳ８３）。

モルフォロジ処理部１４３は、特定されたセグメント画像データに対して、モルフォロジ処理の一種であるエロージョン処理を実行し、処理結果を第２データ格納部１４４に格納する（ステップＳ８５）。なお、モルフォロジ処理については、例えば「L. Vincent, Morphological grayscale reconstruction in image analysis: Applications and efficient algorithms, IEEE Transactions on Image Processing, 2(2):176-201, 1993.」を参照のこと。

図２１を用いて、ステップＳ８５におけるエロージョン処理について説明する。図２１において、曲線２０２で囲まれた領域（以下、内側領域と呼ぶ）はユーザにより指定された領域である。ステップＳ８５におけるエロージョン処理においては、ユーザにより指定された領域を、曲線２０１で囲まれた領域によりマスクする。以下では、曲線２０１と曲線２０２との間の領域を外側領域と呼ぶ。

輝度値抽出部１４５は、ステップＳ８３において特定されたセグメント画像データの内側領域に対応する領域を、第１データ格納部１４２に格納されている当該セグメント画像データに対応する断層画像データから特定する（ステップＳ８７）。

輝度値抽出部１４５は、ステップＳ８７において特定された領域に含まれるボクセルの輝度値を抽出し、輝度値格納部１４９に格納する（ステップＳ８９）。

図２２に、輝度値格納部１４９に格納されているデータの一例を示す。図２２の例では、内側領域及び外側領域について、断層画像データにおけるボクセルの座標と、そのボクセルの輝度値とが格納されている。

輝度値抽出部１４５は、ステップＳ８３において特定されたセグメント画像データの外側領域に対応する領域を、第１データ格納部１４２に格納されている当該セグメント画像データに対応する断層画像データから特定する（ステップＳ９１）。

輝度値抽出部１４５は、ステップＳ９１において特定された領域に含まれるボクセルの輝度値を抽出し、輝度値格納部１４９に格納する（ステップＳ９３）。

輝度値抽出部１４５は、未処理のセグメント画像データがセグメント画像データ格納部１２に有るか判断する（ステップＳ９５）。未処理のセグメント画像データが有る場合（ステップＳ９５：Ｙｅｓルート）、次のセグメント画像データについて処理するため、ステップＳ８３の処理に戻る。一方、未処理のセグメント画像データが無い場合（ステップＳ９５：Ｎｏルート）、処理は端子Ｂを介して図２３のステップＳ９７に移行する。

図２３の説明に移行し、輝度値抽出部１４５は、輝度値格納部１４９に格納されている輝度値を用いて、内側領域についてのヒストグラムＨ_ｉｎ（Ｉ_{ｉ，ｊ，ｋ}）及び外側領域についてのヒストグラムＨ_ｏｕｔ（Ｉ_{ｉ，ｊ，ｋ}）を生成する（ステップＳ９７）。そして、輝度値抽出部１４５は、生成されたヒストグラムをヒストグラムデータ格納部１４６に格納する。

図２４に、ステップＳ９７において生成されるヒストグラムを示す。図２４の縦軸は頻度であり、横軸は輝度値である。内側領域についてのヒストグラムＨ_ｉｎ（Ｉ_{ｉ，ｊ，ｋ}）は、輝度値が３００付近にピークがある。一方、外側領域についてのヒストグラムＨ_ｏｕｔ（Ｉ_{ｉ，ｊ，ｋ}）は、輝度値が７０付近にピークがある。

マッピング部１４７は、ヒストグラムデータ格納部１４６に格納されているヒストグラムを用いて関数を生成する（ステップＳ９８）。また、マッピング部１４７は、断層画像データにおける各ボクセルに対し、生成された関数を用いて確率データをマッピングする（ステップＳ９９）。そして、マッピング部１４７は、マッピングの結果を確率マップ格納部１４８に格納する。そして元の処理に戻る。ステップＳ９９においては、以下の関数を用いてマッピングをする。

数式１７によって求められる確率は、断層画像データにおけるボクセル（ｉ，ｊ，ｋ）が左室流体に含まれる確率を表している。確率マップ生成処理によって生成されたデータの一例を図２５に示す。図２５は、特定のｚ座標における確率データを示している。図２５においては、確率が１に近いほど明るい色になり、確率が０に近いほど暗い色になっている。なお、本実施の形態においては左室流体についての確率マップを生成しているため、図２５において左室流体以外の領域についての確率は意味を持たない。

図１２Ｂの説明に戻り、確率マップ生成処理が終了すると、マージ部１５は、第２形状データ格納部１３８に格納されているデータ及び確率マップ格納部１４８に格納されている確率データをマージすることにより初期値を生成する（ステップＳ１５）。そして、マージ部１５は、生成された初期値を初期値格納部１６に格納し、元の処理に戻る。ステップＳ１５においては、以下の式に従って初期値を生成する。

ここで、λは０＜λ＜１であり、例えば０．５とする。Ｄ１_{ｉ，ｊ，ｋ}は、第２形状データ格納部１３８に格納されている、ボクセル（ｉ，ｊ，ｋ）におけるボクセル値である。Ｄ２_{ｉ，ｊ，ｋ}は、確率マップ格納部１４８に格納されている、ボクセル（ｉ，ｊ，ｋ）の確率データである。ｅｐｓは十分に小さな値であり、例えば０．０１とする。

図２６に、ステップＳ１５において生成された初期値の一例を示す。図２６は、特定のｚ座標における初期値を示している。図２６においては、値が１に近いほど明るい色になり、値が０に近いほど暗い色になっている。また、数式１８においてｅｐｓの条件を利用しているため、左室流体以外であるとみなされる領域は除去されている（すなわち、値が０になっている）。

図１１の説明に戻り、領域抽出部１７は、初期値格納部１６に格納されている初期値を用いて領域抽出処理を実行し（ステップＳ７）、処理結果を抽出結果格納部１８に格納する。

ステップＳ７においては、以下の反応拡散方程式を利用して、初期値格納部１６に格納されている初期値に対してフィルタをかけることにより、左室流体の領域を抽出している。

ここで、ｕはボクセル値である。αは反応係数であり、βは拡散係数である。

第１項は反応項であり、以下の関数を用いる。

反応項は、ｕを代表値に集める効果をもたらす。上記の反応拡散方程式の場合、ｕが０から０．５までの値である場合には、ｕは次第に０に近付く。ｕが０．５から１までの値である場合には、ｕは次第に１に近付く。

第２項は拡散項であり、ガウシアンフィルタと同様の効果（すなわち、平滑化の効果）をもたらす。

数式１９を離散化すると、以下のようになる。

ここで、ｕ_{ｉ，ｊ，ｋ}は、ボクセル（ｉ，ｊ，ｋ）におけるボクセル値である。Ｄは以下のようになる。

また、外縁のボクセルの値ｕ_０は、以下のようにする。

セグメント画像データにおいて指定された領域に含まれるボクセルの値ｕ_１は、ノイマン境界条件として、以下のようにする。

セグメント画像データにおいて指定された領域に含まれないボクセルの値ｕ_２は、ノイマン境界条件として、以下のようにする。

心室と心房の境界に該当するボクセルの値ｕ_３は、ディリクレ境界条件として、以下のようにする。

以上のような条件のもとで領域抽出処理を実行する（すなわち、所定回数計算を実行する）と、例えば図２７に示すようなデータが得られる。図２７は、特定のｚ座標におけるボクセル値を示している。領域抽出処理により、ボクセル値は０又は１に近付くようになっている。そのため、図２６に示した初期値と比較すると、境界付近に存在していたぼかされた部分は少なくなっている。

図１１の説明に戻り、出力部１９は、抽出結果格納部１８に格納されているデータを用いて、抽出した領域を示すデータを生成し、表示装置に表示させる（ステップＳ９）。そして処理を終了する。

図２８乃至図３０は、表示装置に表示される画像の一例を示す図である。図２８においては、図２７と同じｚ座標における断層画像データと、図２７に示した領域（白い部分）を示す太線２８０とが重ねて表示されている。図２９においては、図２８とは異なるｚ座標における断層画像データと、領域抽出処理によって抽出された領域を示す太線２９０とが重ねて表示されている。また、図３０においては、図２８及び図２９とは異なるｚ座標における断層画像データと、領域抽出処理によって抽出された領域を示す太線３００とが重ねて表示されている。

図２８乃至図３０に示すように、領域抽出処理によって抽出された領域には、乳頭筋部分が含まれるようになっている。また、造影剤の影響によるムラが無くなっている。すなわち、輝度値のみによって左室流体に該当する領域を抽出した場合と比較して、シミュレーションに用いる領域を適切に抽出することができるようになる。

以上本技術の一実施の形態を説明したが、本技術はこれに限定されるものではない。例えば、上で説明した形状データ生成装置１の機能ブロック図は必ずしも実際のプログラムモジュール構成に対応するものではない。

また、上で説明した処理フローは、処理結果が変わらなければ処理の順番を入れ替えることも可能である。さらに、並列に実行させるようにしても良い。

例えば、初期値生成処理においては形状データ生成処理を実行した後に確率マップ生成処理を実行しているが、これらの処理を並行で実行してもよい。

また、上では左室流体を例にして説明したが、その他の複雑な形状を有する臓器であっても本実施の形態の処理を適用することができる。

なお、上で述べた形状データ生成装置１は、コンピュータ装置であって、図３１に示すように、メモリ２５０１とＣＰＵ（Central Processing Unit）２５０３とハードディスク・ドライブ（ＨＤＤ：Hard Disk Drive）２５０５と表示装置２５０９に接続される表示制御部２５０７とリムーバブル・ディスク２５１１用のドライブ装置２５１３と入力装置２５１５とネットワークに接続するための通信制御部２５１７とがバス２５１９で接続されている。オペレーティング・システム（ＯＳ：Operating System）及び本実施例における処理を実施するためのアプリケーション・プログラムは、ＨＤＤ２５０５に格納されており、ＣＰＵ２５０３により実行される際にはＨＤＤ２５０５からメモリ２５０１に読み出される。ＣＰＵ２５０３は、アプリケーション・プログラムの処理内容に応じて表示制御部２５０７、通信制御部２５１７、ドライブ装置２５１３を制御して、所定の動作を行わせる。また、処理途中のデータについては、主としてメモリ２５０１に格納されるが、ＨＤＤ２５０５に格納されるようにしてもよい。本発明の実施例では、上で述べた処理を実施するためのアプリケーション・プログラムはコンピュータ読み取り可能なリムーバブル・ディスク２５１１に格納されて頒布され、ドライブ装置２５１３からＨＤＤ２５０５にインストールされる。インターネットなどのネットワーク及び通信制御部２５１７を経由して、ＨＤＤ２５０５にインストールされる場合もある。このようなコンピュータ装置は、上で述べたＣＰＵ２５０３、メモリ２５０１などのハードウエアとＯＳ及びアプリケーション・プログラムなどのプログラムとが有機的に協働することにより、上で述べたような各種機能を実現する。

以上述べた本実施の形態をまとめると以下のようになる。

本実施の形態に係る形状データ生成方法は、（Ａ）複数の断層画像のうち対象物体が占める領域が指定されている断層画像を用いて、対象物体を表す三次元のボクセルデータを生成し、（Ｂ）対象物体が占める領域が指定されている断層画像から、指定されている領域に含まれるボクセルの輝度値を抽出し、抽出された当該輝度値を用いて、対象物体が占める領域に含まれるボクセルである確率を算出するための関数を生成し、（Ｃ）複数の断層画像を包含するボクセル空間における各ボクセルについて、当該ボクセルの輝度値と関数とを用いて確率を算出し、（Ｄ）対象物体の形状を表す三次元のボクセルデータとボクセル空間における各ボクセルについて算出された確率とを用いて、対象物体の形状を表す三次元のボクセルデータを再生成する処理を含む。

上で述べたように、対象物体が占める領域に含まれるボクセルである確率を用いて三次元のボクセルデータを再生成すれば、高精度の三次元形状になる。

また、上で述べた輝度値を抽出する処理において、（ｂ１）指定されている領域以外の領域に含まれるボクセルの輝度値を抽出し、関数を生成する処理において、（ｂ２）指定されている領域に含まれるボクセルの輝度値についての第１の度数分布と、指定されている領域以外の領域に含まれるボクセルの輝度値についての第２の度数分布とを生成し、第１の度数分布の値を第１の度数分布の値と第２の度数分布の値との和で除した値を確率とする関数を生成してもよい。

このような関数とすれば、指定されている領域に含まれるボクセルの輝度値と同程度の輝度値を有するボクセルの確率が高くなる。

また、上で述べた輝度値を抽出する処理において、（ｂ３）対象物体が占める領域が指定されている断層画像のうち指定されている領域に対しエロージョン処理を実行し、当該エロージョン処理により特定される領域のうち指定されている領域以外の領域から輝度値を抽出してもよい。

エロージョン処理により特定される領域は指定されている領域に肉付けされたような領域になるので、このようにすれば、指定されている領域以外の領域の輝度値を適切に抽出できるようになる。

また、上で述べた対象物体の形状を表す三次元のボクセルデータを再生成する処理において、（ｄ１）対象物体の形状を表す三次元のボクセルデータに含まれるボクセルの輝度値と、当該ボクセルに対応するボクセル空間上のボクセルの確率との加重平均を算出してもよい。

このようにすれば、算出された確率が適切に反映された三次元のボクセルデータを生成できるようになる。

また、本形状データ生成方法が、（Ｅ）再生成された対象物体の形状を表す三次元のボクセルデータに含まれるボクセルの輝度値を、当該輝度値についての反応拡散方程式を解くことにより、対象物体の内部であることを示す第１の輝度値又は対象物体の外部であることを示す第２の輝度値のいずれかに近付けてもよい。

このようにすれば、対象物体と対象物体ではない部分との境界がぼやけている場合であっても、境界を明瞭にすることができるようになる。

また、上で述べた対象物体の形状を表す三次元のボクセルデータを生成する処理において、（ａ１）対象物体の標準的な形状を表す三角形メッシュデータを、指定されている領域に合うように変形し、（ａ２）変形後の三角形メッシュデータにより特定される形状の内部が第１の輝度値になり且つ外部が第２の輝度値となるように、対象物体の形状を表す三次元のボクセルデータを生成してもよい。

このように標準的な形状を利用すれば、最終的に生成される三次元のボクセルデータが現実的なものになる。

また、上で述べた対象物体の形状を表す三次元のボクセルデータを生成する処理において、（ａ３）対象物体の形状を表す三次元のボクセルデータに含まれるボクセルのうち、対象物体と対象物体以外の部分との境界に該当するボクセルを特定し、（ａ４）特定された境界に該当するボクセルの輝度値が第１の輝度値と第２の輝度値との平均値になるようにガウシアンフィルタを実行してもよい。

ボクセル値が第１の輝度値又は第２の輝度値のいずれかであるならば、その平均値が境界に該当すると考えるのが妥当であるからである。

また、上で述べた対象物体の形状を表す三次元のボクセルデータを生成する処理において、（ａ５）対象物体の形状を表す三次元のボクセルデータに含まれるボクセルのうち対象物体以外の部分のボクセルをシード点としてリジョングローイング法による処理を実行してもよい。

このようにすれば、対象物体以外の部分を抽出することができるので、対象物体の内部と外部とを分けることができるようになる。

また、上で述べた対象物体は心臓であってもよい。心臓のような複雑な臓器であっても、上で述べた処理であれば高精度の三次元形状になる。

なお、上記方法による処理をコンピュータに行わせるためのプログラムを作成することができ、当該プログラムは、例えばフレキシブルディスク、ＣＤ−ＲＯＭ、光磁気ディスク、半導体メモリ、ハードディスク等のコンピュータ読み取り可能な記憶媒体又は記憶装置に格納される。尚、中間的な処理結果はメインメモリ等の記憶装置に一時保管される。

Claims (9)

1. コンピュータに、
   複数の断層画像のうち対象物体が占める領域が指定されている断層画像を用いて、前記対象物体を表す三次元のボクセルデータを生成させ、
   前記対象物体が占める領域が指定されている断層画像から、指定されている領域に含まれるボクセルの輝度値を抽出させ、抽出された当該輝度値を用いて、前記対象物体が占める領域に含まれるボクセルである確率を算出するための関数を生成させ、
   前記複数の断層画像を包含するボクセル空間における各ボクセルについて、当該ボクセルの輝度値と前記関数とを用いて前記確率を算出させ、
   前記対象物体の形状を表す三次元のボクセルデータと前記ボクセル空間における各ボクセルについて算出された確率とを用いて、前記対象物体の形状を表す三次元のボクセルデータを再生成させる
   ことを特徴とする形状データ生成プログラム。
2. 前記輝度値を抽出させる処理において、
   前記指定されている領域以外の領域に含まれるボクセルの輝度値を抽出させ、
   前記関数を生成させる処理において、
   前記指定されている領域に含まれるボクセルの輝度値についての第１の度数分布と、前記指定されている領域以外の領域に含まれるボクセルの輝度値についての第２の度数分布とを生成させ、前記第１の度数分布の値を前記第１の度数分布の値と前記第２の度数分布の値との和で除した値を前記確率とする関数を生成させる
   ことを特徴とする請求項１記載の形状データ生成プログラム。
3. 前記輝度値を抽出させる処理において、
   前記対象物体が占める領域が指定されている断層画像のうち前記指定されている領域に対しエロージョン処理を実行させ、当該エロージョン処理により特定される領域のうち前記指定されている領域以外の領域から前記輝度値を抽出させる
   ことを特徴とする請求項２記載の形状データ生成プログラム。
4. 前記対象物体の形状を表す三次元のボクセルデータを再生成させる処理において、
   前記対象物体の形状を表す三次元のボクセルデータに含まれるボクセルの輝度値と、当該ボクセルに対応する前記ボクセル空間上のボクセルの確率との加重平均を算出させる
   ことを特徴とする請求項１乃至３のいずれか１つ記載の形状データ生成プログラム。
5. 前記コンピュータに、さらに、
   再生成された前記対象物体の形状を表す三次元のボクセルデータに含まれるボクセルの輝度値を、当該輝度値についての反応拡散方程式を解かせることにより、前記対象物体の内部であることを示す第１の輝度値又は前記対象物体の外部であることを示す第２の輝度値のいずれかに近付けさせる
   ことを特徴とする請求項１乃至４のいずれか１つ記載の形状データ生成プログラム。
6. 前記対象物体の形状を表す三次元のボクセルデータを生成させる処理において、
   前記対象物体の標準的な形状を表す三角形メッシュデータを、前記指定されている領域に合うように変形させ、
   変形後の前記三角形メッシュデータにより特定される形状の内部が前記第１の輝度値になり且つ外部が前記第２の輝度値となるように、前記対象物体の形状を表す三次元のボクセルデータを生成させる
   ことを特徴とする請求項５記載の形状データ生成プログラム。
7. 前記対象物体の形状を表す三次元のボクセルデータを生成させる処理において、
   前記対象物体の形状を表す三次元のボクセルデータに含まれるボクセルのうち、前記対象物体と前記対象物体以外の部分との境界に該当するボクセルを特定させ、
   特定された前記境界に該当するボクセルの輝度値が前記第１の輝度値と前記第２の輝度値との平均値になるようにガウシアンフィルタを実行させる
   ことを特徴とする請求項６記載の形状データ生成プログラム。
8. コンピュータが、
   複数の断層画像のうち対象物体が占める領域が指定されている断層画像を用いて、前記対象物体を表す三次元のボクセルデータを生成し、
   前記対象物体が占める領域が指定されている断層画像から、指定されている領域に含まれるボクセルの輝度値を抽出し、抽出された当該輝度値を用いて、前記対象物体が占める領域に含まれるボクセルである確率を算出するための関数を生成し、
   前記複数の断層画像を包含するボクセル空間における各ボクセルについて、当該ボクセルの輝度値と前記関数とを用いて前記確率を算出し、
   前記対象物体の形状を表す三次元のボクセルデータと前記ボクセル空間における各ボクセルについて算出された確率とを用いて、前記対象物体の形状を表す三次元のボクセルデータを再生成する
   処理を実行する形状データ生成方法。
9. 複数の断層画像のうち対象物体が占める領域が指定されている断層画像を用いて、前記対象物体を表す三次元のボクセルデータを生成する第１生成部と、
   前記対象物体が占める領域が指定されている断層画像から、指定されている領域に含まれるボクセルの輝度値を抽出し、抽出された当該輝度値を用いて、前記対象物体が占める領域に含まれるボクセルである確率を算出するための関数を生成し、前記複数の断層画像を包含するボクセル空間における各ボクセルについて、当該ボクセルの輝度値と前記関数とを用いて前記確率を算出する算出部と、
   前記対象物体の形状を表す三次元のボクセルデータと前記ボクセル空間における各ボクセルについて算出された確率とを用いて、前記対象物体の形状を表す三次元のボクセルデータを再生成する第２生成部と、
   を有する形状データ生成装置。