begin

filenames=(/"BG-1pctCO2","+3.2 K, high sliding", "+3.2 K, no isostasy","+3.2 K","+3.8 K, high sliding", "+3.8 K, no isostasy","+3.8 K"/)

filelist=systemfunc("ls ./data/scalars_runs_sensitivity/*e21.B1850G.f09_g17_gl4.CMIP6*cism.scalars.nc")
nfiles=dimsizes(filelist)
;print(nfiles)

;read first all times
fall=addfiles(filelist,"r")
tmin=min(fall[:]->times)
tmax=max(fall[:]->times)
time=ispan(0,tmax,1)
nt=dimsizes(time)
;print(nt)

;now read all the variables
area=new((/nfiles,nt/),double)
areapct=new((/nfiles,nt/),double)
volGt=new((/nfiles,nt/),double)
volSLE=new((/nfiles,nt/),double)
SLR=new((/nfiles,nt/),double)
SMBf=new((/2*nfiles,nt/),double)
BMBf=new((/2*nfiles,nt/),double)
Cf=new((/2*nfiles,nt/),double)
MBf=new((/2*nfiles,nt/),double)


do f=0,nfiles-1
timecheck=fall[f]->times
dimtimecheck=dimsizes(timecheck)
do t=1,dimtimecheck-1
	if(timecheck(t-1) .eq. timecheck(t)) then
		print("PROBLEM: TIMES ARE NON-MONOTONIC")
		print("Casename: +"+filenames(f))
		print("Year "+timecheck(t-1)+" repeated twice")
	end if
end do
delete(timecheck)
delete(dimtimecheck)
end do

do f=0,nfiles-1
	t0=min(fall[f]->times)
	t1=max(fall[f]->times)
	if(f .eq. 0) then
		bgt0=t0
		bgt1=t1
	end if
	print("Simulation "+filenames(f)+": starting at time "+t0)
	y0=t0-t0
	y1=t1-t0
 	;print(t0)
	;print(t1)
	;print(y0)
	;print(y1)
	if(f .ge. 1) then
	SMBf(f,bgt0:t0)=SMBf(0,bgt0:t0)
	MBf(f,bgt0:t0)=MBf(0,bgt0:t0)
	Cf(f,bgt0:t0)=Cf(0,bgt0:t0)
	end if	
	area(f,t0+1:t1-1)=fall[f]->area(y0+1:y1-1)
	volGt(f,t0+1:t1-1)=fall[f]->volGt(y0+1:y1-1)
	volSLE(f,t0+1:t1-1)=fall[f]->volSLE(y0+1:y1-1)
	SMBf(f,t0+1:t1-1)=fall[f]->SMBf(y0+1:y1-1)
	BMBf(f,t0+1:t1-1)=fall[f]->BMBf(y0+1:y1-1)
	Cf(f,t0+1:t1-1)=fall[f]->Cf(y0+1:y1-1)
	MBf(f,t0+1:t1-1)=fall[f]->MBf(y0+1:y1-1)
	do t=t0+1,t1-1
		if(SMBf(f,t) .eq. 0 .or. ismissing(SMBf(f,t))) then
			SMBf(f,t) = SMBf(f,t+1)
		end if
		if(BMBf(f,t) .eq. 0) then
			BMBf(f,t) = BMBf(f,t+1)
		end if
		if(MBf(f,t) .eq. 0) then
			MBf(f,t) = MBf(f,t+1)
		end if
		if(Cf(f,t) .eq. 0) then
			Cf(f,t) = Cf(f,t+1)
		end if
	end do
	delete(t0)
	delete(t1)
	delete(y0)
	delete(y1)
end do
	
	volGt=where(volGt .gt. 0.,volGt,volGt@_FillValue)
	volSLE=where(volGt .gt. 0.,volSLE,volGt@_FillValue)	
	;area=where(area .gt. 0., area,area@_FillValue)
	SMBf=where(SMBf .ne. 0.,SMBf,SMBf@_FillValue)
	MBf=where(MBf .ne. 0.,MBf,MBf@_FillValue)
	BMBf=where(BMBf .ne. 0.,BMBf,BMBf@_FillValue)
	Cf=where(Cf .ne. 0.,Cf,Cf@_FillValue)
	areapct=area/(area(0,2))*100.
	SLR=3.2347184506371*1.e15*917.*1.e-15/361.-volSLE
	SLR=where(volSLE .gt. 0.,SLR,SLR@_FillValue)

forcinglength=(/20,23,23,23,23,23,23/)
	
