JP7155670B2 - 医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法 - Google Patents

Description

本開示は、医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法に関する。

医用画像診断において、特定の人体組織の観察に適したボリュームレンダリング像を得たい場面がある。例えば、胸部や頭部にある内臓や脳等の所定の人体組織の観察を行いたい場合に、内臓や脳等は骨に囲まれており、骨領域はむしろ診断の妨げになる。ＣＴ画像から生成されるボリュームレンダリング像では一般に骨が鮮明に描写され、内臓や血管は隠れてしまうことがあるため、可視化にあたっては切断を行うなど工夫が必要となることがある。骨領域は比較的高信号であるが、信号値分布が広範で造影血管と被るため、カラーマップを工夫して骨領域だけを透明化したオパシティカーブを設計することは容易ではなかった。

たとえば、下記の特許文献１では、ＣＴ画像の信号値に基づいて領域の抽出とラベリングを行い、ラベリングされた各肋骨領域が適正か否かを、別途構築した統計上の肋骨の体軸方向の高さ幅情報を有する統計データベースを用いた手法が開示されている。これにより、骨領域（肋骨領域）を検出して除去（透明化）させることができるが、特許文献１の方法では、肋骨に先天性の奇形があったり、骨折を伴う場合には判定誤りが発生し、適切に除去できなくなる。

また、領域拡張法（リージョングローイング法）を用いて非観察対象としたい骨領域を抽出した３次元マスクを作成し、３次元マスクを参照しながらマスク処理により非観察対象としたい骨領域を除去したボリュームレンダリング像を生成する方法が開示されている（特許文献２、３参照）。領域拡張法とは、非観察対象領域の画素を指定し、その画素を開始点（拡張開始点）として、近傍画素を次々と抽出する方法である。

特開２００９－２４０５６９号公報 特許４０８７５１７号公報 特許５２５７９５８号公報

しかしながら、領域拡張法を用いた３次元マスクの作成は、ユーザによる複数の拡張開始点の指示が必須で、拡張の打ち切り段階もユーザが指示する必要があるため、自動化が難しく、ユーザごとに結果にバラつきが発生するという問題がある。また領域拡張法に限らず、３次元マスクの作成は、ユーザのスキルや経験が必要であり、その作成に時間や手間がかかる。

これに対し、本発明者は、被写体領域の外郭を単一の幾何形状で近似し、近似した幾何形状の情報に基づいてボクセルの不透明度を制御することで、容易に、所望の観察対象を鮮明に表示させる方法を開発している（特願２０１８－１０３７６７）。

本開示は、この方法を更に拡張するものであり、被写体領域の外郭を複数の幾何形状で近似することで、所望の観察対象を、選択的にかつ鮮明に表示することを目的とする。

本開示の一実施形態によると、複数の断層画像を取得する画像取得手段と、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出手段と、前記断層画像毎に、観察対象や非観察対象の位置関係に応じて、抽出され被写体領域の外郭を近似する複数の幾何形状を算出し、当該複数の幾何形状の情報を取得する幾何情報取得手段と、を備え、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正する医用画像処理装置が提供される。

また、前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する複数の幾何形状の情報に基づいて算出する補正倍率算出手段と、を更に備え、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算するようにしても良い。

また、前記補正倍率算出手段は、算出した補正倍率を格納する補正テーブルを作成するテーブル作成手段と、を更に備え、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正テーブルを参照して得られる補正倍率を乗算するようにしても良い。

また、前記補正倍率算出手段は、前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくするようにしても良い。

また、前記補正倍率算出手段は、前記複数の幾何形状が重複する範囲においては、各々の幾何形状の情報に基づいて算出される補正倍率のうち、最も値が大きい補正倍率を採用するようにしても良い。

また、抽出された被写体領域に複数の個別領域を設定する個別領域設定手段と、を更に備え、前記幾何情報取得手段は、個別領域毎に被写体領域の外郭を幾何形状で近似し、個別領域毎に幾何形状の情報を取得するようにしても良い。

また、前記個別領域設定手段は、抽出された被写体領域を横方向に分割した２つの個別領域を設定し、更に、被写体領域の横方向の中央付近に前記個別領域より小さい個別領域を設定するようにしても良い。

なお本開示において、横方向とは、断層画像の座標系を規定する２つの座標軸方向（Ｘ軸方向、Ｙ軸方向）のうち、該断層画像に映る被写体の左右軸の方向と近いほうの座標軸方向のことを指す。また、縦方向とは、断層画像の座標系を規定する２つの座標軸方向（Ｘ軸方向、Ｙ軸方向）のうち、該断層画像に映る被写体の背腹軸（前後軸）の方向と近いほうの座標軸方向のことを指す。ここで、断層画像のＸ軸およびＹ軸は、例えば、モダリティでの撮影時に適宜設定されるが、本開示では、被写体の左右軸と方向が近いほうの座標軸をＸ軸、被写体の背腹軸（前後軸）と方向が近いほうの座標軸をＹ軸とする。すなわち、本開示では、横方向は断層画像のＸ軸方向に相当し、縦方向は断層画像のＹ軸方向に相当するものとして説明する。

また、前記幾何情報取得手段は、前記中央付近に設定した個別領域に対応する幾何形状の幾何中心に対してオフセットを加えることで、前記中央付近に設定した個別領域に対応する幾何形状の情報を補正するようにしても良い。

また、前記幾何情報取得手段は、各幾何形状の幾何中心に対して、スライス位置が所定の位置から末端の方向に位置するにつれ所定の割合で増大するオフセットを加えることで、各幾何形状の情報を補正するようにしても良い。

また、前記幾何情報取得手段は、断層画像の総数をＳｚ、スライス位置をｚ、前記所定の位置をＺｔｈとすると、ｚ≧Ｚｔｈの場合は（ｚ－Ｚｔｈ）／（Ｓｚ－１－Ｚｔｈ）に比例したオフセットを各幾何形状の幾何中心のＹ座標に加えるようにしても良い。

また、前記幾何形状は楕円であり、前記幾何情報取得手段は、前記幾何形状の情報として、楕円の中心座標、楕円の横方向のサイズ、及び楕円の縦方向のサイズを取得し、前記補正倍率算出手段は、前記中心座標、前記横方向のサイズ、及び前記縦方向のサイズに基づいて、前記補正倍率を算出するようにしても良い。

また、前記補正倍率算出手段は、前記中心座標における前記補正倍率を最大とし、前記中心座標と前記横方向のサイズ及び前記縦方向のサイズの比が同一となる前記楕円内部の各同心楕円上の前記補正倍率を、各同心楕円の横方向のサイズ及び縦方向のサイズが大きくなるにつれ小さくするようにしても良い。

また、前記幾何情報取得手段は、前記中心座標に所定のオフセットを加算して前記中心座標を補正し、前記横方向のサイズ及び前記縦方向のサイズに所定の横方向のサイズに対する倍率及び縦方向のサイズに対する倍率を乗算して前記横方向のサイズ及び前記縦方向のサイズを補正し、前記補正倍率算出手段は、補正された前記中心座標、前記横方向のサイズ、及び前記縦方向のサイズに基づいて、前記補正倍率を算出するようにしても良い。

また、前記オフセット、前記横方向のサイズに対する倍率、及び前記縦方向のサイズに対する倍率をユーザに調整させる調整手段と、を更に備えるようにしても良い。

また、被写体領域の抽出精度を所定の基準で判断し、所定の基準を満たさない前記断層画像の被写体領域を、所定の基準を満たす他の前記断層画像の被写体領域に基づいて補正する被写体領域補正手段と、を更に備えるようにしても良い。

また、前記被写体領域補正手段は、被写体領域の横方向のサイズが画像の横方向のサイズと一致する前記断層画像における被写体領域の横方向のサイズを、被写体領域の横方向のサイズが画像の横方向のサイズより小さい他の前記断層画像における被写体領域の横方向のサイズ及び縦方向のサイズの比と同一となるように、補正するようにしても良い。

また、信号値と所定の閾値との差分値に応じて当該信号値に対応する不透明度を減衰させるように前記カラーマップを調整するカラーマップ調整手段と、を更に備え、前記ボクセル作成手段は、調整された前記カラーマップを参照して、ボクセルデータを作成するようにしても良い。

また、代表的な信号値と色値及び不透明度との対応関係に基づいて、所定範囲の信号値と色値及び不透明度との対応関係を定義するカラーマップを作成するカラーマップ作成手段と、を更に備え、前記カラーマップ取得手段は、作成した前記カラーマップを取得するようにしても良い。

また、前記ボクセル作成手段は、各ボクセルの色値に対して、当該ボクセルの不透明度に基づいて当該ボクセルにおける勾配ベクトルを算出し、勾配ベクトルとあらかじめ定義された光源ベクトルに基づいて拡散反射成分と環境光成分から構成される陰影値を算出し、前記色値に前記陰影値を乗算することにより、前記色値を改変する陰影付加手段と、を更に備えるようにしても良い。

また、前記陰影付加手段は、前記陰影付加手段は、前記勾配ベクトルを算出するにあたり、当該ボクセルのＸ軸方向、Ｙ軸方向、Ｚ軸方向の各々２近傍のボクセルの不透明度の差分値を前記勾配ベクトルのＸ方向成分、Ｙ方向成分、Ｚ方向成分として算出し、あらかじめ定義されたＺ軸方向変倍率に基づいて補正を施したベクトルようにしても良い。

また、前記レンダリング手段は、前記ボクセルデータを生成するボリュームレンダリング像に投影変換した座標系を視点座標系とすると、視点座標系において、前記ボリュームレンダリング像の各画素よりＺ軸方向に沿って、Ｚ軸の上限値より下限値に向けて、視点座標系のボクセル座標毎に座標変換を行って前記ボクセルデータより不透明度を取得しながら、不透明ボクセルを探索し、最初に見つかった不透明ボクセルの視点座標系におけるＺ座標を、前記ボリュームレンダリング像の画素毎に記録した探索制御マスクを作成する探索制御マスク作成手段と、前記ボリュームレンダリング像の画素毎に、前記探索制御マスクからＺ座標を取得し、取得したＺ座標よりＺ軸の下限値に向けてＺ軸方向に沿って、所定の光強度をもつ仮想光線を照射する際、視点座標系のボクセル座標毎に座標変換を行って前記ボクセルデータより不透明度を取得し、不透明ボクセルが見つかった場合、当該ボクセル座標に対して座標変換を行って前記ボクセルデータより色値を取得し、当該ボクセルの不透明度に基づいて前記光強度を減衰させるとともに、当該ボクセルの不透明度及び色値並びに前記減衰させた光強度に基づいて累積輝度値を算出する処理を繰り返し、算出された累積輝度値に基づいて、前記ボリュームレンダリング像の当該画素に対応する画素値として与えるレイキャスティング手段と、を備えるようにしても良い。

また、前記探索制御マスク作成手段及び前記レイキャスティング手段は、前記座標変換を行って前記ボクセルデータより不透明度または色値を取得する際、所定の回転を定義した回転行列、ＸＹＺ軸各方向のオフセット値、ＸＹＺ軸方向の拡大又は縮小倍率、Ｚ軸方向の変倍率、注視点から視点までの距離を含む前記所定の座標変換のパラメータを取得し、前記視点座標系のボクセルの整数値の座標を、前記パラメータに基づいて前記ボクセルデータの座標系に変換を行って、前記ボクセルデータの実数値の座標を算出し、算出した実数値の座標の近傍の複数の整数値の座標に対応する前記ボクセルデータの複数のボクセルを特定し、特定した複数のボクセルの不透明度または色値に基づいて前記ボクセルデータより取得される不透明度または色値として算出するようにしても良い。

また、前記レンダリング手段は、前記ボクセルデータに基づいて３Ｄテクスチャを生成する３Ｄテクスチャ生成手段と、前記３Ｄテクスチャに対して所定の座標変換を行って変換後３Ｄテクスチャを生成する座標変換手段と、３次元空間のＸＹ座標面上の四角形をＺ軸方向に並べた積層四角形を設定する積層四角形設定手段と、所定の視点からＺ軸方向に平行な視線上の前記四角形のＸＹ座標に対応する前記変換後３Ｄテクスチャのボクセルの色値を前記ボクセルの不透明度に基づいて前記視点から遠い四角形の順にアルファブレンディングして取得し、前記ボリュームレンダリング像の画素値として与える画素値算出手段と、を備えるようにしても良い。

また、前記座標変換手段は、所定の回転を定義した回転行列、視野角度、視点位置、クリッピング位置、ＸＹＺ軸各方向のオフセット値、ＸＹＺ軸方向の拡大又は縮小倍率、Ｚ軸方向の変倍率を含む所定の座標変換のパラメータを取得し、
前記３Ｄテクスチャに対して、前記取得したパタメータを用いた前記所定の座標変換を行って前記変換後３Ｄテクスチャを生成するようにしても良い。

また、前記座標変換手段及び前記画素値算出手段は、ビデオカードに搭載されたＧＰＵ及びフレームメモリを用いて実行するようにしても良い。

また、本開示の一実施の形態によると、複数の断層画像を取得する画像取得手段と、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出手段と、前記断層画像毎に、抽出された被写体領域の外郭を複数の幾何形状で近似し、当該複数の幾何形状の情報を取得する幾何情報取得手段と、前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する複数の幾何形状の情報に基づいて算出する補正倍率算出手段と、を備え、前記補正倍率算出手段は、前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、前記複数の幾何形状が重複する範囲においては、各々の幾何形状の情報に基づいて算出される補正倍率のうち、最も値が大きい補正倍率を採用し、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算して補正する医用画像処理装置が提供される。
また、本開示の一実施の形態によると、複数の断層画像を取得する画像取得手段と、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出手段と、前記断層画像毎に、抽出された被写体領域の外郭を複数の幾何形状で近似し、当該複数の幾何形状の情報を取得する幾何情報取得手段と、を備え、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正し、各ボクセルの色値に対して、当該ボクセルの不透明度に基づいて当該ボクセルにおける勾配ベクトルを算出し、勾配ベクトルとあらかじめ定義された光源ベクトルに基づいて拡散反射成分と環境光成分から構成される陰影値を算出し、前記色値に前記陰影値を乗算することにより、前記色値を改変する陰影付加手段を更に備える医用画像処理装置が提供される。
また、前記陰影付加手段は、前記勾配ベクトルを算出するにあたり、当該ボクセルのＸ軸方向、Ｙ軸方向、Ｚ軸方向の各々２近傍のボクセルの不透明度の差分値を前記勾配ベクトルのＸ方向成分、Ｙ方向成分、Ｚ方向成分として算出し、あらかじめ定義されたＺ軸方向変倍率に基づいて補正を施したベクトルを、当該ボクセルにおける勾配ベクトルとして算出しても良い。
また、本開示の一実施形態によると、コンピュータの制御部が、複数の断層画像を取得する画像取得ステップと、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得ステップと、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成ステップと、記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリングステップと、前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出ステップと、前記断層画像毎に、観察対象や非観察対象の位置関係に応じて、抽出された被写体領域の外郭を近似する複数の幾何形状を算出し、当該複数の幾何形状の情報を取得する幾何情報取得ステップと、を含み、前記ボクセル作成ステップは、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正する医用画像処理方法が提供される。

また、本開示の一実施形態によると、コンピュータを、複数の断層画像を取得する画像取得手段、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段、前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出手段、前記断層画像毎に、観察対象や非観察対象の位置関係に応じて、抽出された被写体領域の外郭を近似する複数の幾何形状を算出し、当該複数の幾何形状の情報を取得する幾何情報取得手段、として機能させ、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正するプログラムが提供される。

また、本開示の一実施形態によると、コンピュータの制御部が、複数の断層画像を取得する画像取得ステップと、前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出ステップと、前記断層画像毎に、観察対象や非観察対象の位置関係に応じて、抽出された被写体領域を近似する複数の幾何形状を算出し、当該複数の幾何形状の情報を取得する幾何情報取得ステップと、前記断層画像の各画素の信号値に応じたオパシティカーブにより規定される不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正する不透明度補正ステップと、前記断層画像の各画素に対応する補正後の不透明度を格納したデータを作成するデータ作成ステップと、を含むデータ作成方法が提供される。

本開示によれば、カラーマップを参照して断層画像の各画素の信号値を色及び不透明度に変換してボクセルデータを作成する際、カラーマップ（オパシティカーブ）に定義された不透明度を、断層画像の被写体領域の外郭に近似させた複数の幾何形状の情報に基づいて補正する。このように、被写体領域の外郭を近似した複数の幾何形状に基づいて不透明度をコントロールすることで、所望の観察対象を選択的にかつ鮮明に可視化させたレンダリング像を得ることができる。

医用画像処理装置１のハードウェア構成を示す図 ボリュームレンダリング処理の概要を示す図 医用画像処理装置１の機能構成を示す図 補正テーブル作成部２６の内部構成を示す図 レンダリング部２８の内部構成を示す図 医用画像処理装置１の全体動作を示すフローチャート 調整されたカラーマップＣｍａｐの性質を示す図 補正テーブル作成処理を示すフローチャート 被写体領域の抽出について説明する図 被写体領域の補正について説明する図 被写体領域の外郭を単一の楕円で近似する例を示す図 補正倍率Ｓαの２次元分布パターンを示す図 ＸＹ方向の振幅倍率と補正倍率との関係を示す図 被写体領域の外郭を複数の楕円で近似する例を示す図 補正量ｄｙ（ｚ）の関数と胴体のサジタル像（ＭＰＲ像）を表す図 補正倍率Ｓαの２次元分布パターンを示す図 補正倍率Ｓαを３次元分布表示した図 ボクセル作成処理を示すフローチャート レイキャスティング法の概要を示す図 レイキャスティング法によるレンダリング処理を示すフローチャート 探索制御マスク算出処理を示すフローチャート 探索制御マスク算出処理を説明する図 不透明ボクセル探索処理を示すフローチャート レイキャスティング処理を示すフローチャート ３Ｄテクスチャマッピング法ついて説明する図 色補正処理を示すフローチャート ３Ｄテクスチャマッピング法によるレンダリング処理を示すフローチャート 投影画面設定処理を示すフローチャート レンダリング処理を示すフローチャート パラメータ調整画面４０の表示例を示す図 従来手法により生成した胸部レンダリング像（正面） 提案手法１により生成した胸部レンダリング像（正面） 提案手法１により生成した胸部レンダリング像（背面） 提案手法２により生成した胸部レンダリング像（背面） 提案手法１により生成した胸部レンダリング像（側面） 提案手法２により生成した胸部レンダリング像（側面） 従来手法により生成した腹部レンダリング像（正面） 提案手法１により生成した腹部レンダリング像（正面） 提案手法１により生成した腹部レンダリング像（背面） 提案手法２により生成した腹部レンダリング像（背面） 提案手法１により生成した腹部レンダリング像（側面） 提案手法２により生成した腹部レンダリング像（側面）

以下図面に基づいて、本開示の実施の形態を詳細に説明する。本実施の形態では、Ｘ線ＣＴ装置により撮影されたＣＴ画像に基づいてボリュームレンダリング像を生成する場合について説明するが、本開示はＣＴ画像以外の医用画像（ＭＲＩ画像やＰＥＴ画像等）に基づいてボリュームレンダリング像を生成する場合にも適用可能である。

図１は、本実施の形態における医用画像処理装置１のハードウェア構成を示すブロック図である。図１に示すように、医用画像処理装置１は、制御部１１、記憶部１２、メディア入出力部１３、通信制御部１４、入力部１５、表示部１６、周辺機器Ｉ／Ｆ部１７等が、バス１８を介して接続される汎用のコンピュータで実現される。但し、これに限ることなく、用途、目的に応じて様々な構成を採ることが可能である。

