标准差场，相关系数及相关系数场的程序示例

Aires

把脚本贴出来，供大家参考~有问题多交流~
数据文件说明：
dset  d:\work\oisst\%y4%m2.dat
options template
undef  -999.9
title SST monthly for 1981.12-2003.02 derived with the revised OI version2 analysis
xdef 360 linear  0.5 1
ydef 180 linear  -89.5 1
zdef 1 linear  1 1
tdef 255 linear Dec1981 1mon
vars 1
sst     1   99   clim
endvars

一、画1981年~2003年冬季（12月~2月）标准差场
************提取冬季（12月~2月平均）逐年SST场**************
'reinit'
'open d:\work\oisst\oisstmon.ctl'
'set fwrite d:\work\sst_winter_1981_2002.grd'
'set gxout fwrite'
'set x 1 360'
'set y 1 180'
'set z 1'
'set t 1'
i=1
while(i<255)
'd ave(sst,t='i',t='i+2')'
i=i+12
endwhile
'disable fwrite'
****************画出冬季SST标准差场******************
'reinit'
'open d:\work\sst_winter_1981_2002.ctl'
'set grads off'
'set grid off'
'set map 1 1 10'
'set lon 85 135'
'set lat 0 25'
'set z 1'
'set t 1 22'
'define anom=sst-ave(sst,t=1,t=22)'
'set xlopts 1 4 0.15'
'set ylopts 1 4 0.15'
'set t 1'
'set gxout shaded'
'set clevs 0.6'
'set rbcols 0 15'
'd sqrt(ave(anom*anom,t=1,t=22))'
'set gxout contour'
'set clskip 2'
'set clopts 1 0 0.15'
'd sqrt(ave(anom*anom,t=1,t=22))'
'printim d:\work\work4_rms.png white'
;

二、画Nino3指数（逐年12月）与SCS冷舌指数（冬季）（两个序列标准化的）的时间序列，并求相关系数（程序中略）
************冬季SST场做平均，输出22个SST场************
'reinit'
'open d:\work\oisst\oisstmon.ctl'
'set fwrite d:\work\scs.grd'
'set gxout fwrite'
'set x 1 360'
'set y 1 180'
'set z 1'
'set t 1'
i=1
while(i<=255)
'd ave(sst,t='i',t='i+2')'
i=i+12
endwhile
'disable fwrite'
*************提取逐年12月份SST场*****************
'reinit'
'open d:\work\oisst\oisstmon.ctl'
'set fwrite d:\work\nino.grd'
'set gxout fwrite'
'set x 1 360'
'set y 1 180'
'set z 1'
i=1
while(i<=255)
'set t 'i
'd sst'
i=i+12
endwhile
'disable fwrite'
************求SCS冷舌指数，Nino3指数，求相关系数并绘图***********
'reinit'
'open d:\work\scs.ctl'
'open d:\work\nino.ctl'
'set grads off'
'set grid off'
'set x 1'
'set y 1'
'set z 1'
'set t 1 22'
'define scs=tloop(aave(sst.1,lon=106,lon=111,lat=5,lat=10))'
'define nino=tloop(aave(sst.2,lon=-150,lon=-90,lat=-5,lat=5))'
'set xlopts 1 4 0.15'
'set ylopts 1 4 0.15'
'set ccolor 1'
'set cstyle 1'
'set cmark 2'
'set cthick 6'
'd (scs-ave(scs,t=1,t=22))/0.58'
'set ccolor 1'
'set cstyle 3'
'set cmark 0'
'set cthick 6'
'd (nino-ave(nino,t=1,t=22))/1.27'
'set strsiz 0.15'
'set string 1 c 6'
'draw string 2.54 6.87 Nino3'
'draw string 8.69 1.34 SCS cold tongue'
'set strsiz 0.2'
'draw string 6.08 7.29 R=0.727'
'printim d:\work\corrcoef.png white'
;

三、作逐年12月的Nino3指数与逐年SST场的相关系数分布场
*************提取逐年12月份的SST数据****************
'reinit'
'open d:\work\oisst\oisstmon.ctl'
'set fwrite d:\work\sst_dec.grd'
'set gxout fwrite'
'set x 1 360'
'set y 1 180'
'set z 1'
i=1
while(i<=255)
'set t 'i
'd sst'
i=i+12
endwhile
'disable fwrite'
*********画Nino3指数与SST的相关系数分布场*************
'reinit'
'open d:\work\sst_dec.ctl'
'set grads off'
'set grid off'
'set x 1'
'set y 1'
'set z 1'
'set t 1 22'
'define mean=tloop(aave(sst,lon=-150,lon=-90,lat=-5,lat=5))'
'define nino3=(mean-ave(mean,t=1,t=22))/1.27'
'set lon 50 290'
'set lat -40 40'
'set z 1'
'set t 1'
'define corr=tcorr(nino3,sst,t=1,t=22)'
'set xlopts 1 4 0.15'
'set ylopts 1 4 0.15'
'set gxout shaded'
'set clopts 1 0 0.15'
'set clevs 0.6'
'set ccols 0 15'
'd corr'
'set gxout contour'
'set clevs -0.4 -0.6'
'set cstyle 3'
'set ccolor 1'
'd corr'
'set clevs 0.4 0.6 0.8'
'set cstyle 1'
'set ccolor 1'
'd corr'
'printim d:\work\corrcoef.png white'
;

风儿飘飘

真心不错，感谢楼主的分享，努力学习中

风荷1

虽然不懂，也支持。

river

楼主发飙了

做个霸气的木头

支持下小A

Aires

做个霸气的木头 发表于 2013-3-31 09:04

                               
登录/注册后可看大图

支持下小A

我只是放在这以后备用，呵呵

做个霸气的木头

Aires 发表于 2013-3-31 09:07

                               
登录/注册后可看大图

我只是放在这以后备用，呵呵

顺便方便大家了，是吧

小小人

谢谢楼主，强大。

Aires

做个霸气的木头 发表于 2013-3-31 09:14

                               
登录/注册后可看大图

顺便方便大家了，是吧

主要是方便我们同学了，呵呵，这是作业

liuzf406

i=1
while(i<255)
'd ave(sst,t='i',t='i+2')'
i=i+12

这是指的jfm平均值吧？为何说是冬季（12月~2月）的平均值呢？可以解释一下吗？
谢谢！

Aires

liuzf406 发表于 2013-3-31 11:50

                               
登录/注册后可看大图

i=1
while(i

数据的第一时次就是1981年12月，所以我先把冬季做的平均才做的