利用站点资料画出等值线图的GrADS示例

grads画图坐标设置⼀直听到有⼈抱怨，GrADS的坐标轴怎么那么固定，要设置个起始间隔还不⽀持时间轴，要在经纬度后⾯加个⼩圆圈的度数怎么就那么难，为啥不能四周都弄上坐标轴？好了，你的这些抱怨将会随着你看到这个帖⼦⽽消失，下⾯给出这个脚本的⽤法介绍：注意：在你display变量前需要先设置如下两个命令'set xlpos -20''set ylpos -20'这样能够屏蔽GrADS⾃带的坐标轴，否则会出现坐标轴重叠[code=gs]请先display变量，再运⾏该脚本该脚本主要实现了⾃定义的坐标轴显⽰，您可以免费使⽤该脚本该脚本的⽤法规则如下：#axis -param1 val1 -param2 val2...其中，-param表⽰需要⾃定义的参数类型，val表⽰该参数的具体设定值参数和参数值必须成对出现，例如#axis -type b -position o -sinterval 2上⾯的语句表⽰绘制类型是底部的x轴，刻度位置是朝外，每两个⼤的刻度之间显⽰两个⼩刻度线#注意：-type是必选参数所有参数如下所⽰：-type ：表⽰绘制类型，参数值有：b/l/t/r 分别表⽰底部，左侧，上部，右侧-label：表⽰是否绘制数值标签，参数值有：on/off 分别表⽰绘制和不绘制-position：表⽰刻度的位置，参数值有：o/i/m 分别表⽰外侧，内侧，中间-start：表⽰刻度从该处开始绘制，参数值为该轴范围内的数字当该轴为时间轴时，表⽰开始绘制的时次（1,2,3...）-end：表⽰刻度的结束位置，参数值范围同start，时间轴时表⽰结束的时刻(1,2,3...) -interval：表⽰刻度的间隔，隔多少值绘制⼀个刻度和标签例如，当x轴表⽰经度110-150E，如果interval为10，默认情况下会标出110,120,130,140,150不给定该参数时，默认情况下，会绘制5个刻度，以此确定默认间隔-sinterval：表⽰没两个⼤刻度之间要绘制⼏个⼩的刻度，不给定该参数则不绘制-size：表⽰刻度的长度，单位是英⼨，默认为0.1-color：表⽰坐标轴和刻度标签的颜⾊，默认为1-lfont：表⽰坐标轴标签的字体，默认为当前环境下设置的字体-lsize：表⽰坐标轴标签的字体⼤⼩，默认为0.12-lthick：表⽰坐标轴标签字体的粗细，默认为0.3-langle：表⽰坐标轴标签旋转的⾓度，默认为不旋转-suffix：表⽰需要在每个坐标轴标签后⾯添加的⽂字或者符号后缀，⽐如⽤来添加度数等-asuffix：表⽰当坐标轴为经纬度轴时，是否⾃动添加E/W/N/S/EQ这种标记，默认为添加，参数值为 on/off ，分别表⽰⾃动添加和关闭⾃动添加-hoffset：表⽰坐标轴标签在⽔平⽅向上的偏移量，单位为英⼨，正负均可-voffset：表⽰坐标轴标签在垂直⽅向上的偏移量，单位为英⼨，正负均可-tformat：当所绘制坐标轴为时间轴时产⽣作⽤，表⽰要显⽰的时间格式参数值为y m d h 这四个字母的任意组合，分别表⽰年⽉⽇时如：-tformat my 则会将时间轴标签显⽰为JAN1951这种格式-v：当绘制图形是1-D图（set gxout line等）的时候是必选参数，表⽰当前display的是哪个变量[/code]来具体的看两个例⼦，这两个例⼦只是为了描述脚本的⽤法，所以不⼀定美观。

第一步将文本格式的不规则站点数据转换成grads可识别的二进制格式。

1）首先将文本格式的站点数据按下列方式排列，保存在txt文件中（图1）。

图1其中，第一列是站号（5位字符），第二列是该站点的纬度，第三列是该站点的经度，第四列是变量（降水，温度等）。

2）接下来，用下列Fortran程序将txt文件转换成二进制格式。

图23）然后，为新生成的二进制站点数据station.grd创建站点描述文件station.ctl。

dset D:/data/station/station.grddtype stationstnmap D:/data/station/station.mapundef 999.9title daily tdtdef 1 linear Jan2002 1movars 1td 0 99 tdpoint dataendvars创建好后，进入grads安装目录下的win32文件夹，双击stnmap.exe，按照提示输入该ctl文件名（包括路径），回车。

成功后在D:/data/station下会生成相应的地图文件station.map（即station.ctl中stnmap行的设置）。

