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

医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法 Download PDF

Info

Publication number
JP7180123B2
JP7180123B2 JP2018103767A JP2018103767A JP7180123B2 JP 7180123 B2 JP7180123 B2 JP 7180123B2 JP 2018103767 A JP2018103767 A JP 2018103767A JP 2018103767 A JP2018103767 A JP 2018103767A JP 7180123 B2 JP7180123 B2 JP 7180123B2
Authority
JP
Japan
Prior art keywords
opacity
tomographic image
voxel
pixel
image
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2018103767A
Other languages
English (en)
Other versions
JP2019205791A (ja
Inventor
敏雄 茂出木
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Dai Nippon Printing Co Ltd
Original Assignee
Dai Nippon Printing Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Dai Nippon Printing Co Ltd filed Critical Dai Nippon Printing Co Ltd
Priority to JP2018103767A priority Critical patent/JP7180123B2/ja
Publication of JP2019205791A publication Critical patent/JP2019205791A/ja
Application granted granted Critical
Publication of JP7180123B2 publication Critical patent/JP7180123B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Apparatus For Radiation Diagnosis (AREA)
  • Image Processing (AREA)
  • Image Generation (AREA)

Description

本開示は、医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法に関する。
医用画像診断において、特定の人体組織の観察に適したボリュームレンダリング像を得たい場面がある。例えば、胸部や頭部にある内臓や脳等の所定の人体組織の観察を行いたい場合に、内臓や脳等は骨に囲まれており、骨領域はむしろ診断の妨げになる。CT画像から生成されるボリュームレンダリング像では一般に骨が鮮明に描写され、内臓や血管は隠れてしまうことがあるため、可視化にあたっては切断を行うなど工夫が必要となることがある。骨領域は比較的高信号であるが、信号値分布が広範で造影血管と被るため、カラーマップを工夫して骨領域だけを透明化したオパシティカーブを設計することは容易ではなかった。
たとえば、下記の特許文献1では、CT画像の信号値に基づいて領域の抽出とラベリングを行い、ラベリングされた各肋骨領域が適正か否かを、別途構築した統計上の肋骨の体軸方向の高さ幅情報を有する統計データベースを用いた手法が開示されている。これにより、骨領域(肋骨領域)を検出して除去(透明化)させることができるが、特許文献1の方法では、肋骨に先天性の奇形があったり、骨折を伴う場合には判定誤りが発生し、適切に除去できなくなる。
また、領域拡張法(リージョングローイング法)を用いて非観察対象としたい骨領域を抽出した3次元マスクを作成し、3次元マスクを参照しながらマスク処理により非観察対象としたい骨領域を除去したボリュームレンダリング像を生成する方法が開示されている(特許文献2、3参照)。領域拡張法とは、非観察対象領域の画素を指定し、その画素を開始点(拡張開始点)として、近傍画素を次々と抽出する方法である。
特開2009-240569号公報 特許4087517号公報 特許5257958号公報
しかしながら、領域拡張法を用いた3次元マスクの作成は、ユーザによる複数の拡張開始点の指示が必須で、拡張の打ち切り段階もユーザが指示する必要があるため、自動化が難しく、ユーザごとに結果にバラつきが発生するという問題がある。また領域拡張法に限らず、3次元マスクの作成は、ユーザのスキルや経験が必要であり、その作成に時間や手間がかかる。
本開示は、前述した問題点に鑑みてなされたものであり、その目的とすることは、特定の領域を鮮明に可視化することが可能な、医用画像処理装置等を提供することである。
本開示の一実施例によると、複数の断層画像を取得する画像取得手段と、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得手段と、を備え、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正する医用画像処理装置が提供される。
本開示の一実施形態によれば、カラーマップを参照して断層画像の各画素の信号値を色及び不透明度に変換してボクセルデータを作成する際、カラーマップ(オパシティカーブ)に定義された不透明度を、断層画像の被写体領域の外郭を近似した幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正する。このように、被写体形状(幾何形状)及びスライス位置の情報に基づいて不透明度をコントロールすることで、容易に特定の領域を鮮明に可視化させたレンダリング像を得ることができる。
また、前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報、及び当該断層画像のスライス位置の情報に基づいて算出する補正倍率算出手段と、を更に備え、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算するようにしても良い。これにより、断層画像の各画素の信号値に応じた不透明度に対して、幾何形状の情報及びスライス位置の情報に基づいて算出された補正倍率を乗算することで、不透明度を補正することができる。
また、前記補正倍率算出手段は、算出した補正倍率を格納する補正テーブルを作成するテーブル作成手段と、を更に備え、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正テーブルを参照して得られる補正倍率を乗算するようにしても良い。これにより、不透明度の補正倍率を格納した補正テーブルを参照して、不透明度を補正することができる。
また、前記補正倍率算出手段は、前記補正倍率を、前記幾何形状の幾何中心(重心)から遠ざかるにつれ小さくし、かつ、スライス位置が中央から末端に位置するほど小さくするようにしても良い。これにより、被写体領域の中心から遠ざかるつれ不透明度が減衰するので、被写体領域の中心部から離れた骨領域や体表が透過され、被写体領域の中心部に近い臓器等が鮮明に可視化される。また、スライス位置が中央から遠ざかるつれ不透明度が減衰するので、特にスライス末端近くで被写体領域に占める骨領域(頭頂骨や顎骨)が増大する頭部CTにおいて、骨領域を効果的に透過させることができる。
また、前記幾何形状は楕円であり、前記幾何情報取得手段は、前記幾何形状の情報として、楕円の中心座標、楕円の横方向のサイズ、及び楕円の縦方向のサイズを取得し、前記スライス位置の情報は、スライス総数、及びスライス順位であり、前記補正倍率算出手段は、前記中心座標、前記横方向のサイズ、前記縦方向のサイズ、前記スライス総数、及び前記スライス順位に基づいて、前記補正倍率を算出するようにしても良い。これにより、被写体領域の外郭を楕円で近似し、楕円パラメータ(中心座標、横方向のサイズ、縦方向のサイズ)、及びスライスパラメータ(スライス総数、スライス順位)に基づいて、不透明度の補正倍率を算出することができる。なお本開示において、横方向とは、断層画像の座標系を規定する2つの座標軸方向(X軸方向、Y軸方向)のうち、該断層画像に映る被写体の左右軸の方向と近いほうの座標軸方向のことを指す。また、縦方向とは、断層画像の座標系を規定する2つの座標軸方向(X軸方向、Y軸方向)のうち、該断層画像に映る被写体の背腹軸(前後軸)の方向と近いほうの座標軸方向のことを指す。ここで、断層画像のX軸およびY軸は、例えば、モダリティでの撮影時に適宜設定されるが、本開示では、被写体の左右軸と方向が近いほうの座標軸をX軸、被写体の背腹軸(前後軸)と方向が近いほうの座標軸をY軸とする。すなわち、本開示では、横方向は断層画像のX軸方向に相当し、縦方向は断層画像のY軸方向に相当するものとして説明する。
また、前記補正倍率算出手段は、前記中心座標における前記補正倍率を最大とし、前記中心座標と前記横方向のサイズ及び前記縦方向のサイズの比が同一となる前記楕円内部の各同心楕円上の前記補正倍率を、各同心楕円の横方向のサイズ及び縦方向のサイズが大きくなるにつれ小さくし、かつ、前記スライス総数/2とスライス順位との差が大きくなるにつれ前記補正倍率を小さくするようにしても良い。これにより、楕円内部の不透明度の補正倍率がスライス方向を考慮して好適に算出される。
また、前記補正倍率算出手段は、前記中心座標に所定のオフセットを加算して前記中心座標を補正し、前記横方向のサイズ及び前記縦方向のサイズに所定の横方向のサイズに対する倍率及び縦方向のサイズに対する倍率を乗算して前記横方向のサイズ及び前記縦方向のサイズを補正し、前記スライス総数に所定のスライス総数に対する倍率を乗算して補正し、前記スライス順位を所定のスライスオフセットを加算して補正し、補正された前記中心座標、前記横方向のサイズ前記縦方向のサイズ、前記スライス総数、及び前記スライス順位に基づいて、前記補正倍率を算出するようにしても良い。これにより、楕円パラメータ(中心座標、横方向のサイズ、縦方向のサイズ)、及びスライスパラメータ(スライス総数、スライス順位)を補正したうえで、不透明度の補正倍率を算出することができる。
また、前記オフセット、前記横方向のサイズに対する倍率、及び前記縦方向のサイズに対する倍率、前記スライス総数に対する倍率、及び前記スライスオフセットをユーザに調整させる調整手段と、を更に備えるようにしても良い。これにより、楕円パラメータ(中心座標、横方向のサイズ、縦方向のサイズ)、及びスライスパラメータ(スライス総数、スライス順位)の補正量をユーザが調整することができる。
また、本開示の一実施態様によると、複数の断層画像を取得する画像取得手段と、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得手段と、前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報に基づいて算出する補正倍率算出手段と、を備え、前記補正倍率算出手段は、前記被写体領域の中心部から離れた骨領域や体表を透過させ、中心部に近い臓器を鮮明に可視化するように前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報に基づいて補正し、かつ前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算する医用画像処理装置が提供される。
本開示の一実施形態によれば、カラーマップを参照して断層画像の各画素の信号値を色及び不透明度に変換してボクセルデータを作成する際、カラーマップ(オパシティカーブ)に定義された不透明度を、断層画像の被写体領域の外郭を近似した幾何形状の情報に基づいて補正する。このように、被写体形状(幾何形状)に基づいて不透明度をコントロールすることで、容易に特定の領域を鮮明に可視化させたレンダリング像を得ることができる。また、補正倍率算出手段によって、前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報に基づいて算出する際に、被写体領域の中心部から離れた骨領域や体表を透過させ、中心部に近い臓器を鮮明に可視化するように前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくする。断層画像の各画素の信号値に応じた不透明度に対して、幾何形状の情報に基づいて算出された補正倍率を乗算することで、不透明度を補正することができる。また、被写体領域の中心から遠ざかるつれ不透明度が減衰するので、被写体領域の中心部から離れた骨領域や体表が透過され、被写体領域の中心部に近い臓器等が鮮明に可視化される。
また、前記補正倍率算出手段は、算出した補正倍率を格納する補正テーブルを作成するテーブル作成手段と、を更に備え、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正テーブルを参照して得られる補正倍率を乗算するようにしても良い。これにより、不透明度の補正倍率を格納した補正テーブルを参照して、不透明度を補正することができる。
また、前記幾何形状は楕円であり、前記幾何情報取得手段は、前記幾何形状の情報として、楕円の中心座標、楕円の横方向のサイズ、及び楕円の縦方向のサイズを取得し、前記補正倍率算出手段は、前記中心座標、前記横方向のサイズ、及び前記縦方向のサイズに基づいて、補正倍率を算出するようにしても良い。これにより、楕円パラメータ(中心座標、横方向のサイズ、縦方向のサイズ)に基づいて、不透明度の補正倍率を算出することができる。
また、前記補正倍率算出手段は、前記中心座標における前記補正倍率を最大とし、前記中心座標と前記横方向のサイズ及び前記縦方向のサイズの比が同一となる前記楕円内部の各同心楕円上の前記補正倍率を、各同心楕円の横方向のサイズ及び縦方向のサイズが大きくなるにつれ小さくしても良い。これにより、楕円内部の不透明度の補正倍率が好適に算出される。
また、前記補正倍率算出手段は、前記中心座標に所定のオフセットを加算して前記中心座標を補正し、前記横方向のサイズ及び前記縦方向のサイズに所定の横方向のサイズに対する倍率及び縦方向のサイズに対する倍率を乗算して前記横方向のサイズ及び前記縦方向のサイズを補正し、補正された前記中心座標、前記横方向のサイズ、及び前記縦方向のサイズに基づいて、前記補正倍率を算出するようにしても良い。これにより、楕円パラメータ(中心座標、横方向のサイズ、縦方向のサイズ)を補正したうえで、不透明度の補正倍率を算出することができる。
また、前記オフセット、前記横方向のサイズに対する倍率、前記縦方向のサイズに対する倍率をユーザに調整させる調整手段と、を更に備えるようにしても良い。これにより、楕円パラメータの補正量をユーザが調整することができる。
また、前記断層画像の幾何形状の近似精度を所定の基準で判断し、所定の基準に満たない前記断層画像の幾何形状の情報を、所定の基準を満たす他の前記断層画像の幾何形状の情報に基づいて補正する幾何情報補正手段と、を更に備えるようにしても良い。これにより、近似精度が良好でない断層画像の幾何形状の情報(幾何パラメータ)を補正することができる。
また、前記幾何形状は楕円であり、前記幾何情報取得手段は、前記幾何形状の情報として、楕円の中心座標、楕円の横方向のサイズ及び楕円の縦方向のサイズを取得し、前記幾何情報補正手段は、楕円の横方向のサイズが画像の横方向のサイズと一致する前記断層画像における楕円の横方向のサイズを、楕円の横方向のサイズが画像の横方向のサイズより小さい他の前記断層画像における楕円の横方向のサイズ及び縦方向のサイズの比と同一となるように、補正するようにしても良い。これにより、被写体領域が画像範囲に収まらない断層画像であっても、幾何パラメータ(楕円パラメータ)を良好に算出することができる。
また、信号値と所定の閾値との差分値に応じて当該信号値に対応する不透明度を減衰させるように前記カラーマップを調整するカラーマップ調整手段と、を更に備え、前記ボクセル作成手段は、調整された前記カラーマップを参照して、ボクセルデータを作成するようにしても良い。これにより、カラーマップ(オパシティカーブ)の不透明度を信号値に応じて調整することができる。
また、代表的な信号値と色値及び不透明度との対応関係に基づいて、所定範囲の信号値と色値及び不透明度との対応関係を定義するカラーマップを作成するカラーマップ作成手段と、を更に備え、前記カラーマップ取得手段は、作成した前記カラーマップを取得するようにしても良い。これにより、例えば、代表的な信号値(例えば、10箇所程度の特徴的な信号値)と当該信号値に対する色及び不透明度をユーザが設定するだけで、自動でカラーマップを作成することができる。
また、前記レンダリング手段は、前記ボクセルデータを生成するボリュームレンダリング像に投影変換した座標系を視点座標系とすると、視点座標系において、前記ボリュームレンダリング像の各画素よりZ軸方向に沿って、Z軸の上限値より下限値に向けて、視点座標系のボクセル座標毎に座標変換を行って前記ボクセルデータより不透明度を取得しながら、不透明ボクセルを探索し、最初に見つかった不透明ボクセルの視点座標系におけるZ座標を、前記ボリュームレンダリング像の画素毎に記録した探索制御マスクを作成する探索制御マスク作成手段と、前記ボリュームレンダリング像の画素毎に、前記探索制御マスクからZ座標を取得し、取得したZ座標よりZ軸の下限値に向けてZ軸方向に沿って、所定の光強度をもつ仮想光線を照射する際、視点座標系のボクセル座標毎に座標変換を行って前記ボクセルデータより不透明度を取得し、不透明ボクセルが見つかった場合、当該ボクセル座標に対して座標変換を行って前記ボクセルデータより色値を取得し、当該ボクセルの不透明度に基づいて前記光強度を減衰させるとともに、当該ボクセルの不透明度及び色値並びに前記減衰させた光強度に基づいて累積輝度値を算出する処理を繰り返し、算出された累積輝度値に基づいて、前記ボリュームレンダリング像の当該画素に対応する画素値として与えるレイキャスティング手段と、を備えるようにしても良い。これにより、レイキャスティング法によるレンダリング処理が実行される。
また、前記探索制御マスク作成手段及び前記レイキャスティング手段は、前記座標変換を行って前記ボクセルデータより不透明度または色値を取得する際、
所定の回転を定義した回転行列、XYZ軸各方向のオフセット値、XYZ軸方向の拡大又は縮小倍率、Z軸方向の変倍率、注視点から視点までの距離を含む前記所定の座標変換のパラメータを取得し、前記視点座標系のボクセルの整数値の座標を、前記パラメータに基づいて前記ボクセルデータの座標系に変換を行って、前記ボクセルデータの実数値の座標を算出し、算出した実数値の座標の近傍の複数の整数値の座標に対応する前記ボクセルデータの複数のボクセルを特定し、特定した複数のボクセルの不透明度または色値に基づいて前記ボクセルデータより取得される不透明度または色値として算出するようにしても良い。これにより、レンダリング処理中の座標変換処理が実行される。
また、前記レンダリング手段は、前記ボクセルデータに基づいて3Dテクスチャを生成する3Dテクスチャ生成手段と、前記3Dテクスチャに対して所定の座標変換を行って変換後3Dテクスチャを生成する座標変換手段と、3次元空間のXY座標面上の四角形をZ軸方向に並べた積層四角形を設定する積層四角形設定手段と、所定の視点からZ軸方向に平行な視線上の前記四角形のXY座標に対応する前記変換後3Dテクスチャのボクセルの色値を前記ボクセルの不透明度に基づいて前記視点から遠い四角形の順にアルファブレンディングして取得し、前記ボリュームレンダリング像の画素値として与える画素値算出手段と、を備えるようにしても良い。これにより、3Dテクスチャマッピング法によるレンダリング処理が実行される。
また、前記座標変換手段は、所定の回転を定義した回転行列、視野角度、視点位置、クリッピング位置、XYZ軸各方向のオフセット値、XYZ軸方向の拡大又は縮小倍率、Z軸方向の変倍率を含む所定の座標変換のパラメータを取得し、前記3Dテクスチャに対して、前記取得したパタメータを用いた前記所定の座標変換を行って前記変換後3Dテクスチャを生成するようにしても良い。これにより、レンダリング処理中の座標変換処理が実行される。
また、前記座標変換手段及び前記画素値算出手段は、ビデオカードに搭載されたGPU及びフレームメモリを用いて実行するようにしても良い。これにより、レンダリング処理中の座標変換処理及び画像値算出処理が、高速に実行される。
また、本開示の一実施形態によると、コンピュータの制御部が、複数の断層画像を取得する画像取得ステップと、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得ステップと、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成ステップと、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリングステップと、前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得ステップと、を備え、前記ボクセル作成ステップは、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正する医用画像処理方法が提供される。
本開示の一実施形態によれば、カラーマップを参照して断層画像の各画素の信号値を色及び不透明度に変換してボクセルデータを作成する際、カラーマップ(オパシティカーブ)に定義された不透明度を、断層画像の被写体領域の外郭を近似した幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正する。このように、被写体形状(幾何形状)及び当該断層画像のスライス位置の情報に基づいて不透明度をコントロールすることで、従来のように3次元マスクを作成することなく、特定の領域を鮮明に可視化させたレンダリング像を得ることができる。
また、本開示の一実施形態によると、コンピュータの制御部が、複数の断層画像を取得する画像取得ステップと、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得ステップと、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成ステップと、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリングステップと、前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得ステップと、前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報に基づいて算出する補正倍率算出ステップと、を備え、前記補正倍率算出ステップは、前記被写体領域の中心部から離れた骨領域や体表を透過させ、中心部に近い臓器を鮮明に可視化するように前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、前記ボクセル作成ステップは、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報に基づいて補正し、かつ前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出ステップにより算出された補正倍率を乗算する医用画像処理方法が提供される。
本開示の一実施形態によれば、カラーマップを参照して断層画像の各画素の信号値を色及び不透明度に変換してボクセルデータを作成する際、カラーマップ(オパシティカーブ)に定義された不透明度を、断層画像の被写体領域の外郭を近似した幾何形状の情報に基づいて補正する。このように、被写体形状(幾何形状)
に基づいて不透明度をコントロールすることで、従来のように3次元マスクを作成することなく、特定の領域を鮮明に可視化させたレンダリング像を得ることができる。また、前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報に基づいて算出する際に、被写体領域の中心部から離れた骨領域や体表を透過させ、中心部に近い臓器を鮮明に可視化するように前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくする。断層画像の各画素の信号値に応じた不透明度に対して、幾何形状の情報に基づいて算出された補正倍率を乗算することで、不透明度を補正することができる。また、被写体領域の中心から遠ざかるつれ不透明度が減衰するので、被写体領域の中心部から離れた骨領域や体表が透過され、被写体領域の中心部に近い臓器等が鮮明に可視化される。
また、本開示の一実施形態によると、コンピュータを、複数の断層画像を取得する画像取得手段、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段、前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得手段、として機能させ、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正するプログラムが提供される。
また、本開示の一実施形態によると、コンピュータを、複数の断層画像を取得する画像取得手段、信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段、前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段、前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段、前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得手段、前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報に基づいて算出する補正倍率算出手段、として機能させ、前記補正倍率算出手段は、前記被写体領域の中心部から離れた骨領域や体表を透過させ、中心部に近い臓器を鮮明に可視化するように前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報に基づいて補正し、かつ前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算するプログラムが提供される。
本開示の一実施形態に係るプログラムを汎用のコンピュータにインストールすることによって、本開示の一実施形態に係る医用画像処理装置を得ることができる。
また、本開示の一実施形態によると、コンピュータの制御部が、複数の断層画像を取得する画像取得ステップと、前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得ステップと、前記断層画像の各画素の信号値に応じたオパシティカーブにより規定される不透明度を、当該断層画像に対応する幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正する不透明度補正ステップと、前記断層画像の各画素に対応する補正後の不透明度を格納したデータを作成するデータ作成ステップと、を含むデータ作成方法が提供される。本開示の一実施形態によれば、断層画像の各画素に適用される不透明度を参照用のデータとして作成しておくことができるので、断層画像からボクセルデータを作成する際、不透明度を都度算出する必要がなくなり、レンダリング結果を得るまでの処理時間を短縮できる。
本開示により、特定の領域を鮮明に可視化させることできる。
医用画像処理装置1のハードウェア構成を示す図 ボリュームレンダリング処理の概要を示す図 医用画像処理装置1の機能構成を示す図 補正テーブル作成部26の内部構成を示す図 レンダリング部28の内部構成を示す図 医用画像処理装置1の全体動作を示すフローチャート 調整されたカラーマップCmapの性質を示す図 補正テーブル作成処理を示すフローチャート 楕円パラメータP(z)を取得する処理を説明する図 楕円パラメータP(z)の補正について説明する図 断層画像Do(z)毎に、被写体形状の外郭が楕円形状で近似されていく様子を示す図 補正倍率Sαの分布パターンを示す図 XY方向の振幅倍率と補正倍率との関係を示す図 ボクセル作成処理を示すフローチャート レイキャスティング法の概要を示す図 レイキャスティング法によるレンダリング処理を示すフローチャート 探索制御マスク算出処理を示すフローチャート 探索制御マスク算出処理を説明する図 不透明ボクセル探索処理を示すフローチャート レイキャスティング処理を示すフローチャート 3Dテクスチャマッピング法ついて説明する図 色補正処理を示すフローチャート 3Dテクスチャマッピング法によるレンダリング処理を示すフローチャート 投影画面設定処理を示すフローチャート レンダリング処理を示すフローチャート 従来手法により生成したレンダリング像(正面) 提案手法により生成したレンダリング像(正面) 提案手法により生成したレンダリング像(背面) 提案手法により生成したレンダリング像(側面) 第2の実施の形態における補正倍率Sα(x、y、z)の分布パターンを示す図 XYZ方向の振幅倍率と補正倍率との関係を示す図 第3の実施形態における医用画像処理装置1aの機能構成を示す図 第3の実施形態における医用画像処理装置1aの全体動作を示すフローチャート パラメータ調整画面40の表示例を示す図 第4の実施形態における医用画像処理装置1bの機能構成を示す図 体軸断面、冠状断面、矢状断面について説明する図 第4の実施形態における医用画像処理装置1bの全体動作を示すフローチャート MPR像の表示例を示す図 MPR像の表示例を示す図
以下図面に基づいて、本開示の実施の形態を詳細に説明する。本実施の形態では、X線CT装置により撮影されたCT画像に基づいてボリュームレンダリング像を生成する場合について説明するが、本開示はCT画像以外の医用画像(MRI画像やPET画像等)に基づいてボリュームレンダリング像を生成する場合にも適用可能である。
[第1の実施の形態]
図1は、本実施の形態における医用画像処理装置1のハードウェア構成を示すブロック図である。図1に示すように、医用画像処理装置1は、制御部11、記憶部12、メディア入出力部13、通信制御部14、入力部15、表示部16、周辺機器I/F部17等が、バス18を介して接続される汎用のコンピュータで実現される。但し、これに限ることなく、用途、目的に応じて様々な構成を採ることが可能である。
制御部11は、CPU(Central Processing Unit)、GPU(Graphics Processing Unit)、ROM(Read Only
Memory)、RAM(Random Access Memory)、フレームメモリ(Frame Memory)等によって構成される。CPU、GPUは、記憶部12、ROM、記録媒体等に格納されるプログラムをRAM、フレームメモリ上のワークメモリ領域に呼び出して実行し、バス18を介して接続された各装置を駆動制御し、医用画像処理装置1が行う後述する処理を実現する。
GPU及びフレームメモリは、医用画像処理装置1のビデオカードに搭載されている。ROMは、不揮発性メモリであり、コンピュータのブートプログラムやBIOS等のプログラム、データ等を恒久的に保持している。RAM、フレームメモリは、揮発性メモリであり、記憶部12、ROM、記録媒体等からロードしたプログラム、データ等を一時的に保持するとともに、制御部11が各種処理を行う為に使用するワークエリアを備える。
記憶部12は、HDD(Hard Disk Drive)等であり、制御部11が実行するプログラム、プログラム実行に必要なデータ、OS(Operating System)等が格納される。プログラムに関しては、OSに相当する制御プログラムや、後述する処理をコンピュータに実行させるためのアプリケーションプログラムが格納されている。これらの各プログラムコードは、制御部11により必要に応じて読み出されてRAM、フレームメモリに移され、CPU、GPUに読み出されて各種の手段として実行される。
メディア入出力部13(ドライブ装置)は、データの入出力を行い、例えば、CDドライブ(-ROM、-R、-RW等)、DVDドライブ(-ROM、-R、-RW等)等のメディア入出力装置を有する。通信制御部14は、通信制御装置、通信ポート等を有し、コンピュータとネットワーク間の通信を媒介する通信インタフェースであり、ネットワークを介して、他のコンピュータ間との通信制御を行う。ネットワークは、有線、無線を問わない。
入力部15は、データの入力を行い、例えば、キーボード、マウス等のポインティングデバイス、テンキー等の入力装置を有する。入力部15を介して、コンピュータに対して、操作指示、動作指示、データ入力等を行うことができる。表示部16は、液晶パネル等のディスプレイ装置、ディスプレイ装置と連携してコンピュータのビデオ機能を実現するための論理回路等(ビデオアダプタ等)を有する。なお、入力部15及び表示部16は、タッチパネルディスプレイのように、一体となっていてもよい。
周辺機器I/F(Interface)部17は、コンピュータに周辺機器を接続させるためのポートであり、周辺機器I/F部17を介してコンピュータは周辺機器とのデータの送受信を行う。周辺機器I/F部17は、USB(Universal Serial Bus)やIEEE1394やRS-232C等によって構成されており、通常複数の周辺機器I/Fを有する。周辺機器との接続形態は有線、無線を問わない。バス18は、各装置間の制御信号、データ信号等の授受を媒介する経路である。
医用画像処理装置1は、1台のコンピュータで構成されてもよいし、複数のコンピュータがネットワークを介して構成されてもよい。例えば、医用画像処理装置1が、サーバとクライアント端末で構成される場合、クライアント端末においてデータの入力を受け付けて、サーバが各種の処理を行い、クライアント端末が処理結果を表示するようにしてもよい。以下の説明では、簡素な構成例として、医用画像処理装置1が1台のコンピュータで構成された例を説明する。
図2は、医用画像処理装置1が実行するボリュームレンダリング処理の概要を示す図である。図2に示すように、複数の断層画像(図の例では、512×512ピクセルの370枚の胸部CT画像)に基づいてボリュームレンダリング処理を実行し、ボリュームレンダリング像として、遠位の視点から被写体を観察する全体レンダリング像(図左下)や視点を被写体内に自由に移動させて気管支や大腸等を観察する仮想内視鏡像(図右下)を生成し表示する。本実施の形態では、主に、全体レンダリング像(図左下)を生成する場合について説明する。
図3は、医用画像処理装置1の機能構成を示す図である。図3に示すように、医用画像処理装置1は、画像取得部21、階調圧縮部22、領域指定部23、カラーマップ取得部24、カラーマップ調整部25、補正テーブル作成部26、ボクセル作成部27、レンダリング部28、及びデータ出力部29を備える。
画像取得部21は、CT装置により撮影された断層画像群Doを取得する。断層画像群Doは、被写体(人体)を頭尾軸に沿って所定のスライス間隔で連続的に撮影した複数の断層画像(CT画像)からなる。各断層画像はDICOM形式の2次元の画像データである。DICOM形式は、1ファイルにヘッダ部と画像データ部を含む医療画像で一般的に用いられる画像フォーマットであり、画像撮影時のパラメータや診断情報を保存しておくことができる。
1つの断層画像は、例えば、512×512ピクセルの画像(CT画像)である。断層画像の各画素には、信号値vが付与されており、CT画像の場合、信号値vはCT値である。本実施の形態では、信号値v(CT値)は16ビット(-32768≦v≦32767)のデータとする(但し、信号値vのビット数は特に限定されない)。
CT値は、水を基準として表現した組織のX線減弱係数であり、CT値により組織や病変の種類等が判断できるようになっている(単位はHU(Hounsfield Unit))。CT値は、水と空気のX線減弱係数で標準化されており、水のCT値は0、空気のCT値を-1000である。また、脂肪のCT値は-120~-100程度であり、通常組織のCT値は0~120程度であり、骨のCT値は1000前後を示す。
断層画像群Doは、XYの2次元データである断層画像をZ軸方向に積層したものであり、XYZの3次元データとして表現可能である。例えば断層画像群Doは、以下のように、XYZの3次元データ(Do(x、y、z))として定義される。尚、本開示において、X軸は人体の左右軸、Y軸は人体の背腹軸(前後軸)、Z軸は人体の頭尾軸(上下軸)に相当するものとする。XY軸は、例えば、CT装置での撮影時に設定される。また、本開示において、X軸方向を横方向、Y軸方向を縦方向と呼ぶ場合がある。
(式1)
-32768≦Do(x、y、z)≦32767
0≦x≦Sx-1、0≦y≦Sy-1、0≦z≦Sz-1
解像度:Rxy、解像度Rz
(式1)において、Sxは断層画像の横方向(X軸方向)の画素数、Syは断層画像の縦方向(Y軸方向)の画素数、Szはスライス枚数(断層画像の枚数)を表す。RxyはX軸方向及びY軸方向の断層画像の解像度であり、画素の間隔の逆数、すなわち単位距離あたりの画素数を示す。Rzはスライスの解像度であり、断層画像のスライス間隔(例えば、0.5mmや1mm)の逆数、すなわち単位距離あたりのスライス枚数を表す。
尚、断層画像群Doのうち、所定のz番目(0≦z≦Sz-1)の断層画像を表す場合は、「断層画像Do(z)」と表記することがある。
階調圧縮部22は、画像取得部21により取得した断層画像群Doの信号値を、必要に応じて8ビットに階調圧縮する(後述の(式6)参照)。8ビットに階調圧縮された断層画像Do(x、y、z)は以下のように定義される。
(式2)
0≦Do(x、y、z)≦255
0≦x≦Sx-1、0≦y≦Sy-1、0≦z≦Sz-1
解像度:Rxy、解像度Rz
このように階調圧縮をすることで、断層画像群を保持するためのメモリ容量を半分に抑えることができる。たとえ信号値の階調が16ビットあっても、カラーマップCmapにより、変換される色値(具体的にはRGB値)及び不透明度の階調はディスプレイの階調に合わせて8ビットに制限されるため、階調圧縮に伴う画質劣化は殆ど生じない。
領域指定部23は、ユーザから被写体の関心領域(ROI:Region of
Interest)の指定を受け付ける。例えば、関心領域ROIを直方体で指定する場合は、以下のように、X方向ROI(X方向の開区間)、Y方向ROI(Y方向の開区間)、Z方向ROI(Z方向の開区間)を指定する。
(式3)
X方向ROI=(Xs、Xe)
Y方向ROI=(Ys、Ye)
Z方向ROI=(Zs、Ze)
尚、関心領域を指定しない場合は、Xs=0、Xe=Sx-1、Ys=0、Ye=Sy-1、Zs=0、Ze=Sz-1となる。
カラーマップ取得部24は、断層画像群Doに適用するカラーマップCmapを取得する。カラーマップ取得部24は、取得したカラーマップCmapを更に調整する場合には、カラーマップ調整部25に出力し、カラーマップCmapの調整を行わない場合には、ボクセル作成部27に出力する。
カラーマップCmapは、信号値vと色値(具体的にはRGB値)及び不透明度(α値)との対応関係を定義するものであり、信号値vを24ビットの色(RGB値)及び8ビットの不透明度(α値)に変換する関数(実体的には2次元のデータテーブル)として表現可能である。例えば、16ビットの断層画像群Do((式1)参照)に適用されるカラーマップCmapは次のように定義される。
(式4)
0≦Cmap(v、n)≦255
-32768≦v≦32767
n=0(R)、1(G)、2(B)、3(α)
また、8ビットの断層画像群Do((式2)参照)に適用されるカラーマップCmap(v、n)は次のように定義される。
(式5)
0≦Cmap(v、n)≦255
0≦v≦255
n=0(R)、1(G)、2(B)、3(α)
(式4)、(式5)に示すカラーマップCmap(v、n)(0≦n≦3)のうち、特にカラーマップCmap(v、n)(0≦n≦2)は、信号値vを色値(RGB値)に変換する関数に相当する。本開示では、信号値vを色値に変換するカラーマップCmap(v、n)(0≦n≦2)を「カラーパレット」と称する。
また、(式4)、(式5)に示すカラーマップCmap(v、n)(0≦n≦3)のうち、特にカラーマップCmap(v、3)は、信号値vを不透明度に変換する関数に相当し、一般的に「オパシティカーブ」と呼ばれる。
カラーマップCmap(v、3)(オパシティカーブ)は、例えば、通常組織(16ビットの信号値v=0~120程度)の不透明度が255に設定(不透明度=255は不透明(光がすべて反射)であることを示す)され、骨領域(16ビットの信号値v=1000前後)の不透明度が1~254の中間値に設定(不透明度=1~254は半透明(入射光の一部が反射され、その他は透過)であることを示す)され、その他の組織等の部位を不透明度が0に設定(不透明度=0は透明(入射光の全てが透過)であることを示す)される。
カラーマップ調整部25は、必要に応じて、カラーマップCmap(v、n)のうち、信号値vと不透明度との対応関係を定義するCmap(v、3)(オパシティカーブともいう)を信号値に応じて調整し(後述の(式7)(式8)参照)、調整したカラーマップCmapをボクセル作成部27に出力する。
補正テーブル作成部26は、断層画像群Doの各画素に対応する不透明度の補正倍率を格納した補正テーブルSαを作成し、ボクセル作成部27に出力する。
補正テーブル作成部26は、図3に示すように、幾何情報取得部26a、幾何情報補正部26b、及び補正倍率算出部26cから構成される。
幾何情報取得部26aは、断層画像Do(z)毎に、断層画像Do(z)の被写体領域の外郭を所定の幾何形状で近似し、幾何形状の情報(幾何パラメータP(z)、後述の(式9)参照)を取得する。
幾何情報補正部26bは、幾何形状の近似精度を所定の基準で判断し、近似精度が良好でない断層画像Do(z)の幾何パラメータP(z)を、近似精度が良好な他の断層画像Do(z)の幾何パラメータP(z)に基づいて補正する(後述の(式10)参照)。
補正倍率算出部26cは、断層画像Do(z)毎に、断層画像Do(z)に対応する幾何パラメータP(z)に基づいて、断層画像Do(z)の各画素(x、y、z)に対応する不透明度の補正倍率Sα(x、y、z)を算出し、補正倍率Sα(x、y、z)を格納した3次元の補正テーブルSαを作成する(後述の(式11)参照)。
この際、補正倍率Sα(x、y、z)は、断層画像Do(z)の被写体領域の外郭を近似した幾何形状の幾何中心(後述する実施形態では、楕円の中心)から遠ざかるにつれその倍率が小さくなるように算出される。
尚、補正倍率Sα(x、y、z)は3次元のデータ配列であり、そのままデータテーブル(補正テーブル)として扱えるため、補正倍率と補正テーブルは実質的に同一である。このため、補正倍率と補正テーブルは同一の符号「Sα」で表す。
ボクセル作成部27は、断層画像群Do、カラーマップ取得部24により取得されたカラーマップCmap(又はカラーマップ調整部25により調整されたカラーマップCmap)、補正テーブル作成部26により作成された補正テーブルSα、及び領域指定部23により指定された関心領域ROIに基づいて、ボクセル構造体V(ボクセルデータ)を作成し、レンダリング部28に出力する。
具体的には、ボクセル作成部27は、カラーマップCmapを参照することで、断層画像群Doの各画素(x、y、z)の信号値vを、信号値vに応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセル構造体V(ボクセルデータ)を作成する。
特に本開示では、ボクセル作成部27は、カラーマップCmapを参照して断層画像群Doの各画素(x、y、z)の信号値vを色値及び不透明度に変換してボクセル構造体V(ボクセルデータ)を作成する際、カラーマップCmap(v、3)(オパシティカーブ)に定義された不透明度をそのまま用いるのではなく、各画素(x、y、z)に対応する補正倍率Sα(x、y、z)を乗算することで不透明度を補正する(後述の(式13)参照)。
前述したように、補正倍率Sα(x、y、z)は、被写体領域の外郭を近似した幾何形状の幾何中心(重心)から遠ざかるにつれ倍率が小さくなるので、被写体領域の中心から離れたボクセルほど不透明度が減衰する。これにより、骨領域(胸骨や頭蓋骨)や皮膚、外部の寝台、固定用治具等が透明化され、被写体領域の中心付近にある臓器等が鮮明に可視化されたレンダリング像を得ることができる。
レンダリング部28は、ボクセル作成部27により作成されたボクセル構造体V(ボクセルデータ)に基づいて、ボリュームレンダリング処理を実行し、ボリュームレンダリング像(3次元画像)を生成し表示する。レンダリング部28は、図5に示すように、第1レンダリング部28a及び第2レンダリング部28bを含む。
第1レンダリング部28aは、CPUによりレンダリング処理を実行するレンダリング部であり、第2レンダリング部28bは、コンピュータグラフィックスAPIのオープン標準規格であるOpenGL(Open Graphics Library)を用いてGPUによりレンダリング処理を実行するレンダリング部である。レンダリング処理は、第1レンダリング部28a又は第2レンダリング部28bにより実行される。
データ出力部29は、ユーザからデータの保存操作を受け付けることによって、各種データを出力して記憶部12に保存する。例えば、データ出力部29は、レンダリング部28により生成し表示されたボリュームレンダリング像(3次元画像)を出力して保存する。
また、データ出力部29は、ボクセル構造体V(ボクセルデータ)を作成する際に参照可能なデータ(例えば、補正後の不透明度を格納した不透明度テーブルTα等)を出力しても良い。これにより、同じ断層画像群Doに基づいて再度レンダリングを行う際、いくつかの処理(例えば、補正テーブルSαの作成処理等)を省略することができ、レンダリング像を得るまでの処理時間を短縮できる。
次に、図6~図25を参照しながら、医用画像処理装置1の動作について説明する。
図6は、医用画像処理装置1の全体動作を示すフローチャートである。
まず、医用画像処理装置1の制御部11(画像取得部21)は、レンダリング対象となるDICOM形式の複数の断層画像(断層画像群)Doを取得する(図6のステップS1)。断層画像群Doの取得方法は任意であるが、例えば、断層画像群Doが予め記憶部12に記憶されている場合は、制御部11は、記憶部12から断層画像群Doを読込み取得することができる。
また、断層画像群DoがUSB等の記憶媒体に記憶されている場合は、制御部11は、周辺機器I/F部17に接続した当該記憶媒体から断層画像群Doを読込み取得することができる。或いは、サーバに断層画像群Doが記憶されている場合は、制御部11は、通信制御部14を介して当該サーバから断層画像群Doをダウンロードして取得することができる。
続いて、制御部11(階調圧縮部22)は、必要に応じて、断層画像群Doに対して階調圧縮処理を施す(図6のステップS2)。具体的には、制御部11は、以下の手順で16ビットの断層画像群Do((式1)参照)を8ビットの断層画像群Do((式2)参照)に変換する。
(式6)
(1)断層画像群DoのSz/2番目の中間スライスDo(z/2)における全ての画素の最小値Dmin、最大値Dmaxを算出する。
(2)下限値Lmin=(Dmax-Dmin)・γ+Dmin、上限値Lmax=(Dmax-Dmin)・(1-γ)+Dminを設定する。ここで、γは階調圧縮画像のコントラスト調整幅で、0に近いほどコントラストが増大する(但し、輝度が小さくなる)。通常はγ=0.1に設定する。
(3)断層画像群Doを以下の計算式で8ビットの断層画像群Doに変換する(下限値Lmin~上限値Lmaxの範囲で256段階に圧縮する)。
Do(x、y、z)
=(Do(x、y、z)-Lmin)・255/(Lmax-Lmin)
但し、Do(x、y、z)>255の場合はDo(x、y、z)=255、Do(x、y、z)<0の場合はDo(x、y、z)=0に飽和させる。
続いて、制御部11(領域指定部23)は、必要に応じて、ユーザから被写体の関心領域ROIの指定を受け付ける(図6のステップS3)。例えば、前記した(式3)のように、関心領域ROIを直方体等で指定する。
続いて、制御部11(カラーマップ取得部24)は、断層画像群Doに適用するカラーマップCmapを取得する(ステップS4)。ここで、断層画像群Doの信号値が16ビットの場合((式1)参照)は、16ビット対応のカラーマップCmap((式4)参照)を取得し、断層画像群Doの信号値が8ビットに圧縮されている場合((式2)参照)は、8ビット対応のカラーマップCmap((式5)参照)を取得する。
カラーマップCmapの取得方法は特に限定されず、例えば、制御部11は、予め用意されている複数のカラーマップCmapの中から、断層画像群Doに適用するカラーマップCmapをユーザに選択させることで取得することができる。また、制御部11は、前回レンダリング時に使用したカラーマップCmapを記憶部12から読込み、断層画像群Doに適用するカラーマップCmapとして取得しても良い。
また、制御部11は、カラーマップ作成画面(不図示)においてユーザが作成したカラーマップCmapを取得しても良い。この場合、制御部11は、カラーマップ作成画面において、代表的な信号値v(例えば、10箇所程度の特徴的な信号値)と当該信号値vに対する色値(RGB値)及び不透明度(α値)をユーザに設定させ、ユーザが設定した代表的な信号値vと色値及び不透明度との対応関係に基づいて所定範囲の信号値vを色値及び不透明度に変換するカラーマップCmapを自動生成する。
続いて、制御部11は、ステップS4において取得したカラーマップCmapを更に調整する場合(ステップS5;「Yes」)、ステップS6へ移行し、カラーマップCmapの調整を行う。具体的には、制御部11(カラーマップ調整部25)は、以下のように、カラーマップCmap(v、n)のうち、信号値vと不透明度との対応関係を定義するカラーマップCmap(v、3)(オパシティカーブ)を調整する。
断層画像群Doの信号値が16ビットの場合((式1)参照)は、制御部11は、次の手順で、16ビット対応のカラーマップCmap(v、3)((式4)参照)を調整する。
(式7)
(1)断層画像群Doの中で、Sz/2番目の中間スライスDo(z/2)における全ての画素の最小値Dmin、最大値Dmaxを算出する。
(2)下限値Lmin=(Dmax-Dmin)・γ+Dmin、上限値Lmax=(Dmax-Dmin)・(1-γ)+Dminを設定する。ここで、γは階調圧縮画像のコントラスト調整幅で、0に近いほどコントラストが増大する(但し、輝度が小さくなる)。通常はγ=0.1に設定する。
(3)16ビットの信号値vを以下の計算式で8ビットの信号値v8に変換する(下限値Lmin~上限値Lmaxの範囲で256段階に圧縮する)。
v8=(v-Lmin)・255/(Lmax-Lmin)
但し、v8>255の場合はv8=255、v8<0の場合はv8=0に飽和させる。
(4)次のように、16ビット対応のカラーマップCmap(v、3)((式4)参照)を調整する。
v8>Svの場合:Cmap(v、3)←Cmap(v、3)(v8-Sv)δv8≦Svの場合:調整しない
尚、Cmap(v、n)(0≦Cmap(v、n)≦255、-32768≦v≦32767、n=0(R)、1(G)、2(B)、3(α))であり、8ビット換算の閾値をSvとし、減衰係数をδとする(δは負の実数で、通常はδ=-0.5)。閾値Svは、不透明度を減衰させる信号値の閾値であり、骨領域を減衰させる場合には、Sv=128程度に設定する。
一方、断層画像群Doの信号値が8ビットの場合((式2)参照)は、制御部11は、次のように、8ビット対応のカラーマップCmap((式5)参照)を調整する。
(式8)
v>Svの場合:Cmap(v、3)←Cmap(v、3)(v-Sv)δ
v≦Svの場合:調整しない
尚、Cmap(v、n)(0≦Cmap(v、n)≦255、0≦v≦255、n=0(R)、1(G)、2(B)、3(α))であり、閾値をSvとし、減衰係数をδとする(δは負の実数で、通常はδ=-0.5)。閾値Svは、不透明度を減衰させる信号値の閾値であり、骨領域を減衰させる場合には、Sv=128程度に設定する。
図7(a)は、(式7)においてCmap(v、3)(オパシティカーブ)に乗算される(v8-Sv)δの減衰特性を示す図である。図に示すように、信号値v8>閾値Svの場合、信号値v8と閾値Svとの差分値(v8-Sv)に応じた減衰比率(v8-Sv)δで不透明度が減衰する。例えば、図7(b)に示すように、内臓領域と骨領域の中心からの距離に差が無い場合、両者の不透明度に対して後述のステップS7の補正テーブルに基づく補正を加えても、内臓領域は骨領域に被ってしまう。そこで同図に示すように、内臓領域と骨領域の信号値の差を利用し、閾値Svを内臓領域と骨領域の信号値の境界に設定することで、骨領域の不透明度を効果的に減衰させることができ、内臓領域が鮮明に可視化される。
一方、カラーマップCmapを調整しない場合には(ステップS5;「No」)、制御部11は、上記したカラーマップCmapの調整を省略し、ステップS7へ移行する。
図8~図13を参照して、ステップS7における補正テーブルSαを作成する処理について説明する。
図8は、補正テーブルSαを作成する処理を示すフローチャートである。
まず、制御部11(幾何情報取得部26a)は、断層画像Do(z)(0≦z≦Sz-1)毎に、断層画像Do(z)の被写体形状の外郭を所定の幾何形状で近似し、幾何パラメータP(z)を取得する(図8のステップS21)。本実施の形態では、被写体形状の外郭を楕円で近似し、幾何パラメータP(z)として楕円のパラメータ(楕円パラメータP(z))を取得するものとする。
図9は、断層画像Do(z)の楕円パラメータP(z)を取得する処理を説明する図である。
制御部11は、図9(a)に示すように、断層画像Do(z)に対して、信号値が所定の閾値以上の画素を「白」、閾値未満の画素を「黒」で塗りつぶして2値化し、被写体領域を抽出する。図9(a)の例では、2値化された断層画像Do(z)の「白」の領域が被写体領域である。閾値としては、例えば、空気以外を被写体領域として抽出する場合には、信号値vが16ビットの場合は「-700」、信号値vが8ビットの場合には「10」などに設定する。
続いて、制御部11は、図9(b)に示すように、抽出した被写体領域(白領域)を全て含む最小の矩形領域R(z)を特定し、図9(c)のように、矩形領域R(z)に内接する楕円E(z)を、被写体形状の外郭を近似する楕円として求める。そして、制御部11は、図9(c)に示すように、楕円E(z)の楕円パラメータP(z)として、楕円E(z)の中心座標(Cx(z)、Cy(z))、横方向のサイズW(z)、及び縦方向のサイズH(z)を以下のように算出する。尚、図9(c)に示すように、矩形領域R(z)の最小座標を(Xmin、Ymin)、最大座標を(Xmax、Ymax)とする。
(式9)
中心座標Cx(z)=(Xmax+Xmin)/2
Cy(z)=(Ymax+Ymin)/2
横方向のサイズW(z)=Xmax-Xmin+1
縦方向のサイズH(z)=Ymax-Ymin+1
制御部11は、全ての断層画像Do(z)(0≦z≦Sz-1)について、被写体領域の外郭を近似する楕円E(z)を求め、楕円パラメータP(z)を算出する。
続いて、制御部11(幾何情報補正部26b)は、幾何形状(楕円)の近似精度が良好でない断層画像Do(z)について、幾何パラメータP(z)(楕円パラメータP(z))の補正を行う(図8のステップS22)。
図10は、断層画像Do(z)の楕円パラメータP(z)の補正について説明する図である。図10(a)の断層画像Do(z)は、楕円近似が良好に行えない断層画像の例である。図10(a)の断層画像は、被写体領域(幅方向)の両端が切れており、被写体領域が画像範囲に収まっていない(撮影時の条件によってこのような画像が得られる場合がある)。このような場合、被写体領域を正確に楕円近似できない。
すなわち、図10(b)に示すように、断層画像から得られる楕円E(z)の横方向のサイズW(z)が、断層画像の横方向のサイズSxと同一となり(W(z)=Sx)、実際の被写体領域の幅(Sxより大きい幅)を正確に反映したものでなくなる。このため、本実施の形態では、被写体領域(幅方向)の両端が切れた断層画像Do(z)について横方向のサイズW(z)の補正を行う。
具体的には、制御部11は、まず、ステップS21により取得された楕円の横方向のサイズW(z)がW(z)=Sxとなる1以上の断層画像Do(z)を、被写体領域の両端が切れた断層画像(幾何形状の近似精度が良好でない断層画像)として特定する。
次に、制御部11は、特定された断層画像Do(z)毎に、断層画像Do(z)とスライス位置が最も近く、かつ、ステップS21により取得された楕円の横方向のサイズW(z)がW(z)<Sxを満たす断層画像(「断層画像Do(z’)」と表記)を特定する。
断層画像Do(z’)は、楕円E(z’)の横方向のサイズW(z’)が画像サイズSxより小さいため、被写体領域(幅方向)が画像範囲に収まっている断層画像(幾何形状の近似精度が良好な断層画像)である。そして、制御部11は、被写体領域の両端が切れた断層画像Do(z)から得られる楕円E(z)の横方向のサイズW(z)を、断層画像Do(z’) から得られる楕円E(z’)の横方向のサイズW(z’)及び縦方向のサイズH(z’)の比と同一となるように、以下の式で補正する。
(式10)
W(z)=H(z)・W(z’)/H(z’)
これにより、図10(c)に示すように、横方向のサイズW(z)が、真の被写体領域の外郭を近似するような楕円E(z)の横方向のサイズW(z)に補正され、被写体領域(幅方向)が画像範囲に収まらない場合でも、楕円パラメータP(z)を良好に算出することができる。
以上、図8のステップS21及びS22の処理により、図11に示すように、断層画像Do(z)毎に、被写体形状の外郭が楕円形状で近似され、全ての断層画像Do(z)(0≦z≦Sz-1)について楕円パラメータP(z)(0≦z≦Sz-1)が得られる。
そして、制御部11(補正倍率算出部26c)は、断層画像Do(z)毎に、断層画像Do(z)に対応する幾何パラメータP(z)(楕円パラメータP(z))に基づいて、断層画像Do(z)の各画素(x、y、z)に対応する不透明度の補正倍率Sα(x、y、z)を算出し、補正テーブルSαを作成する(図8のステップS23)。
例えば、z番目の断層画像Do(z)の各画素(x、y、z)の補正倍率Sα(x、y、z)は、以下のように算出される。
(式11)
Sα(x、y、z)=1-k・r(x、y、z)・amp
但し、0≦Sα(x、y、z)≦1
r(x、y、z)={(x-Cx(z)-xoffset)/(a・W(z)/2)}+{(y-Cy(z)-yoffset)/(b・H(z)/2)}
(式11)において、kは減衰係数であり初期値は1.0である。k=1.0の場合、骨及び寝台が不可視となる。aは横方向のサイズW(z)の補正係数、bは縦方向のサイズH(z)の補正係数である。補正しない場合はa=b=1.0に設定する。
xoffset、yoffsetはXY面におけるオフセット値であり、初期値はxoffset=yoffset=0(単位:ボクセル)である。
図12は、(式11)により算出される補正倍率Sα(x、y、z)の分布パターン(同心楕円分布)を示す図である。図に示すように、楕円の中心(Cx(z)、Cy(z))の補正倍率Sαは1(最大)となり、楕円内部は中心から遠ざかるにつれ補正倍率Sαが小さくなる。具体的には、楕円の中心(Cx(z)、Cy(z))と横方向のサイズW(z)及び縦方向のサイズH(z)の比が同一となる各同心楕円上の補正倍率Sαが、各同心楕円の径(横方向のサイズ及び縦方向のサイズ)が大きい程小さくなる(径の2乗に反比例して小さくなる)。楕円外部の補正倍率Sαは一律に0とする。これにより、被写体領域の中心から外側にいくにつれ不透明度が減衰するので、骨領域(胸骨や頭蓋骨)や皮膚、外部の寝台、固定用治具等が透明化され、被写体領域の中心付近にある臓器等が鮮明に可視化される。
尚、(式11)のampはXY平面の振幅倍率であり、初期値はamp=1.0である。図13に示すように、被写体中心部の補正倍率を平坦にし、飽和領域を広くしたい場合は、amp>1.0に設定する(但し、陰影感が弱くなる)。
制御部11は、全ての断層画像Do(z)について、補正倍率Sα(x、y、z)(0≦x≦Sx-1、0≦y≦Sy-1、0≦z≦Sz-1)を算出し、算出された補正倍率Sα(x、y、z)を格納する3次元の補正テーブルSαを作成する。尚、前述したように、補正倍率Sα(x、y、z)は3次元のデータ配列であり、そのままデータテーブル(補正テーブル)として扱えるため、補正倍率と補正テーブルは実質的に同一である。
以上、図8~図13を参照して、補正テーブルSαを作成する処理について説明した。
図6のフローチャートの説明に戻る。
制御部11(ボクセル作成部27)は、ステップS1において取得された断層画像群Do(又はステップS2において階調圧縮された断層画像群Do)、ステップS4において取得されたカラーマップCmap(又はステップS6において調整されたカラーマップCmap)、ステップS7において作成された補正テーブルSα、及びステップS3において指定された関心領域ROIに基づいて、ボクセル構造体V(ボクセルデータ)を作成する(図6のステップS8)。
本実施の形態では、ボクセル構造体V(ボクセルデータ)は、24ビットの色値(RGB値)を保持するボクセル構造体Vcと、8ビットの不透明度(α値)を保持するボクセル構造体Vαとの2つのボクセル構造体から構成されるものとする。このようにボクセル構造体を2つに分けて作成することで、色値及び不透明度に関するボクセル処理を個別に実行できるため、メモリアクセス等の向上が図られる。
図14は、ボクセル構造体V(Vc、Vα)を作成する処理を示すフローチャートである。
まず、制御部11は、カラーマップCmap(v、n)(0≦n≦2)(=カラーパレット)、及び関心領域ROIに基づいて、24ビットのRGB値を保持するボクセル構造体Vc(x、y、z、n)(0≦V(x、y、z、n)≦255、0≦n≦2、0≦x≦Sx-1、0≦y≦Sy-1、0≦z≦Sz-1;解像度:Rxy、Rz)を以下のように作成する(図14のステップS801)。
(式12)
(1)Xs<x<Xe、Ys<y<Ye、及びZs<z<Zeを全て満たす場合(関心領域ROI内の場合)
Vc(x、y、z、n)=Cmap(Do(x、y、z)、n) (0≦n≦2)
(2)x≦Xs、x≧Xe、y≦Ys、y≧Ye、z≦Zs、又はz≧Zeのいずれかを満たす場合(関心領域ROI外の場合)
Vc(x、y、z、n)=0 (0≦n≦2)
次に、制御部11は、カラーマップCmap(v、3)(=オパシティカーブ)、補正テーブルSα、及び関心領域ROIに基づいて、8ビットの不透明度を保持するボクセル構造体Vα(x、y、z)(0≦Vα(x、y、z)≦255、0≦x≦Sx-1、0≦y≦Sy-1、0≦z≦Sz-1;解像度:Rxy、Rz)を以下のように作成する(図14のステップS802)。
(式13)
(1)Xs<x<Xe、Ys<y<Ye、及びZs<z<Zeを全て満たす場合(関心領域ROI内の場合)
Vα(x、y、z)=Cmap(Do(x、y、z)、3)
Vα(x、y、z)=Vα(x、y、z)・Sα(x、y、z)
(2)x≦Xs、x≧Xe、y≦Ys、y≧Ye、z≦Zs、又はz≧Zeのいずれかを満たす場合(関心領域ROI外の場合)
Vα(x、y、z)=0
(式13)のVα(x、y、z)=Vα(x、y、z)・Sα(x、y、z)は、本開示の要旨に係る処理であり、これにより、カラーマップCmapに定義された不透明度に対して補正倍率Sαが乗算され、不透明度が補正される。特に、補正倍率Sα(x、y、z)は、被写体領域の外郭を近似した楕円の中心から遠ざかるにつれその倍率が小さくなるので、被写体領域の中心から離れたボクセルほど不透明が減衰する。これにより、骨領域(胸骨や頭蓋骨)や皮膚、外部の寝台、固定用治具等が透明化され、被写体領域の中心付近にある臓器等が鮮明に可視化されたレンダリング像を得ることができるようになる。
尚、メモリアクセスの向上等を目的に、32ビットのボクセル構造体Vを色値を保持するボクセル構造体Vcと不透明度を保持するボクセル構造体Vαとに分けているが、ボクセル構造体を分けなくてもよい。この場合、制御部11は、以下のように、32ビットの色値及び不透明度を保持する1つのボクセル構造体V(x、y、z、n)(0≦V(x、y、z、n)≦255、0≦n≦3、0≦x≦Sx-1、0≦y≦Sy-1、0≦z≦Sz-1;解像度:Rxy、Rz)を作成する。
(式14)
(1)Xs<x<Xe、Ys<y<Ye、及びZs<z<Zeを全て満たす場合(関心領域ROI内の場合)
V(x、y、z、n)=Cmap(Do(x、y、z)、n) (0≦n≦3)
V(x、y、z、3)=V(x、y、z、3)・Sα(x、y、z)
(2)x≦Xs、x≧Xe、y≦Ys、y≧Ye、z≦Zs、又はz≧Zeのいずれかを満たす場合(関心領域ROI外の場合)
V(x、y、z、n)=0 (0≦n≦3)
上記した(式12)~(式14)に示すボクセル構造体を作成する処理は、断層画像群Doの信号値が16ビット((式1)参照)の場合は、16ビット対応のカラーマップCmap((式4)参照)を用いて実行され、断層画像群Doの信号値が8ビット((式2)参照)に階調圧縮されている場合は、8ビット対応のカラーマップCmap((式5)参照)を用いて実行される。
制御部11は、ボクセル構造体V(Vc、Vα)を作成した後、必要に応じてボクセル構造体V(Vc、Vα)に対して、公知のスムーシング処理を施したり、公知の陰影処理を施しても良い。
続いて、制御部11(第1レンダリング部28a又は第2レンダリング部28b)は、ボクセル構造体Vに基づいて、ボリュームレンダリング像を生成する(図6のステップS9)。以降、CPU(第1レンダリング部28a)により実行されるレイキャスティング法によるレンダリング処理と、GPU(第2レンダリング部28b)により実行される3Dテクスチャマッピング法によるレンダリング処理について順に説明する。
最初に、図15~図20を参照しながら、レイキャスティング法によるレンダリング処理について説明する。
図15は、レイキャスティング法の概要を示す図である。図15(a)に示すように、投影面(レンダリング像)の各点からボクセルデータに対して任意の方向にレイを発射し、最深点(レイ方向で光が届く最も深い点)よりその軌跡を逆にたどり、投影面上の点の輝度値を算出する。
また、図15(b)に示すように、計算を簡単にするため、ボクセルデータに座標変換を施し、投影面の各点よりZ軸負方向の最深点までレイを進め、同時に累積輝度計算を行う。
図16は、レンダリング処理を示すフローチャートである。まず、制御部11は、座標変換パラメータを設定する(ステップS31)。後述する探索制御マスクの算出(ステップS32)、レイキャスティング処理(ステップS33)の各処理で、視点座標系からボクセル座標系への座標変換処理を行う。そのため、制御部11は、以下のように、座標変換処理に共通する座標変換パラメータを事前に設定・算出しておく。
<座標変換パラメータ>
回転パラメータ行列R:
R=[R11 R12 R13;
R21 R22 R23;
R31 R32 R33]
(ボクセル座標系から視点座標系への座標変換を行うための3×3の行列の逆行列、GUI側はボクセル座標系から視点座標系への座標変換を指示するが、レンダリング側は視点座標系からボクセル座標系に座標変換を行う。)
X軸方向のオフセットXoff
Y軸方向のオフセットYoff
Z軸方向のオフセットZoff
X軸方向のROI:Xs-Xe(0≦Xs、Xe≦Sx-1)
Y軸方向のROI:Ys-Ye(0≦Ys、Ye≦Sy-1)
Z軸方向のROI:Zs-Ze(0≦Zs、Ze≦Sz-1)
透視変換パラメータ:原点からZ軸方向の距離Dist(正の実数値、平行投影の場合はDist=0)
拡大縮小倍率Scale(XYZ軸方向で同一)
Z方向変倍率Scz=Rxy/Rz
Z方向変倍後のZ方向サイズSz’=Sz・Scz
座標変換サブサンプル・オフセット:X軸方向dx、Y軸方向dy、Z軸方向dz
制御部11は、上記した回転パラメータ行列Rに対して、初期値として単位行列を設定する。すなわち、R11=1、R12=0、R13=0、R21=0、R22=1、R23=0、R31=0、R32=0、R33=1と設定する。
そして、GUIの指示に従い、X軸中心回転Rx、Y軸中心回転Ry、Z軸中心回転Rz(角度単位:ラジアン)のいずれかを逐次指定し、以下のように、各々回転行列Aを生成して回転パラメータ行列Rに右から乗算して、回転パラメータ行列Rを更新する。これにより、GUIの指示により生成されるボクセル座標系から視点座標系への回転行列の逆行列が算出される。
回転行列Aを
A=[A11 A12 A13;
A21 A22 A23;
A31 A32 A33]
とすると、
X軸中心回転Rxの場合の回転行列Aの各要素は、
A11=1、A12=0、A13=0
A21=0、A22=cosRx、A23=sinRx
A31=0、A32=sinRx、A33=cosRx
Y軸中心回転Ryの場合の回転行列Aの各要素は、
A11=cosRy、A12=0、A13=sinRy
A21=0、A22=1、A23=0
A31=-sinRy、A32=0、A33=cosRy
Z軸中心回転Rzの場合の回転行列Aの各要素は、
A11=cosRz、A12=sinRz、A13=0
A21=-sinRz、A22=cosRz、A23=0
A31=0、A32=0、A33=1
となる。
回転パラメータ行列Rは、R←R×Aと更新される。
以上の座標変換処理は、探索制御マスクの算出(ステップS32)、レイキャスティング処理(ステップS33)の各処理の中で逐次実行される。
次に、図17のフローチャートを参照して、図16のステップS32における探索制御マスクの算出処理について説明する。
探索制御マスクM(x、y、l)(0≦x≦Sx-1、0≦y≦Sy-1、0≦l≦L-1)は、算出するレンダリング像の画素(x、y)毎に先頭の不透明ボクセルのZ座標を記憶したものである。座標変換時にサブサンプリング(ボリュームレンダリング像のジャギーを除去するため、アンチエリアシングともいわれる)を行う場合は、制御部11は、以下の処理をL回(0≦l≦L-1)実行する。
まず、制御部11は、処理回数l=0とし、サブサンプル・オフセット値をdx=dy=dz=0と初期化する(ステップS41)。
続いて、制御部11は、x=2m及びy=2n(m、nは整数)の値をとる画素(x、y)に対して、Zs=Sz’-1から座標変換を行いながら先頭の不透明ボクセルを探索し、そのZ座標をM(x、y、l)に記録する。x=Sx-1又はy=Sy-1のいずれかの値をもつ画素(x、y)に対しても同様の探索を行う(ステップS42)。これにより、図18(a)のように、縦横偶数番目の画素について不透明ボクセルのZ座標が算出される。不透明度ボクセルを探索する処理は、図19にて後述する。
続いて、制御部11は、x=2m+1(x<Sx-1)及びy=2n(m、nは整数)の値をとる画素(x、y)に対して、M(x-1、y、l)=M(x+1、y、l)の場合、M(x、y、l)=M(x-1、y、l)を与え、M(x-1、y、l)≠M(x+1、y、l)の場合、上記同様Zs=Sz’-1から座標変換を行いながら先頭の不透明ボクセルを探索(図19参照)し、そのZ座標をM(x、y、l)に記録する(ステップS43)。これにより、図18(b)のように、X方向に奇数番目の画素について不透明ボクセルのZ座標が算出される。隣接する左右画素の不透明ボクセルのZ座標が同一の場合、そのZ座標で補間する(不透明ボクセルの探索は行わない)。
続いて、制御部11は、全てのx(0≦x≦Sx-1)及びy=2n+1(nは整数)の値をとる画素(x、y)に対して、M(x、y-1、l)=M(x、y+1、l)の場合、M(x、y、l)=M(x、y-1、l)を与え、M(x、y-1、l)≠M(x、y+1、l)の場合、上記同様Zs=Sz’-1から座標変換を行いながら先頭の不透明ボクセルを探索(図19参照)し、そのZ座標をM(x、y、l)に記録する(ステップS44)。これにより、図18(c)のように、Y方向に奇数番目の画素について不透明ボクセルのZ座標が算出される。隣接する上下画素の不透明ボクセルのZ座標が同一の場合、そのZ座標で補間する(不透明ボクセルの探索は行わない)。
続いて、制御部11は、処理回数をl←l+1に更新し、サブサンプル・オフセット値を、dx←dx+1/L、dy←dy+1/L、dz←dz+1/Lに更新する(ステップS45)。
制御部11は、処理回数l=L-1となるまで(ステップS46;Yes)、上記したステップS42~S45の処理を繰り返す。
次に、図19のフローチャートを参照して、探索制御マスク算出処理(図17)において実行される不透明度ボクセルの探索処理について説明する。
まず、制御部11は、zをz=Zsと初期化し、探索対象画素(x、y)を入力する(ステップS51)。続いて、制御部11は、3次元座標(x、y、z)に対して座標変換を行い、ボクセルα値(n=3)を算出する(ステップS52)。座標変換を行い、ボクセルα値を算出する処理(座標変換処理)は後述する。
ステップS52において算出されたボクセルα値がα=0の場合、制御部11は、z←z-m(例えばm=8)に更新する(ステップS53)。そして、z<0の場合、制御部11は、M(x、y、l)=-1に設定し(ステップS54)、処理を終了する。z≧0の場合、制御部11は、ステップS52に戻る。
一方、ステップS52において算出されたボクセルα値がα>0の場合、制御部11は、まず、iをi=0と初期化する(ステップS55)。続いて、制御部11は、i←i+1に更新し(ステップS56)、i<mかつz+i<Zsの場合、3次元座標(x、y、z+i)に対して座標変換を行い、ボクセルα値を算出する(ステップS57)。座標変換を行い、ボクセルα値を算出する処理(座標変換処理)は後述する。
一方、i≧mまたはz+i≧Zsの場合、制御部11は、M(x、y、l)=z+i-1に設定し(ステップS58)、処理を終了する。
また、ステップS57において算出されたボクセルα値がα=0の場合、制御部11は、M(x、y、l)=z+i-1に設定し(ステップS58)、処理を終了する。
一方、ステップS57において算出されたボクセルα値がα>0の場合、ステップS56に戻り、ステップS58において不透明ボクセルのZ座標が設定されるまで、ステップS56~S57の処理を繰り返す。
不透明度ボクセル探索処理(図19)のステップS52及びステップS57において実行される座標変換処理について説明する。
座標変換処理は、視点座標系をボクセル座標系に変換する処理であり、GUI側の変換処理とは逆になる。GUI側では関心領域ROIによるクリッピング、スケーリング、Z方向変倍処理、オフセット(XYZ方向同時)、回転、透視変換の順に行うものと仮定し、制御部11は、与えられた視点座標系の3次元座標値(x、y、z)(整数値)に対応する8ビットのボクセル構造体Vα(x、y、z)又は24ビットのボクセル構造体Vc(x、y、z、n)の実数の座標値(xr、yr、zr)を以下のように算出する。
制御部11は、視点座標系の座標値(x、y、z)(整数値)を次のように実数値(xx、yy、zz)に変換する。
xx=x-Sx/2+dx
yy=y-Sy/2+dy
zz=z-Sz’/2+dz
続いて、原点からZ軸方向の距離Dist>0の場合、制御部11は、以下のように透視変換を行う。
zz’=Dist・zz/(Dist+zz)
xx’=xx・(Dist-zz’)/Dist
yy’=yy・(Dist-zz’)/Dist
透視変換後の(xx’、yy’、zz’)を(xx、yy、zz)とする。
続いて、制御部11は、以下のように、回転パラメータ行列Rを用いて回転処理を行う。
xx’=R11・xx+R12・yy+R13・zz
yy’=R21・xx+R22・yy+R23・zz
zz’=R31・xx+R32・yy+R33・zz
回転処理後の(xx’、yy’、zz’)を(xx、yy、zz)とする。
制御部11は、スケーリング、Z方向変倍処理、及びオフセット(XYZ方向同時)を同時に行い、次のようにボクセル構造体の対応するボクセルの座標値(xr、yr、zr)(実数値)を取得する。
xr=xx/Scale+Sx/2-Xoff
yr=yy/Scale+Sy/2-Yoff
zr=zz/Scale/Scz+Sz/2-Zoff
続いて、制御部11は、算出したボクセルの座標値(xr、yr、zr)(実数値)に対して、小数点以下を切り捨て整数化した座標値を(xi、yi、zi)(整数値)とし、切り捨てた小数点以下の端数を(wx、wy、wz)(0≦wx、wy、wz<1)とする(すなわち、xr=xi+wx、yr=yi+wy、zr=zi+wz)。
そして、制御部11は、X軸・Y軸・Z軸方向のROIを考慮して、視点座標系の3次元座標値(x、y、z)(整数値)に対応するボクセルのα値(Vα)を次のように3通りの場合に分けて決定する。
1)xi<Xs、xi>Xe、yi<Ys、yi>Ye、zi<Zs、又はzi>Zeのいずれかを満たす場合(クリッピング範囲)
Vα=0
2)上記1)の条件を満たさない場合において、xi+1>Xe、yi+1>Ye、又はzi+1>Zeのいずれかを満たす場合(補間しない)
Vα=Vα(xi、yi、zi)
3)上記1)2)の条件を満たさない場合(補間する)
Vα=(1-wz)(1-wy)(1-wx)・Vα(xi、yi、zi)+(1-wz)(1-wy)・wx・Vα(xi+1、yi、zi)+(1-wz)・wy・(1-wx)・Vα(xi、yi+1、zi)+(1-wz)・wy・wx・Vα(xi+1、yi+1、zi)+wz・(1-wy)(1-wx)・Vα(xi、yi、zi+1)+wz・(1-wy)・wx・Vα(xi+1、yi、zi+1)+wz・wy・(1-wx)・Vα(xi、yi+1、zi+1)+wz・wy・wx・Vα(xi+1、yi+1、zi+1)
尚、後述するレイキャスティング処理(図20参照)のようにボクセルのRGB値を取得する場合は、制御部11は、視点座標系の3次元座標値(x、y、z)(整数値)に対応するボクセルのRGB値(Vc(n))(0≦n≦2)を次のように3通りの場合に分けて決定する。
1)xi<Xs、xi>Xe、yi<Ys、yi>Ye、zi<Zs、又はzi>Zeのいずれかを満たす場合(クリッピング範囲)
Vc(n)=0 (0≦n≦2)
2)上記1)の条件を満たさない場合において、xi+1>Xe、yi+1>Ye、又はzi+1>Zeのいずれかを満たす場合(補間しない)
Vc(n)=Vc(xi、yi、zi、n) (0≦n≦2)
3)上記1)2)の条件を満たさない場合(補間する)
Vc(n)=(1-wz)(1-wy)(1-wx)・Vc(xi、yi、zi、n)+(1-wz)(1-wy)・wx・Vc(xi+1、yi、zi、n)+(1-wz)・wy・(1-wx)・Vc(xi、yi+1、zi、n)+(1-wz)・wy・wx・Vc(xi+1、yi+1、zi、n)+wz・(1-wy)(1-wx)・Vc(xi、yi、zi+1、n)+wz・(1-wy)・wx・Vc(xi+1、yi、zi+1、n)+wz・wy・(1-wx)・Vc(xi、yi+1、zi+1、n)+wz・wy・wx・Vc(xi+1、yi+1、zi+1、n) (0≦n≦2)
次に、図20のフローチャートを参照して、図16のステップS33におけるレイキャスティング処理について説明する。
制御部11は、まず、生成する24ビット(RGB)のレンダリング画像Image(x、y、n)の初期値を全て0に設定する(Image(x、y、n)=0、n=0(R)、1(G)、2(B))。そして、サブサンプル回数Lとして、各2次元座標(x、y)(0≦x≦Sx-1、0≦y≦Sy-1)に対して、以下の処理を実行する。
まず、制御部11は、光源Light(n)(n=0(R)、1(G)、2(B)、0≦Light(n)≦255)を設定し、サブサンプル・オフセット値をdx=dy=dz=0、処理回数l=0に初期化する(ステップS61)。
続いて、制御部11は、z=M(x、y、l)とし(ステップS62)、z<0の場合、ステップS70へ移行し、画素(x、y)におけるRGB画素値を決定する。
一方、z≧0の場合、制御部11は、仮想光強度Trans=1.0、累積輝度値Energy(n)=0.0(0≦n≦2)に初期化し(ステップS63)、3次元座標(x、y、z)を座標変換してボクセルα値(Vα)を取得し、Alpha=Vα/255を算出する(ステップS64)。
ステップS64において算出したAlphaがAlpha=0かつz=0の場合、ステップS70へ移行し、画素(x、y)におけるRGB画素値を決定する。
一方、ステップS64において算出したAlphaがAlpha>0の場合、制御部11は、3次元座標(x、y、z)を座標変換してボクセルRGB値Vc(n)を取得する(ステップS67)。ボクセルのRGB値Vc(n)を取得は、前述した座標変換処理により行われる。
続いて、制御部11は、累積輝度をEnergy(n)=Energy(n)+Trans/Alpha・Vc(n)/255、透過光強度をTrans=Trans・(1.0-Alpha)と更新する(ステップS68)。Alpha=1.0又はTrans<0.001の場合、ステップS70へ移行し、画素(x、y)におけるRGB画素値を決定する。それ以外の場合、制御部11は、z=z-1に更新し(ステップS69)、z≧0の場合、ステップS64に戻り、z≦0の場合、ステップS70へ移行し、画素(x、y)におけるRGB画素値を決定する。
一方、ステップS64において算出したAlphaがAlpha=0かつz>0の場合、制御部11は、Zs=z-1から座標変換して先頭の不透明ボクセルのZ座標z’を探索する(ステップS65)。不透明ボクセルの探索は、前述した不透明ボクセル探索処理(図19参照)により実行される。
探索されたz’がz’<0の場合、ステップS70へ移行し、画素(x、y)におけるRGB画素値を決定する。z’≧0の場合、制御部11は、z=z'に更新し(ステップS66)、ステップS64の処理に戻る。
ステップS70において、制御部11は、画素(x、y)におけるRGB画素値を次のように決定する。
Image(x、y、n)=Image(x、y、n)+k・Energy(n)・Light(n)/L
ここで、kは強度倍率であり、初期値はk=1.0に設定される。
続いて、制御部11は、処理回数をl←l+1に更新し、サブサンプル・オフセット値を、dx←dx+1/L、dy←dy+1/L、dz←dz+1/Lに更新する(ステップS71)。制御部11は、処理回数l=L-1となるまで(ステップS72;Yes)、上記したステップS62~S71の処理を繰り返し、最終的なレンダリング画像Image(x、y、n)(n≦0≦2)を得る。
以上、図15~図20を参照しながら、CPU(第1レンダリング部28a)により実行されるレイキャスティング法によるレンダリング処理について説明した。
次に、図21~図25を参照しながら、GPU(第2レンダリング部28b)により実行される3Dテクスチャマッピング法によるレンダリング処理について説明する。
GPUによるレンダリングはOpenGLシステムを用いて実行するため、制御部11は、作成したボクセル構造体Vを3DテクスチャとしてOpenGLシステムへ登録する。具体的には、制御部11は、24ビットのボクセル構造体Vc((式12)参照)と8ビットのボクセル構造体Vα(式(13)参照)を32ビットのボクセル構造体V((式14)参照)に統合し、glTexImage3D関数を用いて3Dテクスチャを生成し、3DテクスチャをOpenGLシステムに登録する。
また、制御部11は、図21に示すように、3Dテクスチャを構成する各スライスを四角形で定義し、3次元空間のXY座標面上のZ軸方向に並べた積層四角形を設定し、テクスチャ座標と四角形座標との対応付けを行い、3Dテクスチャを各四角形に貼り付ける。3Dテクスチャマッピング法では、四角形に貼り付けた3Dテクスチャマップを視線方向に合わせて座標変換することで、任意方向のレンダリング像を生成する。
尚、GPU(第2レンダリング部28b)によるレンダリング処理では、レンダリング時にモアレが発生する場合がある。制御部11は、レンダリング時のモアレの発生を抑えるため、次のように、関心領域ROI境界面の透明ボクセルの色補正処理を事前に行う。
図22は色補正処理のフローチャートである。制御部11は、境界面ボクセル(x、y、z)を抽出する(ステップS81)。境界面ボクセル(x、y、z)とは、x=Xs、x=Xe、y=Ys、y=Ye、z=Zs、又はz=Zeのいずれかを満たすボクセルである。
続いて、制御部11は、抽出した境界面ボクセルの27近傍ボクセルより、不透明ボクセルを抽出し、抽出した不透明ボクセルの平均色(Ra、Ga、Ba)を算出する(ステップS82)。不透明ボクセルとは、Vα(x+i、y+j、z+k)>0(-1≦i≦1、-1≦j≦1、-1≦k≦1)の条件を満たすボクセルである。尚、不透明ボクセルが抽出されない場合は、制御部11は、Ra=Ga=Ba=0とする。
そして、制御部11は、以下のように、境界面ボクセルの不透明度を0に設定し、RGB値を算出した不透明ボクセルの平均色に設定する(ステップS83)。
Vα(x、y、z)=0
Vc(x、y、z、1)=Ra、Vc(x、y、z、2)=Ga、Vc(x、y、z、3)=Ba
制御部11は、全ての境界面ボクセル(x、y、z)について、上記した色補正処理(ステップS81~83)を実行する。
次に、図23のフローチャートを参照して、3Dテクスチャマッピング法によるレンダリング処理の流れを説明する。まず、制御部11は、OpenGLの透視変換パラメータを設定する(ステップS91)。仮想内視鏡モードでレンダリングを行う場合は、制御部11は、OpenGLシステムに対して、次の透視変換パラメータを設定する。
<透視変換パラメータ>
・カメラの視野角度(単位:度、平行投影の場合は0、焦点距離に相当)
・視点座標(カメラの位置、通常はZ軸上に設定)
・注視点座標(被写体の凝視点、通常はZ軸上で原点z=0に固定)
・近方クリッピング位置(視点からの距離、視点より若干距離を置く)
・遠方クリッピング位置(視点からの距離、通常は1000など無限大に固定しクリッピングしない)
一方、通常のレンダリングを行う際は、カメラの視野角度を0に設定し、平行投影モードに設定する。
続いて、制御部11は、OpenGLの投影画面設定を行う(ステップS92)。
図24は、投影画面設定の処理を示すフローチャートである。図24に示すように、制御部11は、レンダリング画像のスクリーンサイズ(縦横画素数、縦横アスペクト比率)を設定し(ステップS101)、平行投影(通常の外観レンダリング)又は透視投影(内視鏡モード)の設定を行う(ステップS102)。そして、透視投影が設定された場合は、制御部11は、透視投影パラメータ(カメラの視野角度(焦点距離)、視点座標、注視点座標、近方クリッピング位置、遠方クリッピング位置)の設定を行う(ステップS103)。
図23のフローチャートの説明に戻る。
続いて、制御部11は、以下のように、OpenGLの座標変換パラメータ(回転、スケール、移動、Z方向変倍パラメータRxy/Rz)を設定する(ステップS93)。
<OpenGLの座標変換パラメータ>
回転パラメータ行列R:
R=[R11 R12 R13 R14;
R21 R22 R23 R24;
R31 R32 R33 R34;
R41 R42 R43 R44]
(ボクセル座標系から視点座標系への座標変換を行うための4×4の行列で、レンダリング側もGUI側と同一方向の座標変換をOpenGLに指示する。実際は3×3の座標変換行列で定義できるが、OpenGLの仕様上、4×4の座標変換行列に拡張)
X軸方向のオフセットXoff
Y軸方向のオフセットYoff
Z軸方向のオフセットZoff
拡大縮小倍率Scale(XYZ軸方向で同一)
Z方向変倍率Scz=Rxy/Rz
Z方向変倍後のZ方向サイズSz’=Sz・Scz
制御部11は、上記した回転パラメータ行列Rに対して、初期値として単位行列を設定する。すなわち、R11=1、R12=0、R13=0、R14=0、R21=0、R22=1、R23=0、R24=0、R31=0、R32=0、R33=1、R41=0、R42=0、R43=0、R44=1と設定する。
そして、GUIの指示に従い、X軸中心回転Rx、Y軸中心回転Ry、Z軸中心回転Rz(角度単位:ラジアン)のいずれかを逐次指定し、以下のように、各々回転行列Aを生成して回転パラメータ行列Rに左から乗算して、回転パラメータ行列Rを更新する。
回転行列Aを
A=[A11 A12 A13 A14;
A21 A22 A23 A24;
A31 A32 A33 A34;
A41 A42 A43 A44]
とすると、
X軸中心回転Rxの場合の回転行列Aの各要素は、
A11=1、A12=0、A13=0、A14=0
A21=0、A22=cosRx、A23=sinRx、A24=0
A31=0、A32=sinRx、A33=cosRx、A34=0
A41=0、A42=0、A43=0、A44=1
Y軸中心回転Ryの場合の回転行列Aの各要素は、
A11=cosRy、A12=0、A13=sinRy、A14=0
A21=0、A22=1、A23=0、A24=0
A31=-sinRy、A32=0、A33=cosRy、A34=0
A41=0、A42=0、A43=0、A44=1
Z軸中心回転Rzの場合の回転行列Aの各要素は、
A11=cosRz、A12=sinRz、A13=0、A14=0
A21=-sinRz、A22=cosRz、A23=0、A24=0
A31=0、A32=0、A33=1、A34=0
A41=0、A42=0、A43=0、A44=1
となる。
回転パラメータ行列Rは、R←A×Rと更新される。
図23のフローチャートの説明に戻る。
制御部11は、OpenGLによるレンダリング処理を実行する(ステップS94)。
図25は、OpenGLによるレンダリング処理を示すフローチャートである。図25に示すように、制御部11は、拡大縮小倍率Scale(XYZ軸方向で同一)、及びZ方向変倍率Sczに基づいて3Dテクスチャマップに対してスケーリング、及びZ方向変倍処理を施し(ステップS111)、回転パラメータ行列Rに基づいて3Dテスチャマップに対して回転処理を施し(ステップS112)、X軸方向のオフセットXoff、Y軸方向のオフセットYoff、及びZ軸方向のオフセットZoffに基づいて3Dテスチャマップに対してオフセット処理を施すことで(ステップS113)、3Dテクスチャマップに対して所定の座標変換を行う。
続いて、制御部11は、3次元空間のXY座標面上の四角形をZ軸方向に並べた積層四角形をOpenGLで設定し、3Dテスクチャマップの対応付けを行う(ステップS114、図21参照)。そして、制御部11は、OpenGLによるレンダリング処理を実行し、積層四角形に3Dテクスチャマップを貼り付けながらスキャンコンバージョンし、フレームメモリ上で奥からZ方向(視点方向)にアルファブレンディング処理を行う(ステップS115)。すなわち、所定の視点からZ軸方向に平行な視線上の四角形のXY座標に対応する座標変換後の3Dテスクチャマップのボクセルの色値(RGB値)を、ボクセルの不透明度(α値)に基づいて視点から遠い四角形の順にアルファブレンディングして取得し、ボリュームレンダリング像の画素値として与える。
以上、図21~図25を参照しながら、GPU(第2レンダリング部28b)により実行される3Dテクスチャマッピング法によるレンダリング処理について説明した。
制御部11は、CPU(第1レンダリング部28a)により実行されるレイキャスティング法によるレンダリング処理(図16)、又はGPU(第2レンダリング部28b)により実行される3Dテクスチャマッピング法によるレンダリング処理(図23)によりボリュームレンダリング像を生成し、生成したボリュームレンダリング像を表示部16に表示する。
図6のフローチャートの説明に戻る。
制御部11(データ出力部29)は、ユーザからデータの保存操作を受け付けることによって、各種データを出力して記憶部12に保存することができる(図6のステップS10)。例えば、制御部11は、ステップS9において生成し表示されたボリュームレンダリング像(3次元画像)を出力して保存する。
また、制御部11は、ボクセル構造体V(ボクセルデータ)を作成する際に参照可能なデータ(参照データ)を出力しても良い。これにより、同じ断層画像群Doに基づいて再度レンダリングを行う際、いくつかの処理を省略することができる。参照データとしては、例えば、以下のものが挙げられる。
<参照データ1>
ボクセル構造体V(Vc、Vα)の作成((式12)、(式13)参照)に使用されたカラーマップCmap(v、n)(0≦n≦3)のデータ
<参照データ2>
補正テーブルSα((式11)参照)の作成に使用された楕円パラメータP(z)((式9)、(式10)参照)のデータ
<参照データ3>
ボクセル構造体Vαの作成((式13)参照)に使用された補正テーブルSα((式11)参照)のデータ
<参照データ4>
断層画像群Doの各画素(x、y、z)に適用される補正後の不透明度を格納した不透明度テーブルTα(x、y、z)(3次元α値分布)
ここで、不透明度テーブルTαは、以下のように定義される。
(式15)
Tα(x、y、z)=Cmap(Do(x、y、z)、3)・Sα(x、y、z)
(0≦Tα(x、y、z)≦255、0≦x≦Sx-1、0≦y≦Sy-1、0≦z≦Sz-1;解像度:Rxy、Rz)
<参照データ5>
ボクセル構造体V(Vc、Vα)の作成((式12)、(式13)参照)に使用されたカラーマップCmap(v、n)(0≦n≦3)と、補正テーブルSα(式11参照)の作成に使用された楕円パラメータP(z)((式9)、(式10)参照)を含むデータ
<参照データ6>
ボクセル構造体V(Vc、Vα)の作成((式12)、(式13)参照)に使用されたカラーマップCmap(v、n)(0≦n≦3)と、ボクセル構造体Vαの作成((式13)参照)に使用された補正テーブルSα((式11)参照)を含むデータ
<参照データ7>
ボクセル構造体Vcの作成((式12)参照)に使用されたカラーマップCmap(v、n)(0≦n≦2)(カラーパレット)と、不透明度テーブルTα(x、y、z)を含むデータ
以上のような参照データを保存しておくことで、同じ断層画像群Doに基づいて再度レンダリングを行う際に、いくつかの処理を省略でき、レンダリング像を得るまでの処理時間を短縮できる。
例えば、参照データ1(カラーマップCmap)を保存しておけば、図6のステップS5~S7(カラーマップの取得・調整処理)を省略できる。
また、参照データ2(楕円パラメータP(z))を保存しておけば、図8のステップS21、S22(補正テーブル作成処理の一部)を省略できる。
また、参照データ3(補正テーブルSα)を保存しておけば、図6のステップS4(補正テーブル作成処理)を省略できる。
また、参照データ4(不透明度テーブルTα)を保存しておけば、図6のステップS4(補正テーブル作成処理)を省略できるとともに、(式13)の不透明度の補正処理(Vα(x、y、z)=Vα(x、y、z)・Sα(x、y、z))を省略できる。すなわち、不透明度テーブルTαに格納されている不透明度をボクセル構造体Vαに設定するだけでよい(Vα(x、y、z)=Tα(x、y、z))。
また、参照データ5(カラーマップCmap+楕円パラメータP(z))を保存しておけば、図6のステップS5~S7(カラーマップの取得・調整処理)を省略できるとともに、図8のステップS21、S22(補正テーブル作成処理の一部)を省略できる。
また、参照データ6(カラーマップCmap+補正テーブルSα)を保存しておけば、図6のステップS4~S7(補正テーブルの作成処理、及びカラーマップの取得・調整処理)を省略できる。
また、参照データ7(カラーマップCmap(カラーパレット)+不透明度テーブルTα)を保存しておけば、図6のステップS4~S7(補正テーブルの作成処理、及びカラーマップの取得・調整処理)を省略できるとともに、(式13)の不透明度の補正処理(Vα(x、y、z)=Vα(x、y、z)・Sα(x、y、z))を省略できる。
以上、図6~図25を参照しながら、医用画像処理装置1の動作について説明した。
[実施例]
最後に、従来手法と提案手法により生成されるレンダリング像の比較を行う。従来手法とは、図6のステップS4において取得したカラーマップCmapを断層画像群Doにそのまま適用してボクセル構造体V(Vc、Vα)を作成し、レンダリング像を生成する手法である。すなわち、従来手法では、以下の(式16)に基づいてボクセル構造体V(Vc、Vα)を作成する。
(式16)
(1)Xs<x<Xe、Ys<y<Ye、及びZs<z<Zeを全て満たす場合(関心領域ROI内の場合)
Vc(x、y、z、n)=Cmap(Do(x、y、z)、n)(0≦n≦2)
Vα(x、y、z)=Cmap(Do(x、y、z)、3)
(2)x≦Xs、x≧Xe、y≦Ys、y≧Ye、z≦Zs、又はz≧Zeのいずれかを満たす場合(関心領域ROI外の場合)
Vc(x、y、z、n)=0 (0≦n≦2)
Vα(x、y、z)=0
提案手法は、本実施の形態に基づいてレンダリング像を生成する方法(図6の処理によりレンダリング像を生成する方法)である。ここで、補正倍率Sα((式11)参照)を算出する際のパラメータを、k=1.1(減衰係数)、amp=1.0(振幅倍率)、xoffset=2(X方向のオフセット)、yoffset=5(Y方向のオフセット)、a=0.58(横方向のサイズの補正係数)、b=0.58(縦方向のサイズの補正係数)と設定した。
また、提案手法、従来手法の双方ともカラーマップCmapの調整(図6のステップS6の処理)は行わないものとし、双方とも第1レンダリング部28aを用いてレンダリングを実行した(第2レンダリング部28bを用いても殆ど同様なレンダリング像が得られる)。
図26は、従来手法により生成したレンダリング像を正面から観察した画像であり、図27は、提案手法により生成したレンダリング像を正面から観察した画像である。図26(従来手法)では、肋骨及び胸骨の透明化が不十分であり、内部の臓器等を観察することが困難であるのに対し、図27(提案手法)では、胸骨は多少残るものの肋骨全体が透明化され、内部の臓器等を良好に観察することができる。図28、図29は、それぞれ提案手法により生成したレンダリング像を背面、側面から観察した画像である。図28、図29に示すように、椎骨は残るものの、図27と同様に肋骨が良好に透明化されていることが分かる。
以上、本開示の第1の実施の形態について説明した。第1の実施の形態によれば、カラーマップCmapを参照して断層画像群Doの各画素(x、y、z)の信号値vを色値及び不透明度に変換してボクセル構造体V(ボクセルデータ)を作成する際、カラーマップCmap(v、3)(オパシティカーブ)に定義された不透明度をそのまま用いるのではなく、被写体領域の外郭を近似した楕円の情報(楕円パラメータP(z))に基づいて算出された補正倍率Sα(x、y、z)を乗算することで不透明度を補正する((式13)参照)。そして、補正された不透明度を用いてボクセル構造体(ボクセルデータ)を作成し、レンダリング像を生成する。これにより、従来のように3次元マスクを作成することなく、特定の領域を鮮明に可視化させたレンダリング像を得ることができる。従来の3次元マスクの代わりに、本開示では不透明度の補正倍率Sα(x、y、z)に関する補正テーブルの作成を必要とするが、これらは前述のように事前に補正テーブルを作成する方法をとらず、ボクセル作成時に各ボクセルごとに(式11)で都度算出する方法をとることもでき、3次元マスクのようにメモリ上に確保する必要はない。
また、補正倍率Sαは、被写体領域の外郭を近似した楕円の中心から遠ざかるにつれその倍率が小さくなるので、被写体領域の中心から外側にいくにつれ不透明度が減衰し、骨領域(胸骨や頭蓋骨)や皮膚、外部の寝台、固定用治具等が透明化され、被写体領域の中心付近にある臓器等が鮮明に可視化される。
また、補正倍率Sαは、被写体領域の外郭の形状を少数の変数で表現可能な楕円パラメータP(z)(中心座標、横方向のサイズ、縦方向のサイズ)に基づいて解析的に算出されるので、複雑な処理を伴わずに、不透明度の制御を行える。
また、楕円の近似精度が良好でない断層画像Do(z)の楕円パラメータP(z)を、近似精度が良好な他の断層画像Do(z)の楕円パラメータP(z)に基づいて補正することで、撮影条件等によって楕円近似が良好に行えない場合でも、楕円パラメータP(z)を良好に算出することができる。
また、カラーマップCmap(v、3)(オパシティカーブ)を信号値に応じて調整することで、骨領域をより効果的に減衰させることができる(図7参照)。
また、従来のマスク処理のように、各ボクセルの可視/不可視の制御(オン/オフの制御)では、単一のボクセルに可視の内臓領域と不可視の骨領域が混在する境界領域を表現することができず、自然なレンダリング像が得られない場合があった。これに対して、本実施の形態では、マスク処理のようなオン/オフの制御でなく、被写体領域の外部に向かって不透明度を徐々に減衰させ、オンとオフの中間状態をもたせながら滑らかな制御を行うことで、自然なレンダリング像を得ることができる。
[第2の実施の形態]
頭部は球形状に近いため、胸部に比べ、被写体領域に占める骨領域の面積比がスライス方向(Z軸方向)に大きく変化する。特にスライス末端近くでは被写体領域に占める骨領域(頭頂骨や顎骨)の面積比が増大するため、骨領域(頭頂骨や顎骨)を透過させて脳の頭頂領域や脳幹領域を精度よく可視化することが困難となる場合がある。例えば、頭頂部では断層像の中心領域も骨領域であるため、同心楕円のパターンでは骨領域を透明にすることはできない。そこで、第2の実施の形態では、断層画像のスライス位置を更に考慮して補正倍率Sαを算出することで、特に頭部CTにおけるレンダリングの精度向上を図る。
第2の実施の形態では、第1の実施形態における図8のステップS23の処理(補正テーブルSαを作成する処理)が次のように変わる。すなわち、図8のステップS23において、制御部11(補正倍率算出部26c)は、断層画像Do(z)毎に、断層画像Do(z)に対応する幾何パラメータP(z)(楕円パラメータP(z))及び断層画像Do(z)のスライス位置の情報(以下、「スライスパラメータPs(z)」と表記)に基づいて、断層画像Do(z)の各画素(x、y)に対応する不透明度の補正倍率Sα(x、y、z)を算出し、算出した補正倍率Sα(x、y、z)が格納された3次元の補正テーブルSαを作成する。
例えば、z番目の断層画像Do(z)の各画素(x、y)の補正倍率Sα(x、y、z)は、以下のように算出される。
まず、図8のステップS23の(式11)と同様に、制御部11は、断層画像Do(z)に対応する幾何パラメータP(z)(楕円パラメータP(z))に基づいて、補正倍率Sα(x、y、z)を、次のように算出する。
(式17)
Sα(x、y、z)=1-k・r(x、y、z)・amp
但し、0≦Sα(x、y、z)≦1
r(x、y、z)={(x-Cx(z)-xoffset)/(a・W(z)/2)}+{(y-Cy(z)-yoffset)/(b・H(z)/2)}2s
更に第2の実施の形態では、制御部11は、次のように、スライスパラメータPs(z)に基づいて、Z方向の補正倍率Sαz(z)を算出し、Sα(x、y、z)に乗算する。ここで、スライスパラメータPs(z)とは、スライス総数Sz及びスライス順位zのことを言う。
(式18)
Sαz(z)=1.0-kz・{(z-Sz/2-zoffset)/Sz・c/2}}・ampz
Sα(x、y、z)=Sα(x、y、z)・Sαz(z)
但し、0≦Sα(x、y、z)、Sαz(z)≦1
(式17)(式18)において、kはXY平面、kzはZ方向の減衰係数である。初期値はk=1.0、kz=0.0であり、この場合、骨及び寝台が不可視となる。頭部CTの場合はkz>0.0に設定する。
a、bはそれぞれ横方向のサイズW(z)、縦方向のサイズH(z)の補正係数であり、cはZ方向幅Sz/2の補正係数(スライス総数倍率)である。補正しない場合はa=b=c=1.0に設定する。
xoffset、yoffsetはXY面、zoffsetはZ方向のオフセット値であり、初期値はxoffset=yoffset=zoffset=0(単位:ボクセル)である。
(式18)によって、断層画像Do(z)のスライス位置が中央から末端に位置するにつれ(スライス総数Sz/2とスライス順位zとの差が大きくなるにつれ)、補正倍率Sα(x、y、z)が一律に減衰し、楕円の中心であってもSα(x、y、z)は1未満の値をとる。これによって、スライス末端近くの骨領域(頭頂骨や顎骨)を効果的に透過させることができる。
図30は、補正倍率Sα(x、y、z)の分布パターンを示す図である。
図30の上図は、(式17)で算出される楕円パラメータP(z)のみに基づいた補正倍率Sα(x、y、z)の分布パターン(同心楕円分布)を示し、図30の下図は、(式18)で算出される楕円パラメータP(z)及びスライスパラメータPs(z)に基づいた補正倍率Sα(x、y、z)の分布パターン(同心楕円分布)を示す。図30に示すように、Z方向の補正倍率Sαz(z)が乗算されることによって、補正倍率Sα(x、y、z)がスライス位置に応じて一律に減衰する。すなわち、補正倍率Sαの分布パターンが、第1の実施の形態では2次元の同心楕円分布であったのに対し、第2の実施の形態では3次元の同心楕円体分布に拡張される。
尚、(式17)のampはXY面、(式18)のampzはZ方向の振幅倍率であり、初期値はamp=ampz=1.0である。図31に示すように、被写体中心部の補正倍率を平坦にし、飽和領域を広くしたい場合は、amp、ampz>1.0に設定する(但し、陰影感が弱くなる)。
制御部11は、全ての断層画像Do(z)について、(式17)(式18)によりZ方向の補正倍率Sαz(z)を乗算した補正倍率Sα(x、y、z)(0≦x≦Sx-1、0≦y≦Sy-1、0≦z≦Sz-1)を算出し、算出された補正倍率Sα(x、y、z)を格納する3次元の補正テーブルSαを作成する。
以上のように第2の実施の形態によれば、断層画像のスライス位置を考慮して補正倍率Sαを算出するので、更に柔軟な不透明度の制御を行うことができる。特に、頭部CT画像をレンダリングする際に有効である。頭部は球形状に近いため、胸部に比べ、被写体領域に占める骨領域の面積比がスライス方向(Z軸方向)に大きく変化し、特に骨領域の面積比が大きくなるスライス末端近くの頭頂骨や顎骨の領域を透過させて脳の頭頂領域や脳幹領域を精度よく可視化することが困難となる場合がある。第2の実施の形態によれば、断層画像のスライス位置が中央から末端に位置するにつれ、当該断層画像の全画素に対応する補正倍率Sαを一律に小さくするので、スライス末端近くで被写体領域に占める骨領域(頭頂骨や顎骨)が増大する頭部においても、骨領域を効果的に透過させることができる。
[第3の実施の形態]
第3の実施の形態は、GUI上で、補正テーブルSαを作成する際の各種パラメータをユーザが調整できるようにしたものである。
図32は、第3の実施の形態における医用画像処理装置1aの機能構成を示す図である。図32に示すように、第3の実施形態では、第1の実施の形態における医用画像処理装置1の機能構成(図3参照)に加え、パラメータ調整部30を更に備える。
パラメータ調整部30は、パラメータ調整画面40を表示し、補正倍率Sαを算出する際の各種パラメータをユーザに調整させる。具体的には、補正倍率Sαが(式11)に示したように楕円パラメータP(z)に基づいて算出される場合(第1の実施の形態の場合)には、楕円の横方向のサイズの補正係数a、縦方向のサイズの補正係数b、X方向のオフトセットxoffset、Y方向のオフセットyoffset、XY平面の減衰係数kを調整させる。
一方、補正倍率Sαが(式17)(式18)に示したように楕円パラメータP(z)及びスライスパラメータPs(z)に基づいて算出される場合(第2の実施の形態の場合)には、上記パラメータに加え、スライス総数倍率c、及びZ方向のオフセットzoffset、Z方向の減衰係数kzを更に調整可能とする。
図33は、第3の実施の形態における医用画像処理装置1aの全体動作を示すフローチャートである。ボリュームレンダリング像を生成するまでの処理(ステップS201~ステップS209)は、図6の第1の実施の形態のステップS1~ステップS9と同様である。第3の実施の形態では、ステップS209においてボリュームレンダリング像を生成し表示部16に表示した後、ユーザがパラメータの調整が必要と判断した場合(ステップS210;「Yes」)、パラメータ調整画面40においてパラメータをユーザに調整させる(ステップS211)。
図34は、パラメータ調整画面40の例を示す。図に示すように、パラメータ調整画面40において、XY平面の減衰係数k(図のXY方向減衰[%])、Z方向の減衰係数kz(図のZ方向減衰[%])、楕円の横方向のサイズの補正係数a(図の楕円倍率X)、楕円の縦方向のサイズの補正係数b(図の楕円倍率Y)、スライス総数倍率c(図の楕円倍率Z)、X方向のオフトセットxoffset(図の中心オフセットX)、Y方向のオフセットyoffset(図の中心オフセットY)、Z方向のオフセットzoffset(図の中心オフセットZ)を調整することができる。
パラメータの調整が終了すると(パラメータ調整画面40の「OK」ボタンを押下すると)、ステップS207に戻り、制御部11は、ステップS207~ステップS209の処理を再度実行する。すなわち、制御部11は、調整されたパラメータに基づいて、再度、補正テーブルSαを作成し直した上で(ステップS207)、ボクセル構造体V(ボクセルデータ)を作成し(ステップS208)、レンダリング処理を行う(ステップS209)。ステップS207~ステップS209の処理は、ステップS209において生成し表示されるレンダリング像を確認しながら、ユーザがパラメータの調整が十分と判断するまで(ステップS210;「No」)、繰り返し実行する。良好なレンダリング結果が得られ、ユーザがパラメータの調整が十分と判断すると(ステップS210;「No」)、制御部11は、必要に応じてレンダリング結果等のデータを出力し(ステップS212)、処理を終了する。ステップS212のデータの出力処理は、図6の第1の実施の形態のステップS10と同様である。
以上、第3の実施の形態によれば、補正テーブルSαを作成する際の各種パラメータをGUI上でユーザが調整できる。特に、被写体領域に撮影治具(寝台、固定具等)が映りこむと、被写体領域の外郭を近似する楕円の径が実際の被写体幅(頭幅や胸幅)より大きめに算出されたり、楕円の中心位置がズレたりし、楕円の近似精度が悪化する場合がある。第3の実施の形態によれば、GUI上で楕円の径や中心位置等を適宜調整することができるので、上記のような場合でも、良好なレンダリング結果を得ることができる。
尚、本実施の形態では、レンダリング結果を見ながらパラメータの調整が行えるよう、レンダリング後にパラメータを調整可能とする例を説明したが、パラメータはレンダリング後だけでなく、補正テーブルSαを作成する前の任意の段階で調整可能である。また、ステップS207における調整されたパラメータに基づく補正テーブルSαの再作成を行わずに、ステップS208におけるボクセル構造体V(ボクセルデータ)の再作成の段階で、(式17)(式18)に基づいて補正倍率をボクセルごとに都度算出する方法をとることもできる。
[第4の実施の形態]
第4の実施の形態は、更にMPR(Multi Plane
Reconstruction)像を生成し表示することを特徴とする。また、第3の実施形態におけるパラメータ調整を、MPR像を確認しながら行うことができる。
図35は、第4の実施の形態における医用画像処理装置1bの機能構成を示す図である。図35に示すように、第4の実施の形態では、第3の実施の形態における医用画像処理装置1aの機能構成(図32参照)に加え、MPR生成部31、MPR表示部32を更に備える。
MPR生成部31は、断層画像群Doに基づいて、任意断面の画像を生成し表示する。例えば、図36に示すように、XY平面に平行な体軸断面(アキシャル(Axial)断面)、XZ平面に平行な冠状断面(コロナル(Coronal)断面)、ZY平面に平行な矢状断面(サジタル(Sagittal)断面)のMPR像を生成する。また、XYZ空間を斜めに切断した斜断面(オブリーク(Oblique)断面)のMPR像を生成しても良い。
特に、MPR生成部31は、断層画像群Doの各画素(x、y、z)の信号値vに対して、各画素(x、y、z)に対応する補正倍率Sα(x、y、z)を乗算することで信号値vを補正し、補正された信号値vに基づいてMPR像を生成する。これにより、被写体形状(幾何形状)に応じて信号値vを減衰させたMPR像が得られる。
MPR表示部32は、MPR生成部31により生成されたMPR像を表示部16に表示する。また、信号値補正前と補正後の結果を比較するために、補正前の断層画像群Do(原画像)に基づいて生成したMPR像を併せて表示しても良い。
図37は、第4の実施の形態における医用画像処理装置1bの全体動作を示すフローチャートである。補正テーブルSαを作成するまでの処理(ステップS301~ステップS307)は、図33の第3の実施の形態のステップS201~ステップS207(図6の第1の実施の形態のステップS1~ステップS7)と同様である。
第4の実施の形態では、制御部11は、ステップS307において補正テーブルSαを作成した後、MPR像を生成し表示部16に表示する(ステップS308)。具体的には、以下のように、断層画像群Doの信号値v(=Do(x、y、z))を補正したうえで、MPR像を生成し表示する。
断層画像群Doの信号値が16ビットの場合((式1)参照)には、制御部11は、次の手順で断層画像群Doの信号値v(=Do(x、y、z))を補正する。
(式19)
(1)断層画像群Doの中で、Sz/2番目の中間スライスDo(z/2)における全ての画素の最小値Dmin、最大値Dmaxを算出する。
(2)下限値Lmin=(Dmax-Dmin)・γ+Dminを設定する。ここで、γは階調圧縮画像のコントラスト調整幅で、0に近いほどコントラストが増大するが輝度が小さくなる。通常はγ=0.1に設定する。
(3)次のように、断層画像群Doの信号値v(=Do(x、y、z))を補正する。
Do(x、y、z)≧Lminの場合:
Do(x、y、z)=(Do(x、y、z)-Lmin)・Sα(x、y、z)+Lmin
Do(x、y、z)≦Lminの場合:
Do(x、y、z)=Lmin
一方、断層画像群Doの信号値が8ビットの場合((式2)参照)には、制御部11は、次のように信号値v(Do(x、y、z))を補正する。
(式20)
Do(x、y、z)=Do(x、y、z)・Sα(x、y、z)
図38、図39は、頭部CT画像に基づいて生成されたMPR像の表示例を示す図である。図38、図39は、補正前の断層画像群Do(原画像)に基づいて生成されたMPR像と、補正後の断層画像Doに基づいて生成されたMPR像と、を体軸断面、冠状断面、矢状断面の各断面について比較表示した例である。図に示すように、補正倍率Sαにより信号値vを補正することによって、頭蓋領域の信号値が減衰していることが分かる。このように、補正倍率Sαの適用効果をMPR像の表示により確認することができる。
制御部11は、MPR像を生成し表示した後、ユーザがパラメータの調整が必要と判断した場合(ステップS309;「Yes」)、パラメータ調整画面40(図34参照)においてユーザにパラメータを調整させる(ステップS310)。
パラメータの調整が終了すると、ステップS307に戻り、制御部11は、ステップS307~ステップS308の処理を再度実行する。すなわち、制御部11は、調整されたパラメータに基づいて、再度、補正テーブルSαを作成し直した上で(ステップS308)、MPR像を生成し表示部16に表示する(ステップS308)。ステップS307~ステップS308の処理は、ステップS308において生成し表示されるMPR像を観察しながら、ユーザがパラメータの調整が十分と判断するまで(ステップS309;「No」)、繰り返し実行する。
良好なMPR像が得られ、ユーザがパラメータの調整が十分と判断すると(ステップS309;「No」)、制御部11は、ボクセル構造体V(ボクセルデータ)を作成し(ステップS311)、レンダリング処理を行い(ステップS312)、必要であればレンダリング結果等のデータを出力し(ステップS312)、処理を終了する。ステップS311~ステップS313の処理は、図6の第1の実施の形態のステップS8~ステップS10と同様である。
以上、第4の実施の形態によれば、断層画像群Doに基づいてMPR像を生成し表示する。特に、断層画像群Doの各画素の信号値vに対して、補正倍率Sαを乗算して信号値vを補正し、補正された信号値vに基づいてMPR像を生成する。これにより、3次元マスクを作成することなく、特定の領域を鮮明に可視化させたMPR像を得ることができる。従来の3次元マスクの代わりに、本開示では不透明度の補正倍率Sα(x、y、z)に関する補正テーブルの作成を必要とするが、これらは前述のように事前にステップ307で補正テーブルを作成する方法をとらず、ステップ308においてMPR像を生成する際に各画素ごとに(式17)(式18)で都度算出する方法をとることもでき、3次元マスクのようにメモリ上に確保する必要はない。
また、補正テーブルSαの適用効果を、ボリュームレンダリング像を生成する前に、MPR像にて確認することができ、また、MPR像の結果を見ながら補正テーブルSαのパラメータ調整を行うことができる。これにより、計算量が多いレンダリング処理の実行前に、補正テーブルSαのパラメータ調整を行うことができるので、所望のレンダリング像を得るまでの全体時間を短縮することができる。
尚、本実施の形態では、MPR像の結果を見ながらパラメータの調整が行えるよう、MPR表示後にパラメータを調整可能とする例を説明したが、パラメータはMPR表示後だけでなく、補正テーブルSαを作成する前の任意の段階で調整可能である。また、第3の実施の形態のように、レンダリング後にパラメータの調整が行えるようにしても良い。
以上、添付図面を参照しながら、本開示に係る医用画像処理装置等の好適な実施形態について説明したが、本開示はかかる例に限定されない。当業者であれば、本願で開示した技術的思想の範疇内において、各種の変更例又は修正例に想到し得ることは明らかであり、それらについても当然に本開示の技術的範囲に属するものと了解される。
1、1a、1b:医用画像処理装置
21:画像取得部
22:階調圧縮部
23:領域指定部
24:カラーマップ取得部
25:カラーマップ調整部
26:補正テーブル作成部
26a:幾何情報取得部
26b:幾何情報補正部
26c:補正倍率算出部
27:ボクセル作成部
28:レンダリング部
28a:第1レンダリング部
28b:第2レンダリング部
29:データ出力部
30:パラメータ調整部
31:MPR生成部
32:MPR表示部
40:パラメータ調整画面
Sα:補正倍率、補正テーブル

