[ncl-talk] problem with function (sum)
Xia Chu
xchu1 at uwyo.edu
Fri Apr 15 16:19:04 MDT 2016
Dear helpers,
I'm trying to do summation of a 3D array using function (sum). When the
data is 60x512x512 in dimension it works fine. However, with a larger
data set of 60x1536x1536, the summation came to 4.29497e+09 which is
definitely wrong. I had to use dim_sum_n to sum two dimension first and
then another sum. I also find summation value varies if I dim_sum_n
different dimensions first. This problem only stands out with variable
QVAPOR (not QCLOUD or other mixing ratios) which have values at every
grids.
The code is attached. I'm working on Yellowstone supercomputer. The
script has my working directory. You should be able to run the script on
Yellowstone. I'm using ncl/6.3.0.
Would you recommend a good way to do summation? I'm trying to calculate
total mass.
I would appreciate any of your help! Thanks in advance.
Best,
Xia
--
Xia Chu
PhD Candidate
Graduate Research Assistant
Department of Atmospheric Science
University of Wyoming
Dept. 3038, 1000 E. University Ave.
Laramie, Wyoming 82071
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160415/9de4d5c1/attachment.html
-------------- next part --------------
;TB is the total mass using function sum
;TB2 is the total mass using demensional summation first
;I'm expection Tb and TB2 to be the same values or at least close. But it is not the case...
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
case = "33.3m_1536"
var = "QVAPOR"
dir = "/glade/scratch/xchu/ideal_case/BLT/"+case+"/joint_wrf/aux_result/"
fileTB = systemfunc("ls "+dir+"auxhist*")
tsize = dimsizes(fileTB)
Firsttime = True
do tt = 0,2;tsize-1
a=addfile(fileTB(tt),"r")
dx = a at DX
dy = a at DY
times = wrf_user_list_times(a)
;get dimension
mdims1 = getfilevardimsizes(a,"QVAPOR")
;process for dz
phb1 = wrf_user_getvar(a, "PHB", 0)
ph1 = wrf_user_getvar(a, "PH", 0)
hgt_m1=(phb1+ph1)/9.81
delete(phb1)
delete(ph1)
dzTB=new((/mdims1(1),mdims1(2),mdims1(3)/),"float")
do k=0,mdims1(1)-1
dzTB(k,:,:)=(hgt_m1(k+1,:,:)-hgt_m1(k,:,:))
end do
delete(hgt_m1)
if (Firsttime .eq. True)then
mdims = getfilevardimsizes(a,"QVAPOR")
TB = new((/tsize/),float)
TB2 = new((/tsize/),float)
Firsttime = False
;calculate mass
varTB=wrf_user_getvar(a,var,0)/wrf_user_getvar(a,"ALT",0)*dzTB*dx*dy
;dimension summation
tempTB = dim_sum_n(varTB,(/1,2/))
;sum
TB(tt) =sum(tempTB)
TB2(tt) =sum(varTB)
else
varTB=wrf_user_getvar(a,var,0)/wrf_user_getvar(a,"ALT",0)*dzTB*dx*dy
tempTB = dim_sum_n(varTB,(/1,2/))
TB(tt) =sum(tempTB)
TB2(tt) =sum(varTB)
end if
end do
print(TB(0:2))
print(TB2(0:2))
More information about the ncl-talk
mailing list