求助，ncl程序运行出现错误

islandowner

在做pop分析时遇到问题，跪求指点
首先贴出错误的提示


选择的sst数据信息如下

; ==============================================================
; POPanalysis: Calculate POPs
;              This is a user donated script.
;             *It is NOT supported software.
; ==============================================================
; Calculate POPs from SST data
; ==============================================================
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"

load "/home/gyy/pop/PrnOscPat_driver.ncl"    ; creates the POPs from EOFs
load "/home/gyy/pop/pltEofPop.ncl"      ; plot EOF/POP
load "/home/gyy/pop/pltClim.ncl"        ; plot climatology

; ==============================================================
;                      MAIN
; ==============================================================

  wcStrt = systemfunc("date")

  netCDF  = True                        ; create nc for POP
  dirNC   = "/home/gyy/pop/"                        ; dir for output nc

  pltDir  = "/home/gyy/pop/"                        ; dir for plots
  pltEOF  = True
  pltPOP  = True
  pltCLIM = True   

  pltType = "ps"
  pltName = "POP"

; ==============================================================
; User defined parameters that specify region of globe and time period
; ==============================================================

  latS   = 29.75                        ; spatial domain
  latN   = 29.75
  lonL   = 30.25
  lonR   = 119.75

  kPOP   = 1                                ; decide which POP we want to look at
  POPphase = (/"precursor","peak"/)        ; (/ 0, 1/)       

  yrStrt = 1958
  yrLast = 2008

  neof   = 10                           ; number of EOFs
  nPOP   = 2                                ; number of POPs

  optEOF = True      
  optETS = False

  REDUCE = False                         ; reduce spatial density
  nskip  =  2                            ; only when REDUCE=True
  mskip  =  2

; ==============================================================
; Open the file
; ==============================================================
  patho  = "/home/gyy/pop/"                                         ; path for output  

  pathi  = "/home/gyy/pop/"                                         ; path to source
  ifile  = "SODA_SST.nc"

  varname= "temp"                                        ; variable name
  f      = addfile (pathi+ifile, "r")                        ; entire file

; ==============================================================
; Read only the user specified period and region
; ==============================================================

  TIME   = f->time                                     ; ALL time on file
  YYYY   = cd_calendar(TIME,-1)/100                        ; ALL years
  iYYYY  = ind(YYYY.ge.yrStrt .and. YYYY.le.yrLast)         ; indices of year subset
  nTIME  = dimsizes(iYYYY)
                                                        ; read full year subset
  if (getfilevartypes(f, varname).eq."short") then      ; clarity
      DATA = short2flt(f->$varname$(iYYYY,{latS:latN},:))         
  else
      DATA = f->$varname$(iYYYY,{latS:latN},:)
  end if

  printVarSummary(DATA)                                    ; variable overview

; ==============================================================
; If dataset longitudes are -180=>180, flip to (nominally) span 0=>360   
; Only if both'lonL' and 'lonR' are positive.
; ==============================================================
  if (lonL.ge.0 .and. lonR.gt.0 .and. DATA&lon(0).lt.0) then
      DATA    = lonFlip(DATA)                                ; flip longitudes
     ;printVarSummary(DATA)  
  end if

; ==============================================================
; Subset to the time period (full years) and region of interest
; Optionally: reduce the grid point density for less memory.
; ==============================================================
  if (REDUCE) then
      data  = DATA(:, ::nskip,{lonL:lonR:mskip},:)
  else
      data  = DATA(:,:,{lonL:lonR},:)
  end if
  delete(DATA)                                         ; reduce memory

  printVarSummary(data)                                ; data to be used

; ==============================================================
; Remove annual cycle from the data: create anomalies
; .  use only full years to find annual cycle
; ==============================================================  
  clim  = clmMonTLL(data)                                ; calculcate climatology
;printVarSummary(clim)

  data         = calcMonAnomTLL(data, clim)                        ; remove annual cycle
;printVarSummary(data)                                      ; variable overview

  ntim  = dimsizes(data&time)                       

; =================================================================
; remove linear trend in time:  data(time,lat,lon): time is dimension # 0)
; =================================================================
  data = dtrend_n(data, False, 0)

; =================================================================
; normalize data at each gridpoint by local standard deviation at each grid pt
; =================================================================
  data = dim_standardize_n(data,1,0)

; =================================================================
; low pass (lp) filter
; remove short time scale variance (everything shorter than 18 months)
; =================================================================
  ihp   = 0                                                ; lowpass filter
  sigma = 1
  nWgt  = 37                                                ; loose 18 months each end
  fca   = 1./18.                                        ; 18 months
  wgt         = filwgts_lanczos(nWgt,ihp,fca,data@_FillValue,sigma)
  datalp= wgt_runave_n_Wrap(data,wgt,0,0)               ; running average (time)
  datalp@long_name = "Low Pass: "+data@long_name
  printVarSummary(datalp)                               ; (time,lat,lon)

; =================================================================
; Reorder (lat,lon,time) the input data so time is rightmost dimension.
; This is the order expected by the 'eofunc'
; =================================================================

  xlp   = datalp(lat|:,lon|:,time|:)   
  delete(datalp)                                        ; no longer needed
; =================================================================
; find EOFs: no weighting is performed here due to limited lat range
; =================================================================

  eof    = eofunc_Wrap(xlp, neof, optEOF)                      ; find first neof EOFs
  eof_ts = eofunc_ts_Wrap (xlp, eof, optETS)                ; time series

  printVarSummary( eof )                                ; examine EOF variables
  printVarSummary( eof_ts )

