【raytracing】レンズに入らなかった光線を除外した光線図を描く方法[Python]

2025年5月26日

URLをコピーしました！

raytracing
アパチャーに入らなかった光線を除外する方法
アパチャーに入らなかった光線を除外した光線図

raytracing

前回、raytracingライブラリでアパチャーのサイズの取得方法とレンズ位置の図示の方法を紹介しました。

あわせて読みたい

【raytracing】アパチャーのサイズの取得方法とレンズ位置の図示の方法[Python] raytracing 前回、焦点距離を変化させ、特定の高さ（サイズ）の光線を得る方法を紹介しました。今回はaytracingライブラリでアパチャーのサイズの取得方法とレンズ位置…

今回はraytracingライブラリでレンズに入らなかった光線を除外した光線図を描く方法を紹介します。

前回のおさらいをするとこんなプログラムを書いて、こんな結果となっていました。

from raytracing import *
import matplotlib.pyplot as plt

path = ImagingPath()
path.append(Space(d=50))
path.append(Lens(f=10, diameter=5))
path.append(Space(d=50))

angle_list = np.linspace(-10, 10, 21)

#----------位置、角度情報取得ここから----------
y_list = []; z_list = []; theta_list = []
for angle in angle_list:
    ray = Ray(y=0, theta=np.radians(angle))
    trace = path.trace(ray)

    y_optics = []; z_optics = []; theta_optics = []
    for data in trace:
        y_optics.append(data.y)
        z_optics.append(data.z)
        theta_optics.append(data.theta)
    y_list.append(y_optics)
    z_list.append(z_optics)
    theta_list.append(theta_optics)
#----------位置、角度情報取得ここまで----------

#----------結像位置取得ここから----------
cross_list = []
for y_ray, z_ray in zip(y_list, z_list):
    z1 = z_ray[-2]; z2 = z_ray[-1]
    y1 = y_ray[-2]; y2 = y_ray[-1]
    cross = -y1*(z2-z1)/(y2-y1)
    cross_point = cross + z1
    cross_list.append(cross_point)
#----------結像位置取得ここまで----------

#----------光線プロットここから----------
fig = plt.figure()
plt.clf()

for y_optics, z_optics in zip(y_list, z_list):
    plt.plot(z_optics, y_optics, color="red", lw=0.5)

#----------結像位置プロットここから----------
for point in cross_list:
    plt.annotate("", xy=[point, 2], xytext=[point, 7], arrowprops=dict())
#----------結像位置プロットここまで----------

ylim_min, ylim_max = plt.ylim()

#----------焦点距離取得＆プロットここから----------
position = 0
for optics in path:
    optics_position = optics.B
    position = position+optics_position
    ffl = optics.frontFocalLength()
    bfl = optics.backFocalLength()
    dia = optics.apertureDiameter
    if ffl != None and bfl != None:
        plt.scatter(position+ffl, 0, color="black", s=5)
        plt.scatter(position-bfl, 0, color="black", s=5)
        if dia == float(inf):
            plt.annotate("", xy=[position, 0.5*ylim_min], xytext=[position, 0.5*ylim_max], arrowprops=dict(arrowstyle='<->'))
        else:
            plt.annotate("", xy=[position, -dia/2], xytext=[position, dia/2], arrowprops=dict(arrowstyle='<->'))
#----------焦点距離取得＆プロットここまで----------

plt.xlabel("Distance")
plt.ylabel("Height")
plt.show()
#----------光線プロットここまで----------

実行結果

　問題点はアパチャーに入らなかった光線が途中まで描かれていることで、これをどうにかして取り除きたいというのが今回の目的です。

それでは始めていきましょう。

アパチャーに入らなかった光線を除外する方法

まずはアパチャーに入らなかった光線を除外する方法を紹介していきます。

やり方としては光線解析をそれぞれの光学系で行なっていき、アパチャーに入らなかった時点でその光線は取り除く方法ととりあえず全部光線解析を行い、後からアパチャーに入っているか評価する方法が考えられます。

今回は全部光線解析をして、後からアパチャーに入っているか評価する方法で行なっていきます。

そのためまずはアパチャーの情報を取得しますが、これは「y」や「z」、「theta」の情報と同時に取れるので、その部分を改変しました。

from raytracing import *
import matplotlib.pyplot as plt

path = ImagingPath()
path.append(Space(d=50))
path.append(Lens(f=10, diameter=5))
path.append(Space(d=50))

