侯平均计算

onenuister

最近毕设要计算急流区域1979-2018年300hPa7月的纬向风侯平均值，所使用的是日平均资料。我先用grads计算了第一侯的数据，程序如下(grads不能一次性读40年数据，我就分成两段时间提取了)：
***1979-2018年300hPa7月第一侯平均uwnd
'reinit'
'set fwrite f:\002\05\PWJ_fir_uwnd1979-2010.grd'
'set gxout fwrite'
year=1979
while(year<=2010)
'sdfopen f:\002\uwnd.'year'.nc'
'set lat -90 90'
'set lon 0 357.5'
'set lev 300'
***闰年
if(((math_fmod(year,4)=0)&(math_fmod(year,100)!=0))|(math_fmod(year,400)=0))
'define u1=ave(uwnd,t=183,t=187)'
'define u2=aave(u1,lon=70,lon=90,lat=60,lat=65)'
'd u2'
else
***平年
'define u1=ave(uwnd,t=182,t=186)'
'define u2=aave(u1,lon=70,lon=90,lat=60,lat=65)'
'd u2'
endif
year=year+1
endwhile
'disable fwrite'
'c'
'reinit'
'set fwrite f:\002\05\PWJ_fir_uwnd2011-2018.grd'
'set gxout fwrite'
year=2011
while(year<=2018)
'sdfopen f:\002\uwnd.'year'.nc'
'set lat -90 90'
'set lon 0 357.5'
'set lev 300'
***闰年
if(((math_fmod(year,4)=0)&(math_fmod(year,100)!=0))|(math_fmod(year,400)=0))
'define u1=ave(uwnd,t=183,t=187)'
'define u2=aave(u1,lon=70,lon=90,lat=60,lat=65)'
'd u2'
else
***平年
'define u1=ave(uwnd,t=182,t=186)'
'define u2=aave(u1,lon=70,lon=90,lat=60,lat=65)'
'd u2'
endif
year=year+1
endwhile
'disable fwrite'
;

然后我又用FORTRAN读取了一下我计算的数据，程序如下：
program ceshi
    implicit none
    integer k
    parameter(k=40)
    integer i,j,it
    real PWJ(k)
    open(10,file='f:\002\05\PWJ_fir_uwnd1979-2010.grd',form='binary')
    open(11,file='f:\002\05\PWJ_fir_uwnd2011-2018.grd',form='binary')
    do it=1,32
        read(10)PWJ(it)
    enddo
    do it=33,k
        read(11)PWJ(it)
    enddo
    print*,PWJ
    end

结果如图所示，数据刚好就以2011年为分界，而且照常理来说每个月数据重复率不会这么高，已经苦思冥想了好几天了，还是不知道哪里有问题，希望得到大家指点

onenuister

我后来又分成1979-2000和2001-2018两段时间，结果就以2001年分界了，真的想不通啊做毕设太难了

lightmoon

楼主你这样，你在你endwhile前边加一个close 1应该就好了

onenuister

lightmoon 发表于 2020-3-17 18:28
楼主你这样，你在你endwhile前边加一个close 1应该就好了

太谢谢你了，数据果然没刚才那么奇怪了

软壳

楼主，想问一下，你的grads用了两个ave，结果是区域平均吗

onenuister

软壳 发表于 2020-4-5 21:19
楼主，想问一下，你的grads用了两个ave，结果是区域平均吗

第一个ave是进行侯平均，第二个aave是进行区域平均

软壳

onenuister 发表于 2020-4-6 10:08
第一个ave是进行侯平均，第二个aave是进行区域平均

楼主，你出来的grd文件，写了ctl吗？如何写的

onenuister

软壳 发表于 2020-4-8 19:57
楼主，你出来的grd文件，写了ctl吗？如何写的

dset f:\002\05\PWJ_fir_uwnd1979-2010.grd
title PWJ_fir
undef -9.96921e+36
xdef 1 linear 0 2.5
ydef 1 linear -90 2.5
zdef 1 levels 300
tdef 32 linear 00z01jul1979 1yr
vars 1
u 1 99
endvars

软壳

onenuister 发表于 2020-4-9 10:22
dset f:\002\05\PWJ_fir_uwnd1979-2010.grd
title PWJ_fir
undef -9.96921e+36

不知道楼主能不能帮我看看我昨晚发的帖子，我的问题差不多也是侯平均多年循环然后区域平均。
是1979年-2019年每年日降水资料（nc文件），需要提取青藏高原每年夏季90天的侯平均降水（18侯），一共738侯。然后再进行区域平均。图已经出来，但是感觉图不对，不知道问题出在哪里，哪位大神帮忙找出问题所在？
提取41年的夏季90天（18侯）的数据，生成dat文件：
'reinit'
'set fwrite E:\biyelunwenziliao\gs\pecip1979_2019.dat'
'set gxout fwrite'
year=1979
while(year<=2019)
'sdfopen E:\biyelunwenziliao\rizongliang\precip.'year'.nc'
'set lat 25 40'
'set lon 74 104'
'set z 1'
tt=152
while(tt<=241)
'set t 'tt''
'define pre=ave(precip,t=152,t=156)'
'd pre'
tt=tt+5
endwhile
'close 1'
year=year+1
endwhile
'disable fwrite'
;
有dat文件出来，我编写了ctl:
dset E:\biyelunwenziliao\gs\pecip1979_2019.dat
title CPC GLOBAL PRCP V1.0
undef -9.96921e+36
options yrev
xdef 61 linear 74.25 0.5
ydef 31 linear 24.75 0.5
zdef 1 linear 0 1
tdef 738 linear 00Z01JUN0001 1440mn
vars 1
pre=>pre  0  t,y,x  Daily total of precipitation
endvars
利用这个ctl画出青藏高原东部（我自己的研究区域）区域平均的图，我写的gs:
'reinit'
'open E:\biyelunwenziliao\gs\1979_2019dat.ctl'
'set grads off'
'set lon 74'
'set lat 25'
'set t 1 738'
'set z 1'
'set parea 1.0 10.5 0.8 7.75'
'set gxout line'
'set ccolor 8'
'set cmark 3'
'define p= tloop(aave(pre,lon=99,lon=103,lat=28,lat=35))'
'd p'
出来的图是738侯的图，可是看图很多数据都是接近0左右，考虑到是夏季，候平均以后的数据在区域平均，降水量不应该这样小，我就单独利用相同的方法处理某一年的数据画图，果然同样的年份，出来的图不一样，我现在不知道问题出在哪里。拜托各位了。

onenuister

软壳 发表于 2020-4-9 10:32
不知道楼主能不能帮我看看我昨晚发的帖子，我的问题差不多也是侯平均多年循环然后区域平均。
是1979年 ...

tt=152
while(tt<=241)
'set t 'tt''
'define pre=ave(precip,t=152,t=156)'
'd pre'
tt=tt+5
endwhile
你这里循环一直计算的都是一个值呀？