%--------------------------------------------------------------------------
% function [] = FTE_AdLIGO_thermal_mirrors()
%
% Matlab example which calculates the thermal distortions on the itm and
% etm in an AdLIGO cavity, and produces finesse map files for use in
% simulations.  It also returns the radii of curvature for the distortions.
%
% Charlotte Bond    13.02.2012
%--------------------------------------------------------------------------
%

function [] = FTE_AdLIGO_thermal_mirrors()

    % Choose coating absorption
    option = 'coating';

    % Names to save plots
    imagename1 = 'etm_thermal_dis';
    imagename2 = 'etm_thermal_res';
    imagename3 = 'etm_thermal_dis_map';
    imagename4 = 'etm_thermal_res_map';
    imagename5 = 'itm_thermal_dis';
    imagename6 = 'itm_thermal_res';
    imagename7 = 'itm_thermal_dis_map';
    imagename8 = 'itm_thermal_res_map';
    
    % Filenames for map files
    file1 = 'etm_thermal_dis_map.txt';
    file2 = 'etm_thermal_res_map.txt';
    file3 = 'itm_thermal_dis_map.txt';
    file4 = 'itm_thermal_res_map.txt';
    filea1 = 'etm_aperture.txt';
    filea2 = 'itm_aperture.txt';
    
    % Number
    num_terms = 100;

    % Beam parameters (AdLIGO parameters)
    % Two beam parameters at the two mirrors
    P = 770e3;
    p = 0;
    l = 0;
    w_itm = 5.3e-2;
    w_etm = 6.2e-2;
    lambda = 1064e-9;
    
    beam_etm = FT_init_thermal_beam(P,p,l,w_etm,lambda);
    beam_itm = FT_init_thermal_beam(P,p,l,w_itm,lambda);
    
    % Mirror parameters -> same parameters for ETM and ITM
    % 1ppm absorption
    Cabs = 1e-6;
    K = 1.38;
    E = 1;
    alpha = 5.4e-7;
    Y = 7.3e10;
    pois = 0.17;
    a = 0.17;
    h = 0.2;
    % Parameters we don't need
    dn_dT = 0;
    Babs = 0;
    
    mirror = FT_init_thermal_mirror(dn_dT,Cabs,Babs,K,E,alpha,Y,pois,a,h);
    
    r = -a:2*a/401:a;
    
    % Thermal distortion - etm
    [U_etm,c_etm,p_etm] = FT_thermal_distortion(beam_etm,mirror,r,option,num_terms);
    Up_etm = c_etm*r.^2 + p_etm;
    Ur_etm = U_etm - Up_etm;
    
    [Umap_etm,rmap_etm,Rc_etm] = FT_mirror_map_from_thermal_distortion(beam_etm,mirror,r,option,num_terms);
 
    % Creating aperture map and writing maps to file
    Amap_etm = FT_create_aperture_map(Umap_etm);
    FT_write_surface_map(file1,Umap_etm);
    FT_write_surface_map(file2,rmap_etm);
    FT_write_surface_map(filea1,Amap_etm);
    
    % Thermal distortion - itm
    [U_itm,c_itm,p_itm] = FT_thermal_distortion(beam_itm,mirror,r,option,num_terms);
    Up_itm = c_itm*r.^2 + p_itm;
    Ur_itm = U_itm - Up_itm;
    
    [Umap_itm,rmap_itm,Rc_itm] = FT_mirror_map_from_thermal_distortion(beam_itm,mirror,r,option,num_terms);
    
    % Creating aperture maps and writing maps
    Amap_itm = FT_create_aperture_map(Umap_itm);
    FT_write_surface_map(file3,Umap_itm);
    FT_write_surface_map(file4,rmap_itm);
    FT_write_surface_map(filea2,Amap_itm);
    
    % Plotting results
    figure
        hold on;
        plot(r,U_etm,'r','LineWidth',3)
        plot(r,Up_etm,'b--','LineWidth',3)
        grid on;
        xlabel('r [m] ')
        ylabel('Thermal distortion [m] ')
        title('Thermal distortion: ETM')
        legend('Thermal distortion','Approximate paraboloid')
        set(gca,'XLim',[-a a])
        
        plotall(imagename1,0,0)
        
    figure
        plot(r,Ur_etm,'r','LineWidth',3)
        grid on;
        xlabel('r [m] ')
        ylabel('Thermal distortion [m] ')
        title('Residual thermal distortion: ETM')
        set(gca,'XLim',[-a a])
        
        plotall(imagename2,0,0)
            
    FT_plot_mirror_map(Umap_etm,0)
    plotall(imagename3,0,0)
    FT_plot_mirror_map(rmap_etm,0)
    plotall(imagename4,0,0)
        
    figure
        hold on;
        plot(r,U_itm,'r','LineWidth',3)
        plot(r,Up_itm,'b--','LineWidth',3)
        grid on;
        xlabel('r [m] ')
        ylabel('Thermal distortion [m] ')
        title('Thermal distortion: ITM')
        legend('Thermal distortion','Approximate paraboloid')
        set(gca,'XLim',[-a a])
        
        plotall(imagename5,0,0)
        
    figure
        plot(r,Ur_itm,'r','LineWidth',3)
        grid on;
        xlabel('r [m] ')
        ylabel('Thermal distortion [m] ')
        title('Residual thermal distortion: ITM')
        set(gca,'XLim',[-a a])
        
        plotall(imagename6,0,0)
        
    FT_plot_mirror_map(Umap_itm,0)
    plotall(imagename7,0,0)
    FT_plot_mirror_map(rmap_itm,0)
    plotall(imagename8,0,0)
    
    % Outputting radii of curvature of distortions
    disp(sprintf('ETM approximate (weighted) radius of curvature: %g [m]',Rc_etm))
    disp(sprintf('ITM approximate (weighted) radius of curvature: %g [m]',Rc_itm))


end

