JP4864909B2 - 画像処理装置 - Google Patents

Description

本発明は、画像処理装置に関し、特に、画像上に存在するノイズの低減を行える同方法及び同装置に関するものである。

現在、画像処理技術は、様々な分野で利用されるようになっている。

画像処理は、ビデオテープレコーダやデジタルカメラ等で取得される画像の劣化やその改質等に対処するために行われたり、構造物が設計通りに製造されているかを検査するために構造物のパターンや構造そのものを明瞭に把握すること等を目的として行われる。

Ｘ線ＣＴ装置、ＳＰＥＣＴ装置、ＭＲＩ装置等の種々の医用画像診断装置においても、様々な画像処理が施されている。血流ないし造影剤流の描出、あるいは病変部抽出や臓器等の輪郭抽出等を施すことについては、その効用が広く認められている。

画像処理技術は、ノイズ抑制技術、特徴抽出技術、パターン認識技術等の各種要素技術からなり、それぞれの技術単独で、又は適宜組み合わせて利用される。なお、このような要素技術の中でも特に、画像に含まれるランダムなノイズを低減する技術は、撮像や再構成等した物体をより鮮明に再現するためには、欠かせない技術である。

しかしながら、従来の画像処理技術、とりわけノイズ低減技術には、更なる改良が求められている。例えば、ノイズ低減技術としては、いわゆる「平滑化」が広く知られている。この平滑化とは、ある画素（ｉ，ｊ）について入力値ｆ（ｉ，ｊ）があったときに、この画素（ｉ，ｊ）近傍の平均濃度を当該画素（ｉ，ｊ）についての出力値ｇ（ｉ，ｊ）とするものである。具体的には、上記画素（ｉ，ｊ）近傍のｎ×ｎ画素を用いるとすれば、出力値ｇ（ｉ，ｊ）は、

として求められる。ただし、上式（１）におけるａ，ｂ，ｃ，ｄは整数である。また、上式（１）における１／（ｂ−ａ＋１）（ｄ−ｃ＋１）は、いわゆる重みと呼ばれるものである。なお、図１３は、ａ，ｂ，ｃ，ｄ＝−１，１，−１，１である場合を示している。

ところで一般に、分散がσ^２であるような母集団の分布から、独立にとられたｎ個のサンプルの平均値を計算すると、当該平均値の分散がσ^２／ｎになることが知られているから、上記式（１）によれば、上記にいう「母集団」及び「その分散σ^２」が、それぞれ、各画素（ｉ，ｊ）の値が含む、ノイズに起因する成分を確率変数とした確率分布、及び、その分散に該当するために、各画素の値ｆ（ｉ，ｊ）におけるノイズの寄与分を低下させることができる。

しかしながら、これを単純に適用するのみでは、いわゆる「エッジぼけ」が発生し、画像の空間分解能が損なわれ、全体がぼやけた感じになってしまう。上述した医用画像を例として言えば、細密な血管構造をなるべく少ないノイズで描写したいという場合にあっても、上記（１）式によるノイズ抑制処理によれば、本来、血管構造を描写していない画素をも含めて平均化（平滑化）が行われるため、ノイズは抑制されるにしても、血管構造を表すコントラストも平滑化によって低下してしまい、細密な血管構造の描写が困難になる場合がある。

本発明は、上記事情に鑑みてなされたものであり、その目的とするところは、画像のぼけを生じることなく、ノイズを十分に抑制することができることをはじめ、その他の画像処理技術、例えばパターン認識技術等にも有効に貢献し得る画像処理装置を提供することにある。

本発明は、画像を構成する第１の画素と第２の画素との類似度を、検定によって判定するものであって、前記画像を含む時間方向に異なる複数枚の画像から判定する判定手段と、この判定手段による判定結果に基づいて、前記第１の画素と前記第２の画素とを重み付き平均する平均手段と、を具備したことを特徴とする画像処理装置を提供する。

本発明によると、画像のぼけを生じることなく、ノイズを十分に抑制することができることをはじめ、その他の画像処理技術、例えばパターン認識技術等にも有効に貢献し得る。

以下では、本発明の実施の形態について図を参照しつつ説明する。
本実施形態はまず、ある画像を構成する画素について、ベクトル値又はスカラー値た画素値を構成し、前記画素値及び該画素値とは別に構成された他の画素値間の適合度を定量化し、前記画素についての新たな画素値を、前記他の画素値を利用して構成するに際し、前記適合度が大きい場合には、当該他の画素の寄与を大きくし、前記適合度が小さい場合には、当該他の画素の寄与を小さくして、当該新たな画素値を構成することを特徴とする画像処理方法である。また、本発明は、前記「別に構成された他の画素値」が、前記画像を構成する他の画素に基づき構成されてよい。すなわちこの場合、画素値は、一の画素及び他の画素につき、それぞれ、ベクトル値又はスカラー値たる「一の画素値」及び「他の画素値」として構成されることになる。さらに、本発明は特に、前記適合度の関数である重み関数を、前記他の画素値のそれぞれについて求められた前記適合度に作用させて当該他の画素値のそれぞれの重みを決定し、次にこの重みを用いた当該他の画素値の重み付き平均を算出することによって、前記新たな画素値を構成するときには、前記適合度が大きい場合には前記重みを大きくすることで、前記他の画素値の前記重み付き平均における寄与を大きくし、前記適合度が小さい場合には前記重みを小さくすることで、前記他の画素値の前記重み付き平均における寄与を小さくして、当該新たな画素値を構成するものでもある。このような処理によれば、一の画素と「類似である」と判定される他の画素が重視されて、新たな画素値が構成されることになるから、従来のように、空間分解能を損なうようなことや、あるいは前記画像が動画像である場合には、時間分解能を損なうようなことがない。（重み関数及び重み）本発明では、前記重み関数を、前記適合度に関する非負の単調増加関数とすると好ましい。また、前記一の画素及び前記他の画素をそれぞれｘ及びｙ、前記一の画素値及び前記他の画素値をそれぞれｖ（ｘ）＝（ｖ_１（ｘ），ｖ_２（ｘ），…，ｖ_Ｋ（ｘ））及びｖ（ｙ）＝（ｖ_１（ｙ），ｖ_２（ｙ），…，ｖ_Ｋ（ｙ））とするとともに、構成すべき前記新たな画素値をｖ′（ｘ）＝（ｖ′_１（ｘ），ｖ′_２（ｘ），…，ｖ′_Ｋ（ｘ））とすれば、前記重み付け平均をとる処理が、

（ただし、Ｎ（ｘ）は前記他の画素が含まれる範囲を表す。）
と表される形態とすると好ましい（以下では、新たな画素値ｖ′（ｘ）を求める上式右辺等のような形式を、「本実施形態に係るコヒーレント・フィルタ」と称す。）。

（本実施形態の適用可能な画像処理技術の例）
本発明は、上記したような処理、すなわち前記新たな画素値を、例えば前記画像の全面に関して構成することにより、前記画像のノイズを低減すること、あるいはまた、前記画像におけるパターン認識を実施すること等の各種応用が可能であり、後述するようにこれら画像処理技術の各要素において、本発明は優れた効果を発揮する。

（適合度）
本実施形態においては、前記適合度が、前記一の画素値及び前記他の画素値をそれぞれ構成するスカラー値に対し、統計的検定法を適用した結果求められる危険率に基づいて定量化し得るようにしてもよい。なお、ここにいう「統計的検定法」とは、後述するように例えば、「χ二乗検定法」等を考えることができる。また一般に、ここで導入された危険率と前記適合度との関係は、特に、一方が増えれば他方が減るという場合を含む。すなわち、前記新たな画素値を構成する際に、危険率が大きく（適合度が小さく）なれば、前記他の画素に対する重みが小さくなり（構成すべき新たな画素値に対する寄与が小さくなり）、危険率が小さく（適合度が大きく）なれば、前記他の画素に対する重みが大きくなる（構成すべき新たな画素値に対する寄与が大きくなる）、という場合である。

（「他の画素」の選択範囲）
また、本実施形態は、前記他の画素が、前記一の画素の周囲における所定の領域から選択され、特に、前記所定の領域が、前記一の画素を中心とした近傍であるようにしてもよい。このように構成すると、一の画素に類似であると判定される他の画素は、通常、当該一の画素の周囲に存在する可能性が高いことにより、上記のような領域の選択を行うと、新たな画素値を有効に求めるに際し、無駄な演算を省くことが可能となる。

（画素値の構成法）
本実施形態では、前記画素値（一の画素値及び他の画素値の双方を指す。）の構成法を種々選択することが可能であり、より具体的に第一には、前記画像が複数枚の画像である場合には、前記画素値は、前記複数枚の画像の各々を通じた同一点の画素が有するスカラー値を並べたベクトル値とすることが可能である。

このとき、前記複数枚の画像は、動画を構成する複数枚の静止画像であってもよく、そのような場合において本発明は、前記新たな画素値を求める際に、空間分解能を損なわないばかりでなく、時間分解能をも損なわない。なお、このような場合における好適な適用例は、例えば、前記動画が、医用画像診断装置により取得されたダイナミックＣＴ像である場合等である。

また、画素値の構成法の第二には、前記画像が、同一被写体に関する複数種類の画像であって、上述と同様に、当該複数種類の画像の各々を通じた同一点の画素が有するスカラー値を並べたベクトル値とすることが可能である。なお、このような場合における好適な適用例は、例えば、前記画像がＳＰＥＣＴ装置又はＰＥＴ装置におけるマルチウインドウ撮影法により取得された画像であり、前記複数種類の画像とは、異なる放射性同位体からそれぞれ発せられたガンマ線に基づいた、それぞれ異なる画像である場合等である。また、別の好適な適用例は、例えば前記画像がカラー画像であり、前記複数種類の画像とは、光の三原色にそれぞれ分解して取得された、それぞれ異なる画像である場合等である。

さらに、画素値の構成法の第三には、前記画像は一枚の画像であり、当該画素値を、その一枚の画像から構成することも可能である。より具体的には例えば、当該画素値を、前記一枚の画像上のある画素が有するスカラー値と該画素の近傍に在る画素が有するスカラー値を並べたベクトル値として構成することが可能である。なお、このような画素値の構成法によれば、あらゆるデジタル画像に対し、本発明を適用することが可能である。

（本実施形態の画像処理を実現する好適な装置構成）
最後に、上記した画像処理方法を実現する装置構成としては、ある画像を構成する画素につき、ベクトル値又はスカラー値として構成された画素値及び別に構成された他の画素値との間の適合度を判定する適合度定量化手段と、前記画素についての新たな画素値を、前記他の画素値を利用して構成するに際し、前記適合度が大きい場合には、当該他の画素の寄与を大きくし、前記適合度が小さい場合には、当該他の画素の寄与を小さくして、当該新たな画素値を構成する画素値演算手段と、を有するものとすることが好ましい。

以下詳細に説明する。以下ではまず、本発明に係る「コヒーレント・フィルタ」の概要に関する、より一般的な説明（項目番号“Ｉ”。以下同様。）をし、その後、該コヒーレント・フィルタを各種画像処理に適用した例（項目番号“ＩＩ”及び“ＩＩＩ”〜“ＶＩＩＩ”）について順次説明することとする。

また、そのような各種適用例の説明に続いてさらに、本実施形態におけるコヒーレント・フィルタの最も一般的な形態及びその応用例（項目番号“ＩＸ”）について説明し、最後に、本実施形態の補足事項（項目番号“Ｘ”）についてまとめて説明することとする。

（ Ｉ 本発明に係るコヒーレント・フィルタの一般的説明 ）
まず、本発明に係る「コヒーレント・フィルタ」を説明するための、準備的な説明を行う。

（ Ｉ−１ 画素値ｖ（ｘ） ）
一般に、カメラ等の撮像手段を介して取得されたデジタル画像は、複数の画素（ｐｉｘｅｌ）から構成されている（あるいは、当該画像をそのような画素の集合として考えることができる。）。以下の説明では、当該画素の位置をベクトルｘ（すなわち座標値のベクトル）として表し、画素ｘが有する値（例えば濃淡を表わす数値）をＫ次元ベクトルとして表す。２次元画像の場合、画素ｘとは画像上における位置を表す座標値（ｘ、ｙ）を示す２次元ベクトルである。ある画素ｘについて定義される「画素値ｖ（ｘ）」を、

と表記する。なお、この（３）式の右辺における、ｖ_１（ｘ），ｖ_２（ｘ），…，ｖ_Ｋ（ｘ）それぞれを、以下では、画素ｘについての「スカラー値」と呼ぶことにする。

例えば、画像が「カラー画像」であるとき、各画素が、それぞれ三原色（赤，緑，青）の明るさ（スカラー値）を有することから、これら各画素の画素値ｖ（ｘ）は、その次元がＫ＝３のベクトルであると考えることができる（上記（３）式の右辺各項で、その添え字が例えば「赤」，「緑」及び「青」である場合を想定されたい。後述の（３．１）式参照。）。また例えば、画像がＫ枚の静止画像から構成される動画像であって、第ｎ番目の画像の各画素はスカラー値ｖ_ｎ（ｘ）を持つという場合には、Ｋ枚の静止画像上、共通の同一点（同一座標）の画素ｘの持つ画素値（スカラー値）を並べて構成される、Ｋ次元ベクトル値ｖ_ｎ（ｘ）＝（ｖ_１（ｘ），ｖ_２（ｘ），…，ｖ_Ｋ（ｘ））が以下で述べるベクトル値としての画素値である。

（ Ｉ−２ 適合度ないし危険率ｐ（ｘ，ｙ）と重みｗ（ｐ（ｘ，ｙ）） ）
上記画素ｘに対して、適当な画素の集合Ｎ（ｘ）を考える（この集合Ｎ（ｘ）は画素ｘを含んでよい。）。次に、Ｎ（ｘ）の要素であるそれぞれの画素ｙと、前記画素ｘとの間で、重みｗ（ｐ（ｘ，ｙ））を考える。この重みｗ（ｐ（ｘ，ｙ））は、次に記す性質を有する。

（ Ｉ−２−１ 適合度、ないし危険率ｐ（ｘ，ｙ） ）
まず、ｗ（ｐ（ｘ，ｙ））の値を左右する関数ｐ（ｘ，ｙ）の意味について述べる。このｐ（ｘ，ｙ）は、本発明にいう「適合度」を定量化する手段であり、一般的にいえば、画素ｘと画素ｙ∈Ｎ（ｘ）とが、何らかの意味でどの程度類似しているか（例えば、両画素ｘ及びｙの上記画素値ｖ（ｘ）及びｖ（ｙ）間に認められる統計的差異の程度）、を示す具体的数値を与える。

より具体的には例えば、ｐ（ｘ，ｙ）が小さな値を与えるときには、画素ｘと画素ｙとが、その画素値ｖ（ｘ）及びｖ（ｙ）間に「統計的に有意な差がなく（＝適合度が大きく）」、類似である可能性が高いと判断され、ｐ（ｘ，ｙ）が大きな値を与えるときには、「統計的に有意な差があり（＝適合度が小さい）」、の如く判断されるということである。

ところで、画素値ｖ（ｘ）及びｖ（ｙ）（ないしスカラー値ｖ_１（ｘ），…，ｖ_Ｋ（ｘ）及びｖ_１（ｙ），…，ｖ_Ｋ（ｙ））には、必ずノイズが含まれていると考えなければならない。例えば、画像がＣＣＤ撮像素子により取得された場合を考えると、それを構成する各画素については、素子内の暗電流や外界から入射する光量の不規則変動に起因するノイズ等が存在する。

このようなノイズは、一般に、全画素についてまちまちな値をとるため、画素ｘと画素ｙとが、仮に（外界における）同一物体を反映したものである場合であっても、実際に観測される画像上では、同一の値を持たないことがある。このことを逆にいえば、いずれも同一物体を反映した画素ｘと画素ｙにおいて、それぞれのノイズを除去した状況を仮に想定すれば、これらは該同一物体を表象するものとして画像上に表示され（＝そのように認識され）るし、また、両者は本来同一の（あるいはごく近い）画素値を有する。

そこで、上述したノイズの性質を踏まえ、上記のｐ（ｘ，ｙ）に関し、統計的検定法でよく知られている「帰無仮説」の概念を用いると、このｐ（ｘ，ｙ）については、具体的に次のように言うことができる。すなわち、帰無仮説Ｈ「画素ｘと画素ｙとはそれぞれのノイズを除去した場合に同一の画素値を有する」言いかえれば「ｖ（ｘ）＝ｖ（ｙ）。ただし、両画素のノイズに起因する差異を除く」を立てる（つまり、このような命題が成立する場合、「両画素ｘ及びｙが類似である（＝適合度が大きい）」と考える。）と、関数ｐ（ｘ，ｙ）は、この仮説Ｈを棄却する場合の危険率（あるいは、有意水準）であるということができる（この場合、ｐ（ｘ，ｙ）は、その値域が［０，１］であるような関数として定義される（ｐ（ｘ，ｙ）∈［０，１］）。）。

したがって、危険率ｐ（ｘ，ｙ）が大きい場合、すなわち棄却が誤りである危険性が大きい場合には上記仮説Ｈを満たす可能性が高いといえ、逆に小さい場合、すなわち棄却が誤りである危険性が小さい場合には仮説Ｈを満たさない可能性が高いということができる（なお、統計的検定における周知事項ではあるが、仮説Ｈが「棄却」されないといっても、それが「真」であることを意味するわけではない。この場合、仮説Ｈが示す命題が、否定し得ないことを意味するに過ぎない。）。

