JP5885405B2 - 撮像装置、干渉縞解析プログラム及び干渉縞解析方法 - Google Patents

Description

本発明は撮像装置に関し、特にシアリング干渉計を用いて被検体の情報を得る撮像装置、撮像装置に用いられるプログラム、及び解析方法に関する。

Ｘ線を含む様々な波長の光の干渉を利用して被検体を撮像、計測する技術が知られている。

この技術について簡単に説明をする。

位相波面の揃った（コヒーレントな）光を被検体に対して照射すると、光は被検体の形状や組成によって波面が変化する。波面に変化が生じた光を何らかの方法を用いて干渉させて干渉パターン（干渉縞）を形成させ、この干渉パターンを解析して光の位相波面を回復すると、被検体の位相、散乱、吸収に関する情報を算出することができる。

シアリング干渉計は上述のような光の干渉を利用して光のシア波面を計測する干渉計であり、シアリング干渉計で検出される干渉パターンは被検体によって生じる光の波面変化を微分した情報を有している。

この技術の代表的な適用例としてはレンズなどの形状測定を行う波面測定技術がある。

また、その他の適用例として、Ｘ線を用いて被検体の微分位相像を得る技術もある。

この技術はＸ線を被写体に照射し、被検体の形状や組成によって生じたＸ線の位相差を計測する技術であり、被検体の内部構造に関する情報を有する微分位相像を算出することが可能である。

干渉によって得られた干渉パターンから被検体によって生じた光の波面変化を計算する方法は位相回復法と呼ばれる。

位相回復法には幾つかの種類があるが、その一つにフーリエ変換法と呼ばれる方法があり、なかでも非特許文献１に記載されているような、干渉パターンに対して窓関数を掛けてからフーリエ変換を行う方法を窓フーリエ変換法と呼ぶ。

窓フーリエ変換法は、窓関数を用いないフーリエ変換法と比較すると一般的にノイズに対して耐性があるといった特徴がある。

"Ｗｉｎｄｏｗｅｄ Ｆｏｕｒｉｅｒ ｔｒａｎｓｆｏｒｍ ｍｅｔｈｏｄ ｆｏｒ ｄｅｍｏｄｕｌａｔｉｏｎ ｏｆ ｃａｒｒｉｅｒ ｆｒｉｎｇｅｓ，¨Ｏｐｔ．Ｅｎｇ．４３（７）１４７２−１４７３（Ｊｕｌｙ ２００４）

窓フーリエ変換法に用いる窓関数の大きさ（目安としては半値全幅がよく用いられる）が小さいほど、干渉パターンを局所的に周波数成分で表現することができる。そのため、空間分解能が向上する。

しかし、窓関数の大きさが小さいほど、波数空間において隣接するスペクトル同士の重なりが大きくなり、周波数の分解能が低下するという課題が生じる。

このように窓フーリエ変換は、空間分解能と周波数分解能の一方を高くすると他方が低くなるという課題を有する。

そこで、本発明は隣接するスペクトルの重なりの影響を低減することができる計算方法、プログラムと、撮像装置を提供することを目的とする。

その目的を達成するために、本発明の一側面としての撮像装置は、シアリング干渉計と、前記シアリング干渉計により得られる干渉パターンから被検体の情報を算出する算出手段と、を備えた撮像装置であって、前記算出手段は、前記干渉パターンを窓フーリエ変換することにより得られる波数空間上の互いに異なる座標におけるフーリエ成分を表す式を３以上用い、前記フーリエ成分を表す３以上の式を連立方程式として解くことにより、前記被検体の情報を算出し、前記被検体の情報は、前記被検体の吸収像、散乱像、微分位相像、位相像、の少なくとも１つに関する情報であることを特徴とする。

本発明のその他の側面については、以下で説明する実施の形態で明らかにする。

本発明によると、窓フーリエ変換を用いた位相回復を行う際に隣接するスペクトルの重なりの影響を低減することができる計算方法、プログラムと、撮像装置を提供することができる。