制御部１１は、ＣＰＵ（Central Processing Unit）、ＧＰＵ（Graphics Processing Unit）、ＲＯＭ（Read Only
Memory）、ＲＡＭ（Random Access Memory）、フレームメモリ(Frame Memory)等によって構成される。ＣＰＵ、ＧＰＵは、記憶部１２、ＲＯＭ、記録媒体等に格納されるプログラムをＲＡＭ、フレームメモリ上のワークメモリ領域に呼び出して実行し、バス１８を介して接続された各装置を駆動制御し、医用画像処理装置１が行う後述する処理を実現する。

ＧＰＵ及びフレームメモリは、医用画像処理装置１のビデオカードに搭載されている。ＲＯＭは、不揮発性メモリであり、コンピュータのブートプログラムやＢＩＯＳ等のプログラム、データ等を恒久的に保持している。ＲＡＭ、フレームメモリは、揮発性メモリであり、記憶部１２、ＲＯＭ、記録媒体等からロードしたプログラム、データ等を一時的に保持するとともに、制御部１１が各種処理を行う為に使用するワークエリアを備える。

記憶部１２は、ＨＤＤ（Hard Disk Drive）等であり、制御部１１が実行するプログラム、プログラム実行に必要なデータ、ＯＳ（Operating System）等が格納される。プログラムに関しては、ＯＳに相当する制御プログラムや、後述する処理をコンピュータに実行させるためのアプリケーションプログラムが格納されている。これらの各プログラムコードは、制御部１１により必要に応じて読み出されてＲＡＭ、フレームメモリに移され、ＣＰＵ、ＧＰＵに読み出されて各種の手段として実行される。

メディア入出力部１３（ドライブ装置）は、データの入出力を行い、例えば、ＣＤドライブ（－ＲＯＭ、－Ｒ、－ＲＷ等）、ＤＶＤドライブ（－ＲＯＭ、－Ｒ、－ＲＷ等）等のメディア入出力装置を有する。通信制御部１４は、通信制御装置、通信ポート等を有し、コンピュータとネットワーク間の通信を媒介する通信インタフェースであり、ネットワークを介して、他のコンピュータ間との通信制御を行う。ネットワークは、有線、無線を問わない。

入力部１５は、データの入力を行い、例えば、キーボード、マウス等のポインティングデバイス、テンキー等の入力装置を有する。入力部１５を介して、コンピュータに対して、操作指示、動作指示、データ入力等を行うことができる。表示部１６は、液晶パネル等のディスプレイ装置、ディスプレイ装置と連携してコンピュータのビデオ機能を実現するための論理回路等（ビデオアダプタ等）を有する。なお、入力部１５及び表示部１６は、タッチパネルディスプレイのように、一体となっていてもよい。

周辺機器Ｉ／Ｆ（Interface）部１７は、コンピュータに周辺機器を接続させるためのポートであり、周辺機器Ｉ／Ｆ部１７を介してコンピュータは周辺機器とのデータの送受信を行う。周辺機器Ｉ／Ｆ部１７は、ＵＳＢ（Universal Serial Bus）やＩＥＥＥ１３９４やＲＳ－２３２Ｃ等によって構成されており、通常複数の周辺機器Ｉ／Ｆを有する。周辺機器との接続形態は有線、無線を問わない。バス１８は、各装置間の制御信号、データ信号等の授受を媒介する経路である。

医用画像処理装置１は、１台のコンピュータで構成されてもよいし、複数のコンピュータがネットワークを介して構成されてもよい。例えば、医用画像処理装置１が、サーバとクライアント端末で構成される場合、クライアント端末においてデータの入力を受け付けて、サーバが各種の処理を行い、クライアント端末が処理結果を表示するようにしてもよい。以下の説明では、簡素な構成例として、医用画像処理装置１が１台のコンピュータで構成された例を説明する。

図２は、医用画像処理装置１が実行するボリュームレンダリング処理の概要を示す図である。図２に示すように、複数の断層画像（図の例では、５１２×５１２ピクセルの３７０枚の胸部ＣＴ画像）に基づいてボリュームレンダリング処理を実行し、ボリュームレンダリング像として、遠位の視点から被写体を観察する全体レンダリング像（図左下）や視点を被写体内に自由に移動させて気管支や大腸等を観察する仮想内視鏡像（図右下）を生成し表示する。本実施の形態では、主に、全体レンダリング像（図左下）を生成する場合について説明する。

図３は、医用画像処理装置１の機能構成を示す図である。図３に示すように、医用画像処理装置１は、画像取得部２１、階調圧縮部２２、領域指定部２３、カラーマップ取得部２４、カラーマップ調整部２５、補正テーブル作成部２６、ボクセル作成部２７、レンダリング部２８、データ出力部２９、及びパラメータ調整部３０を備える。

画像取得部２１は、ＣＴ装置により撮影された断層画像群Ｄｏを取得する。断層画像群Ｄｏは、被写体（人体）を頭尾軸に沿って所定のスライス間隔で連続的に撮影した複数の断層画像（ＣＴ画像）からなる。各断層画像はＤＩＣＯＭ形式の２次元の画像データである。ＤＩＣＯＭ形式は、１ファイルにヘッダ部と画像データ部を含む医療画像で一般的に用いられる画像フォーマットであり、画像撮影時のパラメータや診断情報を保存しておくことができる。

１つの断層画像は、例えば、５１２×５１２ピクセルの画像（ＣＴ画像）である。断層画像の各画素には、信号値ｖが付与されており、ＣＴ画像の場合、信号値ｖはＣＴ値である。本実施の形態では、信号値ｖ（ＣＴ値）は１６ビット（－３２７６８≦ｖ≦３２７６７）のデータとする（但し、信号値ｖのビット数は特に限定されない）。