（ Ｉ−２−２ 重みｗ（ｐ（ｘ，ｙ）） ）
さて、重みｗ（ｐ（ｘ，ｙ））は、その表され方から明らかな通り、上記したような危険率ｐ（ｘ，ｙ）の関数（より一般には、適合度の関数（適合度をρ（ｘ，ｙ）とすれば、ｗ（ρ（ｘ，ｙ））となるように構成できる）であり、また、この重みｗ（ｐ（ｘ，ｙ））を求めるため、ｘ及びｙの組み合わせそれぞれについて求められた危険率ｐ（ｘ，ｙ）に作用させる重み関数ｗは、一般的にいうと、上記「棄却」を具現化する作用を有するものである。具体的には、危険率ｐ（ｘ，ｙ）が大きい場合には重み関数ｗの値、すなわち重みｗ（ｐ（ｘ，ｙ））が大きな正の値をとり、その逆の場合には小さな正の値（又は“０”）をとる、等というように調整されている（重み関数ｗの具体的形式については後述する。）。つまり、重みｗ（ｐ（ｘ，ｙ））は、画素ｘと画素ｙとが、上記仮説Ｈに示される命題を満たすらしい場合には、大きい値をとり、その逆の場合には小さい値をとる。一例として特に、ｗのとりうる値が”０”かまたは”０”でない一定値の２通りしかないように構成してもよい。

なお、以上までに述べた仮説Ｈ、危険率ｐ（ｘ，ｙ）、重みｗ（ｐ（ｘ，ｙ））間の関係を、図１にまとめて示しておく。また、重み関数ｗ（ｔ）は、より一般に、「ｔ∈［０，１］で定義される非負の単調増加関数」ということができ、また、該ｗ（ｔ）の満たすべき性質は、少なくともそのようであればよい。

（ Ｉ−３ コヒーレント・フィルタ ）
以上までの準備的説明により、本発明に係る「コヒーレント・フィルタ」は次のように導かれる。すなわちまず、画像を構成するある画素ｘに対し、集合Ｎ（ｘ）の要素たる画素ｙのすべてについて上記した重みｗ（ｐ（ｘ，ｙ））を計算する。次に、これら複数の重みｗ（ｐ（ｘ，ｙ））を用いて、当該画素ｘを構成する新たなスカラー値ｖ′_ｋ（ｘ）を、以下の式で計算する。すなわち、

ただし、ｋ＝１，２，…，Ｋである。そして、この式で求められたｖ′_ｋ（ｘ）を用いて、当該画素ｘの変換後の画素値（新たな画素値）ｖ′（ｘ）を、

として構成する。

ここに、上記（４）式で表される、画素値ｖ（ｙ）＝（ｖ_１（ｙ），ｖ_２（ｙ），…，ｖ_Ｋ（ｙ））（ｙ＝ｘである場合を含む。）を、ｖ′（ｘ）＝（ｖ′_１（ｘ），ｖ′_２（ｘ），…，ｖ′_Ｋ（ｘ））に変換するフィルタが、本発明に係る「コヒーレント・フィルタ」の、本実施形態における形式である。これはその表式から明らかな通り、画素値を構成するスカラー値ｖ_ｋ（ｙ）の重み付け平均値を表している。

このような処理は、以下のような結果をもたらす。すなわち、画素値ｖ′（ｘ）は、画素ｘとノイズを除いて同一の画素値をとることが確からしい（＝上記仮説Ｈの命題を満たす可能性が高い）画素ｙを重視した重み付け平均値ｖ′_ｋ（ｘ）から構成されたベクトルを表すこととなる。また、このような画素ｙが十分な数存在するならば、画素値ｖ′（ｘ）は、画素ｘが本来有すべきその真値から外れることなく、上記したような平均化の作用によりノイズのみを抑制した値を有することとなる。

なお、危険率ｐ（ｘ，ｙ）が小さく、したがって、帰無仮説Ｈが「棄却」され、重みｗ（ｐ（ｘ，ｙ））が小さくなるような場合であっても、上記記述又は（３）式の表式からもわかる通り、必ずしもこれを完全に「棄却」するとは限らない。このようなことは、後述する重み関数ｗの具体的形式に依存するところであるが、危険率ｐ（ｘ，ｙ）が“０”（＝０％）に近いような場合でも、ｗ（ｐ（ｘ，ｙ））≠０（ただし、ｐ（ｘ，ｙ）が“１”に近い場合に比べて、より小さな正の値ではある。）としてよい（なお、ｐ（ｘ，ｙ）＝１である場合とは、後述するように、ｖ（ｘ）＝ｖ（ｙ）のときである。）。すなわち、完全な棄却ということではなく、小さな寄与は認めてよいということである（なおこのような場合に、ｗ（ｐ（ｘ，ｙ））＝０とするのであれば、完全な棄却を行うのと同義である。後述の（１４）式参照）。

このような処理は、一般的に次のように言える。すなわち、ある画像を構成する（複数の）画素ｘが存在するとき、この画素ｘとある任意の画素ｙ（上記ではｙ∈Ｎ（ｘ）とされた。）との適合度を定量化し（上記では、ｐ（ｘ，ｙ）に基づいていた。）、該適合度が大きい場合には、画素値ｖ（ｙ）を利用した重み付き平均化処理において、当該画素ｙについて大きな寄与を認め、適合度が小さい場合には小さな寄与しか認めないようにすることで、当該画素ｘのノイズを有効に抑制する画像処理方法である、といえる。いわば、画素ｘと画素ｙとが「似たもの同士」のときには、該画素ｙを前記平均化処理に、より貢献させ、「似ていないもの同士」のときには、該画素ｙを殆ど又は全く無視する、と言い換えてもよい。

このような処理を画像全体に施すことにより、画像のぼけを殆ど生じることなく極めて高いノイズ抑制効果を発揮することができる。また、ノイズ抑制という用途に限定せず、例えばパターン認識の分野においても、重み関数、あるいはコヒーレント・フィルタを好適な具体的形式にすることによって、優れた効果を発揮することができる。なお、これらの効果が具体的にどのようなものであるかについては、以下各種画像処理に関する説明中において述べる。

（ ＩＩ 様々な適用例；画素値ｖ（ｘ）の構成法の差異による分類 ）
さて、本発明は、上記したコヒーレント・フィルタの一般的形態に基づき、様々な分野への適用が可能である。そして、これは例えば、図２に示すようにまとめることができる。なお、図２では、様々な分野における各種適用例を、画素値ｖ（ｘ）の構成法の差異に応じて分類したものを示している。すなわちこの図２によれば、上記したベクトル値たる画素値ｖ（ｘ）が、「複数の画像から」、（例えば一度の撮影等により得られる）「複数種類の画像から」、又は「１枚の画像から」の、それぞれの方式により構成することが可能であることがわかり、本実施形態に関する各種適用例が、これらそれぞれに分類されることがわかる。

以下では、この図２に示す分類表に沿って、上記のような一般的形態となるコヒーレント・フィルタを、より具体的な場面で適用した各種適用例について、また、当該図２に示す３種の画素値ｖ（ｘ）構成法が、具体的にどのように実現されるかについて、順次説明していく（ただし、図２に示す順とは異なる）。

（ ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（複数枚の画像からノイズを除去する））
まず、上記コヒーレント・フィルタを、Ｘ線ＣＴ装置における、いわゆる「ダイナミックＣＴ（ｄｙｎａｍｉｃ ＣＴ）」撮影に適用した場合について説明する。ここに、Ｘ線ＣＴ装置１００とは、図３に示すように、Ｘ線管１０１及びＸ線検出器１０２、データ収集部１０３（ＤＡＳ；Ｄａｔａ Ａｑｕｉｓｉｔｉｏｎ Ｓｙｓｔｅｍ）、前処理部１０４、メモリ部１０５、再構成部１０６及び画像表示部１０７、並びにこれら各部を制御する制御部１０８及びＸ線照射条件、撮影モード等その他各種の設定・入力が可能な入力部１０９等からなる。

このＸ線ＣＴ装置１００によれば、Ｘ線管１０１及びＸ線検出器１０２を被検体Ｐ周囲で回転させながらＸ線を被検体Ｐに曝射することにより、当該被検体Ｐの内部の様子を例えば断層像等その他のＣＴ画像として観察することができる。

この際、上記各構成要素は概略次のように作用する。すなわち、高電圧装置１０１ａの電圧印加によってＸ線管１０１より発せられ被検体Ｐを透過したＸ線は、Ｘ線検出器１０２によりアナログ電気信号に変換され、以下、データ収集部１０３によるデジタル変換処理、前処理部１０４による各種補正処理等を受け、投影データとしてメモリ部１０５に蓄えられる。そして、前記Ｘ線管１０１及びＸ線検出器１０２が、被検体Ｐ周囲を例えば１周（＝３６０°）した結果得られる投影データに基づいて、前記再構成部１０６による断層像等その他のＣＴ画像の再構成を行い、画像表示部１０７は当該ＣＴ画像を表示する。なお、再構成されたＣＴ画像は、記憶装置１０Ｍに記憶させることも可能である。

ここで上記した「ダイナミックＣＴ」撮影とは、上記Ｘ線管１０１及びＸ線検出器１０２が被検体Ｐの同一部位を反復撮影（反復スキャン。連続回転型ＣＴ装置では、連続回転による反復撮影がしばしば行われる。）して次々に投影データを取得するとともに、該投影データに基づいて次々に再構成処理を行って時系列的な一連の画像を得る撮影方式のことをいう（この場合、画像表示部１０７における画像表示は、例えば図示しないカウンタ等によって、その画像の元となった投影データ収集に係るスキャン開始点又は終点から一定時間後に行われるように制御される。）。

したがって、このように取得・表示される画像は、映画等と同様に時系列的な複数枚の静止画像からなる、いわゆる動画像となる。なお、このような撮影方式は、典型的には、被検体Ｐに対し造影剤を注入し、その経時変化を観察・解析して、例えば血管における狭窄や閉塞等その他病変部の病態を分析するために用いられる。また、造影剤投与の前後２回だけに限り同一部位のＣＴ撮影を行う方式も、広義のダイナミックＣＴ撮影と考えることができる。

さて、従来においては、上記のような「ダイナミックＣＴ」撮影時、例えばＫ回の撮影を実施する間に被検体Ｐに何らかの変化（例えば、造影剤の濃度変化や呼吸動等が一般的に考えられる。）があった場合、空間解像度を損なわず画像ノイズを抑制するためには、時間方向の平滑化を行うほかなかった。その結果、時間分解能が損なわれるという弊害は避け得なかった。

ところが、ダイナミックＣＴ撮影により取得される画像は、上述したように、動画像であって時間的変化を仔細に観察する目的で行うものであるから、その時間分解能が損なわれるというのは、本来、好ましい状況とは言えない。

本発明に係るコヒーレント・フィルタを利用すれば、時間分解能を損ねず、Ｋ枚の静止画像のすべて（複数枚の画像）につきそのノイズを抑制することが可能な、次のような処理（以後「ダイナミック・コヒーレント・フィルタ処理」と呼ぶ。）を実施することができる。

まず、上記のようにして得られた動画像たるＫ枚の静止画像につき定義される画素ｘについては、既に述べたように、画素値ｖ（ｘ）として、

を構成することができる。ここで右辺各項における添え字１，２，…，Ｋは、Ｋ枚の各静止画像の各々を通して振られた番号である（上述（３）式に関する説明参照）。

次に、この場合における重み関数ｗ１の具体的形式を、例えば次の式により与える。

ただし、ｙ∈Ｎ（ｘ）であって、かつ、この集合Ｎ（ｘ）は、画素ｘにつき任意に設定してよい（＝どのような基準によって設定してもよい。）。しかし実際上は、画素ｘと該画素ｘから遠く離れた位置にある画素ｙとが仮説「ｖ（ｘ）＝ｖ（ｙ）。ただし、両画素のノイズに起因する差異を除く」を満たす可能性は一般に低いといえるから、集合Ｎ（ｘ）をｘに近接している画素の集合という基準で限定することは、演算速度向上等の実用的な意義がある。

したがってここでは、その一例として、集合Ｎ（ｘ）を、当該画素ｘを中心としたその周囲の矩形状エリアに含まれる画素の集合、とする。より具体的に、集合Ｎ（ｘ）としては、例えば、いま注目している静止画像一枚を構成する全画素が１２８×１２８画素であるような場合に、前記画素ｘを中心とした３×３画素分のエリアとしたり、また、５１２×５１２画素であるような場合に、当該画素ｘを中心とした１３×１３画素分のエリア等としてもよい。

また、上記（６）式におけるσ_ｋは、ｋ枚目の静止画像の各画素が、そのどれにも共通な一定の程度で有するものと仮定して推定されたノイズの標準偏差であり、一方Ｃは、重みｗ１（ｐ（ｘ，ｙ））が、上記（４）式に代入された場合における作用の程度を決定調節可能なするパラメータである。

以下、これらσ_ｋ及びＣについての説明を順に行う。

まず、（６）式におけるσ_ｋについて説明する（以下では、分散σ_ｋ ^２として説明する。）。このσ_ｋ ^２は、上述したように、ｋ枚目の静止画像上の各画素のスカラー値が有するノイズ成分の分散である。そしてまた、上記（６）式における分散σ_ｋ ^２は、ｋ枚目の画像の各画素のスカラー値について一定値たる分散σ_ｋ ^２を持つノイズを含んでいるものと仮定して推定したものである。一般に、このような仮定は、次に記すようなことを背景として、十分な正当性を持つ。

まず、被検体Ｐの大きさ、Ｘ線管１０１及びＸ線検出器１０２、再構成部１０６等の構造が一定で、かつ、照射Ｘ線のエネルギを一定にした状態では、ＣＴ画像のノイズは、照射Ｘ線量、すなわちこれと比例関係にあるＸ線管１０１における管電流と照射時間との積（いわゆる管電流時間積（ｍＡ・ｓ））によって決定される。一方、ＣＴ画像のノイズは加法的であり、概ねガウス分布に従うことも知られている。すなわち、ある画素ｘの画素値ｖ（ｘ）を構成する任意のスカラー値ｖ_ｎ（ｘ）（ｎ＝１，２，…，Ｋ）について、その真値（ノイズの寄与分を除去した値）をｖ_ｎ ^０（ｘ）とすると、これらの差の値ｖ_ｎ（ｘ）−ｖ_ｎ ^０（ｘ）は、概ね平均０、分散σ_ｋ ^２のガウス分布に従う（なお、照射Ｘ線量ないし管電流時間積ｍ・Ａｓとノイズの分散σ_ｋ ^２とは、概ね反比例関係にある。）。

また、この分散σ_ｋ ^２は、画素ｘの位置そのもの（上で述べたように、例えば各座標値ｘ＝（ｘ，ｙ））にも依存するが、通常のＸ線ＣＴ装置１００においては、Ｘ線管１０１及びＸ線検出器１０２の間に、Ｘ線照射量を調節する物理的なＸ線フィルタ（例えば銅箔や金属塊等により構成された、いわゆる「ウェッジ」あるいは「Ｘ線フィルタ」と呼称されるもの。）（図示せず）を備えているため、これを無視することができる。なぜならばウェッジは、被検体Ｐが水とほぼ同じ密度を持つ物質から構成されていることを利用して、どのＸ線検出器１０２においても同程度のＸ線量が検出されるよう、照射されるＸ線量の一部又は全部を調節（吸収ないし遮蔽）する作用を有するものであり、従ってこのようなウェッジによれば、結果的に、ノイズの分散σ_ｋ ^２を画素ｘの位置に殆ど依らない一定値にする効果を生じるからである（ちなみに、このウェッジは、一般に、Ｘ線検出器１０２のダイナミックレンジを有効に利用することを本来の目的として設置されるものである。）。

以上のことから、ダイナミックＣＴ撮影により取得されたＫ枚の静止画像上においては、ｋ枚目の静止画像上におけるすべての画素について、一定値たる分散σ_ｋ ^２を推定することは妥当である。むろん、画素ごとに分散が異なる場合について本実施例を拡張することも容易に推考することが可能であろう（詳しくは後述するＸ−３項で説明する）。

さて次に、上記（５）式を具体的に演算するためには、その分散σ_ｋ ^２として、どのような数値をあてるか、が問題となる。このようなことが問題となるのは、通常、ノイズの分布の形は想定できても（上記ではガウス分布）、分散σ_ｋ ^２の具体値は不明であることが多いからである。

更に、一般的に、毎回の撮影毎に照射線量（Ｘ線管電流×照射時間（ｍＡｓ））を変更して撮影を行ってもよい。

さて、ｋ枚目の画像（ｋ＝１，２，…，Ｋ）に於いて各画素のスカラー値が持つノイズの分散をσ_ｋ ^２とし、ｋ枚目の画像の撮影に用いた照射線量をＲ_ｋとするとき、σ_ｋ ^２はＲ_ｋに比例する。従って少なくともひとつのｋ＝ｋ_０についてσｋ_０ ^２が指定できれば、他のｋに関しても

によってσ_ｋ ^２を正確に推定することができる。

本実施形態（このような事情が当てはまる）に於いては少なくともひとつのｋについて、以下のような方法でσ_ｋ ^２の具体的数値の推定を行うことができる。