実施形態の撮像装置の概略図。 １次元のトールボット干渉計に用いる格子と干渉パターンの例の模式図。 ２次元のトールボット干渉計に用いる格子と干渉パターンの例の模式図。 実施形態における位相回復に用いた座標を示した波数空間の模式図。 実施例と比較例のシミュレーションに用いた被検体と検出されるモアレの模式図。 実施例により得られた１２８×１２８画素の微分位相像。 比較例により得られた１２８×１２８画素の微分位相像。 窓フーリエ変換法における逐次位相変換のイメージ図。

周波数分解能を保ちつつ空間分解能を向上させる又は空間分解能を保ちつつ周波数分解能を向上させるためには、波数空間において隣接するスペクトルの重なりの影響を考慮して位相回復を行えばよいことが本発明の発明者が鋭意検討した結果わかった。スペクトルの重なりの影響を考慮して位相回復を行う方法として、例えば、スペクトルのフィッティングにより隣接するスペクトルを分離して位相回復を行う方法が考えられる。

しかし、スペクトルのフィッティングを行うには多大なデータを扱うことになる。例えば、１０００×１０００画素の画像を得たとする。それぞれの画素を中心として窓関数を適用し、フーリエ変換を実施すると、元画像の１画素毎に１０００×１０００画素の波数空間のデータが得られる。この結果、データの合計は１０００の４乗となり、スペクトルの分離と位相回復に多大な時間又はコンピュータのリソースを必要とする。

そこで、上述のようなスペクトルのフィッティングによりスペクトルを分離する方法よりも短時間又は低リソースで行うことができる位相回復方法を行う撮像装置について、添付の図面に基づいて以下で説明をする。なお、各図において、同一の部材については同一の参照番号を付し、重複する説明は省略する。

本実施形態ではシアリング干渉計としてトールボット干渉計を使用した撮像装置について説明する。ただし、本特許が適用可能なシアリング干渉計はこれに限られるわけではなく、さまざまな形態のシアリング干渉計に適用可能である。

図１は本実施形態における撮像装置の構成を示した図である。図１に示した撮像装置１は、トールボット干渉計２と、算出手段として計算機６１０と、を備えている。トールボット干渉計２は、光源としてのＸ線源１１０と、Ｘ線を回折する回折格子３１０、Ｘ線の一部を遮る遮蔽格子４１０、Ｘ線を検出する検出器５１０を有している。更に撮像装置１は計算機６１０による計算結果に基づいた画像を表示する画像表示装置７１０と接続されて撮像システムを構成している。

以下、各構成について説明をする。

Ｘ線源としては、連続Ｘ線を出射するＸ線源を用いても、特性Ｘ線を出射するＸ線源を用いても良く、平行Ｘ線（平行光）を出射するＸ線源を用いても、発散Ｘ線（発散光）を出射するＸ線源を用いても良い。但し、本明細書においてＸ線とはエネルギーが２以上１００ｋｅＶ以下の光を指す。

また、Ｘ線源１１０からのＸ線は回折格子３１０で回折されることにより、干渉パターンを形成する必要があるため、Ｘ線源１１０からのＸ線には干渉パターンを形成できる程度の空間的可干渉性が求められる。

Ｘ線源１１０からのＸ線は回折格子３１０により回折され、トールボット距離と呼ばれる所定の距離をおいて明部と暗部が配列した干渉パターンを形成する。但し、本明細書では、Ｘ線（光）の強度が大きい所を明部、小さい所を暗部とする。

本実施形態に用いられる回折格子３１０は位相型の回折格子である。回折格子として振幅型の回折格子を用いることもできるが、位相型の回折格子の方がＸ線量（光量）の損失が少ないので有利である。

図２（ａ）は１次元の干渉パターンを形成する位相格子３１０ａの平面上の構成の一例を表した図である。３１１は位相の基準となる部分、３１２は３１１に対してπ分Ｘ線の位相が変化する部分である。図２（ｂ）は位相格子３１０ａによって形成される干渉パターン８１０ａの明部８１１と暗部８１２を表している。

図３（ａ）は２次元の干渉パターンを形成する位相格子３１０ｂの構成例である。３１１は位相の基準となる部分、３１２は３１１に対してＸ線の位相をπ変化させる部分である。図３（ｂ）は位相格子３１０ｂによって形成される干渉パターン８１０ｂの明部８１１と暗部８１２を表している。

