JP2014520641A - 時間フィルタリングおよび空間フィルタリングによる一連の蛍光透視画像における雑音を低減させるための方法 - Google Patents

Abstract

本発明は、蛍光透視法を使用する医用イメージングの分野にある。本発明は、Ｘ線検出器によって得られる一連の蛍光透視画像における雑音を低減させるための方法に関する。
本発明による方法は、一連の各画像（ｘ_ｎ）について以下の連続ステップ、すなわち、
− 瞬間ｎに得られる画像（ｘ_ｎ）に時間フィルタリング（２２）を適用するステップであって、前記時間フィルタリングが、蛍光透視雑音低減（ＦＮＲ）処理を含むステップと、
− 瞬間ｎに得られる画像（ｙ_ｎ）に空間フィルタリング（２３）を適用するステップであって、前記空間フィルタリングが、以下のサブステップ、すなわち、
− 瞬間ｎに得られる画像（ｙ_ｎ）を、カーブレット変換を使用して時間領域からカーブレット領域に変換するサブステップ（２３１）であって、変換された各画像が１組の係数によって表されるサブステップ、
− 閾値化関数を使用して画像（ｙ_ｎ）の係数を閾値化するサブステップ（２３２）であって、この閾値化関数が、第３の所定の閾値を下回る係数を無効化し、第３の所定の閾値を上回る係数を保存または調整するサブステップ、
− 係数が閾値化された画像（Ｚ_ｎ）を、逆カーブレット変換カーブレット領域から時間領域に変換するサブステップ（２３３）を含むステップとを含む。

Description

本発明は、Ｘ線イメージングの分野に関し、より詳細には、蛍光透視法を使用する医用イメージングの分野に関する。本発明は、Ｘ線検出器によって得られる一連の蛍光透視画像における雑音を低減させるための方法に関する。

医用イメージングの分野では、Ｘ線蛍光透視イメージングにより、たとえば、肝臓の化学塞栓、椎体形成、動脈瘤のカテーテル法、または血管狭窄の処置など、いわゆる低侵襲外科的処置を施している間、患者の一連の画像を医師に提示することが可能になる。これらの画像は、特に外科用器械を案内する助けとなる。蛍光透視法によって案内される処置は通常、患者の血管網にカテーテルを挿入することによって実行される。従来は、造影剤を血管網に注入して、この血管網を不透明にし、血管の可視性を改善する場合もあり、そうでない場合もあった。一般に、このような処置は相対的に長く、患者が曝されるＸ線の線量を制限して、損傷または組織の焼損の発生を回避しなければならない。使用されるＸ線の線量をこうして制限することから、得られる蛍光透視画像の雑音のレベルは相対的に高く、したがって、コントラスト・ノイズ比が相対的に低くなり、画像が読み取りにくくなる。

これらの画像に存在する量子雑音を低減させ、そのコントラスト・ノイズ比を増大させるために、フィルタリング処理が実施される。

したがって一般に、蛍光透視雑音低減すなわちＦＮＲ処理の名前で知られる処理が、Ｘ線検出器によって得られた画像に実行される。このＦＮＲ処理の目的は、画像に存在する情報におけるコントラストを保存しながら、この画像に存在する雑音をフィルタリングすることである。ＦＮＲ処理は、動きのない画像の領域に時間フィルタを適用することによって実行される。画像の領域において動きが存在するか存在しないかは、別々に考えられている画素それぞれの輝度における個々の変化に基づいて検出される。画素は、２つの画素間でのその輝度変化が、雑音の標準偏差に関連する閾値を超えるときに、動いていると考えられる。動いていると判定されたこれらの画素は、フィルタリングされないか、またはほとんどフィルタリングされない。一方では、固定されていると判定された画素には、巡回型時間フィルタリング処理が加えられる。ＦＮＲ処理は、雑音を空間的に色づけしない利点を示す。すなわち、ＦＮＲ処理により、輝度が似かよった画素のグループ化から生じる仮想対象が現れないようにする。しかし、ＦＮＲ処理は、雑音低減能力を相対的に制限し、目的とする対象を画像から消し去る傾向がある。

時間フィルタリングによるＦＮＲ処理と比較すると、空間フィルタリング処理は有益な特性を示す。空間フィルタリング処理の様々な変形形態の中で、画像をウェーブレット領域に変換することによるフィルタリングが、最も普通に使用される技法のうちの１つである。１９８０年代に紹介されたウェーブレット理論は、フーリエ変換から開発された。多くの物理現象を説明するのにフーリエ変換が使用される。具体的には、フーリエ変換は、多くの物理現象、特に定常状態の信号によってモデル化された物理現象を説明するための有効な数学的ツールを形成する。一方で、非定常状態の信号によってモデル化された物理現象の表現は、一般に満足できるものではない。具体的には、フーリエ変換は、周波数が突然変化する信号の各部分を局在化することができない。これを改善するために、窓フーリエ変換が開発された。窓フーリエ変換は、信号に所与の大きさの窓を乗算することにある。したがって、この分析は局所的になり、窓の大きさによって、得られる時間分解能、または２次元信号の場合の空間分解能が決定される。ウェーブレット理論は、窓フーリエ変換に基づいている。ウェーブレット理論は、変換中に窓の大きさを変更できるようにすることによって窓フーリエ変換とは異なる。こうした理由で、多重解像度分析または多重スケール分析といった用語が使用される。他の違いは、信号がもはや、サインとコサインの和には分解されず、「ウェーブレット」と呼ばれる関数に従うことである。ウェーブレットにより、画像をいわゆるスパース表現に変換する利点が得られる。スパース表現は、少数のパラメータで画像を忠実に表すことを必要とする表現を意味するものと理解される。しかし、ウェーブレットは、物体外形など線タイプの不連続性を表現するのに十分適していないが、それというのも、その形状によりウェーブレットが対象物の形状を使用できないからである。結果として、これら不連続性についての情報が、数多くのスケールの係数に含まれる。したがって、再構成誤差を小さくするには、多数の係数を計算することが必要であり、結果として計算コストが高くつく。画像を「リアルタイムに」、すなわち画像のリフレッシュ期間未満の期間内に処理しなければならない蛍光透視の用途では、この計算コストには著しいハードウェア・リソースを必要とする。