do f=nfiles,2*nfiles-1
	SMBf(f,:)=runave(SMBf(f-nfiles,:),forcinglength(f-nfiles),0)
	BMBf(f,:)=runave(BMBf(f-nfiles,:),forcinglength(f-nfiles),0)
	Cf(f,:)=runave(Cf(f-nfiles,:),forcinglength(f-nfiles),0)
	MBf(f,:)=runave(MBf(f-nfiles,:),forcinglength(f-nfiles),0)
end do

;let's plot area/volume/SMB/BMB/Cf/MB 

wtype          = "png"
wtype@wkWidth  = 2000 
wtype@wkHeight = 4000
plotname="./result/timeseries.sensitivity.long"
wks = gsn_open_wks(wtype,plotname) 

plot=new(5,graphic)
minx=500
maxx=40000 ;20000 ;40000
majorlab = ispan(0,maxx,5000) 
minorlab = ispan(0,maxx,1000)
thin=1.5
medium=3
thick=6
linethick = (/thin,medium,thick,thin,medium,thick/)
dashpattern=(/0,15,0,0,15,0/)
linecolor = (/"darkgreen","darkgreen","darkgreen","orange3","orange3","orange3"/)

linethick1=linethick
linecolor1=linecolor
dashpattern1=dashpattern
;------------------------------------------------------------
;FIRST plot: VOLUME
;------------------------------------------------------------

;left Y axis 

;do not plot, panel size
res4a                       = True          
res4a@gsnDraw               = False
res4a@gsnFrame              = False
res4a@vpWidthF              = 0.70       
res4a@vpHeightF             = 0.175

;line thickness, color, pattern
res4a@xyLineThicknesses = linethick1
res4a@xyLineColors      = linecolor1
res4a@xyDashPatterns    = dashpattern1

;X/Y axes:
res4a@tmEqualizeXYSizes    = True
res4a@tmYLLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmYLLabelFontColor    = "black"
res4a@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmXTLabelFontColor    = "black"

res4a@tmYUseLeft           = True
res4a@tmXUseBottom         = False

res4a@tmXBOn               = False
res4a@tmXBLabelsOn	   = False
res4a@tmXTOn               = True
res4a@tmXTLabelsOn	   = True 

res4a@tmYLOn                 = True
res4a@tmYLLabelsOn           = True
res4a@tmYROn                 = False
res4a@tmYRLabelsOn           = False

res4a@tmYLBorderOn           = True
res4a@tmYRBorderOn           = False
res4a@tmXBBorderOn           = False
res4a@tmXTBorderOn           = True

res4a@tiXAxisString         = "Simulation time (kyr)"
res4a@tiXAxisFontColor       = "black"
res4a@tiXAxisFontHeightF     = 0.015
res4a@tiXAxisSide = "Top"
res4a@tiYAxisString          = "GrIS volume (M Gt)"
res4a@tiYAxisFontColor       = "black"
res4a@tiYAxisFontHeightF     = 0.015
res4a@tiYAxisSide = "Left"

res4a@trXMinF  = minx
res4a@trXMaxF  = maxx 
ymin=0.
ymax=3.
res4a@trYMinF  = ymin-(ymax-ymin)*10./100.               
res4a@trYMaxF  = ymax+(ymax-ymin)*0./100.

res4a@tmYLMode          = "Explicit"   
res4a@tmYLValues        = fspan(0,3,4)
res4a@tmYLLabels        = res4a@tmYLValues  
res4a@tmYLMinorValues   = fspan(0,3,7)

res4a@tmXTMode          = "Explicit"                ; explicit labels
res4a@tmXTValues        = majorlab 
res4a@tmXTLabels        = majorlab/1000.         ; make strings
res4a@tmXTMinorValues   = minorlab

;---For higher quality large images
res4a@tiYAxisFontThicknessF  = 4.
res4a@tmXTMajorThicknessF    = 4.
res4a@tmXTMinorThicknessF    = 4.
res4a@tmXTLabelFontThicknessF= 4.
res4a@tmYLMajorThicknessF    = 4.
res4a@tmYLMinorThicknessF    = 4.
res4a@tmYLLabelFontThicknessF= 4.
res4a@tmBorderThicknessF     = 4.

;right Y axis

;do not plot, panel size
res4b                       = True          
res4b@gsnDraw               = False
res4b@gsnFrame              = False
res4b@vpWidthF              = 0.70       
res4b@vpHeightF             = 0.175

;line thickness, color, pattern
res4b@xyLineThicknesses = res4a@xyLineThicknesses
res4b@xyLineColors      = res4a@xyLineColors 
res4b@xyDashPatterns    = res4a@xyDashPatterns

