<html dir="ltr">
<head>
<meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
<style type="text/css" id="owaParaStyle">
<!--
p
{margin-top:0;
margin-bottom:0}
-->
P {margin-top:0;margin-bottom:0;}</style>
</head>
<body ocsi="0" fpstyle="1">
<div style="direction: ltr;font-family: Tahoma;color: #000000;font-size: 10pt;">Hi All
<br>
<br>
Please ignore previous email below. I have managed to figure the area-averaging out but have some concerns regarding the overall code. I am happy that the code does what I would like for dP but not for the other components?<br>
<br>
My overall aim is to: <br>
<br>
1. mask out all ocean regions<br>
2. mask out all the wetting regions<br>
3. area average these changes over southern africa<br>
4. create boxplots of dP and 6 decomposed precipitation components<br>
<br>
I feel the code does this correctly for dP, but not convinced for the other components.<br>
<br>
Any advice would very much be appreciated regarding this issue. Below is the full code to tackle task 1 -4<br>
<br>
Many thanks in advance!<br>
<br>
Kindest Regards<br>
Melissa<br>
<br>
;*********************************************<br>
; boxplot.ncl<br>
;<br>
; Concepts illustrated:<br>
; - Drawing box plots<br>
; - Setting the color of individual boxes in a box plot<br>
; - Setting the width of individual boxes in a box plot<br>
;<br>
;*********************************************<br>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"<br>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"<br>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"<br>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"<br>
;*********************************************<br>
begin<br>
;**********************************************<br>
; Add data in text format and sort in ascending order<br>
;**********************************************<br>
<br>
latS = -30. <br>
latN = 0. <br>
lonL = 10. <br>
lonR = 40.<br>
<br>
model = (/"bcc-csm1-1-m","BNU-ESM","CanESM2","CCSM4","CESM1-BGC","CESM1-CAM5","CSIRO-Mk3-6-0","FIO-ESM","GFDL-CM3","GFDL-ESM2G","GFDL-ESM2M","GISS-E2-H","HadGEM2-CC","HadGEM2-ES","IPSL-CM5A-LR","IPSL-CM5A-MR","MIROC5","MRI-CGCM3","NorESM1-M","NorESM1-ME"/)<br>
<br>
<br>
do gg = 0,dimsizes(model)-1<br>
<br>
;dP<br>
in1=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadprfutureOND/Delta_pr_"+model(gg)+".nc","r")
<br>
<br>
;dPShift<br>
in2=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPSpat_"+model(gg)+"_safrica_PSpat.nc","r")<br>
;dPT<br>
in3=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPT_"+model(gg)+"_safrica_PT.nc","r")<br>
;dPWeak<br>
in4=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPDiv_"+model(gg)+"_safrica_PDiv.nc","r")<br>
;dPTWRH<br>
in5=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPDiv+PT+PRH_"+model(gg)+".nc","r")<br>
;dPCross<br>
in6=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPNL_"+model(gg)+"_safrica_PNL.nc","r")<br>
;dPRH<br>
in7=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPRH_"+model(gg)+"_safrica_PRH.nc","r")<br>
<br>
<br>
lat = in1->lat <br>
lon = in1->lon <br>
time = in1->time<br>
<br>
if ((gg.eq.0) .or. (gg.eq.9) .or. (gg.eq.10) .or. (gg.eq.11)) then<br>
dP = dble2flt(in1->pr(:,:,:))<br>
dPShift = in2->huss(:,:,:)<br>
dPtwrh = in5->huss(:,:,:)<br>
dPcross = in6->huss(:,:,:)<br>
dPT = in3->es2(:,:,:)<br>
dPweak = in4->huss(:,:,:)<br>
dPrh = in7->huss(:,:,:)<br>
<br>
else<br>
<br>
dP = in1->pr(:,0,:,:)<br>
;dPShift = in2->huss(:,0,:,:)<br>
;dPtwrh = in5->huss(:,0,:,:)<br>
;dPcross = in6->huss(:,0,:,:)<br>
;dPT = in3->es2(:,0,:,:)<br>
;dPweak = in4->huss(:,0,:,:)<br>
;dPrh = in7->huss(:,0,:,:)<br>
<br>
end if<br>
<br>
;Land mask<br>
<br>
b=addfile("$NCARG_ROOT/lib/ncarg/data/cdf/landsea.nc","r")<br>
lsm_dP = landsea_mask(b->LSMASK,dP&lat,dP&lon) ; read in land sea mask<br>
lsm_dPShift = landsea_mask(b->LSMASK,dPShift&lat,dPShift&lon) ; read in land sea mask<br>
lsm_dPtwrh = landsea_mask(b->LSMASK,dPtwrh&lat,dPtwrh&lon) ; read in land sea mask<br>
lsm_dPcross = landsea_mask(b->LSMASK,dPcross&lat,dPcross&lon) ; read in land sea mask<br>
lsm_dPT = landsea_mask(b->LSMASK,dPT&lat,dPT&lon) ; read in land sea mask<br>
lsm_dPweak = landsea_mask(b->LSMASK,dPweak&lat,dPweak&lon) ; read in land sea mask<br>
lsm_dPrh = landsea_mask(b->LSMASK,dPrh&lat,dPrh&lon) ; read in land sea mask<br>
<br>
dP = mask(dP,lsm_dP.eq.0,False) ;mask out all ocean points)<br>
dPShift = mask(dPShift,lsm_dPShift.eq.0,False) ;mask out all ocean points)<br>
dPtwrh = mask(dPtwrh,lsm_dPtwrh.eq.0,False) ;mask out all ocean points)<br>
dPcross = mask(dPcross,lsm_dPcross.eq.0,False) ;mask out all ocean points)<br>
dPT = mask(dPT,lsm_dPT.eq.0,False) ;mask out all ocean points)<br>
dPweak = mask(dPweak,lsm_dPweak.eq.0,False) ;mask out all ocean points)<br>
dPrh = mask(dPrh,lsm_dPrh.eq.0,False) ;mask out all ocean points)<br>
<br>
<br>
;Drymask<br>
dP_dry = new ((/dimsizes(time),dimsizes(lat),dimsizes(lon)/), typeof(dP))<br>
dP_dry!0="time"<br>
dP_dry&time=time<br>
dP_dry!1="lat"<br>
dP_dry&lat=lat<br>
dP_dry!2="lon"<br>
dP_dry&lon=lon<br>
<br>
dP_dry = mask(dP,dP.gt.0, False) ; mask out all the wetting regions i.e. any values above 0 for dP<br>
<br>
dPShift_dry = new ((/dimsizes(time),dimsizes(lat),dimsizes(lon)/), typeof(dPShift))<br>
dPShift_dry!0="time"<br>
dPShift_dry&time=time<br>
dPShift_dry!1="lat"<br>
dPShift_dry&lat=lat<br>
dPShift_dry!2="lon"<br>
dPShift_dry&lon=lon<br>
<br>
dPShift_dry = mask(dPShift,dPShift.gt.0, False) ; mask out all the wetting regions i.e. any values above 0 for dP<br>
<br>
dPtwrh_dry = new ((/dimsizes(time),dimsizes(lat),dimsizes(lon)/), typeof(dPtwrh))<br>
dPtwrh_dry!0="time"<br>
dPtwrh_dry&time=time<br>
dPtwrh_dry!1="lat"<br>
dPtwrh_dry&lat=lat<br>
dPtwrh_dry!2="lon"<br>
dPtwrh_dry&lon=lon<br>
<br>
dPtwrh_dry = mask(dPtwrh,dPtwrh.gt.0, False) ; mask out all the wetting regions i.e. any values above 0 for dP<br>
<br>
dPcross_dry = new ((/dimsizes(time),dimsizes(lat),dimsizes(lon)/), typeof(dPcross))<br>
dPcross_dry!0="time"<br>
dPcross_dry&time=time<br>
dPcross_dry!1="lat"<br>
dPcross_dry&lat=lat<br>
dPcross_dry!2="lon"<br>
dPcross_dry&lon=lon<br>
<br>
dPcross_dry = mask(dPcross,dPcross.gt.0, False) ; mask out all the wetting regions i.e. any values above 0 for dP<br>
<br>
dPT_dry = new ((/dimsizes(time),dimsizes(lat),dimsizes(lon)/), typeof(dPT))<br>
dPT_dry!0="time"<br>
dPT_dry&time=time<br>
dPT_dry!1="lat"<br>
dPT_dry&lat=lat<br>
dPT_dry!2="lon"<br>
dPT_dry&lon=lon<br>
<br>
dPT_dry = mask(dPT,dPT.gt.0, False) ; mask out all the wetting regions i.e. any values above 0 for dP<br>
<br>
dPweak_dry = new ((/dimsizes(time),dimsizes(lat),dimsizes(lon)/), typeof(dPweak))<br>
dPweak_dry!0="time"<br>
dPweak_dry&time=time<br>
dPweak_dry!1="lat"<br>
dPweak_dry&lat=lat<br>
dPweak_dry!2="lon"<br>
dPweak_dry&lon=lon<br>
<br>
dPweak_dry = mask(dPweak,dPweak.gt.0, False) ; mask out all the wetting regions i.e. any values above 0 for dP<br>
<br>
dPrh_dry = new ((/dimsizes(time),dimsizes(lat),dimsizes(lon)/), typeof(dPrh))<br>
dPrh_dry!0="time"<br>
dPrh_dry&time=time<br>
dPrh_dry!1="lat"<br>
dPrh_dry&lat=lat<br>
dPrh_dry!2="lon"<br>
dPrh_dry&lon=lon<br>
<br>
dPrh_dry = mask(dPrh,dPrh.gt.0, False) ; mask out all the wetting regions i.e. any values above 0 for dP<br>
<br>
<br>
;if ((gg.eq.0) .or. (gg.eq.9) .or. (gg.eq.10) .or. (gg.eq.11)) then<br>
<br>
dP_fldmean = wgt_areaave(dP_dry(:,{-30:0},{10:40}),1.0, 1.0, 0)<br>
dPShift_fldmean = wgt_areaave(dPShift_dry(:,{-30:0},{10:40}),1.0, 1.0, 0)<br>
dPtwrh_fldmean = wgt_areaave(dPtwrh_dry(:,{-30:0},{10:40}),1.0, 1.0, 0)<br>
dPcross_fldmean = wgt_areaave(dPcross_dry(:,{-30:0},{10:40}),1.0, 1.0, 0)<br>
dPT_fldmean = wgt_areaave(dPT_dry(:,{-30:0},{10:40}),1.0, 1.0, 0)<br>
dPweak_fldmean = wgt_areaave(dPweak_dry(:,{-30:0},{10:40}),1.0, 1.0, 0)<br>
dPrh_fldmean = wgt_areaave(dPrh_dry(:,{-30:0},{10:40}),1.0, 1.0, 0)<br>
<br>
;else<br>
<br>
;dP_fldmean = wgt_areaave(dP_dry(:,0,{-30:0},{10:40}),1.0, 1.0, 0)<br>
;dPShift_fldmean = wgt_areaave(dPShift_dry(:,0,{-30:0},{10:40}),1.0, 1.0, 0)<br>
;dPtwrh_fldmean = wgt_areaave(dPtwrh_dry(:,0,{-30:0},{10:40}),1.0, 1.0, 0)<br>
;dPcross_fldmean = wgt_areaave(dPcross_dry(:,0,{-30:0},{10:40}),1.0, 1.0, 0)<br>
;dPT_fldmean = wgt_areaave(dPT_dry(:,{-30:0},0,{10:40}),1.0, 1.0, 0)<br>
;dPweak_fldmean = wgt_areaave(dPweak_dry(:,0,{-30:0},{10:40}),1.0, 1.0, 0)<br>
;dPrh_fldmean = wgt_areaave(dPrh_dry(:,0,{-30:0},{10:40}),1.0, 1.0, 0)<br>
<br>
<br>
;end if <br>
<br>
print(dP_fldmean)<br>
print(dPShift_fldmean)<br>
<br>
<br>
;Create new arrays of 20x40x240<br>
<br>
if (gg.eq.0) then<br>
dPfldavg = new ((/20/), typeof(dP_fldmean))<br>
<br>
dPfldavg!0="model"<br>
dPfldavg&model=model<br>
<br>
end if<br>
<br>
dPfldavg(gg) = dP_fldmean(0)<br>
<br>
if (gg.eq.0) then<br>
dPfldavgShift = new ((/20/), typeof(dPShift_fldmean))<br>
<br>
dPfldavgShift!0="model"<br>
dPfldavgShift&model=model<br>
<br>
end if<br>
<br>
dPfldavgShift(gg) = dPShift_fldmean(0)<br>
<br>
if (gg.eq.0) then<br>
dPfldavgtwrh = new ((/20/), typeof(dPtwrh_fldmean))<br>
<br>
dPfldavgtwrh!0="model"<br>
dPfldavgtwrh&model=model<br>
<br>
end if<br>
<br>
dPfldavgtwrh(gg) = dPtwrh_fldmean(0)<br>
<br>
if (gg.eq.0) then<br>
dPfldavgcross = new ((/20/), typeof(dPcross_fldmean))<br>
<br>
dPfldavgcross!0="model"<br>
dPfldavgcross&model=model<br>
<br>
end if<br>
<br>
dPfldavgcross(gg) = dPcross_fldmean(0)<br>
<br>
if (gg.eq.0) then<br>
dPfldavgT = new ((/20/), typeof(dPT_fldmean))<br>
<br>
dPfldavgT!0="model"<br>
dPfldavgT&model=model<br>
<br>
end if<br>
<br>
dPfldavgT(gg) = dPT_fldmean(0)<br>
<br>
if (gg.eq.0) then<br>
dPfldavgweak = new ((/20/), typeof(dPweak_fldmean))<br>
<br>
dPfldavgweak!0="model"<br>
dPfldavgweak&model=model<br>
<br>
end if<br>
<br>
dPfldavgweak(gg) = dPweak_fldmean(0)<br>
<br>
if (gg.eq.0) then<br>
dPfldavgrh = new ((/20/), typeof(dPrh_fldmean))<br>
<br>
dPfldavgrh!0="model"<br>
dPfldavgrh&model=model<br>
<br>
end if<br>
<br>
dPfldavgrh(gg) = dPrh_fldmean(0)<br>
<br>
;end do<br>
<br>
<br>
dimt = dimsizes(dPfldavg) ; should be 20<br>
x25 = round(.25*dimt,3)-1 ; -1 to account for NCL indexing starting<br>
x75 = round(.75*dimt,3)-1 ; at 0 <br>
qsort(dPfldavg)<br>
print(dPfldavg) ; Print the values<br>
<br>
dimt1 = dimsizes(dPfldavgShift) ; should be 20<br>
x25 = round(.25*dimt1,3)-1 ; -1 to account for NCL indexing starting<br>
x75 = round(.75*dimt1,3)-1 ; at 0 <br>
qsort(dPfldavgShift)<br>
print(dPfldavgShift) ; Print the values<br>
<br>
dimt2 = dimsizes(dPfldavgtwrh) ; should be 20<br>
x25 = round(.25*dimt2,3)-1 ; -1 to account for NCL indexing starting<br>
x75 = round(.75*dimt2,3)-1 ; at 0 <br>
qsort(dPfldavgtwrh)<br>
print(dPfldavgtwrh) ; Print the values<br>
<br>
dimt3 = dimsizes(dPfldavgT) ; should be 20<br>
x25 = round(.25*dimt3,3)-1 ; -1 to account for NCL indexing starting<br>
x75 = round(.75*dimt3,3)-1 ; at 0 <br>
qsort(dPfldavgT)<br>
print(dPfldavgT) ; Print the values<br>
<br>
dimt4 = dimsizes(dPfldavgcross) ; should be 20<br>
x25 = round(.25*dimt4,3)-1 ; -1 to account for NCL indexing starting<br>
x75 = round(.75*dimt4,3)-1 ; at 0 <br>
qsort(dPfldavgcross)<br>
print(dPfldavgcross) ; Print the values<br>
<br>
dimt5 = dimsizes(dPfldavgweak) ; should be 20<br>
x25 = round(.25*dimt5,3)-1 ; -1 to account for NCL indexing starting<br>
x75 = round(.75*dimt5,3)-1 ; at 0 <br>
qsort(dPfldavgweak)<br>
print(dPfldavgweak) ; Print the values<br>
<br>
dimt6 = dimsizes(dPfldavgrh) ; should be 20<br>
x25 = round(.25*dimt6,3)-1 ; -1 to account for NCL indexing starting<br>
x75 = round(.75*dimt6,3)-1 ; at 0 <br>
qsort(dPfldavgrh)<br>
print(dPfldavgrh) ; Print the values<br>
<br>
end do<br>
<br>
iarr=new((/7,5/),"float",-999.) ; fill with minimum, 25th absile, median, 75th absile, maximum of each timeseries<br>
iarr(0,:) = (/min(dPfldavg),dPfldavg(x25),dim_avg(dPfldavg),dPfldavg(x75),max(dPfldavg)/)
<br>
iarr(1,:) = (/min(dPfldavgShift),dPfldavgShift(x25),dim_avg(dPfldavgShift),dPfldavgShift(x75),max(dPfldavgShift)/)<br>
iarr(2,:) = (/min(dPfldavgtwrh),dPfldavgtwrh(x25),dim_avg(dPfldavgtwrh),dPfldavgtwrh(x75),max(dPfldavgtwrh)/)
<br>
iarr(3,:) = (/min(dPfldavgcross),dPfldavgcross(x25),dim_avg(dPfldavgcross),dPfldavgcross(x75),max(dPfldavgcross)/)<br>
iarr(4,:) = (/min(dPfldavgT),dPfldavgT(x25),dim_avg(dPfldavgT),dPfldavgT(x75),max(dPfldavgT)/)<br>
iarr(5,:) = (/min(dPfldavgweak),dPfldavgweak(x25),dim_avg(dPfldavgweak),dPfldavgweak(x75),max(dPfldavgweak)/)
<br>
iarr(6,:) = (/min(dPfldavgrh),dPfldavgrh(x25),dim_avg(dPfldavgrh),dPfldavgrh(x75),max(dPfldavgrh)/)<br>
<br>
print(iarr)<br>
<br>
x = (/-9.,-6.,-3.,0.,3.,6.,9./)<br>
<br>
<br>
;**********************************************<br>
; create plot<br>
;**********************************************<br>
wks = gsn_open_wks("X11","NCL_JOC_Boxplot_Land_Drying_abs_20_OND")<br>
;**********************************************<br>
; resources for plot background<br>
;**********************************************<br>
res = True ; plot mods desired<br>
res@tmXBLabels = (/"~F8~D~F21~P","~F8~D~F21~P~B~Shift","~F8~D~F21~P~B~twrh","~F8~D~F21~P~B~Cross","~F8~D~F21~P~B~T","~F8~D~F21~P~B~Weak","~F8~D~F21~P~B~RH"/) ; labels for each box ;<br>
res@tiMainString = "OND Box Plot over SA Drying (Land-Only) (10-40E 0S-30S)"<br>
res@trYMinF = -0.6 ; set minimum Y-axis value<br>
res@trYMaxF = 0.6 ; set maximum Y-axis value<br>
res@tiXAxisString = "Precipitation Components" <br>
res@tiYAxisString = "Absolute Change (mm/day/K)"<br>
res@tmXBLabelFontHeightF = 0.015 <br>
<br>
;**********************************************<br>
; resources for polylines that draws the boxes<br>
;********************************************** <br>
llres = True <br>
llres@gsLineThicknessF = 2.5 ; line thickness <br>
;**********************************************<br>
; resources that control color and width of boxes<br>
;********************************************** <br>
opti = True <br>
opti@boxWidth = 2.0 ; Width of box (x units)<br>
;colors = (/"red","blue","green","yellow","magenta","orange","cyan" /)<br>
opti@boxColors = (/"red","blue","green","yellow","purple","orange","cyan" /) ; Color of box(es);<br>
;gsn_define_colormap(wks, colors) <br>
;***********************************************<br>
plot = boxplot(wks,x,iarr,opti,res,llres) ; All 3 options used...<br>
draw(wks) ; box plot does not call<br>
frame(wks) ; these for you<br>
<br>
;end do<br>
<br>
<br>
end<br>
<br>
<div style="font-family: Times New Roman; color: #000000; font-size: 16px">
<hr tabindex="-1">
<div style="direction: ltr;" id="divRpF244188"><font face="Tahoma" color="#000000" size="2"><b>From:</b> Melissa Lazenby<br>
<b>Sent:</b> 21 March 2017 12:41<br>
<b>To:</b> ncl-talk@ucar.edu<br>
<b>Subject:</b> Masking, area averaging and box whisker plots<br>
</font><br>
</div>
<div></div>
<div>
<div style="direction:ltr; font-family:Tahoma; color:#000000; font-size:10pt">
<div style="direction:ltr; font-family:Tahoma; color:#000000; font-size:10pt">Dear NCL Users<span class="Apple-tab-span" style="white-space:pre">
</span>
<div><br>
</div>
<div>I am currently trying to create box whisker plots of drying regions over land regions.</div>
<div><br>
</div>
<div>I have tried creating the masks and field averages in CDO however there seems to be a bug in the CDO code I am running so I am trying to do all the masking and area averaging in an NCL script.</div>
<div><br>
</div>
<div>Essentially what I would like to do is:</div>
<div><br>
</div>
<div>1. Mask out the land regions (DONE)<br>
</div>
<div>2. Mask out the wetting regions (i.e. only have the drying regions) (DONE)<br>
</div>
<div>3. Area average the change in precipitation over a particular domain over the southern African continent (NOT SURE?? area average or dim_avg?)<br>
</div>
<div>4. Create a box whisker plot of the change in precipitation over that region for 20 different CMIP5 models to determine the robustness and spread.</div>
<div><br>
</div>
<div>Below is my code in NCL attempting to do 1-4<br>
<br>
I have managed to do steps 1 and 2 without problems, however I am struggling to complete step 3 i.e. create an field mean using the area average function.
<br>
Essentially I would like 20 output values from my 20 models of the area averages and currently I am only obtaining 1 value not 20.<br>
<br>
;*********************************************<br>
; boxplot.ncl<br>
;<br>
; ;*********************************************<br>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"<br>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"<br>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"<br>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"<br>
;*********************************************<br>
begin<br>
;**********************************************<br>
; Add data in text format and sort in ascending order<br>
;**********************************************<br>
<br>
<br>
model = (/"bcc-csm1-1-m","BNU-ESM","CanESM2","CCSM4","CESM1-BGC","CESM1-CAM5","CSIRO-Mk3-6-0","FIO-ESM","GFDL-CM3","GFDL-ESM2G","GFDL-ESM2M","GISS-E2-H","HadGEM2-CC","HadGEM2-ES","IPSL-CM5A-LR","IPSL-CM5A-MR","MIROC5","MRI-CGCM3","NorESM1-M","NorESM1-ME"/)<br>
<br>
<br>
do gg = 0,dimsizes(model)-1<br>
;dP<br>
in1=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadprfutureOND/Delta_pr_"+model(gg)+".nc","r")
<br>
<br>
;dPShift<br>
in2=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPSpat_"+model(gg)+"_safrica_PSpat.nc","r")<br>
;dPT<br>
in3=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPT_"+model(gg)+"_safrica_PT.nc","r")<br>
;dPWeak<br>
in4=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPDiv_"+model(gg)+"_safrica_PDiv.nc","r")<br>
;dPTWRH<br>
in5=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPDiv+PT+PRH_"+model(gg)+".nc","r")<br>
;dPCross<br>
in6=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPNL_"+model(gg)+"_safrica_PNL.nc","r")<br>
;dPRH<br>
in7=addfile("/research/geog/data2/DATA/mlazenby/Chadwick_Analysis/ChadDeltaPComponentsOND/DeltaPRH_"+model(gg)+"_safrica_PRH.nc","r")<br>
<br>
<br>
lat = in1->lat <br>
lon = in1->lon <br>
time = in1->time<br>
<br>
if ((gg.eq.0) .or. (gg.eq.9) .or. (gg.eq.10) .or. (gg.eq.11)) then<br>
dP = dble2flt(in1->pr(:,:,:))<br>
dPShift = in2->huss(:,:,:)<br>
dPtwrh = in5->huss(:,:,:)<br>
dPcross = in6->huss(:,:,:)<br>
dPT = in3->es2(:,:,:)<br>
dPweak = in4->huss(:,:,:)<br>
dPrh = in7->huss(:,:,:)<br>
<br>
else<br>
<br>
dP = in1->pr(:,0,:,:)<br>
dPShift = in2->huss(:,0,:,:)<br>
dPtwrh = in5->huss(:,0,:,:)<br>
dPcross = in6->huss(:,0,:,:)<br>
dPT = in3->es2(:,0,:,:)<br>
dPweak = in4->huss(:,0,:,:)<br>
dPrh = in7->huss(:,0,:,:)<br>
<br>
end if<br>
<br>
;Land mask<br>
<br>
b=addfile("$NCARG_ROOT/lib/ncarg/data/cdf/landsea.nc","r")<br>
lsm_dP = landsea_mask(b->LSMASK,dP&lat,dP&lon) ; read in land sea mask<br>
dP = mask(dP,lsm_dP.eq.0,False) ;mask out all ocean points)<br>
<br>
;Drymask<br>
<br>
dP_dry = new ((/dimsizes(time),dimsizes(lat),dimsizes(lon)/), typeof(dP))<br>
dP_dry!0="time"<br>
dP_dry&time=time<br>
dP_dry!1="lat"<br>
dP_dry&lat=lat<br>
dP_dry!2="lon"<br>
dP_dry&lon=lon<br>
<br>
dP_dry = mask(dP,dP.gt.0, False) ; mask out all the wetting regions i.e. any values above 0 for dP<br>
<br>
<br>
printVarSummary(dP)<br>
printVarSummary(dP_dry)<br>
<br>
dP_fldmean = wgt_areaave(dP_dry((gg),{-30:0},{10:40}),1.0, 1.0, 0)<br>
;dP_fldmean = dim_avg(dP_dry(0,{-30:0},{10:40}))<br>
<br>
print(dP_dry)<br>
printVarSummary(dP_fldmean)<br>
print(dP_fldmean)<br>
<br>
end do <br>
<br>
end <br>
</div>
<div><br>
</div>
<div>Any advice regarding this code would very much be appreciated. </div>
<div><br>
Kindest Regards<br>
Melissa<br>
</div>
<div><br>
</div>
</div>
</div>
</div>
</div>
</div>
</body>
</html>