JP2017027357A - 画像処理装置、撮像システム、画像処理方法 - Google Patents

Abstract

【課題】被検体画像（特にＳＮ比が低い画像）に対し非等方フィルタリングを適用する場合に、フィルタ特性を適切に決定し、信頼性の高いフィルタ結果を得るための技術を提供する。【解決手段】画像処理装置は、被検体を透過した電磁波により形成される像を撮像装置によって撮影することによって得られたデータを、前記被検体の測定データとして取得する測定データ取得手段と、前記測定データから、互いに異なる物理量の情報を表す第１画像と第２画像を生成する画像生成手段と、前記第１画像に基づいて、非等方フィルタを適用するときに設定すべきフィルタ特性を画像内の位置ごとに決定するフィルタ特性決定手段と、前記第１画像に基づいて決定されたフィルタ特性にしたがって、画像内の位置ごとに前記非等方フィルタのフィルタ特性を変えながら、前記第２画像に対し前記非等方フィルタを適用するフィルタリング手段と、を有する。【選択図】図４

Description

本発明は、被検体画像の画像処理技術に関する。

従来、精密に物質を測定する手段として位相を用いた測定法が用いられてきた。位相を用いた測定法は、波面のそろった（コヒーレントな）入射光に対して干渉を起こし、その干渉縞を計測することで、波長の数分の一から数十分の一の位相差による入射光波面（位相）の変化を測定する。このような測定法を利用した干渉計は、例えばレンズの表面のわずかな凹凸を測定するのに好適な手段である。

さらに、干渉を用いた波面計測手法の中でも近年注目を集めているのが、数十ナノメール以下の波長の光（電磁波）、すなわちＸ線を利用したＸ線位相イメージングである。Ｘ線位相イメージングは、従来の被検体の吸収によるコントラストを画像化するＸ線吸収像とは異なり、被検体に対するＸ線の透過時に生じる入射光の光路長の変化を位相の干渉によって検出する手法である。

Ｘ線位相イメージングの例として、Ｘ線を用いたトールボット干渉法について述べる。Ｘ線トールボット干渉計では、光源からのＸ線が被検体を透過し、それに伴って光の入射位相が変化する。被検体を透過した光は、回折格子と呼ばれる周期的パターンを持った格子で回折されることによって、回折格子からトールボット距離と呼ばれる所定の距離だけ離れた位置に第一の干渉パターンを形成する。この第一の干渉パターンの変化を被検体のない場合と比較解析することで前述の入射光波面の変化を測定する。

前述の周期的なパターンを持った回折格子のパターン周期は、装置長や入射光の波長等の条件によって変化する。一般的なＸ線の場合、その周期は数μｍのオーダーとなる。また、それによって生じる第一の干渉パターンも同様に数μｍのオーダーの周期となることが知られている。このような場合、一般的に用いられている検出器では、分解能がせいぜい数十μｍであるため、第一の干渉パターンを検出することは不可能である。そのため、干渉パターンが形成される位置に第一の干渉パターンとほぼ同じ周期の遮蔽格子を配置する。遮蔽格子で第一の干渉パターンの一部を遮ることにより周期が数百μｍ程度の第二の干渉パターン、すなわちモアレを形成し、このモアレを検出器で検出することによって干渉パターンの変化を間接的に測定することができる。モアレの形成法には周期を調節した遮蔽格子を第一の格子による第一の干渉パターンと同じ向きにアライメントして出す方法（拡大モアレ）と、格子を回転させてモアレの周期と向きを調節する方法（回転モアレ）がある。

このモアレを数値的に解析することによって、被検体に関するいくつかのパラメータを取得できる。代表的なものが前述した入射光の光路長の違いを微分位相として検出するものである。これによって、被検体の屈折率を計算することができる。このような屈折率による画像は被検体の種類、特に軟組織やプラスチックの様な被検体に対してはＸ線吸収像よりも信号雑音比（ＳＮ比）が高くなる利点がある。

近年ではモアレの振幅の変化を計算することで被検体の画素オーダー以下の小角散乱の情報を取得する手法に関しても研究がおこなわれている。つまり二次元Ｘ線トールボット干渉計においては、従来の吸収像に加えて、微分位相像、振幅像、散乱像など、互いに独立した物理量を表す画像を一度の撮像プロセスで取得可能である。

