JP2010125116A - ３次元ずれ量計測方法 - Google Patents

Abstract

【課題】アーチファクトの影響を受けず、短時間でかつ正確に、レジストレーションを行う。
【解決手段】基準立体画像を構成する２次元の投影画像をデータベースＤＢ１に保存する。比較対象立体画像を構成する２次元の投影画像をデータベースＤＢ２に保存する。データベースＤＢ１，ＤＢ２に保存されている２次元の投影画像から基準立体画像と比較対象立体画像との間のＸ，Ｙ，Ｚの３軸の回転ずれ量（Δθ，Δφ，ΔΨ）および平行ずれ量（Δｘ，Δｙ，Δｚ）を求める。
【選択図】 図１３

Description

この発明は、基準となる３次元の立体画像を基準立体画像、この基準立体画像に対して比較対象とする３次元の立体画像を比較対象立体画像とし、基準立体画像と比較対象立体画像との間のずれ量を計測する３次元ずれ量計測方法に関するものである。

従来より、医学分野などにおいては、多様な画像化方式によって立体画像を得ている。例えば、磁気共鳴断層撮影法（ＭＲＩ）、Ｘ線コンピュータ断層撮影法（ＣＴ）、陽電子放射断層撮影法（ＰＥＴ）、単光子射出断層撮影法（ＳＰＥＣＴ）、超音波断層撮影法などによって、３次元の立体画像を得ている。

このような方法で得られた３次元の立体画像は、ボクセル（voxel）値の配列として、通常、コンピュータのメモリなどにデジタル的に保存される。各ボクセルは、互いに直交するＸ，Ｙ，Ｚの３軸で表される３次元空間の座標位置（ｘ，ｙ，ｚ）と結びついており、ボクセル値として色値（典型的にはグレイスケール）が割り当てられる。

図１４にＸ線コンピュータ断層撮影法による３次元の立体画像の取得例を示す。この例では、スキャナ装置としてコーンビームＣＴを用いている。コーンビームＣＴは、中間点Ｃ０を挾んで対向してＸＹ平面内に設けられたＸ線源１と２次元撮像装置２とを備え、このＸ線源１と２次元撮像装置２とを中間点Ｃ０を通るＺ軸を中心として水平方向に定速度で回転させながら、被写体（人体）３の撮影を行う。この撮影に際して、被写体３は、その頭部の中心を中間点Ｃ０に位置させる。

被写体３の撮影は、Ｘ線源１と２次元撮像装置２とが１回転する間、所定の時間間隔で行われる。この結果、１回転する間に任意の異なる角度（以下、この角度を投影角度と呼ぶ）から撮影されたｎ個（例えば、２５６個）の２次元の投影画像を取得して、このｎ個の投影画像を再構成処理することで３次元の立体画像を得ることができる。この３次元の立体画像は、互いに直交するＸ，Ｙ，Ｚの３軸で表される３次元空間の座標位置（ｘ，ｙ，ｚ）と結びつけたボクセルデータ群として表現される。

この３次元の立体画像は、医療や歯科などの分野において、各種の診断を行うために使用される。しかし、情報量が多いために見落としが発生する虞があり、正確な診断を行うために過大な時間を必要とする問題がある。このため、コンピュータ上で、画像処理を行うことで、この問題を解決しようと試みられている。

この場合の画像処理の１つに、レジストレーションと呼ばれる位置合わせを行った後、２つの立体画像を比較するという方法がある。例えば、前回撮影した被写体３（被写体Ａ）の３次元の立体画像を基準立体画像とし、今回撮影した被写体３（被写体Ｂ）の３次元の立体画像を比較対象立体画像とし、基準立体画像と比較対象立体画像との位置合わせを行い、この位置合わせが行われた２つの立体画像の差分画像をディスプレイ上に表示する。これにより、前回撮影時と今撮影時との間の変化部分のみを抽出し、短時間でかつ正確な診断を行うことができるようになる。

この方法において、基準立体画像と比較対象立体画像とのレジストレーションは、６つのパラメータを求めることで実現される。６つのパラメータは、Ｘ，Ｙ，Ｚの各軸の平行ずれ量Δｘ，Δｙ，Δｚと、Ｘ，Ｙ，Ｚの各軸の回転ずれ量Δθ，Δφ，ΔΨである。基準立体画像と比較対象立体画像との間で、これら６つのパラメータを求めるための手法としては、例えば、基準立体画像に対して特徴のある局所領域を定め、この局所領域が比較対象立体画像のどの部分に存在するかという検索を行い、その局所領域の移動量から全体のずれ量を求め、６つのパラメータを取得する方法が一般的である。（例えば、特許文献１参照）。

特開平９−２２４０６号公報 特開２００１−２１２１２６号公報 特開平１０−１２４６６７号公報

しかしながら、上述した基準立体画像と比較対象立体画像とのレジストレーションでは、局所領域の検索を行うために、繰り返し処理の実行回数が比較的多くなり、処理時間が長くなる傾向にある。また、処理時間が一定ではなく、対象によって時間がかかる場合がある。

また、局所領域は、どの部分でもいいという訳ではなく、できるだけ特徴がある部分を選択しなければ、正しい検索結果が得られない。したがって、対象は特徴が抽出できる立体画像に限られる。

また、３次元のボクセルデータ群は情報量が多く、この情報量が多いボクセルデータ群を使用してレジストレーションを行うため、各種処理に膨大な演算を必要とし、高速処理が困難となる。

また、３次元のボクセルデータ群は、２次元の投影画像から再構成処理を行うことで得られるが、この再構成処理の段階で、各種アーチファクト（偽像）が発生することがある。アーチファクトが発生すると、その影響がボクセルデータ群に現れ、正確なレジストレーションの妨げとなる。特に、金属アーチファクトは激しいストリーク状のノイズが発生するため、レジストレーションに大きな影響を与える。

