JP4538661B2 - 光学結像システムにおける深度弁別改良のための方法 - Google Patents

Description

本発明は光学結像システムにおける深度弁別方法に関する。この方法は、特に光学顕微鏡検査において３次元に広がる対象物の画像品質を高めるために適用することができる。

旧来型の光学顕微鏡検査では、３次元に広がる対象物の場合、すなわち、光軸に沿う方向の広がりが使用対物レンズの焦点深度より大きい対象物の場合には、フォーカス領域からの鮮明な画像の上に、焦点の外れた、したがって不鮮明な結像状態にある画像部分が覆い被さるという問題がある。
この問題の克服には、いわゆるピンホールの使用により焦点外の領域からの光を遮断して結像に関わらせないようにする共焦点結像方法が利用されており公知である。このようにして、いわゆる光学スライスが生成される。この共焦点ドット方式の結像では、画像を得るためには対象物を画像平面内で走査する必要がある。この走査はレーザ走査型顕微鏡のスキャナにより、あるいはニポーディスクにより行なうことができる。
様々なフォーカスポジションにおける複数の光学スライス像の撮影により、対象物の３次元描画を可能にする「ｚスタック」を得ることができる。

光学スライスを生成する別な方法として構造化照明の適用がある。これは、Ｍｅｉｅｒ Ｂｅｎ−ＬｅｖｙおよびＥｙａｌ ＰｅｌｅｃによりＷＯ ９７／６５０９において初めて公表された。この方法の改良および拡充については、Ｔｏｎｙ Ｗｉｌｓｏｎ他がＷＯ ９８／４５７４５に、Ｖｏｌｋｅｒ Ｇｅｒｓｔｎｅｒ他がＷＯ ０２／１２９４５に記述している。この３公報はここに参考文献として挙げておく。

ＷＯ ９７／６５０９によると、対象物は、周期性構造（正弦波格子または長方形格子）により照明され、カメラにより撮像され、デジタル化されて、メモリに保存される。続いて、周期性構造は、構造の位相位置に変化が現われるように、画像平面内でシフトされ、改めて撮像が行われ保存される。この過程（シフト、撮像、保存）は、複数回繰り返すことができる。続いて、それぞれの画像はスライス像生成のため相互に差分計算される。その場合、数式としてフーリエ展開式が使用されるが、それには複雑な公式計算器が必要になる。ＷＯ ９８／４５７４５の場合では、個々の撮影間でそれぞれ同じ位相シフトであるという例について、ＷＯ ９７／６５０９からの式を簡略化して導き出したより簡略な公式をスライス像に適用している。

上記方法の適用に必要な物理的周辺条件を実現するのは、実際には非常に困難であることが証明されている。例えば、それぞれの撮像間における照明光度のブレが、生成スライス像にストライプ状のアーチファクトを形成させる原因になる。蛍光性対象物の場合では、同様に欠陥原因となる蛍光色素の経時的退色により追加的な問題が発生する。また、個々の位相シフト過程は一定でなければならないが、実際にはそれを遵守することができない。

したがってＷＯ ０２／１２９４５では、対象物の照明に用いられる光の部分的分離、強度記録、続いての個別撮影の基準化により、時間によって様々である照明光度の影響を調整することが提案されている。一様でない位相シフト過程を考慮するために、連立方程式（そこでは式２２）が用いられている。退色補整のため、最低限度必要な３撮像の代わりに各スライス像につき６撮像を１−２−３−３−２−１の順序で記録し、それより（１−２−３と３−２−１から成る）２つのスライス像を構成、計算して、それら双方の平均値を算出することが提案されている。
この提案の実現には装置面で相当な費用が必要である。その上、追加的撮像のための撮影時間を要し、それに伴い蛍光励起光による照明時に試料が受ける負担も大きくなる。
本発明では、現状技術の欠点を回避して、光学スライス像測定のための改良方法を提供することを課題にしている。

