气象统计方法实习报告.docx

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

2、源程序program ex_grads implicit none integer,parameter:nx=37,ny=17,nz=4,nt=12 integer i,j,iz,itreal var(nx,ny,nz,nt),cl(nx,ny,nt),sum,jp(nx,ny,nz,nt),jfc(nx,ny,nt)! Opening file open(10,file=g:gradsdatah500.dat) do iz=1,nz do it=1,nt read(10,1000) read(10,3000) (var(i,j,iz,it),i=1,nx),j=1,ny) enddo end

3、do1000 format(2i7)2000 format(37f6.2)3000 format(37f8.1)4000 format(37f7.2) close(10)!Output open(16,file=g:gradsdatah500.grd,form=binary) do iz=1,nz do it=1,nt write(16) (var(i,j,iz,it),i=1,nx),j=1,ny) enddo enddo!Calculating the Climatological Field do it=1,nt do i=1,nx do j=1,ny sum=0 do iz=1,nz

4、sum=sum+var(i,j,iz,it) enddo cl(i,j,it)=sum/4 enddo enddoenddo! Output climate-file open(12,file=g:gradsdataclimate.grd,form=binary) do it=1,nt write(12) (cl(i,j,it),i=1,nx),j=1,ny) enddo !Calculating the Anomaly do iz=1,nz do it=1,nt do i=1,nx do j=1,ny jp(i,j,iz,it)=var(i,j,iz,it)-cl(i,j,it) enddo

5、 enddo enddo enddo open(13,file=g:gradsdataanomaly.grd,form=binary)!Output anomaly-file do iz=1,nz do it=1,nt write(13) (jp(i,j,iz,it),i=1,nx),j=1,ny) enddo enddo!Calculating the Mean-square Deviation do it=1,nt do i=1,nx do j=1,ny sum=0 do iz=1,nz sum=sum+(jp(i,j,iz,it)*2 enddo jfc(i,j,it)=sqrt(sum

6、/4) enddo enddoenddo!Output mean-square deviation-file open(14,file=g:gradsdatadeviation.grd,form=binary) do it=1,nt write(14) (jfc(i,j,it),i=1,nx),j=1,ny) enddo end(2)、grads文件open g:gradsdata*.ctl (*为所求场对应的ctl文件名)set lat 0 40set lon 60 150set lev 500enable print g:gradsdata*.gmf （*为所求场名称）i=1while(i

7、max_y) THEN max_y=rxy_ty(ty) 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_

