(57)【要約】選択された投影を逆投影する（１００）ことにより投影（ｑ1…ｑp）からピクセル画像ｆが形成され、それにより中間画像（ｌ1，ｍ）が生成されるとともに、選択された中間画像に関してデジタル画像座標変換（１０２）及び／又はリサンプリング（図３１、１８６、１９２、１９６）が行われる。デジタル画像座標変換（１０２）は、中間画像の成分投影の視野角およびそれらのフーリエ特性を明らかにするように選択され、それにより、中間画像をまばらなサンプルによって正確に表すことができる。その結果として得られる中間画像が部分集合に集められ（１０４）、このプロセスは、十分な投影および中間画像が処理されて集められることによりピクセル画像ｆが形成されるまで再帰的形式で繰り返される。デジタル画像座標変換としては、回転（図１８、１０２）、剪断（図１０Ｂ、１２０、１２２）、伸張、縮小（１０９）等を挙げることができる。リサンプリングとしては、アップサンプリング（１０１、１０６）やダウンサンプリング（１０９）などを挙げることができる。デジタル画像座標変換（２０２）及び／又はリサンプリング（２０４）及び／又はデシメーション（図３２、２０４；図３３、２１２）、最終的な中間画像の再投影（２０８）を行うことにより、１つのピクセル画像（ｆ）から投影（図３２、ｐθ1…ｐθ18）を形成することができる。
2007504912.tif 000002
【特許請求の範囲】
【請求項１】
投影（ｑ1…ｑp）から１つのピクセル画像
【数１】
2007504912.tif 000205
を形成する方法（図５）であって、
（ａ）選択された投影（ｑ1…ｑp）から中間画像（Ｉ1，Ｐ）を生成するステップ（１００）と、
（ｂ）選択された中間画像（Ｉ1，Ｐ）に関してデジタル画像座標変換を行うステップ（１０２）であって、前記座標変換のパラメータが、前記中間画像を生成した投影の視野角および前記中間画像のフーリエ特性を明らかにするように選択されるステップと、
（ｃ）ステップ（ｂ）において生成された変換された中間画像の一部を集めて、集合中間画像（Ｉ2，Ｐ／２）を生成するステップ（１０４）と、
（ｄ）全ての投影および中間画像が処理されて集められることによりピクセル画像
【数２】
2007504912.tif 000206
が形成されるまで、ステップ（ｂ）およびステップ（ｃ）を再帰的な形式で繰り返すステップと
を含み、
前記座標変換パラメータが、前記中間画像の集合（１０４）をまばらなサンプルにより所望の精度で成し得るように選択される、前記方法。
【請求項２】
前記集合（１０４，１０８）が、デジタル画像を加えることによって行われる、請求項１に記載の方法。
【請求項３】
ステップ（ａ）において、選択された投影（ｑ1…ｑp）を逆投影することにより少なくとも幾つかの中間画像（Ｉn，ｍ）を生成する、請求項１に記載の方法。
【請求項４】
ステップ（ａ）において、選択された２つ以上の投影（ｑ1…ｑp）を逆投影することにより少なくとも幾つかの中間画像（Ｉn，ｍ）をそれぞれ形成する、請求項１に記載の方法。
【請求項５】
ステップ（ｄ）において、選択された３つ以上の変換された中間画像を集めることにより少なくとも幾つかの集合中間画像（Ｉn，ｍ）をそれぞれ形成する、請求項１に記載の方法。
【請求項６】
前記デジタル画像座標変換がデジタルフィルタリングを使用して行われる、請求項１に記載の方法。
【請求項７】
選択された座標変換がデジタル画像回転を含む、請求項１に記載の方法。
【請求項８】
選択された座標変換が、デジタル画像剪断（図１０Ｂ、１２０、１２２）または剪断スケーリングを含む、請求項１に記載の方法。
【請求項９】
選択された座標変換が、前記デジタル画像のアップサンプリング（１０１，１０６）及び／又はダウンサンプリング（１０９）を含む、請求項１に記載の方法（図１５）。
【請求項１０】
前記デジタル画像剪断が一次元線形デジタルフィルタによって行われる、請求項８に記載の方法。
【請求項１１】
前記デジタル画像アップサンプリング及び／又はダウンサプリングが一次元線形デジタルフィルタによって行われる、請求項９に記載の方法。
【請求項１２】
前記デジタルフィルタの少なくとも一部がシフト不変フィルタである、請求項１０に記載の方法。
【請求項１３】
前記デジタルフィルタの少なくとも一部がシフト不変フィルタである、請求項１１に記載の方法。
【請求項１４】
前記デジタルフィルタの少なくとも一部が再帰フィルタである、請求項１２又は１３に記載の方法。
【請求項１５】
前記デジタルフィルタの少なくとも一部が高速フーリエ変換（ＦＦＴ）を使用して実施される、請求項１２又は１３に記載の方法。
【請求項１６】
選択された中間画像及び／又は変換された中間画像及び／又は集合中間画像に対して選択されたオーバーサンプリングが適用される、請求項１に記載の方法。
【請求項１７】
非デカルトサンプリングパターンが使用される、請求項１に記載の方法。
【請求項１８】
選択された座標変換が、階層の１つのレベル内または隣り合うレベルにわたって組み合わせることができる、請求項１に記載の方法。
【請求項１９】
一群の線、曲線または面に沿う投影（ｑ1…ｑp）から１つのピクセル画像
【数３】
2007504912.tif 000207
を形成する方法（図３１）であって、
（ａ）中間画像（Ｉl，ｍ）を生成するステップ（１８４）と、
（ｂ）選択された中間画像に関してデジタル画像リサンプリングを行うステップ（１８６）であって、サンプルの位置が、選択された投影の視野角および前記中間画像のフーリエ特性を明らかにするように選択されるステップと、
（ｃ）リサンプリングされた中間画像の選択された一部を集めて、集合中間画像（Ｉz，ｍ）を生成するステップ（１９０）と、
（ｄ）全ての投影および中間画像が処理されて集められることによりピクセル画像が形成されるまで、中間画像のサンプルの密度を高める反復の各レベルで、ステップ（ｂ）およびステップ（ｃ）を再帰的な形式で繰り返すステップと
を含み、
サンプリング方式が、リサンプリングされた中間画像の集合をまばらなサンプルにより所望の精度で成し得るように選択される、前記方法。
【請求項２０】
ステップ（ａ）において、選択された投影の加重逆投影（図１９、１８０、１８２）により少なくとも幾つかの中間画像を形成する、請求項１９に記載の方法。
【請求項２１】
ステップ（ａ）において、選択された２つ以上の投影の加重逆投影（図１９、１８０、１８２）により少なくとも幾つかの中間画像をそれぞれ形成する、請求項１９に記載の方法。
【請求項２２】
ステップ（ｄ）において、選択された３つ以上の変換された中間画像を集めることにより少なくとも幾つかの集合中間画像をそれぞれ形成する、請求項１９に記載の方法。
【請求項２３】
前記中間画像が、一群の線、曲線または面上に位置するサンプルを有している、請求項１９に記載の方法。
【請求項２４】
前記デジタル画像リサンプリングが、一群の線、曲線または面の交差に伴う中間サンプリング方式を使用することによって、一連の低次元デジタルフィルタリング演算により行われる、請求項１９に記載の方法。
【請求項２５】
オーバーサンプリングの選択された度合いが、選択されリサンプリングされた中間画像および集合中間画像に対して適用される、請求項１９に記載の方法。
【請求項２６】
前記集合が、デジタル画像を加えることにより行われる、請求項１９に記載の方法。
【請求項２７】
前記リサンプリングおよび集合が、連続するレベルにわたって組み合わせることができる、請求項１９に記載の方法。
【請求項２８】
リサンプリングステップ（ｂ）が含められる前および後に、少なくとも１つの中間画像が重み付けられる、請求項３０に記載の方法。
【請求項２９】
サンプリング密度の変更がデジタルフィルタリングによって行われる、請求項１９に記載の方法。
【請求項３０】
投影（ｑ1…ｑp）から１つのピクセル画像
【数４】
2007504912.tif 000208
を形成する方法（図８）であって、
（ａ）複数の中間画像（Ｉl，ｍ）を生成するステップ（９９）であって、少なくとも１つの前記中間画像が非デカルト及び／又は周期的でないサンプリングパターンに対応しているステップと、
（ｂ）選択された中間画像に関してデジタル画像アップサンプリングまたはダウンサンプリングを行うステップ（１０６）と、
（ｃ）アップサンプリングされた／ダウンサンプリングされた中間画像に関してデジタル画像座標変換を行うステップと、
（ｄ）ステップ（ｃ）において生成された変換された中間画像の一部を集めて、集合中間画像を生成するステップ（１１０）と、
（ｅ）全ての投影および中間画像が処理されて集められることによりピクセル画像が形成されるまで、ステップ（ｂ）、（ｃ）及び（ｄ）を再帰的な形式で繰り返すステップと
を含み、
前記デジタル画像座標変換のうちの少なくとも１つが、非デカルト及び／又は周期的でないサンプリングパターンで行われ、前記座標変換パラメータが、前記中間画像の集合をまばらなサンプルにより所望の精度で成し得るように選択される、前記方法。
【請求項３１】
前記集合が、デジタル画像を加えることによって行われる、請求項３０に記載の方法。
【請求項３２】
ステップ（ａ）において、選択された投影を逆投影することにより少なくとも幾つかの中間画像を生成する、請求項３０に記載の方法。
【請求項３３】
ステップ（ａ）において、選択された２つ以上の投影を逆投影することにより少なくとも幾つかの中間画像をそれぞれ形成する、請求項３０に記載の方法。
【請求項３４】
ステップ（ｃ）において、選択された３つ以上の変換された中間画像を集めることにより少なくとも幾つかの集合中間画像をそれぞれ形成する、請求項３０に記載の方法。
【請求項３５】
前記デジタル画像座標変換がデジタルフィルタリングを使用して行われる、請求項３０に記載の方法。
【請求項３６】
選択された座標変換がデジタル画像回転を含む、請求項３０に記載の方法。
【請求項３７】
前記デジタル画像リサンプリングアップサンプリング及び／又はダウンサプリングが一次元線形デジタルフィルタによって行われる、請求項３０に記載の方法。
【請求項３８】
前記デジタルフィルタの少なくとも一部がシフト不変フィルタである、請求項３７に記載の方法。
【請求項３９】
前記デジタルフィルタの少なくとも一部が再帰フィルタである、請求項３７に記載の方法。
【請求項４０】
前記デジタルフィルタの少なくとも一部が高速フーリエ変換（ＦＦＴ）を使用して実施される、請求項３８に記載の方法。
【請求項４１】
選択された中間画像及び／又は変換された中間画像及び／又は集合中間画像に対して選択されたオーバーサンプリングが適用される、請求項３０に記載の方法。
【発明の詳細な説明】
【技術分野】
【０００１】
この出願は、参照することによりその全体が本願に組み込まれる２００３年９月９日に提出された仮出願第６０／５０１，３５０号の一部継続出願である。
【０００２】
この発明は、断層撮影法（トモグラフィー）に関し、特に、投影からピクセル画像を形成するための方法および装置に関する。
【背景技術】
【０００３】
［発明の背景］
断層再構成は、Ｘ線コンピュータ断層撮影（ＣＴ）、ポジション放射断層撮影（ＰＥＴ）、単光子放出コンピュータ断層撮影（ＳＰＥＣＴ）、磁気共鳴断層撮影（ＭＲＩ）用の特定の取得方法を含むほぼ全ての画像診断様相の基礎にある良く知られた技術である。また、断層再構成は、製造における非破壊評価（ＮＤＥ）用、セキュリティ走査用、合成開口レーダー（ＳＡＲ）、電波天文学、地球物理学および他の分野においても用途を見出している。
【０００４】
断層データには幾つかの主な形式が存在する。すなわち、（ｉ）平行線の組（セット）に沿って線積分が行われる平行ビーム；（ii）ファン（扇）またはコーン（円錐）として発散する線の組に沿って線積分が行われる発散ビーム；（iii）円、楕円または他の閉じられた或いは開かれた曲線等の曲線の組に沿って積分が行われる曲線ビームがそれである。断層再構成の１つの問題点は、２Ｄ画像または３Ｄ画像をその線積分投影の組から再構成するという点である。断層再構成の他の問題点は、３Ｄ画像をその面積分投影の組すなわち一群の面上でのその積分から再構成するという点である。例えば、３Ｄラドン変換は、原点からの距離および方向が様々な一群の２Ｄ平面上で画像を積分することを含んでいる。断層再構成の問題の幾つか、および、再構成方法の幾つかは、例えばF．Natterer著「The  Mathematics  of  Computerized  Tomography」（チチェスター：John  Wiley，1986）；F．Netterer及びF．Wubbeling著「Mathematical  Methods  in  Image  Reconstruction」（フィラデルフィア、Society  for  Industrial  and  Applied  Mathematics，2001）；A.C．Kak及びM．Soaney著「Principles  of  Computerized  Tomographic  Imaging」（ニューヨーク：ＩＥＥＥプレス、1998）；S.R．Deans著「The  Radon  Transform  and  Some  of  Its  Applications」（ニューヨーク、Wiley，1983）などの標準的な文献に記載されている。
【０００５】
断層再構成における選択方法は、加重されない（平行ビームまたはラドン変換の場合）或いは加重される（他の殆どの場合）逆投影ステップを使用するフィルタ逆投影（ＦＢＰ）または重畳逆投影（ＣＢＰ）である。このステップは、２ＤのＮ×Ｎピクセル画像に関してはＮ3、３ＤのＮ×Ｎ×Ｎボクセル画像に関しては少なくともＮ4のような計算要件スケーリングを伴い、技術的には計算上の障害となる。したがって、画像の解像度をＮから２Ｎへと２倍にすると、計算が約８倍（あるいは、３Ｄの場合には１６倍）増える。かつてよりも多量のデータをリアルタイムで収集できる新たな技術（例えば、多列検出器を用いた心臓イメージング、インターベンショナルイメージング）の出現および３Ｄ取得ジオメトリの急増に伴ってコンピュータがかなり高速になってきている一方で、高速再構成の必要性は高まりを見せている。高速再構成は、画像形成プロセスをスピードアップすることができ、専用の画像再構成コンピュータのコストを低減することができ、あるいは、この両方を達成できる。
【０００６】
逆投影の双対演算は、電子的に記憶された画像の投影を計算するプロセスである再投影である。このプロセスも、断層再構成においては基本的な役割を果たす。人間のオブジェクト（対象物）の実用的な３Ｄイメージングにとって重要な階層コーンビームジオメトリ（階層コーンビーム構成）におけるロングオブジェクト問題のための高速再構成アルゴリズムを構築するべく、逆投影と再投影とを組み合わせて使用することもできる。また、様々な用途では、単一の画像を再構成（復元）するために逆投影ステップおよび再投影ステップの両方が何回も行われる反復再構成アルゴリズムを使用することが有益であり或いは必要でさえある。逆投影ステップおよび再投影ステップのスピードアップは、そのような反復方法の経済的実現可能性を決定する。
【０００７】
再構成をスピードアップするため、何年にもわたって幾つかの方法が提案されてきた。例えば、Brandtらによる米国特許第５，７７８，０３８号は、マルチレベル分解を使用して各段階で視野全体をカバーする解像度の高い画像を生成する２Ｄ平行ビーム断層撮影のための方法について記載している。Nillsonらによる米国特許第６，１５１，３７７号は他の階層型逆投影方法を開示している。これらのシステムはメリットがあるが、更に正確な画像を生成し且つ精度とスピードとの間に多くの自由度を与える方法および装置が必要である。
【発明の開示】
【発明が解決しようとする課題】
【０００８】
したがって、この発明の１つの目的は、コンピュータ断層撮影（ＣＴ）スキャンのための新規な改良された方法および装置を提供することである。
【０００９】
他の目的は、更に正確な画像を生成し且つ精度とスピードとの間に多くの自由度を与えるＣＴスキャンのための方法および装置を提供することである。
【課題を解決するための手段】
【００１０】
［発明の要約］
これらの目的は、本発明によって達成または突破される。選択された投影を逆投影して、中間画像を生成し、中間画像に関してデジタル画像座標変換及び／又はリサンプリングを行うことにより、投影からピクセル画像が形成される。デジタル画像座標変換は、中間画像の成分投影の視野角およびそれらのフーリエ特性を明らかにするように選択され、それにより、中間画像をまばらなサンプルによって正確に表すことができる。その結果として得られる中間画像は部分集合(subset)に集められ、このプロセスは、十分な投影および中間画像が処理されて集められることによりピクセル画像が形成されるまで再帰的形式で繰り返される。
【００１１】
デジタル画像座標変換としては、回転、剪断、伸張、縮小等を挙げることができる。リサンプリングとしては、アップサンプリンやダウンサンプリングなどを挙げることができる。
【００１２】
デジタル画像座標変換及び／又はリサンプリング及び／又はデシメーションおよび最終的な中間画像の再投影を行うことにより、１つのピクセル画像から投影を形成することができる。
【００１３】
添付図面と併せて本発明の実施態様の以下の説明を参照することにより、この発明の前述した特徴および他の特徴並びにこれらを得る形式が更に明らかになるとともに、本発明自体を最も良く理解できる。
【発明を実施するための最良の形態】
【００１４】
［詳細な説明］
［記号およびフォント］
明瞭さを高めるために、数学的な記号およびフォントの以下の方式が使用される。
【００１５】
空間域における関数は、小文字（例えばｆ（ｘ））によって示され、一方、それらのフーリエ変換は、大文字（Ｆ（ω））によって示される。
【００１６】
２つの可変関数の添字は、利便性に応じて様々に示される。関数ｆの以下の３つの表記、すなわち、
【数１】
2007504912.tif 000003
は等価なものである。
【００１７】
連続域関数および離散域関数はそれぞれ、それらの添字を用いて使用される括弧書きの形式によって区別される。すなわち、ｆ（ｘ1，ｘ2）は、２つの連続変数（すなわち、
【数２】
2007504912.tif 000004
）の関数であり、また、ｆ［ｍ1，ｍ2］は
【数３】
2007504912.tif 000005
のサンプル値バージョンであり、したがって、２Ｄ配列である。
【００１８】
一次演算子およびその対応する行列は、フォント形式によって区別される。
【数４】
2007504912.tif 000006
が行列であるとすると、
【数５】
2007504912.tif 000007
はその関連する一次演算子である。例えば、Ａが座標変換である場合には、
【数６】
2007504912.tif 000008
である。
【００１９】
同じ演算子は、時として、ブロック図の内側と外側とでは異なって示される。演算子は、ブロック図の外側では
【数７】
2007504912.tif 000009
として示される場合があるが、ブロック図の内側では
【数８】
2007504912.tif 000010
として示される。
【００２０】
［ハードウェアの概要］
本発明は、ＣＴスキャナを含む様々な画像装置において用途を有している。一般的な画像装置１０（図１１）は、頭部などの対象物からデータを取得し且つ例えば発散ビームの幾何学的構成（ジオメトリ）を成す線積分投影に対応する未加工データ１４を投影プリプロセッサ１６に対して送るスキャナ１２を有している。投影プリプロセッサ１６は、シフト変化する場合があるデータに対して様々な変換、正規化、訂正を適用するとともに、重み付けやフィルタリングも行う。投影プリプロセッサ１６の出力は、前処理された投影の集合体であり、以下ではこれを単に投影と称し、また、シノグラム1とも呼ぶ。このシノグラム1  １８はシノグラム更新プロセッサ２０に対して供給される。シノグラム更新プロセッサ２０は、シノグラム2  ３４からの情報を使用して入力シノグラム1  １８を修正し、これにより、例えば、ビーム硬化を含む様々なアーチファクトを訂正し、または、これを多段階あるいは反復の再構成（復元）処理の一環として行う。
【００２１】
シノグラム更新プロセッサ２０の出力は、高速逆投影プロセッサ２４に対して入力されるシノグラム3  ２２である。高速逆投影プロセッサ２４は、一般に、ここで説明するアルゴリズムを実行するようにプログラムされ及び／又は配線された任意の適したタイプのコンピュータまたは特定目的のハードウェア或いはこれらの組み合わせである。
【００２２】
高速逆投影プロセッサ２４の出力は、画像調整プロセッサ２８に対して入力される電子画像1  ２６である。画像調整プロセッサ２８は、アーチファクト画像または多段階あるいは反復の再構成プロセスにおける更なる処理のための画像についての識別および抽出を含む、電子画像の必要な後処理を行う。
【００２３】
必要に応じて、画像調整プロセッサ２８は、高速再投影プロセッサ３２に対して入力される電子画像2  ３０を生成することができる。高速再投影プロセッサ３２は、一般に、ここで説明するアルゴリズムを実行するようにプログラムされ及び／又は配線された任意の適したタイプのコンピュータまたは特定目的のハードウェア或いはこれらの組み合わせである。必要に応じて、このプロセッサは、逆投影プロセッサ２４によって用いられる同じコンピュータおよびハードウェアの使用を共有することができる。
【００２４】
高速再投影プロセッサ３２の出力は、シノグラム更新プロセッサ２０にフィードバックされるシノグラム2  ３４である。逆投影／再投影プロセスは、適切な結果が得られるまで続けることができる。再投影が常に必要であるとは限らないが、多くの状況において再投影が役に立つ。
【００２５】
電子画像1  ２６が適切な場合、画像調整プロセッサ２８は、記憶／解析／表示装置３８に対して供給される電子画像3  ３６を生成する。電子画像3  ３６をコンピュータメモリ内に記憶し、及び／又は例えば異常や危険な材料に関して電子的に解析し、及び／又は表示し、及び／又は特定の視認可能な形式で印刷することができると考えられる。
【００２６】
［本発明の逆投影方法および再投影方法の概要］
本発明の逆投影方法は、様々な技術を使用することにより、ピクセル（画素）及び／又はボクセル（３Ｄ画素）から成る画像を形成する。ピクセル及び／又はボクセルは、以下では、総称してピクセルと称し、ここでは一般的な方法で導入される。
【００２７】
この説明は、例えばD．Dudgeon及びR．Mersereau著「Multidimensional  Digital  Signal  Processing」（Englewood  Cliffs：Prentice-Hall，1983）に記載されているような多次元信号処理において一般に用いられる専門用語およびプロセスを使用する。本発明のこの明細書本文中の幾つかの用語は、以下の文脈において使用される。サンプリングパターンなる用語は、座標系に対してその位置が規定される空間内の複数の点から成る組のことを示している。サンプリングパターンの例が図２Ａ〜２Ｃおよび図３Ａ〜３Ｃに示されている。デカルトサンプリングパターンとは、２つの互いに垂直な平行線の組の交わりによって形成される複数の点から成る組のことである。連続画像という用語は、座標系上に規定される関数のことであり、例えばｆ（ｘ，ｙ）およびｆ（ｘ，ｙ，ｚ）はそれぞれ２Ｄ関数および３Ｄ関数である。デジタル画像は、サンプリングパターンにおける連続画像の値の配列である。より広義には、連続画像は、スプライン関数等の特定の基準組に対する級数展開の係数として役立つ数の配列によって表すことができる。なお、スプライン関数のうち、ゼロ次スプライン関数のデカルト積は、連続画像としてデジタル画像を表示するための良く知られた平方ピクセル形式を生じる。以下、数のこの配列をデジタル画像とも称する。デジタル演算装置に記憶された全ての画像はデジタルでなければならない。簡単のため、以下では、しばしば、デジタル画像および連続画像の両方を単に画像と称する。この場合、意味は文脈から推測される。この専門用語を用いると、ピクセル画像は、格子すなわち均一に離間された周期的なパターンであるサンプリングパターンに対応するデジタル画像であり、通常、必ずしもデカルトとは限らない。
【００２８】
１つのサンプリングパターンは、それがより少ない総数の連続画像サンプルを生じる場合には、他よりもまばらであると言われている。一般に、よりまばらなサンプリングパターンは、より低いサンプリング密度を有している。オーバーサンプリングとは、所望の精度で連続画像を表示するために必要以上に多くのサンプルを使用することである。対応するデジタル画像は、オーバーサンプリングされたと言われる。１つのサンプリングパターンに関して連続画像に対応するデジタル画像が与えられるとすると、異なるサンプリングパターンにおける同じ連続画像に対応する新たなデジタル画像を所望の精度で生成するプロセスは、デジタル画像リサンプリングと呼ばれる。アップサンプリングおよびダウンサンプリングはそれぞれ、密度が高い或いはまばらなサンプリングパターンに関するリサンプリングの特定の場合である。また、１ファクタ分のアップサンプリングまたはダウンサンプリングは、デジタル画像の変化を伴わず、アップサンプリングまたはダウンサンプリングの一形式と見なされる。デジタル画像付加は、同じサンプリングパターン対して規定される或いは基準に対する展開の場合には同じ基準に対して規定されるデジタル画像のポイント毎の付加のことである。低次元デジタルフィルタリングとは、その次元の一部のみに沿う多次元配列のデジタルフィルタリング、例えば２Ｄ矩形デジタル画像の各縦列の個々の１Ｄフィルタリングのことである。
【００２９】
連続画像の座標変換とは、回転、剪断、伸張または縮小等の演算のことである。デジタル座標変換を規定するため、座標変換によって関連付けられる２つの連続画像、および、共通のサンプリングパターンに対して連続画像を表す対応するデジタル画像を考慮する。他のデジタル画像から１つのデジタル画像を生成するプロセスは、デジタル画像座標変換と呼ばれる。これは、デジタルフィルタリングによって、すなわち、デジタル画像上での離散的な指数演算によって行うことができる。具体例としては、デジタル画像回転、デジタル画像剪断、デジタル画像伸張、デジタル画像縮小、そのような演算の組み合わせを挙げることができる（しかしながら、これらに限定されない）。デジタル画像座標変換を行うための方法は、例えばM．Unser，P．Thevenaz及びL．Yaroslavsky著「Convolution-based  interpolation  for  fast，high-quality  rotation  of  images」（IEEE  Transactions  Image  Processing，Vol.4，pp.1371-1381，1995）に記載されているように既知である。
【００３０】
幾つかのデジタル画像座標変換が図２Ａ〜２Ｃおよび図３Ａ〜３Ｃに示されている。図２Ａは、連続画像（長方形）の外形およびそれを表すデジタル画像におけるサンプリングパターンを示している。黒点による連続画像の値は、デジタル画像内に含まれている。図２Ｂおよび図２Ｃは、同じサンプリングパターンにおける回転された連続画像および剪断された連続画像をそれぞれ示している。この場合、黒点は、図２Ａのデジタル画像をデジタル的に回転し／剪断した形態に含まれる値を示している。
【００３１】
また、図３Ａは、デジタル画像を規定するサンプリングパターンおよび連続画像を示している。図３Ｂは、ｘ，ｙの寸法の幾つかの一定のファクタによって伸張された連続画像を示している。デジタル的に伸張された画像は、黒点による伸張された連続画像の値によって規定される。図３Ｃは、図３Ａと同じ連続画像を示しているが、特定の伸張ファクタによってサンプリングパターンの密度が高くなっている。図３Ｃにおける対応するデジタル画像は、図３Ｂと同じである。より一般的に言うと、デジタル画像伸張または縮小は、幾つかの規則性をもったサンプリングパターンにおいては、デジタル画像アップサンプリングまたはダウンサンプリングに相当する。
【００３２】
なお、０度分の回転、ゼロの剪断パラメータによる剪断、あるいは、１ファクタ分の伸張または縮小などの特定の座標変換の適用では、デジタル画像は不変のままであり、したがって、そのような特定の座標変換の適用は、結果を変えることなく、プロセスの説明において含まれても良く或いは省かれても良い。
【００３３】
連続画像のフーリエ変換によれば、サンプリング理論により、サンプリングパターンの特性を決定することができ、それにより、例えばD．Dudgeon及びR．Mersereau著「Multidimensional  Digital  Signal  Processing」（Englewood  Cliffs：Prentice-Hall，1983）に記載されているように、対応するデジタル画像が所望の精度の連続画像を表示する。同様に、デジタル画像の離散時間フーリエ変換（ＤＴＦＴ）によれば、どのようなデジタル画像座標変換が、所望の精度レベルで変換連続画像を表すデジタル画像を生成するのかを決定することができる。連続画像のフーリエ変換およびデジタル画像のＤＴＦＴの関連する特性は、総称して、フーリエ特性と称される。
【００３４】
加重逆投影演算は、逆投影により生成されるピクセルの位置によって決まる重みによる各投影の重み付けを必要とする。異なる重みは、A.C．Kak及びM．Slaney著「Principles  of  Computerized  Tomographic  Imaging」（ニューヨーク：ＩＥＥＥプレス、１９９８年）に記載されているように、例えば投影が取得された放射線源の位置に応じて、異なる投影において使用することができる。特定のケースとして、重み係数は１と等しくなり得る。これは、重み付けが無いことに対応しているが、重み係数である。具体的に示されていなければ、加重逆投影および非加重投影を総称して逆投影と称する。
【００３５】
このような基礎的な知識を用いて、本発明の幾つかの実施態様について説明する。
【００３６】
逆投影プロセッサ２８（図１）は、本発明を実施するために使用されるアルゴリズムを実行するようにプログラムされている。このアルゴリズムについて詳細に説明するが、最初に、より一般的な言葉で説明する。逆投影プロセスのステップが参照符号を用いてブロック図で示されている。
【００３７】
図４は、ステップ５０においてゼロ角度で個々の投影ｑ1…ｑpを逆投影することにより形成された回転中間画像の総和として、回転に基づく逆投影を示している。逆投影は、ステップ５２において座標変換に晒され且つステップ５４において１つの画像
【数９】
2007504912.tif 000011
を生成するために集められる複数の画像を生成する。この構造は、それ単独で従来の逆投影に相当しており、演算総数を何ら減らすものではない。しかしながら、この構造は、本発明の高速階層型逆投影方法の一部の導入への足掛かりとして役立つ。
【００３８】
図５Ａは、複数の投影ｑ1…ｑpから１つのピクセル画像
【数１０】
2007504912.tif 000012
を形成するための階層型逆投影方法を示している。各投影ｑmがステップ１００において逆投影されることにより、複数の中間画像Ｉ1,1…Ｉ1,Pが生成される。これは、階層型逆投影のゼロ番目の或いは予備のレベルである。ステップ１０２では、選択された中間画像に関してデジタル画像座標変換が行われる。変換された中間画像の一部（ここでは、対）はステップ１０４において集められ、これにより、集合中間画像Ｉ2,1，…Ｉ2,P/2が生成される。これは階層型逆投影の一番目のレベルである。一番目のレベルの集合中間画像は、階層型逆投影の次のレベルに対して入力される新たな中間画像としての機能を果たす。選択された中間画像に対してデジタル画像座標変換を適用するステップ１０６におけるプロセス、および、選択された中間画像を集めて新たな中間画像を生成するステップ１０８におけるプロセスは、ステップ１１６において全ての中間画像が処理されて集められることにより最終的な画像
【数１１】
2007504912.tif 000013
が生成されるまで続けられる。
【００３９】
必要に応じて、幾つかのレベル内の又は幾つかのレベルにわたる演算を組み合わせることができる。例えば、投影の一部に関して、および、図５Ｂに示されるように２つ以上（図５Ａの実施態様においては、正確に２つ）の選択された投影ｑpの選択された視野角での逆投影１１２によって代わりに生成された組Ｉ2,1，…Ｉ2,P/2からの対応する中間画像に関しては、ゼロ番目と一番目のレベルとを組み合わせることができる。あるいは、最初の中間画像の一部は、F．Natterer及びF．Wubbeling著「Mathematical  Methods  and  Image  Reconstruction」（フィラデルフィア：Society  for  Industrial  and  Applied  Mathematics，2001）に記載されるような直接フーリエ再構成方法等の明示的でない逆投影を含む等価プロセスによって生成することができる。
【００４０】
中間画像の成分投影の視野角および中間画像のフーリエ特性を明らかにし、それにより、後述するようにまばらなサンプルによって中間画像の集合体を正確に表すことができるようにするため、様々なデジタル画像座標変換のパラメータが選択される。これらのフーリエ特性は、中間画像の本質的なスペクトルサポート、すなわち、中間画像のフーリエ変換がゼロとは大きく異なっているフーリエ（周波数）領域の範囲に焦点を合わせている。サンプリング理論の教示内容によれば、連続画像のスペクトルサポートは、基礎を成す連続画像を表すデジタル画像を生成するサンプリングパターンの性質を決定し、それにより、連続画像を確実に再構成することができる。特に、図６は、階層アルゴリズムにおける中間画像
【数１２】
2007504912.tif 000014
の典型的なＶ字形（くさび形）スペクトルサポートを示している。
【００４１】
図７（Ａ）および図７（Ｂ）は、座標変換がデジタル画像回転となるように選択される際の、図５（Ａ）に示される二成分階層逆投影アルゴリズムにおける中間画像のフーリエサポートを示している。図７（Ａ）は、仮想中間画像
【数１３】
2007504912.tif 000015
のフーリエ領域サポートを示している。仮想画像
【数１４】
2007504912.tif 000016
は、対応する画像Ｉl,m中に含まれる投影の逆投影から成る。図７（Ｂ）は、中間画像の成分投影の視野角および中間画像のフーリエ特性を明らかにするように座標変換のパラメータを選択することにより、中間画像の垂直帯域幅（破線の長方形の高さ）を最小限に抑えることができ、それにより、まばらなサンプリングが可能になるとともに、計算要件を減らすことができることを示している。図７（Ｃ）は、Ｉl,m，Ｉl-1,2mおよびＩl-1,2m-1の空間領域サンプリング方式を示している。サンプリング点は空間領域内の点線の交点にある。
【００４２】
図８は、中間画像の三成分集合を行う本発明の他の実施態様を示している。この場合、Ｐ＝２×３L個の投影である。具体的に示すため、Ｐ＝１８個の投影の場合が示されている。図５（Ａ）の場合と同様に、ステップ１００においては、投影ｑθ1，ｑθ2，…，ｑθ18が逆投影されることにより、複数の中間画像Ｉd1,1…Ｉd1,18が生成される。これは階層型逆投影のゼロレベルである。
【００４３】
投影が３つにグループ分けされるとともに、各グループの選択された投影（図８の場合には、３つのうちの２つ）がステップ１０２においてデジタル画像座標変換を受ける。変換された中間画像の一部Ｉd1がステップ１０４において集められ、これにより、集合中間画像Ｉd2,1，…，Ｉd2,6が生成される。これは階層型逆投影の一番目のレベル（ｌ＝１）である。
【００４４】
階層型逆投影の二番目のレベル（ｌ＝２）では、選択された集合中間画像Ｉd2,mがステップ１０６においてｙ座標に沿う伸張またはアップサンプリングから成る座標変換を受け、また、ステップ１０８において、選択された集合中間画像に関して他の座標変換
【数１５】
2007504912.tif 000017
が行われる。ここでも、変換された中間画像がステップ１１０において３つのグループで集められ、それにより、中間画像Ｉd3,1，Ｉd3,2が生成される。レベル２からの選択された中間画像は、ステップ１１２においてデジタル画像座標変換を受けるとともに、ステップ１１４において集められ、それにより、次のレベルの中間画像Ｉd4,1が生成される。このプロセスは、投影の数に応じて、ステップ１１６で画像
【数１６】
2007504912.tif 000018
を生成するために必要な程度まで繰り返される。最後のレベル１１２において
【数１７】
2007504912.tif 000019
により示されるデジタル画像座標変換はデジタル画像回転であり、その前のレベル１０２及び１０８におけるデジタル画像座標変換もデジタル画像回転となるように選択することができる。ここでも、中間画像の成分投影の視野角および中間画像のフーリエ特性を明らかにし、それにより、後述するようにまばらなサンプルによって中間画像の集合体を正確に表すことができるように、様々なデジタル画像座標変換のパラメータが選択される。
【００４５】
図９は、図８の三成分逆投影アルゴリズムにおけるスペクトルサポートの進展を示している。番号（１），…，（９）は、図８中のアルゴリズムのブロック図における対応する点と一致している。図示のスペクトルは、デジタル画像Ｉl,mの離散時間フーリエ変換（ＤＴＦＴ）である。
【００４６】
図１０（Ａ）は、三成分二剪断型の階層型逆投影を使用する本発明の他の実施態様におけるアルゴリズムを示している。図１０（ａ）の実施態様は図８の実施態様と類似しており、そのため、適切な場合には図８の参照番号が使用されている。図８の場合と同様に、Ｐ＝２×３L個の投影であり、Ｌ＝２の場合が示されている。図１０（ａ）の実施態様は、２つの点で図８に示された実施態様とは異なる。第１に、ｘ座標に沿う更なるアップサンプリングステップが一番目のレベルの座標変換のステップ１０１に含まれている。第２に、ステップ１０２及び１０８において
【数１８】
2007504912.tif 000020
により示されるデジタル画像座標変換の少なくとも一部が、図１０（ｂ）に示されるように、一連の２つのデジタル画像剪断演算、すなわち、ｙ座標に沿う第１の演算（１２０）と、ｘ座標に沿う第２の演算（１２２）とから成る。
【００４７】
図１１Ａは剪断スケール逆投影を示しており、図１１Ｂは等価な階層型剪断スケール逆投影を示している。図１１Ａにおいて、複数の投影ｑ1，…，ｑ4は、ステップ１４０で逆投影されるとともに、ステップ１４２で剪断スケール座標変換を受ける。その結果として得られる中間画像は、ステップ１４４において集められる。
【００４８】
図１１Ｂにおいて、投影ｑ1，…，ｑ4は、ステップ１４６で逆投影されるとともに、ステップ１４８で剪断スケール座標変換を受け、ステップ１５０で一部（サブセット）に集められる。この集合によって生成された中間画像は、ステップ１５２において再び剪断スケール変換を受け、その結果として得られる中間画像はステップ１５４において集められる。このプロセスは、全ての投影および中間画像が処理されて集められることにより１つの画像
【数１９】
2007504912.tif 000021
が形成されるまで再帰的な形式で続けられる。ここでも、視野角およびフーリエ特性を明らかにし、それにより、サンプリング要件および必要な計算を減らすことができるようにデジタル剪断変換および中間画像のサンプリングのパラメータが選択される。
【００４９】
図１２は、中間画像のスペクトルサポートに対する画像剪断の影響を示している。図１２Ａは、投影の特定の一部が最終的な画像中に現れる際の当該特定の一部のスペクトルサポートを示している。図１２（Ｂ）は、同じ投影から成る中間画像のスペクトルサポートを示している。この場合、座標変換パラメータは、垂直方向で最も高いラジアン周波数を最小にするように選択されている。
【００５０】
図１３は、本発明の他の実施態様におけるアルゴリズム、すなわち、三成分階層型剪断スケール逆投影を示している。投影ｑθ1，…，ｑθ1は、ステップ１６０において逆投影されるとともに、ステップ１６２において剪断スケールデジタル画像座標変換を受ける。その結果として得られる画像の一部（サブセット）はステップ１６４において集められ、それにより生じる中間画像は、ステップ１６６および１６８においてアップサンプリングを受けるとともに、ステップ１６８においてデジタル剪断座標変換を受ける。これらの画像の一部はステップ１７０において集められ、また、それにより生じる選択された画像はステップ１７２において更なる座標変換を受ける。このプロセスは、全ての投影および中間画像が処理されてステップ１７４で集められることにより画像
【数２０】
2007504912.tif 000022
が生成されるまで続けられる。ここでステップ１７２において示される最終的な座標変換
【数２１】
2007504912.tif 000023
は、ピクセル画像
【数２２】
2007504912.tif 000024
におけるサンプリングパターンがデカルトである場合、ピクセルの再配置のみを含んでいる。
【００５１】
図１４は、図１２に示される中間画像のフーリエ特性の性質を使用して前述した実施態様における座標変換のパラメータを見つけるためのアルゴリズムである。
【００５２】
図１５は、本発明の他の実施態様におけるアルゴリズム、すなわち、オーバーサンプリング三成分二剪断階層型逆投影を示している。図１５の実施態様は図１０Ａの実施態様に類似しているため、適切な場合には図１０Ａの参照符号が使用されている。しかしながら、図１５の実施態様は、図１０Ａのステップに加えて、最後の１つ前にダウンサンプリングステップ１０９を含んでおり、これは図１５においては二番目のレベルである。ここでも、中間画像の成分投影の視野角および中間画像のフーリエ特性を明らかにし、それにより、まばらなサンプルによって中間画像の集合体を正確に表すことができるように、様々なデジタル画像座標変換のパラメータが選択される。しかしながら、後述するようにその後の処理の精度を高めるために特定の度合いのオーバーサンプリングが使用される。必要に応じて計算効率を高めるため、二剪断デジタル画像座標変換１０８（図１０Ｂに示される）を含む第２のｘ座標剪断とダウンサンプリングステップ１０９とを組み合わせて、剪断スケール変換を生成し、それにより、プロセス１０８及び１０９を共に図１５（Ｂ）に示されるプロセスに取って代えるようにすることもできる。
【００５３】
図１６は、本発明の他の実施態様におけるアルゴリズム、すなわち、オーバーサンプリング三成分階層型剪断スケール逆投影を示している。図１６の実施態様は図１３の実施態様と類似しているため、適切な場合には図１３の参照符号が使用されている。しかしながら、図１６は、中間のレベルにアップサンプリングの更なるステップ１６１を含んでいるとともに、最後のレベルの前のレベルにダウンスプリング１６９を含んでおり、これにより、中間画像がアップサンプリングされダウンサンプリングされる。図１５の実施態様と同様に、ｘ座標に沿うダウンサプリングステップ１６９  ↓xＵl,mがｘ剪断ステップ１６８
【数２３】
2007504912.tif 000025
の後に続く図１６の実施態様においても、計算効率のため、２つを組み合わせて１つのデジタル画像剪断スケールにすることができる。
【００５４】
［高速階層逆投影アルゴリズムにおける非デカルトサンプリング方式］
図８、１０および１３に示される階層型逆投影の異なる実施態様における中間画像は、効率的な非デカルトサンプリングに従う特定のスペクトルサポートを有している。特に、ｌ番目のレベルにおいて基礎を成す連続領域画像Ｉl,mは、例えば図６、７、９および１２から分かるように、フーリエ空間内にＶ字形を占めている。多次元サンプリング理論から分かるように、このようなスペクトルサポートを伴う画像においては、デカルトグリッドにおけるサンプリングが、適当な非デカルトグリッドにおけるサンプリングよりも効率性が低い。非デカルトサンプリングは、ベースバンドスペクトルのコピーをフーリエ平面内にしっかりと詰め込むことによりサンプリング効率を高めることができる。２Ｄサンプリングの説明のため［？］を参照されたい。特に、５点形(quincunx)サンプリング方式は、中間画像のサイズを減少させ、したがって、ほぼ２のファクタ分だけアルゴリズムの計算コストを減少させる。周期的な非デカルトサンプリングパターンにおけるデジタル画像座標変換は、一次元シフト不変フィルタを使用して効率的に実行できる。同様に、デジタル画像変換のパラメータの選択のための前述した全ての方法は、そのようなサンプリングパターンの場合にも同様に適用される。したがって、前述した全ての実施態様は、周期的な非デカルトサンプリングパターンを使用する実施態様にまで及ぶ。
【００５５】
既に説明した本発明の実施態様の変形例は、例えば陽電子放出断層撮影において起こるＸ線変換および磁気共鳴映像法において起こる３Ｄラドン変換を含む様々な形式の３Ｄ投影の３Ｄ逆投影に適用することができる。
【００５６】
三次元（３Ｄ）Ｘ線変換は、３Ｄの様々な方向での平行線の組に沿う積分の集まりである。各３Ｄ  Ｘ線投影は、線が方向付けられる３Ｄ角度により特徴付けることができる二次元関数である。３Ｄ  Ｘ線データのための階層型逆投影におけるブロック図は、前述したブロック図、例えば図８、１０、１３、１５および１６と同様である。この場合の中間画像は、三次元であり、前述した例の場合のような二次元ではない。各中間画像は、成分投影の３Ｄ角度の平均となる方向でまばらな３Ｄサンプリングパターンでサンプリングされる。アルゴリズムが進むと、３Ｄ  Ｘ線投影のフーリエ解析により説明されるように、このまばらなサンプリング（遅い）方向に沿うサンプルの密度は、その方向で増大する帯域幅に対応するように増大する。その結果、アルゴリズムの全ての段階で、画像が集められる前に、各画像は、この遅い方向に沿ってアップサンプリングされなければならない。また、３Ｄの実施態様で利用可能な余分な次元は、アルゴリズムの様々なステップでの使用に利用できる更に多くの座標変換、例えば３Ｄ回転を与える。２Ｄの場合と同様に、これらのデジタル画像座標変換のパラメータは、中間画像の成分視野角およびフーリエ特性を明らかにするように選択することができる。これらのデジタル画像座標変換は、前述した剪断や剪断スケール等の一連の一次元演算に分解することができる。２Ｄの場合と同様に、画質を向上させるため、三次元の任意の一部においてオーバーサンプリングが実施されてもよい。
【００５７】
３Ｄラドン変換投影は、一次元関数（原点からの平面の変位によってパラメータ化される平行な２Ｄ平面の組に沿う積分の集まり）である。投影の視野角は、平面の組に対して垂直なベクトルの３Ｄ方位角である。３Ｄラドン投影の階層型逆投影のブロック図は、図８、１０、１３、１５または１６の場合と同様である。一番目のレベルにおいて、投影は３Ｄ画像領域上にラドン逆投影される。これらの画像は、二次元に沿って一定であり、したがって、２つの一定の方向に対して垂直な方向でのみサンプリングすれば済む。３Ｄラドン投影のグループが組み合わされると、集合画像の帯域幅は、これらの成分投影の視野角に応じて、１つ又は２つの次元で増大する場合がある。したがって、中間画像を座標変換し且つ他の中間画像に加える前に、２つの次元で中間画像をアップサンプリングすることが必要である。２Ｄの場合と同様に、座標変換は、別個に行われてもよく、アップサンプリング演算と組み合わされてもよく、また、中間画像に関してオーバーサンプリングが実施されてもよい。
【００５８】
本発明の様々な実施態様において、デジタル画像座標変換およびダウンサンプリング演算またはアップサンプリング演算は、一連の低次元（一次元）線形デジタルフィルタリング演算によって実行されてもよい。また、使用されるサンプリングパターンが特定の周期性を有している場合には、後述するように、これらのデジタルフィルタがシフト不変フィルタであってもよい。計算効率のため、デジタルフィルタは、既知のように再帰フィルタ構造を使用して或いはＦＦＴを使用して実施することができる。デジタルフィルタにとって好ましい値を決定するための１つの方法は、M．Unser，A．Aldroubi及びM．Eden著「Fast  B-spline  transforms  for  continuous  image  representation  and  interpolation」（IEEE  Transactions  on  Patern  Analysis  and  Machine  Intelligence，Vol.13，pp.277-285，1991）；M．Unser，A．Aldroubi及びM．Eden著「B-spline  signal  processing：Part  I  -  theory」（IEEE  Transactions  Signal  Processing，Vol.41，pp.821-832，1993）；M．Unser，A．Aldroubi及びM．Eden著「B-spline  signal  processing：Part  II  -  efficient  design  and  applications」（IEEE  Transactions  Signal  Processing，Vol.41，pp.834-848，1993）において説明されているように、スプライン補間の理論を使用することである。
【００５９】
本発明の階層型逆投影方法は、任意の線源軌道を持つファンビームおよびコーンビームを含む離散ビーム構成等の広範囲の断層イメージング構成に適用できる。円形の走査軌道を持つ一般的なファンビーム構成が図１７に示されている。放射線源は、半径Ｒの画像の周りで、半径Ｄの円形の軌道上を移動する。線源角度βにおけるファンビーム投影は、ファン角度γによってパラメータ化される放射線の組に沿う線積分測定値に対応している。ＴSTは、画像ディスクの最も近い端までの放射線源の放射線に沿う距離であり、ＴENDは、画像ディスクの最も遠い端までの放射線源の距離である。
【００６０】
図１８は、本発明の一実施態様におけるアルゴリズム、すなわち、ファンビーム加重逆投影に適用できる階層三成分回転型投影を示している。図１８の実施態様は図８の実施態様と類似しているため、適切な場合には図８の参照符号が使用される。図１８の実施態様は、幾つかの点で図８に示される実施態様と異なっている。第１に、Ｐ＝４×３L個の投影であり、Ｌ＝２の場合が示されている。第２として、ゼロ番目のレベルでは、ステップ９９において、ここでは
【数２４】
2007504912.tif 000026
で示される加重逆投影により、最初の中間画像Ｉd1,mがファンビームから生成される。第３として、最後の段階では、２つの中間画像ではなく、４つの中間画像が集められる。第４としては、後述するように、階層の初期のレベルの殆どで使用されるサンプリングパターンは、好ましくは非デカルトとなるように選択される。ここでも、最後のレベルの１つ前で中間画像Ｉd3,mに関してデカルトサンプリングパターンが使用される場合、最後のデジタル画像座標変換は、画像ピクセルの並べ換えだけを必要とする。また、図５（Ｂ）の場合と同様に、階層のゼロ番目および一番目のレベルにおける演算を組み合わせることにより、中間画像Ｉd2,mがこれらの３つの成分投影の直接加重ファンビーム逆投影によって或いは他の手段によって生成されるようにすることもできる。
【００６１】
選択された数のレベルにおいては、図１９に示されるように、アップサンプリングステップおよびデジタル画像回転のカスケード（図１８におけるステップ１０６及び１０８）の前後に更なる重み付け（加重）ステップを含ませることにより、図１８の実施態様を修正することが有益である。中間画像Ｉ2,mは、ステップ１０６及び１０８の前後でそれぞれステップ１８０及び１８２において空間的に変化する重みにより加重される。後述するように、この重み付けは、中間画像のサンプリング要件を減らすため、したがって、計算を減らすために使用することができる。
【００６２】
中間画像のサンプリング方式はアルゴリズムの性能に影響を与える。望ましいサンプリング方式は、画像情報を失うことなく最も少ないサンプルを使用する方式である。ここでも、中間画像の成分投影の視野角および中間画像のフーリエ特性を明らかにし、それにより、後述するようにまばらなサンプルによって中間画像の集合体を正確に表すことができるように、様々なデジタル画像座標変換及びリサンプリング演算のパラメータを選択することができる。これらのパラメータを選択するための他の方法は、後述するように、放射線または曲線の特定の組の交わりに基づいている。
【００６３】
図２０Ａ〜図２０Ｃは、線源角度が間隔Δβで均一に離間された場合における、図１８に示された三成分階層アルゴリズムの複数のレベルにわたる中間画像の経過（進行）を示している。アルゴリズムにおいて中間画像を形成する投影のファン（扇）が示されている。ゼロ番目のレベルにおいて、各中間画像Ｉdl,mは、図２０（Ａ）に示されるように、ファン（扇）がβ＝０で方向付けられた１つの投影によって形成されている。再帰の一番目のレベルの後、各中間画像は、図２０（Ｂ）に示されるようなファン（扇）を有する３つの投影から形成される。二番目のレベルの後、各画像は、図２０（Ｃ）に示されるように９個の逆投影されたファン（扇）によって形成される。
【００６４】
図２１Ａおよび図２１Ｂには、図２０（Ｃ）に示される成分ファンを有する中間画像Ｉ3,mに関して、座標変換を選択し且つアルゴリズムでパターンをサンプリング／リサンプリングするための交差に基づく方法が示されている。Ｉ3,mにおけるサンプリングパターンは、中心の成分ファン上に位置する複数の点によって形成されており、これは図２０（Ａ）に示されるパターンと一致している。図２１（Ａ）および図２１（Ｂ）に示されるように、Ｉ3,mにおけるサンプリング点は、中心ファンのそれぞれの側に外側成分ファン（扇）を有する中心ファンの半分の交差によって決定される。
【００６５】
逆投影の精度を高め且つ計算要件を減らすためには、２つの制約を適用することにより、結果として得られるサンプリングパターンを修正することが有益である。第１に、放射線に沿うサンプルの密度は、最終的な画像の所望のサンプリング密度を超えないように制限される。第２に、サンプリングパターンは、それぞれの側の画像ディスクの外側で各放射線上に少なくとも１つのサンプルが含められるようにされる。図２２におけるプロットは、ファン角度γrにおける中心ファンの特定の放射線に沿うサンプリング点の位置を表示している。サンプリング点は一定の間隔γ’で位置している。ここで、γ’は、図２１（ａ）、（ｂ）及び（ｃ）に示される外側ファン（扇）のうちの１つのファン角度である。図２２の連続曲線上に位置する点は、図２１Ａおよび図２１Ｂに示される方法を使用して得られた当初の交点である。破線の曲線上の点は、上記２つの制約を使用して修正された点である。？図２３（Ｄ）？は、この修正された交差方法を使用して得られたサンプリングパターンの一例を示している。図示のファン（扇）は、このサンプリングパターンが生成された中間画像の中心ファンである。
【００６６】
本発明の前述した実施態様の場合と同様に、リサンプリング演算およびデジタル画像回転演算を一連の一次元演算に分解することが有益である。図１８において↑Ｕが記されたブロックは、中心ファン上でサンプリングされた中間画像を、同じファン上のサンプリング点の更に細かい組へとアップサンプリングすることを表している。これは、ファンの各放射線に沿う別個のアップサンプリングによって行うことができる。図１８において↑Ｕおよび
【数２５】
2007504912.tif 000027
によって記されたアップサンプリングおよび回転のカスケードに関しては、図２３（Ａ）〜図２３（Ｄ）に示されるように、１Ｄ演算への分解を一緒に行うことが便利である。４つのパネルのそれぞれにおいて、２つの点線の円は、半径Ｄの画像の境界および半径Ｒの（大きい）円上の線源軌道を表している。サンプリング点は、小さな円によって表されている。図２３（Ａ）に示されるサンプリングパターンおよび中心ファン（Fana）を有するデジタル中間画像は、リサンプリングされるとともに図２３（Ｄ）に示される中心ファン（Fand）の角度まで回転されなければならない。この場合、最終的なサンプリングパターンが図２３（ｄ）に示されている。これは、以下の２つのステップで達成される。
【００６７】
（ｉ）Fanaの各放射線に沿って個別に、画像を、FanaとFandとの交差によって規定される図２３（ｂ）に示されるサンプリングパターンへとアップサンプリングする。したがって、これらの点は図２３（ｃ）に示されるようにFand上にある。
【００６８】
（ii）Fandの各放射線に沿って個別に、画像を、図２３（ｄ）に示されるサンプリングパターンへとアップサンプリングする。
【００６９】
これらのステップは、１Ｄアップサンプリング演算を使用して、リサンプリングと回転とを組み合わせた演算を行う。
【００７０】
周期的なサンプリングパターンでサンプリングされた中間画像の場合において本発明の実施態様に関し説明したフーリエ解析技術は、離散ビームの場合における空間的に変化するサンプリング方式を考え出すことにまで及ぶ。これらの技術は、二次元および三次元の両方で任意の軌道上の投影および逆投影に適用されるのに一般に十分である。線または曲線の周期的でない系における逆投影の場合に関しては、局所スペクトルサポートの概念がスペクトルサポートの概念に取って代わる。これが図２４に示されている。この場合、図２４の左側に示される１つのファン（扇）の逆投影によって１つの中間画像が生成される。後述するように、ｒおよびθによってパラメータ化された図示の点におけるこの連続中間画像の局所スペクトルは、図２４（ｂ）の右側に示されるフーリエ領域における線分である。複数のファンビーム投影の逆投影から成る中間画像の１つの点における局所スペクトルサポートが図２５に示されている。左側には、その点の位置が示されており、また、中間画像における成分投影の視野角の範囲も示されている。右側の蝶ネクタイ形状の領域は、対応する局所スペクトルサポートである。
【００７１】
局所スペクトルサポートのこの解析は、中間画像における局所的なサンプリング要件を決定するために使用される。その結果として得られる空間的に不均一なサンプリングパターンが、２つの典型的な中間ファン逆投影画像に関して、図２６Ａおよび図２６Ｂに点線の円弧で示されている。画像の境界が破線の円によって示されており、また、この境界の外側では、数少ない点がとられるだけで済む。この局所フーリエサンプリング方法は、任意の寸法にわたる線、曲線または面にわたった任意の投影形状のためのサンプリング方式を見出すために直接に適用されてもよい。
【００７２】
階層型ファンビーム逆投影における座標変換およびリサンプリングのためのサンプリングパターンを決定するためのフーリエに基づく方法は、後述するように、中間画像の中心ファン以外の線や曲線上に位置するサンプリングパターンにまで更に及ぶことができる。そのような有利な組の例が図２７（Ａ）および図２７（Ｂ）に示されている。
【００７３】
［一般的な発散ビームアルゴリズム］
ファンビーム逆投影に関して前述した方法は、後述するように現代の診断放射線医学で最も重要なものの１つである３Ｄコーンビームを含む他の発散ビーム形状にまで直接に及ぶ。発散放射線の線源が画像化される対象物の周囲を所定の軌道（すなわち、円）に沿って移動するファンビームの幾何学的な形状と同様に、一般的な発散ビームの幾何学的形状においても、発散放射線の線源は、画像化される対象物の周囲を３Ｄ空間内の所定の軌道（必ずしも円形ではない）に沿って移動する。検出器面は、画像化される対象物を通して線積分を測定する。
【００７４】
一般的な発散ビーム形状における階層逆投影アルゴリズムの１つの実施態様も同様に図１８に類似したブロック図によって説明することができるが、以下の方法で修正が加えられる。第１に、ゼロレベルにおいては、ステップ９９において発散ビーム投影がゼロ逆投影される。この場合、適切な発散ビームシングルビュー(single-view)逆投影
【数２６】
2007504912.tif 000028
は、基準「ゼロ」線源位置に対応している。そして、これにより、最初の中間画像が生成される。第２として、線源軌道は必ずしも円形ではないため、中間画像の成分発散ビームは、ファンビーム形状の場合と同様に単に回転されるだけではなく、互いに対して並進されてもよい。したがって、ステップ１８における座標変換
【数２７】
2007504912.tif 000029
が選択される。第３としては、線源軌道および位置における対称性の存在に応じて、ファンビームアルゴリズムにおけるπ／２回転に取って代わるピクセル再配置などの「自由な」座標変換が存在しても良く或いは存在していなくてもよい。
【００７５】
ファンビームの場合と同様に、最初の中間画像は、アルゴリズムによって階層的に処理される。ファンビームの場合と同様に、位置および方向が互いに近い中間画像が集められ、それにより、集められた集合中間画像がまばらにサンプリングされてもよい。中間画像における３Ｄサンプリングパターンは、互いに対して回転され且つ並進される成分加重逆投影発散ビームの構造を検討することにより決定される。図２８に示される１つのサンプリングパターンは、中間画像の中心成分ビームに対応する発散ビームの放射線に沿ってサンプルが位置される場合のものである。各放射線に沿うサンプル間隔は、その中間画像の全ての成分発散ビームが放射線に沿って十分にサンプリングされるべく選択される。あるいは、より一般的な方法は、ファンビームの場合に関して前述したように、局所フーリエ方法を使用してサンプリングパターンを見つけることである。中間画像がどのようにしてその成分投影から成っているのかが分かれば、全ての中間画像の局所３Ｄフーリエ構造が決定される。ファンビームの場合において説明したように、局所フーリエサポートのためのサンプリング要件に適合する、中間画像の各点における３Ｄ局所サンプリング行列関数が見つけられる。この行列関数は、その後、サンプルの位置を決定するために画像領域にわたって積分される（可能であれば、数に関して）。
【００７６】
回転および並進と組み合わせて新たなサンプリングビームへとアップサンプリングする別個の方法が、ファンビームの場合と同様に行われる。この方法は、３Ｄ座標変換演算およびリサンプリング演算を、一連の１Ｄリサンプリング演算まで減少させる。図２９（ａ）に示されるように、各発散ビームは、方位角で分けられる垂直平面ファンビームの組の交わりとして見なされてもよい。この場合、傾斜平面ファンビームの組は異なる仰角を成している。分離できる座標変換のステップは以下の通りである（図２９（ｂ）に示されている）。
【００７７】
１．最初の発散ビームは、最初の放射線と新たな発散ビームの垂直平面との交点の組へとリサンプリングされる。したがって、これらの点は、最終的なサンプリング点によって共有される平面上に位置される。
２．ファンビームの場合における分離可能なリサンプリングのステップ１，２は、点の最後の組へとリサンプリングするように平面ごとに別個に行われる。
【００７８】
効果的なサンプリング方式を適切に用いると、発散ビームにおける高速階層逆投影アルゴリズムは、所望の精度で大きなスピードアップを達成することが期待できる。ファンビームの場合と同様に、連続的レベルで修正される（例えば、頂点の位置が原点から更に遠く離れるように移動する）擬似ビームを使用することもできる。
【００７９】
当業者であれば分かるように、本発明の方法は、幾何学的形状を画像化する実施例またはここで説明した特定の実施態様に限定されない。本発明の方法は、他の幾何学的構成を伴う逆投影の一般的な問題にも同様に適用することができる。
【００８０】
図３０は、線源位置βpにおける個々の投影の一般化された逆投影によって形成されるアップサンプリングされた中間画像の総和として、リサンプリングに基づく逆投影を示している。この逆投影は、回転に基づく逆投影４に類似しているが、２つの点でそれと異なる。第１に、最初のステップでは、４の場合と同様に、ゼロ位置においてではなく線源位置または方向βpにおいて、ｐ番目の（できれば処理された）投影が加重逆投影１８４を受ける。例えばファンビーム投影の場合では、各投影がその線源角度の方向で逆投影される。第２としては、ステップ１８８での加算（付加）による集合の前に、最初の中間画像のそれぞれがステップ１８６においてアップサンプリング演算を受ける。これは必要なことである。なぜなら、これらのＰ個の単一投影中間画像のそれぞれのサンプリングパターンは、効率的で且つまばらなパターンとなるように選択されるからであり、そのため、一般的には各投影毎に異なっており、しばしばデカルトではない。中間画像が集められる前に、これらの中間画像を共通の密度が高いグリッドへとリサンプリングする必要がある。この構造は、それ単独では、演算総数を何ら減らすものではない。しかしながら、この構造は、本発明の高速階層型逆投影方法の一部の導入への足掛かりとして役立つ。
【００８１】
図３１は、本発明の他の実施態様、すなわち、複数の投影ｑ1…ｑpから１つのピクセル画像
【数２８】
2007504912.tif 000030
を形成するためのリサンプリングに基づく階層型逆投影である。図３１は、図３０の二成分階層版である。最初に、非階層型の場合と同様に、ステップ１８４において、各投影は、その個々の方向で、その方向に適したサンプリングパターンへと逆投影され、それにより１つの中間デジタル画像が生成される。これは階層のゼロ番目のレベルである。階層の一番目のレベルでは、これらの中間デジタル画像が、ステップ１８６において、選択された画像対に共通の高密度サンプリングパターンへとアップサンプリングされる。このアップサンプリングは、通常、非デカルトサンプリングパターンに対してのものであるが、図２３または図２９（ｂ）に示されるように一連の一次元リサンプリング演算によって行うことができる。結果として得られた選択されたアップサンプリング画像は、ステップ１９０において対を成して集められ、これにより新たな中間画像が生成される。二番目のレベルにおいて、新たな中間画像は、ステップ１９２で再びアップサンプリングされるとともにステップ１９４で集められ、これにより新たな中間画像が生成される。このプロセスは、全ての中間画像および投影が処理されて最後の収集ステップ１９８の後に最終的な画像
【数２９】
2007504912.tif 000031
が生成されるまで続けられる。前述した実施態様の場合と同様に、レベル内およびレベルにわたって演算を組み合わせることができる。
【００８２】
階層中の各段階でのサンプリングパターンおよび図３１に示される実施態様におけるアップサンプリング演算のパラメータは、前述した任意の方法によって選択されてもよい。例えば、ファンビーム投影の場合、所定の中間画像における１つの可能なサンプリングパターンは２つのファンの交差部に位置している。すなわち、第１のファンは中心成分線源位置で方向付けられ、第２のファンは外側成分線源位置で方向付けられている。あるいは、アップサンプリングステップのパラメータおよびサンプリングパターンは、フーリエまたは局所フーリエ特性および中間画像中に含まれる投影の視野角に基づいて決定することができる。特に、周期的でないサンプリングパターンの場合には、サンプリングパターンを見つけるために、ファンビーム回転に基づくアルゴリズムに関して前述した局所フーリエ方法を使用することができる。すなわち、投影をどのように組み合わせて中間画像を形成できるかについての知識があれば、局所スペクトルサポートを決定することができ、また、この局所スペクトルサポートを使用して局所サンプリング行列関数を計算することができる。なお、この局所サンプリング行列関数は、画像領域にわたって積分されると、中間画像におけるサンプリングパターンを生成する。
【００８３】
図３２は高速階層型再投影のためのブロック図である。再投影は、所定の電子画像から断層投影を計算するプロセスである。再投影アルゴリズムは、フローグラフ転置(flow  graph  transposition)のプロセスを逆投影アルゴリズムの任意のブロック図に対して適用し、できれば重み付け演算の一部に変更を施すことにより見出される。プロセスにおいては、演算がそれらの随伴行列または双対形に取って代えられる。したがって、再投影のためのブロック図は、逆投影における対応するブロック図の逆バージョンに類似しているように見える。図３２に示される再投影プロセスは、図８に示される逆投影アルゴリズムから得られる１つのそのような再投影の実施態様である。
【００８４】
一番目のレベルにおいては、入力画像ｆのコピーが上側中間画像として線図の上側の枝に保たれ、また、下側の枝においては、ステップ２００において画像ｆが−π／２だけ回転され、それにより下側中間画像が生成される。二番目のレベルでは、線図の上半分において、回転されなかった中間画像がステップ２０２において３つの別個のデジタル画像座標変換を受ける。その一部は画像が不変のままであっても良く、それにより３つの異なる上側中間画像が生成される。同様のプロセスが下側の枝にも適用され、それにより３つの下側中間画像が生成される。上側および下側中間画像（全部で６個）のそれぞれは、その後、ステップ２０４においてデシメーション（ローパスフィルタリングの後、ダウンサンプリングする）のプロセスを受け、それにより新たな中間画像が生成される。図３２に示される実施態様の場合には、２・３L（Ｌ＝２）個の視野角しか存在せず、そのため、三番目のレベルが再帰的階層における最後のレベルとなる。三番目の最後のレベルにおいて、中間画像は、その一部の画像が不変のままであってもよい別個の座標変換（２０６）を受け、それにより１８個の中間画像が生成される。再帰の一部ではない最後のステップは異なっている。すなわち、各中間画像は、ステップ２０８においてゼロ度で再投影を受ける。ゼロ度での再投影は、ピクセルの垂直線を合計して一次元信号を生成することに相当している。これらの一次元信号（５１８）は、アルゴリズムによって生成される出力投影である。
【００８５】
アルゴリズムにおけるデジタル画像座標変換のパラメータは、中間画像のフーリエ特性の知識によって選択される。これらのパラメータは、対応する逆投影アルゴリズムのパラメータに対して単に関連付けられている。容易に分かるように、再投影ブロック図は逆投影ブロック図のフロー転置であることから、再投影ブロック線図の全ての枝は、対応する枝を逆投影ブロック線図中に有しており、また、対応する枝における座標変換は互いの数学的な随伴行列である。図１０の二剪断逆投影アルゴリズムに対応するこの再投影アルゴリズムの形態において、二番目のレベルにおける座標変換（２０２）ではｘ剪断後にｙ剪断が行われ、また、最後のレベルにおける座標変換（２０６）では、ｘ剪断後にｙ剪断が行われた後、ｘにおいて部分デシメーションが行われる（これらの３つの演算は剪断スケールまで減らすことができる）。これらの剪断のパラメータは、逆投影アルゴリズムで使用される対応するパラメータのマイナスである。ｘにおけるデシメーションのパラメータは、逆投影アルゴリズムの一番目のレベルでのｘにおけるアップサンプリングのパラメータと同じである。（図１３に示される剪断スケール逆投影アルゴリズムに対応する）このアルゴリズムの剪断スケール版において、二番目のレベルにおける座標変換（２０２）はｘにおける剪断である（各剪断のパラメータは図１３における対応するパラメータのマイナスである）。最後のレベルにおける座標変換（２０６）は剪断スケールである。再投影アルゴリズムで使用される剪断スケールは、逆投影において使用される対応する剪断スケールの数学的な随伴行列である。また、デシメーションファクタ（係数）のパラメータも逆投影における対応するアップサンプリングファクタ（係数）と同じである。
【００８６】
逆投影アルゴリズムの場合と同様に、アルゴリズムにおける中間画像のオーバーサンプリングは、まずアルゴリズムの最初に画像をアップサンプリングするとともにアルゴリズムの最後に同じファクタ分だけこれらの画像をダウンサプリングすることにより実施することができる。これらの演算は、レベル内およびレベルにわたって組み合わせることもできる。
【００８７】
図３３は、デシメーションに基づく加重再投影のブロック図である。図３１のフローグラフ転置と同様に、図３３は、１８個の異なる線源角度の投影組への画像の再投影を示している。最初に、所定の画像が２つの平行な経路に沿って処理される。各経路において、画像は、３つの異なるまばらなサンプリングパターンへの３つの平行なリサンプリング（２１０）を受ける。このリサンプリングは、一次元デシメーション（ローパスフィルタリングの後、まばらなグリッドへのリサンプリングを行う）を使用して別個に行うことができる。フィルタのパラメータは、所望の投影および中間画像のフーリエ特性によって決定される。所望の投影の局所フーリエ解析は、投影が平行な直線上に位置していない場合に使用される。アルゴリズムの最後のレベルにおいて、各中間画像は、まばらなサンプリングパターンへの３つの平行なリサンプリング（２１２）を再び受ける。最後に、画像に関して加重投影が行われ、それにより投影ｐθが生成される。これは、一次元投影を生成するための画像ピクセルの加重合計を含んでいる。
【実施例】
【００８８】
［実施および実験結果］
Ｃプログラミング言語において本発明の好ましい実施態様を実施し、また、３８４ＭＢ  ＲＡＭを有するＳｕｎ  Ｕｌｔｒａ  ５  ワークステーションによる数値実験で当該実施態様の試験を行った。試験画像はシェップローガンヘッドファントム（断層アルゴリズムの数値評価で使用される標準的な試験画像）であった。アルゴリズムのパラメータを変えることにより、精度と計算コストとの間でトレードオフを作ることができる。精度とは、再構成（復元）された画像の画質のことである。視覚の質は簡単に定量化できないが、本発明者らは、再構成された画像と、ラドン変換がそこから数値的に計算されたオリジナル画像との間の誤差を測定する。使用される誤差の指標は相対２乗平均平方根誤差（ＲＲＭＳＥ）である。ｆ［ｍ2，ｍ1］の断層投影からＮ×Ｎ画像
【数３０】
2007504912.tif 000032
を再構成する際のＲＲＭＳＥは以下のように計算される。
【００８９】
【数３１】
2007504912.tif 000033
【００９０】
平行ビームデータの場合には、試験画像が２５６×２５６のサイズを有し、視野角の数が４８６個であるとともに、個々の投影をフィルタリングするためにシェップローガンフィルタ（正弦関数によって乗じられる理想的なランプフィルタ）が使用された。図３０は、０．７５〜１．０の範囲のオーバーサンプリングパラメータγの様々な値における、２つのアルゴリズムでのＲＲＭＳＥ誤差と実行時間との間の関係を示している。二剪断アルゴリズムが円によって表され、また、剪断スケールアルゴリズムが正方形によって表されている。各アルゴリズムは２種類のフィルタ（ＭＯＭＳ１６と呼ばれる３次（破線）スプラインフィルタおよび５次（実線）スプラインフィルタ）を使用して実行される。アルゴリズムのそれぞれのタイプにおいては、γが減少するにつれて、アルゴリズムの誤差が減少するとともに、実行時間が増大する。他のものに何ら接続されないプロット点は、接続された点によって表されるアルゴリズムの非オーバーサンプリング版である。線形補間を使用する従来のアルゴリズムの実行時間と比べると、実行時間は１４ｓであり、また、その出力画像のＲＲＭＳＥ誤差は０．０４８６である（ここに示された高速アルゴリズムよりも悪い）。
【００９１】
平行ビームデータにおけるアルゴリズムの出力からの幾つかのサンプル画像が図３５に示されている。これらのサンプル画像は、左から右に向かって縦列方向で、従来の逆投影アルゴリズム、二剪断アルゴリズム、剪断スケールアルゴリズムからの出力画像である。γ＝０．８２のオーバーサンプリングが２つの高速アルゴリズムに対して適用された。画像の下側の横列は、上側の横列における対応する画像の部分を詳細に示している。
【００９２】
ファンビームの幾何学的形状における高速階層アルゴリズムは、５１２×５１２  ２Ｄシェップローガンファントムにおいてうまく試験がなされた。考慮された取得形状は、Ｄ＝１．０６×Ｎ＝５４４の線源半径、９７２個の線源角度、１０２５個の等角検出器を有していた。様々な補間方法を使用して、高速アルゴリズムの規則的なオーバーサンプリング版が実施された。その結果として再構成された画像が、線形補間を使用する従来のファンビームアルゴリズムと比較された。全ての実験においては、シェップローガンフィルタを使用して投影がフィルタリングされた。
【００９３】
図３６（Ａ）および図３６（Ｂ）は再構成された画像を示している。図３６（Ａ）は従来の再構成画像であり、図３６（Ｂ）は高速アルゴリズムの再構成画像である。高速アルゴリズムの結果は、従来のアルゴリズムの結果に匹敵している。高速アルゴリズムの点像分布関数は、従来のアルゴリズムにおけるそれに匹敵している。図３７（Ｂ）は従来のアルゴリズムのＰＳＦを示しており、図３７（Ｃ）は線形補間およびγ＝０．４オーバーサンプリングを伴う高速アルゴリズムのＰＳＦを示している。図３７（Ｃ）は、従来のアルゴリズムのＰＳＦのｘ軸を貫くスライスと、非オーバーサンプリングアルゴリズムと、γ＝０．４オーバーサンプリング高速アルゴリズムとを比較している。ＰＳＦの類似性から、高速アルゴリズムにおいて匹敵する画質が確認される。
【００９４】
使用されるサンプリング方式は、後に拡張が示唆されないファン交差に基づく方法であった。前述したこの方式の欠点とは別に、実施の容易さにより、実験で試験（検討）が行われた実施態様では、アルゴリズムの正しい演算のために必要のない多くのサンプル点が使用された。これらの非効率性にもかかわらず、Ｎ＝５１２およびＤ＝１．４１Ｎにおいては、高速（線形）アルゴリズムに対する従来のアルゴリズムにおける（データ依存の）乗算の比率が６．４であり、加算演算の比率が３．０であった。なお、原点に非常に近い線源に関して（Ｄ＝２．０６Ｒ）この実験で使用された幾何学的構成は、線源の頂点に近いサンプルの密度が高いことから、この非拡張型の実施においては特にやりがいがある。殆どの実用的なシステムにおいて、線源は、遠く離れており、これらの効果を低減する。また、前述した他のサンプリング方式を使用すると、スピードアップファクタがかなり高くなる。
【００９５】
［本発明で使用される逆投影アルゴリズムおよび再投影アルゴリズムの詳細な説明］
逆投影は、処理された投影から画像を形成するために使用されるプロセスである。再投影は、所定の画像から投影を計算するために使用される逆のプロセスである。両方の演算は、図１に示されるように投影からの画像再構成（画像復元）において使用される。最初に、従来の逆投影について説明し、その後、本発明の逆投影方法および再投影方法について説明する。また、２Ｄ画像の平行ビーム投影の場合に関して最初に説明を行い、その後、より一般的な２Ｄおよび３Ｄの幾何学的構成に関して説明する。
【００９６】
１つの画像をその投影から評価するために使用される伝統的な好ましい方法は、フィルタ補正逆投影（ＦＢＰ）アルゴリズムである。ＦＢＰは、いわゆるランプフィルタまたは修正されたランプフィルタを用いて投影を最初にフィルタリングした後、これらのフィルタリングされた投影を逆投影することを含んでいる。逆投影の前に更なる前処理ステップが投影に対して適用されてもよい。視野角θp（ｐ＝１，…Ｐは投影インデックスである）で取得された処理済みの投影は、変数ｔまたはθに対する依存がはっきりと示される必要があるかどうかに応じて、ｑ（ｔ，ｐ）またはｑpまたはｑθpにより交互に示される。以下では、処理済の投影を略して投影と称する。
【００９７】
Ｐ個の投影から成る組によるＦＢＰ再構成
【数３２】
2007504912.tif 000034
（長さＰのベクトル
【数３３】
2007504912.tif 000035
の成分によって特定される角度において）は、以下のように表される。
【００９８】
【数３４】
2007504912.tif 000036
ここで、
定義０．１  逆投影演算子は以下によって規定される。
【００９９】
【数３５】
2007504912.tif 000037
【０１００】
本発明のＯ（Ｎ2logＰ）座標変換に基づく階層逆投影アルゴリズムは、座標変換に関して逆投影の階層的分解に基づいている。座標変換の様々な選択とともに、結論的には最も一般的な座標変換の様々な選択とともに、幾つかの好ましい実施態様の詳細について説明する。
【０１０１】
特に、本発明のＯ（Ｎ2logＰ）回転に基づく階層逆投影アルゴリズムは、画像の回転に関して逆投影の分解に基づいている。
【０１０２】
画像
【数３６】
2007504912.tif 000038
をその回転バージョン
【数３７】
2007504912.tif 000039
に対してマッピングするθ回転演算子
【数３８】
2007504912.tif 000040
は、以下によって規定される。
定義０．２
【数３９】
2007504912.tif 000041
【０１０３】
ここで、
定義０．３  平面内での角度θ分の回転の行列は、
【数４０】
2007504912.tif 000042
である。
【０１０４】
単一角度θにおける１Ｄ関数ｑ（ｔ）の逆投影は、
【数４１】
2007504912.tif 000043
として規定される。特に、ゼロ角逆投影（θ＝０の場合）は以下の通りである。
【０１０５】
【数４２】
2007504912.tif 000044
定義０．４  単一θp逆投影中間画像は、
【数４３】
2007504912.tif 000045
である。
【０１０６】
図４から分かるように、式（３）の逆投影は、以下のように書き換えられてもよい。ここで、ｑpは、その視野角がθpである投影を示している。
【０１０７】
【数４４】
2007504912.tif 000046
【０１０８】
階層型逆投影の１つの実施態様は、幾つかの連続回転の累積結果が依然として回転であるという事実から生じている。特に、
【数４５】
2007504912.tif 000047
である。
【０１０９】
式（６）から分かるように、ある整数Ｌに関してＰ＝２Lである場合には、図４のブロック図を、図５に示されるような階層ツリー構造へと再構成することができる。ｌ番目のレベルのｍ番目の枝における中間画像はＩl,mとして示される。一番目のレベルにおいては、
【数４６】
2007504912.tif 000048
である。中間画像は、個々の投影を逆投影することにより生成される。また、その後のレベルにおいて、すなわち、ｌ＝１，２，３，…，log2Ｐに関しては、
【数４７】
2007504912.tif 000049
である。中間画像は、座標変換（この実施態様では、回転）を受けるとともに、加算によって集められる。投影角度θi，ｉ＝１，…，Ｐの任意の組においては、図４および５の構造が等価となり、それにより、図５に示される階層逆投影アルゴリズムが所望の画像
【数４８】
2007504912.tif 000050
を生成するように、中間回転角度δl,mを選択することができる。
【０１１０】
実際には、視野角（制約）よりも多くの中間回転角度（フリーパラメータ）が存在するため、δl,mの選択においては多くの自由度が存在する。また、実際に使用されるデジタル画像においては、その精度および計算コストが実施およびサンプリングパターンの選択によって決まるデジタル画像座標変換が使用される。これらの様々な自由度またはパラメータは、本発明においては、後述するように階層逆投影アルゴリズムの計算要件を減らすために使用される。
【０１１１】
定義０．５  １つの中間画像Ｉl,mにおける成分投影のインデックスの組Ｎl,mは、図５における入力から中間画像Ｉl,mまでの経路が存在する投影のインデックスの組である。
【０１１２】
例えば、図５から分かるように、Ｎ3,P/4＝｛Ｐ−３，Ｐ−２，Ｐ−１，Ｐ｝である。以下の式は、式（８）を使用して或いはブロック図を検討することよって簡単に分かる。
【数４９】
2007504912.tif 000051
【０１１３】
したがって、Ｎl,mは、中間画像Ｉl,mを形成する投影をリストアップしている。組Ｎl,mによってインデックスが付けられるこれらの同じ投影が、図５に示される方法によって処理されずに、それぞれの対応する視野角で一緒に直接逆投影された場合には、これにより、以下のような仮想中間画像
【数５０】
2007504912.tif 000052
が生成される。
【０１１４】
定義０．６
【数５１】
2007504912.tif 000053
【０１１５】
図５のブロック図を検討して分かるように、中間画像中の投影同士の間の相対角度が最終的な再構成画像
【数５２】
2007504912.tif 000054
において保たれる。言い換えると、アルゴリズムにおける全ての中間画像において、すなわち、∀ｌ，ｍにおいて、
【数５３】
2007504912.tif 000055
である。
【０１１６】
階層アルゴリズムの中間回転角度δl,mは、以下のようにαl,mによって完全に決定される。Ｎl,mの定義または式（９）から容易に分かるように、Ｎl+1,m＝Ｎl,2m-1∪Ｎl,2mであり、その結果、定義
【数５４】
2007504912.tif 000056
により
【数５５】
2007504912.tif 000057
である。
【０１１７】
式（８）、（１０）及び（１１）により、以下のようになる。
【数５６】
2007504912.tif 000058
【０１１８】
中間回転角度δl,m及びそれに応じた相対角度αl,mは、中間画像の帯域幅を減少させるように選択される。小さな帯域幅の画像は、少ないサンプルによって表すことができ、その結果、アルゴリズム全体の計算コストを減少させる。これらの中間画像の帯域幅は、フーリエ領域においてアルゴリズムを理解し且つ階層にわたる中間画像のスペクトルサポートの進展を辿ることにより決定される。
【０１１９】
［フーリエ領域における断層投影］
デジタル画像座標変換のパラメータは、前述したように、選択された投影の視野角θPを明らかにするように選択される。また、パラメータは、中間画像のフーリエ特性に対するそれらの影響を狙って選択される。これらの考慮は、いずれの中間画像を選択して一緒に集めるのかを決定するためにも使用される。
【０１２０】
フーリエ領域で逆投影アルゴリズムを解釈する手がかりは、投影
【数５７】
2007504912.tif 000059
の一次元フーリエ変換
【数５８】
2007504912.tif 000060
を画像ｆの二次元フーリエ変換
【数５９】
2007504912.tif 000061
に関連付ける投影スライス定理である。投影スライス定理は以下の通りである。
【数６０】
2007504912.tif 000062
【０１２１】
図４の逆投影アルゴリズムはフーリエ領域内で解釈できる。ゼロ演算子における逆投影は、そのスペクトルサポートが以下のように水平周波数軸ω1に限定される１つの画像を生成する。
【数６１】
2007504912.tif 000063
【０１２２】
空間領域で画像を回転させると、以下のように、そのフーリエ変換も同じ角度分だけ回転することは良く知られている。
【数６２】
2007504912.tif 000064
【０１２３】
以上から分かるように、式（５）における逆投影演算子の機能は、各投影のスペクトル成分を再構成画像中で適切な角度まで回転させてそれらを全て一緒に加えることである。
【０１２４】
片側帯域幅（ｔ変数）がπに等しい投影を想定すると、階層アルゴリズムにおける一般的な仮想連続中間画像
【数６３】
2007504912.tif 000065
は、図６に示されるように、成分投影によって決定される角度範囲にＶ字形のスペクトルサポートを有している。特に、Ｉl,mを図５（Ａ）に示されるアルゴリズムのブロック図において中間画像となるようにする。この場合、Ｎl,m＝｛ｂ，ｂ＋１，…，ｅ｝である。すなわち、この画像の成分投影の視野角は｛θp：ｐ＝ｂ，ｂ＋１，…，ｅ｝である。式（１０）の関係により、Ｉl,mのスペクトルサポートは、回転角αl,mでの
【数６４】
2007504912.tif 000066
のスペクトルサポートの回転バージョンである。図６においてＶ字形で占められる角度範囲は［γ−β，γ＋β］＝［θb−αl,m，θe−αl,m］であることが分かる。第１および第２の座標におけるＩl,mの帯域幅Ωs，Ωfはそれぞれ、図６に破線で示されるＶ字形における境界（有境）長方形上に示されている。したがって、αl,mの選択は、Ｉl,mの帯域幅を制御するための１つの方法を与える。
【０１２５】
サンプリング理論は、そのような連続中間画像Ｉl,mをデカルトサンプリングパターン上の１つのサンプルによってどのように表すことができるのかを決定付ける。特に、第１および第２の座標におけるΩs，Ωfの帯域幅を有する連続画像は、座標軸に合わせられるデカルトパターン上でサンプリングされてもよい。この場合、第１および第２の座標におけるサンプル間隔Δs及びΔfは、Δs＜π／Ωs及びΔf＜π／Ωfのようになる。連続画像を表すために必要なサンプルの数を最小にするためには、境界長方形の面積に比例する量１／（ΔsΔf）＞π2ΩsΩfを最小にしなければならない。この面積を最小にする回転角度、したがって、中間画像Ｉl,mのためのサンプリング要件を最小にする回転角度は以下のようになる。
【数６５】
2007504912.tif 000067
【０１２６】
［座標変換のためのパラメータの選択］
式（１５）は、図５に示される回転に基づく階層型逆投影の実施態様におけるデジタル画像座標回転に関して、最適な回転パラメータδl,mを決定するために式（１２）において使用される。インデックスｂ及びｅは、組Ｎl,mにおける最小のもの及び最大のものにそれぞれ対応するように決定される。この選択を用いると、中間画像
【数６６】
2007504912.tif 000068
となり、第１の座標（遅い帯域幅）および第２の座標（速い帯域幅）におけるその帯域幅はそれぞれ以下のようになる。
【数６７】
2007504912.tif 000069
【０１２７】
帯域幅を更に最小にするためには、差θe−θbを最小にしなければならない。これは、複数の投影からそれらの視野角間の最小最大角度分離をもって生成される中間画像を各ステップで集めることを選択することにより達成される。
【０１２８】
回転角度と、どの中間画像を集めるかの選択とに加え、アルゴリズムにおいては、様々な中間画像のデジタル画像座標変換に適したサンプルパターンを選択することが必要である。これらはサンプリング要件を使用して選択され、また、サンプリング要件は、中間画像の帯域幅およびスペクトルサポートを使用して決定される。
【０１２９】
特に、１つの投影の逆投影によって形成される最初の中間画像
【数６８】
2007504912.tif 000070
においては、図５の場合と同様に、遅い帯域幅Ωs（Ｉl,m）＝０であり、ｘ2座標に沿う１つのサンプルで足りる。図３ｂの場合と同様に複数の投影を逆投影することにより形成される最初の中間画像においては、遅い帯域幅が更に大きく、ｘ2座標に沿って更に多くのサンプルが必要とされる。
【０１３０】
図７（Ａ）〜図７（Ｃ）には、２つの仮想中間画像
【数６９】
2007504912.tif 000071
及びそれらの対応する中間画像の総和が空間領域および周波数領域において示されている。図７（Ａ）および図７（Ｂ）は
【数７０】
2007504912.tif 000072
およびＩのフーリエ領域サポートを示している。図７（Ｃ）は、Ｉl,m，Ｉl-1,2m，Ｉl-1,2m-1の空間領域サンプリング方式を示している。サンプリング点は、空間領域内の水平線と垂直線との交点にある。その各成分よりも大きな遅い帯域幅を有する集合中間画像Ｉl,mは、更に密度が高いサンプリングパターンを必要とする。逆に、階層アルゴリズムのそれよりも前のレベルの中間画像は、更にまばらなサンプリングパターンを必要とし、そのため、計算を節約できる。
【０１３１】
［計算コストおよび節約］
計算コストは簡単に見積もられる。同様に、投影視野角が［０，π）の範囲でＰ個あり、角度間隔がΔθ＝π／Ｐであり、Ｎl,mが式（９）で定められるとすると、式（１６）（ｂ＝２l-1（ｍ−１）＋１およびｅ＝２l-1（ｍ−１）＋２l-1を使用する）により、Ωs（Ｉl,m）＝πsin（π（２l-1−１）／（２Ｐ））及びΩf（Ｉl,m）＝πとなる。したがって、ｌ番目のレベルにおけるデジタル画像Ｉdl,mのサイズ（点の数）は以下のようになる。
【数７１】
2007504912.tif 000073
【０１３２】
Ｓ個のピクセルから成る１つのデジタル画像を回転させるコストがＯ（Ｓ）であり且つサイズ（Ｉdl,m）＝Ｏ（Ｎ2２l／Ｐ）であるとすると、アルゴリズムの複雑度は以下のようになる。
【数７２】
2007504912.tif 000074
【０１３３】
一般的なＰ＝Ｏ（Ｎ）において、これは、Ｏ（Ｎ2logＮ）であり、従来の逆投影のＯ（Ｎ3）よりも十分に有利な、画像のサイズを用いたコストのスケーリングである。
【０１３４】
［改良された三成分回転型階層逆投影］
サンプリングされた画像の一定の特定角度分（０，±π／２およびπ）の回転が計算上安価な演算であることは容易に分かる。なぜなら、これらの演算は、補間を伴っておらず、せいぜい既存のピクセルの再配置を含んでいるにすぎないからである。本発明のアルゴリズムの計算効率は、これらの自由演算をアルゴリズム中に組み入れることにより向上させることができる。二成分（バイナリ）アルゴリズム（図５）は、２つではなく３つの画像が段階毎に加えられるように修正することができる。そのようなトリプルのそれぞれの中心画像は、０ラジアンだけ回転される（自由演算）。−π／２の自由回転を使用するため、それは最後の段階で盛り込まれる。すなわち、成分画像のうちの１つは、−π／２分だけ回転されるとともに、他の回転されていない成分画像に対して加えられる。二成分段階と三成分段階とのこのような特定の組み合わせにより、ある整数Ｌにおいてサイズが２×３Lの視野角の組が得られる。投影の任意の組および数にアルゴリズムを合わせることができるが、アルゴリズムの記述および解析は、それらが正確に２×３Lを数えて以下のように均一に分配されると仮定すれば簡略化することができる。
【数７３】
2007504912.tif 000075
【０１３５】
したがって、視野角は、それぞれが３Lから成る２つの組に分けることができる。この場合、一方の組はθ＝０を中心とし、他方の組はθ＝π／２を中心とする。この三成分アルゴリズムのブロック図は、Ｌ＝２の特定の場合に関して、すなわち、Ｐ＝２×３2＝１８個の投影に関して、図８に示されている。デジタル中間画像Ｉdl,mは、基礎を成す連続中間画像Ｉl,mのサンプリングされたバージョンである。ブロック図において、
【数７４】
2007504912.tif 000076
のマークが付されたブロックは、ゼロ角の逆投影演算子を表している。
【数７５】
2007504912.tif 000077
のマークが付されたブロックは、デジタル画像回転演算子を表している。後述するように、デジタル画像およびデジタル画像回転におけるサンプリング周期および回転角度は、計算コストを最小にするように選択される。
【０１３６】
選択された集合中間デジタル画像を伸張して新たな中間画像を生成することができる。図８において↑yＵl,mの標示が付されたブロックは、第２の（遅い）座標に沿ったアップサンプリングを表している。このアップサンプリングは、アルゴリズムの進行に伴う中間画像の帯域幅の増大に対応するために必要である。式（１８）及び（２２）における視野角の設定に関し、ｌ番目のレベルの中間画像の帯域幅は、Ωs（Ｉl,m）＝πsin（Δθ（３l-1−１）／２）およびΩf（Ｉl,m）＝πである。（１），…，（９）の番号が付されたアルゴリズムの重要な点における正規化されたスペクトルサポートが図９に示されている。アルゴリズムが進行すると（ｌが増大すると）、遅い帯域幅が増大し、必要な遅いサンプリング周期Δs＜π／Ωs（Ｉl,m）＝sin（３l-1−１）／２）-1が減少する。様々なアップサンプリングブロックにおけるアップサンプリングファクタＵl,mは、これらの要件を満たすように選択される。
【０１３７】
セクション０．０．０．１は、２つの剪断を含む離散画像の分離可能な回転がどのようにして逆投影アルゴリズム中に組み入れられるかを示している。逆投影画像の画質を向上させるために、中間画像の系統的なオーバーサンプリングが導入される。セクション０．０．０．１は、このオーバーサンプリングを含むように修正されたアルゴリズムを示している。
【０１３８】
［座標変換における中間回転角度を以下のように選択できる。］
一般的なＬにおいて、階層は（Ｌ＋１）個のレベルから成る。ゼロレベルにおいては、各投影ｑθがゼロ
【数７６】
2007504912.tif 000078
で逆投影され、それにより対応する画像
【数７７】
2007504912.tif 000079
が生成される。この画像は、その後、３つにグループ分けされて組み合わされ、それにより第３の同数の画像が生成される。レベルｌでのグループ分けは以下の関係によって規定される。その後のレベルｌ＝１，…，Ｌにおいては、
【数７８】
2007504912.tif 000080
となる。
【０１３９】
式（１８）における視野角の設定において、最適な中間回転角度は、δl,3m-1＝０（ｌ＝１，２，…，Ｌ），δL+1,1＝０，δL+1,2＝−π／２である。図１０の検討時には、その連続画像Ｉl,mの投影のインデックスＮl,mが以下のようになることは容易に分かる。
【数７９】
2007504912.tif 000081
【０１４０】
したがって、式（２１）により、ｂ＝３l-1（ｍ−１）＋１およびｅ＝３l-1（ｍ−１）＋３l-1であり、式（１８）及び（１５）を使用すると、α*l,m＝−π／４（１−３-L）＋Δθ（３l-1（ｍ−１／２）−１／２）となる。二成分と同様に、ｌ＝１，…，Ｌにおいては以下のようになるのが分かる。
【数８０】
2007504912.tif 000082
【０１４１】
最後の段階において、本発明者らは、π／２ラジアンおよび０ラジアン分の自由回転を使用する。
【０１４２】
［０．０．０．１  座標変換は以下のように剪断に基づくことができる。］
中間画像のデジタル画像回転は、図１０（Ｂ）に示されるように、一連の２つのデジタル画像剪断に置き換えることができる。この場合、ｘ座標およびｙ座標の剪断は、以下のように規定される：
定義０．７  ｘ−剪断およびｙ−剪断演算子は以下によって規定される。
【数８１】
2007504912.tif 000083
ここで、ｘ−剪断行列およびｙ−剪断行列はそれぞれ
【数８２】
2007504912.tif 000084
である。
【０１４３】
この他の実施態様は、回転行列の２剪断因数分解から得られる。すなわち、
【数８３】
2007504912.tif 000085
となる。ここで、
【数８４】
2007504912.tif 000086
である。この分解により分かるように、デジタル画像の場合、完全な帯域制限された補間および十分に帯域制限されたデジタル画像を２剪断デジタル画像座標変換への入力として使用すると、２回剪断された画像は、単なるθ分だけ回転された画像であり、１／cosθのファクタ分だけ効率的にｘにおいてダウンサンプリングされるとともにｙにおいてアップサンプリングされる。サンプリングパターンにおけるこれらの変化を補正するには、２剪断演算が行われる度に非整数リサンプリングが必要となる。その代わり、中間画像を異なってサンプリングすることにより、最終的な画像に関しては累積効果を扱うだけで済む。
【０１４４】
都合のよい選択は、最初のデジタル中間画像をｘにおいてアップサンプリングし且つｙにおいてダウンサンプリングすることにより、階層の初めに１回リサンプリングを補正することである。図１０に示される実施態様の場合と同様にシングルビュー(single-view)逆投影によって最初の中間画像が生成されると、ｙにおけるダウンサンプリングはフリー（自由）になる。なぜなら、このダウンサンプリングは、異なるｙ密度をもったサンプリングパターンへの最初の逆投影だけを含んでいるからである。したがって、それは図１０のブロック図にははっきりと示されていない。最初の段階でのｘにおけるアップサンプリングは、後の段階でエイリアシングの問題を回避する。最初の画像のそれぞれにおいて、最終的な画像に対する累積効果的なダウンサンプリング／アップサンプリングファクタ
【数８５】
2007504912.tif 000087
が、最終的な画像への経路上における全ての２剪断変換にわたって計算され、また、最初の画像においては上記のように逆の演算が適用される。最初の画像Ｉd1,p，ｐ＝１，…，Ｐが異なるｘサンプリング密度およびｙサンプリング密度をもってリサンプリングされるのが分かる。
【０１４５】
そのため、結果として得られるアルゴリズムは、図１０（Ａ）に示されるように、（θpに応じた）適切なファクタによって全ての画像
【数８６】
2007504912.tif 000088
（ｐ＝１，…，Ｐ）を最初にリサンプリングし、その後、図１０（Ｂ）に示されるように２剪断デジタル画像座標変換を行うことを含んでいる。２つの剪断は、各中間画像のサンプリングパターンの変化を考慮して行われる。
【０１４６】
［剪断スケールアルゴリズム］
剪断スケールに基づく階層逆投影アルゴリズムは、スケーリング（伸張／縮小）されて剪断される単一投影画像の総計として逆投影を規定することに基づいている。
【０１４７】
定義０．８  ｘ−剪断スケール演算子
【数８７】
2007504912.tif 000089
は以下によって規定される。
【数８８】
2007504912.tif 000090
ここで、ｘ−剪断スケール行列は
【数８９】
2007504912.tif 000091
である。
【０１４８】
そのため、式（３）は以下のように書き換えられてもよい。
【数９０】
2007504912.tif 000092
【０１４９】
以下は、容易に示すことができる剪断スケールおよび累積された剪断の効果を表す幾つかの結果である。
【０１５０】
【数９１】
2007504912.tif 000093
【０１５１】
式（２４）と同様に、式（２５）は、一連のｎ個の剪断スケールも１つの剪断スケールであることが分かるように誘導的に適用できる。この特性を使用して、式（２３）を階層アルゴリズム構造に書き換えることができる。図１１（Ａ）は式（２３）（Ｐ＝４の場合）のブロック図表示であり、図１１（Ｂ）はその階層的等価物である。
【０１５２】
２つの構造の等価の要件は、その解が常に存在する、所定の投影角度θpと未知の剪断スケールパラメータ
【数９２】
2007504912.tif 000094
との間の方程式の１つの系を定める。
【０１５３】
非整数のスケールファクタ（倍率）を使用せずに整数のスケールファクタを使用することは、計算的に（演算または記憶スペースの点で）有利である。例えば、まさしく純粋なｘ剪断
【数９３】
2007504912.tif 000095
である整数剪断パラメータσ1＝１．０を有するｘ−剪断スケール
【数９４】
2007504912.tif 000096
を使用することが有利である。純粋なｘ剪断だけが使用されるように、また、レベルｌ＝２，３，…，Ｌにおいてｘ剪断スケーリングが避けられるようにアルゴリズムを設計することができる。
【０１５４】
一連の剪断スケールをもって式（２３）において各剪断スケールを置き換える代わりに、式（２４）および（２５）を使用して、各剪断スケールは、一連の剪断がその後に続く１つの剪断スケールに取って代えられる。すなわち、
【数９５】
2007504912.tif 000097
である。ここで、α1，α2，…，αnは自由パラメータである。
【０１５５】
その結果として得られる階層構造が図１３Ａに示されている。都合のよい選択は、レベルｌ＝１においてｘ−剪断スケール
【数９６】
2007504912.tif 000098
を行い、その後のレベルｌ＝２，３，…，Ｌにおいてｘ−剪断を行うことである。したがって、一番目のレベルにおいて基礎を成す連続領域画像における式は以下のようになる。
【数９７】
2007504912.tif 000099
【０１５６】
ｌ＝２，３，…，Ｌにおいて、中間画像Ｉl+1,mは、前のレベルにおける３つの中間画像（｛Ｉl,3m-i｝2i=0）のｘ剪断されたバージョンから成る。すなわち、
【数９８】
2007504912.tif 000100
である。
【０１５７】
また、０ラジアン分および−π／２ラジアン分の自由回転を使用するため、最後のレベルにおいては、
【数９９】
2007504912.tif 000101
である。
【０１５８】
図１３に示されるアルゴリズムにおける座標変換のパラメータ（様々な剪断スケールおよび剪断において）は過少決定される。実際には、一番目のレベルでは、関連する投影の角度によって（非整数）スケールファクタσ1だけが完全に決定される。中間画像を最もまばらにサンプリングできるようにするため、中間剪断ファクタ（｛αl,m｝l=2,3,…,L）は、中間画像の帯域幅を最小にするように選択される。そのような任意の組｛αl,m｝が与えられると、図１３の階層アルゴリズム全体によって生成された逆投影画像
【数１００】
2007504912.tif 000102
が逆投影式（２３）と等しくなるように、レベルｌ＝１においてσ2（剪断ファクタ）を常に選択することができる。特に、階層の最初のレベルにおける剪断スケールパラメータは、以下のように組｛αl,m｝によって決まる。
【数１０１】
2007504912.tif 000103
ここで、μ（ｐ，ｌ）＝［ｐ／３l-1］は、レベル１のｐ番目の投影からルートへと至る階層中の経路を表している。
【０１５９】
最後のレベルから逆方向へ向かう一連の剪断にしたがった、図１３のブロック図の上側半分ｐ≦Ｐ／２）を検討することにより分かるように、階層のこの上側半分において、基礎を成す連続中間画像Ｉl,mは、１つの剪断により、その対応する仮想中間画像
【数１０２】
2007504912.tif 000104
に関連付けられる。すなわち、
【数１０３】
2007504912.tif 000105
となる。
【０１６０】
式（１２）の微分と同様、中間剪断ファクタが以下のように｛βl,m｝によって決まるのが分かる。
【数１０４】
2007504912.tif 000106
【０１６１】
剪断係数αl,mの選択の自由は、βパラメータの選択の自由を与え、これにより、中間投影の帯域幅を最小にするようにβパラメータを使用できる。これにより、それらのサンプリング要件が減り、アルゴリズムの計算効率が向上する。
【０１６２】
中間画像Ｉl,3m-1のスペクトルサポートをＷl,3m-iで表し、また、そのｙ帯域幅（すなわち、遅い帯域幅）をΩy（Ｗl,3m-i）で表すとすると、Ωy（Ｗl,3m-i）を最小にする最適なβl,3m-i（ｉ＝０，１，２）は、以下のようになる。
【数１０５】
2007504912.tif 000107
【０１６３】
視野角が均一に分けられる場合には、
【数１０６】
2007504912.tif 000108
となるのが分かる。有利な選択は
【数１０７】
2007504912.tif 000109
であり、これにより、α*l,3m-1＝０となって、必要な剪断演算のほぼ１／３を削減できる。
【０１６４】
組｛βl,m｝は、階層の最後のレベルから最初のレベルに向かって逆方向にアルゴリズム中の中間画像のスペクトルサポートを辿ることにより決定される。スライス投影定理から分かるように、
【数１０８】
2007504912.tif 000110
のスペクトルサポートが以下のようになることは明らかである。
【数１０９】
2007504912.tif 000111
【０１６５】
サンプリング要件のため、スペクトルサポート
【数１１０】
2007504912.tif 000112
は、図１２にも示されるように、
【数１１１】
2007504912.tif 000113
として示されるその終点によって特徴付けられる。
【数１１２】
2007504912.tif 000114
【０１６６】
この組
【数１１３】
2007504912.tif 000115
の上側および下側の帯域エッジはそれぞれ
【数１１４】
2007504912.tif 000116
で表される。
【数１１５】
2007504912.tif 000117
である場合には、
【数１１６】
2007504912.tif 000118
となる。三成分組み合わせ規則により、これらのスペクトルサポート終点から成る上側、中間、下側の第３の組は、
【数１１７】
2007504912.tif 000119
として示される。すなわち、
【数１１８】
2007504912.tif 000120
となる。また、式（３０）およびアフィン変換に関するフーリエ理論に従えば、
【数１１９】
2007504912.tif 000121
となる。その結果、最適な剪断ファクタは以下のようになる。
【数１２０】
2007504912.tif 000122
【０１６７】
図１４に示されるアルゴリズムFINDALPHASは、剪断スケール逆投影アルゴリズム（図１３）においてレベルＬからレベル２へと階層を下側へ詳しく検討しながら、最適な剪断ファクタα*l,mを見つける。アルゴリズムへの入力は、（Ｌ＋１）番目のレベルの開始における画像のスペクトルサポートの終点から成る組である。すなわち、ＥL+1,m（ｍ＝１，２）である。特に、角度が均一に分けられた場合には、
【数１２１】
2007504912.tif 000123
となる。
【０１６８】
［アップサンプリングファクタの認定］
アルゴリズムにおけるｙ−アップサンプリングファクタ｛Ｕl,m｝および本発明の異なる実施態様は、遅い方向で中間画像のサンプリング密度を決定する。これらのアップサンプリングファクタは、整数値に制限されつつ、中間画像のサンプリング要件を満たすように選択することができる。この整数への制限は、計算上有利であり、同じ理由により分母の値が低い有理数の制限へと広げられてもよい。
【０１６９】
これらのファクタを選択するという問題は、剪断スケールアルゴリズムの場合には、２剪断の場合よりも若干簡単である。前者の場合、図１３を検討すると、中間画像Ｉl,mの遅いサンプリング間隔が
【数１２２】
2007504912.tif 000124
であることが容易に分かる。ここで、μ（ｌ’，ｌ，ｍ）＝［ｍ／３l'-l］は、Ｉl,mからアルゴリズムの最後のノードへと至る経路を表している。サンプリング理論は、遅い方向のサンプリング間隔が以下のようになっていることを要件としている。
【数１２３】
2007504912.tif 000125
ここで、Ωy（Ｗl,m）はＩl,mのｙ帯域幅である。アップサンプリングファクタ
【数１２４】
2007504912.tif 000126
の組が与えられたアルゴリズム全体の計算コストは、以下のようになるように表すことができる。
【数１２５】
2007504912.tif 000127
ここで定数ｃl,mは、中間デジタル画像Ｉdl,mに対して適用されたデジタル座標変換の相対計算コストを示している。これらの定数は、使用されるデジタルフィルタの次数によって決定されるとともに、座標変換の効果的な実施に関する記述で後述するように座標変換の特定の実施によって決定される。制約式（３５）に制約される
【数１２６】
2007504912.tif 000128
を最小にするための最良の組
【数１２７】
2007504912.tif 000129
は、動的なプログラミングまたは包括的なサーチによって解くことができる。視野角の所定の組において、アップサンプリングファクタの最良の組は、予め計算してルックアップテーブルに記憶させることができる。
【０１７０】
２剪断回転アルゴリズムの場合には、２剪断回転の使用により中間画像の事実上の部分アップ・ダウンサンプリングが生じるという事実により、問題が若干複雑になる。特に、アルゴリズムの隣り合うレベルにおける遅いサンプリング周期Δl,mは以下のように関連付けられる。すなわち、
【数１２８】
2007504912.tif 000130
となる。ここで、
【数１２９】
2007504912.tif 000131
である。
【０１７１】
そのため、遅い方向のサンプリング間隔に対する制約は以下のようになる。
【数１３０】
2007504912.tif 000132
ここで、Ωy（Ｗl,m）はＩl,mのｙ−帯域幅である。ここで、式（３６）によって与えられる計算コストを、式（３７）に制約されて最小に抑えることができる。
【０１７２】
［アルゴリズムのオーバーサンプリングバージョン］
オーバーサンプリングをこのアルゴリズムに組み入れて、再構成の精度を高めることができる。そのようにする１つの方法は、アルゴリズムにおいて全ての中間画像を均一にオーバーサンプリングすることである。すなわち、アルゴリズムにおいて補間が行われるときにはいつでも、演算される画像が、ある所定の定数γによって少なくともオーバーサンプリングされる。特に、デジタル座標変換またはリサンプリングに制約される全ての中間画像のサンプリング周波数に対するナイキストレートの比率（遅い方向および速い方向の両方において）は、好ましくはγよりも小さくなければならない。このオーバーサンプリングは、本発明のアルゴリズムの修正を極力少なくしつつ組み入れられることが好ましい。
【０１７３】
［遅い方向でのオーバーサンプリング］
本発明の前述したアルゴリズムにおいて、遅い方向でのサンプリング周波数は、アップサンプリングファクタＵl,mによって制御される。これは、アルゴリズムのオーバーサンプリングバージョンにおいてオーバーサンプリング条件を維持するための有用なツールであり、したがって、これにより、本発明のアルゴリズムの構造が変わらないことが分かる。１／γのファクタ分のアップサンプリングは、遅い方向のサンプリング間隔に対する制約をファクタγ分だけ単に修正することにより達成される。すなわち、
【数１３１】
2007504912.tif 000133
であることを要する。
【０１７４】
［速い方向でのオーバーサンプリング］
速い方向においては、オーバーサンプリングを制御するために計算上安価な整数アップサンプリングが使用されず、そのため、アルゴリズムは、部分アップサンプリングを含むように修正される。２剪断階層逆投影アルゴリズムにおいて、遅い方向でのファクタ｛Ｕl,m：ｍ＝１，２，…，２，３L｝分のそのような部分アップサンプリングは、レベルｌ＝１に既に含まれている。したがって、本発明者らは、アップサンプリングファクタＵl,mを単に増大させることにより、このオーバーサンプリングを組み入れ、その後、最後のデジタル座標変換が行われた後に画像をダウンサンプリングして、画像を所望のサンプリング方式に戻す（この場合、Δf＝Δs＝１．０）。したがって、２剪断アルゴリズムに対するこの修正は、１つの更なるレベルの部分リサンプリングしか伴わない。このｘ座標リサンプリングは、計算効率を向上させるため、Ｌ番目のレベルのデジタルｘ−剪断と組み合わせられる。このアルゴリズムのブロック図が図１５に示されている。
【０１７５】
剪断スケール階層逆投影アルゴリズムにおいては、ｘ−アップサンプリング演算がｘ−剪断スケール（アルゴリズムの初めにおいては、
【数１３２】
2007504912.tif 000134
）と組み合わされることにより、余分なレベルのリサンプリングが回避される。ｘ−アップサンプリングが、単に、剪断ファクタを０に設定した状態におけるｘ−剪断スケーリングの特定の場合であることは容易に分かる。これらの２つの演算の組み合わせもｘ−剪断スケールである。すなわち、
【数１３３】
2007504912.tif 000135
である。アルゴリズムの最後のダウンサンプリングをＬ番目のレベル
【数１３４】
2007504912.tif 000136
でｘ−剪断と組み合わせることにより、それらを有効にデジタル画像の剪断スケールにすることができる。αによるｘ−剪断の後にＵによるｘ−ダウンサンプリングを行うことは、有効なｘ−剪断スケールである。すなわち、
【数１３５】
2007504912.tif 000137
である。
【０１７６】
これらのアップサンプリングファクタおよびダウンサンプリングファクタの正確な値は、パラメータγおよび中間画像のスペクトル構造の両方によって決定される。速い方向において、オーバーサンプリング条件は、速い方向（ｙ座標）のサンプリング間隔がΔl,mf≦γπ／Ωf（Ｉl,m）（ｌ＝２，３，…，Ｌ）を満たすことである。そのため、必要な更なるアップサンプリングは以下の通りである。
【数１３６】
2007504912.tif 000138
【０１７７】
その結果、オーバーサンプリングされたアルゴリズムのアップサンプリングファクタ
【数１３７】
2007504912.tif 000139
は、以下のように、オーバーサンプリングされない２剪断アルゴリズムのアップサンプリングファクタＵl,mから修正される。
【数１３８】
2007504912.tif 000140
【０１７８】
すなわち、一番目のレベルにおけるアップサンプリングファクタは、アルゴリズムにおける全ての中間画像がパラメータγにしたがってオーバーサンプリングされるように修正される。その結果、最後の回転後、画像は１／ηという速いサンプリング間隔を有するようになり、したがって、Ｌ番目のレベルにおいて、η分だけ画像をダウンサンプリングして、これらの画像を単位サンプリング方式に戻さなければならない。
【０１７９】
なお、ここでは、入力投影がナイキストレートで正確にサンプリングされるものと仮定しているため、速い方向でのオーバーサンプリング条件は一番目のレベルの画像
【数１３９】
2007504912.tif 000141
（ｐ＝１，…，Ｐ）においては満たされない。
【０１８０】
この三成分オーバーサンプリング剪断スケ−ル型アルゴリズムのブロック図が図１６に示されている。
【０１８１】
［類似演算の倒れ込みシーケンス］
同じ１つの座標に作用する一連の２つの演算が存在する場合にはいつでも、これらの演算は、計算コストおよびリサンプリング精度を向上させるために組み合わされてもよい。以下は、前述した実施例に加えてｘ−剪断とｘ−アップサンプリング／ダウンサンプリングとを組み合わせた或いはｘ−剪断スケールとｘ−アップサンプリング／ダウンサンプリングとを組み合わせた他の実施例である。
【０１８２】
図１０の非オーバーサンプリング２剪断アルゴリズムは、Ｌ番目の段階において６つの中間画像のうちの４つをｙ−剪断した後にｘ−剪断することを含んでいる。本発明者らは、一連の２つの演算（ｘ−ダウンサンプリングの後に、この段階のｘ−剪断を行う）を１つの演算に倒れ込ませることによって、図１５のオーバーサンプリングバージョンにおいて画像の最後のダウンサウプリングを組み入れる。これにより、補間のカスケードの長さは、オーバーサンプリングされない場合と変わらない。
【０１８３】
［任意の基数ファクタ／分岐ファクタの階層］
視野角の組が２＊３Lの形式を成していない場合においては、全てのアルゴリズムを簡単に修正することができる。好ましい実施態様は、階層の全ての分岐が中間画像のトリプレットの集合体を含んでいる場合或いは±π／２または０ラジアン分の回転を使用する場合であるが、任意の数の中間画像が各段階で組み合わされてもよい。
【０１８４】
任意の数（例えばＭ個）の中間画像がアルゴリズムの特定のレベルで与えられると、それらのうちの３×［Ｍ／３］個は、組み合わされて３個から成るグループに分けられる（ここで、［ｘ］は、ｘ以下の最も大きい整数である）。中間画像の残りの数が２つである場合、これらは、その後、集められて１つの対を成し、次のレベルで１つの画像を生成する。１つの中間画像しか残らない場合、その中間画像は、変更されることなく階層の次のレベルへ送られてもよい。
【０１８５】
階層の分岐ファクタ（階層の１つのノードで集まる枝の数）は、視野角の任意の数に適応させないが階層の深さを減らすように変更されても良く、それにより、画質が向上されてもよい。その場合、画像を２つ（対）或いは３つ（トリプレット）ずつ集めるのではなく、それよりも多い数でグループ化することが有益な場合もある。
【０１８６】
座標変換のパラメータに関する先の規定は、任意の分岐ファクタ（分岐因子）を有するノードへと簡単に拡張させることができる。例えば、
【数１４０】
2007504912.tif 000142
となるような回転に基づくアルゴリズムにおいては、式（１０）及び（１５）に記載される関係が依然として保たれ、したがって、回転角度はδl,mi＝αl+1,m−αl,miによって規定される。
【数１４１】
2007504912.tif 000143
となるような剪断スケールに基づくアルゴリズムの場合には、式（２９）が依然として保たれる。集められる複数の投影から成るＭ個の組のうちのｉ番目の組を表すように標記Ｅil,mを拡張させると、式（３４）の右側と同一のα*l,miに関する式が得られる。
【０１８７】
［０．１  他の画像変換に基づく階層アルゴリズム］
逆投影式（３）は、以下の形式に書き換えられてもよい。
【数１４２】
2007504912.tif 000144
この式は、任意の行列Ａθに関するものであるが、Ａθの第１の横列が［cosθsinθ］に限られる。Ａθの残りの２つの入力（エントリ）に関しては任意の値を自由に選択できることにより、階層アルゴリズムで使用される座標変換の設計に自由度が得られる。行列Ａθは、θに関連する幾つかのパラメータベクトルδlに関して、
【数１４３】
2007504912.tif 000145
のように因数分解できるが、Ａθの２つの下側の入力（エントリ）の自由度に起因して完全に決定されない。この因数分解は、図５などの対応するブロック図を伴う逆投影式（４０）の階層的分解を得るために使用できる。この場合、
【数１４４】
2007504912.tif 000146
によって示される座標変換ステップは、
【数１４５】
2007504912.tif 000147
として規定される画像座標変換演算子を表している。
【０１８８】
ここで説明したこの発明の特定の実施態様は、行列Ａθ、その因数分解、関連する座標変換のこのような更に一般的な選択の特定の場合であることは言うまでもない。中間画像のフーリエスペクトルに対するこれらの座標変換の影響は、前述したケースと同様に解析される。なぜなら、空間領域内の行列Ａによるアフィン変換の結果は、周波数領域内の行列Ａ-TすなわちＡの逆転置によるアフィン変換だからである。したがって、同様の考え方を使用することにより、計算要件を減らす目的で、変換において自由パラメータを選択することができる。そのため、本発明の階層逆投影アルゴリズムで使用されるデジタル画像座標変換のクラスは、特定の実施態様に関して前述した変換に加えて、多くの他の変換を含んでいる。また、行列Ａ（δl,m）を三角行列の積へと因数分解することができるため、必要に応じて、一軸座標変換のカスケードとして座標変換を行うことができる。
【０１８９】
［０．２  デジタル画像座標変換およびリサンプリングの効率的且つ正確な実施］
本発明の階層逆投影・再投影アルゴリズムの精度および速度は、様々なデジタル画像座標変換演算およびリサンプリング演算の実施の仕様によって決まる。精度の向上には、一般に、高次のフィルタまたは補間が必要であるが、これは一般に計算コストを増大させる。高次フィルタリングまたは補間のコストは、演算を更に低次の演算に分解することにより減少させることができる。使用されるフィルタがシフト不変フィルタである場合には、計算及び／又はメモリの要件を更に減らすことができる。低次再帰フィルタを使用して、あるいは、高速フーリエ変換（ＦＦＴ）を使用することにより、高次有限インパルス応答フィルタを効果的に実施することができる。
【０１９０】
特に、連続画像に関して分離可能な表示基準がとられる場合には、デジタル画像配列中の個々の垂直線の部分遅延によってデジタル画像のデジタルｙ−剪断を達成できる。同様に、デジタルｘ−剪断を、１つの横列の同時の部分遅延として表すことができる。また、１Ｄ信号の部分遅延は、シフト不変フィルタを使用して達成することができる。同様の分解は、３Ｄ画像演算および剪断演算において知られている。
【０１９１】
デジタル画像のリサンプリングは、通常、更に低い次元の演算を使用して行うこともできる。この演算は、しばしば、シフト不変であり得る。例えば、整数ファクタＵによる１つの座標に沿う画像アップサンプリングは、Ｕ個の異なる計算上効率的な部分遅延に分解されてもよい。これは、本質的に、デジタル信号処理において良く知られている、いわゆる多相分解である。１つの座標に沿う合理的であるが整数ではないリサンプリングは、整数ダウン・アップサンプリングの１つのカスケードに分解できる。これはそれぞれ効率的に行われる。
【０１９２】
より一般的なデジタル画像リサンプリングも、更に低い次元の演算に分解することができる。ＣＦ1で示される一群の曲線上に位置するサンプルＳ1を有するサンプリングパターンからのデジタル２Ｄ画像ｆを、異なる一群の曲線ＣＦ2上に位置する点Ｓ2を有する他のパターンへとリサンプリングして、デジタル画像ｈを生成することを考える。２つの一群の曲線が十分な密集度で交差する場合には、１つのファン（扇）から他の回転されたファンへのリサンプリングに関して図２３を参照して前述した方法を一般的な曲線において使用できる。あるいは、ＣＦ1およびＣＦ2の両方と所望の密集度で交差する第３の一群の曲線ＣＦ3を導入することができる。この場合、デジタル画像は、ＣＦ1から、そのＣＦ3との交点へ、その後ＣＦ3とＣＦ2との交点へ、最終的にはＣＦ2上へと所望のサンプリングパターンにリサンプリングすることができる。
【０１９３】
このプロセスは、例えば曲線の代わりに面を考慮して、最初にプロセスを面上のリサンプリングのプロセスへと減らし、その後、２Ｄプロセスを使用して面上のリサンプリングを曲線上のリサンプリングへと減らすことにより、３Ｄへと一般化される。
【０１９４】
計算効率を向上させるため、デジタル座標変換およびリサンプリングは、しばしば組み合わせることができる。例えば、図１３に示される剪断スケールアルゴリズムにおいて使用されるデジタルｘ−剪断スケール演算は、個々の水平線上におけるリサンプリング演算へと分解することができる。
【０１９５】
［１  発散ビーム高速階層逆投影アルゴリズム］
［１．１  ファンビーム投影および逆投影］
対象物を中心とする円の上に検出器が等角度間隔で位置する場合を考える。二次元画像ｆ（ｘ，ｙ）の線源角度βでのファンビーム断層投影は、
【数１４６】
2007504912.tif 000148
によって示されるとともに、γによってパラメータ化され且つ原点から半径Ｄの円上の線源位置に中心付けられるファンの放射線に沿う線積分の組として規定される。関数
【数１４７】
2007504912.tif 000149
は、半径Ｒのディスクの外側でゼロの値をとると仮定する。
【０１９６】
図１７に示されるファンビーム放射線Ｖ（β，γ，Ｔ）は、以下によって規定される。
【数１４８】
2007504912.tif 000150
【０１９７】
線源角度βおよびファン角度γでのファンビーム投影は、
【数１４９】
2007504912.tif 000151
となる。ｆは半径Ｒのディスクの外側ではゼロであるため、ディスク内で、すなわち、ＴSTとＴENDとの間で積分を行うだけで済む。
【０１９８】
ファンビームの幾何学的構成を用いて計算された断層撮影においては、別個の線源角度の組｛βp：ｐ＝１，２，…，Ｐ｝で投影を利用でき、また、各ファン内で、放射線の角度には｛γj：ｊ＝１，２，…，Ｊ｝によってインデックスが付される。等角度のファンビーム構成の場合には、線源を中心とする円弧上で検出器が等しく分配され、そのため、ファン角度が均一に離間される。等間隔の検出器の場合、検出器は、線源位置から原点まで延びる線に対して垂直な線上で均等に分配される。ここで説明したファンビーム構成における高速逆投影アルゴリズムは、奇数Ｊおよびγj＝Δγ・（ｊ−（Ｊ＋１）／２）をもった等角度間隔の分配を仮定している。しかしながら、アルゴリズムは、他のファン角度構成へと容易に拡張することができる。
【０１９９】
Ｐ個のファンビーム投影の組｛Ｒβi（γ）：ｉ＝１，２，…，Ｐ｝からの再構成アルゴリズムは、その後に加重逆投影（５）を伴う各ファンビーム投影のスケーリングおよびフィルタリングとして表されてもよい。すなわち、
【数１５０】
2007504912.tif 000152
となる。
ここで、Ｗ（Ｔ）は適当な加重関数であり、
【数１５１】
2007504912.tif 000153
はそれぞれ、線源位置βにおける線源と画像点
【数１５２】
2007504912.tif 000154
との間の放射線に沿う距離、および、その放射線のファン角度であり、また、ｑβ（γ）は加重され且つフィルタリングされた投影である。
【数１５３】
2007504912.tif 000155
ここで、ｇ（γ）は適当なフィルタである。
【０２００】
定義１．１  １つの線源角度βにおける関数ｑの加重逆投影は、以下によって規定される。
【数１５４】
2007504912.tif 000156
また、加重逆投影演算子
【数１５５】
2007504912.tif 000157
は、以下によって規定される。
【０２０１】
【数１５６】
2007504912.tif 000158
したがって、
【数１５７】
2007504912.tif 000159
【０２０２】
また、以下のことが容易に分かる。
【数１５８】
2007504912.tif 000160
ここで、
【数１５９】
2007504912.tif 000161
はβ演算子による回転を示しており、その結果、以下のようになる。
【数１６０】
2007504912.tif 000162
【０２０３】
したがって、平行ビームの場合と同様に、Ｐ個のファンビーム投影の加重逆投影は、
【数１６１】
2007504912.tif 000163
を用いて、図４から分かるように、加重ゼロ逆投影画像の総和として表されてもよい。実際には、式（４３）及び（４４）に規定された逆投影演算子の近い相似形（アナログ）は、適当な関数定義
【数１６２】
2007504912.tif 000164
を用いた一般的な形態の投影からの二次元以上での関数の再構成に対してより一般的に適用される。本発明の方法は、座標変換
【数１６３】
2007504912.tif 000165
の適当な定義を用いて、これらの他の用途にまで及ぶ。
【０２０４】
［１．２  高速階層型逆変換］
ファンビーム構成における高速階層逆変換アルゴリズムは、平行ビームの場合におけるそれと類似している。このアルゴリズムは、投影角度βが互いに近い複数の投影によって形成される中間画像がまばらにサンプリングされるという事実を生かして、三成分階層構造でファンビーム投影を組み合わせる。そのブロック図が図１８に示されている。二成分または他の（可能であれば、混合して）基数構造においても同様のアルゴリズムを得ることができる。
【０２０５】
図１８において基礎を成す連続画像の組み合わせを規定する式は、以下の通りである。
【数１６４】
2007504912.tif 000166
【０２０６】
次に説明する実施態様においては、以下のように線源角度が［０，２π）の範囲で均一に離間されているとする。
【数１６５】
2007504912.tif 000167
【０２０７】
その後、平行ビームの場合と同様に、式（２２）を使用して中間回転角度が選択される。この場合、ΔβがΔθに取って代わる。
【０２０８】
図１８のブロック図に係るアルゴリズムの３つのレベルにおける中間画像の成分ファンが図２０に示されている。
【０２０９】
［１．２．１  ファンビームサンプリング方式］
ファンビームアルゴリズムの１つの実施態様は、次に説明する平行ビームの場合の相似形（類似形）によって得られるアルゴリズムにおける中間画像のサンプリング方式を使用する。変形例および改良例については後述する。
【０２１０】
単一ファンビーム逆投影画像（１つのファンビームの加重逆投影によって生成される１つの画像
【数１６６】
2007504912.tif 000168
）は、まばらなサンプリングに従う構造を有している。式（４３）から分かるように、
【数１６７】
2007504912.tif 000169
となる。したがって、β＝０で方向付けられるファン中のγによってインデックスが付された放射線上の特定のＴ＝Ｔ1における１つのサンプルは、放射線全体に沿う画像の値を特定するのに足りる。すなわち、このサンプルはＷ（Ｔ）に単に比例している。ここで、Ｔは、線源から放射線中の点までの距離である。
【０２１１】
アルゴリズムにおける中間画像は、それぞれが成分投影から生成されるそのような単一ファンビーム逆投影画像の幾つかの回転されたバージョンの総和である。本発明者らは、成分投影によってまたがる角度間隔の中心にある角度におけるファンを中心成分ファンと称する。これは、実際の投影（例えば三成分アルゴリズムの場合）に対応しており、あるいは仮想投影（例えば二成分アルゴリズムの場合）に対応している。例えば、図２０（ａ）のファンは図２０（ｂ）および図２０（ｃ）の両方における中心成分ファンである。
【０２１２】
平行ビームアルゴリズムにおいては、中心成分に位置合わされたデカルトサンプリングパターンで中間画像がサンプリングされることを思い出されたい。回転に基づく平行ビームの場合、中間画像は、平行な垂直放射線に沿ってサンプリングされる。これらの平行な放射線のそれぞれに沿うサンプルの間隔は、画像の他の成分投影がその放射線に沿って十分に（ナイキストサンプリング基準により）サンプリングされるように選択される。必要な間隔は、垂直放射線の交点の間隔に全く等しい。この場合、回転された平行な放射線の組は、極値的成分投影すなわち中心投影から視野角が最も離れた成分投影に対応している。
【０２１３】
平行ビームの直接的な相似形の場合、ここでは、中間画像は、中心成分ファンの放射線に沿ってサンプリングされる。中心ファンの放射線に沿うサンプリング点は、極値的成分ファンと中心ファンとの交点である。これにより、各放射線に沿って均一に離間されないサンプルが得られる。２つのファン
【数１６８】
2007504912.tif 000170
を考慮する。等角ファンビーム構成が奇数Ｊの検出器を有している、すなわち、
【数１６９】
2007504912.tif 000171
であると仮定する。β’−ファンと交差するβ−ファン（ファン角度γγに対応している）のｒ番目の放射線上の点は、以下の解を有する方程式
【数１７０】
2007504912.tif 000172
によって決定される。
【数１７１】
2007504912.tif 000173
【０２１４】
関数
【数１７２】
2007504912.tif 000174
は、ファンＶ（β’，・，・）が、ファンＶ（β，・，・）のｒ番目の放射線に沿ってどのように変化するのかを示している。これは、ファンの放射線に沿って変化するサンプリングレートに関する情報を持っている。傾斜がステップ状になればなるほど、サンプルがまばらになる。任意のγ’における局所的なサンプリングレートは、（ｄＴ／ｄγ’）-1に比例している。典型的なＴ（γ’）が図２２に破線で示されている。
【０２１５】
新たなサンプリング関数
【数１７３】
2007504912.tif 000175
を得るために、以下の２つの修正が
【数１７４】
2007504912.tif 000176
に対して成される。
【０２１６】
１．中間画像は、半径Ｒのディスク内でサンプリングされなければならない。この場合、サンプリングは、ディスクの外側の各放射線上の少なくとも１サンプルの余裕をもって行われる（１つのサンプルにおいてはＴ＜ＴSTであり、１つのサンプルにおいてはＴ＞ＴENDである）。
【数１７５】
2007504912.tif 000177
であると定義する。すなわち、γ’の値はＴENDに一致する。隣り合う線源角度を有するファンは、Ｔ＞ＴENDにおいて交点を欠いていてもよい。それを修正するため、Ｔ＞ＴENDにおける
【数１７６】
2007504912.tif 000178
の傾斜は、
【数１７７】
2007504912.tif 000179
におけるそれにとどめられる。
【０２１７】
２．最後の画像は、単位サンプリング間隔をもってデカルトパターンでサンプリングされる。この単位サンプリング間隔は、中間画像のサンプリングに必要な最も小さいサンプリング間隔である。γ番目の放射線に沿う単位サンプリングレートが達成される点γ’*は、γに関する方程式
【数１７８】
2007504912.tif 000180
を解くことにより決定される。単位サンプリング間隔制約を維持するために、γ’≧γ’*における
【数１７９】
2007504912.tif 000181
の傾斜は１／Δγに設定される。
【０２１８】
図２２は、典型的なβおよびγにおける
【数１８０】
2007504912.tif 000182
を示している。また、この図は、組
【数１８１】
2007504912.tif 000183
である放射線上の実際のサンプリング点のＴ−値を示している。整数ｊは、画像境界のそれぞれの側に少なくとも１サンプルの余裕が存在するように選択される。
【０２１９】
ここで概説した原理によって規定されるサンプル点の位置は、ファンの幾何学的構成、成分投影の組Ｎl,mに関する式（２１）、中間画像に関して選択された回転角度δl,m（等間隔の視野角の場合においては式（２２）で与えられる）、
【数１８２】
2007504912.tif 000184
における選択形式から計算できることは言うまでもない。
【０２２０】
［ファンビームサンプリングにおける分離可能な回転およびアップ補間］
概要で説明したように、アップサンプリング演算および回転演算し組み合わされたアップサンプリングは、計算上効率のよい一次元リサンプリング演算へ分解できる。
【０２２１】
［１．２．２  局所フーリエ構造に基づくサンプリング方式］
画像
【数１８３】
2007504912.tif 000185
が与えられると、本発明者らは、
【数１８４】
2007504912.tif 000186
が小さな欠くことができない帯域幅を有し、そのため、点の組
【数１８５】
2007504912.tif 000187
でサンプリングできるような、サンプリング関数
【数１８６】
2007504912.tif 000188
を見つけたくなる。本発明者らは、このサンプリング関数を以下のようにして見つける。直ぐ後で説明するように、
【数１８７】
2007504912.tif 000189
がその成分投影からどのようにして成っているかを知ることにより、本発明者らは行列関数
【数１８８】
2007504912.tif 000190
を見つける。この行列関数は、関数ｆが点
【数１８９】
2007504912.tif 000191
においてどのようにして局所的にサンプリングされるべきかを示している。その後、本発明者らは、この行列関数を画像領域にわたって積分することにより、サプリング関数
【数１９０】
2007504912.tif 000192
を得る。
【０２２２】
本発明者らのアルゴリズムにおいて、本発明者らは、中間画像ｆが幾つかの角度δp∈［βmin，βmax］において
【数１９１】
2007504912.tif 000193
の形式を有していることを知っている。Ｗ（Ｔ）による重み付けを無視し、投影ｑpがナイキストレートで正確にサンプリングされると仮定すると、本発明者らは、画像ｆ内の各点におけるスペクトルサポートの局所構造が分かる。帯域制限された各投影ｑpのファン逆投影により、１つの投影からファン逆投影された１つの画像は、ファンのスポークの方向に極僅かな空間帯域幅を有している。一方、垂直方向における当該画像の帯域幅は、図２４から分かるように、ファンの頂点からの距離に反比例する。これらのファンの組が回転されて一緒に加えられると、これにより得られる画像の局所的なスペクトルサポートは、個々のファンのスペクトルサポートの結合となる。
【０２２３】
例えば、図２５に示されるように、成分ファンがβmin〜βmaxの範囲の角度を有している場合、画像領域中の点
【数１９２】
2007504912.tif 000194
における局所的なスペクトルサポートは、θmin〜θmaxの範囲で方向付けられた半径を有し且つ角度θにおいてΩＲ／ｒθの径方向帯域幅を有するＶ字形にそらされる。ここで、Ωは、投影がナイキストレートでサンプリングされたと仮定した場合、画像面の中心にある逆投影された投影の空間帯域幅である。等角度の場合において、Ωは、投影のサンプル数を、逆投影された投影がカバーする画像面の原点を通る円弧の長さで割ったものである。数学的に、スペクトルサポートは以下のようになる。
【数１９３】
2007504912.tif 000195
ここで、
【数１９４】
2007504912.tif 000196
となる。
【０２２４】
画像内の１つの点
【数１９５】
2007504912.tif 000197
において局所的な二次元フーリエ構造の知識が与えられると、その点における行列関数
【数１９６】
2007504912.tif 000198
は、その点におけるスペクトルサポートが画像全体にわたって均一であった場合（例えば平行ビームの場合）に画像を効率的にサンプリングするサンプリング行列である。ここで対象となっている中間画像において、これは、２つの別個の小帯域幅サンプリング方向および大帯域幅サンプリング方向を生む。本発明者らは、最初の座標を前記小帯域幅方向となるように定める。
【０２２５】
本発明者らは、このサンプリング行列関数
【数１９７】
2007504912.tif 000199
を
【数１９８】
2007504912.tif 000200
によって概算する。ここで、
【数１９９】
2007504912.tif 000201
【０２２６】
ここで、
【数２００】
2007504912.tif 000202
である。角度θmidは、線源角度（βmin＋βmax）／２および
【数２０１】
2007504912.tif 000203
に対応する投影の角度を示している。スペクトルサポートに関する幾何学的な論拠によって選択されるこれらのパラメータにより、サンプリング行列による変換時のスペクトルサポートは［−π，π］×［−π，π］に制限される。
【０２２７】
画像全体におけるサンプリング関数は、このサンプリング行列関数を画像全体にわたって積分することによって（微分方程式：
【数２０２】
2007504912.tif 000204
（ｉ，ｊ＝１，２）を解くことによって）見出される。これは数に関して解かれてもよい。この正確な所定のパターンが使用されない場合であっても、局所フーリエサポート解析は任意のサンプリングパターンの有効性を評価する。
【０２２８】
結果として得られるサンプリングパターンが幾つかの中間画像に関して図２６（Ａ）および図２６（Ｂ）に示されている。ファンビームの場合において、この局所フーリエに基づく方法は、前述した交差に基づく方法によって得られるサンプリングパターンに類似するサンプリングパターンを生成する。１つのサンプリングパターンから他のサンプリングパターンへリサンプリングするための前述した分離可能な方法は、これらのパターンとともに使用することもできる。この局所フーリエサンプリング方法は、任意の寸法にわたる線、曲線または面にわたった任意の投影形状のためのサンプリング方式を見出すために直接に適用されてもよい。平行ビーム構成の場合、この方法は、前述した方法を軽減する。
【０２２９】
［１．２．３  他のサンプリング方式］
その頂点が線源軌道上にある１つのファン（「サンプリングファン」と呼ばれる）上にサンプルが配置されるサンプリング方式は、幾つかの欠点を有している。サンプリング点が各放射線毎に別個に選択されるため、二次元的には必ずしも最適とはならない方式になってしまう。アルゴリズムの最後の画像は均一な直角グリッド上に単位間隔をもってサンプリングされるが、前述した方式を使用する中間画像は、画像の特定の領域（例えばファンの頂点に近い領域）において必要な密度よりも高い密度でサンプリングされる。これを修正するため、前述した方式は、そのようなオーバーサンプリングされた領域でのサンプリングをまばらにするように修正されてもよい。そのような２つの可能性が図２７に示されている。これらの可能性の両方は、一番目のレベルでは画像がファン上に効率的にサンプリングされるが最後のレベルでは長方形グリッド上に画像をサンプリングする必要があるという事実を組み入れている。これらの方式は、ファン上のサンプリングからグリッド上のサンプリングへの漸進的な移行を組み入れようとしている。
【０２３０】
図２７（Ａ）において、サンプルは、その頂点が線源半径からよりも原点から遠く離れて位置する擬似ファン上で選択される。連続的なレベルにおいて、中間画像は大きな範囲の線源角度からのファンから形成され、原点からの擬似ファンの頂点の距離は増大する。最後のレベルにおいて、頂点は無限にある。すなわち、放射線は、長方形グリッドの平行な線である。図２７（Ｂ）において、サンプルは、ファン上ではなく、あまり発散しない線源位置に近いビーム上に位置されている。連続的なレベルにおいて、ビームは、更に平行になり、あまり発散しなくなる。二次元（または、３Ｄ画像の場合には３Ｄ）的な観点をとるこれらのサンプリング方式または他のサンプリング方式は、高速アルゴリズムに寄与する。
【０２３１】
［１．２．４  他の最適化］
・完全な軌道をもったファンビームの場合、線源角度は、［０，２π）の範囲で分けられる。これにより、補間を伴わずピクセルの再配置だけを伴う−π／２，−π，−３π／２分の回転をうまく利用できる。
【０２３２】
・アルゴリズムの全てのレベルにおけるサンプリングパターン中の点の位置は、予め計算してルックアップテーブルに記憶することができ、あるいは、データの処理とともにオンザフライ式で計算することができる。いずれの場合も、図２３（Ａ）および図２３（Ａ）において見られるような放射線指標（インデックス）に応じたサンプル位置の滑らかな変化をうまく利用することにより計算及び／又は記憶を減らすことができる。
【０２３３】
［１．２．５  オーバーサンプリング高速階層型逆投影］
平行ビームの場合と同様に、ファクタγ＜１．０分のオーバーサンプリングを使用して、精度を高めることができる。これをファンビームの場合において達成するための１つの方法は、ファン角度Δ’β＝Δβ・γを使用して密度の高いサンプリングパターンを決定することである。
【０２３４】
［１．２．６  ショートスキャンおよび他の角度分配のためのアルゴリズム］
この発明の方法は、いわゆるショートスキャンデータ取得フォーマットにも適用できる。この場合、線源角度の範囲は［０，２π）よりも狭い。違いは、逆投影において使用される重み付けにある。すなわち、ショートスキャンを担うために付加的な重み係数（例えばパーカー重み付け）が導入される。同様に、アルゴリズムは、線源角度の不均一な分配、あるいは、完全な画像の再構成を許容しないかもしれない線源角度の僅かな組の不均一な分配に一般化する。また、アルゴリズムにおける演算の順序は、当業者に知られた様々な方法（例えばパイプライン方式）で修正されてもよい。例えば、演算は、１つの投影が利用可能になると直ぐに処理を始めることができるように順序付けられてもよい。
【図面の簡単な説明】
【０２３５】
【図１】本発明において使用される装置のブロック図である。
【図２】図２Ａ、２Ｂ及び２Ｃは本発明の幾つかの実施態様で使用されるサンプリングパターンの図である。
【図３】図３Ａ、３Ｂ及び３Ｃは本発明の幾つかの実施態様で使用される更なるサンプリングパターンの図である。
【図４】逆投影の既知の方法を示す図である。
【図５Ａ】本発明の一実施態様におけるアルゴリズムを示す図である。
【図５Ｂ】図５Ａの実施態様において中間画像が生成される形式を示す図である。
【図６】図５Ａの中間画像を生成するために使用されるフーリエ特性を示す図である。
【図７】図７Ａ、７Ｂ及び７Ｃは、座標変換がデジタル画像回転である場合における、図５Ａに示される逆投影アルゴリズムでの中間画像のフーリエサポートを示す図である。
【図８】本発明の他の実施態様で使用されるアルゴリズムを示す図である。
【図９】図８のアルゴリズムにおけるスペクトルサポートの進展を示す図であり、ブロック（１…９）が図８中の対応する点と一致している図である。
【０２３６】
【図１０Ａ】本発明の実施態様におけるアルゴリズムを示す図である。
【図１０Ｂ】図１０Ａの実施態様で使用される画像座標変換を示している。
【図１１】図１１Ａは剪断スケール逆投影を示す図であり、図１１Ｂは階層型剪断スケール逆投影を示す図である。
【図１２】図１２Ａ及び図１２Ｂは中間画像のスペクトルサポートに対する画像剪断の影響を示す図である。
【図１３】本発明の他の実施態様におけるアルゴリズムを示す図である。
【図１４】最適な剪断ファクタを見つけるためのアルゴリズムである。
【図１５】本発明の他の実施態様におけるアルゴリズムを示す図である。
【図１６】本発明の更に他の実施態様におけるアルゴリズムを示す図である。
【図１７】円形の走査軌道を持つ一般的なファンビーム構成を示す図である。
【図１８】本発明の他の実施態様におけるアルゴリズムを示す図である。
【０２３７】
【図１９】図１８に示されるアルゴリズムにおける中間画像の重み付けを示している。
【図２０】図２０Ａ、２０Ｂ及び２０Ｃは図１８の二番目の階層レベルにおけるサンプリング点を示している。
【図２１】図２１Ａ及び図２１Ｂは図１８のアルゴリズムで使用されるサンプリングパターンの図である。
【図２２】図２０Ａ〜２０Ｃに示される方法を使用して得られる当初の交点を示す図である。
【図２３】図２３Ａ、２３Ｂ、２３Ｃ及び２３Ｄは図１８で使用されるアップサンプリングにおける回転のためのサンプリング点を示している。
【図２４】図１８のアルゴリズムによって生成される中間画像の１つの点における局所スペクトルサポートを示している。
【図２５】図１８のアルゴリズムで使用される不均一なサンプリングパターンの図である。
【図２６】図２６Ａ及び図２６Ｂはリサンプリングおよび座標変換におけるサンプリングパターンを示している。
【図２７】図２７Ａ及び図２７Ｂは本発明において使用できる他のサンプリング方式を示している。
【図２８】本発明で使用される発散ビームの２つの図を示している。
【０２３８】
【図２９】図２９Ａはコーンビームを示す図であり、図２９Ｂはリサンプリングを示す図である。
【図３０】リサンプリング投影において使用されるアルゴリズムの図である。
【図３１】本発明の実施態様におけるリサンプリングのために使用される他のアルゴリズムである。
【図３２】高速階層型再投影において使用されるアルゴリズムの図である。
【図３３】高速階層型再投影のための他のアルゴリズムの図である。
【図３４】本発明を使用する実験の結果を示すグラフである。
【図３５】本発明を用いて生成されるサンプル画像を含んでいる。
【図３６】図３６Ａは従来のアルゴリズムを使用する再構成画像の表示であり、図３６Ｂは本発明の高速アルゴリズムを用いて得られる結果を示している。
【図３７】図３７Ａ及び図３７Ｂは従来のアルゴリズムに匹敵するアルゴリズムの点像分布関数の図であり、図３７Ｃは本発明の高速アルゴリズムの点像分布関数を示している。
【図１】
2007504912.tif 000209
【図２】
2007504912.tif 000210
【図３】
2007504912.tif 000211
【図４】
2007504912.tif 000212
【図５Ａ】
2007504912.tif 000213
【図５Ｂ】
2007504912.tif 000214
【図６】
2007504912.tif 000215
【図７】
2007504912.tif 000216
【図８】
2007504912.tif 000217
【図９】
2007504912.tif 000218
【図１０Ａ】
2007504912.tif 000219
【図１０Ｂ】
2007504912.tif 000220
【図１１】
2007504912.tif 000221
【図１２】
2007504912.tif 000222
【図１３】
2007504912.tif 000223
【図１４】
2007504912.tif 000224
【図１５】
2007504912.tif 000225
【図１６】
2007504912.tif 000226
【図１７】
2007504912.tif 000227
【図１８】
2007504912.tif 000228
【図１９】
2007504912.tif 000229
【図２０】
2007504912.tif 000230
【図２１】
2007504912.tif 000231
【図２２】
2007504912.tif 000232
【図２３】
2007504912.tif 000233
【図２４】
2007504912.tif 000234
【図２５】
2007504912.tif 000235
【図２６】
2007504912.tif 000236
【図２７】
2007504912.tif 000237
【図２８】
2007504912.tif 000238
【図２９】
2007504912.tif 000239
【図３０】
2007504912.tif 000240
【図３１】
2007504912.tif 000241
【図３２】
2007504912.tif 000242
【図３３】
2007504912.tif 000243
【図３４】
2007504912.tif 000244
【図３５】
2007504912.tif 000245
【図３６】
2007504912.tif 000246
【図３７】
2007504912.tif 000247
【国際調査報告】
2007504912.tif 000248