JP4436459B2 - Phase recovery phase contrast image - Google Patents

Abstract

A method of obtaining in image of the phase change introduced by an object in penetrating radiation incident on the object includes irradiating the object with penetrating radiation having high lateral spatial coherence, and receiving at least a portion of the radiation at a detector after the radiation has merged from the object and thereby obtaining and storing at least two intensity records for the received radiation each including intensity values at predetermined intervals. These values are utilized to derive a grid of values defining an image of the phase change introduced by the object in the penetrating radiation. The intensity records are obtained at a uniform finite distance after the radiation has emerged from the object, and are for respective different energy distributions of the detected radiation. Apparatus is also disclosed.

Description

が被写体面ｚ＝０でのスカラー複素数強度であるとする。ここに、φ(x,y)＝ｋΨ(x,y)は、被写体を通る波の透過により現れる位相「収差」であり、Ｉ(x,y)＝ａ²(x,y)は対応する強度である。
光源のない半空間ｚ＞０への強度ｕ(x,y,0)＝ｕ₀(x,y)の伝播は、ヘルムホルツの方程式により決定される。
(∇²+k²)u=0. (2)
ｕ(x,y,0)＝ｕ₀(x,y)などのヘルムホルツ方程式の解は、第１レーリー・ゾンマーフェルト積分により与えられる。この積分は、ｚ>>λについて次のように書くことができる。

解析を単純にするために、被写体が光源からの距離に比べて小さく、また、弱く散乱されると仮定する。すなわち、

ここに、Ｒ’＝Ｒ₁Ｒ₂／(Ｒ₁＋Ｒ₂)であり、ｚ＝Ｒ₂(Ｒ₂>>λ)は、「イメージ」面の位置である。
条件(４)は、被写体に入射した波面が弱い曲率を有することを必要とする。また、条件(５)は、単純な物理的意味を有する。特に、被写体面ｚ＝０とイメージ面ｚ＝Ｒ₂との間に焦点がないことを保証する(｜∇₂Ψ｜<<１／Ｒ’であるので)。条件(５)がないと、イメージ面での位相の分布は、焦点面の可能なずれのため、より複雑に被写体面での強度分布に関連する。条件(６)は、入射光束が一様に被写体を照射することと、所定のｘ線波長で被写体内の屈折率の虚数部分の変化が実数部分の変化より十分弱いことを本質的に要求する。この最後の要求は、「位相」被写体にも、また、位相／吸収が混合される被写体にも適用する。
条件(４から６)は回折角が小さいことを保証する。それゆえ、これらの条件の下では、画像はフレネル積分を用いて計算できる：

ここで、Ｒ＝Ｒ₁＋Ｒ₂であり、

である。
積分(７)に、固定位相式を適用することにより、次式が得られる。

ここで、総和は、画像点(ｘ’，ｙ’)に対応する全ての停留点(stationary points)(ｘ_s，ｙ_s)に対して計算され、sgnＳ″_sは点(ｘ_s，ｙ_s)で評価された位相関数Ｓの二次の偏導関数の行列Ｓ″の負の固有値の数であり、detＳ″_sはその行列式｜Ｏ(ｋ^-1)｜≦const/kである。条件(６)は、(９)において残りの項Ｏ(ｋ^-1)が、先行する項よりもずっと小さいことを保証する。
位相関数(８)の偏導関数を計算し、それらをゼロに等しくすることにより、停留点(ｘ_s，ｙ_s)に対して以下の式が得られる：
ｘ’＝Ｍｘ_s＋Ｒ₂Ψ’_x(ｘ_s，y_s)
ｙ’＝Ｍｙ_s＋Ｒ₂Ψ’_y(ｘ_s，y_s) (１０)
ここで、Ｍ＝(Ｒ₁＋Ｒ₂)／Ｒ₁であり、偏導関数：Ψ’_q＝ｄΨ／ｄｑに対してΨ’_qの形式で表記している。式(１０)は、点(ｘ_s，y_s)での物体平面上から発し、画像の与えられた点(ｘ’，y’，Ｒ₂)上で終わる光線の軌道を記述する；定数Ｍは画像の倍率の係数である。条件(５)は式(１０)が唯一の解(ｘ_s，y_s)を有すること、すなわち、画像の各点とそれに対応する物体の特有の点とを接続する唯一の光があることを保証する。結果として、式(１０)において総和記号は省略できる。
位相関数の二次の偏導関数の行列を計算することは簡単である：

この行列は波面の曲率の変化を記述する。条件(５)によりこの行列は非縮退であり、両固有値は正であり(これにより、sgnＳ″＝０)、その行列式は次式により近似され得る。
detＳ″＝(１＋Ｒ’∇²Ψ)／Ｒ’²， Ｒ’｜∇²Ψ｜≪１ (１２)
今、固有位相の式(９)を単純化できる：

画像平面上での対応する強度分布は次式となる。
Ｉ(x’,y’,Ｒ₂)≒Ｍ²Ｉ(ｘ_s,y_s)[１−Ｒ’∇²Ψ(ｘ_s,y_s)] (１４)
今、停留点の式(１０)を用いてさらに式(１４)を単純化すると、
Ｉ(ｘ_s,y_s)＝Ｉ(Ｍ^-1x’−Ｍ^-1Ｒ₂Ψ’ｘ，Ｍ^-1y’−Ｍ^-1Ｒ₂Ψ’y)
≒Ｉ(Ｍ^-1x’，Ｍ^-1y’)−Ｒ’(∇Ｉ・∇Ψ)(Ｍ^-1x’，Ｍ^-1y’) (１５)
(１５)を(１４)に代入し、ｘ＝Ｍ^-1ｘ’、ｙ＝Ｍ^-1ｙ’と表記すると、次式が得られる。
Ｍ²Ｉ(Ｍ_x,Ｍ_y,Ｒ₂)≒
Ｉ(x,y)[１−Ｒ’∇²Ψ(x,y)−Ｒ’∇logＩ(x,y)・∇Ψ(x,y)] (１６)
この式は、弱い位相の物体(条件(４)−(６)により記述されるような)の画像を、光源から物体及び物体から画像への距離の関数として記述する。吸収が無視できる場合、画像コントラストは焦点距離Ｒ’に、及び、物体を通過する際に入射ビームに導かれる位相シフトのラプラシアンに比例する。なお、∇Ψ＝０及び∇²Ψ＝０の場合、式(１６)は光源と画像平面の間に物体がない場合における強度の正確な値を当然に与える。
物体の材料の吸収端から離れたｘ線波長に対して、線形な減衰係数μは近似的にλ³に比例すること及び次式が知られている。

ここでｒ_eは古典的な電子の半径であり、ρは物体の電子密度である。これらの事実及び式(１６)は、受けた近似において画像コントラストが波長に対して単純な二次式的な依存性を有することを意味する。それゆえ、光源の多色性は、必らずしもこの方法の適用に対して障害とはならず、因子λ2が、それぞれの波長成分の二乗に依存するスペクトルで重み付けられた合計により置換されることは当然に考えられる。
式(１６)は、図１においてＤ₁とＤ₂で示される２つの画像平面における強度測定から位相を明示的に取出すために使用される。これがどの様に実現されるか見るために、その式を次の様に書き変えると便利である：
−∇・[Ｉ(x，y)∇Ψ(x，y)]＝Ｆ(x，y) (１８)
ここで、Ｆ(x，y)＝[Ｍ²Ｉ(Ｍ_x，Ｍ_y，Ｒ₂)−Ｉ(x，y)]／Ｒ’である。式(１８)は、物体平面上の位相の横方向の導関数を、物体及び画像平面上の強度分布に関連付ける。
式(１８)は、Ｒ₂→０のときに形式的に得られる平面波の強さの搬送式(ＴＩＥ)と同様である。しかし、両者にはいくつかの重要な相違点がある。先ず第１に、ＴＩＥと違って、式(１８)は、強さが測定される２つの面Ｄ₁とＤ₂間の距離が無限小であると仮定しないことである。第２に、式(１８)は画像拡大を明確に考慮し、従って点光源を用いた位相画像形成にもっと適している。平面波ＴＩＥを直接的に球面入力波に適用しようとすれば、式(１８)の右辺では異なる関数となることは明らかであり、即ち、関数Ｆ(x,y)の分母、即ち、距離Ｒ’はＲ₂に置き換えられる。その結果起こる誤差は、入力球面波によって発生された特別な波面の曲げひずみに起因し、Ｒ₂＜＜Ｒ₁のときのみ無視できる。最後に、我々の式(１８)を導き出す方法は、定常位相形式のフレネル積分への適用に基づいたものであるが、これにより、対物面上の振幅分布ｕ₀(x,y)に関する有効性条件式(４)−(６)を得ることが可能となる。これらは、光源、サンプルおよび画像形成設計の幾何学的パラメータの光学特性において有効な条件である。これに対して、ヘルムホルツ式から平面波ＴＩＥを導き出す従来の方法は、回折された振幅u(x,y,z)，z>0に必要条件を含み、実際の立証が困難である。
式(１８)を解くユニークな方法を見出すために、そのためのいくつかの境界条件を設定する必要がある。このような境界条件は条件式(４)を用いることにより容易に得られる。Ωが座標の原点を中心として半径ｄの円を含む平面(x,y,0)内での単に接続された領域である場合は、Ωの境界Γ上に
Ψ(x,y)/Γ=0 (19)
を得る。
実際には、ディリクレット(Dirichlet)境界条件(１９)だけでなく、Γ上に任意の均一な境界条件、例えばノイマンまたは周期条件など設定可能である。この状態は、平面波ＴＩＥも位相回復のために使用される天文学上の光学に適用される状態とは大変異なっている。そこでは、湾曲波面は、本実施例のようには１次の非撹乱または均等に撹乱された波面によって囲まれていない。従って、境界条件は常に直接に測定しなければならない。
この１次の非撹乱または均等に撹乱された波面が測定されない場合は、画像フィールドの端部で画像を人為的にかすませる工程、または反復的に位相を回復し、対象物の端部の周囲の領域における強さ分布を再演算して、均等な境界条件を効果的に得る工程が使用されうるであろう。
境界条件式(１９)を用いた式(１８)は、I(x,y)≧C²＞０、ただしＣは式(６)で定義された定数である場合、位相のユニークな解法を有し、一方、ノイマンおよび周期的境界条件式を用いた同様の問題解決法は、任意の付加定数までユニークである。強さはゼロではないという必要条件は条件式(６)により満足されている。
問題式(１８)−(１９)とノイマンおよび周期的境界条件式の問題の解法は、右辺の関数F(x,y)における誤差に関して安定している。しかし、M²I(Mx,My,R₂)とI(x,y)における誤差は、好ましくは、これらの強さに対してだけでなく、式(１６)と条件(５)により強さ自体よりもはるかに小さくあらねばならないそれらの差異に対して比較しても、充分小さくなければならないことに留意すべきである。従って、測定された強さのデータは、位相の量的な回復の信頼性を最適にするために、ノイズからは完全にクリーンであることが特に好ましい。
対象面上において均等な強さ、即ち、I(x,y)≡I₀の場合、２乗(-D,D)x(-D,D)，D>dにおける周期的境界条件を用いた式(１８)の単純な解法が、下記の