; +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
; A driver for the POP matrix operations.      
; +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

  pop_results  = PrnOscPat_driver(eof, eof_ts, kPOP)   ; returns a list
  print(pop_results)                                   ; containing variables

; =================================================================
; For convenience & clarity extract each of the two elements
;     and place them into separate variables reather than acces
;     via 'list' syntax.
; =================================================================

  pop_ts    = pop_results[0]     ; POP time series  (nEOF,time)
  printVarSummary(pop_ts)
;print("pop_ts: min="+min(pop_ts)+"   max="+max(pop_ts))

  pop_pat   = pop_results[1]     ; POP patterns     (nEOF,lat,lon)
  printVarSummary(pop_pat)
;print("pop_pat: min="+min(pop_pat)+"   max="+max(pop_pat))

  delete(pop_results)            ; no longer needed

; ===========================================================
; Miscellaneous for file creation and plots
; ===========================================================
  sfx    = get_file_suffix(ifile, 0)
  froot  = sfx@fBase

  yyyymm = cd_calendar(eof_ts&time,-2)/100  
  yrfrac = yyyymm_to_yyyyfrac(yyyymm, 0.0)   

  yyyymm!0 = "time"
  yrfrac!0 = "time"
  yyyymm&time = eof_ts&time
  yrfrac&time = eof_ts&time

if (netCDF) then
; ===========================================================
; NETCDF: write POP data to file: Use simple method of writing netCDF
; ===========================================================

  pathNC = dirNC+"POP_"+froot+".nc"
  system("/bin/rm -f " + pathNC)    ; remove file if exists
  fout = addfile(pathNC,"c")
  setfileoption(fout,"DefineMode",True)

; create global attributes
  fileAtt                = True
  fileAtt@creation_date = systemfunc("date")
  fileattdef(fout,fileAtt)

  fout->yyyymm      = yyyymm   
  fout->yrfrac      = yrfrac   
  fout->POP_TS      = pop_ts
  fout->POP_PATTERN = pop_pat

  fout->EOF_TS      = eof_ts
  fout->EOF         = eof
end if

;============================================================
; PLOTS
;============================================================

if (pltEOF) then
  pltEofPop(eof_ts, eof, yrfrac, froot, pltType, "EOF", "BlueYellowRed")
end if    ; pltEOF

if (pltPOP) then
  pltEofPop(pop_ts, pop_pat, yrfrac, froot, pltType, "POP", "BlueYellowRed")
end if    ; pltPOP

if (pltCLIM) then
    pltClim(clim, froot, pltType, "CLIM", "BlAqGrYeOrReVi200")
end if    ; pltCLIM


wallClockElapseTime(wcStrt, "POPanalysis", 0)

islandowner

好像是原始数据的维数和脚本中的DATE维数出现了问题，不知如何更改，不知哪位大神知道指导一下，谢谢

lanlaniris

这是你所有的错误信息吗？感觉上面还有的样子……82行上面什么变量是应该有4维的，而你的只有3维。应该是temp吧，有时空4维，但是你读数据只描述了3维

if (getfilevartypes(f, varname).eq."short") then      ; clarity
      DATA = short2flt(f->$varname$(iYYYY,{latS:latN},:))         
  else
      DATA = f->$varname$(iYYYY,{latS:latN},:)
  end if

winters

我也做过pop  印象当中 某个版本的ncl有bag，就是求数组转置的函数 参数有设反的bag，也可以检查下这个

islandowner

lanlaniris 发表于 2015-11-7 12:05
这是你所有的错误信息吗？感觉上面还有的样子……82行上面什么变量是应该有4维的，而你的只有3维。应该是te ...

谢谢回复！大概知道输数组维数的问题，把它的维数也改了很多次，改成四维的，但运行一直还是出错，最近刚学的ncl，大神能帮忙给一下怎么修改的意见么？谢谢

islandowner

winters 发表于 2015-11-7 17:44
我也做过pop  印象当中 某个版本的ncl有bag，就是求数组转置的函数 参数有设反的bag，也可以检查下这个

谢谢回复！楼主的ncl是什么版本的来着？楼主ncl小白，近期刚学的，大神可以介绍一下怎么检查么？谢谢

lanlaniris

islandowner 发表于 2015-11-7 19:04
谢谢回复！大概知道输数组维数的问题，把它的维数也改了很多次，改成四维的，但运行一直还是出错，最近刚 ...

你怎么改成四维的？改完以后的报错信息是什么
你之后的错误应该都是这一条引起的，这个改了，应该一系列都没问题了。

lanlaniris

islandowner 发表于 2015-11-7 19:04
谢谢回复！大概知道输数组维数的问题，把它的维数也改了很多次，改成四维的，但运行一直还是出错，最近刚 ...

建议你把temp变量的大致信息打出来看一下，有可能你的维数顺序也有错。

winters

islandowner 发表于 2015-11-7 19:06
谢谢回复！楼主的ncl是什么版本的来着？楼主ncl小白，近期刚学的，大神可以介绍一下怎么检查么？谢谢

不是大神，也是刚入门小白一只，运行的时候就只敲“ncl”回车，就会看到版本啦，我用的是现在比较新的版本6.3.0。如果某个函数在那个版本有bag 他的网页上通常会说名的，查看函数的时候上点心就行啦~ 之前好像是6.2.1还是6.2.0的那个版本的时候出的问题

islandowner

可以留个联系方式方便请教么？qq或者微信都可以，谢谢！