前回つくった信頼区間つきでスペクトル密度をもとめる関数を、NINO3監視指数に適用する。NINO3の読み込みは以前に作ったread_NINOモジュールを用いる。
下の例ではtaperにはHammingフィルタ、周波数空間でのスムージングに1-2-1フィルタ、信頼度は90%を用いた。
黒線がスペクトル密度、赤点線が信頼区間を示す。縦の黒点線は2年と7年に対応する補助線である。
from read_NINO import read_NINO import matplotlib.mlab as mlab import matplotlib.pyplot as plt import numpy as np from psd_withError import psd_withError #import psd_withError module data=read_NINO() x=data["NINO3 ANOM"] # case1 Pxx,f,upper,lower=psd_withError(x, NFFT=len(x),Fs=12,window=np.hamming(len(x)),smooth=[1,2,1],Confidencelevel=0.90) plt.figure() plt.loglog(f,Pxx,'k-o') plt.loglog(f,upper,'r--') plt.loglog(f,lower,'r--') plt.loglog([1./2.,1./2.],[1.E-5,1.E2],'k--') #vertical line for 2 years plt.loglog([1./7.,1./7.],[1.E-5,1.E2],'k--') #vertical line for 7 years plt.savefig("NINO_psd_withError.png") plt.show()
0 件のコメント:
コメントを投稿