%% Read data
clear all

% Read in 1D data for normalization
oneD_withTD = importdata('/media/jschlup/Storage/Research/H2/1D_Detailed/WSSCI_Meeting/2015/followup/trans_nobindiff_output/MA_with_TD_UpdatedMech/H2_p01_0phi0_4000tu0298.kg');
oneD_noTD = importdata('/media/jschlup/Storage/Research/H2/1D_Detailed/WSSCI_Meeting/2015/followup/trans_nobindiff_output/MA_no_TD_UpdatedMech/H2_p01_0phi0_4000tu0298.kg');
[peaksource_withTD, r] = max(oneD_withTD.data(:,32));
peaksource_withTD = peaksource_withTD./oneD_withTD.data(r,3); % Divide by density to get units to match NGA
[peaksource_noTD, r]= max(oneD_withTD.data(:,32));
peaksource_noTD = peaksource_noTD./oneD_noTD.data(r,3); % Divide by density to get units to match NGA



% Read in Cyl. data
% [cyl_x_t,cyl_y_t,cyl_z_t,cyl_nx_t,cyl_ny_t,cyl_nz_t,cyl_dx,cyl_dy,cyl_dz] = config_reader('/media/jschlup/Storage/Research/H2/2D_Detailed/Other_Cases/cellular_tubular_forcesteadystate_from_DTinside/config_128z');
% cyl_data = NGAdatareader('/media/jschlup/Storage/Research/H2/2D_Detailed/Other_Cases/cellular_tubular_forcesteadystate_from_DTinside/data_withSoret_forcesteady_WBC_128z_1.244E+00');
% cyl_optdata = NGAdatareader('/media/jschlup/Storage/Research/H2/2D_Detailed/Other_Cases/cellular_tubular_forcesteadystate_from_DTinside/optdata_withDT_forcesteady_1.244E+00');

% Read in Cart. data
[cart_x_t,cart_y_t,cart_z_t,cart_nx_t,cart_ny_t,cart_nz_t,cart_dx,cart_dy,cart_dz] = config_reader('/media/jschlup/Storage/Research/H2/2D_Detailed/Papers/thermal_diffusion_updatemech/config_16dxperlF');
cart_TD_data = NGAdatareader('/media/jschlup/Storage/Research/H2/2D_Detailed/Papers/thermal_diffusion_updatemech/data_twomodes_withTD_4.000E-01');
cart_noTD_data = NGAdatareader('/media/jschlup/Storage/Research/H2/2D_Detailed/Papers/no_thermal_diffusion_updatemech/data_twomodes_noTD_4.000E-01');
cart_TD_prodrates = NGAdatareader('/media/jschlup/Storage/Research/H2/2D_Detailed/Papers/thermal_diffusion_updatemech/prodRates_twomodes_withTD_4.000E-01');
cart_noTD_prodrates = NGAdatareader('/media/jschlup/Storage/Research/H2/2D_Detailed/Papers/no_thermal_diffusion_updatemech/prodRates_twomodes_noTD_4.000E-01');

%% Compute stuff
lf = 0.000671608;
for i=1:length(cart_x_t)-1
    xm(i) = (cart_x_t(i+1)+cart_x_t(i))/2;
end
for i=1:length(cart_y_t)-1
    ym(i) = (cart_y_t(i+1)+cart_y_t(i))/2;
end
for i=1:length(cart_z_t)-1
    zm(i) = (cart_z_t(i+1)+cart_z_t(i))/2;
end

%% Plot
close all
% Initial disturbance plot
Width=cart_x_t(end);
Height=cart_y_t(end);
A = 1e-04;
k1 = 20;
B = 1e-04;
k2 = 13;
%b=abs(flamex_t-Width+A+B);
b=0.05;
x_profile=A*cos(k1*2*pi*cart_y_t/Height)-B*cos(k2*2*pi*cart_y_t/Height); %Perturbed
figure()
plot(x_profile/lf,cart_y_t/lf,'k','LineWidth',2); axis image
xlim([-14 14]); ylim([0 30]);
set(gca,'FontSize',16,'LineWidth',2);
xlabel('$x/l_F$'); ylabel('$y/l_F$');
text(3,27,'Products','Interpreter','Latex','FontSize',16);
text(-10,27,'Reactants','Interpreter','Latex','FontSize',16);
annotation('arrow',[0.25 0.35],[1/4 1/4]); annotation('arrow',[0.25 0.35],[0.5 0.5]);
annotation('arrow',[0.25 0.35],[3/4 3/4]); 
text(-8,4,'$S_{eff}$','Interpreter','Latex','FontSize',16);
text(-8,23,'$S_{eff}$','Interpreter','Latex','FontSize',16);
ax = gca;
ax.XMinorTick = 'on'; ax.XAxis.MinorTickValuesMode = 'manual'; ax.XAxis.MinorTickValues = -14:2:14; 
ax.YMinorTick = 'on'; ax.YAxis.MinorTickValuesMode = 'manual'; ax.YAxis.MinorTickValues = 5:5:30;
export_fig('figures/2DFPconfig.png');