なお、局所領域を用いないレジストレーションとして、立体画像に現れている被写体の構造に詳しい専門家が、位置合わせする立体画像のそれぞれに目印をつけ、一方の立体画像を移動回転させることによって位置合わせする方法がある（例えば、特許文献２参照）。

しかしながら、この方法では、位置合わせ精度が選択された目印の数と位置によって制限されてしまう。また、目印の数が少な過ぎると、位置合わせが不正確になる。また、目印が多くても必ずしも位置合わせが正確となるわけではなく、位置合わせに要する計算を複雑にするだけである。また、オペレータが手作業で対応点を選択する必要があり、時間がかかるとともに熟練が必要となる。さらに、すべての立体画像で適当な構造上の目印を見出せるわけでもない。

本発明は、このような課題を解決するためになされたもので、その目的とするところは、アーチファクトの影響を受けず、短時間でかつ正確に、レジストレーションを行うことが可能な３次元ずれ量計測方法を提供することにある。

このような目的を達成するために本発明は、基準となる３次元の立体画像を基準立体画像、この基準立体画像に対して比較対象とする３次元の立体画像を比較対象立体画像とし、基準立体画像と比較対象立体画像との間のずれ量を計測する３次元ずれ量計測方法において、基準立体画像を構成する２次元の投影画像と比較対象立体画像を構成する２次元の投影画像とから、基準立体画像と比較対象立体画像との間の互いに直交するＸ，Ｙ，Ｚの３軸の回転ずれ量および平行ずれ量を計測するようにしたものである。

この発明によれば、再構成処理後の３次元の立体画像（ボクセルデータ群）ではなく、再構成処理前の２次元の投影画像を使用して、基準立体画像と比較対象立体画像との間の回転ずれ量（Δθ，Δφ，ΔΨ）と平行ずれ量（Δｘ，Δｙ，Δｚ）を求めることが可能である。このため、３次元画像処理を必要とせず、局所領域や目印を定める必要もなく、高速処理が実現可能となる。また、再構成処理で発生するアーチファクトの影響も受けない。

本発明において、レジストレーションに必要な６つのパラメータは、例えば次のようにして求める。基準立体画像を構成する２次元の投影画像中の所定の投影角度から見た投影画像を第１の基準投影画像として抽出する（第１ステップ）。比較対象立体画像を構成する２次元の投影画像中の第１の基準投影画像と同じ投影角度付近の複数の投影画像を第１の比較候補投影画像として抽出する（第２ステップ）。第１の基準投影画像と第１の各比較候補投影画像との間の回転角度のずれを求める（第３ステップ）。第１の各比較候補投影画像を第３ステップで求めた回転角度のずれ分回転補正して第１の比較候補補正投影画像とする（第４ステップ）。第１の基準投影画像と第１の各比較候補補正投影画像との相関値を求め、その中の相関値が最も大きい比較候補補正投影画像を第１の比較確定投影画像とする（第５ステップ）。第１の比較確定投影画像の作成元の第１の比較候補投影画像の投影角度と第１の基準投影画像の投影角度との差をＺ軸を中心とする回転ずれ量ΔΨとして求める（第６ステップ）。また、第１の比較確定投影画像の作成元の第１の比較候補投影画像と第１の基準投影画像との間の回転角度のずれをＸ軸を中心とする回転ずれ量Δθとして求める（第７ステップ）。第１の比較確定投影画像の第１の基準投影画像との間の水平方向のずれ量をＹ軸方向への平行ずれ量Δｙ、垂直方向のずれ量をＺ軸方向への平行ずれ量Δｚとして求める（第８ステップ）。基準立体画像を構成する２次元の投影画像中の第１の基準投影画像に対して投影角度が９０゜ずれた投影画像を第２の基準投影画像として抽出する（第９ステップ）。比較対象立体画像を構成する２次元の投影画像中の第１の比較確定投影画像の作成元の第１の比較候補投影画像に対してその投影角度が第２の基準投影画像と同方向に９０゜ずれた投影画像を第２の比較候補投影画像として抽出する（第１０ステップ）。第２の基準投影画像と第２の比較候補投影画像との間の回転角度のずれをＹ軸を中心とする回転ずれ量Δφとして求める（第１１ステップ）。第２の比較候補投影画像を回転ずれ量Δφ分回転補正して第２の比較確定投影画像とする（第１２ステップ）。第２の基準投影画像と第２の比較確定投影画像との間の水平方向のずれ量をＸ軸方向への平行ずれ量Δｘとして求める（第１３ステップ）。

なお、上述した第６ステップに代えて、第１の基準投影画像と第１の各比較候補補正投影画像との相関値から最大相関値を有する比較候補補正投影画像を推測し、この推測した比較候補補正投影画像の投影角度と第１の基準投影画像の投影角度との差をＺ軸を中心とする回転ずれ量ΔΨとして求めるステップを設けるようにしてもよい。

また、上述した第１０ステップに代えて、比較対象立体画像を構成する２次元の投影画像中に第１の比較確定投影画像の作成元の第１の比較候補投影画像に対してその投影角度が第２の基準投影画像と同方向に９０゜ずれた投影画像が実在しない場合、その投影角度が近い実在する複数の２次元の投影画像から補間して第２の比較候補投影画像を作成するステップを設けるようにしてもよい。

本発明によれば、再構成処理後の３次元の立体画像（ボクセルデータ群）ではなく、再構成処理前の２次元の投影画像を使用して、基準立体画像と比較対象立体画像との間の回転ずれ量（Δθ，Δφ，ΔΨ）と平行ずれ量（Δｘ，Δｙ，Δｚ）を求めることが可能となるので、アーチファクトの影響を受けず、短時間でかつ正確に、レジストレーションを行うことができるようになる。

以下、本発明を図面に基づいて詳細に説明する。図１は本発明に係る３次元ずれ量計測方法の実施に用いるＸ線コンピュータ断層撮影システムの一例を示す図である。同図において、図１４と同一符号は図１４を参照して説明した構成要素と同一或いは同等構成要素を示し、その説明は省略する。

