JP6014738B2 - 三次元画像の投影方法 - Google Patents

Description

本発明は、三次元画像の投影方法に関し、より詳しくは、放射線医学影像に用いられる三次元画像の投影方法に関する。

科学技術の進歩に伴って、放射線医学影像も近代発展した有効な診断工具である。その結像原理は物理現象に対して適切な数値モデルを構築し、像再構成アルゴリズムにより撮影時に受信する信号を三次元（three-dimensional, 3D）影像に還元する。ＣＴスキャン(computed tomography、CT)、陽電子放射形コンピューター断層撮影法(positron emission tomography、PET)、単光子放射型ＣＴスキャン法(single-photon emission computed tomography、SPECT)等の高等医療器材にこの技術が使用されている。

像再構成アルゴリズムの歴史は、その初期にはフィルタ補正逆投影法（filtered back projection、FBP）、マトリックス反転トモシンセシス（matrix inversion tomosynthesis、MITS）等の分析式解法が採用され、近年ではコンピューター技術の進歩及びハードウェアの進歩に伴い、高影像品質の逐次技術、例えば、最尤度期待値最大化法（maximum likelihood expectation maximization、MLEM）、サブセット化による期待値最大化法（ordered subset expectation maximization、OSEM）、代数的再構成法（algebraic reconstruction technique、ART）、凸射影法（projection onto convex sets、POCS）等の再構成法が流行した。しかしながら、好ましい影像品質を得るためには、さらに精確なシステムマトリックスモデルにより逐次再構成を行う必要がある。

なお、ＣＴスキャン(computed tomography, CT)を例にすると、撮影器が信号を受信して離散配列に変換すると共にコンピューターに保存し、基本的な影像再構成数学モデルは下述の離散化のマトリックスの表示に用いられる：

仮に、ベクターＣ(ｍ×１)を受信した信号として離散配列方式でパソコンに保存され、Ｇがシステムマトリックス(ｍ×ｎ)とすると、要求される三次元（three-dimensional、3D）影像はｆ(ｎ×１)となる。光子は空間中で直線移動するため、乱射等の他の物理要素を考慮せずともよい。仮に、撮影器が信号を受信させ、１つの体積を有する光源が発射され、何れか１つの検出器により受光されると、Ｇは簡略化されて幾何学的なマトリックスとして描写され、ｇｉｊは第iのボクセル（voxel）の第ｊの検出器に対する幾何学検出効果を表す。異なる位置では各検出器が受信する確率に差異があり、確率データをフィードバックさせて再構成過程で、空間分布の差異の修正を行う必要がある。

上述の方式では、大きなマトリックスの演算には時間が掛かるというボトルネックがあり、仮に、ｆが１つの５１２×５１２×５１２の三次元影像だとすると、ｎ＝１３４２１７７２８個の影像画素(pixel)を有し、３０７２×３０７２を有する１つの検出器により３６０度の撮影を行い、１度隔てる毎に１回スキャンを行い、ｍ＝３０７２×３０７２×３６０となり、システムマトリックスの次元数は３３９７３８６２４０×１３４２１７７２８〜４．５６ｅ＋１７となり、コンピューターの演算ソースには制限があるため、如何に精確にシステムを描写し、疎行列（sparse matrix）の概念によりコンピューターの演算時間を短縮するかが、逐次再構成の応用の課題であった。

現在の方法では、従来の技術としてレイトレーシング(ray tracing)が常用され、具体的なプログラム最適化の実施方法としてSiddon R L(１９８５)が文献にて提出しており、多くの影像処理技術は共にレイトレーシングを青写真として開発された。例えば、光子の空間中の運行軌跡のシミュレーションや水晶体透過効率の検査等である。ＣＴスキャン(computed tomography, CT)を例にすると、仮想光源及び検出器が無体積であり、光源の中心及び検出器の面積の中心を代表とし、点対点のビーム接続により撮影状況を描写し、ビーム接続が空間中で影像中のボクセル（voxel）を通過する交差の長さを確率とし、１本のビーム接続がボクセル（voxel）を通過する長さから影像空間全体の長さを引いて影像再構成に使用するｇｉｊとする。この方法は撮影時の光源の体積及び検出器の体積の影響を無視している。また、高解析度を追求すると、ボクセル（voxel）が単一の検出器より小さい場合、ボクセル（voxel）が無視される現象が発生する。上述の問題を解決するため、ビームの数量を増加させる方式により、検出器の数学モデルをより細く切断し、これにより検出器の面積が数値演算に与える影響を低減させる。例えば、Huesman R H(２０００)の文献を例にすると、各検出ユニットが３×３個に分割され、２組の検出ユニットが本来の９倍に各々切断されるため、光子を受けると１対のペアになり、ビームの総数量が本来の７２９倍に変化し、ビームの傾度及び比例により圧縮及び簡略化を行うが、Ｏ(ｎ４)の対価に対しては、効果が限定された。総括すると、レイトレーシング(ray tracing)方式により、モデルが簡易になり演算しやすくなるが、ただし、日々高まる撮影結果の精確性の要求に応えるには、前記技術では複数本のサブ検出ビーム(sub-ray)に切り分けて再構成処理を行わねばならず、演算時間が大幅に増加し、実用できなかった。

さらに、ボクセル（voxel）を１つの球体に結像する方法があり、この技術は球対称体素（spherically symmetric volume elements(blobs)）と呼ばれ、適切な数学公式と組み合わされ、点が球体の中心に至る距離により、球体内部の密度の調整を行う。この技術は高解析度が必要な影像処理に応用され、例えば、Marabini R L(１９９８)はこの方法を電子顕微鏡に応用した。ただし、前述の方法は計算に非常に時間が掛かり、学術研究分野或いはシミュレーションで常用されるが、商業用途には広がらなかった。

