clearvars -except cmaps
clear global
close all

global mconfig ivar2 ivar1 its ici nikki output_dir case_list_str vnum ...
   bintype var1_str var2_str indvar_name_set ...
   indvar_ename_set indvar_units_set %#ok<*NUSED>

nikki='Dt80';
mconfig='fullmic';
global_var
get_var_comp

load(sprintf('../pfm_summary/%s_%s_pfm_bincomp.mat',nikki,mconfig))
var1_str = pfm.misc.var1_str;
var2_str = pfm.misc.var2_str;

tmpvarname=fieldnames(pfm(1));
fldnms=fieldnames(pfm(1).(tmpvarname{1}));
fldnms=fldnms(1:end-1);

plot_var={'cloud_M1_path','mean_surface_ppt'};

figure('position',[0 0 800 480])
tl=tiledlayout(7,4);
lsty = {':', '-'};
tiles = [];
aux_line_color = [.5 .5 .5];

selected_cases = [3 2];
ivar1 = selected_cases(1); ivar2 = selected_cases(2);
its = 1;
bintau_struct = loadnc('bin', indvar_name_set);
amptau_struct = loadnc('amp', indvar_name_set);
its = 2;
binsbm_struct = loadnc('bin', indvar_name_set);
ampsbm_struct = loadnc('amp', indvar_name_set);

for ipvar=1:length(plot_var)
alb_idx=find(contains(indvar_name_set,plot_var{ipvar}));
ivar=alb_idx;
ifn=1;

nexttile(ipvar*2+15,[3 2])
mr=pfm.(indvar_name_set{ivar}).(fldnms{ifn});
nanimagesc(mr)
cb=colorbar;
rectangle('position',[selected_cases(2)-.5, selected_cases(1)-.5, 1, 1],'EdgeColor',...
   aux_line_color,'LineWidth',2)
xlab_key=extractBefore(var2_str,digitsPattern);
ylab_key=extractBefore(var1_str,digitsPattern);
xlab=[initVarName_dict(xlab_key{1}) initVarUnit_dict(xlab_key{1})];
ylab=[initVarName_dict(ylab_key{1}) initVarUnit_dict(ylab_key{1})];

if ipvar==length(bintype)
   cb.Label.String='TAU/SBM ratio'; 
else
   ylabel(ylab)
end
title(['(',char(ipvar+98),') ',indvar_ename_set{ivar} indvar_units_set{ivar} ' summary'],...
   'FontWeight','normal')
xticks(1:length(var2_str))
yticks(1:length(var1_str))
xticklabels(extractAfter(var2_str,lettersPattern))
yticklabels(extractAfter(var1_str,lettersPattern))
set(gca,'FontSize',11)
tiles = [tiles gca];

if strcmp(fldnms{ifn},'mr')
   colormap(cmaps.BrBG)
   set(gca,'ColorScale','log')
   caxis([.5 2])

   mpath=pfm.(indvar_name_set{ivar}).mpath_sbm;
   [XX,YY]=meshgrid(1:length(var2_str),1:length(var1_str));
   mpath_bin_str=sprintfc('%0.3g',mpath);

   for ivar1=1:length(var1_str)
      for ivar2=1:length(var2_str)

         % ----- get text color -----
         ngrads=size(cmaps.coolwarm_r,1);
         clr_idx=roundfrac(pfm.(indvar_name_set{ivar}).rsq(ivar1,ivar2),1/ngrads)*ngrads;
         clr_idx=round(clr_idx); % in case prev line outputs double
         % ----- got text color -----

         if isnan(clr_idx) continue, end
         if clr_idx==0 clr_idx=1; end

         text(ivar2+0.015,ivar1-0.015,mpath_bin_str{ivar1,ivar2},'FontSize',11,...
            'HorizontalAlignment','center',...
            'Color',cmaps.coolwarm_r(clr_idx,:)*.1)
         text(ivar2,ivar1,mpath_bin_str{ivar1,ivar2},'FontSize',11,...
            'HorizontalAlignment','center',...
            'Color',cmaps.coolwarm_r(clr_idx,:))


      end
   end

end

