summaryrefslogtreecommitdiff
path: root/examples/unstable_cavity.m
diff options
context:
space:
mode:
Diffstat (limited to 'examples/unstable_cavity.m')
-rw-r--r--examples/unstable_cavity.m151
1 files changed, 151 insertions, 0 deletions
diff --git a/examples/unstable_cavity.m b/examples/unstable_cavity.m
new file mode 100644
index 0000000..cf7afb6
--- /dev/null
+++ b/examples/unstable_cavity.m
@@ -0,0 +1,151 @@
+% Illustrates the use of cav.m by modeling a marginally stable optical cavity.
+%------------------------------------------------------------------------------------
+% PROGRAM NAME: unstable_cavity
+% AUTHOR: Andri M. Gretarsson
+% DATE: Aug. 24, 2004
+%
+% SYNTAX: unstable_cavity
+%
+% Uses cav.m to model the fields inside and outside a marginally unstable linear cavity
+% The example is based on the LIGO Livingston recycling cavity in its cold state.
+% In the 1D mode, the output is a set of 9 graphs showing the
+% the field inside the cavity and in reflection as a function of length around a resonance.
+% Due to the marginally stable nature of the cavity in this example, we get interesing
+% structure in the fields. The program runs as a script rather than a function and
+% numerous parameters can be set near the beginning of the script. These are:
+%
+% R1, R2 - the radii of curvature of the input and end cavity mirrors
+% r1, r2 - the amplitude reflection coefficients of the two mirrors
+% l1, l2 - the amplitude loss coefficients of the two mirrors
+% L0 - the cavity length
+% lambda - the optical wavelength
+% l_in, m_in - Hermite-Gaussian mode numbers for the input light (referred to
+% the basis set by qin (see next line).
+% qin - the complex radius of curvature of the input light
+% npts - for a 1D calculation, the number of points in the calculation, for a 2D
+% calculation the number of grid points is npts^2
+% twoD - Whether the calculation should be done in 1D (cross section) or 2D. The 2D
+% form is significantly slower. To improve running speeds it may be a good
+% idea to reduce npts. NOTE: It will help to maximize the figure window and
+% rerun the program to improve visibility in the 1D mode.
+%
+% INPUT ARGUMENTS: none
+% OUTPUT ARGUMENTS: none
+% Last Updated: June 30, 2007 by AMG
+%
+%------------------------------------------------------------------------------------
+% SYNTAX: unstable_cavity
+%------------------------------------------------------------------------------------
+
+% Mirror properties
+Ritmx=-14.76e3; % ITMx radius of curvature [meters]
+Ritmy=-14.52e3; % ITMy radius of curvature [meters]
+Rrm=20.78e3; % RM radius of curvature [meters]
+R1=Rrm;
+R2=(Ritmx+Ritmy)/2/1.45; % 1.45 is due to the index of refraction
+r1=sqrt(0.9725); % Measured reflectance coefficients assuming 150ppm loss in ITMs.
+r2=sqrt(0.9731); % see: http://www.ligo-la.caltech.edu/ilog/pub/ilog.cgi?group=detector&date_to_view=08/19/2003&anchor_to_scroll_to=2003:08:25:10:47:28-rana
+%r1=sqrt(1-0.027); % nominal from COC (Helena Armandula) webpage
+%r2=sqrt(1-0.0288);
+%r1=0.96; % test
+l1=0;%sqrt(150e-6); % ITM loss (nominal)
+l2=0;%sqrt(150e-6); % ITM loss (nominal)
+L0=9.2;
+
+% Incident light
+lambda=1.064e-6;
+l_in=0;
+m_in=0;
+[qin,pin]=prop(q_(0.03875,7.1e3,lambda),free(L0),[l_in,m_in]); % Input beam is matched to arms so use arm cavity beam at ITMx to get q value.
+%[qin,pin]=prop(q_(0.03685,-6.9564e3,lambda),free(L0),[l_in,m_in]);
+win=w_(qin,lambda);
+
+% Plot domain
+twoD=0;
+if twoD
+ npts=9;
+ xseed=[-2*win:win/npts:2*win];
+ yseed=xseed';
+ [x,y]=meshgrid(xseed,yseed);
+ thecolormap='bone';
+else
+ npts=200;
+ xseed=[-2.5*win:win/npts:2.5*win];
+ x=xseed;
+ y=zeros(size(x));
+end
+
+
+zin=HermiteGaussianE([l_in,m_in,qin,lambda(1),pin],x,y);
+
+
+n_bounces=70;
+n_lengths=9;
+lowerdL=-15e-9; upperdL=15e-9;
+deltaL=[lowerdL:(upperdL-lowerdL)/(n_lengths-1):upperdL];
+Pcoherent_out=zeros(n_lengths,1);
+for s=1:n_lengths
+ L=round(L0/lambda)*lambda+deltaL(s); % Cavity length
+
+ [qrefl,qtrans,qcav,prefl,ptrans,pcav]=...
+ cav(qin,l_in,m_in,lambda,L,R1,R2,r1,r2,l1,l2,n_bounces,[1.46,1,1,],[0,0]);
+
+ zrefl=HermiteGaussianE([l_in * ones(size(qrefl)),m_in * ones(size(qrefl)),...
+ qrefl,lambda(1) * ones(size(qrefl)),prefl],x,y);
+ zcav=HermiteGaussianE([l_in * ones(size(qcav)),m_in * ones(size(qcav)),...
+ qcav,lambda(1) * ones(size(qcav)),pcav],x,y);
+ ztrans=HermiteGaussianE([l_in * ones(size(qtrans)),m_in * ones(size(qtrans)),...
+ qtrans,lambda(1) * ones(size(qtrans)),ptrans],x,y);
+ if length(size(zcav))==3
+ zsumrefl=sum(zrefl,3);
+ zsumcav=sum(zcav,3);
+ zsumtrans=sum(ztrans,3);
+ else
+ zsumrefl=sum(zrefl,2);
+ zsumcav=sum(zcav,2);
+ zsumtrans=sum(ztrans,2);
+ end
+
+ if ~twoD
+ interpdom=[min(xseed):(max(xseed)-min(xseed))/999:max(xseed)];
+ dx_interp=interpdom(2)-interpdom(1);
+ if s==1
+ inputintens_interp=interp1(x,abs(zin).^2,interpdom,'spline');
+ inputpower=sum(inputintens_interp*dx_interp.*abs(interpdom)*pi);
+ end
+ cavintens_interp=interp1(x,abs(zsumcav).^2,interpdom,'spline');
+ cavpower=sum(cavintens_interp*dx_interp.*abs(interpdom)*pi);
+ reflintens_interp=interp1(x,abs(zsumrefl).^2,interpdom,'spline');
+ reflpower=sum(reflintens_interp*dx_interp.*abs(interpdom)*pi);
+ transintens_interp=interp1(x,abs(zsumtrans).^2,interpdom,'spline');
+ transpower=sum(transintens_interp*dx_interp.*abs(interpdom)*pi);
+ Pcoherent_out(s)=reflpower+transpower;
+ end
+
+ if twoD
+ figure(1);
+ subplot(3,3,s); h=pcolor(x,y,abs(zsumcav)); % h=pcolor(x,y,abs(zsumcav).*cos(angle(zsumcav)));
+ set(h,'EdgeColor','none'); axis square;
+ colormap(thecolormap); %colorbar;
+ set(h,'EdgeColor','none'); set(h,'FaceColor','interp');
+ set(gca,'Visible','off'); set(gcf,'Color','black');
+ set(gca,'Zlim',[-35 35]);
+ drawnow;
+ h=figtext3D(3,1,10,['\Delta{L}: ',num2str(deltaL(s)*1e9,'%0.3g'),' nm'],10);
+ set(h,'color','white');
+ shg;
+ else
+ figure(1); subplot(3,3,s);
+ set(gcf,'Color','default');
+ plot(x*100,zin.*conj(zin),'m--','linewidth',2); hold on;
+ h=plot(x*100,zsumrefl.*conj(zsumrefl),'-',x*100,zsumcav.*conj(zsumcav),'-');
+ set(h,'linewidth',2);
+ legend('input','refl','cav');
+ drawnow; hold off;
+ xlabel('cm'); ylabel('intensity');
+ figtext(0.5,9,['\Delta{L}: ',num2str(deltaL(s)*1e9,'%0.3g'),' nm'],8);
+ figtext(0.5,7.5,['P_{cav}/P_{in} = ',num2str(cavpower,'%0.2g')],8);
+ shg;
+ end
+
+end \ No newline at end of file