商業分野については、前述のレイトレーシング(ray tracing)方法に内挿技術が組み合わされて改良された線駆動（ray driven）及び画素駆動（pixel driven）技術が常用される。前述の方式は内挿及び距離の関係により重量の分配が行われ、空間中のいかなるボクセル（voxel）のｚ方向データが圧縮されることで少量のデータが損失するのが無視されず、撮影角度が大きい場合、誤差が大きくなった。

再構成の精確性に対する要求に応えるには、演算コストが高くなって逐次再構成に時間が掛かるため、市場では大部分の商品は分析式解法が主流である。このため、新世代のモデルベース（model-based）逐次再構成の発展が将来的に主流になる。例えば、対応する足跡（footprint）グループ化演算法及び平行加速技術が設計され、医療メーカーも逐次再構成法の技術開発を積極的に推し進めており、精確なモデルに逐次再構成技術を組み合わせて、高解析度の３Ｄ影像を獲得し、演算のハードルを突破しようとしている。

本発明者は上記の欠点が改善可能と考え、鋭意検討を重ねた結果、合理的設計で上記の課題を効果的に改善する本発明の提案に到った。そこで、本発明は、以上の従来技術の課題を解決する為になされたものである。即ち、本発明の目的は、計算に費やす時間を短縮させる三次元画像の投影方法を提供する。

上述した課題を解決し、上記目的を達成するための本発明は、本発明に係る三次元画像の投影方法は、光源により放射線場が発生して検出モジュールに投射され、検出モジュールはアレイ状に配列される複数の検出器を有する、光源を提供する工程と、光源と前記検出モジュールとの間には複数の空間ユニットが堆積されて形成される三次元空間配列のボクセルを有する工程と、放射線場の検出器と光源との間のサブ放射線場を各々検出して定義する工程と、各検出器がサブ放射線場で対応する各空間ユニットに基づいて演算し、サブ幾何学的因子を獲得する工程と、各検出器が前記サブ放射線場で対応する各空間ユニットが有するサブ幾何学的因子が組み合わされて幾何学的因子マトリックスが形成される工程とを含む。

また、上述の演算は、検出器、各検出器に対応するサブ放射線場及び複数の空間ユニットが第一検出平面及び第二検出平面にそれぞれ投影され、検出器、各検出器に対応するサブ放射線場及び複数の空間ユニットは第一検出平面及び第二検出平面に対応する平面検出器、平面放射線場及び隣接されて配列されて二次元平面となる複数の一次元配列をそれぞれ有し、各一次元配列は隣接して配列される複数のグリッドを有し、第一検出平面は第一軸及び第二軸で構成され、第二検出平面は前記第一軸及び第三軸で構成される工程と、
平面検出器が第一検出平面の第一軸及び第二検出平面の第一軸までそれぞれ回転し、平面検出器及び複数のグリッドが対応する位置関数変数を発生させる工程と、位置関数変数に基づいて第一検出平面及び第二検出平面上での各グリッドに関して各平面検出器が対応する平面放射線場で占める幾何学的比率値をそれぞれ計算する工程と、第一検出平面及び第二検出平面上の相互に対応するグリッドが有する幾何学的比率値を相乗して各検出器が対応する空間ユニットが有するサブ幾何学的因子を獲得する工程と、各検出器がサブ放射線場で対応する各空間ユニットが有するサブ幾何学的因子が組み合わされて幾何学的因子マトリックスが形成される工程とを含む。

本発明の好適な実施の形態において、上述の位置関数変数に基づいて第一検出平面及び第二検出平面上での各グリッドに関して各平面検出器が対応する平面放射線場で占める幾何学的比率値をそれぞれ計算する工程では、以下の工程をさらに含む。第一検出平面及び第二検出平面で各平面検出器が対応する平面放射線場の経路の計算上の複数の放射線場交点をそれぞれ調べ、経路の計算は第一検出平面及び第二検出平面で第一検出平面の前記第二軸及び第二検出平面の第三軸にそれぞれ平行する工程。位置関数変数に基づいて各グリッドが第一検出平面及び第二検出平面で経路の計算上の対応する複数の特定の交点をそれぞれ調べる工程。計算を行い、複数の放射線場交点の位置関数及び複数の特定の交点の位置関数を比較し、経路の計算上で各グリッドが各平面検出器が対応する平面放射線場で占める幾何学的比率値を決定する工程。

本発明の好適な実施の形態において、上述の第一検出平面及び第二検出平面で各平面検出器が対応する平面放射線場の経路の計算上の複数の放射線場交点をそれぞれ調べる工程の後には、以下の工程をさらに含む。位置関数変数に基づいて、経路の計算が平面放射線場で占有する長さを比例により算出する工程。経路の計算が平面放射線場で占有する長さに基づいて、複数の放射線場交点の位置関数を算出する工程。

本発明の好適な実施の形態において、上述の経路の計算が平面放射線場で占有する長さを比例により算出する工程は、以下の工程をさらに含む。各平面検出器及び対応する平面放射線場により包囲されて形成される第一領域を調べる工程。平面放射線場の経路の計算により包囲されて形成される第二領域を調べる工程。位置関数変数に基づいて、第一領域及び第二領域の比例を比較して経路の計算が平面放射線場で占有する長さを獲得する工程。

本発明の好適な実施の形態において、上述の計算を行い、複数の放射線場交点の位置関数及び複数の特定の交点の位置関数を比較し、経路の計算上の各グリッドが各平面検出器が対応する平面放射線場で占める幾何学的比率値を決定する工程の前には、以下の工程をさらに含む。境界距離が画定される工程。経路の計算で選択された距離を選択し、選択された距離は境界距離に等しいかそれより大きい工程。経路の計算ではグリッドを中心として選択された距離は半径の範囲内になり、平面放射場内に前記特定の交点があるか否かを判断する工程。

本発明の好適な実施の形態において、上述の前記平面放射場内に前記特定の交点を有する工程では、以下の工程をさらに含む。サブ幾何学的因子に加重係数を乗算する工程。