情報の冗長度を減らすために、他の理論が開発されてきた。これらには特に、等方性ウェーブレットおよび異方性ウェーブレット、あるいは双対木、直交、または他のものが含まれる。

２０００年代の初めに、カーブレット理論が出現した。この理論も、多重スケール分析の原理に基づく。カーブレットの幾何形状により、カーブレットは物体外形の規則性を利用できるようになる。したがって、画像は、相対的に少ない数の係数で表現することができる。しかし、カーブレット変換には、画像にアーティファクト、特に雑音の空間色付けが入り込むという欠点がある。アーティファクトの混入は、コントラスト・ノイズ比が低いときに強くなる。しかし、蛍光透視イメージングの分野では、コントラスト・ノイズ比はしばしば１未満である。したがって、アーティファクトの存在が不可避となり、その使用をこれまで、はるかに高いコントラスト・ノイズ比の画像に限定してきた。

本発明の１つの目的は、特に、コントラスト・ノイズ比が相対的に低い蛍光透視画像に適した雑音低減を提案することである。具体的には、雑音低減処理は、高い雑音低減能力を示し、計算コストを低く抑え、画像から対象となる物体の損失を制限し、いかなるアーティファクトも混入させないようにしなければならない。これを目的にして、本発明の主題は、Ｘ線検出器によって得られた一連の蛍光透視画像における雑音を低減させる方法であり、各画像が空間領域で画素の行列によって形成され、各画素が信号レベルを表す値を有し、この方法は、各画像に対して以下の連続したステップを含む。
− 瞬間ｎに得られる画像に時間フィルタリングを適用するステップであって、画像の各画素に対して以下のサブステップを有するステップ。
− 瞬間ｎに得られる画像内にあると考えられる画素の値と、時間フィルタリングが適用された、瞬間ｎ−１に得られる画像内の対応する画素の値との間の差を決定するサブステップ。
− この差が第１の所定の閾値を下回る場合、第１の組のパラメータを用いて、対応する画素の値を訂正するサブステップ。
− この差が第２の所定の閾値を上回る場合、第２の組のパラメータを用いて、対応する画素の値を訂正するサブステップ。
− 瞬間ｎに得られる画像に空間フィルタリングを適用するステップであって、以下のサブステップを含むステップ。
− 瞬間ｎで得られた画像を、カーブレット変換を使用して空間領域からカーブレット領域に変換するサブステップであって、カーブレット領域内の各画像が１組の係数によって表現されるサブステップ。
− 閾値化関数を使用して画像の係数を閾値化するサブステップであって、第３の所定の閾値を下回る係数を無効化し、第３の所定の閾値を上回る係数を保存または調整するサブステップ。
− 係数が閾値化された画像を、逆カーブレット変換を使用して、カーブレット領域から空間領域に変換するサブステップ。

瞬間ｎに得られる画像をカーブレット領域に変換するステップは、有利には、６つもしくは９つのスケール、および／または１６のオリエンテーションを有するカーブレット変換を使用する。

第１の具体的な実施形態では、瞬間ｎに得られる画像をカーブレット領域に変換するステップは、ＵＳＦＦＴを用いた離散カーブレット変換を使用する。

第２の具体的な実施形態では、瞬間ｎに得られる画像をカーブレット領域に変換するステップは、ラッピングを用いた離散カーブレット変換を使用する。

空間フィルタリングの閾値化関数は、ハード閾値化関数、すなわち各係数が第３の所定の閾値を下回る場合に無効化され、閾値を上回る場合に保存される拠り所となる閾値化関数とすることができる。

空間フィルタリングの第３の所定の閾値Ｔ_ｓは、以下の連続したステップによって決定することができる。
− 瞬間ｎに得られる画像と同じ大きさであり、全ての画素にわたって均一な値１である画像Ｙを作成するステップ。
− フーリエ変換を使用して画像Ｙを変換するステップ。
− 次の関係式を使用して、画像Ｙのフーリエ変換でのピークの振幅Ｆを正規化するステップ。

ここで、ＭおよびＰは、画像Ｙの大きさを示す。
− 画像Ｙのフーリエ変換を、カーブレット変換を使用して空間領域からカーブレット領域に変換するステップであって、カーブレット領域内の各画像が１組の係数

によって表現されるステップ。
− 以下の関係式を使用して、係数

についてノルム

を決定するステップ。

ここで、ＭおよびＰは、画像Ｙの大きさを示し、ｌはカーブレットのオリエンテーションを示す。− 以下の関係式を使用して、所定の閾値Ｔ_ｓを決定するステップ。

ここで、σは、瞬間ｎに得られる、時間フィルタリングが適用される画像における雑音の標準偏差であり、ｒは、カーブレットの最も粗いスケールにおいて値０をとり、最も細かいスケールおいて値４をとり、その他のスケールにおいて値１をとる。

具体的な一実施形態では、瞬間ｎに得られる画像に時間フィルタリングを適用するステップは、瞬間ｎに得られる画像内の対象となる画素の値と、時間フィルタリングが適用された、瞬間ｎ−１に得られる画像内の対応する画素の値との間の差を決定するサブステップに続いて、以下のさらなるサブステップを含む。
− この差が、第１の所定の閾値と第２の所定の閾値の間にある場合、第３の組のパラメータを用いて、対応する画素の値を訂正するサブステップ。

要約すれば、本発明による雑音低減法は、引き続いて、カーブレット変換による時間フィルタリングおよび空間フィルタリングを有する。雑音の色付けに関しては、フィルタリングを実行する順序が重要である。これは、時間フィルタリングを先に実行しないと、カーブレット変換による空間フィルタリングがアーティファクトをもたらすからである。さらに、本発明は特に、空間フィルタリングについて直観的にどんな形状をも必要としないという利点を有している。