Ｋ回の撮影のうち、被検体Ｐに殆ど変化がなかったと仮定することのできるＮ回（１＜Ｎ≦Ｋ）の画像を用いて、実測により、分散σ_ｋ ^２に対する期待値Ｅ［σ_ｋ ^２］を求める方法が有効である。以下説明を簡単にするために、これらＮ枚の画像における照射線量は同じであり、従ってｋ＝１，２，…Ｎに関してσ_ｋ ^２は一定（σ^２と書く）と仮定する。これらＮ枚の画像における、ある画素ｘ_ｆの画素値ｖ（ｘ_ｆ）を構成する各スカラー値ｖ_１（ｘ_ｆ），ｖ_２（ｘ_ｆ），…，ｖ_Ｋ（ｘ_ｆ）が含むノイズは、上述したように平均０、分散σ^２のガウス分布に従うと予想されるから、これらの平均値

を用いると、真の分散σ^２に対する期待値Ｅ［σ^２］を、

として求めることができる。そして、この分散の期待値Ｅ［σ^２］は、上述した通り、Ｋ枚すべての静止画像上の全画素ｘにつき妥当するものと考えることができ、真の分散σ^２の代用として用いるのに、一定程度以上確からしさが保証された値である。したがって、上記（６）式の実際の演算においては、このＥ［σ^２］を（６）式のσ^２に代入すればよい。

なお、このようなＥ［σ^２］は、より具体的には、Ｋ枚の静止画像中、例えば１枚目と２枚目の静止画像に基づく実測値により求めてもよい（上記（７）及び（８）式で言えば、Ｎ＝２とすることに該当する。）。また、上記（７）及び（８）式の実際の演算に供される画素ｘ_ｆについては、例えば、空気や骨が撮像されている部分を除いた適当な画素ｘ_ｆのみを選定する（複数選定した場合は得られるＥ［σ^２］すべての平均をとる）等といった工夫を施してもよい。さらに、その他一般的には、被検体Ｐの動きによる影響を抑える工夫等を施すと尚よい（なお、（６）式におけるノイズの分散σ^２の評価については、最後に述べる［本実施形態の補足事項」においても再び触れる。）。

これらＮ枚の画像の撮影において照射線量が一定でない場合においても、σ_ｋ ^２がＲ_ｋに比例することを利用して正しくσ_ｋ ^２を推定することは容易に推考できるであろう。

さて次に、上記（６）式におけるパラメータＣについての説明を行う。まず、（６）式においては、上記一般的形態で述べた危険率ｐ（ｘ，ｙ）の考え方が、以下のようにして含まれている。すなわち、（６）式の右辺分子における根号内の表式は、いわゆるχ二乗分布に従うとされる当該χ^２値に一致するものであり、これを（２σ）^２で除し、括弧の全体をｅの肩に置いた値は、危険率ｐ１（ｘ，ｙ）そのものである。つまり、

たものである。

結局、（６）式においては、上記したような一般的形態で述べた危険率ｐ（ｘ，ｙ）が陽には表示されてはいないが、重みｗ１（ｐ（ｘ，ｙ））の実態は、上述したように、まさしく危険率（＝ｐ１（ｘ，ｙ））の関数であると見ることができ（（１０）式）、すなわち「適合度の関数」である（ただし、危険率と適合度とは、上述したように、一方が増えれば他方も増加する関係にある。）。

そして、上記（１０）式からわかるように、パラメータＣは、重みｗ１（ｐ（ｘ，ｙ））が、危険率ｐ１（ｘ，ｙ）にどの程度敏感に反応するかを決める効果がある。つまり、Ｃを大きくすると、ｐ１（ｘ，ｙ）がわずかに小さくなるだけで、ｗ１（ｐ（ｘ，ｙ））は０に近づく。また、Ｃを小さくするとそのような過敏な反応を抑制することができる。なお、Ｃとして、具体的には１乃至１０程度とすればよく、好適にはＣ＝３とするとよい。

この実施形態においては、両画素ｘ及びｙに関する類似判定、言い換えると、両画素ｘ及びｙに関する上述した帰無仮説Ｈの棄却の判定は、上述したことから明らかなように、上記危険率ｐ１（ｘ，ｙ）に基づいて、いわゆるχ二乗検定法（統計的検定法）によって決定されている。

また、上記（６）式の表式からわかるように、本発明においては、危険率ｐ（ｘ，ｙ）をｘ，ｙの組み合わせそれぞれについて計算した後、重みｗ（ｐ（ｘ，ｙ））を求めるといった手順を踏む必要は必ずしもなく、危険率ｐ（ｘ，ｙ）を具体的に求めずに、合成関数としての（ｗｏｐ）を、直接計算する構成としてもよい。

以上述べたように、分散σ^２の推定をし（例えば、（８）式のＥ［σ^２］）、かつ、パラメータＣを適当に決める（例えば、Ｃ＝３）ことにより、（６）式を用いて、ある画素ｘにつき定義される集合Ｎ（ｘ）（上述したように、例えば画素ｘを中心とした３×３画素分のエリア等）に含まれるすべての画素ｙについて、具体的な重みｗ１（ｐ（ｘ，ｙ））を求めることができる。後は、上記（４）式におけるｗ（ｐ（ｘ，ｙ））に代えて、このｗ１（ｐ（ｘ，ｙ））を用いることにより、コヒーレント・フィルタの具体的な数値演算を実施することが可能となる。そしてその結果、時間分解能は勿論のこと、空間分解能をも損なわずに、ノイズを強く抑制した画素値ｖ′（ｘ）＝（ｖ′_１（ｘ），ｖ′_２（ｘ），…，ｖ′_Ｋ（ｘ））（＝（５）式）、すなわちそのようなＫ枚の静止画像ないし動画像を、得ることができる。

このような画像処理を、概念的に把握しやすいよう図示したものが、図４である。すなわちまず、図４（ａ）においては、１，２，…，Ｋ枚ある静止画像において、ある画素ｘにつき、該画素ｘを中心とした３×３画素分の矩形状エリアＮ^３×３（ｘ）が想定されている。この矩形状エリアＮ^３×３（ｘ）の左角隅における画素を、ｙ_１とすれば、この画素ｙ_１は、図４に併せて示すように、画素値ｖ（ｙ_１）を有している。

そして、この画素値ｖ（ｙ_１）を構成するスカラー値ｖ_１（ｙ_１），ｖ_２（ｙ_１），…，ｖ_Ｋ（ｙ_１）と画素値ｖ（ｘ）におけるスカラー値ｖ_１（ｘ），ｖ_２（ｘ），…，ｖ_Ｋ（ｘ）とのそれぞれにより、上記（６）式によって重みｗ１（ｐ（ｘ，ｙ_１））が計算される（図４（ｂ））。また、矩形状エリアＮ^３×３（ｘ）の残る画素ｙ_２．…，ｙ_８についても同様で、結局図４（ｂ）に示すように、ｗ１（ｐ（ｘ，ｙ_１）），…，ｗ１（ｐ（ｘ，ｙ_８））及び、ｗ１（ｐ（ｘ，ｘ））が得られる。（この場合、（９）式より危険率ｐ（ｘ，ｘ）は、“１”であり、したがって重みｗ１（ｐ（ｘ，ｘ））も、（１０）式より“１”である（＝最大の重み付けがされている））。

次に、このようにして得られた、重みｗ１（ｐ（ｘ，ｙ_１）），…，ｗ１（ｐ（ｘ，ｙ_８）），ｗ１（ｐ（ｘ，ｘ））を、対応する画素の、ｋ枚目の画像におけるスカラー値ｖ_ｋ（ｙ_１），ｖ_ｋ（ｙ_２），…，ｖ_ｋ（ｙ_８），ｖ_ｋ（ｘ）にそれぞれ乗算して総和を取り（上記（４）式における分子に該当する。）、これを矩形状エリアＮ^３×３（ｘ）に関する重みｗ１の総和（同じく（４）式の分母に該当する。）により除せば、当該ｋ枚目の画像における画素ｘについての、ノイズが抑制されたスカラー値ｖ′_ｋ（ｘ）を求めることができる（図４（ｃ））。また、ｋ＝１，２，…，Ｋのすべての画像につき、同じ重みｗ１（ｐ（ｘ，ｙ_１）），…，ｗ１（ｐ（ｘ，ｙ_８）），ｗ１（ｐ（ｘ，ｘ））を用いて、ノイズが抑制されたスカラー値ｖ′_ｋ（ｘ）を求めることによって、画素ｘにおけるノイズが抑制された画素値ｖ′_ｋ（ｘ）＝（ｖ′_１（ｘ），ｖ′_２（ｘ），…，ｖ′_Ｋ（ｘ））が得られる。すべての画素ｘにつき、上記演算を繰り返せば、ノイズを抑制したＫ枚の画像が得られる。

このようにして算出された画素値ｖ′（ｘ）で構成される画像は、例えば、図５に示すようになる。図５（ａ）は、全部でＫ＝２３枚の静止画像のうち、コヒーレント・フィルタを適用する前のオリジナルの画像（ｋ＝１１）を示し、（ｂ）及び（ｃ）は、ｋ＝１１及び２３について、コヒーレント・フィルタをかけた場合の画像を、それぞれ示している。図５（ａ）のオリジナル画像で見られるランダムなノイズは、図５（ｂ）及び（ｃ）においては、十分に抑制されていることがわかる。

なお、以上までに述べた各処理は、例えば図６に示すようなフローチャートに則ってこれを行えばよく、また、当該各処理に係る演算・画像表示等を実際のＸ線ＣＴ装置１００上で実現するためには、例えば、図３に示すように、分散値推定部１１１、重み演算部１１２及び画素値演算部１１３により構成される画像処理部１１０を設けて、これを実施すればよい。

このうち重み演算部１１２は、上述した手順通り、画素値ｖ（ｘ）及びｖ（ｙ）から直接重みｗ１（ｐ（ｘ，ｙ））を求める構成となっている。したがって当該演算部１１２は、危険率ｐ１（ｘ，ｙ）の値を具体的に求めることなく（すなわち、「危険率演算部（本発明にいう「適合度定量化部」を内蔵し）、重みを直接に求める装置である。なお、上記したような構成ではなく、具体的に危険率ｐ１（ｘ，ｙ）の値を求める「危険率演算部（適合度定量化部）」と、その出力に基づいて重みｗ１（ｐ（ｘ，ｙ））を求める「重み演算部」という、二段の手順を踏む構成としてもよい。いずれにせよ、重み演算部１１２は、分散値推定部１１１により推定された分散σ^２と、ｖ（ｘ）及びｖ（ｙ）を用いて重みｗ１（ｐ（ｘ，ｙ））を算出する。

また、画素値演算部１１３は、画素値ｖ（ｘ）及びｖ（ｙ）、並びに重み演算部１１２により数値演算された重みｗ１（ｐ（ｘ，ｙ））を使って、画素値ｖ′（ｘ）を演算する。すなわち当該演算部１１３は、元となる画像のノイズを抑制する処理、すなわちコヒーレント・フィルタの適用を実際に行う（以下、これを「コヒーレント・フィルタをかける」と表現する。）。

上記のようなダイナミック・コヒーレント・フィルタ処理においてＫ枚の静止画像から構成される動画像に、コヒーレント・フィルタをかける場合には、上記画像処理部１１０における処理は、一旦すべての静止画像を再構成した後、これらを上記記憶装置１０Ｍに蓄え、後処理として後にこれらに対してコヒーレント・フィルタをかけるようにしてもよいが、本発明はこのような形態に限定されるものではなく、上述した連続スキャン、連続投影データ収集、連続再構成及び連続表示という流れの中で、コヒーレント・フィルタをかける処理をリアルタイムに実施する（以下、これを「リアルタイム・コヒーレント・フィルタ処理」と呼ぶ。）のでもよい。

リアルタイム・コヒーレント・フィルタ処理の好ましい実施例においては、新しい画像が撮影され再構成されるたびに、以下のような処理を行う。最初に得られた画像（画像番号１）から最新の画像（画像番号Ｍ）までのうち、画像番号Ｍ，Ｍ−１，…，Ｍ−Ｋ＋１を持つＫ枚の静止画像上、共通の同一点（同一座標）の画素ｘの持つ画素値（スカラー値）を並べてＫ次元ベクトル値ｖ（ｘ）＝（ｖ_Ｍ（ｘ），ｖ_Ｍ−１（ｘ），…，ｖ_{Ｍ−Ｋ＋１}（ｘ））を構成する。こうして、上記の「ダイナミック・コヒーレント・フィルタ処理」と全く同様にコヒーレント・フィルタをかけることができる。ただし、画素値演算部１１３は実際には画素値ｖ′（ｘ）の全ての要素を計算するのではなく、最新の画像（画像番号Ｍ）に対応するスカラー値ｖ_Ｍ′（ｘ）だけを計算する。この結果、計算速度が向上するので、リアルタイムでノイズが抑制された最新の画像を表示できる。

この「リアルタイム・コヒーレント・フィルタ処理」の別の好ましい実施例として、最初のＫ枚の画像が得られた時点で、上記と全く同様にコヒーレント・フィルタをかけてｖ_１′（ｘ），…，ｖ_Ｋ′（ｘ）を求めておき、以後は、Ｋ次元ベクトル値を画像番号Ｍ，Ｍ−１，…，Ｍ−Ｋ＋１を持つＫ枚の静止画像を用いてｖ（ｘ）＝（ｖ_Ｍ（ｘ），ｖ_Ｍ−１′（ｘ），…，ｖ_{Ｍ−Ｋ＋１}′（ｘ））によって構成し、これに対して上記のリアルタイム・コヒーレント・フィルタ処理を適用するように構成してもよい。なお、これらのリアルタイム・コヒーレント・フィルタ処理の際に画素値ベクトルｖ（ｘ）の次元Ｋを、マニュアル設定、あるいは自動設定によって、随時変更できるように構成しておくと便利である。

（ ＩＶ Ｘ線ＣＴ撮影に本発明を適用する場合２；（ノイズを低減した１枚の画像を得る））
さて次に、上記ダイナミックＣＴとは異なり、ノイズを低減した１枚の画像を得る具体例について説明する。従来、一般に、ノイズの少ない１枚のＣＴ画像を得たいという場合、それを達成するためには、通常のＸ線管電流を使ってＫ回繰り返して撮影したＫ枚のＣＴ画像に対して時間方向の単純平均（例えば、本実施形態における画素値ｖ（ｘ）の概念を利用すれば、当該単純平均とは｛ｖ_１（ｘ）＋…＋ｖ_Ｋ（ｘ）｝／Ｋである。）をした結果１枚の画像を得る方式と、通常のＫ倍のＸ線管電流を使った一度の撮影により１枚の画像を得る方式との二通りが考えられる。いずれの方式にしても、ノイズの分散は通常のＸ線管電流を使って１回撮影して得られる１枚の画像の１／Ｋ倍になることに変わりはない（前者の方式ではノイズの平均化により、後者の方式では上述したように照射Ｘ線量とノイズの分散が反比例関係にあることにより）。したがって従来においては、操作の効率性から、Ｘ線管の性能が許す範囲では後者の一度に撮影する方式により、また、Ｘ線管の性能が大量のＸ線管電流に耐えられない場合は、前者のＫ回繰り返して撮影する方式により、１枚のＣＴ画像を得ることが一般的に行われていた。

しかしながら、これらの方式では、ノイズを抑制しようとすればする程、照射Ｘ線量が増大する。つまり被検体Ｐの被曝量が増大することになる。

ここで、本発明に係るコヒーレント・フィルタを利用すれば、このような１枚のＣＴ画像に含まれるノイズを、被検体Ｐに対する被曝量を増大させずに抑制することが可能である。

まず、本実施形態では、被検体Ｐについて何らの変化もないと仮定できる時間中に、Ｋ枚のＣＴ画像を、通常より低いＸ線照射量で撮影する。したがって、この場合、上記ダイナミック・コヒーレント・フィルタ処理と同様に、Ｋ枚のＣＴ画像につき定義される画素ｘについて、

というように、画素値ｖ（ｘ）を構成することができる。なお、この場合においては、Ｘ線照射量が少ないので、一枚一枚のＣＴ画像に含まれるノイズの分散σ^２は比較的大きいことになる。

次に、この場合における重み関数として、上記（６）式と同様な具体的形式を与える。

ただし、この場合においては、最終的に得たい画像は、１枚のＣＴ画像であるから、次のような方式をとることが考えられる。

例えば、上記ダイナミックＣＴ撮影と同様にして、Ｋ枚の画像すべてにコヒーレント・フィルタをかけ、その結果として得られるＫ枚の画像を重み付平均あるいは単純平均して、目的とする１枚のＣＴ画像を構成する方式や、また、上記Ｋ枚の画像の中から任意に選択された１枚に対してだけコヒーレント・フィルタをかける方式（＝図４（ｃ）の演算は、選択された当該１枚の画像についてのみ行われる。）、等とすればよい。

ちなみに、前者の方式によれば、Ｋ枚の画像すべてにコヒーレント・フィルタをかけた結果であるＫ枚の画像から、これらを平均した１枚の画像を得るための計算は、具体的には例えば、

である。ここに、（１１）式の左辺にあるｖ^＊′（ｘ）は、単純平均の結果得られる１枚の画像の、画素ｘにおける画素値（スカラー）である。このとき、ｖ^＊′（ｘ）が含むノイズの分散（σ^２）^＊は、概ね、

