生成模拟干涉图

移相算法求解相位，相位解包裹，泽尼克多项式拟合程序

最近在实验室折腾相位测量，发现光干涉数据处理的三个关键环节：移相算法提取相位、相位解包裹操作、泽尼克多项式拟合。这几个步骤环环相扣，实测中经常需要代码实现，咱们今天就手撕几个代码片段来聊聊。

移相算法实战

假设我们拿到四幅移相干涉图（相位差π/2），核心代码其实就五行：

import numpy as np I1 = 100 + 50*np.cos(phase) # 0相位 I2 = 100 + 50*np.cos(phase + np.pi/2) # π/2 I3 = 100 + 50*np.cos(phase + np.pi) # π I4 = 100 + 50*np.cos(phase + 3*np.pi/2) # 3π/2 # 计算包裹相位 calc_phase = np.arctan2(I4 - I2, I1 - I3) # 关键操作

这里有个坑要注意：np.arctan2直接算出[-π, π]的包裹相位，但实际相位可能是连续变化的。比如当真实相位是3π时，这里会被压缩成-π，这时候就得靠相位解包裹来展开。

相位解包裹的暴力美学

最简单的行展开法代码长这样：

def unwrap_row(phase): offset = 0 unwrapped = np.zeros_like(phase) for i in range(1, len(phase)): delta = phase[i] - phase[i-1] if delta > np.pi: offset -= 2*np.pi elif delta < -np.pi: offset += 2*np.pi unwrapped[i] = phase[i] + offset return unwrapped

这方法在噪声大的地方容易翻车，实测中发现用skimage的unwrap_phase更靠谱：

from skimage import restoration unwrapped_phase = restoration.unwrap_phase(wrapped_phase)

但自己实现的版本有助于理解原理——本质上就是在相位跳变超过π时加减2π的补偿操作。

泽尼克多项式拟合实战

测完相位数据后，我们经常要分解波前像差。泽尼克系数计算可以用现成轮子：

from zernike import RZern cart = RZern(6) # 6阶多项式 cart.fit(phase_data) # 拟合数据 coef = cart.coef # 获取系数

想自己造轮子的话，得先实现泽尼克基函数。比如径向多项式部分：

def zernike_radial(n, m, rho): if (n - m) % 2 != 0: return 0 total = 0 for k in range((n - m)//2 + 1): num = (-1)**k * np.math.factorial(n - k) den = (np.math.factorial(k) * np.math.factorial((n + m)//2 - k) * np.math.factorial((n - m)//2 - k)) total += num / den * rho**(n - 2*k) return total

这个递归实现效率不高，但胜在直观。实际项目里建议用矩阵最小二乘法求解系数，避免重复计算基函数。

调试时记得验证拟合效果——把拟合后的波前与原始数据相减，残差应该接近随机噪声。遇到过拟合的情况就把泽尼克阶数调低，实测发现前15项通常能搞定大部分光学系统的像差分析。