添付図面を参照してなされる以下の説明を読めば、本発明がよりよく理解され、他の利点が明らかになろう。

蛍光透視画像での物体のコントラストを決定するために考慮されるポイントの一例を表す図である。 本発明によって雑音を低減するための方法の実現可能なステップを表す図である。 ハード閾値化関数およびソフト閾値化関数を、曲線によってそれぞれ表す図である。 並進運動および回転運動する物体を示す一連の蛍光透視画像における雑音を低減させるためのこの方法の性能を、グラフによって示す図である。 空間フィルタリングの実行を可能にする、雑音を低減させるための方法のサブステップの一例を表す図である。 図５の空間フィルタリング用に使用される閾値化関数の一例を示す図である。

残りの説明において、Ｘ線検出器によって得られる一連の蛍光透視画像ｘ_ｎについて考えてみる。ここで、ｎは一連の画像における取得の様々な瞬間を示す。各画像は、信号レベルを表す値ｘ_ｎ（ｍ，ｐ）をそれぞれが有する画素の行列からなり、ｍおよびｐは、それぞれ行列の行および列を示す。値は、たとえば、０〜２５５の間の数値である。したがって、数値は８ビットで符号化される。数値は、１４ビットで符号化することもできる。コントラスト・ノイズ比（ＣＮＲ）は次式のように定義される。

ここで、雑音は、対象となる物体がない画像の領域にわたって決定されることが好ましく、コントラストは、画像内の対象となる物体に対応するポイントから決定することが好ましく、σは雑音の標準偏差、すなわち対象となる任意の物体を含まない画像の領域に対応する画素値の標準偏差を示す。図１には、「ステント」の名前でよく知られている血管内蔵式人工臓器のコントラストを確定していると考えられる（白色の）ポイントが示してある。

図２には、一連の蛍光透視画像における雑音を低減させるための方法の実現可能なステップが示してある。この方法は、瞬間ｎに得られる画像ｘ_ｎについて記述されている。しかし、この方法は、この一連の画像のそれぞれに連続して適用することができる。

第１のステップ２１で、一連の前処理によって画像ｘ_ｎが処理される。この前処理は、特に、画像内の隣接する画素をグループ分けするステップを含む。たとえば、各画素は、所望の解像度または所望の雑音レベルの関数としてグループ分けすることができる。前処理はまた、Ｘ線検出器の較正の段階を含むことができる。具体的には、検出器の各画素の適用される利得およびオフセットを決定することができる。さらに、この較正はまた、受けた放射に関わらず、ランダム値または一定値をとる画素の識別を含むことができる。これら欠陥のある画素は、隣接する画素の値を使用して、それらの値を補間することによって訂正することができる。前処理はまた、各画素の平均値が一連の全ての画像において一定になるよう、画像の再スケーリング、すなわち各画像の全画素の値の一括調整を含むことができる。最後に、この前処理は、画像を適応させるサブステップを含むことができる。Ｘ線の物理的性質が与えられている場合、各画像ｘ_ｎに存在する雑音はポアソン雑音である。ここで、このような雑音は、雑音低減処理には十分に適していない。さらに、各画像ｘ_ｎは、Ｘ線の減衰を表す画像であるが、放射線濃度の画像が、分析にはより自然な選択である。Ａｓｃｏｍｂｅの公式により、擬似ガウス雑音を用いて放射線濃度の画像を得ることが可能になる。

ここで、ｘ_０，ｎ（ｍ，ｎ）は、画像ｘ_ｎを適応させるサブステップの前の画素（ｍ，ｎ）の値である。

第２のステップ２２で、時間フィルタによって画像ｘ_ｎがフィルタリングされる。時間フィルタリングは、たとえば、「蛍光透視雑音低減」（ＦＮＲ）と呼ばれるフィルタリングである。このフィルタリングは、特に、米国特許第６，３１４，１６０Ｂ１号明細書に記載されている。ｙ_ｎを、時間フィルタリング後に得られる画像ｘ_ｎとする。Δ（ｍ，ｐ）は、瞬間ｎに得られる画像ｘ_ｎにおける値ｘ_ｎ（ｍ，ｐ）と、瞬間ｎ−１での取得に対応するフィルタリングされた画像ｙ_ｎ−１における値ｙ_ｎ−１（ｍ，ｐ）との間の差を示す。
Δ（ｍ，ｐ）＝ｘ_ｎ（ｍ，ｐ）−ｙ_ｎ−１（ｍ，ｐ）
表記を簡略化するために、ｘ_ｎは、単に画像ｘ_ｎの画素の値を示し、ｙ_ｎは、画像ｙ_ｎの画素の値を示し、Δは、対象となる画素の２つの値ｘ_ｎとｙ_ｎ−１の間の差を示す。ＦＮＲフィルタリングは、適応モードおよび非適応モードの２つのモードに分割することができる。適応モードは、物体の動きを考慮するので、そう呼ばれる。ＦＮＲフィルタリングは、画素毎に、設定閾値での差Δを比較することにある。第１の閾値Ｓ_１は、相対的に小さい差Δに対応する。第２の閾値Ｓ_２は、相対的に大きい差Δに対応する。次いで、この差Δおよび対象となるモードの結果に応じて、係数αおよびβによってＦＮＲフィルタリングがパラメータ表示される。αはフィルタの強度を示す。

ここで、σｘは、画像ｘｎでの標準偏差であり、σｙは、画像ｙ_ｎ−１での標準偏差である。
定数ＡＡ、ＢＢ、ＣＣ、およびＤＤは、以下のように定義される。

適応モードにおいては、閾値Ｓ_１およびＳ_２が以下のように定義される。
Ｓ_１＝σ
Ｓ_２＝２σ
係数αおよびβは以下のように定義される。
− ｜Δ｜＜Ｓ_１の場合、