となり、またここで、ｗ（ｐ（ｘ，ｙ））≧０であるから、（１２）式の右辺において、総和されている項は、常に１よりも小さい。したがって、本実施形態によるノイズ抑制効果は、従来の方式により達成されたノイズ抑制効果（σ^２／Ｋ）よりも優れている。逆に言えば、もし本実施形態を用いて従来と同等のノイズ抑制効果を得ようとするならば、従来の方式よりも、必要なＸ線照射量が少なくて済み、従って、被検体Ｐの被曝量を低減させることができる。

本実施形態を用いて得られる効果は、たとえば図７に示すようなものである。図７（ａ）は、従来方式、すなわちＫ回繰り返して撮影した画像を単純に時間方向に平均化して得られる画像を示し、（ｂ）は、本実施形態に従って、Ｋ枚の画像すべてにコヒーレント・フィルタをかけ、その結果として得られるＫ枚の画像を単純平均して構成した画像を示している。図７（ａ）に見られるランダムなノイズは、図７（ｂ）においては、十分に抑制されていることがわかる。（なお、図７は、０．５秒スキャンの撮影を３回行った（Ｋ＝３）場合で、かつ、上記集合Ｎ（ｘ）として１３×１３画素分の矩形状エリアを用い、Ｃ＝３とした条件下の画像である。）
また、図７（ｃ）は、図７（ｂ）におけるノイズ抑制効果があまりに優れているため、従来の画像に比べて見慣れない画像となる点に配慮して構成されたものであり、具体的には、図７（ａ）と図７（ｂ）との平均をとった画像である。すなわち、図７（ｃ）の画像は、その画素値をｖ^＊＊′（ｘ）、図７（ａ）における画素値をｖ（ｘ）（＝｛ｖ_１（ｘ）＋…＋ｖ_Ｋ（ｘ）｝／Ｋ）としたとき、

において、ｕ＝１／２としたものである。なお、本発明では、本実施形態そのもの（図７（ｂ）が得られる）に限らず、例えば上記（１３）式におけるｕの値（ただし、０≦ｕ≦１）を、診断者の好みに応じて調整できるようにした構成も、好ましい。

（Ｖ Ｘ線ＣＴ撮影に本発明を適用する場合３；マルチエネルギ撮影（複数種類の画像のノイズを低減する） ）
さて次に、本発明を、いわゆる「マルチエネルギ撮影」に適用した場合について説明する。ここに、マルチエネルギ撮影とは、Ｘ線管１０１に印加する電圧を変更したり、Ｘ線管１０１と被検体Ｐとの間に、薄い金属板等で構成されたＸ線フィルタを挿入すること等によって、照射するＸ線のエネルギ分布を種々に変え、このようにして得られた複数のエネルギ分布により、同一の被検体Ｐを撮影する方式をいう。この場合には、被検体Ｐは変化していないにもかかわらず、Ｘ線のエネルギ分布に依存して画像が変化し、その結果複数種類の画像が得られる。

ここで、エネルギ分布を変更しながらＫ回の撮影を行ったとすると、当該Ｋ回の撮影で得られたＫ種類の画像の画素ｘにおけるスカラー値を並べれば、画素値ｖ（ｘ）＝（ｖ_１（ｘ），ｖ_２（ｘ），…，ｖ_Ｋ（ｘ））を構成することができる。したがって、この画素値ｖ（ｘ）に対して、（ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（複数枚の画像からノイズを除去する））におけるダイナミック・コヒーレント・フィルタ処理と同様の手段を適用して、コヒーレント・フィルタをかけることができ、その結果、ランダムなノイズを抑制したＫ種類の画像を得ることができる。

（ ＶＩ Ｘ線ＣＴ装置以外の医用画像診断装置に本発明を適用する場合 ）
以下では、本発明に係るコヒーレント・フィルタを、Ｘ線ＣＴ装置以外の医用画像診断装置、例えば、核磁気共鳴イメージング装置（いわゆる「ＭＲＩ装置」）、核医学診断装置（「シンチレーションカメラ、ＳＰＥＣＴ装置、あるいはＰＥＴ装置」）、Ｘ線透視装置等に適用する実施例について説明する。

（ ＶＩ−Ｉ ＭＲＩ装置に本発明を適用する場合）
第一に、ＭＲＩ装置に関する本発明の適用例について説明する。本発明は、このＭＲＩ装置に関し、従来「アベレージング（ａｖｅｒａｇｉｎｇ）」、「高速撮影」、「三次元的撮影」と呼ばれている撮影方式に対して、また「複数のパルス・シーケンスを用いる撮影方式」等に対して、適用することが可能である。

まず、アベレージングとは、その名の通り平均化（ａｖｅｒａｇｉｎｇ）によって画像ノイズを抑制することを目的として、被検体Ｐに変化が生じない間に撮影を繰り返し行い、得られた複数枚（Ｋ枚）の画像を単純平均して１枚の画像を得る撮影方式である。当該Ｋ枚の画像の画素ｘにおけるスカラー値を並べれば、画素値ｖ（ｘ）＝（ｖ_１（ｘ），ｖ_２（ｘ），…，ｖ_Ｋ（ｘ））を構成することができる。また、重み関数ｗについても、例えば上記（５）式を用いてもよい。従って、この画素値ｖ（ｘ）に対して、（ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（複数枚の画像からノイズを除去する））におけるダイナミック・コヒーレント・フィルタ処理と全く同様の手段を適用して、コヒーレント・フィルタをかけることができ、その結果、アベレージングによるものより優れたノイズを抑制した複数枚の画像が得られる。あるいは、（ ＩＶ Ｘ線ＣＴ撮影に本発明を適用する場合２；（ノイズを低減した１枚の画像を得る））と全く同様の処理を適用して、ノイズを抑制した１枚の画像を得ることもできる。

また、高速撮影とは、被検体Ｐに変化が生じている場合に、それを反復して撮影することによって経時変化を表わす動画を得る方式であり、ダイナミックＭＲＩとも呼ばれる。この場合にも、複数枚（Ｋ枚）の静止画像が得られることに変わりはないから、上記アベレージングの場合と同様に（ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（複数枚の画像からノイズを除去する）におけるダイナミック・コヒーレント・フィルタ処理と全く同様に）して、コヒーレント・フィルタを適用し、ノイズを抑制できる。さらに、高速撮影においては、（ＩＩＩ．Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（複数枚の画像からノイズを除去する））におけるリアルタイム・コヒーレント・フィルタ処理と全く同様の処理を適用してもよい。

三次元的撮影とは、一度に三次元的画像を得る撮影法である。この撮影法においても、アベレージング（すなわち、ノイズを抑制することを目的として、反復して撮影を行い、単純平均をとる処理）を行ったり、高速撮影（被検体Ｐに変化が生じている場合に、反復して撮影すること）を行う場合があり、これらの場合には、画素ｘが、ｘ＝（ｘ，ｙ，ｚ）なる三次元ベクトルとして表現されるという点を除けば、上述したアベレージングや、高速撮影におけるコヒーレント・フィルタの適用方法と全く同様である。さらに、三次元的撮影が一度だけ行われる場合にコヒーレント・フィルタを適用するためには、後述する（ＶＩＩＩ １枚の画像から画素値ｖ（ｘ）を構成する）において説明する方式を適用することができる。

一方、ＭＲＩ装置においては、パルス・シーケンスを変更することによって、異なる物理パラメータを画像にすることができる。パルス・シーケンスとは、ＲＦパルス及び傾斜磁場の印加順序、印加する強度と時間、休止時間の組み合わせのことで、これを変化させることによって非常に多様な画像が得られる。典型的に使用されるパルス・シーケンスのごく少数の例を挙げれば、画像のコントラストに最も寄与する物理パラメータの名前を冠して、Ｔ２＊−ｗｅｉｇｈｔｅｄ ｉｍａｇｅ、Ｔ１−ｗｅｉｇｈｔｅｄ ｉｍａｇｅ、ｐｒｏｔｏｎ−ｄｅｎｓｉｔｙ ｉｍａｇｅ等と呼ばれるものが知られている。いずれにせよ、パルス・シーケンスに依存して、被検体は変化していないにもかかわらず、異なる画像が得られ、その結果、実質的に複数種類の画像が得られることになる。従って（ Ｖ Ｘ線ＣＴ撮影に本発明を適用する場合３；複数種類の画像のノイズを低減する）と全く同様にしてコヒーレント・フィルタをかけることができる。

より具体的には、パルス・シーケンスを変えながらＫ回の撮影を行ったとすると、当該Ｋ回の撮影で得られた画像の画素ｘにおけるスカラー値を並べれば、ベクトルｖ（ｘ）＝（ｖ_１（ｘ），ｖ_２（ｘ），…，ｖ_Ｋ（ｘ））を構成することができる。したがって、この場合においても、上記アベレージングの場合と同様に、ダイナミック・コヒーレント・フィルタ処理と全く同様の手段を適用して、コヒーレント・フィルタを適用し、ノイズを抑制できる。

（ ＶＩ−２ シンチレーションカメラ、ＳＰＥＣＴ装置又はＰＥＴ装置に本発明を適用する場合）
第二に、シンチレーションカメラ、ＳＰＥＣＴ装置又はＰＥＴ装置に関する本発明の適用例について説明する。これらの装置は、被検体内に放射性同位元素（以下「ＲＩ」という。）を含む薬剤（放射性薬剤）を投与し、ＲＩの自然崩壊によって生じるガンマ線を検出することにより、被検体内における放射性薬剤の分布を画像化するものであり、特に生体を被検体とした場合には、その生化学代謝・血液循環の程度を表す、いわゆるファンクショナル・マップを得ることを目的として用いられる。シンチレーションカメラは、放射性薬剤の分布を、あたかも被検体を透視したかのような２次元画像として撮影する。また、ＳＰＥＣＴ装置およびＰＥＴ装置は、コンピュータ断層撮影技術を用いて、放射性薬剤の３次元分布を表す３次元画像を撮影する。

シンチレーションカメラ、ＳＰＥＣＴ装置又はＰＥＴ装置では、通常、検出されるガンマ線量が極めて小さいから、ＲＩの量子力学的性質によって生じる光子数のばらつきに起因する、いわゆる「フォトンノイズ」が画像に現れる。

従来、画像の解像度を低下させることなくフォトンノイズを抑制するためには、検出されるガンマ線量を増やす以外の方法がなく、具体的には投与するＲＩの量を増やす方法と、撮影時間を長くする方法があった。

投与しうるＲＩの量が増加できる場合（例えば、医学研究用の実験動物に対してなら投与量増大は許される。）には、それによってフォトンノイズを低減できるが、このような方法は作業者の放射線被曝線量の増加を伴う。また、人体を診断する（＝人体が被検体である）場合には、そもそも上記したようなＲＩ量の増加が自由には許されていない（法的に規制されている）。さらに、自然物に天然に含まれるＲＩの出す放射線を用いた撮影を行う場合には、その量を増やすことは不可能である。

また、撮影時間を長くする方法を用いようとする場合、有効なノイズ低減効果を得るためには非常に長時間（たとえば数時間）に渡って被検体が動かないことが必要である。この条件は、生体（特に人体）を被検体とする場合には実現困難であり、また撮影時間が長いことによって、撮影装置の処理能力（あるいはコスト・パフォーマンス）が低下するという問題があった。従って現実的には、許容可能な時間（たとえば１時間以内）を限度とする、適当な撮影時間を設定するしかなく、画像に強いノイズが現れることは避けられない。そこで、このノイズを抑制するには画像の平滑化によって行うしかなく、従って、平滑化に伴う弊害として画像の解像度が損なわれるのはやむをえないと考えられていた。

このような場合において、本発明によれば、次のような処理を行うことができる。すなわち、前記の適当な撮影時間と同等の時間の間に、より短時間の撮影を複数回行い、ノイズの比較的多い複数の画像を得る。そして、（ ＩＶ Ｘ線ＣＴ撮影に本発明を適用する場合２；（ノイズを低減した１枚の画像を得る））と全く同様の処理を適用して、空間解像度を損なうことなく、ノイズがより抑制された画像を得ることができ、しかも、ＲＩの量を増やすことなく実現できる。

さらに、シンチレーションカメラ、ＳＰＥＣＴ装置又はＰＥＴ装置においても、Ｘ線ＣＴなどと同様に、同一の被検体を繰り返し撮影して、被検体内の放射性薬剤の分布の経時変化を追跡する撮影方式、すなわち「ダイナミックシンチグラム」、「ダイナミックＳＰＥＣＴ」又は「ダイナミックＰＥＴ」と称される撮影方式がある。この場合にも、（ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（複数枚の画像からノイズを除去する））における、ダイナミック・コヒーレント・フィルタ処理もしくはリアルタイム・コヒーレント・フィルタ処理と全く同様の処理を適用することができ、画像のノイズを効果的に抑制することができる。

さらに、ＰＥＴ装置やＳＰＥＣＴ装置では、ＲＩがその種類に固有のエネルギーのガンマ線を放出することを利用して、複数の種類のＲＩを同時に被検体に投与することで、一度に複数の種類のＲＩについて、それぞれの分布画像を得る、いわゆる「マルチウインドウ撮影法」が使用されることがある。これは、検出したガンマ線の持つエネルギを計測し、その値が人為的に設定した数値範囲（エネルギウインドウ）に該当するかどうかによって、該複数の種類のＲＩのうちのどれによって放射されたガンマ線であるかを判別する技術に基づいている。

より具体的には、適当な数のエネルギウインドウ（例えば、ウインドウＷ１及びＷ２）を設定することにより、そのうちの一のウインドウ（例えば、ウインドウＷ１）に該当するガンマ線に基づく画像、他のウインドウ（例えば、ウインドウＷ２）に該当するガンマ線に基づく画像、というように別個に画像を形成することにより、一度の撮影によって、上記した複数のＲＩ分布画像を得る。

すなわち、イメージング装置を用いて同じ対象物を撮影して異なる処理条件で複数種の画像を得る。より一般的に言えば、「Ｋ種類のＲＩを用いて、一度の撮影により、相異なるＫ種類の画像を得る撮影方式である」ということができる。従って（ Ｖ Ｘ線ＣＴ撮影に本発明を適用する場合３；複数種類の画像のノイズを低減する）と全く同様にしてコヒーレント・フィルタをかけることができる。

すなわち、この場合においても、このようなＫ種類の画像上の画素ｘにつき、その画素値ｖ（ｘ）＝（ｖ_１（ｘ），ｖ_２（ｘ），…，ｖ_Ｋ（ｘ））を構成することができる。また、各画素値ｖ（ｘ）が含むノイズの分散は、画素ｘの位置、及び、どの種類の画像であるかによって異なるが、画像の種類ごとに、それぞれ検出されたガンマ線の総量にほぼ比例するものと考えてよく、この性質から、各画素ｘについてのノイズの分散を推定することができる。以上のようにして、「マルチウインドウ撮影」に対しても本発明を適用してコヒーレント・フィルタをかけることができ、その結果、ランダムなノイズを抑制したＫ種類の画像を得ることができる。

なお、核医学診断装置においては、解像度向上に効果のある再構成法（ファンビーム用ＭＵ−ＥＮ再構成法及びファンビーム、パラレルホールコリメータ用超解像再構成法等）があるが、これらと本発明に係るコヒーレント・フィルタを組み合わせることには何の原理的制約もなく、その結果、高解像度かつランダムなノイズが十分抑制された画像を作成することも可能である。

（ ＶＩ−３ Ｘ線透視装置に本発明を適用する場合 ）
第三に、Ｘ線透視装置に関する本発明の適用例について説明する。まず、ここにいう「Ｘ線透視装置」とは、低線量のＸ線を被検体に連続して曝射しながら、一定時間（例えば１／２４秒）ごとに繰り返し撮影を行う装置であって、１回の撮影を「フレーム」と呼び、フレームごとに得られるＸ線画像を次々と連続的に表示することによって、いわゆる透視像（動画像）として観察することができる装置を想定している。したがって、一般的には、いわゆる「Ｃアーム」ないし「Ωアーム」と呼ばれる構成を備えた装置が該当することを始め、上記したＸ線ＣＴ装置１００においても、いま述べたような機能を搭載するものであれば、以下に述べる形態は当然に適用可能である。

さて、このようなＸ線透視装置においては、上記したように低線量であるがゆえに、前記透視像はフォトン・ノイズに起因するランダムなノイズを多く含む。このような問題を解決するためにＸ線照射量を十分に増大させると、連続したＸ線照射を行うのであるから、被検体に対する被曝量が甚だしく増大するのみならず、検査者に対する被曝も著しく増加してしまう。

本発明に係るコヒーレント・フィルタを適用すれば、これらの問題点を解決して、ノイズを十分に抑制した透視像を提供することが可能となる。

具体的には、以下の説明のために、１フレームごとに得られた画像に順次画像番号１，２，…，Ｍを付けることにすると、最も新しいＫ枚、すなわち画像番号Ｍ，Ｍ−１，…，Ｍ−Ｋ＋１の画像をメモリに蓄え、これらＫ枚の画像に対して、（ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（複数枚の画像からノイズを除去する））における、ダイナミック・コヒーレントフィルタ処理もしくはリアルタイム・コヒーレント・フィルタ処理と全く同様の処理を適用することができる。

