; ============================================================== ; eof_1.ncl ; ; Concepts illustrated: ; - Calculating EOFs ; - Using coordinate subscripting to read a specified geographical region ; - Rearranging longitude data to span -180 to 180 ; - Calculating symmetric contour intervals ; - Drawing filled bars above and below a given reference line ; - Drawing subtitles at the top of a plot ; ; ============================================================== ; Calculate EOFs of the Sea Level Pressure over the North Atlantic. ; ============================================================== ; The slp.mon.mean file can be downloaded from: ; http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.surface.html ; ============================================================== 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 ; ============================================================== ; User defined parameters that specify region of globe and ; ============================================================== latS = -10. latN = 10. lonL = -180. lonR = 180. levD = 0. levU = 19. yrStrt = 1998 yrLast = 2012 season = "JJA" ; choose Dec-Jan-Feb seasonal mean neof = 3 ; number of EOFs optEOF = True optEOF@jopt = 0 ; This is the default; most commonly used; no need to specify. ;;optEOF@jopt = 1 ; **only** if the correlation EOF is desired optETS = False ; ============================================================== ; Open the file: Read only the user specified period ; ============================================================== f = addfile ("JJA_test.nc", "r") TIME = f->time printVarSummary(TIME) ; variable overview YYYY = cd_calendar(TIME,-1)/100 ; entire file iYYYY = ind(YYYY.ge.yrStrt .and. YYYY.le.yrLast) ; LH = f->r(iYYYY,:) LH = f->r ; lat = f->lat printVarSummary(LH) ; variable overview ; printVarSummary(lat) ; ============================================================== ; dataset longitudes span 0=>357.5 ; Because EOFs of the North Atlantic Oscillation are desired ; use the "lonFlip" (contributed.ncl) to reorder ; longitudes to span -180 to 177.5: facilitate coordinate subscripting ; ============================================================== ; LH = lonFlip( LH ) ; printVarSummary(sst) ; note the longitude coord ; ============================================================== ; compute desired global seasonal mean: month_to_season (contributed.ncl) ; ============================================================== ; SST = month_to_season (sst, season) ; nyrs = dimsizes(SST&time) ; printVarSummary(SST) ; ================================================================= ; create weights: sqrt(cos(lat)) [or sqrt(gw) ] ; ================================================================= ; rad = 4.*atan(1.)/180. ; lat = f->lat ; 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) ; clat = f->lat ; clat =tofloat(f->lat) ; clat = sqrt( cos(rad*clat) ) ; gw for gaussian grid ; printVarSummary(lat) ; printVarSummary(clat) ; ================================================================= ; weight all observations ; ================================================================= ; wLH = LH ; copy meta data ; wLH = LH*conform(LH, clat, 1) ; wLH@long_name = "Wgt: "+wLH@long_name ; ================================================================= ; Reorder (lat,lon,time) the *weighted* input data ; Access the area of interest via coordinate subscripting ; ================================================================= ; x = wLH({lat|latS:latN},{lon|lonL:lonR},time|:) ; x = wLH({lev|levD:levU},time|:) eof = eofunc_Wrap(LH, neof, optEOF) eof_ts = eofunc_ts_Wrap (LH, eof, optETS) printVarSummary( eof ) ; examine EOF variables printVarSummary( eof_ts ) ; ================================================================= ; Normalize time series: Sum spatial weights over the area of used ; ================================================================= ; dimx = dimsizes( x ) ; mln = dimx(1) ; sumWgt = mln*sum( clat({lat|latS:latN}) ) ; eof_ts = eof_ts/sumWgt ; ================================================================= ; Extract the YYYYMM from the time coordinate ; associated with eof_ts [same as SLP&time] ; ================================================================= ; yyyymm = cd_calendar(eof_ts&time,-2)/100 ;;yrfrac = yyyymm_to_yyyyfrac(yyyymm, 0.0); not used here ;============================================================ ; PLOTS ;============================================================ wks = gsn_open_wks("pdf","JJA_test") gsn_define_colormap(wks,"BlWhRe") ; choose colormap plot = new(neof,graphic) ; create graphic array ; only needed if paneling ; EOF patterns res = True res@gsnDraw = False ; don't draw yet res@gsnFrame = False ; don't advance frame yet ;---This resource not needed in V6.1.0 res@gsnSpreadColors = True ; spread out color table res@gsnAddCyclic = False ; plotted dataa are not cyclic res@mpFillOn = False ; turn off map fill res@mpMinLatF = latS ; zoom in on map res@mpMaxLatF = latN res@mpMinLonF = lonL res@mpMaxLonF = lonR res@cnFillOn = True ; turn on color fill res@cnLinesOn = False ; True is default ;res@cnLineLabelsOn = False ; True is default res@lbLabelBarOn = False ; turn off individual lb's ; set symmetric plot min/max symMinMaxPlt(eof, 16, False, res) ; contributed.ncl ; panel plot only resources resP = True ; modify the panel plot resP@gsnMaximize = True ; large format resP@gsnPanelLabelBar = True ; add common colorbar resP@lbLabelAutoStride = True ; auto stride on labels ; yStrt = yyyymm(0)/100 yStrt = 1998 yLast = 2012 ; yLast = yyyymm(nyrs-1)/100 resP@txString = "SST: "+season+": "+yStrt+"-"+yLast ;******************************************* ; first plot ;******************************************* do n=0,neof-1 res@gsnLeftString = "EOF "+(n+1) res@gsnRightString = sprintf("%5.1f", eof@pcvar(n)) +"%" plot(n)=gsn_csm_contour_map_ce(wks,eof(n,:),res) end do gsn_panel(wks,plot,(/neof,1/),resP) ; now draw as one plot ;******************************************* ; second plot ;******************************************* ; EOF time series [bar form] rts = True rts@gsnDraw = False ; don't draw yet rts@gsnFrame = False ; don't advance frame yet rts@gsnScale = True ; force text scaling ; these four rtsources allow the user to stretch the plot size, and ; decide exactly where on the page to draw it. rts@vpHeightF = 0.40 ; Changes the aspect ratio rts@vpWidthF = 0.85 rts@vpXF = 0.10 ; change start locations rts@vpYF = 0.75 ; the plot rts@tiYAxisString = "Pa" ; y-axis label rts@gsnYRefLine = 0. ; reference line rts@gsnXYBarChart = True ; create bar chart rts@gsnAboveYRefLineColor = "red" ; above ref line fill red rts@gsnBelowYRefLineColor = "blue" ; below ref line fill blue ; panel plot only resources rtsP = True ; modify the panel plot rtsP@gsnMaximize = True ; large format rtsP@txString = "SST: "+season+": "+yStrt+"-"+yLast year = yyyymm/100 ; create individual plots do n=0,neof-1 rts@gsnLeftString = "EOF "+(n+1) rts@gsnRightString = sprintf("%5.1f", eof@pcvar(n)) +"%" plot(n) = gsn_csm_xy (wks,year,eof_ts(n,:),rts) end do gsn_panel(wks,plot,(/neof,1/),rtsP) ; now draw as one plot end