Convert WRF out data to ARL dataΒΆ
ARL meteorological data format is specified using in HYSPLIT model. This is an example script for converting WRF out netCDF data to ARL data.
#---- Set data folder
datadir = 'E:/Temp/Yaqiang'
#---- Set output data file
outfn = os.path.join(datadir, 'test_wrfout1.arl')
#---- Read a netCDF data file
infn = os.path.join(datadir, 'wrfout_d01_1984-05-30_12_00_00-subset')
print infn
inf = addfile(infn)
print 'NetCDF file has been opened...'
#---- Set output ARL data file
arlf = addfile(outfn, 'c', dtype='arl')
#---- Set variable and level list
wvar2d = ['PSFC','HGT','T2','U10','V10','PBLH','RAINNC']
wvar3d = ['P','U','V','W','T','QVAPOR']
avar2d = ['PRSS','SHGT','T02M','U10M','V10M','PBLH','TPPA']
avar3d = ['PRES','UWND','VWND','WWND','TEMP','SPHU']
wv = inf['P']
nx = wv.dimlen(wv.ndim - 1)
ny = wv.dimlen(wv.ndim - 2)
levels = wv.dimvalue(wv.ndim - 3)
nz = len(levels)
arlf.setlevels(levels)
arlf.set2dvar(avar2d)
for l in levels:
arlf.set3dvar(avar3d)
#---- Constant for poisson's equation
cp = 1004.0 # J/kg/K; specific heat
grav = 9.81 # m/s**2; gravity
rdry = 287.04 # J/kg/K; gas constant
rovcp = rdry / cp # constant for poisson's equation
#---- Write ARL data file
arlf.setx(wv.dimvalue(wv.ndim - 1))
arlf.sety(wv.dimvalue(wv.ndim - 2))
tNum = inf.timenum()
fhour = 6
for t in range(0, tNum):
print 'Time index: ' + str(t)
atime = inf.gettime(t)
dhead = arlf.getdatahead(inf.proj, 'AWRF', 1, fhour)
#Pre-write index record without checksum - will be over-write latter
arlf.writeindexrec(atime, dhead)
#Checksum list
ksumlist = []
# Write 2d variables
ksums = []
for avname,wvname in zip(avar2d, wvar2d):
print avname
gdata = inf[wvname][t,:,:]
if avname == 'PRSS':
gdata = gdata * 0.01
elif avname == 'TPPA':
gdata = gdata * 0.001
ksum = arlf.writedatarec(atime, 0, avname, fhour, 99, gdata)
ksums.append(ksum)
ksumlist.append(ksums)
# Write 3d variables
for lidx in range(0, nz):
ksums = []
print lidx
pp = inf['P'][t,lidx,:,:]
pb = inf['PB'][t,lidx,:,:]
pres = pp + pb
uwnd = inf['U'][t,lidx,:,:]
vwnd = inf['V'][t,lidx,:,:]
temp = inf['T'][t,lidx,:,:]
#potential to ambient temperature
temp = (temp + 300.) * (pres / 100000.) ** rovcp
sphu = inf['QVAPOR'][t,lidx,:,:]
#convert mixing ratio to specific humidity
sphu = sphu / (1. + sphu)
wwnd = inf['W'][t,lidx+1,:,:]
#convert vertical velocity from m/s to hPa/s using omega = -W g rho
wwnd = -wwnd * grav * pres * 0.01 / (rdry * temp * (1.0 + 0.6077 * sphu))
pres = pres * 0.01
ksum = arlf.writedatarec(atime, lidx + 1, 'PRES', fhour, 99, pres)
ksums.append(ksum)
ksum = arlf.writedatarec(atime, lidx + 1, 'UWND', fhour, 99, uwnd[:,:nx])
ksums.append(ksum)
ksum = arlf.writedatarec(atime, lidx + 1, 'VWND', fhour, 99, vwnd[:ny,:])
ksums.append(ksum)
ksum = arlf.writedatarec(atime, lidx + 1, 'WWND', fhour, 99, wwnd)
ksums.append(ksum)
ksum = arlf.writedatarec(atime, lidx + 1, 'TEMP', fhour, 99, temp)
ksums.append(ksum)
ksum = arlf.writedatarec(atime, lidx + 1, 'SPHU', fhour, 99, sphu)
ksums.append(ksum)
ksumlist.append(ksums)
#Re-write index record with checksum
arlf.writeindexrec(atime, dhead, ksumlist)
fhour += 1
arlf.close()
print 'Finished!'