;X/Y axes:
res4b@tmEqualizeXYSizes    = True
res4b@tmYRLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmYRLabelFontColor    = "black"
res4b@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmXTLabelFontColor    = "black"

res4b@tmYUseLeft           = False
res4b@tmXUseBottom         = False

res4b@tmXBOn               = False
res4b@tmXBLabelsOn	   = False
res4b@tmXTOn               = True
res4b@tmXTLabelsOn	   = True 

res4b@tmYLOn                 = False
res4b@tmYLLabelsOn           = False
res4b@tmYROn                 = True
res4b@tmYRLabelsOn           = True

res4b@tmYLBorderOn           = False
res4b@tmYRBorderOn           = True
res4b@tmXBBorderOn           = False
res4b@tmXTBorderOn           = True

res4b@tiXAxisString         = "Simulation time (kyr)"
res4b@tiXAxisFontColor       = "black"
res4b@tiXAxisFontHeightF     = 0.015
res4b@tiXAxisSide = "Top"

res4b@tiYAxisString          = "m SLE"
res4b@tiYAxisFontColor       = "black"
res4b@tiYAxisFontHeightF     = 0.015
res4b@tiYAxisSide = "Right"

res4b@trXMinF  = res4a@trXMinF
res4b@trXMaxF  = res4a@trXMaxF
res4b@trYMinF  = res4a@trYMinF*1.e6/1.e3/361                
res4b@trYMaxF  = res4a@trYMaxF*1.e6/1.e3/361

res4b@tmYRMode          = "Explicit"   
res4b@tmYRValues        = fspan(0,8,5)
res4b@tmYRLabels        = res4b@tmYRValues  
res4b@tmYRMinorValues   = fspan(0,8,9)

;res4b@tmYMajorGrid                = True    ; Turn on grid lines
;res4b@tmYMajorGridLineDashPattern = 2       ; Dashed lines
;res4b@tmYMajorGridThicknessF      = 1.

;---For higher quality large images
res4b@tiYAxisFontThicknessF  = 4.
res4b@tmXTMajorThicknessF    = 4.
res4b@tmXTMinorThicknessF    = 4.
res4b@tmXTLabelFontThicknessF= 4.
res4b@tmYRMajorThicknessF    = 4.
res4b@tmYRMinorThicknessF    = 4.
res4b@tmYRLabelFontThicknessF= 4.
res4b@tmBorderThicknessF     = 4.

plot(0) = gsn_csm_xy2(wks,time,volGt(1:,:),volSLE(1:,:),res4a,res4b)

delete(res4a)
delete(res4b)

;------------------------------------------------------------
;Third plot: AREA
;------------------------------------------------------------

;left Y axis 

;do not plot, panel size
res4a                       = True          
res4a@gsnDraw               = False
res4a@gsnFrame              = False
res4a@vpWidthF              = 0.70       
res4a@vpHeightF             = 0.175

;line thickness, color, pattern
res4a@xyLineThicknesses = linethick1
res4a@xyLineColors      = linecolor1
res4a@xyDashPatterns    = dashpattern1

;X/Y axes:
res4a@tmEqualizeXYSizes    = True
res4a@tmYLLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmYLLabelFontColor    = "black"
res4a@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmXTLabelFontColor    = "black"

res4a@tmYUseLeft           = True
res4a@tmXUseBottom         = False 

res4a@tmXBOn               = False 
res4a@tmXBLabelsOn	   = False
res4a@tmXTOn               = True 
res4a@tmXTLabelsOn	   = True

res4a@tmYLOn                 = True
res4a@tmYLLabelsOn           = True
res4a@tmYROn                 = False
res4a@tmYRLabelsOn           = False

res4a@tmYLBorderOn           = True
res4a@tmYRBorderOn           = False
res4a@tmXBBorderOn           = False 
res4a@tmXTBorderOn           = True

res4a@tiXAxisString         = "Simulation time (kyr)"
res4a@tiXAxisFontColor       = "black"
res4a@tiXAxisFontHeightF     = 0.015
res4a@tiXAxisSide = "Top"

res4a@tiYAxisString          = "GrIS area (Mkm~S~2~N~)"
res4a@tiYAxisFontColor       = "black"
res4a@tiYAxisFontHeightF     = 0.015
res4a@tiYAxisSide = "Left"

res4a@trXMinF  = minx 
res4a@trXMaxF  = maxx
ymin=0.
ymax=2.
res4a@trYMinF  = ymin-(ymax-ymin)*0./100.               
res4a@trYMaxF  = ymax+(ymax-ymin)*10./100.