本発明の好適な実施の形態において、上述の第一検出平面及び第二検出平面上で相互に対応するグリッドが有する幾何学的比率値を相乗して各検出器が対応する空間ユニットが有するサブ幾何学的因子を獲得する工程の後には、以下の工程をさらに含む。第一検出平面及び第二検出平面の第一軸方向に平行する複数のグリッドの数をそれぞれ計算する工程。サブ幾何学的因子に均一化パラメータを乗算する工程。ここでは、均一化パラメータは１と前記複数のグリッドの数との比率である。

本発明の好適な実施の形態において、上述の光源は回転運動をさらに行う。

本発明の好適な実施の形態において、上述の検出モジュールはＸ線断層撮影合成法（Tomosynthesis）或いはＣＴスキャン（CT）の検出モジュールである。

上述のように、本発明の三次元画像の投影方法は、まず本来三次元空間の数学モデルが２つの二次元の検出平面に投影され、平面検出器が検出平面の軸まで回転し、煩雑な演算を簡略化させ、計算に掛かる時間を短縮させる。

本発明に係る三次元画像の投影方法を示すフローチャート図である。 本発明の光源及び検出器の関係を示す概念図である。 図１に示すサブ幾何学的因子を計算するフローチャート図である。 図２に示すサブ放射線場及び空間ユニットが第一検出平面に投影させる概念図である。 図２に示すサブ放射線場及び空間ユニットが第二検出平面に投影させる概念図である。 図４Ａに示す平面検出器が第一軸まで回転させる概念図である。 本発明に係る幾何学的比例値を計算するフローチャート図である。 本発明の計算平面放射線場の経路の計算での長さ及び放射線場交点にある位置関数のフローチャート図である。 本発明の経路の計算が平面放射線場で占有する長さを比例により算出するフローチャート図である。 本発明の位置関数変数に基づいてグリッドが第一検出平面で経路の計算上で対応する特定の交点を調べるフローチャート図である。 図５に示す経路の計算で、交点線段を調べて、計算する概念図である。 ボクセルを示す概念図である。 ボクセルが任意の回転なく放射線場と交差を行う概念図である。 ボクセルが回転後に放射線場と交差を行う概念図である。 本発明に係る三次元画像の投影方法の均一化パラメータ工程を示す概念図である。

本発明における好適な実施の形態について、添付図面を参照して説明する。尚、以下に説明する実施の形態は、特許請求の範囲に記載された本発明の内容を限定するものではない。また、以下に説明される構成の全てが、本発明の必須要件であるとは限らない。

まず、本発明の三次元画像の投影方法は影像再構成（image reconstruction）に用いられ、Ｘ線断層撮影合成法（Tomosynthesis）或いはＣＴスキャン（computed tomography、CT）等の逐次再構成に応用される。

図１は本発明に係る三次元画像の投影方法を示すフローチャート図であり、図２は本発明の光源及び検出器の関係を示す概念図である。三次元画像の投影方法Ｓ１００は、以下の工程を含む。まず、工程Ｓ１１０では光源２０を提供する。光源２０により放射線場２１０が発生して検出モジュール４０に投射され、検出モジュール４０はアレイ状に配列される複数の検出器ｊを有する。

本実施形態では、光源２０は、例えばＸ線を発生させる光源であり、検出モジュール４０は光源２０が発生させるＸ線を検出するために使用される。実際上、光源の種類に制限はなく、実際の使用時の要求に合わせて決定される。このほか、実施形態において、光源２０はさらに回転運動を行い、測定対象物及び検出モジュール４０の周りを回転し、ただし本発明はこれに制限されず、実際の使用の要求に基づいて決定される。また、本実施形態に係る検出モジュール４０はＸ線断層撮影合成法(Tomosysnesis)やＣＴスキャン(CT)の検出モジュールに使用される検出モジュールである。本発明は上述の検出モジュールの種類に制限はなく、実際の使用の要求に基づいて決定される。

次は、工程Ｓ１２０を実行し、光源２０と検出モジュール４０との間に複数の空間ユニットｉが堆積されて形成される三次元空間配列のボクセル３０(Voxel)を有する。具体的には、図２に示すように、本実施形態に係る光源２０及び上述のアレイ状に配列される複数の検出器ｊにより第一軸Ｚ、第二軸Ｘ及び第三軸Ｙで構成される三次元空間が定義され、且つ三次元空間内には複数の空間ユニットｉが堆積されて形成される三次元空間配列のボクセル３０ (Voxel)を有する。また、説明が必要な点は、前述の光源２０が具設される位置に特段の制限はなく、本実施形態では、前記光源は三次元空間の原点に具設される。また、これを基準の中心点とすることで、計算の煩雑さを低減させる。続いて、工程Ｓ１３０を実行し、放射線場２１０の検出器ｊと光源２０との間のサブ放射線場２１２を各々検出して定義し、これは前記放射線場２１０の一部に属する。

続いて、工程Ｓ１４０を実行し、各検出器ｊがサブ放射線場２１２で対応する各空間ユニットｉに基づいて演算を行い、サブ幾何学的因子ｇｉｊを獲得する。

図３は図１に示すサブ幾何学的因子を計算するフローチャート図である。上述の演算は工程Ｓ１４２乃至工程Ｓ１４８を含む。工程Ｓ１４２を実行し、検出器ｊ、各検出器ｊに対応するサブ放射線場２１２及び複数の空間ユニットｉが第一検出平面Ｐ１及び第二検出平面Ｐ２にそれぞれ投影され、第一検出平面Ｐ１は第一軸Ｚ及び第二軸Ｘで構成され、第二検出平面Ｐ２は第一軸Ｚ及び第三軸Ｙで構成される（図２及び図３参照）。また、ちなみに、上述の検出平面は前述の第一軸、第二軸及び第三軸の間の何れか二軸が選択されるが、特別な制限はなく、実際の計算時の需要により決定される。