さらに近年、二次元格子を用いた二次元Ｘ線トールボット干渉計の研究も進んでいる。二次元Ｘ線トールボット干渉計は、例えばＸ軸方向の微分位相像とＸ軸に垂直なＹ軸方向の微分位相像を同時に取得できる特徴があるため、より正確な屈折率情報を取得することができる。このような二方向の微分位相像から、被検体に対する明瞭な像を取得する方法が、特許文献１に開示されている。特許文献１の方法は、逆リース変換を利用して被検体のエッジ境界部を明確化することが可能であり、被検体の輪郭に関するＳＮ比のさらなる向上が期待できる。その他にも微分位相からさらに微分してラプラシアンを求めてもよい。

このような技術はＸ線干渉計のようにＸ線を使用するモダリティ、特に被検体が生体のような被曝に弱いものである場合に大いに有効となる。このような被検体の場合、Ｘ線照射量が低ければ低いほど望ましい。当然照射量を減らすことはＳＮ比の悪化とトレードオフの関係になっており、最終的な画像に対してはよりエッジ部を強調する手法が望まれる。

効果的にエッジ部を強調する手法としては、例えば特許文献２の手法がある。特許文献２の手法はノイズリダクションの改良手法で、非等方なノイズリダクションフィルタリングを被検体の場所・形状に合わせて適用するプロセスを繰り返すことで効果的に強調したい情報部分のみを強調することが可能となる。このような手法はＸ線トールボット干渉計にも応用可能であるが、特許文献２の手法は一次情報として得られた像に一定のＳＮ比がすでに得られていることが前提となり、そこから非等方フィルタの方向を検出していく必要がある。一方で、低被曝化のために照射量を可能な限り低下させると、特許文献２の手法の適用が困難な程度（非等方フィルタの方向の検出ができない程度）までＳＮ比が悪化する場合もある。

このような例の一つが生体中の軟骨の撮像である。軟骨の撮像は従来の吸収Ｘ線では撮像できず、Ｘ線トールボット干渉計の応用用途として有望視されている撮像対象の一つである。この軟骨そのものは真空中の屈折率が高いが、これが生体組織中に存在すると他の生体組織との相対屈折率が低下し、撮像が不可能ではないが十分なコントラストが得られない場合がある。

特開２０１４−１０８３５８号公報 特表２０１２−５０４００３号公報

本発明は上記実情に鑑みなされたものであり、被検体画像（特にＳＮ比が低い画像）に対し非等方フィルタリングを適用する場合に、信頼性の高いフィルタ結果を得るための技術を提供することを目的とする。

本発明の第一態様は、被検体を透過した電磁波により形成される像を撮像装置によって撮影することによって得られたデータを、前記被検体の測定データとして取得する測定データ取得手段と、前記測定データから、互いに異なる物理量の情報を表す第１画像と第２画像を生成する画像生成手段と、前記第１画像に基づいて、非等方フィルタを適用するときに設定すべきフィルタ特性を画像内の位置ごとに決定するフィルタ特性決定手段と、前記第１画像に基づいて決定されたフィルタ特性にしたがって、画像内の位置ごとに前記非等方フィルタのフィルタ特性を変えながら、前記第２画像に対し前記非等方フィルタを適
用するフィルタリング手段と、を有することを特徴とする画像処理装置を提供する。

本発明の第二態様は、被検体を透過した電磁波により形成される像を撮影し、前記被検体の測定データを出力する撮像装置と、前記撮像装置により得られた前記測定データに対し画像処理を適用する、本発明の第一態様に係る画像処理装置と、を有することを特徴とする撮像システムを提供する。

本発明の第三態様は、コンピュータが、被検体を透過した電磁波により形成される像を撮像装置によって撮影することによって得られたデータを、前記被検体の測定データとして取得するステップと、コンピュータが、前記測定データから、互いに異なる物理量の情報を表す第１画像と第２画像を生成するステップと、コンピュータが、前記第１画像に基づいて、非等方フィルタを適用するときに設定すべきフィルタ特性を画像内の位置ごとに決定するステップと、コンピュータが、前記第１画像に基づいて決定されたフィルタ特性にしたがって、画像内の位置ごとに前記非等方フィルタのフィルタ特性を変えながら、前記第２画像に対し前記非等方フィルタを適用するステップと、を含むことを特徴とする画像処理方法を提供する。

本発明の第四態様は、本発明の第三態様に係る画像処理方法の各ステップをコンピュータに実行させることを特徴とするプログラムを提供する。