ＣＴ値は、水を基準として表現した組織のＸ線減弱係数であり、ＣＴ値により組織や病変の種類等が判断できるようになっている(単位はＨＵ（Hounsfield Unit））。ＣＴ値は、水と空気のＸ線減弱係数で標準化されており、水のＣＴ値は０、空気のＣＴ値を－１０００である。また、脂肪のＣＴ値は－１２０～－１００程度であり、通常組織のＣＴ値は０～１２０程度であり、骨のＣＴ値は１０００前後を示す。

断層画像群Ｄｏは、ＸＹの２次元データである断層画像をＺ軸方向に積層したものであり、ＸＹＺの３次元データとして表現可能である。例えば断層画像群Ｄｏは、以下のように、ＸＹＺの３次元データ（Ｄｏ（ｘ、ｙ、ｚ））として定義される。尚、本開示において、Ｘ軸は人体の左右軸、Ｙ軸は人体の背腹軸（前後軸）、Ｚ軸は人体の頭尾軸（上下軸）に相当するものとする。ＸＹ軸は、例えば、ＣＴ装置での撮影時に設定される。また、本開示において、Ｘ軸方向を横方向、Ｙ軸方向を縦方向と呼ぶ場合がある。

（式１）
－３２７６８≦Ｄｏ（ｘ、ｙ、ｚ）≦３２７６７
０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｚ≦Ｓｚ－１
解像度：Ｒｘｙ、解像度Ｒｚ

（式１）において、Ｓｘは断層画像の横方向（Ｘ軸方向）の画素数、Ｓｙは断層画像の縦方向（Ｙ軸方向）の画素数、Ｓｚはスライス枚数（断層画像の枚数）を表す。ＲｘｙはＸ軸方向及びＹ軸方向の断層画像の解像度であり、画素の間隔の逆数、すなわち単位距離あたりの画素数を示す。Ｒｚはスライスの解像度であり、断層画像のスライス間隔（例えば、０．５ｍｍや１ｍｍ）の逆数、すなわち単位距離あたりのスライス枚数を表す。

尚、断層画像群Ｄｏのうち、所定のｚ番目（０≦ｚ≦Ｓｚ－１）の断層画像を表す場合は、「断層画像Ｄｏ（ｚ）」と表記することがある。

階調圧縮部２２は、画像取得部２１により取得した断層画像群Ｄｏの信号値を、必要に応じて８ビットに階調圧縮する（（式６）参照）。８ビットに階調圧縮された断層画像Ｄｏ（ｘ、ｙ、ｚ）は以下のように定義される。

（式２）
０≦Ｄｏ（ｘ、ｙ、ｚ）≦２５５
０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｚ≦Ｓｚ－１
解像度：Ｒｘｙ、解像度Ｒｚ

このように階調圧縮をすることで、断層画像群を保持するためのメモリ容量を半分に抑えることができる。たとえ信号値の階調が１６ビットあっても、カラーマップＣｍａｐにより、変換される色値（具体的にはＲＧＢ値）及び不透明度の階調はディスプレイの階調に合わせて８ビットに制限されるため、階調圧縮に伴う画質劣化は殆ど生じない。

領域指定部２３は、ユーザから被写体の関心領域（ＲＯＩ：Region of
Interest）の指定を受け付ける。例えば、関心領域ＲＯＩを直方体で指定する場合は、以下のように、Ｘ方向ＲＯＩ（Ｘ方向の開区間）、Ｙ方向ＲＯＩ（Ｙ方向の開区間）、Ｚ方向ＲＯＩ（Ｚ方向の開区間）を指定する。

（式３）
Ｘ方向ＲＯＩ＝（Ｘｓ、Ｘｅ）
Ｙ方向ＲＯＩ＝（Ｙｓ、Ｙｅ）
Ｚ方向ＲＯＩ＝（Ｚｓ、Ｚｅ）
尚、関心領域を指定しない場合は、Ｘｓ＝０、Ｘｅ＝Ｓｘ－１、Ｙｓ＝０、Ｙｅ＝Ｓｙ－１、Ｚｓ＝０、Ｚｅ＝Ｓｚ－１となる。

カラーマップ取得部２４は、断層画像群Ｄｏに適用するカラーマップＣｍａｐを取得する。カラーマップ取得部２４は、取得したカラーマップＣｍａｐを更に調整する場合には、カラーマップ調整部２５に出力し、カラーマップＣｍａｐの調整を行わない場合には、ボクセル作成部２７に出力する。

カラーマップＣｍａｐは、信号値ｖと色値（具体的にはＲＧＢ値）及び不透明度(α値)との対応関係を定義するものであり、信号値ｖを２４ビットの色（ＲＧＢ値）及び８ビットの不透明度（α値）に変換する関数（実体的には２次元のデータテーブル）として表現可能である。例えば、１６ビットの断層画像群Ｄｏ（（式１）参照）に適用されるカラーマップＣｍａｐは次のように定義される。

（式４）
０≦Ｃｍａｐ（ｖ、ｎ）≦２５５
－３２７６８≦ｖ≦３２７６７
ｎ＝０（Ｒ）、１（Ｇ）、２（Ｂ）、３（α）

また、８ビットの断層画像群Ｄｏ（（式２）参照）に適用されるカラーマップＣｍａｐ（ｖ、ｎ）は次のように定義される。

（式５）
０≦Ｃｍａｐ（ｖ、ｎ）≦２５５
０≦ｖ≦２５５
ｎ＝０（Ｒ）、１（Ｇ）、２（Ｂ）、３（α）

（式４）、（式５）に示すカラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦３）のうち、特にカラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦２）は、信号値ｖを色値（ＲＧＢ値）に変換する関数に相当する。本開示では、信号値ｖを色値に変換するカラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦２）を「カラーパレット」と称する。

また、（式４）、（式５）に示すカラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦３）のうち、特にカラーマップＣｍａｐ（ｖ、３）は、信号値ｖを不透明度に変換する関数に相当し、一般的に「オパシティカーブ」と呼ばれる。

カラーマップＣｍａｐ（ｖ、３）（オパシティカーブ）は、例えば、通常組織（１６ビットの信号値ｖ＝０～１２０程度）の不透明度が２５５に設定（不透明度＝２５５は不透明（光がすべて反射）であることを示す）され、骨領域（１６ビットの信号値ｖ＝１０００前後）の不透明度が１～２５４の中間値に設定（不透明度＝１～２５４は半透明（入射光の一部が反射され、その他は透過）であることを示す）され、その他の組織等の部位を不透明度が０に設定（不透明度＝０は透明（入射光の全てが透過）であることを示す）される。

カラーマップ調整部２５は、必要に応じて、カラーマップＣｍａｐ（ｖ、ｎ）のうち、信号値ｖと不透明度との対応関係を定義するＣｍａｐ（ｖ、３）（オパシティカーブともいう）を信号値に応じて調整し（（式７）（式８）参照）、調整したカラーマップＣｍａｐをボクセル作成部２７に出力する。

補正テーブル作成部２６は、断層画像群Ｄｏの各画素に対応する不透明度の補正倍率を格納した補正テーブルＳαを作成し、ボクセル作成部２７に出力する。

補正テーブル作成部２６は、図４に示すように、被写体領域抽出部２６ａ、被写体領域補正部２６ｂ、個別領域設定部２６ｃ、幾何情報取得部２６ｄ、及び補正倍率算出部２６ｅから構成される。

被写体領域抽出部２６ａは、断層画像Ｄｏ（ｚ）毎に、断層画像Ｄｏ（ｚ）の被写体領域を抽出する（図９、（式９）参照）。

被写体領域補正部２６ｂは、断層画像Ｄｏ（ｚ）毎に、被写体領域抽出部２６ａにより抽出された被写体領域の抽出精度を所定の基準で判断し、所定の基準を満たさない（抽出精度が良好でない）断層画像Ｄｏ（ｚ）の被写体領域を、所定の基準を満たす（抽出精度が良好な）他の断層画像Ｄｏ（ｚ）の被写体領域の情報に基づいて補正する（図１０、（式１０）参照）。

個別領域設定部２６ｃは、必要に応じて、抽出された被写体領域に複数の個別領域を設定する（図１４参照）。

幾何情報取得部２６ｄは、断層画像Ｄｏ（ｚ）毎に、被写体領域の外郭を１以上の幾何形状で近似し、幾何形状の情報を取得する。

具体的には、被写体領域に個別領域が設定されていない場合には、幾何情報取得部２６ｄは、被写体領域の外郭を単一の幾何形状で近似し、当該単一の幾何形状の情報（幾何パラメータＰ（ｚ）、（式１１）参照））を取得する。

被写体領域に複数の個別領域が設定されている場合には、幾何情報取得部２６ｄは、個別領域毎に被写体領域の外郭を幾何形状で近似し、個別領域毎に幾何形状の情報（幾何パラメータＰ１（ｚ）、Ｐ２（ｚ）・・・、（式１３）参照）を取得する。すなわち、被写体領域の外郭が複数の幾何形状で近似され、複数の幾何形状の情報が取得される。

補正倍率算出部２６ｅは、断層画像Ｄｏ（ｚ）毎に、断層画像Ｄｏ（ｚ）に対応する１以上の幾何パラメータＰ（ｚ）に基づいて、断層画像Ｄｏ（ｚ）の各画素（ｘ、ｙ、ｚ）に対応する不透明度の補正倍率Ｓα（ｘ、ｙ、ｚ）を算出し、補正倍率Ｓα（ｘ、ｙ、ｚ）を格納した３次元の補正テーブルＳαを作成する（式１２）（式１５）参照）。

補正倍率Ｓα（ｘ、ｙ、ｚ）は、断層画像Ｄｏ（ｚ）の被写体領域の外郭を近似した幾何形状の幾何中心（後述する実施形態では、楕円の中心）から遠ざかるにつれその倍率が小さくなるように算出される（（式１２）（式１５）参照）。

被写体領域の外郭を複数の幾何形状で近似した場合、複数の幾何形状が互いに重複することがある（図１６参照）。この場合、補正倍率算出部２６ｅは、複数の幾何形状が重複する範囲の補正倍率Ｓα（ｘ、ｙ、ｚ）として、各々の幾何形状の情報に基づいて算出される補正倍率のうち、最も値が大きい補正倍率を割り当てる（（式１５）参照）。

尚、補正倍率Ｓα（ｘ、ｙ、ｚ）は３次元のデータ配列であり、そのままデータテーブル（補正テーブル）として扱えるため、補正倍率と補正テーブルは実質的に同一である。このため、補正倍率と補正テーブルは同一の符号「Ｓα」で表す。

ボクセル作成部２７は、断層画像群Ｄｏ、カラーマップ取得部２４により取得されたカラーマップＣｍａｐ（又はカラーマップ調整部２５により調整されたカラーマップＣｍａｐ）、補正テーブル作成部２６により作成された補正テーブルＳα、及び領域指定部２３により指定された関心領域ＲＯＩに基づいて、ボクセル構造体Ｖ（ボクセルデータ）を作成し、レンダリング部２８に出力する。

具体的には、ボクセル作成部２７は、カラーマップＣｍａｐを参照することで、断層画像群Ｄｏの各画素（ｘ、ｙ、ｚ）の信号値ｖを、信号値ｖに応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセル構造体Ｖ（ボクセルデータ）を作成する。

特に、ボクセル作成部２７は、カラーマップＣｍａｐを参照して断層画像群Ｄｏの各画素（ｘ、ｙ、ｚ）の信号値ｖを色値及び不透明度に変換してボクセル構造体Ｖ（ボクセルデータ）を作成する際、カラーマップＣｍａｐ（ｖ、３）（オパシティカーブ）に定義された不透明度をそのまま用いるのではなく、各画素（ｘ、ｙ、ｚ）に対応する補正倍率Ｓα（ｘ、ｙ、ｚ）を乗算することで不透明度を補正する（（式１７）参照）。

レンダリング部２８は、ボクセル作成部２７により作成されたボクセル構造体Ｖ（ボクセルデータ）に基づいて、ボリュームレンダリング処理を実行し、ボリュームレンダリング像（３次元画像）を生成し表示する。レンダリング部２８は、図５に示すように、第１レンダリング部２８ａ及び第２レンダリング部２８ｂを含む。

第１レンダリング部２８ａは、ＣＰＵによりレンダリング処理を実行するレンダリング部であり、第２レンダリング部２８ｂは、コンピュータグラフィックスＡＰＩのオープン標準規格であるＯｐｅｎＧＬ（Open Graphics Library）を用いてＧＰＵによりレンダリング処理を実行するレンダリング部である。レンダリング処理は、第１レンダリング部２８ａ又は第２レンダリング部２８ｂにより実行される。

データ出力部２９は、ユーザからデータの保存操作を受け付けることによって、各種データを出力して記憶部１２に保存する。例えば、データ出力部２９は、レンダリング部２８により生成し表示されたボリュームレンダリング像（３次元画像）を出力して保存する。

また、データ出力部２９は、ボクセル構造体Ｖ（ボクセルデータ）を作成する際に参照可能なデータ（例えば、補正後の不透明度を格納した不透明度テーブルＴα等）を出力しても良い。これにより、同じ断層画像群Ｄｏに基づいて再度レンダリングを行う際、いくつかの処理（例えば、補正テーブルＳαの作成処理等）を省略することができ、レンダリング像を得るまでの処理時間を短縮できる。

パラメータ調整部３０は、パラメータ調整画面４０（図３０参照）において、幾何パラメータＰ（ｚ）（（式１１）（式１３）参照）や補正倍率Ｓα（（式１２）（式１５）参照）を算出する際の各種パラメータをユーザに調整させる。

次に、図６～図３０を参照しながら、医用画像処理装置１の動作について説明する。
図６は、医用画像処理装置１の全体動作を示すフローチャートである。

まず、医用画像処理装置１の制御部１１（画像取得部２１）は、レンダリング対象となるＤＩＣＯＭ形式の複数の断層画像（断層画像群）Ｄｏを取得する（図６のステップＳ１）。断層画像群Ｄｏの取得方法は任意であるが、例えば、断層画像群Ｄｏが予め記憶部１２に記憶されている場合は、制御部１１は、記憶部１２から断層画像群Ｄｏを読込み取得することができる。

また、断層画像群ＤｏがＵＳＢ等の記憶媒体に記憶されている場合は、制御部１１は、周辺機器Ｉ／Ｆ部１７に接続した当該記憶媒体から断層画像群Ｄｏを読込み取得することができる。或いは、サーバに断層画像群Ｄｏが記憶されている場合は、制御部１１は、通信制御部１４を介して当該サーバから断層画像群Ｄｏをダウンロードして取得することができる。

続いて、制御部１１（階調圧縮部２２）は、必要に応じて、断層画像群Ｄｏに対して階調圧縮処理を施す（図６のステップＳ２）。具体的には、制御部１１は、以下の手順で１６ビットの断層画像群Ｄｏ（（式１）参照）を８ビットの断層画像群Ｄｏ（（式２）参照）に変換する。

（式６）
（１）断層画像群ＤｏのＳｚ／２番目の中間スライスＤｏ（ｚ／２）における全ての画素の最小値Ｄｍｉｎ、最大値Ｄｍａｘを算出する。
（２）下限値Ｌｍｉｎ＝（Ｄｍａｘ－Ｄｍｉｎ）・γ＋Ｄｍｉｎ、上限値Ｌｍａｘ＝（Ｄｍａｘ－Ｄｍｉｎ）・（１－γ）＋Ｄｍｉｎを設定する。ここで、γは階調圧縮画像のコントラスト調整幅で、０に近いほどコントラストが増大する（但し、輝度が小さくなる）。通常はγ＝０．１に設定する。
（３）断層画像群Ｄｏを以下の計算式で８ビットの断層画像群Ｄｏに変換する（下限値Ｌｍｉｎ～上限値Ｌｍａｘの範囲で２５６段階に圧縮する）。
Ｄｏ（ｘ、ｙ、ｚ）
＝（Ｄｏ（ｘ、ｙ、ｚ）－Ｌｍｉｎ）・２５５／（Ｌｍａｘ－Ｌｍｉｎ）
但し、Ｄｏ（ｘ、ｙ、ｚ）＞２５５の場合はＤｏ（ｘ、ｙ、ｚ）＝２５５、Ｄｏ（ｘ、ｙ、ｚ）＜０の場合はＤｏ（ｘ、ｙ、ｚ）＝０に飽和させる。

続いて、制御部１１（領域指定部２３）は、必要に応じて、ユーザから被写体の関心領域ＲＯＩの指定を受け付ける（図６のステップＳ３）。例えば、前記した（式３）のように、関心領域ＲＯＩを直方体等で指定する。

続いて、制御部１１（カラーマップ取得部２４）は、断層画像群Ｄｏに適用するカラーマップＣｍａｐを取得する（ステップＳ４）。ここで、断層画像群Ｄｏの信号値が１６ビットの場合（（式１）参照）は、１６ビット対応のカラーマップＣｍａｐ（（式４）参照）を取得し、断層画像群Ｄｏの信号値が８ビットに圧縮されている場合（（式２）参照）は、８ビット対応のカラーマップＣｍａｐ（（式５）参照）を取得する。

カラーマップＣｍａｐの取得方法は特に限定されず、例えば、制御部１１は、予め用意されている複数のカラーマップＣｍａｐの中から、断層画像群Ｄｏに適用するカラーマップＣｍａｐをユーザに選択させることで取得することができる。また、制御部１１は、前回レンダリング時に使用したカラーマップＣｍａｐを記憶部１２から読込み、断層画像群Ｄｏに適用するカラーマップＣｍａｐとして取得しても良い。

また、制御部１１は、カラーマップ作成画面（不図示）においてユーザが作成したカラーマップＣｍａｐを取得しても良い。この場合、制御部１１は、カラーマップ作成画面において、代表的な信号値ｖ（例えば、１０箇所程度の特徴的な信号値）と当該信号値ｖに対する色値（ＲＧＢ値）及び不透明度（α値）をユーザに設定させ、ユーザが設定した代表的な信号値ｖと色値及び不透明度との対応関係に基づいて所定範囲の信号値ｖを色値及び不透明度に変換するカラーマップＣｍａｐを自動生成する。

続いて、制御部１１は、ステップＳ４において取得したカラーマップＣｍａｐを更に調整する場合（ステップＳ５；「Ｙｅｓ」）、ステップＳ６へ移行し、カラーマップＣｍａｐの調整を行う。具体的には、制御部１１（カラーマップ調整部２５）は、以下のように、カラーマップＣｍａｐ（ｖ、ｎ）のうち、信号値ｖと不透明度との対応関係を定義するカラーマップＣｍａｐ（ｖ、３）（オパシティカーブ）を調整する。

断層画像群Ｄｏの信号値が１６ビットの場合（（式１）参照）は、制御部１１は、次の手順で、１６ビット対応のカラーマップＣｍａｐ（ｖ、３）（（式４）参照）を調整する。

（式７）
（１）断層画像群Ｄｏの中で、Ｓｚ／２番目の中間スライスＤｏ（ｚ／２）における全ての画素の最小値Ｄｍｉｎ、最大値Ｄｍａｘを算出する。
（２）下限値Ｌｍｉｎ＝（Ｄｍａｘ－Ｄｍｉｎ）・γ＋Ｄｍｉｎ、上限値Ｌｍａｘ＝（Ｄｍａｘ－Ｄｍｉｎ）・（１－γ）＋Ｄｍｉｎを設定する。ここで、γは階調圧縮画像のコントラスト調整幅で、０に近いほどコントラストが増大する（但し、輝度が小さくなる）。通常はγ＝０．１に設定する。
（３）１６ビットの信号値ｖを以下の計算式で８ビットの信号値ｖ８に変換する（下限値Ｌｍｉｎ～上限値Ｌｍａｘの範囲で２５６段階に圧縮する）。
ｖ８＝（ｖ－Ｌｍｉｎ）・２５５／（Ｌｍａｘ－Ｌｍｉｎ）
但し、ｖ８＞２５５の場合はｖ８＝２５５、ｖ８＜０の場合はｖ８＝０に飽和させる。
（４）次のように、１６ビット対応のカラーマップＣｍａｐ（ｖ、３）（（式４）参照）を調整する。
ｖ８＞Ｓｖの場合：Ｃｍａｐ（ｖ、３）←Ｃｍａｐ（ｖ、３）（ｖ８－Ｓｖ）^δｖ８≦Ｓｖの場合：調整しない

尚、Ｃｍａｐ（ｖ、ｎ）（０≦Ｃｍａｐ（ｖ、ｎ）≦２５５、－３２７６８≦ｖ≦３２７６７、ｎ＝０（Ｒ）、１（Ｇ）、２（Ｂ）、３（α））であり、８ビット換算の閾値をＳｖとし、減衰係数をδとする（δは負の実数で、通常はδ＝－０．５）。閾値Ｓｖは、不透明度を減衰させる信号値の閾値であり、骨領域を減衰させる場合には、Ｓｖ＝１２８程度に設定する。

一方、断層画像群Ｄｏの信号値が８ビットの場合（（式２）参照）は、制御部１１は、次のように、８ビット対応のカラーマップＣｍａｐ（（式５）参照）を調整する。

（式８）
ｖ＞Ｓｖの場合：Ｃｍａｐ（ｖ、３）←Ｃｍａｐ（ｖ、３）（ｖ－Ｓｖ）^δ
ｖ≦Ｓｖの場合：調整しない

尚、Ｃｍａｐ（ｖ、ｎ）（０≦Ｃｍａｐ（ｖ、ｎ）≦２５５、０≦ｖ≦２５５、ｎ＝０（Ｒ）、１（Ｇ）、２（Ｂ）、３（α））であり、閾値をＳｖとし、減衰係数をδとする（δは負の実数で、通常はδ＝－０．５）。閾値Ｓｖは、不透明度を減衰させる信号値の閾値であり、骨領域を減衰させる場合には、Ｓｖ＝１２８程度に設定する。

図７（ａ）は、（式７）においてＣｍａｐ（ｖ、３）（オパシティカーブ）に乗算される（ｖ８－Ｓｖ）^δの減衰特性を示す図である。図に示すように、信号値ｖ８＞閾値Ｓｖの場合、信号値ｖ８と閾値Ｓｖとの差分値（ｖ８－Ｓｖ）に応じた減衰比率（ｖ８－Ｓｖ）^δで不透明度が減衰する。

例えば、図７（ｂ）に示すように、内臓領域と骨領域の中心からの距離に差が無い場合、両者の不透明度に対して後述のステップＳ７の補正テーブルに基づく補正を加えても、内臓領域は骨領域に被ってしまう。そこで同図に示すように、内臓領域と骨領域の信号値の差を利用し、閾値Ｓｖを内臓領域と骨領域の信号値の境界に設定することで、骨領域の不透明度を効果的に減衰させることができる。

一方、カラーマップＣｍａｐを調整しない場合には（ステップＳ５；「Ｎｏ」）、制御部１１は、上記したカラーマップＣｍａｐの調整を省略し、ステップＳ７へ移行する。

図８～図１７を参照して、ステップＳ７における補正テーブルＳαを作成する処理について説明する。

図８は、補正テーブルＳαを作成する処理を示すフローチャートである。
まず、制御部１１（被写体領域抽出部２６ａ）は、断層画像Ｄｏ（ｚ）（０≦ｚ≦Ｓｚ－１）毎に、断層画像Ｄｏ（ｚ）の被写体領域を抽出する（図８のステップＳ２１）。

図９は、ある断層画像Ｄｏ（ｚ）の被写体領域を抽出する処理を説明する図である。図９（ａ）に示すように、制御部１１は、断層画像Ｄｏ（ｚ）に対して、信号値が所定の閾値以上の画素を「白」、閾値未満の画素を「黒」で塗りつぶして２値化する。閾値としては、例えば、空気以外を被写体領域として抽出する場合には、信号値ｖが１６ビットの場合は「－７００」、信号値ｖが８ビットの場合には「１０」などに設定する。

そして、制御部１１は、２値化された断層画像Ｄｏ（ｚ）の「白」の領域を被写体領域として抽出する。具体的には、図９（ｂ）に示すように、被写体領域（白領域）の中心座標ｃ（ｚ）、被写体領域の横方向のサイズｗ（ｚ）、被写体領域の縦方向のサイズｈ（ｚ）を以下の（式９）で算出する。尚、図９（ｂ）の矩形領域Ｒ（ｚ）は、被写体領域を全て含む最小の矩形であり、その最小座標を（Ｘｍｉｎ、Ｙｍｉｎ）、最大座標を（Ｘｍａｘ、Ｙｍａｘ）とする。

（式９）
被写体領域の中心座標ｃ（ｚ）＝（ｃｘ（ｚ）、ｃｙ（ｚ））
ｃｘ（ｚ）＝（Ｘｍａｘ＋Ｘｍｉｎ）／２
ｃｙ（ｚ）＝（Ｙｍａｘ＋Ｙｍｉｎ）／２
被写体領域の横方向のサイズｗ（ｚ）＝Ｘｍａｘ－Ｘｍｉｎ＋１
被写体領域の縦方向のサイズｈ（ｚ）＝Ｙｍａｘ－Ｙｍｉｎ＋１

続いて、制御部１１（被写体領域補正部２６ｂ）は、断層画像Ｄｏ（ｚ）（０≦ｚ≦Ｓｚ－１）毎に、ステップＳ２１により抽出された被写体領域の抽出精度を所定の基準で判断し、所定の基準を満たさない（抽出精度が良好でない）断層画像Ｄｏ（ｚ）の被写体領域を、所定の基準を満たす（抽出精度が良好である）他の被写体領域の情報に基づいて補正する（図８のステップＳ２２）。

図１０は、断層画像Ｄｏ（ｚ）の被写体領域の補正について説明する図である。図１０（ａ）の断層画像Ｄｏ（ｚ）は、被写体領域の抽出が良好に行えない断層画像の例である。図１０（ａ）の断層画像は、被写体領域（幅方向）の両端が切れており、被写体領域が画像範囲に収まっていない（撮影時の条件によってこのような画像が得られる場合がある）。このような場合、被写体領域を正確に抽出できない。

すなわち、図１０（ｂ）に示すように、断層画像から得られる被写体領域の横方向のサイズｗ（ｚ）が、断層画像の横方向のサイズＳｘと同一となり（ｗ（ｚ）＝Ｓｘ）、実際の被写体領域の幅（Ｓｘより大きい幅）を正確に反映したものでなくなる。このため、本実施の形態では、被写体領域の抽出精度を、被写体領域の横方向のサイズｗ（ｚ）がｗ（ｚ）＜Ｓｘを満たすか否かで判断し、ｗ（ｚ）＜Ｓｘを満たさない（つまり、ｗ（ｚ）＝Ｓｘとなる）場合には、被写体領域の横方向のサイズｗ（ｚ）の補正を行う。

具体的には、制御部１１は、まず、ステップＳ２１により抽出された被写体領域の横方向のサイズｗ（ｚ）がｗ（ｚ）＝Ｓｘとなる１以上の断層画像Ｄｏ（ｚ）を、被写体領域の両端が切れた断層画像（被写体領域の抽出精度が良好でない断層画像）として特定する。

次に、制御部１１は、特定された断層画像Ｄｏ（ｚ）毎に、断層画像Ｄｏ（ｚ）とスライス位置が最も近く、かつ、ステップＳ２１により抽出された被写体領域の横方向のサイズｗ（ｚ）がｗ（ｚ）＜Ｓｘを満たす断層画像（「断層画像Ｄｏ（ｚ’）」と表記）を特定する。

断層画像Ｄｏ（ｚ’）は、抽出された被写体領域の横方向のサイズｗ（ｚ’）が画像サイズＳｘより小さいため、被写体領域（幅方向）が画像範囲に収まっている断層画像（被写体領域の抽出精度が良好な断層画像）である。そして、制御部１１は、被写体領域の両端が切れた断層画像Ｄｏ（ｚ）から抽出された被写体領域の横方向のサイズｗ（ｚ）を、断層画像Ｄｏ（ｚ’） から抽出された被写体領域の横方向のサイズｗ（ｚ’）及び縦方向のサイズｈ（ｚ’）の比と同一となるように、以下の式で補正する。

（式１０）
ｗ（ｚ）＝ｈ（ｚ）・ｗ（ｚ’）／ｈ（ｚ’）

これにより、図１０（ｃ）に示すように、被写体領域の横方向のサイズｗ（ｚ）が、真の被写体領域の横方向のサイズｗ（ｚ）に補正され、被写体領域（幅方向）が画像範囲に収まらない場合でも、被写体領域を精度よく求めることができる。

図８のフローチャートの説明に戻る。
続いて、抽出された被写体領域に複数の個別領域を設定しない場合（ステップＳ２３；Ｎｏ）には、制御部１１(幾何情報取得部２６ｄ)は、断層画像Ｄｏ（ｚ）毎に、断層画像Ｄｏ（ｚ）から抽出された被写体領域の外郭を単一の幾何形状で近似し、幾何パラメータＰ（ｚ）を取得する（図８のステップＳ２４）。本実施の形態では、被写体形状の外郭を楕円Ｅ（ｚ）で近似し、幾何パラメータＰ（ｚ）として楕円Ｅ（ｚ）のパラメータ（楕円パラメータＰ（ｚ））を取得するものとする。

具体的には、図１１に示すように、制御部１１は、ステップＳ２１において抽出した被写体領域（またはステップＳ２２において補正された被写体領域）を囲う矩形領域Ｒ（ｚ）に内接する楕円Ｅ（ｚ）を求め、この楕円Ｅ（ｚ）を被写体形状の外郭を近似する楕円として求める。そして、制御部１１は、楕円Ｅ（ｚ）の楕円パラメータＰ（ｚ）として、以下の（式１１）に示すように、楕円Ｅ（ｚ）の中心座標Ｃ（ｚ）、横方向のサイズＷ（ｚ）、及び縦方向のサイズＨ（ｚ）を算出する。

（式１１）
楕円Ｅ（ｚ）の中心座標Ｃ（ｚ）＝（Ｃｘ（ｚ）、Ｃｙ（ｚ））
Ｃｘ（ｚ）＝ｃｘ（ｚ）＋ｘｏｆｆｓｅｔ
Ｃｙ（ｚ）＝ｃｙ（ｚ）＋ｙｏｆｆｓｅｔ
楕円Ｅ（ｚ）の横方向のサイズＷ（ｚ）＝ａ・ｗ（ｚ）
楕円Ｅ（ｚ）の縦方向のサイズＨ（ｚ）＝ｂ・ｈ（ｚ）

（式１１）において、ａ、ｂは楕円Ｅ（ｚ）の横方向のサイズ、縦方向のサイズの補正係数（倍率）である。補正しない場合はａ＝ｂ＝１．０に設定する。ｘｏｆｆｓｅｔ、ｙｏｆｆｓｅｔはＸＹ面における楕円のオフセット値であり、初期値はｘｏｆｆｓｅｔ＝ｙｏｆｆｓｅｔ＝０である。

制御部１１は、全ての断層画像Ｄｏ（ｚ）（０≦ｚ≦Ｓｚ－１）について、被写体領域の外郭を近似する単一の楕円Ｅ（ｚ）を求め、（式１１）の楕円パラメータＰ（ｚ）を算出する。

そして、制御部１１（補正倍率算出部２６ｅ）は、断層画像Ｄｏ（ｚ）毎に、断層画像Ｄｏ（ｚ）に対応する楕円パラメータＰ（ｚ）に基づいて、断層画像Ｄｏ（ｚ）の各画素（ｘ、ｙ、ｚ）に対応する不透明度の補正倍率Ｓα（ｘ、ｙ、ｚ）を算出し、補正テーブルＳαを作成する（図８のステップＳ２５）。
例えば、ｚ番目の断層画像Ｄｏ（ｚ）の各画素（ｘ、ｙ、ｚ）の補正倍率Ｓα（ｘ、ｙ、ｚ）は、以下のように算出される。

（式１２）
Ｓα（ｘ、ｙ、ｚ）＝１－ｋ・ｒ（ｘ、ｙ、ｚ）・ａｍｐ
ｒ（ｘ、ｙ、ｚ）＝｛（ｘ－Ｃｘ（ｚ））／（Ｗ（ｚ）／２）｝^２＋｛（ｙ－Ｃｙ（ｚ））／（Ｈ（ｚ）／２）｝^２
但し、０≦Ｓα（ｘ、ｙ、ｚ）≦１

図１２は、（式１２）により算出される補正倍率Ｓα（ｘ、ｙ、ｚ）の２次元分布パターン（同心楕円分布）を示す図である。図に示すように、楕円Ｅ（ｚ）の中心Ｃ（ｚ）の補正倍率Ｓαは１（最大）となり、楕円Ｅ（ｚ）の内部は中心Ｃ（ｚ）から遠ざかるにつれ補正倍率Ｓαが小さくなる。具体的には、楕円Ｅ（ｚ）の中心Ｃ（ｚ）と横方向のサイズＷ（ｚ）及び縦方向のサイズＨ（ｚ）の比が同一となる各同心楕円上の補正倍率Ｓαが、各同心楕円の径（横方向のサイズ及び縦方向のサイズ）が大きい程小さくなる（径の２乗に反比例して小さくなる）。楕円Ｅ（ｚ）の外部の補正倍率Ｓαは一律に０とする。

尚、（式１２）において、ｋはＸＹ平面の減衰係数であり、初期値はｋ＝１．０である。また、（式１２）のａｍｐはＸＹ面の振幅倍率であり、初期値はａｍｐ＝１．０である。図１３に示すように、被写体中心部の補正倍率を平坦にし、飽和領域を広くしたい場合は、ａｍｐ＞１．０に設定する（但し、陰影感が弱くなる）。

制御部１１は、全ての断層画像Ｄｏ（ｚ）について、補正倍率Ｓα（ｘ、ｙ、ｚ）（０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｚ≦Ｓｚ－１）を算出し、算出された補正倍率Ｓα（ｘ、ｙ、ｚ）を格納する３次元の補正テーブルＳαを作成する。尚、前述したように、補正倍率Ｓα（ｘ、ｙ、ｚ）は３次元のデータ配列であり、そのままデータテーブル（補正テーブル）として扱えるため、補正倍率と補正テーブルは実質的に同一である。

一方、被写体領域に複数の個別領域を設定する場合（ステップＳ２３；Ｙｅｓ）には、まず、制御部１１（個別領域設定部２６ｃ）は、断層画像Ｄｏ（ｚ）毎に、ステップＳ２１において抽出した被写体領域（またはステップＳ２２において補正された被写体領域）に複数の個別領域を設定する（図８のステップＳ２６）。

本実施の形態では、図１４（ａ）に示すように、制御部１１は、被写体領域を横方向に２分割（本実施の形態では等分割）した個別領域Ｒ１（ｚ）、Ｒ２（ｚ）と、被写体領域の横方向の中央付近に個別領域Ｒ１（ｚ）、Ｒ２（ｚ）より小さい個別領域Ｒ３（ｚ）を設定するものとする。個別領域Ｒ３（ｚ）は、少なくとも被写体領域の中心座標ｃ（ｚ）（（式９）参照）を含む領域である。尚、図１４（ａ）のように、個別領域同士に一部で重なりあう領域があっても良い。

図１４（ａ）は、特に、胴体部を対象とした場合の有効な個別領域の設定方法を示している。一般的な人の胸郭はＸ軸方向に２つの肺野が存在し、ほぼ中央に心臓、胸骨、椎骨が存在することから、図１４（ａ）のように、３つの個別領域を設定することが望ましい。但し、個別領域の設定数および設定範囲は被写体領域における観察対象や非観察対象の位置関係などに応じて適宜設定されるものであり、図１４（ａ）の例に限定されるものではない。

続いて、制御部１１（幾何情報取得部２６ｄ)は、断層画像Ｄｏ（ｚ）毎に、被写体領域の外郭を複数の楕円で近似し、複数の楕円パラメータＰ（ｚ）を取得する。具体的には、図１４（ｂ）に示すように、３つの個別領域Ｒ１（ｚ）、Ｒ２（ｚ）、Ｒ３（ｚ）毎に、被写体領域の外郭を楕円Ｅ１（ｚ）、Ｅ２（ｚ）、Ｅ３（ｚ）で近似し、以下の（式１３）に示すように、楕円パラメータＰ１（ｚ）、Ｐ２（ｚ）、Ｐ３（ｚ）を算出する（図８のステップＳ２７）。（式１３）によれば、楕円Ｅ１（ｚ）、Ｅ２（ｚ）は合同な楕円となり、楕円Ｅ３（ｚ）は楕円Ｅ１（ｚ）、Ｅ２（ｚ）より小さい楕円となる。