res4a@tmYLMode          = "Explicit"   
res4a@tmYLValues        = fspan(0,2,3)
res4a@tmYLLabels        = res4a@tmYLValues  
res4a@tmYLMinorValues   = fspan(0,2,5)

res4a@tmXTMode          = "Explicit"                ; explicit labels
res4a@tmXTValues        = majorlab 
res4a@tmXTLabels        = majorlab/1000.         ; make strings
res4a@tmXTMinorValues   = minorlab

;---For higher quality large images
res4a@tiYAxisFontThicknessF  = 4.
res4a@tmXTMajorThicknessF    = 4.
res4a@tmXTMinorThicknessF    = 4.
res4a@tmXTLabelFontThicknessF= 4.
res4a@tmYLMajorThicknessF    = 4.
res4a@tmYLMinorThicknessF    = 4.
res4a@tmYLLabelFontThicknessF= 4.
res4a@tmBorderThicknessF     = 4.

;right Y axis

;do not plot, panel size
res4b                       = True          
res4b@gsnDraw               = False
res4b@gsnFrame              = False
res4b@vpWidthF              = 0.70       
res4b@vpHeightF             = 0.175

;line thickness, color, pattern
res4b@xyLineThicknesses = res4a@xyLineThicknesses
res4b@xyLineColors      = res4a@xyLineColors 
res4b@xyDashPatterns    = res4a@xyDashPatterns

;X/Y axes:
res4b@tmEqualizeXYSizes    = True
res4b@tmYRLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmYRLabelFontColor    = "black"
res4b@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmXTLabelFontColor    = "black"

res4b@tmYUseLeft           = False
res4b@tmXUseBottom         = False

res4b@tmXBOn               = False
res4b@tmXBLabelsOn	   = False
res4b@tmXTOn               = True
res4b@tmXTLabelsOn	   = True 

res4b@tmYLOn                 = False
res4b@tmYLLabelsOn           = False
res4b@tmYROn                 = True
res4b@tmYRLabelsOn           = True

res4b@tmYLBorderOn           = False
res4b@tmYRBorderOn           = True
res4b@tmXBBorderOn           = False
res4b@tmXTBorderOn           = True

res4b@tiXAxisString         = "Simulation time (kyr)"
res4b@tiXAxisFontColor       = "black"
res4b@tiXAxisFontHeightF     = 0.015
res4b@tiXAxisSide = "Top"

res4b@tiYAxisString          = "GrIS area (%)"
res4b@tiYAxisFontColor       = "black"
res4b@tiYAxisFontHeightF     = 0.015
res4b@tiYAxisSide = "Right"

res4b@trXMinF  = res4a@trXMinF
res4b@trXMaxF  = res4a@trXMaxF
res4b@trYMinF  = res4a@trYMinF/(area(0,2))*100.                
res4b@trYMaxF  = res4a@trYMaxF/(area(0,2))*100.

res4b@tmYRMode          = "Explicit"   
res4b@tmYRValues        = fspan(0,100,6)
res4b@tmYRLabels        = res4b@tmYRValues  
res4b@tmYRMinorValues   = fspan(0,100,11)


;---For higher quality large images
res4b@tiYAxisFontThicknessF  = 4.
res4b@tmXTMajorThicknessF    = 4.
res4b@tmXTMinorThicknessF    = 4.
res4b@tmXTLabelFontThicknessF= 4.
res4b@tmYRMajorThicknessF    = 4.
res4b@tmYRMinorThicknessF    = 4.
res4b@tmYRLabelFontThicknessF= 4.
res4b@tmBorderThicknessF     = 4.

plot(1) = gsn_csm_xy2(wks,time,area(1:,:),areapct(1:,:),res4a,res4b)

delete(res4a)
delete(res4b)

;------------------------------------------------------------
;Fourth plot: SMB
;------------------------------------------------------------

;left Y axis 

;do not plot, panel size
res4a                       = True          
res4a@gsnDraw               = False
res4a@gsnFrame              = False
res4a@vpWidthF              = 0.7       
res4a@vpHeightF             = 0.4

;line thickness, color, pattern
res4a@xyLineThicknesses = linethick
res4a@xyLineColors      = linecolor 
res4a@xyDashPatterns    = dashpattern

;X/Y axes:
res4a@tmEqualizeXYSizes    = True
res4a@tmYLLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmYLLabelFontColor    = "black"
res4a@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmXTLabelFontColor    = "black"

res4a@tmYUseLeft           = True
res4a@tmXUseBottom         = True

res4a@tmXBOn               = True 
res4a@tmXBLabelsOn	   = True
res4a@tmXTOn               = False 
res4a@tmXTLabelsOn	   = False

