气象程序设计及绘图期末复习课件.ppt

1、期期 末末 复复 习习气象程序设计及绘图 学习目标通过学习，了解FORTRAN语言的特点，基本成份及使用方法，具有阅读程序、编制程序、调试程序的基本能力。在此基础上，根据气象数据处理与分析的特点，切实地掌握用GrADS软件对数据进行图形化显示的基本方法，初步具备编程绘图的能力。气象程序设计及绘图专业基础课程序设计FORTRAN语言GrADS绘 图气象数据资料分析 基本程序结构program mainimplicit noneintegerrealcharactercomplexlogicalFORTRAN程序设计主程序开始不采用I-N规则变量（包括数组）申明及初始化 基本程序结构program

2、 mainimplicit noneinteger:real:character:complex:logical:FORTRAN程序设计主程序开始不采用I-N规则变量（包括数组）申明及初始化注意：character类型的使用方法赋值语句data 变量列表/初值列表/read*,open(文件号，file=路径及文件名)read(文件号，格式）。顺序结构选择结构循环结构获取数据程序算法单分支IF(条件条件)THENIF块块END IF双分支IF（条件）（条件）THENIF块块ELSEELSE块块END IF多分支IF（条件（条件1）THEN块块1ELSE IF（条件（条件2）THEN块块2ELS

3、E IF（条件（条件3）THEN块块3 .ELSE IF（条件（条件n）THEN块块nELSE块块n+1END IF逻辑if语句IF（条件）（条件）语句语句SELECT CASE（选择表达式）（选择表达式）CASE（控制表达式（控制表达式1）块块1CASE（控制表达式（控制表达式2）块块2.CASE（控制表达式（控制表达式n）块块nCASE DEFAULT默认块默认块END SELECT可以相互嵌套DO 循环变量循环变量=E1，E2，E3 循环体循环体END DOR=MAX(INT(E2-E1+E3)/E3)，0)DO WHILE(逻辑表达式逻辑表达式)循环体循环体 END DO可以相互嵌套p

