;************************************************************* ; regrid_12.ncl ; ; Concepts illustrated: ; - Interpolating from a fixed grid to a lower resolution using conservative remapping ; - Computing dispersion statistics ; - Drawing a cylindrical equidistant map ;************************************************************* 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 ;********************************************* ; Read global grid [180W to 180E] ; ELEV: [LAT | 5401] x [LON | 10801] *includes* cyclic pt ; We do *not* want to read the cyclic point. ;********************************************* diri = "/Users/shea/Data/netCDF/" fili = "ETOPO2_GLOBAL_2_ELEVATION.nc" ; NGDC f = addfile(diri+fili,"r") Z = short2flt( f->ELEV(:,0:10799) ) ; Do not read the cyclic pt printVarSummary(Z) ; Z: min=-10654 max=8593 ;********************************************* ; areal average conservative remapping ;********************************************* NLATR = 180 MLONR = 360 LATR = latGlobeFo(NLATR, "LATR", "latitude" , "degrees_north") LONR = lonGlobeFo(MLONR, "LONR", "longitude", "degrees_east") ; 0-360 LONR = (/LONR-180/) ; -180 to 180 LONR&LONR = LONR ZR = area_conserve_remap_Wrap (Z&LON,Z&LAT,Z, LONR,LATR, False) ; 5.2.0 ;ZR = area_hi2lores_Wrap (Z&LON,Z&LAT,Z, True, 1.0, LONR,LATR, False) ; 5.1.0 ZR@long_name = ZR@long_name +": 1x1 grid" printVarSummary(ZR) ; ZR (1x1): min=-10439.7 max=5713.99 ;********************************************* ; Calculate assorted statistics and print ; stat_dispersion is a bit slow for *large* variables ;********************************************* opt = True opt@PrintStat = True ;statZ = stat_dispersion(Z , opt ) ; bit slow (5401,10800) statZR = stat_dispersion(ZR, opt ) ;********************************************* ; Plot ;********************************************* wks = gsn_open_wks ("ps", "regrid") gsn_define_colormap(wks, "BlAqGrYeOrReVi200") setvalues NhlGetWorkspaceObjectId() "wsMaximumSize": 500000000 end setvalues res = True res@gsnMaximize = True res@gsnAddCyclic = False ; dataset allready has cyclic pt res@gsnSpreadColors = True res@cnFillMode = "RasterFill" res@cnLinesOn = False ; Turn off contour lines res@cnFillOn = True ; Turn on contour fill res@cnLevelSelectionMode = "ManualLevels" ; set manual contour levels res@cnMinLevelValF = -5000 ; set min contour level res@cnMaxLevelValF = 5000 ; set max contour level res@cnLevelSpacingF = 500 res@lbLabelAutoStride = True res@mpCenterLonF = 0 res@mpOutlineOn = True res@mpFillOn = False plot = gsn_csm_contour_map_ce(wks,Z , res) plot = gsn_csm_contour_map_ce(wks,ZR, res) end