看这个程序有没有问题：程序没运行完时停止，不提示出错

yuguotianqing · 发表于 2014-2-22 09:29:49

本帖最后由 yuguotianqing 于 2014-2-22 13:30 编辑

[Fortran] 纯文本查看 复制代码

Program ForcingPrepare 
implicit none

! 1. 程序变量声明

integer,parameter :: step=15341
integer,parameter :: mmw=3
character(len=4) :: a='data'
character(len=1) :: b='_'
character(len=7) :: c
character(len=1) :: d='_'
character(len=8) :: e

integer :: num=0
integer :: ncols, nrows, i, j, k, l, m, n  ! 列数 ncols，行数nrows
real :: xllcorner, yllcorner !dem最左下角的经度xllcorner，纬度yllcorner
real :: cellsize
real :: LAT, LNG
real :: p, nump,stepp
character(len=26) :: rainname(step)
character(len=21) :: forcename
real, allocatable :: rain(:,:)
real, allocatable :: tw(:,:)
real, allocatable :: dem(:,:)  

! 2. 输入输出数据文件的打开部分：
write(*,"('VIC模型气象强迫数据文件准备程序计算开始！')")   !在屏幕上打出程序标题
write(*,*)
write(*,"('程序正在处理中，由于处理数据量较大，不要中断程序，请耐心等待......')") 
    !打开DEM文件、降水名称文件、温度风速文件
open(7,file="dem.txt",status='old',err=999)
open(8,file="rainfilename.txt",status='old',err=999)
open(9,file="tempwind.txt",status='old',err=999)

open(4, file="test.txt")

!读出dem文件的列数、行数、左下角经度、左下角纬度 
read(7,"(14x,i5)")ncols 
read(7,"(14x,i5)")nrows 
read(7,"(14x,f8.4)")xllcorner
read(7,"(14x,f7.4)")yllcorner
read(7,"(14x,f10.7)")cellsize
read(7,*)


!将DEM文件的高程数据输入到二维数组dem(i,j)中去
allocate(dem(ncols,nrows)) 

do j=1,nrows
    read(7,*)(dem(i,j),i=1,ncols)
end do


do k=1,step
    read(8,"(a26)")rainname(k)
        !write(*,*)rainname(k)
end do


!!!!!打开降水文件

do k=1,step
   open(10+k,file=rainname(k),status='old',err=999)
end do

allocate(tw(mmw,step)) 
do j=1,step
    read(9,*)(tw(i,j),i=1,mmw)
end do


! 7. 土壤、植被参数处理与输出
do j=1,nrows
        do i=1,ncols
                
        num=num+1
                
                if(dem(i,j)/=-9999) then

          LAT=yllcorner+(nrows-j+0.5)*cellsize
                  LNG=xllcorner+(i-0.5)*cellsize
                  write(c,"(f7.4)")LAT
          write(e,"(f8.4)")LNG 
                  forcename=a//b//c//d//e
          open(10+step+num,file=forcename)
                       
                  
                  do k=1,step
             read(10+k,"(/////)")  !把第k个降水文件的头语句跳过去

             allocate(rain(ncols,nrows)) 
             do n=1,nrows
                read(10+k,*)(rain(m,n),m=1,ncols)!把第K个降水文件的降水数据赋值到二维数组中
             end do
             write(10+step+num,"(f6.2)",advance='no')rain(i,j)!把第K个降水文件中当前格网点的降水写到该格网文件中
                         
                         do l=1,mmw
                            write(10+step+num,"(f7.2)",advance='no')tw(l,k)!把最高温、最低温、风速写到降水后面
                         end do
                         
                         write(10+step+num,*) !换行，准备写当前格网的第二天数据 
                         
                         deallocate(rain)
                                                       
                  end do

                  do k=1,step
                      rewind(10+k)   !将所有降水空间分布文件的读数据光标打回文件开头，非常重要！！！
          end do      
                                                      
                endif 

          nump=num
                  stepp=step
                  p=nump/stepp*100
                  !write(*,*)p

                  write(*,"(f5.2)", advance='no')p

                   write(*,"(a1)")"%"
                  
                        
        end do
        
end do

write(*,*)
write(*,*)
write(*,"('         ******************************************')")
write(*,"('         *恭喜您！数据处理完毕，您可以查看结果了! *')")
write(*,"('         ******************************************')")
write(*,*)
stop
!文件找不到就结束程序部分：
999 write(*,"('输入的文件不存在或者有错误，程序无法执行，请检查后再运行程序！')")
stop
end program ForcingPrepare