res4a@tmYLOn                 = True
res4a@tmYLLabelsOn           = True
res4a@tmYROn                 = False
res4a@tmYRLabelsOn           = False

res4a@tmYLBorderOn           = True
res4a@tmYRBorderOn           = True
res4a@tmXBBorderOn           = True
res4a@tmXTBorderOn           = True

res4a@tiXAxisString         = "Simulation time (kyr)"
res4a@tiXAxisFontColor       = "black"
res4a@tiXAxisFontHeightF     = 0.015
res4a@tiXAxisSide = "Bottom"

res4a@tiYAxisString          = "SMB"
res4a@tiYAxisFontColor       = "black"
res4a@tiYAxisFontHeightF     = 0.015
res4a@tiYAxisSide = "Left"

res4a@trXMinF  = minx 
res4a@trXMaxF  = maxx
ymin=-500
ymax=300
res4a@trYMinF  = ymin-(ymax-ymin)*11./100.               
res4a@trYMaxF  = ymax+(ymax-ymin)*11./100.

res4a@tmYLMode          = "Explicit"   
res4a@tmYLValues        = ispan(-500,300,100)
res4a@tmYLLabels        = res4a@tmYLValues 
res4a@tmYLMinorValues   = ispan(-500,300,50)

res4a@gsnYRefLine = 0
res4a@gsnYRefLineThicknessF = 2.
res4a@gsnYRefLineColor = "black"

;---For higher quality large images
res4a@tiYAxisFontThicknessF  = 4.
res4a@tmXTMajorThicknessF    = 4.
res4a@tmXTMinorThicknessF    = 4.
res4a@tmXTLabelFontThicknessF= 4.
res4a@tmYLMajorThicknessF    = 4.
res4a@tmYLMinorThicknessF    = 4.
res4a@tmYLLabelFontThicknessF= 4.
res4a@tmBorderThicknessF     = 4.

;right Y axis

;do not plot, panel size
res4b                       = True          
res4b@gsnDraw               = False
res4b@gsnFrame              = False
res4b@vpWidthF              = 0.70       
res4b@vpHeightF             = 0.4

;line thickness, color, pattern
res4b@xyLineThicknesses = res4a@xyLineThicknesses
res4b@xyLineColors      = res4a@xyLineColors 
res4b@xyDashPatterns    = res4a@xyDashPatterns

;X/Y axes:
res4b@tmEqualizeXYSizes    = True
res4b@tmYRLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmYRLabelFontColor    = "black"
res4b@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmXTLabelFontColor    = "black"

res4b@tmYUseLeft           = False
res4b@tmXUseBottom         = False

res4b@tmXBOn               = False
res4b@tmXBLabelsOn	   = False
res4b@tmXTOn               = True
res4b@tmXTLabelsOn	   = True 

res4b@tmYLOn                 = False
res4b@tmYLLabelsOn           = False
res4b@tmYROn                 = True
res4b@tmYRLabelsOn           = True

res4b@tmYLBorderOn           = False
res4b@tmYRBorderOn           = True
res4b@tmXBBorderOn           = False
res4b@tmXTBorderOn           = True

res4b@tiXAxisString         = "Simulation time (kyr)"
res4b@tiXAxisFontColor       = "black"
res4b@tiXAxisFontHeightF     = 0.015
res4b@tiXAxisSide = "Bottom"

res4b@tiYAxisString          = "mm/yr"
res4b@tiYAxisFontColor       = "black"
res4b@tiYAxisFontHeightF     = 0.015
res4b@tiYAxisSide = "Right"

res4b@trXMinF  = res4a@trXMinF
res4b@trXMaxF  = res4a@trXMaxF
res4b@trYMinF  = res4a@trYMinF/361.                
res4b@trYMaxF  = res4a@trYMaxF/361.

res4b@tmYRMode          = "Explicit"   
res4b@tmYRValues        = ispan(-14,8,2)/10.
res4b@tmYRLabels        = res4b@tmYRValues  
res4b@tmYRMinorValues   = ispan(-14,8,1)/10.


;---For higher quality large images
res4b@tiYAxisFontThicknessF  = 4.
res4b@tmXTMajorThicknessF    = 4.
res4b@tmXTMinorThicknessF    = 4.
res4b@tmXTLabelFontThicknessF= 4.
res4b@tmYRMajorThicknessF    = 4.
res4b@tmYRMinorThicknessF    = 4.
res4b@tmYRLabelFontThicknessF= 4.
res4b@tmBorderThicknessF     = 4.


plot(2) = gsn_csm_xy2(wks,time,SMBf(nfiles+1:,:),SMBf(nfiles+1:,:)/361.,res4a,res4b)