この実施の形態において、Ｘ線源１および２次元撮像装置２は、３次元ずれ量計測装置４に接続されている。３次元ずれ量計測装置４は、ＣＰＵ４Ａと、ＲＡＭ４Ｂと、ＲＯＭ４Ｃと、ハードディスクなどの記憶装置４Ｄと、インタフェース４Ｅ〜Ｈと、ディスプレイ４Ｉと、キーボード４Ｊと、マウス４Ｋとを備えている。

ＣＰＵ４Ａは、ＲＡＭ４Ｂにアクセスしながら、ＲＯＭ４Ｃや記憶装置４Ｄに格納されたプログラムに従って動作する。記憶装置４Ｄには、本実施の形態特有のプログラムとして３次元ずれ量計測プログラムが格納されている。この３次元ずれ量計測プログラムは、例えばＣＤ−ＲＯＭなどの記録媒体に記録された状態で提供され、この記録媒体から読み出されて記憶装置４Ｄにインストールされている。

以下、図２〜図６に示すフローチャートを参照して、記憶装置４Ｄに格納された３次元ずれ量計測プログラムに従ってＣＰＵ４Ａが実行する本実施の形態特有の処理動作について説明する。

なお、この３次元ずれ量計測プログラムに従う処理動作には、基準立体画像と比較対象立体画像との間の６つのパラメータ（Δθ，Δφ，ΔΨ，Δｘ，Δｙ，Δｚ）を求めるずれ量計測処理と、この求めた６つのパラメータを用いて基準立体画像と比較対象立体画像との間の位置合わせを行う位置合わせ処理と、この位置合わせされた基準立体画像と比較対象立体画像との照合を行う照合処理とが含まれる。

〔基準立体画像の取得〕
先ず、基準立体画像を取得するために、被写体Ａの撮影を行う（図２：ステップＳ１０１）。この被写体Ａの撮影は、Ｘ線源１と２次元撮像装置２とが１回転する間、所定の時間間隔で行う。この結果、１回転する間の任意の角度（投影角度）から撮影されたｎ個（例えば、２５６個）の２次元の投影画像が得られる（ステップＳ１０２）。

ＣＰＵ４Ａは、この得られたｎ個の２次元の投影画像を記憶装置４Ｄに保存した後（ステップＳ１０３）、このｎ個の投影画像を再構成処理することで３次元の立体画像（ボクセルデータ群）を得る（ステップＳ１０４）。そして、この再構成処理によって得られた３次元の立体画像を基準立体画像として記憶装置４Ｄに保存する（ステップＳ１０５）。

〔比較対象立体画像の取得〕
次に、比較対象立体画像を取得するために、被写体Ｂの撮影を行う（図３：ステップＳ２０１）。この被写体Ｂの撮影も、Ｘ線源１と２次元撮像装置２とが１回転する間、基準立体画像の取得時と同じ所定の時間間隔で行う。この結果、１回転する間の任意の角度（投影角度）から撮影されたｎ個（例えば、２５６個）の２次元の投影画像が得られる（ステップＳ２０２）。

ＣＰＵ４Ａは、この得られたｎ個の２次元の投影画像を記憶装置４Ｄに保存した後（ステップＳ２０３）、このｎ個の投影画像を再構成処理することで３次元の立体画像（ボクセルデータ群）を得る（ステップＳ２０４）。そして、この再構成処理によって得られた３次元の立体画像を比較対象立体画像として記憶装置４Ｄに保存する（ステップＳ２０５）。

〔基準立体画像と比較対象立体画像との間のずれ量の計測〕
本実施の形態の３次元ずれ量計測装置４では、再構成処理された立体画像ではなく、再構成処理に用いられる前の２次元の投影画像から、基準立体画像と比較対象立体画像との間のずれ量を求める。

先ず、ＣＰＵ４Ａは、記憶装置４Ｄにアクセスし、被写体Ａの２次元の投影画像中、所定の投影角度ｔ（この例では、ｔ＝０゜）から見た投影画像を第１の基準投影画像Ｇ１Ａ（図７参照）として抽出する（図４：ステップＳ３０１）。

また、被写体Ｂの２次元の投影画像中、第１の基準投影画像Ｇ１Ａと同じ投影角度ｔ付近の複数の投影画像を第１の比較候補投影画像Ｇ１Ｂとして抽出する（ステップＳ３０２）。この例では、投影角度ｔ付近の複数の投影画像として、投影角度ｔ１，ｔ２，ｔ３の３つの投影画像を第１の比較候補投影画像Ｇ１Ｂ₁，Ｇ１Ｂ₂，Ｇ１Ｂ₃として抽出する。なお、この場合、投影角度ｔ２は、第１の比較候補投影画像側の投影角度ｔに最も近い投影角度とし、投影角度ｔ１，ｔ３は投影角度ｔ２を中心とする前後の投影角度とする。

そして、ＣＰＵ４Ａは、回転不変型位相限定相関法（ＲＩＰＯＣ：Rotation Invariant Phase Only Correlation）や回転位相限定相関法（ＲＰＯＣ：Rotation Phase Only Correlation）、またはこれに類する２次元の回転角度を算出する方法により、第１の基準投影画像Ｇ１Ａと第１の比較候補投影画像Ｇ１Ｂ₁，Ｇ１Ｂ₂，Ｇ１Ｂ₃それぞれとの間の回転角度のずれΔθ１，Δθ２，Δθ３を求める（ステップＳ３０３）。

なお、ＲＩＰＯＣやＲＰＯＣについては、例えば特許文献３などに示されているので、ここでの詳細な説明は省略する。例えば、ＲＰＯＣでは、片方の画像を回転させながら相関値を求め、最大の相関値を持つ回転角度を求める。