本発明によれば、被検体画像（特にＳＮ比が低い画像）に対し非等方フィルタリングを適用する場合に、信頼性の高いフィルタ結果を得ることができる。

本発明の実施形態の撮像システムの構成を模式的に示す図。 被検体の構成と、微分位相像および逆リース画像の例を示す図。 比較例によるフィルタリング結果を示す図。 実施例のフィルタリング処理に関する機能ブロック図。 実施例のフィルタリング処理の流れを示すフローチャート。 実施例によるフィルタリング結果を示す図。 実施例の結果と比較例の結果と理論値を比較する図。 実施例の結果と比較例の結果と理論値を比較する図。

本発明は、被検体画像に対し非等方フィルタを用いたフィルタリングを行う画像処理装置及び画像処理方法に関する。被検体の撮像装置・撮像法は問わないが、本発明は、一回の撮像プロセスによって得られる測定データを基に、互いに異なる（独立した）物理量の情報を表す複数種類の被検体画像を取得できる撮像装置に対し好ましく適用することができる。この種の撮像装置としては、被検体を透過した電磁波により形成される像を撮影する撮像装置、例えば、電磁波が被検体を透過するときに生じる波面（位相）の変化を強度（コントラスト）として画像化する位相イメージング装置などがある。位相イメージング装置の測定データには、微分位相、吸収率、散乱（振幅変化）などの複数の物理量の情報が含まれている。撮像に用いる電磁波としては、光、Ｘ線、ガンマ線などいかなる波長帯域のものでもよい。ただし、詳しくは後述するが、本実施形態のフィルタリング手法はＳＮ比の低い画像に対して適用可能なため、Ｘ線、ガンマ線、中性子線、電子線などの放射線を利用する撮像装置における低線量撮像に対し特に好ましく適用できる。なお本明細書においてＸ線とはエネルギが２ｋｅＶ以上２００ｋｅＶ以下の電磁波を指す。

以下、本発明の好ましい実施の形態として、二次元微分干渉計の一つである二次元トー
ルボットＸ線位相イメージング装置で得られる被検体画像を例に挙げて説明する。なお、以下に述べるフィルタリング手法は二次元だけでなく一次元の微分干渉計で得られた画像にも応用可能である。

（撮像システム）
図１はトールボットＸ線位相イメージングの撮像システムの構成を示した図である。撮像システムは、Ｘ線トールボット干渉計からなる撮像装置１０と、撮像装置１０で取得したＸ線画像を処理する画像処理装置１１とから構成される。撮像装置１０は、Ｘ線を発生させるＸ線源（光源）１１０と、Ｘ線を回折する回折格子１３０と、Ｘ線の一部を遮る遮蔽格子１４０と、Ｘ線を検出する検出器１５０とを備えている。画像処理装置１１は、撮像装置１０に接続された演算部（コンピュータ）１６０と、演算部１６０による演算結果に基づいた画像を表示する画像表示装置１７０とを有している。

演算部１６０は、ＣＰＵ（中央演算処理装置）、ＲＡＭ（ランダムアクセスメモリ）、補助記憶装置などのハードウェア資源を備えた汎用のコンピュータにより構成できる。後述する画像処理、各種演算、および制御は、補助記憶装置に格納されたプログラムをＣＰＵが読み込み実行することで実現されるものである。なお、演算部１６０の機能のうちの一部又は全部をＡＳＩＣ（Ａｐｐｌｉｃａｔｉｏｎ Ｓｐｅｃｉｆｉｃ Ｉｎｔｅｇｒａｔｅｄ Ｃｉｒｃｕｉｔ）のような回路で構成することもできる。

Ｘ線源１１０からのＸ線は回折格子１３０により回折され、トールボット距離と呼ばれる所定の距離をおいて明部と暗部が配列方向に並んだ干渉パターン１８０を形成する。

