JP6743983B2 - Ｘ線位相差撮像システム - Google Patents

Description

本発明は、Ｘ線位相差撮像システムに関し、特に、複数の格子を用いて撮像を行うＸ線位相差撮像システムに関する。

従来、複数の格子を用いて撮像を行うＸ線位相差撮像システムが知られている。このようなＸ線位相差撮像システムは、たとえば、国際公開第２０１４／０３０１１５号に開示されている。

国際公開第２０１４／０３０１１５号に開示されているＸ線位相差撮像システムは、タルボ・ロー干渉計によってＸ線撮像を行い、縞走査法によって位相コントラスト画像を生成するように構成されている。

ここで、タルボ・ロー干渉計では、マルチスリットと、位相格子と、吸収格子とを用いて撮像が行われる。具体的には、複数の格子のうちいずれか１つを、格子のパターンと直交する方向に並進させながら複数回撮像する。また、縞走査法とは、格子を並進移動させながら複数回撮像されたＸ線画像の各画素の画素値の強度変化に基づいて、位相コントラスト画像を生成する手法である。具体的には、縞走査法は、Ｘ線画像の各画素の画素値の強度変化が、格子の周期の関数の各位相におけるデータであるとして、関数の波形を決定し、波形が決定された関数に基づいて、位相コントラスト画像を生成する手法である。位相コントラスト画像には、吸収像、位相微分像および暗視野像が含まれる。吸収像とは、Ｘ線が被写体を通過した際に生じるＸ線の減衰に基づいて画像化した像である。また、位相微分像とは、Ｘ線が被写体を通過した際に発生するＸ線の位相のずれをもとに画像化した像である。また、暗視野像とは、物体の小角散乱に基づくＶｉｓｉｂｉｌｉｔｙの変化によって得られる、Ｖｉｓｉｂｉｌｉｔｙ像のことである。また、暗視野像は、小角散乱像とも呼ばれる。「Ｖｉｓｉｂｉｌｉｔｙ」とは、鮮明度のことである。

縞走査法を用いて位相コントラスト画像を生成する際、複数の格子のいずれか１つを、少なくとも格子の１周期分、格子のパターンと直交する方向に並進移動させながら撮像を行う。具体的には、縞走査法では、格子を格子のパターンと直交する方向に、格子の周期の１／Ｎずつ並進移動させながら撮像する。したがって、Ｎ回並進移動させた際の格子の総移動量が、移動させる格子の周期と等しくなる。また、格子の総移動量と移動させる格子の周期とが等しくなるため、各画素の画素値の強度変化の周期と移動させる格子の周期とが等しくなる。なお、Ｎは正の整数である。

国際公開第２０１４／０３０１１５号

しかしながら、国際公開第２０１４／０３０１１５号のＸ線位相差撮像システムで用いられている一般的な縞走査法では、Ｘ線画像の各画素の画素値の強度変化が１／Ｎ周期の位相の点のデータであることを前提として位相コントラスト画像を生成しているため、たとえば、熱変動による格子の相対位置の変化や格子移動機構の移動精度に起因して、各ステップにおける格子の移動量が１／Ｎ周期でない場合、関数の波形が実データから乖離する。したがって、格子の総移動量と格子の周期とにずれが生じた場合、得られる画像にアーチファクトが生じるという問題点がある。

この発明は、上記のような課題を解決するためになされたものであり、この発明の１つの目的は、格子の格子パターンと直交する方向における総移動量と格子の周期とにずれが生じている場合でも、得られる位相コントラスト画像にアーチファクトが生じることを抑制することが可能なＸ線位相差撮像システムを提供することである。なお、格子の格子パターンと直交する方向における総移動量と格子の周期とのずれとは、熱変動や格子移動機構の精度等に起因する、極めて小さいずれである。

上記目的を達成するために、この発明の一の局面におけるＸ線位相差撮像システムは、Ｘ線源と、Ｘ線源からＸ線が照射される第１格子と第１格子からのＸ線が照射される第２格子とを含む複数の格子と、Ｘ線源から照射されたＸ線を検出する検出器と、複数の格子の少なくともいずれか１つを移動させる格子移動機構と、検出器により検出された各画素の画素値の変化を表す強度変化に基づいて、位相コントラスト画像を生成する画像処理部とを備え、画像処理部は、強度変化の周期と、周期を変数とするとともに格子の移動に伴う画素値の強度変化を表す関数とにより、位相コントラスト画像を生成するように構成されている。

この発明の一の局面におけるＸ線位相差撮像システムは、上記のように、強度変化の周期と、周期を変数とするとともに格子の移動に伴う画素値の強度変化を表す関数とにより、位相コントラスト画像を生成する画像処理部を備える。これにより、関数の周期が変数であるため、関数の周期の最適解を求めることにより、格子の移動によって実際に取得された強度変化の周期を近似することができる。その結果、各画素における画素値の強度変化の周期と格子の周期とにずれが生じている場合でも、関数の周期を強度変化の周期に合わせることが可能となるので、得られる位相コントラスト画像にアーチファクトが生じることを抑制することができる。

上記一の局面におけるＸ線位相差撮像システムにおいて、好ましくは、強度変化の周期は、格子の周期のずれと格子移動機構による格子の移動量のずれとのうち少なくとも一方のずれを含んだ周期である。これにより、格子の周期および格子の移動量のどちらにずれが生じていた場合（または、両方にずれが生じた場合）でも、ずれを含んだ強度変化の周期に、関数の周期（変数）を合わせることにより、生成する位相コントラスト画像にアーチファクトが生じることを抑制することができる。

上記一の局面におけるＸ線位相差撮像システムにおいて、好ましくは、画像処理部は、複数の格子の少なくともいずれか１つを移動させて取得した強度変化のデータと関数とに基づいて、強度変化の周期を決定するように構成されている。このように構成すれば、実際に計測した強度変化のデータから強度変化の周期を取得することができる。その結果、強度変化の周期と格子の周期とにずれが生じている場合でも、ずれを含んだ強度変化の周期を取得することができる。

