Does the script below for computing 500hPa geopotential heights from CCSM4 data
look ok? The resulting values are slightly higher than I expected, by maybe 100 m,
but perhaps I am wrong in my concern.
Thanks, Michael
nlat=192
nlon=288
mns=(/"JAN","FEB","MAR","APR","MAY","JUN","JUL","AUG","SEP","OCT","NOV","DEC"/)
a=addfile("/data/Share/regcm_inputdata/CCSM4/historical/ccsm4.shumJAN.1971.nc","r")
lat=a->lat
lon=a->lon
lev=a->lev(::-1)
coefa=a->a(::-1)
coefb=a->b(::-1)
p0=a->p0
b=addfile("/data/Share/regcm_inputdata/CCSM4/historical/orog_fx_CCSM4_historical_r0i0p0.nc","r")
topo=b->orog ; height above sea level (m)
data=new((/32,12,nlat,nlon/),float) ; 1973-2004
do iy=0,31
do im=0,11
a=addfile("/data/Share/regcm_inputdata/CCSM4/historical/ccsm4.shum"+mns(im)+"."+(iy+1973)+".nc","r")
q=a->hus(:,::-1,:,:) ; kg/kg spec hum
w=q/(1.-q) ; convert spec hum to mixing ratio
q=(/w/)
time=a->time
b=addfile("/data/Share/regcm_inputdata/CCSM4/historical/ccsm4.pres"+mns(im)+"."+(iy+1973)+".nc","r")
ps=b->ps
ps=ps/100.
copy_VarCoords(b->ps,ps)
c=addfile("/data/Share/regcm_inputdata/CCSM4/historical/ccsm4.air"+mns(im)+"."+(iy+1973)+".nc","r")
t=c->ta(:,::-1,:,:)
copy_VarCoords(t,q)
zsfc=new((/dimsizes(time),nlat,nlon/),float)
do it=0,dimsizes(time)-1
zsfc(it,:,:)=topo
end do
copy_VarCoords(ps,zsfc)
t=t*(1.+0.61*q) ; virtual T
copy_VarCoords(q,t)
p=new((/dimsizes(time),dimsizes(coefa),nlat,nlon/),float)
do i=0,dimsizes(coefa)-1
p(:,i,:,:)=doubletofloat(coefa(i))*p0+doubletofloat(coefb(i))*ps
end do
copy_VarCoords(t,p) ; suppose to be hPa
z=p(time|:,lat|:,lon|:,lev|::-1)
z=hydro(p(time|:,lat|:,lon|:,lev|::-1),t(time|:,lat|:,lon|:,lev|::-1),zsfc)
z=z(:,:,:,::-1)
p1=p(time|:,lat|:,lon|:,lev|:)
h500=int2p(p1,z,500.,2)
h500_less=h500(:,:,:,0)
copy_VarCoords(ps,h500_less)
v=h500_less(lat|:,lon|:,time|:)
data(iy,im,:,:)=dim_avg_Wrap(v)
delete(time)
delete(a)
delete(b)
delete(c)
delete(q)
delete(ps)
delete(t)
delete(zsfc)
delete(p)
delete(z)
delete(p1)
delete(h500)
delete(h500_less)
delete(v)
delete(w)
end do
end do
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Thu Mar 7 21:54:46 2013
This archive was generated by hypermail 2.1.8 : Wed Mar 13 2013 - 14:19:38 MDT