通常、回折格子１３０による第一の干渉パターン１８０の周期は数μｍから十数μｍ程度である。そこで、第一の干渉パターン１８０と同じかわずかに異なる周期をもつ遮蔽格子１４０を、第一の干渉パターン１８０が形成される位置に配置する。そうすると、第一の干渉パターン１８０と遮蔽格子１４０によりモアレが形成され、干渉パターンの周期を数十μｍ以上あるいは無限に拡大することができる。モアレの周期は、用いる位相回復方法と検出器の画素サイズを考慮して適宜決めることができるが、本実施形態においては、モアレの周期が画素サイズの２倍以上であることが好ましい。尚本実施形態においては、周期が大きく、検出器の撮像範囲内に明暗が形成されない場合であっても、第一の干渉パターン１８０と遮蔽格子１４０により形成される強度分布のことをモアレと呼ぶ。このモアレ（空間的な周期性をもつパターン）を二次元イメージセンサである検出器１５０により画像化し、二次元画像を得る。このような仕組みにより、数十μｍ平方程度の分解能の検出器１５０で、数μｍから十数μｍの周期をもつ干渉パターンのイメージングを可能にしている。ただし、検出器１５０の空間分解能が十分に高い場合には、遮蔽格子１４０を省略し、第一の干渉パターン１８０をそのまま撮像してもよい。以下、検出器１５０で得られた周期的なパターンをもつ画像を干渉像と呼び、検出器１５０から出力される干渉像の画像データを測定データと呼ぶ。

回折格子１３０には、周期的にＸ線の位相を変調する位相型の回折格子（位相格子）でも、周期的にＸ線の振幅を変調する振幅型の回折格子（遮蔽格子）でもよいが、Ｘ線の損失が少ないため位相格子が用いられることが多い。遮蔽格子１４０には、Ｘ線透過部とＸ線遮蔽部とが配列された格子が用いられることが多い。

測定時には、被検体１２０を回折格子１３０の前に設置する。Ｘ線は一般に透過性が高いために被検体１２０を透過するが、その際に透過した物質の元素組成と密度に応じた位相の変化が生じる。この位相の変化は第一の干渉パターン１８０の配置に影響を与える。そのため遮蔽格子１４０によるモアレにも歪みを生じさせる。したがって、演算部１６０で、被検体１２０が無いときの測定データと被検体１２０を設置したときの測定データの
間での干渉像の歪みを計算することで、被検体１２０によるＸ線の位相変化の情報（微分位相）を取得できる。また、Ｘ線が被検体１２０を透過する際にＸ線のエネルギの一部は吸収されるが、その吸収率は組織の元素組成によって違いがある。したがって、演算部１６０で、被検体１２０が無いときの測定データと被検体１２０を設置したときの測定データの間での強度比を計算することで、被検体１２０によるＸ線吸収の情報を取得できる。また、被検体１２０の内部に微小な構造が密集した部分（例えば線維組織や微小石灰化）が存在すると、その部分において小角散乱が生じ、測定されるＸ線の振幅が変化する。したがって、演算部１６０で、被検体１２０が無いときの測定データと被検体１２０を設置したときの測定データの間での振幅の変化を解析することで、被検体１２０によるＸ線の小角散乱の情報を取得できる。結果として、検出器１５０で得られる測定データを基に、微分位相、吸収率、小角散乱の３つの独立した物理量の情報を取得できる。各物理量の情報の二次元分布を画像化したものを、それぞれ、微分位相像、吸収像、散乱像（又はビジビリティコントラスト像）と呼ぶ。以上がＸ線位相イメージング装置の概要である。尚、微分位相の情報を積分すると位相の情報を取得できるが、微分位相と位相は独立した物理量ではない。

（比較例）
まず、比較例として、特許文献１に開示されている逆リース変換を利用した画像処理方法（以後、従来方法１と呼ぶ）と特許文献２に開示されている非等方フィルタリング方法（以後、従来方法２と呼ぶ）を単純に組み合わせた例について説明する。

本明細書では、被検体として現実の生体組織を模したファントムを想定し、コンピュータ・シミュレーションによって得られた干渉像（モアレ像）を基に、以下の画像処理を行う。図２Ａはファントムの構造である。このファントムは球形の二重の構造物を有している形になっている。第一の層２０１は吸収率の高い構造物である。これは例えば生体組織でいえば骨にあたる部分を示す。第二の層２０２は吸収率が低く屈折率も低い構造物である。これは生体組織でいえば前述の生体組織中の軟骨を模擬している。物理定数もそれらの状況に準じた吸収率・屈折率に設定し、大きさを２ｃｍとした。また、皮膚や筋肉組織といった他の生体組織中の軟骨像を再現するため、このファントムは水中に設置されていると仮定した。これによって第二の層２０２との屈折差が小さくなり、位相像のコントラストが弱くなる。

