% Code to generate plots for NGeo 2017 paper clear all close all %% Figure 1 %MIPAS SF6 age load('MIPAS_age_on_theta.mat') MIPASgamma=age_interp; figure [C,H]=contourf(1:132,lats,squeeze(MIPASgamma(:,12,:)),0.5:.5:8,'linewidth',0.1);%due to some interpolation,there are values up to 30 years or so. ignore. colormap(pmkmp(15)) hcb=colorbar; set(gca,'Clim',[0.5 8], 'fontsize',18) w=H.LineWidth; H.LineWidth=.1; colorTitleHandle = get(hcb,'Title'); titleString = 'yrs'; h1=ylabel('latitude'); xlabel('year') set(gca,'XTick',[1:12:132],'linewidth',1) set(gca,'XTickLabel',['02';'03';'04';'05';'06';'07';'08';'09';'10';'11';'12']) set(colorTitleHandle ,'String',titleString,'fontsize',14); set(gcf, 'PaperUnits', 'inches'); set(gcf, 'PaperSize', [10 2.8]); set(gcf, 'PaperPosition', [0 0 10 2.8]); %N2O age load( 'Age_from_GOZCARDS_N2O.mat') n2oage=zeros(18,132); n2oage(:,25:end)=gammalag_andrews_th(:,12,:);% line up timing with MIPAS timing n2oage(n2oage==0)=NaN; figure [C,H]=contourf(1:132,lat,squeeze(n2oage),0.5:.5:8,'linewidth',0.1); hcb=colorbar; colormap(pmkmp(15)) set(gca,'Clim',[0.5 8], 'fontsize',18) w=H.LineWidth; H.LineWidth=.1; colorTitleHandle = get(hcb,'Title'); titleString = 'yrs'; h1=ylabel('latitude'); xlabel('year') set(gca,'XTick',[1:12:132],'linewidth',1) set(gca,'XTickLabel',['02';'03';'04';'05';'06';'07';'08';'09';'10';'11';'12']) set(colorTitleHandle ,'String',titleString,'fontsize',14); set(gcf, 'PaperUnits', 'inches'); set(gcf, 'PaperSize', [10 2.8]); set(gcf, 'PaperPosition', [0 0 10 2.8]); %and on any given surface, for 2007-2011 as in supplemental: figure contourf(1:60,lat,squeeze(gammalag_andrews_th(:,12,37:37+60-1)),0:.5:5.5) %37 to start in 2007 colormap(pmkmp(11)) set(gca,'Clim',[0 5.5],'fontsize',20) title('N_2O-Age 500 K, Andrews 2001') hcb=colorbar; colorTitleHandle = get(hcb,'Title'); titleString = 'years'; titleString = 'years'; ylabel('latitude') xlabel('year') set(gca,'XTick',[1:12:60]) set(gca,'XTickLabel',['07'; '08' ; '09'; '10'; '11'; '12']) set(gcf, 'PaperUnits', 'inches'); set(gcf, 'PaperSize', [12 4]); set(colorTitleHandle ,'String',titleString,'fontsize',18); set(gcf, 'PaperPosition', [0 0 12 4]); %WACCM ideal age load('WACCM_free_age.mat') figure [C,H]=contourf(1:132,lat,squeeze(age_interp(:,12,240+24:419-24)),0.5:.5:8,'linewidth',0.1); colormap(pmkmp(15)) hcb=colorbar; set(gca,'Clim',[0.5 8], 'fontsize',18) w=H.LineWidth; H.LineWidth=.1; colorTitleHandle = get(hcb,'Title'); titleString = 'yrs'; h1=ylabel('latitude'); xlabel('year') set(gca,'XTick',[1:12:132],'linewidth',1) set(gca,'XTickLabel',['02';'03';'04';'05';'06';'07';'08';'09';'10';'11';'12']) set(colorTitleHandle ,'String',titleString,'fontsize',14); set(gcf, 'PaperUnits', 'inches'); set(gcf, 'PaperSize', [10 2.8]); set(gcf, 'PaperPosition', [0 0 10 2.8]); %WACCM SF6 age load('WACCM_free_SF6_age.mat') figure [C,H]=contourf(1:132,lats,squeeze(age_interp(:,3,24:end-24)),0.5:.5:8,'linewidth',0.1); colormap(pmkmp(15)) hcb=colorbar; set(gca,'Clim',[0.5 8], 'fontsize',18) w=H.LineWidth; H.LineWidth=.1; colorTitleHandle = get(hcb,'Title'); titleString = 'yrs'; h1=ylabel('latitude'); xlabel('year') set(gca,'XTick',[1:12:132],'linewidth',1) set(gca,'XTickLabel',['02';'03';'04';'05';'06';'07';'08';'09';'10';'11';'12']) set(colorTitleHandle ,'String',titleString,'fontsize',14); set(gcf, 'PaperUnits', 'inches'); set(gcf, 'PaperSize', [10 2.8]); set(gcf, 'PaperPosition', [0 0 10 2.8]); %% Figure 2 load('LinzNGeo2017_Figure2.mat') load('trop_heights_MIPAS_forplots.mat') colors1=parula(12); addpath cbrewer colors2=cbrewer('div','PiYG',11); colors3=cbrewer('div','PRGn',11); colors4=cbrewer('div','RdYlBu',11); figure dx1=dgammaMIPAS_std; y=Thlevels; x=MIPASdgamma; fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],colors1(1,:)+[.6 .6 .4]-.2*[-.15 .05 .05],'linestyle','none'); hold on h3=plot(x,y,'color',colors1(1,:)-[-.15 .05 .05],'linewidth',2); dx1=WACCMfree_SF6dgamma_std; x=WACCMSF6; y=Thlevels; fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],colors3(7,:),'linestyle','none'); hold on h4=plot(x,y,'color',colors3(10,:),'linewidth',2); dx1=WACCMdgammastd; x=WACCMdgamma; y=Thlevels; fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],[.85 .85 1],'linestyle','none'); hold on h5=plot(x,y,'color','b','linewidth',2); dx1=N2Ostd; x=N2Odgamma; y=Thlevels(7:10); fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],[.8 .8 .8],'linestyle','none'); hold on h2=plot(x,y,'color','k','linewidth',2); plot(MIPASdgamma, Thlevels,'color',colors1(1,:)-[-.15 .05 .05],'linewidth',2); plot(WACCMSF6,Thlevels,'color',colors3(10,:),'linewidth',2); plot(WACCMdgamma,Thlevels,'color','b','linewidth',2); dx1=MoverMu_std/(365*60*60*24); x=MoverMu/(365*60*60*24); y=Thlevels; hold on h7=plot(x,y,'color','b','linewidth',2,'linestyle',':'); %intentionally have left off std because it makes plot messy set(gca,'fontsize',22) xlabel('\Delta \Gamma (years)') ylabel('\Theta (K)') set(gca,'Ylim',[400 1200]) set(gca,'Xlim',[0 2.5]) set(gca,'yscale','log') legend([h3 h2 h4 h5 h7],{'MIPAS SF_6','N_2O', 'WACCM SF_6','WACCM age','WACCM theory'},'Position',[0.215 .28 0.2 0.2],'fontsize',22);%,'WACCM M/\mu') legend('boxoff') hold on ax1=gca; set(ax1,'box','off','color','none') hold on ax1_pos = get(ax1,'Position'); % position of first axes hold on ax2 = axes('Position',ax1_pos,... 'XAxisLocation','top',... 'YAxisLocation','right','Ylim',[400 1200],'Yscale','log',... 'Ytick',thlevs_for_ax2,'Yticklabel',[18 21 24 27 30 33 36],... 'Xlim',[0 2.5],'Xticklabel',[''],... 'Color','none'); ylabel(ax2,'Tropical altitude (km)') set(ax2,'fontsize',22) set(ax2,'box','off','color','none') set(gcf, 'PaperUnits', 'inches'); set(gcf, 'PaperSize', [8 6]); set(gcf, 'PaperPosition', [0 0 8 6]); %% Figure 3 load('NGeo2017_Figure3.mat') Thlevels=[390:10:500 525:25:700 750:50:1200]; N2O_averages=[nanmean(totalmass(:,37:36+12)./(extage(:,37:36+12)-tropage(:,37:36+12)),2)... nanmean(totalmass(:,49:36+24)./(extage(:,49:36+24)-tropage(:,49:36+24)),2)... nanmean(totalmass(:,61:36+36)./(extage(:,61:36+36)-tropage(:,61:36+36)),2)... nanmean(totalmass(:,73:36+48)./(extage(:,73:36+48)-tropage(:,73:36+48)),2)... nanmean(totalmass(:,85:36+60)./(extage(:,85:36+60)-tropage(:,85:36+60)),2)]; N2O_stddev=nanstd(N2O_averages,1,2); colors1=parula(12); figure dx1=dev_Moverdgamma/365/24/60/60;%units! x=MIPAS_overturning/365/24/60/60; y=Thlevels; fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],colors1(1,:)+[.6 .6 .4]-.2*[-.15 .05 .05],'linestyle','none'); hold on h1=plot(x,y,'color',colors1(1,:)-[-.15 .05 .05],'linewidth',3); dx1=N2O_stddev(7:10)/365/24/60/60; x=N2O_overturning(7:10)/365/24/60/60; y=Thlevels(7:10); fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],[.8 .8 .8],'linestyle','none'); hold on h2=plot(x,y,'color','k','linewidth',2); dx1=JRAstd/365/24/60/60; x=JRA55_overturning/365/24/60/60; y=Thlevels; fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],colors1(4,:)+[.35 .35 .1],'linestyle','none'); hold on h3=plot(x,y,'color',colors1(4,:),'linewidth',3,'linestyle','--'); dx1=WACCM_free_std;%units are right x=Mu_comp_WACCM_free; y=Thlevels; fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],[.85 .85 1],'linestyle','none'); hold on h8=plot(x,y,'color','b','linewidth',3,'linestyle',':'); dx1=MERRA_std;%units are right x=Mu_comp_MERRA; y=Thlevels; fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],colors1(8,:)+[.1 .1 .1],'linestyle','none'); hold on h4=plot(x,y,'color',colors1(8,:)-[.15 .15 .15],'linewidth',3,'linestyle','--'); dx1=std_ERA*365/365/24/60/60; x=Mu_comp_ERA*365/365/24/60/60; y=Thlevels; fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],colors1(10,:)+[.02 .1 .15],'linestyle','none'); hold on h5=plot(x,y,'color',colors1(10,:)-[.15 .15 .15],'linewidth',3,'linestyle','--'); dx1=N2O_stddev(7:10)/365/24/60/60; x=N2O_overturning(7:10)/365/24/60/60; y=Thlevels(7:10); fill([x-dx1;flipud(x+dx1)],[y';flipud(y')],[.8 .8 .8],'linestyle','none'); hold on h6=plot(x,y,'color','k','linewidth',3); xlabel('overturning (kg/s)') ylabel('\Theta (K)') set(gca,'Xlim',[0 15*10^9]) set(gca,'Ylim',[400 1200]) set(gca,'yscale','log') set(gca,'fontsize',22) set(gca,'XTick',[0:2.5:15]*10^9) legend([h1 h2 h3 h4 h5 h8],{'MIPAS SF_6', 'GOZCARDS N_2O', 'JRA55', 'MERRA', 'ERA-I','WACCM' },'location','northeast','fontsize',22); legend('boxoff') hold on ax1=gca; set(ax1,'box','off','color','none') hold on ax1_pos = get(ax1,'Position'); hold on ax2 = axes('Position',ax1_pos,... 'XAxisLocation','top',... 'YAxisLocation','right','Ylim',[400 1200],'Yscale','log',... 'Ytick',thlevs_for_ax2,'Yticklabel',[18 21 24 27 30 33 36],... 'Xlim',[0 15*10^9],'XTick',[0:2.5:15]*10^9,'Xticklabel',[''],... 'Color','none'); ylabel(ax2,'Tropical altitude (km)') set(ax2,'fontsize',22) set(ax2,'box','off','color','none') set(gcf, 'PaperUnits', 'inches'); set(gcf, 'PaperSize', [8 6]); set(gcf, 'PaperPosition', [0 0 8 6]);