（式１３）
＜楕円パラメータＰ１（ｚ）＞
楕円Ｅ１（ｚ）の中心座標Ｃ１（ｚ）＝（Ｃｘ１（ｚ）、Ｃｙ１（ｚ））
Ｃｘ１（ｚ）＝ｃｘ（ｚ）－ｗ（ｚ）／４＋ｘｏｆｆｓｅｔ
Ｃｙ１（ｚ）＝ｃｙ（ｚ）＋ｙｏｆｆｓｅｔ＋ｄｙ（ｚ）
楕円Ｅ１（ｚ）の横方向のサイズＷ１（ｚ）＝ｗ（ｚ）・ａ／２
楕円Ｅ１（ｚ）の縦方向のサイズＨ１（ｚ）＝ｈ（ｚ）・ｂ／２
＜楕円パラメータＰ２（ｚ）＞
楕円Ｅ２（ｚ）の中心座標Ｃ２（ｚ）＝（Ｃｘ２（ｚ）、Ｃｙ２（ｚ））
Ｃｘ２（ｚ）＝ｃｘ（ｚ）＋ｗ（ｚ）／４＋ｘｏｆｆｓｅｔ
Ｃｙ２（ｚ）＝ｃｙ（ｚ）＋ｙｏｆｆｓｅｔ＋ｄｙ（ｚ）
楕円Ｅ２（ｚ）の横方向のサイズＷ２（ｚ）＝ｗ（ｚ）・ａ／２
楕円Ｅ２（ｚ）の縦方向のサイズＨ２（ｚ）＝ｈ（ｚ）・ｂ／２
＜楕円パラメータＰ３（ｚ）＞
楕円Ｅ３（ｚ）の中心座標Ｃ３（ｚ）＝（Ｃｘ３（ｚ）、Ｃｙ３（ｚ））
Ｃｘ３（ｚ）＝ｃｘ（ｚ）＋ｘｏｆｆｓｅｔ
Ｃｙ３（ｚ）＝ｃｙ（ｚ）＋ｈ（ｚ）・ｂ・ｄ／２＋ｙｏｆｆｓｅｔ＋ｄｙ（ｚ）
楕円Ｅ３（ｚ）の横方向のサイズＷ３（ｚ）＝ｗ（ｚ）・ａ・ｅ／２
楕円Ｅ３（ｚ）の縦方向のサイズＨ３（ｚ）＝ｈ（ｚ）・ｂ・ｅ／２

（式１３）において、ａ、ｂは楕円Ｅ１（ｚ）～Ｅ３（ｚ）の横方向のサイズ、縦方向のサイズの補正係数（倍率）であり、各楕円に共通するパラメータである。補正しない場合はａ＝ｂ＝１．０に設定する。ｘｏｆｆｓｅｔ、ｙｏｆｆｓｅｔはＸＹ面におけるオフセット値であり、各楕円に共通するパラメータである（初期値はｘｏｆｆｓｅｔ＝ｙｏｆｆｓｅｔ＝０）。

また（式１３）において、ｄ、ｅは中央の楕円Ｅ３（ｚ）固有のパラメータであり、それぞれ、楕円Ｅ３（ｚ）のＹ方向のオフセット値と楕円Ｅ３（ｚ）のサイズを決めるパラメータである。ｄには、楕円Ｅ１（ｚ）、Ｅ２（ｚ）の縦方向のサイズ（ｈ（ｚ）・ｂ／２）に対する比率が設定され、ｅには、楕円Ｅ１（ｚ）、Ｅ２（ｚ）の縦方向のサイズ（ｈ（ｚ）・ｂ／２）および横方向のサイズ（ｗ（ｚ）・ａ／２）に対する共通の比率が設定される（即ち、０≦ｄ、ｅ≦１）。例えば、胸部ＣＴ画像の場合には、観察対象となる心臓が良好に観察できるように、パラメータｄ、ｅによって中央の楕円Ｅ３（ｚ）の位置およびサイズをチューニングする。

また（式１３）において、ｄｙ（ｚ）は、スライス位置ｚに応じた各楕円の中心座標（Ｙ座標）の補正量（オフセット値）であり、各楕円に共通するパラメータである。ｄｙ（ｚ）は、以下のように定義される。

（式１４）
（１）Ｚｔｈ＞０かつｚ＜Ｚｔｈの場合
ｄｙ（ｚ）＝０
（２）ｚ≧Ｚｔｈの場合
ｄｙ（ｚ）＝（ｚ－Ｚｔｈ）・ｄｙｍ／（Ｓｚ－１－Ｚｔｈ）
（３）Ｚｔｈ＜０かつｚ＜－Ｚｔｈの場合
ｄｙ（ｚ）＝（ｚ＋Ｚｔｈ）・ｄｙｍ／Ｚｔｈ
（４）ｚ≧－Ｚｔｈの場合
ｄｙ（ｚ）＝０
尚、補正しない場合はＺｔｈ＝ｄｙｍ＝０とする。

図１５は、上記した補正量ｄｙ（ｚ）の関数と、胴体部のサジタル像（ＭＰＲ像）を表示した図である。サジタル像から分かるように、椎骨は通常、頭を支えるためにＺ軸方向（体軸方向）に向かって湾曲している（Ｓ字カーブを描く）。胸部から頭部側（ｚ＜Ｚｔｈ）の椎骨の湾曲については、被写体領域も椎骨の湾曲と連動して変化するため、被写体領域内における椎骨の相対的な位置はさほど変化しない。一方、胸部下（ｚ≧Ｚｔｈ）から腰部側の椎骨の湾曲については、被写体領域の変化以上に椎骨の位置が体表面より内側（Ｙ軸方向）にずれるため、被写体領域内における椎骨の相対的な位置が大きく変わる。

このため、胸部下（ｚ≧Ｚｔｈ）では、椎骨のＹ軸方向への位置ずれに対処するため、（式１３）に示したように、各楕円Ｅ１（ｚ）～Ｅ３（ｚ）の各中心座標Ｃｙ１（ｚ）～Ｃｙ３（ｚ）に対して、（式１４）で定義されるｄｙ（ｚ）のオフセットを加えることが望ましい。これにより、スライス位置ｚが椎骨のＹ軸方向への位置ずれが大きくなる分岐点である閾値Ｚｔｈ以上になると、各楕円の中心座標（Ｙ座標）に対してｚ－Ｚｔｈに比例したオフセットが加えられるため、椎骨のＹ軸方向への位置ずれに対処でき、胸部下で椎骨領域が透過されずに残ってしまうことを防ぐことができる。尚、椎骨の湾曲度合は個人差があるため、ＣＴ画像ごとにＺｔｈ、ｄｙｍを調整する必要があるが、サジタル像（ＭＰＲ像）に基づいてこれらのパラメータを適宜調整することも可能である。

制御部１１は、全ての断層画像Ｄｏ（ｚ）（０≦ｚ≦Ｓｚ－１）について、被写体領域の外郭を近似する複数の楕円Ｅ１（ｚ）～Ｅ３（ｚ）を求め、（式１３）の楕円パラメータＰ１（ｚ）～Ｐ３（ｚ）を算出する。

そして、制御部１１（補正倍率算出部２６ｅ）は、断層画像Ｄｏ（ｚ）毎に、断層画像Ｄｏ（ｚ）に対応する複数の楕円パラメータＰ１（ｚ）～Ｐ３（ｚ）に基づいて、断層画像Ｄｏ（ｚ）の各画素（ｘ、ｙ、ｚ）に対応する不透明度の補正倍率Ｓα（ｘ、ｙ、ｚ）を算出し、補正テーブルＳαを作成する（図８のステップＳ２９）。
例えば、ｚ番目の断層画像Ｄｏ（ｚ）の各画素（ｘ、ｙ、ｚ）の補正倍率Ｓα（ｘ、ｙ、ｚ）は、以下のように算出される。

（式１５）
Ｓα（ｘ、ｙ、ｚ）＝１－ｋ・ｒ（ｘ、ｙ、ｚ）・ａｍｐ
ｒ（ｘ、ｙ、ｚ）＝ｍｉｎ［ｒ１（ｘ、ｙ、ｚ）、ｒ２（ｘ、ｙ、ｚ）、ｒ３（ｘ、ｙ、ｚ）］
ｒ１（ｘ、ｙ、ｚ）＝｛（ｘ－Ｃｘ１（ｚ））／（Ｗ１（ｚ）／２）｝^２＋｛（ｙ－Ｃｙ１（ｚ）））／（Ｈ１（ｚ）／２）｝^２
ｒ２（ｘ、ｙ、ｚ）＝｛（ｘ－Ｃｘ２（ｚ））／（Ｗ２（ｚ）／２）｝^２＋｛（ｙ－Ｃｙ２（ｚ））／（Ｈ２（ｚ）／２）｝^２
ｒ３（ｘ、ｙ、ｚ）＝｛（ｘ－Ｃｘ３（ｚ））／（Ｗ３（ｚ）／２）｝^２＋｛（ｙ－Ｃｙ３（ｚ）））／（Ｈ３（ｚ）／２）｝^２

図１６は、（式１５）により算出される補正倍率Ｓα（ｘ、ｙ、ｚ）の２次元分布パターン（複数の同心楕円分布）を示す図である。図に示すように、各楕円Ｅ１（ｚ）～Ｅ３（ｚ）の中心Ｃ１（ｚ）～Ｃ３（ｚ）の補正倍率Ｓαは１（最大）となり、各楕円Ｅ１（ｚ）～Ｅ３（ｚ）の内部は各楕円の中心Ｃ１（ｚ）～Ｃ３（ｚ）から遠ざかるにつれ補正倍率Ｓαが小さくなる。各楕円Ｅ１（ｚ）～Ｅ３（ｚ）の外部の補正倍率Ｓαは一律に０とする。

また（式１５）に示すように、Ｓα（ｘ、ｙ、ｚ）＝１－ｋ・ｒ（ｘ、ｙ、ｚ）・ａｍｐにおけるｒ（ｘ、ｙ、ｚ）として、各楕円Ｅ１（ｚ）～Ｅ３（ｚ）の楕円方程式ｒ１（ｘ、ｙ、ｚ）～ｒ３（ｘ、ｙ、ｚ）の値の最小値（ｒ（ｘ、ｙ、ｚ）＝ｍｉｎ［ｒ１（ｘ、ｙ、ｚ）、ｒ２（ｘ、ｙ、ｚ）、ｒ３（ｘ、ｙ、ｚ）］）を与える。すなわち、楕円同士が重複する範囲の補正倍率Ｓαには、各楕円から算出される補正倍率のうち、最も値が大きい補正倍率が割り当てられる。

尚、（式１５）において、ｋはＸＹ平面の減衰係数であり、初期値はｋ＝１．０である。
また、（式１５）のａｍｐは、（式１２）と同様、ＸＹ面の振幅倍率であり、初期値はａｍｐ＝１．０である。図１３に示したように、被写体中心部の補正倍率を平坦にし、飽和領域を広くしたい場合は、ａｍｐ＞１．０に設定する（但し、陰影感が弱くなる）。

制御部１１は、全ての断層画像Ｄｏ（ｚ）について、補正倍率Ｓα（ｘ、ｙ、ｚ）（０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｚ≦Ｓｚ－１）を算出し、算出された補正倍率Ｓα（ｘ、ｙ、ｚ）を格納する３次元の補正テーブルＳαを作成する。図１７は、実際の胸腹部ＣＴ画像を対象として、３次元の補正テーブルＳα（補正倍率）を作成し、可視化（３次元描画）したものである。例えば、制御部１１は、図１７のように作成した補正テーブルＳαを表示部１６に表示させることで、補正テーブルＳαが適切に作成されているか、ユーザに確認させてもよい。

以上、図８～図１７を参照して、補正テーブルＳαを作成する処理について説明した。
図６のフローチャートの説明に戻る。

制御部１１（ボクセル作成部２７）は、ステップＳ１において取得された断層画像群Ｄｏ（又はステップＳ２において階調圧縮された断層画像群Ｄｏ）、ステップＳ４において取得されたカラーマップＣｍａｐ（又はステップＳ６において調整されたカラーマップＣｍａｐ）、ステップＳ７において作成された補正テーブルＳα、及びステップＳ３において指定された関心領域ＲＯＩに基づいて、ボクセル構造体Ｖ（ボクセルデータ）を作成する（図６のステップＳ８）。

本実施の形態では、ボクセル構造体Ｖ（ボクセルデータ）は、２４ビットの色値（ＲＧＢ値）を保持するボクセル構造体Ｖｃと、８ビットの不透明度（α値）を保持するボクセル構造体Ｖαとの２つのボクセル構造体から構成されるものとする。このようにボクセル構造体を２つに分けて作成することで、色値及び不透明度に関するボクセル処理を個別に実行できるため、メモリアクセス等の向上が図られる。

図１８は、ボクセル構造体Ｖ（Ｖｃ、Ｖα）を作成する処理を示すフローチャートである。
まず、制御部１１は、カラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦２）（＝カラーパレット）、及び関心領域ＲＯＩに基づいて、２４ビットのＲＧＢ値を保持するボクセル構造体Ｖｃ（ｘ、ｙ、ｚ、ｎ）（０≦Ｖ（ｘ、ｙ、ｚ、ｎ）≦２５５、０≦ｎ≦２、０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｚ≦Ｓｚ－１；解像度：Ｒｘｙ、Ｒｚ）を以下のように作成する（図１８のステップＳ８０１）。

（式１６）
（１）Ｘｓ＜ｘ＜Ｘｅ、Ｙｓ＜ｙ＜Ｙｅ、及びＺｓ＜ｚ＜Ｚｅを全て満たす場合（関心領域ＲＯＩ内の場合）
Ｖｃ（ｘ、ｙ、ｚ、ｎ）＝Ｃｍａｐ（Ｄｏ（ｘ、ｙ、ｚ）、ｎ） （０≦ｎ≦２）
（２）ｘ≦Ｘｓ、ｘ≧Ｘｅ、ｙ≦Ｙｓ、ｙ≧Ｙｅ、ｚ≦Ｚｓ、又はｚ≧Ｚｅのいずれかを満たす場合（関心領域ＲＯＩ外の場合）
Ｖｃ（ｘ、ｙ、ｚ、ｎ）＝０ （０≦ｎ≦２）

次に、制御部１１は、カラーマップＣｍａｐ（ｖ、３）（＝オパシティカーブ）、補正テーブルＳα、及び関心領域ＲＯＩに基づいて、８ビットの不透明度を保持するボクセル構造体Ｖα（ｘ、ｙ、ｚ）（０≦Ｖα（ｘ、ｙ、ｚ）≦２５５、０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｚ≦Ｓｚ－１；解像度：Ｒｘｙ、Ｒｚ）を以下のように作成する（図１８のステップＳ８０２）。

（式１７）
（１）Ｘｓ＜ｘ＜Ｘｅ、Ｙｓ＜ｙ＜Ｙｅ、及びＺｓ＜ｚ＜Ｚｅを全て満たす場合（関心領域ＲＯＩ内の場合）
Ｖα（ｘ、ｙ、ｚ）＝Ｃｍａｐ（Ｄｏ（ｘ、ｙ、ｚ）、３）
Ｖα（ｘ、ｙ、ｚ）＝Ｖα（ｘ、ｙ、ｚ）・Ｓα（ｘ、ｙ、ｚ）
（２）ｘ≦Ｘｓ、ｘ≧Ｘｅ、ｙ≦Ｙｓ、ｙ≧Ｙｅ、ｚ≦Ｚｓ、又はｚ≧Ｚｅのいずれかを満たす場合（関心領域ＲＯＩ外の場合）
Ｖα（ｘ、ｙ、ｚ）＝０

（式１７）のＶα（ｘ、ｙ、ｚ）＝Ｖα（ｘ、ｙ、ｚ）・Ｓα（ｘ、ｙ、ｚ）により、カラーマップＣｍａｐに定義された不透明度に対して補正倍率Ｓαが乗算され、不透明度が補正される。

尚、メモリアクセスの向上等を目的に、３２ビットのボクセル構造体Ｖを色値を保持するボクセル構造体Ｖｃと不透明度を保持するボクセル構造体Ｖαとに分けているが、ボクセル構造体を分けなくてもよい。この場合、制御部１１は、以下のように、３２ビットの色値及び不透明度を保持する１つのボクセル構造体Ｖ（ｘ、ｙ、ｚ、ｎ）（０≦Ｖ（ｘ、ｙ、ｚ、ｎ）≦２５５、０≦ｎ≦３、０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｚ≦Ｓｚ－１；解像度：Ｒｘｙ、Ｒｚ）を作成する。

（式１８）
（１）Ｘｓ＜ｘ＜Ｘｅ、Ｙｓ＜ｙ＜Ｙｅ、及びＺｓ＜ｚ＜Ｚｅを全て満たす場合（関心領域ＲＯＩ内の場合）
Ｖ（ｘ、ｙ、ｚ、ｎ）＝Ｃｍａｐ（Ｄｏ（ｘ、ｙ、ｚ）、ｎ） （０≦ｎ≦３）
Ｖ（ｘ、ｙ、ｚ、３）＝Ｖ（ｘ、ｙ、ｚ、３）・Ｓα（ｘ、ｙ、ｚ）
（２）ｘ≦Ｘｓ、ｘ≧Ｘｅ、ｙ≦Ｙｓ、ｙ≧Ｙｅ、ｚ≦Ｚｓ、又はｚ≧Ｚｅのいずれかを満たす場合（関心領域ＲＯＩ外の場合）
Ｖ（ｘ、ｙ、ｚ、ｎ）＝０ （０≦ｎ≦３）