次に、ＣＰＵ４Ａは、求めた回転角度のずれΔθ１，Δθ２，Δθ３分、第１の比較候補投影画像Ｇ１Ｂ₁，Ｇ１Ｂ₂，Ｇ１Ｂ₃をそれぞれ回転補正して、第１の比較候補補正投影画像Ｇ１Ｂ₁’，Ｇ１Ｂ₂’，Ｇ１Ｂ₃’を得る（ステップＳ３０４）。

そして、位相限定相関法（ＰＯＣ：Phase Only Correlation）などの２次元画像の相関値を算出する方法により、第１の基準投影画像Ｇ１Ａと第１の比較候補補正投影画像Ｇ１Ｂ₁’，Ｇ１Ｂ₂’，Ｇ１Ｂ₃’それぞれとの相関値を求め、その中の相関値が最も大きい比較候補補正投影画像を第１の比較確定投影画像Ｇ１ＢＲとする（ステップＳ３０５）。この例では、第１の比較候補補正投影画像Ｇ１Ｂ₂’が第１の基準投影画像Ｇ１Ａとの相関値が最も大きく、第１の比較確定投影画像Ｇ１ＢＲとされるものとする。

この場合、第１の比較候補補正投影画像Ｇ１Ｂ２’の相関値が最も大きいということは、第１の基準投影画像Ｇ１Ａと第１の比較候補投影画像Ｇ１Ｂ２とがほぼ同じ向きから撮影された画像であることを意味している。

次に、ＣＰＵ４Ａは、この相関値が最も大きい第１の比較候補投影画像Ｇ１Ｂ₂’を第１の比較確定投影画像Ｇ１ＢＲとし、この第１の比較確定投影画像Ｇ１ＢＲの作成元の第１の比較候補投影画像Ｇ１ＢＲ’（Ｇ１Ｂ₂）の投影角度ｔ’（ｔ２）と第１の基準投影画像Ｇ１Ａの投影角度ｔとの差を求め、これを３次元の立体画像のＺ軸を中心とする回転ずれ量ΔΨ（ΔΨ＝ｔ’−ｔ）とする（ステップＳ３０６）。

また、ＣＰＵ４Ａは、第１の比較確定投影画像Ｇ１ＢＲの作成元の第１の比較候補投影画像Ｇ１ＢＲ’と第１の基準投影画像Ｇ１Ａとの間の回転角度のずれを求め、これを３次元の立体画像のＸ軸を中心とする回転ずれ量Δθ（Δθ＝Δθ２）とする（ステップＳ３０７）。

また、ＣＰＵ４Ａは、第１の比較確定投影画像Ｇ１ＢＲと第１の基準投影画像Ｇ１Ａとの間の水平方向のずれ量Δｙを３次元の立体画像のＹ軸方向への平行ずれ量として、垂直方向のずれ量Δｚを３次元の立体画像のＺ軸方向への平行ずれ量として求める（図５：ステップＳ３０８，Ｓ３０９）。なお、この場合の平行ずれ量Δｙ，Δｚは、第１の比較候補投影画像Ｇ１Ｂ２’と第１の基準投影画像Ｇ１Ａとの相関値をＰＯＣなどの相関法によって求める際に、位置ずれ量として得ることができる。

次に、ＣＰＵ４Ａは、記憶装置４Ｄにアクセスし、被写体Ａの２次元の投影画像中、第１の基準投影画像Ｇ１Ａに対して投影角度が９０゜異なる投影画像を第２の基準投影画像Ｇ２Ａ（図８参照）として抽出する（ステップＳ３１０）。この例では、被写体Ａの２次元の投影画像中、投影角度ｔ＋９０゜の投影画像を第２の基準投影画像Ｇ２Ａとして抽出する。

また、被写体Ｂの２次元の投影画像中、第１の比較確定投影画像Ｇ１ＢＲの作成元の第１の比較候補投影画像Ｇ１ＢＲ’に対してその投影角度が第２の基準投影画像Ｇ２Ａと同方向に９０゜異なる投影画像を第２の比較候補投影画像Ｇ２Ｂとして抽出する（ステップＳ３１１）。この例では、被写体Ｂの２次元の投影画像中、投影角度ｔ’＋９０゜の投影画像を第２の基準投影画像Ｇ２Ｂとして抽出する。

そして、ＣＰＵ４Ａは、ＲＩＰＯＣやＲＰＯＣまたはこれに類する２次元の回転角度を算出する方法により、第２の基準投影画像Ｇ２Ａと第２の比較候補投影画像Ｇ２Ｂとの間の回転角度のずれΔφを求め、これを３次元の立体画像のＹ軸を中心とする回転ずれ量とする（ステップＳ３１２）。

また、第２の比較候補投影画像Ｇ２Ｂを回転ずれ量Δφ分回転補正して第２の比較確定投影画像Ｇ２Ｒとし（ステップＳ３１３）、ＰＯＣなどの相関法を用いて第２の基準投影画像Ｇ２Ａと第２の比較確定投影画像Ｇ２Ｒとの間の水平方向の位置ずれ量Δｘを求め、これを３次元の立体画像のＸ軸方向への平行ずれ量とする（ステップＳ３１４）。

〔基準立体画像と比較対象立体画像との間の位置合わせ〕
ＣＰＵ４Ａは、このようにしてレジストレーションに必要な６つのパラメータΔθ，Δφ，ΔΨ，Δｘ，Δｙ，Δｚを求めた後、この求めたパラメータΔθ，Δφ，ΔΨ，Δｘ，Δｙ，Δｚに基づいて、記憶装置４に保存されている基準立体画像（ボクセルデータ群）と比較対象立体画像（ボクセルデータ群）との間の位置合わせを行う（図６：ステップＳ４０１）。

すなわち、比較対象立体画像をＸ軸を中心としてΔθ回転補正し、Ｙ軸を中心としてΔφ回転補正し、Ｚ軸を中心としてΔΨ回転補正する。また、Ｘ軸方向へΔｘ、Ｙ軸方向へΔｙ、Ｚ軸方向へΔｚ平行移動し、基準立体画像と比較対象立体画像とのＸ，Ｙ，Ｚの３次元空間における位置を合わせる。なお、この際、比較対象立体画像ではなく、基準立体画像の方を回転補正したり、平行移動するようにしてもよい。