この課題は、主請求項に基づく光学結像システムに対する深度弁別改良方法によって解決される。従属請求項には有利な実施態様が記載されている。

以下では図に基づき本発明をさらに詳しく説明する。
各図はそれぞれ次のものを表わしている ：
図１は、透過型照明の例で、構造化照明を描いた簡略光学図式である。そこには画像形成光路（ピンホール型光路）が記されている。図に描かれた光学的配置の焦点面に位置する１次元周期性構造（透過型格子）（３）が、光源（１）および後の集光光学系（２）によって照明される。光の進行方向で見て、格子の前には平行平面のガラスプレート（４）が配置されている。平行平面プレートの光軸に対する角度は指定されたとおりに調整することができる。続いて、当構造が対象物からは照明側にある光学系（５、６）（コンデンサ）を通ってプレパラート平面（７）に結像する。プレパラートから出た光は、さらに対レンズ（８、９）（対物レンズおよび鏡筒レンズ）を通って、例えばデジタルカメラのＣＣＤマトリックスを配置させることのできる、その後にある焦点面（１０）に結像する。その場合、平行平面のガラスプレート（４）は、指定どおり傾斜させて、プレパラート平面（７）内の対象物に結像する格子構造（３）のシフトに用いられる。対物レンズ（８）は、落射方式蛍光観察の場合では、同時にコンデンサとして用いるのが好ましい。ＷＯ ０２／１２９４５に記述されているように、対象物に現われる輝度分布の記録はデジタルカメラにより、今やガラスプレート（４）の少なくとも３箇所について行われている。

その場合、正弦波／余弦波格子の例における輝度分布Ｉ_ｉ（ｘ、ｙ）は簡単な次式
によって表わすことができる。ただし、ｉ＝０…Ｎ−１は投射された格子のｉ番目の位相位置、Ｎは記録数、ｍ（ｘ、ｙ）は対象物の変調度（つまり、ポジションｘ、ｙにおける検索画像情報）およびφｉは位相値である。

この方程式には３つの未知数Ｉ_０、ｍおよびφ_０が含まれている。したがって、φｉが定められた変数（ｉ＝１、２、３）とすれば、これら未知数の数値決定は少なくとも３つの測定を通して行うことができる。解は測定値から最小二乗法によって得ることができる。その目的のため、次のように、方程式（１）をよりコンパクトな形式で表わし、余弦関数を加法定理に従って書き換えることができる ：
ただし、
ｆ_１（φｉ）＝ｃｏｓφｉ
ｆ_１（φｉ）＝ｓｉｎφｉ
ａ_０（ｘ、ｙ）＝ Ｉ_０（ｘ、ｙ）
ａ_１（ｘ、ｙ）＝ Ｉ_０（ｘ、ｙ）・ｍ（ｘ、ｙ）・ｃｏｓφ_０（ｘ、ｙ）
ａ_２（ｘ、ｙ）＝ −Ｉ_０（ｘ、ｙ）・ｍ（ｘ、ｙ）・ｓｉｎφ_０（ｘ、ｙ）
である。

以上より、関数ｆ_１およびｆ_２は、原則的には自由に選択可能な位相シフトφｉに依存している。非正弦形態であっても周期性構造の場合では、輝度分布Ｉｉ（ｘ、ｙ）はやはり行列展開によって近似的に求めることができる。その場合の計算原理は基本的には変わりない。最小二乗法の解は、次の行列式で求められる ：

ただし、
Ｎ＝測定数（ここでは位相ステップ）および

である。

この行列方程式の解は行列Ｍの逆数処理により数学的に求められる。

このように、ベクトルｂで表わされる測定および格子の位相シフトから検索未知数を求めることができる。特に、ポジション（ｘ、ｙ）における変調度は次式で表わされる。