図２Ａのファントムを図１の撮像装置１０で撮像した場合に観測される干渉像をシミュレーションにより求め、その干渉像を基に生成した微分位相像の例が図２Ｂおよび図２Ｃである。図２ＢがＸ軸方向の微分位相像であり、図２ＣがＹ軸方向の微分位相像である。またここで被曝量は５ｍＧｙを想定した。これは医療用のＸ線撮像としては標準的かやや多めの照射線量である。中央の骨に対応する箇所のコントラストは付いているが、軟骨の端部に関するコントラストはほぼ不可視になっている。これらの微分位相像に対し従来方法１を適用することで生成した画像（逆リース画像と呼ぶ）が図２Ｄである。図２Ｂおよび図２Ｃの微分位相像に比べて幾分軟骨端部が可視化されているが、十分とは言い難い。例えば医療応用では、軟骨の存在が確認できるだけでは十分でなく、軟骨の厚みが重要な診断基準となることがある。図２Ｄは軟骨の存在を確認できるものの、正確な厚み測長を行うには困難が伴う。このように、被検体と撮像条件の組み合わせによっては、従来方法１で得られる画像が診断に不十分なことがあると考えられる。

次に、図２Ｄの逆リース画像に対してフィルタによるノイズリダクションを適用した結果を紹介する。図３Ａは等方的なフィルタを図２Ｄの逆リース画像に適用した結果である。フィルタの大きさは２５×２５ピクセルを用いた。平滑化手法としてメディアンフィルタを用いた。その結果、軟骨部において若干コントラストの改善が見られるものの、コントラストが不均一であり、また軟骨のエッジを示す箇所の幅が広いため、厚み測長時の定
量性が低く、厚み測長が困難であることが分かる。

図３Ｂは、図２Ｄの逆リース画像に従来方法２を適用した結果を示す。ここでは８１×７ピクセルの大きさの非等方フィルタを使用した。非等方フィルタの方向は、二次元トールボットＸ線位相イメージングで取得される一次情報であるＸ軸方向およびＹ軸方向の微分位相から計算した。二つの微分位相を合わせることで位相勾配のベクトルが設定できるため、そのベクトル方向に非等方フィルタの短辺の向きが合うようにして先ほどと同じメディアンフィルタを適用した。その結果全体としてぼやけてしまい特に軟骨のエッジ部のコントラストが消滅した。

このように、従来方法１と従来方法２の単純な組み合わせをこのサンプルに適用すると、ＳＮ比が十分にとれない、あるいはＳＮ比が悪化することが分かる。

（実施例）
次に、本実施形態の画像処理装置１１によって実行される画像処理の一実施例について説明する。本実施例では、非等方フィルタを適用するときのフィルタ特性（方向、形状など）を、微分位相ではなく、同じ測定データから生成された吸収像に基づいて決定する。骨の部分に関しては吸収像のコントラストは微分位相像および逆リース画像に比べて高いため、吸収像における構造のエッジ方向を補助情報としてフィルタの方向を決定する。このようにフィルタ特性の決定に用いる画像（第１画像）とフィルタ処理を適用する画像（第２画像）とが、同じ測定データから生成された画像のうち異なる物理量の情報を表す画像である点に、本実施例の特徴の一つがある。ここで、第１画像と第２画像の組み合わせは任意であるが、第１画像が第２画像に比べてＳＮ比が高い（観察したい構造のコントラストが高い）画像であることが好ましい。フィルタを適用する第２の画像よりもＳＮ比が高い第１の画像を用いてフィルタ特性を決定することで、第２の画像を用いてフィルタ特性を決定するよりも適切にフィルタ特性を決めることができ、被検体の情報に近い値を得ることができる。

図４は画像処理装置１１によって提供されるフィルタリング処理に関する機能ブロック図である。画像処理装置１１は、その機能として、測定データ取得部５００、画像生成部５０１、フィルタ特性決定部５０２、フィルタリング部５０３を有する。測定データ取得部５００は、撮像装置１０によって得られた被検体１２０の測定データを取得する機能であり、画像生成部５０１は、測定データを解析し、吸収像、微分位相像、散乱像、逆リース画像などの画像を生成する機能である。また、フィルタ特性決定部５０２は、非等方フィルタのフィルタ特性を決定する機能であり、フィルタリング部５０３は、フィルタ特性決定部５０２で決定されたフィルタ特性にしたがって非等方フィルタを画像に適用する機能である。この非等方フィルタは、画像を一方向に沿って平滑化する、方向依存性をもつフィルタである。本実施例では、８１×７ピクセルの大きさをもち、その短辺方向に沿って平滑化効果（ノイズリダクション効果）を有する非等方フィルタを用いる。これらの機能は、演算部１６０が必要なプログラムを実行することで実現されるものである。

