用Python从零实现光线追迹：手把手教你模拟单个折射球面的光路

用Python从零实现光线追迹手把手教你模拟单个折射球面的光路光线追迹是几何光学的核心计算工具它能精确模拟光线在光学系统中的传播路径。本文将带你用Python实现一个完整的折射球面光线追迹器从公式推导到代码实现再到可视化分析完整复现光学设计中的关键计算过程。1. 环境配置与基础理论1.1 Python科学计算三件套我们需要以下核心库import numpy as np import matplotlib.pyplot as plt from dataclasses import dataclassNumPy处理向量化数学运算Matplotlib实现光线路径可视化Dataclass创建光线数据结构1.2 折射球面的符号规则光学计算需要统一的符号约定这是后续所有计算的基础参数类型正方向定义代码表示沿轴距离与光线同向float垂直距离光轴上方float角度值顺时针旋转弧度制注意所有角度计算需转换为弧度制math.radians()和math.degrees()可实现角度与弧度转换2. 核心算法实现2.1 光线数据结构设计使用Python的dataclass封装光线参数dataclass class OpticalRay: L: float # 物方截距 U: float # 物方孔径角(弧度) wavelength: float 589.3 # 钠黄光默认波长(nm)2.2 折射计算核心函数实现公式(1)-(4)的向量化计算def refract_sphere(ray: OpticalRay, r: float, n1: float, n2: float) - OpticalRay: r: 球面曲率半径 n1/n2: 入射/折射介质折射率 返回折射后的光线对象 # 计算入射角I sin_I (ray.L - r) / r * np.sin(ray.U) I np.arcsin(sin_I) # 计算折射角I sin_I_prime n1 / n2 * sin_I I_prime np.arcsin(sin_I_prime) # 计算像方孔径角U U_prime ray.U I - I_prime # 计算像方截距L L_prime r r * sin_I_prime / np.sin(U_prime) return OpticalRay(L_prime, U_prime)2.3 近轴光与远轴光对比通过参数控制实现两种计算模式def trace_ray(ray, r, n1, n2, paraxialFalse): if paraxial: # 近轴近似 u ray.U i (ray.L - r)/r * u i_prime n1/n2 * i u_prime u i - i_prime l_prime r r*i_prime/u_prime return OpticalRay(l_prime, u_prime) else: # 实际光线 return refract_sphere(ray, r, n1, n2)3. 可视化与像差分析3.1 光线路径绘制函数def plot_ray_tracing(rays_in, rays_out, r, axNone): if ax is None: fig, ax plt.subplots(figsize(10,6)) # 绘制折射球面 theta np.linspace(-np.pi/3, np.pi/3, 100) x r * np.cos(theta) - r # 球心在(r,0) y r * np.sin(theta) ax.plot(x, y, b-, linewidth2) # 绘制入射光线 for ray in rays_in: x_in [ray.L * np.cos(ray.U), 0] y_in [ray.L * np.sin(ray.U), 0] ax.plot(x_in, y_in, r-) # 绘制折射光线 for ray in rays_out: x_out [0, ray.L * np.cos(ray.U)] y_out [0, ray.L * np.sin(ray.U)] ax.plot(x_out, y_out, g--) ax.set_aspect(equal) ax.grid(True) return ax3.2 球差现象演示创建不同孔径角的入射光线并观察成像差异# 参数设置 r 100 # 曲率半径(mm) n1, n2 1.0, 1.5 # 空气到玻璃 angles [5, 10, 15, 20] # 不同入射角度(度) # 光线追迹计算 rays_in [OpticalRay(200, np.radians(a)) for a in angles] rays_out [refract_sphere(ray, r, n1, n2) for ray in rays_in] # 可视化 plt.figure(figsize(12,6)) ax plt.gca() plot_ray_tracing(rays_in, rays_out, r, ax) plt.title(球差现象演示(远轴光)) plt.show()4. 完整Jupyter Notebook实现4.1 交互式参数调节使用IPython widgets创建交互界面from ipywidgets import interact, FloatSlider interact( LFloatSlider(min50, max300, step10, value200), UFloatSlider(min0, max30, step1, value10), rFloatSlider(min50, max200, step10, value100), n1FloatSlider(min1.0, max1.5, step0.1, value1.0), n2FloatSlider(min1.0, max2.0, step0.1, value1.5) ) def interactive_tracing(L, U, r, n1, n2): ray_in OpticalRay(L, np.radians(U)) ray_out refract_sphere(ray_in, r, n1, n2) plt.figure(figsize(10,6)) ax plot_ray_tracing([ray_in], [ray_out], r) ax.set_title(f折射结果: L{ray_out.L:.2f}mm, U{np.degrees(ray_out.U):.2f}°) plt.show()4.2 性能优化技巧对于复杂系统可采用以下优化策略向量化计算使用NumPy同时处理多条光线def batch_trace(rays, r, n1, n2): L np.array([ray.L for ray in rays]) U np.array([ray.U for ray in rays]) sin_I (L - r)/r * np.sin(U) I np.arcsin(sin_I) sin_I_prime n1/n2 * sin_I I_prime np.arcsin(sin_I_prime) U_prime U I - I_prime L_prime r r * sin_I_prime / np.sin(U_prime) return [OpticalRay(l,u) for l,u in zip(L_prime, U_prime)]JIT编译使用Numba加速计算from numba import jit jit(nopythonTrue) def refract_kernel(L, U, r, n1, n2): # 与前面相同的计算逻辑 ...5. 进阶应用与扩展5.1 像差系数计算基于追迹结果可计算初级像差def spherical_aberration(rays_in, rays_out): 计算球差 focal_lengths [ray.L / np.tan(ray.U) for ray in rays_out] return max(focal_lengths) - min(focal_lengths)5.2 材料色散处理考虑折射率随波长的变化class Material: def __init__(self, name, dispersion_formula): self.name name self.formula dispersion_formula def n(self, wavelength): # 实现色散公式(如Sellmeier方程) ...在实际项目中这种光线追迹基础模块可以扩展为完整的镜头设计工具。一个常见的坑是符号规则的一致性——曾经因为角度符号处理不当导致整个系统的像差分析结果完全错误。后来通过单元测试验证每个计算步骤才最终定位到问题所在。