〔位置合わせされた基準立体画像と比較対象立体画像との照合〕
ＣＰＵ４Ａは、このようにして基準立体画像と比較対象立体画像との位置合わせを行った後、すなわち基準立体画像と比較対象立体画像とのレジストレーションを行った後、位置合わせされた基準立体画像と比較対象立体画像とを照合する（ステップＳ４０２）。

この例では、位置合わせされた基準立体画像（ボクセルデータ群）と比較対象立体画像（ボクセルデータ群）との相関演算を行い、この相関演算の結果得られた相関値と予め定められた照合用の閾値とを比較することによって、基準立体画像と比較対象立体画像との照合を行う。

なお、３次元の立体画像間の相関演算については特許文献１などに示されているので、ここでの説明は省略する。この特許文献１において、Ｎ＝３として３次元フーリエ変換、３次元逆フーリエ変換を施し、３次元の立体画像間の相関演算を行い、この相関演算の結果得られた相関値と予め定められた照合用の閾値とを比較し、相関値が閾値以上であれば基準立体画像と比較対象立体画像とは一致（照合ＯＫ）との照合結果を得る一方、相関値が閾値未満であれば、基準立体画像と比較対象立体画像とは不一致（照合ＮＧ）との照合結果を得る。

ＣＰＵ４Ａは、ステップＳ４０２での照合の結果がＯＫであれば（ステップＳ４０３のＹＥＳ）、照合ＯＫの表示を行う（ステップＳ４０４）。照合の結果がＮＧであれば（ステップＳ４０３のＮＯ）、照合ＮＧの表示を行うとともに（ステップＳ４０５）、位置合わせされた基準立体画像と比較対象立体画像との差分画像をディスプレイ４Ｉ上に表示する（ステップＳ４０６）。なお、相関演算の結果得られた相関値を照合結果とともに表示してもよい。

〔複数投影画像の相関値によるΔΨ分解能の向上〕
上述した例では、第１の基準投影画像Ｇ１Ａと第１の各比較候補補正投影画像Ｇ１Ｂ’との相関値を求め、その中の相関値が最も大きい比較候補補正投影画像を第１の比較確定投影画像Ｇ１ＢＲとし（ステップＳ３０５）、この第１の比較確定投影画像Ｇ１ＢＲの作成元の第１の比較候補投影画像Ｇ１ＢＲ’の投影角度ｔ’と第１の基準投影画像Ｇ１Ａの投影角度ｔとの差をＺ軸を中心とする回転ずれ量ΔΨ（ΔΨ＝ｔ’−ｔ）として求めるようにしたが（ステップＳ３０６）、投影角度の分解能が少なく、投影画像の個数が少ない場合（投影角度の分解能が低い場合）や高精度を求める場合などには、複数の相関値からその相関値が最大となる比較候補補正投影画像を推測することにより、より正確な回転ずれ量ΔΨを求めることが可能である。

例えば、図９に示すように、第１の基準投影画像Ｇ１Ａと第１の比較候補補正投影画像Ｇ１Ｂ₁’，Ｇ１Ｂ₂’，Ｇ１Ｂ₃’との相関値から関数フィッティング（例えば、２次曲線）により最大相関値を有する比較候補補正投影画像を推測し、この推測した比較候補補正投影画像の投影角度ｔ’と第１の基準投影画像Ｇ１Ａの投影角度ｔとの差をＺ軸を中心とする回転ずれ量ΔΨ（ΔΨ＝ｔ’−ｔ）として求めるようにする。

〔被写体Ｂの２次元の投影画像中に投影角度ｔ’＋９０゜の投影画像が存在しない場合〕
上述した例では、被写体Ｂの２次元の投影画像中、投影角度ｔ’＋９０゜の投影画像を第２の比較候補投影画像Ｇ２Ｂとして抽出するようにしたが（ステップＳ３１１）、正確に９０゜異なる投影画像が実在しない場合もあり得る。実際には、実在しない場合が多い。

この場合、例えば、図１０に示すように、ｔ’＋９０゜の前後の実在する投影画像ａと投影画像ｂとを抽出し、この抽出した投影画像ａと投影画像ｂとから補間して投影角度ｔ’＋９０゜の投影画像を作成し、この作成した投影角度ｔ’＋９０゜の投影画像を用いるようにするとよい。

この例において、投影画像ａ，ｂ間の１ステップの投影角度をθ_ST、投影画像ａの１ステップの投影角度θ_STに対する角度比率をα、投影画像ｂの１ステップの投影角度θ_STに対する角度比率を（１−α）とすると、投影画像ａの投影角度はｔ’＋９０゜−θ_ST・αと表され、投影画像ｂの投影角度はｔ’＋９０゜−θ_ST・（１−α）と表される（図１１参照）。この場合、投影画像ａの各画素に重み（１−α）を乗算し、投影画像ｂの各画素に重みαを乗算し、加算する。これによって、投影角度ｔ’＋９０゜の投影画像を作成し、第２の比較候補投影画像Ｇ２Ｂとする。

なお、図１２に示すように、投影画像ａを使用して第２の基準投影画像Ｇ２Ａとの回転ずれ量Δφａを求め、投影画像ｂを使用して第２の基準投影画像Ｇ２Ａとの回転ずれ量Δφｂを求め、この回転ずれ量Δφａ，Δφｂから関数フィッティングなどの補間によって、ｔ’＋９０゜の投影画像を使用した場合の回転ずれ量Δφを求めるようにしてもよい。この場合、Ｘ軸方向への平行ずれ量Δｘについても、投影画像ａと投影画像ｂを使用して平行ずれ量Δｘａ，Δｘｂを求めることにより、この平行ずれ量Δｘａ，Δｘｂから補間した値として求めることができる。

