wrf输出结果 时间序列

追风的阳光

下面是脚本：

oad "$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/csm/contributed.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRF_contributed.ncl"
begin
   in  = addfile("/home/zssapr/WRF3.5.1/20110710/data/wrf_data/wrfout_d02_2011-07-10_00:00:00","r")

   U =in->U
   ua = wrf_user_unstagger(U,U@stagger)
   V = in->V
   va = wrf_user_unstagger(V,V@stagger)
   lat = in->XLAT(0,:,0)
   lon = in->XLONG(0,0,:)

   Time   = WRF_Times2Udunits_c(in->Times, 0)  ; convert to "hours since"
   scale=1.e05
   dv = new ( dimsizes(ua), typeof(ua), ua@_FillValue )
   dv@long_name  = "divengence"
   dv@units      = "scaled"
   dv = uv2dv_cfd (ua,va,lat,lon, 2)

   vr = new ( dimsizes(ua), typeof(ua), ua@_FillValue )
   vr@long_name  = "vorticity"
   vr@units      = "scaled"
   vr = uv2vr_cfd (ua,va,lat,lon, 2) * scale

   svbit = new ( dimsizes(ua), typeof(ua), ua@_FillValue )
   svbit =dv/vr


   time=new(41,integer)
   bottom_top=new(27,integer)
   south_north=new(81,integer)
   west_east=new(81,integer)
   svbit!0="time"
   svbit!1="bottom_top"
   svbit!2="south_north"
   svbit!3="west_east"
   svbit@time=Time
   svbit@bottom_top=bottom_top
   svbit@south_north=south_north
   svbit@west_east=west_east
   znu     = in->ZNU(0,:)                 ; znu(bottom_top)
   Svbit       =svbit(:,:,6,7 )
   wks  = gsn_open_wks("ps","press_time")
   gsn_define_colormap(wks,"BlueWhiteOrangeRed")
   res = True
   res@cnFillOn = True
   plot = gsn_csm_contour(wks, Svbit(bottom_top|:,time|:),res)
end
NCL学习过程得到了许多同学的帮助，在此谢谢了，希望以后能帮助更多的人

longlivehj

   time=new(41,integer)
   bottom_top=new(27,integer)
   south_north=new(81,integer)
   west_east=new(81,integer)
坐标数组全部是missing，没有赋值。

   svbit@time=Time
   svbit@bottom_top=bottom_top
   svbit@south_north=south_north
   svbit@west_east=west_east

坐标变量的赋值要用“&”，“@”是属性赋值时采用的。

另外，我总觉得最后不需要自己去定义坐标信息，用copy_VarCoords应该可以解决。

追风的阳光

longlivehj 发表于 2014-5-4 19:36
time=new(41,integer)
   bottom_top=new(27,integer)
   south_north=new(81,integer)

谢谢你的热心帮助！

andrewsoong

追风的阳光 发表于 2014-5-4 19:53
谢谢你的热心帮助！

做好了放上来啊~~~~

追风的阳光

andrewsoong 发表于 2014-5-4 20:07
做好了放上来啊~~~~

明天放。。现在连不上校内网。。

追风的阳光

longlivehj 发表于 2014-5-4 19:36
time=new(41,integer)
   bottom_top=new(27,integer)
   south_north=new(81,integer)

你好，我今天做又发现了一个问题，上面脚本中的svbit是四维的（时间，高度，经度，纬度），我想取每一时刻高度上的一个区域内的svbit平均值，于是添加了下面这个脚本： Svbit = new ( dimsizes(FNM), typeof(ua) )
copy_VarCoords(ua, svbit )   
dimp   = dimsizes( ua )
ntime   = dimp(0)
nlev   = dimp(1)
Svbit   =  new ( (/ntime,nlev/), typeof(ua) )
do i=0,120,1
do  j=0,26,1
      Svbit(i,j)=dim_sum_n(svbit(i,j,1:5,1:5),0)
end do
end do

这个脚本我想把svbit转换成二维的（FNM是wrf文件中的一个变量（时间，高度），我就图方便把FNM的维度信息赋给Svbit）。
运行后提示说：
Dimension sizes on right hand side of assignment do not match dimension sizes of left hand side
即Svbit(i,j)=dim_sum_n(svbit(i,j,1:5,1:5),0)俩边维数不一样，请问应该怎样调整呢？谢谢！

longlivehj

追风的阳光 发表于 2014-5-5 15:34
你好，我今天做又发现了一个问题，上面脚本中的svbit是四维的（时间，高度，经度，纬度），我想取每一时 ...

根据你的描述，以下两行代码即可：
copy_VarCoords(ua, svbit)
Svbit = dim_avg_n_Wrap(svbit(:, :, 1 : 5, 1 : 5), (/2, 3/))

你程序出错的原因是：
svbit(i,j,1:5,1:5)是二维数组，dim_sum_n(svbit(i,j,1:5,1:5),0)后是一维数组（5个元素），而等号左边Svbit(i,j)是1个元素，这样是无法赋值的。

另外注意dim_sum_n是求和，不是求均值。

追风的阳光

longlivehj 发表于 2014-5-5 15:45
根据你的描述，以下两行代码即可：
copy_VarCoords(ua, svbit)
Svbit = dim_avg_n_Wrap(svbit(:, :, 1  ...

谢谢！明白了，/2,3/就是固定time，lev，取区域的平均。。但“dim_sum_n(svbit(i,j,1:5,1:5),0)后是一维数组（5个元素）”，为何是5个元素呢？

longlivehj

追风的阳光 发表于 2014-5-5 16:42
谢谢！明白了，/2,3/就是固定time，lev，取区域的平均。。但“dim_sum_n(svbit(i,j,1:5,1:5),0)后是一维 ...

svbit(i,j,1:5,1:5)是(/5, 5/)的数组，再对第一维求和，不就是5个元素么！

追风的阳光

longlivehj 发表于 2014-5-5 16:58
svbit(i,j,1:5,1:5)是(/5, 5/)的数组，再对第一维求和，不就是5个元素么！

谢谢！。。。。。。。