上記した（式１６）～（式１８）に示すボクセル構造体を作成する処理は、断層画像群Ｄｏの信号値が１６ビット（（式１）参照）の場合は、１６ビット対応のカラーマップＣｍａp（（式４）参照）を用いて実行され、断層画像群Ｄｏの信号値が８ビット（（式２）参照）に階調圧縮されている場合は、８ビット対応のカラーマップＣｍａp（（式５）参照）を用いて実行される。

制御部１１は、ボクセル構造体Ｖ（Ｖｃ、Ｖα）を作成した後、必要に応じて以下のように、陰影処理を施しても良い。

まず、制御部１１は、平行光源の光源ベクトル（Ｌｘ、Ｌｙ、Ｌｚ）（単位ベクトル）を設定する。例えば、（Ｌｘ、Ｌｙ、Ｌｚ）＝（０．５７７３５、０．５７７３５、０．５７７３５）と設定する。また、制御部１１は、環境光成分Ａｂ（０≦Ａｂ≦１、例えばＡｂ＝０．２）を設定する。

ここで、後述するステップＳ９のボリュームレンダリング処理を、ＣＰＵ（第１レンダリング部２８ａ）により実行する場合、制御部１１は、Ｘｓ≦ｘ≦Ｘｅ、Ｙｓ≦ｙ≦Ｙｅ、Ｚｓ≦ｚ≦Ｚｅの範囲のボクセル（ｘ、ｙ、ｚ）の勾配ベクトル（Ｇｘ、Ｇｙ、Ｇｚ）を以下の（式１９）で算出する（エッジ部も陰影計算を行う）。但し、ｘ＋１＞ＸｅのときＶα（ｘ＋１、ｙ、ｚ）＝０、ｘ－１＜ＸｓのときＶα（ｘ－１、ｙ、ｚ）＝０、ｙ＋１＞ＹｅのときＶα（ｘ、ｙ＋１、ｚ）＝０、ｙ－１＜ＹｓのときＶα（ｘ、ｙ－１、ｚ）＝０、ｚ＋１＜ＺｓのときＶα（ｘ、ｙ、ｚ＋１）＝０、ｚ－１＜ＺｓのときＶα（ｘ、ｙ、ｚ－１）＝０とする。

一方、後述するステップＳ９のボリュームレンダリング処理を、ＧＰＵ（第２レンダリング部２８ｂ）により実行する場合、制御部１１は、Ｘｓ＜ｘ＜Ｘｅ、Ｙｓ＜ｙ＜Ｙｅ、Ｚｓ＜ｚ＜Ｚｅの範囲のボクセル（ｘ、ｙ、ｚ）の勾配ベクトル（Ｇｘ、Ｇｙ、Ｇｚ）を以下の（式１９）で算出する（エッジ部の陰影計算は行わず、後述の色補正処理（図２６参照）を行う）。

（式１９）
Ｇｘ＝｛Ｖα（ｘ＋１、ｙ、ｚ）－Ｖα（ｘ－１、ｙ、ｚ）｝・Ｒｘｙ／Ｒｚ
Ｇｙ＝｛Ｖα（ｘ、ｙ＋１、ｚ）－Ｖα（ｘ、ｙ－１、ｚ）｝・Ｒｘｙ／Ｒｚ
Ｇｚ＝｛Ｖα（ｘ、ｙ、ｚ＋１）－Ｖα（ｘ、ｙ、ｚ－１）｝
Ｇ＝｛Ｇｘ^２＋Ｇｙ^２＋Ｇｚ^２｝^１／２
尚、Ｒｘｙ／Ｒｚはあらかじめ定義されるＺ軸方向変倍率である。（式１９）のように、Ｇｘ、ＧｙにＺ軸方向変倍率Ｒｘｙ／Ｒｚを乗算する代わりに、ＧｚにＺ軸方向変倍率Ｒｘｙ／Ｒｚの逆数を乗算するのでもよい。

続いて、Ｇ≧１の場合、制御部１１は、輝度値（陰影値）Ｓ（ｘ、ｙ、ｚ）（０≦Ｓ（ｘ、ｙ、ｚ）≦１）として、拡散反射成分のみ算出し、Ｓ（ｘ、ｙ、ｚ）＝（１－Ａｂ）（Ｇｘ・Ｌｘ＋Ｇｙ・Ｌｙ＋Ｇｚ・Ｌｚ）／Ｇ＋Ａｂを与える。Ｇ＜１の場合、制御部１１は、輝度値（陰影値）Ｓ（ｘ、ｙ、ｚ）として、Ｓ（ｘ、ｙ、ｚ）＝０を与える（視線を変えても勾配ベクトルは変わらないため、鏡面反射成分を加えることはできない）。

そして、制御部１１は、以下の（式２０）により、色値を保持するボクセル構造体Ｖｃの成分を改変する。

（式２０）
Ｖｃ（ｘ、ｙ、ｚ、ｎ）＝Ｓ（ｘ、ｙ、ｚ）・Ｖｃ（ｘ、ｙ、ｚ、ｎ）（０≦ｎ≦２）

以上のように、制御部１１は、必要に応じて陰影処理を施しても良い。
図６のフローチャートの説明に戻る。

続いて、制御部１１（第１レンダリング部２８ａ又は第２レンダリング部２８ｂ）は、ボクセル構造体Ｖに基づいて、ボリュームレンダリング像を生成する（図６のステップＳ９）。以降、ＣＰＵ（第１レンダリング部２８ａ）により実行されるレイキャスティング法によるレンダリング処理と、ＧＰＵ（第２レンダリング部２８ｂ）により実行される３Ｄテクスチャマッピング法によるレンダリング処理について順に説明する。

最初に、図１９～図２４を参照しながら、レイキャスティング法によるレンダリング処理について説明する。
図１９は、レイキャスティング法の概要を示す図である。図１９（ａ）に示すように、投影面（レンダリング像）の各点からボクセルデータに対して任意の方向にレイを発射し、最深点（レイ方向で光が届く最も深い点）よりその軌跡を逆にたどり、投影面上の点の輝度値を算出する。
また、図１９（ｂ）に示すように、計算を簡単にするため、ボクセルデータに座標変換を施し、投影面の各点よりＺ軸負方向の最深点までレイを進め、同時に累積輝度計算を行う。

図２０は、レンダリング処理を示すフローチャートである。まず、制御部１１は、座標変換パラメータを設定する（ステップＳ３１）。後述する探索制御マスクの算出（ステップＳ３２）、レイキャスティング処理（ステップＳ３３）の各処理で、視点座標系からボクセル座標系への座標変換処理を行う。そのため、制御部１１は、以下のように、座標変換処理に共通する座標変換パラメータを事前に設定・算出しておく。

＜座標変換パラメータ＞
回転パラメータ行列Ｒ：
Ｒ＝［Ｒ１１ Ｒ１２ Ｒ１３；
Ｒ２１ Ｒ２２ Ｒ２３；
Ｒ３１ Ｒ３２ Ｒ３３］
（ボクセル座標系から視点座標系への座標変換を行うための３×３の行列の逆行列、ＧＵＩ側はボクセル座標系から視点座標系への座標変換を指示するが、レンダリング側は視点座標系からボクセル座標系に座標変換を行う。）
Ｘ軸方向のオフセットＸｏｆｆ
Ｙ軸方向のオフセットＹｏｆｆ
Ｚ軸方向のオフセットＺｏｆｆ
Ｘ軸方向のＲＯＩ：Ｘｓ－Ｘｅ（０≦Ｘｓ、Ｘｅ≦Ｓｘ－１）
Ｙ軸方向のＲＯＩ：Ｙｓ－Ｙｅ（０≦Ｙｓ、Ｙｅ≦Ｓｙ－１）
Ｚ軸方向のＲＯＩ：Ｚｓ－Ｚｅ（０≦Ｚｓ、Ｚｅ≦Ｓｚ－１）
透視変換パラメータ：原点からＺ軸方向の距離Ｄｉｓｔ（正の実数値、平行投影の場合はＤｉｓｔ＝０）
拡大縮小倍率Ｓｃａｌｅ（ＸＹＺ軸方向で同一）
Ｚ方向変倍率Ｓｃｚ＝Ｒｘｙ／Ｒｚ
Ｚ方向変倍後のＺ方向サイズＳｚ’＝Ｓｚ・Ｓｃｚ
座標変換サブサンプル・オフセット：Ｘ軸方向ｄｘ、Ｙ軸方向ｄｙ、Ｚ軸方向ｄｚ

制御部１１は、上記した回転パラメータ行列Ｒに対して、初期値として単位行列を設定する。すなわち、Ｒ１１＝１、Ｒ１２＝０、Ｒ１３＝０、Ｒ２１＝０、Ｒ２２＝１、Ｒ２３＝０、Ｒ３１＝０、Ｒ３２＝０、Ｒ３３＝１と設定する。
そして、ＧＵＩの指示に従い、Ｘ軸中心回転Ｒｘ、Ｙ軸中心回転Ｒｙ、Ｚ軸中心回転Ｒｚ（角度単位：ラジアン）のいずれかを逐次指定し、以下のように、各々回転行列Ａを生成して回転パラメータ行列Ｒに右から乗算して、回転パラメータ行列Ｒを更新する。これにより、ＧＵＩの指示により生成されるボクセル座標系から視点座標系への回転行列の逆行列が算出される。

回転行列Ａを
Ａ＝［Ａ１１ Ａ１２ Ａ１３；
Ａ２１ Ａ２２ Ａ２３；
Ａ３１ Ａ３２ Ａ３３］
とすると、
Ｘ軸中心回転Ｒｘの場合の回転行列Ａの各要素は、
Ａ１１＝１、Ａ１２＝０、Ａ１３＝０
Ａ２１＝０、Ａ２２＝ｃｏｓＲｘ、Ａ２３＝ｓｉｎＲｘ
Ａ３１＝０、Ａ３２＝ｓｉｎＲｘ、Ａ３３＝ｃｏｓＲｘ
Ｙ軸中心回転Ｒｙの場合の回転行列Ａの各要素は、
Ａ１１＝ｃｏｓＲｙ、Ａ１２＝０、Ａ１３＝ｓｉｎＲｙ
Ａ２１＝０、Ａ２２＝１、Ａ２３＝０
Ａ３１＝－ｓｉｎＲｙ、Ａ３２＝０、Ａ３３＝ｃｏｓＲｙ
Ｚ軸中心回転Ｒｚの場合の回転行列Ａの各要素は、
Ａ１１＝ｃｏｓＲｚ、Ａ１２＝ｓｉｎＲｚ、Ａ１３＝０
Ａ２１＝－ｓｉｎＲｚ、Ａ２２＝ｃｏｓＲｚ、Ａ２３＝０
Ａ３１＝０、Ａ３２＝０、Ａ３３＝１
となる。
回転パラメータ行列Ｒは、Ｒ←Ｒ×Ａと更新される。

以上の座標変換処理は、探索制御マスクの算出（ステップＳ３２）、レイキャスティング処理（ステップＳ３３）の各処理の中で逐次実行される。

次に、図２１のフローチャートを参照して、図２０のステップＳ３２における探索制御マスクの算出処理について説明する。
探索制御マスクＭ（ｘ、ｙ、ｌ）（０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｌ≦Ｌ－１）は、算出するレンダリング像の画素（ｘ、ｙ）毎に先頭の不透明ボクセルのＺ座標を記憶したものである。座標変換時にサブサンプリング（ボリュームレンダリング像のジャギーを除去するため、アンチエリアシングともいわれる）を行う場合は、制御部１１は、以下の処理をＬ回（０≦ｌ≦Ｌ－１）実行する。

まず、制御部１１は、処理回数ｌ＝０とし、サブサンプル・オフセット値をｄｘ＝ｄｙ＝ｄｚ＝０と初期化する（ステップＳ４１）。

続いて、制御部１１は、ｘ＝２ｍ及びｙ＝２ｎ（ｍ、ｎは整数）の値をとる画素（ｘ、ｙ）に対して、Ｚｓ＝Ｓｚ’－１から座標変換を行いながら先頭の不透明ボクセルを探索し、そのＺ座標をＭ（ｘ、ｙ、ｌ）に記録する。ｘ＝Ｓｘ－１又はｙ＝Ｓｙ－１のいずれかの値をもつ画素（ｘ、ｙ）に対しても同様の探索を行う（ステップＳ４２）。これにより、図２２（ａ）のように、縦横偶数番目の画素について不透明ボクセルのＺ座標が算出される。不透明度ボクセルを探索する処理は、図２３にて後述する。

続いて、制御部１１は、ｘ＝２ｍ＋１（ｘ＜Ｓｘ－１）及びｙ＝２ｎ（ｍ、ｎは整数）の値をとる画素（ｘ、ｙ）に対して、Ｍ（ｘ－１、ｙ、ｌ）＝Ｍ（ｘ＋１、ｙ、ｌ）の場合、Ｍ（ｘ、ｙ、ｌ）＝Ｍ（ｘ－１、ｙ、ｌ）を与え、Ｍ（ｘ－１、ｙ、ｌ）≠Ｍ（ｘ＋１、ｙ、ｌ）の場合、上記同様Ｚｓ＝Ｓｚ’－１から座標変換を行いながら先頭の不透明ボクセルを探索（図２３参照）し、そのＺ座標をＭ（ｘ、ｙ、ｌ）に記録する（ステップＳ４３）。これにより、図２２（ｂ）のように、Ｘ方向に奇数番目の画素について不透明ボクセルのＺ座標が算出される。隣接する左右画素の不透明ボクセルのＺ座標が同一の場合、そのＺ座標で補間する（不透明ボクセルの探索は行わない）。

続いて、制御部１１は、全てのｘ（０≦ｘ≦Ｓｘ－１）及びｙ＝２ｎ＋１（ｎは整数）の値をとる画素（ｘ、ｙ）に対して、Ｍ（ｘ、ｙ－１、ｌ）＝Ｍ（ｘ、ｙ＋１、ｌ）の場合、Ｍ（ｘ、ｙ、ｌ）＝Ｍ（ｘ、ｙ－１、ｌ）を与え、Ｍ（ｘ、ｙ－１、ｌ）≠Ｍ（ｘ、ｙ＋１、ｌ）の場合、上記同様Ｚｓ＝Ｓｚ’－１から座標変換を行いながら先頭の不透明ボクセルを探索（図２３参照）し、そのＺ座標をＭ（ｘ、ｙ、ｌ）に記録する（ステップＳ４４）。これにより、図２２（ｃ）のように、Ｙ方向に奇数番目の画素について不透明ボクセルのＺ座標が算出される。隣接する上下画素の不透明ボクセルのＺ座標が同一の場合、そのＺ座標で補間する（不透明ボクセルの探索は行わない）。

続いて、制御部１１は、処理回数をｌ←ｌ＋１に更新し、サブサンプル・オフセット値を、ｄｘ←ｄｘ＋１／Ｌ、ｄｙ←ｄｙ＋１／Ｌ、ｄｚ←ｄｚ＋１／Ｌに更新する（ステップＳ４５）。
制御部１１は、処理回数ｌ＝Ｌ－１となるまで（ステップＳ４６；Ｙｅｓ）、上記したステップＳ４２～Ｓ４５の処理を繰り返す。

次に、図２３のフローチャートを参照して、探索制御マスク算出処理（図２１）において実行される不透明度ボクセルの探索処理について説明する。
まず、制御部１１は、ｚをｚ＝Ｚｓと初期化し、探索対象画素（ｘ、ｙ）を入力する（ステップＳ５１）。続いて、制御部１１は、３次元座標（ｘ、ｙ、ｚ）に対して座標変換を行い、ボクセルα値（ｎ＝３）を算出する（ステップＳ５２）。座標変換を行い、ボクセルα値を算出する処理（座標変換処理）は後述する。

ステップＳ５２において算出されたボクセルα値がα＝０の場合、制御部１１は、ｚ←ｚ－ｍ（例えばｍ＝８）に更新する（ステップＳ５３）。そして、ｚ＜０の場合、制御部１１は、Ｍ（ｘ、ｙ、ｌ）＝－１に設定し（ステップＳ５４）、処理を終了する。ｚ≧０の場合、制御部１１は、ステップＳ５２に戻る。

一方、ステップＳ５２において算出されたボクセルα値がα＞０の場合、制御部１１は、まず、ｉをｉ＝０と初期化する（ステップＳ５５）。続いて、制御部１１は、ｉ←ｉ＋１に更新し（ステップＳ５６）、ｉ＜ｍかつｚ＋ｉ＜Ｚｓの場合、３次元座標（ｘ、ｙ、ｚ＋ｉ）に対して座標変換を行い、ボクセルα値を算出する（ステップＳ５７）。座標変換を行い、ボクセルα値を算出する処理（座標変換処理）は後述する。

一方、ｉ≧ｍまたはｚ＋ｉ≧Ｚｓの場合、制御部１１は、Ｍ（ｘ、ｙ、ｌ）＝ｚ＋ｉ－１に設定し（ステップＳ５８）、処理を終了する。

また、ステップＳ５７において算出されたボクセルα値がα＝０の場合、制御部１１は、Ｍ（ｘ、ｙ、ｌ）＝ｚ＋ｉ－１に設定し（ステップＳ５８）、処理を終了する。

一方、ステップＳ５７において算出されたボクセルα値がα＞０の場合、ステップＳ５６に戻り、ステップＳ５８において不透明ボクセルのＺ座標が設定されるまで、ステップＳ５６～Ｓ５７の処理を繰り返す。

不透明度ボクセル探索処理（図２３）のステップＳ５２及びステップＳ５７において実行される座標変換処理について説明する。
座標変換処理は、視点座標系をボクセル座標系に変換する処理であり、ＧＵＩ側の変換処理とは逆になる。ＧＵＩ側では関心領域ＲＯＩによるクリッピング、スケーリング、Ｚ方向変倍処理、オフセット（ＸＹＺ方向同時）、回転、透視変換の順に行うものと仮定し、制御部１１は、与えられた視点座標系の３次元座標値（ｘ、ｙ、ｚ）（整数値）に対応する８ビットのボクセル構造体Ｖα（ｘ、ｙ、ｚ）又は２４ビットのボクセル構造体Ｖｃ（ｘ、ｙ、ｚ、ｎ）の実数の座標値（ｘｒ、ｙｒ、ｚｒ）を以下のように算出する。

制御部１１は、視点座標系の座標値（ｘ、ｙ、ｚ）（整数値）を次のように実数値（ｘｘ、ｙｙ、ｚｚ）に変換する。
ｘｘ＝ｘ－Ｓｘ／２＋ｄｘ
ｙｙ＝ｙ－Ｓｙ／２＋ｄｙ
ｚｚ＝ｚ－Ｓｚ’／２＋ｄｚ

続いて、原点からＺ軸方向の距離Ｄｉｓｔ＞０の場合、制御部１１は、以下のように透視変換を行う。
ｚｚ’＝Ｄｉｓｔ・ｚｚ／（Ｄｉｓｔ＋ｚｚ）
ｘｘ’＝ｘｘ・（Ｄｉｓｔ－ｚｚ’）／Ｄｉｓｔ
ｙｙ’＝ｙｙ・（Ｄｉｓｔ－ｚｚ’）／Ｄｉｓｔ
透視変換後の（ｘｘ’、ｙｙ’、ｚｚ’）を（ｘｘ、ｙｙ、ｚｚ）とする。

