运行combined eofs时为什么eof_cdata前面的数全缺省

银魈

报的错

显示这个错误

错误行

错误行

错误行周围

周围

eof_cdata(0:1,0:359)全是缺省，而其他值正常

发现eof_cdata(0:1,0:359)全是缺省，而其他值正常

错误行以前的程序：
undef("read_rename")
function read_rename(f[1]:file, varName[1]:string       \
                    ,iStrt[1]:integer, iLast[1]:integer \
                    ,latS[1]:numeric , latN[1]:numeric  )
; Utility to force specific named dimensions
; This is done for historical reasons (convenience)
begin
   work    = f->$varName$(iStrt:iLast,{latS:latN},:)   ; (time,lat,lon)
   work!0  = "time"                                    ; CAM model names
   work!1  = "lat"
   work!2  = "lon"
   return(work)
end
; =========================>  MAIN  <==============================
begin
   neof    =  2

   latS    = -15
   latN    =  15

   ymdStrt = 20000101                         ; start yyyymmdd
   ymdLast = 20181231                         ; last  

   yrStrt  = ymdStrt/10000
   yrLast  = ymdLast/10000

   pltDir  = "./"                             ; plot directory
   pltType = "png"
   pltName = "mjoclivar"   

   diri    = "/public/home/liuxc/draw4/data/"                             ; input directory

   filpre  = "precip.day.anomalies.2000-2018.nc"
   filu200 = "uwnd.day.200.anomalies.2000-2018.nc"
   filu850 = "uwnd.day.850.anomalies.2000-2018.nc"

;************************************************
; create BandPass Filter
;************************************************
  ihp      = 2                             ; bpf=>band pass filter
  nWgt     = 201
  sigma    = 1.0                           ; Lanczos sigma
  fca      = 1./100.
  fcb      = 1./20.
  wgt      = filwgts_lanczos (nWgt, ihp, fca, fcb, sigma )

;***********************************************************
; Find the indices corresponding to the start/end times
;***********************************************************

   f       = addfile (diri+filpre , "r")                        
   TIME    = f->time                          ; days since ...

   YMD     = cd_calendar(TIME, -2)            ; entire (time,6)

   iStrt   = ind(YMD.eq.ymdStrt)              ; index start
   iLast   = ind(YMD.eq.ymdLast)              ; index last
   delete([/ TIME, YMD /])

;***********************************************************
; Read anomalies
;***********************************************************

   work    = read_rename(f,"precip",iStrt,iLast,latS,latN) ; (time,lat,lon)
   PRE     = dim_avg_n_Wrap(work, 1)                         ; (time,lon)

   f       = addfile (diri+filu850 , "r")                        
   work    = read_rename(f,"U850",iStrt,iLast,latS,latN) ; (time,lat,lon)
   U850    = dim_avg_n_Wrap(work, 1)          ; (time,lon)

   f       = addfile (diri+filu200 , "r")                        
   work    = read_rename(f,"U200",iStrt,iLast,latS,latN) ; (time,lat,lon)
   U200    = dim_avg_n_Wrap(work, 1)          ; (time,lon)
   PRE@__FillValue = U850@_FillValue
   dimw    = dimsizes( work )
   ntim    = dimw(0)
   nlat    = dimw(1)
   mlon    = dimw(2)
   delete(work)

   lon     = PRE&lon                                          
   time    = PRE&time                        
   date    = cd_calendar(time, -2)            ; yyyymmdd

;************************************************
; Apply the band pass filter to the original anomalies
;************************************************
   pre   = wgt_runave_n_Wrap ( PRE, wgt, 0, 0) ; (time,lon)
   u850  = wgt_runave_n_Wrap (U850, wgt, 0, 0)
   u200  = wgt_runave_n_Wrap (U200, wgt, 0, 0)

;************************************************
; remove temporal means of band pass series: *not* necessary
;************************************************
   pre   = dim_rmvmean_n( pre, 0)              ; (time,lon)
   u850  = dim_rmvmean_n(u850, 0)
   u200  = dim_rmvmean_n(u200, 0)

;************************************************
; Compute the temporal variance at each lon
;************************************************
   var_pre  = dim_variance_n_Wrap( pre, 0)     ; (lon)
   var_u850 = dim_variance_n_Wrap(u850, 0)
   var_u200 = dim_variance_n_Wrap(u200, 0)

