第2章 数据处理

第2章数据处理在使用GrADS绘图前，必须使所用数据文件满足GrADS的数据格式要求。

2.1各类数据格式的特点一般二进制无格式直接或者顺序记录的数据格式：可以是格点数据或站点数据，可以用Fortran或者C语言读写。需要数据描述文件.ctl

GRIB码数据格式：是自定义的数据格式，压缩率高，占用空间小。例如NCEP提供的数据资料。NETCDF（.nc）等通用数据格式：自定义数据格式，精确性好，便于传输。

文件中自带描述文件。2.2数据文件的转换

1.转换方法

转换文件的数据存放格式，一般是通过PowerStation或VisualFortran等软件使用Fortran或者C语言来编程转换。

2.数据存放形式

二进制数据排放顺序从内循环到外循环依次是：

x（经度）：从西到东

y（纬度）：从南到北

z（高度层数）：从低层到高层

vars（各种物理变量）

t（时次）x(lon)y(lat)z(lev)vars(不同变量)time

任何一维可省略。

12345………..x

y………54321Z3.举例（格点资料）

现有ASCII码（十进制存放格式）数据资料文件u.dat、v.dat和sst.dat，其空间范围：60-150°E，0-40°N；层次：u、v为850、200hPa；时段：1982.1-1985.12；分辨率：2.5*2.5。要求编写出将这三个文件转换成二进制无格式直接存取（Grads格式.grd）文件的Fortran程序。

编程时数组大小：X方向：150-60=90，90/2.5+1=37Y方向：40-0=40，40/2.5+1=17Z方向：u、v为850、200hPa

，Z为2sst为海表温度，只有一层，Z为1T时次：月资料，4年，共48个月Fortran程序编写如下：

!定义一个X,Y,Z方向的格点数以及总时次ntparameter(nx=37,ny=17,nz=2,nt=48)!定义数组

dimensionu(nx,ny,nz,nt),v(nx,ny,nz,nt),sst(nx,ny,nt)

!打开原始数据文件,如果数据文件和Fortran程序在一个文件夹下，打开数据文件时可以不写路径。

open(1,file='u.dat')open(2,file='v.dat')open(3,file='sst.dat')!打开目标文件，recl为记录的长度，对于直接存取文件必须指定记录长度。

