- 积分
- 116
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2015-5-25
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
先贴上程序
PROGRAM sub_02
IMPLICIT NONE
INTEGER,PARAMETER :: ns=756,ny=58,nm=12,nd=31,ny2=56
INTEGER year,mon,day,i,j
INTEGER is,iy,im,id
REAL a(12)
CHARACTER*5 sta(ns)
REAL,ALLOCATABLE ::staid(:),lat(:),lon(:),hgt(:),pressure(:,:,:,:),tema(:,:,:,:),temh(:,:,:,:),teml(:,:,:,:),rh(:,:,:,:),rain(:,:,:,:),wind(:,:,:,:),sun(:,:,:,:)
allocate (staID(NS),lat(NS),lon(NS),hgt(NS),pressure(NS,NY,NM,ND),temA(NS,NY,NM,ND),temH(NS,NY,NM,ND),temL(NS,NY,NM,ND),rh(NS,NY,NM,ND),rain(NS,NY,NM,ND),wind(NS,NY,NM,ND),sun(NS,NY,NM,ND))
print *,'ok'
!!!!日值说明,日平均气压、日平均气温、日最高气温、日最低气温、日平均相对湿度、日降水量、日平均风速、日日照时数
!===============READ================================================================
OPEN(11,FILE='C:\FORTRAN\sub_02\station.txt')
DO i=1,ns
READ(11,*) staID(i)
ENDDO
CLOSE(11)
pressure=1.e36
tema=1.e36
temh=1.e36
teml=1.e36
rh =1.e36
wind=1.e36
sun =1.e36
print *,'ok1'
!!!!批处理!!!!
!===========1950-2005=========================================================
DO is=1,ns
WRITE(sta(is),'(I5)') INT(staid(is))
print *,'ok2'
OPEN(22,FILE='C:\FORTRAN\sub_02\1\SURF_CLI_CHN_MUL_DAY-'//sta(is)//'.txt',err=8000)
200 READ(22,*,err=100)(a(j),j=1,4),year,mon,day,(a(j),j=5,12)
lat(is)=a(2)
lon(is)=a(3)
hgt(is)=a(4)
pressure(is,year-1949,mon,day)=a(5)
tema(is,year-1949,mon,day)=a(6)
temh(is,year-1949,mon,day)=a(7)
teml(is,year-1949,mon,day)=a(8)
rh (is,year-1949,mon,day)=a(9)
rain(is,year-1949,mon,day)=a(10)
wind(is,year-1949,mon,day)=a(11)
sun (is,year-1949,mon,day)=a(12)
GOTO 200
CLOSE(22)
100 CONTINUE
8000 CONTINUE
ENDDO
!==========2006&2007==================================================
DO is=1,ns
WRITE(sta(is),'(I5)') INT(staid(is))
OPEN(33,FILE='C:\FORTRAN\sub_02\2\SURF_CLI_CHN_MUL_DAY_0607-'//sta(is)//'.txt',err=400)
300 READ(33,*,err=400)(a(j),j=1,4),year,mon,day,(a(j),j=5,12)
lat(is)=a(2)
lon(is)=a(3)
hgt(is)=a(4)
pressure(is,year-1949,mon,day)=a(5)
tema(is,year-1949,mon,day)=a(6)
temh(is,year-1949,mon,day)=a(7)
teml(is,year-1949,mon,day)=a(8)
rh (is,year-1949,mon,day)=a(9)
rain(is,year-1949,mon,day)=a(10)
wind(is,year-1949,mon,day)=a(11)
sun (is,year-1949,mon,day)=a(12)
GOTO 300
CLOSE(33)
400 CONTINUE
ENDDO
!============READ END================================================
!!!!空白处理!!!!
DO IS=1,NS
lat(is)=INT(lat(is)/100)+MOD(lat(is),100.)/60.
lon(is)=INT(lon(is)/100)+MOD(lon(is),100.)/60.
IF(hgt(is).ge.100000)hgt(is)=1.e36
!!!!当海拔高度为100000+xxxxx时,海拔高度为估测
DO iy=1,NY
DO im=1,nm
DO id=1,nd
IF(abs(pressure(is,year-1949,mon,day)).eq.32744.or.abs(pressure(is,year-1949,mon,day)).eq.32766.or.pressure(is,year-1949,mon,day).eq.32700)then
pressure(is,year-1949,mon,day)=1.e36
ENDIF
IF(abs(temA(is,year-1949,mon,day)).eq.32744.or.abs(temA(is,y&
ear-1949,mon,day)).eq.32766.or. &
temA(is,year-1949,mon,day).eq.32700)then
temA(is,year-1949,mon,day)=1.e36
ENDIF
IF(abs(temH(is,year-1949,mon,day)).eq.32744.or.abs(temH(is,y&
ear-1949,mon,day)).eq.32766.or. &
temH(is,year-1949,mon,day).eq.32700)then
temH(is,year-1949,mon,day)=1.e36
ENDIF
IF(abs(temL(is,year-1949,mon,day)).eq.32744.or.abs(temL(is,y&
ear-1949,mon,day)).eq.32766.or. &
temL(is,year-1949,mon,day).eq.32700)then
temL(is,year-1949,mon,day)=1.e36
!!!!以上是对日平均气压、日平均气温、日最高气温和日最低气温在空白、缺测、微量或结冰情况下的处理均为1.e36,和初始值相同
ENDIF
IF(abs(rh(is,year-1949,mon,day)).eq.32744.or.abs(rh(is,year-&
1949,mon,day)).eq.32766)then
rh(is,year-1949,mon,day)=1.e36
ELSEIF(rh(is,year-1949,mon,day).eq.32700)then
rh(is,year-1949,mon,day)=0
ENDIF
!!!!以上是对日平均相对湿度的处理,空白、缺测情况下为1.e36,微量或结冰情况下为0
IF(abs(rain(is,year-1949,mon,day)).eq.32744.or.abs(rain(is,y&
ear-1949,mon,day)).eq.32766)then
rain(is,year-1949,mon,day)=1.e36
ELSEIF(rain(is,year-1949,mon,day).eq.32700)then
rain(is,year-1949,mon,day)=0
ELSEIF(rain(is,year-1949,mon,day).gt.30000.and.rain(is,year-&
1949,mon,day).lt.31000)then
rain(is,year-1949,mon,day)=rain(is,year-1949,mon,day)-30000
ELSEIF(rain(is,year-1949,mon,day).gt.31000.and.rain(is,year-&
1949,mon,day).lt.32000)then
rain(is,year-1949,mon,day)=rain(is,year-1949,mon,day)-31000
ELSEIF(rain(is,year-1949,mon,day).gt.32000.and.rain(is,year-&
1949,mon,day).lt.33000)then
rain(is,year-1949,mon,day)=rain(is,year-1949,mon,day)-32000
ENDIF
!!!!以上是对日降水量的处理,空白、缺测情况下为1.e36,微量或结冰情况下为0,30xxx时,有雨雪,降
!!!!水量为xxx;31xxx时,有雪(雨夹雪、雪暴),降水量为xxx;32xxx时,有雾露霜,降水量为xxx
IF(abs(wind(is,year-1949,mon,day)).eq.32744.or.abs(wind(is,y&
ear-1949,mon,day)).eq.32766)then
wind(is,year-1949,mon,day)=1.e36
ELSEIF(wind(is,year-1949,mon,day).gt.1000.and.wind(is,year-1&
949,mon,day).lt.2000)then
wind(is,year-1949,mon,day)=wind(is,year-1949,mon,day)-1000
ENDIF
!!!!以上是对日平均风速的处理,空白、缺测情况下为1.e36,1000+xxx(即风速大于等于xxx时,加上1000)时,
!!!!风速为xxx
IF(abs(sun(is,year-1949,mon,day)).eq.32744.or.abs(sun(is,yea&
r-1949,mon,day)).eq.32766.or. &
sun(is,year-1949,mon,day).eq.32700)then
sun(is,year-1949,mon,day)=1.e36
ENDIF
!!!!以上是对日日照时数的处理,空白、缺测、微量或结冰情况下的处理均为1.e36
ENDDO
ENDDO
ENDDO
ENDDO
!==================================================================
OPEN(44,FILE='C:\FORTRAN\sub_02\RAIN.TXT')
DO is=1,ns
DO iy=1,ny
DO im=1,nm
WRITE(44,*)staid(is),lat(is),lon(is),(rain(is,iy,im,id),&
id=1,31)
ENDDO
ENDDO
ENDDO
CLOSE(44)
DEALLOCATE (staid,lat,lon,hgt,pressure,tema,temh,teml,rh,rain,&
wind,sun)
END
出现的错误
现在卡在读取第一个站点50136的数据,请各路大神看一下,3q
|
|