一些常用的grads命令与ncl命令对应关系

游子

本文以grads命令为主线，逐条对应写下去。grads文件分为数据文件、数据描述文件和批处理文件。grads界面在windows系统、linux系统和mac系统中大致相同。而ncl的开发基本上是基于linux系统。那你要反驳了，ncl不是也可以装在windows系统中吗？那我要说你naive了，因为cgwin是模拟linux环境的。
1.grads的数据文件与ncl的什么文件对应？
grads只支持按照一定顺序存储的二进制数据文件，后缀名以.grd或者.dat或者.grb结束。这种数据必须以时间为最外层，然后每个变量按照向量形式存储，每个变量由外向内的存储顺序是高度（或等压面）==》纬度==》经度。这必须注意，否则画图容易出现一堆一堆的乱线条，这就说明你的数据没有按照grads的要求存储。
再看看ncl，可以说ncl支持绝大多数各种数据的读写，包括netcdf，hdf,以及二进制数据甚至ASCII码(如果说是十进制数据或许你会更熟悉)，前两种数据一般都有头文件，不能用C语言或者Fortran读取，都需要插件才可以读取，ncl可以直接读取，matlab中也可以直接读取NETCDF格式(.nc)的数据。grads中可以读取按照说明存取的NC数据，这种数据必须又正确的时间说明，也就是说时间必须是真实的，有些模式模拟出来的数据grads的sdfopen命令是打不开的，因为一般模式都是nonleap run，都是平年，没有设定闰年，造成了时间说明不真实，grads就会报错。那么grads'如何使用NC数据呢？所以建议使用ncl转换数据，将NC数据，hdf数据或者十进制数据转换成grd数据，供grads使用。
这样说明是在是太空洞了，那么下面我举个例子吧。
eg1) 使用grads将netcdf数据转换成grd数据
'reinit'
var.1='air';var.2='hgt';var.3='uwnd'
var.4='vwnd';var.5='omega';var.6='shum'
j=6
while (j<=6)
'set fwrite /disk3/users/Rao_Jian/ERA-Interim-daily/entropy/'var.j'.daily.1979-2010.grd'
'set gxout fwrite'
i=1979
while (i<=2010)
'sdfopen /disk3/users/lbq/ERA-Interim-daily/pressure/'var.j'.interim.'i'.nc'
tt=1
if(i=1980|i=1984|i=1988|i=1992|i=1996|i=2000|i=2004|i=2008)
while (tt<=366)
'set t 'tt
zz=1
while (zz<=37)
'set x 1 240'
'set y 1 121'
'set z 'zz
'd 'var.j''
zz=zz+1
endwhile
tt=tt+1
endwhile
else
while (tt<=365)
'set t 'tt
zz=1

zz=1
while (zz<=37)
'set x 1 240'
'set y 1 121'
'set z 'zz
'd 'var.j''
zz=zz+1
endwhile
tt=tt+1
endwhile
endif
i=i+1
name='/disk3/users/lbq/ERA-Interim-daily/pressure/'var.j'.interim.'i'.nc'
'close 1'
endwhile
'disable fwrite'
j=j+1
endwhile

eg2)使用ncl将数据输出成二进制数据，grads可以使用(截取部分)
patho           = "/disk3/users/Rao_Jian/Hadley/"
system("rm "+patho+"ev_ts.grd")
system("rm "+patho+"ev_ts2.grd")
system("rm "+patho+"ev_ts.txt")
system("rm "+patho+"ev_ts2.txt")
do nt=0,719
fbindirwrite(patho+"ev_ts.grd",ev_ts(:,nt));写成二进制数据，供grads使用
end do
fbindirwrite(patho+"ev_ts2.grd",ev_ts(time|:,evn|:))
asciiwrite(patho+"ev_ts.txt",ev_ts(time|:,evn|:))；写成十进制数据，可以贴到EXCEL中使用
asciiwrite(patho+"ev_ts2.txt",ev_ts)
;print(ev_ts(0,:))

此外ncl中还有其他的读写函数，如
fbinrecread,fbinrecwrite,fbinwrite,fbinread
ncl读取netcdf3/4、hdf4、grib-1/2也是通过文件操作的
eg3)ncl支持直接读取格式文件的读法
fi              = addfile(pathi+"HadISST_sst.nc","r")
sst0          = fi->sst(960:1679,:,:)     ; load 50 year data during from 1950 to 2009