なお、医用画像診断装置を離れるが、「Ｘ線透視装置」として、いわゆる非破壊検査に用いられる当該装置にも、上記の方法が適用できる。

また、本発明は、超音波診断装置により得られる断層像についても、同様に適用することが可能である。この場合も、被写体を連続的に撮影しながら次々と、画像を生成することに変わりはなく、従ってＸ線透視装置におけるのと同様の構成によってコヒーレント・フィルタをかけ、ノイズを抑制することができる。

（ ＶＩＩ より一般的に本発明を適用する場合；カラー画像等 ）
本発明は、上記した医用画像診断装置以外の、画像を扱う「装置一般」、あるいは「画像一般」に対しても適用することが可能である。以下では、そのような一般的な装置、あるいは画像に対する本発明の適用例について述べる。なお、ここでは、特に重要な分野である一般の「カラー画像」を対象とした場合に関する説明を、主に行うこととする。

カラー画像においては、既に述べたように、当該画像を構成する各画素ｘについて、光の三原色の明るさを要素とした３次元ベクトルを構成することができる。すなわち、

（ここでＶ_Ｒ，Ｖ_Ｇ，Ｖ_Ｂは、光の三原色、Ｒ（赤），Ｇ（緑），Ｂ（青）の光を意味する）である。これは、図２に示すとおり、「複数種類の画像から画素値ｖ（ｘ）を構成」する場合に該当することがわかる。

また、例えば当該カラー画像を、ＣＣＤ撮像素子、ＣＭＯＳ撮像素子等により構成されたデジタルカメラで撮影する場合には、各画素の画素値が含むノイズの分散は、当該画素ｘに入射した光量によって決まる。具体的には、スカラー値ｖ_Ｒ（ｘ），ｖ_Ｇ（ｘ），ｖ_Ｂ（ｘ）と、これらに含まれるノイズの分散σ^２ _Ｒ（ｘ），σ^２ _Ｇ（ｘ），σ^２ _Ｂ（ｘ）とは、それぞれほぼ比例する関係にある。この性質を利用して、ノイズの分散の推定値を求めることができる。したがって、上記（６）式と同様な重み関数ｗの具体的形式を与え、これにより、上記（４）式のコヒーレント・フィルタを構成することにより、ノイズを抑制した画素値ｖ′（ｘ）＝（ｖ′_Ｒ（ｘ），ｖ′_Ｇ（ｘ），ｖ′_Ｂ（ｘ））を得ることができる。

「カラー」ということを最広義にとらえ、上記にいう赤、緑及び青以外にも、例えば赤外線及び紫外線等を含め、上記（３．１）式を、より高次のベクトルとして、あるいはＲＧＢにこだわらないベクトルとして、構成することもできる。たとえば、カラー画像をベクトル値の画素値を用いて表現する方法が、必ずしも三原色だけに限らないことは公知であり、そのような三原色を用いないベクトル表示を用いた場合に対して、本発明を適用してもよい。また、紫外光Ｖ_Ｕを含めた４原色，紫色Ｖ_Ｐを分離した５原色（ちなみにこれは蝶類の視覚に相当する）等の場合も、たとえば５原色の場合ならば画素ｘがＶ（ｘ）＝（Ｖ_Ｒ（ｘ），Ｖ_Ｇ（ｘ），Ｖ_Ｂ（ｘ），Ｖ_Ｕ（ｘ），Ｖ_Ｐ（ｘ）という５次元ベクトル値を持つものとして取り扱い、上記の構成と同様にしてコヒーレント・フィルタをかけることができる。

さらに、例えば被写体の光反射スペクトルを計測しようとする場合等において、被写体を照明する照明光の波長を段階的あるいは連続的に切り替え、照明の波長ごとに繰り返しモノクロカメラによって撮影を行い、あるいはカラーフィルタを段階的に切り替えて繰り返し撮影を行うことによって多波長の撮影を行う撮影装置において、これらの多数枚の（各々単波長もしくは単色の）画像をもとに画素値ｖ（ｘ）を構成することができる。すなわち、（ Ｖ Ｘ線ＣＴ撮影に本発明を適用する場合３；（複数種類の画像のノイズを低減する））と全く同様の方法を適用して、コヒーレント・フィルタをかけることができ、その結果、ランダムなノイズを低減した複数種類の画像が得られる。

さらにまた、これらの撮像装置の多くは、一つの画像を得るために、通常１回の撮影を行う。しかしながら、本発明を用いれば、上記「Ｘ線ＣＴ撮影」に関し（ ＩＶ Ｘ線ＣＴ撮影に本発明を適用する場合２；（ノイズを低減した１枚の画像を得る））において述べたように、むしろ同じ撮影時間の間に、短時間の撮影を繰り返し行って、複数枚（Ｋ枚）の画像を得た方が良い。なぜならば、（短時間撮影で光量が少なくなるために、これらＫ枚の画像はランダムなノイズを多く含むにもかかわらず、）これらＫ枚の画像に対してコヒーレント・フィルタをかけて、ランダムなノイズをより強く抑制した１枚の画像を構成することができるからである。さらに、これらＫ枚の画像がそれぞれ例えばカラー画像である場合に、これを「光の三原色によって得られる複数種類（３種類）の画像が複数組（Ｋ組）ある」ものとして扱っても良い。すなわち、それぞれの画像を構成する各画素ｘについて、光の三原色の明るさを要素とした３次元ベクトルを構成でき、従ってＫ枚の画像全体では３Ｋ次元ベクトルの画素値ｖ（ｘ）を構成できる。この３Ｋ次元ベクトルの画素値を対象にして、コヒーレント・フィルタをかけてもよい。

また、ビデオカメラのように、動画像を撮影する撮像装置においては、（ ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（ノイズを低減した複数枚の画像を得る））で述べたダイナミック・コヒーレント・フィルタ処理あるいはリアルタイム・コヒーレント・フィルタ処理を適用することができる。さらに、上記と同様に「光の三原色によって得られる複数種類（３種類）の画像が複数組（Ｋ組）ある」ものとして扱って、コヒーレント・フィルタをかけることもできる。

なお、上記したデジタルカメラの他、ビデオカメラ、暗視カメラ、高速度撮影カメラ等（以下、一般に「撮像装置」という。）により撮影された画像においても、上記と同様にして、本発明を適用することが可能である。特に、暗視カメラ及び高速度撮影カメラについては、本来、１フレームあたりの光量が少ない状況で使用されるカメラであり、画像に含まれるランダムなノイズが多いが、空間解像度はもちろんのこと、ことに高速度撮影カメラにおいては時間分解能を損なうのは好ましくない。従って、本発明を適用するのに最も好適な装置の一つである。

したがって、上記撮像装置においても、このような撮影方法により、よりノイズを抑制した画像を得ることができる。

一方、「撮像」するのではなく、何らかの記録媒体に記録され又は送信される信号に基づく画像の再生が可能な装置、例えばビデオテープレコーダやＤＶＤプレーヤ等において、再生される画像に対して本発明に係るコヒーレント・フィルタを適用すれば、ランダムなノイズを抑制した再生画像を得ることができる。

さらに、これらの装置において、静止画像（いわゆる「ポーズ画像」、あるいは印刷用の画像）を生成・表示する場合にも、コヒーレント・フィルタを適用することによって、ランダムなノイズを抑制した静止画像が得られ、印刷等にも適している。このためには（ ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（複数枚の画像からノイズを除去する））で述べた方法、（ Ｖ Ｘ線ＣＴ撮影に本発明を適用する場合３；（複数種類の画像のノイズを低減する））で述べた方法、後述する（ ＶＩＩＩ １枚の画像から画素値ｖ（ｘ）を構成する ）に示す方法、あるいはこれらを組み合わせて用いることができる。

（ ＶＩＩＩ １枚の画像から画素値ｖ（ｘ）を構成する ）
以上までの適用例は、図２に示す分類のうち、画素値ｖ（ｘ）を、「複数の画像から」及び「複数種類の画像から」、それぞれ構成するものに該当するものであった。本発明では、これら二つの画素値ｖ（ｘ）の構成法の他、図２に併せて示すように、「１枚の画像から画素値ｖ（ｘ）を構成」することも考えられ、これは以下に記すように行うことができる。

単一の２次元画像、例えば上記Ｘ線ＣＴ装置１００やＭＲＩ装置等により取得される単一のＣＴ画像、あるいは上記撮像装置等によるモノクロの撮影による単一の画像においては、画素ｘの値はただ一つのスカラー値により構成される。言い換えれば、ベクトル値たる画素値ｖ（ｘ）は、１次元ベクトルとなる。

このような場合において本発明に係るコヒーレント・フィルタを適用するための一つの方法は、画素値ｖ（ｘ）を、上記「１次元ベクトル」そのものとして考える。すなわち、上記（３）式において、Ｋ＝１とし、

とする。後は、上述した手順により、同様にしてコヒーレント・フィルタの構成、ないし画像ノイズの抑制を実行することができる。本発明は、このように「ベクトル値」としての「画素値」ｖ（ｘ）が、１次元ベクトル、すなわちスカラー値である場合も含む。

また、この場合における重みｗ３については、例えば、次式を用いることができる。

ここで、Ｄは、（６）式におけるパラメータＣと同様な性質を有する。すなわち、Ｄが小さければ、重みｗ３（ｐ（ｘ，ｙ））はより小さく、Ｄが大きければより大きくなる。

このような重みｗ３（ｐ（ｘ，ｙ））及び上記（４）式により求められる画素値ｖ′（ｘ）に基づく画像は、例えば図８に示すようになる。図８（ａ）（ｂ）（ｃ）（ｄ）は、血管の３Ｄ抽出例を示す造影された画像であり、図８（ｅ）は他の抽出例として脳の灰白質を抽出した画像である。図８（ａ）は、コヒーレント・フィルタによる処理の対象となるＣＴ画像（元の画像）を示す。図８（ｂ）（ｄ）は、（６）式の中で用いられる上記集合Ｎ（ｘ）として３×３画素分の矩形状エリアを用いてコヒーレント・フィルタをかけた結果を示し、図８（ｂ）では上記の式（１４）においてＤ＝９とした場合の結果、図８（ｄ）では（ｃ）のコヒーレント・フィルタをかけた画像に対して更に閾値処理を施した結果を示している。また図８（ｃ）では、（ａ）の元の画像に対して閾値処理を施した画像を示している。

図８（ａ）と図８（ｂ）の画像を比較してみると、図８（ｂ）の画像は、コヒーレント・フィルタをかけたことにより、図８（ａ）の画像に比べてノイズが適度に抑制され、診断として十分な品位が保証されていることが分かる。

図８（ｃ）と図８（ｄ）の画像を比較してみると、元の画像に対して閾値処理しか施されていない図８（ｃ）の画像は、閾値による血管抽出が困難であるのに対し、コヒーレント・フィルタをかけて閾値処理した図８（ｄ）の画像は、閾値による血管抽出ができ、細い血管も失われていないことが分かる。

また図８（ｅ）によれば、灰白質の領域抽出が出来ていることが分かる。

ところで、このような方式によると、確かに上記図８に示すような効果が得られるものの、画素ｘと画素ｙの適合度を計測するのに１次元ベクトル値ｖ（ｘ）＝ｖ（ｘ）を用いているため、適合度の判別性能が比較的低く、言い換えれば、統計的に不確かな危険率ｐ（ｘ，ｙ）しか得られない。このため、解像度を損なわずに行えるランダムなノイズの抑制の効果は、さほど高いとは言えない。したがって、このような場合においてより好ましくは、次のような別の方式を採るとよい。すなわち、画素ｘの近傍のいくつかの画素からなる集合Ｚ（ｘ）を用いて、多次元のベクトルを構成し、これを画素ｘの画素値ｖ（ｘ）として用いるという方式が有効である。例えばその一例として、図９に示すように、画素ｘ＝（ｘ，ｙ）自身のスカラー値、すなわちｖ（ｘ）＝ｖ（ｘ，ｙ）に加え、該画素ｘに上下左右で隣接する４個の画素（ｘ＋１，ｙ），（ｘ−１，ｙ），（ｘ，ｙ＋１），（ｘ，ｙ−１）からなる集合Ｚ（ｘ）を用いて（以下、これら画素ｘ周囲の集合Ｚ（ｘ）の要素となる画素を「ｚ_１（ｘ），ｚ_２（ｘ），ｚ_３（ｘ），ｚ_４（ｘ）」と表す。）、これら４個の画素が持つスカラー値を付け加えた５次元ベクトルを考える。すなわち、

である。このようなベクトルは、図９に示す画像Ｇを構成するすべての画素について、同様に構成することができる（なお、図９は、全体の画像の一部のみを拡大して示している。）。

このようにして、１枚の画像から構成された多次元ベクトルの画素値ｖ（ｘ）を用いて、本発明に係るコヒーレント・フィルタを構成することによって、１次元ベクトルの画素値を用いる場合に比べて、高精度で適合度を判定できるようになる。その結果、画像の解像度を損なうことなく、ランダムなノイズをより強く抑制できる。

なお、このようにして画素値を構成した場合に、画素の適合度を定量化するために用いる上記した帰無仮説Ｈ「画素ｘと画素ｙとは、それぞれのノイズを除去した場合に同一の画素値を有する」とは、帰無仮説Ｈ′「画素ｘ及びその周囲の画素ｚ_１（ｘ）〜ｚ_４（ｘ）の合計５個の画素と、画素ｙ及びその周囲の画素ｚ_１（ｙ）〜ｚ_４（ｙ）の合計５個の画素とが、それぞれ、ノイズを除いて対応する画素と同一の画素値を有する」という意味に他ならない。

より具体的にみれば、図９に示すように、ある画素ｘと集合Ｚ（ｘ）（これは、上述したように例えば画素ｘ周囲のエリアである。）内に存在するある画素ｙとが共に、例えば、図に示されるような画像Ｇ上に現された何らかの物体Ｅ１の領域Ｅ１の内側に存在する場合には、例えば（９）式等によって計算される危険率ｐ（ｘ，ｙ）の値は小さく、例えば（６）式等によって計算される重みｗ（ｐ（ｘ，ｙ））の値は大きくなる。なぜならこの場合、画素ｚ_２（ｘ）及びｚ_４（ｘ）並びに画素ｚ_２（ｙ）及びｚ_４（ｙ）とは、当該物体の像Ｅ１の領域にあり、画素ｚ_１（ｘ）及びｚ_３（ｘ）並びに画素ｚ_１（ｙ）及びｚ_３（ｙ）とは、物体Ｅ２の像内にあるから（３．３）式で表わされるベクトル値Ｖ（ｘ），Ｖ（ｙ）は類似しているためである。

さて、上記の実施例においては、特に物体の像Ｅ１の縁にある画素ｘに注目すると、帰無仮説Ｈ′を棄却する場合の危険率ｐ（ｘ，ｙ）が小さいような画素ｙは、一般に、あまり多くは見つからないだろうことが推測される（例えば図９において、いずれも物体Ｅ１の像の縁Ｅ１Ｅの内側に存在する画素ｘと画素ｙ_ｐでも、帰無仮説Ｈ′を棄却する場合の危険率ｐ（ｘ，ｙ）は、同図に示す画素ｘ及びｙの場合に対して比較的大きくなる）。

より一般的に言えば、ある物体の像の縁に存在する画素ｘに関し、ほとんどのｙ∈Ｎ（ｘ）について、例えば（９）式等によって計算される危険率ｐ（ｘ，ｙ）の値は高く、例えば（６）式等によって計算される重みｗ（ｐ（ｘ，ｙ））は小さくなる。すると、重み付け平均をしても平滑化効果が弱いために、例えば（４）式によって計算されるコヒーレント・フィルタによるノイズ抑制効果は乏しい。

そこで、上記方式を更に改良した別の方式を以下説明する。すなわち、本方式では、画素値ｖ（ｘ）を、

と構成するのである。ただし、この（３．４）式におけるｖ_１〜ｖ_４は、上記（３．３）式におけるｖ（ｘ＋１，ｙ），ｖ（ｘ−１，ｙ），ｖ（ｘ，ｙ＋１）及びｖ（ｘ，ｙ−１）を一定の規則に従って（例えば値の大きい順に）並べ直したものである。

このようにすると、上記帰無仮説Ｈ′は、帰無仮説Ｈ′′「画素ｘ及びその周囲の画素ｚ_１（ｘ）〜ｚ_４（ｘ）の合計５個の画素と、対応する画素ｙ及びその周囲の画素ｚ_１（ｙ）〜ｚ_４（ｙ）の合計５個の画素とが、それぞれ、〔画素ｚ_１（ｘ）〜ｚ_４（ｘ）内及びｚ_１（ｙ）〜ｚ_４（ｙ）内でその順番を入れ替えることを許しつつ、〕ノイズを除いて同一の値を持つ」という、より緩やかなものとなる（〔〕内が付加された。）。

そして、この場合、例えば図１０に示す物体の像Ｅ２の縁Ｅ２Ｅの内側に存在する画素ｘと画素ｙとの関係においても、上記帰無仮説Ｈ′′を棄却する場合の危険率ｐ（ｘ，ｙ）は小さくなる（なぜなら画素ｚ_２（ｘ）、ｚ_３（ｘ）及びｚ_４（ｘ）並びに画素ｚ_１（ｙ）、ｚ_２（ｙ）及びｚ_４（ｙ）が上記（３．４）式におけるｖ_１乃至ｖ_３として一致し、画素ｚ_１（ｘ）及び画素ｚ_３（ｙ）が同じくｖ_４として一致する）から、結果、上記（３．３）式による画素値ｖ（ｘ）の構成によるよりも、ノイズ抑制効果が大きくなり、しかも空間解像度は損なわれない。