第二步：创建相应的格点文件。

1）运行Fortran程序生成相应的二进制格点文件。

program wgrid! to get gridsinteger lon,latparameter (lon=25,lat=27) ！具体格点数要根据实际绘图区域及格距计算得到real gtd(lon,lat)! redefine the grids numbers when necessaryopen(10,file='D:\data\station\wgrid.dat',form='binary')do i=1,londo j=1,latgtd(i,j)=0.0end doend dowrite(10)((gtd(i,j),i=1,lon),j=1,lat)end2）创建二进制格点文件wgrid.grd的描述文件grid.ctl，如下:dset D:\data\station\wgrid.dattitle gtdundef -9.99E33xdef 25 linear 60 0.2ydef 27 linear 20 0.2 *需要与上一步中的格点数一致zdef 1 linear 1 1tdef 1 linear jan2002 1movars 1gtd 0 99 grid tdendvars此处注意：station.ctl与grid.ctl中时间设置必须一致第三步：生成等值线图片。

以下技巧总结都是笔者从学习实践过程中总结出来的，基本的问题。

不求全面，希望对读者学习有用，如果有问题，敬请留言指正，以促进交流学习！1、软件综述：grads软件是一款绘图软件除了绘制图形，还可以提取数据，主要应用是在大气科学中，当然只要是数据处理成grads能够读取的数据文件就可以进行相关绘图。

软件版本问题，软件本身不是很大，我接触到1.8、1.9、2.0版本的，1.8版本的安装很多情况还要修改环境变量、1.9版本的不识别‘sdfopen’命令，最稳定的版本是2.0版本，所以笔者推荐学习者安装2.0版本，选择默认安装路径就可以。

2、文件类型简述：grads处理的是网格数据，可以处理的数据类型有：grd、grib、nc（海洋常用的数据），cdf（雷达卫星数据），其中nc、cdf数据都是自带描述文件，不需要ctl，grib数据要通过命令生成ctl、index数据才可以调用，常用的是grd数据，需要ctl。

3、数据文件转换：grads软件识别的数据是二进制无格式数据，文件类型是‘binary’，写入和生成时是不需要格式的如read(20) sst(i,j,iz,it),20为文件号，通常是十进制数据与grd数据间转换，这里给一个grd转换成txt数据的fortran程序：parameter(nx=56,ny=41,nz=1,nt=360)dimension sst(nx,ny,nz,nt)real sstopen(15,file='sst.grd',form='binary') ！固定的用form=‘binary’就是二进制数据open(16,file='sst.txt') ！新建txt文件do it=1,ntdo iz=1,nzread(15) ((sst(i,j,iz,it),i=1,nx),j=1,ny) ！read后只有文件号，数据是无格式的enddoenddodo it=1,ntdo iz=1,nzwrite(16,*) ((sst(i,j,iz,it),i=1,nx),j=1,ny) ！输出时是txt文件可直接看的数据，有格式输出，有*enddoenddoclose(15)close(16)end写程序时：注意格点数要与数据对应，如：上程序对应的数据是经度90~200，纬度-20~60，时间：1971.01~2000.12共360个月的海面温度数据，数据格点精度2*2 ，nx=(200-90)/2+1,ny=(60-(-20))/2+1,nt=360，nz=1，大气的数据要根据数据的层次确定几层。

以下技巧总结都是笔者从学习实践过程中总结出来的，基本的问题。

不求全面，希望对读者学习有用，如果有问题，敬请留言指正，以促进交流学习！1、软件综述：grads软件是一款绘图软件除了绘制图形，还可以提取数据，主要应用是在大气科学中，当然只要是数据处理成grads能够读取的数据文件就可以进行相关绘图。

软件版本问题，软件本身不是很大，我接触到1.8、1.9、2.0版本的，1.8版本的安装很多情况还要修改环境变量、1.9版本的不识别‘sdfopen’命令，最稳定的版本是2.0版本，所以笔者推荐学习者安装2.0版本，选择默认安装路径就可以。

2、文件类型简述：grads处理的是网格数据，可以处理的数据类型有：grd、grib、nc（海洋常用的数据），cdf（雷达卫星数据），其中nc、cdf数据都是自带描述文件，不需要ctl，grib数据要通过命令生成ctl、index数据才可以调用，常用的是grd数据，需要ctl。

