(57)【要約】【課題】  既存の装置に対し改良を必要とせず少なくとも一回の撮影で高い空間分解能を持つ被検体の位相を取得可能な位相撮影方法と装置を提供する。【解決手段】  線源からの量子ビームを、被検体と少なくとも１枚の位相格子を介して検出器に入射し、当該検出器を構成する各ピクセルでのビームの強度分布から前記被検体の位相画像を求める位相画像撮影方法であって、ビームの強度分布は、少なくとも、吸収（a0）、Visibility（V）、そして、位相（φ）の情報を含んでおり、隣接する少なくとも３つのピクセルでの吸収（a0）、Visibility（V）、そして、位相（φ）を画像の変数近似によりほぼ同一であると仮定し、少なくとも１枚の測定画像から、吸収、Visibility、位相を求める位相画像撮影方法と装置。【選択図】  図１
2019240165.tif 000002
【特許請求の範囲】
【請求項１】
線源からの量子ビームを、被検体と少なくとも１枚の位相格子を介して検出器に入射し、当該検出器を構成する各ピクセルでのビームの強度分布から前記被検体の位相画像を求める位相画像撮影方法であって、
ビームの強度分布は、少なくとも、吸収（a0）、Visibility（V）、そして、位相（φ）の情報を含んでおり、
隣接する少なくとも３つのピクセルでの吸収（a0）、Visibility（V）、そして、位相（φ）を画像の変数近似によりほぼ同一であると仮定し、
少なくとも１枚の測定画像から、吸収、Visibility、位相を求めることを特徴とする位相画像撮影方法。
【請求項２】
請求項１に記載した位相画像撮影方法であって、
前記検出器の各ピクセル(m,n)でのビームの強度分布は、以下のように表され、
【数１５】
2019240165.tif 000017
ただし、a0は吸収、VはVisibility、φは位相、ｆは背景による干渉縞を表わしており、かつ、
前記画像の変数近似において、吸収を…、a0'、a0、a0"…、Visibilityを…、V'、V、V"…、位相を…、φ'、φ、φ"…としたとき、以下のようにして、
【数１６】
2019240165.tif 000018
少なくとも１枚の測定画像から、吸収、Visibility、位相を求めることを特徴とする位相画像撮影方法。
【請求項３】
請求項１に記載した位相画像撮影方法において、
前記吸収、Visibility、位相を、測定値と理論値との差のLp(p≧0)ノルム最小化を用いて求めることを特徴とする位相画像撮影方法。
【請求項４】
請求項３に記載した位相画像撮影方法において、前記Lpノルム最小化において、p≧1に対して、変数変換によって非凸最小化問題を凸関数最小化問題に変換して、吸収、Visibility、位相を求めることを特徴とする位相画像撮影方法。
【請求項５】
線源からの量子ビームを、被検体と少なくとも１枚の位相格子を介して検出器に入射し、当該検出器を構成する各ピクセルでのビームの強度分布から前記被検体の位相画像を求める位相画像撮影方法であって、
予め与えられた近似位相像を用いて、各ピクセルに対して求めるピクセルの位相値（又は、構造）との類似度によって重みを与える重み付け処理を施し、その後、
前記の重みを用いて、少なくとも１枚の測定画像から、吸収、Visibility、位相を求めることを特徴とする位相画像撮影方法。
【請求項６】
請求項５に記載した位相画像撮影方法において、
位相値（又は構造）が類似する範囲において、吸収を…、a0'、a0、a0"…、Visibilityを…、V'、V、V"…、位相を…、φ'、φ、φ"…としたとき、以下のようにして、
【数１７】
2019240165.tif 000019
少なくとも１枚の測定画像から、吸収、Visibility、位相を求めることを特徴とする位相画像撮影方法。
【請求項７】
請求項６に記載した位相画像撮影方法において、
測定値と理論計算値との差の重み付けのLp(p≧0)ノルム最小化を用いて、吸収、Visibility、位相を求めることを特徴とする位相画像撮影方法。
【請求項８】
請求項７に記載した位相画像撮影方法において、
前記重み付けLpノルム最小化において、p≧1に対して、変数変換によって非凸最小化問題を凸関数最小化問題に変換して、吸収、Visibility、位相を求めることを特徴とする位相画像撮影方法。
【請求項９】
請求項５に記載した位相画像撮影方法において、前記予め与えられた近似位相像は、前記請求項１から５に記載した方法のいずれかによって得られた位相画像であることを特徴とする位相画像撮影方法。
【請求項１０】
請求項５に記載した位相画像撮影方法において、前記予め与えられた近似位相像は、Fourier変換法により求められる位相画像であることを特徴とする位相画像撮影方法。
【請求項１１】
請求項５に記載した位相画像撮影方法において、前記予め与えられた近似位相像は、別の装置で撮影した同じ物体の画像、或いは、前回の測定データ、或いは、構造が似ていることが既知の別の物体の画像であることを特徴とする位相画像撮影方法。
【請求項１２】
請求項５に記載した位相画像撮影方法において、前記重み付けを、Non-local means を用いた計算で、あるいは、バイラテラルフィルタ(Bilateral filter)などのフィルタで、あるいは、画像が類似しているほど大きい重みが求められる手法で行うことを特徴とする位相画像撮影方法。
【請求項１３】
請求項５に記載した位相画像撮影方法において、前記近似位相像に代えて、前記被検体の不連続な場所の吸収撮影画像を用いることを特徴とする位相画像撮影方法。
【請求項１４】
量子ビームを発生する線源と、
前記線源からのビームを照射する被検体を保持する手段と、
前記被検体からのビームを、少なくとも１枚の位相格子を介して受光する検出器とを備え、更に、
前記検出器を構成する各ピクセルでのビームの強度分布から前記被検体の位相画像を求める処理部とを備えた位相画像撮影装置において、
前記処理部は、前記請求項１〜４のいずれか１項に記載した位相画像撮影方法を実行することにより、吸収、Visibility、位相を求めることを特徴とする位相画像撮影装置。
【請求項１５】
量子ビームを発生する線源と、
前記線源からのビームを照射する被検体を保持する手段と、
前記被検体からのビームを、少なくとも１枚の位相格子を介して、受光する検出器とを備え、更に、
前記検出器を構成する各ピクセルでのビームの強度分布から前記被検体の位相画像を求める処理部とを備えた位相画像撮影装置において、
前記処理部は、前記請求項５〜１３のいずれか１項に記載した位相画像撮影方法を実行することにより、吸収、Visibility、位相を求めることを特徴とする位相画像撮影装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、Ｘ線をはじめとする量子ビームを照射して得られる位相の変化を用いて被検体の画像を得る位相画像撮影方法に関し、更には、それを用いた位相画像撮影装置に関する。
【背景技術】
【０００２】
Ｘ線は、可視光とは異なり、物質を透過することができることから、Ｘ線で物体の内部構造を観察することが可能になっており、今日においては、病院でのＸ線撮影や工業の非破壊検査など、幅広い分野に応用されている。通常のＸ線撮影では、Ｘ線が物体を透過する際に物体の吸収によってＸ線の強度が変化することで、コントラスト画像を形成している。このような吸収率を利用している装置では、物体を構成する元素が重いほどＸ線の吸収は大きくなり、コントラストが高い画像が得られる。しかし、生体の軟組織などの軽元素（原子量の小さい）からなる弱吸収物体を観察する場合には、コントラストがつき難いという問題がある。【０００３】
これに対して、以下の先行技術文献に示すように、Ｘ線位相イメージング法の研究が1990年代から行われている。
【先行技術文献】
【特許文献】
【０００４】
【特許文献１】国際公開ＷＯ２００４０５８０７０号
【特許文献２】国際公開ＷＯ２０１００５０４８３号
【非特許文献】
【０００５】
【非特許文献１】A. Bravin, P. Coan, and  P. Suortti, "X-ray phase-contrast imaging: from pre-clinical applications towards clinics," Phys. Med. Biol. 58, R1 (2013).
【非特許文献２】F. Pfeiffer, C. Grunzweig, O. Bunk, G. Frei, E. Lehmann, and C. David, "Neutron Phase Imaging and Tomography," Phys. Rev. Lett. 96, 215505 - Published 2 June 2006; Erratum Phys. Rev. Lett. 97, 069905 (2006).
【非特許文献３】A. Momose, S. Kawamoto, I. Koyama, Y. Hamaishi, K. Takai, and Y. Suzuki, "Demonstration of X-Ray Talbot interferometry," Jpn. J. Appl. Phys. 42, L866-L868 (2003).
【非特許文献４】T. Weitkamp, A. Diaz, C. David, F. Pfeiffer, M. Stampanoni, P. Cloetens, and E. Ziegler, "X-ray phase imaging with a grating interferometer," Opt. Express 12(16), 6296-6304 (2005).
【非特許文献５】M. Takeda, H. Ina and S. Kobayashi, "Fourier-transform method of fringe-pattern analysis for computer-based topography and interferometry," J. Opt. Soc. Am. 72, 156-160 (1982).
【非特許文献６】Y. Ge, K. Li, J. Garrett, and G. H. Chen, "Grating based x-ray differential phase contrast imaging without mechanical phase stepping,"Opt. Express 22, 14246-14252 (2014).
【非特許文献７】C. Kottler, F. Pfeiffer, O. Bunk, C. Gr¨unzweig, and C. David, "Grating interferometer based scanning setup for hard X-ray phase contrast imaging," Rev. Sci. Instrum. 78, 043710 (2007).
【非特許文献８】C. Arboleda, Z. Wang Z, and M. Stampanoni, "Tilted-grating approach for scanning-mode X-ray phase contrast imaging," Opt. Express 22(13), 15447-15458 (2014).
【発明の概要】
【発明が解決しようとする課題】
【０００６】
これらの手法では、Ｘ線が物体を通過する時生じる位相の変化を測定している。即ち、Ｘ線が軽元素からなる物体を通過する時の位相の変化が、吸収の変化よりはるかに大きいことを利用して、高いコントラストの画像を生成している[非特許文献１]。また、これらの手法は、近年、中性子線などの線源へも応用されている[非特許文献２]。【０００７】
開発された数多くの撮影手法の中でも、特に、回折格子を用いてTalbot干渉を利用する位相イメージング法は、通常のＸ線管を使うことが可能で病院などで使えるコンパクトな装置を作ることが可能になっていることから、注目されている。【０００８】
一方、位相撮影方法においては、撮影された強度分布画像から位相を復元する手法として縞走査法が知られている[非特許文献３、４、特許文献１]。この手法では、回折格子を格子周期方向に、少なくとも、3回移動して、複数枚の画像撮影を行っている。なお、毎回の移動距離はほぼサブミクロンレベルであることから、装置には非常に高い精度／安定性が求められると共に、使用する周囲環境においても温度変化などの厳しい管理が要求される。そのうえ、３回以上の移動が必要なので、動いている（或は変化している）物体の動的測定、例えば、細胞培養の視覚モニタリング などへの応用は難しくなる。【０００９】
そのためFourier変換法が提案されている[非特許文献５、特許文献２]。この方法では、撮影された干渉縞像にFourier変換を用いた処理を行うことによって位相を計算しており、1回の撮影で位相を取得することができる。しかし、Fourier変換法では多くの高周波成分が失われることから空間分解能が低下してしまう。このことは、被検体の細かい構造などが失われることを意味する。【００１０】
この他にも、特殊な回折格子を用いる方法[非特許文献６]や一次元センサーを用いる方法[非特許文献７、８]などが提案されているが、いずれも、新しいハードウェアが必要である。【００１１】
そこで、本発明は、上述した従来技術に鑑みてなされたものであり、より具体的には、既存の装置や設備に対する改良を必要とせず（特に、ハードウェアの改良／変更なしに）、少なくとも一回の撮影により高い空間分解能を持つ被検体の位相を取得することができる位相撮影方法、及び、それを利用した位相撮影装置を提供することをその目的とする。
【課題を解決するための手段】
【００１２】
上記の目的を達成するために、本発明によれば、まず、線源からの量子ビームを、被検体と少なくとも１枚の位相格子を介して検出器に入射し、当該検出器を構成する各ピクセルでのビームの強度分布から前記被検体の位相画像を求める位相画像撮影方法であって、ビームの強度分布は、吸収（a0）、Visibility（V）、そして、位相（φ）の情報を含んでおり、隣接する少なくとも３つのピクセルでの吸収（a0）、Visibility（V）、そして、位相（φ）を画像の変数近似によりほぼ同一であると仮定し、少なくとも１枚の測定画像から、吸収、Visibility、位相を求める位相画像撮影方法が提案される。【００１３】
また、本発明によれば、上記の目的を達成するため、線源からの量子ビームを、被検体と少なくとも１枚の位相格子を介して検出器に入射し、当該検出器を構成する各ピクセルでのビームの強度分布から前記被検体の位相画像を求める位相画像撮影方法であって、予め与えられた近似位相像を用いて、各ピクセルに対して求めるピクセルの位相値（又は、構造）との類似度によって重みを与える重み付け処理を施し、その後、上記の重みを用いて、少なくとも１枚の測定画像から、吸収、Visibility、位相を求める位相画像撮影方法が提案される。【００１４】
加えて、本発明によれば、量子ビームを発生する線源と、前記線源からのビームを照射する被検体を保持する手段と、前記被検体からのビームを、少なくとも１枚の位相格子を介して、受光する検出器とを備えており、更に、前記検出器を構成する各ピクセルでのビームの強度分布から前記被検体の位相画像を求める処理部とを備えた位相画像撮影装置において、前記処理部は、前記前者に記載した位相画像撮影方法である位相画像撮影装置が提案される。
【発明の効果】
【００１５】
上述した本発明によれば、既存の装置や設備に対して改良を必要とせず（特に、ハードウェアの改良／変更なしに）、少なくとも一回の撮影により高い空間分解能を持つ被検体の位相を取得することができる位相撮影方法、及び、位相撮影装置が提供されるという実用的に優れた効果が発揮される。
【図面の簡単な説明】
【００１６】
【図１】本発明の一実施の形態（実施例１）である手法Ａの原理を説明する図である。
【図２】上記手法Ａの原理を説明するための図である。
【図３】上記手法Ａの変形例を説明するための図である。
【図４】上記手法Ａの三つの変数に関する凸関数最小化問題を説明するための図である。
【図５】上記手法Ａが適用可能な設備や装置であるTalbot干渉計の構成を示す図である。
【図６】上記手法Ａが適用可能な設備や装置であるTalbot-Lau干渉計の構成を示す図である。
【図７】上記手法Ａが適用可能な設備や装置であるLau干渉計の構成を示す図である。
【図８】本発明の他の実施の形態（実施例２）である手法Ｂの重みについて説明する図である。
【図９】上記手法Ｂの重みの計算について説明する図である。
【図１０】本発明の手法Ａ、Ｂにより数値ファントムを用いて位相の復元を行った例を、従来技術と比較して示した図である。
【図１１】本発明の手法Ａ、Ｂにより得られたポリマー球の位相画像を、従来技術と比較して示した図である。
【図１２】Talbot干渉計でポリマーを撮影した実データを用いて本発明の手法Ａ、Ｂにより得られた位相画像を、従来技術と比較して示した図である。
【図１３】本発明の位相撮影方法や位相撮影装置が適用可能なTalbot干渉計を利用した高速位相ＣＴ撮影装置の構成の一例を示す図である。
【発明を実施するための形態】
【００１７】
以下、本発明の実施の形態について説明するが、それに先立ち、ここでは、添付の図面を参照しながら、本発明の原理について詳細に説明する。
【実施例１】
【００１８】
＜本発明の原理：手法Ａ＞
Ｘ線、中性子線、可視光などを含め、所謂、量子ビームを発生する線源を用いて干渉縞を形成する装置では、検出器（二次元）により直接測定されるのは、位相ではなく、被検体によって変形したモアレ干渉縞（あるいは、Talbot効果による自己像）の強度分布である。この検出器で測定される強度分布は、一般的に、次のように書くことができる。
【数１】
2019240165.tif 000003
【００１９】
ただし、a0は吸収、VはVisibility、φは微分位相（或いは、差分位相。単に「位相」とも記す）であり、fはモアレ干渉縞の正弦波を表す波形である。位相撮影では、位相を求めるため、それぞれ、被検体がない時の背景の撮影とある時の撮影を行っており、ここでfは背景の撮影から求められる。背景の撮影が不可能な時、撮影される一枚の画像の中で物体の位相の変化が小さい（すなわちモアレ干渉縞の変形が小さい）局部を用いてfを推定することも可能である。【００２０】
このように、式（１）では、検出器の各ピクセル（又は、「画素」とも記す）において三つの変数があり、位相を求めるためには、少なくとも三つの異なるデータが必要である。既述のように、従来の縞走査法では、回折格子をM（M≧3）回移動して撮影することで複数枚の撮影画像を得る。そして、得られた複数枚の撮影データを用い、例えば最小二乗法などで位相を求めている。この手法では高い空間分解能で位相を求めることができるが、しかし、上述したように、回折格子の移動において高い精度が要求される。縞走査法で回折格子の毎回の移動距離は、均等で、かつ、回折格子の周期の1/Mである。回折格子の周期は、一般的に、数μｍであるため、移動はほぼサブミクロンレベルで行う。従って、装置／設備には厳しい安定性などが要求される。【００２１】
本発明では、少なくとも一回の撮影で三つの変数を、即ち、吸収（a0）、Visibility（V）、位相（φ）が求められる方法を提供する。上記で説明したように、位相を求めるためには、少なくとも三つの異なるデータが必要であるが、本発明では、画像の変数近似によりこの問題を解決する。【００２２】
具体的には、図１にも示すように、検出器（二次元）の隣接している三つのピクセルにおいては、上記の式（１）は、次のように、吸収（a0）、Visibility（V）、そして、位相（φ）はほぼ同一であると仮定する。
【数２】
2019240165.tif 000004
【数３】
2019240165.tif 000005
なお、ここでは、以下の説明でも同様であるが、ただ三つのピクセルを用いて説明する。しかし、当然ながらより多くのピクセルを用いて計算することも可能である。【００２３】
すると、上記式（１）から、以下のように書くことができる。これは、ピクセル(m,n)において、一つのデータI-1(m,n)が増えたことを意味する。
【数４】
2019240165.tif 000006
【００２４】
さらに、同じ近似によれば、以下にも示すように、少なくとも三つのデータを得ることができる。
【数５】
2019240165.tif 000007
【００２５】
以上のように、画像の変数近似によれば、三つの変数に対して三つのデータが存在することから、位相（φ）を求めることが可能になる。【００２６】
＜例＞
画像の変数近似、即ち、上記式（２）と（３）の近似によれば、図２にも示すように、一つのピクセルに相当するデータが、一つから三つに増えることとなる。【００２７】
なお、ここで、上述した画像の変数近似は、垂直方向に限られているわけではなく、図３（ａ）〜（ｃ）のように、水平、あるいは、斜めなど、種々の方向についても可能である。ここで説明のために具体的な一例として隣接している三つのピクセルにおいて吸収（a0）、Visibility（V）、そして、位相（φ）はほぼ同一であると仮定している。当然なことながら、図３（ｄ）、（ｅ）のように、周りのもっと多くのピクセルにおいて吸収（a0）、Visibility（V）、そして、位相（φ）はほぼ同一であると仮定して本発明手法で吸収、Visibility、位相を求めることができる。当然なことながら、吸収とVibilityは隣接している五つのピクセルにおいてほぼ同一であると仮定し、位相は隣接している三つのピクセルにおいてほぼ同一と仮定するなどそれぞれ違うピクセル数で変数近似を行うのも可能である。【００２８】
このように、画像の変数近似によって得たデータから最小化問題を解くことで、位相が求められる。具体的には、次のように、上記の式（５）で得られるデータと上記の式（１）の理論値の差の最小化（測定値と理論値との差のLp(p≧0)ノルム最小化）の問題を解くことで求められる。
【数６】
2019240165.tif 000008
ただしIk(k=0,±1)は式(５)と同じ画像の変数近似で得られる複数個のデータである。ここでp=2とすると最小二乗問題になる。【００２９】
この問題の式(６)は、以下の式（７）の初期値を与えて様々な既存の反復手法（例えば、最急降下法など）で解くことが可能である。
【数７】
2019240165.tif 000009
そして、求められる最適解{a0,V,φ}が、それぞれ、求める吸収（a0）、Visibility（V）、そして、位相（φ）になる。【００３０】
しかし、式(６)は変数φの三角関数を含んでいるため、非凸関数の最小化問題となっている。従って、初期値の取り方によっては、誤った解（例えば、局所最小値）へ収束する可能性がある。そこで、本発明では、p≧1の場合について、次のような新しい解析手法を提案する。【００３１】
まず、（a0、V、φ）に対して次のような変数変換を行う。
【数８】
2019240165.tif 000010
【００３２】
この変換により、問題の式(６)は次のように書き変えることができる。
【数９】
2019240165.tif 000011
【００３３】
式(９)は、三つの変数（a0、bc、bs）に関する凸関数の最小化問題になっている。従って、図４に示すように、如何なる初期値から始めても、正しい最適解に収束する保証がある。最適解a0,bc,bsを求めた後、次の変数変換により、吸収（a0）、位相（φ）、Visibility（V）が求められる。
【数１０】
2019240165.tif 000012
【００３４】
従って、上記に詳述した本発明になる位相撮影方法（手法Ａ）によれば、既存の装置や設備を利用して、ハードウェアの大幅な改良を必要とせずに、少なくとも一回の撮影により高い空間分解能を持つ被検体の位相を取得することが可能となる。即ち、従来の縞走査法と比べると：
Ａ．回折格子の移動がないため縞走査法での高い安定性の要求がなくなる。
Ｂ．1回の撮影で済むので装置の操作が簡単になる。
更に、従来のFourier変換法と比べると：
Ｃ．より高い空間分解能を保つことで、より多くの細かい内部構造が見える画像を提供することが可能となる（以下の図１０（ｂ）や図１１（ｃ）を参照）。
Ｄ．また、その解析手法においては、凸関数の最小化問題に変換して解いていることから、与えた具体的な初期値にかかわらず、最適解を提供しており、このことは、特に、実際の応用において非常に重要なことであろう。【００３５】
＜本発明手法が適用可能な設備や装置＞
上述した位相撮影方法（手法Ａ）は、既存の装置や設備にも適用可能であり、例えば図５〜７にもその概略を示すように、Ｘ線源（或いは中性子線源）を用いて干渉縞を形成するTalbot干渉計（図５）、Talbot-Lau干渉計（図６）、あるいは、Lau干渉計（図７）などに適用可能である。【００３６】
Ｘ線Talbot干渉計は、図５に示すように、Ｘ線を発生するＸ線源１０を利用し、Ｘ線の干渉を発生させる位相回折格子２１、被検体Ｓで生じた干渉強度分布をモアレ強度分布に変換する吸収回折格子２２、そして、Ｘ線検出器３０から構成されている。なお、この装置では干渉により位相格子に似ている周期画像（自己像）が格子から一定距離だけ離れたところに生じ、この現象はTalbot効果と呼ばれる。そして自己像が生じる位置に二番目の回折格子をおいてモアレ干渉縞を形成している。従ってTalbot干渉計で直接測定されるのは位相ではなく被検体によって変形したモアレ干渉縞の強度分布であり、画像検出器で測定される強度分布I(m,n)は、例えば、演算素子であるＣＰＵや、記憶手段であるＲＡＭやＲＯＭ等を含んだ画像データ処理部に送られ、所定の演算処理によって吸収（a0）、Visibility（V）、そして、位相（φ）として出力されることとなる。なお、これらの干渉計において干渉縞を検出する二次元の画像検出器３０は、多数の半導体Ｘ線検出素子であるピクセルを平面状に配置して構成されており、例えば、特開２００２−２６３００号などに開示された電磁波検出器及び画像検出器を採用すればよい。そして、各ピクセルからは、上述したＸ線強度（I）が出力される。【００３７】
また、図６に示すTalbot-Lau干渉計では、更に、Ｘ線源１０の前に回折格子２３を備えており、そして、図７に示すLau干渉計では、Ｘ線源と被検体の間に位相格子２４を備えており、これらの干渉計でも、上記と同様に、画像検出器３０で測定される強度分布I(m,n)は、例えば、演算素子であるＣＰＵや、記憶手段であるＲＡＭやＲＯＭ等を含んだ画像データ処理部１００に送られ、所定の演算処理によって吸収（a0）、Visibility（V）、そして、位相（φ）として出力されることとなる。【００３８】
このように、本発明の位相撮影方法は、既存の装置や設備にも容易に適用可能であり、画像データ処理部による当該装置や設備により得られた撮影データ（強度信号(I)）を用いた所定の演算により、少なくとも一回の撮影から、高い空間分解能を持つ被検体の位相を取得することが可能となる。
【実施例２】
【００３９】
＜手法Ｂ＞
続いて、本発明では、予め不完全な位相情報が既知の場合に、少なくとも一回の撮影により、上述した手法Ａよりも精度のより高い解が求められる方法を提供する。ここで、上記の既知の位相情報とは、従来の手法であるFourier変換法、あるいは、上述した本発明手法Ａ等により求められる位相像から取得することができる。【００４０】
次の本発明手法Ｂの説明では、この既知の位相情報を「近似位相像」と呼ぶ。【００４１】
この手法でも、発明手法Ａと同様に、最小化問題を解いて位相を求める。上記発明手法Ａとの違いは、ここでは、重み関数W(m,n)を用いて重み付きエネルギー関数の最小化問題を解いている点である。ここで任意の一つのピクセル(m,n)を固定して手法の説明をする。ほかの場所も同じく計算することができる。この重みは、近似位相像を用いて、ピクセル(m,n)において、広い範囲で位相値（又は構造）が近い点を探し、類似している程より大きい重みを与えている。そして、以下の重み付きの最小化問題を解く。
【数１１】
2019240165.tif 000013
【００４２】
ただし、Il,k(l,k=0,±1,…)は近似式によって（上記の式（５）のように）得られる複数個のデータであり、bCとbSは式（８）で定義された変数である。重みWl,kは、以下のように定義する。
【数１２】
2019240165.tif 000014
ここでW(m+l,n+k)(l,k=0,±1,…)は位相を求める注目画素(m,n)と画素(m+l,n+k)(l,k=0,±1,…)の位相値の類似度を表し、その計算方法は後述する。【００４３】
ここで重みの効果について説明する。例えば、図８のように画素(m,n)が物体と物体の境界であり、隣接する画素(m,n-1)の位相値との間に大きな差があるとする。この時、上記の式（３）の仮定は正しくない。従って、画素(m,n-1)のデータを使って求められる位相像も誤差が大きくなる。【００４４】
本発明の手法Ｂでは、この問題を重みWl,kを用いて緩和する。具体的には、まず、図９に示すように、近似位相像を用いて重みを計算する。ピクセル(m,n)とピクセル(m,n-1)では、その位相値の差が大きく、互いに類似していないことから、重みは小さくなる。したがって、上記の式（１２）では、重みW0,-1は以下のようになる。
【数１３】
2019240165.tif 000015
【００４５】
すると、上記の式（１１）でこの重みを使うことによりピクセル(m,n-1)のデータは計算から排除される。これは、上記の式（３）の近似がより正しいものに変換されたことを意味する。結果として、不連続な境界があってもより正確な位相が求められる。【００４６】
続いて、上述した重み関数W(m,n)の計算方法を述べる。重みの計算方法としては、様々な重みを使うことが可能である。その一例として、以下では、Non-local means を用いた具体的な重み計算方法について説明する。【００４７】
例えば、図９のような近似位相像があるとする。注目ピクセル(A5)における位相値を求める状況を考える。このとき、周りの各ピクセルの重みは、近似位相像の位相値（又は、構造）によって決められる。例えば、図９にも示すように、各ピクセル(B5)の重みは、注目ピクセル(A5)を中心とした対象ブロックと(B5)を中心とした参照ブロックを用いることにより、次のように決められる。
【数１４】
2019240165.tif 000016
ただし、Cは作用する全部のピクセルの重みの総和を１にするための正規化定数である。【００４８】
ここで、h1>0は事前に定めたパラメータである。式（１４）はブロック間の類似度に応じた重みであり、位相値（又は、構造）が類似している程より大きな重みが与えられる。【００４９】
＜実験結果＞
次に、上述した本発明の位相撮影方法及び位相撮影装置による効果を確認するため、数値ファントムとＸ線微分位相像を用いて本発明の手法Ａ及びＢにより処理した。コントラスト画像の作り方について以下に説明する。この実験では、まず、仮定した物体の位相（φ）と吸収（a0）、Visibility（V）を用いて式(1)により１枚の干渉縞データを生成する。その干渉縞画像を使って、Fourier変換法と本発明（手法Ａ又は手法Ｂ）によりそれぞれ求めた位相像を示す。【００５０】
図１０は、人体内部のＣＴ画像を模擬した数値ファントムを用いて位相復元を行った例を示しており、図１０(a)は、従来の手法であるFourier変換法で求めた位相像、図１０(b)は、本発明の手法Ａで得られた位相像、そして、図１０(c)は、本発明の手法Ｂで得られた位相像である。これらの位相像から明らかなように、本発明の手法Ａ又は手法Ｂで得られた位相像は、従来のFourier変換法で求めた位相像よりも明らかに空間分解能が上がっており、更に、手法Ａで得られた位相像に比較し、手法Ｂで得られた位相像のほうが境界と細部が正確に求められていることが分かる。【００５１】
図１１は、ポリマー球（直径5/16インチ）の微分位相像を用いてシミュレーションを行った結果であり、(a)の被検体の位相像（正解位相像）と共に、(b)のFourier変換法で求めた位相像、そして、(c)本発明の手法Ａによる位相像、(d)本発明の手法Ｂで得られた位相像を示している。この結果からも明らかなように、図１１(b)の従来手法（Fourier変換法）では解像度が大幅に下がる。一方、図１１(c)の本発明の手法Ａでは、全体的に解像度が上がっており、更に、図１１(d)の本発明の手法Ｂでは境界までほぼ正しく復元されている。【００５２】
図１２は、Talbot干渉計で撮影したポリマーの実データを用いて位相復元を行った結果を示している。撮影はSpring-8（大型放射光施設）にて行った。撮影した画像のサイズは2048*1024ピクセルであり、モアレ縞の一周期は22ピクセルになっている。図１２(a)は従来の手法であるFourier変換法で求めた位相像、図１２(b)は本発明の手法Ａで得られた位相像、そして、図１２(c)は本発明の手法Ｂで得られた位相像である。これらの位相像から明らかなように、本発明の手法Ａ又は手法Ｂで得られた位相像は、従来のFourier変換法で求めた位相像と比較して雑音の影響も小さくなり、明らかに空間分解能が上がっており、更に、手法Ａで得られた位相像と比較して、手法Ｂで得られた位相像のほうが境界と細部が正確に求められていることが分かる。【００５３】
以上に詳述したように、本発明の位相撮影方法では、少なくとも一回の撮影で被検体の位相を高精度で取得することができる。【００５４】
従来の縞走査法と比べると、
Ａ: 回折格子の移動が必要でないため、縞走査法における高い安定性の要求が解消される。
Ｂ: 1回の撮影で位相像が得られることから、装置や設備の操作がより簡単になる。
Ｃ: 1回の撮影で位相像が得られることから、動いている物体や変化している物体への動的測定がより正しくできる。【００５５】
加えて、従来のFourier変換法と比べると、
Ｄ: より高い空間分解能を保つことで、より多くの細かい内部構造が見える画像を提供することが可能となる。
Ｅ: 本発明の手法の精度は撮影される干渉縞の周期に依存するが、Fourier変換法で普段使われている周期より（一周期あたり5、6ピクセル程度）大きくなっても（例えば、数十ピクセルで一周期）精度が大幅に下がらない位相像を求めることができる。回折格子を用いたＸ線位相撮影において、モアレ縞の周期が大きくなるとVisibilityが大幅に上がるので本手法では従来手法より画質の良い画像が得られる。【００５６】
本発明の手法ＢではNon-local means を用いて具体的な重み計算を行っているが、当然、別の手法、例えば、バイラテラルフィルタ(Bilateral filter)などのフィルタ、あるいは、フィルタ以外でも、画像が類似しているほど大きい重みが求められる他の手法など、様々な手法が適用可能である。【００５７】
上記の本発明の手法Ｂでは、近似位相像を用いて重みを計算している。この近似位相像は、本手法Ａ又はFourier変換法の画像に限られることはなく、他の如何なる手法で求めた位相画像でも重み計算に使えるのは当然である。例えば、別の装置やハードウェアで撮影した同じ物体の画像、或いは、前回の測定データ（もしも存在するなら）など、或いは、構造が似ていることが既知の別の物体の画像などが挙げられる。さらに、近似位相像でなくても、物体の不連続な場所の位置情報があれば可能である。一例として、同じ物体の吸収撮影画像を用いることなどが挙げられる。【００５８】
＜その他の本発明が適用可能な設備や装置＞
上述したことから、本発明になる位相撮影方法や位相撮影装置は、特に、その操作の簡便性（回折格子の移動がない）を利用することにより、図１２にも示すように、Talbot干渉計を利用した高速位相ＣＴ撮影装置等において利用することが好適であろう。なお、ここでは、縞走査は行わない。【００５９】
なお、本発明が応用できる技術は上述したものに止まらず、例えば、Ｘ線を照射した際の位相シフト分布の線積分データから位相シフト分布を画像生成する位相Ｘ線ＣＴ、体内に投与した放射性薬剤分布の画像を生成する核医学イメージング装置であるＰＥＴ（ポジトロンエミッションＣＴ）やＳＰＥＣＴ（単光子放射型ＣＴ）、超音波・マイクロ波・音波・地震波などの波動を用いたＣＴ、電子線ＣＴ、投影データからの画像再構成を利用したＭＲＩ（磁気共鳴イメージング）などにも応用可能である。即ち、本発明において「物体」または「画像」とは、画像化する物理量の空間分布のことを指す。【００６０】
また、位相シフト、量子ビームの位相シフト、回折、又は回折を含む画像の数値は、光学素子の追加あるいはその位置変更により取得した複数のセットから抽出され、当該抽出された前記量子ビームの位相シフト、回折、又は回折を含む前記画像の数値を用いて位相値を復元することも可能である。
【産業上の利用性】
【００６１】
本発明は、Ｘ線をはじめとする量子ビームを照射して得られる位相の変化を用いて被検体の画像を得る位相画像撮影方法、更には、それを用いた位相画像撮影装置を提供する。
【符号の説明】
【００６２】
１０…線源、２１〜２４…回折（位相）格子、Ｓ…被検体、１００…画像データ処理部
【図１】
2019240165.tif 000020
【図２】
2019240165.tif 000021
【図３】
2019240165.tif 000022
【図４】
2019240165.tif 000023
【図５】
2019240165.tif 000024
【図６】
2019240165.tif 000025
【図７】
2019240165.tif 000026
【図８】
2019240165.tif 000027
【図９】
2019240165.tif 000028
【図１０】
2019240165.tif 000029
【図１１】
2019240165.tif 000030
【図１２】
2019240165.tif 000031
【図１３】
2019240165.tif 000032
【手続補正書】
【提出日】令和2年11月18日(2020.11.18)
【手続補正１】
【補正対象書類名】明細書
【補正対象項目名】００２８
【補正方法】変更
【補正の内容】
【００２８】
このように、画像の変数近似によって得たデータから最小化問題を解くことで、位相が求められる。具体的には、次のように、上記の式（５）で得られるデータと上記の式（１）の理論値の差の最小化（測定値と理論値との差のLp(p≧0)ノルム最小化）の問題を解くことで求められる。
【数６】
2019240165.tif 000033
ただしIk(k=0,±1)は式(５)と同じ画像の変数近似で得られる複数個のデータである。ここでp=2とすると最小二乗問題になる。
【手続補正２】
【補正対象書類名】明細書
【補正対象項目名】００３２
【補正方法】変更
【補正の内容】
【００３２】
この変換により、問題の式(６)は次のように書き変えることができる。
【数９】
2019240165.tif 000034
【手続補正３】
【補正対象書類名】明細書
【補正対象項目名】００４１
【補正方法】変更
【補正の内容】
【００４１】
この手法でも、発明手法Ａと同様に、最小化問題を解いて位相を求める。上記発明手法Ａとの違いは、ここでは、重み関数W(m,n)を用いて重み付きエネルギー関数の最小化問題を解いている点である。ここで任意の一つのピクセル(m,n)を固定して手法の説明をする。ほかの場所も同じく計算することができる。この重みは、近似位相像を用いて、ピクセル(m,n)において、広い範囲で位相値（又は構造）が近い点を探し、類似している程より大きい重みを与えている。そして、以下の重み付きの最小化問題を解く。
【数１１】
2019240165.tif 000035
【国際調査報告】
2019240165.tif 000036
2019240165.tif 000037
2019240165.tif 000038
2019240165.tif 000039
2019240165.tif 000040
2019240165.tif 000041
2019240165.tif 000042
2019240165.tif 000043