また、被写体Ｂの２次元の投影画像中にｔ’＋９０゜の投影画像が存在しない場合、単純に、その投影角度が最も近い実在する投影画像を第２の比較候補投影画像Ｇ２Ｂとして抽出するようにしてもよい。

以上の説明から分かるように、本実施の形態によれば、再構成処理後の３次元の立体画像（ボクセルデータ群）ではなく、再構成処理前の２次元の投影画像を使用して、基準立体画像と比較対象立体画像との間の回転ずれ量（Δθ，Δφ，ΔΨ）と平行ずれ量（Δｘ，Δｙ，Δｚ）が求められる。このため、３次元画像処理を必要とせず、局所領域や目印を定める必要もなく、高速処理が実現される。また、再構成処理で発生するアーチファクトの影響も受けない。これにより、短時間でかつ正確に、レジストレーションを行うことができるようになる。また、照合スピードもアップし、照合結果の精度も高まる。

図１３に３次元ずれ量計測装置４の要部の機能ブロック図を示す。３次元ずれ量計測装置４は、基準立体画像を構成する２次元の投影画像を記憶するデータベースＤＢ１と、比較対象立体画像を構成する２次元の投影画像を記憶するデータベースＤＢ２と、第１の基準投影画像抽出部４−１と、第１の比較候補投影画像抽出部４−２と、第１の回転角度ずれ量算出部４−３と、第１の回転補正部４−４と、第１の比較確定投影画像決定部４−５と、回転角度ずれΔΨ算出部４−６と、回転角度ずれΔθ算出部４−７と、平行ずれ量Δｙ，Δｚ算出部４−８と、第２の基準投影画像抽出部４−９と、第２の比較候補投影画像抽出部４−１０と、第２の回転角度ずれ量算出部４−１１と、第２の回転補正部４−１２と、平行ずれ量Δｘ算出部４−１３とを備えている。

第１の基準投影画像抽出部４−１は、データベースＤＢ１から、所定の投影角度ｔの投影画像を第１の基準投影画像Ｇ１Ａとして抽出する。第１の比較候補投影画像抽出部４−２は、データベースＤＢ２から、撮影角度ｔ付近の複数の投影画像を第１の比較候補投影画像Ｇ１Ｂとして抽出する。

第１の回転角度ずれ量算出部４−３は、第１の基準投影画像抽出部４−１からの第１の基準投影画像Ｇ１Ａおよび第１の比較候補投影画像抽出部４−２からの第１の各比較候補投影画像Ｇ１Ｂを入力とし、第１の基準投影画像Ｇ１Ａと第１の各比較候補投影画像Ｇ１Ｂとの間の回転角度のずれを求める。

第１の回転補正部４−４は、第１の比較候補投影画像抽出部４−２からの第１の各比較候補投影画像Ｇ１Ｂ入力とし、第１の回転角度ずれ量算出部４−３で求められた各回転角度のずれ量分、第１の各比較候補投影画像Ｇ１Ｂを回転補正して、第１の比較候補補正投影画像Ｇ１Ｂ’を得る。

第１の比較確定投影画像決定部４−５は、第１の基準投影画像抽出部４−１からの第１の基準投影画像Ｇ１Ａと第１の回転補正部４−４からの第１の各比較候補補正投影画像Ｇ１Ｂ１’とを入力とし、第１の基準投影画像Ｇ１Ａと第１の各比較候補補正投影画像Ｇ１Ｂ１’との相関値を求め、その中の相関値が最も大きい比較候補補正投影画像を第１の比較確定投影画像Ｇ１ＢＲとする。

回転角度ずれΔΨ算出部４−６は、第１の比較確定投影画像決定部４−５からの第１の比較確定投影画像Ｇ１ＢＲを入力とし、第１の比較確定投影画像Ｇ１ＢＲの作成元の第１の比較候補投影画像Ｇ１ＢＲ’の投影角度ｔ’を第１の比較候補投影画像抽出部４−２を介して読み出し、第１の基準投影画像Ｇ１Ａの投影角度ｔを第１の基準投影画像抽出部４−１を介して読み出し、第１の比較候補投影画像Ｇ１ＢＲ’の投影角度ｔ’と第１の基準投影画像Ｇ１Ａの投影角度ｔとの差をＺ軸を中心とする回転ずれ量ΔΨ（ΔΨ＝ｔ’−ｔ）として求める。

回転角度ずれΔθ算出部４−７は、第１の比較確定投影画像決定部４−５からの第１の比較確定投影画像Ｇ１ＢＲを入力とし、第１の比較確定投影画像Ｇ１ＢＲの作成元の第１の比較候補投影画像Ｇ１ＢＲ’を第１の比較候補投影画像抽出部４−２より取得し、第１の基準投影画像抽出部４−１からの第１の基準投影画像Ｇ１Ａと第１の比較候補投影画像Ｇ１ＢＲ’との間の回転角度のずれをＸ軸を中心とする回転ずれ量Δθとして求める。

平行ずれ量Δｙ，Δｚ算出部４−８は、第１の比較確定投影画像決定部４−５からの第１の比較確定投影画像Ｇ１ＢＲを入力とし、第１の基準投影画像抽出部４−１からの第１の基準投影画像Ｇ１Ａと第１の比較確定投影画像Ｇ１ＢＲとの間の水平方向のずれ量をＹ軸方向への平行ずれ量Δｙ、垂直方向のずれ量をＺ軸方向への平行ずれ量Δｚとして求める。

第２の基準投影画像抽出部４−９は、データベースＤＢ１から、第１の基準投影画像Ｇ１Ａに対して投影角度が９０゜ずれた投影画像（ｔ＋９０゜の投影画像）を第２の基準投影画像Ｇ２Ａとして抽出する。

