clear classes
close all
clear all 
clc

set(0,'DefaultLineLineWidth',2)
set(0,'DefaultAxesFontSize',12)
set(0,'DefaultAxesXGrid','on')
set(0,'DefaultAxesYGrid','on')
set(0,'DefaultAxesZGrid','on')
legFontSize = 10;

%%
printFigs = 1;

svnRoot = '/Users/kissel/ligo/svncommon/SeiSVN/seismic/';
resultsDir = [svnRoot 'Common/SPI/Results/FiberCollimators/'];
toolsDir = [svnRoot 'Common/MatlabTools/'];

fcsn = 'S0272502';
measDate = '2025-06-03';
author = ', J. Kissel';
measNote = 'Pre-Bake';
figTag = [measDate '_spifc_' fcsn '_prebake'];

%% Settings and raw measurements
w0d = 1.05e-3; % Desired waist size, T2400413
z0d = 0; % "The best, most collimated beam has its waist at the collimator lens"
         % and I set up the beam profiler distances with the collimator
         % lens at z = 0 (+/- 5 mm).
lambda = 1064e-9;

% Measured Data (Use D4sig, mean values "as is," i.e. as beam *diameter* in [um])
zPosition_m = [0.508 0.991 1.499 2.007 3.251 4.496 5.41];
xRadius_m =   [1907  1802  1785  1849  2293  2956  3492]*1e-6/2;
yRadius_m =   [1970  1853  1806  1870  2293  2986  3546]*1e-6/2;

% Set Initial Guess:
measRange_m = [0 20];   % essentially, the desired xlim of the plot
gRadius_m = 1800e-6/2; % Guess based on the data
zGuess_m = 1.5;        % Guess based on the data

%% Fit the data
goo = beamPath;

% Creat a position vector for the model
zModel_m = measRange_m(1):.01:measRange_m(2);

% Initialize the fitting object
goo.seedWaist(gRadius_m,zGuess_m);

% Fit the data
goox = goo.fitBeamWidth(zPosition_m,xRadius_m);
gooy = goo.fitBeamWidth(zPosition_m,yRadius_m);

% Assume only one waist within the measurement range, and find its
% properties
[z0x_m, q0x] = goox.getWaists(measRange_m(1),measRange_m(2));
[z0y_m, q0y] = gooy.getWaists(measRange_m(1),measRange_m(2));
w0x_m = q0x.waistSize;
w0y_m = q0y.waistSize;
zRx_m = pi*(w0x_m.^2)/lambda;
zRy_m = pi*(w0y_m.^2)/lambda;
wzzRx_m = w0x_m * sqrt(2);  % w0 * sqrt(1 + (zR / zR)^2)
wzzRy_m = w0y_m * sqrt(2);
thetaDx = lambda/(pi*w0x_m);
thetaDy = lambda/(pi*w0y_m);
wffx = (measRange_m(end) - z0x_m)*tan(thetaDx);
wffy = (measRange_m(end) - z0y_m)*tan(thetaDy);

%%
figNum = 1; 

xLimits = {[0 6];
           measRange_m};
yLimits = {[-0.25e-3 2.0e-3];
           [-0.25e-3 7.0e-3]};
figExt = {'Near';
          'Far'};

%%
for iZoom = 1:2

    figure(figNum)
    figNames{figNum} = [resultsDir figTag '_beamprofile_fit_' figExt{iZoom} 'Field.pdf'];
    figNum = figNum + 1;
    clf
    hold on
    plot(zPosition_m,xRadius_m,'bo')
    goox.plotBeamWidth(zModel_m,'b')
    pw0x = plot(z0x_m,w0x_m,'bx');
    pwZx = plot(z0x_m+zRx_m,wzzRx_m,'bs');
    divx = plot([z0x_m measRange_m(end)],[0 wffx],'b--');
    plot(zPosition_m,yRadius_m,'ro')
    gooy.plotBeamWidth(zModel_m,'r')
    pw0y = plot(z0y_m,w0y_m,'rx');
    pwZy = plot(z0y_m+zRy_m,wzzRy_m,'rs');
    divy = plot([z0y_m measRange_m(end)],[0 wffy],'r--');
    hold off
    set(pw0x,'MarkerSize',14)
    set(pw0y,'MarkerSize',14)
    set(pwZx,'MarkerSize',16)
    set(pwZy,'MarkerSize',16)
    leg = legend('Measured Data, X',...
        'Fit, X',...
        ['       (w0x, z0x) = (' num2str(w0x_m*1e3) ',' num2str(z0x_m) ') [mm,m]'],...
        'Fit Waist, w0x',...
        'waist at zRx (+ z0x); w(zR+z0x) = sqrt(2)*w0x',...
        'Divergence, waist at end = (end - z0x)*atan(\lambda / [\pi * w0x])',...
        'Measured Data, Y',...
        'Fit, Y',...
        ['       (w0y, z0y) = (' num2str(w0y_m*1e3) ',' num2str(z0y_m) ') [mm,m]'],...
        'Fit Waist, w0y',...
        'waist at zRy (+ z0y); w(zR+z0y) = sqrt(2)*w0y',...
        'Divergence, waist at end = (end - z0y)*atan(\lambda / [\pi * w0y])',...
        'Location','NorthWest');
    xlabel('Distance From Collimator [m]')
    ylabel('Beam Waist [m]')
    set(gca,'XMinorGrid','on','YMinorGrid','on')
    set(leg,'FontSize',legFontSize)
    title({[fcsn ' SuK Fiber Collimator Characterization'];...
           [measDate ' (' measNote '); ' figExt{iZoom} ' Field']})
    xlim(xLimits{iZoom})
    ylim(yLimits{iZoom})

end

%%
if printFigs
    addpath(toolsDir)
    mergeCommand = [svnRoot 'Common/Misc/pdfmerge ' resultsDir figTag '_beamprofile_fit.pdf'];

    for iFig = 1:length(figNames)
        figure(iFig)
        FillPage('w')
        IDfig(author)
        saveas(gcf,figNames{iFig})

        mergeCommand = [mergeCommand ' ' figNames{iFig}];
    end

    system(mergeCommand);
end