xmax = xm(end)/lf; ymax = ym(end)/lf;

% Temperature contour
figure('units','normalized','outerposition',[0 0 1 1])
contourf(xm/lf,ym/lf,cart_TD_data{7}{10}',50,'EdgeColor','none') ; colormap jet; axis image;  
caxis([300 1600]); 
pos=get(gca,'pos');
% h = colorbar; h.Label.String = '$T [K]$'; h.Label.Interpreter = 'Latex'; h.TickLabelInterpreter = 'Latex'; h.Ticks = [500 1000 1500]; h.LineWidth = 2; %title(h, '$T$ [K]','Interpreter','Latex')
set(gca,'pos',[pos(1) pos(2) pos(3)*0.8 pos(4)*0.8]); 
xlabel('$x/l_F$'); ylabel('$y/l_F$');
xlim([0 xmax]); ylim([0 ymax]);
% title('With Thermal Diffusion');
ax = gca;
set(ax,'XTick',[0 20 40 60 80 100]); set(ax,'YTick',[0 10 20]);
set(ax,'FontSize',36,'LineWidth',2); set(ax,'XMinorTick','on'); ax.XAxis.MinorTickValuesMode = 'manual'; ax.XAxis.MinorTickValues = 0:5:120;
export_fig('figures/2DFP_contourT_withTD.png');

figure('units','normalized','outerposition',[0 0 1 1])
contourf(xm/lf,ym/lf,cart_noTD_data{7}{10}',50,'EdgeColor','none') ; colormap jet; axis image; 
caxis([300 1600]);
pos=get(gca,'pos');
% h = colorbar; h.Label.String = '$T [K]$'; h.Label.Interpreter = 'Latex'; h.TickLabelInterpreter = 'Latex'; h.Ticks = [500 1000 1500]; h.LineWidth = 2; %title(h, '$T$ [K]','Interpreter','Latex')
set(gca,'pos',[pos(1) pos(2) pos(3)*0.8 pos(4)*0.8]);
xlabel('$x/l_F$'); ylabel('$y/l_F$');
xlim([0 xmax]); ylim([0 ymax]);
% title('No Thermal Diffusion');
ax = gca;
set(ax,'XTick',[0 20 40 60 80 100]); set(ax,'YTick',[0 10 20]);
set(ax,'FontSize',36,'LineWidth',2); set(ax,'XMinorTick','on'); ax.XAxis.MinorTickValuesMode = 'manual'; ax.XAxis.MinorTickValues = 0:5:120;
export_fig('figures/2DFP_contourT_noTD.png');

f = figure('units','normalized','outerposition',[0 0 1 1]);
contourf(xm/lf,ym/lf,cart_noTD_data{7}{10}',50,'EdgeColor','none') ; colormap jet; axis image; 
caxis([300 1600]);
pos=get(gca,'pos');
h = colorbar; h.Label.String = '$T [K]$'; h.Label.Interpreter = 'Latex'; h.TickLabelInterpreter = 'Latex'; h.Ticks = [500 1000 1500]; h.LineWidth = 2; %title(h, '$T$ [K]','Interpreter','Latex')
set(gca,'pos',[pos(1) pos(2) pos(3)*0.8 pos(4)*0.8]);
xlabel('$x/l_F$'); ylabel('$y/l_F$');
xlim([0 xmax]); ylim([0 ymax]);
% title('No Thermal Diffusion');
ax = gca;
set(ax,'XTick',[0 20 40 60 80 100]); set(ax,'YTick',[0 10 20]);
set(ax,'FontSize',36,'LineWidth',2); set(ax,'XMinorTick','on'); ax.XAxis.MinorTickValuesMode = 'manual'; ax.XAxis.MinorTickValues = 0:5:120;
c = get(f, 'children');
set(c(2), 'visible', 'off');
f.CurrentAxes.Children.Visible = 'off';
export_fig('figures/2DFP_contourT_colorbar.png');

% Source term contour
figure('units','normalized','outerposition',[0 0 1 1])
contourf(xm/lf,ym/lf,cart_TD_prodrates{7}{7}'./peaksource_withTD,50,'EdgeColor','none') ; colormap jet; axis image;
caxis([0 3.5])
% pos=get(gca,'pos');
% h = colorbar; h.Label.String = '$\dot{\omega}_{H_2O}/\dot{\omega}_{H_2O}^{1D}$'; h.Label.Interpreter = 'Latex'; h.TickLabelInterpreter = 'Latex'; h.Ticks = [1 2 3]; h.LineWidth = 2; %title(h, '$\dot{\omega}_{H_2O}/\dot{\omega}_{H_2O}^{1D}$','Interpreter','Latex')
set(gca,'pos',[pos(1) pos(2) pos(3)*0.8 pos(4)*0.8]);
xlabel('$x/l_F$'); ylabel('$y/l_F$');
xlim([0 xmax]); ylim([0 ymax]);
ax = gca;
set(ax,'XTick',[0 20 40 60 80 100]); set(ax,'YTick',[0 10 20]);
set(ax,'FontSize',36,'LineWidth',2); set(ax,'XMinorTick','on'); ax.XAxis.MinorTickValuesMode = 'manual'; ax.XAxis.MinorTickValues = 0:5:120;
export_fig('figures/2DFP_contoursrc_withTD.png');

figure('units','normalized','outerposition',[0 0 1 1])
contourf(xm/lf,ym/lf,cart_noTD_prodrates{7}{7}'./peaksource_noTD,50,'EdgeColor','none') ; colormap jet; axis image; 
caxis([0 3.5]);
% pos=get(gca,'pos');
% h = colorbar; h.Label.String = '$\dot{\omega}_{H_2O}/\dot{\omega}_{H_2O}^{1D}$'; h.Label.Interpreter = 'Latex'; h.TickLabelInterpreter = 'Latex'; h.Ticks = [1 2 3]; h.LineWidth = 2; %title(h, '$\dot{\omega}_{H_2O}/\dot{\omega}_{H_2O}^{1D}$','Interpreter','Latex')
set(gca,'pos',[pos(1) pos(2) pos(3)*0.8 pos(4)*0.8]);
xlabel('$x/l_F$'); ylabel('$y/l_F$');
xlim([0 xmax]); ylim([0 ymax]);
ax = gca;
set(ax,'XTick',[0 20 40 60 80 100]); set(ax,'YTick',[0 10 20]);
set(ax,'FontSize',36,'LineWidth',2); set(ax,'XMinorTick','on'); ax.XAxis.MinorTickValuesMode = 'manual'; ax.XAxis.MinorTickValues = 0:5:120;
export_fig('figures/2DFP_contoursrc_noTD.png');

f = figure('units','normalized','outerposition',[0 0 1 1]);
contourf(xm/lf,ym/lf,cart_noTD_prodrates{7}{7}'./peaksource_noTD,50,'EdgeColor','none') ; colormap jet; axis image; 
caxis([0 3.5]);
% pos=get(gca,'pos');
h = colorbar; h.Label.String = '$\dot{\omega}_{\rm H_2O}/\dot{\omega}_{\rm H_2O}^{\rm 1D}$'; h.Label.Interpreter = 'Latex'; h.TickLabelInterpreter = 'Latex'; h.Ticks = [1 2 3]; h.LineWidth = 2; %title(h, '$\dot{\omega}_{H_2O}/\dot{\omega}_{H_2O}^{1D}$','Interpreter','Latex')
set(gca,'pos',[pos(1) pos(2) pos(3)*0.8 pos(4)*0.8]);
xlabel('$x/l_F$'); ylabel('$y/l_F$');
xlim([0 xmax]); ylim([0 ymax]);
ax = gca;
set(ax,'XTick',[0 20 40 60 80 100]); set(ax,'YTick',[0 10 20]);
set(ax,'FontSize',36,'LineWidth',2); set(ax,'XMinorTick','on'); ax.XAxis.MinorTickValuesMode = 'manual'; ax.XAxis.MinorTickValues = 0:5:120;
c = get(f, 'children');
set(c(2), 'visible', 'off');
f.CurrentAxes.Children.Visible = 'off';
export_fig('figures/2DFP_contoursrc_colorbar.png');