load("data_FigS2.mat")

figure('Color', 'w', 'Units', 'normalized', 'Position', [0.1 0.1 0.8 0.6])

col_tdse = [my_Green; 0.5 0.9 0.7];
col_smm = [0 0 0; 0.5 0.5 0.5];

for idx = 1:numel(near_field_list)
    Espec = Espec_NF(:,idx);
    
    subplot(3,3,2)
    plot(co_dir_nf(idx)*[1 1], [-10 2], 'Color', my_B)
    hold on
    plot(co_rec_nf(idx)*[1 1], [-10 2], 'Color', my_R)
        
    plot(Espec(Espec<co_dir_nf(idx)), log10(abs(spec_NF((Espec<co_dir_nf(idx)),idx)).^2)-2, 'Color', col_smm(idx,:), 'LineWidth', 1.5)
    plot(Espec(Espec>co_dir_nf(idx)), log10(abs(spec_NF((Espec>co_dir_nf(idx)),idx)).^2)-2, '--', 'Color', col_smm(idx,:), 'LineWidth', 1.5)

    ylabel 'log10(Yield)'
    xlim([0 4])
    ylim([-9 0])
    set(gca, 'FontSize', 12)

    subplot(3,3,3)
    plot(d, near_field(:,idx), 'LineWidth', 1.5, 'Color', col_tdse(idx,:))
    hold on
    xlabel 'distance from surface (Å)'
    ylabel 'enhancement'
    ylim([0 8])
    set(gca, 'FontSize', 12)
end

for idx = 1:numel(decay_list)
    Espec = Espec_decay(:,idx);
    
    subplot(3,3,(idx-1)*3+4)
    plot(co_dir_dec(idx)*[1 1], [-10 2], 'Color', my_B)
    hold on
    plot(co_rec_dec(idx)*[1 1], [-10 2], 'Color', my_R)
    
    % plot(out.E_list/Up, log10(out.spec), 'Color', my_Green, 'LineWidth', 1.5)

    plot(Espec(Espec<co_dir_dec(idx)), log10(abs(spec_decay((Espec<co_dir_dec(idx)),idx)).^2)-2, 'Color', col_smm(1,:), 'LineWidth', 1.5)
    plot(Espec(Espec>co_dir_dec(idx)), log10(abs(spec_decay((Espec>co_dir_dec(idx)),idx)).^2)-2, '--', 'Color', col_smm(1,:), 'LineWidth', 1.5)

    ylabel 'log10(Yield)'
    xlim([0 4])
    ylim([-9 0])
    set(gca, 'FontSize', 12)
end

for idx = 1:numel(intensity_list)
    Espec = Espec_intensity(:,idx);
    
    subplot(3,3,(idx-1)*3+5)
    plot(co_dir_int(idx)*[1 1], [-10 2], 'Color', my_B)
    hold on
    plot(co_rec_int(idx)*[1 1], [-10 2], 'Color', my_R)
    
    % plot(out.E_list/Up, log10(out.spec), 'Color', my_Green, 'LineWidth', 1.5)

    plot(Espec(Espec<co_dir_int(idx)), log10(abs(spec_intensity((Espec<co_dir_int(idx)),idx)).^2)-2, 'Color', col_smm(1,:), 'LineWidth', 1.5)
    plot(Espec(Espec>co_dir_int(idx)), log10(abs(spec_intensity((Espec>co_dir_int(idx)),idx)).^2)-2, '--', 'Color', col_smm(1,:), 'LineWidth', 1.5)

    xlim([0 4])
    ylim([-9 0])
    set(gca, 'FontSize', 12)
end

for idx = 1:numel(wavelength_list)
    Espec = Espec_wavelength(:,idx);
    
    subplot(3,3,(idx-1)*3+6)
    plot(co_dir_wavelength(idx)*[1 1], [-10 2], 'Color', my_B)
    hold on
    plot(co_rec_wavelength(idx)*[1 1], [-10 2], 'Color', my_R)
    
    % plot(out.E_list/Up, log10(out.spec), 'Color', my_Green, 'LineWidth', 1.5)

    plot(Espec(Espec<co_dir_wavelength(idx)), log10(abs(spec_wavelength((Espec<co_dir_wavelength(idx)),idx)).^2)-2, 'Color', col_smm(1,:), 'LineWidth', 1.5)
    plot(Espec(Espec>co_dir_wavelength(idx)), log10(abs(spec_wavelength((Espec>co_dir_wavelength(idx)),idx)).^2)-2, '--', 'Color', col_smm(1,:), 'LineWidth', 1.5)

    xlim([0 4])
    ylim([-9 0])
    set(gca, 'FontSize', 12)
end

subplot(3,3,7)
xlabel 'Energy (Up)'
subplot(3,3,8)
xlabel 'Energy (Up)'
subplot(3,3,9)
xlabel 'Energy (Up)'

% add TDSE comparison

subplot(3,3,2)
plot(tdse.E_exp, log10(tdse.spec_exp), 'LineWidth', 1.5, 'Color', col_tdse(1,:))
plot(tdse.E_nf, log10(tdse.spec_nf), 'LineWidth', 1.5, 'Color', col_tdse(2,:))

subplot(3,3,4)
plot(tdse.E_lnf200, log10(tdse.spec_lnf200), 'LineWidth', 1.5, 'Color', col_tdse(1,:))

subplot(3,3,5)
plot(tdse.E_int12, log10(tdse.spec_int12), 'LineWidth', 1.5, 'Color', col_tdse(1,:))

subplot(3,3,6)
plot(tdse.E_wave1100, log10(tdse.spec_wave1100), 'LineWidth', 1.5, 'Color', col_tdse(1,:))

subplot(3,3,7)
plot(tdse.E_lnf50, log10(tdse.spec_lnf50), 'LineWidth', 1.5, 'Color', col_tdse(1,:))

subplot(3,3,8)
plot(tdse.E_int11, log10(tdse.spec_int11), 'LineWidth', 1.5, 'Color', col_tdse(1,:))

subplot(3,3,9)
plot(tdse.E_wave500, log10(tdse.spec_wave500), 'LineWidth', 1.5, 'Color', col_tdse(1,:))