により得られる。ここで、Ψ_mnおよび^〜Ｆ_mnはそれぞれΨ(x,y)および(経験的に測定可能な)関数^〜Ｆ(x,y)=[M₂I(Mx,My,R₂)-I₀]/(I₀R’)のフーリエ係数である。
対象面上において強さ分布が不均等である一般的な場合は、式(１６)の最後の２項が大きさにおいて比較可能である。しかし、実際の多くの状況では、証明光波は均一な強さを有し、対象物による吸収は穏やかであるが無視できない程度であり、式(１６)の最終項はその他の項よりも小さくできる。このような場合、式(１８)は撹乱法により解くことができる。特に、式(２０)と同様であるが、測定された可変強さ分布I(x,y)により置換された関数^〜Ｆ(x,y)の分子に一定強さI₀を用いた式は、式(１８)の解の１次近似であると見なされることを示すことは容易である。
最後に、式(１６)の最終項が前の項と同じような大きさであるときは、方式(１８)−(１９)は平面でのラプラス式の公知の方法の１つを用いて、数値的に解くことができる。
上記の導出および放射線写真方法および装置は、２次元における１つ以上の検出された強さ分布に含まれる位相コントラスト情報を引き出すための実際的な技法を提供するものである。図１に示す構成は、吸収コントラスト放射線写真で従来から使用されている構成の直接的な進歩した類似である。位相コントラスト画像を引き出す開示の方法は、特に、画素状または格子構造または増加方法で読み出される平面検出器に適し、また、このような検出器とコンピュータとの組合せにも適している。
さらに、従来の方法により吸収コントラストの検出値を処理することにより、２つ以上の２次元平面Ｄ₁とＤ₂放射線写真記録から位相コントラストと吸収コントラスト画像の両方を同時に引き出すことができる。
位相検索を行うために異なった距離Ｄ₁、Ｄ₂での画像強度分布を記録することは常に適当、もしくは実用的ではない。例えば、光軸に沿う異なる位置での二つの平面での強度測定には通常、検出システムＤ(図１において矢印ｄで示す)を正確に動かすことが必要であり、これは実像を対象とした用途では望ましくないことがある。このアプローチでの移送検索についての別の方法では、二つかそれ以上の入射エネルギー、例えば簡単な場合では異なった波長の異なった分布で物体から一定距離隔てたところの像面での二つかそれ以上の強度分布を測定することでなされているか、より望ましくはCd-Mn-Teを利用した二次元エネルギー感応ガタ検出器の如くの検出器を利用して行っている。図２を再度参照するとして、周波数をv=ck/(2π)とし、波数をｋ、真空中での光速をｃとして、

である複振幅の時間及び空間コヒーレンスの任意の状態における自由半空間ｚ≧０での近軸波動場(paraxial wavefield)を考えてみる。今、S(r、v)を」(パワー)スペクトル密度とすると、
S(r_⊥,R,v)-S(r_⊥,v)=-R∇_⊥・[S(r_⊥,v)∇_⊥Ψ(r_⊥,v)], (22)
となり、kψ(r_⊥,v)=argU(r_⊥,v)であり、全ての関数の引数のリストからz=０を省いている。モノクロームの場合では、U(r,v)＝U(r)δ(v-v₀)、S(r,v)=I(r)δ(v-v₀)を得、位相kψ(r,v₀)はkψ(r_⊥)=argU(r_⊥)に一致する。
ＴＩＥは、式(２２)を周波数範囲にわたって積分することにより得られる。

但し、I(r)=∫S(r,v)dvは、(時間平均)強度[１４]を表す。前述したように、この種のＴＩＥは、物体面と像面での強度分布の間の差を物体面における位相微分(phase derivative)に関係付けている。
スペクトル密度の空間変化は、対応する位相変化に比べて非常に弱いので、式(２２)は下記のように簡単にすることができる。

このような状態は位相物体、即ち屈折率の実部が特定波長にあって虚部よりも非常に強い物体には典型的である。これは、例えば高エネルギーＸ線で撮像した薄い生体組織試料にはあてはまることである。そのような物体の場合では、吸収による従来の画像ではコントラストが弱いことがあり、それに対して位相コントラスト像でなら、もっと沢山の情報が得られる[５、８、９]。
試料を透過するＸ線の伝播が近軸であって、直線状であると仮定する。この仮定の下では、物体面での位相とスペクトル密度とは、下記のように表すことができる。

但し、上付文字「in」のある値は試料への入射波に対応し、積分は光軸と平行な線に沿って試料を介するものであり、１-δは屈折率の実部を表し、μは線形吸収係数を表す。スペクトルに見られる全ての波長が、試料素材の吸収縁から充分に遙かに隔てているのであれば、下記式が得られる。

但し、σ=v₀/vであり、v₀はスペクトルから得られるある一定の周波数を表す。式(２５)から(２８)を式(２４)に代入すると、下記式が得られる。

但し、r_⊥は全ての関数の引数のリストから省いてあり、M₀=∫μ(r_⊥,z’,v₀)dz’、ψ₀=∫δ(r_⊥,z’,v₀)dz’を表している。式(２９)における角括弧([])内の最後の二つの項は、例えば近軸近似法が有効であれば１よりもより小さいものでなければならない。
試料に入射する放射線の特性、即ちSⁱⁿ(v)とψⁱⁿ(v)とが既知であり、像面での所定周波数におけるスペクトル密度S(R、v)が(例えば、エネルギー感応型検出器を使って)測定できるのであれば、式(２９)は異なった波長での二つのそのような測定値から得る位相検索に利用することができる。実際、例えばS(R、v)とS(R、v₀)とが測定されているのであれば、例えば式(２９)をvとv₀について書き換え、これらの式からexp{-σ³M₀)を除外すると、下記の式が導き出せる。

位相Ψ₀(ｒ_⊥)はポアソン方程式(３０)を解くことにより再生(又は回復)することができる。もし画像がかき乱されていない１次ビームによってとり込まれているならば、任意の均一の境界条件はΨ₀に割り当てることができることに注意する必要がある。典型的な均一な境界条件を有する式(３０)に対する解は、少なくとも任意の加算定数までは、常に唯一である。位相再生の安定性は、２つの周波数の比σ＝ν₀／ν及び式(３１)におけるスペクトル強度の比によって決定されることが式(３０)及び(３１)から明らかである。実験においては、これらの値は、式(３０)を用いて位相の再生が実験データにおけるエラーに関して安定する場合において、１とは十分に異なるように選択することができる。従って、この技術は、光軸に沿った２つの異なる位置においてビームの交差断面における強度分布から位相再生のための上述の方法に対して利点を有する。実に、この後者の場合においては、交差断面の間のｚ方向の距離は微小である必要があり、かつ２つの強度分布の間の差も微小である必要があり、従って、位相のラプラシアンの計算は、”０を０で除算する”タイプの数値的に不安定な処理である。
本発明の上述の実施例において、小さい多色光源を有するインライン型位相イメージング(画像化)のために最適な式(３０)のアナログ値を得ることができる。入射ビームの準球状(疑似的な球形状)の性質のために波面において表れる付加的な曲線を明確に考慮に入れることができる。次の式(２２)の”準球”のアナログ値について考える。
M²S(Mr_⊥,R,ν)-S(r_⊥,ν)=-R’∇_⊥・[S(r_⊥,ν)∇_⊥Ψ(r_⊥,ν)]}, (32)
ここで、Ｍ＝(Ｒ₁＋Ｒ)／Ｒ₁は大きさの係数である。Ｒ’＝(Ｒ₁ ^-1＋Ｒ^-1)^-1及びＲ₁は光源から物体までの距離である。式(３２)から処理を開始し、上述のように処理を進め、式(３０)に対して同様な位相に対するポアソン方程式を得ることは簡単であり、ここで、Ψⁱⁿ(ｒ_⊥，ν)をΨ’ⁱⁿ(ｒ_⊥，ν)＝Ψⁱⁿ(ｒ_⊥，ν)−ｒ_⊥ ²／(２Ｒ₁)に取り替え、ＲをＲ’に取り替え、γを次式のγ’に取り替える。

エネルギー検出型検出器を用いない若干異なった実験装置においては、周波数ν及びν₀を有する２つの異なった単色の入射ビームを用いて画像平面において強度測定を実行することが可能である。このとき、上記式(３０)及び(３２)は、すべてのスペクトル強度を対応する強度で置きかえることでいまだ満足している。
いま、従来の(エネルギー検出型でない)検出器のみを用いて得られた単色の画像からの位相を検索するための方法について概説することを提案する。ここで、吸収の横方向の変動が次式のように弱いということを仮定する。

ここで、

はν＝ν₀におけるサンプルにおける平均的な吸収である一方、どこでも｜σ³Ｍ₀’(ｒ_⊥)｜＜＜１である。定数

はサンプルに入射する全体の強度の知識及び画像の全体の強度から決定することができる。吸収の弱いローカルな変動Ｍ₀’(ｒ_⊥)は未知であると仮定される。
最初の未加工の近似においては、項Ｍ₀’(ｒ_⊥)を無視するために選択してもよい。そして、位相に対する次式のポアソン方程式は方程式(２９)から表記法

を用いて、容易に得ることができる。

従って、上記位相は式(３５)を解き、入射の放射特性を提供する(すなわちＳⁱⁿ(ν)及びΨⁱⁿ(ν)が先験に基づいて既知である。)ことにより、ただ１つの多色画像Ｉ(Ｒ)≡Ｉ(ｒ_⊥，Ｒ)＝∫Ｓ(ｒ_⊥，Ｒ，ν)ｄνから見つけることができる。
いま、もし吸収の弱いローカルな変動Ｍ₀’(ｒ_⊥)を考慮に入れる必要があるとき、上記位相は、異なるスペクトル要素を有する入射ビームを用いて得られた２つの多色の画像から見つけることができる。式(３４)を式(２９)に代入し、周波数範囲にわたって積分することにより次式を得る。

