I'm currently working with wavefront and psf.
I'm trying to convert a wavefront map extracted from opticstudio to a psf using a rather complicated but accurate numerical method, not a linear approximation method, but I'm not getting the right results.
I am attaching the matlab code I am currently using.
Please let me know if there are any fixes.
Thanks.
CODE:
function psf_out=wavefront2psf(wvf, center_wvf, data_space_wvf, lambda, film_Dist, rect, data_space_psf)
% wvf: wavefront map of n x n size from opticsutdio
% center_wvf: center coordinate of the wavefront map
% data_space_wvf: data space of the wavefront map
% lambda: the wavelength
% film_Dist: distance between film and the last lens
% rect : point spread are estimated within the rectangle
% data_space_psf: data space of the point spread function
dr=1;
min_dA=2;
cx=center_wvf(1);
cy=center_wvf(2);
rr=[];
th=[];
RR=0;
TH=0;
AA=pi;
[h_wf, ~]=size(wvf);
for p=1:dr:((h_wf/2)-3)
r=p;
r1=r;
r2=r+dr;
area_dR=pi*(r2^2-r1^2);
n_th=ceil(area_dR/min_dA);
n_th=ceil(n_th/4)*4;
if mod(n_th,2)==1
n_th=n_th+1;
end
ths=linspace(0, 2*pi, n_th+1);
ths=ths(1:end-1);
dth=ths(2)-ths(1);
dA=dth/2*(r2^2-r1^2);
r_=(r1+r2)/2*ones(size(ths));
a_=dA*ones(size(ths));
RR=[RR r_];
TH=[TH ths];
AA=[AA a_];
end
[xc, yc]=pol2cart(TH, RR);
xc=xc+cx;
yc=yc+cy;
wsf_rth=interp2(wvf, xc, yc);
x_=(xc-cx)*data_space_wvf;
y_=(yc-cy)*data_space_wvf;
hh=[rect(1):data_space_psf:rect(2)];
ww=[rect(3):data_space_psf:rect(4)];
phase_sum=zeros(size(hh,2), size(ww,2));
[h_psf, w_psf]=size(phase_sum);
for p=1:h_psf
for q=1:w_psf
px=(p-h_psf/2)*data_space_psf;
py=(q-w_psf/2)*data_space_psf;
ipx=x_-px;
ipy=y_-py;
path=sqrt(ipx.^2+ipy.^2+(film_Dist).^2);
phase=(wsf_rth)+(path-film_Dist)/(lambda);
phase_sum(p,q)=sum(AA.*exp(i*phase(:)'));
end
end
psf_out=sqrt(abs(phase_sum).^2);
psf_out=psf_out/sum(psf_out(:));
end
[Mod note: moved to more appropriate forum for ZOS-API-related discussions.]