上記一の局面におけるＸ線位相差撮像システムにおいて、好ましくは、画像処理部は、強度変化の波形の形状と関数の波形の形状とをフィッティングすることにより強度変化の周期を決定するように構成されている。このように構成すれば、強度変化の波形の形状と関数の波形の形状との偏差に基づいて強度変化の周期を決定することができる。その結果、複数の周期を用いて取得した偏差を比較することにより強度変化の周期を決定することが可能となるので、たとえば、複数の周期を用いて取得した画像を画像処理することにより強度変化の周期を決定する場合と比較して、強度変化の周期を決定する処理を簡素化することができる。

この場合、好ましくは、画像処理部は、少なくとも格子移動機構により移動される格子の周期の設計値よりも大きい値および小さい値の両方の値を周期として用いて強度変化の波形の形状と関数の波形の形状との偏差を求めるように構成されている。このように構成すれば、格子の周期の設計値よりも大きい値により算出される偏差と格子の周期の設計値よりも小さい値により算出される偏差とを比較することにより、偏差を最小にするための周期が、格子の周期の設計値よりも大きい値か小さい値かを判断することができる。その結果、偏差が最小になる周期を容易に決定することができる。

上記強度変化の波形の形状と関数の波形の形状とをフィッティングすることにより強度変化の周期を決定する構成において、好ましくは、画像処理部は、所定の周期に対応する関数の係数を取得し、係数および所定の周期により強度変化の波形の形状と関数の波形の形状との偏差を取得するように構成されている。このように構成すれば、所定の周期に対応する偏差を取得することができる。その結果、複数の周期に対応した偏差を比較することにより、偏差が最小となる周期を容易に決定することができる。

上記一の局面におけるＸ線位相差撮像システムにおいて、好ましくは、画像処理部は、あらかじめ取得された強度変化の周期から、関数の係数を取得し、取得した係数およびあらかじめ取得された強度変化の周期に基づき、位相コントラスト画像を生成するように構成されている。ここで、格子の周期は非常に小さい値なので、製作誤差によりわずかに周期がずれる場合がある。格子の周期からわずかにずれた周期が既知であっても、一般的な縞走査法では、ｋ／Ｎ周期のデータを取得する必要がある。しかし、格子移動機構の精度が十分でないと、既知の周期をＮ等分した移動量で格子を並進移動させることができないので、ｋ／Ｎ周期のデータを取得することができず、アーチファクトが生じる。格子の周期のずれは計測することによりあらかじめ取得することが可能であるため、強度変化の周期と格子の周期とのずれは、あらかじめ取得することが可能である。したがって、上記のように構成すれば、実際の計測値と、あらかじめ取得された強度変化の周期を代入した関数とをフィッティングすることにより、関数の係数を取得することができる。その結果、格子の移動量がずれている場合でも、実データに精度よく近似した関数を決定することが可能となり、生成される位相コントラスト画像にアーチファクトが生じることを抑制することができる。なお、ｋは１からＮまでの正の整数である。

上記一の局面におけるＸ線位相差撮像システムにおいて、好ましくは、複数の格子は、Ｘ線源と第１格子との間に配置された第３格子をさらに含んでいる。このように構成すれば、第３格子によってＸ線源から照射されるＸ線の可干渉性を高めることができる。その結果、Ｘ線源の焦点径に依存することなく第１格子の自己像を形成させることが可能となるので、Ｘ線源の選択の自由度を向上させることができる。

本発明によれば、上記のように、格子の格子パターンと直交する方向における総移動量と格子の周期とにずれが生じている場合でも、得られる位相コントラスト画像にアーチファクトが生じることを抑制することが可能なＸ線位相差撮像システムを提供することができる。

本発明の第１実施形態によるＸ線位相差撮像システムをＹ方向から見た模式図である。 画素における画素値の強度変化を取得する処理を説明するための模式図である。 各画素における画素値の強度変化の周期のずれを説明するためのグラフ（Ａ）、および生成する位相コントラスト画像に生じるアーチファクトを説明するための模式図（Ｂ）である。 周期を変動させた場合のアーチファクトの変化を説明するための模式図（Ａ）〜（Ｃ）である。 本発明の第１実施形態による位相コントラスト画像の生成方法を説明するためのフローチャートである。 本発明の第２実施形態による位相コントラスト画像の生成方法を説明するためのフローチャートである。 本発明の第１実施形態の変形例によるＸ線位相差撮像システムをＹ方向から見た模式図である。

以下、本発明を具体化した実施形態を図面に基づいて説明する。

［第１実施形態］
図１〜図５を参照して、本発明の第１実施形態によるＸ線位相差撮像システム１００の構成について説明する。

（Ｘ線位相差撮像システムの構成）
まず、図１を参照して、本発明の第１実施形態によるＸ線位相差撮像システム１００の構成について説明する。

図１は、Ｘ線位相差撮像システム１００をＹ方向から見た図である。図１に示すように、Ｘ線位相差撮像システム１００は、Ｘ線源１と、第１格子２と、第２格子３と、検出器４と、画像処理部５と、装置制御部９と、格子移動機構１０とを備えている。なお、本明細書において、Ｘ線源１から第１格子２に向かう方向をＺ２方向、その逆向きの方向をＺ１方向とする。また、Ｚ方向と直交する面内の左右方向をＹ方向とし、紙面の奥に向かう方向をＹ２方向、紙面の手前側に向かう方向をＹ１方向とする。また、Ｚ方向と直交する面内の上下方向をＸ方向とし、上方向をＸ１方向、下方向をＸ２方向とする。

Ｘ線源１は、装置制御部９からの信号に基づいて高電圧が印加されることにより、Ｘ線を発生させるとともに、発生させたＸ線を検出器４（Ｚ２方向）に向けて照射するように構成されている。

