HGEXAMPLE Hermite Gaussian Examples
by Chuck DiMarzio Northeastern University December 2008
!! 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
First, let's do a circular aperture
gamma=1; % % radius=2; xmax=5; dx=0.02; ymax=xmax; dy=dx; mmax=20;nmax=20; xaxis=(-xmax:dx:xmax); yaxis=(-ymax:dy:ymax); [x,y]=meshgrid(xaxis,yaxis); rsq=x.^2+y.^2; u=(rsq<radius^2)*1/pi/radius^2; % Show source pattern fig26=figure;imagesc(xaxis,yaxis,abs(u).^2); axis image;colormap(flipud(gray));colorbar; %print('-deps',[epspath,'9-26-hgex.eps']); % % First, find the best w % wlist=(0.7:0.02:1.2)*radius; c00=zeros(1,length(wlist)); for nw=1:length(wlist); h=hgm(0,0,x,y,wlist(nw),wlist(nw)); c00(nw)=sum(sum(u.*conj(h)))*dx*dy; end; % Show optimization of w fig27=figure;plot(wlist,c00,'k');grid on; xlabel('w, Gaussian radius');ylabel('C_{00}, Coefficient'); w=wlist(find(c00==max(c00))); c=zeros(mmax,nmax); for m=0:2:mmax; for n=0:2:nmax; h=hgm(m,n,x,y,w,w); c(m+1,n+1)=sum(sum(u.*conj(h)))*dx*dy; % m+1 because subscripts % start at 1. c(1,1) is c00 end; end;


Show the coefficients
fig28=figure;imagesc(20*log10(abs(c)));axis image;colormap(flipud(gray));
caxis([-50,0]);colorbar;

Approximation
yslice=floor(length(xaxis)/2); ua=zeros(size(x)); usave=zeros(size(x,1),mmax+1); for m=0:2:mmax; for n=0:2:nmax; h=hgm(m,n,x,y,w,w); ua=ua+c(m+1,n+1)*h; % remember the index starts at 1 end; usave(:,m+1)=ua(:,yslice); end;
Image the original source
fig29=figure;imagesc(xaxis,yaxis,abs(ua).^2);
axis image;colormap(flipud(gray));colorbar;

Plot alternate even orders (1:4:end)
fig30=figure;plot(xaxis,20*log10(abs(usave(:,1:4:end))),'k');grid on; xlabel('x, Radial Distance');ylabel('Irradiance, dB'); moose=axis;axis([moose(1:2),moose(4)+[-60,0]]);

Far field;
newxmax=10; newdx=0.1; newymax=newxmax; newdy=newdx; newxaxis=(-newxmax:newdx:newxmax); newyaxis=(-newymax:newdy:newymax); [newx,newy]=meshgrid(newxaxis,newyaxis); yslice=floor(length(newxaxis)/2); uaf=zeros(size(newx)); usavef=zeros(size(newx,1),mmax+1); for m=0:2:mmax; for n=0:2:nmax; h=hgm(m,n,newx,newy,w,w); uaf=uaf+c(m+1,n+1)*exp(1i*pi/2*(1+m+n))*h; % remember the index starts at 1 end; usavef(:,m+1)=uaf(:,yslice); end;
plot diffraction pattern in db.
fig31=figure;imagesc(newxaxis,newyaxis,20*log10(abs(uaf)));
caxis([-50,-20]);
axis image;colormap(flipud(gray));colorbar;

Plot alternate even orders (1:4:end)
fig32=figure;plot(newxaxis*pi/w^2,20*log10(abs(usavef(:,1:4:end))),'k'); grid on; moose=axis;axis([moose(1:2),moose(4)+[-60,0]]); xlabel('x, Radial Distance');ylabel('Irradiance, dB');
