1、气象统计实习报告实习报告书课程名称:气象统计方法课程实践姓名:学号: 班级: 级气科 班*实习一 求500hPa高度场气候场、距平场和均方差场 实习时间:第9周周三1、2节1. 资料介绍有一500hPa高度场资料,文件名h500.dat,范围:60150E,040N.时段:1982.11985.12共48个月。水平分辨率:2.5*2.5,格点数:37*17。2要求编fortran程序,求500hPa高度场的(1)气候场;(2)距平场;(3)均方差场。并能用Grads做出图形,实习报告中气候场、距平场、均方差场任意给出两张图,图注要清楚,即要注明是哪个时间的图形,并做简单分析。注:h500.Fo
2、r给出了如何用fortran读取ASCII码资料h500.dat. FORTRANprogram sx1implicit noneinteger nx,ny,mo,yrparameter(nx=37,ny=17,mo=12,yr=4)real var(nx,ny,mo,yr)real at(nx,ny,mo),xd(nx,ny,mo,yr),sx(nx,ny,mo)integer i,j,m,t,it,iy,irecopen(5,file=d:studyformshixioneh500.dat) do iy=1,4 do m=1,12 read(5,1000) read(5,3000)(var
3、(i,j,m,iy),i=1,nx),j=1,ny)enddoenddo close(5)!计算气候场atdo t=1,12 do j=1,ny do i=1,nx at(i,j,t)=0 do it=1,4 at(i,j,t)=at(i,j,t)+var(i,j,t,it) enddo at(i,j,t)=at(i,j,t)/4 enddo enddoenddo!求距平场xddo t=1,12 do j=1,ny do i=1,nx xd(i,j,t,1)=0 do it=1,4 xd(i,j,t,it)=var(i,j,t,it)-at(i,j,t) enddo enddo enddoen
4、ddo !求均方差场sxdo t=1,12 do j=1,ny do i=1,nx sx(i,j,t)=0 do it=1,4 sx(i,j,t)=sx(i,j,t)+(var(i,j,t,it)-at(i,j,t)*2 enddo sx(i,j,t)=sqrt(sx(i,j,t)/4) enddo enddoenddo!写入气候场open(10,file=d:studyformshixioneat.grd,form=unformatted,access=direct,recl=nx*ny)irec=0do t=1,12irec=irec+1write(10,rec=irec)(at(i,j,
5、t),i=1,nx),j=1,ny)enddoclose(10)!写入距平场open(11,file=d:studyformshixionexd.grd,form=unformatted,access=direct,recl=nx*ny)irec=0do it=1,4 do t=1,12irec=irec+1write(11,rec=irec)(xd(i,j,t,it),i=1,nx),j=1,ny) enddoenddoclose(11)!写入均方差场open(12,file=d:studyformshixionesx.grd,form=unformatted,access=direct,r
6、ecl=nx*ny)irec=0do t=1,12irec=irec+1write(12,rec=irec)(sx(i,j,t),i=1,nx),j=1,ny)enddoclose(12)1000 format(2i7)2000 format(37F6.2)3000 format(37f8.1)4000 format(37f7.2)end program sx1 运行结果: Grads文件气候场reinitenable print d:studyformshixioneat.gmfopen d:studyformshixioneat.ctlset grid offset grads offse
7、t lat 0 40set lon 60 150set lev 500mon=1while(mon=12)set t mond hdraw title 1982yearmonmonthprintcmon=mon+1endwhiledisable print;距平场reinitenable print d:studyformshixionesx.gmfopen d:studyformshixionesx.ctlset grid offset grads offset lat 0 40set lon 60 150set lev 500year=1982while(year=1984)mon=1wh
8、ile(mon=12)set t mond hdraw title 500hPa yearyearmonmonth anomalyprintcmon=mon+1endwhileyear=year+1endwhiledisable print;均方差reinitenable print d:studyformshixionesx.gmfopen d:studyformshixionesx.ctlset grid offset grads offset lat 0 40set lon 60 150set lev 500mon=1while(monmax_y) THEN max_y=rxy_ty(t
9、y) k=ty END IFEND DOPRINT (全年平均气温绝对值最大自相关系数rxy_ty=,f7.4,/,滞后时间长度 k=,I2),rxy_ty(k),kk=0DO tw=1,N/2 DO i=1,N-tw rtw(i)=(w(i)-avr_w)/sw)*(w(i+tw)-avr_w)/sw) rxy_tw(tw)=rxy_tw(tw)+rtw(i) END DO rxy_tw(tw)=rxy_tw(tw)/(N-tw) rxy_tw(tw)=ABS(rxy_tw(tw) IF(rxy_tw(tw)max_w) THEN max_w=rxy_tw(tw) k=tw END IFEN
10、D DOPRINT (冬季平均气温绝对值最大自相关系数rxy_tw=,f7.4,/,滞后时间长度 k=,I2),rxy_tw(k),kk=0END 运行成果:*实习四 求给定数据的一元线性回归方程实习时间:第12周周三1、2节利用下表数据,以环流指标为预报因子,气温为预报量,计算气温和环流指标之间的一元线性回归方程,并对回归方程进行检验。年份气温T环流指标19510.93219521.22519532.22019542.4261955-0.52719562.5241957-1.128195802419596.21519602.71619613.2241962-1.13019632.522196
11、41.23019651.82419660.63319672.42619682.52019691.2321970-0.835答案:F=20.18F=4.41,回归方程显著program sx4implicit noneinteger ireal: x(20),y(20), f1=0,f2=0,b,avx,avy,b0,sx=0,sy=0,sxy=0,f,rxydata(x(i),i=1,20)/32,25,20,26,27,24,28,24,15,16,24,30,22,30,24,33,26,20,32,35/data(y(i),i=1,20)/0.9,1.2,2.2,2.4,-0.5,2.5
12、,-1.1,0,6.2,2.7,3.2,-1.1,2.5,1.2,1.8,0.6,2.4,2.5,1.2,-0.8/do i=1,20f1=f1+x(i)*y(i)end doavx=sum(x)/20 !求均值avy=sum(y)/20do i=1,20f2=f2+x(i)*2end dob=(20*f1-sum(x)*sum(y)/(20*f2-sum(x)*2) !求bb0=avy-b*avxprint *,b=,bprint *,b0=,b0print*,y=,b0,b,xdo i=1,20 sx=sx+(x(i)-avx)*2sy=sy+(y(i)-avy)*2sxy=sxy+(x(
13、i)-avx)*(y(i)-avy)end do sx=sqrt(sx/20) sy=sqrt(sy/20)sxy=sxy/20rxy=sx*b/syf=rxy*2*18/(1-rxy*2)if(f4.41) print *,F=,f,F=4.41,回归方程显著if(f=4.41) print *,F=,f,F=4.41,回归方程不显著End 运行成果:*实习七 计算给定数据的11年滑动平均和累积距平 实习时间:第15周周三1、2节利用数据ma.dat,编写11点滑动平均的程序,ma.for给出了阅读资料的fortran程序。数据在文件夹中单独给出。要求:实习报告中附出程序,并给出原数据和滑动
14、后数据的图形(1张图)和累积距平数据图形(1张图)program sx7 implicit noneinteger n,ih,nyearparameter(n=85,ih=11,nyear=1922)real:x(n),x1(n-ih+1)=0,x2(n)=0,averx=0integer i,jopen(2,file=d:studyformshiyansevenma.dat)open(3,file=d:studyformshiyansevenmaresult.dat)open(4,file=d:studyformshiyansevenaccresult.dat)read(2,*)(x(i),
15、i=1,n)do j=1,n-ih+1 !moving average do i=1,ih x1(j)=x1(j)+x(i+j-1) end do x1(j)=x1(j)/ihend dodo i=1,n-ih+1write(3,*) x1(i)end dodo i=1,n !average of x averx=averx+x(i)end doaverx=averx/ndo i=1,n !accelarate do j=1,i x2(i)=x2(i)+(x(j)-averx) end doend do do i=1,n write(4,*) x2(i)end do End 图形显示:*实习八
16、 对给定的海温数据进行EOF分析 实习时间:第16周周三1、2节给出海表温度距平数据资料sstpx.grd,以及相应的数据描述文件sstpx.ctl,对其进行EOF分析,资料的时空范围可以根据sstpx.ctl获知。 数据在文件夹中单独给出,距平或者标准化距平处理后再进行EOF。Zhunsst.for给出了如何读取资料,Ssteof.for为对距平或者标准化距平处理后的资料进行EOF分析。要求:实习报告中给出第一特征向量及其时间系数,并分析其时空特征。 program main parameter(mo=12,my=43,nx=18,ny=12,mt=516) dimension avf(mo
17、,nx,ny),df(mo,nx,ny) dimension sst(nx,ny,mt),sst3(nx,ny,mt) open(1,file=d:studyformshixieightsstpx.grd,form=unformatted,access=direct,recl=nx*ny) irec=1 do it=1,mt read(1,rec=irec)(sst(i,j,it),i=1,nx),j=1,ny) irec=irec+1 end do! average do j=1,ny do i=1,nx do k=1,mo do it=k,mt,12 avf(k,i,j)=avf(k,i,
18、j)+sst(i,j,it) end do avf(k,i,j)=avf(k,i,j)/my end do end do end do! variance do j=1,ny do i=1,nx do k=1,mo do it=k,mt,12 df(k,i,j)=df(k,i,j)+(sst(i,j,it)-avf(k,i,j)*2 end do df(k,i,j)=sqrt(df(k,i,j)/my) end do end do end do ! standardizing do j=1,ny do i=1,nx do k=1,mo do it=k,mt,12 if(sst(i,j,it)=
19、-999.0)then sst3(i,j,it)=-999.0 else sst3(i,j,it)=(sst(i,j,it)-avf(k,i,j)/df(k,i,j) end if end do end do end doend do! output file open(2,file=d:studyformshixieightstandard.grd,form=unformatted,access=direct,recl=nx*ny) irec=1 do it=1,mt write(2,rec=irec)(sst3(i,j,it),i=1,nx),j=1,ny) irec=irec+1 end
20、 do close(2) close(1)end 运行程序后:得到zhunsst.grd文件,即图中standard下面为老师给出程序,为运行方便将部分简略。 PROGRAM EOFPW PARAMETER(M=516,N=216,MNH=216,KS=-1,KV=10,KVT=2) PARAMETER(ff=-999.0,nx=18,ny=12) DIMENSION F(N,M),A(MNH,MNH),S(MNH,MNH),ER(mnh,4),DF(N),V(MNH),AVF(N),evf(N,KVT),tCF(M,KVT),data(Nx,ny),nf(N) open(11,file=d:
21、studyformshixieightzhunsst.grd,form=unformatted,access=direct,recl=nx*ny) do 132 it=1,m read(11,rec=it)(data(i,j),i=1,nx),j=1,ny) do 132 jj=1,ny do 132 ii=1,nx kkkk=nx*(jj-1)+ii f(kkkk,it)=data(ii,jj)132 continue close(11) CALL Test1(n,m,ff,f,nf) write(*,*)ok2 CALL TRANSF(N,M,F,nf,AVF,DF,KS) write(*
22、,*)ok3 CALL FORMA(N,M,MNH,F,A) write(*,*)ok4 CALL JCB(MNH,A,S,0.00001) write(*,*)ok5 CALL ARRANG(KV,MNH,A,ER,S) write(*,*)ok6 CALL TCOEFF(KVT,KV,N,M,MNH,S,F,V,evf,tcf,ER) write(*,*)ok7 call test3(N,ff,nf,evf,kvt) write(*,*)ok8 open(21,file=d:studyformshixieightevf.grd,form=unformatted,access=direct,
23、recl=nx*ny) irec=0 do 668 kk=1,kvt irec=irec+1668 write(21,rec=irec)(evf(nx*(j-1)+i,kk),i=1,nx),j=1,ny) close(21) open(21,file=d:studyformshixieighttcf.grd,form=unformatted,access=direct,recl=kvt) irec=0 do 345 it=1,m irec=irec+1345 write(21,rec=irec) (tcf(it,iik),iik=1,kvt) close(21)106 format(10f8
24、.4) open(21,file=d:studyformshixieightdats.dat) write(21,106)(er(iiii,3),iiii=1,kv) close(21) STOP END SUBROUTINE Test1(n1,m,ff,f,nf) DIMENSION F(N1,M) DIMENSION nF(N1) do i=1,n1 nf(i)=0.0 enddo do i=1,n1 do k=1,m if(f(i,k).eq.ff)then f(i,k)=0.0 nf(i)=nf(i)+1 endif enddo enddo return end SUBROUTINE TRANSF(n1,m,f,nf,avf,df,ks) DIMENSION F(N1,M),AVF(N1) DIMENSION DF(N1) DIMENSION nF(N1)
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1