求孔径 的PSF函数。前文给过一个了,这里再给一份代码用于求PSF。
clc; clear all; close all
% 调用第一篇里的函数了。
Znm = Zernike(1,128) * 1e-6; % 单位换算到 微米
% Znm(Znm~=0) = 1 * 1e-6;
figure(10);
imagesc(Znm);colormap(jet);colorbar
lambda = 1.22 * 1e-6;
k = 2*pi / lambda;
D = 1.54; % m
FoV = 10; % 视场大小,单位应该是 :秒
THld = lambda/D * 206265; % 衍射极限,单位rad换算到 "
dth = THld / 3;
Znm = exp(1j * k * Znm); % 复数下的波前
figure(1); subplot(131);imshow(Znm,[]);colormap(jet);colorbar
FFTSize = 1024; % PSF图像的分辨率
halo = padarray(Znm,[FFTSize-size(Znm,1),FFTSize-size(Znm,1)],'post');
subplot(132);imshow(halo,[]);colormap(jet);colorbar % 只能用imshow画图
% 比较这两张图,认识到padarry的用法
xx = 1 - size(Znm,1)/2; % Znm的中心坐标
yy = FFTSize/2 - 1; % PSF的中心坐标
% 可以单独运行 figure(11);imagesc(circshift(halo,[xx xx])), 弄清楚其用法
halo = circshift(fft2(circshift(halo,[xx xx])),[yy yy]);
sz = zeros(1,2);
cen = sz;
sz = size(halo);
cen = (sz+2-mod(sz,2))/2;
dk = 2*pi./(FFTSize .* 0.01); % g.spacing_=0.01
kx = ((1:sz(1))-cen(1))*dk;
ky = ((1:sz(2))-cen(2))*dk;
dTH = dk/k*206265; % k 方向上的波矢?
thx_ = kx/k*206265; % x 方向上的波矢?
thy_ = ky/k*206265; % y 方向上的波矢?
SELx = abs(thx_)<=(FoV+2*dTH);
SELy = abs(thy_)<=(FoV+2*dTH);
thx_ = thx_(SELx);
thy_ = thy_(SELy);
[THX_,THY_] = meshgrid(thx_,thy_);
Nth = ceil(FoV/dth);
thx = (-Nth:Nth)*dth;
thy = thx;
[THX,THY] = meshgrid(thx,thy);
halo = halo(SELx,SELy);
HALO = qinterp2(THX_,THY_,halo,THX,THY);
PSF = abs(HALO).^2;
PSF(isnan(PSF)) = 0;
PSFmax = max(PSF(:)); % figure(104); imagesc(PSF0)
PSF0 = PSF/PSFmax; % make the brightest value =1.
figure(3); imagesc(thx,thy,log10(PSF0),[-4 0]);
daspect([1 1 1]);colormap(jet);axis xy;
运行后的结果如下
对于孔径为(这里没有给出这样孔径的代码):
其PSF为: