爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 10644|回复: 27

[求助] fortran算出来的全都是0,求高手指教。

[复制链接]

新浪微博达人勋

发表于 2012-11-30 16:25:19 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
本人最近写的fortran程序不知为什么能够算出来但是全都是0。很纠结,求大神指教!!!
附上本人的程序两个,水平较差请见谅。
1.空间滤波的程序:
program lbx
implicit none
real,dimension(1:144,1:73,1:17,1:720)::hgt
real,dimension(0:12,1:73,1:6)::a,b
real,dimension(1:144,1:73,1:6)::lb0t3,lb6t12
real,dimension(0:12,1:144,1:73,1:6)::lb
real,dimension(1:144,1:73,1:17)::ahgt
real::w
integer::i,j,k,t,irec,m
!读取数据
open(11,file='hgt.dat',status='old',form='unformatted',access='direct',recl=144*73)
irec=1
do t=1,720
  do k=1,17
  read(11,rec=irec)((hgt(i,j,k,t),i=1,144),j=1,73)
  irec=irec+1
  end do
end do
close(11)
print*,'*'
!这就开始空间滤波了
!首先对位势高度时间平均,这样再固定气压面,位势高度就仅仅是水平空间的函数了
do k=12,17
  do j=1,73
   do i=1,144
     ahgt(i,j,k)=0
     do t=1,720
     ahgt(i,j,k)=ahgt(i,j,k)+hgt(i,j,k,t)
     end do
     ahgt(i,j,k)=ahgt(i,j,k)/720
   end do
  end do
end do
print*,'**'
!第二步,对所有谐波的系数进行初始设置
a=0
b=0
!第三步,计算每个谐波的系数
do k=1,6
do m=0,12
  w=360/144
  do j=1,73
   do i=1,144
    a(m,j,k)=a(m,j,k)+ahgt(i,j,k)*cosd(m*w*i)
    b(m,j,k)=b(m,j,k)+ahgt(i,j,k)*sind(m*w*i)
   end do
  end do
end do
end do
!再除以空间序列的长度
a=a/72
b=b/72
do k=1,6
do j=1,73
  a(0,j,k)=a(0,j,k)/2
end do
end do
print*,'***'
!第四步,累加各个谐波,得到各个波段的空间滤波场
do k=1,6
do m=1,12
  do j=1,73
   do i=1,144
   lb(m,i,j,k)=a(0,j,k)+a(m,j,k)*cosd(m*w*i)+b(m,j,k)*sind(m*w*i)
   end do
  end do
end do
end do
!得到平均场即0波滤波
do k=1,6
do j=1,73
  do i=1,144
   lb(0,i,j,k)=a(0,j,k)
  end do
end do
end do
!得到0-3波和6-12波的空间滤波
do k=1,6
do j=1,73
  do i=1,144
   lb0t3(i,j,k)=lb(0,i,j,k)+lb(1,i,j,k)+lb(2,i,j,k)+lb(3,i,j,k)
   lb6t12(i,j,k)=lb(0,i,j,k)+lb(6,i,j,k)+lb(7,i,j,k)+lb(8,i,j,k)+lb(9,i,j,k)+lb(10,i,j,k)+lb(11,i,j,k)+lb(12,i,j,k)
  end do
end do
end do
print*,'****'
!第五步,结果输出
open(11,file='lb.dat',status='replace',form='unformatted',access='direct',recl=144*73*6)
irec=1
do m=0,12
  write(11,rec=irec)(((lb(m,i,j,k),i=1,144),j=1,73),k=1,6)
  irec=irec+1
end do
close(11)
open(11,file='lb0-3.dat',status='replace',form='unformatted',access='direct',recl=144*73*6)
irec=1
  write(11,rec=irec)(((lb0t3(i,j,k),i=1,144),j=1,73),k=1,6)
close(11)
open(11,file='lb6-12.dat',status='replace',form='unformatted',access='direct',recl=144*73*6)
irec=1
  write(11,rec=irec)(((lb6t12(i,j,k),i=1,144),j=1,73),k=1,6)
close(11)
print*,'ok!'
!完成
end program