8、w=rxy_tw(tw) k=tw END IFEND DOPRINT (冬季平均气温绝对值最大自相关系数rxy_tw=,f7.4,/,滞后时间长度 k=,I2),rxy_tw(k),kk=0!落后交叉相关系数DO tyw=1,N/2 DO i=1,N-tyw rtyw(i)=(y(i)-avr_y)/sy)*(w(i+tyw)-avr_w)/sw) rxy_tyw(tyw)=rxy_tyw(tyw)+rtyw(i) END DO rxy_tyw(tyw)=rxy_tyw(tyw)/(N-tyw) rxy_tyw(tyw)=ABS(rxy_tyw(tyw) IF(rxy_tyw(tyw)max

9、_yw) THEN max_yw=rxy_tyw(tyw) k=tyw END IFEND DOPRINT (全年平均温度与冬季平均气温之间的落后交叉相关系数rxy_tyw=,f7.4,/,滞后时间长度 k=,I2),rxy_tyw(k),kEND（2）、程序运行结果：实习三 分析中国夏季降水线性趋势的分布特征 1.资料介绍及要求：利用数据160zhan-rainfall-summer.txt，编写求1982-2006年中国160站各站夏季降水线性倾向率，给出分布图，并进行简单分析。read-rain.for给出了阅读资料的fortran程序。数据在文件夹中单独给出。2.实习结果（1）.mat

10、lab程序%编写求1982-2006年中国160站各站夏季降水线性倾向率clear allclcfid=fopen(E:/160zhan-rainfall-summer.txt,rt);tline=fgets(fid);data1=fscanf(fid,%f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f,28,160);data2=data1;fclose(fid);for i=1:160; j(i,1:25)=data2(i,4:28); n1=1982:1:2006; pp(

11、i,:)=polyfit(n1,j(i,1:25),1);endb=pp(:,1);jd=data2(:,3);wd=data2(:,2);jdc=75:0.5:135;wdc=18:.5:55;bz=griddata(jd,wd,b,jdc,wdc,cubic);c=contour(jdc,wdc,bz)xlabel(精度);ylabel(纬度);title(1982-2006年中国160站各站夏季降水线性倾向率分布图)（2）.程序运行结果实习四 求给定数据的一元线性回归方程1、资料介绍及要求利用下表数据，以环流指标为预报因子，气温为预报量，计算气温和环流指标之间的一元线性回归方程，并对回归

12、方程进行检验。年份气温T环流指标19510.93219521.22519532.22019542.4261955-0.52719562.5241957-1.128195802419596.21519602.71619613.2241962-1.13019632.52219641.23019651.82419660.63319672.42619682.52019691.2321970-0.835答案：F=20.18F=4.41，回归方程显著2、实习结果（1）、MATLAB程序%实习四 求给定数据的一元线性回归方程ClimateData=xlsread(F:气象统计方法实验数据气象统计实验四数据.

13、xls); %从Excel文件读取数据x=ClimateData(:,3); %提取ClimateData的第三列，即环流指标y=ClimateData(:,2); %提取ClimateData的第三列，即气温Txdata=ones(size(x,1),1),x; %在原始数据x的左边加一列1，即模型包含常数项b,bint,r,rint,s=regress(y,xdata); %调用regress函数作一元线性回归yhat=xdata*b; %计算y的估计值%定义元胞数组，以元胞数组形式显示系数的估计值和估计值得95%置信区间head1=系数的估计值,估计值的95%置信下限,估计值的95%置信

14、上限;head1;num2cell(b,bint)%定义元胞数组，以元胞数组形式显示y的真实值、y的估计值、残差和残差的95%置信区间head2=y的真实值,y的估计值,残差,残差的95%置信下限,残差的95%置信上限;%同时显示y的真实值，y的估计值、残差和残差的95%置信区间head2;num2cell(y,yhat,r,rint)%定义元胞数组，以元胞数组形式显示判定系数、F统计量的观测值、检验的P值和误差方差的估计值head3=判定系数,F统计量的观测值,检验的P值,误差方差的估计值;head3;num2cell(s)（2）、程序运行结果ans = 系数的估计值 估计值的95%置信下限

15、 估计值的95%置信上限 7.5095 4.6554 10.3637 -0.2343 -0.3433 -0.1253ans = y的真实值 y的估计值 残差 残差的95%置信下限 残差的95%置信上限 0.9000 0.0123 0.8877 -1.5633 3.3388 1.2000 1.6523 -0.4523 -3.0331 2.1285 2.2000 2.8237 -0.6237 -3.1171 1.8696 2.4000 1.4180 0.9820 -1.5611 3.5251 -0.5000 1.1837 -1.6837 -4.1250 0.7576 2.5000 1.8866 0

16、.6134 -1.9531 3.1800 -1.1000 0.9494 -2.0494 -4.4072 0.3084 0 1.8866 -1.8866 -4.2867 0.5136 6.2000 3.9952 2.2048 0.1971 4.2125 2.7000 3.7609 -1.0609 -3.3535 1.2317 3.2000 1.8866 1.3134 -1.1840 3.8108 -1.1000 0.4808 -1.5808 -3.9959 0.8342 2.5000 2.3552 0.1448 -2.4137 2.7034 1.2000 0.4808 0.7192 -1.800

17、1 3.2385 1.8000 1.8866 -0.0866 -2.6717 2.4986 0.6000 -0.2220 0.8220 -1.5996 3.2437 2.4000 1.4180 0.9820 -1.5611 3.5251 2.5000 2.8237 -0.3237 -2.8318 2.1844 1.2000 0.0123 1.1877 -1.2301 3.6056 -0.8000 -0.6906 -0.1094 -2.4794 2.2606ans = 判定系数 F统计量的观测值 检验的P值 误差方差的估计值 0.5313 20.4045 2.6673e-004 1.5134（3

18、）、结果分析 从输出的结果看，常数项和回归系数的估计值分别为7.5095和-0.2343，从而可以写出线性回归方程为 回归系数估计值的置信区间为 -0.3433,-0.1253。 对回归直线进行显著性检验，原假设和对立假设分别为 检验P的值为2.667310-4 plot(x,y,k.,Markersize,15) %原始数据散点图hold onplot(x,yhat,linewidth,3) %回归直线图xlabel(环流指标(x)%标注x轴ylabel(气温(y)%标注y轴legend(原始散点,回归直线)%加标注框实习五 对给定的海温数据进行EOF分析 1、资料介绍给出海表温度距平数据资

19、料sstpx.grd，以及相应的数据描述文件sstpx.ctl，对其进行EOF分析，资料的时空范围可以根据sstpx.ctl获知。 数据在文件夹中单独给出，距平或者标准化距平处理后再进行EOF。Zhunsst.for给出了如何读取资料，Ssteof.for为对距平或者标准化距平处理后的资料进行EOF分析。2、要求实习报告中给出第一特征向量及其时间系数，并分析其时空特征。 3、实习结果：（1）、FORTRAN源程序! prepare data for eof analysis! the program is to normalize sea surface temperature(SST)! mt: the length of time series;! mo: the month numbers; my: the year numbers;! sst: sea s