[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