続いて、制御部１１は、以下のように、回転パラメータ行列Ｒを用いて回転処理を行う。
ｘｘ’＝Ｒ１１・ｘｘ＋Ｒ１２・ｙｙ＋Ｒ１３・ｚｚ
ｙｙ’＝Ｒ２１・ｘｘ＋Ｒ２２・ｙｙ＋Ｒ２３・ｚｚ
ｚｚ’＝Ｒ３１・ｘｘ＋Ｒ３２・ｙｙ＋Ｒ３３・ｚｚ
回転処理後の（ｘｘ’、ｙｙ’、ｚｚ’）を（ｘｘ、ｙｙ、ｚｚ）とする。

制御部１１は、スケーリング、Ｚ方向変倍処理、及びオフセット（ＸＹＺ方向同時）を同時に行い、次のようにボクセル構造体の対応するボクセルの座標値（ｘｒ、ｙｒ、ｚｒ）（実数値）を取得する。
ｘｒ＝ｘｘ／Ｓｃａｌｅ＋Ｓｘ／２－Ｘｏｆｆ
ｙｒ＝ｙｙ／Ｓｃａｌｅ＋Ｓｙ／２－Ｙｏｆｆ
ｚｒ＝ｚｚ／Ｓｃａｌｅ／Ｓｃｚ＋Ｓｚ／２－Ｚｏｆｆ

続いて、制御部１１は、算出したボクセルの座標値（ｘｒ、ｙｒ、ｚｒ）（実数値）に対して、小数点以下を切り捨て整数化した座標値を（ｘｉ、ｙｉ、ｚｉ）（整数値）とし、切り捨てた小数点以下の端数を（ｗｘ、ｗｙ、ｗｚ）（０≦ｗｘ、ｗｙ、ｗｚ＜１）とする（すなわち、ｘｒ＝ｘｉ＋ｗｘ、ｙｒ＝ｙｉ＋ｗｙ、ｚｒ＝ｚｉ＋ｗｚ）。

そして、制御部１１は、Ｘ軸・Ｙ軸・Ｚ軸方向のＲＯＩを考慮して、視点座標系の３次元座標値（ｘ、ｙ、ｚ）（整数値）に対応するボクセルのα値（Ｖα）を次のように３通りの場合に分けて決定する。

１）ｘｉ＜Ｘｓ、ｘｉ＞Ｘｅ、ｙｉ＜Ｙｓ、ｙｉ＞Ｙｅ、ｚｉ＜Ｚｓ、又はｚｉ＞Ｚｅのいずれかを満たす場合（クリッピング範囲）
Ｖα＝０
２）上記１）の条件を満たさない場合において、ｘｉ＋１＞Ｘｅ、ｙｉ＋１＞Ｙｅ、又はｚｉ＋１＞Ｚｅのいずれかを満たす場合（補間しない）
Ｖα＝Ｖα（ｘｉ、ｙｉ、ｚｉ）

３）上記１）２）の条件を満たさない場合（補間する）
Ｖα＝（１－ｗｚ）（１－ｗｙ）（１－ｗｘ）・Ｖα（ｘｉ、ｙｉ、ｚｉ）＋（１－ｗｚ）（１－ｗｙ）・ｗｘ・Ｖα（ｘｉ＋１、ｙｉ、ｚｉ）＋（１－ｗｚ）・ｗｙ・（１－ｗｘ）・Ｖα（ｘｉ、ｙｉ＋１、ｚｉ）＋（１－ｗｚ）・ｗｙ・ｗｘ・Ｖα（ｘｉ＋１、ｙｉ＋１、ｚｉ）＋ｗｚ・（１－ｗｙ）（１－ｗｘ）・Ｖα（ｘｉ、ｙｉ、ｚｉ＋１）＋ｗｚ・（１－ｗｙ）・ｗｘ・Ｖα（ｘｉ＋１、ｙｉ、ｚｉ＋１）＋ｗｚ・ｗｙ・（１－ｗｘ）・Ｖα（ｘｉ、ｙｉ＋１、ｚｉ＋１）＋ｗｚ・ｗｙ・ｗｘ・Ｖα（ｘｉ＋１、ｙｉ＋１、ｚｉ＋１）

尚、後述するレイキャスティング処理（図２４参照）のようにボクセルのＲＧＢ値を取得する場合は、制御部１１は、視点座標系の３次元座標値（ｘ、ｙ、ｚ）（整数値）に対応するボクセルのＲＧＢ値（Ｖｃ（ｎ））（０≦ｎ≦２）を次のように３通りの場合に分けて決定する。

１）ｘｉ＜Ｘｓ、ｘｉ＞Ｘｅ、ｙｉ＜Ｙｓ、ｙｉ＞Ｙｅ、ｚｉ＜Ｚｓ、又はｚｉ＞Ｚｅのいずれかを満たす場合（クリッピング範囲）
Ｖｃ（ｎ）＝０ （０≦ｎ≦２）
２）上記１）の条件を満たさない場合において、ｘｉ＋１＞Ｘｅ、ｙｉ＋１＞Ｙｅ、又はｚｉ＋１＞Ｚｅのいずれかを満たす場合（補間しない）
Ｖｃ（ｎ）＝Ｖｃ（ｘｉ、ｙｉ、ｚｉ、ｎ） （０≦ｎ≦２）

３）上記１）２）の条件を満たさない場合（補間する）
Ｖｃ（ｎ）＝（１－ｗｚ）（１－ｗｙ）（１－ｗｘ）・Ｖｃ（ｘｉ、ｙｉ、ｚｉ、ｎ）＋（１－ｗｚ）（１－ｗｙ）・ｗｘ・Ｖｃ（ｘｉ＋１、ｙｉ、ｚｉ、ｎ）＋（１－ｗｚ）・ｗｙ・（１－ｗｘ）・Ｖｃ（ｘｉ、ｙｉ＋１、ｚｉ、ｎ）＋（１－ｗｚ）・ｗｙ・ｗｘ・Ｖｃ（ｘｉ＋１、ｙｉ＋１、ｚｉ、ｎ）＋ｗｚ・（１－ｗｙ）（１－ｗｘ）・Ｖｃ（ｘｉ、ｙｉ、ｚｉ＋１、ｎ）＋ｗｚ・（１－ｗｙ）・ｗｘ・Ｖｃ（ｘｉ＋１、ｙｉ、ｚｉ＋１、ｎ）＋ｗｚ・ｗｙ・（１－ｗｘ）・Ｖｃ（ｘｉ、ｙｉ＋１、ｚｉ＋１、ｎ）＋ｗｚ・ｗｙ・ｗｘ・Ｖｃ（ｘｉ＋１、ｙｉ＋１、ｚｉ＋１、ｎ） （０≦ｎ≦２）

次に、図２４のフローチャートを参照して、図２０のステップＳ３３におけるレイキャスティング処理について説明する。
制御部１１は、まず、生成する２４ビット（ＲＧＢ）のレンダリング画像Ｉｍａｇｅ（ｘ、ｙ、ｎ）の初期値を全て０に設定する（Ｉｍａｇｅ（ｘ、ｙ、ｎ）＝０、ｎ＝０（Ｒ）、１（Ｇ）、２（Ｂ））。そして、サブサンプル回数Ｌとして、各２次元座標（ｘ、ｙ）（０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１）に対して、以下の処理を実行する。

まず、制御部１１は、光源Ｌｉｇｈｔ（ｎ）（ｎ＝０（Ｒ）、１（Ｇ）、２（Ｂ）、０≦Ｌｉｇｈｔ（ｎ）≦２５５）を設定し、サブサンプル・オフセット値をｄｘ＝ｄｙ＝ｄｚ＝０、処理回数ｌ＝０に初期化する（ステップＳ６１）。

続いて、制御部１１は、ｚ＝Ｍ（ｘ、ｙ、ｌ）とし（ステップＳ６２）、ｚ＜０の場合、ステップＳ７０へ移行し、画素（ｘ、ｙ）におけるＲＧＢ画素値を決定する。

一方、ｚ≧０の場合、制御部１１は、仮想光強度Ｔｒａｎｓ＝１．０、累積輝度値Ｅｎｅｒｇｙ（ｎ）＝０．０（０≦ｎ≦２）に初期化し（ステップＳ６３）、３次元座標（ｘ、ｙ、ｚ）を座標変換してボクセルα値（Ｖα）を取得し、Ａｌｐｈａ＝Ｖα／２５５を算出する（ステップＳ６４）。

ステップＳ６４において算出したＡｌｐｈａがＡｌｐｈａ＝０かつｚ＝０の場合、ステップＳ７０へ移行し、画素（ｘ、ｙ）におけるＲＧＢ画素値を決定する。

一方、ステップＳ６４において算出したＡｌｐｈａがＡｌｐｈａ＞０の場合、制御部１１は、３次元座標（ｘ、ｙ、ｚ）を座標変換してボクセルＲＧＢ値Ｖｃ（ｎ）を取得する（ステップＳ６７）。ボクセルのＲＧＢ値Ｖｃ（ｎ）を取得は、前述した座標変換処理により行われる。

続いて、制御部１１は、累積輝度をＥｎｅｒｇｙ（ｎ）＝Ｅｎｅｒｇｙ（ｎ）＋Ｔｒａｎｓ／Ａｌｐｈａ・Ｖｃ（ｎ）／２５５、透過光強度をＴｒａｎｓ＝Ｔｒａｎｓ・（１．０－Ａｌｐｈａ）と更新する（ステップＳ６８）。Ａｌｐｈａ＝１．０又はＴｒａｎｓ＜０．００１の場合、ステップＳ７０へ移行し、画素（ｘ、ｙ）におけるＲＧＢ画素値を決定する。それ以外の場合、制御部１１は、ｚ＝ｚ－１に更新し（ステップＳ６９）、ｚ≧０の場合、ステップＳ６４に戻り、ｚ≦０の場合、ステップＳ７０へ移行し、画素（ｘ、ｙ）におけるＲＧＢ画素値を決定する。

一方、ステップＳ６４において算出したＡｌｐｈａがＡｌｐｈａ＝０かつｚ＞０の場合、制御部１１は、Ｚｓ＝ｚ－１から座標変換して先頭の不透明ボクセルのＺ座標ｚ’を探索する（ステップＳ６５）。不透明ボクセルの探索は、前述した不透明ボクセル探索処理（図２３参照）により実行される。

探索されたｚ’がｚ’＜０の場合、ステップＳ７０へ移行し、画素（ｘ、ｙ）におけるＲＧＢ画素値を決定する。ｚ’≧０の場合、制御部１１は、ｚ＝ｚ'に更新し（ステップＳ６６）、ステップＳ６４の処理に戻る。

ステップＳ７０において、制御部１１は、画素（ｘ、ｙ）におけるＲＧＢ画素値を次のように決定する。
Ｉｍａｇｅ（ｘ、ｙ、ｎ）＝Ｉｍａｇｅ（ｘ、ｙ、ｎ）＋ｋ・Ｅｎｅｒｇｙ（ｎ）・Ｌｉｇｈｔ（ｎ）／Ｌ
ここで、ｋは強度倍率であり、初期値はｋ＝１．０に設定される。

続いて、制御部１１は、処理回数をｌ←ｌ＋１に更新し、サブサンプル・オフセット値を、ｄｘ←ｄｘ＋１／Ｌ、ｄｙ←ｄｙ＋１／Ｌ、ｄｚ←ｄｚ＋１／Ｌに更新する（ステップＳ７１）。制御部１１は、処理回数ｌ＝Ｌ－１となるまで（ステップＳ７２；Ｙｅｓ）、上記したステップＳ６２～Ｓ７１の処理を繰り返し、最終的なレンダリング画像Ｉｍａｇｅ（ｘ、ｙ、ｎ）（ｎ≦０≦２）を得る。

以上、図１９～図２４を参照しながら、ＣＰＵ（第１レンダリング部２８ａ）により実行されるレイキャスティング法によるレンダリング処理について説明した。

次に、図２５～図２９を参照しながら、ＧＰＵ（第２レンダリング部２８ｂ）により実行される３Ｄテクスチャマッピング法によるレンダリング処理について説明する。

ＧＰＵによるレンダリングはＯｐｅｎＧＬシステムを用いて実行するため、制御部１１は、作成したボクセル構造体Ｖを３ＤテクスチャとしてＯｐｅｎＧＬシステムへ登録する。具体的には、制御部１１は、２４ビットのボクセル構造体Ｖｃ（（式１６）参照）と８ビットのボクセル構造体Ｖα（式（１３）参照）を３２ビットのボクセル構造体Ｖ（（式１８）参照）に統合し、ｇｌＴｅｘＩｍａｇｅ３Ｄ関数を用いて３Ｄテクスチャを生成し、３ＤテクスチャをＯｐｅｎＧＬシステムに登録する。

また、制御部１１は、図２５に示すように、３Ｄテクスチャを構成する各スライスを四角形で定義し、３次元空間のＸＹ座標面上のＺ軸方向に並べた積層四角形を設定し、テクスチャ座標と四角形座標との対応付けを行い、３Ｄテクスチャを各四角形に貼り付ける。３Ｄテクスチャマッピング法では、四角形に貼り付けた３Ｄテクスチャマップを視線方向に合わせて座標変換することで、任意方向のレンダリング像を生成する。

尚、ＧＰＵ（第２レンダリング部２８ｂ）によるレンダリング処理では、レンダリング時にモアレが発生する場合がある。制御部１１は、レンダリング時のモアレの発生を抑えるため、次のように、関心領域ＲＯＩ境界面の透明ボクセルの色補正処理を事前に行う。

図２６は色補正処理のフローチャートである。制御部１１は、境界面ボクセル（ｘ、ｙ、ｚ）を抽出する（ステップＳ８１）。境界面ボクセル（ｘ、ｙ、ｚ）とは、ｘ＝Ｘｓ、ｘ＝Ｘｅ、ｙ＝Ｙｓ、ｙ＝Ｙｅ、ｚ＝Ｚｓ、又はｚ＝Ｚｅのいずれかを満たすボクセルである。

続いて、制御部１１は、抽出した境界面ボクセルの２７近傍ボクセルより、不透明ボクセルを抽出し、抽出した不透明ボクセルの平均色（Ｒａ、Ｇａ、Ｂａ）を算出する（ステップＳ８２）。不透明ボクセルとは、Ｖα（ｘ＋ｉ、ｙ＋ｊ、ｚ＋ｋ）＞０（－１≦ｉ≦１、－１≦ｊ≦１、－１≦ｋ≦１）の条件を満たすボクセルである。尚、不透明ボクセルが抽出されない場合は、制御部１１は、Ｒａ＝Ｇａ＝Ｂａ＝０とする。

そして、制御部１１は、以下のように、境界面ボクセルの不透明度を０に設定し、ＲＧＢ値を算出した不透明ボクセルの平均色に設定する（ステップＳ８３）。
Ｖα（ｘ、ｙ、ｚ）＝０
Ｖｃ（ｘ、ｙ、ｚ、１）＝Ｒａ、Ｖｃ（ｘ、ｙ、ｚ、２）＝Ｇａ、Ｖｃ（ｘ、ｙ、ｚ、３）＝Ｂａ

制御部１１は、全ての境界面ボクセル（ｘ、ｙ、ｚ）について、上記した色補正処理（ステップＳ８１～８３）を実行する。

次に、図２７のフローチャートを参照して、３Ｄテクスチャマッピング法によるレンダリング処理の流れを説明する。まず、制御部１１は、ＯｐｅｎＧＬの透視変換パラメータを設定する（ステップＳ９１）。仮想内視鏡モードでレンダリングを行う場合は、制御部１１は、ＯｐｅｎＧＬシステムに対して、次の透視変換パラメータを設定する。

＜透視変換パラメータ＞
・カメラの視野角度（単位：度、平行投影の場合は０、焦点距離に相当）
・視点座標（カメラの位置、通常はＺ軸上に設定）
・注視点座標（被写体の凝視点、通常はＺ軸上で原点ｚ＝０に固定）
・近方クリッピング位置（視点からの距離、視点より若干距離を置く）
・遠方クリッピング位置（視点からの距離、通常は１０００など無限大に固定しクリッピングしない）

一方、通常のレンダリングを行う際は、カメラの視野角度を０に設定し、平行投影モードに設定する。

続いて、制御部１１は、ＯｐｅｎＧＬの投影画面設定を行う（ステップＳ９２）。
図２８は、投影画面設定の処理を示すフローチャートである。図２８に示すように、制御部１１は、レンダリング画像のスクリーンサイズ（縦横画素数、縦横アスペクト比率）を設定し（ステップＳ１０１）、平行投影（通常の外観レンダリング）又は透視投影（内視鏡モード）の設定を行う（ステップＳ１０２）。そして、透視投影が設定された場合は、制御部１１は、透視投影パラメータ（カメラの視野角度（焦点距離）、視点座標、注視点座標、近方クリッピング位置、遠方クリッピング位置）の設定を行う（ステップＳ１０３）。

図２７のフローチャートの説明に戻る。
続いて、制御部１１は、以下のように、ＯｐｅｎＧＬの座標変換パラメータ（回転、スケール、移動、Ｚ方向変倍パラメータＲｘｙ／Ｒｚ）を設定する（ステップＳ９３）。

＜ＯｐｅｎＧＬの座標変換パラメータ＞
回転パラメータ行列Ｒ：
Ｒ＝［Ｒ１１ Ｒ１２ Ｒ１３ Ｒ１４；
Ｒ２１ Ｒ２２ Ｒ２３ Ｒ２４；
Ｒ３１ Ｒ３２ Ｒ３３ Ｒ３４；
Ｒ４１ Ｒ４２ Ｒ４３ Ｒ４４］
（ボクセル座標系から視点座標系への座標変換を行うための４×４の行列で、レンダリング側もＧＵＩ側と同一方向の座標変換をＯｐｅｎＧＬに指示する。実際は３×３の座標変換行列で定義できるが、ＯｐｅｎＧＬの仕様上、４×４の座標変換行列に拡張）
Ｘ軸方向のオフセットＸｏｆｆ
Ｙ軸方向のオフセットＹｏｆｆ
Ｚ軸方向のオフセットＺｏｆｆ
拡大縮小倍率Ｓｃａｌｅ（ＸＹＺ軸方向で同一）
Ｚ方向変倍率Ｓｃｚ＝Ｒｘｙ／Ｒｚ
Ｚ方向変倍後のＺ方向サイズＳｚ’＝Ｓｚ・Ｓｃｚ

制御部１１は、上記した回転パラメータ行列Ｒに対して、初期値として単位行列を設定する。すなわち、Ｒ１１＝１、Ｒ１２＝０、Ｒ１３＝０、Ｒ１４＝０、Ｒ２１＝０、Ｒ２２＝１、Ｒ２３＝０、Ｒ２４＝０、Ｒ３１＝０、Ｒ３２＝０、Ｒ３３＝１、Ｒ４１＝０、Ｒ４２＝０、Ｒ４３＝０、Ｒ４４＝１と設定する。
そして、ＧＵＩの指示に従い、Ｘ軸中心回転Ｒｘ、Ｙ軸中心回転Ｒｙ、Ｚ軸中心回転Ｒｚ（角度単位：ラジアン）のいずれかを逐次指定し、以下のように、各々回転行列Ａを生成して回転パラメータ行列Ｒに左から乗算して、回転パラメータ行列Ｒを更新する。

