begin dirdat = "./analdat/" hinddir = "../../climGCMs/p5degRes/" plotdir = "./epsfiles/" ;type = "x11" type = "png" MONINIT = 12 plotindividual = False ; NN, BN, AN separately RNORM = True ; use regression based forecast RPROB = False ; plot only statistically significant. NOT DONE YET plotmean = True ; plot Ensemble mean plotanom = True ; plot Ensemble anomaly plotcombined = True ; plot probability of BN,NN, and AN in one map plotOneModel = True ; plot individual model results chisquare = False ; chisquare is not convincing; if has to be used ; for it to have more meaning with plotOneModel = True ensmem = (/10,10,10,12,11,28,10/) ; read when running z1prec*seasonal*.ncl ; used for chisquare calculations ; Order of GCM is important. It should follow the order in z1* gcmtype = (/"CMC1-CanCM3","CMC2-CanCM4","COLA-CCSM4", "GFDL-CM2p5","NASA-GMAO","NCEP-CFSv2","NCAR-CESM1" /) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;; Do Monthly and Seasonal Ensemble ;;;;;; The folowing should change if the year changes ;;;;;; Also work only for 1 and 2 month leads ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; yrmonind1 =(/2016, 2016, 2016, 2016, 2016, 2016,\ 2016, 2016, 2016, 2016, 2016, 2017 /) yrmonind2 =(/2016, 2016, 2016, 2016, 2016, 2016,\ 2016, 2016, 2016, 2016, 2017, 2017 /) yrsnind1 =(/2016, 2016, 2016, 2016, 2016, 2016,\ 2016, 2016, 2016, 2017, 2017, 2017 /) yrsnind2 =(/2016, 2016, 2016, 2016, 2016, 2016,\ 2016, 2016, 2017, 2017, 2017, 2017 /) seasonname1 = (/"FMA","MAM","AMJ","MJJ","JJA","JAS", \ "ASO","SON","OND","NDJ","DJF","JFM" /) seasonname2 = (/"MAM","AMJ","MJJ","JJA","JAS","ASO", \ "SON","OND","NDJ","DJF","JFM","FMA" /) monthname1 = (/"Feb","Mar","Apr","May","Jun","Jul", \ "Aug","Sep","Oct","Nov","Dec","Jan" /) monthname2 = (/"Mar","Apr","May","Jun","Jul","Aug", \ "Sep","Oct","Nov","Dec","Jan","Feb" /) ndaymon1 = (/28.,31.,30.,31.,30.,31.,31.,30.,31.,30.,31.,31/) ; Lead+1 ndaymon2 = (/31.,30.,31.,30.,31.,31.,30.,31.,30.,31.,31.,28/) ; Lead+2 valmonind1 = (/1,2,3,4,5,6,7,8,9,10,11,0/) ;to index observation given valmonind2 = (/2,3,4,5,6,7,8,9,10,11,0,1/) ;MONINIT-1 to avoid if then ndayszn1 = (/89.,92.,91.,92.,92.,92.,92.,91.,92.,92.,90.,90./) ndayszn2 = (/92.,91.,92.,92.,92.,92.,91.,92.,92.,90.,90.,89./) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;; Unlikely for a need to change below ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; chiF = (/0.10,0.21,0.45,0.71,1.39,2.41,3.22,4.60,5.99,9.21,13.82/) chiP = (/0.95,0.90,0.80,0.70,0.50,0.30,0.20,0.10,0.05,0.01,0.001/) siglevel = 0.05 yrStart = 1982 ; to get common start year yrEobs = 2015 ; to get common end year yrEhind = 2010 ; to get common end year monJan = 1 ; All tseries start from Jan monDec = 12 ; All tseries end in Dec VARN = "Prec" Panom = "Anomaly" Pmean = "Mean" Tbelow = "BelowTerc" Tabove = "AboveTerc" Tnormal = "NormalTerc" Pbelow = "BelowProb" Pabove = "AboveProb" Pnormal = "NormalProb" Psig = "PValue" LeadSmon1 = 1 ; JAN index for DEC intial time LeadEmon1 = 1 ; JAN index for DEC intial time LNmon1 = LeadSmon1+"monLead" LeadSszn1 = 1 ; JAN index for DEC intial time LeadEszn1 = 3 ; MAR index for DEC intial time LNszn1 = LeadSszn1+"monLead" LeadSmon2 = 2 ; FEB of following yr index for DEC intial time LeadEmon2 = 2 ; FEB of following yr index for DEC intial time LNmon2 = LeadSmon2+"monLead" LeadSszn2 = 2 ; FEB of following yr index for DEC intial time LeadEszn2 = 4 ; APR of following yr index for DEC intial time LNszn2 = LeadSszn2+"monLead" ;;;;;;;;; pmin = 40. ; lowest probability level pmax = 80. plev = 10. raindrymon = 10.0 ;;;;;;;;; latS = -13.5 latN = 24.5 lonW = 20.5 lonE = 53.0 mlat = toint((latN-latS)/0.5 + 1) mlon = toint((lonE-lonW)/0.5 + 1) GHAlat = fspan(latS, latN, mlat) GHAlon = fspan(lonW, lonE, mlon) GHAlat@units = "degrees_north" GHAlon@units = "degrees_east" lat = GHAlat lon = GHAlon ;;; dirSHP = "/data/WRF/data/GHAshapefiles/" f1 = "gha.shp" a = addfile(hinddir+"precChirpsMonthly1982to2015r5km.nc","r") x = a->prec ; jan1982-Dec2015 nyrs = yrEobs-yrStart+1 y = reshape(x, (/nyrs,12,mlat,mlon/)) ;; Use monSind to index variables already shifted by Lead. No need 2 use if/then monSind = MONINIT - 1 Month1 = monthname1(monSind) Month2 = monthname2(monSind) ind1 = valmonind1(monSind) ind2 = valmonind2(monSind) obsA1 = y(:,ind1,:,:) obsA2 = y(:,ind2,:,:) ValidYr1 = yrmonind1(monSind) ValidYr2 = yrmonind2(monSind) print(VARN+" forecast for "+Month1+" and "+Month2) omean1 = dim_avg_n(obsA1,0) omean1!0 = "lat" omean1!1 = "lon" omean1&lat = GHAlat omean1&lon = GHAlon omean1@_FillValue = -9999. omean2 = dim_avg_n(obsA2,0) omean2!0 = "lat" omean2!1 = "lon" omean2&lat = GHAlat omean2&lon = GHAlon omean2@_FillValue = -9999. ; printVarSummary(obsA1) delete([/x,y, a, ind1,ind2,obsA1, obsA2/]) ;;;;;;;;;;;;;;;;;;;;;;;; fensanom1 = dirdat+Panom+VARN+LNmon1+Month1+"Monthly.nc" fensanom2 = dirdat+Panom+VARN+LNmon2+Month2+"Monthly.nc" fensmean1 = dirdat+Pmean+VARN+LNmon1+Month1+"Monthly.nc" fensmean2 = dirdat+Pmean+VARN+LNmon2+Month2+"Monthly.nc" if (RNORM) then fensbelow1 = dirdat+Pbelow+VARN+LNmon1+Month1+"Monthly.nc" fensbelow2 = dirdat+Pbelow+VARN+LNmon2+Month2+"Monthly.nc" fensabove1 = dirdat+Pabove+VARN+LNmon1+Month1+"Monthly.nc" fensabove2 = dirdat+Pabove+VARN+LNmon2+Month2+"Monthly.nc" fensnormal1 = dirdat+Pnormal+VARN+LNmon1+Month1+"Monthly.nc" fensnormal2 = dirdat+Pnormal+VARN+LNmon2+Month2+"Monthly.nc" fenspsig1 = dirdat+Psig+VARN+LNmon1+Month1+"Monthly.nc" fenspsig2 = dirdat+Psig+VARN+LNmon2+Month2+"Monthly.nc" ;figProb = "RanNorm" figProb = "Tercile" else fensbelow1 = dirdat+Tbelow+VARN+LNmon1+Month1+"Monthly.nc" fensbelow2 = dirdat+Tbelow+VARN+LNmon2+Month2+"Monthly.nc" fensabove1 = dirdat+Tabove+VARN+LNmon1+Month1+"Monthly.nc" fensabove2 = dirdat+Tabove+VARN+LNmon2+Month2+"Monthly.nc" fensnormal1 = dirdat+Tnormal+VARN+LNmon1+Month1+"Monthly.nc" fensnormal2 = dirdat+Tnormal+VARN+LNmon2+Month2+"Monthly.nc" figProb = "Tercile" end if ; RNORM ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Plot results ; Set some basic resources for all res = True res@gsnFrame = False res@gsnDraw = False res@gsnLeftString = "" res@gsnRightString = "" res@tfDoNDCOverlay = False res@gsnMaximize = False res@gsnAddCyclic = False res@mpDataBaseVersion = "MediumRes" res@mpDataSetName = "Earth..4" res@mpGeophysicalLineThicknessF = 4 res@mpNationalLineThicknessF = 4 res@pmTickMarkDisplayMode = "Always" res@mpLimitMode = "Corners" res@mpLeftCornerLatF = latS res@mpLeftCornerLonF = lonW res@mpRightCornerLatF = latN res@mpRightCornerLonF = lonE res@cnFillDrawOrder = "Predraw" res@cnLineDrawOrder = "Predraw" ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;Masking fill_specs = (/"ocean","land"/) mask_specs = (/\ "burundi","djibouti","eritrea",\ "ethiopia","kenya","rwanda","somalia",\ "sudan","tanzania","uganda","African Lakes" /) res@mpFillOn = True res@mpOutlineOn = False res@mpFillBoundarySets = "NoBoundaries" res@mpOutlineBoundarySets = "National" res@mpFillAreaSpecifiers = fill_specs res@mpSpecifiedFillColors = (/"white","white"/) res@mpAreaMaskingOn = 1 res@mpMaskAreaSpecifiers = mask_specs ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;; Define Colors tempKlevels = (/278, 280, 282, 284, 286, \ 288, 290, 292, 294, 296, \ 298, 300, 302, 304, 306, \ 308, 310, 312, 314 /) tempClevels = (/ 4, 6, 8, 10, 12, \ 14, 16, 18, 20, 22, \ 24, 26, 28, 30, 32, \ 34, 36, 38, 40 /) tempcolors = read_colormap_file("temp_19lev") colorsunshine = read_colormap_file("sunshine_9lev") ;;;; Create anomay color from blue to white to yellow to red cyananomaly = (/ (/0.95,1,1/), (/0.9, 1, 1/), (/0.8, 1, 1/), \ (/0.7, 1,1/), (/0.6, 1, 1/), (/0.5, 1, 1/), \ (/0.4, 1,1/), (/0.3, 1, 1/), (/0.2, 1, 1/), \ (/0.1, 1,1/), (/0,1,1/) /) x = (/"blue","white" /) opt = True opt@NumColorsInRange = (/ 9/) yb = span_named_colors(x,opt) yblue = yb(2::,:) x = (/"lightyellow","yellow" /) opt = True opt@NumColorsInRange = (/ 8/) yy1 = span_named_colors(x,opt) yy2 = yy1(2::,:) yyellow = yy2(1:2,:) x = (/"yellow","red" /) opt = True opt@NumColorsInRange = (/ 9/) yr1 = span_named_colors(x,opt) yyellow(1,:) = yr1(3,:) yred = yr1(4::,:) dim = dimsizes(yblue) y1d = dim(0) dim = dimsizes(yyellow) y2d = dim(0) dim = dimsizes(yred) y3d = dim(0) y = y1d+y2d+y3d coloranomaly = new((/y,dim(1)/),typeof(yblue)) coloranomaly(0:y1d-1,:) = yblue coloranomaly(y1d:(y1d+y2d-1),:) = yyellow coloranomaly(y1d+y2d:y-1,:) = yred blueanomaly = coloranomaly(0:y1d-2,:) blueanomaly = blueanomaly(::-1,:) redanomaly = coloranomaly(y1d:y-1,:) delete([/yb,yy1,yy2,yr1,y1d,y2d,y3d,dim,yblue,yred,x,opt,y,y1d,y2d/]) ;;;;;;;;;; ;; Temp color colorabove = redanomaly colorbelow = blueanomaly colornormal = cyananomaly ;;;; opt = True opt@fout = "colAnomalyTT2m.txt" fmtx = "3(f7.5,1x) " write_matrix(coloranomaly,fmtx,opt) opt@fout = "colAboveTT2m.txt" fmtx = "3(f7.5,1x) " xabove = colorabove write_matrix(xabove,fmtx,opt) opt@fout = "colBelowTT2m.txt" xbelow = colorbelow ;xbelow = toint(colorbelow*255) write_matrix(xbelow,fmtx,opt) ; print(colorbelow) opt@fout = "colNormalTT2m.txt" ;fmtx = "3(i4,1x) " xnormal = colornormal ;xnormal = toint(colornormal*255) write_matrix(xnormal,fmtx,opt) print(zew) ;;;;;;;;;;;;;;; ;;; Define plotting resources: ensmean ensmean = res ensmean@cnInfoLabelOn = False ensmean@cnConstFLabelOn = False ensmean@cnFillOn = True ensmean@cnLinesOn = False ensmean@lbOrientation = "Vertical" ensmean@tiMainOffsetYF = 0.02 ensmean@cnLevelSelectionMode = "ExplicitLevels" ensmean@cnLevels = wetlevels ensmean@cnFillColors = wetcolors ;;;;ensanomly ensanom = res ensanom@cnInfoLabelOn = False ensanom@cnConstFLabelOn = False ensanom@cnFillOn = True ensanom@cnLinesOn = False ensanom@lbOrientation = "Vertical" ensanom@tiMainOffsetYF = 0.02 ensanom@cnLevelSelectionMode = "ManualLevels" ensanom@cnFillPalette = coloranomaly ;;;;;;; ;;;; individual and combined probabilities ensp = res ensp@cnInfoLabelOn = False ensp@cnConstFLabelOn = False ensp@cnFillOn = True ensp@cnLinesOn = False ensp@lbOrientation = "Vertical" ensp@tiMainOffsetYF = 0.02 ensp@cnLevelSelectionMode = "ManualLevels" ensp@cnMinLevelValF = pmin ensp@cnMaxLevelValF = pmax ensp@cnLevelSpacingF = plev ; set contour spacing ensprob = ensp ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;; plot parameters for probability ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 1. Plot Multimodel Ensemble Average rainfall Monthly if (plotmean) then a = addfile(fensmean1,"r") fprec = a->ensmeanunbiased tmp = dim_avg_n_Wrap(fprec,0) wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon1+"_"+ \ Month1 + ValidYr1 +"_Total") cntr = gsn_csm_contour_map(wks,tmp,ensmean) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec/]) a = addfile(fensmean2,"r") fprec = a->ensmeanunbiased tmp = dim_avg_n_Wrap(fprec,0) wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon2+"_"+ \ Month2 + ValidYr2 +"_Total") cntr = gsn_csm_contour_map(wks,tmp,ensmean) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec/]) end if ;; 2. Plot Ensemble mean anomaly Monthly if (plotanom) then a = addfile(fensanom1,"r") fprec = a->ensmeananom tmp = dim_avg_n_Wrap(fprec,0) wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon1+"_"+ \ Month1 + ValidYr1 +"_Anomaly") print("Monthly anom range for "+ Month1 +" for Ensemble is") printMinMax(tmp,0) cmin = min(tmp) cmax = max(tmp) maxlev = 25 mnmx = nice_mnmxintvl(cmin,cmax,maxlev,False) ensanom@cnMinLevelValF = mnmx(0) ensanom@cnMaxLevelValF = mnmx(1) ensanom@cnLevelSpacingF = mnmx(2) delete(ensanom@cnFillPalette) ensanom@cnFillPalette = coloranomaly symMinMaxPlt (tmp,maxlev,False,ensanom) cntr = gsn_csm_contour_map(wks,tmp,ensanom) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec,mnmx,cmin,cmax,maxlev/]) a = addfile(fensanom2,"r") fprec = a->ensmeananom tmp = dim_avg_n_Wrap(fprec,0) wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon2+"_"+ \ Month2 + ValidYr2 +"_Anomaly") print("Monthly anom range for "+ Month2 +" for Ensemble is") printMinMax(tmp, 0) cmin = min(tmp) cmax = max(tmp) maxlev = 25 mnmx = nice_mnmxintvl(cmin,cmax,maxlev,False) ensanom@cnMinLevelValF = mnmx(0) ensanom@cnMaxLevelValF = mnmx(1) ensanom@cnLevelSpacingF = mnmx(2) delete(ensanom@cnFillPalette) ensanom@cnFillPalette = coloranomaly symMinMaxPlt (tmp,maxlev,False,ensanom) cntr = gsn_csm_contour_map(wks,tmp,ensanom) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec,mnmx,cmin,cmax,maxlev/]) end if ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 3. Plot Compined Probabilities Monthly ;;; if (plotcombined) then abelow = addfile(fensbelow1,"r") aabove = addfile(fensabove1,"r") anormal= addfile(fensnormal1,"r") tbelow = abelow->tercileBelow tabove = aabove->tercileAbove tnormal= anormal->tercileNormal ;; Do we need to compute Chisquare significnce if (chisquare) then t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) tercileInd = t0 tercileInd = t0@_FillValue tercileProb = tercileInd do iy = 0, mlat-1 do ix = 0, mlon-1 tmp = dimsizes(gcmtype) * \ ((t0(iy,ix) - 0.333) + \ (t1(iy,ix) - 0.333) + \ (t2(iy,ix) - 0.333))/0.333 indx = closest_val(tmp,chiF) if (.not.ismissing(indx) .and. indx .lt. dimsizes(chiF)) then if (chiP(indx) .lt. siglevel ) then tmp1 = ((/t0(iy,ix),t1(iy,ix), t2(iy,ix) /)) tercileInd(iy,ix) = maxind(tmp1) tercileProb(iy,ix) = (/tmp1(maxind(tmp1))/) delete(tmp1) else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if delete([/tmp,indx/]) end do end do end if ; chisquare ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Start combine plot from below normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 0., 1.0, 0.) x = t0 * tmpx delete([tmpx]) else x = where(t0 .gt. t1 .and. t0 .gt. t2 .and. t0 .gt. pmin, \ t0, tbelow@_FillValue) end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensp@lbTitleString = "BN" ensp@lbTitleFontHeightF = 0.02 ensp@lbTitleJust = "CenterCenter" ensp@lbTitlePosition = "Top" ensp@lbTitleOffsetF = -0.05 ensp@lbLabelsOn = False ensp@pmLabelBarOrthogonalPosF = -0.06 ensp@lbBottomMarginF = 0.3 ensp@cnFillPalette = colorbelow ; ensp@cnConstFEnableFill = False ensp@cnConstFLabelOn = False ;ensp@cnMissingValFillColor = "gray" ensp@cnNoDataLabelOn = False wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon1+"_"+ \ Month1 + ValidYr1 +"_"+figProb+"Prob") cntrb = gsn_csm_contour_map(wks,x,ensp) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) delete([/t0,t1,t2,x,ensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Near Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 1., 1.0, 0.) x = t1 * tmpx delete(tmpx) else x = where(t1 .gt. t0 .and. t1 .gt. t2 .and. t1 .gt. pmin,\ t1, tnormal@_FillValue) end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp = True bensp@gsnDraw = False bensp@gsnFrame = False bensp@cnInfoLabelOn = False bensp@cnFillOn = True bensp@cnLinesOn = False bensp@lbOrientation = "Vertical" bensp@cnFillPalette = colornormal bensp@lbTitleString = "NN" bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterCenter" bensp@lbTitlePosition = "Top" bensp@lbTitleOffsetF = -0.05 bensp@lbLabelsOn = False bensp@cnLevelSelectionMode = "ManualLevels" bensp@cnMinLevelValF = pmin bensp@cnMaxLevelValF = pmax bensp@cnLevelSpacingF = plev bensp@pmLabelBarOrthogonalPosF = -0.115 bensp@lbBottomMarginF = 0.3 bensp@cnFillDrawOrder = "Predraw" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntrn = gsn_csm_contour(wks,x,bensp) delete([/t0,t1,t2,x,bensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Above Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 2., 1.0, 0.) x = t2 * tmpx delete(tmpx) else x = where(t2 .gt. t0 .and. t2 .gt. t1 .and. t2 .gt. pmin, \ t2, tabove@_FillValue) end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp@lbOrientation = "vertical" bensp@lbTitleString = " AN" bensp@lbTitleOffsetF = -0.05 bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterLeft" bensp@lbTitlePosition = "Top" bensp@lbLabelsOn = True bensp@pmLabelBarOrthogonalPosF = 0.0 bensp@cnFillPalette = colorabove bensp@lbBottomMarginF = 0.3 bensp@lbLeftMarginF = 0.1 ;bensp@cnMissingValFillColor = "gray" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntra = gsn_csm_contour(wks,x,bensp) ;;;;;;;;;;;;;;;;;;;;;;;;;; ;plot gray for precip less than raindrymon x = where(omean1 .ge. raindrymon , omean1@_FillValue, 0.) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon xensp = True xensp@gsnDraw = False xensp@gsnFrame = False xensp@cnInfoLabelOn = False xensp@cnFillOn = True xensp@cnLinesOn = False xensp@lbOrientation = "Vertical" xensp@cnFillPalette = (/"gray","gray"/) xensp@lbLabelsOn = False xensp@lbLabelBarOn = False xensp@cnFillDrawOrder = "Predraw" xensp@cnConstFEnableFill = True xensp@cnConstFLabelOn = False cntrx = gsn_csm_contour(wks,x,xensp) ;;;;;;;;;;;;;;;;;;;;;;; overlay(cntrb,cntrn) overlay(cntrb,cntra) overlay(cntrb,cntrx) draw(cntrb) frame(wks) if (chisquare) then delete([/tercileProb,tercileInd/]) end if delete([/abelow,aabove,anormal, tbelow, tabove,tnormal/]) delete([/t0,t1,t2,x, cntrb,cntra,cntrn,cntrx/]) ;;;;;;; End combined plot for Month 1 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;; Do combined plot for Month2 abelow = addfile(fensbelow2,"r") aabove = addfile(fensabove2,"r") anormal= addfile(fensnormal2,"r") tbelow = abelow->tercileBelow tabove = aabove->tercileAbove tnormal= anormal->tercileNormal ;printMinMax(tnormal,0) ;print(num(tnormal .gt. pmin)) ;; Do we need to compute Chisquare significnce if (chisquare) then t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) tercileInd = t0 tercileInd = t0@_FillValue tercileProb = tercileInd do iy = 0, mlat-1 do ix = 0, mlon-1 tmp = dimsizes(gcmtype) * \ ((t0(iy,ix) - 0.333) + \ (t1(iy,ix) - 0.333) + \ (t2(iy,ix) - 0.333))/0.333 indx = closest_val(tmp,chiF) if (.not.ismissing(indx) .and. indx .lt. dimsizes(chiF)) then if (chiP(indx) .lt. siglevel ) then tmp1 = ((/t0(iy,ix),t1(iy,ix), t2(iy,ix) /)) ; print(chiP(indx)+" "+maxind(tmp1)+" "+tmp1+" "+tmp1(maxind(tmp1))) tercileInd(iy,ix) = maxind(tmp1) tercileProb(iy,ix) = (/tmp1(maxind(tmp1))/) delete(tmp1) else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if delete([/tmp,indx/]) end do end do end if ; chisquare ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Start combine plot from below normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 0., 1.0, 0.) x = t0 * tmpx delete([tmpx]) else x = where(t0 .gt. t1 .and. t0 .gt. t2 .and. t0 .gt. pmin, \ t0, tbelow@_FillValue) end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensp@lbTitleString = "BN" ensp@lbTitleFontHeightF = 0.02 ensp@lbTitleJust = "CenterCenter" ensp@lbTitlePosition = "Top" ensp@lbTitleOffsetF = -0.05 ensp@lbLabelsOn = False ensp@pmLabelBarOrthogonalPosF = -0.06 ensp@lbBottomMarginF = 0.3 ensp@cnFillPalette = colorbelow ; ensp@cnConstFEnableFill = False ensp@cnConstFLabelOn = False ;ensp@cnMissingValFillColor = "gray" ensp@cnNoDataLabelOn = False wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon2+"_"+ \ Month2 + ValidYr2 +"_"+figProb+"Prob") cntrb = gsn_csm_contour_map(wks,x,ensp) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) delete([/t0,t1,t2,x,ensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Near Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 1., 1.0, 0.) x = t1 * tmpx delete(tmpx) else x = where(t1 .gt. t0 .and. t1 .gt. t2 .and. t1 .gt. pmin,\ t1, tnormal@_FillValue) end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp = True bensp@gsnDraw = False bensp@gsnFrame = False bensp@cnInfoLabelOn = False bensp@cnFillOn = True bensp@cnLinesOn = False bensp@lbOrientation = "Vertical" bensp@cnFillPalette = colornormal bensp@lbTitleString = "NN" bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterCenter" bensp@lbTitlePosition = "Top" bensp@lbTitleOffsetF = -0.05 bensp@lbLabelsOn = False bensp@cnLevelSelectionMode = "ManualLevels" bensp@cnMinLevelValF = pmin bensp@cnMaxLevelValF = pmax bensp@cnLevelSpacingF = plev bensp@pmLabelBarOrthogonalPosF = -0.115 bensp@lbBottomMarginF = 0.3 bensp@cnFillDrawOrder = "Predraw" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntrn = gsn_csm_contour(wks,x,bensp) delete([/t0,t1,t2,x,bensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Above Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 2., 1.0, 0.) x = t2 * tmpx delete(tmpx) else x = where(t2 .gt. t0 .and. t2 .gt. t1 .and. t2 .gt. pmin, \ t2, tabove@_FillValue) end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp@lbOrientation = "vertical" bensp@lbTitleString = " AN" bensp@lbTitleOffsetF = -0.05 bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterLeft" bensp@lbTitlePosition = "Top" bensp@lbLabelsOn = True bensp@pmLabelBarOrthogonalPosF = 0.0 bensp@cnFillPalette = colorabove bensp@lbBottomMarginF = 0.3 bensp@lbLeftMarginF = 0.1 ;bensp@cnMissingValFillColor = "gray" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntra = gsn_csm_contour(wks,x,bensp) ;;;;;;;;;;;;;;;;;;;;;;;;;; ;plot gray for precip less than raindrymon x = where(omean2 .ge. raindrymon , omean2@_FillValue, 0.) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon xensp = True xensp@gsnDraw = False xensp@gsnFrame = False xensp@cnInfoLabelOn = False xensp@cnFillOn = True xensp@cnLinesOn = False xensp@lbOrientation = "Vertical" xensp@cnFillPalette = (/"gray","gray"/) xensp@lbLabelsOn = False xensp@lbLabelBarOn = False xensp@cnFillDrawOrder = "Predraw" xensp@cnConstFEnableFill = True xensp@cnConstFLabelOn = False cntrx = gsn_csm_contour(wks,x,xensp) ;;;;;;;;;;;;;;;;;;;;;;; overlay(cntrb,cntrn) overlay(cntrb,cntra) overlay(cntrb,cntrx) draw(cntrb) frame(wks) if (chisquare) then delete([/tercileProb,tercileInd/]) end if delete([/abelow,aabove,anormal, tbelow, tabove,tnormal/]) delete([/t0,t1,t2,x, cntrb,cntra,cntrn,cntrx/]) ;;;;;; End combined plot for Month 2 end if ; plotcombined ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Individual probability plots ;;; Above normal rainfall probability ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 4. Plot Individual probabilities Monthly if (plotindividual) then ;;;; Do Month1 Above a = addfile(fensabove1,"r") tercile = a->tercileAbove x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensprob@cnFillPalette = colorabove wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon1+"_"+ \ Month1 + ValidYr1 +"_"+figProb+"AboveNormal") cntra = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntra,dirSHP+f1,lnres) draw(cntra) frame(wks) delete([/a,x,tercile/]) ;;; Do Month2 Above a = addfile(fensabove2,"r") tercile = a->tercileAbove x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensprob@cnFillPalette = colorabove wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon2+"_"+ \ Month2 + ValidYr2 +"_"+figProb+"AboveNormal") cntra = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntra,dirSHP+f1,lnres) draw(cntra) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Normal probability ;;; Do for Month1 Normal a = addfile(fensnormal1,"r") tercile = a->tercileNormal x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colornormal wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon1+"_"+ \ Month1 + ValidYr1 +"_"+figProb+"NearNormal") cntrn = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrn,dirSHP+f1,lnres) draw(cntrn) frame(wks) delete([/a,x,tercile/]) ;;;; Do for Month2 Normal a = addfile(fensnormal2,"r") tercile = a->tercileNormal x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colornormal wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon2+"_"+ \ Month2 + ValidYr2 +"_"+figProb+"NearNormal") cntrn = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrn,dirSHP+f1,lnres) draw(cntrn) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Below probability ;;;; Do for Month1 Below a = addfile(fensbelow1,"r") tercile = a->tercileBelow x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colorbelow wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon1+"_"+ \ Month1 + ValidYr1 +"_"+figProb+"BelowNormal") cntrb = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) draw(cntrb) frame(wks) delete([/a,x,tercile/]) ;;;; Do for Month2 Below a = addfile(fensbelow2,"r") tercile = a->tercileBelow x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colorbelow wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNmon2+"_"+ \ Month2 + ValidYr2 +"_"+figProb+"BelowNormal") cntrb = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) draw(cntrb) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;; end if ; plotindividual monthly ;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;; Loop through gcmtype for Monthly ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; if (plotOneModel) then do imodel = 0, dimsizes(gcmtype)-1 ;; 1. Plot Multimodel Ensemble Average rainfall Monthly if (plotmean) then a = addfile(fensmean1,"r") fprec = a->ensmeanunbiased tmp = fprec(imodel,:,:) wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon1+"_"+ Month1 + ValidYr1 +"_Total") cntr = gsn_csm_contour_map(wks,tmp,ensmean) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec/]) a = addfile(fensmean2,"r") fprec = a->ensmeanunbiased tmp = fprec(imodel,:,:) wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon2+"_"+ Month2 + ValidYr2 +"_Total") cntr = gsn_csm_contour_map(wks,tmp,ensmean) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec/]) end if ;; 2. Plot Ensemble mean anomaly Monthly if (plotanom) then a = addfile(fensanom1,"r") fprec = a->ensmeananom tmp = fprec(imodel,:,:) wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon1+"_"+ Month1 + ValidYr1 +"_Anomaly") print("Monthly anom range for "+ Month1 +" for "+gcmtype(imodel)) printMinMax(tmp,0) cmin = min(tmp) cmax = max(tmp) maxlev = 25 mnmx = nice_mnmxintvl(cmin,cmax,maxlev,False) ensanom@cnMinLevelValF = mnmx(0) ensanom@cnMaxLevelValF = mnmx(1) ensanom@cnLevelSpacingF = mnmx(2) delete(ensanom@cnFillPalette) ensanom@cnFillPalette = coloranomaly symMinMaxPlt (tmp,maxlev,False,ensanom) cntr = gsn_csm_contour_map(wks,tmp,ensanom) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec,mnmx,cmin,cmax,maxlev/]) a = addfile(fensanom2,"r") fprec = a->ensmeananom tmp = fprec(imodel,:,:) wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon2+"_"+ Month2 + ValidYr2 +"_Anomaly") print("Monthly anom range for "+ Month2 +" for "+gcmtype(imodel)) printMinMax(tmp,0) cmin = min(tmp) cmax = max(tmp) maxlev = 25 mnmx = nice_mnmxintvl(cmin,cmax,maxlev,False) ensanom@cnMinLevelValF = mnmx(0) ensanom@cnMaxLevelValF = mnmx(1) ensanom@cnLevelSpacingF = mnmx(2) delete(ensanom@cnFillPalette) ensanom@cnFillPalette = coloranomaly symMinMaxPlt (tmp,maxlev,False,ensanom) cntr = gsn_csm_contour_map(wks,tmp,ensanom) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec,mnmx,cmin,cmax,maxlev/]) end if ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 3. Plot Compined Probabilities Monthly ;;; if (plotcombined) then abelow = addfile(fensbelow1,"r") aabove = addfile(fensabove1,"r") anormal= addfile(fensnormal1,"r") tbelow = abelow->tercileBelow tabove = aabove->tercileAbove tnormal= anormal->tercileNormal if (RNORM .and. RPROB) then asiglev = addfile(fenspsig1,"r") psigval = asiglev->psig end if ;; Do we need to compute Chisquare significnce if (chisquare) then t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) tercileInd = t0 tercileInd = t0@_FillValue tercileProb = tercileInd do iy = 0, mlat-1 do ix = 0, mlon-1 tmp = ensmem(imodel) * \ ((t0(iy,ix) - 0.333) + \ (t1(iy,ix) - 0.333) + \ (t2(iy,ix) - 0.333))/0.333 indx = closest_val(tmp,chiF) if (.not.ismissing(indx) .and. indx .lt. dimsizes(chiF)) then if (chiP(indx) .lt. siglevel ) then tmp1 = ((/t0(iy,ix),t1(iy,ix), t2(iy,ix) /)) ; print(chiP(indx)+" "+maxind(tmp1)+" "+tmp1+" "+tmp1(maxind(tmp1))) tercileInd(iy,ix) = maxind(tmp1) tercileProb(iy,ix) = (/tmp1(maxind(tmp1))/) delete(tmp1) else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if delete([/tmp,indx/]) end do end do end if ; chisquare ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Start combine plot from below normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 0., 1.0, 0.) x = t0 * tmpx delete([tmpx]) else x = where(t0 .gt. t1 .and. t0 .gt. t2 .and. t0 .gt. pmin, \ t0, tbelow@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensp@lbTitleString = "BN" ensp@lbTitleFontHeightF = 0.02 ensp@lbTitleJust = "CenterCenter" ensp@lbTitlePosition = "Top" ensp@lbTitleOffsetF = -0.05 ensp@lbLabelsOn = False ensp@pmLabelBarOrthogonalPosF = -0.06 ensp@lbBottomMarginF = 0.3 ensp@cnFillPalette = colorbelow ; ensp@cnConstFEnableFill = False ensp@cnConstFLabelOn = False ;ensp@cnMissingValFillColor = "gray" ensp@cnNoDataLabelOn = False wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon1+"_"+ Month1 + ValidYr1 +"_"+figProb+"Prob") cntrb = gsn_csm_contour_map(wks,x,ensp) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) delete([/t0,t1,t2,x,ensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Near Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 1., 1.0, 0.) x = t1 * tmpx delete(tmpx) else x = where(t1 .gt. t0 .and. t1 .gt. t2 .and. t1 .gt. pmin,\ t1, tnormal@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp = True bensp@gsnDraw = False bensp@gsnFrame = False bensp@cnInfoLabelOn = False bensp@cnFillOn = True bensp@cnLinesOn = False bensp@lbOrientation = "Vertical" bensp@cnFillPalette = colornormal bensp@lbTitleString = "NN" bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterCenter" bensp@lbTitlePosition = "Top" bensp@lbTitleOffsetF = -0.05 bensp@lbLabelsOn = False bensp@cnLevelSelectionMode = "ManualLevels" bensp@cnMinLevelValF = pmin bensp@cnMaxLevelValF = pmax bensp@cnLevelSpacingF = plev bensp@pmLabelBarOrthogonalPosF = -0.115 bensp@lbBottomMarginF = 0.3 bensp@cnFillDrawOrder = "Predraw" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntrn = gsn_csm_contour(wks,x,bensp) delete([/t0,t1,t2,x,bensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Above Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 2., 1.0, 0.) x = t2 * tmpx delete(tmpx) else x = where(t2 .gt. t0 .and. t2 .gt. t1 .and. t2 .gt. pmin, \ t2, tabove@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp@lbOrientation = "vertical" bensp@lbTitleString = " AN" bensp@lbTitleOffsetF = -0.05 bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterLeft" bensp@lbTitlePosition = "Top" bensp@lbLabelsOn = True bensp@pmLabelBarOrthogonalPosF = 0.0 bensp@cnFillPalette = colorabove bensp@lbBottomMarginF = 0.3 bensp@lbLeftMarginF = 0.1 ;bensp@cnMissingValFillColor = "gray" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntra = gsn_csm_contour(wks,x,bensp) ;;;;;;;;;;;;;;;;;;;;;;;;;; ;plot gray for precip less than raindrymon x = where(omean1 .ge. raindrymon , omean1@_FillValue, 0.) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon xensp = True xensp@gsnDraw = False xensp@gsnFrame = False xensp@cnInfoLabelOn = False xensp@cnFillOn = True xensp@cnLinesOn = False xensp@lbOrientation = "Vertical" xensp@cnFillPalette = (/"gray","gray"/) xensp@lbLabelsOn = False xensp@lbLabelBarOn = False xensp@cnFillDrawOrder = "Predraw" xensp@cnConstFEnableFill = True xensp@cnConstFLabelOn = False cntrx = gsn_csm_contour(wks,x,xensp) ;;;;;;;;;;;;;;;;;;;;;;; overlay(cntrb,cntrn) overlay(cntrb,cntra) overlay(cntrb,cntrx) draw(cntrb) frame(wks) if (chisquare) then delete([/tercileProb,tercileInd/]) end if delete([/abelow,aabove,anormal, tbelow, tabove,tnormal/]) delete([/t0,t1,t2,x, cntrb,cntra,cntrn,cntrx/]) if (RNORM .and. RPROB) then delete([/asiglev,psigval /]) end if ;;;;;;; End combined plot for Month 1 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;; Do combined plot for Month2 abelow = addfile(fensbelow2,"r") aabove = addfile(fensabove2,"r") anormal= addfile(fensnormal2,"r") tbelow = abelow->tercileBelow tabove = aabove->tercileAbove tnormal= anormal->tercileNormal if (RNORM .and. RPROB) then asiglev = addfile(fenspsig2,"r") psigval = asiglev->psig end if ;; Do we need to compute Chisquare significnce if (chisquare) then t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) tercileInd = t0 tercileInd = t0@_FillValue tercileProb = tercileInd do iy = 0, mlat-1 do ix = 0, mlon-1 tmp = ensmem(imodel) * \ ((t0(iy,ix) - 0.333) + \ (t1(iy,ix) - 0.333) + \ (t2(iy,ix) - 0.333))/0.333 indx = closest_val(tmp,chiF) if (.not.ismissing(indx) .and. indx .lt. dimsizes(chiF)) then if (chiP(indx) .lt. siglevel ) then tmp1 = ((/t0(iy,ix),t1(iy,ix), t2(iy,ix) /)) ; print(chiP(indx)+" "+maxind(tmp1)+" "+tmp1+" "+tmp1(maxind(tmp1))) tercileInd(iy,ix) = maxind(tmp1) tercileProb(iy,ix) = (/tmp1(maxind(tmp1))/) delete(tmp1) else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if delete([/tmp,indx/]) end do end do end if ; chisquare ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Start combine plot from below normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 0., 1.0, 0.) x = t0 * tmpx delete([tmpx]) else x = where(t0 .gt. t1 .and. t0 .gt. t2 .and. t0 .gt. pmin, \ t0, tbelow@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensp@lbTitleString = "BN" ensp@lbTitleFontHeightF = 0.02 ensp@lbTitleJust = "CenterCenter" ensp@lbTitlePosition = "Top" ensp@lbTitleOffsetF = -0.05 ensp@lbLabelsOn = False ensp@pmLabelBarOrthogonalPosF = -0.06 ensp@lbBottomMarginF = 0.3 ensp@cnFillPalette = colorbelow ; ensp@cnConstFEnableFill = False ensp@cnConstFLabelOn = False ;ensp@cnMissingValFillColor = "gray" ensp@cnNoDataLabelOn = False wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon2+"_"+ Month2 + ValidYr2 +"_"+figProb+"Prob") cntrb = gsn_csm_contour_map(wks,x,ensp) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) delete([/t0,t1,t2,x,ensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Near Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 1., 1.0, 0.) x = t1 * tmpx delete(tmpx) else x = where(t1 .gt. t0 .and. t1 .gt. t2 .and. t1 .gt. pmin,\ t1, tnormal@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp = True bensp@gsnDraw = False bensp@gsnFrame = False bensp@cnInfoLabelOn = False bensp@cnFillOn = True bensp@cnLinesOn = False bensp@lbOrientation = "Vertical" bensp@cnFillPalette = colornormal bensp@lbTitleString = "NN" bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterCenter" bensp@lbTitlePosition = "Top" bensp@lbTitleOffsetF = -0.05 bensp@lbLabelsOn = False bensp@cnLevelSelectionMode = "ManualLevels" bensp@cnMinLevelValF = pmin bensp@cnMaxLevelValF = pmax bensp@cnLevelSpacingF = plev bensp@pmLabelBarOrthogonalPosF = -0.115 bensp@lbBottomMarginF = 0.3 bensp@cnFillDrawOrder = "Predraw" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntrn = gsn_csm_contour(wks,x,bensp) delete([/t0,t1,t2,x,bensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Above Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 2., 1.0, 0.) x = t2 * tmpx delete(tmpx) else x = where(t2 .gt. t0 .and. t2 .gt. t1 .and. t2 .gt. pmin, \ t2, tabove@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp@lbOrientation = "vertical" bensp@lbTitleString = " AN" bensp@lbTitleOffsetF = -0.05 bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterLeft" bensp@lbTitlePosition = "Top" bensp@lbLabelsOn = True bensp@pmLabelBarOrthogonalPosF = 0.0 bensp@cnFillPalette = colorabove bensp@lbBottomMarginF = 0.3 bensp@lbLeftMarginF = 0.1 ;bensp@cnMissingValFillColor = "gray" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntra = gsn_csm_contour(wks,x,bensp) ;;;;;;;;;;;;;;;;;;;;;;;;;; ;plot gray for precip less than raindrymon x = where(omean2 .ge. raindrymon , omean2@_FillValue, 0.) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon xensp = True xensp@gsnDraw = False xensp@gsnFrame = False xensp@cnInfoLabelOn = False xensp@cnFillOn = True xensp@cnLinesOn = False xensp@lbOrientation = "Vertical" xensp@cnFillPalette = (/"gray","gray"/) xensp@lbLabelsOn = False xensp@lbLabelBarOn = False xensp@cnFillDrawOrder = "Predraw" xensp@cnConstFEnableFill = True xensp@cnConstFLabelOn = False cntrx = gsn_csm_contour(wks,x,xensp) ;;;;;;;;;;;;;;;;;;;;;;; overlay(cntrb,cntrn) overlay(cntrb,cntra) overlay(cntrb,cntrx) draw(cntrb) frame(wks) if (chisquare) then delete([/tercileProb,tercileInd/]) end if delete([/abelow,aabove,anormal, tbelow, tabove,tnormal/]) delete([/t0,t1,t2,x, cntrb,cntra,cntrn,cntrx/]) if (RNORM .and. RPROB) then delete([/asiglev,psigval /]) end if ;;;;;; End combined plot for Month 2 end if ; plotcombined ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Individual probability plots ;;; Above normal rainfall probability ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 4. Plot Individual probabilities Monthly if (plotindividual) then ;;;; Do Month1 Above a = addfile(fensabove1,"r") tercile = a->tercileAbove x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensprob@cnFillPalette = colorabove wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon1+"_"+ Month1 + ValidYr1 +"_"+figProb+"AboveNormal") cntra = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntra,dirSHP+f1,lnres) draw(cntra) frame(wks) delete([/a,x,tercile/]) ;;; Do Month2 Above a = addfile(fensabove2,"r") tercile = a->tercileAbove x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensprob@cnFillPalette = colorabove wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon2+"_"+ Month2 + ValidYr2 +"_"+figProb+"AboveNormal") cntra = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntra,dirSHP+f1,lnres) draw(cntra) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Normal probability ;;; Do for Month1 Normal a = addfile(fensnormal1,"r") tercile = a->tercileNormal x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colornormal wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon1+"_"+ Month1 + ValidYr1 +"_"+figProb+"NearNormal") cntrn = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrn,dirSHP+f1,lnres) draw(cntrn) frame(wks) delete([/a,x,tercile/]) ;;;; Do for Month2 Normal a = addfile(fensnormal2,"r") tercile = a->tercileNormal x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colornormal wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon2+"_"+ Month2 + ValidYr2 +"_"+figProb+"NearNormal") cntrn = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrn,dirSHP+f1,lnres) draw(cntrn) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Below probability ;;;; Do for Month1 Below a = addfile(fensbelow1,"r") tercile = a->tercileBelow x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colorbelow wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon1+"_"+ Month1 + ValidYr1 +"_"+figProb+"BelowNormal") cntrb = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) draw(cntrb) frame(wks) delete([/a,x,tercile/]) ;;;; Do for Month2 Below a = addfile(fensbelow2,"r") tercile = a->tercileBelow x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colorbelow wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNmon2+"_"+ Month2 + ValidYr2 +"_"+figProb+"BelowNormal") cntrb = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) draw(cntrb) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;; end if ; plotindividual monthly ;;;;;;;;;; end do ; imodel Monthly end if ; plotOneModel delete([/omean1,omean2, ensmean,ensanom,ensp,ensprob/]) ;;;;;;;;;;;;;;;;;;;;; ;;;;; Done all Monthly; Now Start Seasonal ;;;;; Plot-Seasonal ;;;;;;;;;;;;;;; ;;;;;;;;; raindryszn = 30.0 ; gray areas with < raindry ;;;;;;;;; yrAllObs = ispan(yrStart,yrEobs,1) ind10 = get1Dindex(yrAllObs,yrEhind) all29inds = ispan(0,ind10,1) all29yrs = yrAllObs(all29inds) all28inds = ispan(1,ind10,1) all28yrs = yrAllObs(all28inds) ;; Use monSind to index variables already shifted by Lead. No need 2 use if/then monSind = MONINIT - 1 Season1 = seasonname1(monSind) Season2 = seasonname2(monSind) ValidYr1 = yrsnind1(monSind) ValidYr2 = yrsnind2(monSind) print(VARN+" forecast for "+Season1+" and "+Season2) a = addfile(hinddir+"precChirpsMonthly1982to2015r5km.nc","r") x = a->prec ; jan1982-Dec2015 y = month_to_seasonN(x,(/Season1,Season2/)) ;;;;;;;; ; First DJF of Jan-Dec data is only JF so ignore if (Season1 .eq. "DJF") then obsA1 = y(0,all28inds,:,:) * 3.0 else obsA1 = y(0,all29inds,:,:) * 3.0 end if if (Season2 .eq. "DJF") then obsA2 = y(1,all28inds,:,:) * 3.0 else obsA2 = y(1,all29inds,:,:) * 3.0 end if omean1 = dim_avg_n(obsA1,0) omean1!0 = "lat" omean1!1 = "lon" omean1&lat = GHAlat omean1&lon = GHAlon omean1@_FillValue = -9999. omean2 = dim_avg_n(obsA2,0) omean2!0 = "lat" omean2!1 = "lon" omean2&lat = GHAlat omean2&lon = GHAlon omean2@_FillValue = -9999. ; printVarSummary(obsA1) delete([/x,y, a, obsA1, obsA2/]) ;;;;;;;;;;;;;;;;;;;;;;;; fensanom1 = dirdat+Panom+VARN+LNszn1+Season1+"Seasonal.nc" fensanom2 = dirdat+Panom+VARN+LNszn2+Season2+"Seasonal.nc" fensmean1 = dirdat+Pmean+VARN+LNszn1+Season1+"Seasonal.nc" fensmean2 = dirdat+Pmean+VARN+LNszn2+Season2+"Seasonal.nc" if (RNORM) then fensbelow1 = dirdat+Pbelow+VARN+LNszn1+Season1+"Seasonal.nc" fensbelow2 = dirdat+Pbelow+VARN+LNszn2+Season2+"Seasonal.nc" fensabove1 = dirdat+Pabove+VARN+LNszn1+Season1+"Seasonal.nc" fensabove2 = dirdat+Pabove+VARN+LNszn2+Season2+"Seasonal.nc" fensnormal1 = dirdat+Pnormal+VARN+LNszn1+Season1+"Seasonal.nc" fensnormal2 = dirdat+Pnormal+VARN+LNszn2+Season2+"Seasonal.nc" fenspsig1 = dirdat+Psig+VARN+LNszn1+Season1+"Seasonal.nc" fenspsig2 = dirdat+Psig+VARN+LNszn2+Season2+"Seasonal.nc" ;figProb = "RanNorm" figProb = "Tercile" else fensbelow1 = dirdat+Tbelow+VARN+LNszn1+Season1+"Seasonal.nc" fensbelow2 = dirdat+Tbelow+VARN+LNszn2+Season2+"Seasonal.nc" fensabove1 = dirdat+Tabove+VARN+LNszn1+Season1+"Seasonal.nc" fensabove2 = dirdat+Tabove+VARN+LNszn2+Season2+"Seasonal.nc" fensnormal1 = dirdat+Tnormal+VARN+LNszn1+Season1+"Seasonal.nc" fensnormal2 = dirdat+Tnormal+VARN+LNszn2+Season2+"Seasonal.nc" figProb = "Tercile" end if ; RNORM ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Define plotting resources: ensmean ensmean = res ensmean@cnInfoLabelOn = False ensmean@cnConstFLabelOn = False ensmean@cnFillOn = True ensmean@cnLinesOn = False ensmean@lbOrientation = "Vertical" ensmean@tiMainOffsetYF = 0.02 ensmean@cnLevelSelectionMode = "ExplicitLevels" ensmean@cnLevels = wetlevels ensmean@cnFillColors = wetcolors ;;;;ensanomly ensanom = res ensanom@cnInfoLabelOn = False ensanom@cnConstFLabelOn = False ensanom@cnFillOn = True ensanom@cnLinesOn = False ensanom@lbOrientation = "Vertical" ensanom@tiMainOffsetYF = 0.02 ensanom@cnLevelSelectionMode = "ManualLevels" ensanom@cnFillPalette = coloranomaly ;;;;;;; ;;;; individual and combled probabilities ensp = res ensp@cnInfoLabelOn = False ensp@cnConstFLabelOn = False ensp@cnFillOn = True ensp@cnLinesOn = False ensp@lbOrientation = "Vertical" ensp@tiMainOffsetYF = 0.02 ensp@cnLevelSelectionMode = "ManualLevels" ensp@cnMinLevelValF = pmin ensp@cnMaxLevelValF = pmax ensp@cnLevelSpacingF = plev ; set contour spacing ensprob = ensp ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;; plot parameters for probability ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 5. Plot Multimodel Ensemble Average rainfall Seasonal if (plotmean) then a = addfile(fensmean1,"r") fprec = a->ensmeanunbiased tmp = dim_avg_n_Wrap(fprec,0) wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn1+"_"+ \ Season1 + ValidYr1 +"_Total") cntr = gsn_csm_contour_map(wks,tmp,ensmean) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec/]) a = addfile(fensmean2,"r") fprec = a->ensmeanunbiased tmp = dim_avg_n_Wrap(fprec,0) wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn2+"_"+ \ Season2 + ValidYr2 +"_Total") cntr = gsn_csm_contour_map(wks,tmp,ensmean) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec/]) end if ;;;;;;;;;;;;;;;;; ;; 6. Plot Ensemble mean anomaly Seasonal if (plotanom) then a = addfile(fensanom1,"r") fprec = a->ensmeananom tmp = dim_avg_n_Wrap(fprec,0) wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn1+"_"+ \ Season1 + ValidYr1 +"_Anomaly") print("Seasonal anom range for "+ Season1 ) printMinMax(tmp,0) cmin = min(tmp) cmax = max(tmp) maxlev = 25 mnmx = nice_mnmxintvl(cmin,cmax,maxlev,False) ensanom@cnMinLevelValF = mnmx(0) ensanom@cnMaxLevelValF = mnmx(1) ensanom@cnLevelSpacingF = mnmx(2) delete(ensanom@cnFillPalette) ensanom@cnFillPalette = coloranomaly symMinMaxPlt (tmp,maxlev,False,ensanom) cntr = gsn_csm_contour_map(wks,tmp,ensanom) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec,mnmx,cmin,cmax,maxlev/]) a = addfile(fensanom2,"r") fprec = a->ensmeananom tmp = dim_avg_n_Wrap(fprec,0) wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn2+"_"+ \ Season2 + ValidYr2 +"_Anomaly") print("Seasonal anom range for "+ Season2 ) printMinMax(tmp,0) cmin = min(tmp) cmax = max(tmp) maxlev = 25 mnmx = nice_mnmxintvl(cmin,cmax,maxlev,False) ensanom@cnMinLevelValF = mnmx(0) ensanom@cnMaxLevelValF = mnmx(1) ensanom@cnLevelSpacingF = mnmx(2) delete(ensanom@cnFillPalette) ensanom@cnFillPalette = coloranomaly symMinMaxPlt (tmp,maxlev,False,ensanom) cntr = gsn_csm_contour_map(wks,tmp,ensanom) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec,mnmx,cmin,cmax,maxlev/]) end if ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; 7. Plot Combined Probability Seasonal if (plotcombined) then abelow = addfile(fensbelow1,"r") aabove = addfile(fensabove1,"r") anormal= addfile(fensnormal1,"r") tbelow = abelow->tercileBelow tabove = aabove->tercileAbove tnormal= anormal->tercileNormal ;; Do we need to compute Chisquare significnce if (chisquare) then t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) tercileInd = t0 tercileInd = t0@_FillValue tercileProb = tercileInd ; printVarSummary(tercileProb) do iy = 0, mlat-1 do ix = 0, mlon-1 tmp = dimsizes(gcmtype) * \ ((t0(iy,ix) - 0.333) + \ (t1(iy,ix) - 0.333) + \ (t2(iy,ix) - 0.333))/0.333 indx = closest_val(tmp,chiF) if (.not.ismissing(indx) .and. indx .lt. dimsizes(chiF)) then if (chiP(indx) .lt. siglevel ) then tmp1 = ((/t0(iy,ix),t1(iy,ix), t2(iy,ix) /)) ; print(chiP(indx)+" "+maxind(tmp1)+" "+tmp1+" "+tmp1(maxind(tmp1))) tercileInd(iy,ix) = maxind(tmp1) tercileProb(iy,ix) = (/tmp1(maxind(tmp1))/) delete(tmp1) else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if delete([/tmp,indx/]) end do end do end if ; chisquare ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Start combine plot from below normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 0., 1.0, 0.) x = t0 * tmpx delete([tmpx]) else x = where(t0 .gt. t1 .and. t0 .gt. t2 .and. t0 .gt. pmin, \ t0, tbelow@_FillValue) end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensp@lbTitleString = "BN" ensp@lbTitleFontHeightF = 0.02 ensp@lbTitleJust = "CenterCenter" ensp@lbTitlePosition = "Top" ensp@lbTitleOffsetF = -0.05 ensp@lbLabelsOn = False ensp@pmLabelBarOrthogonalPosF = -0.06 ensp@lbBottomMarginF = 0.3 ensp@cnFillPalette = colorbelow ; ensp@cnConstFEnableFill = False ensp@cnConstFLabelOn = False ;ensp@cnMissingValFillColor = "gray" ensp@cnNoDataLabelOn = False wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn1+"_"+ \ Season1 + ValidYr1 +"_"+figProb+"Prob") cntrb = gsn_csm_contour_map(wks,x,ensp) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) delete([/t0,t1,t2,x,ensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Near Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 1., 1.0, 0.) x = t1 * tmpx delete(tmpx) else x = where(t1 .gt. t0 .and. t1 .gt. t2 .and. t1 .gt. pmin,\ t1, tnormal@_FillValue) end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp = True bensp@gsnDraw = False bensp@gsnFrame = False bensp@cnInfoLabelOn = False bensp@cnFillOn = True bensp@cnLinesOn = False bensp@lbOrientation = "Vertical" bensp@cnFillPalette = colornormal bensp@lbTitleString = "NN" bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterCenter" bensp@lbTitlePosition = "Top" bensp@lbTitleOffsetF = -0.05 bensp@lbLabelsOn = False bensp@cnLevelSelectionMode = "ManualLevels" bensp@cnMinLevelValF = pmin bensp@cnMaxLevelValF = pmax bensp@cnLevelSpacingF = plev bensp@pmLabelBarOrthogonalPosF = -0.115 bensp@lbBottomMarginF = 0.3 bensp@cnFillDrawOrder = "Predraw" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntrn = gsn_csm_contour(wks,x,bensp) delete([/t0,t1,t2,x,bensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Above Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 2., 1.0, 0.) x = t2 * tmpx delete(tmpx) else x = where(t2 .gt. t0 .and. t2 .gt. t1 .and. t2 .gt. pmin, \ t2, tabove@_FillValue) end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp@lbOrientation = "vertical" bensp@lbTitleString = " AN" bensp@lbTitleOffsetF = -0.05 bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterLeft" bensp@lbTitlePosition = "Top" bensp@lbLabelsOn = True bensp@pmLabelBarOrthogonalPosF = 0.0 bensp@cnFillPalette = colorabove bensp@lbBottomMarginF = 0.3 bensp@lbLeftMarginF = 0.1 ;bensp@cnMissingValFillColor = "gray" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntra = gsn_csm_contour(wks,x,bensp) ;;;;;;;;;;;;;;;;;;;;;;;;;; ;plot gray for precip less than raindry x = where(omean1 .ge. raindryszn , omean1@_FillValue, 0.) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon xensp = True xensp@gsnDraw = False xensp@gsnFrame = False xensp@cnInfoLabelOn = False xensp@cnFillOn = True xensp@cnLinesOn = False xensp@lbOrientation = "Vertical" xensp@cnFillPalette = (/"gray","gray"/) xensp@lbLabelsOn = False xensp@lbLabelBarOn = False xensp@cnFillDrawOrder = "Predraw" xensp@cnConstFEnableFill = True xensp@cnConstFLabelOn = False cntrx = gsn_csm_contour(wks,x,xensp) ;;;;;;;;;;;;;;;;;;;;;;; overlay(cntrb,cntrn) overlay(cntrb,cntra) overlay(cntrb,cntrx) draw(cntrb) frame(wks) if (chisquare) then delete([/tercileProb,tercileInd/]) end if delete([/abelow,aabove,anormal, tbelow, tabove,tnormal/]) delete([/t0,t1,t2,x, cntrb,cntra,cntrn,cntrx/]) ;;;;;;; End combined plot for Season2 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;; Do combined plot for Season2 abelow = addfile(fensbelow2,"r") aabove = addfile(fensabove2,"r") anormal= addfile(fensnormal2,"r") tbelow = abelow->tercileBelow tabove = aabove->tercileAbove tnormal= anormal->tercileNormal ;printMinMax(tnormal,0) ;print(num(tnormal .gt. pmin)) ;; Do we need to compute Chisquare significnce if (chisquare) then t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) tercileInd = t0 tercileInd = t0@_FillValue tercileProb = tercileInd ; printVarSummary(tercileProb) do iy = 0, mlat-1 do ix = 0, mlon-1 tmp = dimsizes(gcmtype) * \ ((t0(iy,ix) - 0.333) + \ (t1(iy,ix) - 0.333) + \ (t2(iy,ix) - 0.333))/0.333 indx = closest_val(tmp,chiF) if (.not.ismissing(indx) .and. indx .lt. dimsizes(chiF)) then if (chiP(indx) .lt. siglevel ) then tmp1 = ((/t0(iy,ix),t1(iy,ix), t2(iy,ix) /)) ; print(chiP(indx)+" "+maxind(tmp1)+" "+tmp1+" "+tmp1(maxind(tmp1))) tercileInd(iy,ix) = maxind(tmp1) tercileProb(iy,ix) = (/tmp1(maxind(tmp1))/) delete(tmp1) else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if delete([/tmp,indx/]) end do end do end if ; chisquare ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Start combine plot from below normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 0., 1.0, 0.) x = t0 * tmpx delete([tmpx]) else x = where(t0 .gt. t1 .and. t0 .gt. t2 .and. t0 .gt. pmin, \ t0, tbelow@_FillValue) end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensp@lbTitleString = "BN" ensp@lbTitleFontHeightF = 0.02 ensp@lbTitleJust = "CenterCenter" ensp@lbTitlePosition = "Top" ensp@lbTitleOffsetF = -0.05 ensp@lbLabelsOn = False ensp@pmLabelBarOrthogonalPosF = -0.06 ensp@lbBottomMarginF = 0.3 ensp@cnFillPalette = colorbelow ; ensp@cnConstFEnableFill = False ensp@cnConstFLabelOn = False ;ensp@cnMissingValFillColor = "gray" ensp@cnNoDataLabelOn = False wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn2+"_"+ \ Season2 + ValidYr2 +"_"+figProb+"Prob") cntrb = gsn_csm_contour_map(wks,x,ensp) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) delete([/t0,t1,t2,x,ensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Near Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 1., 1.0, 0.) x = t1 * tmpx delete(tmpx) else x = where(t1 .gt. t0 .and. t1 .gt. t2 .and. t1 .gt. pmin,\ t1, tnormal@_FillValue) end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp = True bensp@gsnDraw = False bensp@gsnFrame = False bensp@cnInfoLabelOn = False bensp@cnFillOn = True bensp@cnLinesOn = False bensp@lbOrientation = "Vertical" bensp@cnFillPalette = colornormal bensp@lbTitleString = "NN" bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterCenter" bensp@lbTitlePosition = "Top" bensp@lbTitleOffsetF = -0.05 bensp@lbLabelsOn = False bensp@cnLevelSelectionMode = "ManualLevels" bensp@cnMinLevelValF = pmin bensp@cnMaxLevelValF = pmax bensp@cnLevelSpacingF = plev bensp@pmLabelBarOrthogonalPosF = -0.115 bensp@lbBottomMarginF = 0.3 bensp@cnFillDrawOrder = "Predraw" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntrn = gsn_csm_contour(wks,x,bensp) delete([/t0,t1,t2,x,bensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Above Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = dim_avg_n_Wrap(tbelow,0) t1 = dim_avg_n_Wrap(tnormal,0) t2 = dim_avg_n_Wrap(tabove,0) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 2., 1.0, 0.) x = t2 * tmpx delete(tmpx) else x = where(t2 .gt. t0 .and. t2 .gt. t1 .and. t2 .gt. pmin, \ t2, tabove@_FillValue) end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp@lbOrientation = "vertical" bensp@lbTitleString = " AN" bensp@lbTitleOffsetF = -0.05 bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterLeft" bensp@lbTitlePosition = "Top" bensp@lbLabelsOn = True bensp@pmLabelBarOrthogonalPosF = 0.0 bensp@cnFillPalette = colorabove bensp@lbBottomMarginF = 0.3 bensp@lbLeftMarginF = 0.1 ;bensp@cnMissingValFillColor = "gray" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntra = gsn_csm_contour(wks,x,bensp) ;;;;;;;;;;;;;;;;;;;;;;;;;; ;plot gray for precip less than raindry x = where(omean2 .ge. raindryszn , omean2@_FillValue, 0.) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon xensp = True xensp@gsnDraw = False xensp@gsnFrame = False xensp@cnInfoLabelOn = False xensp@cnFillOn = True xensp@cnLinesOn = False xensp@lbOrientation = "Vertical" xensp@cnFillPalette = (/"gray","gray"/) xensp@lbLabelsOn = False xensp@lbLabelBarOn = False xensp@cnFillDrawOrder = "Predraw" xensp@cnConstFEnableFill = True xensp@cnConstFLabelOn = False cntrx = gsn_csm_contour(wks,x,xensp) ;;;;;;;;;;;;;;;;;;;;;;; overlay(cntrb,cntrn) overlay(cntrb,cntra) overlay(cntrb,cntrx) draw(cntrb) frame(wks) if (chisquare) then delete([/tercileProb,tercileInd/]) end if delete([/abelow,aabove,anormal, tbelow, tabove,tnormal/]) delete([/t0,t1,t2,x, cntrb,cntra,cntrn,cntrx/]) ;;;;;; End combined plot for Season 2 end if ; plotcombined ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Individual probability plots ;;; Above normal rainfall probability ;;; 8. Plot Individual probabilities Seasonal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; if (plotindividual) then ;;;; Do Season1 Above a = addfile(fensabove1,"r") tercile = a->tercileAbove x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensprob@cnFillPalette = colorabove wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn1+"_"+ \ Season1 + ValidYr1 +"_"+figProb+"AboveNormal") cntra = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntra,dirSHP+f1,lnres) draw(cntra) frame(wks) delete([/a,x,tercile/]) ;;; Do Season2 Above a = addfile(fensabove2,"r") tercile = a->tercileAbove x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensprob@cnFillPalette = colorabove wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn2+"_"+ \ Season2 + ValidYr2 +"_"+figProb+"AboveNormal") cntra = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntra,dirSHP+f1,lnres) draw(cntra) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Normal probability ;;; Do for Season1 Normal a = addfile(fensnormal1,"r") tercile = a->tercileNormal x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colornormal wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn1+"_"+ \ Season1 + ValidYr1 +"_"+figProb+"NearNormal") cntrn = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrn,dirSHP+f1,lnres) draw(cntrn) frame(wks) delete([/a,x,tercile/]) ;;;; Do for Season2 Normal a = addfile(fensnormal2,"r") tercile = a->tercileNormal x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colornormal wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn2+"_"+ \ Season2 + ValidYr2 +"_"+figProb+"NearNormal") cntrn = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrn,dirSHP+f1,lnres) draw(cntrn) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Below probability ;;;; Do for Season1 Below a = addfile(fensbelow1,"r") tercile = a->tercileBelow x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colorbelow wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn1+"_"+ \ Season1 + ValidYr1 +"_"+figProb+"BelowNormal") cntrb = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) draw(cntrb) frame(wks) delete([/a,x,tercile/]) ;;;; Do for Season2 Below a = addfile(fensbelow2,"r") tercile = a->tercileBelow x = dim_avg_n(tercile,0) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colorbelow wks = gsn_open_wks(type,plotdir+"Ensemble_"+VARN+"_"+LNszn2+"_"+ \ Season2 + ValidYr2 +"_"+figProb+"BelowNormal") cntrb = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) draw(cntrb) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;; end if ; plotindividual seasonal ;;;;;;;;;;;;; ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;; Loop through gcmtype for Seasonal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; if (plotOneModel) then do imodel = 0, dimsizes(gcmtype)-1 ;; 1. Plot Multimodel Ensemble Average rainfall Seasonal if (plotmean) then a = addfile(fensmean1,"r") fprec = a->ensmeanunbiased tmp = fprec(imodel,:,:) wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn1+"_"+ Season1 + ValidYr1 +"_Total") cntr = gsn_csm_contour_map(wks,tmp,ensmean) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec/]) a = addfile(fensmean2,"r") fprec = a->ensmeanunbiased tmp = fprec(imodel,:,:) wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn2+"_"+ Season2 + ValidYr2 +"_Total") cntr = gsn_csm_contour_map(wks,tmp,ensmean) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec/]) end if ;; 2. Plot Ensemble mean anomaly Seasonal if (plotanom) then a = addfile(fensanom1,"r") fprec = a->ensmeananom tmp = fprec(imodel,:,:) wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn1+"_"+ Season1 + ValidYr1 +"_Anomaly") print("Seasonal anom range for "+ Season1 +" for "+gcmtype(imodel)) printMinMax(tmp,0) cmin = min(tmp) cmax = max(tmp) maxlev = 25 mnmx = nice_mnmxintvl(cmin,cmax,maxlev,False) ensanom@cnMinLevelValF = mnmx(0) ensanom@cnMaxLevelValF = mnmx(1) ensanom@cnLevelSpacingF = mnmx(2) delete(ensanom@cnFillPalette) ensanom@cnFillPalette = coloranomaly symMinMaxPlt (tmp,maxlev,False,ensanom) cntr = gsn_csm_contour_map(wks,tmp,ensanom) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec,mnmx,cmin,cmax,maxlev/]) a = addfile(fensanom2,"r") fprec = a->ensmeananom tmp = fprec(imodel,:,:) wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn2+"_"+ Season2 + ValidYr2 +"_Anomaly") print("Seasonal anom range for "+ Season2 +" for "+gcmtype(imodel)) printMinMax(tmp,0) cmin = min(tmp) cmax = max(tmp) maxlev = 25 mnmx = nice_mnmxintvl(cmin,cmax,maxlev,False) ensanom@cnMinLevelValF = mnmx(0) ensanom@cnMaxLevelValF = mnmx(1) ensanom@cnLevelSpacingF = mnmx(2) delete(ensanom@cnFillPalette) ensanom@cnFillPalette = coloranomaly symMinMaxPlt (tmp,maxlev,False,ensanom) cntr = gsn_csm_contour_map(wks,tmp,ensanom) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntr,dirSHP+f1,lnres) draw(cntr) frame(wks) delete([/wks,cntr,a,tmp,fprec,mnmx,cmin,cmax,maxlev/]) end if ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 3. Plot Compined Probabilities Seasonal ;;; if (plotcombined) then abelow = addfile(fensbelow1,"r") aabove = addfile(fensabove1,"r") anormal= addfile(fensnormal1,"r") tbelow = abelow->tercileBelow tabove = aabove->tercileAbove tnormal= anormal->tercileNormal if (RNORM .and. RPROB) then asiglev = addfile(fenspsig1,"r") psigval = asiglev->psig end if ;; Do we need to compute Chisquare significnce if (chisquare) then t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) tercileInd = t0 tercileInd = t0@_FillValue tercileProb = tercileInd do iy = 0, mlat-1 do ix = 0, mlon-1 tmp = ensmem(imodel) * \ ((t0(iy,ix) - 0.333) + \ (t1(iy,ix) - 0.333) + \ (t2(iy,ix) - 0.333))/0.333 indx = closest_val(tmp,chiF) if (.not.ismissing(indx) .and. indx .lt. dimsizes(chiF)) then if (chiP(indx) .lt. siglevel ) then tmp1 = ((/t0(iy,ix),t1(iy,ix), t2(iy,ix) /)) ; print(chiP(indx)+" "+maxind(tmp1)+" "+tmp1+" "+tmp1(maxind(tmp1))) tercileInd(iy,ix) = maxind(tmp1) tercileProb(iy,ix) = (/tmp1(maxind(tmp1))/) delete(tmp1) else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if delete([/tmp,indx/]) end do end do end if ; chisquare ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Start combine plot from below normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 0., 1.0, 0.) x = t0 * tmpx delete([tmpx]) else x = where(t0 .gt. t1 .and. t0 .gt. t2 .and. t0 .gt. pmin, \ t0, tbelow@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensp@lbTitleString = "BN" ensp@lbTitleFontHeightF = 0.02 ensp@lbTitleJust = "CenterCenter" ensp@lbTitlePosition = "Top" ensp@lbTitleOffsetF = -0.05 ensp@lbLabelsOn = False ensp@pmLabelBarOrthogonalPosF = -0.06 ensp@lbBottomMarginF = 0.3 ensp@cnFillPalette = colorbelow ; ensp@cnConstFEnableFill = False ensp@cnConstFLabelOn = False ;ensp@cnMissingValFillColor = "gray" ensp@cnNoDataLabelOn = False wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn1+"_"+ Season1 + ValidYr1 +"_"+figProb+"Prob") cntrb = gsn_csm_contour_map(wks,x,ensp) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) delete([/t0,t1,t2,x,ensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Near Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 1., 1.0, 0.) x = t1 * tmpx delete(tmpx) else x = where(t1 .gt. t0 .and. t1 .gt. t2 .and. t1 .gt. pmin,\ t1, tnormal@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp = True bensp@gsnDraw = False bensp@gsnFrame = False bensp@cnInfoLabelOn = False bensp@cnFillOn = True bensp@cnLinesOn = False bensp@lbOrientation = "Vertical" bensp@cnFillPalette = colornormal bensp@lbTitleString = "NN" bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterCenter" bensp@lbTitlePosition = "Top" bensp@lbTitleOffsetF = -0.05 bensp@lbLabelsOn = False bensp@cnLevelSelectionMode = "ManualLevels" bensp@cnMinLevelValF = pmin bensp@cnMaxLevelValF = pmax bensp@cnLevelSpacingF = plev bensp@pmLabelBarOrthogonalPosF = -0.115 bensp@lbBottomMarginF = 0.3 bensp@cnFillDrawOrder = "Predraw" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntrn = gsn_csm_contour(wks,x,bensp) delete([/t0,t1,t2,x,bensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Above Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 2., 1.0, 0.) x = t2 * tmpx delete(tmpx) else x = where(t2 .gt. t0 .and. t2 .gt. t1 .and. t2 .gt. pmin, \ t2, tabove@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp@lbOrientation = "vertical" bensp@lbTitleString = " AN" bensp@lbTitleOffsetF = -0.05 bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterLeft" bensp@lbTitlePosition = "Top" bensp@lbLabelsOn = True bensp@pmLabelBarOrthogonalPosF = 0.0 bensp@cnFillPalette = colorabove bensp@lbBottomMarginF = 0.3 bensp@lbLeftMarginF = 0.1 ;bensp@cnMissingValFillColor = "gray" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntra = gsn_csm_contour(wks,x,bensp) ;;;;;;;;;;;;;;;;;;;;;;;;;; ;plot gray for precip less than raindrymon x = where(omean1 .ge. raindrymon , omean1@_FillValue, 0.) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon xensp = True xensp@gsnDraw = False xensp@gsnFrame = False xensp@cnInfoLabelOn = False xensp@cnFillOn = True xensp@cnLinesOn = False xensp@lbOrientation = "Vertical" xensp@cnFillPalette = (/"gray","gray"/) xensp@lbLabelsOn = False xensp@lbLabelBarOn = False xensp@cnFillDrawOrder = "Predraw" xensp@cnConstFEnableFill = True xensp@cnConstFLabelOn = False cntrx = gsn_csm_contour(wks,x,xensp) ;;;;;;;;;;;;;;;;;;;;;;; overlay(cntrb,cntrn) overlay(cntrb,cntra) overlay(cntrb,cntrx) draw(cntrb) frame(wks) if (chisquare) then delete([/tercileProb,tercileInd/]) end if delete([/abelow,aabove,anormal, tbelow, tabove,tnormal/]) delete([/t0,t1,t2,x, cntrb,cntra,cntrn,cntrx/]) if (RNORM .and. RPROB) then delete([/asiglev,psigval /]) end if ;;;;;;; End combined plot for Season 1 ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;;;;;;;; Do combined plot for Season2 abelow = addfile(fensbelow2,"r") aabove = addfile(fensabove2,"r") anormal= addfile(fensnormal2,"r") tbelow = abelow->tercileBelow tabove = aabove->tercileAbove tnormal= anormal->tercileNormal if (RNORM .and. RPROB) then asiglev = addfile(fenspsig2,"r") psigval = asiglev->psig end if ;; Do we need to compute Chisquare significnce if (chisquare) then t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) tercileInd = t0 tercileInd = t0@_FillValue tercileProb = tercileInd do iy = 0, mlat-1 do ix = 0, mlon-1 tmp = ensmem(imodel) * \ ((t0(iy,ix) - 0.333) + \ (t1(iy,ix) - 0.333) + \ (t2(iy,ix) - 0.333))/0.333 indx = closest_val(tmp,chiF) if (.not.ismissing(indx) .and. indx .lt. dimsizes(chiF)) then if (chiP(indx) .lt. siglevel ) then tmp1 = ((/t0(iy,ix),t1(iy,ix), t2(iy,ix) /)) ; print(chiP(indx)+" "+maxind(tmp1)+" "+tmp1+" "+tmp1(maxind(tmp1))) tercileInd(iy,ix) = maxind(tmp1) tercileProb(iy,ix) = (/tmp1(maxind(tmp1))/) delete(tmp1) else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if else tercileInd(iy,ix) = tercileInd@_FillValue tercileProb(iy,ix) = tercileProb@_FillValue end if delete([/tmp,indx/]) end do end do end if ; chisquare ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Start combine plot from below normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 0., 1.0, 0.) x = t0 * tmpx delete([tmpx]) else x = where(t0 .gt. t1 .and. t0 .gt. t2 .and. t0 .gt. pmin, \ t0, tbelow@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if ; chisquare x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensp@lbTitleString = "BN" ensp@lbTitleFontHeightF = 0.02 ensp@lbTitleJust = "CenterCenter" ensp@lbTitlePosition = "Top" ensp@lbTitleOffsetF = -0.05 ensp@lbLabelsOn = False ensp@pmLabelBarOrthogonalPosF = -0.06 ensp@lbBottomMarginF = 0.3 ensp@cnFillPalette = colorbelow ; ensp@cnConstFEnableFill = False ensp@cnConstFLabelOn = False ;ensp@cnMissingValFillColor = "gray" ensp@cnNoDataLabelOn = False wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn2+"_"+ Season2 + ValidYr2 +"_"+figProb+"Prob") cntrb = gsn_csm_contour_map(wks,x,ensp) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) delete([/t0,t1,t2,x,ensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Near Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 1., 1.0, 0.) x = t1 * tmpx delete(tmpx) else x = where(t1 .gt. t0 .and. t1 .gt. t2 .and. t1 .gt. pmin,\ t1, tnormal@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp = True bensp@gsnDraw = False bensp@gsnFrame = False bensp@cnInfoLabelOn = False bensp@cnFillOn = True bensp@cnLinesOn = False bensp@lbOrientation = "Vertical" bensp@cnFillPalette = colornormal bensp@lbTitleString = "NN" bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterCenter" bensp@lbTitlePosition = "Top" bensp@lbTitleOffsetF = -0.05 bensp@lbLabelsOn = False bensp@cnLevelSelectionMode = "ManualLevels" bensp@cnMinLevelValF = pmin bensp@cnMaxLevelValF = pmax bensp@cnLevelSpacingF = plev bensp@pmLabelBarOrthogonalPosF = -0.115 bensp@lbBottomMarginF = 0.3 bensp@cnFillDrawOrder = "Predraw" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntrn = gsn_csm_contour(wks,x,bensp) delete([/t0,t1,t2,x,bensp@cnFillPalette/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; Now do Above Normal ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; t0 = tbelow(imodel,:,:) t1 = tnormal(imodel,:,:) t2 = tabove(imodel,:,:) if (chisquare) then tmpx = tercileInd tmpx = where(tmpx .eq. 2., 1.0, 0.) x = t2 * tmpx delete(tmpx) else x = where(t2 .gt. t0 .and. t2 .gt. t1 .and. t2 .gt. pmin, \ t2, tabove@_FillValue) if (RNORM .and. RPROB) then pval = psigval(imodel,:,:) pval = where(pval .le. 0.05, 1.0, 0.0) x = x * pval delete(pval) end if end if x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon bensp@lbOrientation = "vertical" bensp@lbTitleString = " AN" bensp@lbTitleOffsetF = -0.05 bensp@lbTitleFontHeightF = 0.02 bensp@lbTitleJust = "CenterLeft" bensp@lbTitlePosition = "Top" bensp@lbLabelsOn = True bensp@pmLabelBarOrthogonalPosF = 0.0 bensp@cnFillPalette = colorabove bensp@lbBottomMarginF = 0.3 bensp@lbLeftMarginF = 0.1 ;bensp@cnMissingValFillColor = "gray" bensp@cnNoDataLabelOn = False bensp@cnConstFEnableFill = False bensp@cnConstFLabelOn = False cntra = gsn_csm_contour(wks,x,bensp) ;;;;;;;;;;;;;;;;;;;;;;;;;; ;plot gray for precip less than raindrymon x = where(omean2 .ge. raindrymon , omean2@_FillValue, 0.) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon xensp = True xensp@gsnDraw = False xensp@gsnFrame = False xensp@cnInfoLabelOn = False xensp@cnFillOn = True xensp@cnLinesOn = False xensp@lbOrientation = "Vertical" xensp@cnFillPalette = (/"gray","gray"/) xensp@lbLabelsOn = False xensp@lbLabelBarOn = False xensp@cnFillDrawOrder = "Predraw" xensp@cnConstFEnableFill = True xensp@cnConstFLabelOn = False cntrx = gsn_csm_contour(wks,x,xensp) ;;;;;;;;;;;;;;;;;;;;;;; overlay(cntrb,cntrn) overlay(cntrb,cntra) overlay(cntrb,cntrx) draw(cntrb) frame(wks) if (chisquare) then delete([/tercileProb,tercileInd/]) end if delete([/abelow,aabove,anormal, tbelow, tabove,tnormal/]) delete([/t0,t1,t2,x, cntrb,cntra,cntrn,cntrx/]) if (RNORM .and. RPROB) then delete([/asiglev,psigval /]) end if ;;;;;; End combined plot for Season 2 end if ; plotcombined ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Individual probability plots ;;; Above normal rainfall probability ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 4. Plot Individual probabilities Seasonal if (plotindividual) then ;;;; Do Season1 Above a = addfile(fensabove1,"r") tercile = a->tercileAbove x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensprob@cnFillPalette = colorabove wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn1+"_"+ Season1 + ValidYr1 +"_"+figProb+"AboveNormal") cntra = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntra,dirSHP+f1,lnres) draw(cntra) frame(wks) delete([/a,x,tercile/]) ;;; Do Season2 Above a = addfile(fensabove2,"r") tercile = a->tercileAbove x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon ensprob@cnFillPalette = colorabove wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn2+"_"+ Season2 + ValidYr2 +"_"+figProb+"AboveNormal") cntra = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntra,dirSHP+f1,lnres) draw(cntra) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Normal probability ;;; Do for Season1 Normal a = addfile(fensnormal1,"r") tercile = a->tercileNormal x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colornormal wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn1+"_"+ Season1 + ValidYr1 +"_"+figProb+"NearNormal") cntrn = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrn,dirSHP+f1,lnres) draw(cntrn) frame(wks) delete([/a,x,tercile/]) ;;;; Do for Season2 Normal a = addfile(fensnormal2,"r") tercile = a->tercileNormal x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colornormal wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn2+"_"+ Season2 + ValidYr2 +"_"+figProb+"NearNormal") cntrn = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrn,dirSHP+f1,lnres) draw(cntrn) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;;; Below probability ;;;; Do for Season1 Below a = addfile(fensbelow1,"r") tercile = a->tercileBelow x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colorbelow wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn1+"_"+ Season1 + ValidYr1 +"_"+figProb+"BelowNormal") cntrb = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) draw(cntrb) frame(wks) delete([/a,x,tercile/]) ;;;; Do for Season2 Below a = addfile(fensbelow2,"r") tercile = a->tercileBelow x = tercile(imodel,:,:) x!0 = "lat" x!1 = "lon" x&lat = GHAlat x&lon = GHAlon delete(ensprob@cnFillPalette) ensprob@cnFillPalette = colorbelow wks = gsn_open_wks(type,plotdir+gcmtype(imodel)+"_"+VARN+"_"+ \ LNszn2+"_"+ Season2 + ValidYr2 +"_"+figProb+"BelowNormal") cntrb = gsn_csm_contour_map(wks,x,ensprob) lnres = True lnres@gsLineColor = "black" lnres@gsLineThicknessF = 2.0 poly1 = gsn_add_shapefile_polylines(wks,cntrb,dirSHP+f1,lnres) draw(cntrb) frame(wks) delete([/a,x,tercile/]) ;;;;;;;;;;;;;;;;;;; end if ; plotindividual monthly ;;;;;;;;;; end do ; imodel Seasonal end if ; plotOneModel delete([/omean1,omean2, ensmean,ensanom,ensp,ensprob/]) end