第１格子２は、Ｘ方向に所定の周期（ピッチ）ｐ_１で配列される複数のスリット２ａおよび、Ｘ線位相変化部２ｂを有している。各スリット２ａおよびＸ線位相変化部２ｂはそれぞれ、直線状に延びるように形成されている。また、各スリット２ａおよびＸ線位相変化部２ｂはそれぞれ、平行に延びるように形成されている。第１格子２は、いわゆる位相格子である。

第１格子２は、Ｘ線源１と、第２格子３との間に配置されており、Ｘ線源１からＸ線が照射される。第１格子２は、タルボ効果により、第１格子２の自己像２０（図２参照）を形成するために設けられている。可干渉性を有するＸ線が、スリットが形成された格子を通過すると、格子から所定の距離（タルボ距離）離れた位置に、格子の像（自己像２０）が形成される。これをタルボ効果という。

第２格子３は、Ｘ方向に所定の周期（ピッチ）ｐ_２で配列される複数のＸ線透過部３ａおよびＸ線吸収部３ｂを有する。各Ｘ線透過部３ａおよびＸ線吸収部３ｂはそれぞれ、直線状に延びるように形成されている。また、各Ｘ線透過部３ａおよびＸ線吸収部３ｂはそれぞれ、平行に延びるように形成されている。第２格子３は、いわゆる、吸収格子である。第１格子２、第２格子３はそれぞれ異なる役割を持つ格子であるが、スリット２ａおよびＸ線透過部３ａはそれぞれＸ線を透過させる。また、Ｘ線吸収部３ｂはＸ線を遮蔽する役割を担っており、Ｘ線位相変化部２ｂはスリット２ａとの屈折率の違いによってＸ線の位相を変化させる。

第２格子３は、第１格子２と検出器４との間に配置されており、第１格子２を通過したＸ線が照射される。また、第２格子３は、第１格子２からタルボ距離離れた位置に配置される。第２格子３は、第１格子２の自己像２０と干渉して、検出器４の検出表面上にモアレ縞（図示せず）を形成する。

検出器４は、Ｘ線を検出するとともに、検出されたＸ線を電気信号に変換し、変換された電気信号を画像信号として読み取るように構成されている。検出器４は、たとえば、ＦＰＤ（Ｆｌａｔ Ｐａｎｅｌ Ｄｅｔｅｃｔｏｒ）である。検出器４は、複数の変換素子（図示せず）と複数の変換素子上に配置された画素電極（図示せず）とにより構成されている。複数の変換素子および画素電極は、所定の周期（画素ピッチ）で、Ｘ方向およびＹ方向にアレイ状に配列されている。また、検出器４は、取得した画像信号を、画像処理部５に出力するように構成されている。

画像処理部５は、画像生成部６、制御部７および記憶部８を含んでいる。画像生成部６は、検出器４から出力された画像信号に基づいて、Ｘ線画像（図示せず）を生成するように構成されている。また、画像生成部６は、生成したＸ線画像に基づいて、位相コントラスト画像１１（図３（Ｂ）参照）を生成するように構成されている。画像生成部６は、たとえば、ＧＰＵ（Ｇｒａｐｈｉｃｓ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）または画像処理用に構成されたＦＰＧＡ（Ｆｉｅｌｄ−Ｐｒｏｇｒａｍｍａｂｌｅ Ｇａｔｅ Ａｒｒａｙ）などのプロセッサを含む。

また、制御部７は、検出器４で検出したＸ線の強度変化３０（図２参照）の周期ｐなどを取得するように構成されている。制御部７が強度変化３０の周期ｐを取得する詳細な構成については後述する。また、制御部７は、たとえば、ＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）などのプロセッサを含む。

記憶部８は、画像生成部６が生成したＸ線画像や位相コントラスト画像１１を生成するためのプラグラムなどを保存するように構成されている。記憶部８は、たとえば、ＨＤＤ（ハードディスクドライブ）や不揮発性のメモリなどを含む。

装置制御部９は、格子移動機構１０を介して、第２格子３を格子面内において格子パターンの方向（Ｙ方向）と直交する方向（Ｘ方向）にステップ移動させるように構成されている。また、装置制御部９は、たとえば、ＣＰＵなどのプロセッサを含む。

格子移動機構１０は、装置制御部９からの信号に基づいて、第２格子３を格子面内（ＸＹ面内）において格子の格子パターンの方向（Ｙ方向）と直交する方向（Ｘ方向）にステップ移動させるように構成されている。具体的には、格子移動機構１０は、第２格子３の周期ｐ_２をＮ分割し、ｐ_２／Ｎずつ第２格子３をステップ移動させる。なお、Ｎは正の整数であり、第１実施形態では、たとえば、Ｎ＝４である。また、格子移動機構１０は、たとえば、ステッピングモータやピエゾアクチュエータなどを含む。

（各画素における画素値の強度変化を取得する構成）
次に、図２を参照して、画像処理部５が各画素における画素値の強度変化３０を取得する構成について説明する。

図２は、第１実施形態による画像生成部６が縞走査法によって位相コントラスト画像１１を生成する際に格子移動機構１０によって第２格子３を並進移動させた際の第１格子２の自己像２０と第２格子３との位置関係を示す模式図（Ａ）〜（Ｄ）および、第２格子３を並進移動させて撮影したＸ線画像の各画素における画素値の強度変化３０を示したグラフである。

第１実施形態では、装置制御部９は、格子移動機構１０を介して第２格子３をｐ_２／４ずつ４回並進移動させて撮影を行うように構成されている。すなわち、第１実施形態では、図２（Ａ）〜図２（Ｄ）に示すように、Ｘ線位相差撮像システム１００は、第２格子３をｐ_２／４ずつＸ２方向に移動させつつ撮像を行うように構成されている。また、制御部７は、第２格子３をｐ_２／４ずつＸ２方向に移動させつつ撮像を行った際の各画素の画素値の強度変化として、図２のグラフに示すような強度変化３０を取得する。

