- 积分
- 36
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2013-3-1
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
各位大侠,
我最近写了一个fortran 程序(gfortran 编译器), 貌似一切都没有问题, 可是最后写成的grads站点数据,画图时, 400hp以上的数据都没有, 1000到400的数据貌似都是对的。 请各位大侠帮忙啊。现贴出我的程序和CTL文件。 现在的情况是: 程序正常运行,stnmap 也显示正常, 就是从300hpa以上就啥也显示不了了。
感觉就是300语句open或者write的不对呢
请大侠看看我open(300,然后的写法有什么错误么? gfortran 写高空站点数据用什么格式好啊? gfortran 不认binary
high_test.f90
!!提取高空的数据---一次处理一个月------------------
!!!!-------------------------------------------------------------------------------------------------
program main
implicit none
integer,parameter:: nn=800,lev=11,tt=2
!integer,parameter:: nn=800,lev=11, lev1=11,tt=2
real,parameter:: lo1=70,lo2=150,la1=15,la2=65
real,parameter::undef=9999
integer:: i,j,k,l,g,gg,m,it,j_s
integer:: nline(tt),num(tt)
integer:: nflag,nlev
integer:: ierr,irec
integer::kind1,kind2
real::tim,press(lev)
real:: lon(nn,tt),lat(nn,tt),ls(nn,tt),hgt(nn,lev,tt),wd(nn,lev,tt),&
ttd(nn,lev,tt),fx(nn,lev,tt),fs(nn,lev,tt)
real:: lonh(nn,tt),lath(nn,tt),hgth(nn,lev,tt),wdh(nn,lev,tt),&
ttdh(nn,lev,tt),tdh(nn,lev,tt),fxh(nn,lev,tt),fsh(nn,lev,tt)
real:: est(nn,lev,tt), estd(nn,lev,tt), qvs(nn,lev,tt),qes(nn,lev,tt)
real :: rh(nn,lev,tt)
integer:: zj(nn,tt),sta(nn,tt),stah(nn,tt),st(nn)
character(len=2)::yy,hour
character(len=4)::yr
character(len=2)::mon
character(len=100)::dirname,data_in,data_out
character(len=200)::filename, openname
character:: a(4)*1,b(10)*1
character(len=4):: hp(lev)
character(len=2):: sh(2)
character*8 stid(nn)
logical alive,live
character(len=10):: temp(5)
data a/'0','1','2','3'/
data b/'0','1','2','3','4','5','6','7','8','9'/
data hp/'1000','925','850','700','500','400','300','250','200','150', &
'100'/
data sh/'08','20'/
REAL , PARAMETER :: SVP1=0.6112
REAL , PARAMETER :: SVP2=17.67
REAL , PARAMETER :: SVP3=29.65
REAL , PARAMETER :: SVPT0=273.15
REAL :: ALIQ, BLIQ, CCLIQ, DLIQ
integer:: int_yr,start_date, end_date
integer:: dates(12)=(/ 31,28,31,30,31,30,31,31,30,31,30,31/)
integer:: reck
!-----------------------------屏幕读入处理时间-------------------------------
write(*,*) "enter year month:******(201401)"
read(*,"(a4,a2)")yr,mon
yy=yr(3:4)
read(yr, "(i4)") int_yr
write(*,*) "enter start_date:(01)"
read(*,"(i2)")start_date
write(*,*) "enter end_date:(31)"
read(*,"(i2)")end_date
write(*,*) "enter finish"
write(*,*) "yr,mon,start_date,end_date", yr,mon,start_date,end_date
!--------------------------------------判断是否闰年----------------------------------------
if((mod(int_yr,4).eq.0.and.mod(int_yr,100).ne.0) .or. mod(int_yr,400).eq.0) then
dates(2)=29
else
dates(2)=28
endif
!------------------------------------prepare for rh calculating -----------------------------
ALIQ = SVP1*1000.
BLIQ = SVP2
CCLIQ = SVP2*SVPT0
DLIQ = SVP3
!!!---------------------输入、输出文件夹---------------------------
data_in="D:\obs_pro\data\micaps\"
data_out="D:\obs_pro\output\met_data\high\high_stn\"
do l=1,lev
read(hp(l),*) press(l)
end do
write(*,*) press(:)
!!!!----
do j=start_date, end_dateme,'(3a,i2.2,2a)') trim(data_out),trim(yr),trim(mon),j,sh(k),'.grd'
write(*,*) 'output name', filename
! OPEN(300,FILE=trim(filename),FORM='unformatted', recordtype=' stream')
OPEN(300,FILE=trim(filename),FORM='unformatted',access='sequential')
do l=1,lev
write(openname,'(8a,i2.2,2a)') trim(data_in),trim(yr),trim(mon),'\high\plot\',trim(adjustl(hp(l))), &
'\',trim(yy),trim(mon),j,sh(k),".000"
print*,'read file', openname
open(100,file=trim(openname),status='old')
read(100,*)
read(100,*) (temp(it),it=1,5),nline(k)
print*, 'nline(k)', nline(k)
do g=1,nline(k)
read(100,*)sta(g,k),lon(g,k),lat(g,k),ls(g,k),zj(g,k),hgt(g,l,k),wd(g,l,k),&
ttd(g,l,k),fx(g,l,k),fs(g,l,k)
end do !g
close(100)
end do
do m=1,nline(k)
write(stid(m),'(i8)') stah(m,k) ! convert integer to character
enddo !m
print*, 'before write 300 lev='
tim=0.0
nlev=lev
nflag=0
irec=1
do m=1,nline(k)
write(300) stid(m), lat(m,k),lon(m,k),tim,nlev,nflag
! write(300) stid(m), lat(m,k),lon(m,k),tim,5,nflag
do l=1,lev
! do l=1,5
print*, 'before write 300, st,lev,press', m,l
write(300) press(l),hgt(m,l,k),wd(m,l,k),ttd(m,l,k),fx(m,l,k),fs(m,l,k)
print*, press(l),hgt(m,l,k),wd(m,l,k),ttd(m,l,k),fx(m,l,k),fs(m,l,k)
print*, 'press,hgt,wd(m,l,k),ttd(m,l,k),fx(m,l,k),fs(m,l,k)'
enddo !l1
enddo !m
nlev=0
write(300) stid(m-1), lat(m-1,k),lon(m-1,k),tim,nlev,nflag
close(300)
end do
end do
end program
下面是CTL文件high_test.ctl
dset D:\obs_pro\output\met_data\high\high_stn\2015110120.grd
dtype station
options sequential
stnmap D:\obs_pro\output\met_data\high\high_stn\2015110120.map
undef 9999
title high
TDEF 1 linear 12Z01NOV2015 12hr
vars 5
hgt 11 99 geopotential height
temp 11 99 temperature
ttd 11 99 T-Td
windd 11 99 wind direction
winds 11 99 wind speed
endvars
|
|