図５のフローチャートを参照して、画像処理装置１１によって提供されるフィルタリング処理の流れを説明する。まず測定データ取得部５００が、撮像装置１０によって得られた測定データを取得する（ステップＳ４００）。このとき、撮像装置１０からリアルタイムで測定データを取得してもよいし、過去に撮像されストレージ等に保存されている測定データを読み込んでもよい。次に画像生成部５０１が、測定データを解析し、吸収像、Ｘ軸方向の微分位相像、Ｙ軸方向の微分位相像を生成する（ステップＳ４０１）。次に画像生成部５０１は、逆リース変換により、微分位相像から逆リース画像を取得する（ステップＳ４０２）。吸収像、微分位相像、逆リース画像などが既に生成され、ストレージ等に保存されている場合には、ステップＳ４０１、Ｓ４０２ではその画像データを読み込んで
もよい。なお、逆リース変換については特許文献１に開示された方法を利用するため、ここでは詳しい説明は省略する。

次に、フィルタ特性決定部５０２は、吸収像を用いて、画像中の構造（例えば骨）のエッジ方向（境界の接線方向）を抽出し、そのエッジ方向に直交する法線方向を求める。そのためには、画像中の画素値の勾配の向きを画像内の位置ごとに計算すればよい。本実施例では、フィルタ特性決定部５０２がまず吸収像を二値化する（ステップＳ４０３）。二値化とは、画像内の各画素の値を１（Ｔｒｕｅ）か０（Ｆａｌｓｅ）に変換する操作である。本実施例では、吸収像の画素値が閾値より大きいか否かを判定することで二値化を行う。このときの閾値は固定値でもよいし、動的（適応的）に決めてもよい。さらに必要に応じて二値画像に対してメディアンフィルタを適用することが好ましい。これは骨の領域とそれ以外の領域をより明確に分けるための操作である。このようにして取得した二値画像が図６Ａである。閾値とメディアンフィルタのパラメータを調整することで骨の領域を抽出可能である（図６Ａの白色の領域が骨の領域である。）。

次に、フィルタ特性決定部５０２は、二値画像に対してコンボリューションによるぼけを与える（ステップＳ４０４）。本実施例では、二値画像に対しガウシアンフィルタをコンボリューションする。その結果得られた画像（勾配画像と呼ぶ）が図６Ｂである。ステップＳ４０４の処理は、画像中の構造のエッジ（例えば骨の境界）を稜線とする画素値の連続的な勾配を形成するためである。このような勾配画像を用いることで、画像内の各位置におけるフィルタ特性を適切に設定することが容易になる。このとき、画像の隅までぼけを与えることが好ましく、このぼけは等方的であることが好ましい。

次に、フィルタ特性決定部５０２は、勾配画像から画像内の各位置の勾配を計算する（ステップＳ４０５）。本実施例では、勾配の向きだけを取得するが、勾配の大きさをも取得してフィルタ特性の決定に利用してもよい。図６Ｂの勾配画像から計算された勾配の向きを示したのが図６Ｃである。矢印が勾配の向きを表している。中央から外向きに勾配ベクトルが向いていることが分かる。これによって画像全領域に対して骨の境界を判別し、その境界の法線方向を得ることができたことになる。今回の被検体は球状になっているが、本手法を用いれば複雑な境界線を持つ骨の像に関してもその境界の法線方向を求めることが可能となる。フィルタ特性決定部５０２は、勾配の向きと非等方フィルタの方向（短辺方向）とが平行になるように、画像内の位置ごとのフィルタ特性を決定する（ステップＳ４０６）。

次に、フィルタリング部５０３は、ステップＳ４０６で決定されたフィルタ特性にしたがって、画像内の位置ごとにフィルタ特性（本実施例ではフィルタの方向）を変えながら、逆リース画像の全ての画素に対し非等方フィルタを適用する（ステップＳ４０７）。これにより、勾配の向きに沿ったノイズリダクションが実施される。