angle_list = np.linspace(-10, 10, 21)

#----------位置、角度情報取得ここから----------
y_list = []; z_list = []; theta_list = []; aperture_list = []
for angle in angle_list:
    ray = Ray(y=0, theta=np.radians(angle))
    trace = path.trace(ray)

    y_optics = []; z_optics = []; theta_optics = []; aperture_optics = []
    for data in trace:
        y_optics.append(data.y)
        z_optics.append(data.z)
        theta_optics.append(data.theta)
        aperture_optics.append(data.apertureDiameter)
    y_list.append(y_optics)
    z_list.append(z_optics)
    theta_list.append(theta_optics)
    aperture_list.append(aperture_optics)
#----------位置、角度情報取得ここまで----------

print(aperture_list)

実行結果
[[inf, inf, inf, 5, 5, 5], [inf, inf, inf, 5, 5, 5], 
[inf, inf, inf, 5, 5, 5], [inf, inf, inf, 5, 5, 5], 
[inf, inf, inf, 5, 5, 5], [inf, inf, inf, 5, 5, 5], 
[inf, inf, inf, 5, 5, 5], [inf, inf, inf, 5, 5, 5], 
[inf, inf, inf, 5, 5, inf], [inf, inf, inf, 5, 5, inf], 
[inf, inf, inf, 5, 5, inf], [inf, inf, inf, 5, 5, inf], 
[inf, inf, inf, 5, 5, inf], [inf, inf, inf, 5, 5, 5], 
[inf, inf, inf, 5, 5, 5], [inf, inf, inf, 5, 5, 5], 
[inf, inf, inf, 5, 5, 5], [inf, inf, inf, 5, 5, 5], 
[inf, inf, inf, 5, 5, 5], [inf, inf, inf, 5, 5, 5], 
[inf, inf, inf, 5, 5, 5]]

次にそれぞれの光線、そしてそれぞれの光学系において光線の高さである「y」とアパチャーのサイズを比較・評価します。

そして全ての光線の高さが全ての光学系におけるアパチャーサイズ（実際にはその半分）より小さければその光線の情報を新たなリストに格納するという流れで行います。

from raytracing import *
import matplotlib.pyplot as plt

path = ImagingPath()
path.append(Space(d=50))
path.append(Lens(f=10, diameter=5))
path.append(Space(d=50))

angle_list = np.linspace(-10, 10, 21)

#----------位置、角度情報取得ここから----------
y_list = []; z_list = []; theta_list = []; aperture_list = []
for angle in angle_list:
    ray = Ray(y=0, theta=np.radians(angle))
    trace = path.trace(ray)

    y_optics = []; z_optics = []; theta_optics = []; aperture_optics = []
    for data in trace:
        y_optics.append(data.y)
        z_optics.append(data.z)
        theta_optics.append(data.theta)
        aperture_optics.append(data.apertureDiameter)
    y_list.append(y_optics)
    z_list.append(z_optics)
    theta_list.append(theta_optics)
    aperture_list.append(aperture_optics)
#----------位置、角度情報取得ここまで----------

#----------aperture size確認ここから----------
y_list_aperture = []; z_list_aperture = []; theta_list_aperture = []; aperture_list_aperture = []
for y_ray, z_ray, theta_ray, aperture_ray in zip(y_list, z_list, theta_list, aperture_list):
    aperture_test = []
    for y_optics, aperture_optics in zip(y_ray, aperture_ray):
        if abs(y_optics) <= aperture_optics/2:
            aperture_test.append(0)
        else:
            aperture_test.append(1)
    if not 1 in aperture_test:
        y_list_aperture.append(y_ray)
        z_list_aperture.append(z_ray)
        theta_list_aperture.append(theta_ray)
        aperture_list_aperture.append(aperture_ray)
#----------aperture size確認ここまで----------

print(y_list_aperture)
print(aperture_list_aperture)

