% beam propagation along REFL path

RM1_Rc = 1.7 ;

% given as distances from RM1
zs = [-1162 -314 1377 511 3460+94 3460+466.5 3460+788.5 3460+1092.5] * 1e-3 ;
wx = [4038.9 3950.6 2304.8 1650.9 670.26 793.73 1484.15 2253.78] * 1e-6 / 2 ;
wy = [4206.2 4315.0 2335.9 1805.1 778.95 711.29 1387.24 2119.24] * 1e-6 / 2 ;

reflMeas = beamPath ;

reflMeas.addComponent(component.curvedMirror(RM1_Rc, 0, 'RM1'))
reflMeas.addComponent(component.curvedMirror(-0.6, 856e-3, 'RM2'))
reflMeas.addComponent(component.curvedMirror(RM1_Rc, 1751e-3, 'M5'))
reflMeas.addComponent(component.lens(333.6e-3, 2977e-3, 'L101'))
reflMeas.addComponent(component.lens(-166.8e-3, 3180e-3, 'L102'))

% qin = beamq ;
% qin.q = -1 + 1i ;
% reflMeas.seedq = qin ;
% reflMeas.seedz = -1500e-3 ;

qH1X_RM1_pred = 2.8071 + 13.3724i; % Lisa/Paul prediction 
qH1Y_RM1_pred = 2.5071 + 13.0988i; % Lisa/Paul prediction 

reflMeasX = reflMeas.duplicate ;
reflMeasY = reflMeas.duplicate ;

reflMeasX.seedq.q = qH1X_RM1_pred ;
reflMeasX.seedz = -1e-6 ;
reflMeasY.seedq.q = qH1Y_RM1_pred ;
reflMeasY.seedz = -1e-6 ;

reflXfit = reflMeasX.fitBeamWidth(zs,wx) ;
reflYfit = reflMeasY.fitBeamWidth(zs,wy) ;

zdomain = -10:0.01:5 ;
figure(2)
subplot(2,1,1)
hold on
reflXfit.plotBeamWidth(zdomain, 'r')
plot(zs,wx,'ro')
reflXfit.plotWaists(zdomain, 'r')
reflYfit.plotBeamWidth(zdomain, 'b')
reflYfit.plotWaists(zdomain, 'b')
plot(zs,wy,'bx')
ylabel('Beam Width [m]')
%xlim([-0.05,1.2])
title(['refl Path, red=horiz, blue=vert, RM1 =', num2str(RM1_Rc), 'm'])
reflMeasX.plotComponents(zdomain)
hold off
subplot(2,1,2)
hold on
reflXfit.plotGouyPhase(zdomain, 'r')
reflYfit.plotGouyPhase(zdomain, 'b')
reflMeasX.plotComponents(zdomain)
ylabel('Gouy Phase [deg]')
xlabel('axial dimension, z [m]')
%xlim([-0.05,1.2])
hold off


% figure(3)
% subplot(2,1,1)
% hold on
% reflYfit.plotBeamWidth(zdomain, 'r')
% reflYfit.plotWaists(zdomain)
% plot(zs,wy,'ko')
% ylabel('Beam Width [m]')
% title('refl Path, vertical')
% %xlim([-0.05,1.2])
% reflMeasY.plotComponents(zdomain)
% hold off
% subplot(2,1,2)
% hold on
% reflYfit.plotGouyPhase(zdomain, 'r')
% reflMeasY.plotComponents(zdomain)
% ylabel('Gouy Phase [deg]')
% xlabel('axial dimension, z [m]')
% %xlim([-0.05,1.2])
% hold off


%%


% beam size at PRM