図４Ａは図２に示すサブ放射線場及び空間ユニットが第一検出平面に投影される概念図であり、図４Ｂは図２に示すサブ放射線場及び空間ユニットが第二検出平面に投影される概念図である。上述の工程Ｓ１４２に示すように、第一検出平面Ｐ１(即ち、ＸＺ平面)に投影された後、図４Ａに示すように、検出器ｊ、各検出器ｊに対応するサブ放射線場２１２及び複数の空間ユニットｉは第一検出平面Ｐ１に対応する平面検出器４０ｊ、平面放射線場２１２ａ及び隣接して配列される二次元平面の複数の一次元配列３２を有し、各一次元配列３２は隣接して配列される複数のグリッド３２０を有する。

同様に、第二検出平面Ｐ２(ＹＺ平面)に投影された後、図４Ｂに示すように、検出器ｊ、各検出器ｊに対応するサブ放射線場２１２及び複数の空間ユニットｉは第二検出平面Ｐ２に対応する平面検出器４０ｊ、平面放射線場２１２ｂ及び隣接して配列される二次元平面の複数の一次元配列３４を有し、各一次元配列３４は隣接して配列される複数のグリッド３４０を有する。

上述のように、本実施形態は上述の投影モデル(即ち検出器ｊ、各検出器ｊに対応するサブ放射線場及び複数の空間ユニットｉの投影)が第一検出平面Ｐ１(即ちＸＺ平面)及び第二検出平面Ｐ２ (即ちＹＺ平面)に投影される。対照については、図２の各空間ユニットｉは図４Ａの第一検出平面Ｐ１のグリッド３２０に対応する。同様に、図２の各空間ユニットｉは図４Ｂの第二検出平面Ｐ２のグリッド３４０に対応する。

次に、工程Ｓ１４４を実行し、平面検出器４０ｊが前記第一検出平面Ｐ１の第一軸Ｚ及び第二検出平面Ｐ２の第一軸Ｚまでそれぞれ回転する（図３参照）。換言すれば、第一検出平面Ｐ１で、平面検出器４０ｊが第一軸Ｚ(即ちＺ軸)まで回転し、同様に、第二検出平面Ｐ２で、平面検出器４０ｊが第一軸Ｚ(即ちＺ軸)まで回転する。

上述の工程では、本発明の三次元画像の投影方法Ｓ１００は、まず、本来の三次元空間の数学モデルが２つの二次元の検出平面に投影され、平面検出器４０ｊが対応する第一軸Ｚまで回転し、煩雑な演算を簡略化させ、計算に掛かる時間を短縮させる。他の実施形態では、第一検出平面Ｐ１で、平面検出器４０ｊがＸ軸(即ち第二軸Ｘ)まで回転し、第二検出平面Ｐ２で、平面検出器４０ｊがＹ軸(即ち第三軸Ｙ)まで回転する。同様に、上述の演算の簡略化の目的及び効果を達成し、本発明はこれに対して制限はなく、実際の演算時の要求によって決定される。

以下、第一検出平面Ｐ１で、平面検出器４０ｊが第一軸Ｚ(即ちＺ軸)まで回転するのを例に説明する。
図５は図４Ａに示す平面検出器が第一軸に回転される概念図である（図４Ａ及び図５参照）。ちなみに、説明しやすくするため、図５は図４Ａの第一検出平面Ｐ１の内の何れか１つのグリッド３２０を図示し、光源２０の位置は原点(０，０)である。

本実施形態では、第一検出平面Ｐ１で、平面検出器４０ｊが第一軸Ｚ(即ちＺ軸)まで回転し、平面検出器４０ｊ及び複数のグリッド３２０が対応する位置関数変数を発生させる。例えば、平面検出器４０ｊの表面の中心点ＰＣ１の位置関数は(ｕ，ｗ)であり、グリッド３２０の表面の中心点ＰＣ２の位置関数は(ｘ'，ｚ')であり、Ｚ'軸及び第一軸Ｚは夾角θを有する（図４Ａ参照）。

図５に示すように、平面検出器４０ｊが第一軸Ｚ(即ちＺ軸)まで回転した後、平面検出器４０ｊの表面の中心点ＰＣ１の位置関数(ｕ，ｗ)が位置関数変数(０，ｈ)に変換され、グリッド３２０の表面の中心点ＰＣ２の位置関数(ｘ'，ｚ')が位置関数変数(ｘ＿ｒｏｔ，ｚ＿ｒｏｔ)に変換される。距離ｈは光源２０から平面検出器４０ｊの表面の中心点までの距離である。

具体的には、上述の平面検出器４０ｊが第一軸Ｚ(即ちＺ軸)まで回転し、これはｓｉｎ或いはｃｏｓの三角関数を乗算する方式により第一軸Ｚまで回転する。このさい、グリッド３２０の表面の中心点ＰＣ２の位置関数(ｘ'，ｚ')が位置関数変数(ｘ＿ｒｏｔ，ｚ＿ｒｏｔ)に変換され、関連する数式は下述する：

当然ながら、平面検出器４０ｊが第二軸Ｘ或いは第三軸Ｙまで回転すると、同様の原理から類似する数式(１)(２)が導き出される。

続いて、工程Ｓ１４６を実行し（図３参照）、位置関数変数に基づいて第一検出平面Ｐ１及び第二検出平面Ｐ２での各グリッドに関して各平面検出器が対応する平面放射線場が占める幾何学的比率値がそれぞれ計算される。本工程では、工程Ｓ１４６が同時に実行され、即ち、第一検出平面Ｐ１及び第二検出平面Ｐ２が同時に処理されるが、ただし、第一検出平面Ｐ１が先に処理された後、第二検出平面Ｐ２が処理されてもよく、或いは、第二検出平面Ｐ２が先に処理された後、第一検出平面Ｐ１が処理されてもよい。実際の使用の要求に基づいて決定されるため、本発明では制限されない。

