(57)【要約】【課題】微細構造を持つ物体であっても、忠実性の高い３次元画像を取得すること【解決手段】データ処理装置１３０は、第１波長λ１の光から得られる試料１０５の第１ホログラムを表す第１位相分布と、第２波長λ２の光から得られる試料１０５の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行う処理手段１３２と、前記第１位相分布と前記第２位相分布の夫々のうねりが波長毎に異なることに起因して、前記位相アンラッピングによって得られた位相分布において発生する低周波エラー４２１を補正する補正手段１３４と、を有する。【選択図】図１
2019078635.tif 000002
【特許請求の範囲】
【請求項１】
光を物体光と参照光に分割し、物体を透過した前記物体光と、前記参照光とを干渉させることによって生成されたホログラムを表すデータを処理するデータ処理装置であって、
第１波長の光から得られる前記物体の第１ホログラムを表す第１位相分布と、第２波長の光から得られる前記物体の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行う処理手段と、
前記位相アンラッピングによって得られた位相分布における、前記第１位相分布と前記第２位相分布の夫々のうねりが波長毎に異なることに起因して発生するエラー成分を補正する補正手段と、
を有することを特徴とするデータ処理装置。
【請求項２】
前記補正手段は、２πの位相飛びを有する閉領域を前記位相アンラッピングによって得られた位相分布から抽出し、前記閉領域で発生する前記エラー成分を補正することを特徴とする請求項１に記載のデータ処理装置。
【請求項３】
前記補正手段は、第２波長の４倍の２乗よりも小さい面積の閉領域を抽出し、前記閉領域で発生する前記エラー成分を補正することを特徴とする請求項１または２のうちいずれか１項に記載のデータ処理装置。
【請求項４】
前記補正手段は、前記エラーの位相に２πを加算または減算することを特徴とする請求項１乃至３のうちいずれか１項に記載のデータ処理装置。
【請求項５】
前記補正手段は、前記加算または前記減算を繰り返すことを特徴とする請求項４に記載のデータ処理装置。
【請求項６】
光を物体光と参照光に分割し、物体を透過した前記物体光と、前記参照光とを干渉させることによって生成されたホログラムを表すデータを処理するデータ処理装置であって、
第１波長の光から得られる前記物体の第１ホログラムを表す第１位相分布と、第２波長の光から得られる前記物体の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行う処理手段と、
前記第１位相分布と前記第２位相分布の一方の前後の位相が融合することに起因して、前記位相アンラッピングによって得られた位相分布において発生するエラー成分を補正する補正手段と、
を有することを特徴とするデータ処理装置。
【請求項７】
前記補正手段は、ローパスフィルタ演算処理を適用した分布との差分を用いて前記エラーの領域を抽出することを特徴とする請求項６に記載のデータ処理装置。
【請求項８】
前記補正手段は、微分法を用いて前記エラーの領域を抽出することを特徴とする請求項６に記載のデータ処理装置。
【請求項９】
前記補正手段は、抽出された前記エラーの領域を拡大し、拡大された前記領域のエラー成分を補正することを特徴とする請求項７または８に記載のデータ処理装置。
【請求項１０】
前記補正手段は、前記エラーを削除して補間することによって補正することを特徴とする請求項６乃至９のうちいずれか１項に記載のデータ処理装置。
【請求項１１】
前記補正手段は、平面補間を用いることを特徴とする請求項１０に記載のデータ処理装置。
【請求項１２】
前記補正手段によって補正された前記データから前記物体の位相分布、二次元像、屈折率分布、及び三次元像の少なくとも一つを再生する再生部を更に有することを特徴とする請求項１乃至１１のうちいずれか１項に記載のデータ処理装置。
【請求項１３】
光を物体光と参照光に分割し、物体を透過した前記物体光と前記参照光とを干渉させ、ホログラムを生成することによって、前記物体の屈折率分布または位相分布を測定する測定装置であって、
前記ホログラムを光電変換する撮像素子と、
前記撮像素子から得られたデータを処理する請求項１乃至１２のうちいずれか１項に記載のデータ処理装置と、
を有することを特徴とする測定装置。
【請求項１４】
光を物体光と参照光に分割し、物体を透過した前記物体光と、前記参照光とを干渉させることによって生成されたホログラムを表すデータを処理するデータ処理方法であって、
第１波長の光から得られる前記物体の第１ホログラムを表す第１位相分布と、第２波長の光から得られる前記物体の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行うステップと、
前記位相アンラッピングによって得られた位相分布における、前記第１位相分布と前記第２位相分布の夫々のうねりが波長毎に異なることに起因して発生するエラー成分を補正するステップと、
を有することを特徴とするデータ処理方法。
【請求項１５】
光を物体光と参照光に分割し、物体を透過した前記物体光と、前記参照光とを干渉させることによって生成されたホログラムを表すデータを処理するデータ処理方法であって、
第１波長の光から得られる前記物体の第１ホログラムを表す第１位相分布と、第２波長の光から得られる前記物体の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行うステップと、
前記第１位相分布と前記第２位相分布の一方の前後の位相が融合することに起因して、前記位相アンラッピングによって得られた位相分布において発生するエラー成分を補正するステップと、
を有することを特徴とするデータ処理方法。
【請求項１６】
コンピュータに、光を物体光と参照光に分割し、物体を透過した前記物体光と、前記参照光とを干渉させることによって生成されたホログラムを表すデータを処理するデータ処理方法を実行させるためのプログラムであって、
前記データ処理方法は、
第１波長の光から得られる前記物体の第１ホログラムを表す第１位相分布と、第２波長の光から得られる前記物体の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行うステップと、
前記位相アンラッピングによって得られた位相分布における、前記第１位相分布と前記第２位相分布の夫々のうねりが波長毎に異なることに起因して発生するエラー成分を補正するステップと、
を有することを特徴とするプログラム。
【請求項１７】
コンピュータに、光を物体光と参照光に分割し、物体を透過した前記物体光と、前記参照光とを干渉させることによって生成されたホログラムを表すデータを処理するデータ処理方法を実行させるためのプログラムであって、
前記データ処理方法は、
第１波長の光から得られる前記物体の第１ホログラムを表す第１位相分布と、第２波長の光から得られる前記物体の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行うステップと、
前記第１位相分布と前記第２位相分布の一方の前後の位相が融合することに起因して、前記位相アンラッピングによって得られた位相分布において発生するエラー成分を補正するステップと、
を有することを特徴とするプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、測定装置、データ処理装置、データ処理方法およびプログラムに関する。
【背景技術】
【０００２】
非特許文献１は、二波長を用いて試料を３次元撮像する方法を開示している。非特許文献１は、第１波長で取得した試料の位相分布と、第２波長で取得した試料の位相分布の差分値を利用して、試料による位相変化量が２πよりも大きい場合でも位相分布を計測可能にする特許文献１の二波長位相アンラッピング法を、３次元撮像に適用している。非特許文献２は、Ｔｅｎｅｎｂａｕｍ微分法を開示している。
【先行技術文献】
【特許文献】
【０００３】
【特許文献１】米国特許第６８０９８４５号
【非特許文献】
【０００４】
【非特許文献１】Ｍ．  Ｒｉｎｅｈａｒｔ，  Ｍ．  Ｇｉａｃｏｍｅｌｌｉ，  Ｋ．  Ｃｈａｌｕｔ，  ａｎｄ  Ａ．  Ｗａｘ，  ”Ｑｕａｎｔｉｔａｔｉｖｅ  Ｐｈａｓｅ  Ｍｉｃｒｏｓｃｏｐｙ  ｗｉｔｈ  Ｍｕｌｔｉ−Ｗａｖｅｌｅｎｇｔｈ  Ｕｎｗｒａｐｐｉｎｇ  ａｎｄ  Ｔｏｍｏｇｒａｐｈｉｃ  ３Ｄ  Ｒｅｃｏｎｓｔｒｕｃｔｉｏｎ”  ＯＳＡ  Ｄｉｇｉｔａｌ  Ｈｏｌｏｇｒａｐｈｙ  ａｎｄ  Ｔｈｒｅｅ−Ｄｉｍｅｎｓｉｏｎａｌ  Ｉｍａｇｉｎｇ  （ＤＨ）  ２００８  ｐａｐｅｒ：ＪＭＡ７
【非特許文献２】Ｙ．Ｓｕｎ，Ｓ．Ｄｕｔｈａｌｅｒ，ａｎｄ Ｂ．Ｊ．Ｎｅｌｓｏｎ，”Ａｕｔｏｆｏｃｕｓｉｎｇ ａｌｇｏｒｉｔｈｍ ｓｅｌｅｃｔｉｏｎ ｉｎ ｃｏｍｐｕｔｅｒ ｍｉｃｒｏｓｃｏｐｙ，”ｉｎ Ｐｒｏｃ．ＩＥＥＥ／ＲＳＪＩｎｔｅｒｎａｔｉｏｎａｌ Ｃｏｎｆｅｒｅｎｃｅｏｎ ＩｎｔｅｌｌｉｇｅｎｔＲｏｂｏｔｓ ａｎｄ Ｓｙｓｔｅｍｓ，ｐｐ．７０−７６（２００５）
【発明の概要】
【発明が解決しようとする課題】
【０００５】
しかしながら、微細構造を持つ試料では、回折の影響が各波長で異なるため、特許文献１に開示されている方法を用いると、位相分布に複雑なアンラッピングエラーが生じてしまう。そして、この位相分布から非特許文献１を用いて３次元画像を取得すると、忠実性の高い再生像は得られない。
【０００６】
本発明は、微細構造を持つ物体であっても、忠実性の高い３次元画像を取得することが可能な測定装置、データ処理装置、データ処理方法およびプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００７】
本発明の一側面としてのデータ処理装置は、光を物体光と参照光に分割し、物体を透過した前記物体光と、前記参照光とを干渉させることによって生成されたホログラムを表すデータを処理するデータ処理装置であって、第１波長の光から得られる前記物体の第１ホログラムを表す第１位相分布と、第２波長の光から得られる前記物体の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行う処理手段と、前記位相アンラッピングによって得られた位相分布における、前記第１位相分布と前記第２位相分布の夫々のうねりが波長毎に異なることに起因して発生するエラー成分を補正する補正手段と、を有することを特徴とする。
【０００８】
また、本発明の別の側面としてのデータ処理装置は、光を物体光と参照光に分割し、物体を透過した前記物体光と、前記参照光とを干渉させることによって生成されたホログラムを表すデータを処理するデータ処理装置であって、第１波長の光から得られる前記物体の第１ホログラムを表す第１位相分布と、第２波長の光から得られる前記物体の第２ホログラムを表す第２位相分布と、の差分を用いた位相アンラッピングを行う処理手段と、前記第１位相分布と前記第２位相分布の一方の前後の位相が融合することに起因して、前記位相アンラッピングによって得られた位相分布において発生するエラー成分を補正する補正手段と、を有することを特徴とする。
【０００９】
上記データ処理装置を有する測定装置やデータ処理装置によって実行されるデータ処理方法等も本発明の別の側面を構成する。
【発明の効果】
【００１０】
本発明によれば、微細構造を持つ物体であっても、忠実性の高い３次元画像を取得することが可能な測定装置、データ処理装置、データ処理方法およびプログラムを提供することができる。
【図面の簡単な説明】
【００１１】
【図１】本発明の測定装置の構成図である。（実施例１、２）
【図２】方位角を説明するための図である。（実施例１、２）
【図３】３次元屈折率分布を求めるフローチャートである。（実施例１、２）
【図４】従来手法で求めた位相分布である。
【図５】図４（ｂ）の部分拡大図である。
【図６】実施例１の閉空間の探索方法を説明するための図である。
【図７】実施例１の低周波エラーを補正した位相分布である。
【図８】実施例１の高周波エラーを補正した位相分布である。
【図９】実施例１と比較例において生成された３次元屈折率分布の再生像である。
【図１０】実施例２の高周波エラーを補正した位相分布である。
【図１１】実施例２において生成された３次元屈折率分布の再生像である。
【発明を実施するための形態】
【００１２】
以下、添付図面を参照して、本発明の好適な実施の形態について説明する。なお、各図において、同一の部材については同一の参照番号を付し、重複する説明は省略する。
【００１３】
まず、本発明を実施するための測定装置の構成について説明する。測定装置および測定装置によって実行される測定方法は、試料の屈折率分布（三次元）や位相分布（二次元）の測定に使用することができ、試料の二次元形状または三次元形状を算出及び再生（再構成）する用途に適用することができる。
【００１４】
また、本発明の測定装置は、デジタルホログラムを生成するデジタルホログラム生成装置としても機能する。デジタルホログラムは、試料１０５の透過光と参照光による干渉縞を撮像素子１０８でデジタル化し、コンピュータ（画像処理手段、データ処理手段）１３０による後処理で透過光の電場を算出する際に使用される。
【００１５】
デジタルホログラフィを用いると、電場の位相変化を定量的に求めることができるため、無色透明な試料の構造を観察することが可能となる。医学、生物学などの分野において、生体試料を観察することは病気の診断や生体構造を理解する上で欠くことのできない。生体試料は通常、無色透明あるいは半透明であるため、像のコントラストを上げるために染色を行うと、試料を生きたまま観察することができない場合がある。デジタルホログラフィは、顕微鏡に干渉計測を応用することで非染色の試料を観察する方法の一つである。
【００１６】
病気の診断をする際に、医師は細胞の核などの形状を根拠の１つとしている。そのため試料の３次元情報を取得することが必要とされる。そこで、試料を異なる入射角で照明し、得られた複数枚の透過位相分布に対して再構成処理を行うことで試料の３次元的な屈折率分布を取得する方法が提案されている。しかしながら、大きい位相差を発生させる試料に適用した場合、適切な値が得られないことがある。本発明はこの問題を解決することができる。
【００１７】
図１は、測定装置の構成図である。本実施例の測定装置は、光を物体光と参照光に分割し、試料を透過した物体光と参照光とを干渉させ、干渉縞（ホログラム）を生成することによって、試料の屈折率分布または位相分布を測定する。
【００１８】
図１に示す測定装置は、軸外し型ホログラムを使用する。１０１は第１波長の光を射出する光源（第１光源）であり、レーザーなどのコヒーレンシーが高いものを用いる。１０２は一部の光を透過し、一部の光を反射するハーフミラーである。１１１は試料（物体）を透過する光である物体光を、１１２は試料を透過せずに物体光１１１と干渉させる光である参照光を示している。物体光１１１と参照光１１２は、例えば、ハーフミラー１０２により分割された光源１０１からの光である。
【００１９】
１０３は２軸のガルバノミラーを、１０４はコンデンサーレンズを、１０５は試料を、１０６は対物レンズを、それぞれ示している。１０７は物体光１１１と参照光１１２を統合するハーフミラーである。
【００２０】
１０８は撮像素子（またはデジタルセンサー）である。ここでは参照光１１２を撮像素子１０８に垂直に入射させ、物体光１１１を参照光１１２に対して傾斜させる。撮像素子１０８に対する入射角が、物体光１１１と参照光１１２で異なるため、撮像素子１０８では、干渉縞が形成される。このようにして得られた干渉縞は、物体光１１１の電場の情報を含んでいる。撮像素子１０８は干渉縞を光電変換し、不図示のＡ／Ｄ変換器によってデジタル化されたデータが生成される。
【００２１】
試料１０５に対する物体光１１１の角度（方位角）は、試料１０５の共役面に配置されたガルバノミラー１０３の角度で変化させることが可能である。物体光１１１の角度を変化させ、複数方向からの物体光１１１の電場を取得することで試料１０５の３次元像を再生することが可能となる。ここで、物体光１１１と光軸ＯＡのなす角度をθ、方位角をφと定義する。
【００２２】
図２（ａ）、（ｂ）は、方位角を説明するための図である。方位角φは、光軸ＯＡに垂直な平面に投射された物体光１１１と前記平面上の光軸ＯＡから所定方向に延びる直線とのなす角度である。光軸方向をＺ軸方向、光軸方向に垂直な平面（試料面）をＸＹ平面とし、光軸ＯＡから所定方向に延びる直線をＸ軸する。
【００２３】
本実施例では、θを一定とし、φを変化させて、試料１０５に対する物体光１１１の角度を変化させる。具体的にはθ＝５０度、φ＝０，３，６，・・・，３５７度とする。各φに対して撮像を行うため、撮像回数は１２０回となる。θが一定のため、物体光１１１が参照光１１２と成す角度が一定となり、得られる干渉縞のピッチは常に一定となる。撮像素子１０８上の干渉縞のピッチが、撮像素子１０８の画素サイズに対して不適切な場合、撮像素子１０８の前に倍率を変化させる光学系を挿入してもよい。
【００２４】
１２１は第１波長とは異なる第２波長の光を射出する光源（第２光源）であり、光源１０１と同様に、レーザーなどのコヒーレンシーが高いものを用いる。１２２はハーフミラーであり、第１波長の光路と第２波長の光路を共通化する役割を果たす。通常、第１波長の光源１０１か第２波長の光源１２１のどちらか一方を遮光して、各波長のホログラムを撮像する。ただし、フーリエ空間にて第１波長の成分と第２波長の成分が分離できる場合は、２種類の波長のレーザーを同時に照射してもよい。
【００２５】
第１波長、第２波長はいずれも限定されない。試料１０５の厚さがｄである場合、後述する合成波長Λは、ｄΔｎを少し上回る程度、例えば、ｄΔｎ×２程度が好ましい。例えば、ｄΔｎ＜Λ＜ｄΔｎ×３であってもよい。ここで、Δｎは試料１０５の屈折率差、例えば、試料１０５の厚さが５μｍで、屈折率差が０．２の場合、合成波長Λは、Λ＝５×０．２×２＝２［μｍ］となる。第２波長λ２を０．５μｍとすると、合成波長Λは、後述する式（２）で与えられるので、第１波長λ１は、次式に示すように、０．６６７μｍ程度となることが好ましい。
【００２６】
λ１＝Λ・λ２／（Λ−λ２）＝１／１．５＝０．６６７［μｍ］
１３０は、データ処理装置（画像処理装置（画像処理手段））であり、図２に示すデータ処理方法（画像処理方法）を実行する。本実施例では、データ処理装置１３０は、パーソナルコンピュータ（プロセッサ）として構成されている。なお、光源１０１、１２１を不図示のシャッターを駆動することにより遮光、導光したり、電源のオンオフしたりする不図示の光源駆動手段、ガルバノミラー１０３の角度を駆動する不図示のミラー駆動手段が設けられている。光源駆動手段、ミラー駆動手段は、不図示の制御手段（コントローラ）に接続され、制御手段が、光源駆動手段、ミラー駆動手段による駆動を制御している。この制御手段とデータ処理装置１３０が一体に構成され、データ処理装置１３０が制御手段の機能を有してもよい。
【００２７】
第１波長の光および撮像素子１０８から、試料１０５の第１ホログラムを表す第１位相分布が生成され、データ処理装置１３０に入力される。第２波長の光および撮像素子１０８から、試料１０５の第２ホログラムを表す第２位相分布が生成され、データ処理装置１３０に入力される。
【００２８】
データ処理装置（データ処理手段）１３０は、記憶手段１３８を有し、表示部（表示手段）１４０とネットワーク１５０に接続されている。
【００２９】
データ処理装置１３０は、位相アンラッピング手段１３２、補正手段１３４および再生部（再構成部）１３６の機能を有する。
【００３０】
位相アンラッピング手段１３２は、第１位相分布と第２位相分布の差分と次数を用いた二波長位相アンラッピングを行う。位相アンラッピング（位相接続）は、２π周期で折り畳まれている試料（物体）の位相分布を実際の位相分布に展開することをいう。
【００３１】
補正手段１３４は、第１補正手段と第２補正手段の機能を有しているが、いずれか一方の機能のみを有していても本発明を構成する。第１補正手段は、第１位相分布と第２位相分布の夫々のうねりが波長毎に異なることに起因して、位相アンラッピングによって得られた位相分布において発生するエラー成分（後述する低周波エラー）を補正する。第２補正手段は、第１位相分布と第２位相分布の一方の前後の位相が融合することに起因して、位相アンラッピングによって得られた位相分布において発生するエラー成分（後述する高周波エラー）を補正する。なお、補正手段１３４は低周波エラーと高周波エラーを同時に補正してもよい。
【００３２】
再生部１３６は、補正手段１３４によって補正されたデータから試料（物体）１０５の位相分布、二次元像、屈折率分布、及び三次元像の少なくとも一つを再生する。
【００３３】
記憶手段１３８は、撮像素子１０８によって光電変換され、ホログラムを表すデータ、データ処理方法（画像処理方法）のプログラム（データ処理プログラム）及びそれに使用される情報、データ処理方法によって処理されたデータを記憶する。
【００３４】
表示手段１４０は、試料１０５の屈折率分布、位相分布、それらの断面図、三次元再生像、二次元再生像、各種の制御情報等を表示する。ネットワーク１５０は、ＬＡＮ、インターネット等の有線または無線のネットワークであり、ネットワーク１５０には、不図示のコンピュータ、出力装置（プリンタ、ディスプレイ）、記憶装置等が接続されている。このため、データ処理装置１３０は、ネットワーク１５０を介して測定装置に接続されたコンピュータであってもよいし、測定装置に接続されていないコンピュータであってもよい。
【００３５】
図３は、データ処理装置１３０が実行する３次元再生像計算方法（データ処理方法）を示すフローチャートであり、Ｓはステップを表している。図３に示すフローチャートは、コンピュータに各ステップを実行させるためのプログラムとして実現が可能であり、かかるプログラムは、例えば、コンピュータ可読媒体に記憶可能である。
【００３６】
まず、Ｓ３０１で、ガルバノミラー１０３の角度を不図示のミラー駆動手段を制御することによって方位角φを指定する。ここでは第１方位角として、０度を指定する。Ｓ３０２で、第１波長のホログラム（第１ホログラム）のデータを読み込む。Ｓ３０３で、物体光１１１の電場の算出を行う。物体光１１１の電場の算出方法は種々存在するが、ここではフーリエ変換を用いた方法を用いる。具体的には、ホログラムをフーリエ変換してスペクトルを求め、物体情報を含む＋１次のスペクトルを計算機上の絞りで抽出し、抽出したスペクトルを座標の原点にシフトし、逆フーリエ変換で物体光の電場を算出する。３次元屈折率の分布の算出方法も複数存在するが、ここではＲｙｔｏｖ近似に基づいた逆伝搬法を用いた手法を使用する。Ｒｙｔｏｖ近似を用いるには、物体光１１１の電場の対数を取り、位相分布を算出する。そこで、Ｓ３０４で物体光１１１の電場の対数を取る。この対数の虚部が位相分布に相当する。第２波長に対しても同等の処理をＳ３０５〜３０７で行う。ここでは第１波長の処理の後に第２波長の処理を行ったが、この順序は逆でもよい。
【００３７】
次に、位相アンラッピング手段１３２が行う二波長位相アンラッピング法（Ｓ３１１〜Ｓ３１４）について説明する。忠実性の高い３次元屈折率分布を求めるためには、物体光１１１の電場の位相変化量が２πを超える場合、位相アンラッピングをする必要がある。隣り合う画素でπ以上変化する箇所がない位相分布では、単一波長で得られた位相分布からアンラッピングをすることが可能である。これに対し、隣り合う画素でπ以上変化する箇所がある位相分布では、単一波長の位相分布からアンラッピングをすることは難しい。この問題を克服するアンラッピング手法として、特許文献１に記載の二波長位相アンラッピング法が存在する。
【００３８】
ここで、二波長位相アンラッピング法の具体的な処理を説明する。まず、Ｓ３１１で、位相の差分ψ１（ｘ，ｙ）−ψ２（ｘ，ｙ）求める。ψ１（ｘ，ｙ）は第１位相分布、ψ２（ｘ，ｙ）は第２位相分布である。次に、Ｓ３１２で、第１位相分布と第２位相分布の差分に、合成波長Λを第１波長λ１で割った係数Λ／λ１を掛けた以下の式で位相分布ψ（ｘ，ｙ）を算出する。
【００３９】
【数１】
2019078635.tif 000003
【００４０】
なお合成波長はΛ以下の式で求められる。
【００４１】
【数２】
2019078635.tif 000004
【００４２】
波長がλ１からΛに増加することに伴って、アンラッピングできる範囲もΛ／λ１倍拡大される。Λを用いずに表現するとλ２／｜λ１−λ２｜倍拡大される。但し、拡大倍率に比例してノイズの量も増大される。そこで、以下の２つのＳによりノイズの増大を防ぐ。先ずＳ３１３で、以下の式によりアンラッピングの次数分布ｎ（ｘ，ｙ）を求める。次数ｎとは、‐π〜πに折り畳まれた位相を元に戻す（位相アンラップ）するために２πに掛けられる係数（整数）である。例えば、‐２π〜２πの位相分布を計測した場合、その値は−π〜πに折り畳まれている。この場合、−２π〜−πの領域では、ｎ＝−１、−π〜πの領域では、ｎ＝０、π〜２πの領域では、ｎ＝１となる。
【００４３】
【数３】
2019078635.tif 000005
【００４４】
続いて、Ｓ３１４で、求めた次数を以下の式に代入し、ノイズ量が増加していない位相分布ψ’（ｘ，ｙ）を得る。ψ’（ｘ，ｙ）が二波長位相アンラッピング法で求めた位相分布である。
【００４５】
【数４】
2019078635.tif 000006
【００４６】
図４に、屈折率１．３４の液体に浸した屈折率１．４６、直径５μｍの球を試料１０５とした場合の、物体光１１１の位相分布の計算例を示す。球の中心を（ｘ，ｙ，ｚ）＝（０，０，０）とし、物体光１１１の角度を（θ，φ）＝（５０，０）とした。第１波長λ１は０．５４３μｍ、第２波長λ２は０．５９３μｍである。計算に使用した画素サイズは試料面換算で０．０５３４μｍ、画素数は１０２４×１０２４、対物レンズのＮＡは１．４、計算機絞りのＮＡは０．７である。ｚ軸は光軸、ｘｙは光軸に垂直な平面（試料面）である。
【００４７】
図４（ａ）はＳ３１１で求めた位相の差分であり、横軸はｘ軸、左側の縦軸はｙ軸、右側の縦軸は位相を表している。図４（ｂ）は図４（ａ）のｙ＝０での断面図であり、横軸はｘ軸、縦軸は位相、λ１は第１波長、λ２は第２波長、diffはそれらの位相の差分を表している。
【００４８】
図４（ｃ）はＳ３１２で求めた合成波長Λの位相分布であり、横軸はｘ軸、左側の縦軸はｙ軸、右側の縦軸は位相を表している。図４（ｄ）は図４（ｃ）のｙ＝０での断面図であり、横軸はｘ軸、縦軸は位相を表している。合成波長は式（２）で求められ、６．４４μｍとなる。
【００４９】
図４（ｅ）はＳ３１４で得られた位相分布であり、横軸はｘ軸、左側の縦軸はｙ軸、右側の縦軸は位相を表している。図４（ｆ）は図４（ｅ）のｙ＝０での断面図であり、横軸はｘ軸、縦軸は位相を表している。図（４）（ｆ）では、図４（ｄ）に示す位相分布がある程度滑らかにされていることが分かる。
【００５０】
しかしながら、球の形状に対応した滑らかな位相分布が期待されるのに対し、実際に算出された図４（ｆ）の位相分布には、例えば、４２１や４２２にアンラッピングエラーが発生している。４２１を低周波エラー、４２２を高周波エラーと呼ぶことにする。特許文献１で発生する１〜数画素の範囲の２π飛びエラーであれば、ローパスフィルタ演算処理などで補正することが可能である。しかしながら、微細な構造を持つ試料１０５を測定した場合、波長毎の位相分布に異なる回折波が生じ、複雑なエラーが発生するため、ローパスフィルタ演算処理などの単純は処理でエラー成分を補正することは難しい。複雑なエラーの例としては、数画素を超える範囲に及ぶエラーや、特異点を含むエラーが挙げられる。
【００５１】
ここで、低周波エラーおよび高周波エラーについて説明する。低周波エラーおよび高周波エラーは、そもそも位相の差分の分布を算出した時点で、図４（ｂ）に示すように、低周波エラー４０１、高周波エラー４０２として発生している。これらが、図４（ｄ）に示すように、低周波エラー４１１、高周波エラー４１２として残り、最終的に図４（ｆ）に示すように、低周波エラー４２１、高周波エラー４２２となっている。そこで、エラーが最初に発生した位相の差分を用いて、各エラーについて説明する。
【００５２】
まず、低周波エラーの原因について説明する。図５（ａ）は、図４（ｂ）の低周波エラー４０１が発生している領域の拡大図である。第１波長λ１の位相分布５０１および第２波長λ２の位相分布４０２において、回折の効果によりうねりが生じている。このうねりは波長ごとに異なる。一方の波長のうねりの山と、他方の波長のうねりの谷が重なる領域では、差分にもうねりが生じる。この差分のうねりが大きくなると、Ｓ３１３の次数算出時に、エラーが発生する。エラーが±１の場合、Ｓ３１４で±２πのエラーとなる。このように、低周波エラーは、第１位相分布と第２位相分布の夫々のうねりが波長毎に異なることに起因して、二波長位相アンラッピングによって得られた位相分布（の閉領域）において発生するエラー成分である。
【００５３】
次に、高周波エラーの原因について説明する。図５（ｂ）は、図４（ｂ）の高周波エラー４０２が発生している領域の拡大図である。点線５１１は第１波長λ１の位相分布を示している。一点鎖線４１２は第２波長λ２の位相分布を示している。図５（ｃ）は、図５（ｂ）をアンラッピングした位相分布である。この位相分布は位相の遅れを表しているため、相対的に、球の中心に近いｘ＝−３．５μｍ付近では大きい値、球の中心から遠いｘ＝−４．５μｍ付近では小さい値となることが期待される。よって、第２波長λ２の位相分布が本来の分布である。これに対し、第１波長λ１の位相分布は原点に近い方が小さい値となっている。これは領域５２２において位相が急激に変化し、前後の位相が融合したためである。具体的には、領域５２１の位相（前の位相）５３１と、領域５２３の位相−２πの位相（後の位相）５３２とが領域５２２において融合したためである。図５（ｃ）のように、一方の波長のみで前後の位相の融合が起こると、位相の差分に大きなエラーが発生する。このエラーは領域５２２のように非常に狭い領域で発生し、高周波エラーとなる。このように、高周波エラーは、第１位相分布と第２位相分布の一方の前後の位相が融合することに起因して、位相アンラッピングによって得られた位相分布において発生する（傾きが閾値以上の）エラー成分である。高周波エラーの周波数は低周波エラーの周波数よりも高い。
【００５４】
次に、補正手段１３４が行う低周波エラーおよび高周波エラーの補正について説明する。この処理は、低周波エラー補正の処理（Ｓ３２１〜Ｓ３２２）と高周波エラー補正の処理（Ｓ３２３〜Ｓ３２５）で構成される。
【００５５】
まず、Ｓ３２１で、低周波エラーの領域を抽出する。低周波エラーは、２種類の位相分布のうねりが異なるために発生する。よって、２π飛びに特異点は無く、閉領域で発生する。そこで、２π飛びをしている領域のうち閉領域を探索する。そして、探索した閉領域のうち面積が小さいものを低周波エラーとする。面積が小さいものに限定するのは、試料１０５が発生させるπ以上の位相差を誤って補正しないためである。このπ以上の位相差とは、λ１あるいはλ２の単一の波長でも発生する位相差である。
【００５６】
図６は、閉空間の探索方法を説明するための図である。図６（ａ）は、図４（ｃ）に対応する位相分布（１）であり、ｘｙ平面上の３２グリッド（画素）×３２グリッド（画素）の大きさの領域を表している。データ処理装置１３０は、これを上下左右に１グリッドずつシフトさせて和を取って図６（ｂ）に示す位相分布（２）を取得する。但し、シフト量は例示である。次に、図６（ｃ）に示すように、データ処理装置１３０は、（２）−（１）×４によって輪郭を抽出する。４倍しているのは図６（ｂ）において、４回位相分布を足し合わせているからである。次に、図６（ｄ）に示すように、（３）＞πでπ飛びの領域を囲い、ライブラリーを用いて閉領域を抽出する。πよりも大きい領域は２π位相飛びの領域を表している。非閉領域の典型は渦領域（Ｖｏｒｔｅｘ）である。２π飛びはあるが、途切れているため閉領域とはならず、渦領域は高周波エラーとして補正される。
【００５７】
次に、抽出された閉領域を所定の閾値と比較して、所定の面積よりも小さい閉領域を抽出する。閉領域として決定する領域の面積が大き過ぎると測定したい試料１０５を消してしまうため、面積の上限を規定する。面積の直径は、測定したい試料１０５の直径の半分以下程度とすればよい。例えば、５μｍの試料１０５を測定するのであれば数μｍとすればよい。具体的には、直径を２μｍ程度、面積を２２［μｍ２］程度とすればよい。波長が長いほど測定できる位相差も大きくなるため、規定の面積も大きくしてもよい。そこで、面積を長さの絶対量で決めるのではなく、波長の関数で表してもよい。例えば、第２波長λ２を０．５μｍとすると、面積は（４λ２）２となるので、閉領域を、（４×第２波長λ２）（第２波長λ２の４倍）の２乗以下の面積に設定してもよい。（４×第１波長λ１）の２乗以下の面積に設定してもよい。
【００５８】
次に、Ｓ３２２で、低周波エラーを補正するべく、小領域の位相アンラッピングを行う。Ｓ３２１で抽出した低周波エラー領域に対して、位相飛びが発生した領域の位相が周りの位相よりも小さい場合は２πを加え、大きい場合は２πを差し引く。この処理で低周波エラーが補正される。なお、２πを加える、または２πを差し引く処理で十分な補正が得られない場合は、２πを加算または減算する処理を繰り返し行ってもよい（即ち、２π×Ｎ（Ｎは整数）を加算または減算してもよい）。
【００５９】
Ｓ３２３で、高周波エラーがある領域を抽出する。高周波エラーの領域を抽出する方法は種々あるが、ここではローパスフィルタ演算処理を用いて得た画像と比較して差分が閾値を超える領域を高周波エラー領域とする。ローパスフィルタの演算は、位相分布をフーリエ変換し、高周波成分を指定したＮＡの絞りで削除して、逆フーリエ変換することにより行う。Ｓ３２３で得られた領域が高周波エラー領域をカバーできていない場合に備え、Ｓ３２４で領域を数画素分拡大する。Ｓ３２５で、抽出した高周波エラー領域の位相を取り除き、その領域の位相を補間で求める。補間方法は複数存在するが（スプライン補間等）、ここでは平面補間を用いる。
【００６０】
Ｓ３２３の高周波エラーの領域の別の抽出方法として、Ｔｅｎｅｎｂａｕｍ微分法などの微分法を用いてもよい。Ｔｅｎｅｎｂａｕｍ微分法は、非特許文献２に記載され、ある程度の広がりを持った領域毎の微分値を、以下の式で算出する。
【００６１】
【数５】
2019078635.tif 000007
【００６２】
＊はコンボリューションを表わす。ψは位相分布である。この手法で位相分布の微分値を求め、微分値が閾値を超える領域を、高周波エラー領域としてもよい。
【００６３】
最後に、再生部１３６が行う３次元像の算出について説明する（Ｓ３３１〜Ｓ３３３）。Ｓ３３１で、Ｓ３２５で求めた位相分布を虚部とした物体光の電場の対数をフーリエ変換してスペクトルを算出する。次に、Ｓ３３２で、試料１０５への入射角を考慮してスペクトルを、周波数空間において半径が１／λ、中心位置が（−ｆｘ０、−ｆｙ０、−ｆｚ０）の球（Ｅｗａｌｄ球）面上に貼り付ける。ここで、ｆｘ０、ｆｙ０、ｆｚ０は、物体光のベクトルをフーリエ空間で表した際の周波数である。上述の処理（Ｓ３０２〜Ｓ３３２）を物体光１１１の全方位角に対して繰り返し行い（Ｓ３３３、Ｓ３３５）、全てのＥｗａｌｄ球の平均値を取ることで３次元的なスペクトルが得られる。最後に、Ｓ３３４で、得られた３次元のスペクトルを逆フーリエ変換して３次元屈折率の分布を算出する。
【００６４】
ここでは低周波エラーの補正の後に高周波エラーの補正を行ったが、順序はこれに限ったものではなく、高周波エラーの補正の後に低周波エラーの補正を行ってもよい。また位相アンラップは電場の対数の虚部に適用したが、これに限定するものではなく、いかなる逆伝搬法のアンラッピング処理部に適用してもよい。また軸外し型のホログラムを使用する構成で説明したが、位相シフト型のホログラムを使用する構成でもよい。
【実施例１】
【００６５】
実施例１は、屈折率１．３４の液体に浸した屈折率１．４６、直径５μｍの球を試料１０５とした例について説明する。物体光１１１の角度（θ，φ）＝（５０，０）に対して、特許文献１に記載の二波長位相アンラッピング法で求めた位相分布は、図４（ｅ）、（ｆ）に示した通りである。
【００６６】
図７（ａ）は、Ｓ３２１で抽出した２π飛びが発生している小領域（閉領域）である。この小領域は、図４（ｅ）の低周波エラーの領域に対応している。ここでは、面積が波長の４倍を２乗したものよりも小さい領域のみ抽出した。面積の代わりに幅を用い、領域の幅が指定した幅より狭い領域のみを抽出してもよい。図７（ｂ）は、Ｓ３２２で抽出した領域を位相アンラッピングした位相分布である。図７（ｃ）は、図７（ｂ）のｙ＝０での断面図である。図４（ｂ）に対して、低周波エラー４０１が正しく補正されている。
【００６７】
図８（ａ）は、Ｓ３２３で抽出した高周波エラーの領域である。この領域は、図４（ｅ）の高周波エラーの領域に対応している。領域抽出用のローパスフィルタ演算処理のＮＡは０．５、ローパスフィルタを用いて得た画像と比較して差分が閾値３／４πを超える領域を高周波エラー領域とした。図８（ｂ）は、Ｓ３２４で図８（ａ）を拡大した領域である。ここでは３画素分拡大した。図８（ｃ）はＳ３２５で拡大した領域の位相分布を平面補間により求めた位相分布である。図８（ｄ）は図８（ｃ）のｙ＝０での断面図である。図７（ｃ）に対して、高周波エラーが正しく補正されている。
【００６８】
この様にして求めた位相分布を、Ｓ３３１にてスペクトルに変換し、Ｓ３３２で３次元化する。物体光の角度（θ，φ）＝（５０，３），（５０，３），・・・，（５０，３５７）に対しても同様の処理を行い、Ｓ３３３にて３次元逆フーリエ変換で屈折率分布を算出する。２次元の計算では画素数は１０２４×１０２４としたが、メモリの制約上３次元スペクトルの領域の画素数は２５６×２５６×２５６とした。その結果、実空間での画素サイズは、０．０５３４μｍから４倍の０．２１４μｍに拡大されている。
【００６９】
図９（ａ）が実施例１で得られた３次元屈折率分布のｙ＝０でのｘｚ断面である。なお光軸はｚ軸である。図９（ｂ）は図９（ａ）のｚ＝０での断面図であり、図９（ｃ）は図９（ａ）のｘ＝０での断面図である。図９（ｄ）は、比較のために本発明の技術を用いずに特許文献１の技術で求めた３次元屈折率分布のｙ＝０でのｘｚ断面である。図９（ｅ）は図７（ｄ）のｚ＝０での断面図であり、図９（ｆ）は図９（ｄ）のｘ＝０での断面図である。図９（ａ）と図９（ｄ）を比較すると、本発明の技術を用いて得られた再生像である図９（ａ）の方が、球の形状をより忠実に再現している。また図９（ｆ）で発生している振動成分が、本発明を用いて得られた再生像である図９（ｃ）では低減されている。
【００７０】
本実施例によれば、アンラッピングエラーが低減され、大きな回折を発生させる微細構造を持つ試料に対しても忠実性の高い３次元再生像が得られる。
【実施例２】
【００７１】
実施例２では、高周波エラーがある領域を微分法で求める。試料および計算条件は実施例１と同じである。実施例２の説明において、実施例１と同じ内容については、説明を省略する。図１０（ａ）は式（２）による微分法を用いて求めた、高周波エラーの領域である。図１０（ｂ）は、Ｓ３２４で図１０（ａ）を拡大した領域である。図１０（ｃ）は拡大した領域の位相分布を平面補間により求めた位相分布である。図１０（ｄ）は図１０（ｃ）のｙ＝０での断面図である。図７（ｃ）に対して、高周波エラーが正しく補正されている。
【００７２】
以下、実施例１と同様に、アンラッピングエラーで発生した乱れが補正された再生像を求める。図１１（ａ）が得られた３次元屈折率分布のｙ＝０でのｘｚ断面である。図１１（ｂ）は図１１（ａ）のｚ＝０での断面図であり、図１１（ｃ）は図１１（ａ）のｘ＝０での断面図である。図１１（ａ）と図９（ｄ）を比較すると、本発明の技術を用いて得られた再生像である図１１（ａ）の方が、球の形状をより忠実に再現している。また、図９（ｆ）で発生している振動成分が、本発明を用いて得られた再生像である図１１（ｃ）では低減されている。
【００７３】
本実施例によれば、アンラッピングエラーが低減され、大きな回折を発生させる微細な構造を持つ試料に対しても忠実性の高い３次元再生像が得られる。
【００７４】
以上、本発明の好ましい実施例について説明したが、本発明はこれらの実施例に限定されず、その要旨の範囲内で種々の変形及び変更が可能である。例えば、高周波エラーを補正する際に補正対象とする周波数（または傾き）を低周波エラーを含むように調節すれば、高周波エラーのみならず低周波エラーも取り除くことができる。この場合は、ステップＳ３２１，Ｓ３２２は不要となる。
【符号の説明】
【００７５】
１０５…試料、１０８…撮像素子、１１１…物体光、１１２…参照光、１３０…データ処理装置（データ処理手段）、４２１…低周波エラー、４２２…高周波エラー、λ１…第１波長、λ２…第２波長
【図１】
2019078635.tif 000008
【図２】
2019078635.tif 000009
【図３】
2019078635.tif 000010
【図４】
2019078635.tif 000011
【図５】
2019078635.tif 000012
【図６】
2019078635.tif 000013
【図７】
2019078635.tif 000014
【図８】
2019078635.tif 000015
【図９】
2019078635.tif 000016
【図１０】
2019078635.tif 000017
【図１１】
2019078635.tif 000018