delete(res4a)
delete(res4b)

;------------------------------------------------------------
;Fifth plot: Calving
;------------------------------------------------------------

;left Y axis 

;do not plot, panel size
res4a                       = True          
res4a@gsnDraw               = False
res4a@gsnFrame              = False
res4a@vpWidthF              = 0.70       
res4a@vpHeightF             = 0.175

;line thickness, color, pattern
res4a@xyLineThicknesses = linethick
res4a@xyLineColors      = linecolor 
res4a@xyDashPatterns    = dashpattern

;X/Y axes:
res4a@tmEqualizeXYSizes    = True
res4a@tmYLLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmYLLabelFontColor    = "black"
res4a@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmXTLabelFontColor    = "black"

res4a@tmYUseLeft           = True
res4a@tmXUseBottom         = True

res4a@tmXBOn               = True 
res4a@tmXBLabelsOn	   = True
res4a@tmXTOn               = False 
res4a@tmXTLabelsOn	   = False

res4a@tmYLOn                 = True
res4a@tmYLLabelsOn           = True
res4a@tmYROn                 = False
res4a@tmYRLabelsOn           = False

res4a@tmYLBorderOn           = True
res4a@tmYRBorderOn           = True
res4a@tmXBBorderOn           = True
res4a@tmXTBorderOn           = True

res4a@tiXAxisString         = "Simulation time (kyr)"
res4a@tiXAxisFontColor       = "black"
res4a@tiXAxisFontHeightF     = 0.015
res4a@tiXAxisSide = "Bottom"

res4a@tiYAxisString          = "Ice discharge (Gt/yr)"
res4a@tiYAxisFontColor       = "black"
res4a@tiYAxisFontHeightF     = 0.015
res4a@tiYAxisSide = "Left"

res4a@trXMinF  = minx 
res4a@trXMaxF  = maxx
ymin=-400
ymax=0
res4a@trYMinF  = ymin-(ymax-ymin)*12./100.               
res4a@trYMaxF  = ymax+(ymax-ymin)*12./100.

res4a@tmYLMode          = "Explicit"   
res4a@tmYLValues        = ispan(-400,0,100)
res4a@tmYLLabels        = res4a@tmYLValues  
res4a@tmYLMinorValues   = ispan(-400,0,50)

res4a@gsnYRefLine = 0
res4a@gsnYRefLineThicknessF = 2.
res4a@gsnYRefLineColor = "black"

;---For higher quality large images
res4a@tiYAxisFontThicknessF  = 4.
res4a@tmXTMajorThicknessF    = 4.
res4a@tmXTMinorThicknessF    = 4.
res4a@tmXTLabelFontThicknessF= 4.
res4a@tmYLMajorThicknessF    = 4.
res4a@tmYLMinorThicknessF    = 4.
res4a@tmYLLabelFontThicknessF= 4.
res4a@tmBorderThicknessF     = 4.

;right Y axis

;do not plot, panel size
res4b                       = True          
res4b@gsnDraw               = False
res4b@gsnFrame              = False
res4b@vpWidthF              = 0.70       
res4b@vpHeightF             = 0.175

;line thickness, color, pattern
res4b@xyLineThicknesses = res4a@xyLineThicknesses
res4b@xyLineColors      = res4a@xyLineColors 
res4b@xyDashPatterns    = res4a@xyDashPatterns

;X/Y axes:
res4b@tmEqualizeXYSizes    = True
res4b@tmYRLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmYRLabelFontColor    = "black"
res4b@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmXTLabelFontColor    = "black"

res4b@tmYUseLeft           = False
res4b@tmXUseBottom         = False

res4b@tmXBOn               = False
res4b@tmXBLabelsOn	   = False
res4b@tmXTOn               = True
res4b@tmXTLabelsOn	   = True 

res4b@tmYLOn                 = False
res4b@tmYLLabelsOn           = False
res4b@tmYROn                 = True
res4b@tmYRLabelsOn           = True

res4b@tmYLBorderOn           = False
res4b@tmYRBorderOn           = True
res4b@tmXBBorderOn           = False
res4b@tmXTBorderOn           = True

res4b@tiXAxisString         = "Simulation time (kyr)"
res4b@tiXAxisFontColor       = "black"
res4b@tiXAxisFontHeightF     = 0.015
res4b@tiXAxisSide = "Bottom"

res4b@tiYAxisString          = "mm/yr"
res4b@tiYAxisFontColor       = "black"
res4b@tiYAxisFontHeightF     = 0.015
res4b@tiYAxisSide = "Right"