ここで、ａ＝Ｉ(Ｒ)−∫Ｓ^out(ν)[１−Ｒ∇_⊥ ²Ψⁱⁿ(ν)]ｄν、
ｂ＝−∫Ｓ^out(ν)[１−Ｒ∇_⊥ ²Ψⁱⁿ(ν)]σ³ｄν、
ｃ＝Ｒ∫Ｓ^out(ν)σ²ｄν、及びｄ＝−∫Ｓ^out(ν)σ⁵ｄνである。
(先験で既知の分布Ｓ_j ⁱⁿ(ν)及びΨ_j ⁱⁿ(ν)，ｊ＝１，２によって記述された)異なるスペクトル要素を有する入射ビームを用いて測定された２つの画像Ｉ_j(Ｒ)，ｊ＝１，２が与えられたときに、Ｍ₀’は式(３６)の対応する対から除去することができ、次式を得る。

ここで、上述のように、正確に、入射ビームのそれぞれに対して定義されたすべての量を用いる。例えば、
ａ_j＝Ｉ_j(Ｒ)−∫Ｓ_j ^out(ν)[１−Ｒ∇_⊥ ²Ψ_j ⁱⁿ(ν)]ｄν，ｊ＝１，２ほかであり、同様に、有効的な純粋の吸収−コントラスト画像を得るために、Ｍ₀を解くことができる。
これらの式は、Ｘ線位相が、強度輸送方程式(Transport of Intensity equation)を用いて、光軸に沿った固定位置で得られる多色画像から再生されることができることを示している。それは、擬似平面又は擬似球面の近軸入射波でもって実行されることができ、そのスペクトル密度及び単色成分はプリオリ(priori)として知られている。このアプローチは、好ましくは、位相回復(phase retrieval)実験に先立つ光源(source)の実質的に完全な特性記述を含むであろう。本発明の上記実施の形態は、下記のいずれか１つを含む種々の手法で、有効な位相差(phase−contrast)及び吸収差(absorption−contrast)の画像が得られることを示している。
(１)２つの異なる波長(光子エネルギ)を伴った単色の入射ビームでもって得られる２つの画像から。
(２)エネルギに感応する検出器(例えば、カドミウム・マンガン・テルル化物を基材とする)を用いたＸ線エネルギの関数である画像強度の測定及び多色の入射ビームでもって、少なくとも２つのエネルギのレンジ(ビン：貯蔵場)中において十分な画像強度を与えるべきエネルギ値の貯蔵でもって。
(３)例えば、箔状吸収材の使用に基づくエネルギフィルタの使用による、検出器におけるエネルギ選択性を用いた画像強度の測定及び多色入射ビームでもって。実用上は、これは、これらの間に適当な箔状フィルタを備えた１対のＸ線フィルムを使用して、同時に適切な強度データを記録するといったことを含むであろう。
(４)異なるスペクトル組成を有する入射ビームでもって得られる２つの多色画像から、及びエネルギに感応するのではない従来の検出器から。
上記のこれらの位置のそれぞれにおいて、もし吸収の横方向の空間的ばらつきが無視小であれば、位相回復には画像の測定はただ１つで十分である。あるいは、画像強度データから個別の有効な位相差及び吸収差の画像を得ることが可能である。
図３は、それぞれの異なるエネルギ分布に対する複数の強度記録を得るための装置の変形例を図示している。このような測定は、例えば、エネルギに感応する検出器Ｄ_Eでもって、あるいは調整可能な制御器Ｃにより波長形状を変化させてチュニング(同調)することができる光源でもって、あるいは多色光源及び交換可能な単色光源、又は物体(Ｆ₁)の前か物体の後のどちらかに配置されるフィルタ、例えばエネルギに感応しないＣＣＤ検出器Ｄの前に配置されたシャッタＦ₂でもって、実行されることができる。もし、サンプル中における吸収が無視小であれば、位相はただ１つの多色画像の強度から見出されることができ、得られる入射放射線の特性はプリオリ(priori)として知られている。
２つ又はこれより多い異なるスペクトルの分布の使用に対する他のアプローチ(すなわち、強度の輸送方程式に基づくアプローチに対する代替案)は、データからのスペクトル分布のデコンボルーション(deconvolution)、あるいはこれに代わる∇²Ψを抽出するための最小二乗法によるフィッティング手法を含むことができる。バイエシアン(Bayesian)解析に基づくさらなる精巧な方法もまた用いられることができるであろう。ここで、フーリエ光学のアプローチに基づくテクニックによる位相及び吸収の回復を簡潔に概説することが提案される。
球面波の場合における位相の回復処理は、平面波及び完全直進に対する処理に類似している。かくして、このセクションでは、後者についての説明にとどめられる。
平面波の場合においては、Ｘ線波関数に対するキルヒホッフの式のカウリ(Cowley)形が、ある小さな角の近似でもって、次のように簡略化されることが証明できる。

ここにおいて、我々は、実空間関数はより低い場合において存在し、そして対応する周波数−空間関数(フーリエ変換実施後)はより高い場合、例えばＱ(ｆ_xｆ_y)＝Ｆ[ｑ(ｘ,ｙ)]において存在する、といった約束事を採用する。この誘導に含まれるステップのさらなる論議については、１９７５年のＪ．Ｍカウリ(Cowley)による回折物理(北オランダ：アムステルダム)を参照されたい。関連する理論的な論議はまた、１９９７年のポガニら(Pogany et al)による、Rev. Sci. Inst. ６８、２７７４−２７８２にも見出される。
伝達関数は、φ(ｘ,ｙ)ｔ(ｘ,ｙ)及びμ(ｘ,ｙ)ｔ(ｘ,ｙ)が十分に小さいと仮定すれば、ｑ(ｘ,ｙ)≒１＋ｉφ’(ｘ,ｙ)−μ’(ｘ,ｙ)と記述することができる。ここで、φ’(ｘ,ｙ)＝−φ(ｘ,ｙ)ｔ(ｘ,ｙ)であり、かつμ’(ｘ,ｙ)＝μ(ｘ,ｙ)ｔ(ｘ,ｙ)／2である。伝達関数のフーリエ変換は、この後、Ｑ(ｆ_x，ｆ_y)≒δ(ｆ_x，ｆ_y)＋ｉΦ’(ｆ_x，ｆ_y)−Ｍ’(ｆ_x，ｆ_y)とすることにより与えられる。式(３９)は、この後、次のように記述することができる。

ここにおいて、χ(ｆ_x，ｆ_y)＝πλＲ₂(ｆ_x ²＋ｆ_y ²)である。式(４０)を展開し、両側の逆フーリエ変換を実行することにより、φ’及びμ’中の１次項に対して、関数依存性を除去することができる。

その結果、純粋なお位相の物体の場合は、単一の画像からφｔ−分布を回復することができる。

そして、純粋な吸収の物体の場合は、次のようになる。

位相及び吸収の両効果が重要である物体でのより通常の場合においては、２つの画像Ｉ₁及びＩ₂が用いられることができ、これらは異なる値のＲ₂及び／又はＸ線エネルギでもって収集される(我々はここでは単色の場合について論議しているが、多色の場合は、Ｘ線源に２つの異なる電圧をセットして、２つの異なるスペクトル分布を与えるのが適切であろう)。構成要素に対してどのような吸収エッジもなく、μ₂＝(λ₂／λ₁)³μ₁であり、かつφ₂＝(λ₂／λ₁)φ₁であると仮定すれば、その結果同時に方程式を解くことができ、次の式が得られる。