遮蔽格子４１０はＸ線を透過する透過部とＸ線を遮蔽する遮蔽部が配列した構造を有しており、回折格子３１０からトールボット距離はなれて配置される。これにより、遮蔽格子４１０を透過したＸ線はモアレを形成する。尚、遮蔽部は遮蔽格子を透過したＸ線がモアレを形成する程度にＸ線を遮ればよいため、Ｘ線を完全に遮らなくても良い。

光としてＸ線を用いたトールボット干渉計の場合、回折格子により形成される干渉パターンの周期は通常数μｍから大きくとも数十μ程度であるのに対し、一般的に用いられるＸ線の検出器の解像度は数十μｍから数百μｍ程度である。そのためこの干渉パターンを直接検出することは難しい。そこで、本実施形態のように遮蔽格子４１０を用いてモアレを形成し、このモアレを検出する方法が用いられることが多い。このようにモアレを形成する場合、遮蔽格子のピッチは干渉パターンのピッチと同じか、僅かに異なる値をとることができ、形成したいモアレのピッチによって決めることができる。尚、モアレのピッチは遮蔽格子の遮蔽部と透過部が配列した方向と、干渉パターンの明部と暗部が配列した方向のなす角度によっても変化する。モアレの周期は様々な値をとることができるが、一般的に検出器５１０の検出素子の３画素分以上が望ましいとされている。

図２（ｃ）は図２（ｂ）の干渉パターンに用いられる遮蔽格子４１０ａの平面上の構成を表した図である。また、図３（ｃ）は図３（ｂ）の干渉パターンに用いられる遮蔽格子４１０ｂの平面上の構成を表した図である。図２（ｃ）の遮蔽格子４１０ａ、図３（ｃ）の遮蔽格子４１０ｂ共に、透過部４１１と遮蔽部４１２が周期的に配列している。

図２と図３に示した回折格子と遮蔽格子の組み合わせは一例であり、他の組み合わせを用いることもできる。尚、本発明はこれらの格子の構成に依存しない。また、干渉パターンを直接検出する場合は遮蔽格子は不要である。

検出器５１０は、Ｘ線を検出することのできる検出素子（例えばＣＣＤ）を有しており、遮蔽格子を透過して形成されたモアレの強度分布を検出する。本実施形態の撮像装置ではモアレの強度分布を検出するが、干渉パターンの強度分布を直接検出し、解析しても良い。尚、本実施形態の場合は説明のために干渉パターンとモアレを区別して説明したが、モアレも干渉パターンの一種である。つまり、本実施形態ではモアレを検出し、検出したモアレを解析するのでモアレを用いて説明をするが、干渉パターンを直接検出する場合もモアレを検出する場合と同様に解析を行うことができる。

計算機６１０はトールボット干渉計２の検出器５１０の検出結果に基づいて被検体の微分位相像の情報を計算する。

計算機６１０によって行われる計算方法（位相回復方法）について説明するために、まず比較例としてスペクトルのフィッティングによりスペクトルを分離して微分位相像の情報を計算する位相回復方法について説明をする。

２次元の窓フーリエ変換は次の式によって定義される。

ここで、ｆ（ｘ、ｙ）は元関数、ｇ（ｘ、ｙ）は窓関数、（ｘ、ｙ）は座標、（ｕ，ｖ）は窓関数の中心、（ｋ_ｘ，ｋ_ｙ）は波数を表す。ＷＦ［・・・］は括弧内の関数に対して窓フーリエ変換を行うことを示す演算子である。あるモアレの強度分布Ｉ（ｘ，ｙ）に対して窓フーリエ変換を行うと、それぞれの窓関数の中心位置（ｕ，ｖ）毎に波数空間（ｋ_ｘ，ｋ_ｙ）が得られる。

例えば、離散化した１０００×１０００画素の画像（本実施形態ではモアレの強度分布）を窓フーリエ変換すると、画像の夫々の画素を窓関数の中心位置とした窓毎に１０００×１０００画素の波数空間が得られる。つまり、１０００×１０００画素の波数空間が１０００×１０００枚得られる。これは、１０００×１０００画素の画像の情報が、１０００の４乗の情報に変換されることを意味する。

