[ncl-talk] EOF PC timeseries

[Sri nandini]

Hello ncl-users,

I have an EOF script, which im applying to my data over the 
time+ensemble level. The EOF spatial pattern comes out correctly, 
however, the timeseries does not get plotted. Could someone help me out 
with my code below?

   begin

   latS   =  20.
   latN   =  68.
   lonL   = -78.
   lonR   =  15.

   neof   = 2       ; number of EOFs
   optEOF = True
   optEOF at jopt = 0   ;uses correlation matrix to compute EOFs. The 
default is to use a covariance matrix (optEOF at jopt = 0).
   optETS = False
   yrStrt                                   = 1986
   yrLast                                   = 2006

   f     = addfile ("data.nc", "r")
    time   = f->time
    YYYY   = cd_calendar(time,-1)/100
    iYYYY  = ind(YYYY.ge.yrStrt .and. YYYY.le.yrLast)
    slp=f->hist_trend(iYYYY,:,{latS:latN},{lonL:lonR})
    printVarSummary(slp)    ;[240,100,13,24]

   slp = dtrend_msg_n(ispan(0,dimsizes(slp&time)-1,1),slp,True,False,0) 
; detrend the data
; =================================================================
; create weights:
; =================================================================
   rad    = 4.*atan(1.)/180.
   lat = f->lat({latS:latN})
    if (typeof(lat).eq."double") then
        clat = sqrt( cos(rad*tofloat(lat)) )
    else
        clat = sqrt( cos(rad*lat) )
    end if
    copy_VarCoords(lat, clat) ; contributed
    printVarSummary(clat)
; =================================================================
; weight all observations
; =================================================================
   wSLP   = slp  ; type float
   wSLP   = slp*conform(slp, clat, 2)
   printVarSummary(wSLP)
; =================================================================
; Reorder (lat,lon,time) the *weighted* input data Access the area of 
interest via coordinate subscripting
; =================================================================

   work := reshape(wSLP ,(/ntim*nens,nlat,mlon/))    ;concatanate the 
time+ensemeble dimension together to get total variance
   copy_VarCoords(wSLP(0,0,:,:), work(0,:,:))      ; trick to copy coord 
information
   printVarSummary( work )

   eof    = eofunc_n_Wrap(work,neof,optEOF,0 )
   eof_ts = eofunc_ts_n_Wrap(work, eof, optETS,0)
   printVarSummary( eof )                         ; should be [3,13,24]
   printVarSummary( eof_ts )                      ; should be [3,2400]

    eof_ts = dim_standardize_n( eof_ts, 0, 1)
   ;printVarSummary(eof_ts)

;====================================================================
;Regressnig data onto standardized PC1,2
;======================================================================

   eof_regres=eof ;create an array with meta data
   printVarSummary(eof_regres)

   do ne=0,neof-1
    eof_regres(ne,:)=(/ regCoef_n(eof_ts(ne,:), work,0,0) /)
     end do
   printVarSummary(eof_regres)

   plotfile="EOF_figure"
   wks_type = "png"
   wks = gsn_open_wks(wks_type,plotfile)
   plot                    = new(neof,graphic) ; create graphic array; 
only needed if paneling
   res                     = True
   res at mpProjection        = "LambertConformal"; choose projection
   res at gsnDraw             = False           ; don't draw
   res at gsnFrame            = False           ; don't advance frame

   res at cnFillOn            = True            ; turn on color
   res at cnLinesOn   = False                       ; turn off the contour 
lines
   res at cnLineThicknessF = 0.9
   res at cnLineColor = 0
   res at mpDataBaseVersion   = "MediumRes"
   res at cnLineLabelsOn      = False    ; turn off contour line labels

   res at cnFillDrawOrder     = "PreDraw"       ; draw contours before 
