白光成像仿真：波前差、采样和插值

转载：白光成像仿真：波前差、采样和插值
代码在最后


这是matlab自带的图片

% 初始化环境  
clc; clear; close all;  
  
% 定义参数  
Ns = 513; % 采样点数  
EFL = 100; % 焦距  
D = 15; % 口径  
R = D / 2; % 半口径  
  
% 入瞳相位直角坐标  
x = linspace(-20, 20, Ns);  
[x, y] = meshgrid(x, x);  
[t, r] = cart2pol(x, y); % 转换为极坐标  
stop = zeros(size(x));  
stop(r <= D / 2) = 1; % 定义光透过率函数  
  
% 初始化最终像面坐标  
ximg_all = linspace(-0.1, 0.1, 1024);  
[ximg_all, yimg_all] = meshgrid(ximg_all, ximg_all);  
Zimg_all = zeros(1024);  
Opth = 63.28E-6 * (8^(1/2) * (3 * (r / R).^3 - 2 * (r / R))) .* sin(t); % 加载慧差光程差  
  
% 波长循环  
for WL = 400E-6:50E-6:700E-6 % 波长范围[mm]  
    WF = stop .* exp(1i * 2 * pi * Opth / WL); % 计算波前函数  
  
    % 补零  
    Nzero = 4;  
    xobi = linspace(-20 * Nzero, 20 * Nzero, Ns * Nzero);  
    [xobi, yobi] = meshgrid(xobi, xobi);  
    Eobi_ori = WF;  
    Eobi = zeros(Nzero * Ns);  
    a = (Nzero * Ns / 2 - Ns / 2 + 0.5);  
    b = (Nzero * Ns / 2 + Ns / 2 - 0.5);  
    Eobi(a:b, a:b) = Eobi_ori; % 将原始波前函数嵌入到补零后的波前函数中  
  
    % FFT2  
    PixelSize = WL * EFL / (Nzero * 20);  
    x_img = [-Nzero * Ns / 2:Nzero * Ns / 2 - 1] * PixelSize;  
    [x_img, y_img] = meshgrid(x_img, x_img);  
    E_img = fftshift(fft2(fftshift(Eobi)));  
    I_img = abs(E_img).^2; % 计算光强  
    I_img = I_img / max(I_img(:)); % 归一化光强  
  
    % 插值到最终像面坐标  
    Iimg_all = interp2(x_img, y_img, I_img, ximg_all, yimg_all);  
    Zimg_all = Zimg_all + Iimg_all;  
end  
  
% 归一化总光强  
Zimg_all = Zimg_all / max(Zimg_all(:));  
  
% 显示像面光强分布  
figure;  
surf(ximg_all, yimg_all, Zimg_all);  
shading interp;  
colormap('jet');  
colorbar;  
view(0, 90);  
xlabel('x [mm]');  
ylabel('y [mm]');  
title('像面光强分布（多波长叠加）');  
drawnow;  
  
% 读取并显示图像  
A = imread('football.jpg');  
A = double(A);  
A = imresize(A, 4);  
Agray = flipud((A(:, :, 1) + A(:, :, 2) + A(:, :, 3)) / 3);  
Agray = Agray / max(Agray(:));  
[aA, bA] = size(Agray);  
yA = [1:aA] / 15000;  
xA = [1:bA] / 15000;  
[xA, yA] = meshgrid(xA, yA);  
figure;  
S = surf(xA, yA, Agray);  
S.EdgeColor = 'none';  
colormap gray;  
colorbar;  
view(0, 90);  
xlabel('x [mm]');  
ylabel('y [mm]');  
title('football');  
  
% 图像反转和插值  
Agray = (Agray * (-1)) + 1;  
z_img_ioe = interp2(xA, yA, Agray, ximg_all, yimg_all);  
z_img_ioe(isnan(z_img_ioe)) = 0;  
figure;  
S = surf(ximg_all, yimg_all, z_img_ioe);  
S.EdgeColor = 'none';  
colormap gray;  
colorbar;  
view(0, 90);  
xlabel("x[mm]");  
ylabel("y[mm]");  
  
% 非相干成像  
B1 = fftshift(fft2(fftshift(Zimg_all)));  
B2 = fftshift(fft2(fftshift(z_img_ioe)));  
C = B1 .* B2;  
Img = fftshift(ifft2(fftshift(C)));  
Img = Img / max(Img(:));  
figure;  
surf(ximg_all, yimg_all, Img);  
shading interp;  
colormap gray;  
colorbar;  
view(0, 90);  
xlabel('x [mm]');  
ylabel('y [mm]');  
title('Non-coherent Imaging Result');