実行結果
[[0, 0, np.float64(-1.7453292519943295), np.float64(-1.7453292519943295), 
np.float64(-1.7453292519943295), np.float64(5.235987755982991)], 
[0, 0, np.float64(-0.8726646259971648), np.float64(-0.8726646259971648), 
np.float64(-0.8726646259971648), np.float64(2.6179938779914953)], [0, 0, np.float64(0.0), 
np.float64(0.0), np.float64(0.0), np.float64(0.0)], [0, 0, np.float64(0.8726646259971648), 
np.float64(0.8726646259971648), np.float64(0.8726646259971648), 
np.float64(-2.6179938779914953)], [0, 0, np.float64(1.7453292519943295), 
np.float64(1.7453292519943295), np.float64(1.7453292519943295), 
np.float64(-5.235987755982991)]]
[[inf, inf, inf, 5, 5, inf], [inf, inf, inf, 5, 5, inf], [inf, inf, inf, 5, 5, inf], 
[inf, inf, inf, 5, 5, inf], [inf, inf, inf, 5, 5, inf]]

これでアパチャーに入った光線のみ抽出することができました。

アパチャーに入らなかった光線を除外した光線図

あとは前回までの光線図を作成するプログラムと結合するだけです。

ということでこんな感じです。

from raytracing import *
import matplotlib.pyplot as plt

path = ImagingPath()
path.append(Space(d=50))
path.append(Lens(f=10, diameter=5))
path.append(Space(d=50))

angle_list = np.linspace(-10, 10, 21)

#----------位置、角度情報取得ここから----------
y_list = []; z_list = []; theta_list = []; aperture_list = []
for angle in angle_list:
    ray = Ray(y=0, theta=np.radians(angle))
    trace = path.trace(ray)

    y_optics = []; z_optics = []; theta_optics = []; aperture_optics = []
    for data in trace:
        y_optics.append(data.y)
        z_optics.append(data.z)
        theta_optics.append(data.theta)
        aperture_optics.append(data.apertureDiameter)
    y_list.append(y_optics)
    z_list.append(z_optics)
    theta_list.append(theta_optics)
    aperture_list.append(aperture_optics)
#----------位置、角度情報取得ここまで----------

#----------aperture size確認ここから----------
y_list_aperture = []; z_list_aperture = []; theta_list_aperture = []; aperture_list_aperture = []
for y_ray, z_ray, theta_ray, aperture_ray in zip(y_list, z_list, theta_list, aperture_list):
    aperture_test = []
    for y_optics, aperture_optics in zip(y_ray, aperture_ray):
        if abs(y_optics) <= aperture_optics/2:
            aperture_test.append(0)
        else:
            aperture_test.append(1)
    if not 1 in aperture_test:
        y_list_aperture.append(y_ray)
        z_list_aperture.append(z_ray)
        theta_list_aperture.append(theta_ray)
        aperture_list_aperture.append(aperture_ray)
#----------aperture size確認ここまで----------

y_list = y_list_aperture
z_list = z_list_aperture
theta_list = theta_list_aperture
aperture_list = aperture_list_aperture

#----------結像位置取得ここから----------
cross_list = []
for y_ray, z_ray in zip(y_list, z_list):
    z1 = z_ray[-2]; z2 = z_ray[-1]
    y1 = y_ray[-2]; y2 = y_ray[-1]
    cross = -y1*(z2-z1)/(y2-y1)
    cross_point = cross + z1
    cross_list.append(cross_point)
#----------結像位置取得ここまで----------

#----------光線プロットここから----------
fig = plt.figure()
plt.clf()

for y_optics, z_optics in zip(y_list, z_list):
    plt.plot(z_optics, y_optics, color="red", lw=0.5)

#----------結像位置プロットここから----------
for point in cross_list:
    plt.annotate("", xy=[point, 1], xytext=[point, 3], arrowprops=dict())
#----------結像位置プロットここまで----------

ylim_min, ylim_max = plt.ylim()

#----------焦点距離取得＆プロットここから----------
position = 0
for optics in path:
    optics_position = optics.B
    position = position+optics_position
    ffl = optics.frontFocalLength()
    bfl = optics.backFocalLength()
    dia = optics.apertureDiameter
    if ffl != None and bfl != None:
        plt.scatter(position+ffl, 0, color="black", s=5)
        plt.scatter(position-bfl, 0, color="black", s=5)
        if dia == float(inf):
            plt.annotate("", xy=[position, 0.5*ylim_min], xytext=[position, 0.5*ylim_max], arrowprops=dict(arrowstyle='<->'))
        else:
            plt.annotate("", xy=[position, -dia/2], xytext=[position, dia/2], arrowprops=dict(arrowstyle='<->'))
#----------焦点距離取得＆プロットここまで----------

plt.xlabel("Distance")
plt.ylabel("Height")
plt.show()
#----------光線プロットここまで----------

実行結果