図６Ｄは、図２Ｄの逆リース画像に対して図６Ｃに示す方向に沿ったノイズリダクション（ここではメディアンフィルタを用いた）を実施した結果である。図３Ａ、図３Ｂの比較例と比べると、図６Ｄの画像では軟骨の形状がはっきりしていることが分かる。また、軟骨のエッジを示す個所の幅も狭いため、骨表面（境界）から軟骨のエッジ位置がどれだけ離れているか、すなわち軟骨の厚さ、を測長することに適していることが分かる。

図７Ａ〜図７Ｄは本実施例の効果を比較例と定量的に比較するための図である。図７Ａ〜図７Ｄは、軟骨の近傍領域をドーナツ形状に切り出し、それを長方形状に写像した画像である。図７Ａが理論値、図７Ｂ〜図７Ｄはそれぞれ図３Ａ（比較例１）、図３Ｂ（比較例２）、図６Ｄ（実施例）から切り出した画像である。本実施例で得られた画像（図７Ｄ）では、理論値（図７Ａ）とほぼ同じ位置に、軟骨のエッジ部が明りょうに現れているこ
とが分かる。

また、図８は、図７Ａ〜図７Ｄの画像の横方向の平均値をプロットした図である。縦方向の位置は、図７Ａ〜図７Ｄの紙面における上端を０とした。比較例１および比較例２のプロファイルは理論値のものとずれており、プロファイルの極小値の位置（軟骨のエッジと推定される位置）も理論値から大きくずれている。これに対し、本実施例では、理論値に近いプロファイルが取得できており、比較例よりも定量性の高い画像処理が可能であることが分かる。

（応用例）
本発明に係るフィルタリング手法は様々な形態で応用が可能である。例えば上記実施例は吸収像を参照して微分位相から取得される逆リース画像のノイズリダクションを行った。しかし、非等方フィルタを適用する画像は、微分位相像、その積分像である位相像、微分位相像をさらに微分したラプラシアン像でもよいし、それ以外の画像でもよい。また、吸収コントラストよりも位相コントラストの方が強い場合には、本実施例とは逆に、微分位相像やそれを加工した画像（逆リース画像、位相像など）からフィルタ特性を決定し、吸収像に対し非等方フィルタを適用してもよい。また、これらの技術を位相イメージングで取得できるもう一つの画像である散乱像（モアレ振幅像）のノイズリダクションに用いることも可能である。Ｘ線位相イメージングでは、一次像として吸収像、微分位相像、散乱像が取得できるため、こられの画像、もしくは画像から作成される別の画像の一つを画像処理するために、他の作成された画像を利用する様々な様態が考えられる。また本手法は二次元Ｘ線位相イメージングに限らず、一次元Ｘ線位相イメージングにも応用可能である。また、トールボット干渉計以外の撮像装置にも本手法を応用可能である。上記の実施例では、第１画像の勾配情報を取得するために、二値化及びコンボリューションによるぼけの付与を行っているが、勾配情報の取得方法はこれに限らず、既知の方法を用いることもできる。第１画像で取得した勾配情報を用いて第２画像に適用する非等方フィルタを決定すれば、上記実施例と同様の結果が得られると考えられる。

また、上記実施例では非等方フィルタは形状が８１×７ピクセルで固定されていたが、フィルタ形状はこれに限らないし、また、条件もしくは画像内の場所によって非等方フィルタの形状を変えたり、正方のフィルタを使用したりしてもよい。例えば、本実施例ではフィルタ特性としてフィルタの方向のみを変更したが、勾配の大きさに応じて非等方フィルタの形状、サイズ、係数などを変えてもよい。また、フィルタ手法もメディアンフィルタではなく、普通のガウシアンによる平均化フィルタやその他様々なフィルタを用いることができる。

上述した画像処理装置の具体的な実装は、ソフトウェア（プログラム）による実装と、ハードウェアによる実装のいずれも可能である。例えば、画像処理装置に内蔵されたコンピュータ（マイコン、ＣＰＵ、ＭＰＵ、ＦＰＧＡ等）のメモリにコンピュータプログラムを格納し、当該コンピュータプログラムをコンピュータに実行させて、各処理を実現させてもよい。また、本発明の全部または一部の処理を論理回路により実現するＡＳＩＣ等の専用プロセッサを設けることも好ましい。また、本発明は、クラウド環境におけるサーバーにも適用可能である。