回転行列Ａを
Ａ＝［Ａ１１ Ａ１２ Ａ１３ Ａ１４；
Ａ２１ Ａ２２ Ａ２３ Ａ２４；
Ａ３１ Ａ３２ Ａ３３ Ａ３４；
Ａ４１ Ａ４２ Ａ４３ Ａ４４］
とすると、
Ｘ軸中心回転Ｒｘの場合の回転行列Ａの各要素は、
Ａ１１＝１、Ａ１２＝０、Ａ１３＝０、Ａ１４＝０
Ａ２１＝０、Ａ２２＝ｃｏｓＲｘ、Ａ２３＝ｓｉｎＲｘ、Ａ２４＝０
Ａ３１＝０、Ａ３２＝ｓｉｎＲｘ、Ａ３３＝ｃｏｓＲｘ、Ａ３４＝０
Ａ４１＝０、Ａ４２＝０、Ａ４３＝０、Ａ４４＝１
Ｙ軸中心回転Ｒｙの場合の回転行列Ａの各要素は、
Ａ１１＝ｃｏｓＲｙ、Ａ１２＝０、Ａ１３＝ｓｉｎＲｙ、Ａ１４＝０
Ａ２１＝０、Ａ２２＝１、Ａ２３＝０、Ａ２４＝０
Ａ３１＝－ｓｉｎＲｙ、Ａ３２＝０、Ａ３３＝ｃｏｓＲｙ、Ａ３４＝０
Ａ４１＝０、Ａ４２＝０、Ａ４３＝０、Ａ４４＝１
Ｚ軸中心回転Ｒｚの場合の回転行列Ａの各要素は、
Ａ１１＝ｃｏｓＲｚ、Ａ１２＝ｓｉｎＲｚ、Ａ１３＝０、Ａ１４＝０
Ａ２１＝－ｓｉｎＲｚ、Ａ２２＝ｃｏｓＲｚ、Ａ２３＝０、Ａ２４＝０
Ａ３１＝０、Ａ３２＝０、Ａ３３＝１、Ａ３４＝０
Ａ４１＝０、Ａ４２＝０、Ａ４３＝０、Ａ４４＝１
となる。
回転パラメータ行列Ｒは、Ｒ←Ａ×Ｒと更新される。

図２７のフローチャートの説明に戻る。
制御部１１は、ＯｐｅｎＧＬによるレンダリング処理を実行する（ステップＳ９４）。
図２９は、ＯｐｅｎＧＬによるレンダリング処理を示すフローチャートである。図２９に示すように、制御部１１は、拡大縮小倍率Ｓｃａｌｅ（ＸＹＺ軸方向で同一）、及びＺ方向変倍率Ｓｃｚに基づいて３Ｄテクスチャマップに対してスケーリング、及びＺ方向変倍処理を施し（ステップＳ１１１）、回転パラメータ行列Ｒに基づいて３Ｄテスチャマップに対して回転処理を施し（ステップＳ１１２）、Ｘ軸方向のオフセットＸｏｆｆ、Ｙ軸方向のオフセットＹｏｆｆ、及びＺ軸方向のオフセットＺｏｆｆに基づいて３Ｄテスチャマップに対してオフセット処理を施すことで（ステップＳ１１３）、３Ｄテクスチャマップに対して所定の座標変換を行う。

続いて、制御部１１は、３次元空間のＸＹ座標面上の四角形をＺ軸方向に並べた積層四角形をＯｐｅｎＧＬで設定し、３Ｄテスクチャマップの対応付けを行う（ステップＳ１１４、図２５参照）。そして、制御部１１は、ＯｐｅｎＧＬによるレンダリング処理を実行し、積層四角形に３Ｄテクスチャマップを貼り付けながらスキャンコンバージョンし、フレームメモリ上で奥からＺ方向（視点方向）にアルファブレンディング処理を行う（ステップＳ１１５）。すなわち、所定の視点からＺ軸方向に平行な視線上の四角形のＸＹ座標に対応する座標変換後の３Ｄテスクチャマップのボクセルの色値（ＲＧＢ値）を、ボクセルの不透明度（α値）に基づいて視点から遠い四角形の順にアルファブレンディングして取得し、ボリュームレンダリング像の画素値として与える。

以上、図２５～図２９を参照しながら、ＧＰＵ（第２レンダリング部２８ｂ）により実行される３Ｄテクスチャマッピング法によるレンダリング処理について説明した。

制御部１１は、ＣＰＵ（第１レンダリング部２８ａ）により実行されるレイキャスティング法によるレンダリング処理（図２０）、又はＧＰＵ（第２レンダリング部２８ｂ）により実行される３Ｄテクスチャマッピング法によるレンダリング処理（図２７）によりボリュームレンダリング像を生成し、生成したボリュームレンダリング像を表示部１６に表示する。

図６のフローチャートの説明に戻る。
ステップＳ９においてボリュームレンダリング像を生成し表示部１６に表示した後、ユーザがパラメータの調整が必要と判断した場合（ステップＳ１０；「Ｙｅｓ」）、パラメータ調整画面４０において、楕円パラメータＰ（ｚ）（（式１１）（式１３）参照）や補正倍率Ｓα（（式１２）（式１５）参照）を算出する際の各種パラメータをユーザに調整させる（ステップＳ１１）。

図３０は、パラメータ調整画面４０の例を示す。図に示すように、パラメータ調整画面４０において、（式１１）（式１３）における楕円の横方向のサイズの補正係数ａ（図の楕円倍率Ｘ）、楕円の縦方向のサイズの補正係数ｂ（図の楕円倍率Ｙ）、Ｘ方向のオフトセットｘｏｆｆｓｅｔ（図の中心オフセットＸ）、Ｙ方向のオフセットｙｏｆｆｓｅｔ（図の中心オフセットＹ）、ｅ（図の複合楕円・中心径［％］）、ｄ（図の中心オフセット［％］）、（式１４）におけるｄｙｍ（図のＹ補正量）、（式１２）（式１５）におけるＸＹ平面の減衰係数、ＸＹ面の振幅倍率ａｍｐ、等を調整することができる。

パラメータの調整が終了すると（パラメータ調整画面４０の「ＯＫ」ボタンを押下すると）、ステップＳ７に戻り、制御部１１は、ステップＳ７～ステップＳ９の処理を再度実行する。すなわち、制御部１１は、調整されたパラメータに基づいて、再度、補正テーブルＳαを作成し直した上で（ステップＳ７）、ボクセル構造体Ｖ（ボクセルデータ）を作成し（ステップＳ８）、レンダリング処理を行う（ステップＳ９）。ステップＳ７～ステップＳ９の処理は、ステップＳ９において生成し表示されるレンダリング像を確認しながら、良好なレンダリング結果が得られユーザがパラメータの調整が十分と判断するまで（ステップＳ１０；「Ｎｏ」）、繰り返し実行する。

そして、制御部１１（データ出力部２９）は、ユーザからデータの保存操作を受け付けることによって、各種データを出力して記憶部１２に保存することができる（図６のステップＳ１２）。例えば、制御部１１は、ステップＳ９において生成し表示されたボリュームレンダリング像（３次元画像）を出力して保存する。

また、制御部１１は、ボクセル構造体Ｖ（ボクセルデータ）を作成する際に参照可能なデータ（参照データ）を出力しても良い。これにより、同じ断層画像群Ｄｏに基づいて再度レンダリングを行う際、いくつかの処理を省略することができる。参照データとしては、例えば、以下のものが挙げられる。

＜参照データ１＞
ボクセル構造体Ｖ（Ｖｃ、Ｖα）の作成（（式１６）、（式１７）参照）に使用されたカラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦３）のデータ

＜参照データ２＞
補正テーブルＳα（（式１２）（式１５）参照）の作成に使用された楕円パラメータＰ（ｚ）（（式１１）、（式１３）参照）のデータ

＜参照データ３＞
ボクセル構造体Ｖαの作成（（式１７）参照）に使用された補正テーブルＳα（（式１２）（式１５）参照）のデータ

＜参照データ４＞
断層画像群Ｄｏの各画素（ｘ、ｙ、ｚ）に適用される補正後の不透明度を格納した不透明度テーブルＴα（ｘ、ｙ、ｚ）（３次元α値分布）

ここで、不透明度テーブルＴαは、以下のように定義される。

（式２１）
Ｔα（ｘ、ｙ、ｚ）＝Ｃｍａｐ（Ｄｏ（ｘ、ｙ、ｚ）、３）・Ｓα（ｘ、ｙ、ｚ）
（０≦Ｔα（ｘ、ｙ、ｚ）≦２５５、０≦ｘ≦Ｓｘ－１、０≦ｙ≦Ｓｙ－１、０≦ｚ≦Ｓｚ－１；解像度：Ｒｘｙ、Ｒｚ）

＜参照データ５＞
ボクセル構造体Ｖ（Ｖｃ、Ｖα）の作成（（式１６）（式１７）参照）に使用されたカラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦３）と、補正テーブルＳα（（式１２）（式１５）参照）の作成に使用された楕円パラメータＰ（ｚ）（（式１１）（式１３）参照）を含むデータ

＜参照データ６＞
ボクセル構造体Ｖ（Ｖｃ、Ｖα）の作成（（式１６）（式１７）参照）に使用されたカラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦３）と、ボクセル構造体Ｖαの作成（（式１７）参照）に使用された補正テーブルＳα（（式１２）（式１５）参照）を含むデータ

＜参照データ７＞
ボクセル構造体Ｖｃの作成（（式１６）参照）に使用されたカラーマップＣｍａｐ（ｖ、ｎ）（０≦ｎ≦２）（カラーパレット）と、不透明度テーブルＴα（ｘ、ｙ、ｚ）を含むデータ

以上のような参照データを保存しておくことで、同じ断層画像群Ｄｏに基づいて再度レンダリングを行う際に、いくつかの処理を省略でき、レンダリング像を得るまでの処理時間を短縮できる。

例えば、参照データ１（カラーマップＣｍａｐ）を保存しておけば、図６のステップＳ４～Ｓ６（カラーマップの取得・調整処理）を省略できる。
また、参照データ２（楕円パラメータＰ（ｚ））を保存しておけば、図８のステップＳ２１～Ｓ２４、ステップＳ２１～Ｓ２７の処理（ステップＳ２５、Ｓ２８以外の処理）を省略できる。

また、参照データ３（補正テーブルＳα）を保存しておけば、図６のステップＳ７（補正テーブル作成処理）を省略できる。
また、参照データ４（不透明度テーブルＴα）を保存しておけば、図６のステップＳ７（補正テーブル作成処理）を省略できるとともに、（式１７）の不透明度の補正処理（Ｖα（ｘ、ｙ、ｚ）＝Ｖα（ｘ、ｙ、ｚ）・Ｓα（ｘ、ｙ、ｚ））を省略できる。すなわち、不透明度テーブルＴαに格納されている不透明度をボクセル構造体Ｖαに設定するだけでよい（Ｖα（ｘ、ｙ、ｚ）＝Ｔα（ｘ、ｙ、ｚ））。

また、参照データ５（カラーマップＣｍａｐ＋楕円パラメータＰ（ｚ））を保存しておけば、図６のステップＳ４～Ｓ６（カラーマップの取得・調整処理）を省略できるとともに、図８のステップＳ２１～Ｓ２４、ステップＳ２１～Ｓ２７の処理（ステップＳ２５、Ｓ２８以外の処理）を省略できる。

また、参照データ６（カラーマップＣｍａｐ＋補正テーブルＳα）を保存しておけば、図６のステップＳ４～Ｓ７（補正テーブルの作成処理、及びカラーマップの取得・調整処理）を省略できる。

また、参照データ７（カラーマップＣｍａｐ（カラーパレット）＋不透明度テーブルＴα）を保存しておけば、図６のステップＳ４～Ｓ７（補正テーブルの作成処理、及びカラーマップの取得・調整処理）を省略できるとともに、（式１７）の不透明度の補正処理（Ｖα（ｘ、ｙ、ｚ）＝Ｖα（ｘ、ｙ、ｚ）・Ｓα（ｘ、ｙ、ｚ））を省略できる。

以上、図６～図３０を参照しながら、医用画像処理装置１の動作について説明した。

[実施例]
最後に、従来手法、提案手法１、提案手法２により生成されるレンダリング像の比較を行う。従来手法とは、図６のステップＳ４において取得したカラーマップＣｍａｐを断層画像群Ｄｏにそのまま適用してボクセル構造体Ｖ（Ｖｃ、Ｖα）を作成し、レンダリング像を生成する手法である。すなわち、従来手法では、以下の（式２２）に基づいてボクセル構造体Ｖ（Ｖｃ、Ｖα）を作成する。

（式２２）
（１）Ｘｓ＜ｘ＜Ｘｅ、Ｙｓ＜ｙ＜Ｙｅ、及びＺｓ＜ｚ＜Ｚｅを全て満たす場合（関心領域ＲＯＩ内の場合）
Ｖｃ（ｘ、ｙ、ｚ、ｎ）＝Ｃｍａｐ（Ｄｏ（ｘ、ｙ、ｚ）、ｎ）（０≦ｎ≦２）
Ｖα（ｘ、ｙ、ｚ）＝Ｃｍａｐ（Ｄｏ（ｘ、ｙ、ｚ）、３）
（２）ｘ≦Ｘｓ、ｘ≧Ｘｅ、ｙ≦Ｙｓ、ｙ≧Ｙｅ、ｚ≦Ｚｓ、又はｚ≧Ｚｅのいずれかを満たす場合（関心領域ＲＯＩ外の場合）
Ｖｃ（ｘ、ｙ、ｚ、ｎ）＝０ （０≦ｎ≦２）
Ｖα（ｘ、ｙ、ｚ）＝０

提案手法１とは、本実施の形態に基づいてレンダリング像を生成する方法（図６の処理によりレンダリング像を生成する方法）であり、特に、図８のステップＳ２４～Ｓ２５により補正倍率Ｓαを算出する方法である。すなわち、被写体領域の外郭を単一の楕円で近似し、補正倍率Ｓαを算出する方法である。

提案手法２とは、本実施の形態に基づいてレンダリング像を生成する方法（図６の処理によりレンダリング像を生成する方法）であり、特に、図８のステップＳ２６～Ｓ２８により補正倍率Ｓαを算出する方法である。すなわち、被写体領域の外郭を複数の楕円で近似し、補正倍率Ｓαを算出する方法である。

なお、従来手法、提案手法１、提案手法２のいずれについても、カラーマップＣｍａｐの調整（図６のステップＳ６の処理）は行わないものとし、また、第１レンダリング部２８ａを用いてレンダリング処理を実行するものとした（第２レンダリング部２８ｂを用いても殆ど同様なレンダリング像が得られる）。

（実施例１）
まず、胸部ＣＴ画像を対象としてレンダリング像の比較を行った。
ここで、提案手法１における楕円パラメータＰ（ｚ）（（式１１）参照）および補正倍率Ｓα（（式１２）参照）の各種パラメータを、ｘｏｆｆｓｅｔ＝２（Ｘ方向のオフセット）、ｙｏｆｆｓｅｔ＝５（Ｙ方向のオフセット）、ａ＝０．５８（横方向のサイズの補正係数）、ｂ＝０．５８（縦方向のサイズの補正係数）、ｋ＝１．１（減衰係数）、ａｍｐ＝１．０（振幅倍率）と設定した。

また、提案手法２における楕円パラメータＰ１（ｚ）、Ｐ２（ｚ）、Ｐ３（ｚ）（（式１３）（式１４）参照）および補正倍率Ｓα（（式１５）参照）の各種パラメータを、ｘｏｆｆｓｅｔ＝２（Ｘ方向のオフセット）、ｙｏｆｆｓｅｔ＝５（Ｙ方向のオフセット）、ａ＝０．５８（横方向のサイズの補正係数）、ｂ＝０．５８（縦方向のサイズの補正係数）、ｋ＝１．１（減衰係数）、ａｍｐ＝１．０（振幅倍率）、ｄｙｍ＝９０、Ｚｔｈ＝０．５８×Ｓｚ、ｄ＝０．２２、ｅ＝０．５５と設定した。

図３１は、従来手法により生成したレンダリング像を正面から観察した画像であり、図３２は、提案手法１により生成したレンダリング像を正面から観察した画像である。図３１（従来手法）では、肋骨及び胸骨の透明化が不十分であり、内部の臓器等を観察することが困難であるのに対し、図３２（提案手法１）では、胸骨は多少残るものの肋骨全体が除去され、内部の臓器等を良好に観察することができる。

図３３、図３４は、提案手法１、提案手法２により生成したレンダリング像を背面から観察した画像である。提案手法１（図３３）では、肋骨は良好に除去されるが、椎骨の除去が不十分である。これに対し、提案手法２（図３４）では、椎骨も良好に除去されている。

図３５、図３６は、提案手法１、提案手法２により生成したレンダリング像を側面から観察した画像である。図３３、図３４の場合と同様に、提案手法１（図３５）では、椎骨の除去が不十分であるのに対し、提案手法２（図３６）では、椎骨も良好に除去されている。

（実施例２）
次に、腹部ＣＴ画像を対象としてレンダリング像の比較を行った。
ここで、提案手法１における楕円パラメータＰ（ｚ）（（式１１）参照）および補正倍率Ｓα（（式１２）参照）の各種パラメータを、ｘｏｆｆｓｅｔ＝－５（Ｘ方向のオフセット）、ｙｏｆｆｓｅｔ＝－４０（Ｙ方向のオフセット）、ａ＝０．７５（横方向のサイズの補正係数）、ｂ＝０．６２（縦方向のサイズの補正係数）、ｋ＝１．０５（減衰係数）、ａｍｐ＝１．０（振幅倍率）と設定した。

また、提案手法２における楕円パラメータＰ１（ｚ）、Ｐ２（ｚ）、Ｐ３（ｚ）（（式１３）（式１４）参照）および補正倍率Ｓα（（式１５）参照）の各種パラメータを、ｘｏｆｆｓｅｔ＝－５（Ｘ方向のオフセット）、ｙｏｆｆｓｅｔ＝－４０（Ｙ方向のオフセット）、ａ＝０．７５（横方向のサイズの補正係数）、ｂ＝０．６２（縦方向のサイズの補正係数）、ｋ＝１．０５（減衰係数）、ａｍｐ＝１．０（振幅倍率）、ｄｙｍ＝０．０、Ｚｔｈ＝０．０、ｄ＝０．４、ｅ＝０．２と設定した。

図３７は、従来手法により生成したレンダリング像を正面から観察した画像であり、図３８は、提案手法１により生成したレンダリング像を正面から観察した画像である。図３７（従来手法）では、肋骨の除去が不十分であり、内部の臓器等を観察することが困難であるのに対し、図３８（提案手法１）では、肋骨全体が除去され、内部の臓器等を良好に観察することができる。

図３９、図４０は、提案手法１、提案手法２により生成したレンダリング像を背面から観察した画像である。提案手法１（図３９）では、実施例１と同様に、椎骨の除去が不十分であるのに対し、提案手法２（図４０）では、椎骨も良好に除去されている。

図４１、図４２は、提案手法１、提案手法２により生成したレンダリング像を側面から観察した画像である。図３９、図４０の場合と同様に、提案手法１（図４１）では、椎骨の除去が不十分であるのに対し、提案手法２（図４２）では、椎骨も良好に除去されている。

以上のように、従来手法に比べ、提案手法１、２によれば、非観察対象である骨領域を除去し、観察対象である内部の臓器等を精度よく可視化することができる。

尚、提案手法１の方法では、体表に近い肋骨については効果的に除去することができる。しかし、胸骨や椎骨は、肋骨に比べ信号値が高く、また、体表より内側に配置されている。このため、胸骨や椎骨領域の補正倍率Ｓαは肋骨より高めに算出されるため、提案手法１のように被写体領域の外郭を単一の楕円で近似する方法では胸骨や椎骨を除去することが難しく、胸骨や椎骨に隠れた心臓や血管を良好に観察できない場合がある。一方で、胸骨や椎骨領域の補正倍率Ｓαが小さくなるように楕円の径を縮小することも考えられるが、観察対象の肺、心臓、腎臓などの内臓領域における補正倍率Ｓαも小さくなり、観察対象が部分的に除去されてしまう虞がある。