3、数据文件转换：grads软件识别的数据是二进制无格式数据，文件类型是‘binary’，写入和生成时是不需要格式的如read(20) sst(i,j,iz,it),20为文件号，通常是十进制数据与grd数据间转换，这里给一个grd转换成txt数据的fortran程序：parameter(nx=56,ny=41,nz=1,nt=360)dimension sst(nx,ny,nz,nt)real sstopen(15,file='sst.grd',form='binary') ！固定的用form=‘binary’就是二进制数据open(16,file='sst.txt') ！新建txt文件do it=1,ntdo iz=1,nzread(15) ((sst(i,j,iz,it),i=1,nx),j=1,ny) ！read后只有文件号，数据是无格式的enddoenddodo it=1,ntdo iz=1,nzwrite(16,*) ((sst(i,j,iz,it),i=1,nx),j=1,ny) ！输出时是txt文件可直接看的数据，有格式输出，有*enddoenddoclose(15)close(16)end写程序时：注意格点数要与数据对应，如：上程序对应的数据是经度90~200，纬度-20~60，时间：1971.01~2000.12共360个月的海面温度数据，数据格点精度2*2 ，nx=(200-90)/2+1,ny=(60-(-20))/2+1,nt=360，nz=1，大气的数据要根据数据的层次确定几层。

GrADS绘图软件实用手册2002年1月目录第一章 GrADS绘图软件概述1.GrADS绘图软件简介2.Internet上的GrADS资源2.1GrADS在Internet上的主页2.2 windows环境下GrADS资源3.GrADS绘图软件的安装（windows环境）3.1在windows环境下安装GrADS软件包3.2X server 的安装第二章 GrADS绘图模板1.GrADS示例演示1.1 启动GrADS1.2 退出GrADS1.3 示例演示GrADS命令的使用2.GrADS绘图模板3.GrADS模板的高级应用3.1GrADS描述语言3.2GrADS高级模板的应用第三章 GrADS数据格式1.格点数据描述文件1.1 数据描述文件各项解释1.2 生成model.le.dat和model.le.ctl文件的程序代码片段2.站点数据的格式附录1.如何精确控制图形输出的尺寸—Landscape纸型2.台站资料的显示3.Linux环境下的安装第一章 GrADS绘图软件概述1GrADS绘图软件简介The Grid Analysis and Display System(GrADS) 是一套应用广泛、使用方便的科学数据绘图软件包。

其主要特点：z GrADS属于自由软件，可以从Internet上免费获得。

z可运行于各种Windows 和Unix工作平台。

z GrADS可用于4D数据的分析。

既经度、纬度、层（气压层、高度层等）和时间/xyzt 4维。

数据可以是格点化的数据或离散点数据。

GrADS特别适用于气象类数据的分析。

但也完全可以用于更广泛类型的数据分析。

z GrADS有多种显示方式：等值线、流线、矢量图、风矢量图、站点填图、折线图、直方图等多种两维图形。

z可处理多种数据格式的数据。

GRIB、NetCDF、HDF-SDS等通用数据格式和系统自定义的一种二进制数据格式。

z采用命令行输入的方式交互式地显示图形。

Grads站点资料的多层次绘图看了很多关于如何使用grads和fortran转换站点资料的帖子，大多数都集中在地面或者单层的，对于高空的多个层次很少有提到，而且大多数有些模糊，不好理解，我把自己做的一些整理传上来，供大家参考参考。

欢迎大家指出理解和程序的错误。

一，首先是数据的说明选用micaps高空观测资料850hpa,700hpa,600 hpa,500hpa，200hpa，100 hpa 六个层次作为例子。

为了方便，时间层次就选择1个时间。

micaps高空观测资料的说明如下：diamond 2 10年07月28日20时200百帕高空观测10 07 28 20 200 419因此定义mt=419,nz=4数据每行的排列分别是：站号，经度，纬度，高度……..具体详见micaps 的使用手册，这里不再赘述。

根据这些定义变量。

二，fortran读取站点资料注：1.根据micaps的高空资料定义相应的变量，由于850-100这6个层次的站点是相同的，站号，经纬度等只需定义成一维的，如lat(mt)，而lev(nz)表示有6个层次，也为一维；高度场等又有mt个站的资料又有nz个层次，因此定义为2维的。

2.由于micaps资料分别在不同的文件内，文件名重复，我修改为100_2820,200_2820,500_2820…….然后使用trim函数顺序读取资料。

3.为了更好的说明地面和高空资料的读取，将地形高度high(mt)设置为地面的量，nlev定义为7，即高空六层，地面一层4.以输出地面一个变量，高空三个变量为例。

先进行站点的循环，再进行高空层次的循环。

5.前面我们定义了lev(nz)，现在给lev(nz)赋值，使用select case.,分别赋值为850,700,600,500,200,100。

在使用是发现如果不进行赋值，后续工作也无法进行，会出现错误。

