[Fortran] 纯文本查看 复制代码
program PWV_classify
implicit none
integer :: stat,year,mon,day,dayofyear,hr,stat2,k,last_day,last_mon,last_year,last_dayofyear
integer :: stat3,stat4,x
integer :: stat5,stat6,y
character :: filename1*13,filename_new1*12
character :: filename2*12,filename_new2*12
character :: filename3*12,filename_new3*12
real*8 :: time_all,pwv_mm,temp_deg,pr_hpa,sum,pwv_mean_day
real*8 :: pwv_mean_mon
real*8 :: pwv_mean_year
open (11,file='day.txt')
close (11,status='delete')
call system ('for %f in (*_time.pwv) do echo %f >> day.txt')
open (11,file='day.txt')
!--------------------------------------------------> calculate pwv's day mean
stat=0
do while (stat==0)
read (11,*,iostat=stat) filename1
if (stat/=0) exit
open (22,file=filename1)
filename_new1=filename1(1:4)//'_day.pwv'
open (33,file=filename_new1)
stat2=0
k=0
last_day=0
last_mon=0
do while (stat2==0)
read(22,*,iostat=stat2) time_all,pwv_mm,temp_deg,pr_hpa,year,dayofyear,hr
if (stat2/=0) exit
call convertday2(year,mon,day,hr,time_all,dayofyear)
k=k+1
if (k==1) then
sum = pwv_mm
else
sum = sum + pwv_mm
end if
if (day/=last_day .and. k==1) then
last_year=year
last_mon=mon
last_day=day
last_dayofyear=dayofyear
else if (day/=last_day .and. k/=1) then
pwv_mean_day = sum /(k*1D0)
write (33,'(i4.4,1x,i2.2,1x,i2.2,3x,i3.3,f12.7)') year,last_mon,last_day,last_dayofyear,pwv_mean_day
last_year=year
last_mon=mon
last_day=day
last_dayofyear=dayofyear
k=0
end if
end do
!-------------------------------------> calcuate the last day
pwv_mean_day = sum /(k*1D0)
write (33,'(i4.4,1x,i2.2,1x,i2.2,3x,i3.3,f12.7)') year,mon,day,dayofyear,pwv_mean_day
k=0
!-------------------------------------> calculate the last day
close (22)
close (33)
end do
close (11)
!--------------------------------------------------> calculate pwv's day mean
!--------------------------------------------------> calculate pwv's month mean
open (44,file='month.txt')
close (44,status='delete')
call system ('for %f in (*_day.pwv) do echo %f >> month.txt')
open (44,file='month.txt')
stat3=0
do while (stat3==0)
read (44,*,iostat=stat3) filename2
if (stat3/=0) exit
open (55,file=filename2)
filename_new2=filename2(1:4)//'_mon.pwv'
open (66,file=filename_new2)
stat4=0
x=0
last_mon=0
last_year=0
do while (stat4==0)
read (55,*,iostat=stat4) year,mon,day,dayofyear,pwv_mean_day
if (stat4/=0) exit
x=x+1
if (x==1) then
sum=pwv_mean_day
else
sum=sum+pwv_mean_day
end if
if (last_mon/=mon .and. x==1 ) then
last_mon=mon
last_year=year
else if (last_mon/=mon .and. x/=1) then
pwv_mean_mon=sum/(x*1D0)
write (66,'(i4.4,1x,i2.2,3x,f12.7)') year,last_mon,pwv_mean_mon
last_mon=mon
last_year=year
x=0
end if
end do
!-------------------------------------> calcuate the last mon
pwv_mean_mon=sum/(x*1D0)
write (66,'(i4.4,1x,i2.2,3x,f12.7)') year,mon,pwv_mean_mon
x=0
!-------------------------------------> calcuate the last mon
close (55)
close (66)
end do
close (44)
!--------------------------------------------------> calculate pwv's month mean
!--------------------------------------------------> calculate pwv's year mean
open (77,file='year.txt')
close (77,status='delete')
call system ('for %f in (*_mon.pwv) do echo %f >> year.txt')
open (77,file='year.txt')
stat5=0
do while (stat5==0)
read (77,*,iostat=stat5) filename3
if (stat5/=0) exit
open (88,file=filename3)
filename_new3=filename3(1:4)//'_mon.pwv'
open (99,file=filename_new3)
stat6=0
y=0
last_year=0
do while (stat6==0)
read (88,*,iostat=stat6) year,mon,pwv_mean_mon
if (stat6/=0) exit
y=y+1
if (y==1) then
sum=pwv_mean_mon
else
sum=sum+pwv_mean_mon
end if
if (last_year/=year .and. y==1 ) then
last_year=year
else if (last_year/=year .and. y/=1) then
pwv_mean_year=sum/(y*1D0)
write (99,'(i4.4,3x,f12.7)') last_year,pwv_mean_year
last_year=year
y=0
end if
end do
!-------------------------------------> calcuate the last year
pwv_mean_year=sum/(y*1D0)
write (99,'(i4.4,1x,i2.2,3x,f12.7)') last_year,pwv_mean_year
y=0
!-------------------------------------> calcuate the last year
close (88)
close (99)
end do
close (77)
!--------------------------------------------------> calculate pwv's year mean
end program
!-----------------------------------------------------------------------------
subroutine convertday2 (yr,mon,day,hr,timeall,dayofyear)
integer :: yr,mon,day,hr,dayofyear,k
real*8 :: timeall
integer :: DAY_OF_MONTH(12) , Days_Of_Year
!---------------------------------------------------------> determine year
DAY_OF_MONTH = (/31,28,31,30,31,30,31,31,30,31,30,31/)
if (MOD(yr,4)==0 .and. MOD(yr,100)/=0) then
DAY_OF_MONTH(2) = 29
Days_Of_Year = 366
else if(MOD(yr,4)==0 .and. MOD(yr,100)==0 .and. MOD(yr,400)==0) then
DAY_OF_MONTH(2) = 29
Days_Of_Year = 366
else !(MOD(yr,4)==0 .and. MOD(yr,100)==0 .and. MOD(yr,400)/=0) then
Days_Of_Year = 365
end if
!---------------------------------------------------------> determine year
!---------------------------------------------------------> calaulate dayofyear
! if ( mon > 1 ) then
! dayofyear = day + sum( DAY_OF_MONTH(1:mon-1) )
! else
! dayofyear = day
! end if
!--------------------------------------------------------> calaulate dayofyear
!--------------------------------------------------------> calculate timeall
timeall = yr*1.0D0+(dayofyear-0.99D0)/Days_Of_Year + 1.0D0*hr/(Days_Of_Year*24.0D0)
!--------------------------------------------------------> calculate timeall
!--------------------------------------------------------> dayofyear convert to mon and day
k=1
day=dayofyear
do while (.true.)
day=day-DAY_OF_MONTH(k)
if (day<=0) exit
k=k+1
end do
mon = k
day = day + DAY_OF_MONTH(k)
!--------------------------------------------------------> dayofyear convert to mon and day
end subroutine
F 币
贡献
权杖
窥视卡
发表于 2014-12-6 09:48:04
微信
收藏
点赞
点踩
提升卡
置顶卡
变色卡
千斤顶
显身卡
发表于 2014-12-6 10:55:39
楼主