Claims (28)

  1. 複数の断層画像を取得する画像取得手段と、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、
    前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得手段と、を備え、
    前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正する
    医用画像処理装置。
  2. 前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報、及び当該断層画像のスライス位置の情報に基づいて算出する補正倍率算出手段と、を更に備え、
    前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算する
    請求項に記載の医用画像処理装置。
  3. 前記補正倍率算出手段は、算出した補正倍率を格納する補正テーブルを作成するテーブル作成手段と、を更に備え、
    前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正テーブルを参照して得られる補正倍率を乗算する
    請求項に記載の医用画像処理装置。
  4. 前記補正倍率算出手段は、前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、かつ、スライス位置が中央から末端に位置するほど小さくする
    請求項または請求項に記載の医用画像処理装置。
  5. 前記幾何形状は楕円であり、
    前記幾何情報取得手段は、前記幾何形状の情報として、楕円の中心座標、楕円の横方向のサイズ、及び楕円の縦方向のサイズを取得し、
    前記スライス位置の情報は、スライス総数、及びスライス順位であり、
    前記補正倍率算出手段は、前記中心座標、前記横方向のサイズ、前記縦方向のサイズ、前記スライス総数、及び前記スライス順位に基づいて、補正倍率を算出する
    請求項から請求項のいずれかに記載の医用画像処理装置。
  6. 前記補正倍率算出手段は、前記中心座標における前記補正倍率を最大とし、前記中心座標と前記横方向のサイズ及び前記縦方向のサイズの比が同一となる前記楕円内部の各同心楕円上の前記補正倍率を、各同心楕円の横方向のサイズ及び縦方向のサイズが大きくなるにつれ小さくし、かつ、前記スライス総数/2とスライス順位との差が大きくなるにつれ前記補正倍率を小さくする
    請求項に記載の医用画像処理装置。
  7. 前記補正倍率算出手段は、前記中心座標に所定のオフセットを加算して前記中心座標を補正し、前記横方向のサイズ及び前記縦方向のサイズに所定の横方向のサイズに対する倍率及び縦方向のサイズに対する倍率を乗算して前記横方向のサイズ及び前記縦方向のサイズを補正し、前記スライス総数に所定のスライス総数に対する倍率を乗算して補正し、前記スライス順位を所定のスライスオフセットを加算して補正し、補正された前記中心座標、前記横方向のサイズ、前記縦方向のサイズ、前記スライス総数、及び前記スライス順位に基づいて、前記補正倍率を算出する
    請求項または請求項に記載の医用画像処理装置。
  8. 前記オフセット、前記横方向のサイズに対する倍率、前記縦方向のサイズに対する倍率、前記スライス総数に対する倍率、及び前記スライスオフセットをユーザに調整させる調整手段と、を更に備える
    請求項に記載の医用画像処理装置。
  9. 複数の断層画像を取得する画像取得手段と、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段と、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段と、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段と、
    前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得手段と
    前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報に基づいて算出する補正倍率算出手段と、を備え、
    前記補正倍率算出手段は、前記被写体領域の中心部から離れた骨領域や体表を透過させ、中心部に近い臓器を鮮明に可視化するように前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、
    前記ボクセル作成手段は、
    前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報に基づいて補正し、かつ前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算する
    医用画像処理装置。
  10. 前記補正倍率算出手段は、算出した補正倍率を格納する補正テーブルを作成するテーブル作成手段と、を更に備え、
    前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正テーブルを参照して得られる補正倍率を乗算する
    請求項に記載の医用画像処理装置。
  11. 前記幾何形状は楕円であり、
    前記幾何情報取得手段は、前記幾何形状の情報として、楕円の中心座標、楕円の横方向のサイズ、及び楕円の縦方向のサイズを取得し、
    前記補正倍率算出手段は、前記中心座標、前記横方向のサイズ、及び前記縦方向のサイズに基づいて、前記補正倍率を算出する
    請求項9または請求項10に記載の医用画像処理装置。
  12. 前記補正倍率算出手段は、前記中心座標における前記補正倍率を最大とし、前記中心座標と前記横方向のサイズ及び前記縦方向のサイズの比が同一となる前記楕円内部の各同心楕円上の前記補正倍率を、各同心楕円の横方向のサイズ及び縦方向のサイズが大きくなるにつれ小さくする
    請求項11に記載の医用画像処理装置。
  13. 前記補正倍率算出手段は、前記中心座標に所定のオフセットを加算して前記中心座標を補正し、前記横方向のサイズ及び前記縦方向のサイズに所定の横方向のサイズに対する倍率及び縦方向のサイズに対する倍率を乗算して前記横方向のサイズ及び前記縦方向のサイズを補正し、補正された前記中心座標、前記横方向のサイズ及び前記縦方向のサイズに基づいて、前記補正倍率を算出する
    請求項11または請求項12に記載の医用画像処理装置。
  14. 前記オフセット、前記横方向のサイズに対する倍率、及び前記縦方向のサイズに対する倍率をユーザに調整させる調整手段と、を更に備える
    請求項13に記載の医用画像処理装置。
  15. 前記断層画像の幾何形状の近似精度を所定の基準で判断し、所定の基準に満たない前記断層画像の幾何形状の情報を、所定の基準を満たす他の前記断層画像の幾何形状の情報に基づいて補正する幾何情報補正手段と、を更に備える
    請求項1から請求項14のいずれかに記載の医用画像処理装置。
  16. 前記幾何形状は楕円であり、
    前記幾何情報取得手段は、前記幾何形状の情報として、楕円の中心座標、楕円の横方向のサイズ及び楕円の縦方向のサイズを取得し、
    前記幾何情報補正手段は、楕円の横方向のサイズが画像の横方向のサイズと一致する前記断層画像における楕円の横方向のサイズを、楕円の横方向のサイズが画像の横方向のサイズより小さい他の前記断層画像における楕円の横方向のサイズ及び縦方向のサイズの比と同一となるように、補正する
    請求項15記載の医用画像処理装置。
  17. 信号値と所定の閾値との差分値に応じて当該信号値に対応する不透明度を減衰させるように前記カラーマップを調整するカラーマップ調整手段と、を更に備え、
    前記ボクセル作成手段は、調整された前記カラーマップを参照して、ボクセルデータを作成する
    請求項1から請求項16のいずれかに記載の医用画像処理装置。
  18. 代表的な信号値と色値及び不透明度との対応関係に基づいて、所定範囲の信号値と色値及び不透明度との対応関係を定義するカラーマップを作成するカラーマップ作成手段と、を更に備え、
    前記カラーマップ取得手段は、作成した前記カラーマップを取得する
    請求項1から請求項17のいずれかに記載の医用画像処理装置。
  19. 前記レンダリング手段は、
    前記ボクセルデータを生成するボリュームレンダリング像に投影変換した座標系を視点座標系とすると、視点座標系において、前記ボリュームレンダリング像の各画素よりZ軸方向に沿って、Z軸の上限値より下限値に向けて、視点座標系のボクセル座標毎に座標変換を行って前記ボクセルデータより不透明度を取得しながら、不透明ボクセルを探索し、最初に見つかった不透明ボクセルの視点座標系におけるZ座標を、前記ボリュームレンダリング像の画素毎に記録した探索制御マスクを作成する探索制御マスク作成手段と、
    前記ボリュームレンダリング像の画素毎に、前記探索制御マスクからZ座標を取得し、取得したZ座標よりZ軸の下限値に向けてZ軸方向に沿って、所定の光強度をもつ仮想光線を照射する際、視点座標系のボクセル座標毎に座標変換を行って前記ボクセルデータより不透明度を取得し、不透明ボクセルが見つかった場合、当該ボクセル座標に対して座標変換を行って前記ボクセルデータより色値を取得し、当該ボクセルの不透明度に基づいて前記光強度を減衰させるとともに、当該ボクセルの不透明度及び色値並びに前記減衰させた光強度に基づいて累積輝度値を算出する処理を繰り返し、算出された累積輝度値に基づいて、前記ボリュームレンダリング像の当該画素に対応する画素値として与えるレイキャスティング手段と、
    を備える請求項1から請求項18のいずれかに記載の医用画像処理装置。
  20. 前記探索制御マスク作成手段及び前記レイキャスティング手段は、前記座標変換を行って前記ボクセルデータより不透明度または色値を取得する際、
    所定の回転を定義した回転行列、XYZ軸各方向のオフセット値、XYZ軸方向の拡大又は縮小倍率、Z軸方向の変倍率、注視点から視点までの距離を含む前記所定の座標変換のパラメータを取得し、
    前記視点座標系のボクセルの整数値の座標を、前記パラメータに基づいて前記ボクセルデータの座標系に変換を行って、前記ボクセルデータの実数値の座標を算出し、
    算出した実数値の座標の近傍の複数の整数値の座標に対応する前記ボクセルデータの複数のボクセルを特定し、
    特定した複数のボクセルの不透明度または色値に基づいて前記ボクセルデータより取得される不透明度または色値として算出する
    請求項19に記載の医用画像処理装置。
  21. 前記レンダリング手段は、
    前記ボクセルデータに基づいて3Dテクスチャを生成する3Dテクスチャ生成手段と、
    前記3Dテクスチャに対して所定の座標変換を行って変換後3Dテクスチャを生成する座標変換手段と、
    3次元空間のXY座標面上の四角形をZ軸方向に並べた積層四角形を設定する積層四角形設定手段と、
    所定の視点からZ軸方向に平行な視線上の前記四角形のXY座標に対応する前記変換後3Dテクスチャのボクセルの色値を前記ボクセルの不透明度に基づいて前記視点から遠い四角形の順にアルファブレンディングして取得し、前記ボリュームレンダリング像の画素値として与える画素値算出手段と、
    を備える請求項1から請求項18のいずれかに記載の医用画像処理装置。
  22. 前記座標変換手段は、
    所定の回転を定義した回転行列、視野角度、視点位置、クリッピング位置、XYZ軸各方向のオフセット値、XYZ軸方向の拡大又は縮小倍率、Z軸方向の変倍率を含む所定の座標変換のパラメータを取得し、
    前記3Dテクスチャに対して、前記取得したパタメータを用いた前記所定の座標変換を行って前記変換後3Dテクスチャを生成する
    請求項21に記載の医用画像処理装置。
  23. 前記座標変換手段及び前記画素値算出手段は、ビデオカードに搭載されたGPU及びフレームメモリを用いて実行する
    請求項22に記載の医用画像処理装置。
  24. コンピュータの制御部が、
    複数の断層画像を取得する画像取得ステップと、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得ステップと、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成ステップと、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリングステップと、
    前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得ステップと、を備え、
    前記ボクセル作成ステップは、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正する
    医用画像処理方法。
  25. コンピュータの制御部が、
    複数の断層画像を取得する画像取得ステップと、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得ステップと、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成ステップと、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリングステップと、
    前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得ステップと、
    前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報に基づいて算出する補正倍率算出ステップと、を備え、
    前記補正倍率算出ステップは、前記被写体領域の中心部から離れた骨領域や体表を透過させ、中心部に近い臓器を鮮明に可視化するように前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、
    前記ボクセル作成ステップは、
    前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報に基づいて補正し、かつ前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出ステップにより算出された補正倍率を乗算する
    医用画像処理方法。
  26. コンピュータを、
    複数の断層画像を取得する画像取得手段、
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段、
    前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得手段、として機能させ、
    前記ボクセル作成手段は、前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正するプログラム。
  27. コンピュータを、
    複数の断層画像を取得する画像取得手段と
    信号値と色値及び不透明度との対応関係を定義するカラーマップを取得するカラーマップ取得手段、
    前記カラーマップを参照することで、前記断層画像の各画素の信号値を、信号値に応じた色値及び不透明度に変換し、色値及び不透明度を保持するボクセルの集合であるボクセルデータを作成するボクセル作成手段、
    前記ボクセルデータに基づいてボリュームレンダリング像を生成するレンダリング手段、
    前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得手段、
    前記断層画像毎に、断層画像の各画素に対応する不透明度の補正倍率を、当該断層画像に対応する幾何形状の情報に基づいて算出する補正倍率算出手段、として機能させ、
    前記補正倍率算出手段は、前記被写体領域の中心部から離れた骨領域や体表を透過させ、中心部に近い臓器を鮮明に可視化するように前記補正倍率を、前記幾何形状の幾何中心から遠ざかるにつれ小さくし、
    前記ボクセル作成手段は、
    前記断層画像の各画素の信号値に応じた前記不透明度を、当該断層画像に対応する幾何形状の情報に基づいて補正し、かつ前記断層画像の各画素の信号値に応じた前記不透明度に対して、前記補正倍率算出手段により算出された補正倍率を乗算するプログラム。
  28. コンピュータの制御部が、
    複数の断層画像を取得する画像取得ステップと、
    前記断層画像毎に、断層画像の被写体領域の外郭を所定の幾何形状で近似し、当該幾何形状の情報を取得する幾何情報取得ステップと、
    前記断層画像の各画素の信号値に応じたオパシティカーブにより規定される不透明度を、当該断層画像に対応する幾何形状の情報及び当該断層画像のスライス位置の情報に基づいて補正する不透明度補正ステップと、
    前記断層画像の各画素に対応する補正後の不透明度を格納したデータを作成するデータ作成ステップと、
    を含むデータ作成方法。