open(12,file='mhy.grd',form='unformatted',#access='direct',recl=nx*ny*4)

!把数据文件读入do100it=1,ntdo50iz=1,nzread(1,*)((u(i,j,iz,it),i=1,nx),j=1,ny)50continuedo51iz=1,nzread(2,*)((v(i,j,iz,it),i=1,nx),j=1,ny)51continueread(3,*)((sst(i,j,it),i=1,nx),j=1,ny)100continue!将数据资料写入新文件中，irec重新进行记录累加irec=0do200it=1,ntdo30iz=1,nzirec=irec+1write(12,rec=irec)((u(i,j,iz,it),i=1,nx),j=1,ny)30Continuedo31iz=1,nzirec=irec+1write(12,rec=irec)((v(i,j,iz,it),i=1,nx),j=1,ny)31continueirec=irec+1write(12,rec=irec)((sst(i,j,it),i=1,nx),j=1,ny)200continueendirec=0do200it=1,ntdo30iz=1,nzirec=irec+1write(12,rec=irec)((u(i,j,iz,it),i=1,nx),j=1,ny)30Continuedo31iz=1,nzirec=irec+1write(12,rec=irec)((v(i,j,iz,it),i=1,nx),j=1,ny)31continueirec=irec+1write(12,rec=irec)((sst(i,j,it),i=1,nx),j=1,ny)200continue如果不考虑采用无格式直接方式存取文件的话，可以使用下面的编程形式，较为简单。parameter(nx=37,ny=17,nz=2,nt=48)dimensionu(nx,ny,nz,nt),v(nx,ny,nz,nt),sst(nx,ny,nt)open(1,file='u.dat')open(2,file='v.dat')open(3,file='sst.dat')!打开目标文件open(12,file='mhy.grd',form='binary')

!打开原始数据文件,如果数据文件和Fortran程序在一个文件!夹下，打开数据文件时可以不写路径do100it=1,ntdo50iz=1,nzread(1,*)((u(i,j,iz,it),i=1,nx),j=1,ny)50continuedo51iz=1,nzread(2,*)((v(i,j,iz,it),i=1,nx),j=1,ny)51continueread(3,*)((sst(i,j,it),i=1,nx),j=1,ny)100continue!把数据文件读入do200it=1,ntdo30iz=1,nzwrite(12)((u(i,j,iz,it),i=1,nx),j=1,ny)30Continuedo31iz=1,nzwrite(12)((v(i,j,iz,it),i=1,nx),j=1,ny)31continuewrite(12)((sst(i,j,it),i=1,nx),j=1,ny)200continueend!将数据资料写入新文件中注意事项：二进制写入的数据必须是实型的(real)若数据为整型，写入时要转换成实型的real(iu(i,j,k,it))2.2数据描述文件（.ctl）

GrADS绘图不能直接使用“数据文件”，而是通过“数据描述文件”间接使用“数据文件”。1.数据描述文件的概念

数据描述文件是原始数据文件的描述文件。用以描述原始数据集的基本信息，包括数据集文件名、数据类型、数据结构、变量描述等等。2.数据描述文件的形式

OPTIONS<keywords>定义数据格式选项，keywords有：<yrev><zrev><sequential><byteswapped><template><big-endian><little-endian>

缺省时为direct直接记录格式（form='binary'）<yrev>自北向南的数据<zrev>自高层到低层的数据<sequential>顺序无格式存贮的数据<byteswapped>二进制数据的位存放顺序取反序

<big-endian>在SUN,SGI等机器上生成的数据<little-endian>在PC,DEC等机器上生成的数据服务器个人电脑big-endian、little-endian用于自动改变二进制位存放顺序<template>

多个时间序列原始数据文件用一个数据描述文件统一地描述这些原始数据时采用的选项，这些数据文件的文件名形式由dset定义的形式命名文件名，提示所含数据的时次。例如：一个逐日的数据集每天数据放到一个文件中，每个文件名形式为：

1may92.dat，2may92.dat，3may92.datdset%d1%mc%y2.datoptionstemplatetdef30linear0Z1may19921dy

定义时间范围和增量例如：正确的替换为：%y2两位数年%y44位数年

%m11或2位数月%m22位数月（用0补齐1位数）

%mc3字符月份缩写%d11或2位数天%d22位数天

%h11或2位数小时%h22位小时

设置X方向格点与经度的对应关系XDEFnumber<LINEARStartincrement>或XDEFnumber<LEVELSvalue-list>number（>=1）给定格点数，整形；

LINEAR指明是等间隔分布格点，

Start起点坐标，负数表示西经；

increment网格间距。

LEVELS参数指明是不等间隔分布格点，列出具体每个格点的坐标值（以空格分开）。

设置Y方向格点与纬度的对应关系YDEFnumber

<LINEARStartincrement>或YDEFnumber

<LEVELSvalue-list>或YDEFnumber<mapping

Start>mapping表示映射方式，有以下几种：GAUSR15--高斯R15纬度(Gaussian)GAUSR20--高斯R20纬度GAUSR30--高斯R30纬度GAUSR40--高斯R40纬度在linear中，Start起点坐标负数表示南纬；在高斯映射中，Start起点坐标表示高斯纬度。方向:自西向东,自南向北(永远不变！)单位:个纬度/经度正负值:以东为正/以北为正xdef144linear02.5ydef73linear-902.5方向维数定义格点数坐标值线性增加起始坐标坐标间隔YDEF20GAUSR4015

表示共有20个Y方向网格点，起始点为高斯R40网格下的高斯纬度15（即64.100S）,实际这20个网格点对应的纬度值为：-64.10,-62.34,-60.58,-58.83,-57.07,-55.32,-53.56,-51.80,-50.05,-48.29,-46.54,-44.78,-43.02,-41.27,-39.51,-37.76,-36.00,-34.24,-32.49,-30.73设置垂直网格点与气压面的对应关系ZDEFnumber<LINEARStartincrement>或ZDEFnumber<LEVELSvalue-list>

zdef7levels1000850700500300200100

列出具体坐标值方向:自下而上单位:层hPa/其他

TDEFnumberLINEARstart-timeincrement

设置网格值与时间的映射关系。时间表示法：hh:mmZddmmmyyyyhourminutedaymonyeartdef12linear00:00z21jan20001mo

时间步长只有linear关键字起始时间必须按照上述格式写年份用4位数写21日00:00时按情况可省略2000年1月(不可省略)3mn

分钟6hr

小时1dy

天1mo

月1yr

年例如：12Z1JAN199014:20Z22JAN1987JUN1960

Z：世界时12Z代表北京时间20点VARSnumber

变量描述开始，同时给出变量个数number。

abrevlevsunitsdescription………..abrev，变量名称，为1到8个字符组成的该变量的缩写名，画图或处理数据将用到，以字母开头levs，数字，代表变量层数，0表示只有一层units，单位，为GRIB预留，给99description，对变量的文字描述，最多40个字符用ENDVARS表示数据描述文件结束。

dtype用于说明数据格式的类型当数据为格点binary格式时省略dtype的可选参数

–grib–station站点数据

–hdfsds–netcdfv1.9版本新增•ctl书写注意事项♥♥♥每行顶格写;虽然大小写没关系,推荐在一个关键字中都用大写或都用小写;如写成dset,而不写成DsEt；每行记录的各项以空格分开；每行记录不超过80个字符；注释行请用*开始,并顶格写;不能出现在变量列表中，即：在vars~endvars之间不要写注释行;不要在GrADS中使用制表(tab)键!4.数据描述文件的建立

数据描述文件为纯文本格式，可用一般的编辑器产生（如记事本，写字板等）。描述文件编辑完成后，保存时为.ctl。

注意事项：如果电脑中常用文件类型隐藏时，请修改该设置。打开“我的电脑”，打开在窗口第一行第五个“工具”，选“文件夹选项”，点“查看”，将滚动条下拉至“隐藏已知文件类型的扩展名”，将其前面对号去掉，点“确定”即可。鼠标右键新建文本文档，打开，写入ctl内容，然后保存为file.ctl。练习：为mhy.grd文件（刚才进行数据转换得到的新文件）编写数据描述文件mhy.ctl。dsetmhy.grdundef-9.99E+33titleNCEP/NCARREANALYSISPROJECTxdef37linear60.02.5ydef17linear0.02.5zdef2levels850200tdef48linearJAN19821movars3u299uwind(m/s)v299vwind(m/s)sst099TSFCDATAendvars打开ctlopenh:\grads\sst.win.ctl5、netCDF格式数据–ECMWF/NCEP再分析资料等–常用后缀.nc–自带数据内容说明(自描述)一般分辨率为2.5°×2.5°或者1°×1°

标准NC格式,GrADS可直接处理，fortranga->sdfopenfilename.nc↵特别提醒：虽然GrADS支持这种格式的数据，但是涉及到要用程序对数据进行处理时，Fortran程序并不支持这种格式的数据。如果只是使用这类数据进行画图和GrADS的命令处理，则可以直接操作。如果要用于程序使用的话，应该利用命令“fwrite”将所用数据提取为“.grd”的二进制数据文件，以后再用于编程。6、GRIB格式数据T213,NECP/NCAR等常用后缀.grb需.ctl

和.idx

文件一起使用用grib2ctl.exe生成数据说明(描述)文件.ctl用gribmap–itest.ctl生成映射文件.idx1.8版本需另外下载放到win32e下！用grib2ctl.pl生成.ctl的方法：1、从“开始”菜单的“运行”选项里输入“cmd”进入DOS状态，然后cdc:\PCGrADS\win32e。2、在该路径下输入命令：

grib2ctl路径/model.grb

>

路径/model.ctl

这样就可以在指定路径下生成数据文件model.grb的数据描述文件model.ctl.注意：命令中间的“>”符号不能缺少，生成的数据描述文件的路径和文件名都是用户自己指定的。Grib2ctl.exe所在文件夹下用gribmap–itest.ctl生成映射文件.idx：继续在C:\Program\PCGrADS\win32e下输入命令：

gribmap-i路径/model.ctl

这样在原来存放数据文件model.grd的目录下自动生成数据文件的映射文件“model.idx”。在GrADS的命令窗口ga->提示符后输入命令：

open路径/model.ctl

即可使用model.grd数据.GRIB格式数据ctl写法特别提醒：虽然GrADS支持grib格式的数据，但是涉及到要用程序对数据进行处理时，Fortran程序并不支持这种格式的数据。如果只是使用这类数据进行画图和GrADS的命令处理，则可以直接操作。如果要用于程序使用的话，应该利用命令“fwrite”将所用数据提取为“.grd”的二进制数据文件，以后再用于编程。7、站点数据文本格式，点开可见数据值50136台站06-07年逐日要素值，有说明文件

年份月份日期台站号纬度经度海拔高度日平均气压(0.1hpa)日平均气温(0.1℃)……..1951-2005年7月160站月降水量160站站点号和经纬度文件：

idlonlat50353 51.72 126.6550632 48.77 121.9250527 49.22 119.7550434 50.45 121.7050557 49.17 125.2350745 47.38 123.9250756 47.43 126.9750788 47.23 131.9850873 46.49 130.1750978 45.28 130.9550953 45.75 126.7754094 44.57 129.60 dimensionrain(160,55) dimensionhp(160),ae(160,55) dimensionalat(160),alon(160)integeridx(160) character*8id(160)！读入站点号和经纬度资料

open(21,file='id160.dat')doist=1,160read(21,'(i5,2f7.2)')idx(ist),alat(ist),alon(ist)enddoclose(21)！读入降水资料

open(41,file='r1607')doiy=1,55read(41,*)(rain(ii,iy),ii=1,160)enddoclose(41)!求55年平均值

doidx=1,160hp(idx)=0.0doit=1,55hp(idx)=hp(idx)+rain(idx,it)enddohp(idx)=hp(idx)/55enddo!求距平值

doidx=1,160doit=1,55ae(idx,it)=rain(idx,it)-hp(idx)enddoenddo对数据进行处理！距平值写成grads可以识别的站点数据！将整数数组idx转化成字符型数组id，id是8个字符，这个是固定的

doi=1,160write(id,‘(i5)’)idx(i)

!或者id(i)=char(idx(i))enddo！打开新文件，存储站点数据open(41,file=‘arain7.tu.grd',form='binary')doit=1,55iflag=1！iflag”取1表示有地面变量，取0表示记录中没有地面观测变量。ilev=1！ilev取1表示一个时次记录的开始，取0表示一个时次记录的结束time=0.0！time为本报的时次，不是一个确定的时刻值，只是该时次的！相对值

doist=1,160write(41)id(ist),alat(ist),alon(ist),time,ilev,iflag,ae(ist,it)enddo！给出一个特别记录，表示这个时间组的记录结束

ilev=0！取0表示一个时次记录的结束

write(41)id(160),alat(160),alon(160),time,ilev,iflagenddoclose(41)enddsetd:\data\arain7.tu.grddtypestationstnmapd:\data\rain.mapundef-999.0titleraintdef55linearJUl19511yrvars1rain099anomalyrainfalldataENDVARS站点数据Ctl:站点映射文件名字stnmap用于将观测站点数据与图形显示时的站点图建立链接以加速显示。语法为：

ga->stnmap-i站点数据描述文件路径/**.ctl也可直接进入stnmap交互式操作：stnmap.exe位于pcgrads目录中win