Hi All
I have created 44 plots of CMIP5 model output and I would like to put all my plots obtained into a panel plot of 4 rows by 4 columns but I am not so sure how to go about it. Is there a quick way of doing this? Should I save the output and use that for the panel plot or recode but using 44 input files?
Below is the code of one plot I want in a panel.
;*************************************************
; regline.ncl
;
;*************************************************
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
begin
;************************************************
; Read in Precip Data.
;************************************************
fili = "/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_ACCESS1-3_historical_safrica_climDJF1.nc" ; data
f = addfile (fili , "r") ; add file
lat = f->lat ; get lat
lon = f->lon ; get lon
time = f->time ; get time
;level = f->z ; get level
pr = f->pr ; get precip
printVarSummary(pr)
pr2 = pr(time|0,z|0, {lat|-30:0}, {lon|10:50})
printVarSummary(pr2)
;*************************************************************
;Calculations of max precip for lat and lon values
;**************************************************************
dimpr2 = dimsizes(pr2)
nlat = dimpr2(0)
mlon = dimpr2(1)
pr2MaxLon = new ( mlon, typeof(pr2), pr2@_FillValue)
do ml=0,mlon-1
imax = maxind(pr2(:,ml))
pr2MaxLon(ml) = dble2flt(lat(imax))
end do
print(pr2MaxLon)
print(pr2&lon)
print("-------------------------------")
print("pr2MaxLon: "+pr2&lon+" "+pr2MaxLon)
;Regression Line
rcMaxLon = regline(pr2&lon,pr2MaxLon)
print(rcMaxLon)
print(rcMaxLon@yave)
bMaxLon = rcMaxLon@yintercept
print(bMaxLon)
xMaxLon = pr2&lon
print(xMaxLon)
yMaxLon = rcMaxLon*pr2&lon + bMaxLon
print(yMaxLon)
print("-------------------------------")
print(xMaxLon+" "+yMaxLon)
;************************************************
; create an array to hold both the original data
; and the calculated regression line
;************************************************
data = new ( (/2,dimsizes(pr2MaxLon)/), typeof(pr2MaxLon))
;data(0,:) = pr2MaxLon
; y = mx+b
; m is the slope: rc returned from regline
; b is the y intercept: rc@yave attribute of rc returned from regline
data(1,:) = dble2flt(rcMaxLon)*(dble2flt(xMaxLon)-dble2flt(rcMaxLon@xave)) + dble2flt(rcMaxLon@yave)
;************************************************
; Read in Precip Data with larger domain.
;************************************************
fili2 = "/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_ACCESS1-3_historical_safrica_climDJF1.nc" ; data
f = addfile (fili2 , "r") ; add file
lat1 = f->lat ; get lat
lon1 = f->lon ; get lon
time1 = f->time ; get time
;level1 = f->z ; get level
pr1 = f->pr ; get precip
pr4 = pr1(time|0,z|0,{lat|-40:10},{lon|-10:70})
printVarSummary(pr4)
;************************************************
; Read in Precip Data.
;************************************************
fili = "/mnt/nfs2/geog/ml382/melphd/regressionline/cmap61.nc" ; data
f = addfile (fili , "r") ; add file
lat = f->lat ; get lat
lon = f->lon ; get lon
time = f->time ; get time
;level = f->z ; get level
pr6 = f->precip ; get precip
printVarSummary(pr6)
pr5 = pr6(time|0, {lat|-30:0}, {lon|10:50})
printVarSummary(pr5)
;*************************************************************
;Calculations of max precip for lat and lon values
;**************************************************************
dimpr5 = dimsizes(pr5)
nlat = dimpr2(0)
mlon = dimpr2(1)
pr5MaxLon = new ( mlon, typeof(pr5), pr5@_FillValue)
do ml=0,mlon-1
imax = maxind(pr5(:,ml))
pr5MaxLon(ml) = dble2flt(lat(imax))
end do
print(pr5MaxLon)
print(pr5&lon)
print("-------------------------------")
print("pr5MaxLon: "+pr5&lon+" "+pr5MaxLon)
;Regression Line
regcMaxLon = regline(pr5&lon,pr5MaxLon)
print(regcMaxLon)
print(regcMaxLon@yave)
bMaxLon = regcMaxLon@yintercept
print(bMaxLon)
xMaxLon = pr5&lon
print(xMaxLon)
yMaxLon = regcMaxLon*pr5&lon + bMaxLon
print(yMaxLon)
print("-------------------------------")
print(xMaxLon+" "+yMaxLon)
;************************************************
; create an array to hold both the original data
; and the calculated regression line
;************************************************
datac = new ( (/2,dimsizes(pr5MaxLon)/), typeof(pr5MaxLon))
;datac(0,:) = pr5MaxLon
; y = mx+b
; m is the slope: rc returned from regline
; b is the y intercept: rc@yave attribute of rc returned from regline
datac(1,:) = dble2flt(regcMaxLon)*(dble2flt(xMaxLon)-dble2flt(regcMaxLon@xave)) + dble2flt(regcMaxLon@yave)
;************************************************
; plotting parameters
;************************************************
wks = gsn_open_wks("X11","ACCESS1-3_regline11") ; specifies a plot
gsn_define_colormap(wks,"gsltod") ; choose color map
res = True ; plot mods desired
res@gsnAddCyclic = False
;res@mpOutlineBoundarySets = "national"
res@cnFillOn = True ; color on
res@lbLabelStride = 2 ; every other label
res@gsnSpreadColors = True ; use full range of color map
res@cnLineLabelsOn = False ; no contour line labels
res@lbLabelBarOn = False
res@gsnDraw = False ; do not draw the plot
res@gsnFrame = False
pr4 = pr1(0,0,{-40:10},{-10:70})
dimpr4 = dimsizes(pr4)
latn = dimpr4(0)
lonm = dimpr4(1)
res@mpLimitMode = "Corners" ;
res@mpLeftCornerLonF = lon1(0)
res@mpRightCornerLonF = lon1(lonm-1)
res@mpLeftCornerLatF = lat1(0)
res@mpRightCornerLatF = lat1(latn-1)
res@tmXBMode = "Explicit" ; Define own tick mark labels.
res@tmXBValues = (/ -9.,0.,10.,20.,30.,40.,50.,60.,69 /)
res@tmXBLabels = (/ "10W","0","10E","20E","30E","40E","50E","60E","70E" /)
res@tmYLMode = "Explicit" ; Define own tick mark labels.
res@tmYLValues = (/ -39.,-30.,-20.,-10.,-0.,9 /)
res@tmYLLabels = (/"40S","30S","20S","10S","0","10N" /)
res@vpXF = 0.12 ; default is 0.2 change aspect
res@vpYF = 0.8 ; default is 0.8 ration
res@vpHeightF = 0.4 ; default is 0.6
res@vpWidthF = 0.8 ; default is 0.6
res@cnLevelSelectionMode = "ManualLevels" ; manual levels
res@cnMinLevelValF = 0 ; min level
res@cnMaxLevelValF = 400 ; max level
res@cnLevelSpacingF = 50 ; contour spacing
res@tmXBLabelFontHeightF = 0.02 ; resize tick labels
res@tmYLLabelFontHeightF = 0.02
res@pmLabelBarOrthogonalPosF = .10 ; move whole thing down
;************************************************
; plotting parameters
;************************************************
sres = True ; plot mods desired
sres@gsnMaximize = True ; maximize plot in frame
sres@xyMarkLineModes = (/"Markers","Lines"/) ; choose which have markers
sres@xyMarkers = 16 ; choose type of marker
sres@xyLineColors = (/"blue","black"/)
sres@xyMonoDashPattern = True
sres@gsLineDashPattern = 1
sres@xyDashPattern = 1 ; solid line
sres@xyLineThicknesses = (/1,4/) ; set second line to 2
sres@gsnDraw = False ; do not draw the plot
sres@gsnFrame = False
pr4@long_name = "ACCESS1-3"
pr4@units = "s=-0.17 lat=-11.47"
; reverse the first two colors
setvalues wks
"wkColorMap" : "gsltod"
"wkForegroundColor" : (/0.,0.,0./)
"wkBackgroundColor" : (/1.,1.,1./)
end setvalues
;************************************************
; plotting parameters
;************************************************
cres = True ; plot mods desired
cres@gsnMaximize = True ; maximize plot in frame
cres@xyMarkLineModes = "Lines" ; choose which have markers
;cres@xyMarkers = 16 ; choose type of marker
;cres@xyMarkerSizeF = 0.01 ; Marker size (default 0.01)
cres@xyDashPatterns = 0 ; solid line
cres@xyLineThicknesses = (/1,4/) ; set second line to 2
cres@gsnDraw = False ; do not draw the plot
cres@gsnFrame = False
; reverse the first two colors
setvalues wks
"wkColorMap" : "gsltod"
"wkForegroundColor" : (/0.,0.,0./)
"wkBackgroundColor" : (/1.,1.,1./)
end setvalues
print(rcMaxLon)
plot = gsn_csm_contour_map(wks,pr4({-40:10},{-10:70}),res) ; create plot
plot2 = gsn_csm_xy(wks,pr2&lon,data,sres) ; create plot
plot3 = gsn_csm_xy(wks,pr5&lon,datac,cres) ; create plot
overlay(plot,plot2)
overlay(plot,plot3)
draw(plot)
frame(wks)
end
Below is an attempt of how to go about doing all 44 at once.
;*************************************************
; regline_panel plots.ncl
;
;*************************************************
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
begin
;************************************************
; Read in Precip Data.
;************************************************
;Obs
in = addfile("/mnt/nfs2/geog/ml382/melphd/regressionline/cmap61.nc","r")
prcmap = in->precip(:,:,:)
prcmap = prcmap(time|0, {lat|-30:0}, {lon|10:50})
;Models
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_ACCESS1-0_historical_safrica_climDJF1.nc","r")
pr1 = in->pr(:,0,:,:)
pr1 = pr1(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_ACCESS1-3_historical_safrica_climDJF1.nc","r")
pr2 = in->pr(:,0,:,:)
pr2 = pr2(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_bcc-csm1-1_historical_safrica_climDJF1.nc","r")
pr3 = in->pr(:,0,:,:)
pr3 = pr3(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_bcc-csm1-1-m_historical_safrica_climDJF1.nc,"r")
pr4 = in->pr(:,0,:,:)
pr4 = pr4(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_BNU-ESM_historical_safrica_climDJF1.nc,"r")
pr5 = in->pr(:,0,:,:)
pr5 = pr5(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CanESM2_historical_safrica_climDJF1.nc","r")
pr6 = in->pr(:,0,:,:)
pr6 = pr6(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CCSM4_historical_safrica_climDJF1.nc","r")
pr7 = in->pr(:,0,:,:)
pr7 = pr7(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CESM1-BGC_historical_safrica_climDJF1.nc","r")
pr8 = in->pr(:,0,:,:)
pr8 = pr8(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CESM1-CAM5_historical_safrica_climDJF1.nc","r")
pr9 = in->pr(:,0,:,:)
pr9 = pr9(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CESM1-FASTCHEM_historical_safrica_climDJF1.nc","r")
pr10 = in->pr(:,0,:,:)
pr10 = pr10(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CESM1-WACCM_historical_safrica_climDJF1.nc","r")
pr11 = in->pr(:,0,:,:)
pr11 = pr11(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CMCC-CESM_historical_safrica_climDJF1.nc","r")
pr12 = in->pr(:,0,:,:)
pr12 = pr12(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CMCC-CM_historical_safrica_climDJF1.nc","r")
pr13 = in->pr(:,0,:,:)
pr13 = pr13(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CMCC-CMS_historical_safrica_climDJF1.nc","r")
pr14 = in->pr(:,0,:,:)
pr14 = pr14(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CNRM-CM5_historical_safrica_climDJF1.nc","r")
pr15 = in->pr(:,0,:,:)
pr15 = pr15(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_CSIRO-Mk3-6-0_historical_safrica_climDJF1.nc","r")
pr16 = in->pr(:,0,:,:)
pr16 = pr16(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_EC-EARTH_historical_safrica_climDJF1.nc","r")
pr17 = in->pr(:,0,:,:)
pr17 = pr17(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/FGOALS-g2_climDJFsmall.nc,"r")
pr18 = in->pr(:,:,:)
pr18 = pr18(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_FGOALS-s2_historical_safrica_climDJF1.nc","r")
pr19 = in->pr(:,0,:,:)
pr19 = pr19(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_FIO-ESM_historical_safrica_climDJF1.nc","r")
pr20 = in->pr(:,0,:,:)
pr20 = pr20(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_GFDL-CM3_historical_safrica_climDJF1.nc","r")
pr21 = in->pr(:,0,:,:)
pr21 = pr21(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprreglineGFDLsmalldomain/pr_Amon_GFDL-ESM2G_historical_safrica_climDJF1.nc","r")
pr22 = in->pr(:,:,:)
pr22 = pr22(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprreglineGFDLsmalldomain/pr_Amon_GFDL-ESM2M_historical_safrica_climDJF1.nc","r")
pr23 = in->pr(:,:,:)
pr23 = pr23(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_GISS-E2-H_p1_historical_safrica_climDJF1.nc","r")
pr24 = in->pr(:,0,:,:)
pr24 = pr24(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_GISS-E2-H-CC_p1_historical_safrica_climDJF1.nc","r")
pr25 = in->pr(:,0,:,:)
pr25 = pr25(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_GISS-E2-R_p1_historical_safrica_climDJF1.nc","r")
pr26 = in->pr(:,0,:,:)
pr26 = pr26(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_GISS-E2-R-CC_p1_historical_safrica_climDJF1.nc","r")
pr27 = in->pr(:,0,:,:)
pr27 = pr27(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_HadGEM2-AO_historical_safrica_climDJF1.nc","r")
pr28 = in->pr(:,0,:,:)
pr28 = pr28(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_HadGEM2-CC_historical_safrica_climDJF1.nc","r")
pr29 = in->pr(:,0,:,:)
pr29 = pr29(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_HadGEM2-ES_historical_safrica_climDJF1.nc","r")
pr30 = in->pr(:,0,:,:)
pr30 = pr30(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_inmcm4_historical_safrica_climDJF1.nc","r")
pr31 = in->pr(:,0,:,:)
pr31 = pr31(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_IPSL-CM5A-LR_historical_safrica_climDJF1.nc","r")
pr32 = in->pr(:,0,:,:)
pr32 = pr32(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_IPSL-CM5A-MR_historical_safrica_climDJF1.nc","r")
pr33 = in->pr(:,0,:,:)
pr33 = pr33(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_IPSL-CM5B-LR_historical_safrica_climDJF1.nc","r")
pr34 = in->pr(:,0,:,:)
pr34 = pr34(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_MIROC4h_historical_safrica_climDJF1.nc","r")
pr35 = in->pr(:,0,:,:)
pr35 = pr35(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_MIROC5_historical_safrica_climDJF1.nc","r")
pr36 = in->pr(:,0,:,:)
pr36 = pr36(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_MIROC-ESM_historical_safrica_climDJF1.nc","r")
pr37 = in->pr(:,0,:,:)
pr37 = pr37(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_MIROC-ESM-CHEM_historical_safrica_climDJF1.nc","r")
pr38 = in->pr(:,0,:,:)
pr38 = pr38(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_MPI-ESM-LR_historical_safrica_climDJF1.nc","r")
pr39 = in->pr(:,0,:,:)
pr39 = pr39(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_MPI-ESM-MR_historical_safrica_climDJF1.nc","r")
pr40 = in->pr(:,0,:,:)
pr40 = pr40(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_MPI-ESM-P_historical_safrica_climDJF1.nc","r")
pr41 = in->pr(:,0,:,:)
pr41 = pr41(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_MRI-CGCM3_historical_safrica_climDJF1.nc","r")
pr42 = in->pr(:,0,:,:)
pr42 = pr42(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_NorESM1-M_historical_safrica_climDJF1.nc","r")
pr43 = in->pr(:,0,:,:)
pr43 = pr43(time|0,z|0, {lat|-30:0}, {lon|10:50})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline99/pr_Amon_NorESM1-ME_historical_safrica_climDJF1.nc","r")
pr44 = in->pr(:,0,:,:)
pr44 = pr44(time|0,z|0, {lat|-30:0}, {lon|10:50})
lat =in->lat
time =in->time
;*************************************************************
;Calculations of max precip for lat and lon values
;**************************************************************
number = new(44,float,0)
number = (/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,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44/)
do i=1,number-1
dimpr(i) = dimsizes(pr(i))
nlat = dimpr(i)(0)
mlon = dimpr(i)(1)
pr(i)MaxLon = new ( mlon, typeof(pr(i)), pr(i)@_FillValue)
do ml=0,mlon-1
imax = maxind(pr(i)(:,ml))
pr(i)MaxLon(ml) = dble2flt(lat(imax))
end do
print(pr(i)MaxLon)
print(pr(i)&lon)
print("-------------------------------")
print("pr(i)MaxLon: "+pr(i)&lon+" "+pr(i)MaxLon)
;Regression Line
rcMaxLon(i) = regline(pr(i)&lon,pr(i)MaxLon)
print(rcMaxLon(i))
print(rcMaxLon(i)@yave)
bMaxLon(i) = rcMaxLon(i)@yintercept
print(bMaxLon(i))
xMaxLon(i) = pr(i)&lon
print(xMaxLon(i))
yMaxLon(i) = rcMaxLon(i)*pr(i)&lon + bMaxLon(i)
print(yMaxLon(i)
print("-------------------------------")
print(xMaxLon(i)+" "+yMaxLon(i))
;************************************************
; create an array to hold both the original data
; and the calculated regression line
;************************************************
data = new ( (/2,dimsizes(pr(i)MaxLon)/), typeof(pr(i)MaxLon))
;data(0,:) = pr(i)MaxLon
; y = mx+b
; m is the slope: rc returned from regline
; b is the y intercept: rc@yave attribute of rc returned from regline
data(1,:) = dble2flt(rcMaxLon(i))*(dble2flt(xMaxLon(i))-dble2flt(rcMaxLon(i)@xave)) + dble2flt(rcMaxLon(i)@yave)
end do
;************************************************
; Read in Precip Data with larger domain.
;************************************************
in = addfile("/mnt/nfs2/geog/ml382/melphd/regressionline/cmap5.nc","r")
prcmap1 = in->precip(:,:,:)
prcmap1 = prcmap1(time|0, {lat|-40:10},{lon|-10:70})
;Models
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_ACCESS1-0_historical_safrica_climDJF1.nc","r")
pr1 = in->pr(:,0,:,:)
pr1 = pr1(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_ACCESS1-3_historical_safrica_climDJF1.nc","r")
pr2 = in->pr(:,0,:,:)
pr2 = pr2(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_bcc-csm1-1_historical_safrica_climDJF1.nc","r")
pr3 = in->pr(:,0,:,:)
pr3 = pr3(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_bcc-csm1-1-m_historical_safrica_climDJF1.nc,"r")
pr4 = in->pr(:,0,:,:)
pr4 = pr4(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_BNU-ESM_historical_safrica_climDJF1.nc,"r")
pr5 = in->pr(:,0,:,:)
pr5 = pr5(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CanESM2_historical_safrica_climDJF1.nc","r")
pr6 = in->pr(:,0,:,:)
pr6 = pr6(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CCSM4_historical_safrica_climDJF1.nc","r")
pr7 = in->pr(:,0,:,:)
pr7 = pr7(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CESM1-BGC_historical_safrica_climDJF1.nc","r")
pr8 = in->pr(:,0,:,:)
pr8 = pr8(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CESM1-CAM5_historical_safrica_climDJF1.nc","r")
pr9 = in->pr(:,0,:,:)
pr9 = pr9(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CESM1-FASTCHEM_historical_safrica_climDJF1.nc","r")
pr10 = in->pr(:,0,:,:)
pr10 = pr10(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CESM1-WACCM_historical_safrica_climDJF1.nc","r")
pr11 = in->pr(:,0,:,:)
pr11 = pr11(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CMCC-CESM_historical_safrica_climDJF1.nc","r")
pr12 = in->pr(:,0,:,:)
pr12 = pr12(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CMCC-CM_historical_safrica_climDJF1.nc","r")
pr13 = in->pr(:,0,:,:)
pr13 = pr13(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CMCC-CMS_historical_safrica_climDJF1.nc","r")
pr14 = in->pr(:,0,:,:)
pr14 = pr14(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CNRM-CM5_historical_safrica_climDJF1.nc","r")
pr15 = in->pr(:,0,:,:)
pr15 = pr15(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_CSIRO-Mk3-6-0_historical_safrica_climDJF1.nc","r")
pr16 = in->pr(:,0,:,:)
pr16 = pr16(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_EC-EARTH_historical_safrica_climDJF1.nc","r")
pr17 = in->pr(:,0,:,:)
pr17 = pr17(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/FGOALS-g2_climDJF1.nc,"r")
pr18 = in->pr(:,:,:)
pr18 = pr18(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_FGOALS-s2_historical_safrica_climDJF1.nc","r")
pr19 = in->pr(:,0,:,:)
pr19 = pr19(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_FIO-ESM_historical_safrica_climDJF1.nc","r")
pr20 = in->pr(:,0,:,:)
pr20 = pr20(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_GFDL-CM3_historical_safrica_climDJF1.nc","r")
pr21 = in->pr(:,0,:,:)
pr21 = pr21(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprreglineGFDLlargedomain/pr_Amon_GFDL-ESM2G_historical_safrica_climDJF1.nc","r")
pr22 = in->pr(:,:,:)
pr22 = pr22(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprreglineGFDLlargedomain/pr_Amon_GFDL-ESM2M_historical_safrica_climDJF1.nc","r")
pr23 = in->pr(:,:,:)
pr23 = pr23(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_GISS-E2-H_p1_historical_safrica_climDJF1.nc","r")
pr24 = in->pr(:,0,:,:)
pr24 = pr24(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_GISS-E2-H-CC_p1_historical_safrica_climDJF1.nc","r")
pr25 = in->pr(:,0,:,:)
pr25 = pr25(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_GISS-E2-R_p1_historical_safrica_climDJF1.nc","r")
pr26 = in->pr(:,0,:,:)
pr26 = pr26(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_GISS-E2-R-CC_p1_historical_safrica_climDJF1.nc","r")
pr27 = in->pr(:,0,:,:)
pr27 = pr27(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_HadGEM2-AO_historical_safrica_climDJF1.nc","r")
pr28 = in->pr(:,0,:,:)
pr28 = pr28(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_HadGEM2-CC_historical_safrica_climDJF1.nc","r")
pr29 = in->pr(:,0,:,:)
pr29 = pr29(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_HadGEM2-ES_historical_safrica_climDJF1.nc","r")
pr30 = in->pr(:,0,:,:)
pr30 = pr30(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_inmcm4_historical_safrica_climDJF1.nc","r")
pr31 = in->pr(:,0,:,:)
pr31 = pr31(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_IPSL-CM5A-LR_historical_safrica_climDJF1.nc","r")
pr32 = in->pr(:,0,:,:)
pr32 = pr32(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_IPSL-CM5A-MR_historical_safrica_climDJF1.nc","r")
pr33 = in->pr(:,0,:,:)
pr33 = pr33(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_IPSL-CM5B-LR_historical_safrica_climDJF1.nc","r")
pr34 = in->pr(:,0,:,:)
pr34 = pr34(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_MIROC4h_historical_safrica_climDJF1.nc","r")
pr35 = in->pr(:,0,:,:)
pr35 = pr35(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_MIROC5_historical_safrica_climDJF1.nc","r")
pr36 = in->pr(:,0,:,:)
pr36 = pr36(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_MIROC-ESM_historical_safrica_climDJF1.nc","r")
pr37 = in->pr(:,0,:,:)
pr37 = pr37(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_MIROC-ESM-CHEM_historical_safrica_climDJF1.nc","r")
pr38 = in->pr(:,0,:,:)
pr38 = pr38(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_MPI-ESM-LR_historical_safrica_climDJF1.nc","r")
pr39 = in->pr(:,0,:,:)
pr39 = pr39(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_MPI-ESM-MR_historical_safrica_climDJF1.nc","r")
pr40 = in->pr(:,0,:,:)
pr40 = pr40(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_MPI-ESM-P_historical_safrica_climDJF1.nc","r")
pr41 = in->pr(:,0,:,:)
pr41 = pr41(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_MRI-CGCM3_historical_safrica_climDJF1.nc","r")
pr42 = in->pr(:,0,:,:)
pr42 = pr42(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_NorESM1-M_historical_safrica_climDJF1.nc","r")
pr43 = in->pr(:,0,:,:)
pr43 = pr43(time|0,z|0, {lat|-40:10}, {lon|-10:70})
in = addfile("/mnt/geog/ml382/melphd/regressionline/siczoutputprregline2/pr_Amon_NorESM1-ME_historical_safrica_climDJF1.nc","r")
pr44 = in->pr(:,0,:,:)
pr44 = pr44(time|0,z|0, {lat|-40:10}, {lon|-10:70})
lat =in->lat
time =in->time
;*************************************************************
;Calculations of max precip for lat and lon values
;**************************************************************
dimpr5 = dimsizes(pr5)
nlat = dimpr2(0)
mlon = dimpr2(1)
pr5MaxLon = new ( mlon, typeof(pr5), pr5@_FillValue)
do ml=0,mlon-1
imax = maxind(pr5(:,ml))
pr5MaxLon(ml) = dble2flt(lat(imax))
end do
print(pr5MaxLon)
print(pr5&lon)
print("-------------------------------")
print("pr5MaxLon: "+pr5&lon+" "+pr5MaxLon)
;Regression Line
regcMaxLon = regline(pr5&lon,pr5MaxLon)
print(regcMaxLon)
print(regcMaxLon@yave)
bMaxLon = regcMaxLon@yintercept
print(bMaxLon)
xMaxLon = pr5&lon
print(xMaxLon)
yMaxLon = regcMaxLon*pr5&lon + bMaxLon
print(yMaxLon)
print("-------------------------------")
print(xMaxLon+" "+yMaxLon)
;************************************************
; create an array to hold both the original data
; and the calculated regression line
;************************************************
datac = new ( (/2,dimsizes(pr5MaxLon)/), typeof(pr5MaxLon))
;datac(0,:) = pr5MaxLon
; y = mx+b
; m is the slope: rc returned from regline
; b is the y intercept: rc@yave attribute of rc returned from regline
datac(1,:) = dble2flt(regcMaxLon)*(dble2flt(xMaxLon)-dble2flt(regcMaxLon@xave)) + dble2flt(regcMaxLon@yave)
;************************************************
; plotting parameters
;************************************************
wks = gsn_open_wks("X11","Model_regline_panel") ; specifies a plot
gsn_define_colormap(wks,"gsltod") ; choose color map
res = True ; plot mods desired
res@gsnAddCyclic = False
res@cnFillOn = True ; color on
res@lbLabelStride = 2 ; every other label
res@gsnSpreadColors = True ; use full range of color map
res@cnLineLabelsOn = False ; no contour line labels
res@lbLabelBarOn = False
res@gsnDraw = False ; do not draw the plot
res@gsnFrame = False
pr4 = pr1(0,{-40:10},{-10:70})
dimpr4 = dimsizes(pr4)
latn = dimpr4(0)
lonm = dimpr4(1)
res@mpLimitMode = "Corners" ;
res@mpLeftCornerLonF = lon1(0)
res@mpRightCornerLonF = lon1(lonm-1)
res@mpLeftCornerLatF = lat1(0)
res@mpRightCornerLatF = lat1(latn-1)
res@tmXBMode = "Explicit" ; Define own tick mark labels.
res@tmXBValues = (/ -9.,0.,10.,20.,30.,40.,50.,60.,69 /)
res@tmXBLabels = (/ "10W","0","10E","20E","30E","40E","50E","60E","70E" /)
res@tmYLMode = "Explicit" ; Define own tick mark labels.
res@tmYLValues = (/ -39.,-30.,-20.,-10.,-0.,9 /)
res@tmYLLabels = (/"40S","30S","20S","10S","0","10N" /)
res@vpXF = 0.12 ; default is 0.2 change aspect
res@vpYF = 0.8 ; default is 0.8 ration
res@vpHeightF = 0.4 ; default is 0.6
res@vpWidthF = 0.8 ; default is 0.6
res@cnLevelSelectionMode = "ManualLevels" ; manual levels
res@cnMinLevelValF = 0 ; min level
res@cnMaxLevelValF = 400 ; max level
res@cnLevelSpacingF = 50 ; contour spacing
res@tmXBLabelFontHeightF = 0.02 ; resize tick labels
res@tmYLLabelFontHeightF = 0.02
res@pmLabelBarOrthogonalPosF = .10 ; move whole thing down
;************************************************
; plotting parameters
;************************************************
sres = True ; plot mods desired
sres@gsnMaximize = True ; maximize plot in frame
sres@xyMarkLineModes = (/"Markers","Lines"/) ; choose which have markers
sres@xyMarkers = 16 ; choose type of marker
sres@xyLineColors = (/"blue","black"/)
sres@xyMonoDashPattern = True
sres@gsLineDashPattern = 1
sres@xyDashPattern = 1 ; solid line
sres@xyLineThicknesses = (/1,4/) ; set second line to 2
sres@gsnDraw = False ; do not draw the plot
sres@gsnFrame = False
pr4@long_name = "FGOALS-g2"
pr4@units = "s=-0.17 lat=-11.47"
; reverse the first two colors
setvalues wks
"wkColorMap" : "gsltod"
"wkForegroundColor" : (/0.,0.,0./)
"wkBackgroundColor" : (/1.,1.,1./)
end setvalues
;************************************************
; plotting parameters
;************************************************
cres = True ; plot mods desired
cres@gsnMaximize = True ; maximize plot in frame
cres@xyMarkLineModes = "Lines" ; choose which have markers
;cres@xyMarkers = 16 ; choose type of marker
;cres@xyMarkerSizeF = 0.01 ; Marker size (default 0.01)
cres@xyDashPatterns = 0 ; solid line
cres@xyLineThicknesses = (/1,4/) ; set second line to 2
cres@gsnDraw = False ; do not draw the plot
cres@gsnFrame = False
; reverse the first two colors
setvalues wks
"wkColorMap" : "gsltod"
"wkForegroundColor" : (/0.,0.,0./)
"wkBackgroundColor" : (/1.,1.,1./)
end setvalues
print(rcMaxLon)
Panel the plots, first using rows x columns, then using
; number of plots per row
;
panel_res = True
panel_res@txString = "rows x columns"
panel_res@gsnPanelFigureStrings = (/"A","B","C"/)
gsn_panel(wks,(/xy1,xy2,xy3/),(/3,1/),panel_res)
panel_res@txString = "plots per row"
panel_res@gsnPanelRowSpec = True
gsn_panel(wks,(/xy1,xy2,xy3/),(/2,1/),panel_res)
plot = gsn_csm_contour_map(wks,pr4({-40:10},{-10:70}),res) ; create plot
plot2(i) = gsn_csm_xy(wks,pr(i)&lon,data,sres) ; create plot
plot3(j) = gsn_csm_xy(wks,pr(j)&lon,datac,cres) ; create plot
overlay(plot,plot2)
overlay(plot,plot3)
draw(plot)
frame(wks)
end
Many thanks!
Kind Regards
Melissa
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Mon Mar 31 11:15:20 2014
This archive was generated by hypermail 2.1.8 : Thu Apr 03 2014 - 13:36:27 MDT