(57)【特許請求の範囲】
【請求項１】
Ｘ線の照射軸方向と直交する方向に複数のスリットが配列されて構成された複数の格子がＸ線の照射軸方向に並んで設けられ、前記複数の格子の何れかを移動させて縞走査を行うことにより複数のモアレ縞画像を取得するタルボ干渉計又はタルボ・ロー干渉計を用いたＸ線撮影装置と、
前記Ｘ線撮影装置により取得された複数のモアレ縞画像に基づいて、微分位相画像、Ｘ線吸収画像、小角散乱画像のうちの少なくとも一つの再構成画像を生成する再構成手段を備える画像処理装置と、
を備えるＸ線撮影システムであって、
前記画像処理装置は、
前記Ｘ線撮影装置において格子配置を異ならせて複数回の縞走査を行うことにより得られた複数のモアレ縞画像に基づき前記再構成手段により生成された再構成画像に基づいて、縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を推定する推定手段を備えるＸ線撮影システム。
【請求項２】
前記複数の格子は、周期パターンを形成する第１格子と、前記第１格子の周期パターンをモアレ縞に変換する第２格子とを含み、
前記推定手段は、前記Ｘ線撮影装置において前記第１格子と前記第２格子の相対角度を異ならせて複数回の縞走査を行うことにより得られた複数のモアレ縞画像に基づき前記再構成手段により生成された再構成画像に基づいて、縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を推定する請求項１に記載のＸ線撮影システム。
【請求項３】
前記画像処理装置は、
縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を示すパラメーターを設定する設定手段を備え、
前記再構成手段は、前記Ｘ線撮影装置により取得された複数のモアレ縞画像と、前記設定手段により設定されたパラメーターとに基づいて、再構成画像を生成し、
前記推定手段は、前記生成された再構成画像のモアレ縞の振幅の大きさを表す特徴量が予め定められた基準以下となるときの前記パラメーターを求め、求めたパラメーターを縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置として推定する請求項１又は２に記載のＸ線撮影システム。
【請求項４】
前記生成された再構成画像のモアレ縞の振幅の大きさを表す特徴量は、前記生成された再構成画像の各画素の画素値と予め定められた目標値との差分の二乗和である請求項３に記載のＸ線撮影システム。
【請求項５】
前記設定手段は、前記推定手段により推定された縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を前記パラメーターとして設定し、
前記再構成手段は、前記設定手段により設定されたパラメーターと、前記Ｘ線撮影装置においてＸ線経路上の被写体配置位置に被写体を載置して縞走査を行うことにより取得された被写体有りの複数のモアレ縞画像と、Ｘ線経路上の前記被写体配置位置に前記被写体を載置せずに縞走査を行うことにより取得された被写体無しの複数のモアレ縞画像と、に基づいて被写体の再構成画像を生成する請求項３又は４に記載のＸ線撮影システム。
【請求項６】
前記Ｘ線撮影装置は、前記推定手段により推定された縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置に基づいて、縞走査時における格子の移動量を調整する調整手段を備える請求項１〜４の何れか一項に記載のＸ線撮影システム。
【請求項７】
前記推定手段は、前記再構成手段により生成された再構成画像のうち小角散乱画像に基づいて、縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を推定する請求項１〜６の何れか一項に記載のＸ線撮影システム。
【請求項８】
Ｘ線の照射軸方向と直交する方向に複数のスリットが配列されて構成された複数の格子がＸ線の照射軸方向に並んで設けられ、前記複数の格子の何れかを移動させて縞走査を行うことにより複数のモアレ縞画像を取得するタルボ干渉計又はタルボ・ロー干渉計を用いたＸ線撮影装置により取得された複数のモアレ縞画像に基づいて、微分位相画像、Ｘ線吸収画像、小角散乱画像のうちの少なくとも一つの再構成画像を生成する再構成手段を備える画像処理装置であって、
前記Ｘ線撮影装置において格子配置を異ならせて複数回の縞走査を行うことにより得られた複数のモアレ縞画像に基づき前記再構成手段により生成された再構成画像に基づいて、縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を推定する推定手段を備える画像処理装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、Ｘ線撮影システム及び画像処理装置に関する。
【背景技術】
【０００２】
タルボ干渉計又はタルボ・ロー干渉計を用いたＸ線撮影装置により高精細の画像を得るには、格子のうちのひとつを格子のスリット周期の１／Ｍ（Ｍは正の整数、Ｍ＞２)ずつスリット周期方向に移動させＭ回撮影した画像（モアレ縞画像）を用いて再構成を行う縞走査方式が用いられる。しかし、格子のスリット周期は大きくても２０（μｍ）程度であり、縞走査を行う際の格子の移動量は１回あたり１（μｍ）のオーダーとなり、Ｍ回撮影を行う際の各撮影における格子の位置は０．０１（μｍ）のオーダー（スリット周期に対して１％以下）の精度で決められる必要がある。
【０００３】
格子の位置に誤差があるまま撮影を行うと、測定される物理量に誤差を生じるだけでなく、撮影されたモアレ縞画像に基づいて生成される再構成画像に画像ムラが残存することになり好ましくない。
【０００４】
そこで、例えば、特許文献１には、センサーなどを用いて格子の位置を計測して再構成処理時のパラメーターに反映させることにより、格子の位置の誤差による画質の劣化を防止する技術が記載されている。
【先行技術文献】
【特許文献】
【０００５】
【特許文献１】国際公開第２０１２／０５７０２２号公報
【発明の概要】
【発明が解決しようとする課題】
【０００６】
しかしながら、特許文献１に記載の技術においては、０．０１（μｍ）のオーダーで格子位置を計測できる特殊なセンサーをＸ線撮影装置に付加する必要があり、コストがかかる。
【０００７】
本発明の課題は、縞走査時に移動させた格子の相対位置を検出するための特殊なセンサーをＸ線撮影装置に付加することなく、安価に格子の相対位置の誤差に起因する画像劣化を防止できるようにすることである。
【課題を解決するための手段】
【０００８】
上記課題を解決するため、請求項１に記載の発明は、
Ｘ線の照射軸方向と直交する方向に複数のスリットが配列されて構成された複数の格子がＸ線の照射軸方向に並んで設けられ、前記複数の格子の何れかを移動させて縞走査を行うことにより複数のモアレ縞画像を取得するタルボ干渉計又はタルボ・ロー干渉計を用いたＸ線撮影装置と、
前記Ｘ線撮影装置により取得された複数のモアレ縞画像に基づいて、微分位相画像、Ｘ線吸収画像、小角散乱画像のうちの少なくとも一つの再構成画像を生成する再構成手段を備える画像処理装置と、
を備えるＸ線撮影システムであって、
前記画像処理装置は、
前記Ｘ線撮影装置において格子配置を異ならせて複数回の縞走査を行うことにより得られた複数のモアレ縞画像に基づき前記再構成手段により生成された再構成画像に基づいて、縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を推定する推定手段を備える。
【０００９】
請求項２に記載の発明は、請求項１に記載の発明において、
前記複数の格子は、周期パターンを形成する第１格子と、前記第１格子の周期パターンをモアレ縞に変換する第２格子とを含み、
前記推定手段は、前記Ｘ線撮影装置において前記第１格子と前記第２格子の相対角度を異ならせて複数回の縞走査を行うことにより得られた複数のモアレ縞画像に基づき前記再構成手段により生成された再構成画像に基づいて、縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を推定する。
【００１０】
請求項３に記載の発明は、請求項１又は２に記載の発明において、
前記画像処理装置は、
縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を示すパラメーターを設定する設定手段を備え、
前記再構成手段は、前記Ｘ線撮影装置により取得された複数のモアレ縞画像と、前記設定手段により設定されたパラメーターとに基づいて、再構成画像を生成し、
前記推定手段は、前記生成された再構成画像のモアレ縞の振幅の大きさを表す特徴量が予め定められた基準以下となるときの前記パラメーターを求め、求めたパラメーターを縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置として推定する。
【００１１】
請求項４に記載の発明は、請求項３に記載の発明において、
前記生成された再構成画像のモアレ縞の振幅の大きさを表す特徴量は、前記生成された再構成画像の各画素の画素値と予め定められた目標値との差分の二乗和である。
【００１２】
請求項５に記載の発明は、請求項１〜４の何れか一項に記載の発明において、
前記設定手段は、前記推定手段により推定された縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を前記パラメーターとして設定し、
前記再構成手段は、前記設定手段により設定されたパラメーターと、前記Ｘ線撮影装置においてＸ線経路上の被写体配置位置に被写体を載置して縞走査を行うことにより取得された被写体有りの複数のモアレ縞画像と、Ｘ線経路上の前記被写体配置位置に前記被写体を載置せずに縞走査を行うことにより取得された被写体無しの複数のモアレ縞画像と、に基づいて被写体の再構成画像を生成する。
【００１３】
請求項６に記載の発明は、請求項１〜４の何れか一項に記載の発明において、
前記Ｘ線撮影装置は、前記推定手段により推定された縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置に基づいて、縞走査時における格子の移動量を調整する調整手段を備える。
【００１４】
請求項７に記載の発明は、請求項１〜６の何れか一項に記載の発明において、
前記推定手段は、前記再構成手段により生成された再構成画像のうち小角散乱画像に基づいて、縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を推定する。
【００１５】
請求項８に記載の発明は、
Ｘ線の照射軸方向と直交する方向に複数のスリットが配列されて構成された複数の格子がＸ線の照射軸方向に並んで設けられ、前記複数の格子の何れかを移動させて縞走査を行うことにより複数のモアレ縞画像を取得するタルボ干渉計又はタルボ・ロー干渉計を用いたＸ線撮影装置により取得された複数のモアレ縞画像に基づいて、微分位相画像、Ｘ線吸収画像、小角散乱画像のうちの少なくとも一つの再構成画像を生成する再構成手段を備える画像処理装置であって、
前記Ｘ線撮影装置において格子配置を異ならせて複数回の縞走査を行うことにより得られた複数のモアレ縞画像に基づき前記再構成手段により生成された再構成画像に基づいて、縞走査時の各撮影における前記移動させた格子の基準位置からの相対位置を推定する推定手段を備える。
【発明の効果】
【００１６】
本発明によれば、縞走査時に移動させた格子の相対位置を検出するための特殊なセンサーをＸ線撮影装置に付加することなく、安価に格子の相対位置の誤差に起因する画像劣化を防止することができる。
【図面の簡単な説明】
【００１７】
【図１】本実施形態に係るＸ線撮影システムの全体構成を示す図である。
【図２】マルチスリットの平面図である。
【図３】図１の本体部の機能的構成を示すブロック図である。
【図４】図１のコントローラーの機能的構成を示すブロック図である。
【図５】タルボ干渉計の原理を説明する図である。
【図６】縞走査時の第２格子の相対位置を推定するための撮影〜画像処理（格子相対位置推定処理）の流れを示すフローチャートである。
【図７】シミュレーションにより作成した周期の異なるモアレ縞画像セット１とモアレ縞画像セット２の例を示す図である。
【図８】パラメーターχkに［表２］に示す初期値を設定し、図７のモアレ縞画像セット１、２を用いて再構成処理することにより得られた微分位相画像、Ｘ線吸収画像、小角散乱画像を示す図である。
【図９】図８に示す微分位相画像に補正処理を行った画像を示す図である。
【図１０】縞走査時の各撮影時の第２格子の相対位置とパラメーターχkを一致させて再構成処理を行ったときのＸ線吸収画像、小角散乱画像、微分位相画像である。
【発明を実施するための形態】
【００１８】
［Ｘ線撮影システムの構成］
以下、図面を参照して本発明の実施形態について説明する。
図１に、本実施形態に係るＸ線撮影システムの構成例を示す。Ｘ線撮影システムは、Ｘ線撮影装置１とコントローラー５を備える。Ｘ線撮影装置１はタルボ・ロー干渉計によるＸ線撮影（縞走査）を行い、コントローラー５は当該Ｘ線撮影により得られた複数のモアレ縞画像を用いて被写体の再構成画像を生成する。
【００１９】
Ｘ線撮影装置１は、図１に示すように、Ｘ線源１１、マルチスリット１２、被写体台１３、第１格子１４、第２格子１５、Ｘ線検出器１６、保持部１７、本体部１８等を備える。
Ｘ線撮影装置１は縦型であり、Ｘ線源１１、マルチスリット１２、被写体台１３、第１格子１４、第２格子１５、Ｘ線検出器１６は、この順序に重力方向であるｚ方向に配置される。Ｘ線源１１の焦点とマルチスリット１２間の距離をＤ１（ｍｍ）、Ｘ線源１１の焦点とＸ線検出器１６間の距離をＤ２（ｍｍ）、マルチスリット１２と第１格子１４間の距離をＤ３（ｍｍ）、第１格子１４と第２格子１５間の距離をＤ４（ｍｍ）で表す。なお、被写体台１３の位置は、第１格子１４と第２格子１５との間に設けられていてもよい。
【００２０】
距離Ｄ１は好ましくは５〜５００（ｍｍ）であり、さらに好ましくは５〜３００（ｍｍ）である。
距離Ｄ２は、一般的に撮影室の高さは３（ｍ）程度又はそれ以下であることから、少なくとも３０００（ｍｍ）以下であることが好ましい。なかでも、距離Ｄ２は４００〜３０００（ｍｍ）が好ましく、さらに好ましくは５００〜２０００（ｍｍ）である。
Ｘ線源１１の焦点と第１格子１４間の距離（Ｄ１＋Ｄ３）は、好ましくは３００〜３０００（ｍｍ）であり、さらに好ましくは４００〜１８００（ｍｍ）である。
Ｘ線源１１の焦点と第２格子１５間の距離（Ｄ１＋Ｄ３＋Ｄ４）は、好ましくは４００〜３０００（ｍｍ）であり、さらに好ましくは５００〜２０００（ｍｍ）である。
それぞれの距離は、Ｘ線源１１から照射されるＸ線の波長から、第２格子１５上に第１格子１４による格子像（自己像）が重なる最適な距離を算出し、設定すればよい。
【００２１】
Ｘ線源１１、マルチスリット１２、被写体台１３、第１格子１４、第２格子１５、Ｘ線検出器１６は、同一の保持部１７に一体的に保持され、ｚ方向における位置関係が固定されている。保持部１７はアーム状に形成され、本体部１８に設けられた駆動部１８ａによりｚ方向に移動可能に本体部１８に取り付けられている。
Ｘ線源１１は、緩衝部材１７ａを介して保持されている。緩衝部材１７ａは、衝撃や振動を吸収できる材料であれば何れの材料を用いてもよいが、例えばエラストマー等が挙げられる。Ｘ線源１１はＸ線の照射によって発熱するため、Ｘ線源１１側の緩衝部材１７ａは加えて断熱素材であることが好ましい。
【００２２】
Ｘ線源１１はＸ線管を備え、当該Ｘ線管によりＸ線を発生させて重力方向（ｚ方向）にＸ線を照射する。Ｘ線管としては、例えば医療現場で広く一般に用いられているクーリッジＸ線管や回転陽極Ｘ線管を用いることができる。陽極としては、タングステンやモリブデンを用いることができる。
Ｘ線の焦点径は、０．０３〜３（mm）が好ましく、さらに好ましくは０．１〜１（mm）である。
Ｘ線源１１のＸ線照射方向には、Ｘ線の照射範囲を狭めるための図示しない照射野絞りが設けられている。
【００２３】
マルチスリット１２は回折格子であり、図２に示すようにｘ方向に複数のスリットが所定間隔で設けられている。マルチスリット１２はシリコンやガラスといったＸ線の吸収率が低い材質の基板上に、タングステン、鉛、金といったＸ線の遮蔽力が大きい、つまりＸ線の吸収率が高い材質により形成される。例えば、フォトリソグラフィーによりレジスト層がスリット状にマスクされ、ＵＶが照射されてスリットのパターンがレジスト層に転写される。露光によって当該パターンと同じ形状のスリット構造が得られ、電鋳法によりスリット構造間に金属が埋め込まれて、マルチスリット１２が形成される。
【００２４】
マルチスリット１２のスリット周期は１〜６０（μｍ）である。スリット周期は、図２に示すように隣接するスリット間の距離を１周期とする。スリットの幅（ｘ方向の長さ）はスリット周期の１〜６０（％）の長さであり、さらに好ましくは１０〜４０（％）である。スリットの高さ（ｚ方向の長さ）は１〜５００（μｍ）であり、好ましくは１〜１５０（μｍ）である。
マルチスリット１２のスリット周期をｗ０（μｍ）、第２格子１５のスリット周期をｗ2（μｍ）とすると、スリット周期ｗ０は下記式により求めることができる。
ｗ０＝ｗ2・Ｄ３／Ｄ４
当該式を満たすように周期ｗ０を決定することにより、マルチスリット１２及び第１格子１４の各スリットを通過したＸ線により形成される自己像が、それぞれ第２格子１５上で重なり合い、いわばピントが合った状態とすることができる。
【００２５】
被写体台１３は、Ｘ線源１１からのＸ線照射経路上の被写体配置位置に設けられた、被写体を載置するための台である。
第１格子１４は、マルチスリット１２と同様にｘ方向に所定の周期の複数のスリットが設けられた回折格子である（図２参照）。第１格子１４は、マルチスリット１２と同様にＵＶを用いたフォトリソグラフィーによって形成することもできるし、いわゆるＩＣＰ法によりシリコン基板に微細細線で深掘加工を行い、シリコンのみで格子構造を形成することとしてもよい。第１格子１４のスリット周期は１〜２０（μｍ）である。スリットの幅はスリット周期の２０〜７０（％）であり、好ましくは３５〜６０（％）である。スリットの高さは１〜１００（μｍ）である。
【００２６】
第１格子１４として位相型を用いる場合、スリットの高さはスリット周期を形成する２種の素材、つまりＸ線透過部とＸ線遮蔽部の素材による位相差がπ／８〜１５×π／８となる高さとする。好ましくは、π／２又はπとなる高さである。第１格子１４として吸収型を用いる場合、スリットの高さはＸ線遮蔽部によりＸ線が十分吸収される高さとする。
【００２７】
第１格子１４が位相型である場合、第１格子１４と第２格子１５間の距離Ｄ４は、第１格子のスリット周期をｗ１（μｍ）として、次の条件をほぼ満たすことが必要である。
Ｄ４＝ｐ・（ｗ１2／λ）・｛Ｄ３／（Ｄ３−ｐ・ｗ１2／λ）｝
なお、ｐはタルボ次数で、第１格子１４の高さがπ／２に対応する場合はｐ＝ｎ／２、πに対応する場合はｐ＝ｎ／８(ｎは正の奇数)となる。また、λはＸ線の波長である。
【００２８】
第２格子１５は、マルチスリット１２と同様にｘ方向に所定の周期の複数のスリットが設けられた回折格子である（図２参照）。第２格子１５もフォトリソグラフィーにより形成することができる。第２格子１５のスリット周期は１〜２０（μｍ）である。スリットの幅はスリット周期の３０〜７０（％）であり、好ましくは３５〜６０（％）である。スリットの高さは１〜１００（μｍ）である。第１格子１４及び第２格子１５は、それぞれの格子面がｚ方向に対し垂直（ｘ−ｙ平面内で平行）である。
【００２９】
図１に示すように、第２格子１５に隣接して、第２格子１５をスリット周期方向であるｘ方向に移動させる駆動部１５ａが設けられている。駆動部１５ａとしては、例えばウォーム減速機等の比較的大きな減速比系の駆動機構を単体で又は組合せて用いることができる。
【００３０】
また、本実施形態では、第２格子１５近傍に、第１格子１４に対する第２格子１５の相対角度θ（ｘ−ｙ平面内の相対角度）を調整するための相対角度調整部１５ｂが設けられている。
【００３１】
上記マルチスリット１２、第１格子１４、第２格子１５は、例えば下記のように構成することができる。
Ｘ線源１１の焦点径；３００（μｍ）、管電圧：４０（ｋＶｐ）、付加フィルター：アルミ１．６（ｍｍ）
Ｘ線源１１の焦点からマルチスリット１２までの距離Ｄ１  ：  ２４０（ｍｍ）
マルチスリット１２から第１格子１４までの距離Ｄ３      ：１１０７（ｍｍ）
マルチスリット１２から第２格子１５までの距離Ｄ３＋Ｄ４：１３６４（ｍｍ）
マルチスリット１２のサイズ：１０（ｍｍ四方）、スリット周期：２２．８（μｍ）
第１格子１４のサイズ：５０（ｍｍ四方）、スリット周期：４．３（μｍ）
第２格子１５のサイズ：５０（ｍｍ四方）、スリット周期：５．３（μｍ）
【００３２】
Ｘ線検出器１６は、照射されたＸ線に応じて電気信号を生成する変換素子が２次元状に配置され、当該変換素子により生成された電気信号を画像信号として読み取る。
Ｘ線検出器１６の画素サイズは１０〜３００（μｍ）であり、さらに好ましくは５０〜２００（μｍ）である。
【００３３】
Ｘ線検出器１６は第２格子１５に当接するように保持部１７に位置を固定することが好ましい。第２格子１５とＸ線検出器１６間の距離が大きくなるほど、Ｘ線検出器１６により得られるモアレ縞画像がボケるからである。
Ｘ線検出器１６としては、ＦＰＤ（Flat Panel Detector）を用いることができる。ＦＰＤには、検出されたＸ線を光電変換素子を介して電気信号に変換する間接変換型、検出されたＸ線を直接的に電気信号に変換する直接変換型があるが、何れを用いてもよい。
【００３４】
間接変換型は、ＣｓＩやＧｄ２Ｏ２Ｓ等のシンチレータプレートの下に、光電変換素子がＴＦＴ（薄膜トランジスタ）とともに２次元状に配置されて各画素を構成する。Ｘ線検出器１６に入射したＸ線がシンチレータプレートに吸収されると、シンチレータプレートが発光する。この発光した光により、各光電変換素子に電荷が蓄積され、蓄積された電荷は画像信号として読み出される。
【００３５】
直接変換型は、アモルファスセレンの熱蒸着により、１００〜１０００（μｍ）の膜厚のアモルファスセレン膜がガラス上に形成され、２次元状に配置されたＴＦＴのアレイ上にアモルファスセレン膜と電極が蒸着される。アモルファスセレン膜がＸ線を吸収するとき、電子正孔対の形で物質内に電圧が遊離され、電極間の電圧信号がＴＦＴにより読み取られる。
なお、ＣＣＤ（Charge Coupled Device）、Ｘ線カメラ等の撮影手段をＸ線検出器１６として用いてもよい。
【００３６】
本体部１８は、図３に示すように、制御部１８１、操作部１８２、表示部１８３、通信部１８４、記憶部１８５等を備えて構成されている。
制御部１８１は、ＣＰＵ（Central Processing Unit）やＲＡＭ（Random Access Memory）等から構成され、記憶部１８５に記憶されているプログラムとの協働により撮影制御処理を行う。具体的に、制御部１８１は、Ｘ線源１１、駆動部１５ａ、駆動部１８ａ、Ｘ線検出器１６等の各部に接続されており、例えば、第２格子１５の移動、Ｘ線源１１からのＸ線照射のタイミングやＸ線照射条件、Ｘ線検出器１６による画像信号の読取タイミング等を制御する。
【００３７】
操作部１８２は、曝射スイッチ等を備え、これらの操作に応じた操作信号を生成して制御部１８１に出力する。
表示部１８３は制御部１８１の表示制御に従って、ディスプレイに操作画面やＸ線撮影装置１の動作状況等を表示する。
【００３８】
通信部１８４は通信インターフェイスを備え、ネットワーク上のコントローラー５と通信する。例えば、通信部１８４はＸ線検出器１６によって読み取られ、記憶部１８５に記憶されたモアレ縞画像をコントローラー５に送信する。
【００３９】
記憶部１８５は、制御部１８１により実行されるプログラム、プログラムの実行に必要なデータ、プログラムの実行による処理結果等を記憶している。また、記憶部１８５はＸ線検出器１６によって得られたモアレ縞画像を記憶する。
【００４０】
コントローラー５は、Ｘ線撮影装置１により得られた一連のモアレ縞画像を用いて被写体の再構成画像を生成し、生成した再構成画像に各種処理を施して表示する画像処理装置である。また、コントローラー５は、後述する格子相対位置推定処理を実行し、Ｘ線撮影装置１において取得されたモアレ縞画像から縞走査時の各撮影における第２格子１５の位置を推定する。
【００４１】
コントローラー５は、図４に示すように、制御部５１、操作部５２、表示部５３、通信部５４、記憶部５５を備えて構成されている。
制御部５１は、ＣＰＵ（Central Processing Unit）やＲＡＭ（Random Access Memory）等から構成され、記憶部５５に記憶されているプログラムとの協働により、後述する再構成処理、格子相対位置推定処理をはじめとする各種処理を実行する。制御部５１は、再構成手段、推定手段として機能する。
【００４２】
操作部５２は、カーソルキー、数字入力キー、及び各種機能キー等を備えたキーボードと、マウス等のポインティングデバイスを備えて構成され、キーボードで押下操作されたキーの押下信号とマウスによる操作信号とを、入力信号として制御部５１に出力する。表示部５３のディスプレイと一体に構成されたタッチパネルを備え、これらの操作に応じた操作信号を生成して制御部５１に出力する構成としてもよい。
【００４３】
表示部５３は、例えば、ＣＲＴ（Cathode Ray Tube）やＬＣＤ（Liquid Crystal Display）等のモニターを備えて構成されており、制御部５１の表示制御に従って、操作画面、生成された再構成画像等を表示する。
【００４４】
通信部５４は、通信インターフェイスを備え、ネットワーク上のＸ線撮影装置１やＸ線検出器１６と有線又は無線により通信する。例えば、通信部５４は、Ｘ線撮影装置１に撮影条件や制御信号を送信したり、Ｘ線撮影装置１又はＸ線検出器１６からモアレ縞画像を受信したりする。
【００４５】
記憶部５５は、制御部５１により実行されるプログラム、プログラムの実行に必要なデータ、プログラムの実行による処理結果等を記憶している。
また、記憶部５５は、Ｘ線撮影装置１で取得されたモアレ縞画像、モアレ縞画像に基づき生成された再構成画像等を患者情報に対応付けて記憶する。
【００４６】
［タルボ干渉計、タルボ・ロー干渉計の原理］
ここで、上記Ｘ線撮影装置１のタルボ・ロー干渉計によるＸ線撮影方法を説明する。
図５に示すように、Ｘ線源１１から照射されたＸ線が第１格子１４を透過すると、透過したＸ線がｚ方向に一定の間隔で像を結ぶ。この像を自己像といい、自己像が形成される現象をタルボ効果という。自己像を結ぶ位置に第２格子１５が自己像と概ね平行に配置され、第２格子１５を透過したＸ線によりモアレ縞画像（図５においてＭｏで示す）が得られる。即ち、第１格子１４は、周期パターンを形成し、第２格子１５は周期パターンをモアレ縞に変換する。Ｘ線源１１と第１格子１４との間に被写体（図５においてＨで示す）が存在すると、被写体によってＸ線の位相がずれるため、図５に示すようにモアレ縞画像上のモアレ縞は被写体の辺縁を境界に乱れる。このモアレ縞の乱れを、モアレ縞画像を処理することによって検出し、被写体像を画像化することができる。これがタルボ干渉計の原理である。
【００４７】
Ｘ線撮影装置１では、Ｘ線源１１と第１格子１４との間のＸ線源１１に近い位置に、マルチスリット１２が配置され、タルボ・ロー干渉計によるＸ線撮影が行われる。タルボ干渉計はＸ線源１１が理想的な点線源であることを前提としているが、実際の撮影にはある程度焦点径が大きい焦点が用いられるため、マルチスリット１２によってあたかも点線源が複数連なってＸ線が照射されているかのように多光源化する。これがタルボ・ロー干渉計によるＸ線撮影法であり、焦点径がある程度大きい場合にも、タルボ干渉計と同様のタルボ効果を得ることができる。
【００４８】
［縞走査による撮影］
本実施形態におけるＸ線撮影システムでは、Ｘ線撮影装置１において、縞走査による撮影を行う。縞走査とは、格子（マルチスリット１２、第１格子１４、第２格子１５）のうちの何れか１枚をスリット周期方向（ｘ方向）に動かしてＭ回（Ｍは正の整数、Ｍ＞２）の撮影を行い、Ｍ枚のモアレ縞画像を取得することをいう。理想的には、移動させる格子のスリット周期をｄ（μｍ）とすると、ｄ／Ｍ（μｍ）ずつ格子をスリット周期方向に動かして撮影を行うことを繰り返し、Ｍ枚のモアレ縞画像を取得する。
本実施形態においては、Ｘ線撮影装置１は、第２格子１５を移動させて縞走査による撮影を行う。
【００４９】
具体的には、Ｘ線撮影装置１において、オペレーターにより操作部１８２の曝射スイッチが押下されると、制御部１８１が撮影制御処理を実行することにより、以下の縞走査動作が行われる。なお、例えば、前回の縞走査終了後に、第２格子１５は駆動部１５ａにより原点位置に移動され、撮影制御処理の開始時には第２格子１５は略原点位置に位置しているものとする。
まず、第２格子１５が停止した状態でＸ線源１１によるＸ線の照射が開始される。Ｘ線検出器１６では前回の撮影により残存する不要な電荷を取り除くリセット後、Ｘ線照射のタイミングに合わせて電荷が蓄積され、Ｘ線の照射停止のタイミングに合わせて蓄積された電荷が画像信号として読み取られる。これが１回目の撮影である。１回目の撮影が終了するタイミングで駆動部１５ａによる第２格子１５の移動が開始され、ｄ／Ｍ（実際には、誤差を含む場合がある）（μｍ）移動すると停止され、次の撮影が行われる。このようにして、第２格子１５の移動と停止による（Ｍ−１）回分の撮影が繰り返され、第２格子１５が停止したときにＸ線の照射と画像信号の読み取りが行われる。第１回目の撮影を含めてＭ回の撮影が終了したときに、１枚の再構成画像を生成するのに必要な、モアレ縞の位相が異なる複数のモアレ縞画像を取得するための一連の撮影が終了する。
【００５０】
縞走査による撮影が終了すると、制御部１８１は、通信部１８４によりコントローラー５にモアレ縞画像を送信させる。通信部１８４からコントローラー５に対しては１回の撮影が終了する毎に１枚ずつ送信することとしてもよいし、縞走査によるＭ回分の全てのモアレ縞画像が得られた後、まとめて送信することとしてもよい。
ここで、制御部１８１は、通信部１８４によりコントローラー５にモアレ縞画像を送信する際、モアレ縞画像に対し、撮影条件（縞走査における撮影回数Ｍ、当該モアレ縞画像の撮影順番ｋ、移動させた第２格子１５のスリット周期ｄ、Ｘ線検出器１６の画素サイズＤ、相対角度θ、画像種別等）を付帯させて（例えば、ヘッダ情報等に書き込んで）コントローラー５に送信する。
【００５１】
なお、被写体の再構成画像を取得するには、被写体台１３に被写体を載置した状態での縞走査によるＸ線撮影（被写体有りの縞走査）と被写体台１３に被写体を載置しない状態での縞走査によるＸ線撮影（被写体無しの縞走査）が行われ、被写体有りの一連のモアレ縞画像及び被写体無しの一連のモアレ縞画像が生成される。これを本撮影という。また、モアレ縞画像のうち、被写体有りのモアレ縞画像を被写体モアレ縞画像と呼び、被写体無しのモアレ縞画像をＢＧモアレ縞画像と呼ぶ。また、縞走査により得られた一連のモアレ縞画像をモアレ縞画像セットと呼ぶ。
【００５２】
ここで、本撮影において、被写体有りの縞走査と被写体無しの縞走査とで第１格子１４と第２格子１５の相対角度θを変化させてしまうと、格子が均一に製造できていない場合には、格子の高周波のムラが再構成画像に残存し、被写体由来の信号を損なう恐れがある。そのため、本撮影の被写体有りの縞走査と被写体無しの縞走査は格子のモアレ縞の変化がほとんどない状態で撮影される。
【００５３】
［再構成処理］
コントローラー５においては、通信部５４により本体部１８からの被写体モアレ縞画像セット及びＢＧモアレ縞画像セットが受信されると、制御部５１と記憶部５５に記憶されているプログラムとの協働により、再構成処理が実行され、受信した被写体モアレ縞画像セット及びＢＧモアレ縞画像セットに基づいて微分位相画像、小角散乱画像、Ｘ線吸収画像等の再構成画像が生成される。以下に、再構成処理について述べる。
【００５４】
上述のように、縞走査では、理想的には格子のうちのひとつ（マルチスリット１２、第１格子１４、第２格子１５のうちのひとつ。ここでは第２格子１５。）をそのスリット周期ｄの１／Ｍ（Ｍ＞２）ずつスリット周期方向に移動させ、Ｍ枚の画像を得る。その際、移動させる第２格子１５以外の他の格子は静止しているものと仮定する。ただ、環境条件が変わると縞走査ごとに第２格子１５とその他の格子との位置関係は変化する。そこで、静止した他の格子はそれぞれ原点に静止しているものと考え、静止した格子と移動させる第２格子１５の原点位置との相対位置（スリット周期方向の相対位置）をχ0とおく。
【００５５】
縞走査の各撮影における第２格子１５の原点位置からの相対位置（スリット周期方向の相対位置）をχk（ｋ＝１,…,Ｍ。以下、第２格子１５（格子）の相対位置χkと呼ぶ。）とすると、（式１）のように表される。
【数１】
0006369206.tif 000002
ただし、ｄは縞走査時に移動させる第２格子１５のスリット周期で、αkはｋ枚目の撮影時の相対位置の誤差を表す。
【００５６】
上述のように、タルボ干渉計やタルボ・ロー干渉計では、第１格子１４の自己像と第２格子１５との間でモアレ縞を発生させる。一般に、第１格子１４の自己像と第２格子１５のスリット周期が一致している場合、Ｘ線検出器１６の画素(x,y)で検出されるモアレ縞強度Ik(x,y)は（式２）のように表される。
【数２】
0006369206.tif 000003
ここで、ＤはＸ線検出器１６の画素サイズ、d2は第１格子１４の自己像及び第２格子１５のスリット周期、θは第１格子１４のスリットがｙ軸に平行に配置されている場合の第２格子１５の第１格子１４に対する相対角度（rad）を表している。an（ｎは０以上の整数、ｎ＝±１、±２、・・・）は、モアレ縞強度Ik（x,y）を離散フーリエ変換した際の周期d2／nθ（μm）のモアレ縞成分の振幅の１／２を表す。a0は直流成分で、モアレ縞強度Ik（x,y）の平均値である。また、（式２）は被写体がない場合のモアレ縞の強度を表している。
【００５７】
実際のタルボ干渉計やタルボ・ロー干渉計では、適当な縞走査回数を設定することによりｎ＞２の高次の成分が無視できる。つまり、モアレ縞の強度分布を表す（式２）は（式３）のように近似することが出来る。
【数３】
0006369206.tif 000004
ここで、
【数４】
0006369206.tif 000005
【数５】
0006369206.tif 000006
δｋは、格子のうちの１枚を相対位置χkへ移動させたことによるモアレ縞の位相変化を表すものである。
再構成処理では、この（式３）を形成するパラメーターb0(x,y)、b1(x,y)、φ(x,y)をＭ枚のモアレ縞画像より抽出し、微分位相画像、Ｘ線吸収画像、小角散乱画像を再構成する。
【００５８】
ところで、（式３）は（式６）のように表すことも出来る。
【数６】
0006369206.tif 000007
【００５９】
δkが既知であれば、モアレ縞強度の実測値Ik´(x,y)とモアレ縞強度の理論値Ik(x,y)の差を最小にするようなc0(x,y)、c1(x,y)、c2(x,y)を求めることが出来る。（式６）のパラメーターc0(x,y)、c1(x,y)、c2(x,y)の導出方法は、谷田貝豊彦著、「応用工学光計測入門」、丸善株式会社、1988年発行、p.131-135に記載されており、ここではその結果のみを示す。
行列Ｃを（式７）のように定義すると、パラメーターc0(x,y)、c1(x,y)、c2(x,y)は（式８）のように求めることが出来る。
【数７】
0006369206.tif 000008
【数８】
0006369206.tif 000009
ただし、
【数９】
0006369206.tif 000010
【数１０】
0006369206.tif 000011
【００６０】
（式８）で求められたパラメーターc0(x,y)、c1(x,y)、c2(x,y)より、以下のようにパラメーターb0(x,y)、b1(x,y)、φ(x,y)を求めることができる。
【数１１】
0006369206.tif 000012
【数１２】
0006369206.tif 000013
【数１３】
0006369206.tif 000014
【００６１】
Ｘ線吸収画像IATT(x,y)、小角散乱画像ISAS(x,y)、微分位相画像IDPC(x,y)は、被写体モアレ縞画像セットより求めたパラメーターをb0,sample(x,y)、b1,sample(x,y)、φsample(x,y)、c1,sample(x,y)、c2,sample(x,y)とし、ＢＧモアレ縞画像セットより求めたパラメーターをb0,BG(x,y)、b1,BG(x,y)、φBG(x,y)、c1,BG(x,y)、c2,BG(x,y)として、これらを用いてそれぞれ（式１４）、（式１５）、（式１６（あるいは１７））により算出することができる。
【数１４】
0006369206.tif 000015
【数１５】
0006369206.tif 000016
微分位相画像IDPC(x,y)は、
【数１６】
0006369206.tif 000017
で計算することができるが、
【数１７】
0006369206.tif 000018
として計算してもよい。
【００６２】
即ち、再構成処理では、被写体モアレ縞画像及びＢＧモアレ縞画像のモアレ縞強度Ik´(x,y)及び縞走査のｋ回目の撮影における第２格子１５の相対位置χkに基づいてＸ線吸収画像IATT(x,y)、小角散乱画像ISAS(x,y)、微分位相画像IDPC(x,y)を算出する。そのため、制御部５１は、再構成処理を行うにあたり、予め縞走査のｋ回目の撮影における第２格子１５の相対位置χkをパラメーターとして設定しておく必要がある。
【００６３】
なお、本発明ではデータの扱いやすさから微分位相画像の各画素の値を（式１７）で定義し、微分位相画像の各画素のとりうる値を−１から１（または０から２）としている。微分位相画像の各画素の値を被写体の物理量に応じた値で定義する際には（式１７）にＸ線エネルギーやスリット周期、被写体位置などによって決まる係数を積算すればよい。
【００６４】
［相対位置の推定］
上述のように、再構成処理を行うには、縞走査時の各撮影における第２格子１５の相対位置χkをパラメーター（パラメーターχkとする）として設定しておく必要がある。しかし、縞走査時の実際の格子の相対位置χkと設定したパラメーターχkが一致していない場合には、再構成画像に撮影時のモアレ縞に起因するアーチファクトが発生してしまう。そのため、センサーを用いて縞走査時の相対位置χkを検出してこれをパラメーターχkとして設定することが提案されているが、高精度なセンサーを必要とし、コストがかかる。
【００６５】
上述の再構成処理には、（式１４）〜（式１６）に示すように、被写体モアレ縞画像の各画素の値をＢＧモアレ縞画像の同じ画素の値で除算する、若しくは減算するオフセット処理が含まれている。そのため、仮に、被写体モアレ縞画像の撮影時とＢＧモアレ縞画像の撮影時の被写体の有無および撮影条件が全く同じであれば縞走査時の相対位置χkに誤差があってもこのオフセット処理により減算或いは除算されるため、オフセット処理後の再構成画像ではその誤差は検出できないことになる。そこで、本願発明者は、意図的に撮影条件、特にモアレ縞の周期が異なる条件で撮影したモアレ縞画像セット間で再構成処理を行い、アーチファクトが明瞭に残存した再構成画像を作成し、このアーチファクトが最も小さくなるような再構成処理時のパラメーターχkを求めることで、縞走査時の各撮影における格子の相対位置χkを推定することが可能なことを見出した。
以下、図６を参照して、縞走査時の第２格子１５の相対位置χkを推定するための撮影及び格子相対位置推定処理の流れについて説明する。
【００６６】
まず、Ｘ線撮影装置１において、或るモアレ縞の周期で縞走査による撮影を行い、１セット分（Ｍ枚）のモアレ縞画像セット１を取得する（ステップＳ１）。取得されたモアレ縞画像セット１の画像データは、通信部１８４によりコントローラー５に送信される。
【００６７】
次いで、Ｘ線撮影装置１において、相対角度調整部１５ｂを調整して第１格子１４と第２格子１５との相対角度θを変更することにより、ステップＳ１とは異なる相対角度θで縞走査による撮影を行い、１セット分（Ｍ枚）のモアレ縞画像セット２を取得する（ステップＳ２）。取得されたモアレ縞画像セット２の画像データは、通信部１８４によりコントローラー５に送信される。
【００６８】
ここで、モアレ縞の周期は、（式３）に示すように、第１格子１４と第２格子１５の相対角度θを変化させることで変更することができる。ステップＳ１とステップＳ２の撮影時には、モアレ縞の周期以外は変えないことが好ましい。また、ステップＳ１、ステップＳ２の撮影時には被写体Ｈを置いても置かなくてもどちらでもよいが、被写体Ｈを置かないほうが望ましい。
【００６９】
図７に、シミュレーションにより作成したモアレ縞の周期の異なるモアレ縞画像セット１とモアレ縞画像セット２の例を示す。
図７に示すモアレ縞画像セット１、２は、以下の条件でシミュレーションを行って計算した画像である。
マルチスリット１２と第１格子１４間の距離Ｄ３：１１０７（ｍｍ）
マルチスリット１２と第１格子１４間の距離Ｄ３＋Ｄ４：１３６４（ｍｍ）
マルチスリット１２のスリット周期：２２．８（μｍ）、スリット幅：６（μｍ）
第１格子１４のスリット周期：４．３（μｍ）、スリット幅：２．１５（μｍ）
第２格子１５のスリット周期ｄ：５．３（μｍ）、スリット幅：２．６５（μｍ）
画素サイズ：８５（μｍ）
第１格子１４と第２格子１５の相対角度θ：０．０２°、０．００５°
画像サイズ：２００画素×２００画素
第２格子１５を動かして４回の撮影を行うこととし、第２格子１５の相対位置χkは、モアレ縞画像セット１、２で共通で、χ1＝０、χ2＝（０．２５＋０．００５）ｄ、χ3＝（０．５−０．０１）ｄ、χ4＝（０．７５＋０．００５）ｄとした。
【００７０】
次いで、コントローラー５において、制御部５１によりモアレ縞画像セット１、２を用いて格子相対位置推定処理を実行し、縞走査時の各撮影における第２格子１５の相対位置χkを推定する。具体的には、以下のステップＳ３〜ステップＳ８の処理を実行する。
【００７１】
まず、制御部５１は、縞走査において移動させた格子（第２格子１５）の相対位置χkを示すパラメーターχk（再構成処理において相対位置χkに設定するパラメーター）の初期値を設定する（ステップＳ３）。
【００７２】
ここで、第２格子１５の相対位置χk（ｋ＝１，・・・，Ｍ）は、（式１）に示すとおりであり、これは、（式５）のδkの分子である。１枚目の相対位置χ1はα1となるがこれを定数項とみなし、
【数１８】
0006369206.tif 000019
【数１９】
0006369206.tif 000020
【数２０】
0006369206.tif 000021
として、縞走査における各撮影時の第２格子１５の相対位置χkを１枚目の撮影時の格子位置χ1に対する相対位置と定義し直しても問題ない。
【００７３】
縞走査における各撮影時の第２格子１５の相対位置χkを１枚目の撮影時の格子位置に対する相対位置として定義した場合、ｋ＝１の場合のパラメーターχ1は０とすることができ、Ｍ枚撮影の場合、(Ｍ−１)回の撮影のパラメーターを設定すればよいということになる。
【００７４】
［表１］に、Ｍ＝４の場合の相対位置のパラメーターχkの一例を示す。
【表１】
0006369206.tif 000022
ステップＳ１、ステップＳ２に対応する撮影が同様の環境条件で撮影されている場合、例えば、モアレ縞画像セット１、２が続けて撮影された場合、縞走査時の格子移動による相対位置の変化は２セットの撮影で共通とみなすことが出来る。このような場合、求めるべきパラメーターはモアレ縞画像セット１、２で共通のパラメーターセット１とすることができる。また、ステップＳ１、ステップＳ２に対応する撮影が同様の環境条件で撮影されていない場合、例えば、撮影時間が離れて温度条件が変化したような場合には、モアレ縞画像セット１およびモアレ縞画像セット２は、２〜４回目の撮影時の相対位置が異なるパラメーターセット２とする必要がある。以下の説明では、［表１］のパラメーターセット１の場合を例にとり説明する。
【００７５】
ところで、ｋ＝１の場合の格子の相対位置χ1を０とすることにより、（式２０）に示すようにφ(x,y)に1回目の撮影時の第２格子１５の位置の誤差α1が反映されることとなる。誤差α1がモアレ縞画像セット１、２で同じ場合には（式１６）で算出される微分位相画像にα1は影響しないが、誤差α1がモアレ縞画像セット１、２で異なる場合であっても、後述する微分位相画像への補正処理によりこの誤差は除去あるいは減弱されるため、考慮する必要はない。
【００７６】
Ｘ線撮影装置１においては、縞走査時に１／Ｍ周期分ずつ第２格子１５を移動させる設定になっているため、Ｍ＝４の場合、ステップＳ３においては、パラメーターχkの初期値はモアレ縞画像セット１、２で共通で［表２］のように与えるのでよい。
【表２】
0006369206.tif 000023
【００７７】
次いで、制御部５１は、モアレ縞画像セット１、２のそれぞれに設定されたパラメーターχkに基づいて再構成処理を実行し、再構成画像（微分位相画像、Ｘ線吸収画像、小角散乱画像）を生成する（ステップＳ４）。
ステップＳ４においては、ステップＳ３で設定したパラメーターχkと、第２格子１５のスリット周期ｄを（式１９）に代入してδkを算出し、モアレ縞画像セット１又は２の一方を被写体モアレ縞画像セット、他方をＢＧモアレ縞画像セットとして（式１４）、（式１５）、（式１７）により微分位相画像、Ｘ線吸収画像、小角散乱画像を生成する。
【００７８】
図８は、パラメーターχkに［表２］のパラメーターセット１に示す初期値を設定し、図７のモアレ縞画像セット１、２を用いて再構成処理することにより得られた微分位相画像、Ｘ線吸収画像、小角散乱画像である。図８に示すように、Ｘ線吸収画像、小角散乱画像では、ｙ方向に縞状のアーチファクトが見える。一方、微分位相画像では細かい縞状アーチファクトが見られないがｙ方向に沿って信号が徐々に変化しているのがわかる。
【００７９】
微分位相画像におけるｙ方向に沿ったアーチファクトは、特開２０１２−１７０６１８号公報（特許文献２）に記載のように、第１格子１４と第２格子１５の相対角度θが２つのモアレ縞画像セットで異なることに起因するものである。この相対角度θが異なることに起因する微分位相画像におけるアーチファクトは、特許文献２に記載されているように、画素(x,y)に対してf(x,y)=ax+by+cで表される位相誤差が加わったものである。そこで、微分位相画像に対しては、微分位相画像IDPC(x,y)からf(x,y)=ax+by+cを減算して、上記アーチファクトを補正する処理を行う。
【００８０】
位相誤差を表す式のa、b、cの値は、微分位相画像において、少なくとも３点の測定点を設定し、その位置座標（x1、y1）〜（x3、y3）及び信号値I１〜I3を取得し、以下の（式２１）に各値を代入することにより求めることができる。
【数２１】
0006369206.tif 000024
ここで、相対角度θの違いによるアーチファクトは、図８に示すように、主にｙ方向に生じている。そこで、微分位相画像IDPC(x,y)から減算する式をf(x,y)=by+c(x)で近似しても問題ない。相対角度θが既知であればｙ方向の傾きbはあらかじめ計算することが可能なため、微分位相画像IDPC(x,y)からf(x,y)=byを減算し、その後、補正した微分位相画像IDPC(x,y)のｘ座標ごとにｙ方向の画素値の平均c(x)を計算し減算するとしてもよい。
【００８１】
また、相対角度θの大きさによっては、微分位相画像IDPC(x,y)の面内にほぼｘ軸に平行に不連続が発生することがある。これは、微分位相画像IDPC(x,y)が−１から１の範囲の値しかとれないため、この範囲外の値では値が−１から１に折り返されるためである。このような場合には、ｙ方向に連続になるように微分位相画像IDPC(x,y)に２×整数の値を加算しておいてから、上記の補正処理を行えばよい。あるいは、先に、f(x,y) = byを減算し、その後、ｙ方向に連続になるように２×整数の値を加算し、ｘ座標ごとにｙ方向の画素値の平均c(x)を計算し減算するとしてもよい。
【００８２】
図９に、上記の補正処理を施して算出された微分位相画像を示す。図９に示すように、相対角度θが異なることに起因するアーチファクトの補正処理を施すことにより、微分位相画像にも縞状アーチファクトが発生したことがわかる。ステップＳ４においては、図８に示すＸ線吸収画像、小角散乱画像、図９に示す微分位相画像を再構成画像として生成する。
【００８３】
次いで、制御部５１は、再構成画像におけるモアレ縞の振幅の大きさを示す特徴量Ｅ（（式２２）参照）を算出する（ステップＳ５）。
【００８４】
図１０に、縞走査時の各撮影における第２格子１５の相対位置をχ2＝（０．２５＋０．００５）ｄ、χ3＝（０．５−０．０１）ｄ、χ4＝（０．７５＋０．００５）ｄと設定してモアレ縞画像セット１、２を取得し、パラメーターを縞走査時の第２格子１５の相対位置と同じχ2＝（０．２５＋０．００５）ｄ、χ3＝（０．５−０．０１）ｄ、χ4＝（０．７５＋０．００５）ｄと設定してステップＳ４の再構成処理を行ったＸ線吸収画像、小角散乱画像、微分位相画像である。
【００８５】
図１０に示す画像から、実測（縞走査）の際の第２格子１５の相対位置χkを正しく見積もってパラメーターχkとして設定して再構成処理が行われていれば、図８、図９に示すような再構成画像上の縞状のアーチファクトが発生しないことがわかる。よって、図８、９に見られる縞状のアーチファクトが十分小さくなるような、即ち、再構成画像に残存するモアレ縞の振幅が十分小さくなるような再構成処理時のパラメーターχkを求めれば、それが略撮影時の格子の相対位置χkであるということになる。つまり、Ｘ線吸収画像の画素値IATT(x,y)、小角散乱画像の画素値ISAS(x,y)、微分位相画像の画素値IDPC(x,y)が、評価を行う画像データ範囲内で、それぞれ一定値で与えられる目標値IATT0、ISAS0、IDPC0になるようなパラメーターχkを求めれば、これを撮影時の相対位置χkとして推定することができる。具体的には、以下の（式２２）に示す、再構成画像のモアレ縞の振幅の大きさを表す特徴量Ｅが最も小さくなるときのパラメーターχkを求めることで、縞走査時の各撮影における第２格子１５の相対位置χkを推定することができる。なお、特徴量Ｅの算出に用いる画像データ範囲は、Ｘ線撮影装置１の有効視野全てを用いてもよいし、視野の一部、例えば、中央部のみを用いてもよい。
【数２２】
0006369206.tif 000025
【００８６】
ここで、WATT、WSAS、WDPCは特徴量Ｅに対するＸ線吸収画像、小角散乱画像、微分位相画像の重み付けを表す。また、IATT0、ISAS0には、第1格子１４と第２格子１５の間の相対角度θやモアレ縞画像セット１、２の線量条件に応じた定数を、IDPC0には、第１格子１４と第２格子１５の間の相対角度θを考慮する必要がないため（相対角度θに起因するアーチファクトは補正されているため）０を設定すればよい。あるいは、IATT0、ISAS0、IDPC0には、それぞれ評価に用いる画像データ範囲内のIATT(x,y)、ISAS(x,y)、IDPC(x,y)の平均値を設定してもよい。
ステップＳ５においては、（式２２）により特徴量Ｅを算出する。
【００８７】
次いで、制御部５１は、特徴量Ｅが予め定められた基準値以下であるか否かを判断する（ステップＳ６）。予め定められた基準値とは、特徴量Ｅが取り得る最小値、又は十分小さいとして予め定められた値である。
特徴量Ｅが予め定められた基準値以下ではないと判断した場合（ステップＳ６；ＮＯ）、制御部５１は、再構成処理時のパラメーターχ2〜χMを変化させ（ステップＳ７）、ステップＳ４に戻り、再度再構成画像の生成を行う。
特徴量Ｅが予め定められた基準値以下であると判断した場合（ステップＳ６；ＹＥＳ）、制御部５１は、ステップＳ４において再構成処理で用いたパラメーターχ2〜χMを縞走査時の各撮影における第２格子１５の相対位置χ2〜χMであると推定する（ステップＳ８）。
【００８８】
なお、特徴量Ｅの最小化については、共役勾配法などの一般的な最適化手法を用いて最小値を求めることとしてもよい。
【００８９】
図７に示すモアレ縞画像セット１、２に基づいて、（式２２）のWATT＝WSAS＝１、WDPC＝１として共役勾配法を用いて特徴量Ｅを最小化するパラメーターχ2、χ3、χ4を計算したところ、χ2＝（０．２５＋０．００５＋０．０００２)ｄ、χ3＝（０．５−０．０１＋０．００３)ｄ、χ4＝（０．７５＋０．００５＋０．００００)ｄという値が得られた。実際の縞走査時の相対位置χ2、χ3、χ4に対する誤差が周期の０．１％以下で求められており、この手法を用いることで、縞走査時の各撮影の第２格子１５の相対位置をよく推定できているのがわかる。ここで、吸収画像、小角散乱画像にはモアレ縞画像セット１とモアレ縞画像セット２で差異を与えないシミュレーションを行ったため、IATT0＝ISAS0＝１、IDPC0＝０と設定した。
【００９０】
このように、Ｘ線撮影システムにおいては、図６に示す撮影〜格子相対位置推定処理を実行することで、縞走査時の各撮影における第２格子１５の相対位置を再構成画像から制度良く推定することができるので、縞走査時に移動させた格子の相対位置を検出するための特殊なセンサーをＸ線撮影装置に付加することなく、安価に格子の相対位置の誤差に起因する画像劣化を防止することができる。
【００９１】
ここまでは、Ｘ線吸収画像、小角散乱画像、微分位相画像の３画像をすべて用いて格子の相対位置を推定する手法を述べたが、図８、９に示すように、再構成画像上に残存する縞状アーチファクトの周期が再構成画像により異なる。Ｘ線吸収画像は、微分位相画像、小角散乱画像の縞状アーチファクトの倍の周期の縞状アーチファクトが残存している。これは、モアレ縞の平均値成分と振幅成分に由来する物理量への格子の相対位置誤差の影響が異なることを示している。しかし、実際、小角散乱画像は、（式１５）に示すように、モアレ縞の平均値成分、振幅成分の両方が画像に影響している。よって、格子の相対位置の推定時に、小角散乱画像のみを用いてもよい。つまり、（式２２）のWATT＝WDPC＝０、WSAS＝１として特徴量Ｅを計算してもよい。このようにして、図７に示すモアレ縞画像セット１、２に基づき、特徴量Ｅを最小にするパラメーターχ2、χ3、χ4を共役勾配法を用いて計算すると、それぞれχ2＝（０．２５＋０．００５＋０．０００２)ｄ、χ3＝（０．５０−０．０１＋０．０００７)ｄ、χ4＝（０．７５＋０．００５＋０．０００４）ｄとなり、縞走査時の相対位置χ2、χ3、χ4をよく再現出来ていることがわかる。この手法によれば、格子の相対位置の推定に係る処理時間を短縮することができる。
【００９２】
このように、モアレ縞周期の異なるモアレ縞画像セットに基づく再構成画像上に残存する縞状アーチファクトの振幅を最小にするようなパラメーターχkを求めることにより、実際の縞走査時の第２格子１５の相対位置を精度よく把握することできる。よって、詳細を後述するが、求めたパラメーターχkを用いて、相対位置の誤差（αk−α1）が０になるようにＸ線撮影システムにおける格子の移動量の調整を行ったり、あるいは、求めたパラメーターχkを用いて再構成処理を行ったりすることで、本撮影を精度よく行うことができるようになる。
【００９３】
また、ここまでは、モアレ縞画像セット１、２でパラメーターχkを同じとした場合について説明を行ったが、モアレ縞画像セット１、２でパラメーターが異なる場合、例えば、［表１］のパラメーターセット２の場合でも、処理は同様である。ただし、実際に被写体の再構成画像を生成する際の再構成処理においては、被写体の信号を得るための本撮影と環境がより近い状態で撮影された方のモアレ縞画像セット１又は２のパラメーターを用いて再構成画像の生成を行う。例えば、モアレ縞画像セット１を朝撮影し、モアレ縞画像セット２を本撮影直前に撮影した場合、本撮影はモアレ縞画像セット２と環境が近い状態で行われるため、モアレ縞画像セット２に対応するχkを用いて再構成画像の生成を行う。なお、モアレ縞画像セット１、２でパラメーターが異なる場合には、モアレ縞の周期が広いと格子の相対位置の算出精度が低下するため、モアレ縞画像セット１、２において、格子の相対位置を算出する際に評価に用いる画像領域内で、モアレ縞の位相がπ程度（半周期程度）変化する領域を含むことが望ましい。また、モアレ縞画像セット１、２でパラメーターが異なる場合には、モアレ縞画像セット１、２のうち、モアレ縞の周期の小さなモアレ縞画像セットに対応するχkのほうが精度よく求められる。そのため、例えば、本撮影がモアレ縞画像セット２と環境が近い状態となる場合には、モアレ縞画像セット２のモアレ縞の周期のほうがモアレ縞画像セット１に比べて小さくなるように第１格子１４と第２格子１５の相対角度θを調整するのが望ましい。
【００９４】
また、ここでは、モアレ縞画像セット間のモアレ縞の周期を変化させるのに、第１格子１４と第２格子１５の相対角度θを変化させたが、何れかの格子をｚ方向に移動させることにより、ｘ方向に平行なモアレ縞を発現させて同様に格子の相対位置χkを推定することも可能である。
【００９５】
［推定した相対位置のフィードバック］
以下、推定した相対位置χkを被写体の再構成画像の生成時にフィードバックする手法について説明する。
【００９６】
（Ｘ線撮影装置１での格子の移動量を調整）
本撮影前に、図６を用いて説明した格子の相対位置χkの推定（２つのモアレ縞画像セットの撮影及び格子相対位置推定処理）を行い、コントローラー５の制御部５１は、通信部５４により推定された相対位置χkをＸ線撮影装置１に送信する。推定された相対位置χkを通信部１８４により受信すると、Ｘ線撮影装置１の制御部１８１は、推定されたχkから第２格子１５の相対位置の誤差（αk−α1）を算出し、この（αk−α1）が０となるよう、縞走査時に第２格子１５を移動させる距離（移動量）を調整する。調整後、Ｘ線撮影装置１において、駆動部１５ａにより調整した移動量で第２格子を移動させて本撮影を行い、コントローラー５において、得られた被写体モアレ縞画像セットとＢＧモアレ縞画像に基づく再構成処理を行う。このとき、コントローラー５の制御部５１は、χk＝ｄ（ｋ−１）／Ｍとして再構成画像の生成を行う。これにより、縞走査時の各撮影における第２格子１５の相対位置χkと再構成時に設定するパラメーターχkが一致しないことに起因するアーチファクトを抑制することができる。
【００９７】
環境変化による（αk−α1）が小さいと見積もられるのであれば、Ｘ線撮影装置１において所定のタイミングで実行される装置調整時に、或いは、Ｘ線撮影装置１の電源を入れる毎に、図６に示す相対位置χkの推定〜推定結果に基づく第２格子１５を移動させる距離の調整までを行って、その値（調整した移動量）を記憶部１８５に記憶しておき、本撮影時には、記憶部１８５に記憶されている移動量に基づいて、駆動部１５ａにより第２格子１５を移動させて撮影を行うこととしてもよい。
【００９８】
（再構成処理時の格子位置の設定に反映）
本撮影前に、図６を用いて説明した格子の相対位置χkの推定（２つのモアレ縞画像セットの撮影及び格子相対位置推定処理）を行って、コントローラー５の記憶部５５に推定された相対位置χkを記憶しておき、本撮影の被写体モアレ縞画像セットとＢＧモアレ縞画像セットに基づく再構成処理時に、制御部５１は、記憶部５５に記憶されている相対位置χkをパラメーターとして設定して再構成画像の生成を行う。これにより、縞走査時の各撮影における第２格子１５の相対位置χkと再構成時にパラメーターとして設定するχkが一致しないことに起因するアーチファクトを抑制することができる。
【００９９】
ところで、Ｘ線撮影装置１を立ち上げた後、装置による発熱により装置内部で温度変化が発生する。これにより、第２格子１５の相対位置χkが変化する場合もある。そこで、本撮影の被写体モアレ縞画像セット、ＢＧモアレ縞画像セットの撮影とともに、ＢＧモアレ縞画像セットに対して第１格子１４と第２格子１５の相対角度θを変化させてＢＧモアレ縞画像セット（ＢＧモアレ縞画像セット２と呼ぶ）の撮影を行う。そして、コントローラー５において、制御部５１は、本撮影用のＢＧモアレ縞画像セットとＢＧモアレ縞画像セット２に基づいて、共通のパラメーター（［表１］のパラメーター１）を用いて格子相対位置推定処理を実行し、推定された相対位置χkを本撮影の被写体モアレ縞画像セットとＢＧモアレ縞画像セットによる再構成処理時のパラメーターχkに反映させて被写体の再構成画像を生成する。この際、相対角度調整部１５ｂにより、本撮影時のモアレ縞画像セットに対してＢＧモアレ縞画像セット２のモアレ縞の周期が小さくなるように相対角度θを変化させることが望ましい。
【０１００】
例えば、本撮影ではモアレ縞の周期が最も広がるよう第１格子１４と第２格子１５の相対角度θを調整した状態で撮影を行う。これに対して、ＢＧモアレ縞画像セット２では第１格子１４と第２格子１５の相対角度θを本撮影に対して角度を０．０１〜０．０３°程度変化させて撮影を行えばよい。撮影視野内の４５ｍｍ角の画像データ範囲を特徴量Ｅの算出に用いるとし、第２格子１５のスリット周期ｄを５．３（μｍ）とすると、ＢＧモアレ縞画像セットのモアレ縞の周期が４５ｍｍに対して十分大きい場合には、ＢＧモアレ縞画像セット２は４５ｍｍ角内に１〜５周期分のモアレ縞が発現した状態で撮影することとなる。モアレ縞の位相は、０〜２πで変化するため、特徴量Ｅの算出に用いる画像データ範囲内に1周期分以上のモアレ縞を発現させることにより、とりうるあらゆるモアレ縞の位相に対して再構成画像の誤差が最小となるパラメーターχkを算出することができる。
【０１０１】
なお、ＢＧモアレ縞画像セット２の撮影は、必ずしも本撮影の撮影とあわせて行う必要はなく、予めＢＧモアレ縞画像セット２の撮影を行ってコントローラー５の記憶部５５に記憶しておいてもよい。そして、本撮影終了後、本撮影のＢＧモアレ縞画像セットと記憶部５５に記憶されているＢＧモアレ縞画像セット２に基づいて、図６の格子相対位置推定処理により［表１］のパラメーターセット２（本撮影のＢＧモアレ縞画像セットと記憶部５５に記憶されているＢＧモアレ縞画像セット２で異なるパラメーターχk）を算出し、本撮影のＢＧモアレ縞画像セットについて算出したパラメーターχkを格子の相対位置χkとして本撮影のＢＧモアレ縞画像セットの再構成処理を行う。なお、ＢＧモアレ縞画像セット２は、本撮影に対して、線量を低くして撮影しても構わない。
【０１０２】
以上、本発明に係るＸ線撮影システムの好ましい実施形態について説明したが、上述した本実施形態における記述は、本発明に係る好適な一例であり、これに限定されるものではない。
【０１０３】
例えば、上記実施形態では、撮影時に第２格子１５をマルチスリット１２及び第１格子１４に対して相対移動させる方式のタルボ・ロー干渉計を用いたＸ線撮影装置を例にとり説明したが、本発明は、マルチスリット１２又は第１格子１４の何れかを移動させる方式のタルボ・ロー干渉計を用いたＸ線撮影装置に適用してもよい。また、本発明は、第１格子１４又は第２格子１５の何れかを他の格子に対して移動させる方式のタルボ干渉計を用いたＸ線撮影装置に適用してもよい。
【０１０４】
また、上記実施形態においては、３種の再構成画像を生成する場合について説明したが、微分位相画像又は小角散乱画像の少なくとも一つ生成するＸ線撮影システムであれば本発明を適用することが可能ある。
【０１０５】
その他、Ｘ線撮影システムを構成する各装置の細部構成及び細部動作に関しても、発明の趣旨を逸脱することのない範囲で適宜変更可能である。
【符号の説明】
【０１０６】
１  Ｘ線撮影装置
１１  Ｘ線源
１２  マルチスリット
１３  被写体台
１４  第１格子
１５  第２格子
１５ａ  駆動部
１５ｂ  相対角度調整部
１６  Ｘ線検出器
１７  保持部
１７ａ  緩衝部材
１８  本体部
１８１  制御部
１８２  操作部
１８３  表示部
１８４  通信部
１８５  記憶部
１８ａ  駆動部
５  コントローラー
５１  制御部
５２  操作部
５３  表示部
５４  通信部
５５  記憶部
【図１】
0006369206.tif 000026
【図２】
0006369206.tif 000027
【図３】
0006369206.tif 000028
【図４】
0006369206.tif 000029
【図５】
0006369206.tif 000030
【図６】
0006369206.tif 000031
【図７】
0006369206.tif 000032
【図８】
0006369206.tif 000033
【図９】
0006369206.tif 000034
【図１０】
0006369206.tif 000035