2.一个筛选数据的程序
   program count
character(60)::filename,filename1
    character(2)::p1,p2,p3,p4,p5
    character(79)::a
integer::i,j,e,f,g,k,t,n,m
integer,dimension(1:40)::x
character(15),dimension(1:12)::w
    open(13,file='readfile.txt',status='unknown',form='formatted')
do i=08,10
write(p1,'(i2.2)') i
     do j=1,31
write(p2,'(i2.2)') j
   do e=0,24
write(p3,'(i2.2)') e
       do f=0,60
write(p4,'(i2.2)') f
     do g=0,60
write(p5,'(i2.2)') g
filename='00PCASP_X22012'//p1//p2//p3//p4//p5//'.csv'
filename1=p1//p2//p3//p4//p5//".dat"
    open(11,file=trim(filename),status='old',form='formatted',iostat=it)
    if(it.ne.0) cycle
open(12,file=trim(filename1),status='unknown',form='unformatted',access='direct',recl=40)
     do k=1,36
  read(11,*)a
     end do
  m=0
do k=1,200000
  read(11,*,iostat=n)(w(t),t=1,12),(x(t),t=1,40)
  if (n<0) exit
  print*,x(5)
     write(12,rec=k)(x(t),t=1,40)
  m=m+1
end do
write(13,*)filename1,m
close(11)
close(12)
end do
    end do
end do
end do
    end do
close(13)
end program

以上问题严重困扰着我,请各位不吝赐教。


密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2012-11-30 16:54:09 | 显示全部楼层
使用print屏幕输出一下,看看到哪一步计算变成零的,这样就能确定错误的大概位置了。你这一长串放在这儿,也一下看不出来什么东西
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-11-30 17:26:28 | 显示全部楼层

非常感谢!!!第一个这个我再检查下,但是第二个就是个特别简单的读取再写入,而且我看过读出来的没问题,写进去的用grads打开就全是零了。
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2012-11-30 18:52:40 | 显示全部楼层
lqouc 发表于 2012-11-30 17:26
非常感谢!!!第一个这个我再检查下,但是第二个就是个特别简单的读取再写入,而且我看过读出来的没问题 ...

第二个程序是用来写成二进制文件的?这个程序太高深的了,我没看出来是用来把资料转成二进制文件的
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2012-11-30 18:56:02 | 显示全部楼层
可以把第二个程序索要读取的文件发上来看看
就像river所说,除了这样的问题,就一步步的print出来结果,看哪里出问题了
密码修改失败请联系微信:mofangbao

新浪微博达人勋

0
早起挑战累计收入
发表于 2012-11-30 18:57:24 | 显示全部楼层
又是贴了一大坨。。。无从看起
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-11-30 19:02:06 | 显示全部楼层
river 发表于 2012-11-30 18:52
第二个程序是用来写成二进制文件的?这个程序太高深的了,我没看出来是用来把资料转成二进制文件的

额,这个应该没有错吧,我的程序意思就是把一堆csv格式的文件名很混乱的数据都读出来,找到我要用的部分,然后再写成无格式的dat文件。我在读取的时候全都print了,然后没有做任何处理就写入dat,于是,就变成一个全0场了。
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-11-30 19:05:47 | 显示全部楼层
topmad 发表于 2012-11-30 18:56
可以把第二个程序索要读取的文件发上来看看
就像river所说,除了这样的问题,就一步步的print出来结果,看 ...

00PCASP_X220120825085214.csv (44.86 KB, 下载次数: 3)
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-11-30 19:06:46 | 显示全部楼层
topmad 发表于 2012-11-30 18:56
可以把第二个程序索要读取的文件发上来看看
就像river所说,除了这样的问题,就一步步的print出来结果,看 ...

这个是刚才那个第二个程序的数据,我就先发一个上来。
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-11-30 19:08:21 | 显示全部楼层
mofangbao 发表于 2012-11-30 18:57
又是贴了一大坨。。。无从看起

额,不好意思啊。第一次发帖求助。是在走投无路了,第一个的问题就在第三步求谐波的系数a和b。那之后就全是0了。
密码修改失败请联系微信:mofangbao
您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表