（比較例の説明）
ここで、一般的な縞走査法では、被写体Ｑを配置せずに撮像したＸ線画像における強度変化３０と、被写体Ｑを配置して撮像したＸ線画像における強度変化３０とに基づいて、位相コントラスト画像１１を生成する。たとえば、第２格子３をＮ回並進させて撮影した場合を考える。

被写体Ｑを配置して撮影した際の検出器４で検出されるＸ線の強度をＩ_ｋ（ｘ、ｙ）、被写体Ｑを配置せずに撮影した際のＸ線の強度をＩ_０ｋ（ｘ、ｙ）とし、以下のＳ（ｘ，ｙ）およびＳ_０（ｘ、ｙ）を定義する。

ここで、ｋは１からＮまでの正の整数である。また、ｘおよびｙは、各画素のｘ座標およびｙ座標である。

位相微分像１１ａを示す各画素の画素値は、上記式（１）および式（２）を用いて、以下に示す式（３）によって算出することができる。なお、位相微分像１１ａの各画素の画素値φ（ｘ、ｙ）は、被写体Ｑを配置した場合と配置しない場合における位相値の差分である。

ここで、Ｚ_Ｔは、第１格子２と第２格子３との間の距離である。また、吸収像１１ｂを示す各画素の画素値Ｔ（ｘ、ｙ）は、以下の式（４）によって表される。

被写体Ｑを配置した際のＶｉｓｉｂｉｌｉｔｙであるＶ（ｘ、ｙ）および、被写体Ｑを配置しない場合のＶｉｓｉｂｉｌｉｔｙであるＶ_０（ｘ、ｙ）は、以下の式（５）および式（６）により表される。

暗視野像１１ｃは、被写体Ｑを配置した場合と配置しない場合とのＶｉｓｉｂｉｌｉｔｙの比率であるため、以下の式（７）で表される。

上記式（１）および式（２）に示す通り、一般的な縞走査法では、1ステップにつき、第２格子３の周期ｐ_２をステップ数で割った値だけ移動することを前提として計算された値に基づいて、位相コントラスト画像１１を画像化している。言い換えると、上記式（１）および式（２）は、各ステップで取得される画像が、第２格子３の周期ｐ_２の１／４、２／４、３／４、４／４の各位相における画像であることを前提としているため、第２格子３を並進移動させながら撮影することによって、１／４〜４／４の各位相の画像を確実に取得する必要がある。したがって、取得される強度変化３０の周期ｐと、第２格子３の周期ｐ_２とにずれが生じた場合、１／４〜４／４の各位相からずれた位相における画像（画素値）を用いて演算を行うことになるため、生成される位相コントラスト画像１１にアーチファクト１２が形成される。

（強度変化の周期のずれ）
次に、図３を参照して、強度変化３０の周期ｐと第２格子３の周期ｐ_２とにずれが生じた際に、位相コントラスト画像１１に生じるアーチファクト１２について説明する。

図３（Ａ）に示すグラフは、強度変化３０の周期ｐと第２格子３の周期ｐ_２とが一致している場合の強度変化３０ａ（直線グラフ）および強度変化３０の周期ｐと第２格子３の周期ｐ_２とにずれが生じている場合の強度変化３０ｂ（破線グラフ）を示す模式図である。なお、強度変化３０ｂは、実際にはＮ点のデータであるが、図３（Ａ）に示す例では、Ｎ点のデータのそれぞれをｋ／Ｎ周期におけるデータとみなしている。ここで、Ｘ線位相差撮像システム１００では、Ｘ線源１から生じる熱に起因して、第２格子３が熱変形を起こす場合がある。第２格子３が熱変形を起こした場合、第２格子３の周期ｐ_２が変化する。また、格子移動機構１０は、第２格子３を第１格子２と平行に移動させるように構成されているが、設置の際の位置ずれなどに起因して、第２格子３の移動方向が第１格子２とは平行にならない場合がある。このように、第２格子３の周期ｐ_２の変化および／または第２格子３の移動方向にずれが生じた場合、第２格子３の格子パターンの方向（Ｙ方向）と直交する方向（Ｘ方向）における総移動量と、第２格子３の周期ｐ_２とにずれが生じる。第２格子３の格子パターンの方向（Ｙ方向）と直交する方向（Ｘ方向）における総移動量と、第２格子３の周期ｐ_２とにずれが生じた場合、制御部７が取得する強度変化３０の周期ｐが第２格子３の周期ｐ_２とずれる。すなわち、強度変化３０の周期ｐは、第２格子３の周期ｐ_２のずれと格子移動機構１０による第２格子３の移動量のずれとのうち少なくとも一方のずれを含んだ周期ｐである。

一般的な縞走査法は、上記式（１）〜式（７）に示したように、得られた強度変化３０ｂを、第２格子３の周期ｐ_２と同じ周期である強度変化３０ａを表す関数に当てはめることにより位相コントラスト画像１１を生成している。したがって、強度変化３０の周期ｐと第２格子３の周期ｐ_２とにずれが生じた状態で、上記式（１）〜式（８）を用いた一般的な縞走査法によって位相コントラスト画像１１を生成した場合、図３（Ｂ）に示すように、各画像に縞状のアーチファクト１２が形成される。

そこで、第１実施形態では、画像処理部５は、強度変化３０ｂの周期ｐと、周期を変数とするとともに格子の移動に伴う画素値の強度変化３０ｂを表す関数とにより、位相コントラスト画像１１を生成するように構成されている。具体的には、画像生成部６は、強度変化３０の周期ｐと第２格子３の周期ｐ_２とにずれが生じた状態であっても、位相コントラスト画像１１にアーチファクト１２が生じることを抑制できるように、以下に示す処理によって位相コントラスト画像１１を生成する。

（位相コントラスト画像の生成処理）
第１実施形態では、制御部７は、第２格子３を移動させて取得した強度変化３０ｂのデータと関数とに基づいて、強度変化３０ｂの周期ｐを決定するように構成されている。具体的には、制御部７は、強度変化３０ｂの波形の形状と関数の波形の形状とをフィッティングすることにより強度変化３０ｂの周期ｐを決定するように構成されている。制御部７は、取得した強度変化３０ｂ上のＮ個のデータ点に対して、最小二乗フィッティングを行うことにより、関数の係数ａ〜ｃを決定し、関数の偏差の二乗和Ｚが最小となる周期ｐを取得するように構成されている。