これを図にすると図８のようになる。図８（ａ）はモアレＩ（ｘ、ｙ）の模式図である。窓関数ｇ（ｕ，ｖ）で切り取られる範囲９００は任意の座標（ｕ，ｖ）を中心座標としている。この範囲９００内にフーリエ変換を実行すると、図８（ｂ）に示したような波数空間９０００が得られる。この波数空間は、０次スペクトルの９１１、１次スペクトルの９１２、９１３、９１４、９１５といったスペクトルを有しており、これらのスペクトルから被検体によるＸ線の波面の位相変化やＸ線吸収量、Ｘ線散乱の情報を計算することができる。尚、１次スペクトルはモアレの周期に由来したスペクトルである。

このような波数空間は通常、窓関数の中心座標（ｕ，ｖ）毎に計算される。

つまり、範囲９００から窓関数の中心座標を変化させた範囲９０１内にフーリエ変換を実行すると波数空間９００１が得られる。同様に範囲９０２内にフーリエ変換を実行すると波数空間９００２が、範囲９０３内にフーリエ変換を実行すると波数空間９００３が、範囲９０４内にフーリエ変換を実行すると９００４が得られる。

上述のように窓フーリエ変換において窓関数で切り取られる範囲９００の半径を縮小すると、隣接するスペクトルが重なり合ってしまう可能性がある。

そこで、隣接するスペクトルを分離する。スペクトル９１１〜９１４は窓関数の形状でスペクトルをフィッティングできると考えられるため、本比較例ではこの方法によるフィッティングを用いてスペクトルを分離する。

例えば、窓関数としてガウシアン窓を使用した場合、そのフーリエ変換は同じガウシアン窓となるため、波数空間上のスペクトルもガウシアン窓でフィッティングすればよい。

しかし、ガウシアン窓でフィッティングするにはすべての（ｕ，ｖ）の組み合わせについて窓フーリエ変換を実行し、図８（ｂ）のように波数空間の全座標におけるフーリエ成分を計算する必要がある。そのため、１０００×１０００画素の画像を用いて位相回復を行う場合、上述のようにすべての画素を窓関数の中心として窓フーリエ変換を実行し、得られた波数空間毎に「１０００×１０００画素に相当する波数空間」を用いて位相回復を行う必要がある。そのため位相回復には多大な時間がかかる。特に画像が巨大化すると指数関数的に位相回復に必要な計算時間または計算機のリソースが増す。

そこで、本実施形態では波数空間（ｋ_ｘ，ｋ_ｙ）のマップを作成せず、（式１）における（ｋ_ｘ，ｋ_ｙ）の組み合わせのうち数点のみのフーリエ成分を計算して位相回復を行うことで、計算量を減らす。

本実施形態の計算機６１０によって行われる位相回復方法について説明する。
まず、２次元の位相イメージングを仮定し、モアレが近似的に次のような形で記述できると仮定する。

ここで、ａ（ｘ，ｙ）は被検体による吸収の量を表し、ｂ（ｘ，ｙ）はモアレの振幅を表す。また、Ｐ_１（ｘ，ｙ）およびＰ_２（ｘ，ｙ）は測定する位相を表す。これらは場所によって異なる値をとりうる。ω_１およびω_２はそれぞれｘ方向、ｙ方向のモアレの周期を表す。モアレの形状は（式２）で表される形状に限定されるわけではなく、あくまで一例であり、本実施形態はさまざまなモアレ（干渉パターン）に適用可能である。例えば、モアレの方向が画面のｘ軸方向、ｙ軸方向に沿っていないモアレは（式２）よりも複雑な式で表される。詳しく説明はしないが、このような場合は回転移動変換などを行うことで最終的に（式２）に帰結することができる。

また、（式２）の第３項を０にすると（式２）は１次元のモアレを示す。以下の説明は１次元のモアレにおいても適用できる。

