光学成像的MATLAB仿真

clc,clear close all %initial lamda = 500e-9; % 光波长 R = 10e-2; % 光瞳半径 di = 50e-2; % 像距 trans = 6.5e-6; %相机的像元大小 p0 = (R/di/lamda)*trans; %这是根据实验参数计算的截止频率 %read the image I = imread('细胞.jpg'); figure(1),imshow(I); I = double(I); [Nx,Ny] = size(I); p0 = max(fix([Nx,Ny]/2)); % 这里模拟截止频率为图像最大频率的1/2 fX = [0:fix(Nx/2),ceil(Nx/2)-1:-1:1];%linspace(-10e5,10e5,2000); fY = [0:fix(Ny/2),ceil(Ny/2)-1:-1:1];%linspace(-10e5,10e5,2000); [fx,fy] = meshgrid(fY,fX); %OTF p = sqrt(fx.^2+fy.^2); f =@(fx,fy) 2/pi*(acos(p./p0)-p./p0.*sqrt(1-(p/p0).^2)).*circ(fx,fy,p0); H0 = f(fx,fy); mesh(fftshift(H0)) % 显示OTF三维曲面 %simulation F_I = fft2(I); %物强度频谱 G0 = F_I(1,1); R_I = F_I/G0; %归一化 R_i = R_I.*H0*G0; i = ifft2(R_i); figure(2),imshow(uint8(i));