以上述べたように、本発明においては、１枚の画像からベクトル値である画素値ｖ（ｘ）を構成することができ、従って、コヒーレント・フィルタを用いて空間解像度を損なうことなくノイズを抑制できる。このような方式は、１枚の画像のみあればよいから、基本的にどのような画像に対しても一般的に適用可能である。すなわち、上記した各種適用例（Ｘ線ＣＴ装置、医用画像診断装置、あるいはカラー画像等）は勿論のこと、あらゆる画像に対して、本方式を適用することができる。また、画像が動画像である場合にも（上で何度か述べたように、それは多数枚の静止画像とみることができるから）、その静止画像の一枚一枚につき、本方式（「１枚の画像から画素値ｖ（ｘ）を構成する」）を適用してよいことも勿論である。さらに、本方式を３次元あるいは多次元の画像に用いる場合Ｚ（ｘ）は例えばχ_１に接する６個のピクセルとχ_１自身から構成することができるなど、拡張する方法もごく容易に推考できるので説明するまでもない。

なお、本発明においては、図２に示した「複数の画像から」、「複数種類の画像から」又は「１枚の画像から」の、それぞれ単独の方式により構成された画素値ｖ（ｘ）で、（既に述べたように）コヒーレント・フィルタを構成してよいことは勿論であるが、さらに、これら各方式を適宜組み合わせて、より高い次元のベクトル値として表される画素値ｖ（ｘ）を構成し、これに基づいてコヒーレント・フィルタを構成することも可能である。

既にＶＩＩで述べたように、例えば、同一の被写体に関する２枚のカラー画像があるとすれば、ある画素ｘに対して、上記（２．２）式右辺に示されるスカラー値ｖ_Ｒ（ｘ），ｖ_Ｇ（ｘ）及びｖ_Ｂ（ｘ）が、二枚の画像それぞれから、画素ｘの画素値ｖ（ｘ）として、

という６次元のベクトルを得ることができる。このように、より次元が高いベクトルで表される画素値ｖ（ｘ）を用いれば、例えば（９）式及び（６）式等によって計算されるような危険率ｐ（ｘ，ｙ）として、より確実な値を算出できる。その結果、画像の解像度を損なわずに、ランダムなノイズを一層強く抑制することが可能になる。

（ ＩＸ 本実施形態におけるコヒーレント・フィルタの最も一般的な形態 ）
さて以下では、上記図２に沿った説明から離れ、上記（４）式のように導入されたコヒーレント・フィルタに関する、より一般的な形態について説明する。

上記コヒーレント・フィルタの一般的形態（（４）式）、あるいは各種適用例においては、その最初に述べた帰無仮説Ｈ「画素ｘと画素ｙとはそれぞれのノイズを除去した場合に同一の値を有する」ないし「ｖ（ｘ）＝ｖ（ｙ）。ただし、両画素のノイズに起因する差異を除く」が用いられているが、これは、次に記すような、より一般的な表現とすることができる。

すなわち、Ｍ次元の未知のパラメータａ＝（ａ_１，ａ_２，…，ａ_Ｍ）（ただし、Ｍ＜Ｋ）を含むモデル関数ｆ（ａ，ｖ（ｙ））を導入することにより、画素の位置を表すある二つのベクトルｘ及びｙにおいて、帰無仮説Ｈ_０「ｖ（ｘ）＝ｆ（ａ，ｖ（ｙ））＋ξ（ただし、ξは（既知の）確率分布に従う。）」を立てる。このとき、パラメータａの推定値ａ^〜を、ｖ（ｘ）及びｖ（ｙ）から適当な当てはめ（ｆｉｔｔｉｎｇ）法（例えば、最小二乗法等、後述）を用いて決定することによって、上記帰無仮説を、それと等価な命題であるＨ_０「ｖ（ｘ）＝ｆ（ａ^〜，ｖ（ｙ））＋ξ、（ただしξは（既知の）確率分布に従う。）」に言い換えることができる。ここで、上記ｆ（ａ^〜，ｖ（ｙ））とは、すなわち、画素ｙの画素値ｖ（ｙ）に対して、関数ｆの記述するモデルが許す自由度（すなわちパラメータａ）を最適に調節して、画素ｘの画素値ｖ（ｘ）と最も高い適合度を持つように変換したものである。そして、そのような最適なパラメータがａ^〜である。

その結果、（もし帰無仮説が正しいならばξが従うと想定されるところの）既知の確率分布を用いて、危険率ｐ（ｘ，ｙ）が具体的に計算可能となり、したがって重みｗ（ｐ（ｘ，ｙ））が計算できるから、重み付け平均ｖ′_ｋ（ｘ）を

によって計算することができる。これが、本実施形態における最も一般的なコヒーレント・フィルタの形式である。なお、上記（４）式とこの（１５）式とを比較すれば、前者は後者において、実際に

（ ＩＸ−１ 上記（１５）式に基づく応用例；ダイナミック像に基づく時間濃度曲線の構成 ）
コヒーレント・フィルタの（１５）式に示した形式を用いた、本発明の別の適用例について述べる。ここでは、当該適用例として、上記医用画像診断装置により得られたダイナミックＣＴ像を対象とする時間濃度曲線（ｔｉｍｅ−ｄｅｎｓｉｔｙ ｃｕｒｖｅ）の定量測定を行う際に、本発明を適用することによって、ランダムなノイズを十分に抑制した前記時間濃度曲線を得ることを目的とする実施例について説明する。

まず、ここにいう「ダイナミック像」とは、同一の被写体を反復して撮影した結果得られる、一連の画像（画像番号１，２，・・・，Ｋ）と、それぞれの撮影時刻（ｔ_１，ｔ_２，・・・，ｔ_Ｋ）の組であり、例えば、上記医用画像診断装置、すなわち上記Ｘ線ＣＴ装置１００により得られたダイナミックＣＴ画像、あるいはＸ線透視装置、シンチレーション・カメラ、ＭＲＩ装置、ＳＰＥＣＴ装置又はＰＥＴ装置、超音波診断装置等を用いて撮影される。また、ダイナミック像は、上記一般の撮影装置、すなわちデジタルカメラ、ビデオカメラ、暗視カメラ、高速度撮影カメラ等によって繰り返し撮影を行った結果得られる、一連の画像を用いて構成したものでも良い。

また、ここにいう「時間濃度曲線」とは、上記ダイナミック像中の特定の部位における像の濃度値の経時的変化を表す曲線である。ことに、上記医用画像診断装置においては、人体組織等における血流動態や代謝機能等の詳細を調べる事を目的として、人体の特定組織内の造影剤濃度等の経時的変化を時間濃度曲線として計測することが行われている。また、天体観測等においては、特定の天体の光度変化等を解析する目的で、時間濃度曲線が用いられる。より形式的に明示すると、すなわち、時間濃度曲線とは、時刻ｔ_ｋにおけるある部位の濃度値をｄ_ｋとするとき、対の列｛＜ｔ_ｋ，ｄ_ｋ＞（ｋ＝１，２，・・・，Ｋ）｝として表現される。また、時間濃度曲線の多くの用途においては、必ずしもｄ_ｋの絶対的な値が必要なのではなく、むしろ最初の画像１を基準とする増分（ｄ_ｋ−ｄ_１）だけが得られれば十分である。さらにそのような用途のうちの多くでは、単に（ｄ_ｋ−ｄ_１）に比例するデータＡ（ｄ_ｋ−ｄ_１）（ここにＡは未知の比例係数）だけが得られれば十分である。この場合には、従って、対の列｛＜ｔ_ｋ，Ａ（ｄ_ｋ−ｄ_１）＞（ｋ＝１，２，・・・，Ｋ）｝が、求める時間濃度曲線である。

このような時間濃度曲線を求めるためには、原理的は、上記ダイナミック像を構成する各画像ｋ（ｋ＝１，２，・・・，Ｋ）における、該時間濃度曲線を測定しようとする部位に含まれる画素ｘのスカラー値ｖ_ｋ（ｘ）を用いて、対の列｛＜ｔ_ｋ，ｖ_ｋ，（ｘ）＞｝あるいは、｛＜ｔ_ｋ，Ａ（ｖ_ｋ（ｘ）−ｖ_１（ｘ））＞｝を構成すればよい。

しかし、実用においては、上記医用画像診断装置等によって撮影されたダイナミック像にランダムなノイズが含まれているために、本来測定しようとする時間濃度曲線を正確に求められないという問題がある。

さらに、実用においては、これらのダイナミック像においては、いわゆる「パーシャル・ボリウム効果」が生じる。パーシャル・ボリウム効果とは、すなわち、被検体内の微小な物体の像は、画像上では少数個の画素によって表現されるが、これら少数個の画素には、被検体内の隣接する物体の像も影響を与えるため、これら少数個の画素の画素値は（本来計測しようとする濃度値の変動に比例するものの）比較的小さな変動しか示さない、という現象である。言い換えれば、これら少数個の画素の画素値は僅かな信号しか含まない。従って、パーシャル・ボリウム効果が生じている場合には、どの画素ｘを取っても対の列｛＜ｔ_ｋ，ｖ_ｋ（ｘ）＞（ｋ＝１，２，・・・，Ｋ）｝は非常に信号レベルが低く、本来計測しようとしているのではない組織における濃度値の変化の影響を受け、さらにランダムなノイズが存在するために、本来測定しようとする時間濃度曲線｛＜ｔ_ｋ，ｄ_ｋ＞｝を正確に求められないという問題がある。

そこで、従来は、ランダムなノイズを抑制するために、まず計測しようとする部位に含まれる画素ｘに関する、対の列｛＜ｔ_ｋ，ｖ_ｋ（ｘ）＞（ｋ＝１，２，・・・，Ｋ）｝を構成し、次に時間平均値ｖ_ｋ＃（ｘ）を例えば

によって計算し、またこれに対応する時刻の時間平均値を

によって計算して（ここにｍ＜ｎは適当な整数、Ｇ_ｊは適当な重み係数）、時間濃度曲線｛＜ｔ_ｋ＃，ｖ_ｋ＃（ｘ＞（ｋ＝１，２，・・・，Ｋ）｝を構成し、これを時間濃度曲線として用いる、という時間平均による方法が用いられている。

あるいは、従来は、ランダムなノイズを抑制するために、計測しようとする部位に概ね相当する画素の集合Ｒ（すなわち、画像上の関心領域（ＲＯＩ；Ｒｅｇｉｏｎ ｏｆ Ｉｎｔｅｒｅｓｔ））を設定し、画像番号ｋにおいて、この集合に含まれる全画素ｘ∈Ｒのスカラー値ｖ_ｋ（ｘ）の空間平均値ｖ_ｋ（Ｒ）を下式（１９）によって求めて、対の列｛＜ｔ_ｋ，ｖ_ｋ（Ｒ）＞（ｋ＝１，２，・・・，Ｋ）｝を構成し、これを時間濃度曲線として用いる、という空間平均による方法が用いられている。

制する方式も用いられている。

しかし、これら従来の方式によると、時間平均を行うと時間分解能が損なわれてしまい、また、空間平均を行うと、本来の該測定しようとする部位以外の部位の濃度の経時変化が計測値に混入するという問題点があった。

このような問題点を解決し、より正確な時間濃度曲線を得るために、本発明に係るコヒーレント・フィルタを、上記（ ＩＸ 本実施形態におけるコヒーレント・フィルタの最も一般的な形態 ）にもとづいて適用することができる。

まず、本実施形態のコヒーレント・フィルタにおいて用いるべき、帰無仮説について説明する。計測しようとする部位における真の時間濃度曲線を｛＜ｔ_ｋ，ｄ_ｋ＞（ｋ＝１，２，・・・，Ｋ）｝であると仮定するとき、その一次変換である｛＜ｔ_ｋ，Ａ（ｄ_ｋ−ｄ_１）＞（ｋ＝１，２，・・・，Ｋ）｝（だたしＡは未知の係数）を計測することを目的とする場合において、計測しようとする部位に概ね相当する画素の集合Ｒを設定する。この集合Ｒの要素である任意の画素ｘ∈Ｒについて、条件Ｑ：「もし、この画素ｘが上記真の時間濃度曲線を良く反映し、しかも他の部位の経時的濃度変化の影響をほとんど受けていない」のであれば、（ベクトル値としての）画素値ｖ（ｘ）＝（ｖ_１（ｘ），ｖ_２（ｘ），．．．，ｖ_ｋ（ｘ））について、パーシャル・ボリウム効果およびランダム・ノイズの影響を考慮することによって、

が成り立つと仮定することができる。ここに、ｐ（ｘ）およびｑ（ｘ）は、画素ｘごとに異なるが画像番号ｋ（すなわち撮影時刻ｔ_ｋ）によっては変化しない未知の係数であり、パーシャル・ボリウム効果をモデル化したものである。またγ_ｋ（ｘ）はランダムなノイズをモデル化したものであって、画素ｘごとに、しかも画像番号ｋごとに値が異なるが、その期待値は０であり、またその統計分布は画素ｘにも画像番号ｋにも依存しない（なお以下では、説明のため、該統計分布が平均０、分散σ^２の正規分布である場合の例を用いることにするが、この統計分布が任意の既知の分布によって概ね近似されるのであれば良く、その場合には、以下の実施形態を適合するように改変する方法は自明である。）。

以上の仮定によれば、該集合Ｒの要素である任意の２個の画素ｘ，ｙに関して、もし「画素ｘ，ｙが共に（上記の）条件Ｑを満たす。」という命題が成り立つのであれば、次式の関係が成り立つことが証明できる。

ここに、ａ_１およびａ_２は、画素の組ｘ，ｙごとに異なるが画像番号ｋ（すなわち撮影時刻ｔ_ｋ）によっては変化しない未知の係数である。またξ_ｋはランダムなノイズであって、画素の組ｘ，ｙごとに、しかも画像番号ｋごとに値が異なるが、その期待値は０である。

（２１）式は以下のようにして導かれる。すなわち、（１９）式においてｘにｙを代入して得られる式

を得る。従って、

とおくことによって、（２１）式が導かれる。ここで、（２４）式のａ_１とａ_２は、パーシャル・ボリウム効果を現すパラメータであり、また（２４）式のξ_ｋは、ランダムなノイズを表す。

以上から、「画素ｘ，ｙが共に条件Ｑを満たす。」という命題は、帰無仮説Ｈ_０′「ｖ_ｋ（ｘ）＝ａ_１ｖ_ｋ（ｙ）＋ａ_２＋ξ_ｋ （ｋ＝１，…，Ｋ）である。」と等価であることが示された。

次に帰無仮説Ｈ_０′「ｖ_ｋ（ｘ）＝ａ_１ｖ_ｋ（ｙ）＋ａ_２＋ξ_ｋ （ｋ＝１，…，Ｋ）である。」を、実質的に等価であり、かつ実際に検定できる形式の命題に変換する方法について述べる。この帰無仮説を改めて数学的に厳密な表現で述べると、帰無仮説Ｈ_０′「ある定数ａ_１およびａ_２が存在して、ξ_ｋ＝ｖ_ｋ（ｘ）−ａ_１ｖ_ｋ，（ｙ）−ａ_２（ｋ＝１，…，Ｋ）は平均０、分散（σ^２ｈ（ａ_１））の正規分布に従う。」となる。ここに係数ｈ（ａ_１）は、

である。（（２５）式はａ_１とξ_ｋの定義である（２４）式、および、ランダム変数に関する分散の持つ一般的な性質から直ちに導かれる。）また、上記の分散σ^２の値は、前述の（ ＩＩＩ Ｘ線ＣＴ撮影に本発明を適用する場合１；ダイナミックＣＴ撮影（ノイズを低減した複数枚の画像を得る））等において述べた方法で、簡便かつ実用上十分正確に推定できる。

以上から、もし、上記の定数ａ_１およびａ_２を決定することができれば、上記の帰無仮説Ｈ_０′を検定することが可能である。そして実際上は、これらの定数の最適な推定値ａ_１ ^〜およびａ_２ ^〜が得られれば十分である。

このような、定数ａ_１およびａ_２の最適な推定値の算出には、公知の当てはめ法（ｆｉｔｔｉｎｇ）がそのまま利用できる。そこで、以下では、そのような当てはめ法の典型的な具体例として、線形最小二乗法を用いる場合における概要を説明する。線形最小二乗法を本実施例に適用するには、単に、上記の帰無仮説のξ_ｋの二乗和をＳ（ａ）として、すなわち

を定義する。Ｓ（ａ）の値は定数ベクトルａ＝（ａ_１，ａ_２）、すなわち上記の定数ａ_１およびａ_２、の値に依存する。このＳ（ａ）が最小の値を取るような定数ベクトルａを算出すれば、定数ａ_１およびａ_２に関する、不偏推定の意味での最適な推定値ａ_１ ^〜およびａ_２ ^〜が得られる。なお、線形最小二乗法の具体的な計算方法としては、様々な公知の方法を利用することができ、しかも、これら公知の計算方法はいずれも非常に簡単であり、必要な計算時間はごく僅かである。

このようにして、上記の定数ａ_１，ａ_２の最適な推定値ａ_１ ^〜，ａ_２ ^〜を算出した結果、次式で定義される残差