JP2018103767A 2018-05-30 2018-05-30 医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法 Active JP7180123B2 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2018103767A JP7180123B2 (ja) 2018-05-30 2018-05-30 医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2018103767A JP7180123B2 (ja) 2018-05-30 2018-05-30 医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法

Publications (2)

Publication Number Publication Date
JP2019205791A JP2019205791A (ja) 2019-12-05
JP7180123B2 true JP7180123B2 (ja) 2022-11-30

Family

ID=68768208

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2018103767A Active JP7180123B2 (ja) 2018-05-30 2018-05-30 医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法

Country Status (1)

Country Link
JP (1) JP7180123B2 (ja)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP7496247B2 (ja) 2020-06-15 2024-06-06 株式会社 ディー・エヌ・エー 画像描画処理装置、画像描画処理方法及び画像描画処理プログラム並びに電子ゲーム提供装置
CN112508858B (zh) * 2020-11-17 2024-04-23 杭州依图医疗技术有限公司 一种医学影像处理方法及装置、计算机设备

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2006136619A (ja) 2004-11-15 2006-06-01 Ziosoft Inc 画像処理方法
JP2011019633A (ja) 2009-07-14 2011-02-03 Toshiba Corp X線診断装置及び被曝線量低減用制御プログラム
WO2012073769A1 (ja) 2010-11-29 2012-06-07 株式会社 日立メディコ 画像処理装置及び画像処理方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2006136619A (ja) 2004-11-15 2006-06-01 Ziosoft Inc 画像処理方法
JP2011019633A (ja) 2009-07-14 2011-02-03 Toshiba Corp X線診断装置及び被曝線量低減用制御プログラム
WO2012073769A1 (ja) 2010-11-29 2012-06-07 株式会社 日立メディコ 画像処理装置及び画像処理方法