図６は本発明に係る幾何学的比例値を計算するフローチャート図である。具体的には、工程Ｓ１４６２を実行し、第一検出平面及び第二検出平面で各平面検出器が対応する平面放射線場の経路の計算の複数の放射線場交点をそれぞれ調べ、経路の計算は第一検出平面及び第二検出平面で第一検出平面の第二軸及び第二検出平面の第三軸にそれぞれ平行する。

ここで、図５を参照しながら、第一検出平面Ｐ１を説明する。第一検出平面Ｐ１で平面検出器４０ｊが対応する平面放射線場２１２ａの経路の計算ＰＡ上での第一放射線場交点ＰＡ１及び第二放射交点ＰＡ２を調べ、経路の計算ＰＡは第一検出平面Ｐ１上で第一検出平面Ｐ１の第二軸Ｘに平行する。同様に、工程Ｓ１４４で、平面検出器４０ｊが第二軸Ｘ(即ちＸ軸)まで回転すると、上述の経路の計算は第一検出平面で第一検出平面の第一軸Ｚに平行する。

他の実施形態では、第二検出平面Ｐ２を例にすると、平面検出器４０ｊが第一軸Ｚ(即ちＺ軸)まで回転し、経路の計算は第二検出平面Ｐ２上で第二検出平面Ｐ２の第三軸Ｙに平行する。同様に、平面検出器４０ｊが第三軸Ｙ(即ちＹ軸)まで回転すると、経路の計算は第二検出平面Ｐ２上で第二検出平面Ｐ２の第一軸Ｚに平行する。上述の状況では、本発明には制限がなく、実際の使用の要求に基づいて決定される。

このほか、上述の工程Ｓ１４６２を実行した後、三次元画像の投影方法Ｓ１００は工程Ｓ１４６２ａ乃至工程Ｓ１４６２ｂをさらに含む。図７は本発明の計算平面放射線場の経路の計算での長さ及び放射線場交点にある位置関数のフローチャート図である。

工程Ｓ１４６２ａを実行し、位置関数変数に基づいて、経路の計算ＰＡが平面放射線場２１２ａで占有する長さＬを比例により算出する。例えば、図８は本発明の経路の計算が平面放射線場で占有する長さを比例により算出するフローチャート図である。

図８及び図５を参照して、まず工程Ｓｔ１を実行し、各平面検出器４０ｊ及び対応する平面放射線場２１２ａにより包囲されて形成される第一領域Ｚ１を調べる。換言すれば、第一領域Ｚ１は光源２０、平面放射線場２１２ａ及び第一軸Ｚまで回転した後の平面検出器４０ｊにより包囲される三角形領域である。

次に、工程Ｓｔ２を実行し、平面放射線場２１２ａが経路の計算ＰＡで包囲して形成される第二領域Ｚ２を調べる。換言すれば、第二領域Ｚ２は光源２０、平面放射線場２１２ａ及び経路の計算ＰＡにより包囲されて形成される三角形領域である。

その後、工程Ｓｔ３を実行し、位置関数変数に基づいて、第一領域Ｚ１及び第二領域Ｚ２の比例を比較して経路の計算ＰＡが平面放射線場２１２ａで占有する長さＬを獲得する。図５を例にすると、平面検出器４０ｊの表面の中心点ＰＣ１の位置関数変数は(０，ｈ)であり、グリッド３２０の表面の中心点ＰＣ２の位置関数変数は(ｘ'＿ｒｏｔ，ｚ'＿ｒｏｔ)である。このほか、検出器ｊの長さは例えば２ｄであり、平面検出器４０ｊが第一軸Ｚまで回転した後、平面検出器４０ｊの長さは２ｄ＿ｒｏｔに修正され、修正差異値に関する数式は下述のとおりである：

上述の数式(３)において、距離が遠くなるほど、修正差異値が小さくなる。実際の運用では、光源２０と検出器ｊとの間の距離は非常に遠く、ゆえに、上述の数式(３)で表される修正差異値は小さくなり、このため計算の差異に影響しない。

上述のように、第一領域Ｚ１及び第二領域Ｚ２の比例を比較し、経路の計算ＰＡが平面放射線場２１２ａで占有する長さＬを獲得し、関連する数式は下述のとおりである：

図７を参照しながら、工程Ｓ１４６２ａを実行し、上述の経路の計算ＰＡが平面放射線場２１２ａで占有する長さＬを獲得した後、工程Ｓ１４６２ｂを実行し、経路の計算が平面放射線場で占有する長さに基づいて、複数の放射線場交点の位置関数を算出する。

なお、図５を例にすると、前述の放射線場交点は第一放射線場交点ＰＡ１及び第二放射線場交点ＰＡ２である。第一放射線場交点ＰＡ１と第二放射交点ＰＡ２との間の距離は上述の経路の計算ＰＡが平面放射線場２１２ａで占有する長さＬである。したがって、第一放射線場交点ＰＡ１の位置関数は（−Ｌ／２，ｚ＿ｒｏｔ）であり、第二放射交点ＰＡ２の位置関数は（Ｌ／２，ｚ＿ｒｏｔ）である。

図６を参照しながら、工程Ｓ１４６２が完了し、即ち、平面検出器４０ｊが対応する平面放射線場２１２ａの経路の計算ＰＡ上での第一放射線場交点ＰＡ１及び第二放射交点ＰＡ２を調べた後、続いて工程Ｓ１４６４を実行し、位置関数変数に基づいて各グリッドが第一検出平面及び第二検出平面で経路の計算上で対応する複数の特定の交点をそれぞれ調べる。

第一検出平面Ｐ１を例にすると、図９は本発明の位置関数変数に基づいてグリッドが第一検出平面で経路の計算上で対応する特定の交点を調べるフローチャート図である。また、図１０は図５の経路の計算上で交点線段を調べて計算する概念図である。ちなみに、図１０は光源２０、一部の平面放射線場２１２ａ及び経路の計算ＰＡのみを図示し、即ち、第二領域Ｚ２は以下で説明する。

