(57)【要約】前処理された発散ビームシノグラム（２２）から電子画像を生成するための、発散ビーム逆投影のための高速な方法が提案されている。ここで、シノグラムは、発散ビーム投影の集まりである。本方法は、以下のステップ、すなわち、前記シノグラム（２２）を複数の部分シノグラム（ｇ１、ｇ２、ｇ３、ｇ４）に細分するステップ（１０１０、１０１２、１０１４、１０１６）と、前記部分シノグラムの重み付けされた逆投影を行って（１０１８、１０２０、１０２２、１０２４）、複数の対応する部分画像（ｆ１、ｆ２、ｆ３、ｆ４）を生成するステップと、前記部分画像（ｆ１、ｆ２、ｆ３、ｆ４）を集合させて（１０２６）、前記電子画像（ｆ、２６）を生成するステップと、を含んでいる。部分シノグラムへのシノグラムの細分は、反復手法で行うことが可能である。提案された方法の双対が、電子画像を再投影するための  −すなわち画像から発散ビームシノグラムを生成するための−  高速手段を提供する。
2005518892.tif 000002
【特許請求の範囲】
【請求項１】
逆投影を施せる発散ビームシノグラム（２２）から電子画像を生成する方法であって、
前記シノグラム（２２）を複数の部分シノグラム（図１０中のｇ１、ｇ２、ｇ３、ｇ４）に細分するステップ（１０１０、１０１２、１０１４、１０１６）と、
全体座標系内で前記部分シノグラムの重み付けされた逆投影を行って（１０１８、１０２０、１０２２、１０２４）、前記全体座標系内の正しい位置に複数の対応する部分画像（図１０中のｆ１、ｆ２、ｆ３、ｆ４）を生成するステップと、
前記部分画像（図１０中のｆ１、ｆ２、ｆ３、ｆ４）を集合させて（１０２６）、前記電子画像（図１０中のｆ、２６）を生成するステップと、
を備える方法。
【請求項２】
前記シノグラムの細分は近似（１２１０、１２１２、１２１４、１２１６）である請求項１記載の方法。
【請求項３】
前記シノグラムは、各部分シノグラムが所望サイズの部分画像を示すまで、反復手法で複数の部分シノグラムに細分され、前記細分するステップは、所望数の正確な細分（１５１０、１５１２、１５１４、１５１６）と、所望数の近似の細分（１５１８、１５２０、１５２２、１５２４）とを含んでいる請求項１記載の方法。
【請求項４】
前記正確な細分のステップは、シフティングなしのトランケーションを含んでいる請求項３記載の方法。
【請求項５】
前記部分シノグラムは、サイズが１ピクセル又は１ボクセルと同程度に小さい部分画像に相当する請求項３記載の方法。
【請求項６】
前記集合のステップは反復手法で行われる請求項１記載の方法。
【請求項７】
前記電子画像は断層撮影画像である請求項１記載の方法。
【請求項８】
データ座標内でのオーバーサンプリングを使用して前記電子画像の精度を向上させる前処理ステップ（１６）を更に備える請求項１記載の方法。
【請求項９】
前記シノグラムは前処理された発散ビーム投影を含む請求項１記載の方法。
【請求項１０】
前記シノグラムの細分ステップは、近似であり、前記シノグラムのトランケーション（１１１０）、シフティング（１１１２、１１２０）、及びデシメーション（１１１６）を含むことで、それぞれの生成された部分シノグラムが全体座標系内の部分画像に相当するようにする請求項１記載の方法。
【請求項１１】
前記近似の細分ステップは重み付け（１１１４、１１１８）をも含んでいる請求項１０記載の方法。
【請求項１２】
前記細分は反復的に行われ、該反復の連続段でのシフティング演算が結合されることで計算コストを低減する請求項１０記載の方法。
【請求項１３】
細分が空間において非等方的かつ非均一に行われることで、反復の或るレベルでの部分画像が同一形状又は同一サイズであることを必要としない請求項１記載の方法。
【請求項１４】
区切りがソース軌道の対称特性に従って選択されることで、計算コストと画像精度との間で所望のトレードオフを達成する請求項１３記載の方法。
【請求項１５】
前記方法が、複数のソース軌道から得られた複数のシノグラムに対して個々に適用されることで、最終画像を生成するのに結合される画像を生成する請求項１記載の方法。
【請求項１６】
前記正確な細分ステップは、部分画像用のシノグラム又は部分シノグラムにおける投影のトランケーション演算が、１以上の投影が該部分画像に利用可能になった後に開始されるようになされている請求項１記載の方法。
【請求項１７】
前記細分ステップは、部分画像用のシノグラムの細分が、１以上の投影が該部分画像に利用可能になった後に開始されるようになされている請求項１記載の方法。
【請求項１８】
複数の発散ビームシノグラムが利用可能になるにつれて、これらに対し時間的に連続して適用される請求項１記載の方法。
【請求項１９】
前記発散ビームシノグラムは収束ビームジオメトリ内で得られる請求項１記載の方法。
【請求項２０】
オブジェクトの電子画像（３６）を生成する装置であって、
前記オブジェクトを走査して前記オブジェクトの画像を表すデータを生成する手段（１２）と、
前記データを操作して、処理されたシノグラム（２２）を生成する手段（１６、２０）と、
前記シノグラム（２２）を複数の部分シノグラム（１０１０、１０１２、１０１４、１０１６）に細分する手段（２４）と、
全体座標系内で前記部分シノグラムのそれぞれを逆投影して（１０１８、１０２０、１０２２、１０２４）、前記全体座標系内の正しい位置に複数の対応する部分画像（図１０中のｆ１、ｆ２、ｆ３、ｆ４）を生成する手段（２４）と、
前記部分画像を集合させて前記電子画像を生成する手段（２４）と、
前記電子画像を格納、及び／又は表示、及び／又は解析する手段（３８）と、
を備える装置。
【請求項２１】
前記シノグラムの細分は近似（１２１０、１２１２、１２１４、１２１６）である請求項２０記載の装置。
【請求項２２】
前記シノグラムは、各部分シノグラムが所望サイズの部分画像に相当するまで、反復手法で複数の部分シノグラムに細分され、前記細分する手段は、所望数の正確な細分（１５１０、１５１２、１５１４、１５１６）と、所望数の近似の細分（１５１８、１５２０、１５２２、１５２４）とを行う請求項２０記載の装置。
【請求項２３】
前記正確な細分は、シフティングなしのトランケーションを含んでいる請求項２２記載の装置。
【請求項２４】
前記部分シノグラムは、サイズが１ピクセル又は１ボクセルと同程度に小さい部分画像に相当する請求項２２記載の装置。
【請求項２５】
前記集合させる手段は反復手法で演算する請求項２０記載の装置。
【請求項２６】
前記電子画像は断層撮影画像である請求項２０記載の装置。
【請求項２７】
前記操作する手段は、データ座標内でオーバーサンプリングを行って（１６）、前記電子画像の精度を向上させる請求項２０記載の装置。
【請求項２８】
前記シノグラムの細分ステップは、近似であり、前記シノグラムのトランケーション（１１１０）、シフティング（１１１２、１１２０）、及びデシメーション（１１１６）を含むことで、それぞれの生成された部分シノグラムが前記全体座標系内の部分画像に相当するようにする請求項２０記載の装置。
【請求項２９】
前記近似の細分ステップは重み付け（１１１４、１１１８）をも含んでいる請求項２８記載の装置。
【請求項３０】
電子画像から発散ビームシノグラム（３４）を生成する方法であって、
前記画像を複数の部分画像（図１７及び図９中のｆ１、ｆ２、ｆ３、ｆ４）に分割するステップ（１７１０、１７１２、１７１４、１７１６）と、
全体座標系内の前記部分画像それぞれの部分シノグラム（ｇ１、ｇ２、ｇ３、ｇ４）を計算するステップ（１７１８、１７２０、１７２２、１７２４）と、
前記部分シノグラムを集合させて前記シノグラム（３４）を生成するステップ（１７２６）と、
を備える方法。
【請求項３１】
前記部分シノグラムの計算は近似（１９２６、１９２８、１９３０、１９３２）である請求項３０記載の方法。
【請求項３２】
前記細分は、前記部分画像のそれぞれが所望のサイズを有するまで、反復手法で行われ（例えば図１３）、前記部分シノグラムの計算は、前記反復における所望数のレベルにおいては近似であり、前記反復の残りのレベルにおいては正確である請求項３０記載の方法。
【請求項３３】
前記集合は反復手法で行われる請求項３０記載の方法。
【請求項３４】
前記正確な集合は、シフティングなしで部分シノグラムを加算することによって行われる請求項３０記載の方法。
【請求項３５】
前記部分シノグラムの計算は、近似であり、部分シノグラムのシフティング（１８１０、１８１８）、補間（１８１４）、及び加算によって行われ、シフティング及び補間された各部分シノグラムは前記全体座標系内の部分画像に相当する請求項３０記載の方法。
【請求項３６】
前記部分シノグラムの近似計算ステップは重み付け（１８１２、１８１６）をも含んでいる請求項３５記載の方法。
【請求項３７】
オブジェクトの画像（３８）を生成する装置であって、
前記オブジェクトからデータを生成するスキャナ（１２）と、
前記画像の少なくとも１つの発散ビーム投影を生成するプロセッサ（１６）と、
少なくとも１つの投影から画像を再構成する、或いは反復再構成装置においては逆投影する、手段（２４）と、
前記再構成手段によって生成された画像中のエラーを検出する手段（２８）と、
エラーの訂正後に画像を再投影し、前記シノグラムに更なる訂正を加え、該訂正されたシノグラムを前記再構成又は逆投影手段（２４）に与える（２０）手段（３２）と、
前記エラーが訂正された後に前記再構成手段によって生成された画像を表示する手段（３８）と、を備え、
前記再投影手段は、前記画像を複数の部分画像（図１７中のｆ１、ｆ２、ｆ３、ｆ４）に分割し、全体座標系内で前記部分画像の部分シノグラム（図１７中のｇ１、ｇ２、ｇ３、ｇ４）を計算し、該部分シノグラムを集合させて（１７２６）前記画像シノグラム（図１７中のｇ）を生成する装置。
【請求項３８】
前記部分シノグラムの計算は近似である請求項３７記載の装置。
【請求項３９】
前記細分は、前記部分画像のそれぞれが所望のサイズを有するまで、反復手法で行われ、前記部分シノグラムの計算は、前記反復における所望数のレベルにおいては近似であり、前記反復の残りのレベルにおいては正確である請求項３７記載の装置。
【請求項４０】
前記集合は反復手法で行われる請求項３７記載の装置。
【請求項４１】
前記正確な集合は、シフティングなしで部分シノグラムを加算することによって行われる請求項３７記載の装置。
【請求項４２】
前記部分シノグラムの計算は、近似であり、部分シノグラムのシフティング、補間、及び加算によって行われ、シフティング及び補間された各部分シノグラムは前記全体座標系内の部分画像に相当する請求項３７記載の装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、発散ビーム断層撮影のための方法及び装置に係り、特には、視認可能なほどの劣化や重大な数値上の誤りなしに計算上の複雑さを低減する、発散ビーム断層撮影のための方法及び装置に関する。
【背景技術】
【０００２】
断層撮影の再構成（tomographic reconstruction）は、Ｘ線コンピュータ断層撮影法（ＣＴ）、ポジトロン放射断層撮影法（ＰＥＴ）、単一フォトン放射断層撮影法（ＳＰＥＣＴ）、及び磁気共鳴イメージング法（ＭＲＩ）用の或る捕捉(acquisition)法を含む診断用イメージング法（diagnostic imaging modalities）のほとんど全ての基礎となる周知の技術である。それはまた、非破壊評価（ＮＤＥ）や、より最近では空港での荷物検査にも広く利用されている。現在のＣＴスキャナにおける主データ捕捉モードは、ファンビームジオメトリであり、しばしばヘリカルスキャン又はスパイラルスキャンと組み合わされる。
【０００３】
より高速の画像取得のために、ＣＴスキャナ技術は、面検出器（area detectors）を用いたコーンビーム取得へと移行しつつある。生産されたＣＴスキャナの多くは、すでに、多列（multirow）検出器（典型的には１６列）を有することで、（狭角の）コーンビーム捕捉が得られている。広角のコーンビーム捕捉用の面検出器は、今後数年に渡って医用ＣＴスキャナに導入されるであろう。現在、隠された武器や爆破装置を検出するための高速の３Ｄ断層撮影画像法への必要性が存在する。ＰＥＴやＳＰＥＣＴのようなその他の方法においては、コーンビーム捕捉がすでに主フォーマットとなっている。
【０００４】
ファンビームＣＴ及びコーンビームＣＴにおける再構成上の問題は、特定の捕捉軌跡又は捕捉軌道上を進むソース（source）からそれぞれ扇形又は円錐として発散する線に沿った１組の線積分投影から、オブジェクトの２Ｄ又は３Ｄ画像を復元することである。従って、ファンビーム及びコーンビームのジオメトリは、発散ビームジオメトリの一例である。
【０００５】
断層撮影の再構成のための選択方法は、フィルタされた（filtered）逆投影（ＦＢＰ）又は回旋（convolution）逆投影（ＣＢＰ）であり、このどちらも重み付け（weighted）逆投影ステップを使用する。このステップは、２ＤにおけるＮ×Ｎピクセル画像にとってのＮ３のような尺度化（scaling）や、３ＤにおけるＮ×Ｎ×Ｎボクセル画像にとっての少なくともＮ４のような尺度化を行う計算上の必要性を伴う技術において、計算上のネックとなっている。従って、画像分解能を２倍にしてＮから２Ｎにしたとすると、計算がおおよそ８倍（３Ｄでは１６倍）増加する。これまで以上の大量データをリアルタイムで収集可能な新技術（例えば、多列検出器を備えた心臓画像法、干渉画像法）の出現や、３Ｄコーンビームジオメトリへの移行により、コンピュータがより一層高速になる一方で、高速再構成技術の必要性がまだ高まりつつある。高速再構成は、画像形成プロセスをスピードアップし、特殊目的の画像再構成コンピュータ（以下を参照）のコストを低減し、或いはそれら両方を可能にする。
【０００６】
逆投影の双対（dual）演算は、再投影（reprojection）であり、これは電子的に格納された画像の投影を計算するプロセスである。このプロセスも、断層撮影の再構成において基本的な役割を演じる。逆投影と再投影の組み合わせは、被験体の実際の３Ｄ画像法へのキーとなる、ヘリカルコーンビームジオメトリにおける長いオブジェクトの問題のための高速再構成アルゴリズムを構築するのにも利用可能である。更に、各種応用において、逆投影ステップと再投影ステップの両方が１つの画像の再構成のために数回行われる、反復の再構成アルゴリズムを使用することは、有利であり、或いは必要でさえもある。逆投影ステップと再投影ステップをスピードアップすることは、そのような反復方法の経済的な実行可能性を決定することになろう。
【０００７】
高速方法について議論する場合、断層撮影データ用の２つの主フォーマット、すなわち（i）平行ビームと（ii）発散ビーム及びその特別な場合、を互いに区別することが重要である。平行ビームフォーマットは理論上の操作と計算上の操作の両方にとって最も便利であるのに対し、発散ビームは市販の医用ＣＴスキャナや産業用ＣＴスキャナにおいて最も一般的に見出されるフォーマットである。元々、２Ｄ法のファンビーム再構成は、ヘリカル及びマルチスライスヘリカル３Ｄ（体積測定の）再構成のための技術状態法（stage-of-the-art）におけるキー要素であるが、新たな設計及び将来の設計は、コーンビームジオメトリへと移行することが予測される。従って、２Ｄにせよ３Ｄにせよ、発散ビームジオメトリは主たる画像フォーマットであり、おそらく近い将来においてもそのままであろう。
【０００８】
オブジェクトからの大きなソース距離の限界では、発散ビームジオメトリは平行ビームジオメトリにまで減少する。従って、（本発明の方法を含む）発散ビームジオメトリのための全てのデータ処理方法は、平行ビームジオメトリにも応用可能である。しかし、その逆は真ではない。すなわち、実際に見出された小さな又は適度のソース距離では、２つのジオメトリは十分に相違しており、平行ビーム処理方法は、これを発散ビームデータに直接応用した場合は劣った結果又は不満足な結果を生じる。従って、発散ビームジオメトリは、平行ビームジオメトリとは異なる処理を必要とする。
【０００９】
本発明は、発散ビーム断層撮影法を扱う。再構成方法の幾つかは、再ビニング（rebinning）と呼ばれるプロセスを利用することで、発散ビームデータを平行ビームデータに再配置（又は補間）し、これは次に平行ビームアルゴリズムによって処理される。３Ｄ発散ビームジオメトリのための他の方法は、コーンビームデータを３Ｄランダム変換データ（又はそれの派生データ）に変換することを利用する。その代わりに、限定による元々の発散ビーム法（本発明はその高速バージョンである）が、先立って平行ビームデータを再ビニングすることや、そのデータを３Ｄランダム変換領域に変換することなしに、発散ビームデータを直接処理する。元々の発散ビーム法における再構成プロセスは、投影の前処理（例えば重み付けやフィルタリング）、それに続く重み付けされた発散ビームの逆投影演算、及びおそらくは発散ビームの再投影演算、又は一連のそのような演算によって構成されている。従って、そのような方法は、発散ビームのフィルタされた逆投影（ＤＢ−ＦＢＰ）アルゴリズムと呼ばれる。それらの計算は、平行ビームの場合には、通常、逆投影演算及び再投影演算によって支配される。
【００１０】
発散ビームの再投影は、まず平行ビームの再投影を行い、次にその結果を再ビニングすることによってなし得るが、ここでも、再ビニングを使用せずに元々の発散ビーム再投影アルゴリズムを使用するのが有利である。
【００１１】
再ビニングの欠点は、スピードアップを制限する、アーチファクト（artifacts）及び追加の計算を含むことである。再ビニングはまた、処理が開始可能になる前にデータの多くの部分の捕捉及び操作を必要とし、これはスピードを再び制限することになる。３Ｄランダム領域への変換を利用する方法は、同様な欠点を有している。従って、発散ビーム断層撮影のための方法及び装置にとって、これらの欠点を解消し、高度に柔軟性があり、従来の発散ビーム逆投影を使用する方法と比べて大きな性能向上を提供する必要がある。
【００１２】
特殊目的のハードウェアが、逆投影プロセスをスピードアップするための伝統的な手段である。カスタムチップ、マルチプルプロセッサ、又はそれらの組み合わせを使用する特別なタイプのコンピュータが、必要な計算を行うために設計されている。
【００１３】
このアプローチの欠点は、マルチプルプロセッサやカスタムハードウェアのコストを含み、また、一般目的のコンピュータの性能が急速に向上しているので、要求される特殊目的の構造がすぐに時代後れになってしまうという事実を含んでいることである。従って、発散ビームジオメトリにおける高速処理は、特殊目的のハードウェアを必要とせず、かつ、標準のシリアル又はパラレル構造において容易に実行されることで、より高いコスト効率にする必要がある。
【発明の開示】
【００１４】
従って、本発明の目的は、発散ビーム断層撮影のための新規で改良された方法及び装置を提供することにある。
他の目的は、視認可能なほどの劣化や重大な数値上の誤りなしに計算上の複雑さを低減する、発散ビーム断層撮影のための新規で改良された方法及び装置を提供することにある。
【００１５】
本発明の１つのアスペクトによれば、逆投影を施せる、前処理された発散ビームシノグラムから電子画像を生成する方法が提供され、上記シノグラムは発散ビーム投影の集合である。この方法は、上記シノグラムを複数の部分シノグラムに細分するステップと、全体座標系内で上記部分シノグラムの重み付けされた逆投影を行って、全体座標系内の正しい位置に複数の対応する部分画像を生成するステップと、上記部分画像を集合させて上記電子画像を生成するステップと、を含んでいる。上記シノグラムの部分シノグラムへの細分は、各部分シノグラムが所望サイズの部分画像を示すまで反復手法で行われ、この細分ステップは、所望数の正確な細分と、所望数の近似の細分とを含んでいる。
【００１６】
本発明の他のアスペクトによれば、電子画像を再投影する方法、すなわち該画像から発散ビームシノグラムを生成する方法が提供される。この方法は、
上記画像を複数の部分画像に分割するステップと、全体座標系内の各部分画像の部分シノグラムを計算するステップと、上記部分シノグラムを集合させて上記シノグラムを生成するステップと、を含んでいる。上記画像の細分は、各部分画像が所望のサイズを有するまで反復手法で行われ、上記部分シノグラムの計算は、上記反復における所望数のレベルにおいては近似であり、上記反復の残りのレベルにおいては正確である。
【００１７】
これらの方法は、従来の逆投影アルゴリズムや再投影アルゴリズムと比べて、ＦＦＴと同じ程のスピードアップを提供する。認識可能なほどの劣化や重大な数値上の誤りなしに、実質的な計算上の節約が得られる。
【００１８】
添付の図面を参照しつつ、本発明の実施形態についての以下の記載を参照することにより、本発明における上述した構成及びその他の構成と、それらを得るための方法が、一層明らかになり、そして、本発明が最も良く理解されることになろう。
【発明を実施するための最良の形態】
【００１９】
本発明は、ＣＴスキャナを含む様々なイメージング装置への応用を有している。典型的なイメージング装置１０（図１）は、頭のようなオブジェクトからデータを取得して、発散ビーム投影（ＤＢ投影）１４に相当する生（raw）データを投影プリプロセッサ（前処理装置）１６へ送るスキャナ１２を含んでいる。投影プリプロセッサ１６は、重み付けや、シフト変更であるフィルタリングだけでなく、データに対して様々な変換、正規化、及び訂正を加える。投影プリプロセッサ１６の出力は、発散ビームシノグラム（ＤＢシノグラム１）１８であり、これはシノグラム更新プロセッサ２０に与えられる。このシノグラム更新プロセッサ２０は、発散ビームシノグラム（ＤＢシノグラム２）３４からの情報を用いて、入力されたシノグラム１８を変更可能であり、これは例えば、ビーム硬化（beam-hardening）を含む様々なアーチファクトのための訂正を行いつつ、又はマルチステップの中の一部として、又は反復する再構成手順として行われる。
【００２０】
シノグラム更新プロセッサ２０の出力は、発散ビームシノグラム（ＤＢシノグラム３）２２であり、これは高速逆投影プロセッサ２４への入力となる。この高速逆投影プロセッサ２４は、一般に、コンピュータ又は特殊目的のハードウェア、又はそれらの何らかの組み合わせであり、この中で記述されている正確な又は近似の逆投影アルゴリズムを行うようプログラムされた何らかの適切なタイプのものである。
【００２１】
高速逆投影プロセッサ２４の出力は、電子画像（電子画像１）２６であり、これは画像コンディショニングプロセッサ２８に入力される。この画像コンディショニングプロセッサ２８は、電子画像の必要な後処理を行うものであり、アーチファクト画像や、マルチステップ又は反復再構成プロセスにおける更なる処理のための画像を確認及び抽出することを含んでいてもよい。
【００２２】
必要ならば、画像コンディショニングプロセッサ２８は、高速再投影プロセッサ３２に与えられる電子画像（電子画像２）３０を生成可能である。高速再投影プロセッサ３２は、一般に、コンピュータ又は特殊目的のハードウェア、又はそれらの何らかの組み合わせであり、この中で記述されている正確な又は近似の再投影アルゴリズムを行うようプログラムされた何らかの適切なタイプのものである。必要ならば、このプロセッサは、逆投影プロセッサ２４に採用されたのと同じコンピュータ及びハードウェアの使用を共有することも可能である。
【００２３】
高速再投影プロセッサ３２の出力は、発散ビームシノグラム（ＤＢシノグラム２）３４であり、これはシノグラム更新プロセッサ２０にフィードバックされる。逆投影／再投影のプロセスは、適切な結果が得られるまで継続可能である。再投影は必ずしも必要とされないが、多くの状況において役立つものである。
【００２４】
電子画像２６が適切である場合、画像コンディショニングプロセッサ２８が電子画像（電子画像３）３６を生成し、これは格納／解析／表示装置３８に与えられる。これにより、電子画像３６は、コンピュータメモリに格納可能であり、及び／又は、例えば異常や危険物質のために電子的に解析可能であり、及び／又は、表示可能であり、及び／又は、幾つかの視覚可能な形態で印刷可能である。
【００２５】
本発明の高速発散ビーム逆投影アルゴリズムは、米国特許第６，２８７，２５７号に記載の平行ビーム断層撮影用の高速階層逆投影アルゴリズム（fast hierarchical backprojection algorithm；ＦＨＢＰ）に基づいている。同様に、本発明の高速発散ビーム再投影アルゴリズムは、米国特許第６，２６３，０９６号及び同第６，３５１，５４８号に記載の平行ビーム断層撮影用の高速階層再投影アルゴリズム（fast hierarchical reprojection algorithm；ＦＨＲＰ）に基づいている。逆投影と再投影は、密接に関連しており、それらの高速アルゴリズムは同一原理から得られる。より一般的に使用される高速逆投影は、より詳細に記述することとする。
【００２６】
平行ビームＦＨＢＰを理解するために、画像の中心に原点を有する全体座標系の付随したＮ×Ｎ画像について考える。この大きなＮ×Ｎ画像を、部分画像と呼ばれるＮ／２×Ｎ／２のサイズの小さな４つの画像に分解することを考える。部分画像の各々は、全体座標系の別々の象限に配置される。ラドン変換（Radon transform）の２つの特性が、ＦＨＢＰアルゴリズムの派生に使用される。（i）ボウ−タイ特性（bow-tie property）については、全体座標系の原点に中心を置いた半分サイズの画像にとって、角度変数に関する投影のスペクトル支持体（spectral support）も半分だけ減少される、と言われている。これは、中心に置かれた半分サイズの画像は投影の数の半分を備えるシノグラムから再構成されることが可能であることを意味している。（ii）シフト特性（shift property）については、シフトされた画像の平行投影が原画像のシフトされた平行投影に相当する、と言われている。
【００２７】
今、Ｐ個のフィルタされた投影を備えるシノグラムが、全体画像ｆを再構成するのに利用可能であるとする。部分画像の１つｆｉを再構成するために、Ｐ個の投影がまず部分画像の投影の支持体にまでトランケートされ（truncated）、次に部分画像の中心に置かれたバージョンに相当するようシフトされ（shifted）、最後に角度方向に関してＰ／２個の投影にまでデシメートされる（decimated）。部分画像の中心に置かれたバージョンｆｉが、続いてこれら新たなＰ／２個の投影から再構成され、この部分画像が全体座標系中の正しい位置に戻るようシフトされる。これは、４つの部分画像のそれぞれについて行われ、この４つの部分画像が一緒に集合されると全体画像ｆの再構成が提供される。合計で４ｃ(Ｐ／２)(Ｎ／２)２＝ｃＰＮ２／２の演算が、４つの部分画像を再構成するのに必要であり、ここでｃは定数である。これは、元々の再構成コストを１／２だけ減少させることになる。反復して分解することで、その各ステップで画像サイズを１／２だけ減少させ、１／２の数の投影を用いて半分のサイズの部分画像の形成を許容するので、計算コストの合計をＯ(Ｎ２logＮ)にまで低減させることが可能である。
【００２８】
この分解及び集合のアプローチの双対（dual）は、米国特許第６，２６３，０９６号及び第６，３５１，５４８号に記載の高速階層平行ビーム再投影アルゴリズム（ＦＨＲＰ）で使用されている。
【００２９】
発散ビームの逆投影及び再投影のための本発明は、画像を階層的に分解し、次に、より少ない発散ビーム投影を備える発散ビーム部分シノグラムを用いて部分画像を再構成する（又は再投影する）、という考えに基づいている。しかし、平行ビームＦＨＢＰ又はＦＨＲＰの考えを発散ビームに適用しようとする場合に生じる直接の困難性は、上述のシフト特性（ii）がもはや当てはまらないということである。すなわち、シフトされた部分画像の発散ビーム投影は、元の部分画像のシフトされた発散ビーム投影に相当しないということである。平行ビームＦＨＢＰ及びＦＨＲＰはこの特性を使用するので、平行ビームアプローチは発散ビーム用に変更しなければならない。
【００３０】
本発明のアルゴリズムは、分解用の異なる手順を使用するものであり、この手順は、オブジェクトに固定された、すなわちフルサイズ画像に固定された全体座標系の中心に対する部分画像のシフティングを含まない。平行ビームＦＨＢＰと同様、新規な発散ビーム逆投影アルゴリズムは、２つのタイプの分解、すなわち１つは正確だが低速のもの、もう１つは近似だが高速のものを有している。この新規なアルゴリズムの正確な発散ビーム分解は、データのシフティングを全く必要とせず、部分画像に相当する部分シノグラムへのシノグラムのトランケーションを必要とするのみである。これらの部分シノグラムは、次に、全体座標系内に逆投影されることで、対応する部分画像を全体座標系内の正しい位置に生成する。近似の発散ビーム分解は、平行ビームＦＨＢＰに使用されるものとは異なる手順を使用することで、より小さなサイズの部分画像を再構成するのに使用される投影の数を減少させる。この新規な手順は、投影のみのトランケーション及びシフティングを含み、部分画像のものは含まない。この新規な手順は、ＦＨＢＰには使用されていなかった重み付けを含んでいてもよい。同様なことが、発散ビーム再投影のための新規な手順に適用される。
【００３１】
半分サイズの部分画像を再構成するために、適切にトランケート（重み付け及びフィルタ）された発散ビーム投影が、全体座標系における部分画像の中心の投影位置によって決定される距離だけシフトされ、次に、ソース位置パラメータにおいて２だけデシメートされ、そして、同一量だけシフトバックされる。この部分画像は、続いて、この減少された投影群から良好な近似へと再構成されることが可能である。部分画像に対する正確な重み付け逆投影と近似の重み付け逆投影の両方において、重み付けを含む、画像座標を必要とする全ての計算は、元々のフルサイズ画像を意味する絶対画像座標又は全体画像座標の中でなされる。正確な分解と近似の分解とを反復形式で組み合わせることで、発散ビームジオメトリのための高速階層逆投影アルゴリズムが得られる。その平行ビームＦＨＢＰにおける場合と同様、これらのアルゴリズムは、近似の分解に対する正確な分解の数を選択し、かつ、各種のアルゴリズムパラメータを調整することにより、計算と正確度（精度）との間のトレードオフを提供する。
〔発散ビームイメージングジオメトリ〕
３Ｄ断層撮影のための一般的な発散ビーム捕捉ジオメトリは、画像化されるオブジェクトの周囲の３Ｄ空間内の軌跡上を移動する発散放射線の頂点（通常はソース）と、検出器面とを含んでおり、この面上において、画像化されるオブジェクトを通るソース放射線に沿った線積分が計測される。発散ビームジオメトリではなく収束ビームジオメトリを有する、放射断層撮影法（ＳＰＥＣＴ及びＰＥＴ）に見られるような他のイメージング法は、このジオメトリに変換可能である。
【００３２】
上記の捕捉ジオメトリは、通常、検出器の形状、それに関連する放射線サンプリングスキーム、及びソースの軌跡又は軌道によって分類される。２つの最も一般的な検出器の種類は、（i）等間隔の検出器を有する平面状検出器面（planar detector surface）と、（ii）方位角が等間隔にサンプルされ、かつ高さが直線的等間隔にサンプルされる放射線を有する円柱状検出器面（cylindrical detector surface）である。面検出器を有する３Ｄイメージングにおいて、発散ビームジオメトリはコーンビームジオメトリとしても知られている。
【００３３】
発散ビーム２Ｄ捕捉の特別な場合において、検出器面は、単一線状又は弧状の検出器に減少され、それらはそれぞれ、いわゆる同一直線等間隔ファンビームジオメトリと等角ファンビームジオメトリ（colinear equispaced, and equiangular fan-beam geometries）となる。これらのファンビームジオメトリは、現在市販されているシステムにおいて広く配備されている。
【００３４】
提案されている各種ソース軌道の中で最も一般的なのは、単一円形軌道と、ヘリカル（又はスパイラル）軌道である。しかし、２つの円や、円と線からなるような他の軌道も、実際的な関心を有するものであり、これら各種の軌道を任意に組み合わせたバージョンも実際には生じ得る。
【００３５】
本発明のプロセス、アルゴリズム、及び装置は、一般的な発散ビームジオメトリに適用可能であり、また、先に述べた特別なジオメトリのうちのいずれかや、発散ビームジオメトリの他の例に専門に使用されることも可能である。同様に、本発明の方法は、収束ビームジオメトリを含む（これには限定されない）何らかの数学的に同等なデータ群にも適用可能である。前記の目的のために、この方法は、以下の場合について少々詳細に説明する。
１．単一円のソース軌道
（ａ）同一直線等間隔ファンビームジオメトリ
（ｂ）等角ファンビームジオメトリ
（ｃ）平面コーンビームジオメトリ
２．平面検出器を有する一般のコーンビームジオメトリ
３．ヘリカルソース軌道
一般のプロセスがこれらの場合に対してどのように専門化されるのかを、実証する。このことから、或る１つの検出器形状及びソース軌道のためのアルゴリズムはそれとは異なる検出器形状及びソース軌道のために容易に変更されるものであるということが明らかになる。
〔任意の軌道を有する発散ビームイメージングジオメトリ〕
平面等間隔コーンビームジオメトリは、検出器が１つの平表面上に均一に配置されるものであり、これは図２に示されている。この場合、検出器平面は、各ソース位置毎に任意の位置及び方向を有することが可能となっている。
【００３６】
発散放射線のソースＳは、３Ｄオブジェクトの周囲の軌道上を進行し、オブジェクトに固定された全体座標系の中心Ｏに関するソース位置ベクトル
2005518892.tif 000003
【００３７】
によって定義される。このベクトルは、パラメータθによって特定される軌道上のソース位置を有しており、θｍｉｎ<θ<θｍａｘである。原点Ｏθが検出器平面上への軌道点
2005518892.tif 000004
【００３８】
の正投影（orthogonal projection）となるよう、各θ毎に、検出器に付随した３Ｄ座標系が定義される。正規直交ベクトル
2005518892.tif 000005
【００３９】
が検出器平面内で選択され、ここで、
2005518892.tif 000006
【００４０】
は検出器平面に直交する単位ベクトルであり、Ｏθからソースを向いている。検出器平面内の点
2005518892.tif 000007
【００４１】
が、その座標（ｔ１，ｔ２）を用いて定義される。この系において、ソース軌道点は
2005518892.tif 000008
【００４２】
と表現され、ここでＤ（θ）はソースと検出器との間の距離である。（ＳＰＥＣＴでは、Ｄ（θ）はコーンビームコリメータの焦点距離である。）
2005518892.tif 000009
【００４３】
は、３Ｄ画像（すなわちオブジェクト）内の位置を示すものとし、
2005518892.tif 000010
【００４４】
は、点
2005518892.tif 000011
【００４５】
（図２を参照）を通過するソース放射線が検出器平面と交差するｔ１，ｔ２位置を示すものとする。ソース軌道位置θ及び検出器位置（ｔ１，ｔ２）でのオブジェクトｆの発散ビーム投影（Ｐｆ）(θ，ｔ１，ｔ２)は、（θ，ｔ１，ｔ２）によってパラメータ化されたソース放射線に沿った線積分であり、
【００４６】
【数１】
2005518892.tif 000012
【００４７】
によって与えられる。ここで、
2005518892.tif 000013
【００４８】
はディラックのデルタ関数である。ソース放射線がソースから発散するのではなくソースへ向けて収束するような状況を記述するのにも、同じ捕捉ジオメトリを使用可能である。従って、本発明のプロセスは、収束ビームジオメトリにも等しく適用される。
【００４９】
投影は、ソース軌道に沿ったパラメータによって定義されるＰ個の別々のソース位置θｐ＝ｐΔθ，ｐ＝０，・・・Ｐ−１、で得られるが、通常、必ずしも均一な間隔Δθ＝（θｍａｘ−θｍｉｎ）／Ｐである必要はない。（Ｐｆ）(θＰ，・)は、
2005518892.tif 000014
【００５０】
の全ての値において、ソース位置θＰでの発散ビーム投影と呼ばれる。検出器平面は、通常、ｔ１及びｔ２軸上で異なる間隔Ｔ１及びＴ２を有する、均一な直角グリッド上でサンプルされる。
【００５１】
異なるソース位置での発散ビーム投影の集合は、発散ビームシノグラムと呼ばれる。これと同一の用語が、一組の前処理された発散ビーム投影や、その他の数学的に同等なデータ群においても使用される。シノグラム中の投影の集合は、不完全な場合もあり、単一の投影を備えていてもよい。
【００５２】
一般の発散ビームジオメトリの幾つかの特別な場合について、これから記述する。
〔平面検出器を有するコーンビームヘリカル断層撮影〕
この場合のソース軌道は、図３に示されているように、軸回りの螺旋又は渦巻きであり、
2005518892.tif 000015
【００５３】
で定義される。ここで、半径Ｄは一定であり、ｈは螺旋のピッチとして知られている。
この場合のパラメータθは、ｘ軸に対する（ｘ，ｙ）平面内のソース角度を示している。検出器平面は、一般性を損なうことなしに、ｘ軸を含み、ソースの動きと共に回転すると考えることで、それは各θ毎に（ｘ，ｙ）平面上へのソースベクトル
2005518892.tif 000016
【００５４】
の投影に垂直となる。検出器平面内の原点Ｏθは、この場合、ｚ軸上に落ちる。ｔ２軸はｚ軸と一致し、ｔ１軸は（ｘ，ｙ）平面に平行である。（その他の一般的な選択、例えば、
2005518892.tif 000017
【００５５】
が位置θで軌道の接線に平行であり、
2005518892.tif 000018
【００５６】
がそれに直角であるといったことは、ここで座標回転によって記述される場合まで減少可能である。）
実際には、検出器平面は異って配置されてもよいが、単純な直線座標変換が、実際の検出器平面上で得られる投影を、（ｔ１，ｔ２）座標に沿って均一な間隔を有する上述の仮想の検出器上で得られる投影に変換する。
〔単一円軌道及び平面検出器を有するコーンビームイメージング〕
この場合のソース軌道は、図４に示されているように、ヘリカル軌道の特殊な場合であり、ピッチがゼロ ｈ＝０、すなわち、軌道がｘ，ｙ平面内において原点Ｏでの半径Ｄの単一円
2005518892.tif 000019
【００５７】
の上にある。検出器平面は、一般性を損なうことなしに、ｚ軸を含み、ソースから中心への線
2005518892.tif 000020
【００５８】
に垂直であると考える。
オブジェクト内の点
2005518892.tif 000021
【００５９】
について考える。検出器角度θでの、検出器上へのそのコーンビーム投影の位置
2005518892.tif 000022
【００６０】
は、
【００６１】
【数２】
2005518892.tif 000023
【００６２】
によって与えられる。
これらの明白な表現は、後に、高速アルゴリズムの説明において有用である。
〔円形軌道及び同一直線等間隔検出器を有するファンビーム断層撮影〕
この場合は、図５に示されるように、オブジェクトの単一断面がイメージされ、２Ｄ画像が再構成される。ソース軌道はまた、（ｘ，ｙ）平面内の半径Ｄの円であり、原点Ｏに中心が置かれ、ソース角度θで示されるソース位置を有している。しかし、検出器は（ｘ，ｙ）平面内の単一直線（ｔ軸）に限定され、これは一般性を損なうことなしに、ソースの回転中心Ｏを通過し、ソースから中心への線
2005518892.tif 000024
【００６３】
に垂直である。検出器要素はｔ軸上に均一に配置されている。
2005518892.tif 000025
【００６４】
は２Ｄ画像中の位置を示すものとし、
2005518892.tif 000026
【００６５】
はｔ軸と
2005518892.tif 000027
【００６６】
（図５参照）を通過するソース放射線との交点のｔ位置を示すものとする。ソース角度θ及び検出器位置ｔでのオブジェクトｆのファンビーム投影（Ｐｆ）(θ，ｔ)は、（θ，ｔ）によってパラメータ化されたソース放射線に沿った線積分であり、
【００６７】
【数３】
2005518892.tif 000028
【００６８】
によって与えられる。
投影は、均一な間隔Δθ＝θｍａｘ／Ｐを有するＰ個の別々のソース角度θｐ＝ｐΔθ，ｐ＝０，・・・Ｐ−１、で得られる。（Ｐｆ）(θＰ，・)は、ｔの全ての値において、ソース角度θＰでの投影であり、異なるＰでの投影の集合がシノグラムである。最大のソース角度θｍａｘは、フルスキャンの場合には２πに等しいが、ショートスキャンやオーバースキャンの場合には他の値をとることも可能である。部分的なデータの集合は、角度の何らかの部分集合での投影からなっていてもよい。
【００６９】
発散ビームジオメトリのそれぞれの場合において、平面ではなく曲面の検出器面が、変更された放射線サンプリングスキームと共に、しばしば使用されてもよい。例えば、ヘリカル又は単一円形軌道のコーンビームイメージングでは、検出器面はしばしば、ソース上に中心の置かれた円柱の表面となるよう選択され、放射線サンプリングは方位が等角で、高さｚが等間隔となる。これは、図６に単一円形ソース軌道で記載されており、ここでは、曲面検出器がｚ軸を含むように選択されている。
【００７０】
円柱検出器形状と平面検出器形状との関係は、一般に、ファンビームジオメトリを研究することによって得られる。図７に示されているような、等角放射線センプリングを有するファンビームジオメトリにおいて、検出器ｄは、ソースＳに中心の置かれた弧上で、角度がΔｔの間隔で等間隔に配置されている。従って、このイメージングジオメトリにおいては、同一直線等間隔ファンビームジオメトリに見られる検出器位置ｔの代わりに、ファン（扇形）角度ｔが使用される。今、図７に示されるように、
2005518892.tif 000029
【００７１】
が点
2005518892.tif 000030
【００７２】
を通過する放射線のファン角度位置を示すものとする。その時、投影は再び式（４）によって与えられる。
同様に、本発明の適用される、発散ビームジオメトリの特別な場合であるファンビームジオメトリの他の例が存在するということを考えてみる。
〔発散ビームのフィルタされた逆投影再構成及び再投影〕
発散ビームデータの本来の（直接の）逆転は、近似であっても正確であっても、重み付けされフィルタされた逆投影（weighted filtered backprojection）として明確化可能である。先ず、発散ビームシノグラムが前処理される。例えば、シノグラムを構成する投影が重み付け及びフィルタされて、ソース位置θＰに相当する投影
2005518892.tif 000031
【００７３】
を構成する、変更された発散ビームシノグラムを生成する。このフィルタリングは、ファンビーム又は近似のコーンビーム再構成の場合のように、単純なシフト不変ランプフィルタリング（shift-invariant ramp filtering）であってもよく、又は、正確なコーンビーム再構成の場合のように、もっと複雑なシフト変化フィルタリング（shift varying filtering）であってもよい。重み付け及びフィルタリングは、部分スキャンされた長いオブジェクトを有するヘリカルコーンビームデータからの再構成のための正確又はほぼ正確な方法の場合のように、２以上の投影群を生成してそれらの分担が加算される、もっと一般的なものであってもよい。簡潔性及び一般性のために、逆投影されるデータは、それを元々の発散ビーム投影データから得るのに使用される前処理がどんなものであろうと、「発散ビーム投影データ」又は「発散ビームシノグラム」と呼ぶこととする。固定されたｐでは、それは（ソース位置θｐでの）「発散ビーム投影」と呼ばれ、関数ｇ（ｐ，・）で示される。
【００７４】
その場合、３Ｄ画像は、従来の発散ビームの離散θでの逆投影（divergent-beam discrete-θ backprojection）
【００７５】
【数４】
2005518892.tif 000032
【００７６】
によって再構成される。ここで、
2005518892.tif 000033
【００７７】
は適当な重み関数である。この離散（discrete）逆投影式は、全てのθで計測される投影を使用する逆投影のための積分表現を近似する。離散検出器を使用するために、実際には、
2005518892.tif 000034
【００７８】
はまた
2005518892.tif 000035
【００７９】
でサンプルされる。これは、
2005518892.tif 000036
【００８０】
が通常は利用可能なサンプル位置に相当しないために、逆投影を実行するには、変数
2005518892.tif 000037
【００８１】
でｇの補間を必要とする。２Ｄファンビームの場合における逆投影式は、τがスカラーで
2005518892.tif 000038
【００８２】
が２次元であることを除き、同一の形式を有する。
Ｎ×Ｎ×Ｎ画像での３Ｄコーンビーム逆投影の計算コストは、ｃＮ３Ｐである。ここで、定数ｃは補間の複雑さのような実行詳細に依存する。これに対して、重み付け及びランプフィルタリングの計算コストは、フィルタリングがシフト不変であって、コンボルーション（convolution）がＦＦＴを用いて行われる場合、ほんのＯ(ＰＮ２logＮ)にすぎない。フィルタリングのもっと手の込んだ形式も、しばしば同様なコストで行うことが可能である。従って、逆投影のコストは、多くの場合のようにＰ＝Ｏ（Ｎ）である時においてＯ(ＰＮ３)又はＯ(ＰＮ４)のコストを有する従来のコーンビーム再構成のコストを卓越している。この状況は、フィルタリング及び逆投影のステップの複雑さがそれぞれＯ(Ｎ２logＮ)及びＯ(Ｎ３)である２Ｄファンビーム再構成においても、同様である。
【００８３】
電子的に格納された画像ｆの本来の（直接の）発散ビーム再投影は、式（１）の実行である。実際には、ｆは離散形式でのみ利用可能であり、また、
2005518892.tif 000039
【００８４】
の
2005518892.tif 000040
【００８５】
でのサンプルのみが計算される。例えば、ｆは次の級数
【００８６】
【数５】
2005518892.tif 000041
【００８７】
によって表すことが可能である。ここで、係数βｉｊｋはｆの離散表現であり、
2005518892.tif 000042
【００８８】
は、スプラインのテンソル積のような局所化された基本関数、又は標準ピクセル基本関数、すなわち第ijk番目（ijk-th）のボクセルでの指示関数であり、これらは、２Ｄについては例えば米国特許第６，２６３，０９６号及び第６，２８７，２５７号に記載されており、３Ｄについては米国特許第６，３３２，０３５号に記載されている。その場合の投影は、
【００８９】
【数６】
2005518892.tif 000043
【００９０】
によって計算される。ここで、
2005518892.tif 000044
【００９１】
は基本関数
2005518892.tif 000045
【００９２】
の発散ビーム投影であり、式（１）によって定義される。これら基本関数の投影は、通常、解析表現を用いて計算されるか、又はルックアップテーブルに格納されている。逆投影の場合と同様、本来の再投影のコストは、３ＤではＯ(Ｎ４)であり、２ＤではＯ(Ｎ３)である。
【００９３】
ｔ１，ｔ２が検出器上の適当な座標である曲面検出器の場合、今、
2005518892.tif 000046
【００９４】
が点
2005518892.tif 000047
【００９５】
を通過するソース放射線の検出器平面上の（ｔ１，ｔ２）位置を示すものとする。投影は、再び式（１）によって与えられ、曲面検出器上の基本関数
2005518892.tif 000048
【００９６】
の発散ビーム投影を示す
2005518892.tif 000049
【００９７】
と共に、式（７）を使って計算可能である。同様に、再構成は、適当に変更された投影
2005518892.tif 000050
【００９８】
の重み付け逆投影として再び表すことが可能であり、ここで、ｔ１，ｔ２は曲面検出器上の適当な座標である。逆投影は前記と同一の式（５）によって与えられるが、異なる重み付け関数
2005518892.tif 000051
【００９９】
を伴っている。例えば、２Ｄ等角ファンビームジオメトリにおいては、等間隔同一直線ジオメトリの線検出器上の変位ｔの代わりに、ファン角度ｔが使用される。
2005518892.tif 000052
【０１００】
の離散性に由来する考察は、平面検出器アルゴリズムにおけるものと同じであり、計算の複雑さについての考察である。
既に論じたように、２Ｄであろうが３Ｄであろうが、逆投影ステップ（及び、もし使用されるのであれば再投影）は、通常、計算上は最も高くつくものであり、従って、それは本発明の目標でもある。あらゆる発散ビームジオメトリに応用可能な、高速逆投影の新規な一般プロセスについて記述する。例を単純化するために、まず、２Ｄファンビームジオメトリのための高速逆投影及び再投影について述べる。
〔同一直線等間隔ファンビームジオメトリのための高速の本来の逆投影及び再投影アルゴリズム〕
図８に示されるような、ｆの部分画像
2005518892.tif 000053
【０１０１】
のための逆投影演算について考える。
2005518892.tif 000054
【０１０２】
が、
【０１０３】
【数７】
2005518892.tif 000055
【０１０４】
によって定義される画像トランケーション演算子であるとする。ここで、
2005518892.tif 000056
【０１０５】
である。従って、
2005518892.tif 000057
【０１０６】
が、ｆの部分画像であり、そのサイズはＭ×Ｍで、
2005518892.tif 000058
【０１０７】
に中心が置かれている。ソース角度ｐΔθ，ｐ＝０・・・Ｑ−１、でのＱ個の投影を用いた、全体座標系内での（全体座標系におけるその正しい位置に配置された）部分画像ｆ′上への逆投影は、
【０１０８】
【数８】
2005518892.tif 000059
【０１０９】
によって与えられる。
2005518892.tif 000060
【０１１０】
は、関係する逆投影演算子を示している。従って、特に、
2005518892.tif 000061
【０１１１】
である。
逆投影の局所性のために、投影ｇ（ｐ，・）の一部分のみがｆ′上への逆投影に寄与する。この部分が
2005518892.tif 000062
【０１１２】
によって示される。ここで、
2005518892.tif 000063
【０１１３】
は、各θｐ毎に変数ｔでのｇ（ｐ，ｔ）を部分画像
2005518892.tif 000064
【０１１４】
（図８参照）の支持体（support）の投影の支持体
2005518892.tif 000065
【０１１５】
にトランケートする演算子である。
2005518892.tif 000066
【０１１６】
を決定するトランケーション間隔は、全てのソース角度θｐと、対象の部分画像のサイズＭ及び位置
2005518892.tif 000067
【０１１７】
とで、前もって計算可能である。この代わりに、より単純で控えめな近似を使用可能であるが、シフティングとデシメーションを必要とする投影の必要な支持体を過剰評価することになるため、階層アルゴリズムにおいて計算コストが幾分上昇する。
【０１１８】
逆投影の局所化から、もし
2005518892.tif 000068
【０１１９】
である場合には
【０１２０】
【数９】
2005518892.tif 000069
【０１２１】
となる。すなわち、ｆ′上への逆投影は、近似的に適切にトランケートされた投影のサイズ（Ｍ，Ｐ）の逆投影によって得ることができる。
今、画像ｆを、それぞれのサイズがＭ×ＭのオーバーラップなしのＪ＝（Ｎ／Ｍ）２個の部分画像へと区切ることを考える。
【０１２２】
【数１０】
2005518892.tif 000070
【０１２３】
ここで、ベクトル
2005518892.tif 000071
【０１２４】
は、図９に示されるように、全体座標系における部分画像の中心である。式（１０）を適用すると、以下の、逆投影のための画像上への、逆投影への正確な分解が得られる。
【０１２５】
【数１１】
2005518892.tif 000072
【０１２６】
この分解は、Ｎ／Ｍ＝２におけるものが図１０に示されている。式（１２）中のｇによっても示される発散ビームシノグラム２２は、今述べた部分画像に対応する複数の部分シノグラムに分割される。この部分シノグラムは、１ピクセル又は１ボクセルと同じ程度に小さな、何らかの所望のサイズであってもよい。トランケーション演算子１０１０、１０１２、１０１４及び１０１６によって生成される部分シノグラムは、個々に、それぞれ１０１８、１０２０、１０２２及び１０２４において全体座標系内に逆投影される。逆投影されたシノグラムは部分画像ｆ１、ｆ２、ｆ３及びｆ４を生成し、これらは１０２６で集合されて、図１中の電子画像２６を生成する。このプロセスは、ファンビーム投影の正確な分解を生成する。Ｎ／Ｍ＝２では、図１０に示されるように、４つの部分シノグラム及び部分画像へと分解される。
【０１２７】
図１０を再び参照すると、線２２の下のＰは、シノグラム又は部分シノグラム中の投影の数を示す。この分解は、それ自身では、単一ステップの逆投影ＢＮ，Ｐと比べてスピードアップを提供するものではない。事実、単一の部分画像逆投影ＢＭ，Ｐの計算コストｃＭ２Ｐは、全画像での単一ステップＢＮ，Ｐのものよりも（Ｎ／Ｍ）２＝Ｊ倍だけ小さいのに対し、必要とされるＪ個の部分画像逆投影における合計コストは、依然として同じである（簿記のための僅かな追加の諸経費を伴う可能性もある）。正確な分解は、Ｊ個の並列プロセッサ間で逆投影を区切ることによるスピードアップを提供したり、或いは、Ｊ倍小さなメモリを必要とする単一プロセッサの繰り返し使用を許容するのに、まだ使用可能である。
【０１２８】
高速逆投影アルゴリズムは、以下の追加の特性を伴う分解の考え方を使用する。固定された画像分解能（帯域幅）では、画像を再構成するのに必要とされる投影の数が、画像のサイズに比例する。すなわち、Ｎ×Ｎ画像を再構成するのにＰ個の投影が必要とされる場合、それと同一分解能でＭ×Ｍ画像を再構成するのにはＰ′＝（Ｍ／Ｎ）Ｐ個の投影で十分である。
【０１２９】
従って、高速逆投影アルゴリズムは、減少された数の投影からの部分画像
2005518892.tif 000073
【０１３０】
の再構成に基づいている。この減少プロセスは、演算子
2005518892.tif 000074
【０１３１】
によって行われ、この演算子自体は、今から定義する幾つかの演算子から構成されている。
投影シフト演算子
2005518892.tif 000075
【０１３２】
及び
2005518892.tif 000076
【０１３３】
は、それぞれｔ軸に沿って
2005518892.tif 000077
【０１３４】
だけ投影をシフトさせるものであり、
【０１３５】
【数１２】
2005518892.tif 000078
【０１３６】
として定義される。
Ｌ倍（Ｌ-fold）角度デシメーション演算子
2005518892.tif 000079
【０１３７】
は、角度フィルタリングの後に角度サブサンプリングを行うことにより、１組のＰ個の投影における視界（view）の数をＰ／Ｌ個の投影に減少させる。このフィルタリングは、隣接する角度で投影を平均化するのと同じくらいに単純であってもよく、或いは、もっと手の込んだものであってもよく、そしてこれは、
【０１３８】
【数１３】
2005518892.tif 000080
【０１３９】
と表される。ここで、「＊」はコンボルーション（convolution）を示し、ｈＤはフィルタ零空間（filter kernel）である。このフィルタリングは、それぞれの固定されたｔにおけるｐ（即ち、離散したθ）にのみ発生してもよく、或いは、θ及びｔでの組み合わされたフィルタリングであってもよく、それは分離可能でもそうでなくてもよい。角度サブサンプリングは、全Ｌ個のフィルタされた投影
2005518892.tif 000081
【０１４０】
の中から１つを保持する。
演算子
2005518892.tif 000082
【０１４１】
は、トランケーション演算子、シフト演算子、及び角度デシメーション演算子の合成であり、
【０１４２】
【数１４】
2005518892.tif 000083
【０１４３】
によって定義される。ここで、Ｐ′＝Ｐ／Ｌである。演算子
2005518892.tif 000084
【０１４４】
は、デシメーションの前後に、重み付け関数Ｗｋ（θ，ｔ），ｋ＝１，２ による乗法重み付け（multiplicative weighting）を含んでいる。１つの可能な選択は、
2005518892.tif 000085
【０１４５】
及び
2005518892.tif 000086
【０１４６】
である。もっと一般的に言えば、重み付は、アルゴリズムの正確度を最適にするよう選択可能であり、又は全く省略することも可能である。
図１１は、部分シノグラムの処理され得る方法を一層詳細に示している。例えば、トランケーション
2005518892.tif 000087
【０１４７】
がまず初めにステップ１１１０にて起こり、その後にステップ１１１２でシフティング演算子
2005518892.tif 000088
【０１４８】
が続き、そして重み付けステップ１１１４、その後にデシメーションステップ１１１６、もう１つの重み付けステップ１１１８、そして第２のシフティングステップ１１２０となる。このプロセスは、シノグラム分解演算子
2005518892.tif 000089
【０１４９】
を記述している。
これらの定義により、部分画像ｆ′上への逆投影のための正確な式（１０）は、次の近似によって置き換えられる。
【０１５０】
【数１５】
2005518892.tif 000090
【０１５１】
ここで、ＢＭ，Ｐ／Ｌは、式（９）によって定義されるＰ／Ｌ個の投影を用いた、Ｍ×Ｍ部分画像上への全体座標系内での逆投影である。これにより、式（１２）に類似の区切られた画像のための逆投影演算子の近似分解が得られる。
【０１５２】
【数１６】
2005518892.tif 000091
【０１５３】
この分解は図１２に示されており、この図は、Ｌ＝Ｎ／Ｍ＝２での投影数の減少を伴うファンビーム逆投影のための近似分解を示している。この近似分解を実行するために、高速逆投影プロセッサ２４は、図１２中にｇとして示された、Ｐ個の投影からなる発散ビームシノグラムを分割し、１２１０、１２１２、１２１４及び１２１６によって、ｇ１、ｇ２、ｇ３及びｇ４で示されるそれぞれＰ／２個の投影からなる４つの部分シノグラムを生成するように、プログラムされる。部分シノグラムｇ１、ｇ２、ｇ３及びｇ４が、１２１８、１２２０、１２２２及び１２２４によって全体座標系内に逆投影されて、全体座標系内での正しい位置に部分画像ｆ１、ｆ２、ｆ３及びｆ４が生成される。部分画像が１２２６で集合されて、図１２中にｆとしても示されている電子画像２６が生成される。
【０１５４】
再び図１２を参照して、Ｊ＝（Ｎ／Ｍ）２の部分画像のそれぞれへの逆投影のための投影数はｃＭ２Ｐ／Ｌであり、全ての逆投影では合計ｃＮ２Ｐ／Ｌである。これは、単一ステップの逆投影ＢＮ，ＰにおけるよりもＬ倍少ない。画像サイズでの投影の数の尺度化についての先の観測に基づけば、ＬはＬ＝（Ｎ／Ｍ）と同じ大きさに選択可能である。
2005518892.tif 000092
【０１５５】
の計算コストが含まれる場合、部分画像の最適な数は、米国特許第６，２６３，０９６号に示されるような単一レベル分解のために決定可能である。
高速逆投影アルゴリズムの好ましい実施の形態においては、式（１２）及び（１８）中の分解は反復して適用され、各レベルで式（１１）に従い部分画像をより小さな部分画像へと更に区切る。画像のダイアディック（Ｎ／Ｍ＝２）２レベル区切りの一例が図１３に示されている。異なるレベルでの幾つかの部分画像の中心位置も示されている。それらは全て、ゼロ番目のレベルの画像原点Ｏ、すなわち全体座標系に言及している。それに対応する、２レベル近似分解の一例と、１つの正確な分解及び１つの近似分解の段階からなる２レベル分解の一例が、それぞれ図１４及び１５に示されている。図１４において、Ｐ個の投影を備える発散ビーム部分シノグラム２２は、高速逆投影プロセッサ２４で処理されて、それぞれがＰ／２個の投影を備える部分シノグラム１４１０、１４１２、１４１４及び１４１６に分割される。この近似の細分は、図１１のフローチャートによって記述された近似分解演算子によって行われる。部分シノグラム１４１０の更なる近似分解が、１４１８、１４２０、１４２２及び１４２４で、図１１のフローチャートによって記述される演算子を再び用いて生成されるが、今回は図１３中のより小さな部分画像に対応する異なるパラメータで行われる。Ｐ／４個の投影を備える上記部分画像のそれぞれは、１４２６、１４２８、１４３０及び１４３２で全体座標系内で逆投影されて、図１３に示された全体座標系内における正しい位置に複数の部分画像ｆ１，１、ｆ１，２、ｆ１，３、ｆ１，４が生成される。これら部分画像は１４３４で集合されて、やはり全体座標系内の正しい位置に部分画像ｆ１が生成される。
【０１５６】
１４１６で表される部分シノグラムは、１４３６、１４３８、１４４０及び１４４２で同様に細分される。これらの部分シノグラムは、それぞれ１４４４、１４４６、１４４８及び１４５０で逆投影され、そして１４５２で集合されて部分画像ｆ４が生成される。点で示されたブロック１４５４及び１４５６は、それぞれｆ２及びｆ３のための対応するレベル２の分解を表している。図１４において、分解の全てが、Ｌ＝Ｎ／Ｍ＝２におけるファンビーム投影の近似分解である。
【０１５７】
図１５は、Ｌ＝Ｎ／Ｍ＝２におけるファンビーム投影の２レベル混合（１つは正確な、もう１つは近似の）分解を示している。図１５において、Ｐ個の投影を備えるシノグラム２２は、やはりそれぞれがＰ個の投影を備える正確な部分シノグラム１５１０、１５１２、１５１４及び１５１６に分解される。これら正確な部分シノグラムは、それぞれｇ１、ｇ２、ｇ３及びｇ４で示され、（適当なパラメータと共に）図１１に示された近似分解演算子を用いて、それぞれＰ／２個の投影を備える近似部分シノグラムに分解される。例えば、ｇ１は、１５１８、１５２０、１５２２及び１５２４で分解されて、近似部分シノグラムｇ′１，１、ｇ′１，２、ｇ′１，３及びｇ′１，４が生成される。それらの部分シノグラムは、１５２６、１５２８、１５３０及び１５３２で全体座標系内に逆投影されて、全体座標系内の正しい位置に部分画像ｆ１，１、ｆ１，２、ｆ１，３及びｆ１，４が生成される。それらの部分画像は、１５３４で集合されて部分画像ｆ１が生成される。
【０１５８】
１５１６で生成された部分シノグラムは、図１５中でｇ４で示されており、これはまた更に、１５３６、１５３８、１５４０及び１５４２で４つの部分シノグラムに分解されて、部分シノグラムｇ′４，１、ｇ′４，２、ｇ′４，３及びｇ′４，４が生成される。１５４４、１５４６、１５４８及び１５５０での逆投影は部分画像ｆ４，１、ｆ４，２、ｆ４，３及びｆ４，４を生成し、これらが１５５２で集合されて部分画像ｆ４が生成される。点で示されたブロック１５５４及び１５５６は、それぞれｇ２及びｇ３のための対応する２レベル混合分解を表している。それらのプロセスは部分画像ｆ２及びｆ３
を生成し、部分画像ｆ１、ｆ２、ｆ３及びｆ４が１５５８で集合されて、図１５中にｆで示される画像２６が生成される。
【０１５９】
部分画像が或る所望の最小サイズＭｍｉｎとなるまで反復を継続可能であり、その後、逆投影ＢＰｍｉｎ，Ｍｍｉｎが式（５）を用いて全体座標系内で行われる。最適なＭｍｉｎは、通常、インプリメンテーションに依存する（implementation-dependent）。それは、最小の部分画像が１ピクセル又は１ボクセルのサイズとなるよう、Ｍｍｉｎ＝１と同じくらい小さく選択することが可能である。Ｎは２の累乗であり、分解の全てのレベルでＮ／Ｍ＝２であるとすると、ｌｏｇ２Ｎのレベルが存在することになる。更にＬ＝（Ｎ／Ｍ）＝２であり、かつＰ＝ｋＮであり、ｋを小さな整数であるとすると、最後のレベルはＮ２個の単一ピクセル逆投影ＢＰ／Ｎ，１を含むこととなり、合計コストがｃＰＮになる。
2005518892.tif 000093
【０１６０】
の実行におけるデシメーション及びシフトの演算子のために固定長の補間が使用されるものとすると、各レベルでの計算コストがｃ１ＮＰ（ここでｃ１は定数）となることが示される。従って、ｌｏｇ２Ｎのレベルの合計コストはＯ(ＰＮｌｏｇＮ)であり、そしてＰ＝Ｏ（Ｎ）なので、反復分解アルゴリズムの合計コストはＯ(Ｎ２ｌｏｇＮ)になる。
【０１６１】
実際には、離散した検出器を使用するので、投影データは間隔Δｔを有するｔでサンプルされる。従って、投影シフト演算子
2005518892.tif 000094
【０１６２】
及び
2005518892.tif 000095
【０１６３】
は、離散形式で実行されることになるが、Δｔの非整数倍によるシフトが必要となる。投影は逆投影の前にフィルタされているので、それらは帯域制限され、そして、例えば米国特許第６，２８７，２５７号に記載されているように、補間及び再サンプリングの組み合わせを用いてシフトされることが可能である。特に、Δｔの非整数倍によるシフティングは、Δｔの分数によるシフティングを伴って段階的に行われる、Δｔの整数倍によるシフティングへと分離されることが可能である。分数シフトはｔ変数内の補間によって行われるが、整数シフトは単純な再インデクシングに相当する。更に、この補間は、式（１４）内の一般化されたデシメーションフィルタリングステップにおけるフィルタｈＤ（ｐ，ｔ）を有する放射次元内でのコンボルーション（convolution）に対応する放射フィルタリングによって実行可能である。
【０１６４】
重要なことには、ここで提案される特定の投影シフティングスキームによれば、投影のサンプリング間隔が均一かつ一定のままとなり、その結果、補間及び再サンプリングは離散対離散として行われ、単純なディジタルフィルタとして効率的に実行可能である。そのような分数ｔシフティングに含まれるエラーは、補間の長さで指数関数的に衰微させることが可能であり、その結果、通常は短い補間で十分となり、単純な線形補間でさえ十分となる。ｔ補間の正確度はまた、ｔで投影をオーバーサンプリングすることによって向上し、これは必要ならば、投影のディジタル補間（アップサンプリングの後にフィルタリング）によって効率的に行われ得る。
【０１６５】
近似分解を用いて得られる逆投影の正確度は、Ｐ，Ｎ，Ｍ，Ｌ、画像及び投影の帯域幅中の分解能、オブジェクトと交差するソース放射線の最大ファン角度、検出器のレスポンス、及びアルゴリズム中で使用される補間を含む、プロセス中の各種パラメータに依存する。これらのパラメータは、発散ビームアルゴリズム中で適当に最適化されてもよい。正確度に影響を与える重要なパラメータは、逆投影に使用される投影の数と、部分画像のサイズとの間の比（Ｐ／Ｌ）／Ｍである。角度オーバーサンプリング（angular oversampling）とも呼ばれる、この比を増大させることにより、近似分解の正確度が向上する。
【０１６６】
以下の(i)〜(iii)を含む角度オーバーサンプリングを増大するための幾つかの方法がある：(i)分解の開始前のＰを増加させるための、投影ｇのディジタル角度補間、(ii)近似分解の各段階でのＬ＜（Ｎ／Ｍ）の選択、及び(iii)例えば図１５中に例として示されているような正確な分解ステップの使用。方法(iii)を理解するために、式（１２）における正確な分解は部分画像のサイズを減少させるが、部分シノグラム中の投影の数を減少させることはない、ということを思い出すこと。従って、これにより、因子（Ｎ／Ｍ）による角度オーバーサンプリングの増大と、改善された正確度とが得られる。
【０１６７】
事実、上記方法の中の１つである(ii)又は(iii)に他と比べての小さな計算上の利点を与える実行詳細は別として、
2005518892.tif 000096
【０１６８】
が
2005518892.tif 000097
【０１６９】
に中心の置かれた部分画像の部分画像の中心である時にはいつでも成立する恒等式
2005518892.tif 000098
【０１７０】
のために、角度オーバーサンプリングのための方法(iii)は方法(ii)に等しい。従って、ダイアディック区切りを有する２段階の分解、-１つは正確、もう１つは近似- は、可能な最大値Ｌ＝４の代わりにデシメーションＬ＝２での、それぞれサイズがＮ／４の１６個の部分画像への単一近似分解に等しい。従って、両方法は、同一因子２の角度オーバーサンプリングを提供する。一方、Ｌ＝３を選択することにより、１６個の部分画像への単一レベル近似分解は、因子４／３の角度オーバーサンプリングを提供可能であり、従って、動作点の選択に追加の柔軟性を提供する。一般には、角度オーバーサンプリングを増大させるための方法により、一定因子によるアルゴリズムの計算コストが増加する。従って、それらは計算と正確度との間のトレードオフを調整する方法を提供する。
【０１７１】
本アルゴリズムの好ましい実施の形態は、幾つかの正確な細分、及び幾つかの近似の細分、又は細分因子よりも小さなデシメーション因子を有する近似の細分を含む。正確な細分及び／又はデシメーションの因子の数はオーバーサンプリングの因子を調整するものであり、また、計算と正確度との間の所望のトレードオフが得られるよう選択される。
【０１７２】
最も単純な実施形態は、近似の細分レベルではＬ＝２で、各レベルでＮ／Ｍ＝２の一定因子による細分を使用するものであるが、もちろんその他の選択も可能である。
先に示したように、本アルゴリズムは、全てのレベルにおいてデシメーション因子Ｌの積によるＰ個への分割が必要である。これが、デシメーション因子の選択によって好都合又は効率的に達成され得ない場合は、角度補間及び投影の再サンプリングによってＰを変更することが可能である。
【０１７３】
本アルゴリズムの好ましい実施形態の反復構成を更に示すために、Ｍ／Ｎ＝Ｌ＝２における反復関数として書かれたアルゴリズムのための翻訳言語（pseudo-code）の一例が、図１６に示されている。この翻訳言語は自明ではあるが、簡単に説明する。
【０１７４】
関数ＦＡＮ＿ＦＨＢＰは、全シノグラムＧ上の第１トップレベルの援助コール（invocation call）における反復関数演算子である。それ自身への次のコールにおいて、この関数は画像を部分画像に分解し、最小画像サイズＮｍｉｎに達するまでシノグラム分解を行う。これらの分解のＱは正確であり、その残りは近似である。画像サイズＮｍｉｎに達すると、全体座標系内での正確な逆投影が、部分画像に対応する部分シノグラム上で行われる。Ｎｍｉｎよりも大きな部分画像では、部分画像に相当するよう投影がトランケートされる。これは、正確なシノグラム細分を構成する。正確な分解が完了すると（Ｑ＜１の時）、近似分解のための追加ステップが行われる。図１６において、投影はシフトされ、重みデシメートされ（weight-decimated）、そして再度シフトされる。この部分画像は
2005518892.tif 000099
【０１７５】
で集合される。
部分画像の中心はまた、
2005518892.tif 000100
【０１７６】
として表現される、
2005518892.tif 000101
【０１７７】
に中心の置かれたより大きな部分画像の部分画像ｉの中心
2005518892.tif 000102
【０１７８】
を用いて反復的に計算される。ここで、
2005518892.tif 000103
【０１７９】
は、
2005518892.tif 000104
【０１８０】
に対する部分画像ｆｉの近似の位置ベクトルである。フルサイズ画像の中心は、定義により、全体座標系の原点にある。従って、その後の全ての部分画像の中心は全体座標系において表現される。
【０１８１】
高速逆投影アルゴリズムの好ましい実施の形態においては、反復分解が以下の観点により単純化される。ｔでの投影シフト、投影トランケーション、及び投影重み付けの演算が、（適当な座標変換を用いて）代わりに使用できる。従って、１つの段階の演算
2005518892.tif 000105
【０１８２】
がその後の段階の
2005518892.tif 000106
【０１８３】
と組み合わされて、１つのシフト及び１つの重み付けの演算になる。同様に、最後の部分画像の逆投影ＢＭｍｉｎ，Ｐｍｉｎの直前に生じる、このシフティング及び重み付けの演算
2005518892.tif 000107
【０１８４】
は、逆投影自身における重み付け及び補間の演算と組み合わせ可能である。従って、反復実行において、
2005518892.tif 000108
【０１８５】
は、１つの投影シフトステップと、各段階毎の（多くとも）１つの重み付けのみが必要となる。
信号処理アルゴリズムの設計及び実行における当業者にとっては明らかなことであろうが、本発明のプロセスを定義する各種演算は、関数又はプロセスの原理を変化させなければ、添付の図によって記載された順番とは異なる各種順番で再配置することや、並列に実行することも可能である。幾つかのそのような再配置は、実行、演算、又はコストの観点から有利である場合がある。例えば、部分画像のための投影の数をデシメートする演算は、Ｐ個の全ての投影がトランケート及びシフトされるのを待つことなしに、２以上の投影がこの部分画像のためにトランケート及びシフトされるとすぐに開始されるようにすることが可能である。そのような再配置は、再構成アルゴリズムを配備したり、データの捕捉でオーバラップしたりするのに使用可能である。
【０１８６】
更に、投影は適時に連続的に得られるので、本発明の方法は連続的に更新される断層撮影画像に適用可能である。これは、最も最近に得られた投影とそれよりも古い投影群との間の相違から構成されるシノグラムを逆投影するための提案された方法を用いることにより達成可能である。その逆投影の結果が、続いて現在の画像に加算されて、更新画像が生成される。そのような更新は、連続的手法で続行可能である。
〔高速再投影〕
逆投影と再投影は、双対（dual）の演算、又は互いの随伴（adjoint）である。よって、その実行のための高速の本来のファンビームアルゴリズムは密接に関連している。従って、高速の本来のファンビームアルゴリズムについては、高速逆投影アルゴリズムの全ての記述が再投影アルゴリズムにおいても同様であるとの理解のもとに、それほど詳しくは記述しない。
【０１８７】
図８に示された、ｆのＭ×Ｍ部分画像
2005518892.tif 000109
【０１８８】
について考える。そのソース角度θｑ＝ｑΔθ，ｑ＝０・・・Ｑ−１でのファンビーム投影は、
【０１８９】
【数１７】
2005518892.tif 000110
【０１９０】
である。
位置
2005518892.tif 000111
【０１９１】
に中心を有するＭ×Ｍ部分画像のための全体座標系内の、その関係する再投影演算子を
2005518892.tif 000112
【０１９２】
として示す。特に、
2005518892.tif 000113
【０１９３】
は、全体画像のファンビーム投影である。ＰＭ，Ｑの計算コストは、或る定数ｃにおいてｃＭ２Ｑである。
次に、図９に示された式（１１）の画像区切りについて考える。投影演算子の直線性のために、以下の、同一数Ｐのソース角度でより小さいＭ×Ｍ画像のＪ個の再投影の合計への、再投影の正確な分解は、次のように得られる。
【０１９４】
【数１８】
2005518892.tif 000114
【０１９５】
再投影のためのこの正確な分解は、Ｎ／Ｍ＝２のものが図１７に示されている。図１７には、図１７中にｆでも示されている電子画像３０は、トランケーション演算子１７１０、１７１２、１７１４及び１７１６によって、それぞれ部分画像ｆ１、ｆ２、ｆ３及びｆ４に細分される。この部分画像が１７１８、１７２０、１７２２及び１７２４によって全体座標系内に再投影されて、部分シノグラムｇ１、ｇ２、ｇ３及びｇ４が生成される。この部分シノグラムが１７２６で集合されて、図１７中にｇでも示されているシノグラム３４が生成される。
【０１９６】
逆投影の式（１２）の正確な分解の場合と同じ理由により、再投影の分解は単一ステップの逆投影
2005518892.tif 000115
【０１９７】
と比べてスピードアップを提供するものではないが、特に近似分解との組み合わせにより、その他の可能な適用を有している。
高速再投影アルゴリズムは、先に述べた以下の追加の特性を有する分解の考えを使用する。すなわち、固定された画像分解能（帯域幅）においては、画像を特徴づけるファンビーム投影の数、又は同様に、θ変数における投影の帯域幅は、画像のサイズに比例する、ということである。従って、もしＰ個の投影がＮ×Ｎ画像ｆを特徴づける場合は、Ｐ′＝Ｐ／（Ｎ／Ｍ）個の投影が部分画像ｆ′を特徴づけ、その投影の残りはこの小さい方の組から得られる。
【０１９８】
従って、高速再投影アルゴリズムは、Ｍ×Ｍの部分画像ｆ′の１組のＰ′個の投影を計算することと、より大きい組のＰ個の投影を得るのに「投影拡大（projection augmentation）プロセス」を使用することに基づいている。この投影拡大プロセスは演算子
2005518892.tif 000116
【０１９９】
によって行われ、これはそれ自体が幾つかの演算子から構成されている。投影シフト演算子
2005518892.tif 000117
【０２００】
及び
2005518892.tif 000118
【０２０１】
は、先に式（１３）で定義されている。Ｌ倍角度の補間演算子
2005518892.tif 000119
【０２０２】
は、Ｌ倍角度のアップサンプリング（upsampling）とその後の角度フィルタリングにより、一組のＱ個の投影における視界（views）の数をＱＬ個の投影に増加させる。角度アップサンプリングは、Ｌ−１個の、全ての２投影間のゼロ値（zero-valued）の投影を導入する。アップサンプルされた投影をｇ（ｐ，ｔ）で示すと、フィルタリングは式（１４）及びそれに関係する議論によって記述される。その結果得られる角度補間演算子
2005518892.tif 000120
【０２０３】
は、高速逆投影アルゴリズムにおいて生じる角度デシメーション演算子
2005518892.tif 000121
【０２０４】
の随伴（adjoint）である。投影拡大演算子
2005518892.tif 000122
【０２０５】
は、シフト及び角度補間演算子の組み合わせであり、
【０２０６】
【数１９】
2005518892.tif 000123
【０２０７】
によって定義される。ここで、Ｐ′＝ＱＬである。演算子
2005518892.tif 000124
【０２０８】
は、補間の前後に、重み付け関数Ｗｋ（θ，ｔ），ｋ＝１，２ による乗法重み付けを含んでおり、これは同様な関数及び
2005518892.tif 000125
【０２０９】
におけるような選択基準を有している。例えば、図１８に示されるように、第１のシフティング関数１８１０の後に、第１の重み付け関数１８１２、補間関数１８１４、第２の重み付け関数１８１６、及び第２のシフティング関数１８１８が続く。もしその方がよければ、投影拡大演算子の他の組み合わせも使用可能である。
【０２１０】
これらの定義により、部分画像ｆ′の近似再投影は、
【０２１１】
【数２０】
2005518892.tif 000126
【０２１２】
である。ここで、ＰＭ，Ｐ／Ｌは、Ｐ／Ｌ個の投影のみを計算する、全体座標系内の再投影である。これにより、式（２０）に類似の区切られた画像のための再投影演算の近似分解が得られる。
【０２１３】
【数２１】
2005518892.tif 000127
【０２１４】
この近似分解は、Ｌ＝Ｎ／Ｍ＝２におけるものが図１９に示されている。図１９において、図１中の高速再投影プロセッサ３２が、図１９中にｆでも示された電子画像３０を受け入れて、１９１０、１９１２、１９１４及び１９１６においてそれぞれ部分画像ｆ１、ｆ２、ｆ３及びｆ４に分割するようプログラムされている。全体座標系内での再投影は、１９１８、１９２０、１９２２及び１９２４で行われて、それぞれＰ／２個の投影を有する部分シノグラムｇ１′、ｇ２′、ｇ３′及びｇ４′が生成される。この部分シノグラムは、前述したように、１９２６、１９２８、１９３０及び１９３２で拡大され、それによって得られるＰ個の投影を有する部分シノグラムｇ１、ｇ２、ｇ３及びｇ４のそれぞれが１９３４で集合されて、図１９中にｇでも示されるシノグラム３４が生成される。
【０２１５】
逆投影の近似分解の場合のように、計算コストは単一ステップの再投影ＰＮ，ＰよりもＬ倍小さく、同じ議論が適用される。
高速再投影アルゴリズムの好ましい実施の形態においては、式（２０）及び（２３）における分解が反復して適用され、式（１１）に従って、各レベルで部分画像がもっと小さな部分画像へと区切られる。これは図１３に示されている。これに対応する、２レベル近似分解の例と、１つの正確な分解及び１つの近似分解の段階から構成された２レベル分解の例が、それぞれ図２０及び２１に示されている。図２０においては、Ｌ＝Ｎ／Ｍ＝２でのファンビーム再投影の２レベル近似分解が、図２０中にｆでも示された電子画像３０で開始される。この画像ｆは２０１０、２０１２、２０１４及び２０１６でそれぞれ部分画像ｆ１、ｆ２、ｆ３及びｆ４に分割され、部分画像ｆ１は更に２０１８、２０２０、２０２２及び２０２４で部分画像ｆ１，１、ｆ１，２、ｆ１，３及びｆ１，４に更に細分される。それらの部分画像はそれぞれ２０２６、２０２８、２０３０及び２０３２によって全体座標系内に再投影されて、それぞれＰ／４個の投影を備える部分シノグラムｇ″１，１、ｇ″１，２、ｇ″１，３及びｇ″１，４が生成される。これらの部分シノグラムは２０３４、２０３６、２０３８及び２０４０によって拡大されて、それぞれＰ／２個の投影を備える部分シノグラムｇ′１，１、ｇ′１，２、ｇ′１，３及びｇ′１，４が生成される。これらは２０４２で集合されて、Ｐ／２個の投影を備える部分シノグラムｇ′１が生成される。２０４４での更なる拡大の後、Ｐ個の投影を備える部分シノグラムｇ１が２０４６で部分シノグラムｇ２、ｇ３及びｇ４と集合される。
【０２１６】
部分画像ｆ４は２０４８、２０５０、２０５２及び２０５４で更に分解され、それら部分画像が２０５６、２０５８、２０６０及び２０６２で同様に再投影されて、それぞれ部分シノグラムｇ″４，１、ｇ″４，２、ｇ″４，３及びｇ″４，４が生成される。２０６４、２０６６、２０６８及び２０７０での拡大の後、それによって得られた部分シノグラムｇ′４，１、ｇ′４，２、ｇ′４，３及びｇ′４，４が２０７２で集合され、そして、それによって得られた部分シノグラムが２０７４で更に拡大される。その結果得られたｇ４が２０４６で他の部分シノグラムと集合されて、図２０中ではｇでも示された、図１中の部分シノグラム３４が生成される。点で示されたブロック２０７６及び２０７８は、それぞれｆ２及びｆ３のための、対応するレベル２の分解を表している。
【０２１７】
図２１は、Ｌ＝Ｎ／Ｍ＝２でのファンビーム再投影の２レベル混合（１つは正確な、もう１つは近似の）分解を示している。図２１において、図中にｆでも示された、図１中の電子画像３０が、２１１０、２１１２、２１１４及び２１１６によって、それぞれ部分画像ｆ１、ｆ２、ｆ３及びｆ４に分解される。部分画像ｆ１は、２１１８、２１２０、２１２２及び２１２４によって部分画像ｆ１，１、ｆ１，２、ｆ１，３及びｆ１，４に細分される。これらの部分画像は、２１２６、２１２８、２１３０及び２１３２によって全体座標系内に再投影されて、それぞれＰ／２個の投影を備える部分シノグラムｇ′１，１、ｇ′１，２、ｇ′１，３及びｇ′１，４が生成される。２１３４、２１３６、２１３８及び２１４０による拡大の後に、それぞれＰ個の投影を備える部分シノグラムｇ１，１、ｇ１，２、ｇ１，３及びｇ１，４が生成される。２１４２での集合により、Ｐ個の投影を備える部分シノグラムｇ１が生成され、これが２１４４で部分シノグラムｇ２、ｇ３及びｇ４と集合されて、Ｐ個の投影を備えるシノグラムｇが生成される。同様に、部分画像ｆ４は、２１４６、２１４８、２１５０及び２１５２で更に細分されて、部分画像ｆ４，１、ｆ４，２、ｆ４，３及びｆ４，４が生成される。２１５４、２１５６、２１５８及び２１６０による再投影によって、部分シノグラムｇ′４，１、ｇ′４，２、ｇ′４，３及びｇ′４，４が生成され、これらが２１６２、２１６４、２１６６及び２１６８によって拡大されることで、部分シノグラムｇ４，１、ｇ４，２、ｇ４，３及びｇ４，４が生成される。これらの部分シノグラムは、２１７０で集合されて部分シノグラムｇ４が生成される。点で示されたブロック２１７２及び２１７４は、それぞれｇ２及びｇ３のための、対応するレベル２の分解を表している。
【０２１８】
部分画像が何らかの所望の最小サイズＭｍｉｎになるまで反復を続け、その後に、再投影ＰＰｍｉｎ，Ｍｍｉｎを式（７）を用いて計算することが可能である。
高速再投影アルゴリズムのための以下の考察は、逆投影の場合に議論したものと類似している。
１．Ｐ、Ｎ、Ｍ、Ｌ及びＭｍｉｎを含む、高速再投影アルゴリズムにおける各種パラメータの選択。
【０２１９】
２．ｔにおける離散の操作。
３．正確度とスピードとの間のトレードオフを制御するための、角度オーバーサンプリングと、正確及び近似の分解ステップの組み合わせの使用。
【０２２０】
４．連続する段階間でのシフティング及び重み付けの演算を組み合わせることによる演算子
2005518892.tif 000128
【０２２１】
及び
2005518892.tif 000129
【０２２２】
の単純化を行い、その結果としての、各段階で単一のシフティング及び（多くとも）単一の重み付け演算。
５．演算の順番の可能な再配置。
【０２２３】
６．Ｏ（Ｎ２ｌｏｇＮ）の、その結果としての階層再投影アルゴリズムの計算の複雑性。
〔等角ファンビームジオメトリのための高速の本来のアルゴリズム〕
等角検出器の場合にも、画像をもっと小さな部分画像に階層的に分解するために、同様なスキームを適用可能である。図２２は、この場合における部分画像ｆ′のための投影ジオメトリを示している。
【０２２４】
式（８）と同様な、ＫＭの定義を用いて、逆投影演算子
2005518892.tif 000130
【０２２５】
は式（９）によって与えられる。ここで、
2005518892.tif 000131
【０２２６】
は、等角ジオメトリ（図６参照）のために先に定義されている。ｆ′上への逆投影に寄与するｇ（ｐ，・）の一部分は、再び
2005518892.tif 000132
【０２２７】
である。ここで今、
2005518892.tif 000133
【０２２８】
は、部分画像
2005518892.tif 000134
【０２２９】
の支持体の投影の角度ｔＡ及びｔＢによって決定される角度支持体へ、ファン角度ｔでｇ（ｐ，・）をトランケートする。同一直線等間隔ファンビームジオメトリにおけるのと同じ議論及び討議が適用され、これにより、式（１０）−（１２）と、図９に示された正確な分解とが得られる。
【０２３０】
近似分解はまた、以下の手順を用いて、部分画像ｆ′を再構成するために使用される投影の数を減少させる。ｔＯＯ′が、ｆ及び部分画像ｆ′のそれぞれの中心Ｏ及びＯ′を通過する２つの放射線間の角度であるとする。同一直線等間隔検出器の場合におけるｆ′の中心Ｏ′の投影に対応する間隔
2005518892.tif 000135
【０２３１】
によるトランケートされた投影をｔシフティングするための直接の類比において、部分画像のためのトランケートされた投影は、ｔＯＯ′によってファン角度でｔシフトされる。同一の角度間隔Δｔで均一に配置されたサンプル点が得られ、これはｔ＝０で中心にある。従って、等角の場合は、角度変数ｔが、等間隔ジオメトリに見られるｔの役割を担い、
2005518892.tif 000136
【０２３２】
が図７のように定義される。或いは、前述したのと同じ反復分解プロセスが、ここで適用される。特に、部分画像に対応するファンビーム投影が一度トランケートされると、それらはソース角度θに関してＬによりデシメートされ、その残りの手順に続く。各種演算子の定義は、既にリストされた変化を除き、ほとんど変更がない。近似分解は、図１１によって再び記述される。
【０２３３】
同様な議論及び変更が、高速の本来の再投影アルゴリズムに適用される。
〔平面等間隔検出器を有するコーンビームジオメトリのための高速の本来の逆投影及び再投影アルゴリズム〕
一般の３Ｄ発散ビームジオメトリのためのプロセスは、２Ｄの場合において述べたのと同じ要素を有しており、以下の主な相違を備えている。
【０２３４】
１．画像は、全ての３座標に沿って３Ｄで再構成される。
２．投影が２次元であるため、個々の投影のトランケーション、シフティング及び補間が２次元である。
【０２３５】
３．投影は、円上のソース角度ではなく、軌道に沿ったソース位置θでデシメート（又は補間）される。
４．ソース軌道はしばしば３Ｄにおいて等方的でないため、好ましい分解は等方的でなく、空間的に変化している場合がある。
【０２３６】
図２３に示された、ｆの部分画像
2005518892.tif 000137
【０２３７】
のための逆投影演算について考える。
2005518892.tif 000138
【０２３８】
が式（８）によって定義された画像トランケーション演算子であるとし、これにより、
2005518892.tif 000139
【０２３９】
は、サイズＭ×Ｍ×Ｍで、
2005518892.tif 000140
【０２４０】
に中心の置かれた、ｆの部分画像となるようにする。式（９）及び（１０）がまだ適用され、この場合、τが２Ｄベクトル
2005518892.tif 000141
【０２４１】
によって置換され、また、部分画像
2005518892.tif 000142
【０２４２】
のＭ×Ｍ×Ｍの支持体の投影の２Ｄ支持体Ωへのオブジェクト全体の投影をトランケートするトランケーション演算子
2005518892.tif 000143
【０２４３】
を用いる（図２３を参照）。
今、画像ｆを、それぞれのサイズがＭ×Ｍ×Ｍである、Ｊ＝（Ｎ／Ｍ）３個のオーバーラップの無い部分画像に区切ることについて考える。
【０２４４】
【数２２】
2005518892.tif 000144
【０２４５】
ここで、ベクトル
2005518892.tif 000145
【０２４６】
は、図２４に示された、Ｍ／Ｎ＝２の場合における部分画像の中心である。
部分画像上への逆投影への逆投影のための正確な分解を提供する式（１２）が適用される。それは、４つの「経路（channels）」の代わりに８つの経路を備えていることを除き、図９に類似した図によって記述されている。計算コストの議論は、ＢＮ，Ｐ及びＢＭ，Ｐの計算コストがそれぞれｃＮ３Ｐ及びｃＭ３Ｐであることを除き、２Ｄの場合におけるものと同様であり、２Ｄの場合とちょうど同じ比Ｊを有している。
【０２４７】
２Ｄファンビームの場合もそうであるが、高速３Ｄコーンビーム逆投影アルゴリズムは、以下の更なる特性を有する分解の考えを使用する。固定された画像分解能（帯域幅）においては、画像を再構成するのに必要とされる投影の数は画像のサイズに比例する。すなわち、Ｎ×Ｎ×Ｎの画像を再構成するのにＰ個の投影が必要とされる場合、同一分解能のＭ×Ｍ×Ｍの画像を再構成するにはＰ′＝（Ｍ／Ｎ）Ｐ個の投影で十分である。
【０２４８】
減少を行うための演算子
2005518892.tif 000146
【０２４９】
を用いて、部分画像
2005518892.tif 000147
【０２５０】
が、減少された数の投影から再び再構成される。それは、以下の変更を伴う、以前のものと同様な演算子の構成によって定義される。
１．今、投影トランケーション演算子
2005518892.tif 000148
【０２５１】
は、先に説明したように、２Ｄ領域Ωにトランケートする。
２．投影シフト演算子
2005518892.tif 000149
【０２５２】
は、ファンビームジオメトリの場合における１Ｄシフトの代わりに、ｔ１、ｔ２座標内でベクトル
2005518892.tif 000150
【０２５３】
だけ２Ｄシフトを行う。このシフトは、例えばまずｔ１でシフティングし、次にｔ２で、といったような分離可能な方法で行われてもよい。
３．演算子
2005518892.tif 000151
【０２５４】
は、角度ではなく、ソース位置でデシメートする。軌道が、直角をなす２つの円のような、２つ以上の別々の部分の結合から構成されている場合、ソース位置でのデシメーションは、各部分毎に分離して行われる。
【０２５５】
４．式（１４）当たりのデシメーションに含まれる
2005518892.tif 000152
【０２５６】
変数におけるフィルタリングは、１次元ではなく２次元である。それも、
2005518892.tif 000153
【０２５７】
フィルタリングはθでのフィルタリングから分離可能であってもなくてもよいが、ｔ１及びｔ２での分離可能な形式で行われる場合がある。
近似分解は、再び式（１８）によって与えられる。今ここで、ＢＭ，Ｐ／Ｌは、式（９）によって定義されたＰ／Ｌ個の投影を用いた、Ｍ×Ｍ×Ｍの部分画像上への全体座標系内での発散ビーム逆投影である。Ｌ＝Ｎ／Ｍ＝２において、この分解は、４つの「経路（channels）」の代わりに８つの経路を備えていることを除き、図１０に似た図によって記載される。計算コストについての議論は、ＢＮ，Ｐ及びＢＭ，Ｐ／Ｌの計算コストがそれぞれｃＮ３Ｐ及びｃＭ３Ｐ／Ｌであって、Ｊ＝（Ｎ／Ｍ）３であることを除き、２Ｄの場合のものと同様である。Ｊ個の部分画像の逆投影の合計コストは、従来の全体画像の逆投影ＢＮ，Ｐにおけるものよりも再びＬ倍小さい。
【０２５８】
高速逆投影アルゴリズムの好ましい実施の形態においては、式（１２）及び（１８）の分解が反復して適用され、各レベルで式（２４）に従って部分画像をもっと小さな部分画像へと更に区切る。ファンビームジオメトリにおいて先に示された、プロセス、各種の特徴、及び改善についての記載は、正確又は近似の分解それ自身においてリストされたのと同様な変更を伴って、同様に一般のコーンビームジオメトリに適用される。特に、反復分解アルゴリズムを用いたコーンビーム逆投影の合計コストはＯ（Ｎ３ｌｏｇＮ）となり、これは従来のコーンビーム逆投影におけるＯ（Ｎ４）と比較される。
【０２５９】
ここまで、我々は、画像が３軸の全てに等しく分解されるという、３Ｄ画像分解の最も単純な実施の形態を記述してきた。しかし、コーンビームジオメトリの幾つかの例においては、ソース軌道が３Ｄに等方的でないために、好ましい分解も等方的ではなく、更には空間的に変化している場合もある。重要な例は、単一円及びヘリカルコーンビーム軌道である（図４、６及び３）。これらの場合は、ｚ軸回りのソース軌道の円柱対称（又はほぼ対称）を示し、ｚ軸で選択された画像区切りはｘ軸及びｙ軸でのものと異なる場合がある。軌道の全体が対称軸を持たない時でさえ、それは時々、直角をなす２つの円から構成されるソース軌道の場合のように、部分の結合へと分解される場合がある。それで、逆投影演算は、それぞれそれ自身の分離した軸対称軌道を有する２つの逆投影の加算として行われ得る。
【０２６０】
そのような例においては、画像の部分をもっと少なく区切るか、又はｚ軸においてもっと大きなサイズの部分画像へと区切るだけで十分である。極端な例では、図２５に示されるように、区切りは全体的に（ｘ，ｙ）方向に限定されてもよい。もっと一般的な、非等方的で空間的に不均一な区切りが、図２６に示されている。これは、ｔ２軸に沿った投影への必要とされるシフト、補間又はフィルタリングの演算の数を減少させることにより、計算の節約へとつながる。
【０２６１】
そのように区切り作業を減少させる理由は、単一円形のソース軌道の場合における式（２）及び（３）について検討することにより明らかになる。これらの式は、検出器上への点の投影の位置がθに依存することを示している。小さなｚ及び／又は小さな
2005518892.tif 000154
【０２６２】
を有する点を考える。すなわち、ソース軌道の平面に近い点、及び／又はｚ軸に近い点である。なお、そのような点では、
2005518892.tif 000155
【０２６３】
が、
2005518892.tif 000156
【０２６４】
よりも、θでもっとゆっくりと変化する。これは、たとえそのｚ次元がその（ｘ，ｙ）次元と同じ比率だけ減少されたとしても、中心
2005518892.tif 000157
【０２６５】
がそのように配置された部分画像が、減少された数の視界（views）から正確に再構成され得る、ということを示すために使用可能である。結果として、（ｘ，ｙ）平面内におけるもっと細かい区切りが、減少された数の投影からの正確な再構成を確実にするために十分である。もし、ｔ２でのｈθ／（２π）だけの全体的なシフトが投影にまず適用されて、θの変化に伴う、ｚ軸に沿った検出器上での座標の中心の動きが補償されるのであれば（図３を参照）、同様な議論がヘリカル軌道の場合にも適用される。
【０２６６】
そのような非等方的で不均一な区切りでは、それに応じて反復の分解が変更される。もし図２５の区切りが反復的に使用されるのであれば、それに対応する分解は、反復の各レベルで４つの経路を有し、それはｔ１変数でのみシフトを含むことになり、また、これは図１０−１５によって記載されることになる。一方、もし図２６の区切りが２レベルの近似分解に使用されるのであれば、分解の第１レベルは、（立方体の粗い区切りからなる８つの八分体（octants）に対応する）８つの経路を有することになるが、次のレベルでのそれぞれの分解は、それぞれの粗い画像八分体を区切って得られた７つの画像に対応する７つの経路を有するのみとなる。
【０２６７】
先に議論した発散ビームジオメトリの特定の例を含む、一般の発散ビームジオメトリのための高速再投影の設計は、それに対応する高速逆投影アルゴリズムの双対（dual）である。それは、高速ファンビーム再投影が高速ファンビーム逆投影に対してもたらす関係と同様な関係を、高速発散ビームジオメトリ逆投影に対してもたらす。後者の関係は既に少々詳細に記載してあるので、ここではそれを繰り返し述べない。
〔実行及び実験〕
新規なアルゴリズムは、ＭＡＴＬＡＢＴＭ及びＣプログラミング言語で実行され、２Ｄの場合にはサイズＮ×Ｎ＝５１２×５１２で、３Ｄの場合にはサイズＮ×Ｎ×Ｎ＝１２８×１２８×１２８のシェップローガンヘッドファントム（Shepp-Logan head phantom：断層撮影アルゴリズムの数値評価に使用される標準テスト画像）の再構成を含む幾つかの数値実験でのテストに成功した。我々は、頭蓋骨と脳との間のコントラストが明確な（及び高い）、いわゆる未修正ファントム（unmodified phantom）を使用した。これは、このアルゴリズムのための魅力的なテスト（challenging test）を提示するものであり、それは頭蓋骨の再構成によって生成された僅かなアーティファクトでさえ脳の背景で目立つためである。
【０２６８】
まず、２Ｄファンビーム再構成の結果について記述する。等間隔と等角度のファンビームイメージングジオメトリの両方において、ソースと画像中心間の距離を画像サイズの１．２５倍に、すなわちＤ＝１．２５Ｎに設定し、ファン上に１０２５個の検出器を有する場合における合計のファン角度を２α０＝１．１７ラジアンに設定した。Ｐ＝１５３６個の均一間隔で配置されたソース角度、１０楕円（10-ellipse）のファントムにおけるθＰ∈[０，２π]、を有するフルスキャンのファンビーム投影データを、楕円を通る線積分のための分析表現を用いて生成した。
【０２６９】
我々のシミュレーションにおいては、各レベルで画像を半分サイズの４つの部分画像に分解した。近似分解では、デシメーション因子Ｌ＝２を用いた。画像は、単一ピクセルの部分画像、すなわちＭｍｉｎ＝１、に繰り返し分解された。角度オーバーサンプリングを増大させるために、分解の最初のＱレベルは、角度デシメーションなしの正確な分解とした。Ｑ＝２及びＱ＝３を用いた。ここで、Ｑは「ホールドオフパラメータ（holdoff parameter）」と言うことにする。投影トランケーション演算子
2005518892.tif 000158
【０２７０】
のためのパラメータ（「足跡（footprint）」）は、計算が少々増えるコストで、単純な控えめな近似を用いて大急ぎで計算された。
図２７及び２８は、正確な逆投影アルゴリズムと、２つのファンビームイメージングジオメトリのための我々の高速ファンビーム逆投影アルゴリズムとを比較している。参考及び目安として、図２７（ａ）及び２８（ａ）は元々のファントムを示し、図２７（ｂ）及び２８（ｂ）は、等間隔及び等角度の検出器の両方のための正確なファンビーム逆投影を用いて再構成された画像である。図２７（ｃ）、２８（ｃ）、２７（ｄ）及び２８（ｄ）は、それぞれホールドオフパラメータＱ＝３及びＱ＝２を有する我々の高速アルゴリズムを用いた再構成である。高速アルゴリズムのスピードアップは、ホールドオフがＱ＝３及びＱ＝２においてそれぞれ３０×及び６０×である。正確かつ高速のアルゴリズムを用いた再構成は、知覚可能な品質相違について、Ｑ＝２においてさえもほとんど示していない。図２８中の再構成された画像の列及びコラム（縦行）を通してのプロットが、図２９（ａ）及び２９（ｂ）に示されている。図２９（ａ）は第４１０番目の列を通るスライスを比較しており、図２９（ｂ）は図２８の第３５６番目のコラムを通るスライスを比較している。これらのプロットにおいて、実線は元々のファントムを表し、破線は正確なファンビーム逆投影を表し、鎖線はホールドオフ３の高速アルゴリズムを表し、点線はホールドオフ２の高速アルゴリズムを表している。
【０２７１】
次に、等間隔検出器ジオメトリのための正確及び高速のファンビーム逆投影アルゴリズムの点分布関数（ＰＳＦ’ｓ：point-spread functions）を比較する。位置（１０，１０）、（１２８，１２８）、（４９，３７）及び（２０１，９３）に配置された４つの点ターゲットを有するサイズ２５６×２５６の画像上で、正確及び高速のアルゴリズムを実行した。先に使用したものと同様な以下のパラメータを用いた。すなわち、ソースと中心間の距離Ｄ＝１．２５Ｎ、ソース角度の数Ｐ＝２Ｎ、及び７３７個の等間隔検出器での合計ファン角度が１．２０ラジアンである。点応答（point responses）は、４つの位置の全てで非常に類似しており、本質的にシフト不変の応答を示唆した。（Ｑ＝３の）正確及び高速のアルゴリズムで得られた（１２８，１２８）での点ターゲットにおける点分布関数（ＰＳＦ’ｓ）が、図３０（ａ）及び３０（ｂ）において比較されている。図３０（ａ）は正確なアルゴリズムにおける中心ピクセルＰＳＦを示し、図３０（ｂ）はホールドオフが２での高速アルゴリズムにおける中心ピクセルＰＳＦを示している。これらＰＳＦ’ｓの２つの軸を通るスライスが図３１（ａ）及び３１（ｂ）に示されている。図３１（ａ）は、正確なアルゴリズム（実線）及びＱ＝３の高速アルゴリズム（破線）における第１２８番目の列を通るスライスを表している。図３１（ｂ）は、正確なアルゴリズム（実線）及びＱ＝３の高速アルゴリズム（破線）における第１２８番目のコラムを通るスライスを表している。２つのアルゴリズムのＰＳＦ’ｓはほとんど同一であり、これは高速アルゴリズムの正確さを裏付けている。
【０２７２】
次に、（ｘ，ｙ）平面内で単一円形ソース軌道を有するコーンビーム再構成を用いた実験の結果について議論する。ソースからボリュームの中心までの距離が、側面部の長さの１．２５倍である。コーンビーム投影は、[０，２π]内に均一に配置された５１２個のソース角度で行われた。（方位角における）ファン角度は１．２０ラジアンであり、コーン角度は１．５ラジアンであった。それぞれ図４及び図６に示された平面等間隔及び円柱状の検出器ジオメトリの両方について考える。それぞれの場合において、検出器面は、ｔ１座標において等間隔又は等角度のいずれかの間隔で、垂直軸（ｔ２）に沿って等間隔配置された、３７５×３７５の検出器から構成された。
【０２７３】
平面検出器を有する単一円コーンビーム軌道のための従来の再構成アルゴリズムは、よく知られたフェルドカンプ、デービス及びクレス（ＦＤＫ：Feldkamp, Davis and Kress）アルゴリズムであり、これは、投影に重み付けをし、ｔ１座標に沿ってフィルタリングし、そして次に、重み付けされたコーンビーム逆投影を行うことから構成されている。この重み付けとフィルタリングは、平面検出器を有するファンビームジオメトリにおけるものと同じであり、相違するのは、ファンビーム逆投影ではなくコーンビーム逆投影を用いる点である。これは近似アルゴリズムであり、大きなコーン角度では幾つかの固有のアーティファクトを受ける。それでも、それはコーンビームジオメトリのための標準的な目安である。円柱状検出器ジオメトリにおいては、重み付けされたコーンビーム逆投影を円柱状検出器に組み込むことにより、等角検出器ジオメトリのためのファンビームアルゴリズムをコーンビームジオメトリに拡張することによって得られた、類似のアルゴリズムを用いた。
【０２７４】
デモンストレーション目的のために、対応する高速アルゴリズムは、前述した反復分解を有する逆投影を用いており、これは図２５に示された単純な区切りスキームで行われる。我々のシミュレーションでは、アルゴリズムの特別な実行中に含まれるオーバーヘッドが小ボリュームのための逆投影演算子のスピードアップを相殺するので、部分画像サイズが４×４×１２８にまで減少された時に分解プロセスを終了した。高速アルゴリズムは、従来のＦＤＫアルゴリズムと同じフィルタリング及び重み付けを用いた。
【０２７５】
従来のＦＤＫ再構成の結果と、本発明のプロセスを用いる高速再構成の結果が、それぞれ平面及び円柱状検出器ジオメトリにおいて図３２及び図３３で比較されている。図３２は、単一円形軌道及び平面等間隔検出器ジオメトリにおける３Ｄコーンビーム再構成を比較している。最上位の列はｚ＝−０．２５でのｘｙスライスであり、中央の列はｘ＝−０．１３でのｙｚスライスであり、最下位の列はｙ＝−０．０４でのｘｚスライスである。左のコラムは従来のＦＤＫアルゴリズムを用いた再構成の結果であり、中央及び右のコラムはそれぞれホールドオフ２及び１での本発明の高速アルゴリズムを用いている。図３３は、単一ソース軌道及び円柱状検出器ジオメトリにおける３Ｄコーンビーム再構成の比較である。最上位の列はｚ＝−０．２５でのｘｙスライスであり、中央の列はｘ＝−０．１３でのｙｚスライスであり、最下位の列はｙ＝−０．０４でのｘｚスライスである。左のコラムは従来のＦＤＫアルゴリズムを用いた再構成結果であり、中央及び右のコラムはそれぞれホールドオフ２及び１での本発明の高速アルゴリズムを用いている。図３４は、図３２中の対応する列における画像を通るスライスを表している。最上位のプロットは列１０２を示し、中央のプロットは列８０を示し、最下位のプロットはコラム８１を示している。実線は従来のＦＤＫ再構成のものである。ホールドオフ２での高速アルゴリズムが鎖線で示されており、ホールドオフ１での高速アルゴリズムが点線で示されている。図３５は、図３３中の対応する列における画像を通るスライスを示している。最上位のプロットは列１０２を示し、中央のプロットは列８０を示し、最下位のプロットはコラム８１を示している。実線は従来のＦＤＫ再構成のものであり、鎖線はホールドオフ２での高速アルゴリズムを用いた結果を表し、点線はホールドオフ１での高速アルゴリズムを用いた結果を表している。更なる詳細が、図３４及び３５によって提供され、これらは図３２及び３３の画像を通る或る線に沿った、再構成された密度のプロットを示している。ホールドオフ＝１及びホールドオフ＝２における結果が与えられている。
【０２７６】
ホールドオフ因子Ｑ＝２での本発明を用いたＦＤＫアルゴリズムの高速バージョンは、従来のフェルドカンプ（Feldkamp）アルゴリズムよりも約３倍速く、画像品質における劣化はほとんど見られない。ホールドオフ因子Ｑ＝１でさえも、画像品質は幾つかの応用において容認できるものであり、７倍のスピードアップが図れる。５１２×５１２×５１２の画像ボリューム、図２４又は２６におけるように十分になされた３Ｄ区切り、及び対応する分解においては、スピードアップは１０倍の向上が見込める。（事実、円形又はヘリカル軌道のような、ほぼ円柱対称を有するコーンビームジオメトリにおいては、何らかのｚサイズＨにおける５１２×５１２×Ｈでは１０倍のスピードアップが得られる。）
〔実験の結論〕
２Ｄにおける我々の実験において、新規なファンビーム逆投影アルゴリズムは、５１２×５１２の画像において３０×−６０×の実際のスピードアップを提供し、正確度における認識可能な損失はなかった。３Ｄにおいて、単一円軌道及び平面検出器のコーンビームジオメトリで、分解の最も単純な形式を用いた場合、スピードアップは１２８×１２８×１２８において３×−７×であった。本発明で提案された分解を十分に実行することで、５１２×５１２×５１２の画像において、１０倍のスピードアップ因子が可能である。
【０２７７】
この実験により、異なる検出器ジオメトリにおける本発明のプロセスの有効性が実証された。単純な単一円形軌道ジオメトリにおけるこれらの結果は、ヘリカルコーンビームのような、もっと一般的な発散ビームジオメトリに拡張される。再投影アルゴリズムと逆投影アルゴリズムとの間の、原理及び構造の類似性のために、本発明の高速発散ビーム再投影アルゴリズムにおいても同様な結果が維持される。再投影アルゴリズムと逆投影アルゴリズムの両方が、個々にも組み合わせても有用であり、それらは、例えば米国特許第６，２６３，０９６号及び第６，３５１，５４８号において説明したように、断層撮影再構成への各種の反復アプローチや、断層撮影におけるアーティファクト修正に使用可能である。逆投影と再投影の組み合わせは、ヘリカルコーンビームジオメトリにおける長いオブジェクトの問題のための高速再構成アルゴリズムを構築するためにも使用可能である。
【０２７８】
以上のように、本発明の原理は特定の装置及び応用と関係付けて記述されているが、この記述はほんの一例であって、本発明の範囲を限定するものではないことを理解すべきである。
【図面の簡単な説明】
【０２７９】
【図１】本発明で使用する装置のブロック図である。
【図２】平面等間隔コーンビームジオメトリ（a planar eqispaced cone-beam geometry）の図である。
【図３】平面等間隔ヘリカルコーンビームジオメトリ（a planar eqispaced helical cone-beam geometry）の図である。
【図４】円形ソース軌道及び平面等間隔検出器を有するコーンビーム断層撮影を示す図である。
【図５】同一直線等間隔ファンビームジオメトリ（colinear eqispaced fan-beam geometry）の図である。
【図６】円形ソース軌道及び円柱状検出器を有するコーンビーム断層撮影を示す図である。
【図７】等角コーンビームジオメトリ（an equiangular fan-beam geometry）の図である。
【図８】ＭｘＭ部分画像における同一直線等間隔ファンビーム投影ジオメトリの図である。
【図９】画像の部分画像への区切りの図である。
【図１０】ファンビーム逆投影の正確な分解を示す図である。
【図１１】シノグラム分解演算子Ｏ［Ｌ,Ｍ,δ］によって行われる処理のフローチャートである。
【図１２】ファンビーム逆投影の近似分解の図である。
【図１３】２レベルのダイアディック（dyadic）画像の区切りの図である。
【図１４】ファンビーム逆投影の２レベル近似分解の図である。
【図１５】ファンビーム逆投影の２レベル混合（１つは正確な、１つは近似の）分解の図である。
【図１６】高速ファンビーム階層逆投影のための翻訳言語である。
【図１７】ファンビーム再投影の正確な分解を示す図である。
【図１８】投影拡大演算子(projection augmentation operator)Ｖ［Ｌ,Ｍ,δ］によって行われる処理のフローチャートである。
【図１９】ファンビーム再投影の近似分解の図である。
【図２０】ファンビーム再投影の２レベル近似分解の図である。
【図２１】ファンビーム再投影の２レベル混合（１つは正確な、１つは近似の）分解の図である。
【図２２】部分画像のための等角ファンビーム投影ジオメトリの図である。
【図２３】部分画像のための平面等間隔コーンビーム投影ジオメトリの図である。
【図２４】３Ｄ画像の部分画像への区切りの図である。
【図２５】３Ｄ画像の部分画像への非等方区切りの図である。
【図２６】３Ｄ画像の２レベル非一様反復区切りの図である。
【図２７】同一直線等間隔検出器における正確及び高速のファンビーム逆投影を比較する写真である。
【図２８】等角検出器における正確及び高速のファンビーム逆投影を比較する写真である。
【図２９】図２８（ａ）〜２８（ｄ）の逆投影のスライスを通って得られたプロットである。
【図３０】点分布関数（point spread function）の比較である。
【図３１】点分布関数を通るスライスのプロットの比較である。
【図３２】３Ｄコーンビーム再構成の写真の集まりである。
【図３３】３Ｄコーンビーム再構成の写真の他の集まりである。
【図３４】図３２中の画像を通るスライスのプロットの集まりである。
【図３５】図３３中の画像を通るスライスのプロットの集まりである。
【図１】
2005518892.tif 000159
【図５】
2005518892.tif 000160
【図７】
2005518892.tif 000161
【図８】
2005518892.tif 000162
【図９】
2005518892.tif 000163
【図１０】
2005518892.tif 000164
【図１１】
2005518892.tif 000165
【図１２】
2005518892.tif 000166
【図１３】
2005518892.tif 000167
【図１４】
2005518892.tif 000168
【図１５】
2005518892.tif 000169
【図１６】
2005518892.tif 000170
【図１７】
2005518892.tif 000171
【図１８】
2005518892.tif 000172
【図１９】
2005518892.tif 000173
【図２０】
2005518892.tif 000174
【図２１】
2005518892.tif 000175
【図２２】
2005518892.tif 000176
【図２４】
2005518892.tif 000177
【図２５】
2005518892.tif 000178
【図２６】
2005518892.tif 000179
【図２７（ａ）】
2005518892.tif 000180
【図２７（ｂ）】
2005518892.tif 000181
【図２７（ｃ）】
2005518892.tif 000182
【図２７（ｄ）】
2005518892.tif 000183
【図２８（ａ）】
2005518892.tif 000184
【図２８（ｂ）】
2005518892.tif 000185
【図２８（ｃ）】
2005518892.tif 000186
【図２８（ｄ）】
2005518892.tif 000187
【図２９】
2005518892.tif 000188
【図３０（ａ）】
2005518892.tif 000189
【図３０（ｂ）】
2005518892.tif 000190
【図３１】
2005518892.tif 000191
【図３２】
2005518892.tif 000192
【図３３】
2005518892.tif 000193
【図３４】
2005518892.tif 000194
【図３５】
2005518892.tif 000195
【図２】
2005518892.tif 000196
【図３】
2005518892.tif 000197
【図４】
2005518892.tif 000198
【図６】
2005518892.tif 000199
【図２３】
2005518892.tif 000200
【国際調査報告】
2005518892.tif 000201
2005518892.tif 000202