4、rint*，输出列表write(*,*)输出列表open(文件号，file=路径及文件名)write(文件号，格式）。输出数据注意：带格式的输入输出方式二进制数据：open(文件号，file=路径及文件名，form=binary)read（文件号，格式）。write(文件号，格式）。End program main主程序结束一一个个主主程程序序 (p pr ro og gr ra am m)程程序序若若干干子子程程序序程序入口程序入口需要被调用需要被调用子程序包括：子程序包括：函数子程序函数子程序(function)子例行程序子例行程序(subroutine)子程序又分为：子程序又分为：外部子

5、程序外部子程序和内部子程序和内部子程序外部函数子程序外部函数子程序 由一个实现由一个实现某种特定功能某种特定功能的子程序组成，调用程的子程序组成，调用程序单元调用它得到序单元调用它得到一个函数值一个函数值定义形式：定义形式：类型说明类型说明 FUNCTION 函数名函数名(虚参虚参1,虚参虚参2,)说明语句说明语句 执行语句执行语句 。END FUNCTION 函数名函数名表示函数子程序定义开始表示函数子程序定义开始表示函数子程序定义结束表示函数子程序定义结束外部子例行程序外部子例行程序 子例行程序不仅可求一个值，还可求多个值或子例行程序不仅可求一个值，还可求多个值或不求值而执行某种操作，因此

6、具有更广泛的用途不求值而执行某种操作，因此具有更广泛的用途定义形式：定义形式：SUBROUTINE 子例行程序名子例行程序名(虚参虚参1,虚参虚参2,)说明语句说明语句 执行语句执行语句 END SUBROUTINE 子例行程序名子例行程序名表示子例行程序定义开始表示子例行程序定义开始表示子例行程序定义结束表示子例行程序定义结束调试FORTRAN程序 语法错误（syntax error）。如拼写错或标点符号错，在编译时被编译器检测出来；运行错误（run-time error）。当程序执行时企图做非法运算操作，则发生运行错误（例如，除以0）。这些错误导致程序执行时异常中断。逻辑错误（logica

7、l error）。程序能够正常编译和运行，但结果是错误的。例题：利用1951-2010年1月蒙古高压强度指数（data/h-p.dat）(备注：需要备注：需要标准化以后使用标准化以后使用)和1951-2010年1月中国160站气温资料(data/t1601.dat)（备注：实际存放的资料是备注：实际存放的资料是61年的，只读年的，只读60即可，即可，需要标准化以后使用需要标准化以后使用），根据3.5节计算1月蒙古高压强度与中国160站气温的相关关系。要求以“*.dat”和“*.grd”两种格式保存1951-2010年1月蒙古高压强度与我国气温的同期相关数据。12211122111211()()

8、11().(y)()()111()(y)()11(y)ntttx ynnttttnttnntttttntntttxxyynrxxynnxxyynxxynnyynyn已经标准化的蒙古高压指数integer,parameter:n=60,start=1951!character*10 char1integer i!integer(4)station,nstationreal a(n),b(n),c(160,n),coopen(1,file=d:aatest1h-p.dat,form=formatted)open(2,file=d:aatest1t1601.txt,form=formatted)op

9、en(3,file=d:aatest1mh-t-1.dat,form=formatted)open(4,file=d:aatest1mh-t-1.grd,form=binary)read(1,*)(a(i),i=1,n)close(1)read(2,*)(c(i,j),i=1,160),j=1,n)close(2)do i=1,160do j=1,60b(j)=c(i,j)enddocall COR2(n,a,b,co)write(3,(f12.5)cowrite(4)co !,colev!returns co=0.9548129;colev=0.1952enddoclose(3)close(

10、4)end主程序变量申明获取数据循环嵌套调用子程序输出数据主程序结束注释注释!求相关系数rsubroutine COR2(n,a,b,co)integer nreal coreal a(n),b(n)!标准化a(i),b(i)suma=sum(a(1:n)/nsumb=sum(b(1:n)/n ab=0a2=0b2=0do i=1,nab=ab+(a(i)-suma)*(b(i)-sumb)a2=a2+(a(i)-suma)*2b2=b2+(b(i)-sumb)*2enddo co=ab/sqrt(a2*b2)print*,coreturnend子程序开始子程序变量申明子程序算法子程序结束准备

11、数据准备数据文件文件(包括包括前期数据前期数据处理处理)编写及检编写及检查数据描查数据描述文件述文件*.ctl设定图设定图样进行样进行绘图绘图(可可编写脚编写脚本文件本文件*.gs)保保存存结结果果GrADS绘图流程在“写字板”或者“记事本”中编写 现有现有ASCII码（十进制数据格式）数据资料码（十进制数据格式）数据资料文件文件u850.dat和和v850.dat，其空间范围：，其空间范围：60-150E，0-40N；层次：；层次：u、v为为850Pa；时段；时段：1982.1-1985.12；分辨率：；分辨率：2.5*2.5；数据排放；数据排放顺序满足顺序满足GrADS要求。要求。要求编写

12、出将这要求编写出将这2个文件转换成个文件转换成1个二进制（个二进制（binary）文件的）文件的Fortran程序。程序。例题6-6GrADS要求数据如何排放？v850.dat数据内容：时间说明 PROGRAM PROCESSING IMPLICIT none integer,parameter:m0=37,n0=17,tim=48 real v(m0,n0,tim),u(m0,n0,tim)integer i,j,k open(1,file=d:cvfv850.dat)open(2,file=d:cvfu850.dat)do k=1,tim read(1,1000)read(1,2000)(

13、v(i,j,k),i=1,m0),j=1,n0)read(2,1000)read(2,2000)(u(i,j,k),i=1,m0),j=1,n0)end do1000 format(2i7)2000 format(37F6.2)读时间说明准备数据文件：利用fortran将十进制文件转换为二进制文件十进制数据循环顺序：经度、纬度、层次、时间 write(*,2000)(v(i,j,1),i=1,m0),j=1,n0)open(3,file=d:cvfwind850.grd,form=binary)do k=1,tim write(3)(v(i,j,k),i=1,m0),j=1,n0)write(

14、3)(u(i,j,k),i=1,m0),j=1,n0)end doend二进制数据针对针对wind850.grd二进制数据文件，如果给出如二进制数据文件，如果给出如下数据说明文件，即可用下数据说明文件，即可用GrADS绘图绘图dset d:cvfwind850.grdundef-9.99E+33title NCEP/NCAR REANALYSIS PROJECTxdef 37 linear 60.000 2.500ydef 17 linear 0.000 2.500zdef 1 levels 850tdef 48 linear JAN1982 1movars 2v 1 99 v wind(m/

15、s)u 1 99 u wind(m/s)endvarsCTL文件文件编写ctl文件编写gs文件reinitopen d:/cvf/wind850.ctlenable print d:/cvf/wind850-1.gmfset grads offset grid offset t 13set parea 1 10 1.5 7.5set xlopts 1 4 0.16set ylopts 1 4 0.16set xlint 10set ylint 10define uave1=ave(u,t=1,t=48,12)define vave1=ave(v,t=1,t=48,12)set gxout ve

16、ctorset arrscl 0.5 10.0d uave1;vave1printdisable print；编写gs文件reinitopen d:/cvf/wind850.ctl 打开二进制数据文件的文件enable print d:/cvf/wind850-1.gmf 打开存图文件set grads off 去掉画图角标set grid off 去掉网格线set t 13 设置固定于某个时次，以便绘制二维图形set parea 1 10 1.5 7.5 设置绘图区域set xlopts 1 4 0.16 设置x轴标记set ylopts 1 4 0.16 设置y轴标记set xlint 10 设置x轴标记间隔单位set ylint 10 设置x轴标记间隔单位define uave1=ave(u,t=1,t=48,12)定义新变量，计算1月份平均u风场define vave1=ave(v,t=1,t=48,12)定义新变量，计算1月份平均v风场set gxout vector 设置出图类型set arrscl 0.5 10.0 设置风场矢量单位d uave1;vave1 绘制风场print 将图形窗口的图像存入图片文件disable print 结束存图；绘图程序结束图形窗口图形窗口命令窗口命令窗口观察图片文件