intel visual Fortran 处理站点资料遇到的问题

夏天很美

本人很菜，虚心求教！下面是Fortran处理站点资料源代码，红色加粗部分是我写入的站点资料，虽然能得出数据，但是我用grads怎么也写ctl怎么也得不到.map映射文件，我写的求映射文件的ctl在Fortran代码的后面：
program rain
implicit none
integer i,j,nlev,nflag
real tim
realr6(160,60),r7(160,60),r8(160,60)
real lat(160),lon(160)
realr678(160,60),s(160),a(160),a3(160),s3(160)
realrp6(160,60),rp3(160,60)
reals678(160),rxb(160,60)
reald6(160,60),d3(160,60)
character*8 id(160)


!打开需要读入数据的文件
open(12,file='f:\js\r1606.txt')
open(20,file='f:\js\r1607.txt')
open(30,file='f:\js\r1608.txt')
open(40,file='f:\js\lat_lon.txt')
write(*,*)'ok1'

!打开需要写入数据的文件
open(45,file='f:\js\r678.txt')
open(50,file='f:\js\r678.grd',form='binary')
!open(60,file='f:\js\rp6.grd',form='binary')
!open(70,file='f:\js\rp3.grd',form='binary')
!open(80,file='f:\js\rxb.grd',form='binary')
!open(90,file='f:\js\d6.grd',form='binary')
!open(11,file='f:\js\d3.grd',form='binary')
write(*,*)'ok2'

!读入数据

read(12,*)((r6(j,i),j=1,160),i=1,60)
close(12)
write(*,*)'open hello'
read(20,*)((r7(j,i),j=1,160),i=1,60)
close(20)
write(*,*)'open hello'
read(30,*)((r8(j,i),j=1,160),i=1,60)
close(30)
write(*,*)'open hello'

do j=1,160
read(40,*)id(j),lat(j),lon(j)
enddo
close(40)
write(*,*)id(1),lat(1),lon(1)

write(*,*)'ok3'

!计算每年夏季的降水量
do i=1,60
do j=1,160
r678(j,i)=r6(j,i)+r7(j,i)+r8(j,i)
enddo
enddo
write(45,*)((r678(j,i),j=1,160),i=1,60)
close(45)

!计算每个站点60年的夏季平均降水a（160）
do j=1,160
s(j)=0.0
do i=1,60
s(j)=s(j)+r678(j,i)
enddo
a(j)=s(j)/60.0
enddo
write(*,*)'a(1)='
print*,a(1)
write(*,*)'ok4'

!将降雨量标准化处理，s678(160)标准差，rxb(160,60)标准化数组
do j=1,160
callcal_s(r678(j,:),60,s678(j))
enddo
write(*,*)'标准差子程序运行没有问题'

do i=1,60
do j=1,160
rxb(j,i)=(r678(j,i)-a(j))/s678(j)
enddo
enddo


!计算每个站点1981-2010年30年的夏季平均降水a3(160)
do j=1,160
s3(j)=0.0
do i=1,30
s3(j)=s3(j)+r678(j,i+30)!s3(j)每个站点夏季降水30年的总和
enddo
a3(j)=s3(j)/30.0
enddo
write(*,*)'a3(1)='
print*,a3(1)

write(*,*)'ok5'

!计算每年夏季降水对60年距平d6(160,60),和对30年距平d3(160,60)
do i=1,60
do j=1,160
d6(j,i)=r678(j,i)-a(j)
rp6(j,i)=d6(j,i)/a(j)*100.0!对应60年的降水距平百分率,单位%
d3(j,i)=r678(j,i)-a3(j)
rp3(j,i)=d3(j,i)/a3(j)*100.0!对应30年的降水距平百分率,单位%
enddo
enddo
write(*,*)'d3(1,1)='
print*,d3(1,1)
write(*,*)'d6(1,1)='
print*,d6(1,1)
write(*,*)'ok6'

