(57)【要約】本発明は、リアルタイムでX線画像のグラディエント適応フィルタリング（MGRAF）を伴う多重解像の方法に関する。２Kの隣接行の画像ストリップに対して、ラプラシアン・ピラミッド（L0,…L3）とガウシアン・ピラミッド(G0,…G3)への解像がK番目の段まで行われる。処理動作をそのようなストリップに限定することによって、全ての適切なデータを局所メモリに高速アクセス（キャッシュ）によって確保することが可能である。従来のアルゴリズムと比較した別の高速化は、グラディエント（D）をガウシアン・ピラミッド表現から計算することによって達成される。これらの最適化尺度及び別の最適化尺度のおかげで、グラディエント適応フィルタリングを伴う多重解像を毎秒30を超える（768×564）画像の処理速度に増加させることが可能である。
2006510411.tif 000002
【特許請求の範囲】
【請求項１】
N行の像点を備える入力画像を処理する方法であって：
a)該入力画像のn＜N隣接行を備える画像ストリップはK詳細画像のシーケンスに解像され、該シーケンスは各々の場合に該入力画像の空間周波数の一部だけを有し；
b)該詳細画像のうちの少なくとも１つが修正され；
c)出力画像ストリップが場合によっては該修正されている詳細画像から再構築され；
d)工程a)、b）及びc)が該入力画像の別の画像ストリップについて繰り返され；
e)出力画像が計算出力画像ストリップから再構築されることを特徴とする方法。
【請求項２】
請求項1記載の方法であって、各画像ストリップは、K段を備えているラプラシアン・ピラミッド及びガウシアン・ピラミッドに解像されることを特徴とする方法。
【請求項３】
請求項1記載の方法であって、該画像ストリップ各々は2K行の幅を有することを特徴とする方法。
【請求項４】
請求項1記載の方法であって、解像段j<Kの詳細画像の該修正がフィルタを用いて行われ、該フィルタの係数が画像ストリップから計算される少なくとも１つのグラディエントによって変わってくることを特徴とする方法。
【請求項５】
請求項2及び4記載の方法であって、グラディエントがj番目の解像段のガウシアン・ピラミッド表現から計算されることを特徴とする方法。
【請求項６】
請求項4記載の方法であって、
該フィルタ係数
【数１】
2006510411.tif 000032
が所定のフィルタの係数
【数２】
2006510411.tif 000033
から、
【数３】
2006510411.tif 000034
によって計算され、
【数４】
2006510411.tif 000035
はフィルタ処理によって処理される該像点であり、
【数５】
2006510411.tif 000036
は、該フィルタの中心に対する該係数の位置であり、
【数６】
2006510411.tif 000037
は、画像位置
【数７】
2006510411.tif 000038
における該グラディエントであり、0≦r<1であることを特徴とする方法。
【請求項７】
請求項6記載の方法であって、
【数８】
2006510411.tif 000039
であり、
【数９】
2006510411.tif 000040
は、グラディエント・フィールドと、該グラディエント・フィールドの分散とによって好ましくは変わってくる正の因子であることを特徴とする方法。
【請求項８】
N行の像点を備えるディジタル入力画像を処理するデータ処理装置であって：
データ処理装置は、システム・メモリ及びキャッシュ・メモリを有し、以下の処理工程を行うことが意図されており、該処理工程は：
a)該入力画像のn＜N隣接行を備える画像ストリップをK詳細画像のシーケンスに解像する工程を備え、該シーケンスは各々の場合に該入力画像の空間周波数の一部だけを有し；
更に、b)該詳細画像のうちの少なくとも１つを修正する工程；
c)出力画像ストリップを、場合によっては該修正されている詳細画像から再構築する工程；
d)工程a)、b）及びc)を該入力画像の別の画像ストリップについて繰り返す工程；及び
e)出力画像を該計算出力画像ストリップから再構築する工程を備え；
工程a)乃至c)では、処理データは全て、各場合には該キャッシュ・メモリにあることを特徴とするデータ処理装置。
【請求項９】
請求項8記載のデータ処理装置であって、並列プロセッサ及び/又はベクタ・プロセッサを有することを特徴とするデータ処理装置。
【請求項１０】
X線システムであって：
X線源；
X線検知器；及び
該X線検知器に結合されて、該X線検出器によって送信される該X線入力画像を処理する請求項8又は9記載のデータ処理装置を備えることを特徴とするX線システム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、一般的に、入力画像を処理する方法及びデータ処理装置に関し、特に、リアルタイムでX線画像の多重解像とグラディエント適応フィルタリングとを行う方法及びデータ処理装置に関する。
【背景技術】
【０００２】
画像の自動評価は多数の種々の適用分野において行われている。以下に更に詳細に検討する透視X線画像の処理はよって、一例に過ぎないものとして解されることとする。患者と医員が受けるX線の量を最小にするために、X線画像は、できるだけ低い放射線量を伴って撮影される。しかし、重要な画像細部が画像雑音において喪失されることになるというリスクが存在する。このことを妨げるために、X線の画像又は画像シーケンスに対して空間フィルタ及び時間フィルタを用いて雑音を、当該過程で適切な画像情報を損なうことなく抑制する試みが行われている。
【０００３】
そのような画像処理の意味合いにおいて、入力画像の多重解像として知られているものが多くの場合、行われている。入力画像はこの場合、詳細画像シーケンスに解像され、各々の詳細画像は（空間）周波数で関連した領域すなわちストリップからの画像情報を有する。更に、詳細画像はその解像度、すなわち画像のコンテンツの表現に対する像点の数によってその当該周波数範囲に適応される。詳細画像を修正することによって、特定の周波数範囲に対する影響を当該方法で及ぼすことが可能である。この修正後、詳細画像はもう一度合成させて出力画像を構成することが可能である。
【０００４】
この点で、X線画像の後処理を行う効果的な方法で公知のものがある（特許文献１及び特許文献２参照。）。これらの文献では、多重解像が行われ、それによって得られる詳細画像がフィルタを用いて修正され、その係数は画像グラディエントに基づいて適応されている。この場合のグラディエントは多重解像のより粗い解像段から生じている。MRGAF（多重解像グラディエント適応フィルタリング）として呼ばれているこの方法では、低域通過フィルタリングが、線又は端などの画像構造に垂直な方向には当該構造の方向よりも少ない程度で行われるので、情報を得ることを可能にする、雑音の抑制が行われる。しかし、要求される計算の労力が多大であるという理由で、この方法は現在まで、記憶された画像又は画像シーケンスに対してオフラインでしか実施することが可能でない。
【特許文献１】国際公開特許第98/55916号
【特許文献２】欧州公開特許第996090号
【発明の開示】
【発明が解決しようとする課題】
【０００５】
この背景に鑑みて、本発明の目的は、多重解像を用いて入力画像をより効率的に処理する手段を備えることにあり、好ましくは、画像解析をリアルタイムで行うことが可能であることとする。
【課題を解決するための手段】
【０００６】
この目的は請求項1記載の特徴を有する方法によって達成され、請求項8記載の特徴を有するデータ処理装置によって達成され、請求項10記載の特徴を有するX線システムによって達成される。効果的な精緻化は従属クレイムに記載されている。
【０００７】
本発明による方法は、N行の像点を備える入力画像を処理するのに用いられる。通常、像点は行に垂直な列を有する矩形グリッドに配置されるが、例えば六角形グリッドなどの行構造を有する別の配置も可能である。入力画像は特に、ディジタル化透視X線画像であり得るが、この方法はこれに制限されず、画像の多重解像が行われる同等な適用場面において効果的に用いることが可能である。この方法は以下の工程を備える。
【０００８】
a)入力画像のn<Nの隣接する行を備える画像ストリップがKの詳細画像シーケンスに解像され、各場合における詳細画像は画像ストリップの空間周波数の部分範囲のみを有する。多重解像はよって、全体入力画像のストリップ状の部分を用いて行われる。
【０００９】
b)工程a)において得られる詳細画像の少なくとも１つが、例えば、所定のフィルタ又は画像ストリップから計算されるフィルタを用いて修正される。好ましくは、修正に必要な情報全てが画像ストリップにおいて利用可能である。
【００１０】
c)出力画像ストリップが詳細画像又は（修正詳細画像が存在する場合）修正詳細画像から再構築される。
【００１１】
d)上記工程a)、b)及びc)が、入力画像の別の画像ストリップについて繰り返される、すなわち相当する出力画像ストリップの計算と類似した方法で行われる。ストリップ幅n及び/又は解像度の深さKに対する別の値も、適宜、みなし得る。
【００１２】
e)出力画像が計算出力画像ストリップから再構築される。
【００１３】
その結果、上記方法では、出力画像が、よって、入力画像から生成され（その出力画像は同じサイズのものか異なるサイズのものであり）、所望の種類の修正は入力画像の一部又は全部の空間周波数範囲において行われている。詳細画像の修正を伴う従来の多重解像と比較すれば、この方法との違いは、多重解像が各場合においてn行の画像ストリップ上にて部分単位で行われるという点である。この場合では、各画像ストリップが段Kに解像され、更にもう一度合成されて、出力画像ストリップを得る。この手順の効果は、画像ストリップの処理に対するメモリ要件が完全な画像の処理よりも相応に小さいため、データ処理システム上での効率的な実施に特に適しているので、この方法はワーキング・メモリを用いて高速アクセスによって行うことが可能である。その結果、多くの場合、多重解像をリアルタイムで行うことが初めて実用的となっているほど大きな、速度における利得を達成することが可能である。
【００１４】
この方法の１つの特定の精緻化によれば、工程a)の多重解像では、各画像ストリップはその各々の場合K段を備えているそのラプラシアン・ピラミッド及びガウシアン・ピラミッドに解像される。ガウシアン・ピラミッドの段jでは、段入力画像は、先行段(j−1)の出力画像であり、（「段jのガウシアン・ピラミッド表現」として以下に表す）出力画像は低域通過フィルタリング及び後続解像リダクションによって修正された段入力画像である。（以下に「段jのラプラシアン・ピラミッド表現」として表す）段jでのラプラシアン・ピラミッドの出力画像は同じ画像jのガウシアン・ピラミッド表現を減算することによって得られ、その解像度は、先行段(j−1)のガウシアン・ピラミッド表現から、もう一度増加されており、低域通過フィルタリングされている。入力画像をラプラシアン・ピラミッド又はガウシアン・ピラミッドに解像することは、医療画像処理において頻繁に用いられており、画像ストリップに対して用いることは特に適切である。
【００１５】
好ましくは、工程a)乃至d）では、多重解像にかけられる画像ストリップは各々の場合、幅は2K行であり、Kは多重解像の解像段の数である。2Kの幅を有する画像ストリップは、解像の各段で、各場合において行と列とを2分の１に削減することが行われる場合に、段Kまでラプラシアン・ピラミッド又はガウシアン・ピラミッドに解像するのに必要な最小幅を有する。最も粗い段の詳細画像は、そのような画像ストリップの場合、１行の最小幅を有する。更に、画像ストリップは任意的には、お互いに各場合において、(2K−1)行オフセットされる、すなわち、お互いに各場合において、1行重なり合う。そのように重なり合うことは、好ましくは、画像ストリップの解像段全てに存在し、新たな、重なり合っていない領域の端で行われるフィルタリング処理に関して必要な情報を備える。用いられるフィルタの幅に応じて、画像ストリップ間の幅が1行を超えるように重なり合う場合もある。
【００１６】
詳細画像によって行われる修正の種類は適用場面に応じて異なり得る。好ましくは、解像段j＜Kの詳細画像の修正は、フィルタを利用することによるものであり、このフィルタの係数は、画像ストリップから計算される少なくとも１つのグラディエントによって変わってくる。画像のグラディエントは、画像における局所構造の位置を反映するので、異方性フィルタを定義するのに用いることが可能であり、これを用いることによって当該構造は変わらない状態のままとなるか、増幅させることもあり、当該構造に沿った何れかの雑音を抑制する。
【００１７】
好ましくは、上記方法は、ガウシアン・ピラミッドとラプラシアン・ピラミッドとに解像することと組み合わせられ、グラディエントは解像段jのガウシアン・ピラミッド表現から計算され、同じ段jのラプラシアン・ピラミッド表現のフィルタリングに用いられる。このことは、修正に必要な情報全てを解像段jのデータから得ることが可能であるので、修正をこの段の計算中に直接行うことが可能である。
【００１８】
詳細画像の上記グラディエント適応フィルタリングの特別な設計によれば、フィルタ係数
【００１９】
【数１０】
2006510411.tif 000003
は、例えば二項フィルタなどの所定のフィルタの係数
【００２０】
【数１１】
2006510411.tif 000004
から計算され、
【００２１】
【数１２】
2006510411.tif 000005
はフィルタによって処理される像点であり、
【００２２】
【数１３】
2006510411.tif 000006
はフィルタの中心に対する当該係数の位置であり、
【００２３】
【数１４】
2006510411.tif 000007
の式があてはまる。
【００２４】
上記式では、
【００２５】
【数１５】
2006510411.tif 000008
画像位置
【００２６】
【数１６】
2006510411.tif 000009
におけるグラディエントであり、0≦r<1である。重み付け関数が
【００２７】
【数１７】
2006510411.tif 000010
の場合、相当するフィルタ係数βは減少し、フィルタリングの結果に対するその寄与は減少する。このようにして、雑音の寄与は画像の相当する位置で抑制される。
【００２８】
重み付け関数rは好ましくは
【００２９】
【数１８】
2006510411.tif 000011
として定義され、
【００３０】
【数１９】
2006510411.tif 000012
は好ましくは、グラディエント・フィールド
【００３１】
【数２０】
2006510411.tif 000013
とその分散によって変わってくる因子である。因子rの上記定義は、
【００３２】
【数２１】
2006510411.tif 000014
である、グラディエント
【００３３】
【数２２】
2006510411.tif 000015
に垂直な方向においてr=1であり、かつ
【００３４】
【数２３】
2006510411.tif 000016
である、
【００３５】
【数２４】
2006510411.tif 000017
に平行な方向においてrが最小であるという所望の特性を有する。αとrとの計算を行ううえでの定義は、上記特許文献1記載の定義のものよりも計算労力の点でかなり容易であり、その結果はほぼ同じである。
【００３６】
本発明は更に、N行の像点を備えるディジタル入力画像を処理するデータ処理装置に関し、データ処理装置はシステム・メモリ及びキャッシュ・メモリを有する。データ処理装置は：
a)その各々の場合において入力画像の空間周波数の一部のみを有する、入力画像のn<Nの隣接行を備える画像ストリップをKの詳細画像シーケンスに解像する工程；
b)詳細画像のうちの少なくとも１つを修正する工程；
c)場合によっては修正されている詳細画像から出力画像ストリップを再構築する工程；
d)入力画像の別の画像ストリップについて工程a)、b)及びc)を繰り返す工程；及び
e)出力画像を計算出力画像ストリップから再構築する工程；
の処理工程を行うことが意図されており、工程a)乃至c)では、全ての処理データ（画像ストリップのデータ、画像ストリップの多重解像の関連詳細画像のデータ）は各々の場合、キャッシュ・メモリにある。
【００３７】
そのようなデータ処理装置を用いて、上記方法を非常に効率的にかつすばやく行うことが可能であるが、それは、必要なデータ全てをキャッシュ・メモリに収容し、よって、すばやくアクセスすることが可能であるからである。対照的に、従来の多重解像度では、各場合において、完全な画像が解析され、その結果、中間結果を記憶させるよう、システム・メモリ（ワーキング・メモリ、ハード・ディスクなど）を行わなければならない。計算時間の大部分はよって、システム・メモリとの間のデータの読み書きに費やされる。これらの、時間のかかる処理が割愛されるので、上記のデータ処理装置を用いて画像処理をリアルタイムでも行うことが可能である。
【００３８】
好ましくは、データ処理装置には、並列プロセッサ及び/又はベクタ・プロセッサが装備される。この場合、必要な処理は、並列化によってなおもっと高速化することが可能である。
【００３９】
更に、データ処理装置は、好ましくは、前述の方法の変形も行うことが可能であるように企図される。
【００４０】
本発明は更に：
X線源；
X線検知器；
X線検知器に結合されて、X線検知器によって生成されるX線入力画像を処理するデータ処理装置を備え、該データ処理装置は上記方法で企図される。
【００４１】
そのようなX線システムの効果は、効果的な画像処理はリアルタイムで、すなわち、医療介入中に、行うことが可能であり、その画像処理は当該構造を損なうことなく雑音を抑制する。特に、MRGAF方法はリアルタイムで行うことが可能である。情報を得ることを可能にする、雑音の抑制が理由で、相当する低照射量によってX線写真を撮影し、よって、患者と医員が受ける照射量を最小にすることが可能である。
【発明を実施するための最良の形態】
【００４２】
本発明は、添付図面に示す実施例の例を参照しながら更に説明することとするが、本発明はその例に制限されるものでない。
【実施例】
【００４３】
図1に略図を示すMRGAFアルゴリズムは上記特許文献1及び特許文献2に詳細に記載されており、したがって、以下の概要によってのみ説明することとする。MRGAFアルゴリズムの目的は、画像細部と画像鮮鋭度とを維持すると同時に、入力画像Iにおいて雑音をかなり削減するというものである。アルゴリズムの基本的な考え方は、多重解像と、局所画像グラディエントの関数としての結果詳細画像の異方性低域通過フィルタリングとを有する。
【００４４】
図1に示す例では、512x512の像点（画素）を備える入力画像Iの解像が、K=4の解像段にて行われる。各解像段j=0、1、2、3では、ラプラシアン・ピラミッド表現Λj及びガウシアン・ピラミッド表現Γjとして知られているものが詳細画像として生成される。段入力表現は各場合において、先行段(j−1)のガウシアン・ピラミッド表現Γj−1又は元の入力表現Iである。ガウシアン・ピラミッド表現Γjは当該段入力表現に対するリダクション処理Rを用いて生成され、「リダクション」は、低域通過フィルタリング（平滑化）と、サイズが半分の画像をもたらす2分の１の、後続する解像リダクション（サブサンプリング）を意味する。ラプラシアン・ピラミッド表現Λjは、段入力表現と、リダクションRブロックとエクスパンションEブロックを通過した後のその複製との間の差異として定義される。「エクスパンション」Eは、本明細書及び特許請求の範囲では、（ゼロを挿入することによって）2倍の解像増加と、後続する低域通過フィルタリング（補間）とを有する。この場合、3x3の二項フィルタがリダクションR及びエクスパンションEにおける低域通過フィルタリング処理に用いられる。ラプラシアン・ピラミッド表現Λjはよって、解像段jの高域通過部分を有し、ガウシアン・ピラミッド表現Γjは解像段jの低域通過部分を有する（B.Jahne,  Digitale  Bild-verarbeitung〔Digital  Image  Processing〕,  5th  edition,  Springer  Verlag  Berlin  Heidelberg,  2002,  Section  11.4,5.3参照。）
グラディエント・フィルタリングに備えて、隣接画素間の単純な差異の形成によって、グラディエントΔがラプラシアン・ピラミッド表現Λjから更に計算される。この場合における当該差異は差異形成に用いる画素間の、中心における位置に属する。更に、グラディエントは解像段jで計算されるが、先行する、より細かい解像段(j−1)でのフィルタリングに用いられる。これらの理由で、グラディエントは適切に補間されなければならない。最後に、結果がもう一度、２で除算され、より細かいサンプリングを補償する。帯域通過画像のモジュラスは元の画像における非連続性がある場合に最大値を呈するのみならず、その付近においても同様であるので、より粗い解像段j’＞jのグラディエントがブロックEにおいてエクスパンションが行われ、解像段jのグラディエントに付加される。
【００４５】
最も粗い解像段以外は、全てのラプラシアン・ピラミッド表現Λ0乃至Λ2は、上記のように計算されるグラディエントに適応的に反応するフィルタGAFを用いてフィルタリングされる。フィルタ合成の始点は、この場合、3×3の二項フィルタであり、このフィルタ係数
【００４６】
【数２５】
2006510411.tif 000018
は表現構造の主方向に沿って維持される一方、これらの構造に対するグラディエントの方向における係数は
【００４７】
【数２６】
2006510411.tif 000019
の式によって単純化され、
【００４８】
【数２７】
2006510411.tif 000020
はフィルタの新係数であり、
【００４９】
【数２８】
2006510411.tif 000021
はフィルタリングする対象の画像位置であり、
【００５０】
【数２９】
2006510411.tif 000022
はフィルタ・コアの中心から当該係数に指し向けられるベクトルであり、
【００５１】
【数３０】
2006510411.tif 000023
は元のフィルタ係数であり、rは重み付け関数であり、
【００５２】
【数３１】
2006510411.tif 000024
は像点
【００５３】
【数３２】
2006510411.tif 000025
におけるグラディエントである。重み付け関数rは、グラディエントと、係数方向
【００５４】
【数３３】
2006510411.tif 000026
とのスカラ積とともに指数関数的に、
【００５５】
【数３４】
2006510411.tif 000027
に示すように降下し、c、t及びLはユーザによって定義することが可能なパラメータであり、
【００５６】
【数３５】
2006510411.tif 000028
は、グラディエント・フィールドの雑音の推定分散である。
【００５７】
図1において見られるように、グラディエント適応フィルタGAFの計算は、表現Γjを備えている、既に処理され、再構築されているガウシアン・ピラミッドを前提としている。
【００５８】
図1の図の右側部分は、連続する加算及びエクスパンションEによる（未修正状態にあるか、フィルタリングGAFによって修正された）詳細画像Ajから出力画像Aを合成することを反映している。詳細画像が何らフィルタリングされなかったならば、出力画像Aは入力画像Iに同一となる。
【００５９】
上記MRGAFアルゴリズムの欠点は、今日まで、このアルゴリズムが、計算労力が多いという理由で、記憶画像又は記憶画像シーケンスに対してオフラインでしか行うことが可能でないという点である。しかし、このアルゴリズムによって達成することが可能なかなりの画像の向上が理由で、これを、例えば、実施中の医療介入中にリアルタイムでも行うことができることが望ましいものである。この目的は、種々の最適化を用いて下記の方法で達成されるが、特に、「スーパー行」として知られているものを処理する手法によって達成される。この処理原理は、例として本明細書で検討されるMRGAFアルゴリズムの場合に用いることが可能であるのみならず、原則として、全ての種類の多重解像において、かつ、例えば「サブバンド符号化」などの別の同等なアルゴリズムとともに、用いることが可能である。
【００６０】
上記の元のMRGAFアルゴリズムはレベル単位でデータを処理する。まず、入力画像Iが低域通過フィルタリングされる。画像は通常、プロセッサ（キャッシュ）による高速アクセスによって（バッファ）メモリに入力されるには大き過ぎるので、入力データの一部と処理データの一部を主メモリすなわちシステム・メモリ（ワーキング・メモリRAM及び/又は、例えばハード・ディスクなどのバルク・メモリ）との間での読み書きを行わなければならない。しかし、これは２つの理由で欠点がある。第１は、システム・メモリに対するアクセスは比較的遅いという点であり、第２は、技術的な点で、速度における増加に関して、データ処理ハードウェアほど速く進展しないという点である。したがって、メモリ・アクセスの数が削減されるというような方法でMRGAFアルゴリズムを変更することが図られている。
【００６１】
システム・メモリに対する読み書き処理の数を削減するために、各解像段jでの、ガウシアン・ピラミッド表現の計算、ラプラシアン・ピラミッド表現の計算、更にはグラディエント・フィールドjの計算はお互いに組み合わせられるので、これらの値は段入力画像に及ぶ単一パスにおいて局所で計算される。この目的で、次に粗いガウシアン・ピラミッド表現Γj+1の低域通過フィリタリング値及び解像リダクション値をまず計算しなければならない。このことを最も速く行う方法の略図を図2に示す。3x3の二項低域通過フィルタ(1,2,1;2,4,2;1,2,1)・1/16を用いるかわりに、乗算演算と加算演算を、x方向とy方向とにおいて（すなわち、行方向と列方向とにおいて）連続して1次元の低域通過フィルタ(1,2,1)・1/4を用いることによって節減することが可能であり、更に、中間値biをバッファリングすることによって節減することが可能である。特に、Γj+1（図2における点）から値を計算する、図2に示すアルゴリズムは以下のように進む。
【００６２】
x0,0      x0,1      x0,2
x1,0      x1,1      x1,2
x2,0      x2,1      x2,2
を現行位置x1,1の周りの3x3の近接とし、
ｂ1＝x0,0+  x0,1+  x0,1+  x0,2
を先行する行からバッファ記憶される値とする。
【００６３】
更に、以下の工程が行われ、v1、v2は一時変数であり、ｂ0、ｂ1、ｂ2、…はバッファ変数である。
【００６４】
1.  v1＝x1,0+  x1,1+  x1,1+  x1,2
2.  v2＝x2,0+  x2,1+  x2,1+  x2,2
3.  結果=1/16*(ｂ1+  v1+  v1+  v2)  →  Γj+1に入力する
4.  ｂ1＝ｂ2
（など）
である。
【００６５】
更に、図3に示すように、ガウシアン・ピラミッド表現Γj+1の結果値は、各場合において、ラプラシアン・ピラミッド表現Λjと、x方向におけるグラディエント・フィールドΔxと、y方向におけるグラディエント・フィールドΔy（図3においてマーキングされている点を参照。）における4つの値を計算するのに直接用いられる。ラプラシアン値はこの場合、
ガウシアン・ピラミッド表現Γj+1の現行値の減算と、Γj+1における現行位置の左と上との既に計算された近接値に対してそれと補間される値の減算と、
ガウシアン・ピラミッド表現Γjの相当する値の減算から得られる。
【００６６】
グラディエント値は、ガウシアン・ピラミッド表現Γj+1の現行位置（図3における短破線）の左と上との既に計算された値と、グラディエント・フィールドにおける既に計算された値を用いた補間値との差異である。同時に行われる解像度増加のおかげで、計算差異を正しい「中間位置」に設定することが可能である。
【００６７】
GAFルーチンに必要なデータの上記メモリ最適化計算を用いれば、MRGAFアルゴリズムを既に、約15%速く行うことが可能である。更にかなりの、効率における増加が最小可能データ量によって全体の解像が行われるという革新によって達成されるので、この場合に必要であるデータをメモリに高速アクセス（キャッシュ）によってバッファリングすることが可能である。その場合、入力画像及び出力画像の読み書きの処理しか、なお必要な低速システム・メモリに対してアクセスしない。
【００６８】
メモリ・アドレスでの読み書きコマンドは常に、全体の後続データ・ブロックの、キャッシュとの間の読み書きにつながるので、完全な行の処理が確保される。しかし、1回に全体の入力画像Iは処理されず、できるだけ少ない行が処理される。これは、図4に表すように、最も粗い解像段のガウシアン・ピラミッド表現Γ3は、補間と差異計算に必要な、先行する行とともに単一の行を備えるだけである。したがって、ピラミッド構造が理由で、2K行の「スーパー行」は同時に処理されなければならず、Kは最大解像段である（例えば、図4、5ではK=3である。）。
【００６９】
図3は、最も粗い解像段での、ラプラシアン・ピラミッド・ブロックの計算とグラディエント・ブロックの計算との表現としてみることが可能である。ブロックは、2行に加えて、補間するよう別の先行する行を備える。図3にみられるように、変位はリダクション、再エクスパンション及び補間が理由で生じる。相対行0及び1が入力として与えられる場合、グラディエント適応フィルタリングGAFを、ラプラシアン・ピラミッド・ブロックの行−１及び−２に対してのみ行うことが可能である。これは、yグラディエントの結果データの位置と、3X3のGAFフィルタ・コアによってフィルタリングすることは、フィルタ位置の両側で別のデータの画素を必要とするということとの理由からである。この別のデータは、ラプラシアン・ピラミッド・ブロックの行0並びに−3（図3に図示せず）及び第１列並びに最終列である。
【００７０】
図4は、変位効果が別の解像段でどのように応答するかを示す。フィルタリング領域（濃灰色）は常に、現行位置よりも2行上の行に位置し、ラプラシアン・ピラミッド表現Λj（淡灰色）は現行位置よりも1行上の行に位置するということが分かり、後者は２つの先行する行によって最上部でエンスパンションされて3x3のフィルタリング処理を可能にする。
【００７１】
当該方法の最後の工程では、画像ブロックはフィルタリングされたラプラシアン・ピラミッド表現から再構築される。図5に示すように、フィルタリング・データの2行分の変位は合成工程中に加算され、これによって、再構築画像ブロックの変位は比較的大きくなる。しかし、このことはデータが誤った点で再書き込みされるということを意味するものでない。すなわち、全ての値は、それらが生じた場所に戻される。実際に、これは位置による変位ではなく、むしろ、因果関係条件が理由での時間による変位であり、それは、既に計算された値によってしか補間を行うことが可能でないことを意味する。図5における淡灰色領域はよって、合成まで継続しなければならない先行フィルタリング・データを区別している。しかし、この点で、かなりの量のバッファ・データを、フィルタリング工程と合成工程とを単に入れ替えることによって節減することが可能である。まず、合成になお必要な3つの行のみ（解像段2からの2行と段１からの1行、図5における濃灰色）がフィルタリングされる。更に、合成が行われるので、再構築バッファにおけるデータを残りのフィルタ値（濃灰色）によって上書きすることが可能である。この入れ替えによって、例えば段0の再構築バッファは９つの別の先行する行を確保しなくてよく、むしろ１つだけ確保すればよい。この数は、更に多くの解像段が用いられる場合には、3,5,7,…となる。
【００７２】
以下の計算は、上記方法におけるバッファリング・データの全体メモリ要件を推定する。この計算では、画像幅は512画素であり、3段ピラミッドが用いられ、全ての計算について、4バイトの浮動小数値が用いられる。
ガウシアン・ピラミッド：4*〔512*8+256*(4+1)＋128*(2+1)＋64*(1+1)〕＝23552バイト
ラプラシアン・ピラミッド：4*〔512*（8+2）+256*(4+2)＋128*(2+2)〕＝28672バイト
グラディエント・フィールド：2*4*〔512*（8+1）+256*(4+1)＋128*(2+1)〕＝50176バイト
合成バッファ：4*〔512*（8+1+1）+256*(4+1)＋128*(2+1)〕＝27136バイト
----------------
Σ129536バイト
である。
【００７３】
上記計算は、全ての計算を第２レベルのキャッシュを用いて、それが２５６ｋB=262144バイトの通常サイズを有する場合に行うことが可能である。元の画像の別の19行分の空間もなお、結果がそれに再書き込みを行なう対象である場合に存在する。
【００７４】
4*19*512＝38912バイト
上記方法を比較する場合、MRGAFアルゴリズムの元のバージョンを以下のように概説することが可能である。
【００７５】
ピラミッドの全レベルについて：
1.  x方向における段入力画像のリダクション  →  バッファ
2.  y方向におけるバッファのリダクション  →  ガウシアン・ピラミッド表現
3.  y方向におけるガウシアン・ピラミッド表現のエクスパンション  →  バッファ
4.  x方向におけるバッファのエクスパンション  →  第２バッファ
5.  段入力画像からの第２バッファの減算  →  ラプラシアン・ピラミッド表現
これに対して、新たなアルゴリズムの場合、全てが第２レベル・キャッシュからのデータを用いて画像の単一パスだけの間に行われる。ピラミッド解像、グラディエント計算、適応フィルタリング及び画像合成は、できるだけ狭い画像ストリップに対して行われる。
【００７６】
一時バッファ・メモリのサイズは、ガウシアン・ピラミッドの最も粗い段が、補間する先行する行とともに一行のみを備えるという要件から導き出される。この状況は、既に処理された行を用いてのみ可能である補間による再エクスパンション全てが理由で、読み取り画像ブロックの下限より2行上の行までに対してしか行うことが可能でない（図4参照、yグラディエントは最後の2行については計算することが可能でない。）ということによって複雑になる。2画素の高さのラプラシアン・ピラミッドの最も粗い段の場合、このことは、先行ブロックしかフィルタリングすることが可能でないということを意味する。再構築中には。この変位は段毎に増える。このことは、大量のデータ（少なくとも1ブロックのサイズ）しか、各工程におけるバッファにおいて確保しなくてよく、変位させなくてよい。しかし、「フィルタリング」と「再構築」とを入れ替えることによって、複製処理の数をありがたいことに削減することが可能である。アルゴリズムは図6に略図を示し、このアルゴリズムは以下の工程を備える。
【００７７】
サイズが2＾（解像段数）×（画像幅）の画像ストリップ全てについて：
1.  解像段全てについて：
第２工程におけるx方向とy方向とにおける段画像ストリップをパスする工程で、以下が各位置で行われる：
x方向とy方向とにおける低域通過フィルタリング  →  ガウシアン・ピラミッド表現の１画素
書き込み画素をもう一度、（補間するよう、その近接画素を用いて）４つの画素に対してエクスパンションし、段入力表現の画素からのその直接減算を行う工程  →  ラプラシアン・ピラミッド表現の４画素
ガウシアン・ピラミッド表現の計算画素とその近接画素との間の差異を計算する工程、先行して計算されたグラディエントによる結果を補間する工程  →  xグラディエント・フィールド及びyグラディエント・フィールドにおける４画素
2.  ラプラシアン・ピラミッドの最も粗い段の最後の２つの行と、2番目に粗い段からの第１行とをフィルタリングする工程（これらの行は再構築になお必要であり、別のものは既に利用可能である）  →  再構築バッファ
3.  再構築ピラミッド・バッファから画像ストリップを再構築する工程
4.  最も粗いもの以外の段全てについて：
再構築バッファにおける最下部から最上部までの次工程に必要なデータを複製する工程
現行のラプラシアン・ピラミッド・ストリップをフィルタリングする工程  →  再構築バッファ
以下の記載では、更なる高速化に寄与するアルゴリズムの一連の精緻化を記載する。
【００７８】
図1と、本発明によるMRGAFアルゴリズムを示す図6との比較から、重要な違いが、
新アルゴリズムにおいて、グラディエント・フィールドが各段で低域通過フィルタリング段入力表現から直接計算されるということにおいて分かる。次に粗いピラミッド段がフィルタリングされるまでアルゴリズムはもう待たなくてよいので、フィルタリングはこの場合、全解像レベルで並列に行われる。
【００７９】
式(4)に示すように、元のMRGAFアルゴリズムでは、指数関数が、各フィルタ位置でのフィルタ係数毎に、かつ各解像段で計算される。この点で、関数exp（−x）を1/（1+x）の近似によって置き換えることを提案するが、これによって、かなり低い計算労力に対して同様なプロファイルを備える。このようにして、新フィルタ係数の式(2)に帰結する。
【００８０】
元のMRGAFアルゴリズムでは、式(3)に示すように、各フィルタ係数は２つの因子rによって重み付けされ、そのうちの一方は現行画像位置
【００８１】
【数３６】
2006510411.tif 000029
におけるグラディエントによって計算され、他方は、係数位置
【００８２】
【数３７】
2006510411.tif 000030
におけるグラディエントによって計算される。よって、3x3のフィルタ・コアの場合、フィルタリング画素毎に９つの異なるグラディエント因子を計算しなければならない。全てのフィルタ位置でのグラディエント計算のおかげで、曲線の処理を向上させなければならない。しかし、3x3のフィルタを用いる場合、この手順は必要でないことが明らかであるがそれは、粗いピラミッド段から補間される隣接グラディエントはお互いに非常に類似しているからである。よって、式(3)の代わりに、式(1)に示す単純化された式を提案する。フィルタ・ルーチンの計算はこのことが理由でかなり単純化されるが、それはその場合、反対のフィルタ係数が同じ値を有し、式(2)が、9回ではなく、1回だけしか計算しなくてよいからである。
【００８３】
式(2)の分母におけるグラディエント・フィールドの雑音の分散
【００８４】
【数３８】
2006510411.tif 000031
は、粗いピラミッド段の相当する画素の雑音の分散によっておおよそ置き換えられる。フィルタ結果の品質はこれによって損なわれない。
【００８５】
図6から、グラディエント計算にガウシアン・ピラミッド表現を用いる場合、（Γ3及びΛ3を備えている）最も粗いピラミッド段は不必要なものであるということが分かる。この段の計算はよって、行わなくて済ますことが可能である。3段フィルタリング処理はよって、4段フィルタリング処理によって上記で達成されたものと同じ結果をもたらす。
【００８６】
1.7GHzを備えており、インテル(Intel)社のC++コンパイラを用いたコンパイルを備えているデュアル・ゼオン(Xeon)・ペンティアム（登録商標）（Pentium）4上で、最も好適な場合には、上記の単純化を考慮したプログラムの実施によって、毎秒43.6画像（512×512）（すなわち、毎秒、30を超える(768×564の)画像）に相当する、1つの画像について0.0229秒の実行時間をもたらした。当該方法はよって、リアル・タイム・アプリケーションに適切であるような程度に達している。
【図面の簡単な説明】
【００８７】
【図１】先行技術による、MRGAFアルゴリズムのシーケンスを示す図である。
【図２】次に高い解像段のガウシアン・ピラミッド表現の生成中での低域通過フィルタリング及び解像リダクションにおける変数の利用を示す図である。
【図３】ラプラシアン・ピラミッド表現の計算と、２つの連続するガウシアン・ピラミッド表現からのx方向及びy方向におけるグラディエント・フィールドとを示す図である。
【図４】種々の解像段における像点の位置を示す図である。
【図５】種々の合成段における像点の位置を示す図である。
【図６】本発明による、MRGAFアルゴリズムのシーケンスを示す図である。
【図１】
2006510411.tif 000041
【図２】
2006510411.tif 000042
【図３】
2006510411.tif 000043
【図４】
2006510411.tif 000044
【図５】
2006510411.tif 000045
【図６】
2006510411.tif 000046
【国際調査報告】
2006510411.tif 000047
2006510411.tif 000048
2006510411.tif 000049