− Ｓ_１＜｜Δ｜＜Ｓ_２の場合、

− それ以外の場合、

非適応モードにおいては、閾値Ｓ_１およびＳ_２は以下のように定義される。

係数αおよびβは以下のように定義される。
− ｜Δ｜＜Ｓ_１の場合、

− Ｓ_１＜｜Δ｜＜Ｓ_２の場合、

− それ以外の場合、

実際には、ＦＮＲフィルタリングは、画素の値ｘ_ｎならびに係数αおよびβに基づいて、画像の各画素についての新規の値ｙ_ｎを決定することにある。
ｙ_ｎ＝αｘ_ｎ＋（１−α）ｙ_ｎ−１＋β
すなわち、各画素の値ｘ_ｎは、係数αおよびβに基づいて訂正される。図２では、ＦＮＲフィルタリング・ステップ２２が適応モードにおいて表されている。一例として、第１のサブステップ２２１では、差Δの絶対値が、閾値Ｓ_１、すなわち標準偏差σで各画素ｘ_ｎについて比較される。差Δの絶対値が閾値Ｓ_１を下回る場合、対象となる画素が雑音レベルをほぼ表し、対応するパラメータを用いて、いわゆる強いフィルタリングがサブステップ２２２で適用される。そうでない場合、サブステップ２２３で、差Δの絶対値が、閾値Ｓ_２、すなわち２σで各画素ｘ_ｎについて比較される。差Δの絶対値が閾値Ｓ_２を上回る場合、対象となる画素は、移動した物体のポイントにほぼ対応し、対応するパラメータを用いて、いわゆる弱いフィルタリングがサブステップ２２４で適用される。最後に、差Δの絶対値が閾値Ｓ_１とＳ_２の間にある場合、不確実性が存在する。強いフィルタリングと弱いフィルタリングの間で連続関数を得るために、対応するパラメータを用いて、平均フィルタリングがサブステップ２２５で適用される。

第３のステップ２３で、時間フィルタリングの後で得られた画像ｙ_ｎが、カーブレット変換に基づく空間フィルタによってフィルタリングされる。カーブレット理論は、多重解像度空間または多重スケール空間の原理を使用し、スケールは空間の分割レベルに対応する。この空間は、包絡線の大きさが２^ｊ×２^ｊ／２であるカーブレットに分割され、ここで、ｊはスケールを示す正の整数である。各カーブレットは、対象となるスケールの数、画像におけるその位置、およびそのオリエンテーションに依存する。各スケールについて、画像の実現可能な最良近似が計算される。あるスケールからその他のスケールまで、得られた改善のみが表される。スケールの数が、画像の再構成の質を決定する。カーブレットを用いると、６〜１０の間の複数のスケールにより、コントラスト・ノイズ比に関して非常に満足する結果をもたらすことが可能になる。

空間フィルタリング・ステップ２３の第１のサブステップ２３１で、カーブレット変換により、画像ｙ_ｎが時間領域からカーブレット領域に変換される。離散カーブレット変換を使用することが好ましい。このような変換は、特に、Ｅ．Ｃａｎｄｅｓ、Ｄｅｍａｎｅｔ Ｌ．、Ｄｏｎｏｈｏ Ｄ．、およびＹｉｎｇ Ｌ．著「Ｆａｓｔ Ｄｉｓｃｒｅｔｅ Ｃｕｒｖｅｌｅｔ Ｔｒａｎｓｆｏｒｍｓ」、Ｍｕｌｔｉｓｃａｌｅ Ｍｏｄｅｌ．Ｓｉｍｕｌ．、ｖｏｌ．５、ｎｏ．３、２００６年、８６１〜８９９頁に記載されている。こうした離散変換では、フーリエ空間のグリッドを用いてカーブレットをサンプリングすることができない。したがって、Ｄｏｎｏｈｏらによる文書では、所与のスケールおよびオリエンテーションのカーブレットをこの空間グリッドに変換するための２つのソリューションが提案されている。画像を離散カーブレット領域に変換するための第１のソリューションは、ＵＳＦＦＴを用いた離散カーブレット変換と呼ばれ、ここでＵＳＦＦＴは、「不等間隔の高速フーリエ変換」を意味する。このソリューションでは、所与のスケールおよび角度でカーブレットの方向と向きを合わせるよう、この変換グリッドは傾斜している。第２のソリューションは、ラッピングを用いる離散カーブレット変換と呼ばれる。このソリューションでは、変換グリッドは、所与のスケールでの所与の象限の全ての角度について同一である。このことは、変換に使用される角度の数が４の倍数であることを意味する。ＵＳＦＦＴを用いるカーブレット変換は、ラッピングを用いるカーブレット変換よりも、連続カーブレット理論に忠実である。しかし、その計算コストはより高い。両方のソリューションにおいて、画像ｙ_ｎのカーブレット変換は１組の係数

で表すことができ、ここで指数Ｄは、離散変換を指し、ｊ、ｌ、およびｋはそれぞれ対象となるカーブレットのスケール、オリエンテーション、および位置を示す。空間フィルタの好ましい一実施形態では、使用されるスケールの数は９つに等しい。別の実施形態では、使用されるスケールの数は６つに等しい。後者の実施形態は、計算するための係数の数がより少なく、したがって処理時間がより短いことを意味する。さらに、使用されるオリエンテーションの数は、有利には１６に等しい。この数は、特に多数のオリエンテーションが接近することをその曲線が必要とする実像において、オリエンテーションにおける解像度とカーブレットの計算コストとの間の良好な妥協点を表す。さらに、数１６は４の倍数なので、ラッピングを用いるカーブレット変換を容易にする。

空間フィルタリング・ステップ２３の第２のサブステップ２３２では、画像ｙ_ｎの係数

が閾値化関数によって閾値化される。閾値化ステップ２３２は、ある種の係数