res4b@trXMinF  = res4a@trXMinF
res4b@trXMaxF  = res4a@trXMaxF
res4b@trYMinF  = res4a@trYMinF/361.                
res4b@trYMaxF  = res4a@trYMaxF/361.

res4b@tmYRMode          = "Explicit"   
res4b@tmYRValues        = ispan(-12,0,2)/10.
res4b@tmYRLabels        = res4b@tmYRValues  
res4b@tmYRMinorValues   = ispan(-12,0,1)/10.


;---For higher quality large images
res4b@tiYAxisFontThicknessF  = 4.
res4b@tmXTMajorThicknessF    = 4.
res4b@tmXTMinorThicknessF    = 4.
res4b@tmXTLabelFontThicknessF= 4.
res4b@tmYRMajorThicknessF    = 4.
res4b@tmYRMinorThicknessF    = 4.
res4b@tmYRLabelFontThicknessF= 4.
res4b@tmBorderThicknessF     = 4.


plot(3) = gsn_csm_xy2(wks,time,Cf(nfiles+1:,:),Cf(nfiles+1:,:)/361.,res4a,res4b)

delete(res4a)
delete(res4b)

;------------------------------------------------------------
;Seventh plot: Mass balance
;------------------------------------------------------------

;left Y axis 

;do not plot, panel size
res4a                       = True          
res4a@gsnDraw               = False
res4a@gsnFrame              = False
res4a@vpWidthF              = 0.70       
res4a@vpHeightF             = 0.4

;line thickness, color, pattern
res4a@xyLineThicknesses = linethick
res4a@xyLineColors      = linecolor 
res4a@xyDashPatterns    = dashpattern

;X/Y axes:
res4a@tmEqualizeXYSizes    = True
res4a@tmYLLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmYLLabelFontColor    = "black"
res4a@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4a@tmXTLabelFontColor    = "black"

res4a@tmYUseLeft           = True
res4a@tmXUseBottom         = True

res4a@tmXBOn               = True 
res4a@tmXBLabelsOn	   = True
res4a@tmXTOn               = False 
res4a@tmXTLabelsOn	   = False

res4a@tmYLOn                 = True
res4a@tmYLLabelsOn           = True
res4a@tmYROn                 = False
res4a@tmYRLabelsOn           = False

res4a@tmYLBorderOn           = True
res4a@tmYRBorderOn           = True 
res4a@tmXBBorderOn           = True
res4a@tmXTBorderOn           = True

res4a@tiXAxisString         = "Simulation time (kyr)"
res4a@tiXAxisFontColor       = "black"
res4a@tiXAxisFontHeightF     = 0.015
res4a@tiXAxisSide = "Bottom"

res4a@tiYAxisString          = "Mass balance"
res4a@tiYAxisFontColor       = "black"
res4a@tiYAxisFontHeightF     = 0.015
res4a@tiYAxisSide = "Left"

res4a@trXMinF  = minx 
res4a@trXMaxF  = maxx
ymin=-600
ymax=100
res4a@trYMinF  = ymin-(ymax-ymin)*0./100.               
res4a@trYMaxF  = ymax+(ymax-ymin)*11./100.

res4a@tmYLMode          = "Explicit"   
res4a@tmYLValues        = ispan(-600,100,100)
res4a@tmYLLabels        = res4a@tmYLValues  
res4a@tmYLMinorValues   = ispan(-600,100,50)

res4a@tmXBMode          = "Explicit"                ; explicit labels
res4a@tmXBValues        = majorlab 
res4a@tmXBLabels        = majorlab/1000.        ; make strings
res4a@tmXBMinorValues   = minorlab

res4a@tmGridDrawOrder = "PreDraw"
res4a@xyCurveDrawOrder = "PreDraw"

res4a@gsnYRefLine = 0
res4a@gsnYRefLineThicknessF = 2.
res4a@gsnYRefLineColor = "black"

;---For higher quality large images
res4a@tiYAxisFontThicknessF  = 4.
res4a@tmXBMajorThicknessF    = 4.
res4a@tmXBMinorThicknessF    = 4.
res4a@tmXBLabelFontThicknessF= 4.
res4a@tmYLMajorThicknessF    = 4.
res4a@tmYLMinorThicknessF    = 4.
res4a@tmYLLabelFontThicknessF= 4.
res4a@tmBorderThicknessF     = 4.

;right Y axis

;do not plot, panel size
res4b                       = True          
res4b@gsnDraw               = False
res4b@gsnFrame              = False
res4b@vpWidthF              = 0.70       
res4b@vpHeightF             = 0.175