（式２）を（式１）に代入する。ただし、窓関数ｇ（ｘ，ｙ）の幅は十分に小さく、上記のａ（ｘ，ｙ）、ｂ（ｘ，ｙ）、Ｐ_１（ｘ，ｙ）およびＰ_２（ｘ，ｙ）はその範囲内で一定値であると近似できるとする。そのため、以下の議論ではそれぞれａ，ｂ，Ｐ_１，Ｐ_２と略記する。また、窓関数ｇ（ｘ，ｙ）のフーリエ変換をＧ（ｘ，ｙ）と記述する。

（式２）を（式１）に代入すると下記（式３）が得られる。

図４は２次元の位相イメージングの場合に、ある（ｕ，ｖ）を中心座標として窓フーリエ変換を施した時に得られる波数空間（ｋ_ｘ，ｋ_ｙ）のマップを示した図である。上述のように本実施形態ではこのようなマップを作成せず、波数空間における数点のみのフーリエ成分を算出するが、ここでは本実施形態の説明のためにこのようなマップを用いる。

ここで、（０，０）は原点であり０次スペクトルのピークの位置を示し（ω_１，０）、（−ω_１，０）、（０，ω_２）及び（０，−ω_２）は２次元モアレの１次スペクトルのピークを示す。以下図４（ａ）に示すように（０，０）、（ω_１，０）と（−ω_１，０）を用いて位相回復を行う方法を説明する。

この３つの座標におけるフーリエ成分を表す式を用いて位相を回復する。各々の座標におけるフーリエ成分の値は（式３）から以下のよう表すことができる。

尚、（ω_１，０）と（−ω_１，０）は原点に対して点対称であり、複素共役の関係にあるので、（式６）は（式５）から求めることもできる。

（式４）（式５）（式６）を連立方程式として解く。
（式４）にｅｘｐ（−ｉω_１ｕ）Ｇ（ω_１，０）をかけて（式５）との差をとると、

となる。

同様に、（式４）にｅｘｐ（＋ｉω_１ｕ）Ｇ（−ω_１，０）をかけて（式６）との差をとると、

となる。ただし、ここで上式を導出するためには窓関数のフーリエ変換の特徴として以下の特性を考慮する必要がある。
Ｇ（−ω_ａ，ω_ｂ）＝Ｇ（ω_ａ，ω_ｂ）・・・（ｙ成分も同様）
Ｇ（ω_ａ，ω_ｃ）Ｇ（ω_ｂ，ω_ｃ）＝Ｇ（ω_ａ＋ω_ｂ，ω_ｃ）・・・（ｙ成分も同様）
このため（式４）の第４項と第５項は（式５）および（式６）の第４項と第５項と互いに打ち消しあう。

このように、（０，０）、（ω_１，０）、（−ω_１，０）の３つの座標におけるフーリエ成分を表わす式から（式７）と（式８）が導出できる。

（式７）と（式８）に、検出結果から算出したフーリエ成分の値（｛ＷＦ［Ｉ（ｘ，ｙ）］（ｕ，ｖ，０，０）｝、｛ＷＦ［Ｉ（ｘ，ｙ）］（ｕ，ｖ，ω_１，０）｝、｛ＷＦ［Ｉ（ｘ，ｙ）］（ｕ，ｖ，−ω_１，０）｝と窓関数のフーリエ変換の値（Ｇ（ω_１，０）、Ｇ（０，０）、Ｇ（−ω_１，０）、Ｇ（２ω_１，０）、Ｇ（−２ω_１，０））を代入すれば連立方程式の形でｂとＰ_１が計算できる。

ここで、（ω_１，０）と（−ω_１，０）は原点に対して点対称であり、複素共役の関係にあるため、２つの座標におけるフーリエ成分の値は等しい。よって、（０，０）と（ω_１，０）又は（０，０）と（−ω_１，０）の２つの座標のみのフーリエ成分の値を用いることでｂとＰ_１を計算することが可能である。つまり、本実施形態では、波数空間上の座標におけるフーリエ成分を表す式を３つ用い、それらのフーリエ成分を表す式から導かれる連立方程式を、波数空間上の２つの座標におけるフーリエ成分の値を用いて解く。ここで、波数空間上の２つの座標とは、第１の座標（ここでは原点）と、第１の座標と異なり且つ、第１の座標と原点に対して点対称の関係にない第２の座標（ここでは（ω_１，０）又は（−ω_１，０））のことを指す。