は雑音を表しているが、他は有用な情報を含んでいるという仮定に基づいている。したがって、閾値化ステップ２３２は、有用な情報を含まない係数を識別し、これを取り除くことにある。閾値化関数は、一方ではその形状によって特徴づけられ、他方ではその閾値Ｔで特徴づけられる。２つの主な閾値化関数、いわゆるソフト閾値化関数およびいわゆるハード閾値化関数が存在する。係数の表記をｃ_ｎに簡略化し、閾値化された係数

を呼び出すことにより、ハード閾値化関数を以下のようにモデル化することができる。

ソフト閾値化関数は以下のようにモデル化することができる。

図３Ａは、曲線３１によってハード閾値化関数を表しており、図３Ｂは、曲線３２によってソフト閾値化関数を表している。ソフト閾値化関数は、ハード閾値化関数よりも雑音低減能力が高い。一方で、ソフト閾値化関数は、選択された閾値Ｔに応じて画素の値を修正する。この結果、画像のある種のディテールにおける精度が失われる。閾値Ｔの決定は重要なステップであるが、それというのも、この決定により、アーティファクトの混入と同様に雑音低減能力においても、空間フィルタリングの性能が決定されるからである。

第１の実施形態では、閾値Ｔは、Ｄ．Ｌ．ＤｏｎｏｈｏおよびＪｏｈｎｓｔｏｎｅ，Ｉ．Ｍ．による「Ｉｄｅａｌ ｓｐａｔｉａｌ ａｄａｐｔａｔｉｏｎ ｖｉａ ｗａｖｅｌｅｔ ｓｈｒｉｎｋａｇｅ」、Ｂｉｏｍｅｔｒｉｃａ、ｖｏｌ．８１、１９９４年、４２５〜４５５頁に特に記載されているユニバーサル閾値Ｔ_Ｕであると決定される。この閾値Ｔ_Ｕは、以下の関係式によって定義される。

ここで、Ｎは、サンプルの数、すなわちここでは画像ｙ_ｎにおける画素の数である。ユニバーサル閾値Ｔ_Ｕは最適ではなく、アーティファクトを混入させることがある。それにもかかわらず、ユニバーサル閾値Ｔ_Ｕに比例する閾値は、比例定数の選択がよければ、良好な結果を得ることができる。

第２の実施形態では、閾値Ｔは、Ｊ．Ｌ．Ｓｔａｒｃｋ、Ｃａｎｄｅｓ Ｉ．Ｊ．およびＤｏｎｏｈｏ Ｌ．による「Ｔｈｅ Ｃｕｒｖｅｌｅｔ Ｔｒａｎｓｆｏｒｍ ｆｏｒ Ｉｍａｇｅ Ｄｅｎｏｉｓｉｎｇ」、ＩＥＥＥ Ｔｒａｎｓａｃｔｉｏｎｓ ｏｎ Ｉｍａｇｅ Ｐｒｏｃｅｓｓｉｎｇ、ｖｏｌ．１１、ｎｏ．６、２００２年、６７０〜６８４頁に記載の方法によって決定される。この方法は、全ての画素にわたって均一な値１である、雑音除去される画像ｙ_ｎと同じ大きさの、Ｙで示される画像を作成することにある。画像Ｙは、フーリエ変換によって変換される。画像Ｙにおける座標（０，０）の周波数が中心に位置し、その振幅が画像の大きさで正規化されている。Ｆを中心周波数の振幅とする。その正規化された振幅は｜Ｆ｜で示され、たとえば以下の関係式を使用して得られる。

ここでＭおよびＰは画像Ｙの大きさ、すなわちの画像ｘ_ｎおよび画像ｙｎの大きさを等しく示す。

したがって、画像は、画像の中心に位置しているディラック・ピークを用いて得られ、その振幅は画像の大きさの関数である。この方法の閾値はＴ_ｓで示され、スケールｊおよびオリエンテーションｌに依存する。これは、所与のスケールｊおよびオリエンテーションｌにおける画像Ｙのカーブレットの係数

のノルムの関数として決定される。このノルム

は、以下の関係式によって決定される。

閾値化関数は、以下で定義されるハード閾値化関数である。

乗法因子ｒは、スケールｊの関数として決定される。最も細かいスケール（たとえばｊ＝９）では、この因子ｒは４に等しい。その他のスケールでは、この因子ｒは３に等しい。これらの因子は、コントラスト・ノイズ比の高い画像において良好な雑音低減能力を示す。しかし、蛍光透視画像においては、アーティファクトが混入する。因子ｒを適応させることによって、これらのアーティファクトを回避することができる。具体的には、最も粗いスケール（ｊ＝１）では、因子ｒは０に等しくすることができる。したがって、閾値化が存在しない。最も細かいスケールでは、やはり因子ｒは４に等しくすることができる。また、その他のスケールでは、因子ｒは１に等しくすることができる。

空間フィルタリング・ステップ２３の第３のサブステップ２３３で、画像ｙ_ｎが、逆カーブレット変換によって、カーブレット領域から空間（行列）領域に変換される。このために、閾値化された係数が使用される。時間フィルタリングおよび空間フィルタリングの後に得られる画像をｚ_ｎで示す。

図４には、並進運動および回転運動を有するステントを表す一連の画像について、本発明による雑音を低減させるための方法の性能がグラフによって示してある。グラフにおいて、ｘ軸は、一連の画像ｘ_ｎの平均ＣＮＲを表し、ｙ軸は、雑音低減処理の後に得られる一連の画像の平均ＣＮＲを表す。第１の曲線４１は、本発明による雑音を低減させるための方法、すなわち時間フィルタリングとそれに続くカーブレット変換の空間フィルタリングによって得られる一連の画像ｚ_ｎの平均ＣＮＲを表す。第２の曲線４２は、同じ時間フィルタリングとそれに続くウェーブレット変換の空間フィルタリングによって処理された後の、同じ一連の画像の平均ＣＮＲを表す。グラフには、ＣＮＲが相対的に低い（０．５〜０．７〜１）処理すべき一連の画像について、本発明による雑音を低減させるための方法により、ＣＮＲに関して顕著な改善が得られることが示してある。ＣＮＲがより高い（１．５〜２〜２．５）一連の画像については、雑音低減がはるかに良好である。さらに、本発明による雑音を低減させるための方法により、アーティファクトの混入を著しく低減させることが可能になる。

