- 积分
- 923
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2016-11-2
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
参考ncl官网小波分析示例2,程序如下结果运行出来显示Segmentation fault (核心已转储),求助各位大大,这是读取文件有问题呢,还是建立坐标部分有问题呢?
谢谢大家!!
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
begin
ninoseas = asciiread("F:/NCLxuexi/shuju/chjjl.txt",706,"float")
ninoseas!0 = "time"
ntime = dimsizes(ninoseas)
timeo = fspan(1300.705,2005.,ntime)
ninoseas&time = timeo
time = timeo
N = dimsizes(time)
;************************************
; compute wavelet
;************************************
mother = 0
param = 6.0
dt = 1 ;timesteps in units of years
s0 = 0.25
dj = 0.25
jtot = 1+floattointeger(((log10(N*dt/s0))/dj)/log10(2.))
npad = N
nadof = new(2,float)
noise = 1
siglvl = .05
isigtest= 0
w = wavelet(time,mother,dt,param,s0,dj,jtot,npad,noise,isigtest,siglvl,nadof)
;************************************
; create coodinate arrays for plot
;************************************
power = onedtond(w@power,(/jtot,N/))
power!0 = "period" ; Y axis
power&period = w@period ; convert period to units of years
power!1 = "time" ; X axis
power&time = time
power@long_name = "Power Spectrum"
power@units = "1/unit-freq"
gws = w@gws
; compute significance ( >= 1 is significant)
SIG = power ; transfer meta data
SIG = power/conform (power,w@signif,0)
SIG@long_name = "Significance"
SIG@units = " "
;plot=ShadeCOI(wks,plot(0),w,time,False);加检验线?
;********************************************************************************
; initial resource settings
;********************************************************************************
wks = gsn_open_wks("png","wavelet") ; send graphics to PNG file
res = True ; plot mods desired
res@gsnDraw = False ; Do not draw plot
res@gsnFrame = False ; Do not advance frome
res@cnFillOn = True ; turn on color
res@cnFillPalette = "BlAqGrYeOrReVi200" ; set color map
res@cnFillMode = "RasterFill" ; turn on raster mode
res@cnRasterSmoothingOn = True ; turn on raster smoothing
res@cnLinesOn = False ; turn off contour lines
res@cnLineLabelsOn = False
res@cnInfoLabelOn = False
res@trYReverse = True ; reverse y-axis
res@tmYLMode = "Explicit"
res@tmYLValues = (/1,2,4,8,16,32,64,128/)
res@tmYLLabels = (/"1","2","4","8","16","32","64","128"/)
res@tmLabelAutoStride = True
res@vpHeightF = .4 ;
res@vpWidthF = .7
res@cnLevelSelectionMode = "ExplicitLevels" ; set manual contour levels
res@cnLevels = (/0.5,1.,2.,4./)
res@gsnRightString = "Wavelet Power"
res@gsnLeftString = "NINO3: GISST"
res@tiYAxisString = "Period"
res2 = True ; res2 probability plots
res2@trYReverse = True
res2@tmYLMode = "Explicit"
res2@tmYLValues = (/1,2,4,8,16,32,64,128/)
res2@tmYLLabels = (/"1","2","4","8","16","32","64","128"/)
res2@gsnDraw = False ; Do not draw plot
res2@gsnFrame = False ; Do not advance frome
res2@cnLevelSelectionMode = "ManualLevels" ; set manual contour levels
res2@cnMinLevelValF = 0.00 ; set min contour level
res2@cnMaxLevelValF = 2.00 ; set max contour level
res2@cnLevelSpacingF = 1.00 ; set contour spacing
res2@cnInfoLabelOn = False
res2@cnLinesOn = False ; do not draw contour lines
res2@cnLineLabelsOn = False ; do not draw contour labels
res2@cnFillScaleF = 0.5 ; add extra density
res2@gsnLeftString = ""
res2@gsnRightString = ""
plot = new(2,graphic)
plot(0) = gsn_csm_contour(wks,power,res)
iplot = gsn_csm_contour(wks,SIG,res2)
opt = True
opt@gsnShadeFillType = "pattern"
opt@gsnShadeHigh = 17
iplot = gsn_contour_shade(iplot,0, 0.8, opt)
overlay(plot(0),iplot) ; overlay probability plot onto power plot
scale = w@scale
Cdelta = w@cdelta
powernorm = power
powernorm = power/conform(power,scale,0)
scaleavg = dj*dt/Cdelta*dim_sum_Wrap(powernorm(time|:,{period|2.:8.}))
gws = w@gws
resl = True
resl@gsnFrame = False
resl@gsnDraw = False
resl@trYAxisType = "LogAxis"
resl@trYReverse = True ; reverse y-axis
resl@tmYLMode = "Explicit"
resl@tmYLValues = (/1,2,4,8,16,32,64,128/)
resl@tmYLLabels = (/"1","2","4","8","16","32","64","128"/)
plotg = gsn_csm_xy(wks,gws,power&period,resl)
plotc = gsn_attach_plots(plot(0),plotg,res,resl)
ress = True
ress@xyDashPatterns = (/0,1/)
ress@pmLegendDisplayMode = "Always"
ress@pmLegendOrthogonalPosF = -1.1
ress@pmLegendParallelPosF = 0.2
ress@pmLegendWidthF = 0.15
ress@pmLegendHeightF = 0.1
ress@lgPerimOn = False ; turn off box around
ress@gsnFrame = False
ress@gsnDraw = False
ress@vpHeightF = .4
ress@vpWidthF = .7
ress@xyExplicitLegendLabels = (/"2-8 yr"/)
ress@tiYAxisString = "variance"
plot(1) = gsn_csm_xy(wks,power&time,scaleavg,ress)
pres = True
pres@gsnMaximize = True
pres@gsnPaperOrientation = "portrait"
gsn_panel(wks,plot,(/2,1/),pres)
end
|
|