continents
   res at cnFillPalette = "BlRe"
   res at mpMinLatF            = latS           ; zoom in on map
   res at mpMaxLatF            = latN
   res at mpMinLonF            = lonL
   res at mpMaxLonF            = lonR

   res at mpOutlineOn = True                        ; turn the map outline on

   res at cnLevelSelectionMode = "ManualLevels"     ; manual set levels
   res at cnMinLevelValF       = -0.4
   res at cnMaxLevelValF       = .4
   res at cnLevelSpacingF      = 0.05                 ; 20 contour levels
res at lbLabelBarOn         = False              ; turn off individual lb's
   res at lbBoxEndCapStyle     = "TriangleBothEnds" ; Added in NCL V6.4.0
   res at lbLabelAutoStride    = True               ; Control labelbar spacing
   res at gsnMaximize          = True               ; large format in landscape
   res at gsnAddCyclic         = False              ; plotted dataa are not 
cyclic
   res at gsnMaskLambertConformal = True            ; turn on lc masking

;=============================panel plot only resources
   resP                     = True         ; modify the panel plot
   resP at gsnMaximize         = True         ; large format
   resP at gsnPanelLabelBar    = True         ; add common colorbar
   resP at lbLabelFontHeightF = 0.02
   res at gsnCenterStringFontHeightF = 0.02              ; ditto
;====================Create (but don't draw) both plots

      do n=0,neof-1
      res at gsnStringFontHeightF = 0.02
      res at gsnLeftString  = "EOF "+(n+1)
      res at gsnCenterString= "a) Hist SSH"
      res at gsnRightString = sprintf("%5.1f", eof at pcvar(n)) +"%"
      plot(n)=gsn_csm_contour_map(wks,eof_regres(n,:,:),res)
    end do
   gsn_panel(wks,plot,(/neof,1/),resP)     ; now draw as one plot

;============================================================
; time series (principal component) plot
;============================================================
    eof_ts_ens_avg at long_name = "Standardized"

    yyyymm = cd_calendar(slp&time, -1)         ; create a 
yyyy.fraction_of_year
    eof_ts_ens_avg&time = yyyymm_to_yyyyfrac(yyyymm, 0.0) ; better time axis
    ;eof_ts&time at long_name = "YEAR"
    ntim   = dimsizes(yyyymm)

    ;printVarSummary(frac)
    rts           = True
    rts at gsnDraw   = False       ; don't draw yet
    rts at gsnFrame  = False       ; don't advance frame yet
    rts at tmXBFormat= "f"         ; no unnessary decimal pts

    rts at vpHeightF = 0.40        ; Changes the aspect ratio
    rts at vpWidthF  = 0.85
    rts at vpXF      = 0.10        ; change start locations
    rts at vpYF      = 0.75        ; the plo
    rts at tiYAxisString = "Standardized"          ; y-axis label
    rts at trXMinF = 1986
    rts at trXMaxF = 2006
    rts at trYMinF = -3
    rts at trYMaxF =  3
    rts at gsnYRefLine           = 0.              ; reference line
    rts at gsnXYBarChart         = True            ; create bar chart
    rts at gsnAboveYRefLineColor = "red"           ; above ref line fill red
    rts at gsnBelowYRefLineColor = "blue"          ; below ref line fill blue

  ; panel plot only resources
    rtsP                      = True            ; modify the panel plot
    rtsP at gsnMaximize          = True            ; large format
  ; create individual plots
   do n=0,neof-1
      rts at gsnLeftString  = "EOF "+(n+1)
      rts at gsnRightString = sprintf("%5.1f", eof_regres at pcvar(n)) +"%"
      plot(n) = gsn_csm_xy (wks,work(:,0,0),eof_ts(n,:),rts)
   end do
   gsn_panel(wks,plot,(/neof,1/),rtsP)     ; now draw as one plot
end


-- 
Dr. Sri, Nandini-Weiss
Research Associate

Universität Hamburg
Center for Earth System Research and Sustainability (CEN)
Cluster of Excellence 'Climate, Climatic Change, and Society' (CLICCS)

Bundesstrasse 53, 20146 Hamburg