前述の空間フィルタリング・ステップ２３により、画像内の雑音レベルを効果的に低減させることができる。しかし、このステップには、カーブレット、すなわち細長く波状の物体の形で現れるアーティファクトが混入するという欠点がある。特に、これらのアーティファクトは、振幅の大きい係数に起因するが、ただ雑音を表すだけである。これらのアーティファクトは、雑音レベルが高いときにのみ存在し、蛍光透視画像の場合がそうである。あるソリューションは、閾値を相対的に低く設定することにあり、これにより、アーティファクトは低減するが、雑音低減に関する性能が制限される。この説明の残りでは、フィルタリングの性能を落とすことなくアーティファクトの存在を著しく低減させることを可能にする、空間フィルタリング・ステップを実行するための別の方法を説明する。この方法は、所与の係数とそれに関係する係数、すなわちその親係数、その隣接係数およびそのいとこ係数との間に存在する相関性に依存する。具体的には、ある係数が有用な信号を表している場合、この係数に関係する各係数が、やはり有用な信号を表している可能性が高いことが確かである。逆に、アーティファクトが存在することによって係数の絶対値が高い場合、それに関係する係数もこのアーティファクトによって影響を受ける可能性はたいしてない。

図５には、空間フィルタリング・ステップ２３を実行することを可能にするサブステップの別の例が表してある。第１のサブステップ５１で、時間フィルタリングの後に得られる画像ｙ_ｎが、時間領域からカーブレット領域に変換される。画像ｙ_ｎのカーブレット変換はまた、

または簡単な形ｃ_ｎで示される係数によって表される。６つのスケールおよび１２のオリエンテーションでカーブレット変換がなされることが好ましい。具体的な実施形態では、画像ｙ_ｎは、最も細かいスケール、たとえばスケールｊ＝５でのウェーブレットによって表される。

第２のサブステップ５２で、Ｔ_１、Ｔ_２、およびＴ_３で示される３つの閾値が、Ｔ_１＜Ｔ_２＜Ｔ_３になるように決定される。これらの閾値は、係数のサブバンドの標準偏差σ_ｓ，θの関数として決定される。所与の係数

のサブバンドが、同じスケールｊおよび同じオリエンテーションｌを有する係数のセットをグループ分けする。したがって、これが、対象となるスケールおよびオリエンテーションについての位置のセットである。閾値Ｔ_１、Ｔ_２、およびＴ_３は、閾値Ｔ_１を下回る第１の範囲の値、閾値Ｔ_１とＴ_２の間の第２の範囲の値、閾値Ｔ_２とＴ_３の間の第３の範囲の値、および閾値Ｔ_３を上回る第４の範囲の値を定義する。一例として、対象となる以下の閾値の値が利用される。すなわち、Ｔ_１＝７．５、Ｔ_２＝１１．２５、およびＴ_３＝１５である。

第３のサブステップ５３で、各係数が、その親係数およびその隣接係数に応じて訂正される。相関関係のない隣接係数のみを考慮することが好ましい。簡略化された形で

または

が訂正された係数を示す。所与の係数

について、その親係数は、その前の最も粗いスケール、すなわち係数

に対応する係数である。その隣接係数は、隣接位置ｋ’、すなわち同じスケールおよび同じオリエンテーションを有する隣接した位置に対応する係数

である。対（ｘ，ｙ）が参照する位置において、したがってこれらは、位置（ｘ＋１，ｙ）、（ｘ−１，ｙ）、（ｘ，ｙ＋１）、（ｘ，ｙ−１）、（ｘ＋１，ｙ＋１）、（ｘ＋１，ｙ−１）、（ｘ−１，ｙ＋１）、および（ｘ−１，ｙ−１）である。サブバンドの方向に位置している隣接係数は「相関関係がある」と見なされる。各係数

について、その親係数の絶対値およびその相関関係のない隣接係数の絶対値の平均が、閾値Ｔ_１、Ｔ_２、およびＴ_３と比較される。３つの場合が考えられる。第１の場合は、対象となる係数

の絶対値が、その親係数

の絶対値と同じ範囲に位置するか、それとも、その隣接係数

の絶対値の平均と同じ範囲に位置するかの場合である。この第１の場合では、対象となる係数の絶対値が変更されないままである。

第２の場合は、親係数

の絶対値および隣接係数

の絶対値の平均が同じ範囲に位置しており、対象となる係数

の絶対値とは異なる範囲にある場合である。次いで、訂正された係数

の絶対値は、対象となる係数

の絶対値が低い方の値を有する範囲に位置している場合は、親係数の絶対値および隣接係数の絶対値の平均が位置している範囲の最小値をとり、そうでない場合には、この範囲の最大値をとる。第３の場合は、対象となる係数

の絶対値、その親係数

の絶対値、およびその隣接係数

の絶対値の平均が３つの別々の範囲に位置している場合である。次いで、訂正された係数

の絶対値は、対象となる係数

の絶対値が低い方の値を有する範囲に位置している場合は、隣接係数の絶対値の平均が位置している範囲の最小値をとり、そうでない場合には、この範囲の最大値をとる。第３の場合は、最も粗いスケールにおいて各係数が親をもたないという事実により、親係数の範囲を犠牲にして、相関関係のない隣接係数の範囲に有利に働く。それにもかかわらず、それらの親係数の絶対値が位置している範囲に応じて、その他のスケールの係数を訂正することを考えることも可能である。

第４のサブステップ５４で、訂正された係数

が、閾値化関数によって閾値化される。簡略化された形

または

で、閾値化され訂正された係数を示す。図６には、対象となる閾値の値Ｔ_１、Ｔ_２、およびＴ_３に対応する閾値化関数が表してある。この関数は以下のようにモデル化される。