Ｎ個のデータ点のＸ座標ｘ_ｋは以下の式（８）により表すことができる。

フィッティングには、以下の式（９）に示す関数を用いる。

ここで、ｐは強度変化３０ｂの周期である。また、ａ、ｂおよびｃは、周期ｐが決定されることにより決定される関数の係数である。

画像処理部５は、所定の周期に対応する上記式（９）の係数ａ〜ｃを取得するように構成されている。また、制御部７は、取得した係数ａ〜ｃおよび所定の周期ｐにより強度変化３０ｂの波形の形状と、上記式（９）の波形の形状との偏差の二乗和Ｚを取得するように構成されている。具体的には、制御部７は、上記式（９）の関数とＮ個のデータ点との偏差の二乗和Ｚを計算することにより、式（９）の周期ｐを取得するように構成されている。偏差の二乗和Ｚは、以下の式（１０）のように示すことができる。

ここで、ｍは、強度変化３０ｂの周波数であり、ｍ＝２π／ｐである。

第１実施形態では、制御部７は、最小二乗法を用いてフィッティングを行うため、上記式（１０）の係数ａ、ｂ、およびｃの偏微分をとる。係数ａ、ｂ、およびｃの偏微分は、それぞれ以下の式（１１）〜式（１３）で表される。

ここで、ｙ_ｋは、各画素における画素値である。

最小二乗法の定義により、上記式（１１）〜式（１３）のそれぞれが０になるようなａ、ｂ、ｃを求めればよいため、上記式（１１）〜式（１３）を以下に示す式（１４）〜式（１６）のように変形する。

上記式（１４）〜式（１６）を用いて連立方程式を解くことにより、ａ、ｂ、ｃは以下の式（１７）〜式（１９）によって表すことができる。

このとき、ａ_１、ｂ_１、ｃ_１およびｄは、それぞれ、以下の式（２０）〜式（２３）によって表すことができる。

上記式（１７）〜式（１９）で取得した係数ａ、ｂ、ｃ、および上記式（９）に含まれる周期ｐ（周波数ｍ）は未知数である。したがって、第１実施形態では、画像処理部５は、偏差の二乗和Ｚが最小となる周期ｐを取得するように構成されている。図３（Ａ）に示す例では、強度変化３０ｂの周期ｐは、第２格子３の周期ｐ_２よりも小さくなるようにずれている。そのため、制御部７は、周期ｐの値を第２格子３の周期ｐ_２よりも小さくした複数の周期ｐを用いて偏差の二乗和Ｚを取得すればよい。しかし、実際には強度変化３０ｂの周期ｐは未知数である。そのため、制御部７は、少なくとも格子移動機構１０により移動される第２格子３の周期ｐ_２の設計値よりも大きい値および小さい値の両方の値を周期として用いて強度変化３０ｂの波形の形状と式（９）の波形の形状との偏差の二乗和Ｚを求めるように構成されている。

たとえば、制御部７は、第２格子３の周期ｐ_２よりも大きい値を周期として用いた第１偏差と、第２格子３の周期ｐ_２よりも小さい値を周期として用いた第２偏差とを取得する。制御部７は、取得した第１偏差と第２偏差とを比較し、第１偏差の方が小さければ、第２格子３の周期ｐ_２よりも大きい複数の周期ｐを用いて複数の偏差の二乗和Ｚを取得する。一方、第２偏差の方が小さければ、制御部７は、第２格子３の周期ｐ_２よりも小さい複数の周期ｐを用いて複数の偏差の二乗和Ｚを取得する。制御部７は、取得した複数の偏差の二乗和Ｚが最も小さくなる際の周期ｐを強度変化１０ｂの周期ｐとして決定するように構成されている。

図４は、強度変化３０ｂの周期ｐを第２格子３の周期ｐ_２よりも小さくした複数の周期を用いて位相微分像１１ａ、吸収像１１ｂおよび暗視野像１１ｃを生成した場合の例を示している。図４における解析周期とは、第２格子３の周期ｐ_２を所定の係数で乗算することにより得られる周期の事である。たとえば、解析周期０．９５とは、第２格子３の周期ｐ_２に０．９５をかけた周期（０．９５×ｐ_２）のことである。図４に示すように、解析周期を１．００から０．９５へと小さくするにつれて、各画像に生じているアーチファクト１２が軽減されていることがわかる。また、解析周期を０．９５からさらに小さくするとアーチファクト１２が増加していることもわかる。したがって、図４に示す例では、偏差の二乗和Ｚが最も小さくなる周期は、解析周期０．９５（０．９５×ｐ_２）となる。なお、図４に示す例は、強度変化３０ｂの周期ｐを０．９５×ｐ_２としてシミュレーションを行った例である。したがって、実際に強度変化３０ｂの周期ｐを決定するために、制御部７は、偏差の二乗和Ｚの大きさが小さくなる周期ｐから大きくなる周期ｐの変曲点を含むような複数の周期ｐを設定し、偏差の二乗和Ｚが極小値となるような周期ｐを強度変化３０ｂの周期ｐと決定するように構成されている。

第１実施形態では、画像生成部６は、上記式（１７）〜式（２３）により算出した係数ａ、ｂ、およびｃと決定した周期ｐとを用いて、位相コントラスト画像１１の各画素値を算出する。位相コントラスト画像１１の各画素値を算出する際には、上記式（９）を以下に示す式（２４）に書き換えた式を用いる。

ここで、係数Ａ、Ｂ、Ｃと係数ａ、ｂ、ｃとの間には、以下の式（２５）〜式（２７）に示す関係がある。なお、以下に示す式（２６）の単位はｒａｄである。