位相回復の多重エネルギ法を用いる場合は、次の事項を含む点で、多重距離法よりも有利である。
・画像条件の機械的変化ではなく電子的変化によりデータを記録することができ(例えば、エネルギに感応する検出器又はチューブ(Ｘ線管)電圧の急速なスイッチングにより)、非常に急速な位相回復が可能となる。
・かつ、異なるＸ線エネルギ範囲に対して多数の２次元画像を記録するために、入手可能であれば、２次元のエネルギ分散検出器を用いることができる(例えば、ＣｄＭｎＴｅを基材とするもの)。
実施例１
強度の輸送方程式を利用する本発明にかかる方法の有効性をテストするために、シミュレーションデータでもって数値例が準備され、そして図４〜９に、連続的なステップが示されている。
波長λ＝０．１５４ｎｍ(ＣｕＫ_a放射線)であり、ソース・物体間距離及び物体・画像間距離Ｒ₁＝Ｒ₂＝０.２ｍであり、したがってＲ’＝０.１ｍかつＭ＝２であり、完全な単色Ｘ線の点光源であることが仮定された(図２)。また、６４０×６４０μｍ²のサイズであり、図４中に示された物体平面(ｘ,ｙ,０)上に位相分布を生成する、純粋な位相の物体であることも仮定された。位相値の範囲は、[−０.８６６７,０.６７０６]ラジアンであった。
画像平面上の強度分布Ｉ(ｘ,ｙ,Ｒ₂)は、この後、１２８×１２８画素のグリッド(grid)上で、フレネル(Fresnel)積分(７上)を演算することにより計算された。得られた強度分布Ｉ(ｘ,ｙ,Ｒ₂)のスケール化された画像Ｍ²Ｉ(ｘ／Ｍ,ｙ／Ｍ,Ｒ₂)が、図５中に示されている。比較のため、位相ラプラシアン(Laplacian)の計算された分布、−▽²ψ(ｘ,ｙ,())がプロットされ、これが図６中に示されている。明らかに、画像平面上のスケール化された強度分布と、物体平面上の位相ラプラシアンの分布とはほぼ等しくなっており、これは吸収が無視小の場合においては式(１６)が正当であることを裏付けている。
次のステージでは、画像平面上の演算された強度分布と式(２０)とが用いられ、位相が再構築された。この再構築の結果は、図７に示されている。ここでもまた、これは、図４に示す元の位相と非常によく一致し、この再構築の相対平均平方誤差は４．６５％であった。
この方法の実施はまた、均一でない強度の場合、すなわち式(１６)中の最後の項が重要である場合についても、テストされた。物体平面上の強度分布は、次の式に従ってシミュレーションが行われた。
I(x,y)＝−1.647＋3×exp{−(x²＋y²)／(2＊1280²)}
強度値は、[１,１.１２４３]の範囲を有していた。画像平面上の強度分布Ｉ(ｘ,ｙ,Ｒ₂)は、フレネル積分を評価することにより演算された。強度分布の非均一性は、対応する画像に「ビネット(vignetting)」を惹起し(図８)、これは明らかに物体平面上で強度の非均一性を反映している。しかしながら、物体平面上の強度分布を画像平面上のそれから減算し、スケール化された平面間距離Ｒ’＝０.１ｍで距離を除算し、かつ式(２)を適用することにより、我々は、６．６７％のルート平均平方誤差でもって、元の位相を再構築することができた(図９)。この結果は、物体平面上の強度分布の非均一性が普通の程度である場合は(例えば、入射ビームの非均一性及び／又は吸収効果が弱い場合)、式(２０)の上記の単純な変形が、物体によってつくられる位相変化の定量的な回復に有効に用いることができることを確認している。
実施例２
図１０〜１５は、異なる値のＲ₂における２つの画像からの位相及び吸収の回復の具体例を示している。吸収(μｔ)及びフェーズ(φｔ)の分布としては、ＣＳＩＲＯのロゴ(用語：logo)用いられている(それぞれ、図１０、１１参照)。前者については、値は０.０(黒)又は０.１(白)のいずれかであり、後者(逆のロゴ)については−０.１(黒)又は０.０(白)のいずれかである。画像(強度分布)は、１Åの放射線、０.２×０.２ｍｍ²(５１２×５１２画素)の領域、そしてＲ₂＝８０ｃｍ(図１２)及びＲ₂＝１６０ｃｍ(図１３)でもって、平面波の場合について計算された。これらの条件は、画像が即時には解読不可能であるといったものであった。論議された回復アルゴリズム(４４、４５)を適用し、我々は、完全ではないものの、明白に認識できる、図１４、１５に示された結果を得ている。Technical field
The present invention relates to observing structural features of an object using transparent radiation such as X-rays. In particular, an image is extracted from a record of a two-dimensional intensity distribution of radiation penetrating through an object using a phase change that occurs when radiation incident on the object passes through the object. The present invention further extends from single X-ray imaging to acquisition of individual phase and absorption data.
Background art
In the international publications WO95 / 05725 (PCT / AU / 94/00480) and WO96 / 31098 (PCT / AU96 / 00178) of the international patent applications of the applicant's earlier application, a phase difference using X-rays having a high transmission power Various configurations and conditions suitable for contrast imaging are disclosed. Other documents published earlier include the Soviet Union Patent 1402871 and US Patent 5,319,694. In the conventional absorption contrast X-ray projection method, since the absorption contrast is very weak, there is a limit in seeing the internal structure of the object, but the phase contrast imaging method looks at the internal structure of the object. Very reliable on the top. For example, this is the case when looking at soft fibers of the human body.
Substantial problems for optimal and efficient extraction of phase contrast images of objects from actual records at the detector are introduced in the following two related papers. Phys. Rev. Lett. 77, Nugent et al. 2961-2964 (1994); Optical Society (J. Opt. Soc. Am) A13, 1670-82 (1996). In these papers, X-rays having configurations shown in International Publication WO95 / 05725 and US Pat. No. 5,319,694 are X-rays having monotone waves composed of particle beams (monochrome) in a very narrow energy region. As such, the collection of phase information obtained from intensity propagation measurements is based on a modified radiation field propagation process brought about by the phase modification effect by the object. Radiation that penetrates the object due to a change in the thickness of the object, a sudden change in the electron density of the object, or a local change in the refractive index at the boundary, resulting in a change in the local propagation direction of the radiation. The two-dimensional recording of the intensity distribution is performed. According to the above-mentioned paper by Nguyent et al., Based on Maxwell's equation, the transport-of-intensity equation is derived from the propagation characteristics of monochromatic electromagnetic waves consisting of particle beams (monochromatic) in a very narrow energy region. The phase contrast image is extracted from the intensity record with the solution of this equation. If the solution of the intensity transfer equation is used, the phase is expanded to an orthogonal coordinate function. The choice of function depends on the shape of the sample. Therefore, Zernike polynomials are suitable for circular ones, and Fourier expansion equations are suitable for rectangular ones.
In the above-mentioned International Publication WO96 / 31098, an in-line phase contrast imaging method using a substantial point source and a two-dimensional X-ray image detector at a position away from an object is disclosed. Yes. Compared to the phase contrast imaging method described above, a point source is used, the radiation from the source is multicolor, has a high degree of coherence in the lateral direction of space, and the maximum source diameter (s) is , Depending on the distance (R1) between the source and the object. The greater the spacing between source objects, or the smaller the size of the source, the higher the degree of coherence in the lateral direction of the space (see Nature 384 335-8 Wilkins et al (1996)). The technique proposed in WO96 / 31098 is related to a method used for conventional absorption contrast X-ray imaging, and proposes a technique that can be used more easily than previously proposed techniques. This phase contrast imaging method is effective when X-rays having a large transmittance and a small number of lens members are used. If the number of lens members is small, the technique used for a microscope using visible light or X-rays with low transmission ability cannot be used.
An object of the present invention is to propose a method for obtaining a phase contrast image by recording penetrating radiation from a point source and recording a two-dimensional intensity distribution.
Another object of the present invention is to propose a method suitable for extracting phase and absorption contrast information from an X-ray captured image recorded using a microfocus source that does not require the use of pure monochromatic lines. It is to be.
Disclosure of the invention
(Technical problem to be solved by the invention)
According to one aspect of the present invention, radiation is emitted from an object and recorded at two or more intensity distributions at positions separated by a certain distance to obtain an energy distribution of the radiation, respectively.
According to another aspect of the present invention, using a point source shows a value close to the solution of the intensity transfer differential equation (different from that used in the case of the conventional monotone wave).
The present invention is suitable for separating and retrieving (recovering) phase information from a general intensity record (obtained using a fine focus radiation source) with both phase contrast and absorption contrast. .
According to a first aspect, there is provided a method of obtaining an image formed by using a phase change caused by an object by radiation penetrating the object and penetrating the object,
Has a high coherence in the lateral direction of the space, emits radiation that penetrates the object,
The detection means receives a part of the radiation emitted from the object, records at least two intensities, and each intensity has an intensity value at a predetermined interval,
Using these values, we can extract an array of values and form an image with the phase change caused by radiation penetrating the object,
It is a method characterized in that radiation emitted from an object is captured at a certain distance from the object, the intensity is recorded, and different energy distributions of the detected radiation are generated.
In the first aspect described above, an apparatus for obtaining an image formed by using a phase change caused by an object by radiation penetrating the object and penetrating the object,
Means for providing a source having a high coherence in the transverse direction of the space and emitting radiation penetrating the object;
Detecting means for receiving a part of radiation emitted from the object, and recording at least two intensities by the radiation received by the detecting means, each intensity record having an intensity value in a predetermined interval,
The detection means captures radiation emitted from the object at a certain distance from the object, records the intensity, and further provides energy measurement means, which detects different energy distributions of the detected radiation with the two intensities. It is an apparatus characterized in that it is generated corresponding to recording.
In one embodiment, different energy distributions can be obtained by changing the energy spectrum of radiation irradiating the object. This could be achieved, for example, by changing the output of the radiation source somewhat or by prefilter means. In another embodiment, different energy distributions are obtained by providing detection means for providing intensity as a function of energy in a certain energy band or group of bands. For this purpose, the two-dimensional detection means may be variably wavelength sensitive, or a variable filter opening / closing means may be placed in front of it. Further alternatively, image intensity may be recorded as a function of photon energy for each pixel over multiple ranges of x-ray energy. In order to obtain a high resolution, a composite intensity record may be conveniently obtained for each composite energy distribution of radiation.
In the simplest case, each energy distribution may be at a specific wavelength or photon energy level. The aforementioned derivation uses one or more predetermined boundary conditions to relate one or more phases in the plane of the object to the evolution of the intensity distribution along the propagation direction. It may include solving differential transport-of-intensity equations. An alternative derivation involves solving a Fourier optical equation. Others are of course possible for a particular form or situation.
Preferably, the intensity values also reflect the absorption contrast in the object, and the method further uses these values to derive a grid value that defines the effective net absorption-contrast image of the object. Contains.
The apparatus according to the first aspect of the present invention defines an image of a phase change induced by an object in penetrating radiation when installed on a computer having a suitable operating system and storage means. Further included is a computer program product having a set of machine readable instructions that configures the computer to be executable to utilize the values to derive grid values to perform.
In a second aspect, the present invention describes the phase induced by an object at the incidence of transmitted radiation on the object from one or more two-dimensional intensity records of the transmitted radiation after it has passed the object. A method is provided for obtaining an image of the change in which the radiation has a high lateral spatial coherence upon incidence on the object, and the recording or each recording has a uniform phase shift. It can be obtained at a finite distance in the surrounding field after radiation is generated from an object incorporating a phase-varying element, whether it is or not phase-shifted. The method is
Storing intensity values from that record or each record at predetermined intervals;
By solving the differential equation of intensity transfer that relates the phase at the exit surface of the object to the gradual change of the intensity distribution along the propagation direction, the grid value that defines the image of the phase change induced by the object in transmitted radiation is obtained. To derive the above values and whatever predetermined boundary conditions are to be derived.
The present invention also provides, in its second aspect, a method for obtaining an image of a phase change induced by an object upon incidence of transmitted radiation on the object, the method comprising:
Illuminating the object with transmitted radiation having a high lateral spatial coherence;
One or more finite distance detectors in the surrounding field, whether uniformly phased or not phased, after radiation has been generated from an object incorporating a phase-varying element Receiving at least a portion of the radiation at and thereby obtaining and storing intensity values for the radiation received at predetermined intervals;
By solving the differential equation of intensity transfer that relates the phase at the exit surface of the object to the gradual change of the intensity distribution along the propagation direction, the above values and any predetermined constant boundary condition can be transmitted. Using to derive a grid value defining an image of a phase change induced by an object at the incidence of radiation.
Preferably, the same intensity value is further used to derive a value that defines an effective net absorption-contrast image for the object.
The invention is further induced, in its second aspect, by an object at the incidence of transmitted radiation on the object from one or more two-dimensional intensity records of the transmitted radiation after it has passed through the object. A device for obtaining an image of the phase change, wherein the radiation has a high spatial coherence in the lateral direction when incident on the object, and the recording or each recording has a uniform phase variation. Whether or not phase variation is obtained at a finite distance in the surrounding field after radiation is generated from an object incorporating a phase variation element. The device is
(a) means for storing intensity values from the records or from each record at predetermined intervals;
(b) defining an image of the phase change induced by the object in transmitted radiation by solving a differential equation of intensity transfer relating the phase at the exit surface of the object to a gradual change in intensity distribution along the propagation direction Means for using the above values and any predetermined predetermined boundary conditions to derive the grid values, preferably with means having a program storing machine-readable instructions.
The present invention also provides a set of computer means for execution when installed on computer means having suitable operating system software and storage means, such as the apparatus according to the preceding paragraph. Result in mechanically readable instructions. The present invention still further provides storage means, such as a magnetic disk, CD-ROM or optical storage disk or Internet service, in which the mechanically readable instructions are stored.
The invention also provides, in its second aspect, an apparatus for obtaining an imaging of a phase change induced by an object in transmitted radiation incident on the object, the apparatus comprising:
A point source that irradiates an object with transmitted radiation having a high degree of lateral spatial coherence;
To mix the phase-shifted component with the surrounding field of radiation that is not phase-shifted or that is uniformly phase-shifted, thereby generating intensity values for the radiation acquired at a given interval A detector that receives at least a portion of the radiation at a finite distance after the radiation is emitted from the object;
Includes a machine-readable instruction storage program that solves the transport-of-intensity differential equation related to the evolution of the intensity distribution along the direction of propagation of the phase at the exit surface of the object Thus, in order to derive a grid-like value that defines the imaging of the phase change caused by the object in the transmitted radiation, there is a calculation means capable of using the above values and any predetermined invariant boundary conditions. is doing.
Preferably, the calculation means can also use the same intensity value to derive a value that defines the effective net absorption-contrast imaging for the object.
In an embodiment of the method according to the second aspect of the invention, two of the intensity records are acquired at different distances after the radiation has been emitted from the object, for example at two separate imaging planes. In yet another embodiment, two or more intensity records are obtained at a uniform distance, for example on a single imaging plane. In particular, in the latter embodiment mode, one or more of the intensity records are recorded such that the intensity of the imaging as a function of photon energy is recorded for each pixel over multiple ranges of x-ray energy. The
The intensity transport equation related to the point source case may be the following equation (16), or an equation (18) that is a modification thereof, and its solution is by a variation method. Instead, the equation can be solved with absolute values, especially if the two terms behind equation (16) have the same magnitude.
As used herein, the term “transmission radiation” includes X-rays and neutrons, although the present invention is particularly useful for X-ray radiation, and is more typically polychromatic. However, it may be substantially monochromatic. A particularly effective application is in the range of 0.5 keV to 1 MeV, for example, the range of reliable X-rays is 1 keV to 1 MeV. The phase variation due to the object may be considered as a refraction effect or, more strictly, as a Fresnel diffraction effect. For example, for a finite thickness object in the surrounding medium with different refractive indices, the phase variation will also depend on the thickness of the object in the direction of the locally limited wave vector.
The object may be, for example, a boundary, typically a boundary that exhibits a sharp change in refractive index with respect to transmitted radiation. The present invention is particularly useful for radiation where there is a weak or negligible absorption contrast between the intensities of the radiation passing through each side of the boundary.
The boundary conditions typically do not need to be measured and may have, for example, uniform Dirichlet, Neumann or periodic boundary conditions. The boundary conditions are selected to achieve a unique solution of the phase equation up to at least any constant component.
Preferably, the above solution also utilizes one or more optical conditions. Such conditions include, for example, a slight wavefront curvature of incident radiation, a lack of focus between the object and the image, and uniform illumination of the object.
Incident transmitted radiation is not limited to monochromatic ones. For polychromatic radiation, the above equation has a factor that depends on the spectrum weighted term or the square of each wavelength component.
Brief description of the drawings
The invention will be further described, by way of example only, with reference to the accompanying drawings.
FIG. 1 is a diagram showing a configuration of an X-ray optical system according to an embodiment of the present invention in which two intensity recordings are performed on different detection surfaces.
FIG. 2 is a diagram showing a related coordinate system for mathematical examination described later.
FIG. 3 is a diagram showing a configuration of another X-ray optical system in which two intensity recordings are performed on a common detection surface with respect to different energy distributions of the respective radiations.
4-9 show successive intensity records in a mathematical test case of the method using one approach of the present invention.
10-15 are photomicrographs showing another approach test.
Preferred embodiment
The configuration shown in FIG. 1 consists of a high lateral spatial coherence microfocus source S and eg a film, a photosensitive phosphor plate (eg a Fuji imaging plate). Or an X-ray two-dimensional imaging detector D such as a two-dimensional electron detector such as a charge coupled device (CCD) array.
Here, the term “lateral spatial coherence” refers to the complex wave amplitude correlation between different points transverse to the wave propagation direction. Lateral spatial coherence is assumed to occur when each point on the wavefront has a propagation direction that does not change over time. In practice, a high degree of lateral spatial coherence can be achieved, for example, by using an efficient small sized point source or by observing the beam at a wide distance from the point source. Generally, the lateral coherence length d_⊥= ΛR₁/ S. Here, λ is the wavelength of the X-ray, R1 is the distance from the point source to the object, and s is the maximum point source diameter. For example, for a 20 keV x-ray and a distance from the 200 mm point source to the object, a point source size of about 20 μm in diameter or less is typically appropriate. For the purposes of the present invention, the smaller the point source size, the better, assuming that the total flux from the point source is sufficient. The lateral spatial coherence needs to be preserved by careful selection of the point source x-ray window, for example, so that it has a very uniform thickness and homogeneity. The effect of subspace and temporal coherence on image contrast and resolution is described by Pogany et al. [Revised. Science. Instrum. (Rev.Sci.Instrum.)68  2774-82 (1997)].
Regions of reflectance change transverse to the local direction of propagation, or thickness or density change in the direction of propagation, can lead to significant changes in the local direction of propagation of wavefronts through those regions.
Referring to FIG. 1 again, the spherical wavefront W1 radiated from the point source S is distorted to W2 when passing through the subject O. By recording the intensity of the wavefront sufficiently away from the sample, a change in intensity due to a rapid change in refractive index and thickness or density in the sample can be detected and its position is recorded in the image. This corresponds to one form of differential phase and contrast imaging. The position of the image detector is chosen so that the spatial resolution of the detector is sufficient to resolve the intensity differences resulting from severe distortion of the wavefront, and practical considerations make the contrast most effective. The The numerical value recorded by the detector reflects the absorption in the subject, and the intensity value recorded by the detector also includes extractable absorption / contrast information.
For reasons explained later, it is noted here that the radiation emerging from the subject is out of phase in the surrounding radiation field which is not out of phase or is uniformly out of phase. The subject is arranged to accept the component.
Typically, a steep gradient in refractive index or thickness is imaged as a sudden loss or change in intensity at a corresponding point in the image. This feature of abrupt loss or change in intensity at a given point in the image is essentially independent of wavelength in the first approximation, and therefore even in the case where a polychromatic light source is used. It leads to a very sudden contrast change.
A feature of this arrangement is that for a circular light source distribution, the spatial resolution in the image is the same in all directions and is essentially determined by the size of the light source. Another feature is that the image can be significantly magnified, which allows the use of recording media such as photosensitive phosphor imaging plates with many desirable properties such as wide dynamic range, high sensitivity, and low spatial resolution. It is.
The detector D is typically connected to a suitable computing means M such as a personal computer using a 486 CPU operating at 66 MHz and using a suitable memory and application software.
The calculation means M stores a set of intensity values in the memory at intervals predetermined in two dimensions from the intensity record in the detector D. In this simple case, the detector is a flat planar detector and the spacing is on a uniform square grid. The method for deriving these values depends on the type of detector. For example, the detector comprises a pixel structure, each pixel is measured in a two-dimensional scan, and the numerical value is continuously sent to a computer memory. In another method, the recording is made and stored at the detector and scanned or data is obtained as required by the computer.
In addition, in the control program included in the computer M, the transport-of-intensity differential equation is calculated using the stored intensity values from the detector records and the predetermined uniform boundary conditions. Solving derives another grid of values that define an effective phase and contrast image on the selected surface of the subject.
Because of the simplicity of formal development, we discuss the theoretical point of view in the case of a monochromatic point source, where we define a simple differential equation and weak phase as described by conditions (4)-(6) A solution that can be applied to the subject will be described. The case of a wide pleochroic light source will be considered later.
The point z on the optical axis Z = −R will be described with reference to the coordinate system shown in FIG.₁(R₁>> Consider a monochromatic point light source with wavelength λ = 2π / k located at λ). The light source illuminates a subject located between the light source and a plane (x, y, 0) perpendicular to the optical axis at the point z = 0. The subject has a boundary on the right side with plane z = 0 (FIG. 2). The subject is thin, that is, its z-direction dimension is R₁It is assumed that it is much smaller.

