%--------------------------------------------------------------------------
% FTE_aLIGO_thermal_mirrors.m
%
% Matlab example which calculates the thermal distortions on the input test 
% mass (itm) and end test mass (etm) in an Advanced LIGO arm cavity with 
% 1ppm absorption in the mirrors.  This example also produces finesse map 
% files of the thermal distortions for use in simulations.  It also returns
% the radii of curvature for the distortions.
%
% In this example the overall thermal distortion is calculated and 
% approximated by a parabalid (i.e. curved surface).  The residual
% distortion (overal - approximate parabaloid) is also calculated.
%
% Part of the Simtools package, http://www.gwoptics.org/simtools
% Charlotte Bond    13.02.2012
%--------------------------------------------------------------------------
%

clear all;

% Choose coating absorption
option = 'coating';
    
% 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 (aLIGO parameters)
% Two beam parameters at the two cavity mirrors
P = 770e3;
p = 0;
l = 0;
w_itm = 5.3e-2;
w_etm = 6.2e-2;
lambda = 1064e-9;
    
% Initialising thermal beam structures
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 material parameters
% Case for 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 to calculate thermal distortions
dn_dT = 0;
Babs = 0;
    
% Intiate thermal mirror
mirror = FT_init_thermal_mirror(dn_dT,Cabs,Babs,K,E,alpha,Y,pois,a,h);
    
% Radial coordinate
r = -a:2*a/256:a;
    
% Thermal distortion: etm
[U_etm,c_etm,p_etm] = FT_thermal_distortion(beam_etm,mirror,r,option,num_terms);

% Distortion described by approximate parabaloid: etm 
Up_etm = c_etm*r.^2 + p_etm;

% Residual distortion: etm
Ur_etm = U_etm - Up_etm;

% Creating mirror maps representing overall and residual distortions: 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: etm
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);

% Distortion described by approximate parabaloid: itm
Up_itm = c_itm*r.^2 + p_itm;

% Residual distortion: itm
Ur_itm = U_itm - Up_itm;
    
% Creating mirror maps representing overall and residual distortions: 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 to file: itm
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
    % Plotting etm thermal distortion and approximate parabaloid
    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])
                
figure
    % Plotting etm residual distortion
    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])

% Plotting thermal distortion as map: etm
FT_plot_mirror_map(Umap_etm,0)

% Plotting residual distortion as map: etm
FT_plot_mirror_map(rmap_etm,0)
        
figure
    % Plotting itm thermal distortion and approximate parabaloid
    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])
                
figure
    % Plotting itm residual distortion
    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])

% Plotting thermal distortion as map: itm
FT_plot_mirror_map(Umap_itm,0)
    
% Plotting residual distortion as map: itm
FT_plot_mirror_map(rmap_itm,0)
    
% Outputting radii of curvature of distortions (of approximate parabaloid)
disp(sprintf('ETM approximate (gaussian weighted) radius of curvature: %g [m]',Rc_etm))
disp(sprintf('ITM approximate (gaussian weighted) radius of curvature: %g [m]',Rc_itm))



