爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

搜索
查看: 3041|回复: 1

脚本:夏季Nino3.4指数与同期降水相关分布

[复制链接]
发表于 2018-4-21 23:15:06 | 显示全部楼层 |阅读模式

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

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

x
这是原值相关,判断缺省值如何写:
#打开nino3.4数据文件
table = readtable('f:/nino34.csv', delimiter=',', format='%s%f')
yy = table['yyJJA']
nino34 = table['nino34']
#print nino34
#打开降水数据文件
f = addfile('f:/precip.mon.mean.2.5x2.5.nc')
#prep = f['precip'][0::12,'0:60','60:150']
var  = f['precip'][0,'0:60','60:150']
prepJJA = []
for i in range(1981,2016):
    tidx = (i-1948)*12
    t = f.gettime(tidx)
    data = f['precip'][tidx+5,'0:60','60:150']+\
           f['precip'][tidx+6,'0:60','60:150']+\
           f['precip'][tidx+7,'0:60','60:150']   
    prepJJA.append(data)

years = arange(1981,2016)
nx = var.dimlen(1)
ny = var.dimlen(0)
cor = zeros([ny,nx])
cor = dim_array(cor,var.dims)
for i in arange(0,nx):
    for j in arange(0,ny):
        x = []
        y = []
        for t in syears:
            tidx = t -1981
            if(prepJJA[tidx][j,i]>-9.96921E36):
                x.append(nino34[t-1981])
                y.append(prepJJA[tidx][j,i])
        #scatter(xx, yy, s=2, color='k', label='PM')
        if(len(x)>1):
            slope, intercept, r, p, std,number = linregress(x, y)
            cor[j,i] = r
        else:
            cor[j,i] = nan
#print cor
axesm()
mlayer = shaperead('F:\MeteoInfo\MeteoInfo\Map\country1.shp')
clevs =[-0.53,-0.43,-0.33,-0.28,0,0.28,0.33,0.43,0.53]
cols=[(86,49,5),(138,80,0),(186,145,42),(216,179,107),(245,236,212),\
    (217,238,235),(117,197,186),(46,145,137),(0,99,91),(0,59,47)]
#layer = contourfm(cor,20)  
layer = contourfm(cor,clevs,colors=cols)
colorbar(layer,orientation='horizontal')
xlim(60,150)
ylim(0,60)
title('correlation:nino34(JJA)-prep(JJA) 1981-2015')      

cor.png
密码修改失败请联系微信:mofangbao
 楼主| 发表于 2018-4-21 23:22:55 | 显示全部楼层
如果去掉了它们的线性趋势,遇到了问题,还是基础知识没搞懂。原图网址http://cmdp.ncc-cma.net/pred/cn_enso.php?product=cn_enso_corr&season=JJA#corr#打开nino3.4数据文件
table = readtable('f:/nino34.csv', delimiter=',', format='%s%f')
yy = table['yyJJA']
nino34 = table['nino34']
#print nino34
#打开降水数据文件
f = addfile('f:/precip.mon.mean.2.5x2.5.nc')
prep = f['precip'][0::12,'0:60','60:150']
var  = f['precip'][0,'0:60','60:150']
#print prep.dims
nx = prep.dimlen(2)
ny = prep.dimlen(1)
nt = prep.dimlen(0)
t = prep.dimvalue(0)
#print t
prepJJA = zeros([nt,ny,nx])
prepJJA = dim_array(prepJJA,prep.dims)
#prepJJA_ = zeros([nt,ny,nx])

#prepJJA_ = dim_array(prepJJA,prep.dims)
years = arange(1948,2018)
for i in years:
    j = (i-1948)
    tidx = j*12
    t = f.gettime(tidx)
    data = f['precip'][tidx+5,:,:]+f['precip'][tidx+6,:,:]+f['precip'][tidx+7,:,:]
    prepJJA[j,:,:]  = data

syears = arange(1981,2016)
for i in arange(0,nx):
    for j in arange(0,ny):
#for i in arange(0,1):
#    for j in arange(0,1):
        xx = []
        yy = []
        for t in syears:
            tidx = t -1948
            if(prepJJA[tidx,j,i]>-9.96921E36):
                xx.append(t)
                yy.append(prepJJA[tidx,j,i])
        #print xx
        #print yy
        #scatter(xx, yy, s=2, color='k', label='PM')
        if(len(xx)>1):
            slope, intercept, r, p, std,number = linregress(xx, yy)
            #print slope, intercept, r
            y = polyval([slope, intercept], syears)
        else:
            y = -9.96921E36
        for t in arange(0,2):
            prepJJA_[t,j,i] = prepJJA[t,j,i]-y[t]#减掉拟合值

最后报错
Traceback (most recent call last):
  File "F:\MeteoInfo\correlation.py", line 50, in <module>
    prepJJA_[t,j,i] = prepJJA[t,j,i]-y[t]
TypeError: 'float' object is unsubscriptable


密码修改失败请联系微信:mofangbao
您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

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

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

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