;line thickness, color, pattern
res4b@xyLineThicknesses = res4a@xyLineThicknesses
res4b@xyLineColors      = res4a@xyLineColors 
res4b@xyDashPatterns    = res4a@xyDashPatterns

;X/Y axes:
res4b@tmEqualizeXYSizes    = True
res4b@tmYRLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmYRLabelFontColor    = "black"
res4b@tmXTLabelFontHeightF  = 0.015      ; resize tick labels
res4b@tmXTLabelFontColor    = "black"

res4b@tmYUseLeft           = False
res4b@tmXUseBottom         = False

res4b@tmXBOn               = False
res4b@tmXBLabelsOn	   = False
res4b@tmXTOn               = True
res4b@tmXTLabelsOn	   = True 

res4b@tmYLOn                 = False
res4b@tmYLLabelsOn           = False
res4b@tmYROn                 = True
res4b@tmYRLabelsOn           = True

res4b@tmYLBorderOn           = False
res4b@tmYRBorderOn           = True
res4b@tmXBBorderOn           = False
res4b@tmXTBorderOn           = True

res4b@tiXAxisString         = "Simulation time (kyr)"
res4b@tiXAxisFontColor       = "black"
res4b@tiXAxisFontHeightF     = 0.015
res4b@tiXAxisSide = "Bottom"

res4b@tiYAxisString          = "mm/yr"
res4b@tiYAxisFontColor       = "black"
res4b@tiYAxisFontHeightF     = 0.015
res4b@tiYAxisSide = "Right"

res4b@trXMinF  = res4a@trXMinF
res4b@trXMaxF  = res4a@trXMaxF
res4b@trYMinF  = res4a@trYMinF/361.                
res4b@trYMaxF  = res4a@trYMaxF/361.

res4b@tmYRMode          = "Explicit"   
res4b@tmYRValues        = ispan(-15,3,3)/10.
res4b@tmYRLabels        = res4b@tmYRValues  
res4b@tmYRMinorValues   = fspan(-1.5,0.3,19)


;---For higher quality large images
res4b@tiYAxisFontThicknessF  = 4.
res4b@tmXTMajorThicknessF    = 4.
res4b@tmXTMinorThicknessF    = 4.
res4b@tmXTLabelFontThicknessF= 4.
res4b@tmYRMajorThicknessF    = 4.
res4b@tmYRMinorThicknessF    = 4.
res4b@tmYRLabelFontThicknessF= 4.
res4b@tmBorderThicknessF     = 4.


plot(4) = gsn_csm_xy2(wks,time,MBf(nfiles+1:,:),MBf(nfiles+1:,:)/361.,res4a,res4b)

;*******************************
;LEGEND
;*******************************

;first add a white box

x0=20000
x1=40000
x2=40000
x3=20000
y0=-2800
y1=-3800
y2=-800
y3=-800
x=(/x0,x1,x2,x3/)
y=(/y0,y1,y2,y3/)

pgres=True
pgres@gsFillColor="White"
pgres@tfPolyDrawOrder = "Postdraw"
pgres@gsFillOpacityF=1.0
box=gsn_add_polygon(wks,plot(4),x,y,pgres)

;now add line and text

res_text                    = True                  ; text mods desired
res_text@txFontHeightF      = 0.014              ; change text size
res_text@txJust             = "CenterLeft"          ; text justification
res_text@txPerimOn             = False
res_text@txBackgroundFillColor = "White"

line=new(nfiles,graphic)
text=new(nfiles,graphic)
do i=1,nfiles-1
res_line1 = True
res_line1@gsLineColor = linecolor1(i-1)
res_line1@gsLineThicknessF = linethick1(i-1)
res_line1@gsLineDashPattern = dashpattern1(i-1)
xl1 = (/23000,26000/)
yl1 = (/-550+70*(i-1),-550+70*(i-1)/)
line(i) = gsn_add_polyline(wks,plot(4),xl1,yl1,res_line1)              ; add polyline
text(i) = gsn_add_text(wks,plot(4),(filenames(i))+" run",xl1(1)+1000,yl1(1),res_text); add text
end do

delete(res4a)

;*******************************
;ATTACH PLOTS
;*******************************

resA                     = True
resA@gsnAttachPlotsXAxis = True 
;resA@gsnAttachBorderOn   = True 
resB                     = True
resB@gsnAttachPlotsXAxis = True
;resB@gsnAttachBorderOn   = True 

attach= gsn_attach_plots(plot(0),plot(1:4),resA,resB)

pres = False
maximize_output(wks,pres)

system("convert -density 300 "+plotname+".png -trim +repage "+plotname+".png")



end