(57)【要約】      （修正有）【課題】Ｘ線タルボット干渉計により生成された単一フリンジ模様からの復調において軸外フリンジを抑制する方法を提供する。【解決手段】撮像画像は、交差項位相成分を含む位相変調フリンジ模様を含む。交差項位相成分を減衰するために少なくとも第１の非線形変換（ｈ）のパラメータを決定する。第１の非線形変換のパラメータは、撮影システムが交差項位相成分により導入された復調画像におけるアーチファクトを減少するために決定される。交差項位相成分を減衰するために、決定されたパラメータによる少なくとも第１の非線形変換を撮像画像に適用すること４２０により、撮像画像を復調する４３０。【選択図】図４
2016087473.tif 000002
【特許請求の範囲】
【請求項１】
交差項位相成分を含む位相変調フリンジ模様を含む撮像画像を復調する方法であって、
前記交差項位相成分を減衰するために少なくとも第１の非線形変換のパラメータを決定する工程であって、前記第１の非線形変換の前記パラメータが、前記交差項位相成分により導入された復調画像におけるアーチファクトを減少するために決定される工程と、
前記交差項位相成分を減衰するために、前記決定されたパラメータによる少なくとも前記第１の非線形変換を前記撮像画像に適用することにより、前記撮像画像を復調する工程と、
を有することを特徴とする方法。
【請求項２】
前記第１の非線形変換のパラメータを前記決定する工程は、
前記撮像画像を撮影した撮影システムを使用して撮像された校正画像に対して、前記交差項位相成分の大きさと軸上位相成分の大きさとの比を決定する工程と、
前記決定された比を使用して前記第１の非線形変換のパラメータを決定する工程と、を有することを特徴とする請求項１に記載の方法。
【請求項３】
パラメータを前記決定する工程は、
更なる復調のために前記交差項位相成分の影響を軽減するようにパラメータ化非線形変換を実行する工程と、
交差項位相成分を含む位相変調フリンジ模様の前記撮像画像を撮影した撮影システムを用いて撮像された複数の参照画像を受信する工程と、
参照復調画像を生成するためにマルチショット復調方法を使用して前記参照復調画像に対して第１の復調を行う工程と、
初期の非線形変換を形成するように前記実行されたパラメータ化非線形変換に対する初期のパラメータを識別する工程と、
前記初期の非線形変換に従って変換された前記参照復調画像に適用されたシングルショット復調方法を使用して前記参照復調画像のうちの１つに対して第２の復調を行うことにより、候補復調画像を形成する工程と、
前記非線形変換のパラメータを調整するために、前記候補復調画像を前記参照復調画像と比較する工程と、
を更に有することを特徴とする請求項１に記載の方法。
【請求項４】
前記候補復調画像と前記参照復調画像との間の誤差を減少するように変更された候補復調画像を形成するために、前記調整されたパラメータを用いた前記第２の復調を繰り返す工程を更に有することを特徴とする請求項３に記載の方法。
【請求項５】
前記第２の復調は、前記誤差が所定の値を下回るまで繰り返されることを特徴とする請求項４に記載の方法。
【請求項６】
前記調整されたパラメータは、前記撮影システムにより撮像された更なる画像の復調における交差項位相成分の影響を軽減するために使用されることを特徴とする請求項３記載の方法。
【請求項７】
前記撮像画像を復調する工程は、
前記減衰された交差項位相成分を含む第１の画像を生成するために前記第１の非線形変換を前記撮像画像に適用する工程と、
フーリエ領域において前記第１の画像のｙ項及びｘ項のそれぞれを前記第１の画像からフィルタリングすることにより、ｘ位相画像及び位相画像を生成する工程と、を更に含むことを特徴とする請求項１に記載の方法。
【請求項８】
前記第１の非線形変換を前記撮像画像に適用することにより導入された高調波を減衰するために構成された第２の非線形変換を実行する工程であって、前記第２の非線形変換が前記第１の非線形変換を反転することにより実行される工程と、
前記第１の非線形変換の動作により導入された高調波を減衰するために、前記第２の非線形変換をｘ位相画像及びｙ位相画像の各々に適用する工程と、を更に有することを特徴とする請求項７に記載の方法。
【請求項９】
前記フィルタリングは、
前記第１の画像を前記フーリエ領域に変換する工程と、
少なくとも１つのマスク画像を前記変換された第１の画像に適用する工程と、
前記マスク画像を空間領域に変換する工程と、を有することを特徴とする請求項７に記載の方法。
【請求項１０】
前記適用する工程は、前記変換された第１の画像に、前記ｘ位相画像を生成するために第１のマスクを乗算し、ｙ位相画像を生成するために第２のマスクを乗算する工程を有することを特徴とする請求項９に記載の方法。
【請求項１１】
前記少なくとも１つのマスク画像は、画像から所定の周波数スペクトルを除去するように構成されることを特徴とする請求項１０に記載の方法。
【請求項１２】
分離されたフリンジ模様に対して、前記位相画像から低周波数スペクトルを減算する工程を更に有することを特徴とする請求項８に記載の方法。
【請求項１３】
前記第１の非線形変換は、対数関数による変換を含むことを特徴とする請求項１に記載の方法。
【請求項１４】
前記第２の非線形変換は、指数関数による変換を含むことを特徴とする請求項８に記載の方法。
【請求項１５】
前記第１の非線形変換及び第２の非線形変換は、冪関数による変換を含むことを特徴とする請求項８に記載の方法。
【請求項１６】
交差項位相成分を含む位相変調フリンジ模様を含む撮像画像を復調する方法の各工程をコンピュータに実行させることが可能なコンピュータプログラムを記憶したコンピュータ可読の記憶媒体であって、前記コンピュータプログラムは、
前記交差項位相成分を減衰するために少なくとも第１の非線形変換のパラメータを決定するコードであって、前記第１の非線形変換の前記パラメータが、前記交差項位相成分により導入された復調画像におけるアーチファクトを減少するために決定されるコードと、
前記交差項位相成分を減衰するために、前記決定されたパラメータによる少なくとも前記第１の非線形変換を前記撮像画像に適用することにより、前記撮像画像を復調するコードと、
を含むことを特徴とする記憶媒体。
【請求項１７】
撮影システムであって、
少なくとも１つの２次元格子を有し、且つ前記少なくとも１つの２次元格子がある結果として得られる交差項位相成分を含む位相変調フリンジ模様を含む画像を撮像するように構成されたＸ線タルボット撮影装置と、
前記Ｘ線タルボット撮影装置に結合され、且つ少なくともプロセッサ及び前記プロセッサに結合されたメモリを含むコンピュータプロセッサ装置であって、前記メモリが、前記画像を復調するように前記プロセッサにより実行可能なコンピュータプログラムを記憶する、コンピュータプロセッサ装置と、を備え、前記コンピュータプログラムが、
前記交差項位相成分を減衰するために少なくとも第１の非線形変換のパラメータを決定するコードであって、前記第１の非線形変換の前記パラメータが、前記撮影システムが前記交差項位相成分により導入された復調画像におけるアーチファクトを減少するために決定されるコードと、
前記交差項位相成分を減衰するために、前記決定されたパラメータによる少なくとも前記第１の非線形変換を前記撮像された撮像画像に適用することにより、前記撮像画像を復調するコードと、
を含むことを特徴とする撮影システム。
【請求項１８】
前記撮影システムは、Ｘ線干渉計システムを含むことを特徴とする請求項１７に記載の撮影システム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、一般に、画像内で２つ以上のフリンジ周波数又はフリンジ配向が使用されるフリンジ解析のための空間復調の分野に関する。特に、望ましくない軸外周波数が信号内容を含む所望の周波数と共に現れるＸ線モアレ装置と関連付けられたフリンジの２次元集合の復調を支援する方法、記憶媒体、および撮影システムを説明する。
【背景技術】
【０００２】
現在、Ｘ線タルボットモアレインタフェログラムの復調に対して、主に２つの方法がある。第１の方法は、（時間的）位相シフトアルゴリズム（マルチショットアルゴリズムとしても知られている）による方法である。第２の方法は、ワンショット復調アルゴリズムによる方法である。本説明では、図１に概略的に示されるようなＸ線タルボット干渉計により生成された単一（ワンショット）のフリンジ模様からの復調を考える。２次元（２Ｄ）交差格子によりフリンジ模様が生成される結果、網目状の模様が複雑な調波構造を示す。図２Ａは、合成模様２１０（この場合、骨関節の撮影から導出された）の一例を示し、図２Ｂは、対応するフーリエ変換２２０を示す。フーリエ変換２２０は、９個のサイドローブの調波構造を示す。
【０００３】
フリンジ解析の分野は、フリンジ模様の１つ以上の画像の復調に関する。各画像は、各々が各画像において異なる位相を有する少なくとも１つの正弦波フリンジの合計を含む。そのような画像は以下の形式をとる。
【０００４】
【数１】
2016087473.tif 000003
【０００５】
式中、Ｉｉは、正弦波フリンジの合計から取得された強度値を含む座標ベクトルｒ＝（ｘ，ｙ）により索引を付けられた２次元撮像画像であり、ｉが画像番号１≦ｉ≦ｎを示し、Ｓｉは非正弦波ノイズ画像を示し、Ａは正の強度画像を示し、ｊはフリンジ番号１≦ｊ≦ｋを示す。各フリンジは、通常は正の振幅変調画像ｍｊ、位相画像ωｉ及び実数φｉｊにより示され、画像ｉにおけるフリンジｊに対する位相オフセットである。
【０００６】
説明のために画像Ｉｉは２次元であると仮定するが、次元数を限定する理由はなく、本明細書は２次元以上に対して適用可能である。
【０００７】
フリンジは、投影された正弦曲線、撮影されたテストパターン、コヒーレント源の干渉パターン又はモアレ模様などの種々の手段により作成される。画像は、可視光、Ｘ線、オーディオ波又は２次元以上で撮影される他の何らかの情報源から構成される。
【０００８】
復調処理は、振幅Ａ、振幅変調ｍｊ及び位相ωｊの画像Ｉｉの集合からの回復である。振幅画像は照明の強度を示す単位を有し、ｍｊは、通常は正で実数の無次元振幅関数（ｕｎｉｔｌｅｓｓ  ａｍｐｌｉｔｕｄｅ  ｆｕｎｃｔｉｏｎ）であり、１．０より小さい値を含み、ωｊは、複素数値の位相である（ｍｊは、状況によっては負であってよい）。関数ωｊは、通常は波長の総数の小数部、即ちモジュロ２πラジアンとしてのみ決定され、位相アンラッピングアルゴリズムは、この範囲を越えた偏位を有する総位相を連続関数として再構成するために使用される。
【０００９】
位相ステッピングアリゴリズムを使用することにより、フリンジが種々の位相オフセットφｉｊを含む複数の画像を撮像することと、第１にｋ個の異なるフリンジ模様を互いに分離し、第２に実数値フリンジ模様を復調してすべてのｋ個のフリンジ及びすべての画像において存在する２次元座標値に対する振幅画像ｍｊ及び位相画像ωｉを取得するために種々の技術を使用することとにより、これらの数量を直接算出できる。
【００１０】
フリンジを分離及び復調するために必要な数の画像を撮像するには、いくつかの適応例において過度の処理時間を必要とするため、実用的ではない。また、多くの画像を撮像することにより、撮影システムは被検体の動きなどの問題により影響されやすくなる。従って、より少ない数の画像又は単一の画像を使用して復調を実行することが望ましい。これらの場合、単一の画像のみ又はフリンジを直接分離するのに必要とされるよりも少ない画像を取得し、振幅画像及び位相画像を抽出するために別の手段を使用することが可能である。
【００１１】
画像取得時に位相ステップが認識されていないか、あるいはフリンジの特定の周波数が認識されていない場合、これらのパラメータは、フーリエ変換を使用して画像を解析することにより又は他の技術により測定される。フリンジ模様が適切に帯域制限される、即ちフーリエ変換における各フリンジ模様の画像がフーリエ変換の重なり合わない領域にローカライズされると仮定すると、フリンジ模様は、周波数に依存してフリンジを隔離する多数の技術のうちのいずれかにより分離される。しかし、解決方法に曖昧さを生じさせ且つそのような方法の有効性を制限するいくつかの影響がある。
【００１２】
第１に、一般に信号は、振幅及び位相の変調によりキャリア周波数に変調される。この変調によりフーリエ領域において信号が拡散し、通常は結果として得られるフリンジ模様の帯域幅を制限できない。信号の帯域幅がフーリエ領域におけるその信号と隣接信号との間の空間より大きい場合、信号は他の隣接信号と重なり合う。２つの信号が空間領域及びフーリエ領域において重なり合う場合、それらの和のみが検出可能であり、信号生成処理に関連した更なる制約を認識することなく、フーリエ技術又は空間技術を使用してフリンジを分離することは通常は不可能である。
【００１３】
第２に、信号生成処理自体により、所望されるより多くの周波数を生成させる恐れがある。例えば、撮影システムが２つの交差格子を組み合わせた結果として得られるモアレ模様によりフリンジを生成する場合、結果として得られる画像のスペクトルは、理想の正弦波フリンジよりかなり広い格子のスペクトルに一致する。
【００１４】
第３に、変調周波数が撮影センサの領域にわたり厳密な数の波長を生成しない場合、対応するフーリエ変換は、相対的に低速に中心（キャリア）周波数からドロップするｓｉｎｃ関数を生成するため、フルスペクトルをある程度範囲に含む。
【００１５】
第４に、実数の余弦関数は、原点について対称であり且つエルミート対称性と釣り合うフーリエ領域において２つの画像（「エルミートの対」と呼ばれる）を生成する。しかし、ここでも、位相関数の帯域幅が広すぎる場合、２つのエルミートの対が重なり合う結果、更なる曖昧さが生じる。
【００１６】
これらの影響をすべて組み合わせる結果、クロストークであるフーリエ領域において重なり合う種々の信号が得られる。これらの重なり合う信号を適切に分離できない結果、信号は復調の間に混合されるようになり、歪みを引き起こす。
【００１７】
単純な復調方法は、スペクトル領域及び空間領域におけるフィルタリングと、空間領域における複素周波数を含む乗算によるキャリア周波数の復調又は同等にキャリア周波数から原点への空間領域における窓のシフトとを組み合わせることにより、周波数スペクトルにおいて種々の信号を分離する。
【００１８】
しかし、これらの方法は、重なり合う信号を考慮できないため、クロストークによる歪みを除去しない。
【発明の概要】
【発明が解決しようとする課題】
【００１９】
本発明において、軸上フリンジを分離できるようにするために、軸外項を抑制するためにフリンジ画像上の非線形関数を使用することに関する方法が開示される。別の非線形関数は、第１の非線形関数の結果得られる高調波を除去するために使用される。
【課題を解決するための手段】
【００２０】
本発明の一態様によれば、交差項位相成分を含む位相変調フリンジ模様を含む撮像画像を復調する方法であって、
前記交差項位相成分を減衰するために少なくとも第１の非線形変換のパラメータを決定する工程であって、前記第１の非線形変換の前記パラメータが、前記交差項位相成分により導入された復調画像におけるアーチファクトを減少するために決定される工程と、
前記交差項位相成分を減衰するために、前記決定されたパラメータによる少なくとも前記第１の非線形変換を前記撮像画像に適用することにより、前記撮像画像を復調する工程と、を有することを特徴とする方法が提供される。
【００２１】
前記第１の非線形変換のパラメータを前記決定する工程は、
前記撮像画像を撮影した撮影システムを使用して撮像された校正画像に対して、前記交差項位相成分の大きさと軸上位相成分の大きさとの比を決定する工程と、
前記決定された比を使用して前記第１の非線形変換のパラメータを決定する工程と、を有する。
【００２２】
パラメータを前記決定する工程は、
更なる復調のために前記交差項位相成分の影響を軽減するようにパラメータ化非線形変換を実行する工程と、
交差項位相成分を含む位相変調フリンジ模様の前記撮像画像を撮影した撮影システムを用いて撮像された複数の参照画像を受信する工程と、
参照復調画像を生成するためにマルチショット復調方法を使用して前記参照復調画像に対して第１の復調を行う工程と、
初期の非線形変換を形成するように前記実行されたパラメータ化非線形変換に対する初期のパラメータを識別する工程と、
前記初期の非線形変換に従って変換された前記参照復調画像に適用されたシングルショット復調方法を使用して前記参照復調画像のうちの１つに対して第２の復調を行うことにより、候補復調画像を形成する工程と、
前記非線形変換のパラメータを調整するために、前記候補復調画像を前記参照復調画像と比較する工程と、を更に有する。
【００２３】
方法は、前記候補復調画像と前記参照復調画像との間の誤差を減少するように変更された候補復調画像を形成するために、前記調整されたパラメータを用いた前記第２の復調を繰り返す工程を更に有する。
【００２４】
前記調整されたパラメータは、前記撮影システムにより撮像された更なる画像の復調における交差項位相成分の影響を軽減するために使用される。
【００２５】
前記撮像画像を復調する工程は、
前記減衰された交差項位相成分を含む第１の画像を生成するために前記第１の非線形変換を前記撮像画像に適用する工程と、
フーリエ領域において前記第１の画像のｙ項及びｘ項のそれぞれを前記第１の画像からフィルタリングすることにより、ｘ位相画像及び位相画像を生成する工程と、を更に含む。
【００２６】
前記撮像画像を復調する工程は、
前記第１の非線形変換を前記撮像画像に適用することにより導入された高調波を減衰するために構成された第２の非線形変換を実行する工程であって、前記第２の非線形変換が前記第１の非線形変換を反転することにより実行される工程と、
前記第１の非線形変換の動作により導入された高調波を減衰するために、前記第２の非線形変換をｘ位相画像及びｙ位相画像の各々に適用する工程と、を更に有する。
【００２７】
フィルタリングは、
前記第１の画像を前記フーリエ領域に変換する工程と、
少なくとも１つのマスク画像を前記変換された第１の画像に適用する工程と、
前記マスク画像を空間領域に変換する工程と、を有する。
【００２８】
前記適用する工程は、前記変換された第１の画像に、前記ｘ位相画像を生成するために第１のマスクを乗算し、ｙ位相画像を生成するために第２のマスクを乗算する工程を有する。
【００２９】
前記少なくとも１つのマスクは、前記画像から所定の周波数スペクトルを除去するように構成される。
【００３０】
方法は、分離されたフリンジ模様に対して、前記位相画像から低周波数スペクトルを減算する工程を更に有する。
【００３１】
前記第１の非線形変換は、対数関数による変換を含む。前記第２の非線形変換は、指数関数による変換を含む。
【００３２】
あるいは、前記第１の非線形変換及び第２の非線形変換は、冪関数による変換を含む。
【００３３】
通常、前記撮影システムは、Ｘ線干渉計システムを含む。
【００３４】
他の態様も開示される。
【図面の簡単な説明】
【００３５】
次に、以下の図面を参照して、本発明の少なくとも１つの実施形態を説明する。
【図１】Ｘ線タルボット干渉計に対する実験設定を示す図である。
【図２Ａ】タルボットモアレ画像及びそのフーリエ変換の例示的な一例を示す図である。
【図２Ｂ】タルボットモアレ画像及びそのフーリエ変換の例示的な一例を示す図である。
【図３】撮像及び復調処理を示す概略的なフローチャートである。
【図４】非線形関数を使用する復調の間に軸外項を減少する方法を概略的に示すフローチャートである。
【図５】非線形変換を画像における画素に適用する方法を概略的に示すフローチャートである。
【図６】撮像画像から分離されたフリンジ模様を生成している間に軸外項及び高調波の双方を減少する方法を概略的に示すフローチャートである。
【図７】図７Ａ〜図７Ｆは、図６において説明される処理に関わるフーリエスペクトルの位置を示す図である。
【図８Ａ】図６において説明される処理に関わる分離マスクを示す図である。
【図８Ｂ】図６において説明される処理に関わる分離マスクを示す図である。
【図９】分離されたフリンジ模様を復調する方法を概略的に示すフローチャートである。
【図１０】復調処理中に使用される非線形関数に対するパラメータを推定する方法を概略的に示すフローチャートである。
【図１１】復調処理中に使用される非線形関数に対する最適なパラメータを識別する方法を概略的に示すフローチャートである。
【図１２】図１２Ａ〜図１２Ｈは、図９において説明される処理に関わるスペクトルマスクを示す図である。
【図１３】復調処理中に使用される非線形関数に対するパラメータを推定する方法を概略的に示すフローチャートである。
【図１４】ログオフセット変換を使用してタルボットモアレ画像を復調する方法を概略的に示すフローチャートである。
【図１５】パラメータ化冪関数を使用してタルボットモアレ画像を復調する方法を概略的に示すフローチャートである。
【図１６】図１４において説明される処理に関わる吸収画像を推定する方法を概略的に示すフローチャートである。
【図１７】最初にフリンジ模様を軸上フリンジに分離し、次に各模様を個々に復調することによる復調の方法を概略的に示すフローチャートである。
【図１８Ａ】説明される配列が実行される汎用コンピュータシステムを集合的に形成する概略ブロック図である。
【図１８Ｂ】説明される配列が実行される汎用コンピュータシステムを集合的に形成する概略ブロック図である。
【発明を実施するための形態】
【００３６】
構成
本発明は、１つ以上の画像として提供されたフリンジ（縞）の２つ以上のグループを含むフリンジの復調に関する。撮像された画像の数が別個のフリンジを適切に分離するのに十分でない状況において及びフリンジの周波数帯域が重なり合う場合、復調の問題は曖昧になる。従って、復調画像はクロストークにより歪む。
【００３７】
多くの撮影システムはフリンジ画像を生成する。一例は、２次元Ｘ線タルボット（ＸＴ）画像システムである。図１は、Ｘ線タルボット（Ｘ−ｒａｙ  Ｔａｌｂｏｔ）干渉計１００の一例を示す概略図である。Ｘ線タルボット干渉計１００において、２次元線源格子Ｇ０  １０１、被検体１０２及び２次元位相格子Ｇ１  １１０は、Ｘ線源１０４により照明される。自己像１２０は、タルボット効果（Ｔａｌｂｏｔ  ｅｆｆｅｃｔ）により位相格子Ｇ１  １１０の後ろに形成される。Ｘ線波の位相が被検体１０２によりシフトされるため、自己像は変形する。変形した自己像を解析することにより、被検体１０２の特徴を推定できる。この種のＸ線撮影では、対比を生じさせるために吸収ではなく位相差を使用するため、ＸＴシステムの感度は、従来の吸収Ｘ線よりかなり高い。図１において、自己像の位置における第２の２次元吸収格子Ｇ２  １３０は、明確なモアレフリンジ模様を生成する。一般に、第２の吸収格子Ｇ２  １３０の周期は、自己像の周期に類似するように設計されるため、変形した自己像と第２の吸収格子Ｇ２  １３０とを重ね合わせることにより生成されたモアレ模様は、画像センサ１４０がモアレ模様を解像して撮像するために、元の自己像を拡大したバージョンのパターンを与える。シングルショット方法によると、処理は、画像センサ１４０により撮像された単一の画像に対して実行される。マルチショット方法において、位相変化は、格子１０１、１１０及び１３０（例えば、シフト及び／又は回転）、並びに変化ごとに撮像された画像に課される。
【００３８】
図２Ａは、Ｘ線タルボット干渉計から取得された画像２１０の一例を示す。画像２１０は、動物を撮影した結果として画像センサから取得されたモアレ画像である。動物の骨は画像２１０の右において可視であり、動物の毛は画像２１０の左において識別可能であり、画像にわたり可視であり且つモアレ模様である点の模様は、符号化された振幅及び位相の情報を含む。
【００３９】
タルボット効果の複雑さを無視すると、図２Ａにおける点のこの模様は、Ｇ１格子１１０及びＧ２格子１３０により生成された格子パターンの積から構成されたかのように解釈される。格子パターンは、従来のＸ線撮影センサに対して可視になるように、モアレ効果により十分に低い周波数にシフトされた高周波数パターンである。
【００４０】
格子状の格子の一般のスペクトルは、それ自体が格子であるが、各格子要素の周りのフィルタリング効果により、モアレスペクトルは、実質的に直流成分及び基本波のみから構成され、後続の高調波は非常に低い振幅である。
【００４１】
従って、モアレ画像は、実質的に４つの位相変調余弦強度関数及び振幅変調余弦強度関数の合計から構成される。４つの位相変調余弦強度関数及び振幅変調余弦強度関数のうちの１つはｘ軸に沿って変動し、１つはｙ軸に沿って変動し、各々は線（ｘ＋ｙ）及び（ｘ−ｙ）のそれぞれに沿って変動し、即ち０°、９０°、４５°及び１３５°のそれぞれの角度で変動する。これらの余弦関数の各々は、吸収関数により減衰される。
【００４２】
図２Ｂは、モアレ模様の２３０、２４０、２５０、２６０、２７０などの最強のスペクトル成分が明確に可視である画像２１０のフーリエ変換２２０を示す。図２Ｂにおいて、スペクトル２２０は９個の成分を有する。ラベル付けされた成分は、直流成分２３０、ｘ成分２４０、ｙ成分２５０、（ｘ＋ｙ）交差成分２７０、（ｘ−ｙ）交差成分２６０、（ｙ−ｘ）交差成分２８０及び−（ｘ＋ｙ）交差成分２９０を含む。図２Ｂにおいて、成分２４０及び２５０は軸上項の例であり、成分２６０〜２９０は軸外項の例である。
【００４３】
フーリエ変換は、関数又はベクトルを複素指数関数の積分又は和として示す基底変換である。フリンジは、複素関数とは対照的な正実関数の合計から構成される。しかし、各フリンジは、２つの複素指数関数の和として示されてもよく、従って以下の通りである（式中、ｉ＝√（−１））。
【００４４】
【数２】
2016087473.tif 000004
【００４５】
フーリエ変換は、他のあらゆる次元にわたり別個に変換を計算するだけで、１次元を超える次元に容易に拡張される。即ち、２Ｄ画像のフーリエ変換は、単に列のフーリエ変換又は行のフーリエ変換である。
【００４６】
実数２Ｄフリンジ模様は、フーリエ領域に変換される場合、２つの複素指数関数の合計として示されるため、エルミート対称性の２つのローブとして現れる。２つのローブが重なり合わない場合、複素指数関数のいずれか１つを示す周波数は、他方のローブをゼロにすることにより抽出され、この抽出されたローブの逆フーリエ変換により、複素対数を使用して抽出された位相を用いた求積法において復調されたフリンジが得られる。
【００４７】
フリンジの位相ωが実質的に線形関数であり、且つ振幅がおおよそ一定であって、ｗが画像の幅であり、ｈが画像の高さであり且つｕ及びｖがｘ方向、ｙ方向それぞれの波長の数である場合、フリンジの２次元（２Ｄ）周波数を以下のように規定できる。
【００４８】
【数３】
2016087473.tif 000005
【００４９】
【数４】
2016087473.tif 000006
【００５０】
上記から、位相関数φは、線形であるために単一周波数に近接することが理解されるだろう。上記の関数は、おおよそ以下の通り示される。
【００５１】
【数５】
2016087473.tif 000007
【００５２】
フリンジの集合におけるフリンジ周波数の各々が一定である場合、フリンジの和のフーリエ変換は、φ０及び−φ０の位相をそれぞれ含む各位置（ｕ，ｖ）及び（−ｕ，−ｖ）における２つのデルタ関数を含む。
【００５３】
フーリエ変換の各ピークの位置を測定することにより、各フリンジのキャリア周波数が推定され、各ピークの位相を測定することにより、各フリンジの位相オフセットが推定される。各ピークの振幅のフーリエ変換の直流位置の振幅に対する比率を比較することにより、ｍが推定される。即ち、ｍ≒｜ａ｜／｜Ａ｜である。
【００５４】
図２Ｂに戻ると、フーリエ変換２２０の中心の領域２３０が吸収スペクトルであるのに対し、スペクトル成分２４０は、ｘ位相導関数を符号化するｙ方向にフリンジを含むスペクトルを示し、スペクトル成分２５０は、ｙ位相導関数を符号化するｘ方向にフリンジを含むスペクトルを示し、スペクトル成分２６０は、ｘ位相導関数とｙ位相導関数との差分を含む望ましくない軸外項を含み、スペクトル成分２７０は、ｘ位相導関数及びｙ位相導関数の和を含む望ましくない軸外項を含む。
【００５５】
位相導関数が測定される軸及びこの導関数を符号化するフリンジが変動する軸は、干渉計システム１００において格子の位置を変更することにより個々に変動する。図示されるような格子の配向は、システム１００のスペクトル特性に関して対称である正の導関数及び負の導関数の双方を維持するように選択される。
【００５６】
図２Ｂにおけるラベル付けされていないスペクトル成分は、ラベル付けされたスペクトル成分のエルミートの対であり、モアレ画像が実数であるために現れる。尚、図２Ｂは例示するためだけのものであり、スペクトル成分の重なり合いは、理解を容易にするために示されていない。実質的なクロストークが生じている場合、各キャリア周波数の周りのスペクトルコンテンツはより広く現れ、図２Ｂの黒色の領域は重なり合う。
【００５７】
尚、「軸上（ｏｎ−ａｘｉｓ）」及び「軸外（ｏｆｆ−ａｘｉｓ）」という用語は、センサ１４０ではなく、格子１０１、１１０、１３０の軸をそれぞれ示すため、画像生成処理は、「軸上」項（ｏｎ−ａｘｉｓ  ｔｅｒｍｓ）がフーリエ変換のｕ軸及びｖ軸上に現れることを実際には必要としない。これらの項がｕ軸及びｖ軸に割り当てられない場合、それらの位置は、検出されたこれらの周波数の和及び差分に対応する位置において望ましくない項がある状態で、被検体を含まない校正モアレ画像のフーリエ変換において２つのスペクトルピークを取得することにより取得される。しかし、簡潔にするために、以下においてこれをそのように仮定する。
【００５８】
画像のスペクトルは、上述のフリンジモデルの特殊なケースとして説明され、以下の式をとる。
【００５９】
【数６】
2016087473.tif 000008
【００６０】
式中、Ｓｉは２Ｄノイズ信号であり、Ａは被検体の２Ｄ吸収画像であり、ｍｘ及びｍｙは、軸上キャリアの振幅変調を含む２Ｄ画像であり、ｍｕ及びｍｖは、軸外キャリアの振幅変調を含む２Ｄ画像であり、ωｘ及びωｙは、軸上キャリアの周波数であり、Φｘ及びΦｙは位相画像である。各Φは、原点における位相を規定する単一の任意の位相オフセット及び光路長の導関数を示す位相変調関数を含む。種々の材料が種々の屈折率を有するため、光路長は、Ｘ線が被検体と交差する場合に変動する。
【００６１】
シングルショット方法における破損の主な原因は、格子の交差項間におけるスペクトルの重なり合い又はクロストークである。これは、高帯域信号に関する問題であり、軸外（対角線状に変調された）信号２６０及び２７０がｘ変調信号及びｙ変調信号と重なり合うフーリエ空間画像２２０において見られる。これらの軸外信号は、クロストークによる軸上信号の歪みを引き起こし、且つ軸上信号２４０及び２５０においてそれを超える余分な位相情報を全く提供しないため、望ましくない。軸外信号の振幅は、格子要素におけるセルのスペクトルの結果である。軸外信号が軸上フリンジの直積（ｄｉｒｅｃｔ  ｐｒｏｄｕｃｔ）に対応することになっていた場合、対数関数を使用する準同型処理（ｈｏｍｏｍｏｒｐｈｉｃ  ｐｒｏｃｅｓｓｉｎｇ）は、これらの軸外フリンジ  （ｏｆｆ−ａｘｉｓ  ｆｒｉｎｇｅｓ）を除去する。これは、交差項アーチファクトを減少し、且つ帯域幅及び品質の双方を向上させる効果を有する。
【００６２】
準同型処理の主な目的は、乗法処理を加法処理に分離し、線形フィルタリングを適用して、処理済み画像を元の領域に再度変換することである。線形フィルタリングは、通常、フーリエ変換領域において実行される。
【００６３】
３つの成分、即ちＸ線吸収因子Ａ、ｘ方向フリンジ項及びｙ方向フリンジ項の積としてモデル化された２Ｄ強度モアレ模様の乗法モデルを仮定すると、スペクトル画像は、式２の特定の形式に分解され、以下の通りである。
【００６４】
【数７】
2016087473.tif 000009
【００６５】
ここで、式２の信号又は画像の強度は、対数的にマッピングされて以下の式が得られる。
【００６６】
【数８】
2016087473.tif 000010
【００６７】
次に、３つの成分の元の積は、３つの成分の加法に変換されている。ｘフリンジ項及びｙフリンジ項の積がないため、次に、軸外交差項は、除去され、フーリエ領域において観察不可能である。軸外交差項を除去することにより、ｘフリンジ模様及びｙフリンジ模様は、フーリエ領域においてきれいに分離され、簡素なスペクトルフィルタにより隔離される。
【００６８】
フィルタリング後、２つの画像は、空間領域に再度変換され、別個に解析される。空間領域において結果として得られるフィルタリングされた画像の逆対数、即ち指数関数を計算することにより、対数によりキャリア周波数において生成された高調波が厳密に除去されることは、本明細書において後で考慮される重要なポイントである。元のキャリア周波数のフーリエ領域において重なり合いがない場合、この処理は、ｘキャリア及びｙキャリアに対するフリンジを別個に正常に再構成する。
【００６９】
高調波がべき乗により除去された後、対応する解析信号を算出することにより、変調項及び位相項を分離できる。実数の余弦関数の解析信号は、複素指数関数である。このことから、この複素信号の大きさ及び位相項を算出し、ｘフリンジ模様及びｙフリンジ模様それぞれの変調項ｍｘ及びｍｙ、並びに位相項Φｘ及びΦｙの推定値を与えることができる。
【００７０】
しかし、積モデルは、Ｘ線タルボットシステムの物理的特性に一致せず、実数データに対処する際には有用でない。
【００７１】
物理的な構成例
図１８Ａ及び図１８Ｂは、説明される種々の構成を実施する汎用コンピュータシステム１８００を示す。
【００７２】
図１８Ａに示されるように、コンピュータシステム１８００は、コンピュータモジュール１８０１、キーボード１８０２、マウスポインタデバイス１８０３、スキャナ１８２６、カメラ及びマイク１８８０などの入力デバイス、並びにプリンタ１８１５、ディスプレイデバイス１８１４及びスピーカ１８１７を含む出力デバイスを含む。接続回線１８２１を介して通信ネットワーク１８２０との間で通信を実行するために、コンピュータモジュール１８０１により外部変復調器（モデム）トランシーバデバイス１８１６が使用される。通信ネットワーク１８２０は、インターネットなどのワイドエリアネットワーク（ＷＡＮ）、セルラ通信ネットワーク又はプライベートＷＡＮである。接続回線１８２１が電話回線である場合、モデム１８１６は従来の「ダイヤルアップ」モデムである。あるいは、接続回線１８２１が大容量（例えば、ケーブル）接続回線である場合、モデム１８１６はブロードバンドモデムである。通信ネットワーク１８２０への無線接続のために、無線モデムが更に使用される。
【００７３】
コンピュータモジュール１８０１は、通常、少なくとも１つのプロセッサユニット１８０５及びメモリユニット１８０６を含む。例えばメモリユニット１８０６は、半導体ランダムアクセスメモリ（ＲＡＭ）及び半導体読み取り専用メモリ（ＲＯＭ）を有する。コンピュータモジュール１８０１は、ビデオディスプレイ１８１４、スピーカ１８１７及びマイク１８８０に結合するオーディオビデオインタフェース１８０７、キーボード１８０２、マウス１８０３及びスキャナ１８２６に結合し、任意にジョイスティック又は他のヒューマンインタフェースデバイス（図示せず）に結合するＩ／Ｏインタフェース１８１３、並びに外部モデム１８１６及びプリンタ１８１５に対応するインタフェース１８０８を含む複数の入出力（Ｉ／Ｏ）インタフェースを更に含む。いくつかの実現形態において、モデム１８１６は、コンピュータモジュール１８０１の中に、例えばインタフェース１８０８の中に組み込まれる。コンピュータモジュール１８０１は、コンピュータシステム１８００が接続回路１８２３を介して図１に示されたものに類似する撮影構成を実現するＸ線タルボット撮影装置１８９９に結合できるようにするローカルネットワークインタフェース１８１１を更に有する。インタフェース１８１１は、生体試料などの撮影されている被検体１０２を表す画像センサ１４０により撮像された画像データを受信するために、コンピュータモジュール１８０１により使用される。ネットワークインタフェース１８１１は、とりわけＥｔｈｅｒｎｅｔ（登録商標）回路カード、Ｂｌｕｅｔｏｏｔｈ（登録商標）無線構成又はＩＥＥＥ８０２．１１無線構成などを含むネットワークインタフェースである。あるいは、インタフェース１８１１は、接続回路１８２３を介して画像データストリームを受信するように構成されてもよい。撮影装置１８９９は、顕微鏡又は他の撮影装置として構成され、Ｘ線源などの撮影ソース１８９８を有する。センサ１４０により撮像された生画像データは、通常、本発明に係る処理を含む後続の処理のためにハードディスクドライブ（ＨＤＤ）１８１０に記憶される。
【００７４】
Ｉ／Ｏインタフェース１８０８及び１８１３は、直列接続性及び並列接続性のうちいずれか一方又は双方を供与し、前者は、通常ＵＳＢ（Ｕｎｉｖｅｒｓａｌ  Ｓｅｒｉａｌ  Ｂｕｓ）規格に従って実現され、対応するＵＳＢコネクタ（図示せず）を有する。記憶デバイス１８０９が設けられ、通常、ＨＤＤ１８１０を含む。フロッピー（登録商標）ディスクドライブ及び磁気テープドライブ（図示せず）などの他の記憶デバイスが使用されてもよい。光ディスクドライブ１８１２は、通常、不揮発性データ源として機能するように設けられる。システム１８００への適切なデータ源として、例えば光ディスク（例えば、ＣＤ−ＲＯＭ、ＤＶＤ、Ｂｌｕ−ｒａｙ（登録商標）Ｄｉｓｃ）、ＵＳＢ−ＲＡＭ、ポータブルな外部ハードドライブ及びフロッピー（登録商標）ディスクなどのポータブルメモリデバイスが使用される。
【００７５】
コンピュータモジュール１８０１の構成要素１８０５〜１８１３は、通常、相互接続バス１８０４を介して、及び当業者に知られている従来のコンピュータシステム１８００の動作モードとなる方式で通信する。例えばプロセッサ１８０５は、接続回線１８１８を使用してシステムバス１８０４に結合される。同様に、メモリ１８０６及び光ディスクドライブ１８１２は、接続回線１８１９によりシステムバス１８０４に結合される。説明した構成を実施可能なコンピュータの例には、ＩＭＢ−ＰＣ及びそのコンパチブル、Ｓｕｎ  Ｓｐａｒｃｓｔａｔｉｏｎｓ、Ａｐｐｌｅ  ＭａｃＴＭ又は同様のコンピュータシステムがある。
【００７６】
本明細書において説明される画像復調の方法は、コンピュータシステム１８００を使用して実現され、以下に説明される図３〜図１７の処理は、コンピュータシステム１８００内で実行可能な１つ以上のソフトウェアアプリケーションプログラム１８３３として実現される。特に、画像復調ステップは、コンピュータシステム１８００内で実行されるソフトウェア１８３３中の命令１８３１（図１８Ｂを参照）により実行される。ソフトウェア命令１８３１は、各々が１つ以上の特定のタスクを実行する１つ以上のコードモジュールとして形成される。ソフトウェアは２つの個別の部分に分割されてもよく、その場合、第１の部分及び対応するコードモジュールは、復調方法を実行し、第２の部分及び対応するコードモジュールは、第１の部分とユーザとの間のユーザインタフェースを管理する。
【００７７】
ソフトウェアは、例えば以下に説明される記憶デバイスを含むコンピュータ可読媒体に記憶される。ソフトウェアは、コンピュータ可読媒体からコンピュータシステム１８００にロードされ、コンピュータシステム１８００、より具体的にはプロセッサ１８０５により実行される。そのようなソフトウェア又はコンピュータプログラムが記録されたコンピュータ可読媒体は、コンピュータプログラムである。コンピュータシステム１８００においてコンピュータプログラムを使用することにより、画像復調のために有利な装置が生み出されるのが好ましい。
【００７８】
ソフトウェア１８３３は、通常ＨＤＤ１８１０又はメモリ１８０６に記憶される。ソフトウェアは、コンピュータ可読媒体からコンピュータシステム１８００にロードされ、コンピュータシステム１８００により実行される。従って、例えばソフトウェア１８３３は、光ディスクドライブ１８１２により読み取られる光学的可読ディスク記憶媒体（例えば、ＣＤ−ＲＯＭ）１８２５に記憶される。そのようなソフトウェア又はコンピュータプログラムが記録されたコンピュータ可読媒体は、コンピュータプログラムである。
【００７９】
場合によっては、アプリケーションプログラム１８３３は、１つ以上のＣＤ−ＲＯＭ１８２５に符号化された形態でユーザに供給され、対応するドライブ１８１２を介して読み取られるか、あるいはユーザによりネットワーク１８２０から読み取られる。更に、ソフトウェアは、他のコンピュータ可読媒体からコンピュータシステム１８００にロードされることも可能である。コンピュータ可読媒体とは、記録された命令及び／又はデータを実行及び／又は処理のためにコンピュータシステム１８００に提供する何らかの有形の非一時的記憶媒体である。そのような記憶媒体の例には、デバイスがコンピュータモジュール１８０１の内部にあるか又は外部にあるかに関わらず、フロッピー（登録商標）ディスク、磁気テープ、ＣＤ−ＲＯＭ、ＤＶＤ、Ｂｌｕ−ｒａｙ（登録商標）  Ｄｉｓｃ、ハードディスクドライブ、ＲＯＭ又は集積回路、ＵＳＢメモリ、光磁気ディスク、あるいはＰＣＭＰＩＡカードなどのコンピュータ可読カードがある。コンピュータモジュール１８０１へのソフトウェア、アプリケーションプログラム、命令及び／又はデータの提供に関与してもよい一時的又は非一時的コンピュータ可読送信媒体の例には、無線送信チャネル又は赤外線送信チャネル、並びに別のコンピュータ又はネットワーク化デバイスへのネットワーク接続、Ｅメール送信を含むインターネット又はイントラネット及びウェブサイトなどに記録された情報などがある。
【００８０】
先に述べたアプリケーションプログラム１８３３の第２の部分及び対応するコードモジュールは、ディスプレイ１８１４にレンダリングされるか又は他の方法により表現されるべき１つ以上のグラフィカルユーザインタフェース（ＧＵＩ）を実現するために実行される。コンピュータシステム１８００及びアプリケーションのユーザは、通常はキーボード１８０２及びマウス１８０３の操作によって、ＧＵＩと関連するアプリケーションに制御コマンド及び／又は入力を提供するために機能適応方式でインタフェースを操作する。スピーカ１８１７を介して出力される音声プロンプト及びマイク１８８０を介して入力される音声コマンドを利用するオーディオインタフェースなどの他の形態の機能適応型ユーザインタフェースも実現されてよい。
【００８１】
図１８Ｂは、プロセッサ１８０５及び「メモリ」１８３４の詳細な概略ブロック図である。メモリ１８３４は、図１８Ａのコンピュータモジュール１８０１によりアクセス可能なすべてのメモリモジュール（ＨＤＤ１８０９及び半導体メモリ１８０６を含む）の論理的集合体を表す。
【００８２】
コンピュータモジュール１８０１が初めて起動されると、電源オン自己診断テスト（ＰＯＳＴ）プログラム１８５０が実行される。ＰＯＳＴプログラム１８５０は、通常、図１８Ａの半導体メモリ１８０６のＲＯＭ１８４９に記憶される。ソフトウェアを記憶するＲＯＭ１８４９などのハードウェアデバイスはファームウェアと呼ばれる場合もある。ＰＯＳＴプログラム１８５０は、適正な機能を確保するためにコンピュータモジュール１８０１内部のハードウェアを検査し、通常、プロセッサ１８０５、メモリ１８３４（１８０９、１８０６）及び通常はＯＭ１８４９に記憶されている基本入出力システムソフトウェア（ＢＩＯＳ）モジュール１８５１が正しく動作しているか否かを検査する。ＰＯＳＴプログラム１８５０の実行後に問題がなければ、ＢＩＯＳ１８５１は、図１８Ａのハードディスクドライブ１８１０を起動する。ハードディスクドライブ１８１０の起動によって、ハードディスクドライブ１８１０に常駐するブートストラップローダプログラム１８５２がプロセッサ１８０５を介して実行される。このプログラムは、オペレーティングシステム１８５３をＲＡＭメモリ１８０６にロードし、この時点で、オペレーティングシステム１８５３は動作を開始する。オペレーティングシステム１８５３は、プロセッサ管理、メモリ管理、デバイス管理、記憶装置管理、ソフトウェアアプリケーションインタフェース及びユーザインタフェース全般を含む種々の高レベル機能を遂行するためにプロセッサ１８０５により実行可能なシステムレベルアプリケーションである。
【００８３】
オペレーティングシステム１８５３は、コンピュータモジュール１８０１で実行中の各処理又は各アプリケーションが別の処理に割り当てられたメモリと衝突せずに実行するための十分なメモリを確保するために、メモリ１８３４（１８０９、１８０６）を管理する。更に、各処理を効果的に実行できるように、図１８Ａのシステム１８００で利用可能な異なる種類のメモリは適正に使用されなければならない。従って、メモリ集合体１８３４は、メモリの特定のセグメントがどのように割り当てられるかを例示するのではなく（特に指示のない限り）、コンピュータシステム１８００によりアクセス可能なメモリの全体図及びそのメモリがどのように使用されるかを提示することを意図する。
【００８４】
図１８Ｂに示されるように、プロセッサ１８０５は、制御ユニット１８３９、演算論理ユニット（ＡＬＵ）１８４０及びキャッシュメモリと呼ばれる場合もあるローカルメモリ又は内部メモリ１８４８を含む複数の機能モジュールを含む。キャッシュメモリ１８４８は、通常、レジスタ部分に複数の記憶レジスタ１８４４〜１８４６を含む。１つ以上の内部バス１８４１は、これらの機能モジュールを機能的に相互接続する。通常、プロセッサ１８０５は、接続回線１８１８を使用して、システムバス１８０４を介して外部デバイスと通信するための１つ以上のインタフェース１８４２を含む。メモリ１８３４は、接続回線１８１９を使用してバス１８０４に結合される。
【００８５】
アプリケーションプログラム１８３３は、条件付き分岐命令及びループ命令を含む一連の命令１８３１を含む。プログラム１８３３は、プログラム１８３３の実行時に使用されるデータ１８３２を更に含んでもよい。命令１８３１及びデータ１８３２は、メモリの記憶場所１８２８、１８２９、１８３０及び１８３５、１８３６、１８３７にそれぞれ記憶される。記憶場所１８３０に示される命令により表されるように、命令１８３１と記憶場所１８２８〜１８３０の相対的大きさに応じて、特定の１つの命令が１つの記憶場所に記憶される。あるいは、記憶場所１８２８及び１８２９に示される命令セグメントにより表されるように、命令は複数の部分に分割されてもよく、各部分はそれぞれ別の記憶場所に記憶される。
【００８６】
一般に、プロセッサ１８０５は、プロセッサ１８０５で実行される命令のセットを提供される。プロセッサ１８０５は後続入力を待機し、プロセッサ１８０５は、別の命令セットを実行することにより、その入力に反応する。各入力は、入力デバイス１８０２、１８０３のうち１つ以上により生成されるデータ、ネットワーク１８２０、１８０２のうち一方を介して外部データ源から受信されるデータ、記憶デバイス１８０６、１８０９のうち一方から検索されるデータ、又は対応する読み取り装置１８１２に挿入された記憶媒体１８２５から検索されるデータを含む多数のデータ源のうち１つ以上から提供され、これらはすべて図１８Ａに示されている。場合によっては、１つの命令セットの実行の結果、データが出力されてもよい。実行は、メモリ１８３４にデータ又は変数を記憶することを含んでもよい。
【００８７】
開示される画像復調構成は入力変数１８５４を使用し、それらの入力変数は、メモリ１８３４の対応する記憶場所１８５５、１８５６、１８５７に記憶される。構成は出力変数１８６１を生成し、それらの出力変数は、メモリ１８３４の対応する記憶場所１８６２、１８６３、１８６４に記憶される。中間変数１８５８は、記憶場所１８５９、１８６０、１８６６及び１８６７に記憶される。
【００８８】
図１８Ｂのプロセッサ１８０５を参照すると、レジスタ１８４４、１８４５、１８４６、演算論理ユニット（ＡＬＵ）１８４０及び制御ユニット１８３９は、プログラム１８３３を構成する命令セットの命令ごとに「フェッチ、デコード、及び実行」サイクルを実行するために必要とされるマイクロ操作のシーケンスを実行するように協働する。各フェッチ、デコード、及び実行サイクルは、
（ｉ）  記憶場所１８２８、１８２９、１８３０から命令１８３１を取り出す又は読み取るフェッチ動作、
（ｉｉ）どの命令が取り出されたかを制御ユニット１８３９が決定する復号動作、
（ｉｉｉ）制御ユニット１８３９及び／又はＡＬＵ１８４０が命令を実行する実行動作、を含む。
【００８９】
その後、次の命令に対して別のフェッチ、デコード、及び実行サイクルが実行される。同様に、制御ユニット１８３９が値を記憶場所１８３２に記憶する又は書き込む記憶サイクルが実行される。
【００９０】
図３〜図１７の処理における各ステップ又は部分処理は、プログラム１８３３の１つ以上のセグメントと関連し、プログラム１８３３の前述のセグメントに関する命令セットの命令ごとにフェッチ、デコード、及び実行サイクルを実行するために協働するプロセッサ１８０５のレジスタ部分１８４４、１８４５、１８４６、ＡＬＵ１８４０及び制御ユニット１８３９により実行される。
【００９１】
あるいは又は更に、画像復調の方法は、画像復調の機能又は部分機能を実行する１つ以上の集積回路のような専用ハードウェアで実現されてもよい。例えばそのような専用ハードウェアは、例えばフーリエ変換などの特殊な復調処理関数を実行するように特に構成されたデジタルシグナルプロセッサ又は１つ以上のＡＳＩＣ、ＦＰＧＡ、マイクロプロセッサ及び関連メモリを含む。
【００９２】
図３は、撮影装置１８９９を含むシステム１８００を使用して実現される撮像及び変調処理３００を示す。処理３００は、ＨＤＤ１８１０に記憶され且つプロセッサ１８０５により実行可能なソフトウェアアプリケーションにより実現される。初期ステップ３１０において、試料１０２の１つ以上のフリンジ画像は、センサ１４０を介して装置１８９９により撮像され、対応する生画像データは、コンピュータモジュール１８０１、通常はＨＤＤ１８１０に記憶される。次に処理３１０は、ステップ３２０に従って、フリンジ画像のうちの１つ以上を復調する。これは、通常、復調されている画像をメモリ１８０６にロードすることを含む。その結果、プロセッサ１８０５は、以下において詳述される種々の復調処理を生画像データに対して実行し、対応する復調画像を形成する。復調画像及び必要に応じて復調処理３２０中に形成された１つ以上の中間画像は、ＨＤＤ１８１０に記憶される。ステップ３３０において、復調結果が出力される。これは、記憶された画像のうちの１つ以上をディスプレイ１８１４上に表示すること、ドライブ１８１２などを介してＣＤ又はＤＶＤなどのポータブル記憶媒体に画像を記憶すること、プリンタ１８１５を介して画像を印刷すること、あるいは遠隔再現及び／又は更なる処理のためにネットワーク１８２０を介して画像を通信することを含む。
【００９３】
概観構成
スペクトルに軸外項が存在することにより復調中にクロストークが発生するため、復調画像の品質は劣化する。しかし、実際のＸ線タルボットシステムにおいて、システムの出力は、式２により与えられた変調モデルに適切に一致しない。この結果、上述の対数マッピングは実際には適切に機能しない。この節で開示されるのは、Ｘ線タルボットシステムからの実験データとより適切に適合する変調処理のより一般的なモデルである。
【００９４】
ｘ変調フリンジ及びｙ変調フリンジと交差項との間の重なり合いからの破損なしでそれらを分離できるようにしつつ、フリンジ模様から交差項を除去する変換を取得することが望ましい。これを達成するために、フーリエ変換の上部象限における単一の交差項を除去するより単純な問題を最初に調査することは妥当である。第１の象限の単純化モデルを考慮する。
【００９５】
【数９】
2016087473.tif 000011
【００９６】
次に、画像Ｑから交差項ｅｉ（Φｘ＋Φｙ）を除去する一般的な非線形関数を考慮する。尚、交差項がフーリエ領域又は空間領域において分離可能でない場合、非線形関数のみが交差項を除去する特性を有する。Ｑ＝Ｑ０において解析的な一般的な非線形関数ｈ（Ｑ）のテーラー級数は、以下により与えられる。
【００９７】
【数１０】
2016087473.tif 000012
【００９８】
式中、係数ｈｊは、関数を完全に規定する。
【００９９】
Ｑ＝Ａに関してｈ（Ｑ）を展開することにより、位相画像Φｘ及びΦｙの高調波における以下の展開が取得される。
【０１００】
【数１１】
2016087473.tif 000013
【０１０１】
式中、Ｃｋ，ｌは、位置（ｋ，ｌ）における交差係数行列Ｃの成分である。一般性を失うことなく、Ａ＝１と仮定すると、最大ｋ≦３及びｌ≦３のエントリを含む式５において規定された関数に対する係数行列は、以下により与えられる。
【０１０２】
【数１２】
2016087473.tif 000014
【０１０３】
この場合、各交差項をゼロに設定することにより、関数ｈの係数を繰り返し求められることが分かる。また、尚、行列は、主対角線に関して以下のような対称関係を有する。
【０１０４】
【数１３】
2016087473.tif 000015
【０１０５】
従って、式中、Ｃｋ，ｌ＝０、そしてＣｌ，ｋ＝０であり、軸外項をゼロにするためには、行列の上部三角形における係数を求めるだけでよい。しかし、係数よりも多くの項が行列のあらゆる有限上部三角形にあるため、すべての軸外項を一意に求められることは明らかではない。しかし、式５の特殊なケースにおいて、すべての軸外項の係数をゼロにするように設定する関数ｈ（Ｑ）に対する級数の係数は、以下のように取得される。
【０１０６】
【数１４】
2016087473.tif 000016
【０１０７】
式中、定義によると、ｈ０＝０である。これは、後でより直接的に導出される予測された対数関係に対応する。
【０１０８】
【数１５】
2016087473.tif 000017
【０１０９】
普通は軸外項の振幅が認識されていないため、通常、式９の対数関係を式２に直接適用できない。更に、式２において、２つの軸外項は異なる振幅である。この場合、式のｘ及びｙの従属部分を厳密に因数分解できない。
【０１１０】
しかし、軸外スペクトルを減少する非線形関数の存在は、復調画像の品質を向上させるためにこの関数を使用するシングルショット画像を復調するアルゴリズムを提案する。復調ステップ３２０の一実施形態として使用される図４の復調処理４００において、アルゴリズムの一例を示す。
【０１１１】
図４において、処理４００への入力４１０は、軸上項及び軸外項を含むフリンジを含む画像、並びに非線形変換のためのパラメータである。処理４００の第１のステップ４２０は、軸外スペクトルを減衰するように入力画像４１０において画素の非線形変換を実行するためのものである。入力パラメータにより、非線形変換を説明する。この非線形変換は対数関数に制限されない。それにもかかわらず、非線形変換は、軸上の高調波の冪を制限することと同様に軸外の相対的冪も減少させる特性を有すべきである。
【０１１２】
軸外項がステップ４２０により減衰された後、復調ステップ４３０において画像を復調するために、あらゆる適切な方法が使用されてもよい。一例は窓フーリエ変換方法である。尚、複数のフリンジの集合を含む画像の復調の結果、複数の復調画像が得られ、所望の軸上フリンジを除いて全く復調する必要はなくなる。
【０１１３】
処理４００の出力は復調画像４４０である。
【０１１４】
ステップ４２０の非線形変換を適用する処理５００を図５に示す。処理５００への入力５１０は、非線形変換ｆ（ｖ，Ｓ）に対する値ｖ及びパラメータＳを含む画素ｐとして示された、軸上項及び軸外項を含むフリンジを含む画像である。
【０１１５】
ステップ５２０において、画像からの未処理画素ｐｉは、値ｖｉで取得される。
【０１１６】
ステップ５３０において、画素ｐｉは、非線形変換（ｖｉ，Ｓ）の評価により置換される。
【０１１７】
ステップ５４０において、処理対象の画素が更にある場合、ステップ５２０に戻る。さもなければ、処理５００は、非線形変換により処理された画像を残してステップ５５０において完了する。
【０１１８】
方法５００は、軸外項を減衰するためにクロストークを減少する。しかし、非線形変換を使用することにより、軸上項の高調波が更に生成される結果、画質を更に低下させる恐れのある外生信号が生じる。
【０１１９】
図４の処理の変形により、復調処理中の軸外項及び高調波の双方を抑制できる。
【０１２０】
図１７は、復調処理を２つの部分に分割する、図４に示されたものとは僅かに異なる復調処理１７００を示す。図４のように、処理１７００への入力１７１０は、復調される軸上項及び軸外項を含むフリンジ画像である。
【０１２１】
ステップ１７２０において、入力画像１７００は、各々が分離されたフリンジ模様を含む複数の実画像を得るように処理される。尚、軸上フリンジの画像以外の画像を生成する必要はない。ステップ１７２０は、軸外スペクトル及び生成された高調波を抑制する処理を組み込むべきであることが望ましい。この処理を実行する方法を後で説明する。
【０１２２】
ステップ１７３０において、分離されたフリンジ画像の各々は、所望の軸上周波数ごとに復調画像を得るように復調される。ここでも、適切な既知の復調方法が使用される。ステップ１７４０において、復調画像が出力される。
【０１２３】
図７Ａ〜図７Ｆ、並びに図８Ａ及び図８Ｂの説明図を参照して、軸上項を分離し、且つステップ１７２０において使用されたような軸外項及び高調波の双方を抑制するようなシングルショット処理６００を図６において説明する。
【０１２４】
処理６００のいくつかのステップが空間領域において適用されるが、図７Ａ〜図７Ｆは、処理が進行する時のフーリエ変換における画像の連続した様相を示す。
【０１２５】
図８Ａ及び図８Ｂは、フーリエ領域においてｘキャリア周波数及びｙキャリア周波数を分離するために使用された２つの２値マスク８１０、８２０を示し、黒色は「１」の値を示し、白色は「０」の値を示す。各マスク８１０、８２０は、直流に近接する低周波数項を保存するために、平面における４つの象限（ｑｕａｄｒａｎｔｓ）及び「１」を含む画像幅の約５％のディスクで構成される。「１」を含む２つの象限は、保存されている軸上項の中心にあり、「０」を含む２つの象限は、除去されている軸上項の中心にある。マスク８１０は、ｙ軸の中心にあるｘ位相導関数に対するキャリアを選択するために使用されたマスクＳｘであり、マスク８２０は、ｙキャリアを選択するために使用されたマスクＳｙである。これらのマスクは、フーリエ領域における鋭角の結果生じるあらゆるリンギングアーチファクトを減少するために使用される前に、ガウシアンフィルタで畳み込まれる。
【０１２６】
処理６００への入力６１０は、単一のフリンジ画像、並びに２つのパラメータ化非線形変換体、即ち第１の変換体ｈ及び第２の変換体ｇである。入力画像のフーリエ変換は、図７Ａの図７１０に類似するように見える。図７１０において、９個のスペクトル領域は、中心の吸収スペクトル（直流項に関して低周波数）に対応し、軸上項及び対応するエルミートの対は、直流項に水平及び垂直に隣接し、軸外項及び対応するエルミートの対は、直流項に対角線状に隣接する。
【０１２７】
ステップ６２０において、第１の非線形変換体ｈは、図５に示された処理を使用して入力画像に適用される。ステップ６２０は、軸外項を減衰し、軸上項の高調波を生成するように動作する。この結果、スペクトルは図７Ｂに示されたような関数７２０に変化する。図７Ｂにおいて、白丸は減衰された軸外項を示し、軸上項に隣接した新しい黒丸は高調波を示す。
【０１２８】
ステップ６３０において、実数高速フーリエ変換（ＦＦＴ）は、空間領域からフーリエ（周波数）領域に画像を変換するように画像に適用される。その結果、エルミート対称性の複素画像が得られる。
【０１２９】
ステップ６４０において、高速フーリエ変換の結果に２つのスペクトルマスクＳｘ（８１０）及びｉＳｙ（８２０）の和を乗算し、式中、ｉ＝√（−１）である。この乗算は、フーリエ変換のスペクトルを実質的にｘキャリア及びｙキャリアを含む領域に分離する。
【０１３０】
ステップ６５０において、複素逆高速フーリエ変換は、画像の実数部におけるｘキャリア及び画像の虚数部におけるｙキャリアで分離された２つのフリンジ模様を生成するように適用される。これは、図７Ｃのスペクトル７３０が画像の実数部を示し、吸収スペクトルが中心円であり、所望の軸上ｘキャリアのエルミートの対が中心円に隣接し、望ましくない高周波が端点に示されて図示される。同様に、図７Ｄのスペクトル７４０は、ｙキャリアを含む画像の虚数部のスペクトルを示す。
【０１３１】
ステップ６６０において、第２の非線形変換体ｇは、高調波を減衰するように各フリンジ模様７３０、７４０に個々に適用される。ステップ６７０において、低周波数スペクトルは、各フリンジ模様から算出及び減算され、ゼロ平均で分離されたフリンジを生成する。そのような低周波数スペクトルは、５次１次元ルジャンドル多項式から別個に算出された２５個の項を適合することにより算出される。その後、この適合は、各フリンジ画像７３０、７４０から減算される。
【０１３２】
結果として得られる２つのフリンジ７５０及び７６０のスペクトルは、図７Ｅ及び図７Ｆにそれぞれ示される。図７Ｅ及び図７Ｆにおいて、軸上キャリアのエルミートの対は依然として存在し、中央の吸収スペクトル（直流成分）は除去され、高周波数高調波は減衰される。
【０１３３】
軸外項、高調波及び吸収スペクトルが除去され、これらの２つのフリンジ画像は、復調に備えてステップ６８０において出力される。
【０１３４】
図１２Ａ〜図１２Ｈは、復調アルゴリズムにおいて使用される種々のスペクトルマスクを示す。図１２Ａ〜図１２Ｈのマスクは、いくつかのスペクトルが除去された図８Ａ及び図８Ｂのマスクの実施形態を示す。ｘキャリア及びｙキャリアに対する２つのスペクトルマスクは、復調方法において使用される。各マスクは、図１２Ａのローパスフィルタ部分１２１０、図１２Ｂの残りの軸外項を除去するためのノッチフィルタ１２２０及び図１２Ｃに示された最も低い周波数を除去するためのハイパスフィルタ１２３０から構成された共通部分を有する。ローパスフィルタと、ノッチフィルタと、ハイパスフィルタとの組み合わせは、マスク１２４０として図１２Ｄに示される。図１２Ｅ及び図１２Ｆのマスク１２５０及び１２６０は、解析信号を生成するためにヒルベルト変換をそれぞれ使用するため、位相がｘ位相導関数及びｙ位相導関数をそれぞれ搬送した状態で、信号を複素象限に復元する。フィルタマスク１２４０をヒルベルトマスク１２５０及び１２６０と合成することにより、ｘキャリアに対する図１２Ｇのマスク１２７０及びｙキャリアに対する図１２Ｈのマスク１２８０である２つの特注のスペクトルマスクが得られる。
【０１３５】
図９は、図１２Ｇ及び図１２Ｈにおいて規定されたマスクを使用して分離されたフリンジ画像を復調する方法９００を説明し、それにより、図１７のステップ１７３０の一実施形態を示す。方法９００は、図６の手順により生成されるような分離されたｘフリンジ画像及びｙフリンジ画像を入力９１０とみなす。次に、ｘフリンジ画像及びｙフリンジ画像の各々に対するフーリエ変換をステップ９２０において計算する。
【０１３６】
ステップ９３０において、各フーリエ変換に対応するヒルベルトマスク及びフィルタマスク１２７０及び１２８０を乗算する。フィルタは、最も低い周波数、残りの軸外スペクトルを除去し、復調処理からノイズ及びエイリアシングを除去するためのローパス効果を有する。ヒルベルトマスクは、フーリエ変換においてエルミートの対を除去して正の周波数のみを残すため、フリンジの象限位相部分を復元する。
【０１３７】
ステップ９４０において、逆フーリエ変換がマスクされたフリンジに適用されることにより、復調画像を空間領域に戻す。
【０１３８】
一対の復調画像は、処理９００のステップ９５０における出力である。これらの画像は、振幅及びラップ位相が直接算出される象限において、実数部及び虚数部の双方を含む。
【０１３９】
撮影システム１００により生成された画像に対するモデルが与えられたとすると、図６のステップ６２０及び６６０において示されたような一致した一対のパラメータ化非線形関数は、第１に軸外項を除去し、第２に第１の非線形関数により生成された高調波を除去するために導出されなければならない。次に、この一致した一対の非線形関数を作成するための２つの方法を説明する。
【０１４０】
第１の方法は、係数の値ｍｘ、ｍｙ、ｍｕ、ｍｖが与えられたとすると、軸外項及び高調波の振幅を最小限にするように一対の非線形関数を生成することである。これらの係数が入力画像にわたり一定値を有する場合、一対の非線形関数は、図６の処理６００の結果として得られる軸外項及び高調波の振幅を最小限にするように導出される。
【０１４１】
第２の方法は、測定されたグラウンドトルースと共に復調するために画像試料を収集する方法である。これにより、変調アーチファクトの非常に正確な推定値Ａ、ｍｘ、ｍｙ、ｍｕ、ｍｖ、ωｘ及びωｙが得られる。次に、最適化処理は、図６及び図９における方法を使用して復調された画像における誤差を最小限にする一対の非線形関数を導出するために使用される。
【０１４２】
図１０は、第１の方法に従って生成された非線形関数を使用して画像を復調する方法１０００を示す。方法１０００は、ＨＤＤ１８１０などに記憶され且つプロセッサ１８０５により実行されるソフトウェアで実現される。
【０１４３】
ステップ１０１０において、４つのフリンジ模様を含む入力フリンジ画像は、撮像されるか、あるいは記憶装置などから取得される。ステップ１０１０において、被検体１０２が存在する必要はない。ステップ１０２０において、フリンジ画像は、全体として全体的に又は各画素位置において局所的にフリンジ画像に対する変調係数に対する値ｍｘ、ｍｙ、ｍｕ、ｍｖを推定するように解析される。これらの値は、格子１０１、１１０、１３０のうちの適切なものを含むが、被検体１０２を含まない撮影システム１００を使用して撮像された校正画像のフーリエ変換において直流及び対応するキャリア周波数ピークの振幅Ａ、Ｆｘ、Ｆｙ、Ｆｕ、Ｆｖを算出し、そしてｍｘ＝Ｆｘ／Ａ、ｍｙ＝Ｆｙ／Ａなどを計算することにより、全体的に推定される。推定された変調係数ｍｘ、ｍｙ、ｍｕ、ｍｖは、Ｘ線タルボット干渉計１００の特性を示す。これらの値は、別の復調アルゴリズムを使用して画像の試験的な復調を実行するプロセッサ１８０５により局所的に算出される。そのようなアルゴリズムの例は、従来のフーリエ方法及び窓フーリエ変換（ＷＦＴ）を含む。ステップ１０３０において、非線形変換のためのパラメータは、軸外項及び高調波を減少するようにステップ１０２０において推定された変調係数ｍｘ、ｍｙ、ｍｕ、ｍｖから算出される。ステップ１０４０において、画像は、好ましくは図６及び図９において説明された方法を使用して復調される。
【０１４４】
実際のＸ線タルボット撮影システムにおいて、変調処理は、複雑であり、単純な理論モデルに一致しない。従って、軸外項を減少するために使用された非線形関数のパラメータを選択するために、測定された変調パラメータ及び振幅パラメータを使用するのは困難である。特に、変調係数の振幅は一定ではなく、値は、Ｘ線が逸れず且つスペクトルの帯域幅が狭い場合に最大となり、Ｘ線が例えば種々の材料間のインタフェースにおいて被検体により偏向される場合にゼロに向かって減少する。これらの位置において、変調キャリアの帯域幅は最大である。極端な場合、変調係数は負になる。
【０１４５】
ｍｘ、ｍｙ、ｍｕ、ｍｖの値は、急な位相勾配がある状態で局所的に低下する。例えば、ｘ方向の位相勾配の大きさが大きい場合、ｍｘの値は、低下し、極端に急な勾配の場合にはゼロに向かって低下する。他の変調値についても同じであり、軸外項ｍｕ及びｍｖは、１３５°及び４５°のそれぞれにおいて位相勾配と関連付けられる。
【０１４６】
一対の非線形関数は、全体的に又は局所的に選択される。関数は、全体的に選択される場合、画像におけるクロストーク歪みに最も影響を及ぼすと予想される変調係数に対する軸外項を減少するように選択されるべきである。位相勾配と関連付けられない信号の帯域幅が極端に狭いが、急な位相勾配と関連付けられた高帯域信号の振幅も極端に小さいため、これらの項の選択は明らかではない。従って、このような理解に基づいて、非線形関数は、中程度の変調値に対する軸外項を減少するように選択されるべきであることが予想されるだろう。
【０１４７】
非線形関数が局所的な変調値に基づいて選択される場合、これらの変調値は、復調前に使用可能でないため、より単純な方法により推定されなければならない。また、局所的に変化する非線形関数を使用する場合、フリンジを回復するためのこれらの関数の反転は、非自明であり、慎重に考慮されなければならない。
【０１４８】
変調値を正確に測定するのは困難であるため、最適化処理を使用して復調のために使用される非線形関数を選択する別の方法を説明する。
【０１４９】
そのような方法における第１のステップは、既知のグラウンドトルースで１つ以上の変調画像を生成するステップである。複数の位相シフト画像、例えばすべての倍数９０°でｘ方向及びｙ方向に位相シフトした結果得られる１６個の画像をとり、その後参照復調画像を生成するために既知の復調方法を使用することにより、候補画像のすべての関連パラメータを取得する。復調パラメータは、画素ごとに個々に測定され、Ａ、ｍｘ、ｍｙ、ｍｕ、ｍｖ、ωｘ及びωｙである。この解析ではノイズが無視されるため、画像は、ノイズが信号の顕著な成分ではないように十分な露出時間を有するべきである。
【０１５０】
シングルショット復調処理の評価のためのシングルショット画像を取得するために、２つの方法が使用される。第１に、単一の入力画像は、複数の位相シフト入力画像から取得される。第２に、単一の入力画像は、参照復調画像の再変調によりシミュレートされる。いずれかの方法を使用して復調されたこのシングルショット復調画像は、候補復調画像と呼ばれ、復調品質を評価するために参照復調画像と比較される。
【０１５１】
候補復調画像を取得するためのこの第１の方法は、本物のデータの復調がテストされるが、この復調における誤差の源が、復調誤差＋信号ノイズ及び変調モデルにおけるあらゆる誤差であるという欠点を有する。第２の方法により、ノイズがなく且つ復調モデルにより完全に説明される画像を作成できるため、復調における誤差の源は、完全に復調アルゴリズムの結果とみなされる。
【０１５２】
取得された画像が非線形関数に対する候補パラメータが使用されているシングルショットアルゴリズムを使用して復調される場合、候補復調画像の品質、即ちパラメータの品質は、位相の二乗平均平方根（ＲＭＳ）誤差などの画質計測値を利用して評価される。
【０１５３】
この処理を１回実行した後、パラメータは変形し、処理は繰り返された。このように、最小の二乗平均平方根（ＲＭＳ）誤差を提供するパラメータ化変換の最適なパラメータが取得される。
【０１５４】
そのような処理の一例１１００を図１１に示す。処理１１００への入力１１１０は、位相ステップ参照フリンジ画像の集合である。通常、入力画像１１１０は、工場設定及び撮影装置１００の校正中などの制御された設定環境において撮像される。画像１１１０は、３つの軸又は回転に沿って線形オフセット（シフト）を含む格子１０１、１１０及び１３０の複数の位相変動において撮像される。ステップ１１２０において、入力画像１１１０の集合は、参照復調画像１１３０を形成するために高品質マルチショット又は位相ステップアルゴリズムを使用して復調される。ステップ１１３５において、単一の画像１１１０５は、シングルショット復調処理への入力を形成するように入力画像１１１０の集合から選択される。集合１１１０からのあらゆる単一の画像は、シングルショット処理のために選択される。
【０１５５】
ステップ１１４０において、プロセッサ１８０５は、シングルショット復調処理のための候補パラメータを識別するように動作する。これらのパラメータは、復調処理において使用される非線形関数を説明するパラメータを含むが、それらに限定されない。通常、パラメータは、ＨＤＤ１８１０などのメモリに記録される初期の推定値を含む。初期の推定値は、工場設定として予め決定されるか、あるいは又は更に処理１１００の前のインスタンス化時にメモリに保持される。
【０１５６】
ステップ１１６０において、選択された１つの画像１１０５は、図６において説明されたようなシングルショット復調処理６００に従ってステップ１１４０において識別されたパラメータを使用して復調され、その結果、候補復調画像１１５０が得られる。
【０１５７】
ステップ１１７０において、候補復調画像１１５０は、参照復調画像１１３０と比較され、誤差推定値１１７５を形成する。
【０１５８】
処理１１００の最後のステップ１１８０において、プロセッサ１８０５は、ステップ１１４０において識別された非線形変換を含む復調処理のパラメータを調整するために誤差推定値１１７５を使用し、それにより、調整された復調パラメータ及び非線形変換が提供される。
【０１５９】
シンプレックス法、力まかせ探索又は勾配下降方法などの多数の最適化手順のうちのいずれかを使用することにより、図１１における処理は、折り返し破線１１９０により示されるように、入力となるステップ１１８０の調整されたパラメータをシングルショット復調処理１１６０に戻すことにより繰り返されることが望ましい。選択された最適化手順の結果、誤差推定値１１７５が減少又は改善されるため、復調アルゴリズム１１６０のために改善されたパラメータが得られる。誤差が所定の数の反復に対する所定の値を下回るまで又は誤差の向上が停止するまで、繰返しを継続する。
【０１６０】
第１の実施形態
軸外項を減少するために非線形関数を使用することは、変調パラメータｍｘ、ｍｙ、ｍｕ及びｍｖの相対振幅に決定的に依存する。
【０１６１】
【数１６】
2016087473.tif 000018
【０１６２】
軸外項が、積、即ち、式１０から形成される場合に軸上項振幅の積の振幅と同一になる振幅を有する特殊なケースにおいて、式２からのノイズＳｉが無視されている。この場合、軸外項は、乗法交差項に直接対応し、第１の非線形関数がｆ（Ｉ）＝ｌｏｇＩの状態で除去され、高調波は、第２の非線形関数がｇ（Ｉ）＝ｅｘｐＩの状態で除去される。これにより、以下のような第１の非線形関数の後に理論応答が得られる。
【０１６３】
【数１７】
2016087473.tif 000019
【０１６４】
この特殊なケースは、準同型フィルタリングの使用に対応する。しかし、ｘ変調及びｙ変調の完全な積である信号は、実際に発生する可能性は低い。これは、軸外項を誘発する物理処理が２つの１次元余弦関数の積に対応しないためである。
【０１６５】
Ｘ線タルボット撮影システムは、このモデルに適合する画像を生成しない。代わりに、軸外変調は、通常、式２により与えられるような軸上因子の積とは異なる。変調の主軸が既知であるか、あるいは決定されると仮定すると、主軸がｘ方向及びｙ方向にあるように座標をどこにマッピングするかの問題を考慮できる。また、下付き文字ｕ及びｖで示された２つの更なる対角線の交差項が考慮され、交差項の周波数及び位相は、以下のように単に軸上位相項φｘ、φｙの和及び差分である。
【０１６６】
【数１８】
2016087473.tif 000020
【０１６７】
この場合、画像は以下の通り因数分解される。
【０１６８】
【数１９】
2016087473.tif 000021
【０１６９】
式中、表記φｘ＝ωｘｘ＋Φｘ（ｒ）及びφｙ＝ωｙｙ＋Φｙ（ｒ）は便宜上使用され、因数分解係数は以下の関係に従って推定される。
【０１７０】
【数２０】
2016087473.tif 000022
【０１７１】
式１３の計算には、空間的に変化する吸収Ａ（ｒ）及びパラメータγの推定値を認識する必要がある。軸外項が軸上項の厳密な積である理想のケースにおいて、γ＝１である。しかし、通常、交差項が予想より大きいため、γ＜１である。更にこの場合、Ａ（１−γ）が画像Ｉから減算されているため、画像ノイズ、並びにＡ及びγの推定値により、この項が負になる可能性がある。実ログ関数がゼロ値又は負の値に対して規定されないため、関連項は、小さな正の値の下に制限されるべきである。また、交差項が予想より小さい場合、γ＞１であり、正の値が信号Ｉに追加されるため、結果は常にゼロより大きい。
【０１７２】
吸収は、画像から変調キャリア周波数を除去するストップバンドを有するローパスフィルタを使用して推定される。殆どの場合、より高い周波数が、通常、変調が低く且つｍｕ＝ｍｖという仮定が崩れるエッジにあるため、吸収に対する帯域制限された近似値でさえ有用である。図１６は、吸収画像Ａを推定するための処理１６００を示す。
【０１７３】
また、因数分解係数γは、オフセットログ関数が計算される前に推定される必要がある。これには、すべての振幅変調ｍｘ、ｍｙ、ｍｕ及びｍｖの推定値が必要である。これらの因子の一定の推定値は、各キャリアのそれぞれの中心周波数に近接するフーリエ領域においてピーク信号振幅を算出することにより形成される。あるいは、これらの量の低周波数推定値は、フーリエ領域において元の信号をフィルタリングし、且つ結果として得られる信号を復調することにより算出される。
【０１７４】
実験システムの場合、交差項振幅変調は、通常、軸上変調の積の大きさの２倍である。これは、格子１０１、１１０及び１３０が構成されるセルの形状の結果である。実際に、交差項と軸上項との間には乗算関係はない。
【０１７５】
変調関数ｍｘ、ｍｙ、ｍｕ及びｍｖが画像にわたり一定値を有さないことも重要である。被検体１０２の吸収は、画像センサ１４０に到達する光（特に、Ｘ線電磁エネルギー）の量を単に減少し、式４のすべての項を乗算した吸収項Ａにより考慮される。しかし、振幅変調は、吸収より大きな小数だけ低下する。更に、変調値も撮影されている構造のエッジにおいて低下し、紙及び木などのいくつかの材料により、それらの相対的に小さな吸収値に対して不適切な変調の低下が生じる。これは、Ｇ０格子１０１におけるピンホールの集合によるタルボット画像の畳み込みのためであり、このソース畳み込みの初期の１次元モデリングは、実数データにおいて観察される同一の効果のうちの多くを再現している。この畳み込みの結果、光路長の２次導関数が高い場合、変調は、場合によっては９０％も低下する。これは、材料間の境界、並びに高散乱構造を有する木及び紙などの構造において発生する。この特性は、有用であり、画像における鋭角の位置を決定するために使用される。
【０１７６】
また、変調は、光路長における勾配の方向と変調方向との間の内積に比例して低下する。これは、光路長勾配がπ／４の角度の場合に交差項変調ｍｕが高く、光路長勾配が−π／４の角度の場合にｍｖが高いことを意味する。従って、これらの項は、ｍｘ及びｍｙの積であることから大幅に外れる。
【０１７７】
実際のシステムにおけるこれらのパラメータのいくつかの通常の平均値は、以下の通りである。
【０１７８】
【数２１】
2016087473.tif 000023
【０１７９】
式中、Ｘ−≡平均（ｘ）である。
【０１８０】
この形式において、軸外項は、オフセットログ関数を使用して除去される。
【０１８１】
【数２２】
2016087473.tif 000024
【０１８２】
式１５のオフセットログ関数が適用された後、軸外交差項が大幅に減少されていることが分かるため、軸上項はフーリエ領域において分離される。次に、フーリエ領域において変調されたｘ項及びｙ項は、１つの軸に沿う軸上項を優先的に受け入れるフィルタを使用して分離される。フーリエ領域におけるそのようなフィルタの一例を図８Ａ及び図８Ｂに図示する。
【０１８３】
このフィルタが適用された後、残りの信号は、優先的に１つの軸の軸上項及び総信号の直流成分から構成される。尚、信号の直流成分は、フーリエ領域において除外されない他の軸変調信号からの成分を含む。この他の軸成分の寄与度は、γに関する対数項の級数展開をとることにより、以下の関係に従って推定される。
【０１８４】
【数２３】
2016087473.tif 000025
【０１８５】
式中、ｗ＝（ｘ，ｙ）であり、項の係数は以下により与えられる。
【０１８６】
【数２４】
2016087473.tif 000026
【０１８７】
従って、ｘ変調信号及び総信号（撮像画像）の直流成分を抽出するために図８Ｂのｘ分離フィルタ８２０を適用した後、そのような信号は、処理済みフリンジ信号として示される。
【０１８８】
【数２５】
2016087473.tif 000027
【０１８９】
尚、ｍｘ／γ＜１の場合、この和は収束する。交差項が軸上項に対して小さいと仮定すると、展開における第１の項のみを考慮することが可能である。フィルタリングされた対数の指数を取得することにより、以下が与えられる。
【０１９０】
【数２６】
2016087473.tif 000028
【０１９１】
式１９は、吸収及びｙ振幅変調の影響がある状態を除いて、ｘ変調信号をｙ変調信号から正常に分離できることを示す。
【０１９２】
吸収を除去するために、式は、事前に算出されたような調整された推定吸収により分割される。ｍｙ／γが小さいと仮定することにより、ｘ変調信号を推定できるため、以下の通りである。
【０１９３】
【数２７】
2016087473.tif 000029
【０１９４】
ｙ変調信号に対する同様の式は、以下の通りである。
【０１９５】
【数２８】
2016087473.tif 000030
【０１９６】
変調信号の位相成分及び振幅成分を検索するために、負の周波数成分を除去することにより解析信号を形成できる。解析信号は、ＦＦＴ又は１次元（１Ｄ）空間フィルタを使用することに対して算出される。ローパスフィルタは、交差項周波数において残りのフーリエ成分を除去するように更に適用される。これは、以下のように図１２Ａ〜図１２Ｈに示されるフーリエ領域フィルタを使用することにより達成される。
【０１９７】
【数２９】
2016087473.tif 000031
【０１９８】
これらのフィルタリングされた信号Ｔｘ及びＴｙから、振幅変調ｍｘ及びｍｙの推定は、複素信号の実数成分及び虚数成分を含むアークタンジェントを使用して係数、ラップ位相（モジュロ２π）から決定される。
【０１９９】
式３の真の微分位相Φｘ及びΦｙの推定値を取得するために、画像Ｔｘ及びＴｙにおけるラップ位相は、ある種の位相アンラッピングアルゴリズムを使用してアンラップされなければならない。
【０２００】
Ｘ線位相コントラスト撮影において、種々の構造間の境界において発生する極端に大きな相転移がある。領域における明らかな位相勾配が画素ごとのπラジアンのナイキストサンプリング限界より大きい画像の領域において、正の位相勾配と負の位相勾配との間で曖昧さが発生する。従って、位相アンラッピングはこれらの境界領域において困難である。この問題のため、真の大きさ及び測定された位相の符号をラップ位相値から推定するのは非常に困難である。この問題を軽減するために、振幅変調又は吸収値（例えば、吸収画像Ａ）などの追加情報が、位相アンラッピングを誘導するために使用される。
【０２０１】
例えば、システムにおいて使用された格子に関連した光学効果のため、振幅変調は、位相勾配が大きい画像の領域において大幅に低下する。従って、この振幅変調は、位相が画素ごとのπラジアンのナイキストサンプリング限界より大きなレートで変化する画像の領域を決定するために、ガイドとして使用される。変調値及び周囲の位相値を調査することにより、位相勾配が最大値又は最低値に到達する厳密な位置、並びに値の推定を可能にするため、大きな位相勾配を含むエッジ領域における正確な位相アンラッピングが可能になる。
【０２０２】
他の実施形態において、吸収画像Ａが使用される。例えば、大きな相転移を生じさせる種々の構造間の境界がある場合、吸収画像においてエッジもあることが多い。位相導関数がアンラップされる方向に吸収画像の導関数を調査することにより、位相勾配の符号が推定されるため、位相アンラッピング処理を支援する。
【０２０３】
振幅変調の推定を向上させるために、（ｊ＋１）番目の推定値に対する分離されたフリンジがｘ振幅変調ｍ（ｊ）ｘ及びｙ振幅変調ｍ（ｊ）ｙの前の推定値を使用して算出される式１９における対立軸振幅変調の影響を補正するために反復手順が使用されることにより、式２０及び２１を以下で置換する。
【０２０４】
【数３０】
2016087473.tif 000032
【０２０５】
次に、更新された振幅変調は、上述したような復調処理を使用して計算される。
【０２０６】
図１３は、復調処理中に使用される非線形関数に対するパラメータを推定する処理１３００を示す。処理１３００は、ＨＤＤ１８１０に記憶され且つプロセッサ１８０５により実行可能なソフトウェアで実現されることが好ましい。処理１３００は、被検体１０２が存在せずにシステム１００を使用して撮像された画像である１つ以上の校正画像１３１０で動作する。処理１３００は、ステップ１３２０において最初に入力画像１３１０を復調し、校正復調画像１３２５を形成する。ステップ１３２０の復調は、参照復調画像を生成する目的のために知られている上述したようなマルチショット復調方法を含むあらゆる適切な方法であってもよい。ステップ１３３０において、プロセッサ１８０５は、復調画像１３２５において軸上成分及び軸外成分（即ち、交差項）の代表的な振幅を測定するように動作する。これは、スペクトル領域における画像のピークを測定することにより実行される。プロセッサ１８０５が、ステップ１３３０から取得された測定値を用いて、処理１３００の出力１３５０を形成する非線形変換に対するパラメータを算出するように動作するステップ１３４０が後続する。これらの計算は、式１３及び式１４を参照して上述したような振幅ｍｘ、ｍｙ、ｍｕ、ｍｖに関して、例えば式１５のログオフセット関数を実現する。
【０２０７】
次に、図１６を参照すると、処理１６００に対する入力はフリンジ画像１６１０である。吸収画像を推定するために、変調キャリアが除去されなければならず、これはローパスフィルタにより実行される。すべての変調キャリア周波数を完全に除去するために、ローパスフィルタは、これらの周波数の振幅をゼロに設定できるフーリエ領域において実現される。ステップ１６２０において、フーリエ変換は、ローパスフィルタリングを実行するように、プロセッサ１８０５により画像１６１０に適用される。ステップ１６３０において、フーリエ変換に窓を乗算し、キャリア周波数を除去する。ステップ１６３０において使用される窓の厳密な幅及び関数形式は、特に重要ではない。軸上キャリア周波数間の距離の半分の幅の円形対称ハニング窓（ｃｉｒｃｕｌａｒ−ｓｙｍｍｅｔｒｉｃ  Ｈａｎｎｉｎｇ  Ｗｉｎｄｏｗ）は、本発明者により発明されており、良好な結果をもたらす。推定吸収画像１６５０として出力される空間領域におけるローパスフィルタリング画像を作成するように、ステップ１６４０において逆フーリエ変換が適用される。
【０２０８】
図１４は、ログオフセット変換を使用してフリンジ画像から軸外項が削除されることにより、一対の分離されたフリンジ模様を生成する方法１４００を示す。方法１４００は、図６において説明された処理６００一実施形態である。
【０２０９】
図１４において、処理１４００への入力１４１０は、パラメータｍｘ、ｍｙ、ｍｕ及びｍｖに対する推定値と共に、軸上項及び軸外項の双方を含むフリンジ画像である。入力画像１４１０の単位は、画像センサ１４０からの画素カウントの何らかの倍数であり、整数であり、正又はゼロである。画像１４１０の高品質な領域は、１００カウントより大きな値を有する。ｍｘ、ｍｙの対は同様の値を有し、ｍｕ、ｍｖの対は同様の値を有するが、互いに異なると仮定する。上述したように、これらのパラメータは、フーリエ変換において信号ピークの振幅を測定することにより推定される。
【０２１０】
ステップ１４１１において、吸収画像Ａはフリンジ画像１４１０に対して推定される。これは、図１６において説明されたような方法１６００に従って実行されることが好ましい。
【０２１１】
吸収画像Ａは、範囲の下端において０．０（検出された光を示さない）から範囲の上端において数百又は数千のカウントの範囲の画素値を有する。ステップ１４１２において、吸収画像Ａは、後の分割中のゼロ除算誤差を回避するために１．０を下回る値が１．０により置換されるようにクリップされる。
【０２１２】
ステップ１４１２において、フリンジ画像Ｉは、１．０の平均値及び算出されたような変調値、例えばｍ−ｘ＝ｍ−ｙ＝０．１；  ｍ−ｕ＝ｍ−ｖ＝０．０３を含む画像を生成するように、クリップされた吸収画像Ａにより更に分割される。
【０２１３】
これらの変調値から、ステップ１４１３においてオフセット値が算出され、γ＝ｍｘｘｍｙ／（Ｍｕ＋ｍｖ）である。
【０２１４】
オフセット値γにより、入力画像１４１０のログオフセット変換画像１４２５をステップ１４２０においてプロセッサ１８０５により以下のように算出できる。
【０２１５】
【数３１】
2016087473.tif 000033
【０２１６】
ログ関数を計算する前に、ログ演算の計算の結果として得られる極端な負の値を防止するために、ログへの引数は、０．０１などの何らかの最小値にクリップされることが望ましい。
【０２１７】
ステップ１４３０において、フーリエ変換は、ログオフセット変換画像１４２５をスペクトル領域に変換するようにプロセッサ１８０５により計算される。結果として得られるスペクトル画像１４３５において、軸外項は実質的に減少されるが、高調波は、ログ関数の非線形性により軸上項に導入される。
【０２１８】
ステップ１４４０において、図８において上述したように、フーリエ変換画像１４３５に複素分離マスクを乗算して画像１４４５を得る。これにより、画像１４２５に存在するｙ項及びｘ項を除外する。
【０２１９】
ステップ１４５０において、逆フーリエ変換は、画像１４４５を再度空間領域に変換することにより、実数部及び虚数部を有する空間画像１４５５を得るようにプロセッサ１８０５により計算される。部は、別個の画像として処理されるように分離される。特に、逆変換の際、フリンジ模様は、以下の形式になるように、画像１４４５の実数部からのｘキャリアフリンジ（ｘ位相）及び画像１４４５の虚数部からのｙキャリアフリンジ（ｙ位相）に分離される。
【０２２０】
【数３２】
2016087473.tif 000034
【０２２１】
【数３３】
2016087473.tif 000035
【０２２２】
ステップ１４６０において、これらの位相画像の各々は、以下のように第２の非線形関数により変換される。
【０２２３】
【数３４】
2016087473.tif 000036
【０２２４】
この結果、各々が軸上キャリア及び一定のオフセット項のみを含む２つの画像１４６５が得られる。
【０２２５】
ステップ１４７０において、各画像１４６５の低周波数スペクトルは、各画像の低度の多項フィット、例えば５次２Ｄ多項式を減算することにより除去され、それにより信号をゼロ平均にする。
【０２２６】
それにより、処理１４００は、復調に備えてゼロ平均で分離された一対のフリンジ模様を含む出力１４８０を生成する。
【０２２７】
復調は、図９において上述したように実行される。６画素の期間を含むフリンジ模様に対して、本発明者は、シンプレックス最適化を通して復調動作のために以下の幅を選択している。ローパスフィルタ１２１０の場合、ローパスフィルタの幅はフーリエ変換の幅の５２％であり、平滑化領域は、ローパスフィルタディスクの幅の３０％だけ更に拡張するハニング窓から構成される。同様に、軸外スペクトルカットフィルタ１２２０の幅は、画像幅の３．８％であり、平滑化領域を有さない。直流カットフィルタ１２３０は総カット領域を有さず、平滑化領域の幅は画像幅の１０％である。
【０２２８】
第２の実施形態
ログオフセット変換により、原則的に軸外項によるクロストークを削除できるが、本発明者などによれば、復調品質を制限するこの方法によるいくつかの制限がある。第１に、ログ関数がゼロに向かう傾向のある正規化された画像値として解かれるため、ログ関数を計算する前に以下において画像データをクリップする必要がある。クリッピングにより、実質的に誤差が発生する。第２に、軸外項の変調値が画像を通して同一でないことが多いため、変調値の推定により、通常、誤差が発生する。第３に、吸収を厳密に推定できない。これらの問題の結果、復調画像において望ましくないアーチファクトが発生する。
【０２２９】
次に、交差項を除去するために対数関数を使用しないより一般的な処理方式を考慮する。この処理は、それぞれ元の信号の軸外項を減少し且つ軸上項を復元する２つの一般的な非線形関数Ｈ及びＧを使用する。
【０２３０】
使用される第１の非線形変換は、交差項を最小値に減少できる。この関数の他の望ましい特性は、広範囲の潜在的な軸外振幅及び軸上振幅を受け入れられること、全体的な乗法スケーリングに対して不変であること、及び明白に反転することである。オフセットログ変換以外のすべての軸外項を厳密に消去する既知の非線形関数はないが、これらの項を実質的に減少する非線形関数を導出できる。
【０２３１】
特に、式５のモデルケースが考慮され、以下の通り示される。
【０２３２】
【数３５】
2016087473.tif 000037
【０２３３】
ｖ＝１において有限多項級数解析により規定された一般的な非線形関数ｈ（ｖ）は、以下により与えられる。
【０２３４】
【数３６】
2016087473.tif 000038
【０２３５】
式中、係数ｈｊは関数を完全に規定する。
【０２３６】
ｖ＝１に関してｈ（ｖ）を展開することにより、位相画像Φｘ及びΦｙの高調波における以下の展開が取得される。
【０２３７】
【数３７】
2016087473.tif 000039
【０２３８】
式中、Ｃｋ，ｌは、位置（ｋ，ｌ）における交差係数行列Ｃの成分である。この場合、式６の級数展開におけるすべての係数が係数行列の各項において累乗されるため、無限冪展開の項を繰り返し求めることはできない。尚、Ｑ（ｘ，ｙ）が実数であるため、行列は対称性Ｃｋ，ｌ＝Ｃ−ｋ，−ｌを有する。
【０２３９】
次に、一般的な非線形関数の交差項を最小限にする問題を考慮する。Φｘ＋Φｙ及びΦｘ−Φｙにおける基本的な軸外項を考慮すると、振幅Ｃ２０，１及びＣ２１，０の二乗振幅により測定された双方の軸上項における冪に対して、振幅Ｃ２１，１及びＣ２−１，１の二乗振幅により測定されたこれらの軸外項における冪を測定するＲｏａを最小限にすることが望ましい。そのような測定は、以下により与えられる。
【０２４０】
【数３８】
2016087473.tif 000040
【０２４１】
式２７の冪比を最小限にする非線形関数ｈｋのパラメータをシークすることにより、軸外項を抑制するために使用する潜在的な非線形関数が得られる。あるいは、高次軸外項を最小限にする関数が適していてもよい。
【０２４２】
【数３９】
2016087473.tif 000041
【０２４３】
例えば、２次非線形関数Ｎ＝２の最も単純なケースを考慮すると、パラメータｈ１及びｈ２の関数として式２７の最小値を取得し、定数ｈ０＝０を任意に設定できる。説明を簡略化するために、ａ＝ｂと仮定すると、以下のパラメータが取得される。
【０２４４】
【数４０】
2016087473.tif 000042
【０２４５】
ｃ＝ｄの場合、ｈ２＝−ｃｈ１／（２ａ２）を含む交差項を厳密に消去できる。三次関数を考慮すると、第１の軸外高調波をゼロに設定できる。しかし、これが実行される場合、高次軸外項は、二次非線形関数より大きくなる。従って、式２８を使用して総軸外冪を減少することが好ましい。更に、解析は三次より大きい次数に拡張され、交差項を最適なレベルに順次減少する。
【０２４６】
先行する解析により、軸外交差項を減少するための前方変換を構成できる。第２の後方変換を構成するために、第１の変換を直接反転させること又は画質に関連する異なる計測値に対して別個に最適化される関数を構成することが可能である。一例として、後方変換の計算は、軸上項のみの反転として導出される。尚、双方のケースにおいて、後方関数ｇ（ｖ）は、前方関数ｈ（ｘ）の反転に等しくなくてもよい。
【０２４７】
軸上項の後方関数ｇ（ｖ）は、以下により規定された形式的冪級数として計算される。
【０２４８】
【数４１】
2016087473.tif 000043
【０２４９】
式中、係数ｇｍは、多くの周知のそのようなアルゴリズムがある中の反復的アルゴリズムを使用して、軸上係数Ｃｍ，０及びＣ０，ｍのそれぞれから算出される。ｙ軸上項を復元するための関数も同様に導出される。
【０２５０】
二次変換の場合、逆非線形関数は、以下のようにこのインスタンスにおいて直接計算される。
【０２５１】
【数４２】
2016087473.tif 000044
【０２５２】
しかし、この２次関数には多数の欠点がある。逆関数は、平方根関数の分岐が正確に決定されるようにｈ２が十分に小さいことを必要とする。これには、元のモデルのパラメータ、特に軸上冪及び軸外冪が、モデルが要求する限界を超えないことを保証するようにそれらを監視することが必要であった。
【０２５３】
より望ましい関数は、広範な潜在的な軸外振幅及び軸上振幅を受け入れる特性と、全体的な乗法スケーリングに対して不変である特性と、明白に反転する特性とを有する。このため、これらの特性を有する非線形関数の他の類を考慮するのが妥当である。一例は、以下の冪関数である。
【０２５４】
【数４３】
2016087473.tif 000045
【０２５５】
式中、αは整数ではない。使用に適した逆関数は、以下の通りである。
【０２５６】
【数４４】
2016087473.tif 000046
【０２５７】
式中、βは１／αと異なってもよい。
【０２５８】
第１の非線形変換のパラメータαは、実画像上で関数を使用する前に決定される。第１の線形関数で処理した後、ログオフセット変換に対して上述したように、軸上項は、フーリエ領域においてフィルタを使用してフィルタリングされる。
【０２５９】
第１の非線形変換は、軸外項を減少するが、高調波を軸上項に導入する。従って、第２の非線形関数の目的は、最小のアーチファクトで画像を復調できるように、これらの高調波を減少することである。第２の非線形変換のパラメータβは、処理の再構成誤差を最小限にすることに基づいて、実際の復調前に再度決定される。最後に、画像はフィルタリング及び復調される。
【０２６０】
受信したフリンジ模様を復調するために、元の信号の軸外項を減少し且つＲＭＳ再構成誤差を最小限にする２つの一般的な非線形関数ｈ及びｇを使用することを図１５に示す。図１５は、図６において説明された処理６００の例示的な一実施形態を示す。
【０２６１】
図１５は、ソフトウェアで実現され、ＨＤＤ１８１０に記憶され、且つプロセッサ１８０５により実行可能であることが好ましいフリンジ模様を復調する方法１５００を示すフローチャートである。方法１５００への入力は、フリンジ画像１５１０、並びに２つのパラメータ化非線形変換体１５１５、即ち第１の変換体ｈ及び第２の変換体ｇである。画像１５１０の単位は、画像センサ１４０からのフォトンカウントの何らかの倍数であり、整数であり、正又はゼロである。一般に画質は、フォトンカウントと共に向上し、画像の高品質領域は、１００カウントより大きい値を有する。変換体１５１５は、通常、ＨＤＤ１８１０に記憶された所定の変換体から入手される。画像１５１０は、センサ１４０から直接取り込まれてもよく、あるいはＨＤＤ１８１０における記憶装置から又はネットワーク１８２０を介して検索されてもよい。
【０２６２】
ステップ１５２０において、第１の非線形冪変換は、入力画像１５１０に適用される。画像への非線形変換体（ｈ）の適用は、式３２において与えられた非線形関数ｆ＝ｈに対して図５に示される。この第１の非線形変換のパラメータαは、図１１に示されるように、軸外項及び／又は参照画像と復調画像との間の誤差を最小限にするパラメータを決定する最適化処理を使用して決定されている。ステップ１５２０の結果、第１の変換画像１５２５が得られる。
【０２６３】
ステップ１５３０において、画像１５２５の離散（高速）フーリエ変換はプロセッサ１８０５により計算され、スペクトル画像１５３５を形成する。スペクトル画像１５３５において、軸外項が実質的に減少されるが、高調波は、第１の非線形関数の非線形性により軸上項に導入される。
【０２６４】
ステップ１５４０において、図８Ａ及び図８Ｂにおいて上述したように、フーリエ変換画像１５３５に複素分離マスクを乗算して、マスク複素画像１５４５を得る。
【０２６５】
ステップ１５５０において、プロセッサ１８０５は、逆フーリエ変換を計算し、マスク画像１５４５を再度空間領域に変換して画像１５５５を形成する。画像１５５５において、フリンジ模様は、ステップ１４５０と同様に、画像の実部においてｘキャリアフリンジ及び画像の虚部においてｙキャリアフリンジに分離される。
【０２６６】
ステップ１５６０において、画像１５５５の部分である分離された画像の各々は、第２の非線形冪関数（ｇ）により変換される。画像への非線形変換の適用は、式３３において与えられた非線形関数ｆ＝ｇに対して図５において示される。この第２の非線形変換のパラメータβは、図１１に示されるように、参照画像と復調画像との間の誤差を最小限にするパラメータを決定する最適化処理を使用して事前に決定されている。ステップ１５６０は、画像から高調波を除去するように動作する。
【０２６７】
ステップ１５７０において、低周波数スペクトルは、各部分画像の低度の多項フィット、例えば５次２Ｄ多項式を減算することにより画像１５５５から除去され、それにより信号をゼロ平均にする。
【０２６８】
方法１５００の出力は、復調に備えてゼロ平均で分離された一対のフリンジ模様１５８０である。
【０２６９】
復調は、図９において上述されたように実行されることが好ましい。６画素の期間を含むフリンジ模様に対して、シンプレックス最適化を通して復調動作のために以下の幅が選択されている。ローパスフィルタ１２１０の場合、ローパスフィルタの幅はフーリエ変換の幅の５２％であり、平滑化領域は、ローパスフィルタディスクの幅の３０％だけ更に拡張するハニング窓から構成される。同様に、軸外スペクトルカットフィルタ１２２０の幅は、画像幅の３．８％であり、平滑化領域を有さない。直流カットフィルタ１２３０は総カット領域を有さず、平滑化領域の幅は画像幅の１０％である。
【０２７０】
２つの非線形関数により、事前に最適化されているパラメータが得られる。図１１は、この最適化処理を示し、上述されている。この第２の実施形態において、非線形関数のために推定される必要があるパラメータは、第１の非線形関数に対してαであり、第２の非線形関数に対してβである。２つのパラメータは、候補復調画像１１５０と参照復調画像１１３０との間の平均二乗誤差（ＭＳＥ）又は他の誤差測度を最小限にするように選択される。
【産業上の利用可能性】
【０２７１】
説明した構成は、生体構造及び他の構造から詳細を得るために、コンピュータ及びデータ処理の分野、特にフリンジ画像の復調の分野に適用可能である。
【０２７２】
以上、本発明の他のいくつかの実施形態のみを説明したが、本発明の範囲及び趣旨から逸脱することなく、それらの実施形態に対して変形及び／又は変更を実施可能であり、実施形態は単なる例示であって、限定する意図を持たない。
【０２７３】
本明細書に関連して、「備える」という言葉は、「主に含むが、必ずそれのみを含むとは限らない」又は「有する」又は「含む」という意味であり、「それのみから構成される」ことを意味しない。「備える」を変形した言葉は、対応して変形された意味を有する。
【図２Ａ】
2016087473.tif 000047
【図２Ｂ】
2016087473.tif 000048
【図３】
2016087473.tif 000049
【図４】
2016087473.tif 000050
【図５】
2016087473.tif 000051
【図６】
2016087473.tif 000052
【図７】
2016087473.tif 000053
【図８Ａ】
2016087473.tif 000054
【図８Ｂ】
2016087473.tif 000055
【図９】
2016087473.tif 000056
【図１０】
2016087473.tif 000057
【図１１】
2016087473.tif 000058
【図１２】
2016087473.tif 000059
【図１３】
2016087473.tif 000060
【図１４】
2016087473.tif 000061
【図１５】
2016087473.tif 000062
【図１６】
2016087473.tif 000063
【図１７】
2016087473.tif 000064
【図１８Ａ】
2016087473.tif 000065
【図１８Ｂ】
2016087473.tif 000066
【図１】
2016087473.tif 000067