!写入站点数据
do j=1,160
tim=0.0
nlev=1
nflag=1
write(50)id(j),lat(j),lon(j),tim,nlev,nflag,(r678(j,i),i=1,60)
enddo
tim=0.0
nlev=1
nflag=1
write(50)id(j-1),lat(j-1),lon(j-1),tim,nlev,nflag
close(50)

!do j=1,160
!tim=0.0
!nlev=1
!nflag=1
!write(60)id(j),lat(j),lon(j),tim,nlev,nflag,(rp6(j,i),i=1,60)
!enddo
!tim=0.0
!nlev=1
!nflag=1
!write(60)id(j-1),lat(j-1),lon(j-1),tim,nlev,nflag
!
!do j=1,160
!tim=0.0
!nlev=1
!nflag=1
!write(70)id(j),lat(j),lon(j),tim,nlev,nflag,(rp3(j,i),i=1,60)
!enddo
!tim=0.0
!nlev=0
!nflag=1
!write(70)id(j-1),lat(j-1),lon(j-1),tim,nlev,nflag
!
!do j=1,160
!tim=0.0
!nlev=1
!nflag=1
!write(80)id(j),lat(j),lon(j),tim,nlev,nflag,(rxb(j,i),i=1,60)
!enddo
!tim=0.0
!nlev=0
!nflag=1
!write(80)id(j-1),lat(j-1),lon(j-1),tim,nlev,nflag
!
!do j=1,160
!tim=0.0
!nlev=1
!nflag=1
!write(90)id(j),lat(j),lon(j),tim,nlev,nflag,(d6(j,i),i=1,60)
!enddo
!tim=0.0
!nlev=0
!nflag=1
!write(90)id(j-1),lat(j-1),lon(j-1),tim,nlev,nflag
!
!do j=1,160
!tim=0.0
!nlev=1
!nflag=1
!write(11)id(j),lat(j),lon(j),tim,nlev,nflag,(d3(j,i),i=1,60)
!enddo
!tim=0.0
!nlev=0
!nflag=1
!write(11)id(j-1),lat(j-1),lon(j-1),tim,nlev,nflag

close(1)
close(2)
close(3)
close(4)

end

!子程序部分
subroutine cal_s(a,na,s)!计算标准差子程序，a数组，na数组量，s标准差
implicit none
integer na,i
real a(na),s,ave
s=0.0
call sal_ave(a,na,ave)
do i=1,na
s=s+(a(i)-ave)**2
enddo
s=s/real(na)
s=sqrt(s)
end subroutine cal_s

subroutinesal_ave(a,na,ave) !求平均值子程序，a数组，na数组量,ave平均值
implicit none
integer i,na
real a(na),ave,asum
asum=0.0
do i=1,na
asum=asum+a(i)
enddo
ave=asum/real(na)
end subroutine sal_ave


求映射文件的ctl文件
dset f:\js\r678.grd
dtype station
stnmap f:\js\r678.map
undef -999.0
title rain678
tdef 60 linear jan1951 1yr
vars 1
r678 0 99 summer rainfall data
endvars

夏天很美

1. !写入站点数据
   
2. do i=1,60
   
3. tim=0.0
   
4. nlev=1
   
5. nflag=1
   
6. do j=1,160
   
7. write(50)id(j),lat(j),lon(j),tim,nlev,nflag,r678(j,i)
   
8. write(66,*)id(j),lat(j),lon(j),tim,nlev,nflag,r678(j,i)
   
9. enddo
   
10. nlev=0
    
11. write(50)id(j-1),lon(j-1),lat(j-1),tim,nlev,nflag
    
12. enddo
    
13. close(50)

复制代码

夏天很美

本人已解决了

wangxue910402

,正在看读入站点数据的东西

虎式啸傲

已经按照楼主的方法成功地解决了，谢谢

dyhtt29

{:lxm_24:}{:lxm_24:}

shengzhaoqiu

可以可以。。。。

Hatter

学习了