程序如下：program stationparameter (mt=419,nz=6) ! mt 总站数，nz总层次character stid(mt)*8integer flag,nlevreal lat(mt),lon(mt),tim,high(mt),number(mt),lev(nz) realheight(mt,nz),temp(mt,nz),ttd(mt,nz),ud(mt,nz),us(mt,nz),u(mt,nz),v(m t,nz),td(mt,nz)！顺序读取资料character stationfile(6)*8data stationfile/'850','700','600','500','200','100'/do k=1,nzopen(k,file=''//trim(stationfile(k))//'_2820.000')read(k,*)read(k,*)！跳过micaps开头两行的说明文字do i=1,mtread(k,*)stid(i),lon(i),lat(i),high(i),number(i),height(i,k),temp(i,k ),ttd(i,k),ud(i,k),us(i,k)end doend do !结束k的循环tim=0.0flag=1nlev=7 ！高空+地面open(40,file='20.grd',form='binary')do i=1,mt ！站点循环write(40)stid(i),lat(i),lon(i),tim,nlev,flag,high(i)do k=1,nz ！高空循环！对lev(nz)进行赋值select case(k)case(1)lev(k)=850.0case(2)lev(k)=700.0case(3)lev(k)=600.0case(4)lev(k)=500.0case(5)lev(k)=200.0case defaultlev(k)=100.0end select！在lev(k)后将要输出的所有高空变量一次性输出write(40)lev(k),height(i,k),temp(i,k),ttd(i,k)end doend donlev=0write(40)stid(mt-1),lat(mt-1),lon(mt-1),tim,nlev,flag close(40)end三，ctl及map文件的生成注：看到相关帖子说在站点资料的ctl文件中也要设置zdef，我试了下，设置和不设置都是可以的Ctl文件如下：dset f:\20.grddtype stationstnmap f:\20.mapundef 9999.0title 20zdef 6 levels 850 700 600 500 200 100tdef 1 linear 28jul2010 1DYvars 4high 0 99 Terrain Heighth 6 99 heightt 6 99 temperaturettd 6 99 ttdendvarsCtl文件编写完成后，就可以在grads中生成map文件了，相关命令如下：！ stnmap -i f:/20.ctl注意，斜杠一定是“/”,方向不能反第三步结束后，站点资料的读取就结束了，可以通过grads打开Ctl文件看到相关站点的值了。

以下技巧总结都是笔者从学习实践过程中总结出来的，基本的问题。

不求全面，希望对读者学习有用，如果有问题，敬请留言指正，以促进交流学习！（笔者：阿木）1、软件综述：grads软件是一款绘图软件除了绘制图形，还可以提取数据，主要应用是在大气科学中，当然只要是数据处理成grads能够读取的数据文件就可以进行相关绘图。

软件版本问题，软件本身不是很大，我接触到1.8、1.9、2.0版本的，1.8版本的安装很多情况还要修改环境变量、1.9版本的不识别‘sdfopen’命令，最稳定的版本是2.0版本，所以笔者推荐学习者安装2.0版本，选择默认安装路径就可以。

2、文件类型简述：grads处理的是网格数据，可以处理的数据类型有：grd、grib、nc（海洋常用的数据），cdf（雷达卫星数据），其中nc、cdf数据都是自带描述文件，不需要ctl，grib数据要通过命令生成ctl、index数据才可以调用，常用的是grd数据，需要ctl。

3、数据文件转换：grads软件识别的数据是二进制无格式数据，文件类型是‘binary’，写入和生成时是不需要格式的如read(20) sst(i,j,iz,it),20为文件号，通常是十进制数据与grd数据间转换，这里给一个grd转换成txt数据的fortran程序：parameter(nx=56,ny=41,nz=1,nt=360)dimension sst(nx,ny,nz,nt)real sstopen(15,file='sst.grd',form='binary') ！固定的用form=‘binary’就是二进制数据open(16,file='sst.txt') ！新建txt文件do it=1,ntdo iz=1,nzread(15) ((sst(i,j,iz,it),i=1,nx),j=1,ny) ！read后只有文件号，数据是无格式的enddoenddodo it=1,ntdo iz=1,nzwrite(16,*) ((sst(i,j,iz,it),i=1,nx),j=1,ny) ！输出时是txt文件可直接看的数据，有格式输出，有*enddoenddoclose(15)close(16)end写程序时：注意格点数要与数据对应，如：上程序对应的数据是经度90~200，纬度-20~60，时间：1971.01~2000.12共360个月的海面温度数据，数据格点精度2*2 ，nx=(200-90)/2+1,ny=(60-(-20))/2+1,nt=360，nz=1，大气的数据要根据数据的层次确定几层。