上記式（２４）の係数Ａ、Ｂ、Ｃを用いると、位相微分像１１ａ、吸収像１１ｂおよび暗視野像１１ｃはそれぞれ以下の式（２８）〜式（３０）で表すことができる。

ここで、φは位相微分像である。また、Ｔは吸収像である。また、Ｄは暗視野像である。また、下付き文字のｓおよびｒは、それぞれ、被写体Ｑありで撮像した場合の係数および被写体Ｑなしで撮像した場合の係数を示している。

以上のように、第１実施形態では、制御部７は、最小二乗法を用いて、実際に撮影されたＸ線画像の強度変化３０ｂの波形の形状と、式（９）の波形の形状とのフィッティングを行うように構成されている。したがって、強度変化３０ｂの周期ｐが第２格子３の周期ｐ_２とずれている場合でも、強度変化３０ｂの周期ｐに対応した関数（式（９））を所得することにより、位相コントラスト画像１１を生成することができる。

（位相コントラスト画像の生成方法）
次に、図５を参照して、第１実施形態による位相コントラスト画像１１の生成方法について説明する。

ステップＳ１において、画像生成部６は、被写体Ｑを撮影したＸ線画像を取得する。その後、制御部７は被写体Ｑを撮影したＸ線画像の各画素における画素値の強度変化３０ｂのデータを取得する。次に、ステップＳ２において、制御部７は、第２格子３の周期ｐ_２よりも大きい値の第１周期を関数に代入して第１偏差の二乗和を取得する。また、制御部７は、第２格子３の周期ｐ_２よりも小さい値の第２周期を関数に代入し、第２偏差の二乗和を取得する。なお、ステップＳ２では、第１偏差の二乗和および第２偏差の二乗和を取得するために、上記式（９）、式（１７）〜式（２３）およびｍ＝２π／ｐの関係式を用いて、式（９）に含まれる変数を周期ｐのみにした式を用いる。その後、処理はステップＳ３へ進む。

ステップＳ３において、制御部７は、第１偏差の二乗和および第２偏差の二乗和の値の大きさを比較する。第１偏差の二乗和が第２偏差の二乗和よりも小さい場合、処理はステップＳ４へ進む。第１偏差の二乗和が第２偏差の二乗和よりも大きい場合、処理はステップＳ５へ進む。

ステップＳ４において、制御部７は、第２格子３の周期ｐ_２よりも小さい値の複数の周期ｐを設定する。また、ステップＳ５では、制御部７は、第２格子３の周期ｐ_２よりも大きい値の複数の周期ｐを設定する。その後、処理はステップＳ６へ進む。

ステップＳ６において、制御部７は、ステップＳ４（ステップＳ５）において設定した複数の周期ｐを用いて複数の偏差の二乗和Ｚを取得する。その後、制御部７は、複数の偏差の二乗和Ｚの値が最も小さくなる周期ｐを取得し、関数（式（９））の周期ｐとする。その後、ステップＳ７において、制御部７は、関数（式（９））の係数ａ〜ｃを取得する。その後、ステップＳ８において、画像生成部６は、ステップＳ６において取得した周期ｐが決定された関数（式（９））およびステップＳ７において取得された係数ａ〜ｃを用いて、位相コントラスト画像１１を生成し、処理を終了する。

（第１実施形態の効果）
第１実施形態では、以下のような効果を得ることができる。

第１実施形態では、上記のように、Ｘ線位相差撮像システム１００は、Ｘ線源１と、Ｘ線源１からＸ線が照射される第１格子２と第１格子２からのＸ線が照射される第２格子３とを含む複数の格子と、Ｘ線源１から照射されたＸ線を検出する検出器４と、第２格子３を移動させる格子移動機構１０と、検出器４により検出された各画素の画素値の変化を表す強度変化３０ｂに基づいて、位相コントラスト画像１１を生成する画像処理部５とを備え、画像処理部５は、強度変化３０ｂの周期ｐと、周期を変数とするとともに格子の移動に伴う画素値の強度変化３０ｂを表す関数（式（９））とにより、位相コントラスト画像１１を生成するように構成されている。これにより、関数（式（９））の周期ｐが変数であるため、関数（式（９））の周期ｐの最適解を求めることにより、第２格子３の移動によって実際に取得された強度変化３０ｂの周期ｐを近似することができる。その結果、各画素における画素値の強度変化３０ｂの周期ｐと第２格子３の周期ｐ_２とにずれが生じている場合でも、関数（式（９））の周期ｐを強度変化３０ｂの周期ｐに合わせることが可能となるので、得られる位相コントラスト画像１１にアーチファクト１２が生じることを抑制することができる。

また、第１実施形態では、上記のように、強度変化３０ｂの周期ｐは、第２格子３の周期ｐ_２のずれと格子移動機構１０による第２格子３の移動量のずれとのうち少なくとも一方のずれを含んだ周期である。これにより、第２格子３の周期ｐ_２および第２格子３の移動量のどちらにずれが生じていた場合（または、両方にずれが生じた場合）でも、ずれを含んだ強度変化３０ｂの周期ｐに、関数（式（９））の周期ｐ（変数）を合わせることにより、生成する位相コントラスト画像１１にアーチファクト１２が生じることを抑制することができる。

また、第１実施形態では、上記のように、画像処理部５は、第２格子３を移動させて取得した強度変化３０ｂのデータと関数（式（９））とに基づいて、強度変化３０ｂの周期ｐを決定するように構成されている。これにより、実際に計測した強度変化３０ｂのデータから強度変化３０ｂの周期ｐを取得することができる。その結果、強度変化３０ｂの周期ｐと第２格子３の周期ｐ_２とにずれが生じている場合でも、ずれを含んだ強度変化３０ｂの周期ｐを取得することができる。

