气象统计实习报告.docx

1、气象统计实习报告?实?习?报?告?书课程名称：气象统计方法课程实践姓名：学号： 班级： 级气科 班*实习一 求500hPa高度场气候场、距平场和均方差场 实习时间：第9周周三1、2节1. 资料介绍有一500hPa高度场资料，文件名，范围：60150E，040N.时段：共48个月。水平分辨率：*，格点数：37*17。2要求编fortran程序，求500hPa高度场的（1）气候场；（2）距平场；（3）均方差场。并能用Grads做出图形,实习报告中气候场、距平场、均方差场任意给出两张图，图注要清楚，即要注明是哪个时间的图形，并做简单分析。注：给出了如何用fortran读取ASCII码资料. FORT

2、RANprogram 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:studyformshixione) do iy=1,4 do m=1,12 read(5,1000) read(5,3000)(var(i,j,m,iy),i=1,nx),j=1,ny)enddoenddo close(5)!

3、计算气候场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 enddoenddo !求均方差场sxdo t=1,12 do j=1,ny do i=1,nx sx(i

4、,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:studyformshixione,form=unformatted,access=direct,recl=nx*ny)irec=0do t=1,12irec=irec+1write(10,rec=irec)(at(i,j,t),i=1,nx),j=1,ny)enddoclose(10)!写入距平场open(11,file=d

5、:studyformshixione,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:studyformshixione,form=unformatted,access=direct,recl=nx*ny)irec=0do t=1,12irec=irec+1write(12,rec=irec)(sx(i,j,t)

6、,i=1,nx),j=1,ny)enddoclose(12)1000 format(2i7)2000 format3000 format4000 formatend program sx1 运行结果： Grads文件气候场reinitenable print d:studyformshixioneopen d:studyformshixioneset grid offset grads offset lat 0 40set lon 60 150set lev 500mon=1while(mon=12)set t mond hdraw title 1982yearmonmonthprintcmo

7、n=mon+1endwhiledisable print;距平场reinitenable print d:studyformshixioneopen d:studyformshixioneset grid offset grads offset lat 0 40set lon 60 150set lev 500year=1982while(year=1984)mon=1while(mon=12)set t mond hdraw title 500hPa yearyearmonmonth anomalyprintcmon=mon+1endwhileyear=year+1endwhiledisab

8、le print;均方差reinitenable print d:studyformshixioneopen d:studyformshixioneset grid offset grads offset lat 0 40set lon 60 150set lev 500mon=1while(monmax_y) THEN max_y=rxy_ty(ty) k=ty END IFEND DOPRINT (全年平均气温绝对值最大自相关系数rxy_ty=,/,滞后时间长度 k=,I2),rxy_ty(k),kk=0DO tw=1,N/2 DO i=1,N-tw rtw(i)=(w(i)-avr_w)

9、/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 IFEND DOPRINT (冬季平均气温绝对值最大自相关系数rxy_tw=,/,滞后时间长度 k=,I2),rxy_tw(k),kk=0END 运行成果：*实习四 求给定数据的一元线性回归方程实习时间：第12周周三1、2节利用下表数据，以环流指标为预报因子，气温为

10、预报量，计算气温和环流指标之间的一元线性回归方程，并对回归方程进行检验。年份气温T环流指标1951321952251953201954261955271956241957281958024195915196016196124196230196322196430196524196633196726196820196932197035答案：F=F=，回归方程显着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,

11、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,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

12、)-avy)*2sxy=sxy+(x(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(f print *,F=,f,F=，回归方程显着if(f= print *,F=,f,F=,回归方程不显着End 运行成果：*实习七 计算给定数据的11年滑动平均和累积距平 实习时间：第15周周三1、2节利用数据，编写11点滑动平均的程序，给出了阅读资料的fortran程序。数据在文件夹中单独给出。要求：实习报告中附出程序，并给出原数据和滑动后数据的图形（1张图）

13、和累积距平数据图形（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:studyformshiyanseven)open(3,file=d:studyformshiyanseven)open(4,file=d:studyformshiyanseven)read(2,*)(x(i),i=1,n)do j=1,n-ih+1 !moving average do i=1

14、,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 图形显示：*实习八 对给定的海温数据进行EOF分析 实习时间：第16周周三1、2节给出海表温度距平数据

15、资料，以及相应的数据描述文件，对其进行EOF分析，资料的时空范围可以根据获知。 数据在文件夹中单独给出，距平或者标准化距平处理后再进行EOF。给出了如何读取资料，为对距平或者标准化距平处理后的资料进行EOF分析。要求：实习报告中给出第一特征向量及其时间系数，并分析其时空特征。 program main parameter(mo=12,my=43,nx=18,ny=12,mt=516) dimension avf(mo,nx,ny),df(mo,nx,ny) dimension sst(nx,ny,mt),sst3(nx,ny,mt) open(1,file=d:studyformshixiei

16、ght,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,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 d

17、o 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)=then sst3(i,j,it)= else sst3(i,j,it)=(sst(i,j,it)-avf(k,i,j)/df(k,i,j) end if end do end do end doe

18、nd do! output file open(2,file=d:studyformshixieight,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 do close(2) close(1)end 运行程序后：得到文件，即图中standard下面为老师给出程序，为运行方便将部分简略。 PROGRAM EOFPW PARAMETER(M=516,N=216,MNH=216,KS=-1,KV=10

19、,KVT=2) PARAMETER(ff=,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:studyformshixieight,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

20、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(*,*)ok3 CALL FORMA(N,M,MNH,F,A) write(*,*)ok4 CALL JCB(MNH,A,S, write(*,*)ok5 CALL ARRANG(KV,MNH,A,ER,S) write(*,*)ok6 CALL TCOEFF(KVT,KV,N,M,MNH,S,F,V,e

21、vf,tcf,ER) write(*,*)ok7 call test3(N,ff,nf,evf,kvt) write(*,*)ok8 open(21,file=d:studyformshixieight,form=unformatted,access=direct,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:studyformshixieight,form=unformatted,a

22、ccess=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 open(21,file=d:studyformshixieight) 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)= endd

23、o do i=1,n1 do k=1,m if(f(i,k).then f(i,k)= 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) if goto 30 endif do i=1,n1 avf(i)= enddo if goto 5 endif do i=1,n1 df(i)= enddo5 continue DO 141 I=1,N1 if(nf(i).

24、goto 141 do 12 j=1,m 12 AVF(I)=AVF(I)+F(I,J) AVF(I)=AVF(I)/M DO 14 J=1,M 14 F(I,J)=F(I,J)-AVF(I) 141 CONTINUE IF THEN RETURN ELSE DO 241 I=1,N1 if(nf(i). goto 241 DO 22 J=1,M 22 DF(I)=DF(I)+F(I,J)*F(I,J) DF(I)=SQRT(DF(I)/M) DO 24 J=1,M 24 F(I,J)=F(I,J)/DF(I) 241 CONTINUE ENDIF 30 CONTINUE RETURN END SUBROUTINE FORMA(N,M,MNH,F,A) DIMENSION F(N,M),A(MNH,MNH) IF(M-N) 40,50,50 40 DO 44 I=1,MNH DO 44 J=I,MNH A(I,J)= DO 42 IS=1,N 42 A(I,J)=A(I,J)+F