まず、工程Ｓ１４６４ａを実行し、境界距離ｒ１を画定させる。続いて、工程Ｓ１４６４ｂを実行し、経路の計算で選択された距離を選択し、選択された距離は境界距離に等しいかそれより大きい。次に、工程Ｓ１４６４ｃを実行し、経路の計算でグリッドを中心として選択された距離を半径の範囲内とし、平面放射場内に特定の交点があるか否かを判断する。図１０に示すように、経路の計算ＰＡで、グリッド３２０の表面の中心点ＰＣ２を中心とし、選択された距離ｒ２は半径であり、第一特定の交点ＰＤ１及び第二特定の交点ＰＤ２をさらに描写する。また、図１０によると、前述のグリッド３２０の表面の中心点ＰＣ２の位置関数変数は(ｘ＿ｒｏｔ，ｚ＿ｒｏｔ)であり、第一特定の交点ＰＤ１の位置関数は(ｘ＿ｒｏｔ−ｒ２，ｚ＿ｒｏｔ)であり、第二特定の交点ＰＤ２の位置関数は(ｘ＿ｒｏｔ＋ｒ２，ｚ＿ｒｏｔ)である。

ちなみに、本実施形態は前述のボクセル(voxel)のグリッドの中心点を例にする。このほか、他の実施形態では、ボクセル(voxel)のグリッドの境界を基準点として画定させてもよく、実際の計算時の要求によって決定される。

上述の工程Ｓ１４６４ｃ及び図１０によれば、平面放射場２１２ａ内に第二特定の交点ＰＤ２を有する。続いて、図６に示すように、工程Ｓ１４６６を実行し、計算を行い、複数の放射線場交点の位置関数及び複数の特定の交点の位置関数を比較し、経路の計算で各グリッドが各平面検出器が対応する平面放射線場で占める幾何学的比率値を決定する。

図１０を例にすると、第一検出平面ＰＡで、第一特定の交点ＰＤ１の位置関数は(ｘ＿ｒｏｔ−ｒ２，ｚ＿ｒｏｔ)であり、第二特定の交点ＰＤ２の位置関数は(ｘ＿ｒｏｔ＋ｒ２，ｚ＿ｒｏｔ)であり、第一放射線場交点ＰＡ１の位置関数は（−Ｌ／２，ｚ＿ｒｏｔ）であり、第二放射交点ＰＡ２の位置関数は（Ｌ／２，ｚ＿ｒｏｔ）であり、経路の計算ＰＡが平面放射線場２１２ａで長さＬを占有する。ここでは、上述の特定の交点及び放射線場交点は共に経路の計算ＰＡに位置され、ゆえに、計算上第二軸Ｘ上の位置関数を考慮し、関連する数式は下述する：

上述の数式(６)からわかるように、第一放射線場交点ＰＡ１及び第一特定の交点ＰＤ１の両者の位置関数を比較して両者の位置関数の最大値を獲得し、第二放射交点ＰＡ２及び第二特定の交点ＰＤ２の両者の位置関数を比較して両者の位置関数の最小値を獲得する。前述の最大値及び最小値の位置関数は互いに減少すると共にその値を正数に確保し、最後に経路の計算ＰＡが平面放射線場２１２ａで占有する長さＬを除算し、第一検出平面Ｐ１でのグリッド３２０に関して平面検出器４０ｊが対応する平面放射線場２１２ａで占める幾何学的比率値ｇｘｚを獲得する。

同様の原理及び工程により、第二検出平面Ｐ２に処理を施し、幾何学的比率値ｇｙｚを獲得する。

図３を参照しながら、工程Ｓ１４８を実行する。第一検出平面Ｐ１及び第二検出平面Ｐ２が相互に対応するグリッドが有する幾何学的比率値を相乗して各検出器が対応する空間ユニットが有するサブ幾何学的因子ｇｉｊを獲得する。

図６に示すように、上述のサブ幾何学的因子ｇｉｊを獲得した後、サブ幾何学的因子ｇｉｊに加重係数Ｗを乗算する工程Ｓ１４６８を実行し、加重係数Ｗの計算の説明は図１１Ａ及び図１１Ｂのとおりである。図１１Ａはボクセルを示す概念図であり、図１１Ｂは１つのボクセルが任意の回転がなく放射線場と交差する概念図である。

図１１Ａに示すように、前述の工程Ｓ１４６４ａに示すように、境界距離ｒ１を画定する工程では、本実施形態では、ボクセル(voxel)３０の表面の中心点Ｇ１からボクセル(voxel)３０の境界までの垂直距離を境界距離ｒ１とする。図１１Ｂは１つのボクセル(voxel)３０が任意の回転がなく放射線場２１０と交差する状況を示し、第一交点ＰＥ１の位置関数は(ｘ＿ｒｏｔ-ｒ１，ｚ＿ｒｏｔ)であり、第二交点ＰＥ２の位置関数は(ｘ＿ｒｏｔ＋ｒ１，ｚ＿ｒｏｔ)であり、第一放射線場交点ＰＡ１の位置関数は（−Ｌ／２，ｚ＿ｒｏｔ）であり、第二放射交点ＰＡ２の位置関数は（Ｌ／２，ｚ＿ｒｏｔ）である。このため、放射線場２１０のボクセル３０が交差するＸ座標(即ち第二軸)部分は第一交点ａ及び第二交点ｂであり、且つ第一交点ａの位置関数は（−Ｌ／２，ｚ＿ｒｏｔ）であり、第二交点ｂの位置関数は（ｘ＿ｒｏｔ＋ｒ１，ｚ＿ｒｏｔ）である。計算では、使用者により定められた重量方程式ｆ１(ｘ)に基づいて計算し、この重量方程式ｆ１(ｘ)の重量(weighting)は１である。図１１Ｂの下方について説明すると、計算の面積Ｗはボクセル３０の確率であり、計算の数式は以下のとおりである：