対象物のデフォーカシングと共に格子コントラスト、延いては対象物における変調度も必然的に低下するので、方程式（８）は深度の弁別された画像（光学的スライス）を得るための算出方程式の意味合いを持っている。方程式（２）および（３）は、位相位置φｉをそれぞれの実地検査所与条件へ適合化させることのできる汎用的な解決方法を提供する。
例えば、ステップ数の選択および対象物空間におけるポジション（結像格子ラインの位置）の選択が自由に行える。

各光学システムは、それぞれの光学伝達関数に従って、光学系の開口数ＮＡおよび光波長λにより決定される限界値

ｋ_ｍａｘ＝４πＮＡ／λ

に到るまでの空間周波数を伝達する。

その結果、対象物空間で結像した純正弦波格子でない格子により、まさにこの限界周波数ｋ_ｍａｘに到るまでの格子周期高調波が伝達され、対象物のフーリエ変換画像に局部最大値として再現される。上記の方法によって得られたスリット像にもこの格子成分が含まれていて、ストライプアーチファクトを招来する。このアーチファクトは、投射された格子の調波がフーリエ変換により空間周波数スペクトル内で局在化され、フィルタ（帯域通過フィルタ、帯域消去フィルタ）の照準設置で除去されることによって解消することができる。

続いての逆変換でアーチファクトのない画像が得られる。
この方法は図２に図示されている。フーリエ平面（１１）内では、格子構造（３）の高調波は（１２）で表わされており、（１３）は零位点で、照明の定常成分、すなわち、均一な非構造化照明に相当する。
当該帯域通過フィルタの調波（１２）を光路内に持ち込むことによって、この領域は結像に関与しなくなるのでデジタルカメラには記録されない。

蛍光性対象物の場合では、ストライプアーチファクトの発現するまた別な原因、つまり蛍光漂白（フォトブリーチング、蛍光退色）が発生する。この場合、（例えば、格子）構造の対象物内への結像によって、様々な程度に漂白された蛍光漂白領域が生じて、最終的には、連続撮影された各画像にストライプアーチファクトが発現する。そのような画像スペクトルにおけるアーチファクトの周波数成分は、投射構造の基本周波数において特に強く現われる。フーリエ変換画像における投射構造の空間周波数成分を然るべき帯域通過フィルタでの除去、続いての逆変換によって、ここでも原則としてアーチファクトのない画像を生成することができる。

蛍光の照明時に発生するこのアーチファクトに対する防止策は、同一対象物のシーケンシャル撮像において、蛍光性である対象物の漂白特性を測定し、その結果を利用することによって行うこともできる。この場合、照明による撮像毎に対象物の照明箇所に退色が現われる。漂白強度は照明強度のほかに、蛍光性対象物の光化学特性および光物理特性にも依存している。対象物全体または選択された対象物ポジション（ｘ、ｙ）に対する強度積分法による対象物の蛍光強度測定から、撮像後における画像シーケンスの再階層化が可能である。

そのためには、アナログカメラおよびアナログ／デジタル変換器により、またはＣＣＤカメラによりコンピュータ内に保存されている画像情報から、対ごとに画像蛍光強度の商を求めて画像ポジション（ｘ、ｙ）における強度を、最高値を基準として定める。そのようにして、蛍光漂白による蛍光損失を考慮することができる。その目的のため、画像撮影時に双方同じように照明され、その結果蛍光を発する観察対象領域（ＲＯＩ）を、図３に基づき定義付けすることができる。この領域における蛍光低下の度合は漂白率の尺度になるので、様々な輝度分布記録の基準化に利用できる。このＲＯＩの定義付けは操作者が然るべき入力手段を用いて、あるいは（結像構造のそれぞれの位置から）自動的に行うことができる。

図３の（１４）は、第１記録における対象物上の輝度分布状況を、（１５）は、結像構造（１６）の位相位置が変化した第２記録における状況を表わしている。観察対象領域（ＲＯＩ）は、両記録間で一様な照明がなされるように定義付けするのが得策である。