また、第１実施形態では、上記のように、画像処理部５は、強度変化３０ｂの波形の形状と関数（式（９））の波形の形状とをフィッティングすることにより強度変化３０ｂの周期ｐを決定するように構成されている。これにより、強度変化３０ｂの波形の形状と関数（式（９））の波形の形状との偏差の二乗和Ｚに基づいて強度変化３０ｂの周期ｐを決定することができる。その結果、複数の周期を用いて取得した偏差の二乗和Ｚを比較することにより強度変化３０ｂの周期ｐを決定することが可能となるので、たとえば、複数の周期ｐを用いて取得した画像を画像処理することにより強度変化３０ｂの周期ｐを決定する場合と比較して、強度変化３０ｂの周期ｐを決定する処理を簡素化することができる。

また、第１実施形態では、上記のように、少なくとも格子移動機構１０により移動される第２格子３の周期ｐ_２の設計値よりも大きい値および小さい値の両方の値を周期ｐとして用いて強度変化３０ｂの波形の形状と関数（式（９））の波形の形状との偏差の二乗和Ｚを求めるように構成されている。これにより、第２格子３の周期ｐ_２の設計値よりも大きい値により算出される第１偏差と第２格子３の周期ｐ_２の設計値よりも小さい値により算出される第２偏差とを比較することにより、偏差の二乗和Ｚを最小にするための周期ｐが、第２格子３の周期ｐ_２の設計値よりも大きい値か小さい値かを判断することができる。その結果、偏差の二乗和Ｚが最小になる周期ｐを容易に決定することができる。

また、第１実施形態では、上記のように、画像処理部５は、所定の周期ｐに対応する関数（式（９））の係数ａ〜ｃを取得し、係数ａ〜ｃおよび所定の周期ｐにより強度変化３０ｂの波形の形状と関数（式（９））の波形の形状との偏差の二乗和Ｚを取得するように構成されている。これにより、所定の周期ｐに対応する偏差の二乗和Ｚを取得することができる。その結果、複数の周期ｐに対応した偏差の二乗和Ｚを比較することにより、偏差の二乗和Ｚが最小となる周期ｐを容易に決定することができる。

［第２実施形態］
次に、図１および図６を参照して、本発明の第２実施形態によるＸ線位相差撮像システム２００について説明する。複数の周期ｐを用いて取得した複数の偏差の二乗和Ｚを比較することにより関数の周期ｐを決定する第１実施形態とは異なり、第２実施形態では、画像処理部５は、あらかじめ取得された強度変化３０ｂの周期ｐから、関数（式（９））の係数ａ〜ｃを取得し、取得した係数ａ〜ｃおよびあらかじめ取得された強度変化３０ｂの周期ｐに基づき、位相コントラスト画像１１を生成するように構成されている。なお、上記第１実施形態と同様の構成については同様の符号を付し、説明を省略する。

ここで、第２格子３の周期ｐ_２は非常に小さい値なので、製作誤差によりわずかに周期ｐ_２がずれる場合がある。その場合、第２格子３の周期ｐ_２からわずかにずれた周期（以下、ｐ_ｘとする）が既知であっても、上記式（１）および式（２）を用いる一般的な縞走査法では、ｋ／Ｎ周期のデータを取得する必要がある。しかし、格子移動機構１０の精度が十分でないと、既知の周期ｐ_ｘをＮ等分した移動量（ｐ_ｘ×（ｋ／Ｎ））で第２格子３を並進移動させることができないので、ｋ／Ｎ周期のデータを取得することができず、アーチファクト１２が生じる。そこで、第２実施形態では、制御部７は、あらかじめ第２格子３の周期ｐ_２のずれを計測するように構成されている。また、制御部７は、実際の測定値と、あらかじめ取得された強度変化３０ｂの周期ｐ_２を代入した上記式（９）とをフィッティングさせることにより、係数ａ〜ｃを決定するように構成されている。

（位相コントラスト画像の生成方法）
次に、図６を参照して、第２実施形態による位相コントラスト画像１１の生成方法について説明する。

ステップＳ１において、制御部７は、被写体Ｑを撮像したＸ線画像の各画素における画素値の強度変化３０ｂのデータを取得する。次に、ステップＳ１０において、制御部７は、強度変化３０ｂの周期ｐを取得する。その後、ステップＳ１１において、制御部７は、あらかじめ取得された強度変化３０ｂの周期ｐから、式（９）の係数ａ〜ｃを取得する。その後、処理はステップＳ８へ進む。

ステップＳ８において、画像生成部６は、取得した係数ａ〜ｃおよびあらかじめ取得された強度変化３０ｂの周期ｐに基づき、位相コントラスト画像１１を生成し、処理を終了する。

なお、第２実施形態のその他の構成は、上記第１実施形態と同様である。

（第２実施形態の効果）
第２実施形態では、以下のような効果を得ることができる。

第２実施形態では、上記のように、画像処理部５は、あらかじめ取得された強度変化３０ｂの周期ｐから、関数（式（９））の係数ａ〜ｃを取得し、取得した係数ａ〜ｃおよびあらかじめ取得された強度変化３０ｂの周期ｐに基づき、位相コントラスト画像１１を生成するように構成されている。これにより、実際の計測値と、あらかじめ取得された強度変化３０ｂの周期ｐを代入した関数（式（９））とをフィッティングすることにより、関数（式（９））の係数ａ〜ｃを取得することができる。その結果、第２格子３の移動量がずれている場合でも、実データに精度よく近似した関数（式（９））を決定することが可能となり、生成される位相コントラスト画像１１にアーチファクト１２が生じることを抑制することができる。

なお、第２実施形態のその他の効果は、上記第１実施形態と同様である。

（変形例）
なお、今回開示された実施形態は、すべての点で例示であって制限的なものではないと考えられるべきである。本発明の範囲は、上記した実施形態の説明ではなく、特許請求の範囲によって示され、さらに特許請求の範囲と均等の意味および範囲内でのすべての変更（変形例）が含まれる。

