#!/usr/bin/env python """ ==================================================================================================== Script to plot figures for SAI impact on cloud brightness paper. Created by Jake Gristey (jake.j.gristey@noaa.gov) Aug 2024 ==================================================================================================== """ #-------------------------------------------ADD IMPORTS--------------------------------------------# import numpy as np import matplotlib from matplotlib.ticker import AutoMinorLocator,MaxNLocator import matplotlib.pyplot as plt import matplotlib.image as mpimg from matplotlib.offsetbox import OffsetImage, AnnotationBbox #-----------------------------------------DEFINE FUNCTIONS-----------------------------------------# def main(): """ Direct execution code. """ #*****************************************USER INPUT*******************************************# figs=[1,2,3,4] # list of figure numbers to be plotted letters=True # True to include figure letters on plots, False otherwise showPlots=False # True to show, False to save figDir='./' # dir to save figures dpi=400 # resolution of saved figures [dots per square inch] ftype='.png' # file format string (eg. '.png' or '.pdf') #**********************************************************************************************# ################################################################################################ # FIGURE 1 # ################################################################################################ if 1 in figs: # data (hard coded for now), comes from running Copsey MC code in trial-and-error style # using 100,000 photons every 10 deg to find a case where albedo=0.500 sza=[0,10,20,30,40,50,60,70,80] z=[800,790,760,695,605,496,360,209,69.1] tau=2*100e6*np.pi*(10e-6)**2*np.array(z) # set plot defaults matplotlib.rcParams['axes.linewidth']=1.5 matplotlib.rcParams['font.size']=14 # create figure plt.figure(figsize=(7,8)) # plot schematic (created in powerpoint) plt.subplot(2,1,1) ax=plt.gca() img=mpimg.imread('Fig_01a.png') plt.imshow(img) ax.set_aspect('equal') ax.tick_params(axis='x', which='both', bottom=False, top=False, labelbottom=False) ax.tick_params(axis='y', which='both', left=False, right=False, labelleft=False) #ax.set_axis_off() #imagebox = OffsetImage(img, zoom=0.31) #ab = AnnotationBbox(imagebox, (0.001, 0.001), xycoords='axes fraction', \ #box_alignment=(0, 0),bboxprops=dict(facecolor='none', edgecolor='none')) #ax.add_artist(ab) if letters: ob=matplotlib.offsetbox.AnchoredText('a',2,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) # plot photon penetration depth data plt.subplot(2,1,2) ax=plt.gca() ax.yaxis.set_major_locator(MaxNLocator(nbins=6)) ax.xaxis.set_major_locator(MaxNLocator(nbins=9)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=8, width=1.5) ax.tick_params(which='minor',direction='out', length=4, width=1.5) ax.set_xlabel('Incidence angle [$^\circ$]') ax.set_ylabel('Optical depth') ax.grid() ax2=ax.twinx() ax2.yaxis.set_major_locator(MaxNLocator(nbins=5)) ax2.xaxis.set_major_locator(MaxNLocator(nbins=9)) ax2.xaxis.set_minor_locator(AutoMinorLocator()) ax2.yaxis.set_minor_locator(AutoMinorLocator()) #ax2.xaxis.set_ticks_position('right') ax2.tick_params(which='major',direction='out', length=8, width=1.5) ax2.tick_params(which='minor',direction='out', length=4, width=1.5) ax2.set_ylabel('Geometric depth [m]') ax2.set_title('Median photon penetration depth',y=1.02) #'Median photon penetration depth:\n${\omega}$=1, $N_C$=100$cm^{-3}$, $r_e$=10${\mu}m$, \ #$g$=0.85?',y=1.02) ax2.grid() ax.plot(sza,tau,'-o',lw=2.,zorder=999) ax2.plot(sza,z,'-o',zorder=999) if letters: ob=matplotlib.offsetbox.AnchoredText('b',1,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax2.add_artist(ob) plt.tight_layout() if not showPlots: if letters: plt.savefig(figDir+'Fig_01'+ftype,dpi=dpi,bbox_inches='tight', \ pad_inches = 0.1) if not letters: plt.savefig(figDir+'withoutLetters/Fig_01'+ftype,dpi=dpi, \ bbox_inches='tight',pad_inches = 0.1) ################################################################################################ # FIGURE 2 # ################################################################################################ if 2 in figs: # read SBDART output data fname='../sbdart_output/cloud_albedo_sai.txt' f=open(fname,'r') data=f.readlines() f.close() del data[::34] # remove headers # define conditions used in SBDART aod=[0,0.5,1,1.5,2] tau=[0,1,3,5,10,20,50] sza=[0,20,40,60,80] naod=len(aod) ntau=len(tau) nsza=len(sza) nlev=33 nvar=7 # read output data into array (naod,ntau,nsza,nlev,nvar) data=np.array([data[i].rstrip('\n').split() for i in \ range(len(data))]).astype(np.float).reshape(naod,ntau,nsza,nlev,nvar) # set plot defaults matplotlib.rcParams['axes.linewidth']=1.5 matplotlib.rcParams['font.size']=14 # create figure plt.figure(figsize=(10,8)) letter_list=['a','b','c','d','e','f','g','h','i'] # plot US64 atmospheric profile of temperature plt.subplot(2,3,1) z=[0.,1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,11.,12.,13.,14.,15.,16.,17.,18.,19.,20.,21.,22., \ 23.,24.,25.,30.,35.,40.,45,50.,70.,100.] T=[2.881e+02,2.816e+02,2.751e+02,2.687e+02,2.622e+02,2.557e+02,2.492e+02,2.427e+02, \ 2.362e+02,2.297e+02,2.232e+02,2.168e+02,2.166e+02,2.166e+02,2.166e+02,2.166e+02, \ 2.166e+02,2.166e+02,2.166e+02,2.166e+02,2.166e+02,2.176e+02,2.186e+02,2.196e+02, \ 2.206e+02,2.216e+02,2.265e+02,2.365e+02,2.534e+02,2.642e+02,2.706e+02,2.197e+02, \ 2.100e+02] ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=5)) ax.yaxis.set_major_locator(MaxNLocator(nbins=6)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.plot(T,z,c='k',lw=3) ax.set_xlabel('Temperature [K]') ax.set_ylabel('Altitude [km]') ax.set_ylim(0,30) ax.set_xlim(206,292) ax.fill_between([0,1000],[20,20],[21,21],color='brown', alpha=0.5) ax.fill_between([0,1000],[1,1],[2,2],color='gray', alpha=0.5) ax.text(215, 2, 'Cloud layer',fontsize=10, va='bottom', ha='left',color='gray') ax.text(240, 21, 'SAI aerosol layer',fontsize=10, va='bottom', ha='left',color='brown') ax.text(230, 0, 'Ocean surface',fontsize=10, va='bottom', ha='left',color='blue') ax.text(0.5, 0.5, 'US62 Standard\nAtmosphere',fontsize=10, va='center', ha='center', \ color='k',transform=ax.transAxes) ax.spines['bottom'].set_color('blue') ax.set_title('Radiative Transfer Setup') #img=mpimg.imread('/export/home/jgristey/plots/RT_setup.png') #imagebox = OffsetImage(img, zoom=0.31) #ab = AnnotationBbox(imagebox, (0.001, 0.001), xycoords='axes fraction', \ #box_alignment=(0, 0),bboxprops=dict(facecolor='none', edgecolor='none')) #ax.add_artist(ab) if letters: ob=matplotlib.offsetbox.AnchoredText(letter_list[0],2,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) plt.subplot(2,3,2) ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=4)) ax.yaxis.set_major_locator(MaxNLocator(nbins=6)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.plot(data[0,4,0,:,4],data[0,4,0,:,0],lw=3, label='$F_{\downarrow,direct}$') ax.plot(data[0,4,0,:,2]-data[0,4,0,:,4],data[2,4,0,:,0],'--',lw=3, \ label='$F_{\downarrow,diffuse}$') #ax.plot(data[0,4,0,:,2],data[2,4,0,:,0],lw=3, label='$F_{\downarrow}$') ax.text(0.4, 0.15, 'Diffuse\nfraction\nat 2km'+ \ '\n={:.1f}'.format(((data[0,4,0,-3,2]-data[0,4,0,-3,4])/data[0,4,0,-3,2])*100)+'%', \ transform=ax.transAxes, fontsize=12, verticalalignment='bottom',ha='center', \ bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) ax.legend(loc='upper center') ax.set_xlabel('Irradiance [W m$^{-2}$]') ax.set_ylabel('Altitude [km]') ax.set_ylim(0,30) ax.set_xlim(0,1400) ax.set_title('Without SAI') if letters: ob=matplotlib.offsetbox.AnchoredText(letter_list[1],2,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) plt.subplot(2,3,3) ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=4)) ax.yaxis.set_major_locator(MaxNLocator(nbins=6)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.plot(data[2,4,0,:,4],data[2,4,0,:,0],lw=3, label='$F_{\downarrow,direct}$') ax.plot(data[2,4,0,:,2]-data[2,4,0,:,4],data[2,4,0,:,0],'--',lw=3, \ label='$F_{\downarrow,diffuse}$') #ax.plot(data[2,4,0,:,2],data[2,4,0,:,0],lw=3, label='$F_{\downarrow}$') ax.legend(loc='upper center') ax.set_xlabel('Irradiance [W m$^{-2}$]') ax.set_ylabel('Altitude [km]') ax.set_ylim(0,30) ax.set_xlim(0,1400) ax.set_title('With SAI') ax.text(0.7, 0.15, 'Diffuse\nfraction\nat 2km'+ \ '\n={:.1f}'.format(((data[2,4,0,-3,2]-data[2,4,0,-3,4])/data[2,4,0,-3,2])*100)+'%', \ transform=ax.transAxes, fontsize=12, verticalalignment='bottom', ha='center', \ bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) if letters: ob=matplotlib.offsetbox.AnchoredText(letter_list[2],2,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) plt.subplot(2,3,4, polar=True) ax=plt.gca() theta = [0.0000E+00, 1.0000E+01, 2.0000E+01, 3.0000E+01 , 4.0000E+01 , 5.0000E+01 , \ 6.0000E+01 , 7.0000E+01, 8.0000E+01, 1.0000E+02,1.1000E+02 , 1.2000E+02 , \ 1.3000E+02 , 1.4000E+02 , 1.5000E+02, 1.6000E+02 , 1.7000E+02, 1.8000E+02] r1 = np.genfromtxt( \ '../sbdart_output/cloud_albedo_sai_withRadiances.txt', \ skip_header=21,skip_footer=1) r1=np.append(r1,[2.8665E+02,2.8665E+02,2.8665E+02,2.8665E+02,2.8665E+02,2.8665E+02]) r1=r1.reshape(33,18,19) r2 = np.genfromtxt( \ '../sbdart_output/cloud_albedo_saiOFF_withRadiances.txt', \ skip_header=21,skip_footer=1) r2=np.append(r2,[2.8665E+02,2.8665E+02,2.8665E+02,2.8665E+02,2.8665E+02,2.8665E+02]) r2=r2.reshape(33,18,19) ax.plot(np.deg2rad(theta)+np.pi/2, r2[-3,:,12],c='C8',lw=3,label='Without SAI') ax.plot(np.deg2rad(theta)+np.pi/2, r1[-3,:,12],'--',c='C4',lw=3,label='With SAI') ax.plot(np.deg2rad(theta)+np.pi/2+np.pi/2+np.pi/2, np.flip(r1[-3,:,12]),'--',lw=3,c='C4') ax.plot(np.deg2rad(theta)+np.pi/2+np.pi/2+np.pi/2, np.flip(r2[-3,:,12]),lw=3,c='C8') #ax.tick_params(axis='y',labelrotation=90.) ax.set_thetamin(180) ax.set_thetamax(360) ax.set_rmax(900) ax.set_yticks([0,300,600,900]) ax.set_yticklabels([' 0',' 300',' 600',' 900']) ax.yaxis.set_tick_params(labelsize=13,rotation=20) ax.set_thetagrids([180, 225, 270, 315, 360], labels=['-90', '-45', '0', '45','90'], \ rotation=0,size=13) ax.legend(loc='upper center',bbox_to_anchor=(0.5,0)) ax.text(0,1.16,'d',fontsize=18,transform=ax.transAxes,ha='left',va='top') ax.text(0.2,1.02,'Diffuse Radiance [W m$^{-2}$ sr$^{-1}$]',fontsize=13, \ transform=ax.transAxes,ha='center',va='top') ax.text(0.5,0.1,'Incidence angle [$^\circ$]',fontsize=13,transform=ax.transAxes, \ ha='center',va='top') #if letters: # ob=matplotlib.offsetbox.AnchoredText(letter_list[3],2,borderpad=0., pad=0.2, \ # prop=dict(size=18)) # ob.patch.set(linewidth=1.5) # ax.add_artist(ob) plt.subplot(2,3,5) ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=4)) ax.yaxis.set_major_locator(MaxNLocator(nbins=6)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.plot(data[0,4,0,:,2],data[0,4,0,:,0],lw=3, c='C2',label='$F_{\downarrow}$') ax.plot(data[0,4,0,:,3],data[0,4,0,:,0],'--',lw=3, c='C3',label=r'$F_{\uparrow}$') print('Cloud albedo:',data[0,4,0,-3,3]/data[0,4,0,-3,2]) ax.legend(bbox_to_anchor=(0.6, 0.715)) ax.set_xlabel('Irradiance [W m$^{-2}$]') ax.set_ylabel('Altitude [km]') ax.set_ylim(0,30) ax.set_xlim(0,1400) ax.set_title('Without SAI') ax.text(0.6, 0.15, 'Cloud\nalbedo'+ \ '\n={:.1f}'.format((data[0,4,0,-3,3]/data[0,4,0,-3,2])*100)+'%', \ transform=ax.transAxes, fontsize=12, verticalalignment='bottom', ha='center', \ bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) if letters: ob=matplotlib.offsetbox.AnchoredText(letter_list[4],2,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) plt.subplot(2,3,6) ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=4)) ax.yaxis.set_major_locator(MaxNLocator(nbins=6)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.plot(data[2,4,0,:,2],data[2,4,0,:,0],lw=3, c='C2',label='$F_{\downarrow}$') ax.plot(data[2,4,0,:,3],data[2,4,0,:,0],'--',lw=3, c='C3',label=r'$F_{\uparrow}$') print('Cloud albedo SAI:',data[2,4,0,-3,3]/data[2,4,0,-3,2]) ax.legend(bbox_to_anchor=(0.6, 0.715)) ax.set_xlabel('Irradiance [W m$^{-2}$]') ax.set_ylabel('Altitude [km]') ax.set_ylim(0,30) ax.set_xlim(0,1400) ax.set_title('With SAI') ax.text(0.54, 0.15, 'Cloud\nalbedo'+ \ '\n={:.1f}'.format((data[2,4,0,-3,3]/data[2,4,0,-3,2])*100)+'%', \ transform=ax.transAxes, fontsize=12, verticalalignment='bottom', ha='center', \ bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5)) if letters: ob=matplotlib.offsetbox.AnchoredText(letter_list[5],2,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) plt.tight_layout() if not showPlots: if letters: plt.savefig(figDir+'Fig_02'+ftype,dpi=dpi,bbox_inches='tight', \ pad_inches = 0.1) if not letters: plt.savefig(figDir+'withoutLetters/Fig_02'+ftype,dpi=dpi, \ bbox_inches='tight',pad_inches = 0.1) ################################################################################################ # FIGURE 3 # ################################################################################################ if 3 in figs: # read SBDART output data fname='../sbdart_output/cloud_albedo_sai_noAbs.txt' f=open(fname,'r') data=f.readlines() f.close() del data[::34] # remove headers # define conditions used in SBDART aod=[0,0.5,1,1.5,2] tau=[0,1,3,5,10,20,50] sza=[0,20,40,60,80] naod=len(aod) ntau=len(tau) nsza=len(sza) nlev=33 nvar=7 # read output data into array (naod,ntau,nsza,nlev,nvar) data=np.array([data[i].rstrip('\n').split() for i in \ range(len(data))]).astype(np.float).reshape(naod,ntau,nsza,nlev,nvar) # read SBDART output data fname2= \ '../sbdart_output/cloud_albedo_sai_noAbs_scaledcod.txt' # scaled cod values: 0.619 2.23 3.99 8.51 17.52 43.95 f=open(fname2,'r') data2=f.readlines() f.close() del data2[::34] # remove headers # read output data into array (naod,ntau,nsza,nlev,nvar) data2=np.array([data2[i].rstrip('\n').split() for i in \ range(len(data2))]).astype(np.float).reshape(ntau,nlev,nvar) # calculate TOA radiative effect of aerosols and diffusion brightning for all cloud optical # depths toa_re_tot=data[0,:,0,0,3]-data[2,:,0,0,3] toa_re_db=data2[:,0,3]-data[2,:,0,0,3] toa_re_aer=toa_re_tot-toa_re_db # set plot defaults matplotlib.rcParams['axes.linewidth']=1.5 matplotlib.rcParams['font.size']=14 plt.figure(figsize=(6,4)) ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=6)) ax.yaxis.set_major_locator(MaxNLocator(nbins=8)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.set_xlabel(r'Cloud optical depth') ax.set_ylabel('TOA radiative effect [W m$^{-2}$]') ax.plot(tau,toa_re_aer,'--',c='C0',lw=3,label='Stratospheric aerosol reflection') ax.plot(tau,toa_re_db,':',c='C1',lw=3,label='Diffusion-brightening') ax.plot(tau,toa_re_tot,c='C2',lw=3,label='Total') ax.legend(fontsize=13) ax.axhline(0,c='k') ax.grid() plt.tight_layout() if not showPlots: if letters: plt.savefig(figDir+'Fig_03'+ftype,dpi=dpi,bbox_inches='tight', \ pad_inches = 0.1) if not letters: plt.savefig(figDir+'withoutLetters/Fig_03'+ftype,dpi=dpi, \ bbox_inches='tight',pad_inches = 0.1) ################################################################################################ # FIGURE 4 # ################################################################################################ if 4 in figs: # read SBDART output data fname='../sbdart_output/cloud_albedo_sai.txt' f=open(fname,'r') data=f.readlines() f.close() del data[::34] # remove headers # define conditions used in SBDART aod=[0,0.5,1,1.5,2] tau=[0,1,3,5,10,20,50] sza=[0,20,40,60,80] naod=len(aod) ntau=len(tau) nsza=len(sza) nlev=33 nvar=7 # read output data into array (naod,ntau,nsza,nlev,nvar) data=np.array([data[i].rstrip('\n').split() for i in \ range(len(data))]).astype(np.float).reshape(naod,ntau,nsza,nlev,nvar) # read SBDART output data - 2 STREAM fname='../sbdart_output/cloud_albedo_sai_2stream.txt' f=open(fname,'r') data_2s=f.readlines() f.close() del data_2s[::34] # remove headers # read output data into array (naod,ntau,nsza,nlev,nvar) data_2s=np.array([data_2s[i].rstrip('\n').split() for i in \ range(len(data_2s))]).astype(np.float).reshape(naod,ntau,nsza,nlev,nvar) # set plot defaults matplotlib.rcParams['axes.linewidth']=1.5 matplotlib.rcParams['font.size']=14 plt.figure(figsize=(5,10)) plt.subplot(3,1,1) ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=5)) ax.yaxis.set_major_locator(MaxNLocator(nbins=8)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.plot(aod,100/(data[0,4,0,-3,3]/data[0,4,0,-3,2])*(data[:,4,0,-3,3]/ \ data[:,4,0,-3,2])-100,lw=3,label='20 streams') ax.plot(aod,100/(data_2s[0,4,0,-3,3]/data_2s[0,4,0,-3,2])*(data_2s[:,4,0,-3,3]/ \ data_2s[:,4,0,-3,2])-100,'--',lw=3,label='2 streams') ax.set_xlabel(r'Stratospheric aerosol optical depth') ax.set_ylabel('Cloud albedo\nenhancement [%]') ax.plot(aod[2],100/(data[0,4,0,-3,3]/data[0,4,0,-3,2])*(data[2,4,0,-3,3]/ \ data[2,4,0,-3,2])-100,marker='*',ms=10,c='C0',zorder=999,markeredgecolor='k') ax.legend() #ax.set_ylim(0,30) #ax.set_xlim(0,1400) ax.grid() if letters: ob=matplotlib.offsetbox.AnchoredText('a',2,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) plt.subplot(3,1,2) ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=6)) ax.yaxis.set_major_locator(MaxNLocator(nbins=6)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.plot(tau[1:],100/(data[0,1:,0,-3,3]/data[0,1:,0,-3,2])*(data[2,1:,0,-3,3]/ \ data[2,1:,0,-3,2])-100,lw=3,label='20 streams') ax.plot(tau[1:],100/(data_2s[0,1:,0,-3,3]/data_2s[0,1:,0,-3,2])*(data_2s[2,1:,0,-3,3]/ \ data_2s[2,1:,0,-3,2])-100,'--',lw=3,label='2 streams') ax.set_xlabel(r'Cloud optical depth') ax.set_ylabel('Cloud albedo\nenhancement [%]') ax.plot(tau[4],100/(data[0,4,0,-3,3]/data[0,4,0,-3,2])*(data[2,4,0,-3,3]/ \ data[2,4,0,-3,2])-100,marker='*',ms=10,c='C0',zorder=999,markeredgecolor='k') ax.legend(loc='upper center') ax.set_ylim(0,30) #ax.set_xlim(0,1400) ax.grid() if letters: ob=matplotlib.offsetbox.AnchoredText('b',1,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) plt.subplot(3,1,3) ax=plt.gca() ax.xaxis.set_major_locator(MaxNLocator(nbins=7)) ax.yaxis.set_major_locator(MaxNLocator(nbins=8)) ax.xaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_minor_locator(AutoMinorLocator()) ax.yaxis.set_ticks_position('both') ax.xaxis.set_ticks_position('both') ax.tick_params(which='major',direction='out', length=4, width=1.5) ax.tick_params(which='minor',direction='out', length=2, width=1.5) ax.plot(sza[:-1],100/(data[0,4,:-1,-3,3]/data[0,4,:-1,-3,2])*(data[2,4,:-1,-3,3]/ \ data[2,4,:-1,-3,2])-100,lw=3,label='20 streams') ax.plot(sza[:-1],100/(data_2s[0,4,:-1,-3,3]/data_2s[0,4,:-1,-3,2])*(data_2s[2,4,:-1,-3,3]/ \ data_2s[2,4,:-1,-3,2])-100,'--',lw=3,label='2 streams') ax.set_xlabel('Solar zenith angle [$^{\circ}$]') ax.set_ylabel('Cloud albedo\nenhancement [%]') ax.plot(sza[0],100/(data[0,4,0,-3,3]/data[0,4,0,-3,2])*(data[2,4,0,-3,3]/ \ data[2,4,0,-3,2])-100,marker='*',ms=10,c='C0',zorder=999,markeredgecolor='k') ax.legend() #ax.set_ylim(0,30) #ax.set_xlim(0,1400) ax.axhline(0,c='k') ax.grid() if letters: ob=matplotlib.offsetbox.AnchoredText('c',1,borderpad=0., pad=0.2, \ prop=dict(size=18)) ob.patch.set(linewidth=1.5) ax.add_artist(ob) plt.tight_layout() if not showPlots: if letters: plt.savefig(figDir+'Fig_04'+ftype,dpi=dpi,bbox_inches='tight', \ pad_inches = 0.1) if not letters: plt.savefig(figDir+'withoutLetters/Fig_04'+ftype,dpi=dpi, \ bbox_inches='tight',pad_inches = 0.1) if showPlots: plt.show() if __name__ == '__main__': main()