ここで、

は、以下のように定義されるＳ字状である。

閾値化関数が連続であるように、パラメータｍおよびλが決定される。したがって、これらは以下のように定義される。

したがって、閾値化関数は、絶対値が第１の範囲６１にある係数を無効化し、絶対値が第２の範囲６２にある係数にソフト閾値化を適用し、絶対値が第４の範囲６４にある係数にハード閾値化を適用し、絶対値が第３の範囲６３にある係数に遷移閾値化を適用する。ソフト閾値化により、雑音レベルを低減させることが可能になるが、輪郭が広がる。ハード閾値化によって輪郭は良好に維持されるが、雑音低減に関しては効果的ではない。これら２つのタイプの閾値化を共用することにより、範囲に応じてそれぞれの利点を組み合わせることが可能になる。

画像ｙ_ｎが最も細かいスケール用のウェーブレットによって表されるとき、係数

は、以下のようにモデル化されたハード閾値化関数によって閾値化される。

ここで、閾値Ｔ_４は、対象となる係数ｃ_ｎのサブバンドの標準偏差σ_ｓ，θに基づいて定義される。
Ｔ_４＝４×σ_ｓ，θ

第５のサブステップ５５で、画像ｙ_ｎは、カーブレット領域から逆カーブレット変換による時間領域に変換される。このために、訂正された閾値化係数が使用される。時間フィルタリングおよび空間フィルタリングの後に得られる画像が、

で示してある。

図５を参照しながら説明した空間フィルタリングは、親係数と相関関係のない隣接係数との両方を使用する。空間フィルタリングはまた、相関関係のない隣接係数、または全ての隣接係数を使用することができる。空間フィルタリングはまた、親係数または隣接係数のみを使用することもでき、あるいは、これらの係数をいとこ係数と組み合わせることもできる。これらの係数は、スケールが同じで、位置が基準係数と同じ係数である。空間フィルタリングはまた、係数を訂正するどんなサブステップも含まないことがあり、この場合、閾値化サブステップ５４はカーブレット変換の係数ｃ_ｎに直接実行される。

Claims (13)

1. Ｘ線検出器によって得られる一連の蛍光透視画像における雑音を低減させるための方法であって、各画像（ｘ_ｎ）が画素の行列によって時間領域において形成され、各画素が信号レベルを表す値を有し、前記方法が各画像（ｘ_ｎ）について以下の連続ステップ、すなわち、
   − 瞬間ｎに得られる前記画像（ｘ_ｎ）に時間フィルタリング（２２）を適用するステップであって、前記画像の各画素について以下のサブステップ、すなわち、
   − 前記瞬間ｎに得られる前記画像（ｘ_ｎ）において考慮される前記画素の前記値と、前記時間フィルタリング（２２）が適用された瞬間ｎ−１に得られる画像（ｙ_ｎ−１）における対応する画素の値との間の差を決定するサブステップ（２２１、２２３）、
   − 前記差が第１の所定の閾値Ｓ_１を下回る場合に、第１の組のパラメータを用いて、前記対応する画素の前記値を訂正するサブステップ（２２２）、
   − 前記差が第２の所定の閾値Ｓ_２を上回る場合に、第２の組のパラメータを用いて、前記対応する画素の前記値を訂正するサブステップ（２２４）
   を含むステップと、
   − 前記瞬間ｎに得られる前記画像（ｙ_ｎ）に空間フィルタリング（２３）を適用するステップであって、以下のサブステップ、すなわち、
   − 前記瞬間ｎに得られる画像（ｙ_ｎ）を、カーブレット変換を使用して前記時間領域からカーブレット領域に変換するサブステップ（２３１、５１）であって、前記カーブレット領域内の各画像が１組の係数によって表されるサブステップ、
   − 閾値化関数を使用して、前記画像（ｙ_ｎ）の前記係数を閾値化するサブステップ（２３２、５４）であって、前記閾値化関数が、第３の所定の閾値を下回る係数を無効化し、前記第３の所定の閾値を上回る係数を保存または調整するサブステップ、
   − 係数が閾値化された前記画像（ｙ_ｎ）を、逆カーブレット変換を使用して前記カーブレット領域から前記時間領域に変換するサブステップ（２３３、５５）
   を含むステップと
   を含む、方法。
2. 前記瞬間ｎで得られる前記画像（ｙ_ｎ）を前記カーブレット領域に変換する前記ステップ（２３１）が、６つのスケールを有するカーブレット変換を使用する、請求項１に記載の方法。
3. 前記瞬間ｎで得られる前記画像（ｙ_ｎ）を前記カーブレット領域に変換する前記ステップ（２３１）が、９つのスケールを有するカーブレット変換を使用する、請求項１に記載の方法。
4. 前記瞬間ｎで得られる前記画像（ｙ_ｎ）を前記カーブレット領域に変換する前記ステップ（２３１）が、１６のオリエンテーションを含むカーブレット変換を使用する、請求項１〜３のいずれか一項に記載の方法。
5. 前記瞬間ｎで得られる前記画像（ｙ_ｎ）を前記カーブレット領域に変換する前記ステップ（２３１）が、ＵＳＦＦＴを用いる離散カーブレット変換を使用する、請求項１〜４のいずれか一項に記載の方法。
6. 前記瞬間ｎで得られる前記画像（ｙ_ｎ）を前記カーブレット領域に変換する前記ステップ（２３１）が、ラッピングを用いる離散カーブレット変換を使用する、請求項１〜４のいずれか一項に記載の方法。
7. 前記閾値化関数がハード閾値化関数であり、すなわち、前記係数が、前記第３の所定の閾値を下回る場合に無効化され、そうでない場合には保存される、請求項１〜６のいずれか一項に記載の方法。
8. 前記第３の所定の閾値Ｔ_ｓが、以下の連続ステップ、すなわち、
   − 前記瞬間ｎに得られる前記画像（ｙ_ｎ）と同じ大きさであり、全ての画素にわたって均一な値１である画像Ｙを作成するステップと、
   − フーリエ変換を使用して前記画像Ｙを変換するステップと、
   − 前記画像Ｙの前記フーリエ変換のピークの振幅Ｆを関係式
   

   

   を使用して正規化するステップであって、
   ＭおよびＰが前記画像Ｙの大きさを示すステップと、
   − 前記画像Ｙの前記フーリエ変換を、カーブレット変換を使用して前記時間領域から前記カーブレット領域に変換するステップであって、前記カーブレット領域内の各画像が、１組の係数
   

   

   によって表されるステップと、
   − 前記係数
   

   

   についてのノルム
   

   

   を関係式
   

   

   を使用して決定するステップであって、
   ＭおよびＰが前記画像Ｙの前記大きさを示し、ｌが前記カーブレットの前記オリエンテーションを示すステップと、
   − 前記所定の閾値Ｔ_ｓを決定するステップを、関係式
   

   

   を使用して決定するステップであって、
   σが、前記瞬間ｎに得られる、前記時間フィルタリング（２２）が適用された前記画像（ｙ_ｎ）における雑音の標準偏差であり、因子ｒが、前記カーブレットの最も粗いスケールにおいて値０をとり、最も細かいスケールにおいて値４をとり、その他のスケールにおいて値１をとるステップと
   によって決定される、請求項７に記載の方法。