また、例えば、記憶装置に記録されたプログラムを読み込み実行することで前述した実施形態の機能を実現するシステムや装置のコンピュータによって実行されるステップからなる方法によっても、本発明を実施することができる。この目的のために、上記プログラムは、例えば、ネットワークを通じて、又は、上記記憶装置となり得る様々なタイプの記録媒体（つまり、非一時的にデータを保持するコンピュータ読取可能な記録媒体）から、上記コンピュータに提供される。よって、上記コンピュータ（ＣＰＵ、ＭＰＵ等のデバイ
スを含む）、上記方法、上記プログラム（プログラムコード、プログラムプロダクトを含む）、上記プログラムを非一時的に保持するコンピュータ読取可能な記録媒体は、いずれも本発明の範疇に含まれる。

１０：撮像装置、１１：画像処理装置、１６０：演算部、５００：測定データ取得部、５０１：画像生成部、５０２：フィルタ特性決定部、５０３：フィルタリング部

Claims (13)

1. 被検体を透過した電磁波により形成される像を撮像装置によって撮影することによって得られたデータを、前記被検体の測定データとして取得する測定データ取得手段と、
   前記測定データから、互いに異なる物理量の情報を表す第１画像と第２画像を生成する画像生成手段と、
   前記第１画像に基づいて、非等方フィルタを適用するときに設定すべきフィルタ特性を画像内の位置ごとに決定するフィルタ特性決定手段と、
   前記第１画像に基づいて決定されたフィルタ特性にしたがって、画像内の位置ごとに前記非等方フィルタのフィルタ特性を変えながら、前記第２画像に対し前記非等方フィルタを適用するフィルタリング手段と、
   を有することを特徴とする画像処理装置。
2. 前記フィルタ特性は、前記非等方フィルタの方向を含む
   ことを特徴とする請求項１に記載の画像処理装置。
3. 前記非等方フィルタは、前記非等方フィルタの方向に沿って画像を平滑化するフィルタである
   ことを特徴とする請求項２に記載の画像処理装置。
4. 前記フィルタ特性決定手段は、前記第１画像に含まれる構造のエッジ方向と前記非等方フィルタの方向とが直交するように、画像内の位置ごとの前記非等方フィルタのフィルタ特性を決定する
   ことを特徴とする請求項２又は３に記載の画像処理装置。
5. 前記フィルタ特性決定手段は、前記第１画像の画素値の勾配の向きを画像内の位置ごとに計算し、前記勾配の向きと前記非等方フィルタの方向とが平行になるように、画像内の位置ごとの前記非等方フィルタのフィルタ特性を決定する
   ことを特徴とする請求項２〜４のうちいずれか１項に記載の画像処理装置。
6. 前記第１画像は、前記第２画像に比べてＳＮ比が高い画像である
   ことを特徴とする請求項１〜５のうちいずれか１項に記載の画像処理装置。
7. 前記撮像装置は、位相イメージング装置である
   ことを特徴とする請求項１〜６のうちいずれか１項に記載の画像処理装置。
8. 前記電磁波は、放射線である
   ことを特徴とする請求項７に記載の画像処理装置。
9. 前記第１画像は、前記被検体による前記電磁波の吸収の情報を表す画像である
   ことを特徴とする請求項７又は８に記載の画像処理装置。
10. 前記第２画像は、前記被検体による前記電磁波の位相変化の情報を表す画像である
    ことを特徴とする請求項７〜９のうちいずれか１項に記載の画像処理装置。
11. 被検体を透過した電磁波により形成される像を撮影し、前記被検体の測定データを出力する撮像装置と、
    前記撮像装置により得られた前記測定データに対し画像処理を適用する、請求項１〜１０のうちいずれか１項に記載の画像処理装置と、
    を有することを特徴とする撮像システム。
12. コンピュータが、被検体を透過した電磁波により形成される像を撮像装置によって撮影することによって得られたデータを、前記被検体の測定データとして取得するステップと、
    コンピュータが、前記測定データから、互いに異なる物理量の情報を表す第１画像と第２画像を生成するステップと、
    コンピュータが、前記第１画像に基づいて、非等方フィルタを適用するときに設定すべきフィルタ特性を画像内の位置ごとに決定するステップと、
    コンピュータが、前記第１画像に基づいて決定されたフィルタ特性にしたがって、画像内の位置ごとに前記非等方フィルタのフィルタ特性を変えながら、前記第２画像に対し前記非等方フィルタを適用するステップと、
    を含むことを特徴とする画像処理方法。
13. 請求項１２に記載の画像処理方法の各ステップをコンピュータに実行させることを特徴とするプログラム。

Images