このように、計算は第一交点ａから第二交点ｂに積分する。操作の実施では、ボクセル３０は回転する。図１１Ｃはボクセルが回転した後に放射線場と交差する概念図である。このさい、前述の境界距離ｒ１を計算とすると、放射線場２１０と交差せず、計算の誤差が発生する。

このため、本実施形態では１つの選択された距離ｒ２が定義され、前述の選択された距離ｒ２は境界距離ｒ１より大きい。例えば、第一特定の交点ＰＤ１の位置関数が(ｘ＿ｒｏｔ−ｒ２，ｚ＿ｒｏｔ)であり、第二特定の交点ＰＤ２の位置関数が(ｘ＿ｒｏｔ＋ｒ２，ｚ＿ｒｏｔ)である場合、計算では、使用者により定められた重量方程式f２(ｘ)に基づいて計算し、重量方程式f２(ｘ)は回転角度に関係し、計算結果は図１１Ｃの下方に示す。第一交点ａを第二交点ｂに積分し、積分面積を求める。

前述の重量方程式ｆ１(ｘ)及び重量方程式f２(ｘ)の総積分結果は同じにならなければならず、図１１Ｂ及び図１１Ｃを例にすると、総積分結果は共にである。

実施では、異なる第一交点ａ及び第二交点ｂに基づいてセルが構成され、セルを直接調べて加重係数Ｗを獲得する。この方式により、ボクセル３０が任意で回転せず、境界距離ｒ１が選択されると、重量方程式ｆ１(ｘ)は常に１になり、ゆえに、数式(７)の(ｂ−ａ)を採用して加重係数Ｗを獲得する。このほか、前述の選択された距離ｒ２が境界距離ｒ１より大きい方式を採用すると、上述のセルを調べる方式により加重係数Ｗを直接獲得し、これにより選択された距離ｒ２の部分の計算量を減少させる。

このほか、工程Ｓ１４８を実行した後、本実施形態に係る三次元画像の投影方法は以下の工程をさらに含む。図１２は本発明に係る三次元画像の投影方法の均一化パラメータ工程を示す概念図である。図１２に示すように、工程Ｓ１４８２を実行し、第一検出平面及び第二検出平面の第一軸方向に平行する前記複数のグリッドの数をそれぞれ計算する。本実施形態では、第一軸Ｚ方向に沿うグリッドの数は５つである（図２参照）。

続いて、サブ幾何学的因子に均一化パラメータを乗算する工程Ｓ１４８４を実行し、前記均一化パラメータは１と前記複数のグリッドの数との比率である。本実施形態では、均一化パラメータＮは１／５である。

なお、図１を参照しながら、工程Ｓ１４０及び関連する工程を実行した後、サブ幾何学的因子ｇｉｊを獲得する。最後に工程Ｓ１５０を実行し、各検出器ｊがサブ放射線場２１２で対応する各空間ユニットｉが有するサブ幾何学的因子ｇｉｊが組み合わされて幾何学的因子マトリックスＧｉｊが形成される。

結論として、本発明に係る三次元画像の投影方法は、撮影領域の医学影像再構成に必要な幾何因子が構築される。本発明は、まず、本来の三次元空間の数学モデルが２つの二次元の検出平面に投影され、平面検出器が検出平面の軸まで回転し、これにより煩雑な演算を簡略化させ、計算に掛かる時間を短縮させる。

上述の実施形態は本発明の技術思想及び特徴を説明するためのものにすぎず、当該技術分野を熟知する者に本発明の内容を理解させると共にこれをもって実施させることを目的とし、本発明の特許請求の範囲を限定するものではない。従って、本発明の精神を逸脱せずに行う各種の同様の効果をもつ改良又は変更は、後述の請求項に含まれるものとする。

２０ 光源
２１０ 放射線場
２１２ サブ放射線場
２１２ａ 平面放射線場
２１２ｂ 平面放射線場
３０ ボクセル
３２ 一次元配列
３４ 一次元配列
３２０ グリッド
３４０ グリッド
４０ 検出モジュール
４０ｊ 平面検出器
ａ 第一交点
ｂ 第二交点
２ｄ 検出器の長さ
２ｄ＿ｒｏｔ 平面検出器の長さ
ｈ 距離
Ｇ１ ボクセル表面の中心点
Ｇｉｊ 幾何学的因子マトリックス
ｇｉｊ サブ幾何学的因子
ｇｘｚ、ｇｙｚ 幾何学的比率値
ｉ 空間ユニット
ｊ 検出器
Ｌ 長さ
Ｐ１ 第一検出平面
Ｐ２ 第二検出平面
ＰＡ 経路の計算
ＰＡ１ 第一放射線場交点
ＰＡ２ 第二放射線場交点
ＰＣ１ 平面検出器表面の中心点
ＰＣ２ グリッド表面の中心点
ＰＤ１ 第一特定の交点
ＰＤ２ 第二特定の交点
ＰＥ１ 第一交点
ＰＥ２ 第二交点
ｒ１ 境界距離
ｒ２ 選択された距離
Ｓ１００ 三次元画像の投影方法
Ｓ１１０〜Ｓ１５０ 工程
Ｓ１４２〜１４８ 工程
Ｓ１４６２〜１４６８ 工程
Ｓ１４６２ａ〜Ｓ１４６２ｂ 工程
Ｓ１４８２〜Ｓ１４８４ 工程
Ｓｔ１〜Ｓｔ３ 工程
ｔ 時間
Ｗ 加重係数
Ｙ 第三軸
Ｘ 第二軸
Ｚ 第一軸
Ｚ１ 第一領域
Ｚ２ 第二領域
θ 夾角

Claims (9)