Is the scalar complex intensity at the subject plane z = 0. Here, φ (x, y) = kΨ (x, y) is a phase “aberration” that appears due to transmission of waves through the subject, and I (x, y) = a²(x, y) is the corresponding intensity.
Intensity u (x, y, 0) = u to half space z> 0 without light source₀The propagation of (x, y) is determined by the Helmholtz equation.
(∇²+ k²) u = 0. (2)
u (x, y, 0) = u₀The solution of the Helmholtz equation such as (x, y) is given by the first Rayleigh-Sommerfeld integral. This integral can be written for z >> λ as

To simplify the analysis, it is assumed that the subject is small compared to the distance from the light source and is scattered weakly. That is,

Where R '= R₁R₂/ (R₁+ R₂) And z = R₂(R₂>> λ) is the position of the “image” plane.
Condition (4) requires that the wavefront incident on the subject has a weak curvature. Condition (5) has a simple physical meaning. In particular, the subject plane z = 0 and the image plane z = R.₂Guarantee that there is no focus between₂Ψ | << 1 / R ′). Without condition (5), the phase distribution on the image plane is more complexly related to the intensity distribution on the subject plane due to possible focal plane shifts. Condition (6) essentially requires that the incident light beam illuminate the subject uniformly and that the change in the imaginary part of the refractive index in the subject is sufficiently weaker than the change in the real part at a predetermined x-ray wavelength. . This last requirement applies to "phase" subjects as well as subjects with mixed phase / absorption.
Conditions (4 to 6) ensure that the diffraction angle is small. Therefore, under these conditions, the image can be calculated using Fresnel integration:

Where R = R₁+ R₂And

It is.
By applying a fixed phase equation to the integral (7), the following equation is obtained.