9. 前記閾値化関数が、前記第３の所定の閾値Ｔ_１と第４の所定の閾値Ｔ_２の間に絶対値がある係数にソフト閾値化を適用し、絶対値が第５の所定の閾値Ｔ_３を上回る係数にハード閾値化を適用し、前記第４の所定の閾値Ｔ_２と前記第５の所定の閾値Ｔ_３の間に絶対値がある係数に前記ソフト閾値化と前記ハード閾値化の間の遷移閾値化を適用し、前記所定の閾値Ｔ_１、Ｔ_２、およびＴ_３がＴ_１＜Ｔ_２＜Ｔ_３のようになる、請求項１〜６のいずれか一項に記載の方法。
10. 前記遷移閾値化が、次式
    

    

    のように定義されるＳ字状によって定義され、
    ここで、Ｃ_ｎが閾値化する必要のある係数であり、前記パラメータｍおよびλが次式
    

    

    のように決定される、請求項９に記載の方法。
11. 空間フィルタリングを適用する前記ステップ（２３）がさらに、前記閾値化サブステップ（５４）に先立って、
    − 前記係数を、その親係数、その隣接係数および／またはそのいとこ係数に応じて訂正するサブステップ（５３）であって、所与の係数
    

    

    の親係数が、以前の最も粗いスケールｊ−１における所与の係数と同じオリエンテーションおよび同じ位置を有する係数
    

    

    であり、その隣接係数が、前記所与の係数の位置についての隣接した位置に対応し、同じスケールおよび同じオリエンテーションを有する係数
    

    

    であり、前記いとこ係数が、前記所与の係数と同じスケールおよび同じ位置を有する係数であるサブステップを含む、請求項９または１０のいずれかに記載の方法。
12. 第１の範囲の値（６１）が、前記第３の所定の閾値Ｔ_１を下回る値について定義され、第２の範囲の値（６２）が、前記第３の所定の閾値Ｔ_１と前記第４の所定の閾値Ｔ_２との間で定義され、第３の範囲の値（６３）が、前記第４の所定の閾値Ｔ_２と前記第５の所定の閾値Ｔ_３との間で定義され、第４の範囲の値（６４）が、前記第５の所定の閾値Ｔ_３を上回る値について定義され、前記係数を訂正する前記サブステップ（５４）が以下の訂正を適用し、
    − 対象となる前記係数
    

    

    の絶対値が、その親係数
    

    

    の絶対値と同じ範囲、またはその隣接係数
    

    

    の絶対値の平均と同じ範囲に位置しているとき、訂正された係数
    

    

    の絶対値が、前記係数
    

    

    の前記絶対値の値をとり、
    − 前記親係数
    

    

    の前記絶対値および前記隣接係数
    

    

    の前記絶対値の前記平均が、同じ範囲に位置しており、対象となる前記係数
    

    

    の前記絶対値の範囲とは異なるとき、前記訂正された係数
    

    

    の前記絶対値は、対象となる前記係数
    

    

    の前記絶対値が低い方の値を有する範囲に位置している場合は、前記親係数の前記絶対値および前記隣接係数の前記絶対値の前記平均が位置している範囲から最小値をとり、そうでない場合には、この範囲から最大値をとり、
    − 対象となる前記係数
    

    

    の前記絶対値、その親係数
    

    

    の前記絶対値、およびその隣接係数
    

    

    の前記絶対値の平均が、３つの別々の範囲に位置しているとき、前記訂正された係数
    

    

    の前記絶対値は、対象となる前記係数
    

    

    の前記絶対値が低い方の値を有する範囲に位置している場合は、前記隣接係数の前記絶対値の前記平均が位置している範囲から最小値をとり、そうでない場合には、この範囲から最大値をとる、請求項１１に記載の方法。
13. 前記瞬間ｎに得られる前記画像（ｘ_ｎ）に前記時間フィルタリングを適用するステップ（２２）が、前記瞬間ｎに得られる前記画像（ｘ_ｎ）における対象となる前記画素の前記値と、前記時間フィルタリング（２２）が適用された前記瞬間ｎ−１に得られる前記画像（ｙ_ｎ−１）における前記対応する画素の前記値との間の差を決定するサブステップ（２２１、２２３）に続いて、以下のさらなるサブステップ、すなわち、
    − 前記第１の所定の閾値Ｓ_１と前記第２の所定の閾値Ｓ_２との間に差がある場合、第３の組のパラメータを用いて、前記対応する画素の前記値を訂正するサブステップ（２２５）を含む、請求項１〜１２のいずれか一項に記載の方法。

Images