[ncl-talk] regridding

[igor akkerman]

Hi,

I need to regrid NCEP-NCAR reanalysis to a rotated regional model grid that
is set up as follows:

Sketch of the rotated grid (rot.lat, rot.lon)
(24.5,-27.5) ------ (24.5,27.0)
           |                         |
           |                         |
(-25.0,-27.5) ------ (-25.0,27.0)
This is a regular grid with uniform resol. of 0.5.

and in geograf. grid (geograf. lat, lon):
(53.8,-134.6) ------ (54.2,135.1)
           |                         |
           |                         |
(53.5,-44.7) ------ (53.85,44.2)

I also have an ascii file with rotated grid cell coordinates and equivalent
geographical coordinates called 'hirham_rot_geo_koor.dat".

I tried making a script but can't correctly assign the conversion between
rotated and geographical grid coordinates. Could you please take a look at
the script and suggest how to fix it?

The error message is that SLP_regrid variable has only one coordinate
instead of three.


  WRITE_NETCDF = False

; TEXT FILE CONVERTING ROTATED AND SPHERICAL LATS AND LONS
filename = "hirham_rot_geo_koor.dat"
coords = asciiread(filename,(/100,110,4/),"float")
lat2d = coords(:,:,0) ; reading rotated latitude from ascii file
lon2d = coords(:,:,1) ; reading rotated longitude from ascii file
printMinMax(lon2d, True)

;---NCEP-NCAR reanalysis
    srcFileName = "slp.1948.nc"
    sfile = addfile(srcFileName,"r")

    Opt                = True
    Opt at SrcTitle       = "NCEP-NCAR reanalysis"   ; optional

    Opt at WgtFileName    = "NCEP_to_Rect.nc"

    Opt at ForceOverwrite = True

    SLP = sfile->slp
    SLP at lat2d = lat2d
    SLP at lon2d = lon2d

printVarSummary(SLP)

    dims  = dimsizes(lat2d)
    nlat  = dims(0)
    nlon  = dims(1)

    Opt at SrcFileName     = "Rectilinear.nc"      ; Name of source and
    Opt at DstFileName     = "WRF_SCRIP.nc"    ; destination files

;---Create the destination lat/lon grid
    printVarSummary(lat2d)
    printVarSummary(lon2d)
    printMinMax(lat2d, True)
    printMinMax(lon2d, True)
    lat = fspan( -25, 24.5,nlat)
    lon = fspan( -27.5,27,nlon)
    lat at units = "degrees_north"
    lon at units = "degrees_east"
    lat!0     = "lat"
    lon!0     = "lon"
    lat&lat   = lat
    lon&lon   = lon

    Opt at DstGridType          = "rectilinear"
    Opt at DstGridLat           = lat
    Opt at DstGridLon           = lon

    Opt at InterpMethod         = "bilinear"
    Opt at SrcRegional          = True
    Opt at DstRegional          = True
    Opt at Debug                = True
    SLP_regrid = ESMF_regrid(SLP,Opt)     ; Do the regridding for TMP
;
; The source and destination grid description files and
; weight file will be the same for the next call to
; ESMF_grid, so no need to regenerate them.
;
    Opt at SkipSrcGrid   = True
    Opt at SkipDstGrid   = True
    Opt at SkipWgtGen    = True

;---Reset 0 values to missing values.
    SLP_regrid at _FillValue = default_fillvalue(typeof(SLP_regrid))

    SLP_regrid            = where(SLP_regrid.eq.0.0,SLP_regrid at _FillValue,\
                                  SLP_regrid)

    printVarSummary(SLP_regrid)

;----------------------------------------------------------------------
; Plotting section
;
; This section creates filled contour plots of both the original
; data and the regridded data, and panels them.
;----------------------------------------------------------------------
    wks_slp = gsn_open_wks("png","interpolate_slp")

    res                       = True

    res at gsnDraw               = False
    res at gsnFrame              = False

    res at cnFillOn              = True
    res at cnLinesOn             = False
    res at cnLineLabelsOn        = False
    res at lbLabelBarOn          = False
    res at cnLevelSelectionMode  = "ManualLevels"

    res at gsnPolar              = "NH"
    res at mpMinLatF             = min(lat)
    res at gsnAddCyclic          = False

    nrec      = 0
    dims_orig = tostring(dimsizes(SLP(nrec,:,:)))

    mnmxint_slp = nice_mnmxintvl( min(SLP), max(SLP), 18, False)

;---SLP
    res at cnMinLevelValF  = mnmxint_slp(0)
    res at cnMaxLevelValF  = mnmxint_slp(1)
    res at cnLevelSpacingF = mnmxint_slp(2)/2.   ; Create more levels
    res at tiMainFontHeightF = 0.015

;---SLP regridded
    res at tiMainString  = "rectilinear grid (" + Opt at InterpMethod + "
interpolation)"
    slp_regrid = gsn_csm_contour_map(wks_slp,SLP_regrid(nrec,:,:),res)

;---SLP original
    res at tiMainString = "SLP on original curvilinear grid (" +  \
                        str_join(dims_orig," x ") + ")"

    slp_orig = gsn_csm_contour_map(wks_slp,SLP(nrec,:,:),res)

;---Compare the plots in a panel
    pres                   = True
    pres at gsnMaximize       = True
    pres at gsnPanelLabelBar  = True

    gsn_panel(wks_slp,(/slp_orig,slp_regrid/),(/1,2/),pres)

;---Trim
    system("convert -trim interpolate_slp.png interpolate.png")
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150205/31c5b748/attachment.html