本発明に基づく方法によれば、任意に選択した位相シフトφｉに対して実際に光学スライス像の測定が可能になる。誤差防止のため、位相シフトφｉは高い精度で測定しなければならない。これは、システム全体の影響が考慮される下記の較正方法によれば簡単に確保することができる。
その目的のため、まず参照画像（格子の重なった対象物）が撮影され記録される。それに続き、当初は未知である位相位置への格子の僅かな位相シフト（この「制御信号毎のシフト」は予め較正されていなければならない）およびそこでの新たな画像撮影が行われる。得られた２つの画像は計算によって比較される。これは、一般妥当性の制約を受けることなく、差分形成および合計算出によって、またはメリット関数に通ずるその他各種方法によって行うことができる。当初の結果は一般にはストライプを伴う画像である。

作業ステップは次のとおりである ：極僅かな幅の格子シフト、撮像、参照画像との比較、これらは、メリット関数が極値に達するまで、つまり、例えば差分画像の構成成分が基礎ノイズだけになるまで、または合算画像が最大値に達するまでのあいだ繰り返し行われる。後者が満たされれば較正は終了し、そのようにして得られた制御値は後の使用時に呼び出せるように記憶媒体に保存することができる。したがって、この値は１周期分の格子シフトに相当する。それに代わり、合算画像の一様な輝度分布性（すなわち、ストライプ構造の消失）を評価することも可能である。この場合では、その値は半周期分の格子シフトに相当する。この操作は対象物としてミラーを使用して行うのが好ましい。記述の操作過程は手動、自動のいずれでも行なえる。

以下では図４にしたがって、光源輝度の変動、蛍光照明における対象物の退色および非正弦波格子の影響を測定するための代替方法について説明する。
その目的のため、連続するそれぞれの撮影時に考慮される強度比例スカラ係数θｉ＞０による照明の経時変化について記述する（値ｇ_ｉ、ｏ_ｉ等のポジション（ｘ、ｙ）依存性は、以下では簡易化のため明示的には表わしていない）：

モデル観察量ｇ_ｉは堅固な理想的位相画像ｏ_ｉとの積から得られる。この場合、損われた系を完全記述するには、（Ｎ−１）・係数θ_ｉを決定するだけで足りる。

対象物の蛍光照明における退色は、作用ビーム強度に依存して経時的に現われる。ビーム強度は、投射格子があるためにポジションによって異なっているので、格子関数を考慮した処理が必要である。退色度が照明強度の１次関数であると想定すれば、次のように簡単に記述することが可能である ：

ただし、ポジションによって異なる係数から成るベクトル０＜ｋ_ｉ≤１は、その成分毎に、それぞれまだ退色していない観察対象ｏ_ｉと掛け合わせた場合に、ビーム強度の経過および選択可能な退色指数１＞ｄ＞０に依存した減衰率を表わすというものである ：

式（１０）から明かになるように、各観察対象は常に、先行する全位相撮影の係数と掛け合わされる。これらの先行位相撮影は既に一部退色しており、現位相とオーバラップしている。ｋ_ｉはこの式に代わり、冪指数で、すなわち特定蛍光色の減衰曲線に対応して表わすこともできる。

非正弦波状の格子経過が、特に格子周波数の整数倍毎にストライプアーチファクトを惹き起こす原因になっている。個々の撮影における格子の周波数、位相位置および曲線形態について詳しいデータがあれば、正弦波状の経過を生成することが可能である。修正係数Ｉ_Ｋｏｒｒは式（１０）に類似するが、ここでは格子の周期性に沿って方程式（１２）により考察しなければならない。

観察対象の格子関数は記号ｆ_{Ｇｉｔｔｅｒ}で表示される。正弦波状の理想的な強度経過はｆ_{Ｓｉｎｕｓ}で表わされる。分母の定数ｓはｆ_{Ｇｉｔｔｅｒ}を非常に小さい数で除した割り算を示している。観察対象の格子関数は様々な方法で表現することができる。例えば、台形関数用フーリエ級数の使用によって、次のような簡単な近似法の可能性が提供される。