を具体的に計算することができる。従って、この残差ｒ_ｋ ^〜を用いて、上記の帰無仮説Ｈ_０′を、実質的に等価な帰無仮説Ｈ_０”「ｒ_ｋ ^〜（ｘ，ｙ）（ｋ＝１，…，Ｋ）は平均０、分散（１＋（ａ_１ ^〜）^２）σ^２の正規分布に従う。」と言い換えることができる。これは、実際に検定の計算を実行可能な具体的命題である。

なお、さらに、ベクトルによる表現

（ただし、ベクトルａ及びξは画素の組ｘ，ｙに依存する。）を導入し、また、次式

を言い換えると、帰無仮説Ｈ_０”は「ｖ（ｘ）＝ｆ（ａ^〜，ｖ（ｙ））＋ξ、（ただし、ξは平均０、分散（１＋（ａ_１ ^〜）^２）σ^２の正規分布に従う。）」となり、これは（ ＩＸ 本実施形態におけるコヒーレント・フィルタの最も一般的な形態 ）において述べた帰無仮説Ｈ_０と全く同じ形式である。すなわち、本実施形態は本発明に係るコヒーレント・フィルタの一変形例であることは明らかである。なお、ここで、上記ｆ（ａ^〜，ｖ（ｙ））とは、すなわち、画素ｙの画素値ｖ（ｙ）に対して、パーシャル・ボリウム効果を現すパラメータａを最適に調節して、画素ｘの画素値ｖ（ｘ）と最も高い適合度を持つように変換したものを意味する。

次に、本実施形態において、上記の帰無仮説Ｈ_０”を用いて、コヒーレント・フィルタによって時間濃度曲線を求める方法について説明する。計測しようとする部位に概ね相当する画素の集合Ｒについて、この集合Ｒに含まれるあるひとつの画素ｘ∈Ｒについて、集合Ｒの要素である全ての画素ｙ∈Ｒに対して、以下の計算を行う。すなわち、上記の方法を用いて実際に残差ｒ_ｋ ^〜（ｘ，ｙ）（ｋ＝１，…，Ｋ）を算出し、次に、上記の帰無仮説Ｈ_０”「ｒ_ｋ ^〜（ｘ，ｙ）（ｋ＝１，…，Ｋ）は平均０、分散（１＋（ａ_１ ^〜）^２）σ^２の正規分布に従う。」を棄却する場合の危険率ｐ（ｘ，ｙ）ないし重みｗ（ｐ（ｘ，ｙ））を（６）式により具体的に計算する。そして、重み付き平均ｖ′_ｋ（ｘ）を下式（１５’）によって計算し、画素ｘにおける時間濃度曲線｛＜ｔ_ｋ，ｖ′_ｋ（ｘ）−ｖ′_１（ｘ）＞（ｋ＝１，２，・・・，Ｋ）｝を構成する。

こうして得られた時間濃度曲線は、画素ｘにおける真の時間濃度曲線｛＜ｔ_ｋ，ｄ_ｋ＞｝の一次変換である｛＜ｔ_ｋ，Ａ（ｄ_ｋ−ｄ_１）＞｝（だたしＡは未知の係数）を近似している計測値であり、しかも、（１５）式による重み付き平均の効果によって、ランダムなノイズが抑制されている。また、（１５）式による計算に用いる他の画素ｙの画素値ベクトルに対しては、式から明らかなように、パーシャルボリウム効果の影響を補正したものが用いられている。さらに、本実施形態はコヒーレント・フィルタの共通の特徴である、「時間平均を全く使用せず、また空間平均を画素ｘとの適合度に基づく重みを使って計算する」という性質を有する。従って、本実施形態によって、時間分解能を損なわず、パーシャル・ボリウム効果の影響を抑制し、しかもランダムなノイズが抑制された時間濃度曲線を得ることができる。なお、このようにして時間濃度曲線を求める方式を、特に「コヒーレント・レグレッション法」と称す。

次に、具体的に、医療用のＸ線ＣＴにおけるダイナミックＣＴ撮影等で得られたダイナミック像における、時間濃度曲線の臨床的利用の一例を説明する。この応用例では、造影剤を血管に急速に注入しながら、ダイナミックＣＴ等の撮影を行い、人体組織中に存在する動脈の像の濃度変化を時間濃度曲線として計測することによって、当該組織における血流動態を診断しようとするものである。

この応用例において、多くの場合、人体組織中の動脈は一般に非常に細いために、ＣＴによる断層画像上に現れる動脈の像は、パーシャル・ボリウム効果を生じる。さらに、像にはランダムなノイズが含まれていることは言うまでもない。このため、従来の方法では、動脈に関する十分に正確な時間濃度曲線を得ることは困難であり、強いて計測を行えば、動脈に関する真の時間濃度曲線＜ｔ_ｋ，Ｄ_ｋ＞の一次変換である＜ｔ_ｋ，Ａ（Ｄ_ｋ−Ｄ_１）＞（ここにＤ_ｋは動脈の像に相当する一群の画素の、時刻ｔ_ｋにおける（スカラー値である）画素値を表す。また、ｋ＝１，２，・・・，Ｋ）をある程度近似する測定値＜ｔ_ｋ，（ｖ_ｋ（ｘ）−ｖ_１（ｘ））＞しか得られなかった。この測定値はランダムなノイズを含む。また、パーシャル・ボリウム効果の影響のために、係数Ａは未知のままである。

そこで、本発明に係る上記の方式を適用すれば、＜ｔ_ｋ，Ａ（Ｄ_ｋ−Ｄ_１）＞を十分に近似する測定値＜ｔ_ｋ，（ｖ′_ｋ−３（ｘ）−ｖ′_１（ｘ））＞（ｋ＝１，２，・・・，Ｋ）を得ることができる。一方、同じ断層画像上で観察できる静脈の中には、相当に太いものが存在し、従ってそれらの静脈に関しては、従来の方法で、時間濃度曲線の十分に良い近似値＜ｔ_ｋ，（Ｊ_ｋ−Ｊ_１）＞（ｋ＝１，２，・・・，Ｋ）を得ることができる。ここにＪ_ｋは静脈の像に相当する一群の画素の、時刻ｔ_ｋにおける画素値を表す。

ところで、血液循環に関する時間濃度曲線においては、命題Ｓ：「もし、時刻ｔ_１における血中の造影剤濃度が０であるならば、どの血管ｄに関する時間濃度曲線＜ｔ_ｋ，（ｄ_ｋ−ｄ_１）＞も、その曲線下面積（ＡＵＣ：Ａｒｅａ Ｕｎｄｅｒ Ｃｕｒｖｅ）が一致する」という性質が成り立つことが知られている。ここで言う曲線下面積とは、時間濃度曲線＜ｔ_ｋ，（ｄ_ｋ−ｄ_１）＞の時間ｔに関する積分を意味する。

従って、ある血管ｄに関する時間濃度曲線＜ｔ_ｋ，（ｄ_ｋ−ｄ_１）＞の曲線下面積ＡＵＣ（ｄ）は、例えば次式によって近似的に計算することができる。

従って、静脈に関して従来の方法で得られた時間濃度曲線｛＜ｔ_ｋ，（Ｊ_ｋ−Ｊ_１）＞｝に関する曲線下面積ＡＵＣ（Ｊ）を（３０）式を用いて計算することができる。（ｄにＪを代入すればよい。）また、動脈に関して、仮に、時間濃度曲線｛＜ｔ_ｋ，（Ｄ_ｋ−Ｄ_１）＞｝が知られていれば、曲線下面積ＡＵＣ（Ｄ）を（２６）式を用いて同様に計算することができ、しかも上記命題Ｓに従って

であるため、ＡＵＣ（Ｄ）は計算できない。

一方、本発明に係る方式で得られた時間濃度曲線＜ｔ_ｋ，（ｖ′_ｋ（ｘ）−ｖ′_１（ｘ））＞は、＜ｔ_ｋ，Ａ（Ｄ_ｋ−Ｄ_１）＞を近似するものであり、後者は未知の係数Ａを含んでいる。このため、｛＜ｔ_ｋ，（ｖ′_ｋ（ｘ）−ｖ′_１（ｘ））＞｝から（３０）式を用いて具体的に計算できる曲線下面積ＡＵＣ（ｖ′）は、ＡＵＣ（Ｄ）のちょうどＡ倍でなくてはならない。すなわち、

未知であった係数Ａの値が具体的に決定できる。そこで、この係数Ａの値を用いて時間濃度曲線＜ｔ_ｋ，（ｖ′_ｋ（ｘ）−ｖ′_１（ｘ））／Ａ＞を構成すれば、これは、動脈の時間濃度曲線＜ｔ_ｋ，（Ｄ_ｋ−Ｄ_１）＞を近似するものに他ならない。このように、曲線下面積を用いて、未知であった比例係数Ａの値を決定した時間濃度曲線を構成する方法を「ＡＵＣ法」と呼ぶ。

以上から、ダイナミックＣＴ撮影等で得られたダイナミック像における、時間濃度曲線の臨床的利用において、上記コヒーレント・レグレッション法に、さらに上記ＡＵＣ法を組み合わせることによって、従来の方法では計測が困難あるいは不可能であった、細い動脈の時間濃度曲線に関しても、パーシャル・ボリウム効果およびランダムなノイズの影響を排除し、しかも、未知の比例係数Ａを含まない測定値が得られる。

なお、もちろんＡＵＣ法は、単独で従来の方法で計測された動脈に関する時間濃度曲線＜ｔ_ｋ，（ｖ′_ｋ（ｘ）−ｖ′_１（ｘ））＞に対しても適用でき、（ランダムなノイズやパーシャル・ボリウム効果の影響は排除できないものの、）未知であった比例係数Ａの値を決定した時間濃度曲線を構成できる。

上記コヒーレント・レグレッション法に、さらに上記ＡＵＣ法を組み合わせて得られた時間濃度曲線は、例えばある画素ｘについて、図１１において記号（Ａ）で示すようなものとなる。なお、この図では、従来の方法で構成した時間濃度曲線に上記ＡＵＣ法を単独で適用したもの（Ｂ）も併せて示している。（Ｂ）ではランダムなノイズの影響が明らかに見られるのに対し、（Ａ）ではノイズが十分に抑制されており、しかも時間分解能が全く損なわれていないことが分かる。

なお、時間濃度曲線の濃度値として、上記の説明ではスカラー値ｄ_ｋを用いたが、本実施形態はこの場合に限定されるものではなく、たとえばダイナミック像を構成する個々の画像が、カラー画像であったり、より一般には多種類の画像の組であって、時間濃度曲線の濃度値がベクトル値として表現される場合においても、容易に拡張して適用できることは、言うまでもない。

（ Ｘ 本実施形態の補足事項 ）
以下では、上記までに述べた実施形態の補足事項につき説明する。
（ Ｘ−１ 一般的補足事項 ）
まず、上記実施形態においては、本発明にいう「適合度」を定量化する手段として、帰無仮説Ｈを棄却した場合の「危険率」ｐ（ｘ，ｙ）が想定されているが、本発明は、このような形態に限定されるものではない。最初に述べたように、「適合度」とは、画素ｘ及びｙが何らかの意味で類似であるか否かを表す数値的な『指標』であればよい。

また、これに関連して、上記では、上記各種の帰無仮説を棄却した場合の危険率を求めるのに、χ二乗検定法が用いられているが（（６）式、（９）式参照）、本発明はこれにも限定されない。換言すれば、本発明は、「統計的検定法」、あるいはその具体的形態の一種たる「χ二乗検定法」を利用する形態のみに限定されることはない。

さらに、上記実施形態では、重みｗ（ｐ（ｘ，ｙ））の具体的形式として、都合二種の重み関数ｗ１（（６）式、ないしは（９）式が代入された（１０）式）及びｗ３（（１４）式）が挙げられているが、本発明は、これにも限定されるものではない。既に述べたように、危険率ｐ（ｘ，ｙ）∈［０，１］の関数としての重み関数ｗが満たすべき性質は、ｗ（ｔ）が、定義域ｔ∈［０，１］で定義される非負の単調増加関数である、ということだけである。そして、最も一般的にいえば、重み関数ｗは、「適合度に関する非負の単調増加関数」、であればよい。

（ Ｘ−２ 適合度に関するより一般的な形態）
本発明にいう「適合度」は、上記したような、帰無仮説を棄却する場合の危険率ｐ（ｘ，ｙ）によって数値化を行う方式に限定されるものではない。従って「重み」も、前記危険率ｐ（ｘ，ｙ）に重み関数ｗを作用させて算出する方式に限定されるものではない。以下では、このような適合度等に関するより一般的な形態に関し、具体的な例に沿った説明を行うこととする。

まず、一般的に要約すると、本発明では、画像処理の目的に応じて適切に設定された、或る命題に関し、その命題が処理対象である画像において成り立つということの確からしさを、目的に応じた適切な方式で数値化し、その数値を該画像処理に適用するという構成をとる。例えば、一対の画素ｘ及びｙの適合度を数値化するために、本発明では、より緩やかな尺度を構成してもよい。そのような具体例として、曖昧論理（ｆｕｚｚｙ ｌｏｇｉｃ）等における特性関数（ｍｅｍｂｅｒｓｈｉｐ ｆｕｎｃｔｉｏｎ）を利用して、上記適合度を数値化する、等の構成が考えられる。

（ Ｘ−２−１ パターン認識に本発明を適用する例）
ここでは、適合度等に関するより一般的な形態の具体例として、画像上の特徴的なパターンを識別し抽出する、いわゆる「パターン認識」に、本発明を適用する場合の一つの例について説明する。

この場合には、あらかじめ、識別しようとするパターンの特徴をもとに、「画素ｘは当該パターンを構成する画素である」という命題の確からしさを表す関数ｍ（ｘ）を定義しておく。ただし、ｍ（ｘ）∈［０，１］とする。ｍ（ｘ）は、より詳細には、画素ｘの画素値であるベクトル値ｖ（ｘ）および、画素ｘの近傍にある画素の集合Ｎ（ｘ）に含まれる各画素ｙの持つ画素値ｖ（ｙ）、並びに、画素の位置を表すベクトルｘそのもの、などから、上記命題の確からしさを表す数値を算出するような関数である。従って、ｍ（ｘ）は抽出しようとするパターンと、抽出すべきでないパターンとの差異に応じて、適切に設計される必要がある。また、重み関数を

で定義する。そして、このｍ（ｘ）が、ある「閾値」以上であるような画素の集合Ｘを構成すれば、パターンに該当する画素ｘの領域が得られる。

さて、処理すべき画像Ｐが与えられると、画像の全ての画素ｘについて、具体的な閾値Ｔと前記ｍ（ｘ）を上記重み関数に代入したｗ（Ｔ，ｍ（ｘ））を計算する（これによって、ｗ（Ｔ，ｍ（ｘ））を特性関数とする集合Ｘ（すなわちｗ（Ｔ，ｍ（ｘ））＝１であるような画素ｘの集合）が定義されたことになり、このＸが画像Ｐ上で該パターンの占める領域に他ならない。）。さらに、新しい画像Ｐ’を以下のようにして生成する。すなわち、Ｐ’の画素ｘは、ｗ（Ｔ，ｍ（ｘ））＝１のときにはｖ（ｘ）を画素値とし、さもなければ０（ゼロベクトル）を画素値とする。こうして作られた画像Ｐ’は、画像Ｐから該識別しようとするパターンに該当する以外の画素の画素値を０によって消去したものになっている。

もちろん、ｍ（ｘ）の計算、ｗ（Ｔ，ｍ（ｘ））の計算、Ｐ’の構成の各処理をそれぞれ逐次行っても良いし、あるいは、これらの処理を一体化したソフトウエアによって計算を行うことも可能である。

この例においては、関数ｍ（ｘ）は、パターンの特徴と、画素ｘの画素値であるベクトル値ｖ（ｘ）およびｘそのものの値とを照合して、パターンと画素ｘの「適合度」を数値化する。そして（３４）式で定義される重み関数によって、重みｗ（Ｔ，ｍ（ｘ））が算出され、この重みと画像Ｐを使って、処理結果である新しい画像Ｐ’が構成される。

これに関する更に具体的な例として、航空写真や金属表面の顕微鏡写真等に写る細い線を抽出するための画像処理の構成を説明する。

初めに、「画素ｘは細い線を構成する画素である」という命題の確からしさを数値的に表す指標を構成する方法を説明する。

あらかじめ、画素ｘの近傍にある画素の集合Ｎ（ｘ）を定義する。典型的には、Ｎ（ｘ）は画素ｘを中心とする矩形エリアとして定義すればよい。そして、Ｎ（ｘ）の各画素における画素値（スカラー値あるいはベクトル値）を並べた多次元ベクトル値ｖ（ｘ）を定義する。このベクトルの次元をＫとする。典型的にはＫは２５〜１６９程度とするのが好適である。

次に、後の計算に補助的に用いるＫ次元ベクトル関数である、正規化関数ｎｒｍを次のように定義する。すなわち、ｖ≠０である任意のＫ次元ベクトルｖに対して正規化関数ｎｒｍは、

ただし、ｋ＝１，２，…，Ｋである。この正規化関数によって、任意のＫ次元ベクトルｖについて（ｖがゼロベクトルでない限り）、

となることが保証される。