図４（ｂ）に示した、（０、０）、（０、ω_２）と（０、−ω_２）を用いて同様に計算を行うことでＰ_２が求められる。また、ａに関しても、求めたｂとＰ_１、Ｐ_２を（式４）、（式５）（式６）のいずれかに代入することで得られる。このように、複素共役の関係を利用することで実質３つの座標におけるフーリエ成分の値を用いて（式３）で仮定したａ、ｂ、Ｐ_１、Ｐ_２の４つの値が求められる。

これら、ａ、ｂ、Ｐ_１、Ｐ_２から、被検体の吸収像、散乱像、微分位相像を得ることができ、さらに微分位相像を積分すると位相像を得ることができる。

以上より、本実施形態を用いると、比較例のようにスペクトルのフィッティングを行うよりも高速又は低リソースで窓フーリエ変換を利用した位相回復法を行うことができる。

以上の例では０次スペクトルのピークと１次スペクトルのピークを用いた例について説明した。しかし、ａ、ｂ、Ｐ_１、Ｐ_２の計算に用いる波数空間上の座標（ｋ_ｘ，ｋ_ｙ）の組み合わせや用いる座標の数はこれに限定されない。本実施形態ではフーリエ成分を表わす式を５つ用い、これらの式から導かれる連立方程式を解くことでａ、ｂ、Ｐ_１、Ｐ_２を計算したが、用いる座標の数が増えるほどより正確にａ、ｂ、Ｐ_１、Ｐ_２を計算することができる。例えば、フーリエ成分を表わす式を用いて連立方程式を複数たててＰ１の値を複数求め、最小二乗法を用いて最終的なＰ１としても良い。但し、検出器の画素を単位とする窓関数の区間をＲとすると、Ｒ^２を超える数の式を用いても、算出できるａ、ｂ、Ｐ_１、Ｐ_２の正確さは向上しない。窓フーリエ変換によって得られる波数空間は元の窓関数の範囲内の画素の情報のみを含むからである。一方、用いる座標の数が増えるほど計算量は多くなる。そこで、計算に用いる式の数は、５以上Ｒ^２以下が好ましい。尚、１次元のモアレの位相回復を行う場合や、１次元の微分位相像を得たい場合など、Ｐ_２を求める必要がない場合はフーリエ成分を表わす式を３つ以上用いればよい。この場合も上記のように複素共役の関係を用いれば、２つの座標におけるフーリエ成分の値からａ，ｂ，Ｐ１の値を計算することができる。本明細書において１次元の微分位相像とは、位相像を１方向に微分して得られる像である。Ｐ_２を求める必要がないときは、Ｒを超える数の座標における式を用いても、算出できるａ、ｂ、Ｐ_１、の正確さは向上しないため、座標におけるフーリエ成分を表わす式を３以上Ｒ以下用いることが好ましい。尚、ガウシアン窓を用いる場合は、ガウシアン窓の分散をσとしたとき、９９％の情報が存在する範囲である±３σに含まれる画素を窓関数の区間とする。

また、モアレによっては０次や１次のスペクトルのみならず、さらに高次のスペクトルも存在する。この高次スペクトルのピークを用いた場合でも同様に連立方程式を立てて計算することが可能である。例えば図８（ｂ）で示されたスペクトル９１６，９１７，９１８，９１９といった２次スペクトルなどを利用する方法も考えられる。

また、用いる座標はスペクトルのピークでなくても良い。但し、フーリエ成分の絶対値が大きい座標を用いる方がノイズの影響を受けにくいので好ましい。

また、用いる座標はＸ軸またはＹ軸上にある方が、ＸまたはＹ軸上にない座標を用いる場合と比べて計算が簡潔にできるため好ましい。

本実施形態では、計算機６１０による計算を簡潔にするために複素共役の関係を用いたが、複素共役の関係を用いなくても位相を回復することはできる。この場合、連立方程式に代入するフーリエ成分の値も３つ以上必要になる。