ｂは台形側方の勾配を変更させる変数、ｘは位置座標およびＭはフーリエ係数の使用数である。代替法としてより精確なモデルを使用することもできる。それには顕微鏡ＰＳＦ（点像分布関数）の横成分を持つ真の格子関数の畳み込みが必要である。決定的に重要なのは修正ベクトルＩ_Ｋｏｒｒの位相同期性および周波数同期性の合成である。これが保証されれば、式（１１）によるのと同様に、成分毎の簡単な乗法により観察対象について
次のように修正結果を得ることができる。

したがって、修正された位相画像は次式より得られる ：

式（１５）を方程式（３）および（８）に代入することにより、それに加えてパラメータθ_ｉ、ｄおよびｂも含む、再構築のための一般式が得られる。

提起された問題の解は、例えば極端な課題への変形によって行うことができる。パラメータθ_ｉ、ｄおよびｂはメリット関数の利用下、数値最適化の補助手段により求めることができる。

ここで、Ｆ｛ａ｝ωは方程式（５）から得られる結果ベクトルａのフーリエ変換成分である。この場合、フーリエ変換のほかに別な関数変換で適しているものもある。

上に挙げたメリット関数（方程式（１６））は例に過ぎず、その他のメリット関数も適用可能である。
目的は、式（１６）が最小になるようにθ_ｉ、ｄおよびｂを変更することである。それには、多数の方法があって、その対象は決して勾配検索法または直線検索法に限定されることはない。また、各変換係数は係数α_ｉによって重み付けされている。それにより、アルゴリズムを様々な信号／雑音間距離または優先周波数に適合させることができる。これに関して有利な値は、定常成分に重み付けをした次の関係で得られる ：

ここで説明した修正装置の一般的な機能態様は図４の経過図式から見て取れる。パラメータθ_ｉ、ｄおよびｂを効率よく決定するには、投射の利用が有利である。しかし、パラメータの決定作業はこの１次元に限定されることなく、平面の場合も同様に機能する。

対象物に結像した格子構造において異なった位相位置ｉ＝１、２、３、…Ｎを持つ、少なくとも３つの輝度分布（１８）、（１９）、（２０）、場合によってはそれ以上の輝度分布（２１）…がカメラにより固定フォーカス位置で記録される。これら２次元の輝度分布は、投射によりそれぞれ１次元の分布（２２）、（２３）、（２４）および場合によっては（２５）に転換される。パラメータθ_ｉ、ｄおよびｂの初期値を用いて、公式（１５）に従って修正された位相画像が算出され、それより、公式（３）からの式により第１光学スライス像が算出される（２６）。その場合、パラメータθ_ｉ、ｄおよびｂの初期値は操作者によるインプット、あるいはまた自動見積が可能である。この見積にはメリット関数の最小化／最大化のための公知ツールが利用できる（例えば、シンプレックス、マルカイドなど）。

結果から（２７）では公式（１６）に基づくメリット関数が決定され、（２８）で最小の達成がテストされる。まだ最小に達していなければ、パラメータθ_ｉ、ｄおよびｂの新たな見積が行われ（２９）、（２６）により継続される。メリット関数の最小が達成された時点、すなわち、パラメータθ_ｉ、ｄおよびｂがそれ以上変化するとその値が再び上昇に転じる時点で最適パラメータが決定され、（３０）において公式（１５）、（３）および（８）により光学スライス像が決定され、（３１）でアウトプットされる。それに続き、フォーカス位置を変化させ、上記の方法に従ってまた別な光学スライス像を決定すること等々…ができる。対象物のｚ方向でのこの走査により、スライス像のｚスタックが、したがって対象物の３次元表示が可能になる。

