%
%
% Optical Imaging and Spectroscopy
%
% David J. Brady
% Duke University
% www.opticalimaging.org
%
% Figures 7.3, 7.4 and 7.5
%
%
% interesting otf's
%
% units are  inverse Delta
%
apertureRange=linspace(-2,2,128);
[x,y]=meshgrid(apertureRange,apertureRange);
figure(1);set(gcf,'color','white');
colormap(1-gray(256));
fNumber=1;
Delta=2;
field=sqrt(x.^2+y.^2);
h=P(field*fNumber);
H=conv2(h,h,'same');
H=H/max(max(H));
H1=sinc(Delta*x).*sinc(Delta*y).*H;
subplot(3,1,1);
plot(Delta*apertureRange(65:128),abs(H1(64,65:128)),'-k');
axis([0 3 0 1]);
xlabel('1/\Delta');title('(a) \Delta= 2 \lambda f/#');
Delta=0.5;
H1=sinc(Delta*x).*sinc(Delta*y).*H;
subplot(3,1,3);
plot(Delta*apertureRange(65:128),abs(H1(64,65:128)),'-k');
axis([0 3 0 1]);
xlabel('1/\Delta');title('(c) \Delta = 0.5 \lambda f/#');
Delta=1;
H1=sinc(Delta*x).*sinc(Delta*y).*H;
subplot(3,1,2);
plot(Delta*apertureRange(65:128),abs(H1(64,65:128)),'-k');
axis([0 3 0 1]);
xlabel('1/\Delta');title('(b) \Delta =  \lambda f/#');
%
%
% same stuff in units of inverse wavelengths
%
figure(2);set(gcf,'color','white');
colormap(1-gray(256));
fNumber=1;
Delta=2;
field=sqrt(x.^2+y.^2);
h=P(field*fNumber);
H=conv2(h,h,'same');
H=H/max(max(H));
H1=sinc(Delta*x).*sinc(Delta*y).*H;
subplot(4,1,1);
plot(apertureRange(65:128),abs(H1(64,65:128)),'-k');
axis([0 2 0 1]);
xlabel('1/\lambda f/#');title('(a) \Delta= 2 \lambda f/#');
Delta=0.5;
H1=sinc(Delta*x).*sinc(Delta*y).*H;
subplot(4,1,3);
plot(apertureRange(65:128),abs(H1(64,65:128)),'-k');
axis([0 2 0 1]);
xlabel('1/\lambda f/#');title('(c) \Delta = 0.5 \lambda f/#');
Delta=1;
H2=sinc(Delta*x).*sinc(Delta*y).*H;
subplot(4,1,2);
plot(apertureRange(65:128),abs(H2(64,65:128)),'-k');
axis([0 2 0 1]);
xlabel('1/\lambda f/#');title('(b) \Delta =  \lambda f/#');
subplot(4,1,4);
plot(apertureRange(65:128),abs(H1(64,65:128))-abs(H2(64,65:128)),'-k');
axis([0 2 0 .2]);
xlabel('1/\lambda f/#');title('(d)=(c)-(b)');
figure(3);set(gcf,'color','white');
colormap 'gray';
subplot(2,2,1);hp=surf(apertureRange,apertureRange,H);
shading interp;
light('Position',[-2,2,20]);
lighting phong;
material([0.4,0.6,0.5,30]);
set(hp,'FaceColor',[.5,.5,.5],'BackFaceLighting','lit');
xlabel('1/\lambda f/#');ylabel('1/\lambda f/#');
title('(a) MTF');
subplot(2,2,2);plot(apertureRange,H(:,64),'-k');
xlabel('1/\lambda f/#');
title('(b) MTF cross section');
ptfScale=2;
PTF=sinc(ptfScale*x).*sinc(ptfScale*y);
subplot(2,2,3);hp=surf(ptfScale*apertureRange,ptfScale*apertureRange,PTF);
shading interp;
light('Position',[-2,2,20]);
lighting phong;
material([0.4,0.6,0.5,30]);
set(hp,'FaceColor',[.5,.5,.5],'BackFaceLighting','lit');
xlabel('1/\Delta');ylabel('1/\Delta');
title('(a) PTF');
subplot(2,2,4);plot(ptfScale*apertureRange,PTF(:,64),'-k');
xlabel('1/\Delta');
title('(b) PTF cross section');