尚、窓フーリエ変換は畳み込みの定理により、以下のように表すことも可能である。窓関数が原点で対称な奇関数、すなわちｇ（ｘ）＝ｇ（−ｘ）であると仮定すれば、

となる。ここでＦ［・・・］は通常のフーリエ変換、Ｆ^−１［・・・］は逆フーリエ変換を表す。（式９）から、元関数のフーリエ変換Ｆ［ｆ（ｘ，ｙ）］に、波数空間上の窓関数Ｆ［ｇ（ｕ−ｘ，ｖ−ｙ）ｅｘｐ［ｉｋ_ｘ（ｕ−ｘ）＋ｉｋ_ｙ（ｖ−ｙ）］］を乗じ、その逆フーリエ変換を求めることが、元関数に窓関数を乗じてからフーリエ変換を実行することと同じであることがわかる。本実施形態では（式１）を用いて位相回復を行ったが、（式９）を用いて位相回復を行っても良い。

以上、計算機６１０による位相回復方法について説明をした。上述のような計算を計算機６１０で行うには、上述のような計算を実行させるプログラムを計算機６１０に組み込めばよい。

（実施例）
実施形態に記載した撮像装置を用いて、位相回復を行ったシミュレーション結果を実施例として示す。

シミュレーションは、回折格子として図３（ａ）に示した位相格子３１０ｂ、遮蔽格子として図３（ｃ）に示した遮蔽格子４１０ｂ、検出器として１２８×１２８画素の検出器を備えた撮像装置１を用いた。また、被検体は図５（ａ）に示すような球状の被検体１００１を用い、この被検体１００１が検出器の検出範囲の中央に設置されているとしてシミュレーションを行った。

図５（ｂ）は図５（ａ）の被検体と１２８×１２８画素の検出器を用いた時に検出器によって検出されるモアレである。この検出結果を用いて、上述した位相回復方法により得られた微分位相像を図６の（ａ）と（ｂ）に示す。図６（ａ）がＸ方向の微分位相像で、図６（ｂ）がＹ方向の微分位相像である。

同様のシミュレーションを２５６×２５６画素と５１２×５１２画素の検出結果を用いて行い、実際に計算にかかった時間について、表１にまとめた。計算時間はそれぞれ、０．５秒、０．７秒、１．５秒であった。

（比較例）
実施例と同様に、上述した比較例と同じ方法を用いて位相回復を行ったシミュレーション結果を示す。本比較例の撮像装置は計算機で行われる位相回復方法のみが実施例と異なり、そのほかの構成は実施例と同じである。

用いた計算機のリソース上、一度に窓フーリエ空間全体を計算することは不可能であるため、窓フーリエ変換をすべての（ｕ，ｖ）の組み合わせについて計算を行って逐次（ｕ，ｖ）おける微分位相を決定する方法で実行した。実施例同様、図５（ｂ）に示した検出結果を用いて得られた微分位相像を図７（ａ）と（ｂ）に示す。図７（ａ）がＸ方向の微分位相像で、図７（ｂ）がＹ方向の微分位相像である。

図６と図７を比較すると、同じような微分位相像が得られていることがわかる。

また、実施例同様に同様のシミュレーションを２５６×２５６画素、５１２×５１２画素の検出結果を用いて行い、実際に計算にかかった時間について、表１にまとめた。計算時間はそれぞれ、１０５秒、１０４４秒、２０９３８秒と画素数に従って指数的に長い時間を要することがわかる。
また、実際に計算にかかった時間を実施例と比較すると、実施例の方が計算時間が短縮されていることが分かる。

本発明の骨子は、窓フーリエ変換により得られるフーリエ成分の値を表わす式を用いて位相を回復することによって、波数空間の全座標におけるフーリエ成分を計算せずに位相回復を行うことである。これにより、窓フーリエ変換による位相回復法を短時間又は低リソースで実行することができる。したがって、その理念に基づく手法であるのであれば本発明の実施の形態は上記のみに限定されない。

１ シアリング干渉計
２ 撮像装置
６１０ 計算機

Claims (12)