対象物の蛍光照明における退色現象を考慮する別な式として、局部修正関数の決定に基づいたものがある。それには、その局部の値が格子周期τ＝２π／ωを通じて形成される修正関数ｋ_ｉ（ｘ）として
が想定される。この関数式はｙ、すなわち、ここでは格子周期に対して垂直方向を取る座標による影響も考慮に入れることができるが、しかしここでは見やすくする意味からも省いてある。

ポジション可変の退色関数θ_ｉ（ｘ）使用のもと、近似式ｇ_ｉ（ｘ）＝θ_ｉ（ｘ）Ｉ_０（ｘ）［１＋ｍｃｏｓ（ωｘ＋ψ_１）］ （公式（１）および（９）参照）による然るべき変形により
が得られる。退色関数が一様である場合、この修正関数は１／θ_ｉ（ｘ）に比例する。

１格子周期以上に亘って平均値を取れば有利であることが明かになっている。
さらに、この修正関数は退色関数の端近くではオーバシュートする傾向があり、そのため結果に追加誤差の生じる場合のあることが判明した。それの解消には、オーバシュートの発現を平均値との比較によって認識し、その上で誤った値を取り止めて、それに代わり公式（９）および（１６）から予備設定したグローバル値１／θ_ｉを使用するのが有効であると実証された。

当該式は次のようになる ：

この場合、閾値δは画像内における修正関数局部変数の変動パーセンテージとして、例えば５％と置くのが有利である。１／θ_ｉの代わりに、２次元空間で閾値以下に算定された、補間法によるｋ_ｉ（ｘ）値も使用することができる。
以上のように、本方法では、例えば蛍光団の違いまたはその特性偏差を原因としてポジション毎に変動する退色現象が考慮に入れられている。
本発明の実現は上記の実施例に拘束されるものではない。専門家により展開されたさらなる開発、あるいは光源の輝度変動または対象物の蛍光照明における退色を考慮に入れるための別な方程式は本発明の領域を逸脱するものではない。

構造化照明装置を持つ顕微鏡の原理構造 本発明に基づく帯域通過フィルタのフーリエ平面における模式図 周期性構造の結像位置 線形最適化による修正値測定のための経過図

符号の説明

１ 光源
２ 集光光学系
３ 透過型格子
４ ガラスプレート
５ 光学系
６ 光学系
７ プレパラート平面
８ 対物レンズ
９ 鏡筒レンズ
１０ 焦点面
１１ フーリエ平面
１２ 高調波
１３ 零位点
１８〜２１ 輝度分布
２２〜２５ 1次元分布

Claims (1)

1. ａ）対象物が周期性構造によって照明されるステップと、
   ｂ）そのようにして生じた輝度分布が記録されるステップと、
   ｃ）周期性構造の位相位置がシフトされるステップと、
   ｄ）ステップａ）〜ｃ）が、少なくとも３つの輝度分布が記録されるまで複数回繰り返されるステップと、
   ｅ）対象物の輝度分布を得るため、記録された輝度分布が相互に計算されるステップと、
   ｆ）ステップｃ）からの位相位置シフトが記録されるステップと、
   ｇ）照明の輝度変動が捕捉されるステップと、
   ｈ）対象物の蛍光照明における退色現象が測定され、そしてステップｆ）〜ｈ）で得られた結果に基づいてオフセットが行なわれるステップとからなり、
   ステップｆ）〜ｈ)で得られた結果から次式のメリット関数の援助で線形最適化により決定されることを特徴とする光学システムにおける深度弁別改良のための方法。
   ただし、式中のＦ｛ａ｝がベクトルａの関数変換形態、θ _ｉ が照明の輝度変動を特徴付けるためのスカラ係数、ｄが対象物の蛍光照明における退色度、ｂが対象物への照明の非正弦波状分布を特徴付けるための係数、およびα _ｉ が特に、様々な信号／雑音間距離または優先周波数などの撮影条件への適合化のための重み付け係数である。