次に、任意にひとつの画素ｘを決め、背景が０である画像であって、画素ｘを通る１本の細い直線の像だけがあるような画像の典型例を適当な個数Ｊだけ集め、それらに画像番号１，２，．．．，Ｊを付ける。画像番号ｊの画像におけるｖ（ｘ）の値を計算して、ベクトルｒ（ｊ）を、次式

によって具体的に構成する。この計算をｊ＝１，２，．．．，Ｊについて行う。（こうして構成されるＪ個のベクトルｒ（ｊ）を以下「パターンベクトル」と呼ぶことにする。なお、このパターンベクトルは、本発明にいう「別に構成された画素値」に該当する。）
また、関数ｍ（ｘ）を次式で定義する。

値を取り出す関数である。このように構成した関数ｍ（ｘ）は、画素ｘがパターンベクトルｒ（ｊ）（ｊ＝１，２，．．．，Ｊ）のいずれかに一致もしくは類似する場合に、大きい値をとり、さもなければ小さい値を取るので、「画素ｘは細い線を構成する画素である」という命題の確からしさを数値的に表す指標、すなわち「適合度」を表す関数となっている。さらに、重み関数を

で定義する。ここにＴは閾値となる定数であって、適当な値を設定する。

さて、処理すべき画像Ｐが与えられると、画像の全ての画素ｘについて、具体的な閾値Ｔと前記ｍ（ｘ）を上記重み関数に代入したｗ（Ｔ，ｍ（ｘ））を計算する。（これによって、ｗ（Ｔ，ｍ（ｘ））を特性関数とする集合Ｘ（すなわちｗ（Ｔ，ｍ（ｘ））＝１であるような画素ｘの集合）が定義されたことになり、このＸが画像Ｐ上で該パターンの占める領域に他ならない。）
さらに、新しい画像Ｐ’を以下のようにして生成する。すなわち、Ｐ’の画素ｘは、ｗ（Ｔ，ｍ（ｘ））＝１のときにはｖ（ｘ）を画素値とし、さもなければ０（ゼロベクトル）を画素値とする。こうして作られた画像Ｐ’は、画像Ｐから該識別しようとするパターンに該当する以外の画素の画素値を０によって消去したものになっている。

また、ｗ（Ｔ，ｔ）を曖昧な判別関数、例えば

（ここにＣは正の定数）
のように構成して、Ｐ’の画素ｘはｗ（Ｔ，ｍ（ｘ）Ｖ（ｘ）を画素値とするように構成してもよい。これは、曖昧な論理によって「細い線」にあまり該当しないと判断される画素のコントラストを低下させることによって「細い線」の部分が浮き立つように画像Ｐ’を構成する効果を生じる。

（ Ｘ−３ 分散σ^２の推定法に関する補足事項 ）
次に、上記（６）式、あるいは（１４）式に用いられる分散σ^２（ないし標準偏差σ）の推定方法に関する補足説明を行う。

上記「Ｘ線ＣＴ撮影に適用した場合」においては、上記（６）及び（７）式により得られた分散の期待値Ｅ［σ^２］が、Ｋ枚すべての画像上の全画素ｘについて妥当するという背景から、当該全画素ｘにつき一律に、当該Ｅ［σ^２］を使用する形態となっていた。しかしながら、このような前提が満たされない場合、すなわち各画素ｘ（ここでは、ｘ_１，…，ｘ_Ｆとする。）につき、（分布の形が同一（例えばガウス分布）であるとしても）固有な分散σ^２（ｘ_１），…，σ^２（ｘ_Ｆ）が一般に想定される場合もある。

このようなときは、上記（６）及び（７）式により、当該各画素ｘにつき、個別に、分散の期待値Ｅ［σ^２（ｘ_１）］，…，Ｅ［σ^２（ｘ_Ｆ）］を求めるようにすればよい。後は、このように個別に求められた期待値Ｅ［σ^２（ｘ_１）］，…，Ｅ［σ^２（ｘ_Ｆ）］を、それぞれの場合に応じて（＝着目している画素ｘ_１，…，ｘ_Ｆに応じて）、（３）式に代入し、それぞれ、ｖ′（ｘ）＝（ｖ′_１（ｘ），ｖ′_２（ｘ），…，ｖ′_Ｋ（ｘ））を求めればよい。

また、全画素の画素値に含まれるノイズが共通の分散σ^２を持つと想定する場合には、分散σ^２の推定を行うために次のような方式を採ることも可能である。まず、すべての画素ｘの集合Ω（＝｛ｘ_１，…，ｘ_Ｆ｝）について、上記（７）及び（８）式により求められた分散の期待値Ｅ［σ^２（ｘ）］の平均（σ^２）^＋を求め、これを（６）式におけるσ^２の推定値の第一近似とする。すなわち、

である。ただし、上式において、｜Ω｜は、集合Ωの要素の数（ここではすなわち、｜Ω｜＝Ｆである。）を表している。次に、Ｅ［σ^２（ｘ）］が上記平均値（σ^２）^＋の、例えば数倍以内であるような画素ｘの集合をＭとする。Ｍの要素であるような画素ｘだけを対象にして改めて平均をとったものを（σ^２）^＋＋とすると、

となる。そしてこれは、よりもっともらしいσの推定値として利用することができる。

これは例えば、具体的に次に記すような場合、すなわち図１２（ａ）に示すように、動画像を構成する二つの静止画像Ｇ１及びＧ２が存在しているが、図１２（ｂ）に示すように、これらの画像Ｇ１及びＧ２の差をとると、像Ｚ１及びＺ２が互いに若干ずれている（つまり、像Ｚ１が、像Ｚ２のように動いた）場合等に有効である。

上記（３５）式の分散（σ^２）^＋は、この図において、図１２（ｂ）に示されるような画像Ｇ３を構成するすべての画素ｘ（＝Ω）に基づき求められた分散に該当する。しかしながら、この場合においては、図１２（ｃ）に示すように、当該画像Ｇ３に関する確率分布図において、像Ｇ１及びＧ２とが互いに重ならない部分ζ_１及びζ_２とに起因する、例えば二つの山ζ_１′及びζ_２′を生じさせることになるから、当該分散（σ^２）^＋は過大に見積もられていることになり、従ってこのまま全静止画像につき使用するのは不適切である。一方、図１２（ｂ）に示す部分ζ_Ｍについては、図１２（ｃ）における中央の山ζ_Ｍ′が該当していると考えられるから、当該山ζ_Ｍ′を残し上記二つの山ζ_１′及びｇ_２′を取り除いた状態の確立分布図に関する分散を求める方が好ましい。そして、そのような分散こそが、上記（３２）式における分散（σ^２）^＋＋に該当する。

つまり、上記で「Ｅ［σ^２（ｘ）］が（σ^２）^＋の、数倍以内であるような画素ｘの」集合Ｍのみを対象として平均をとったのは、図１２（ｃ）における山ζ_Ｍ′のみを抽出して平均をとる処理に他ならず、その結果、（σ^２）^＋＋は図１２（ｃ）において、山ζ_１′及びζ_２′を取り除いた状態の確立分布図に関する分散を求めていることになる。

したがって、この場合においては、この（σ^２）^＋＋を、上記（５）式等に代入することにより、各画素ｘ_１，…，ｘ_Ｆについての、変換後の画素値ｖ′（ｘ）＝（ｖ′_１（ｘ），ｖ′_２（ｘ），…，ｖ′_Ｋ（ｘ））を、よりもっともらしい値を得ることができる。

なおまた、ノイズの推定に関し、一般的には、画像の隅に写る位置などに一様の明るさの被写体を置き、この部分に該当する画素の集合の値の分布を計測することによって、ノイズの分布を推定するという方法も有効である。

本発明において、ノイズの分布ないし分散の推定法は、基本的に如何なる方式によるものでもよく、ここに述べたいくつかの方式に限定されるものではない。実際的には、それぞれ各種の実施例において利用可能な先験情報、撮像過程理論、及び計測値ないし実測値等を用いて適宜、最も好適に求められる方式を採用すべきである。なお、この際、当該推定法は、それぞれの実施例において必要とされる実用的精度に合わせて、可能な限り簡易に構成することが、装置の構成を簡易化し、かつ、処理速度を向上する上で好ましい。

（ Ｘ−４ 本発明の更なる別の実施例）
最後に、上記で触れなかった、本発明の更なる別の実施例について補足説明する。
（ Ｘ−４−１ データ圧縮の前処理としての実施例 ）
一般に、ある画像Ｑにランダムなノイズを加えた画像Ｑ′を作ると、画像Ｑ′は画像Ｑに比べて遙かに大きい情報量を持つ（例外は、元の画像Ｑが非常に大きいノイズを含んでいる場合である。）。したがって、もしこのような画像Ｑ′を通信や記録媒体に記録しようとすれば、実質的な意味を持たないノイズを記述するための情報量が、その通信や記録におけるビット数を浪費することになる。逆に言えば、画像Ｑ′からノイズを抑制した画像Ｑ′′を作り、これを通信や記録の対象とすれば、小さな情報量で、同じ意味を持つ画像（しかもそれは鮮明である。）を表現できるようになる。

一般に画像を通信や記録の対象とする場合、当該画像は、いわゆる「圧縮処理」を受けることが多い。これは、可能な限り伝送又は記録すべき情報量を小さくして、高速通信又は（一つの記録媒体に対する）大量記録を可能とすることを目的としている。実用上、特に情報量を小さくしたいのは「動画像」の場合である。これらは非可逆画像圧縮技術において広く知られている事実である。

これらの事項を踏まえれば、ノイズを含んだ画像Ｑ′を圧縮するよりも、当該ノイズを抑制した画像Ｑ′′を圧縮した方が、通信や記録に必要なビット数（すなわち情報量）を小さくすることができる。

そこで、通信又は記録しようとする画像Ｑ′に対し、本発明に係るコヒーレント・フィルタをかけてそのノイズを抑制する前処理（具体的には、「前処理装置」をハードウェアで構成するか、あるいはソフトウェアで構成する。以下同じ。）を実施すれば、時間分解能・空間解像度を損なわずにノイズを抑制した画像Ｑ′′が得られる上、当該画像Ｑ′′を高圧縮率で圧縮することが可能であり、通信・記録に必要な時間・メモリが大幅に節約できる。

（ Ｘ−４−２ ２次元又は３次元領域抽出処理の前処理としての実施例 ）
航空写真などの２次元画像や、上記ＭＲＩ装置等による３次元分布画像から、特定の領域を抽出する処理がしばしば行われる。例えば、前者の航空写真では道路のみの領域を抽出したり、後者のＭＲＩ装置等による人体頭部の３次元分布画像では血管のみの領域を抽出し、これをレンダリングしてコンピュータグラフィックスとして表示するような処理である。特に、後者の血管に関する影像では、これを回転させたり、様々な方向から観察すること等が可能となり、脳血管障害の診断等において重要な技術として認識されている。

このような領域抽出を行う場合、最も簡単かつ基本的に行われている処理は、いわゆる「閾値処理」である。すなわち、閾値処理とは、上記２次元画像や３次元分布画像等を構成する各画素ｘのスカラー値が、ある一定の閾値以上であるかどうかによって、各画素ｘを分類する処理のことをいう。そして、特定の分類に入る画素ｘのみからなる画素の集合をレンダリングすれば、上記領域抽出を経た画像、すなわち例えば「血管のみの立体形状像」を表すことができる（ただし、従来技術としても、「閾値処理」以外にも、様々な方式が使われている。）。

ところで、このような領域抽出処理においては、元となる画像（上記にいう２次元画像や３次元分布画像等）にノイズが少ないことが望ましい。なぜなら、画像にノイズがある場合、例えば、本来血管ではない場所にある画素が、あたかも血管であるかのような画素の値を持ち、逆に、本来血管である場所に存在する画素が、血管ではないかのうような画素の値を持つ、という状態が生じるからである。従ってこの場合には、複雑な後処理によって、そのような誤って分類された画素を取り除く工程が必要になる。また、現状においては、上記元となる画像上のノイズを除くために、従来の技術の項で述べたような「平滑化処理」を行うことが一般的に行われているが、この場合、空間解像度が損なわれ、抽出したい領域の細かい構造（例えば、細い血管など）が失われてしまうという問題がある。

そこで、上記元となる画像、すなわち２次減画像や３次元分布画像に対し、本発明に係るコヒーレント・フィルタをかけてそのノイズを抑制する前処理を実施すれば、空間解像度を損なわずに画像のノイズを抑制することができるから、結果、正確な領域抽出処理を実施することができる。

なお、これに類似する実施例については、既に、上記「１枚の画像から画素値ｖ（ｘ）を構成する」場合において述べた（図８（ｃ）を参照）。

図１は、危険率ｐ（ｘ，ｙ）、重みｗ（ｐ（ｘ，ｙ））及び帰無仮説の帰趨の関係を表した表図である。 図２は、本発明に係るコヒーレント・フィルタの本実施形態における各種適用例を、画素値ｖ（ｘ）の構成法の差異に応じて分類して示す表図である。 図３は、コヒーレント・フィルタを組み込んだＸ線ＣＴ装置の概略構成例を示す図である。 図４は、本発明に係るコヒーレント・フィルタによる画像処理を概念的に説明する説明図である。 図５（ａ）は、ダイナミックＣＴ像を構成する複数の画像につき、本発明に係るコヒーレント・フィルタによるノイズ抑制処理を実施した例を示し、特に未処理のオリジナル画像を示すものである。 図５（ｂ）は、ダイナミックＣＴ像を構成する複数の画像につき、本発明に係るコヒーレント・フィルタによるノイズ抑制処理を実施した例を示し、特に処理済の画像を示すものである。 図５（ｃ）は、ダイナミックＣＴ像を構成する複数の画像につき、本発明に係るコヒーレント・フィルタによるノイズ抑制処理を実施した例を示し、特に処理済の画像を示すものである。 図６は、本実施形態におけるコヒーレント・フィルタによるノイズ抑制処理の流れを示すフローチャートである。 図７（ａ）は、１枚のＣＴ画像につき、本発明に係るコヒーレント・フィルタによるノイズ抑制処理を実施した例を示し、特に未処理のオリジナル画像をしている。 図７（ｂ）は、１枚のＣＴ画像につき、本発明に係るコヒーレント・フィルタによるノイズ抑制処理を実施した例を示し、特に処理済の画像を示している。 図７（ｃ）は、１枚のＣＴ画像につき、本発明に係るコヒーレント・フィルタによるノイズ抑制処理を実施した例を示し、図５（ａ）に示す画像と図５（ｂ）に示す画像との平均をとった画像を示している。 図８（ａ）は、血管の３Ｄ抽出例を示す造影された画像である。 図８（ｂ）は、血管の３Ｄ抽出例を示す造影された画像である。 図８（ｃ）は、血管の３Ｄ抽出例を示す造影された画像である。 図８（ｄ）は、血管の３Ｄ抽出例を示す造影された画像である。 図８（ｅ）は他の抽出例として脳の灰白質を抽出した画像である。 図９は、１枚の画像から画素値ｖ（ｘ）を構成する方法の一例を示す説明図である。 図１０は、１枚の画像から画素値ｖ（ｘ）を構成する方法の他の例を示す説明図である。 図１１は、コヒーレントレグレション法とＡＵＣ法を組み合わせたときの時間濃度曲線と、従来法とＡＵＣ法とを組み合わせたときの同曲線とを示す図である。 図１２は、明細書本文中式（４）等における分散σ２の別の推定法に関し、当該別の推定法による分散σ２＋＋の算出方法を説明する説明図である。 図１３は、従来のノイズ抑制方法である「平滑化処理」を説明する説明図である。

符号の説明

１０１…Ｘ線管、１０２…Ｘ線検出器、１０３…データ収集部、１０４…前処理部、１０５…メモリ部、１０６…再構成部、１０７…画像表示部、１０８…制御部、１０９…入力部。

Claims (5)

1. １つの画像を構成する第１の画素各々と前記画像内であって前記第１の画素各々の近傍の第２の画素との類似度に応じた重みを、前記画像を含む時間方向に異なる複数枚の画像にわたる前記第１の画素に対応する複数のスカラー値で構成されるベクトル値と前記複数枚の画像にわたる前記第２の画素に対応する複数のスカラー値で構成されるベクトル値との比較結果に基づいて計算する計算手段と、
   前記計算された重みによる前記第１の画素各々と前記第２の画素との重み付き平均を前記第１の画素各々の画素値として計算する平均手段と、
   を具備したことを特徴とする画像処理装置。
2. １つの画像を構成する第１の画素各々と前記画像内であって前記第１の画素各々の近傍の第２の画素の集合との重み付け平均を、前記第１の画素各々の画素値として計算する平均手段と、
   前記重み付け平均の処理に用いる重みを求める手段とを備え、
   前記重みは、前記第１の画素各々の座標における異なる時刻の画像の該当する画素値と、前記第２の画素各々の座標における異なる時刻の画像の該当する画素値とに基づいて求められるものであることを特徴とする画像処理装置。
3. 前記重みは、画像のノイズの分散値に基づいて求められるものである請求項１または２記載の画像処理装置。
4. 前記画像は、ダイナミックＣＴ画像である請求項１乃至３記載の画像処理装置。
5. 前記重みは、第１の画素の座標における異なる時刻の画像の該当する画素値と、第２の画素の座標における異なる時刻の画像の該当する画素値とを用いた検定の結果に基づくものであることを特徴とする請求項２記載の画像処理装置。

Images