はじめに
散乱光分布を計算していて、ふとこれレンズで結像するとどう見えるんだろう？と思ったので波動光学に基づいてPythonで計算してみたシリーズの第5回です。
前回までで結像計算を一通り扱い満足していたのですが、
散乱光を結像することの特徴は何か？という点についての結論が出ていませんでした。
散乱光を発する光源自体の大きさは無視できるとすると、集光サイズはPSF(点像分布関数,Point Spread Function)で決まります。
PSFに散乱光が及ぼす影響をデルタ関数、フーリエ変換を通して説明します。
前回記事
δ(デルタ)関数とフーリエ変換
点光源からの光を集光した際は、より小さく集光出来た方が結像性能が良いと言えます。
集光スポットの大きさを0にしたいのです。
“大きさがゼロ”を数式で表すにはδ(デルタ)関数を使います。δ関数は以下の性質を満たします。
$$
\delta(x)=\begin{cases}\infty & (x = 0)\\
0 &(x \neq0)
\end{cases}
$$
そして、デルタ関数はフーリエ変換を使った下式で表せることが知られています。
$$
\delta(x) = \frac{1}{2\pi} \int_{-\infty}^{\infty} e^{ikx} dk
$$
$e^{ikx}$は波数が$k$の波を表しています。
つまり、デルタ関数はすべての波数$k$の波を同じ重み1で足し合わせたものなのです。
PSFがδ関数にならない理由
対象の光学系が集光した時の像の大きさの最小値はPSFが表しているとしてよいです。
そして、ASF(振幅分布関数,Amplitude Spread Function)の絶対値を取って2乗したものがPSFなので
ASFがδ関数になってほしいのですが、これは不可能です。
以下で理由を見ていきましょう。
まず、ASF、PSFの関係を整理します。
瞳関数
↓ フーリエ変換
ASF
↓ 絶対値の2乗
PSF
ASFがδ関数となるには、瞳関数がすべての波数を同じ重みでもつ必要があります。
瞳関数は、瞳面における光の複素振幅の分布を表す関数です。
もし仮に光学系がすべての光を拾えたとしても(NA=1)、
これは空間周波数が0~$\frac{1}{\lambda}$までの光しか通せません。
よって、ASFはδ関数にはなれず有限の広がりを持ち、PSFも同様に広がりを持つのです。
散乱光はPSFを広げる
瞳面上での光の強度が一定でない場合、対応する空間周波数成分に対して重みづけがされることになります。
これはδ関数の式で$e^{ikx}$の前に強度分布を表す関数が入ることになり、
方向としてはPSFは更に広がります。
散乱光の空間分布を計算結果を示します。
この表示はそのまま瞳面での強度分布と考えられ、強度は一定ではないことが分かります。
また位相も一定ではありません。
散乱光分布の求め方や計算条件については下記の記事を参照ください。
https://yohaku-labs.com/pygdm-scat/
次に瞳面で強度一定としたときと、散乱光の強度分布を使ったPSFを示します。
散乱光の方がボケて左右に広がっているのが特徴的ですが、散乱光分布を反映しています。
計算条件詳細は以下の記事を参照ください。
https://yohaku-labs.com/psf-zernike/
散乱光は点光源とみなせるものの、瞳面上での強度分布と位相分布が一定でないため、ASF,PSFを広げると考えられ、結像性能の悪化の一因となります。
PSFには光学系が通した周波数成分が詰まっている
光学系が空間周波数に制限を加えていると考えると、
PSFの中には光学系が通した空間周波数が含まれていると言えます。
PSFはその名の通り点像で、これをFFTすれば光学系が透過する周波数の範囲とどれくらい透過するかを計算できるのではないでしょうか？
Pythonで計算してみましょう。
瞳関数は円形で強度一定とし、NAの値を変えながら計算してみます。
import numpy as np
import matplotlib.pyplot as plt
N=512
NAs=[0.01,0.05,0.1,0.3,0.6,1]
mask_ratio =0.4
f_obj =10#f :対物レンズ焦点距離 mm
lam = 0.532 *10**-3 #: 波長
k_lam = lam**-1
res_mtf=[]
res_freq=[]
for NA in NAs:
#kを等間隔で作成
kx_im = np.linspace(-1/mask_ratio*NA,1/mask_ratio*NA,N)
ky_im = np.linspace(-1/mask_ratio*NA,1/mask_ratio*NA,N)
KX_im,KY_im =np.meshgrid(kx_im,ky_im)
#瞳関数
mask = (np.sqrt(KX_im**2+KY_im**2) < NA).astype(int)
ASF = np.fft.fft2(np.fft.fftshift(mask))
# FFTは左上を原点として扱うため、
# 中心を原点とした瞳関数をfftshiftで並べ替える
ASF = np.fft.fftshift(ASF)
PSF = np.abs(ASF) ** 2
#瞳面での1ピクセルの物理的大きさ
dx = f_obj * NA / (N / 2) *mask_ratio
k = np.fft.fftfreq(N, d=dx)
x_real = lam * f_obj * k
x_real = np.fft.fftshift(x_real)
y_real = np.copy(x_real)
X,Y = np.meshgrid(x_real,y_real)
delt = x_real[0]-x_real[1]
freq = np.fft.fftfreq(N, d=delt)
freq /= k_lam
otf=np.fft.fft2(np.fft.fftshift(PSF))
MTF =np.abs(otf)
res_mtf.append(MTF[0])
res_freq.append(freq)
計算した結果は下記です。
横軸のfrequencyの1は、空間周波数$\frac{1}{\lambda}$に対応します。
縦軸は周波数成分をどの程度、通すかを表しています。
NAが大きくなるほど、より高周波成分を通していること、またすべてのNAでも空間周波数が高くになるにつれて通しづらくなっていることが分かります。
実は今計算したものはOTF(光学伝達関数,Optical Transfer Function)やMTF(変調伝達関数,Modulation Transfer Function)と呼ばれるものです。
カメラで撮った画像の像の品質を定義するのに使われているのをよく見ます。
関係をまとめておきます。
瞳関数
↓フーリエ変換
ASF
↓絶対値の2乗
PSF
↓フーリエ変換
OTF
↓絶対値をとる
MTF
今回はNAだけでしたが、収差を考慮すると瞳関数に複素数成分が追加されることなり、フーリエ変換時に位相ずれが入るため、ASF,PSFは広がり、MTF,OTFの周波数帯が削られることになります。
まとめ
散乱光を結像するとどうなるかという点について考察をしました。
・光学系は透過する空間周波数に制限を加えるため、光を1点に集光できない
・散乱光は瞳面上で強度、位相が一定でないため、PSFを広げる一要因となる
・PSFをフーリエ変換をすることで、OTF、OTFの絶対値をとることでMTFを求められる
・NAが大きいほど、光学系が透過する空間周波数が高周波まで透過する
“光は波長の半分にしか集光できなくてエアリーディスクになる”というのはよく言われていて、私もそのような表面的な理解でした。
今回はPSFってどういうものかという点に着目して考えたことで、理解が深まりました。
結像させようと思うとピント合わせ等色々頑張るのですが、結局のところはすべてフーリエ変換で記述できるのは理論として非常に美しいと思います。
参考文献
フーリエ光学を日本語で学びたい方におすすめの一冊です。 回折やフーリエ変換、ホログラムも分かりやすく解説しており、波動光学における結像理論を体系的に理解できます。
光学計算を進めるには今ではPythonの活用が必須と言えます。
体系的に学びたい方は、オンライン講座を活用するのも一つの方法です。
コメント