第２の比較候補投影画像抽出部４−１０は、第１の比較確定投影画像決定部４−５からの第１の比較確定投影画像Ｇ１ＢＲを入力とし、データベースＤＢ２から、第１の比較確定投影画像Ｇ１ＢＲの作成元の第１の比較候補投影画像Ｇ１ＢＲ’に対してその投影角度が第２の基準投影画像Ｇ２Ａと同方向に９０゜ずれた投影画像（ｔ’＋９０゜の投影画像）を第２の比較候補投影画像Ｇ２Ｂとして抽出する。

第２の回転角度ずれ算出部４−１１は、第２の基準投影画像抽出部４−９からの第２の基準投影画像Ｇ２Ａおよび第２の比較候補投影画像抽出部４−１０からの第２の比較候補投影画像Ｇ２Ｂを入力とし、第２の基準投影画像Ｇ２Ａと第２の比較候補投影画像Ｇ２Ｂとの間の回転角度のずれをＹ軸を中心とする回転ずれ量Δφとして求める。

第２の回転補正部４−１２は、第２の比較候補投影画像抽出部４−１０からの第２の比較候補投影画像Ｇ２Ｂを入力とし、この第２の比較候補投影画像Ｇ２Ｂを第２の回転角度ずれ算出部４−１１で求められた回転角度のずれΔφ分回転補正して、第２の比較確定投影画像Ｇ２Ｒを得る。

平行ずれ量Δｘ算出部４−１３は、第２の回転補正部４−１２からの第２の比較確定投影画像Ｇ２ＢＲを入力とし、第２の基準投影画像抽出部４−９からの第２の基準投影画像Ｇ２Ａと第２の比較確定投影画像Ｇ２ＢＲとの間の水平方向のずれ量をＸ軸方向への平行ずれ量Δｘとして求める。

なお、この３次元ずれ量計測装置４におけるＸ，Ｙ，Ｚの３軸の回転ずれ量（Δθ，Δφ，ΔΨ）および平行ずれ量（Δｘ，Δｙ，Δｚ）の計測処理は、ｎ個の２次元の投影画像を再構成処理することによって比較対象立体画像を作成した後に行うようにしてもよいが、比較対象立体画像を再構成処理によって作成する間に行うようにすれば、比較対象立体画像の作成後に即座に基準立体画像との位置合わせを行って、照合結果を得ることができるようになる。

また、上述した実施の形態では、最初に撮影した被写体Ａの３次元の立体画像を基準立体画像とし、次に撮影した被写体Ｂの３次元の立体画像を比較対象立体画像としたが、最初に撮影した被写体Ａの３次元の立体画像を比較対象立体画像とし、次に撮影した被写体Ｂの３次元の立体画像を基準立体画像としてもよい。

また、基準立体画像や比較対象立体画像はボクセルデータ群として上位装置から与えられるものとしてもよく、この上位装置から与えられた基準立体画像や比較対象立体画像からその立体画像を構成する投影画像を作成し、撮影によって得た投影画像と同様にして、Ｘ，Ｙ，Ｚの３軸の回転ずれ量（Δθ，Δφ，ΔΨ）および平行ずれ量（Δｘ，Δｙ，Δｚ）を計測するようにしてもよい。

また、上述した実施の形態では、被写体を人体としたが、被写体は人体に限られるものではない。例えば、製薬などの分野において、マウスを使用した実験などにも利用することができる。また、産業分野において、金属などの非破壊検査に利用することもできる。このように、本発明は、医療や歯科などの分野に限らず、各種の分野において有効利用することが可能である。

また、上述した３次元ずれ量計測処理装置４において、図１３に機能ブロックとして示した各部の処理動作は、プログラム（３次元ずれ量計測プログラム）に従うＣＰＵの処理動作として実現されるが、同様の処理をハードウェアの回路構成で実現するようにしてもよい。

本発明に係る３次元ずれ量計測方法の実施に用いるＸ線コンピュータ断層撮影システムの一例を示す図である。 このシステムの３次元ずれ量計測装置におけるＣＰＵの処理動作（基準立体画像の取得処理）を示すフローチャートである。 このシステムの３次元ずれ量計測装置におけるＣＰＵの処理動作（比較対象立体画像の取得処理）を示すフローチャートである。 このシステムの３次元ずれ量計測装置におけるＣＰＵの処理動作（基準立体画像と比較対象立体画像との間のずれ量の計測処理）を示すフローチャートである。 図４に続くフローチャートである。 このシステムの３次元ずれ量計測装置におけるＣＰＵの処理動作（位置合わせ処理＋照合処理）を示すフローチャートである。 基準立体画像と比較対象立体画像との間のずれ量の計測処理によってＺ軸を中心とする回転ずれ量ΔΨ、Ｘ軸を中心とする回転ずれ量Δθ、Ｙ軸方向への平行ずれ量ΔｙおよびＺ軸方向への平行ずれ量Δｚが求められる様子を説明する図である。 基準立体画像と比較対象立体画像との間のずれ量の計測処理によってＹ軸を中心とする回転ずれ量ΔφおよびＸ軸方向への平行ずれ量ΔＸが求められる様子を説明する図である。 複数投影画像の相関値によるΔΨ分解の向上を説明する図である。 投影画像ａ，ｂからの補間による第２の比較候補投影画像の生成について説明する図である。 投影画像ａ，ｂからの補間による第２の比較候補投影画像の生成に際して用いる重みについて説明する図である。 投影画像ａ，ｂの回転ずれ量Δφａ，φｂからの補間による回転ずれ量Δφの算出について説明する図である。 このシステムにおける３次元ずれ量計測装置の要部の機能ブロック図である。 Ｘ線コンピュータ断層撮影法による３次元の立体画像の取得例を説明する図である。

符号の説明