nexttile(2,[1,2])
set(gca,'Color','none')
set(gca,'XColor','none')
set(gca,'YColor','none')
colormap(gca,cmaps.coolwarm_r)
cb=colorbar('southoutside');
cb.Label.String='R^2';
cb.Label.Position=[0.5000 3.3 0];
set(gca,'FontSize',11)

nexttile(ipvar*2+3,[3 2])

time = bintau_struct.time;
z = bintau_struct.z;

var_comp_raw_amptau = amptau_struct.(indvar_name_set{ivar});
[var_comp_amptau, linORlog, range] = var2phys(var_comp_raw_amptau,ivar,1);
var_comp_raw_ampsbm = ampsbm_struct.(indvar_name_set{ivar});
var_comp_ampsbm = var2phys(var_comp_raw_ampsbm,ivar,1);

if ~contains(indvar_name_set{ivar}, amp_only_var)
   var_comp_raw_bintau = bintau_struct.(indvar_name_set{ivar});
   var_comp_bintau = var2phys(var_comp_raw_bintau,ivar,1);
   var_comp_raw_binsbm = binsbm_struct.(indvar_name_set{ivar});
   var_comp_binsbm = var2phys(var_comp_raw_binsbm,ivar,1);
end

hold on
plot(time, var_comp_bintau, 'LineWidth', 2, 'DisplayName', 'bin-TAU', 'Color', color_order{1})
plot(time, var_comp_binsbm, 'LineWidth', 2, 'DisplayName', 'bin-SBM', 'Color', color_order{2})
plot(time, var_comp_amptau, 'LineWidth', 1, 'DisplayName', 'AMP-TAU', 'LineStyle', ':', 'Color', color_order{1})
plot(time, var_comp_ampsbm, 'LineWidth', 1, 'DisplayName', 'AMP-SBM', 'LineStyle', ':', 'Color', color_order{2})
xlabel('Time [s]')
ylabel([indvar_ename_set{ivar} indvar_units_set{ivar}])
tiles = [tiles gca];
title(['(' char(ipvar+96) ') time series of ' indvar_ename_set{ivar}], 'FontWeight', 'normal')
hold off
grid
box on
set(gca,'FontSize',11)
end

xlabel(tl,xlab,'fontsize',11)
title(tl,'      Full microphysics - TAU vs SBM','fontsize',18,'fontweight','bold')

axes(tiles(1))
[xl_start1, yl_start1] = coord2norm(gca, selected_cases(2)-.5, selected_cases(1)+.5);
[xr_start1, yr_start1] = coord2norm(gca, selected_cases(2)+.5, selected_cases(1)+.5);
axes(tiles(3))
[xl_start2, yl_start2] = coord2norm(gca, selected_cases(2)-.5, selected_cases(1)+.5);
[xr_start2, yr_start2] = coord2norm(gca, selected_cases(2)+.5, selected_cases(1)+.5);
axes(tiles(2))
lower_ax_pos = get(gca,'position');
xl_end1 = lower_ax_pos(1);
xr_end1 = lower_ax_pos(1) + lower_ax_pos(3);
y_end1 = lower_ax_pos(2);
legend % legends would disappear if I put this statement within the for loop above
axes(tiles(4))
lower_ax_pos = get(gca,'position');
xl_end2 = lower_ax_pos(1);
xr_end2 = lower_ax_pos(1) + lower_ax_pos(3);
y_end2 = lower_ax_pos(2);
legend

annotation('line', [xl_start1 xl_end1], [yl_start1 y_end1], ...
   'LineStyle', '--', 'Color', aux_line_color, 'LineWidth', 1.5)
annotation('line', [xr_start1 xr_end1], [yr_start1 y_end1], ...
   'LineStyle', '--', 'Color', aux_line_color, 'LineWidth', 1.5)
annotation('line', [xl_start2 xl_end2], [yl_start2 y_end2], ...
   'LineStyle', '--', 'Color', aux_line_color, 'LineWidth', 1.5)
annotation('line', [xr_start2 xr_end2], [yr_start2 y_end2], ...
   'LineStyle', '--', 'Color', aux_line_color, 'LineWidth', 1.5)

exportgraphics(gcf,'plots/pngs/6 fullmp_sbm_vs_tau.png', 'Resolution', 300)
exportgraphics(gcf,'plots/6 fullmp_sbm_vs_tau.pdf')