注意：这类文件的后缀名一般为.nc /.hdf/ .grb /.hdfeos/.ccm
2.grads中的描述文件与ncl中的什么对应
描述文件（.ctl)是一个纯文本文件，我们有数据grads还是不能出图，需要一个描述文件来指定他的存储数据个数，维度（时间长度、层数和经纬度信息）。ncl是没有这种文件的，那么ncl在绘图时是通过什么信息来获得维度信息的呢？
其实，在ncl里又一个概念叫做coordinate,其实就是维度信息。在ctl中，是通过xdef、ydef、zdef、tdef来声明数据的经纬度信息的。也就是说，在一个ctl文件中声明的所有数据，其相同维度的信息必须相同(维数可以不同）。那么我们还是举个例子好了
eg4.1)ctl描述文件
dset  /disk3/users/Rao_Jian/ERA-DATA/entropy/Var3D.grd
undef -999
title monthly zonally mean uwnd,vwnd,temperature in THETA-PVLAT coordinate
options yrev
xdef 1   linear 0.000 2.500
ydef 200 linear 0.000 0.75
zdef 16  levels 270 280 290 300 315 330 350 370 400 450 500 550 600 650 700 750
tdef 389 linear jan1979 1mo
vars 3
ulat2 16 -99 means zonal mean u wind in theta-pv coordinate
vlat2 16 -99 means zonal mean v wind in theta-pv coordinate
tlat2 16 -99 means zonal mean temperature in theta-pv coordinate
endvars
;

上面的例子中，声明了三个变量，每个变量都是16层，说明都是空间三维的（可以不是的，此时16改为0），黄色的信息是维度信息，也就是说这三个变量ulat2、vlat2、tlat2相同的经度、维度、层数、时间。你要说了，我的数据哪那么巧都是一样的维度，要是不一样怎么办？
对不起，必须是相同维度的信息完全一样的数据才可以放到一个描述文件中，否则要分开在不同的文件中描述，这也就是为什么在一些使用者的批处理脚本中(.gs)中会打开多个文件的原因之一，因为变量的信息可能不同。
eg4.2)ncl中维度声明形式举例
;下面再来看看ncl中是如何声明维度信息的？
dudP0!0="time"
dudP0!1="lat"
dudP0!2="lon"
dudP0!3="level"
dudP0&time=u&time
dudP0&level=u&levelist
dudP0&lat=u&latitude
dudP0&lon=u&longitude
copy_VarCoords(dudP0,dvdP0)
copy_VarCoords(dudP0,dPTdP0)

在上面的例子中，！是用来给坐标维度命名的，&是用来引用某一维度数值的，copy_VarCoords是ncl的贡献函数库中的子程序，将一个变量的坐标信息原封不动的传给另一个变量。
必须注意的是，绘制地图底图时，还要涉及到属性的概念，用@来引用或者声明。所以在给变量声明经纬度之名时，还要声明经纬度之属性，保证ncl能找到地图坐标，保证等值线或矢量与地图坐标对应，不发生错位。
eg4.3)ncl中声明坐标经纬度属性
pvlvl&lev  = lvl                       ; isentropic levels
  pvlvl&lat  = u&latitude
  pvlvl&lon  = u&longitude
  pvlvl&lat@units  = "degrees_north"
  pvlvl&lon@units  = "degrees_east"
  pvlvl&lev@units  = "K"
  pvlvl&time@units = "months since 1900-01-01"

注意： pvlvl&lat@units  = "degrees_north"  @是接在pvlvl&lat后面的，也就是说，degrees_north属性是变量pvlvl的坐标lat的属性
3.grads中的批处理命令文件与ncl中的什么对应
这个不需要解释了吧，gs文件自然就相当与ncl脚本文件了，好好理解一下吧。说的更具体一点，是相当于ncl中graphics部分，因为ncl不仅具有绘图功能，还具有强大的计算功能，尤其是巨大的函数库。ncl的缺点是循环比较慢，优点就是有许多优化函数，不需要使用者自己编写子程序（如果你愿意，也可以自己编写，将脚本存在够得着的地方，用的时候load)。（再发一贴）

下一篇：http://bbs.06climate.com/forum.php?mod=viewthread&tid=10578 （管理员注）

mofangbao

随缘 发表于 2012-10-18 11:35

                               
登录/注册后可看大图

还带减贡献的啊！！！@mofangbao

哈哈哈，rp问题啊，我这个回复也许还能加一个呢

mofangbao

。。。。。。。。
整理了多少了啊。。

小蜜蜂

期待中

随缘

好期待啊！！！

随缘

还带减贡献的啊！！！@mofangbao

随缘

mofangbao 发表于 2012-10-18 12:09

                               
登录/注册后可看大图

哈哈哈，rp问题啊，我这个回复也许还能加一个呢

减了，哈哈，你RP也不行啊。。。哈哈

lsd1984

标题好吸引人啊！

天目神眉

标题党麽

游子

今天先写一点