１…Ｘ線源、２…２次元撮像装置、３（Ａ，Ｂ）…被写体、４Ａ…ＣＰＵ、４Ｂ…ＲＡＭ、４Ｃ…ＲＯＭ、４Ｄ…記憶装置、４Ｅ〜Ｈ…インタフェース、４Ｉ…ディスプレイ、４Ｊ…キーボード、４Ｋ…マウス、ＤＢ１，ＤＢ２…データベース、４−１…第１の基準投影画像抽出部、４−２…第１の比較候補投影画像抽出部、４−３…第１の回転角度ずれ量算出部、４−４…第１の回転補正部、４−５…第１の比較確定投影画像決定部、４−６…回転角度ずれΔΨ算出部、４−７…回転角度ずれΔθ算出部、４−８…平行ずれ量Δｙ，Δｚ算出部、４−９…第２の基準投影画像抽出部、４−１０…第２の比較候補投影画像抽出部、４−１１…第２の回転角度ずれ量算出部、４−１２…第２の回転補正部、４−１３…平行ずれ量Δｘ算出部。

Claims (6)

1. 基準となる３次元の立体画像を基準立体画像、この基準立体画像に対して比較対象とする３次元の立体画像を比較対象立体画像とし、前記基準立体画像と前記比較対象立体画像との間のずれ量を計測する３次元ずれ量計測方法において、
   前記基準立体画像を構成する２次元の投影画像と前記比較対象立体画像を構成する２次元の投影画像とから前記基準立体画像と前記比較対象立体画像との間の互いに直交するＸ，Ｙ，Ｚの３軸の回転ずれ量および平行ずれ量を計測するずれ量計測ステップ
   を備えることを特徴とする３次元ずれ量計測方法。
2. 請求項１に記載された３次元ずれ量計測方法において、
   前記ずれ量計測ステップは、
   前記基準立体画像を構成する２次元の投影画像中の所定の投影角度から見た投影画像を第１の基準投影画像として抽出する第１ステップと、
   前記比較対象立体画像を構成する２次元の投影画像中の前記第１の基準投影画像と同じ投影角度付近の複数の投影画像を第１の比較候補投影画像として抽出する第２ステップと、
   前記第１の基準投影画像と前記第１の各比較候補投影画像との間の回転角度のずれを求める第３ステップと、
   前記第１の各比較候補投影画像を前記第３ステップで求めた回転角度のずれ分回転補正して第１の比較候補補正投影画像とする第４ステップと、
   前記第１の基準投影画像と前記第１の各比較候補補正投影画像との相関値を求め、その中の相関値が最も大きい比較候補補正投影画像を第１の比較確定投影画像とする第５ステップと、
   前記第１の比較確定投影画像の作成元の前記第１の比較候補投影画像の投影角度と前記第１の基準投影画像の投影角度との差を前記Ｚ軸を中心とする回転ずれ量ΔΨとして求める第６ステップと、
   前記第１の比較確定投影画像の作成元の前記第１の比較候補投影画像と前記第１の基準投影画像との間の回転角度のずれを前記Ｘ軸を中心とする回転ずれ量Δθとして求める第７ステップと、
   前記第１の比較確定投影画像の前記第１の基準投影画像との間の水平方向のずれ量を前記Ｙ軸方向への平行ずれ量Δｙ、垂直方向のずれ量を前記Ｚ軸方向への平行ずれ量Δｚとして求める第８ステップと、
   前記基準立体画像を構成する２次元の投影画像中の前記第１の基準投影画像に対して投影角度が９０゜ずれた投影画像を第２の基準投影画像として抽出する第９ステップと、
   前記比較対象立体画像を構成する２次元の投影画像中の前記第１の比較確定投影画像の作成元の前記第１の比較候補投影画像に対してその投影角度が前記第２の基準投影画像と同方向に９０゜ずれた投影画像を第２の比較候補投影画像として抽出する第１０ステップと、
   前記第２の基準投影画像と前記第２の比較候補投影画像との間の回転角度のずれを前記Ｙ軸を中心とする回転ずれ量Δφとして求める第１１ステップと、
   前記第２の比較候補投影画像を前記回転ずれ量Δφ分回転補正して第２の比較確定投影画像とする第１２ステップと、
   前記第２の基準投影画像と前記第２の比較確定投影画像との間の水平方向のずれ量を前記Ｘ軸方向への平行ずれ量Δｘとして求める第１３ステップと
   を備えることを特徴とする３次元ずれ量計測方法。
3. 請求項１又は２に記載された３次元ずれ量計測方法において、
   前記ずれ量計測ステップで求められた回転ずれ量および平行ずれ量に基づいて前記基準立体画像と前記比較対象立体画像との間の位置合わせを行う位置合わせ処理ステップ
   を備えることを特徴とする３次元ずれ量計測方法。
4. 請求項３に記載された３次元ずれ量計測方法において、
   前記位置合わせ処理ステップによって位置合わせされた前記基準立体画像と前記比較対象立体画像との相関値に基づいて前記基準立体画像と前記比較対象立体画像とを照合する照合ステップ
   を備えることを特徴とする３次元ずれ量計測方法。
5. 請求項２に記載された３次元ずれ量計測方法において、
   前記第６ステップに代えて、
   前記第１の基準投影画像と前記第１の各比較候補補正投影画像との相関値から最大相関値を有する比較候補補正投影画像を推測し、この推測した比較候補補正投影画像の投影角度と前記第１の基準投影画像の投影角度との差を前記Ｚ軸を中心とする回転ずれ量ΔΨとして求めるステップ
   を備えることを特徴とする３次元ずれ量計測方法。
6. 請求項２に記載された３次元ずれ量計測方法において、
   前記第１０ステップに代えて、
   前記比較対象立体画像を構成する２次元の投影画像中に前記第１の比較確定投影画像の作成元の前記第１の比較候補投影画像に対してその投影角度が前記第２の基準投影画像と同方向に９０゜ずれた投影画像が実在しない場合、その投影角度が近い実在する複数の投影画像から補間して第２の比較候補投影画像を作成するステップ
   を備えることを特徴とする３次元ずれ量計測方法。

Images