1. シアリング干渉計と、前記シアリング干渉計により得られる干渉パターンから被検体の情報を算出する算出手段と、を備えた撮像装置であって、
   前記算出手段は、
   前記干渉パターンを窓フーリエ変換することにより得られる波数空間上の互いに異なる座標におけるフーリエ成分を表す式を３以上用い、前記フーリエ成分を表す３以上の式を連立方程式として解くことにより、前記被検体の情報を算出し、
   前記被検体の情報は、
   前記被検体の吸収像、散乱像、微分位相像、位相像、の少なくとも１つに関する情報であることを特徴とする撮像装置。
2. 前記算出手段は、
   前記フーリエ成分を表す式を３以上Ｒ以下用い、
   前記被検体の１次元の、吸収像、散乱像、微分位相像、位相像、の少なくとも１つに関する情報を算出することを特徴とする請求項１に記載の撮像装置。
   但し、Ｒとは検出器の画素を単位とする窓関数の区間を示す。
3. 前記算出手段は、
   前記フーリエ成分を表す式を５以上Ｒ^２以下の座標における前記フーリエ成分を表わす式を用い、
   前記算出手段は、
   前記被検体の２次元の、吸収像、散乱像、微分位相像、位相像、の少なくとも１つに関する情報を算出することを特徴とする請求項１に記載の撮像装置。
   但し、Ｒとは検出器の画素を単位とする窓関数の区間を示す。
4. 前記フーリエ成分を表す３以上の式のうち少なくとも１つは、
   前記干渉パターンの周期に由来する前記フーリエ成分を表す式であることを特徴とする請求項１乃至３のいずれか１項に記載の撮像装置。
5. 前記フーリエ成分を表す３以上の式の１つが、
   前記波数空間の原点における前記フーリエ成分を表す式であることを特徴とする請求項１乃至４のいずれか１項に記載の撮像装置。
6. 前記算出手段は、
   前記波数空間の２以上の座標における前記フーリエ成分を算出した値を用いて前記連立方程式を解くことを特徴とする請求項１乃至５に記載の撮像装置。
7. 前記フーリエ成分を表す式のうち少なくとも２つは、
   前記波数空間の原点に対して点対称の関係にある２つの座標における前記フーリエ成分を表す式であることを特徴とする請求項１乃至６のいずれか１項に記載の撮像装置。
8. 前記フーリエ成分を表す式は、
   前記波数空間のＸ軸又はＹ軸上の座標における前記フーリエ成分を表す式であることを特徴とする請求項１乃至７のいずれか１項に記載の撮像装置。
9. 前記フーリエ成分を表す３以上の式は、それぞれ、下記式に示す窓フーリエ変換を定義する式に前記干渉パターンを表す式を代入して得られる式に、前記波数空間上の座標を代入して得られる式であることを特徴とする請求項１乃至８のいずれか１項に記載の撮像装置。
   
10. 前記算出手段は、
    前記座標におけるフーリエ成分の値と、
    前記座標における窓関数のフーリエ変換の値とを前記フーリエ成分を表す式に代入することにより、前記連立方程式を解くことを特徴とする請求項１乃至９のいずれか１項に記載の撮像装置。
11. シアリング干渉計により得られる干渉パターンから被検体の情報を算出するプログラムであって、
    前記干渉パターンを窓フーリエ変換することにより得られる波数空間上の互いに異なる座標におけるフーリエ成分を表す式を３以上用い、前記フーリエ成分を表す３以上の式を連立方程式として解くことにより、前記被検体の情報を算出し、
    前記被検体の情報は、
    前記被検体の吸収像、散乱像、微分位相像、位相像、の少なくとも１つに関する情報であることを特徴とするプログラム。
12. シアリング干渉計により得られる干渉パターンから被検体の情報を得る方法であって、前記干渉パターンを窓フーリエ変換することにより得られる波数空間上の互いに異なる座標におけるフーリエ成分を表す式を３以上用い、前記フーリエ成分を表す３以上の式を連立方程式として解くことにより、前記被検体の情報を算出し、
    前記被検体の情報は、
    前記被検体の吸収像、散乱像、微分位相像、位相像、の少なくとも１つに関する情報であることを特徴とする方法。