Also Published As

Publication number Publication date
JP2019205791A (ja) 2019-12-05

Similar Documents

Publication Publication Date Title
US7529396B2 (en) Method, computer program product, and apparatus for designating region of interest
JP2020512129A (ja) オブジェクトグリッド増強を用いて高次元画像データから低次元画像データを合成するためのシステムおよび方法
US7424140B2 (en) Method, computer program product, and apparatus for performing rendering
JP2007537770A (ja) 内視鏡画像内の管腔状構造のディスプレイ最適化のための動的なクロップボックス決定方法
JP5562339B2 (ja) 医用画像表示装置及び医用画像表示方法
JP3836097B2 (ja) 医用画像生成装置および方法、ならびに、プログラム
JP7180123B2 (ja) 医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法
JP7275669B2 (ja) マスク生成装置、3次元再構成像生成装置、マスク生成方法、及びプログラム
JP4122314B2 (ja) 投影画像処理方法、投影画像処理プログラム、投影画像処理装置
JP7019594B2 (ja) 撮像システムおよび方法
JP7155670B2 (ja) 医用画像処理装置、医用画像処理方法、プログラム、及びデータ作成方法
JP7003635B2 (ja) コンピュータプログラム、画像処理装置及び画像処理方法
JP7013849B2 (ja) コンピュータプログラム、画像処理装置及び画像処理方法
JP7131080B2 (ja) ボリュームレンダリング装置
JP7095409B2 (ja) 医用画像処理装置、医用画像処理方法、プログラム、及びmpr像生成方法
JP6418344B1 (ja) コンピュータプログラム、画像処理装置及び画像処理方法
CN111210898A (zh) 一种对dicom数据进行处理的方法和装置
JP6436258B1 (ja) コンピュータプログラム、画像処理装置及び画像処理方法
JP7223312B2 (ja) ボリュームレンダリング装置
JP2006223894A (ja) 医用画像生成装置および方法、ならびに、プログラム
JP7283603B2 (ja) コンピュータプログラム、画像処理装置及び画像処理方法
JPH0728976A (ja) 画像表示装置
JP6443574B1 (ja) レイキャスティングプログラム、探索制御データ、探索制御データ生成方法、レイキャスティング装置
JP7206846B2 (ja) ボリュームレンダリング装置
Harris Display of multidimensional biomedical image information

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20210323

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20220222

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20220301

A601 Written request for extension of time

Free format text: JAPANESE INTERMEDIATE CODE: A601

Effective date: 20220428

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20220622

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20221018

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20221031

R150 Certificate of patent or registration of utility model

Ref document number: 7180123

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150