Here, the sum is the sum of all stationary points (x's) corresponding to the image points (x ', y')._s, Y_s) And sgnS ″_sIs a point (x_s, Y_s) Is the number of negative eigenvalues of the matrix S ″ of second order partial derivatives of the phase function S evaluated in_sIs the determinant | O (k^-1) | ≦ const / k. Condition (6) defines the remaining term O (k^-1) Is much smaller than the preceding term.
By calculating the partial derivatives of the phase function (8) and making them equal to zero, the stationary point (x_s, Y_s) Gives the following formula:
x ′ = Mx_s+ R₂Ψ ’_x(x_s, Y_s)
y ’= My_s+ R₂Ψ ’_y(x_s, Y_s(10)
Where M = (R₁+ R₂) / R₁And the partial derivative: Ψ ′_q= Ψ 'for dΨ / dq_qIt is written in the form of Equation (10) is the point (x_s, Y_s) From the object plane at the given point (x ', y', R) of the image₂) Describes the ray trajectory ending above; the constant M is a factor of the magnification of the image. Condition (5) is the only solution (x_s, Y_s), I.e., there is only one light connecting each point of the image with the corresponding unique point of the object. As a result, the summation symbol can be omitted in equation (10).
It is easy to compute a matrix of second partial derivatives of the phase function:

This matrix describes the change in curvature of the wavefront. By condition (5) this matrix is non-degenerate, both eigenvalues are positive (thus sgnS ″ = 0), and the determinant can be approximated by:
detS "= (1 + R'∇²Ψ) / R ′², R ’| ∇²Ψ | << 1 (12)
Now we can simplify the eigenphase equation (9):

The corresponding intensity distribution on the image plane is:
I (x ’, y’, R₂) ≒ M²I (x_s, y_s) [1-R'∇²Ψ (x_s, y_s)] (14)
Now, further simplifying equation (14) using the equation (10) for the stopping point,
I (x_s, y_s) = I (M^-1x'-M^-1R₂Ψ′x, M^-1y'-M^-1R₂Ψ’y)
≒ I (M^-1x ’, M^-1y ′) − R ′ (∇I · ∇Ψ) (M^-1x ’, M^-1y ’) (15)
Substituting (15) into (14), x = M^-1x ', y = M^-1When expressed as y ', the following equation is obtained.
M²I (M_x, M_y, R₂) ≒
I (x, y) [1-R′∇²Ψ (x, y) −R′∇log I (x, y) · ∇Ψ (x, y)] (16)
This equation describes an image of a weak phase object (as described by conditions (4)-(6)) as a function of the distance from the light source to the object and from the object to the image. If the absorption is negligible, the image contrast is proportional to the focal length R 'and to the Laplacian of the phase shift that is directed to the incident beam as it passes through the object. Note that ∇Ψ = 0 and ∇²When ψ = 0, equation (16) naturally gives an accurate value of the intensity when there is no object between the light source and the image plane.
For x-ray wavelengths away from the absorption edge of the material of the object, the linear attenuation coefficient μ is approximately λ^ThreeAnd the following equation is known.

Where r_eIs the classical electron radius and ρ is the electron density of the object. These facts and equation (16) mean that the image contrast has a simple quadratic dependence on wavelength in the approximation received. Therefore, the pleochroism of the light source is not necessarily a hindrance to the application of this method, and the factor λ2 is replaced by a summation weighted with a spectrum that depends on the square of each wavelength component. Naturally it is possible to think.
Equation (16) is expressed as D in FIG.₁And D₂Is used to explicitly extract the phase from the intensity measurements in the two image planes shown in FIG. To see how this is achieved, it is convenient to rewrite the expression as follows:
−∇ [[I (x, y) ∇Ψ (x, y)] = F (x, y) (18)
Where F (x, y) = [M²I (M_x, M_y, R₂) -I (x, y)] / R '. Equation (18) relates the lateral derivative of the phase on the object plane to the intensity distribution on the object and the image plane.
Formula (18) is R₂→ This is the same as the plane wave intensity transfer type (TIE) obtained formally when 0. However, there are some important differences between the two. First of all, unlike TIE, equation (18) gives the two planes D whose strength is measured D₁And D₂Do not assume that the distance between them is infinitesimal. Secondly, equation (18) explicitly considers image magnification and is therefore more suitable for phase imaging using point sources. If the plane wave TIE is applied directly to the spherical input wave, it is clear that the right side of the equation (18) has a different function, that is, the denominator of the function F (x, y), that is, the distance R ′. Is R₂Is replaced by The resulting error is due to the special wavefront bending strain generated by the input spherical wave and R₂<< R₁Can be ignored only when Finally, the method of deriving our equation (18) is based on the application to the stationary phase form of Fresnel integration, but this results in an amplitude distribution u on the object plane.₀It becomes possible to obtain the validity conditional expressions (4) to (6) regarding (x, y). These are effective conditions in the optical properties of the geometric parameters of the light source, sample and imaging design. On the other hand, the conventional method for deriving the plane wave TIE from the Helmholtz equation includes a necessary condition for the diffracted amplitude u (x, y, z), z> 0, and is difficult to prove in practice.
In order to find a unique way to solve equation (18), it is necessary to set some boundary conditions for it. Such a boundary condition can be easily obtained by using conditional expression (4). If Ω is simply a connected region in a plane (x, y, 0) that includes a circle with a radius d centered on the origin of the coordinates, on Ω boundary Γ
Ψ (x, y) / Γ = 0 (19)
Get.
In practice, not only the Dirichlet boundary condition (19) but also any uniform boundary condition such as Neumann or periodic condition can be set on Γ. This state is very different from the state applied to astronomical optics where plane wave TIE is also used for phase recovery. There, the curved wavefront is not surrounded by a primary undisturbed or evenly disturbed wavefront as in this embodiment. Therefore, boundary conditions must always be measured directly.
If this first order undisturbed or evenly disturbed wavefront is not measured, the process of artificially blurring the image at the end of the image field, or iteratively recovering the phase around the end of the object A process of recalculating the intensity distribution in the region to effectively obtain uniform boundary conditions could be used.
Expression (18) using the boundary condition expression (19) is expressed as I (x, y) ≧ C²> 0, where C is a constant defined in equation (6), has a unique solution to the phase, while a similar solution using Neumann and the periodic boundary condition is Unique up to a constant. The requirement that the strength is not zero is satisfied by conditional expression (6).
The problem equations (18)-(19) and the solution of the Neumann and periodic boundary condition problem are stable with respect to the error in the right-hand side function F (x, y). But M²I (Mx, My, R₂) And I (x, y) are preferably not only for these strengths, but also for their differences that must be much smaller than the strengths themselves by equation (16) and condition (5). It should be noted that it must be sufficiently small when compared to. Therefore, it is particularly preferred that the measured intensity data is completely clean from noise in order to optimize the reliability of the quantitative recovery of the phase.
Equal strength on the target surface, i.e. I (x, y) ≡I₀, The simple solution of equation (18) using the periodic boundary condition for square (-D, D) x (-D, D), D> d is

Is obtained. Where Ψ_mnand^~F_mnAre Ψ (x, y) and (empirically measurable) functions, respectively.^~F (x, y) = [M₂I (Mx, My, R₂) -I₀] / (I₀R ′) Fourier coefficient.
In the general case where the intensity distribution is uneven on the target surface, the last two terms of equation (16) can be compared in size. However, in many practical situations, the proof light wave has a uniform intensity, the absorption by the object is moderate but not negligible, and the final term in equation (16) can be smaller than the other terms. . In such a case, equation (18) can be solved by the disturbance method. In particular, a function similar to equation (20) but replaced by the measured variable intensity distribution I (x, y)^~F (x, y) molecule has constant strength I₀It is easy to show that the equation using is considered to be a first-order approximation of the solution of equation (18).
Finally, when the final term of equation (16) is as large as the previous term, schemes (18)-(19) use one of the known methods of Laplace's equation in the plane: It can be solved numerically.
The above derivation and radiographic methods and apparatus provide practical techniques for deriving phase contrast information contained in one or more detected intensity distributions in two dimensions. The configuration shown in FIG. 1 is a direct and advanced analog of the configuration conventionally used in absorption contrast radiography. The disclosed method for deriving phase contrast images is particularly suitable for planar detectors that are read out in a pixel-like or grating structure or augmentation method, and also suitable for the combination of such detectors and computers.
Further, two or more two-dimensional planes D can be obtained by processing the absorption contrast detection value by a conventional method.₁And D₂Both phase contrast and absorption contrast images can be extracted simultaneously from radiographic records.
Different distances D to perform phase search₁, D₂It is not always appropriate or practical to record the image intensity distribution at. For example, intensity measurements in two planes at different positions along the optical axis usually require precise movement of the detection system D (indicated by the arrow d in FIG. 1), which is intended for real images. It may not be desirable for the application. Another approach to transport retrieval with this approach is to use two or more incident energies, for example two or more in the image plane at a certain distance from the object with different distributions of different wavelengths in the simple case. This is done by measuring the intensity distribution of Cd-Mn-Te, or more preferably by using a detector such as a two-dimensional energy-sensitive backlash detector using Cd-Mn-Te. Referring again to FIG. 2, the frequency is v = ck / (2π), the wave number is k, the speed of light in vacuum is c,

Consider a paraxial wavefield in a free half space z ≧ 0 in any state of double amplitude time and spatial coherence. Now, let S (r, v) be the (power) spectral density.
S (r_⊥, R, v) -S (r_⊥, v) =-R∇_⊥・ [S (r_⊥, v) ∇_⊥Ψ (r_⊥, v)], (22)
And kψ (r_⊥, v) = argU (r_⊥, v), omitting z = 0 from the list of arguments of all functions. In the case of monochrome, U (r, v) = U (r) δ (v-v₀), S (r, v) = I (r) δ (v-v₀) And phase kψ (r, v₀) Is kψ (r_⊥) = argU (r_⊥).
TIE is obtained by integrating equation (22) over the frequency range.

However, I (r) = ∫S (r, v) dv represents (time average) intensity [14]. As described above, this type of TIE relates the difference between the intensity distributions at the object plane and the image plane to the phase derivative at the object plane.
Since the spatial change in spectral density is very weak compared to the corresponding phase change, equation (22) can be simplified as follows.

Such a state is typical for a phase object, ie an object whose real part of the refractive index is at a particular wavelength and much stronger than the imaginary part. This is true, for example, for thin biological tissue samples imaged with high energy X-rays. In the case of such an object, contrast may be weak in a conventional image due to absorption, whereas a phase contrast image provides more information [5, 8, 9].
Assume that the X-ray propagation through the sample is paraxial and linear. Under this assumption, the phase and spectral density at the object plane can be expressed as follows:

However, the value with the superscript “in” corresponds to the incident wave to the sample, the integration is through the sample along a line parallel to the optical axis, and 1−δ represents the real part of the refractive index. , Μ represents a linear absorption coefficient. If all wavelengths found in the spectrum are sufficiently far away from the absorption edge of the sample material, the following equation is obtained.

However, σ = v₀/ v and v₀Represents a certain frequency obtained from the spectrum. Substituting equations (25) through (28) into equation (24) yields the following equation:

Where r_⊥Is omitted from the list of all function arguments, and M₀= ∫μ (r_⊥, z ’, v₀) dz ’, ψ₀= ∫δ (r_⊥, z ’, v₀) dz ′. The last two terms in square brackets ([]) in equation (29) must be less than 1, for example if the paraxial approximation method is valid.
Characteristics of the radiation incident on the sample, ie Sⁱⁿ(v) and ψⁱⁿIf (v) is known and the spectral density S (R, v) at a given frequency in the image plane can be measured (eg, using an energy sensitive detector), equation (29) is different. It can be used for phase retrieval from two such measurements at wavelength. In fact, for example, S (R, v) and S (R, v₀) Is measured, for example, the equation (29) is converted to v and v₀From these expressions, exp {-σ^ThreeM₀Excluding), the following equation can be derived.

Phase Ψ₀(r_⊥) Can be regenerated (or recovered) by solving the Poisson equation (30). If the image is captured by an undisturbed primary beam, any uniform boundary condition is Ψ₀Note that can be assigned to. The solution to equation (30) with typical uniform boundary conditions is always unique, at least up to any addition constant. The stability of the phase reproduction is the ratio of the two frequencies σ = ν₀It is clear from equations (30) and (31) that is determined by the ratio of / ν and the spectral intensity in equation (31). In an experiment, these values can be selected to be sufficiently different from 1 when phase recovery is stable with respect to errors in the experimental data using equation (30). This technique therefore has advantages over the above method for phase reconstruction from the intensity distribution at the cross-section of the beam at two different positions along the optical axis. Indeed, in this latter case, the distance in the z direction between the cross-sections needs to be small and the difference between the two intensity distributions also needs to be small, so the calculation of the phase Laplacian Is a numerically unstable process of the type “divide 0 by 0”.
In the above embodiment of the present invention, an analog value of equation (30) can be obtained that is optimal for in-line phase imaging (imaging) with a small multicolor light source. Due to the quasi-spherical (pseudo-spherical) nature of the incident beam, additional curves appearing in the wavefront can be clearly taken into account. Consider the analog value of “quasisphere” in the following equation (22).
M²S (Mr_⊥, R, ν) -S (r_⊥, ν) =-R’∇_⊥・ [S (r_⊥, ν) ∇_⊥Ψ (r_⊥, ν)]}, (32)
Where M = (R₁+ R) / R₁Is a coefficient of magnitude. R ′ = (R₁ ^-1+ R^-1)^-1And R₁Is the distance from the light source to the object. Starting from equation (32) and proceeding as described above, it is easy to obtain a Poisson equation for a similar phase for equation (30), where Ψⁱⁿ(r_⊥, Ν) to Ψ ′ⁱⁿ(r_⊥, Ν) = Ψⁱⁿ(r_⊥, Ν) −r_⊥ ²/ (2R₁), R is replaced with R ', and γ is replaced with γ' in the following equation.

In a slightly different experimental setup that does not use an energy detection type detector, the frequencies ν and ν₀It is possible to perform intensity measurements in the image plane using two different monochromatic incident beams having At this time, the above equations (30) and (32) are still satisfied by replacing all spectrum intensities with corresponding intensities.
We propose to outline a method for retrieving the phase from a monochromatic image obtained using only a conventional (non-energy detection type) detector. Here, it is assumed that the lateral fluctuation of the absorption is weak as follows:

here,

Is ν = ν₀│σ everywhere while mean absorption in samples at^ThreeM₀’(R_⊥) | << 1. constant

Can be determined from knowledge of the overall intensity incident on the sample and the overall intensity of the image. Local fluctuation M with weak absorption₀’(R_⊥) Is assumed to be unknown.
In the first raw approximation, the term M₀’(R_⊥) May be selected to ignore. And the following Poisson equation for the phase is expressed from the equation (29)

Can be easily obtained.

Thus, the phase solves equation (35) and provides the incident radiation characteristics (ie, Sⁱⁿ(ν) and Ψⁱⁿ(ν) is known based on prior experience. ) So that only one multicolor image I (R) ≡I (r_⊥, R) = ∫S (r_⊥, R, ν) dν.
Now, if the local fluctuation M is weakly absorbed₀’(R_⊥) Can be found from two multicolor images obtained using incident beams with different spectral elements. Substituting equation (34) into equation (29) and integrating over the frequency range yields:

Where a = I (R) −∫S^out(ν) [1-R∇_⊥ ²Ψⁱⁿ(ν)] dν,
b = -∫S^out(ν) [1-R∇_⊥ ²Ψⁱⁿ(ν)] σ^Threedν,
c = R∫S^out(ν) σ²dν, and d = −∫S^out(ν) σ^Fivedν.
(Distribution S known in the prior experiment_j ⁱⁿ(ν) and Ψ_j ⁱⁿtwo images I measured with an incident beam having different spectral elements (denoted by (ν), j = 1, 2)_jWhen (R), j = 1,2 is given, M₀'Can be removed from the corresponding pair in equation (36) to obtain:

Here, exactly as described above, all quantities defined for each of the incident beams are used. For example,
a_j= I_j(R) -∫S_j ^out(ν) [1-R∇_⊥ ²Ψ_j ⁱⁿ(ν)] dν, j = 1, 2, etc. Similarly, to obtain an effective pure absorption-contrast image, M₀Can be solved.
These equations show that the X-ray phase can be reconstructed from a multicolor image obtained at a fixed position along the optical axis using the Transport of Intensity equation. It can be performed with a quasi-planar or quasi-spherical paraxial incident wave, whose spectral density and monochromatic component is known as priori. This approach will preferably include a substantially complete characterization of the source prior to the phase retrieval experiment. The above embodiment of the present invention shows that an image of an effective phase-contrast and absorption-contrast can be obtained by various techniques including any one of the following.
(1) From two images obtained with a monochromatic incident beam with two different wavelengths (photon energy).
(2) Measurement of image intensity as a function of X-ray energy using an energy sensitive detector (eg, based on cadmium, manganese, telluride) and a multicolor incident beam, at least two With storage of energy values that should give sufficient image intensity in the energy range (bin: storage).
(3) For example, measurement of image intensity using energy selectivity in a detector and use of a multicolor incident beam by use of an energy filter based on the use of a foil-like absorber. In practice, this would include using a pair of x-ray films with a suitable foil filter between them and simultaneously recording the appropriate intensity data.
(4) From two multicolor images obtained with incident beams having different spectral compositions and from conventional detectors that are not energy sensitive.
At each of these positions described above, if the lateral spatial variation in absorption is negligible, only one image measurement is sufficient for phase recovery. Alternatively, individual effective phase difference and absorption difference images can be obtained from the image intensity data.
FIG. 3 illustrates a modification of the device for obtaining a plurality of intensity records for different energy distributions. Such a measurement can be performed, for example, by an energy sensitive detector D._EWith a light source that can be tuned by tuning the wavelength shape with an adjustable controller C, or with a multicolor light source and a replaceable monochromatic light source or object (F₁) Or a filter placed either after the object, for example a shutter F placed in front of an energy insensitive CCD detector D₂So it can be executed. If the absorption in the sample is negligible, the phase can be found from the intensity of only one multicolor image and the resulting incident radiation characteristic is known as priori.
Another approach to the use of two or more different spectral distributions (ie, an alternative to the approach based on the intensity transport equation) is the deconvolution of the spectral distribution from the data, or an alternative²A least square fitting method for extracting Ψ may be included. More sophisticated methods based on Bayesian analysis could also be used. Here it is proposed to briefly outline the recovery of phase and absorption by techniques based on the Fourier optics approach.
The phase recovery process in the case of a spherical wave is similar to the process for plane waves and perfect straight travel. Thus, in this section, only the latter is described.
In the case of plane waves, it can be proved that the Cowley form of Kirchhoff's equation for the X-ray wave function is simplified with an approximation of a small angle as follows.

Here we have real space functions present in the lower case and the corresponding frequency-space functions (after performing the Fourier transform) are higher, eg Q (f_xf_y) = F [q (x, y)] is adopted. For further discussion of the steps involved in this guidance, see J. See diffraction physics by North Cowley (North Holland: Amsterdam). A related theoretical argument is also found in 1997, Pogany et al, Rev. Sci. Inst. 68, 2774-2782.
Assuming that φ (x, y) t (x, y) and μ (x, y) t (x, y) are sufficiently small, the transfer function is q (x, y) ≈1 + iφ ′ (x, y) −μ ′ (x, y). Here, φ ′ (x, y) = − φ (x, y) t (x, y) and μ ′ (x, y) = μ (x, y) t (x, y) / 2 It is. The Fourier transform of the transfer function is then followed by Q (f_x, F_y) ≒ δ (f_x, F_y) + IΦ ’(f_x, F_y) -M '(f_x, F_y). Equation (39) can then be written as:

Where χ (f_x, F_y) = ΠλR₂(f_x ²+ F_y ²). By expanding the equation (40) and performing the inverse Fourier transform on both sides, the function dependence can be removed for the first order terms in φ ′ and μ ′.

As a result, in the case of a pure phase object, the φt-distribution can be recovered from a single image.

And in the case of a pure absorbing object:

In the more usual case with objects where both phase and absorption effects are important, two images I₁And I₂Can be used, which are different values of R₂And / or collected with x-ray energy (we are discussing here the monochromatic case, but in the multi-color case, the x-ray source is set to two different voltages and two different spectral distributions Would be appropriate). Without any absorbing edge against the component, μ₂= (Λ₂/ Λ₁)^Threeμ₁And φ₂= (Λ₂/ Λ₁) φ₁As a result, the equation can be solved simultaneously, and the following equation is obtained.

The use of the multiple energy method of phase recovery is advantageous over the multiple distance method in that it includes the following matters.
Data can be recorded by electronic changes rather than mechanical changes in image conditions (eg by rapid switching of energy sensitive detector or tube (X-ray tube) voltage) and very rapid phase recovery Is possible.
And if available, a two-dimensional energy dispersion detector can be used to record multiple two-dimensional images for different X-ray energy ranges (eg, based on CdMnTe) ).
Example 1
In order to test the effectiveness of the method according to the invention using the intensity transport equation, numerical examples are prepared with simulation data and the successive steps are shown in FIGS.
Wavelength λ = 0.154 nm (CuK_aRadiation) and the distance between the source and the object and the distance between the object and the image R₁= R₂= 0.2 m, thus R '= 0.1 m and M = 2, and was assumed to be a perfect monochromatic X-ray point source (FIG. 2). Also, 640 × 640μm²It was also assumed that this is a pure phase object that produces a phase distribution on the object plane (x, y, 0) shown in FIG. The range of phase values was [−0.8667, 0.6706] radians.
Intensity distribution I (x, y, R on the image plane₂) Was then calculated by computing the Fresnel integral (on 7) on a 128 × 128 pixel grid. The obtained intensity distribution I (x, y, R₂) Scaled image M²I (x / M, y / M, R₂) Is shown in FIG. For comparison, the calculated distribution of the phase Laplacian, − ▽²ψ (x, y, ()) is plotted and is shown in FIG. Obviously, the scaled intensity distribution on the image plane and the phase Laplacian distribution on the object plane are almost equal, and this is valid when the absorption is negligible. Is backed up.
In the next stage, the calculated intensity distribution on the image plane and Equation (20) were used to reconstruct the phase. The result of this reconstruction is shown in FIG. Again, this is in good agreement with the original phase shown in FIG. 4 and the relative mean square error of this reconstruction was 4.65%.
The implementation of this method was also tested for non-uniform intensity cases, ie where the last term in equation (16) is important. The intensity distribution on the object plane was simulated according to the following equation.
I (x, y) = − 1.647 + 3 × exp {− (x²+ Y²) / (2 * 1280²)}
The intensity value had a range of [1,1.1243]. Intensity distribution I (x, y, R on the image plane₂) Was computed by evaluating the Fresnel integral. The non-uniformity of the intensity distribution causes “vignetting” in the corresponding image (FIG. 8), which clearly reflects the non-uniformity of intensity on the object plane. However, by subtracting the intensity distribution on the object plane from that on the image plane, dividing the distance by the scaled interplane distance R ′ = 0.1 m, and applying equation (2), we have The original phase could be reconstructed with a root mean square error of 6.67% (FIG. 9). The result is that if the intensity distribution non-uniformity on the object plane is of a normal degree (for example, if the incident beam non-uniformity and / or the absorption effect is weak), the above simple expression of equation (20) It has been confirmed that the deformation can be used effectively for quantitative recovery of the phase change created by the object.
Example 2
Figures 10-15 show different values of R₂3 shows a specific example of phase and absorption recovery from two images. As the distribution of absorption (μt) and phase (φt), the CSIRO logo (term: logo) is used (see FIGS. 10 and 11, respectively). For the former, the value is either 0.0 (black) or 0.1 (white), and for the latter (reverse logo), either -0.1 (black) or 0.0 (white) It is. The image (intensity distribution) is 1 放射線 radiation, 0.2 x 0.2 mm²(512 x 512 pixels) area, and R₂= 80 cm (FIG. 12) and R₂= 160 cm (FIG. 13) and calculated for the plane wave case. These conditions were such that the image was not immediately decipherable. Applying the discussed recovery algorithms (44, 45), we have obtained the results shown in FIGS. 14 and 15, which are not complete, but are clearly recognizable.

Claims (33)

A method for obtaining a phase change image caused by an object in transparent radiation incident on an object that is a sample material,
Irradiating said object with transparent radiation having a high spatial coherence transverse to the direction of radiation propagation;
After radiation is emitted from the object, at least a part of the radiation is received by the detection means, and at least two intensity records are obtained and stored for the received radiation, each intensity record being in the plane of the detection means Having intensity values at predetermined spatial intervals in two dimensions ;
Using the intensity value, deriving values arranged in a grid in a vertical and horizontal direction defining a phase change image caused by the object in the penetrating radiation;
After the radiation is emitted from the object, the intensity record is generated for each different energy distribution of the detected radiation at a fixed finite distance from the object in the radiation propagation direction . The method of claim 1, wherein each of the different energy distributions is obtained by changing an energy spectrum of radiation irradiating the object. The method according to claim 1, wherein each of the different energy distributions is obtained by providing the sensing means for providing an intensity value as a function of energy in one or more energy bands. The derivation includes solving one or more intensity transfer differential equations that indicate the relationship between the phase in one plane of the object and the gradual change of the radiation intensity distribution along the radiation propagation direction, the equation comprising: The method according to claim 1, 2 or 3, which can be solved using a predetermined uniform boundary condition. 4. The derivation according to claim 1, 2 or 3, wherein the derivation comprises solving a Fourier optical equation showing the relationship between the phase in one plane of the object and the gradual change of the radiation intensity distribution along the radiation propagation direction. Method. The intensity value indicates the absorption contrast in the object;
6. The method of claim 1, further comprising utilizing the intensity value to derive a longitudinal and horizontal grid of values defining an effective and pure absorption contrast image of the object. The method according to any one of the above. 7. A method according to any one of claims 1 to 6, wherein the penetrating radiation comprises x-ray radiation. The method according to claim 7, wherein the X-ray radiation is in the range of 0.5 KeV to 1 MeV. 9. A method according to claim 7 or 8, wherein the radiation is substantially monochromatic. 9. A method according to claim 7 or 8, wherein the radiation is polychromatic. 11. A method according to any one of claims 7 to 10, wherein the radiation is radiation from a substantially point radiation source having a full width at half maximum of 40 [mu] m or less. An apparatus for obtaining a phase change image produced by said object in transparent radiation incident on said object, comprising:
A radiation source for irradiating said object with transparent radiation having a high spatial coherence transverse to the direction of radiation propagation;
Detecting means for receiving at least a part of the radiation after the radiation is emitted from the object;
The detection means generates at least two intensity records for the received radiation, each intensity record having an intensity value at a predetermined two-dimensional spatial interval in the plane of the detection means ;
The detection means performs the intensity recording at a certain finite distance from the object in the radiation propagation direction after radiation is emitted from the object ,
The apparatus further comprises energy characteristic setting means, and the energy characteristic setting means generates the intensity record for each different energy distribution of the detected radiation. The apparatus according to claim 12, wherein the energy characteristic setting unit changes an energy spectrum of radiation that irradiates the object. 14. The apparatus of claim 13, wherein the energy characteristic setting means includes means for allowing the sensing means to provide intensity values as a function of energy in one or more energy bands. Further comprising a computer program having a plurality of machine-readable instructions,
When the computer program is installed in a computer having an operating system and storage means, the computer program uses the intensity value to define a phase change image that defines a phase change image generated by the object in penetrating radiation. 15. An apparatus according to claim 12, 13 or 14, wherein the apparatus is operated to derive values arranged in One or more intensity transfer differential equations wherein the derivation indicates a relationship between a phase in one plane of the object and a gradual change in intensity distribution along the radiation propagation direction using a predetermined uniform boundary condition 16. The apparatus of claim 15, comprising solving. 16. The apparatus of claim 15, wherein the derivation comprises solving a Fourier optical equation that indicates the relationship between the phase in one plane of the object and a gradual change in intensity distribution along the radiation propagation direction. 18. Apparatus according to any one of claims 12 to 17, further comprising an x-ray radiation source as the radiation source for illuminating the object. The apparatus of claim 18, wherein the x-ray radiation is in the range of 0.5 KeV to 1 MeV. 20. Apparatus according to claim 18 or 19, wherein the radiation is substantially monochromatic. 20. An apparatus according to claim 18 or 19, wherein the radiation is multicolored. The apparatus according to any one of claims 18 to 21, wherein the radiation source is a point-like radiation source having a substantially full width at half maximum of 40 µm or less. After the transmitted radiation passes through the object that is the sample material, a phase change image caused by the object in the transmitted radiation incident on the object is acquired from one two-dimensional intensity record of the transmitted radiation. A method,
When incident on the object, the radiation has a high spatial coherence transverse to the direction of radiation propagation;
Obtaining said intensity Symbol recording, after the radiation from the object conjugated phase variation element comes uniformly in a field of radiation of the surroundings without phase fluctuating or phase change, in the finite distance from the object in the radiation propagation direction And
The method
At predetermined spatial intervals of two-dimensional in the plane of the sensing means, and storing the intensity values of the intensity Symbol recording or al,
By using the intensity value and a predetermined uniform boundary condition by solving an intensity transfer differential equation showing the relationship between the phase of one plane in the object and the gradual change of the intensity distribution along the radiation propagation direction. Deriving from the intensity values values arranged in a grid in a horizontal and vertical direction defining a phase change image produced by the object in penetrating radiation. A method for acquiring a phase change image caused by an object in a transparent radiation incident on an object that is a sample material,
The method
Irradiating the object with transparent radiation having a high spatial coherence transverse to the direction of radiation propagation;
One or more finite distances in the radiation propagation direction from the object in the surrounding radiation field with uniform or non-phase variation after the radiation exits the object conjugated with a phase variation element Receiving at least a portion of the radiation at a detector, thereby obtaining and storing intensity values for the received radiation at a predetermined two-dimensional spatial spacing in a plane of the detector ;
By using the intensity value and a predetermined uniform boundary condition, solving an intensity transfer differential equation showing the relationship between the phase of one plane in the object and the gradual change of the intensity distribution along the radiation propagation direction Deriving from said intensity values values arranged in a grid in a vertical and horizontal direction defining a phase change image produced by said object in penetrating radiation. The intensity value also reflects the absorption contrast in the object,
25. The method of claim 23 or 24, further comprising the step of deriving values arranged in a grid pattern to define an effective and pure absorption contrast image of the object using the intensity value. Method. 26. A method according to claim 23, 24 or 25, wherein the penetrating radiation comprises x-ray radiation. 27. The method of claim 26, wherein the x-ray radiation is in the range of 0.5 keV to 1 MeV. 28. A method according to claim 26 or 27, wherein the radiation is substantially monochromatic. 28. A method according to claim 26 or 27, wherein the radiation is polychromatic. 30. The method of claim 29, wherein the intensity transfer differential equation includes a spectrum weighted term or a spectrum weighted factor that depends on the square of each wavelength component. The boundary condition includes a uniform Dirichlet boundary condition, a uniform Neumann boundary condition, or a uniform periodic boundary condition so as to obtain a unique solution of the equation for a phase leading to at least an arbitrary constant component. 31. The method according to any one of claims 23 to 30, wherein The solution may further include one or more optical elements selected from the group consisting of uniform radiation of the object, no focus between the object and the image, and a small wavefront curvature for the incident radiation. 32. The method of claim 31, wherein a condition is utilized. An apparatus for obtaining a phase change image caused by said object in penetrating radiation incident on the object, comprising:
A radiation source for irradiating said object with transparent radiation having a high spatial coherence transverse to the direction of radiation propagation;
At least a portion of the radiation at a finite distance from the object in a radiation propagation direction in a surrounding radiation field that is not phase-shifted or phase-shifted uniformly after the radiation is emitted from the object conjugated with a phase-varying element A detector for receiving and thereby generating an intensity value for said received radiation at a predetermined two-dimensional spatial spacing in the plane of the detector;
A computer including a program stored with mechanically readable instructions, using the intensity value and a predetermined uniform boundary condition, along a phase of one plane and a radiation propagation direction in the object A computer for deriving values arranged in a grid pattern to define a phase change image provided by the object in the penetrating radiation by solving an intensity transfer differential equation showing a relationship with a gradual change of the intensity distribution And a device comprising:

Images