INCPSF Coherent and Incoherent Transfer Functions and Point-Spread
Functions
SYNTAX: incpsf;
by Chuck DiMarzio Northeastern University December 2009
Requires: FFTAXISSHIFT,FFTAXIS
!! This file may be copied, used, or modified for educational and !! research purposes provided that this header information is not !! removed or altered, and provided that the book is cited in !! publications, as DiMarzio, Charles A., Optics for Engineers, !! CRC Press, Boca Raton, FL, 2011. !! http://www.crcpress.com !! Other distribution is prohibited without permission.
Contents
Setup
xax=fftaxisshift(0:1023); % set up axes with zero center [x,y]=meshgrid(xax,xax); % generate coordinates rsq=x.^2+y.^2; % Compute radius squared fax=fftaxisshift(fftaxis(xax)); % generate frequency axis [fx,fy]=meshgrid(fax,fax); % generate frequency coordinates fsq=fx.^2+fy.^2; % Magnitude of frequency squared ftop=max(max(fx));fmax=ftop/sqrt(length(fx)); % for mask size xmax=length(xax)/2/sqrt(length(fx));
Compute the transfer functions
mask1=single(fsq<(fmax^2)); % Low Pass filter for coherent transfer % function. %mask1a=fftshift(ifft2(abs((fft2(fftshift(mask1)))).^2)); % The incoherent transfer function is the autocorrelation function of the % coherent transfer function.
Compute the point--spread functions
psf1=fftshift(ifft2(ifftshift(mask1)));
%psf1a=fftshift(ifft2(ifftshift(mask1a)));
psf1a=abs(psf1).^2;
mask1a=fftshift(ifft2(fftshift(psf1a)));
Show a slice through the Transfer Functions
fig1=figure;plot(fax,real(mask1(:,1+floor(length(fx)/2))/max(max(mask1))),... fax,real(mask1a(:,1+floor(length(fx)/2)))/max(max(mask1a))); grid on; moose=axis;axis([fmax*4*[-1,1],moose([3,4])]); xlabel('f_x, Spatial Freq., pixel^{-1}'); ylabel('Transfer Function(f_x,0)'); % Check that the imaginary part is zero figure;plot(fax,imag(mask1(:,1+floor(length(fax)/2))/max(max(mask1a))),... fax,imag(mask1a(:,1+floor(length(fax)/2)))/max(max(mask1a))); moose=axis;axis([fmax*4*[-1,1],moose([3,4])]); grid on; xlabel('f_x, Spatial Freq., pixel^{-1}'); ylabel('Transfer Function'); title('Sanity Check: Imaginary parts zero');


Show the transfer functions with imagesc
fig2=figure;imagesc(fax,fax,mask1);axis image;colormap(flipud(gray)); axis(fmax*4*[-1,1,-1,1]); xlabel('f_x, Spatial Freq., pixel^{-1}'); ylabel('f_y, Spatial Freq., pixel^{-1}'); fig3=figure;imagesc(fax,fax,real(mask1a));axis image;colormap(flipud(gray)); axis(fmax*4*[-1,1,-1,1]); xlabel('f_x, Spatial Freq., pixel^{-1}'); ylabel('f_y, Spatial Freq., pixel^{-1}');


Show the point--Spread functions with imagesc
fig4=figure;imagesc(xax,xax,20*log10(abs(psf1)));axis image; axis(xmax*8*[-1,1,-1,1]); moose=caxis;caxis(moose(2)+[-30,0]);colormap(flipud(gray));colorbar; xlabel('x, Position, Pixels'); ylabel('y, Position, Pixels'); % Check angle to make sure psf is real figure;imagesc(xax,xax,angle(psf1));colormap('hsv');axis image;colorbar; axis(xmax*8*[-1,1,-1,1]); xlabel('x, Position, Pixels'); ylabel('y, Position, Pixels'); title('Sanity Check: Angle zero or 180'); % Now the incoherent one fig5=figure;imagesc(xax,xax,10*log10(abs(psf1a)));axis image; axis(xmax*8*[-1,1,-1,1]); xlabel('x, Position, Pixels'); ylabel('y, Position, Pixels'); moose=caxis;caxis(moose(2)+[-30,0]);colormap(flipud(gray));colorbar; % Check angle to make sure psf is real figure;imagesc(xax,xax,angle(psf1a));colormap('hsv');axis image;colorbar; axis(xmax*8*[-1,1,-1,1]); xlabel('x, Position, Pixels'); ylabel('y, Position, Pixels'); title('Sanity Check: Angle zero or 180');




Plot slices of the point--spread functions
fig6=figure;plot(xax,real(psf1(:,1+floor(length(fax)/2))/max(max(psf1))),... xax,real(psf1a(:,1+floor(length(fax)/2)))/max(max(psf1a))); moose=axis;axis([xmax*8*[-1,1],moose([3,4])]); grid on; xlabel('x, Position, Pixels'); ylabel('PSF(x,0)');