たとえば、上記第１〜第３実施形態では、複数の格子として、第１格子２および第２格子３を設ける例を示したが、本発明はこれに限られない。たとえば、図７に示すＸ線位相差撮像システム３００のように、Ｘ線源１と第１格子２との間に、第３格子４０を設ける構成でもよい。第３格子４０は、Ｘ方向に所定の周期（ピッチ）ｐ_３で配列される複数のスリット４０ａおよびＸ線吸収部４０ｂを有している。各スリット４０ａおよびＸ線吸収部４０ｂはそれぞれ、直線状に延びるように形成されている。また、各スリット４０ａおよびＸ線吸収部４０ｂはそれぞれ、平行に延びるように形成されている。また、第３格子４０は、Ｘ線源１と第１格子２との間に配置されており、Ｘ線源１からＸ線が照射される。第３格子４０は、各スリット４０ａを通過したＸ線を、各スリット４０ａの位置に対応する線光源とするように構成されている。これにより、第３格子４０によってＸ線源１から照射されるＸ線の可干渉性を高めることができる。その結果、Ｘ線源１の焦点径に依存することなく第１格子２の自己像２０を形成させることが可能となるので、Ｘ線源１の選択の自由度を向上させることができる。

また、上記第１および第２実施形態では、第１格子２として、位相格子を用いる例を示したが、本発明はこれに限られない。たとえば、第１格子２として、吸収格子を用いてもよい。第１格子２として吸収格子を用いた場合、画像処理部５（画像生成部６）は、第１格子２を透過したＸ線の縞模様と第２格子３とにより位相コントラスト画像１１を生成する。したがって、第１格子２の自己像２０を用いることなく位相コントラスト画像１１を取得することが可能となるので、第１格子２の配置位置の自由度を向上させることができる。しかし、第１格子２として吸収格子を用いる場合、得られる位相コントラスト画像１１の画質が低下するため、高画質の位相コントラスト画像１１を得たい場合は、第１格子２として位相格子を用いる方が好ましい。

また、上記第１実施形態では、偏差の二乗和Ｚを取得する際に用いる複数の周期ｐとして、第２格子３の周期ｐ_２よりも小さい値を用いる例を示したが、本発明はこれに限られない。たとえば、強度変化３０ｂの周期ｐが第２格子３の周期ｐ_２よりも大きい場合、複数の周期ｐとして、第２格子３の周期ｐ_２よりも大きい値を用いてもよい。

また、上記第１実施形態では、強度変化３０ｂの波形の形状と関数の波形の形状のフィッティングとして最小二乗法を用いる例を示したが、本発明はこれに限られない。関数（式（９））の係数ａ〜ｃを取得することが可能であれば、どのようなフィッティング方法を用いてもよい。

また、上記第１および第２に実施形態では、格子移動機構１０が第２格子３を並進移動させる例を示したが、本発明はこれに限られない。第１格子２を並進移動させてもよい。移動させる格子はいずれの格子でもよい。

また、上記第１および第２実施形態では、制御部７によって関数（式（９））の偏差の二乗和Ｚを決定する例を示したが、本発明はこれに限られない。たとえば、図４に示すような周期ｐを変化させた複数の位相コントラスト画像１１を表示する表示部を備え、表示部に表示された位相コントラスト画像１１をユーザが選択することにより、関数（式（９））の偏差の二乗和Ｚを決定するように構成されていてもよい。

１ Ｘ線源
２ 第１格子
３ 第２格子
４ 検出器
５ 画像処理部
１０ 格子移動機構
１２ 位相コントラスト画像
３０、３０ａ、３０ｂ 強度変化
４０ 第３格子
１００、２００、３００ Ｘ線位相差撮像システム
Ｑ 被写体
ｐ 強度変化の周期
ｐ_２ 格子の周期
Ｚ 偏差の二乗和

Claims (8)

1. Ｘ線源と、
   前記Ｘ線源からＸ線が照射される第１格子と、前記第１格子からのＸ線が照射される第２格子とを含む複数の格子と、
   前記Ｘ線源から照射されたＸ線を検出する検出器と、
   前記複数の格子の少なくともいずれか１つを移動させる格子移動機構と、
   前記検出器により検出された各画素の画素値の変化を表す強度変化に基づいて、位相コントラスト画像を生成する画像処理部とを備え、
   前記画像処理部は、前記強度変化の周期と、周期を変数とするとともに格子の移動に伴う画素値の強度変化を表す関数とにより、前記位相コントラスト画像を生成するように構成されている、Ｘ線位相差撮像システム。
2. 前記強度変化の周期は、格子の周期のずれと前記格子移動機構による格子の移動量のずれとのうち少なくとも一方のずれを含んだ周期である、請求項１に記載のＸ線位相差撮像システム。
3. 前記画像処理部は、前記複数の格子の少なくともいずれか１つを移動させて取得した前記強度変化のデータと前記関数とに基づいて、前記強度変化の周期を決定するように構成されている、請求項１に記載のＸ線位相差撮像システム。
4. 前記画像処理部は、前記強度変化の波形の形状と前記関数の波形の形状とをフィッティングすることにより前記強度変化の周期を決定するように構成されている、請求項１に記載のＸ線位相差撮像システム。
5. 前記画像処理部は、少なくとも前記格子移動機構により移動される格子の周期の設計値よりも大きい値および小さい値の両方の値を周期として用いて、前記強度変化の波形の形状と前記関数の波形の形状との偏差を求めるように構成されている、請求項４に記載のＸ線位相差撮像システム。
6. 前記画像処理部は、所定の周期に対応する前記関数の係数を取得し、前記係数および前記所定の周期により前記強度変化の波形の形状と前記関数の波形の形状との偏差を取得するように構成されている。請求項４に記載のＸ線位相差撮像システム。
7. 前記画像処理部は、あらかじめ取得された前記強度変化の周期から、前記関数の係数を取得し、取得した前記係数およびあらかじめ取得された前記強度変化の周期に基づき、前記位相コントラスト画像を生成するように構成されている、請求項１に記載のＸ線位相差撮像システム。
8. 前記複数の格子は、前記Ｘ線源と前記第１格子との間に配置された第３格子をさらに含んでいる、請求項１に記載のＸ線位相差撮像システム。

Images