この点、提案手法２によれば、被写体領域の外郭を、２つの楕円と中央に小さな楕円を重ねた３重の楕円で近似する。前述したように、一般的な人の胸郭はＸ軸方向に２つの肺野が存在し、ほぼ中央に心臓、胸骨、椎骨が存在することから、断層画像においては提案手法１のような単一の楕円よりも提案手法２のような３つの楕円で近似するほうが適合する。このように、複数楕円を用いた提案手法２によれば、提案主法１では除去が難しい、体表の内側に位置する胸骨や椎骨も好適に除去することができる（図３４、３６、４０、４２参照）。例えば、図３５、図３６は、提案手法１、提案手法２により生成したレンダリング像を側面から観察した画像であるが、図３５では左上部（矢印で明示）に胸骨の一部が残っているのに対し、図３６の同箇所ではそれが除去できている。

以上、本実施の形態について説明した。本開示の実施の形態によれば、カラーマップＣｍａｐを参照して断層画像群Ｄｏの各画素（ｘ、ｙ、ｚ）の信号値ｖを色値及び不透明度に変換してボクセル構造体Ｖ（ボクセルデータ）を作成する際、カラーマップＣｍａｐ（ｖ、３）（オパシティカーブ）に定義された不透明度をそのまま用いるのではなく、被写体領域の外郭を近似した幾何形状（楕円）の情報に基づいて不透明度を補正する。そして、補正された不透明度を用いてボクセル構造体（ボクセルデータ）を作成し、レンダリング像を生成する。

特に、図８のステップＳ２６～Ｓ２８の処理に示したように、被写体領域の外郭を近似した複数の楕円の楕円パラメータＰ１（ｚ）、Ｐ２（ｚ）、Ｐ３（ｚ）（（式１３）参照）に基づいて算出された補正倍率Ｓα（ｘ、ｙ、ｚ）（（式１５）参照）を用いて不透明度を補正することで、所望の観察対象を選択的にかつ鮮明に可視化させたレンダリング像を得ることができる。

また、補正倍率Ｓαは、被写体領域の外郭を近似した各楕円の中心から遠ざかるにつれその倍率が小さくなるので、各楕円の中心から離れた位置にある骨領域や皮膚、外部の寝台、固定用治具等が透明化され、各楕円の中心付近に位置する臓器等が鮮明に可視化される。

また、図１４（ａ）に示したように、断層画像の被写体領域を２分割して２つの合同な楕円を配置し、中央に小さな楕円を配置した、３重の楕円形状で近似するようにすれば、体表近くに位置する肋骨だけでなく、体表内部に位置する胸骨や椎骨も効果的に除去することができる（図３４、３６、４０、４２参照）。

また、被写体領域の外郭を複数の楕円で近似する場合、楕円同士が重複する場合があるが（図１６参照）、楕円同士が重複する範囲では、各楕円の楕円方程式の値の最小値に基づいて不透明度の補正倍率Ｓαを算出する。すなわち、楕円同士が重複する範囲の補正倍率Ｓαには、各楕円から算出される補正倍率のうち、最も値が大きい補正倍率が割り当てられる。このため、重複する他の楕円の影響によって観察対象が誤って透過されてしまう虞がない。

また、各楕円の中心座標にスライス位置に応じたオフセットを加えることで（（式１３）（式１４）、図１５参照）、椎骨のＹ軸方向への位置ずれ（湾曲）に対処することができ、特に、Ｙ軸方向への位置ずれが大きくなる胸部下側の椎骨領域を良好に除去することができる。

また、中央の楕円に対しては、固有のパラメータｄ、ｅ（（式１３）参照）によって、楕円の中心位置およびサイズを調整できる。これにより、例えば、中央の楕円を観察対象となる心臓の位置やサイズに応じて細かく調整することができる。

また、補正倍率Ｓαは、被写体領域の外郭の形状を少数の変数で表現可能な楕円パラメータＰ（ｚ）（中心座標、横方向のサイズ、縦方向のサイズ）に基づいて解析的に算出されるので、複雑な処理を伴わずに、不透明度の制御を行える。

また、図１０に示したように、被写体領域の抽出精度が良好でない断層画像（例えば、被写体領域の両端が切れた断層画像）については、被写体領域を補正することで、被写体領域を精度よく求めることができる。

また、カラーマップＣｍａｐ（ｖ、３）（オパシティカーブ）を信号値に応じて調整することで、観察対象と非観察対象の領域が近接しているような場合でも、非観察対象を効果的に透過させ、観察対象を鮮明に可視化することができる（図７参照）。

また、従来のマスク処理のように、各ボクセルの可視／不可視の制御（オン／オフの制御）では、単一のボクセルに可視の内臓領域と不可視の骨領域が混在する境界領域を表現することができず、自然なレンダリング像が得られない場合があった。これに対して、本実施の形態では、マスク処理のようなオン／オフの制御でなく、被写体領域の外部に向かって不透明度を徐々に減衰させ、オンとオフの中間状態をもたせながら滑らかな制御を行うことで、自然なレンダリング像を得ることができる。

なお、本開示では不透明度の補正倍率Ｓα（ｘ、ｙ、ｚ）を格納した補正テーブルを作成するものとしたが、補正テーブルを作成する方法をとらず、ボクセル作成時に各ボクセルごとに（式１５）で都度算出する方法をとることもできる。この場合、メモリ上に補正テーブルの保存領域を確保する必要はない。

以上、添付図面を参照しながら、本開示に係る医用画像処理装置等の好適な実施形態について説明したが、本開示はかかる例に限定されない。当業者であれば、本願で開示した技術的思想の範疇内において、各種の変更例又は修正例に想到し得ることは明らかであり、それらについても当然に本開示の技術的範囲に属するものと了解される。

１： 医用画像処理装置
２１： 画像取得部
２２： 階調圧縮部
２３： 領域指定部
２４： カラーマップ取得部
２５： カラーマップ調整部
２６： 補正テーブル作成部
２６ａ： 被写体領域抽出部
２６ｂ： 被写体領域補正部
２６ｃ： 個別領域設定部
２６ｄ： 幾何情報取得部
２６ｅ： 補正倍率算出部
２７： ボクセル作成部
２８： レンダリング部
２８ａ： 第１レンダリング部
２８ｂ： 第２レンダリング部
２９： データ出力部
３０： パラメータ調整部
４０： パラメータ調整画面
Ｓα： 補正倍率、補正テーブル

Claims (31)

1. 複数の断層画像を取得する画像取得手段と、
   信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、
   前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、
   前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、
   前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出手段と、
   前記断層画像毎に、観察対象や非観察対象の位置関係に応じて、抽出された被写体領域の外郭を近似する複数の幾何形状を算出し、当該複数の幾何形状の情報を取得する幾何情報取得手段と、を備え、
   前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正する医用画像処理装置。
2. 前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する複数の幾何形状の情報に基づいて算出する補正倍率算出手段と、を更に備え、
   前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算する請求項１に記載の医用画像処理装置。
3. 前記補正倍率算出手段は、算出した補正倍率を格納する補正テーブルを作成するテーブル作成手段と、を更に備え、
   前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正テーブルを参照して得られる補正倍率を乗算する請求項２に記載の医用画像処理装置。
4. 前記補正倍率算出手段は、前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくする請求項２または請求項３に記載の医用画像処理装置。
5. 前記補正倍率算出手段は、前記複数の幾何形状が重複する範囲においては、各々の幾何形状の情報に基づいて算出される補正倍率のうち、最も値が大きい補正倍率を採用する請求項４に記載の医用画像処理装置。
6. 抽出された被写体領域に複数の個別領域を設定する個別領域設定手段と、を更に備え、
   前記幾何情報取得手段は、個別領域毎に被写体領域の外郭を幾何形状で近似し、個別領域毎に幾何形状の情報を取得する請求項１から請求項５のいずれかに記載の医用画像処理装置。
7. 前記個別領域設定手段は、抽出された被写体領域を横方向に分割した２つの個別領域を設定し、更に、被写体領域の横方向の中央付近に前記個別領域より小さい個別領域を設定する請求項６に記載の医用画像処理装置。
8. 前記幾何情報取得手段は、前記中央付近に設定した個別領域に対応する幾何形状の幾何中心に対してオフセットを加えることで、前記中央付近に設定した個別領域に対応する幾何形状の情報を補正する請求項７に記載の医用画像処理装置。
9. 前記幾何情報取得手段は、各幾何形状の幾何中心に対して、スライス位置が所定の位置から末端の方向に位置するにつれ所定の割合で増大するオフセットを加えることで、各幾何形状の情報を補正する請求項１から請求項８のいずれかに記載の医用画像処理装置。
10. 前記幾何情報取得手段は、断層画像の総数をＳｚ、スライス位置をｚ、前記所定の位置をＺｔｈとすると、ｚ≧Ｚｔｈの場合は（ｚ－Ｚｔｈ）／（Ｓｚ－１－Ｚｔｈ）に比例したオフセットを各幾何形状の幾何中心のＹ座標に加える請求項９に記載の医用画像処理装置。
11. 前記幾何形状は楕円であり、
    前記幾何情報取得手段は、前記幾何形状の情報として、楕円の中心座標、楕円の横方向のサイズ、及び楕円の縦方向のサイズを取得し、
    前記補正倍率算出手段は、前記中心座標、前記横方向のサイズ、及び前記縦方向のサイズに基づいて、前記補正倍率を算出する請求項２に記載の医用画像処理装置。
12. 前記補正倍率算出手段は、前記中心座標における前記補正倍率を最大とし、前記中心座標と前記横方向のサイズ及び前記縦方向のサイズの比が同一となる前記楕円内部の各同心楕円上の前記補正倍率を、各同心楕円の横方向のサイズ及び縦方向のサイズが大きくなるにつれ小さくする請求項１１に記載の医用画像処理装置。
13. 前記幾何情報取得手段は、前記中心座標に所定のオフセットを加算して前記中心座標を補正し、前記横方向のサイズ及び前記縦方向のサイズに所定の横方向のサイズに対する倍率及び縦方向のサイズに対する倍率を乗算して前記横方向のサイズ及び前記縦方向のサイズを補正し、
    前記補正倍率算出手段は、補正された前記中心座標、前記横方向のサイズ、及び前記縦方向のサイズに基づいて、前記補正倍率を算出する請求項１１または請求項１２に記載の医用画像処理装置。
14. 前記オフセット、前記横方向のサイズに対する倍率、及び前記縦方向のサイズに対する倍率をユーザに調整させる調整手段と、を更に備える請求項１３に記載の医用画像処理装置。
15. 被写体領域の抽出精度を所定の基準で判断し、所定の基準を満たさない前記断層画像の被写体領域を、所定の基準を満たす他の前記断層画像の被写体領域に基づいて補正する被写体領域補正手段と、を更に備える請求項１から請求項１４のいずれかに記載の医用画像処理装置。
16. 前記被写体領域補正手段は、被写体領域の横方向のサイズが画像の横方向のサイズと一致する前記断層画像における被写体領域の横方向のサイズを、被写体領域の横方向のサイズが画像の横方向のサイズより小さい他の前記断層画像における被写体領域の横方向のサイズ及び縦方向のサイズの比と同一となるように、補正する請求項１５に記載の医用画像処理装置。
17. 信号値と所定の閾値との差分値に応じて当該信号値に対応する不透明度を減衰させるように前記カラーマップを調整するカラーマップ調整手段と、を更に備え、
    前記ボクセル作成手段は、調整された前記カラーマップを参照して、ボクセルデータを作成する請求項１から請求項１６のいずれかに記載の医用画像処理装置。
18. 代表的な信号値と色値及び不透明度との対応関係に基づいて、所定範囲の信号値と色値及び不透明度との対応関係を定義するカラーマップを作成するカラーマップ作成手段と、を更に備え、
    前記カラーマップ取得手段は、作成した前記カラーマップを取得する請求項１から請求項１７のいずれかに記載の医用画像処理装置。
19. 前記ボクセル作成手段は、各ボクセルの色値に対して、当該ボクセルの不透明度に基づいて当該ボクセルにおける勾配ベクトルを算出し、勾配ベクトルとあらかじめ定義された光源ベクトルに基づいて拡散反射成分と環境光成分から構成される陰影値を算出し、前記色値に前記陰影値を乗算することにより、前記色値を改変する陰影付加手段と、を更に備える請求項１から請求項１８のいずれかに記載の医用画像処理装置。
20. 前記陰影付加手段は、前記勾配ベクトルを算出するにあたり、当該ボクセルのＸ軸方向、Ｙ軸方向、Ｚ軸方向の各々２近傍のボクセルの不透明度の差分値を前記勾配ベクトルのＸ方向成分、Ｙ方向成分、Ｚ方向成分として算出し、あらかじめ定義されたＺ軸方向変倍率に基づいて補正を施したベクトルを、当該ボクセルにおける勾配ベクトルとして算出する請求項１９に記載の医用画像処理装置。
21. 前記レンダリング手段は、
    前記ボクセルデータを生成するボリュームレンダリング像に投影変換した座標系を視点座標系とすると、視点座標系において、前記ボリュームレンダリング像の各画素よりＺ軸方向に沿って、Ｚ軸の上限値より下限値に向けて、視点座標系のボクセル座標毎に座標変換を行って前記ボクセルデータより不透明度を取得しながら、不透明ボクセルを探索し、最初に見つかった不透明ボクセルの視点座標系におけるＺ座標を、前記ボリュームレンダリング像の画素毎に記録した探索制御マスクを作成する探索制御マスク作成手段と、
    前記ボリュームレンダリング像の画素毎に、前記探索制御マスクからＺ座標を取得し、取得したＺ座標よりＺ軸の下限値に向けてＺ軸方向に沿って、所定の光強度をもつ仮想光線を照射する際、視点座標系のボクセル座標毎に座標変換を行って前記ボクセルデータより不透明度を取得し、不透明ボクセルが見つかった場合、当該ボクセル座標に対して座標変換を行って前記ボクセルデータより色値を取得し、当該ボクセルの不透明度に基づいて前記光強度を減衰させるとともに、当該ボクセルの不透明度及び色値並びに前記減衰させた光強度に基づいて累積輝度値を算出する処理を繰り返し、算出された累積輝度値に基づいて、前記ボリュームレンダリング像の当該画素に対応する画素値として与えるレイキャスティング手段と、
    を備える請求項１から請求項２０のいずれかに記載の医用画像処理装置。
22. 前記探索制御マスク作成手段及び前記レイキャスティング手段は、前記座標変換を行って前記ボクセルデータより不透明度または色値を取得する際、
    所定の回転を定義した回転行列、ＸＹＺ軸各方向のオフセット値、ＸＹＺ軸方向の拡大又は縮小倍率、Ｚ軸方向の変倍率、注視点から視点までの距離を含む前記所定の座標変換のパラメータを取得し、
    前記視点座標系のボクセルの整数値の座標を、前記パラメータに基づいて前記ボクセルデータの座標系に変換を行って、前記ボクセルデータの実数値の座標を算出し、
    算出した実数値の座標の近傍の複数の整数値の座標に対応する前記ボクセルデータの複数のボクセルを特定し、
    特定した複数のボクセルの不透明度または色値に基づいて前記ボクセルデータより取得される不透明度または色値として算出する請求項２１に記載の医用画像処理装置。
23. 前記レンダリング手段は、
    前記ボクセルデータに基づいて３Ｄテクスチャを生成する３Ｄテクスチャ生成手段と、
    前記３Ｄテクスチャに対して所定の座標変換を行って変換後３Ｄテクスチャを生成する座標変換手段と、
    ３次元空間のＸＹ座標面上の四角形をＺ軸方向に並べた積層四角形を設定する積層四角形設定手段と、
    所定の視点からＺ軸方向に平行な視線上の前記四角形のＸＹ座標に対応する前記変換後３Ｄテクスチャのボクセルの色値を前記ボクセルの不透明度に基づいて前記視点から遠い四角形の順にアルファブレンディングして取得し、前記ボリュームレンダリング像の画素値として与える画素値算出手段と、
    を備える請求項１から請求項２０のいずれかに記載の医用画像処理装置。
24. 前記座標変換手段は、
    所定の回転を定義した回転行列、視野角度、視点位置、クリッピング位置、ＸＹＺ軸各方向のオフセット値、ＸＹＺ軸方向の拡大又は縮小倍率、Ｚ軸方向の変倍率を含む所定の座標変換のパラメータを取得し、
    前記３Ｄテクスチャに対して、前記取得したパタメータを用いた前記所定の座標変換を行って前記変換後３Ｄテクスチャを生成する請求項２３に記載の医用画像処理装置。
25. 前記座標変換手段及び前記画素値算出手段は、ビデオカードに搭載されたＧＰＵ及びフレームメモリを用いて実行する請求項２４に記載の医用画像処理装置。
26. 複数の断層画像を取得する画像取得手段と、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、
    前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出手段と、
    前記断層画像毎に、抽出された被写体領域の外郭を複数の幾何形状で近似し、当該複数の幾何形状の情報を取得する幾何情報取得手段と、
    前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する複数の幾何形状の情報に基づいて算出する補正倍率算出手段と、を備え、
    前記補正倍率算出手段は、前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、前記複数の幾何形状が重複する範囲においては、各々の幾何形状の情報に基づいて算出される補正倍率のうち、最も値が大きい補正倍率を採用し、
    前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算して補正する医用画像処理装置。
27. 複数の断層画像を取得する画像取得手段と、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、
    前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出手段と、
    前記断層画像毎に、抽出された被写体領域の外郭を複数の幾何形状で近似し、当該複数の幾何形状の情報を取得する幾何情報取得手段と、を備え、
    前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正し、各ボクセルの色値に対して、当該ボクセルの不透明度に基づいて当該ボクセルにおける勾配ベクトルを算出し、勾配ベクトルとあらかじめ定義された光源ベクトルに基づいて拡散反射成分と環境光成分から構成される陰影値を算出し、前記色値に前記陰影値を乗算することにより、前記色値を改変する陰影付加手段を更に備える医用画像処理装置。
28. 前記陰影付加手段は、前記勾配ベクトルを算出するにあたり、当該ボクセルのＸ軸方向、Ｙ軸方向、Ｚ軸方向の各々２近傍のボクセルの不透明度の差分値を前記勾配ベクトルのＸ方向成分、Ｙ方向成分、Ｚ方向成分として算出し、あらかじめ定義されたＺ軸方向変倍率に基づいて補正を施したベクトルを、当該ボクセルにおける勾配ベクトルとして算出する請求項２７に記載の医用画像処理装置。
29. コンピュータの制御部が、
    複数の断層画像を取得する画像取得ステップと、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得ステップと、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成ステップと、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリングステップと、
    前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出ステップと、
    前記断層画像毎に、観察対象や非観察対象の位置関係に応じて、抽出された被写体領域の外郭を近似する複数の幾何形状を算出し、当該複数の幾何形状の情報を取得する幾何情報取得ステップと、を含み、
    前記ボクセル作成ステップは、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正する医用画像処理方法。
30. コンピュータを、
    複数の断層画像を取得する画像取得手段、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段、
    前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出手段、
    前記断層画像毎に、観察対象や非観察対象の位置関係に応じて、抽出された被写体領域の外郭を近似する複数の幾何形状を算出し、当該複数の幾何形状の情報を取得する幾何情報取得手段、として機能させ、
    前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正するプログラム。
31. コンピュータの制御部が、
    複数の断層画像を取得する画像取得ステップと、
    前記断層画像毎に、断層画像の被写体領域を抽出する被写体領域抽出ステップと、
    前記断層画像毎に、観察対象や非観察対象の位置関係に応じて、抽出された被写体領域の外郭を近似する複数の幾何形状を算出し、当該複数の幾何形状の情報を取得する幾何情報取得ステップと、
    前記断層画像の各画素の信号値に応じたオパシティカーブにより規定される不透明度を、当該断層画像に対応する複数の幾何形状の情報に基づいて補正する不透明度補正ステップと、
    前記断層画像の各画素に対応する補正後の不透明度を格納したデータを作成するデータ作成ステップと、
    を含むデータ作成方法。

Images