varname= "wstar"

  nmos   = 12
  yrStrt = 1980
  yrLast = 2010
  nyrs   = yrLast-yrStrt+1
  year   = ispan(yrStrt,yrLast,1)
  year!0 = "year"

  LATN   =  80.0     ; avoid regions 'loaded' with missing values
  LATS   = -60.0

  ntStrt = 0
  ntLast = nmos*nyrs-1

 ;diri   = "~/Dropbox/PhD/ccmi_wdir/ccmi/pp_wstar/refC1SD/" ;input directory
  diri   = "./refC1SD/" ;input directory
  fili   = systemfunc("cd "+diri+" ; ls *refC1SD_*.nc")
  nfili  = dimsizes(fili)
  print(fili)
  
  do nf = 0, nfili-1
     ff1    = addfile(diri+fili(nf), "r")
     wstar  = ff1->$varname$(ntStrt:ntLast,:,{LATN:LATS},0) ; 'lon' is degenerate; ignore
                                                            ; ff1->$varname$(ntStrt:ntLast,:,:,:)
     wstar  = wstar*1000    ;convert to mm/s
     wstar@units = "mm/s"

     print("-------------------")
     print("nf="+nf+"  "+fili(nf))
     printMinMax(wstar,0)
     print("") 

     ANN = month_to_annual(wstar,1) ; opt=1 to get the annual mean value(s)..
    ;printVarSummary(ANN)   ; [year | nyrs] x [lev | 31] x [lat | 64]
    ;printMinMax(ANN,0  ) 
     
     if (nf.eq.0) then
         dimANN = dimsizes(ANN)
         klev   = dimANN(1)
         nlat   = dimANN(2)
         WS     = new((/nfili,nyrs,klev,nlat/),typeof(wstar), getVarFillValue(wstar))
     end if
     
     WS(nf,:,:,:) = ANN     ; meta data transferred          
  end do   ; nf

; print locations of grid points with all missing values.

;;do kl=0,klev-1
;;  do nl=0,nlat-1
;;     if (all(ismissing(WS(0,:,kl,nl)))) then
;;         print(kl+" "+WS&lat(nl))
;;     end if
;;  end do
;;end do

  nmsg = num(ismissing(WS))
  print("nmsg="+nmsg)
  print("-------------------")

  WS!0    = "model_id"
  WS&year = year
  WS@file_names = fili
  printVarSummary(WS)         ;  [file_id | 6] x [year | 31] x [lev | 31] x [lat | 64] 
  printMinMax(WS, 0)
  print("-------------------")
  
  WSS = dim_standardize_n_Wrap(WS, 0, 1)   ; opt=0 (sample std. dev; 1= year [time])
  WSS@long_name = WS@long_name+": standardized"
  printVarSummary(WSS)                     ;  [file_id | 6] x [year | 31] x [lev | 31] x [lat | 64] 
  printMinMax(WSS, 0)
  print("-------------------")

;***************************************************************
; Read from the regressors ascii files 
;***************************************************************
  
  nmonths = nmos
  year0   = 1960
  year1   = 2010
  nyear   = year1-year0+1
  nind    = 7

  YEAR    = ispan(year0, year1, 1)
  YEAR!0  = "YEAR"

;;diri2 = "~/Dropbox/PhD/ccmi_wdir/ccmi/MLR_basis/regressors/refC1SD/" ;input directory
  diri2 = "./regressors/" ;input directory
  fili2 = systemfunc("cd "+diri2+" ; ls * ")
  nfili2 = dimsizes(fili2)
  print("nfili2="+nfili2)
  print(fili2)
  print("-------------------")
  
  reg = new((/nfili2,nyear,nind/),float)

  do nf2 = 0, nfili2 -1
     data = asciiread(diri2+fili2(nf2), -1, "float")
     y1   = onedtond(data,(/nind,nyear*nmonths/))
     
     y1!0 = "r_index"
     y1!1 = "time"
     y    = y1(time|:,r_index|:) 
     
     ann   = month_to_annual(y,1)
     alpha = dim_standardize_n_Wrap(ann, 0, 0)  ; [year | 51] x [r_index | 7]
     
     reg(nf2,:,:) = alpha
  end do

  reg!0 = "model"
  reg!1 = "YEAR"
  reg!2 = "r_index"
  reg&YEAR = YEAR
  
  printVarSummary(reg)
  printMinMax(reg,0)
  print("-------------------")
  
  NIES1      = new ( (/5,nyrs/), "float" )  
  NIES1(0,:) = reg(0,{yrStrt:yrLast},0)
  NIES1(1,:) = reg(0,{yrStrt:yrLast},1)
  NIES1(2,:) = reg(0,{yrStrt:yrLast},2)
  NIES1(3,:) = reg(0,{yrStrt:yrLast},3)
  NIES1(4,:) = reg(0,{yrStrt:yrLast},4)
  
  NIES1!0 = "reg_index"
  NIES2   = NIES1(YEAR|:,reg_index|:)
  
  NIES2@_FillValue = WSS@_FillValue
  printVarSummary(NIES2)
  print("-------------------")
  
  nies = WSS(0,:,:,:)
  printVarSummary(nies)          ; [year | 31] x [lev | 31] x [lat | 51]
  print("-------------------")   ;       0            1            2       
  
  nmsg_nies = num(ismissing(nies))
  print("nmsg_nies="+nmsg_nies)
  print("-------------------")
  
  if (nmsg_nies.eq.0) then
      print("No missing data: No linear interpolations performed")
      print("-------------------")
  end if
  
  printVarSummary(NIES2)      ;  [YEAR | 31] x [reg_index | 5]
  
  nf = 0
  do gg = 0,klev-1
      do hh = 0,nlat-1
         ;test = reg_multlin_stats(nies3(:,gg,hh),NIES2,0)
          test = reg_multlin_stats(WSS(nf,:,gg,hh),NIES2,0)
          printVarSummary(test)
      end do
  end do