[ncl-talk] regridding
Mary Haley
haley at ucar.edu
Fri Feb 6 12:23:23 MST 2015
Hi Alex,
Thanks for providing the files and a script that I could run without any
changes. :-)
The issue appears to be a misunderstanding of the type of grid you are
starting from (the "source" grid) and the grid you want to regrid to (the
"destination" grid).
Here are the issues:
- The "SLP" variable you read off the file has coordinate arrays already
attached to it, so this is a rectilinear grid. These two lines in your
code are therefore unnecessary and actually wrong:
SLP at lat2d = lat2d
SLP at lon2d = lon2d
- You created "lat" and "lon" arrays using "fspan" and then used these
as the destination grid. This is unnecessary too, because you already have
the destination grid in "lat2d" and lon2d".
Here are the fixes:
- Leave "SLP" untouched when you do the regridding. "ESMF_regrid" will
see the lat/lon coordinate arrays associated with this SLP, and will
automatically use this for the source grid.
- Remove the generation of "lat" and "lon" as they are not needed.
- Use lat2d and lon2d for DstGridLat and DstGridLon, and don't bother
setting DstGridType (ESMF_regrid can figure out the type from DstGridLat
and DstGridLon):
Opt at DstGridLat = lat2d ; curvilinear grid
Opt at DstGridLon = lon2d
- I don't think this "where" code is needed:
;---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)
- For the plotting section, you want to set res at gsnAddCyclic to True for
the original grid, because it's global and needs to have the cyclic point
added. Set it to False for the regional grid.
- I recommend not doing a polar stereographic plot, because the area you
are regridding to is not over the north or south pole.
- I added a MAP_ZOOM logical variable so you could have the option of
zooming in on your regional grid.
See the attached new version of your script and two sample images. One
image is the full global map, and the second one is zoomed in on the
regional grid you read off the ASCII file.
--Mary
-
On Fri, Feb 6, 2015 at 1:43 AM, igor akkerman <igorakkerman89 at gmail.com>
wrote:
> Hi Mary,
>
> Please see the files on your ftp server:
>
> 1) Code (regrid.ncl)
> 2) Original NCEP-NCAR reanalysis file (slp.1948.nc)
> 3) ASCII conversion file from rotated grid to geographical grid.
> (hirham_rot_geo_koor.dat)
>
> Thanks for your help,
>
> Alex
>
> 2015-02-05 20:00 GMT+04:00 Mary Haley <haley at ucar.edu>:
>
>> Dear Igor,
>>
>> We need a little more information.
>>
>> First, it is helpful if you copy-and-paste the actual error message in
>> your email.
>>
>> Second, please include the output from "printVarSummary(SLP)".
>>
>> Really, though, it helps if you can provide the data files so we can run
>> the script and see what the issue is.
>>
>> Thanks,
>>
>> --Mary
>>
>>
>> On Thu, Feb 5, 2015 at 2:52 AM, igor akkerman <igorakkerman89 at gmail.com>
>> wrote:
>>
>>> 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@
>>> _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")
>>>
>>> _______________________________________________
>>> ncl-talk mailing list
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>
>>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150206/d9222bce/attachment-0001.html
-------------- next part --------------
A non-text attachment was scrubbed...
Name: regrid_mod.ncl
Type: application/octet-stream
Size: 4686 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150206/d9222bce/attachment-0001.obj
-------------- next part --------------
A non-text attachment was scrubbed...
Name: interpolate_slp_zoom.png
Type: image/png
Size: 117668 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150206/d9222bce/attachment-0002.png
-------------- next part --------------
A non-text attachment was scrubbed...
Name: interpolate_slp.png
Type: image/png
Size: 361681 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20150206/d9222bce/attachment-0003.png
More information about the ncl-talk
mailing list