;************************************************
; Compute the zonal mean of the temporal variance
;************************************************
  zavg_var_pre  = dim_avg_n_Wrap( var_pre , 0)     
  zavg_var_u850 = dim_avg_n_Wrap( var_u850, 0)
  zavg_var_u200 = dim_avg_n_Wrap( var_u200, 0)

;************************************************
; Normalize by sqrt(avg_var*)
;************************************************
  pre   =  pre/sqrt(zavg_var_pre )          ; (time,lon)
  u850  = u850/sqrt(zavg_var_u850)
  u200  = u200/sqrt(zavg_var_u200)

;************************************************
; Combine the normalized data into one variable
;************************************************
  cdata     = new ( (/3*mlon,ntim/), typeof(pre), getFillValue(pre))
  do ml=0,mlon-1
     cdata(ml       ,:) = (/  pre(:,ml) /)
     cdata(ml+  mlon,:) = (/ u850(:,ml) /)
     cdata(ml+2*mlon,:) = (/ u200(:,ml) /)
  end do

;************************************************
; Compute **combined** EOF; Sign of EOF is arbitrary
;************************************************
  eof_cdata    = eofunc(cdata   , neof, False)      ; (neof,3*mlon)
  print("==============")
  printVarSummary(eof_cdata)
  printMinMax(eof_cdata, True)
  print(eof_cdata(0,0:400))
  eof_ts_cdata = eofunc_ts(cdata,eof_cdata,False)   ; (neof,3*ntim)
  print("==============")                                 
  printVarSummary(eof_ts_cdata)
  printMinMax(eof_ts_cdata, True)

;************************************************
; For clarity, explicitly extract each variable. Create time series
;************************************************

  nvar = 3  ; "olr", "u850", "u200"
  ceof = new( (/nvar,neof,mlon/), typeof(cdata), getFillValue(cdata))

  do n=0,neof-1
     ceof(0,n,:) = eof_cdata(n,0:mlon-1)      ; olr
     ceof(1,n,:) = eof_cdata(n,mlon:2*mlon-1) ; u850
     ceof(2,n,:) = eof_cdata(n,2*mlon:)       ; u200
  end do

  ceof!0   = "var"
  ceof!1   = "eof"
  ceof!2   = "lon"   
  ceof&lon =  pre&lon

  ceof_ts        = new( (/nvar,neof,ntim/), typeof(cdata), getFillValue(cdata))
  ceof_ts(0,:,:) = eofunc_ts_Wrap( pre(lon|:,time|:),ceof(0,:,:),False)   ; (0,neof,ntim)
  ceof_ts(1,:,:) = eofunc_ts_Wrap(u850(lon|:,time|:),ceof(1,:,:),False)   ; (1,neof,ntim)
  ceof_ts(2,:,:) = eofunc_ts_Wrap(u200(lon|:,time|:),ceof(2,:,:),False)   ; (2,neof,ntim)

;**********************************************t*
; Add code contributed by Marcus N. Morgan, Florida Institute of Technology; Feb 2015
; Calculate % variance (pcv_ )accounted for by OLR, U850 and U200
;************************************************

    pcv_eof_pre  = new(neof,typeof(ceof))
    pcv_eof_u850 = new(neof,typeof(ceof))
    pcv_eof_u200 = new(neof,typeof(ceof))

    do n=0,neof-1
       pcv_eof_pre(n)  = avg((ceof(0,n,:)*sqrt(ceof@eval(n)))^2)*100
       pcv_eof_u850(n) = avg((ceof(1,n,:)*sqrt(ceof@eval(n)))^2)*100
       pcv_eof_u200(n) = avg((ceof(2,n,:)*sqrt(ceof@eval(n)))^2)*100
     ;;print("pcv: neof="+(n+1)+":  "+pcv_eof_olr(n)+"  "+pcv_eof_u850(n)+"  "+pcv_eof_u200(n))
    end do

;************************************************
; Change sign of EOFs for spatial structures of PC1 and PC2
; to represent convection over the tropical Indian Ocean and the tropical western Pacific Ocean, respectively
; (Ad hoc approach)
;************************************************

  imax_pre_eof1   = maxind(ceof(0,0,:)) ;括号里俩数组全员缺省        
  imax_pre_eof2   = maxind(ceof(0,1,:))   
  ;print(ceof(0,1,:))
  lonmax_eof1 = ceof&lon(imax_pre_eof1)      ; longitude of max value (i.e. suppressed convection)

银魈

我造哪儿有问题了，降水的原始数据里面有几天缺测值和别的时间不一样，导致这些缺测值被当成数据读进去了