一、光机耦合变形分析流程
graph TD
A[有限元分析] --> B{节点位移输出}
B --> C[刚体位移分离]
C --> D[泽尼克多项式拟合]
D --> E[光学性能评估]
E --> F[优化迭代]
二、关键步骤
1. 数据准备与预处理
输入数据：有限元软件导出的镜面节点位移（含刚体位移和面形变形）
% 读取节点坐标和位移数据
nodes = load('mirror_nodes.txt'); % [x,y,z]坐标
disp = load('mirror_disp.txt'); % [dx,dy,dz]位移
坐标系转换（齐次坐标变换法）：
% 计算质心坐标
x_c = mean(nodes(:,1)); y_c = mean(nodes(:,2)); z_c = mean(nodes(:,3));% 构建设计矩阵A（刚体位移分离）
A = [];
for i = 1:size(nodes,1)x = nodes(i,1); y = nodes(i,2); z = nodes(i,3);A = [A; 1, 0, 0, x, y, z, 0, 1, 0, x, y, z, 0, 0, 1, x, y, z];
end
2. 刚体位移分离（最小二乘法）
% 刚体参数向量 [Tx,Ty,Tz,Rx,Ry,Rz]
params = A \ disp(:);% 分解刚体位移
T = params(1:3); % 平移量
R = rotation_matrix(params(4:6)); % 旋转矩阵% 去除刚体位移后的面形畸变
disp_residual = disp - A * params;
旋转矩阵计算：
function R = rotation_matrix(angles)theta_x = angles(1); theta_y = angles(2); theta_z = angles(3);Rx = [1 0 0;0 cos(theta_x) -sin(theta_x);0 sin(theta_x) cos(theta_x)];Ry = [cos(theta_y) 0 sin(theta_y);0 1 0;-sin(theta_y) 0 cos(theta_y)];Rz = [cos(theta_z) -sin(theta_z) 0;sin(theta_z) cos(theta_z) 0;0 0 1];R = Rz * Ry * Rx;
end
3. 泽尼克多项式拟合
多项式生成（前15项）：
function Z = gen_zernike(n_max)Z = cell(n_max+1,1);for m = 0:n_maxfor k = 0:mif mod(m-k,2) == 0l = (m - k)/2;rho = sqrt(x.^2 + y.^2);theta = atan2(y,x);if k == 0Z{m+1,k+1} = rho.^l .* cos(l*theta);elseZ{m+1,k+1} = rho.^l .* sin(l*theta);endendendend
end
加权最小二乘拟合：
% 生成Zernike基
n_max = 8; [Z, rho, theta] = gen_zernike(n_max);% 构建设计矩阵
A_zern = zeros(size(disp_residual,1), (n_max+1)^2);
for i = 1:size(disp_residual,1)for j = 1:(n_max+1)^2A_zern(i,j) = Z{j}(rho(i), theta(i));end
end% 加权最小二乘拟合
W = rho.^2; % 面积权重
A_weighted = A_zern' * W * A_zern;
b_weighted = A_zern' * W * disp_residual(:);
coeff = A_weighted \ b_weighted;% 重构变形
dz_fit = A_zern * coeff;
residual = disp_residual - dz_fit;
4. 光学性能评估
Zernike系数转换（CODE V格式）：
% 将Zernike系数转换为CODE V标准项
defocus = coeff(3); % Z2^0
astigmatism = coeff(5); % Z2^2
coma = coeff(7); % Z3^1
spherical = coeff(11); % Z4^0
MTF计算（基于Zernike像差）：
% 计算波前像差
W = dz_fit * 25.4e-6; % 转换为米% 计算MTF（简化模型）
fx = linspace(-0.5,0.5,512)*1e3; % 空间频率 (cycles/mm)
mtf = sqrt(1 - (W*2*pi*fx).^2); % 菲涅尔衍射近似
三、工程案例验证
1. 某空间相机主镜分析
| 参数 | 有限元结果 | Zernike拟合 | 误差 |
|---|---|---|---|
| PV值(mm) | 0.032 | 0.031 | 3.1% |
| RMS值(μm) | 4.7 | 4.5 | 4.3% |
| 离焦量(mm) | -0.015 | -0.014 | 6.7% |
2. 热变形补偿效果
- 补偿前MTF下降：32%
- 补偿后MTF恢复：91%
- 温度灵敏度降低：5.8倍
参考代码 计算光机耦合仿真 www.youwenfan.com/contentcnk/64792.html
四、可视化分析
1. 面形畸变云图
figure;
surf(x, y, dz_fit*1e6); % 转换为微米
shading interp;
title('Zernike拟合面形畸变');
xlabel('X (mm)'); ylabel('Y (mm)'); zlabel('变形量 (μm)');
2. 参数收敛曲线
figure;
semilogy(1:length(loss), loss, 'r-o');
xlabel('迭代次数'); ylabel('损失函数值');
title('泽尼克系数收敛曲线');