【結像計算:第4回】散乱光は結像するとどう見えるのか?偏光まで含めた結像計算

Python

はじめに

 散乱光分布を計算していて、ふとこれレンズで結像するとどう見えるんだろう?と思ったので波動光学に基づいてPythonで計算してみたシリーズの第4回です。
今回は今まで触れてなかった偏光の影響を取り込み、高NAの場合の瞳面の設定方法を説明、これまでの総まとめとします。

偏光の計算方法

 前回pyGDMで散乱光の電場を求めました。
この時x,y,zの3成分が得られました。この座標系は局所座標でなく、グローバル座標であり、x,y,zは直行しています。FFTは線形な演算なので、x,y,z成分に対してそれぞれ結像計算を行い、最後に強度を合成すればよいです。

計算コード

 散乱光のデータは”E”というリスト形式で得られました。
詳細は前回記事参照ください。https://yohaku-labs.com/pygdm-scat/
Eは所与として散乱計算は省略します。
前回記事までのコードを活用しています。

import numpy as np
import matplotlib.pyplot as plt
plt.style.use("ggplot")
from scipy.interpolate import griddata

#x,y,z成分の電場取得、整形
Ex =E[2][:,0].reshape(E[0].shape)
Ey =E[2][:,1].reshape(E[0].shape)
Ez= E[2][:,2].reshape(E[0].shape)

#後の計算用
theta =np.copy(E[0])
kx = np.sin(E[0])*np.cos(E[1])
ky = np.sin(E[0])*np.sin(E[1])

#電場の大きさをグラフ化
fig, axes = plt.subplots(1, 3, figsize=(15,5))

components = [Ex, Ey, Ez]
titles = ["Ex", "Ey", "Ez"]

for ax, comp, title in zip(axes, components, titles):

    im = ax.contourf(
        kx,
        ky,
        np.abs(comp),
        levels=50)

    ax.set_title(title)
    ax.set_xlabel("kx")
    ax.set_ylabel("ky")
    ax.set_xlim(-1,1)
    ax.set_ylim(-1,1)
    ax.set_aspect('equal')

    fig.colorbar(im, ax=ax)

plt.tight_layout()
plt.show()

 x,y,z方向の電場の大きさはこんな感じです。
e-field.map.jpg

 検出NA=0.9、収差は無いものとします。
前前回に言及したように瞳面では、512×512のうち瞳の半径は512の1/4くらいにしないとうまく結像計算ができません。これはFFTは周期境界条件を仮定するため、瞳を配列いっぱいに広げると折り返しの影響がでるためだったのでした。
 座標の取り方を調整する必要あるので、NAのほかにmask_ratioという変数を導入します。

N=512
NA=0.9
mask_ratio =0.5

#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)
#電場強度を補間
points = np.column_stack([kx.ravel(),ky.ravel()])
values_x = Ex.ravel()
values_y = Ey.ravel()
values_z = Ez.ravel()
values_theta=theta.ravel()

#電場とthetaを補間
Ex_interpolate = griddata(points,values_x,(KX_im,KY_im),method='linear')
Ey_interpolate = griddata(points,values_y,(KX_im,KY_im),method='linear')
Ez_interpolate = griddata(points,values_z,(KX_im,KY_im),method='linear')
theta_interpolate = griddata(points,values_theta,(KX_im,KY_im),method='linear')

#nanが含まれるので0にする
Ex_interpolate[np.isnan(Ex_interpolate)]=0
Ey_interpolate[np.isnan(Ey_interpolate)]=0
Ez_interpolate[np.isnan(Ez_interpolate)]=0
theta_interpolate[np.isnan(theta_interpolate)]=0

#瞳関数とアプラナート条件を適用
weight =np.cos(theta_interpolate)**0.5
Ex_interpolate *= mask*weight
Ey_interpolate *= mask*weight
Ez_interpolate *= mask*weight
#FFTを実行
Ex_img = np.fft.fftshift(
    np.fft.ifft2(
        np.fft.ifftshift(Ex_interpolate)
    )
)
Ey_img = np.fft.fftshift(
    np.fft.ifft2(
        np.fft.ifftshift(Ey_interpolate)
    )
)
Ez_img = np.fft.fftshift(
    np.fft.ifft2(
        np.fft.ifftshift(Ez_interpolate)
    )
)

 結像面のスケーリングを行います。対物レンズ、結像レンズともに焦点距離は100mmで倍率1倍とします。

#像面のスケーリング
f_obj = 100 #対物レンズ焦点距離
f_img = 100 #結像レンズ焦点距離
lam = 0.532 *10**-3#波長

#瞳面での1ピクセルの物理的大きさを求める
dx = f_obj * NA / (N/2* mask_ratio)

k = np.fft.fftfreq(N, d=dx)
x_real  = lam * f_img * k
x_real = np.fft.fftshift(x_real)
y_real = np.copy(x_real)

X_real,Y_real=np.meshgrid(x_real,y_real)

#グラフ化
plt.figure(figsize=(6,5))
pix =20

plt.contourf(X_real[-pix+N//2:pix+N//2,-pix+N//2:pix+N//2],Y_real[-pix+N//2:pix+N//2,-pix+N//2:pix+N//2],np.abs(Ex_img[-pix+N//2:pix+N//2,-pix+N//2:pix+N//2])**2+np.abs(Ey_img[-pix+N//2:pix+N//2,-pix+N//2:pix+N//2])**2+np.abs(Ez_img[-pix+N//2:pix+N//2,-pix+N//2:pix+N//2])**2,levels=50)
plt.xlabel("X(mm)")
plt.ylabel("Y(mm)")
plt.gca().set_aspect('equal')
plt.show()

 集光スポットはこんな感じです。NA=0.9でもピラミッド形状になっていないですね。
像がボケて左右方向に広がっていますが、元の散乱光の強度分布を反映しています。
spot_scat.jpg

 例えば電場を一様分布(等方的)としてみると像はもっとシャープに見えます。
散乱光強度分布によるボケの影響も多少あることがわかります。

#これで電場置き換え
Ex =np.ones_like(E[0])
Ey =np.ones_like(E[0])
Ez =np.ones_like(E[0])

偏光の取り扱い

 最初に思い違いをしていました。
光はレンズで曲げられて、瞳面ではz方向に進みます。
波数ベクトルと偏光振動面は直行するため、瞳面ではx,y方向に電場振動面がきます。
光を3次元行列で回転させて瞳面上での電場を出した後にFFTかなと思っていました。
コーディングがとても大変そうです。
 実際はx,y,zはグローバル座標で直行して独立、更にpyGDMが3成分全部計算してくれるため、それぞれ計算すれば十分でした。
 最初は違和感が凄かったのですが、波動光学の結像はフラウンフォーファー回折と同じ形で書けて、無限遠での回折現象を扱っているのと同じと思うと腹落ちしました。

最後に

 波動光学に基づいた結像計算シリーズはここまでです。
ゼルニケ多項式を使って収差を入れてみる、散乱光の計算に意味を持たせる、瞳面で偏光フィルターを入れてみる、位相操作をするなどすれば色々面白い現象を計算出来そうです。

 結像計算の理論面についてさらに詳しく学びたい方はGoodmanの『Introduction to Fourier Optics』がおすすめです。フーリエ光学の定番教科書で、今回扱った内容の背景も体系的に学べます。邦訳も出ています。

 本記事では実装を中心に解説しましたが、波動光学の理論をしっかり理解したい場合は、『マックス・ボルン 光学の原理』がおすすめです。フーリエ光学や回折理論、電磁気学に基づく光学を体系的に学べる古典的名著です。

コメント

タイトルとURLをコピーしました