1. 三次元画像の投影方法であって、
   光源により放射線場が発生して検出モジュールに投射され、前記検出モジュールはアレイ状に配列される複数の検出器を有する、光源を提供する工程と、
   前記光源と前記検出モジュールとの間には複数の空間ユニットが堆積されて形成される三次元空間配列のボクセルを有する工程と、
   前記放射線場の前記検出器と前記光源との間のサブ放射線場を各々検出して定義する工程と、
   各検出器が前記サブ放射線場で対応する各空間ユニットに基づいて演算し、サブ幾何学的因子を獲得する工程とを含み、
   ここでは、前記演算は、
   前記検出器、各検出器に対応する前記サブ放射線場及び前記複数の空間ユニットが第一検出平面及び第二検出平面にそれぞれ投影され、前記検出器、各検出器に対応する前記サブ放射線場及び前記複数の空間ユニットは前記第一検出平面及び前記第二検出平面に対応する平面検出器、平面放射線場及び隣接されて配列されて二次元平面となる複数の一次元配列をそれぞれ有し、各一次元配列は隣接して配列される複数のグリッドを有し、前記第一検出平面は第一軸及び第二軸で構成され、前記第二検出平面は前記第一軸及び第三軸で構成される工程と、
   前記平面検出器が前記第一検出平面の前記第一軸及び前記第二検出平面の前記第一軸までそれぞれ回転し、前記平面検出器及び前記複数のグリッドが対応する位置関数変数を発生させる工程と、
   前記位置関数変数に基づいて前記第一検出平面及び前記第二検出平面上での各グリッドに関して各平面検出器が対応する前記平面放射線場で占める幾何学的比率値をそれぞれ計算する工程と、
   前記第一検出平面及び前記第二検出平面上の相互に対応するグリッドが有する幾何学的比率値を相乗して各検出器が対応する空間ユニットが有する前記サブ幾何学的因子を獲得する工程と、
   各検出器が前記サブ放射線場で対応する各空間ユニットが有する前記サブ幾何学的因子が組み合わされて幾何学的因子マトリックスが形成される工程とをさらに含むことを特徴とする三次元画像の投影方法。
2. 前記位置関数変数に基づいて前記第一検出平面及び前記第二検出平面上での各グリッドに関して各平面検出器が対応する前記平面放射線場で占める幾何学的比率値をそれぞれ計算する工程は、
   前記第一検出平面及び前記第二検出平面で各平面検出器が対応する前記平面放射線場の経路の計算上の複数の放射線場交点をそれぞれ調べ、前記経路の計算は前記第一検出平面及び前記第二検出平面で前記第一検出平面の前記第二軸及び前記第二検出平面の前記第三軸にそれぞれ平行する工程と、
   前記位置関数変数に基づいて各グリッドが前記第一検出平面及び前記第二検出平面で前記経路の計算上の対応する複数の特定の交点をそれぞれ調べる工程と、
   計算を行い、前記複数の放射線場交点の位置関数及び前記複数の特定の交点の位置関数を比較し、前記経路の計算上で前記各グリッドが各平面検出器が対応する前記平面放射線場で占める幾何学的比率値を決定する工程とを含むことを特徴とする、請求項１に記載の三次元画像の投影方法。
3. 前記第一検出平面及び前記第二検出平面で各平面検出器が対応する前記平面放射線場の前記経路の計算上の前記複数の放射線場交点をそれぞれ調べる工程の後には、
   前記位置関数変数に基づいて、前記経路の計算が前記平面放射線場で占有する長さを比例により算出する工程と、
   前記経路の計算が前記平面放射線場で占有する長さに基づいて、前記複数の放射線場交点の位置関数を算出する工程とをさらに含むことを特徴とする、請求項２に記載の三次元画像の投影方法。
4. 前記経路の計算が前記平面放射線場で占有する長さを比例により算出する工程は、
   各平面検出器及び対応する前記平面放射線場により包囲されて形成される第一領域を調べる工程と、
   前記平面放射線場の前記経路の計算により包囲されて形成される第二領域を調べる工程と、
   前記位置関数変数に基づいて、前記第一領域及び前記第二領域の比例を比較して前記経路の計算が前記平面放射線場で占有する長さを獲得する工程とをさらに含むことを特徴とする、請求項３に記載の三次元画像の投影方法。
5. 前記計算を行い、前記複数の放射線場交点の位置関数及び前記複数の特定の交点の位置関数を比較し、前記経路の計算上の前記各グリッドが各平面検出器が対応する前記平面放射線場で占める幾何学的比率値を決定する工程の前には、
   境界距離が画定される工程と、
   前記経路の計算で選択された距離を選択し、前記選択された距離は前記境界距離に等しいかそれより大きい工程と、
   前記経路の計算では前記グリッドを中心として前記選択された距離は半径の範囲内になり、前記平面放射場内に前記特定の交点があるか否かを判断する工程とをさらに含むことを特徴とする、請求項２に記載の三次元画像の投影方法。
6. 前記サブ幾何学的因子に加重係数を乗算する工程をさらに含むことを特徴とする、請求項２に記載の三次元画像の投影方法。
7. 前記第一検出平面及び前記第二検出平面上で相互に対応するグリッドが有する幾何学的比率値を相乗して各検出器が対応する空間ユニットが有する前記サブ幾何学的因子を獲得する工程の後には、
   前記第一検出平面及び前記第二検出平面の前記第一軸方向に平行する前記複数のグリッドの数をそれぞれ計算する工程と、
   前記サブ幾何学的因子に均一化パラメータを乗算する工程とをさらに含み、前記均一化パラメータは１と前記複数のグリッドの数との比率であることを特徴とする、請求項１に記載の三次元画像の投影方法。
8. 前記光源は回転運動をさらに行うことを特徴とする、請求項１に記載の三次元画像の投影方法。
9. 前記検出モジュールはＸ線断層撮影合成法或いはＣＴスキャンの検出モジュールであることを特徴とする、請求項１に記載の三次元画像の投影方法。

Images