Re: rcm2grid_Wrap returns missing values

Hi Kat

I see the "problem".

The rcm2rgrid is northern hemisphere centric.
The underlying algorithm expects the latitudes to
increase [eg, lat2d(j+1,0) > lat2d(j,0) . The grid
you are using does just the opposite.

You can see this by using

    print(mm2dlatS(:,0))

The following shows the expected order

    print(mm2dlatS(::-1,0))

The solution is to use NCL's syntax to reverse the
lat (also, the corrsponding lon and data) order.

Replace:
    newmadi = rcm2rgrid_Wrap(mm2dlatS, mm2dlonS, masked_madiS,
newlat500, newlon500,1)

With:
    newmadi = rcm2rgrid_Wrap(mm2dlatS(::-1,:), mm2dlonS(::-1,:),
masked_madiS(::-1,:), newlat500, newlon500,1)

I will send you a modified script offline.

Good luck
D

> Hi,
> I am having problems with regridding some data non-continuous (i.e. snow
> cover data) from a 2d curvilinear grid (Lambert conformal) to a
> rectilinear grid. Whilst there are no _FillValues in the original data,
> the data returned from the rcm2grid_Wrap function has been replaced with
> almost all _FillValues after regridding.
>
> I was hoping someone might have some ideas on what is going on here,
> since the regrid function appears to run successfully and the output
> variable is the correct size and has the correct attributes - just the
> data is not what I expected.
>
>
> Note: most of the values in the original data are zero, except for the
> patches of snow which are defined by floating point numbers ranging from
> 0.001-120.0. There are clear patches of snow in the data that sum to
> over 1500km2.
>
> Also note: I am not changing the resolution of the data (just the grid),
> so there shouldn't be problems with too much information lost as a
> result of changing resolution.
>
> Thanks in advance.
> Kat Bormann
>
> The relevant parts of my script (and reports back to screen) are shown
> below:
>
> ;READ DATA FROM BINARY FILe
> ;(floating point, 2D array, direct data order)
> data = fbindirread(filename.img,0,(/2300,1400/),"float")
> printVarSummary(data)
> print("# original data missing values="+num(ismissing(data)))
>
> ;read spatial information for binary data and assign to data variable
> latdat = addfile(latfilename.nc,"r")
> londat = addfile(lonfilename.nc,"r")
> mm2dlat = mmlatdat->mmlat2d
> mm2dlon = mmlondat->mmlon2d
> data@lat2d = mm2dlat
> data@lon2d = mm2dlon
>
> ;confirm spatial data assignment
> ;set some resources
> res=True
> res@gsnAddCyclic = False
> res@cnFillMode="CellFill"
> res@cnRasterSmoothingOn = False
> res@cnLinesOn = False
> res@mpDataBaseVersion = "HighRes"
>
> ;plot data on map (this looks fine and as expected)
> wks = gsn_open_wks("ps","plotfile")
> gsn_csm_contour_map(wks,data,res)
>
> ;SUBSET DATA EXTENT TO LOCAL REGION
> ;to reduce regridding interpolations and runtimes
> ;set local spatial boundaries
> latS = -36.5
> latN = -35.5
> lonW = 148.0
> lonE = 149.0
>
> ;find array positions for specified local spatial boundary
> ji = region_ind (mm2dlat,mm2dlon, latS, latN, lonW, lonE)
> jStrt = ji(0) ;lat start
> jLast = ji(1) ;lat last
> iStrt = ji(2) ;lon start
> iLast = ji(3) ;lon last
>
> ;use array positions to subset lat, lon and data variables
> mm2dlatS= mm2dlat(jStrt:jLast, iStrt:iLast)
> mm2dlonS= mm2dlon(jStrt:jLast, iStrt:iLast)
> dataS = (/data(jStrt:jLast, iStrt:iLast)/) ;value only assignment
>
> ;reassign lat,lon coordinates for local spatial region
> dataS@lat2d = mm2dlatS
> dataS@lon2d = mm2dlonS
> dataS!0 = "mm2dlatS"
> dataS!1 = "mm2dlonS"
> printVarSummary(dataS)
> print("# spatial subset data missing values="+num(ismissing(dataS)))
>
> ;REGRID DATA
> ;regrid data within the local boundary to a rectilinear grid
> ;find actual coordinates of new local region after subsetting the data
> mdsize = dimsizes(dataS)
> TLlat = mm2dlatS(1,1)
> TLlon = mm2dlonS(1,1)
> BRlat = mm2dlatS(mdsize(0)-2, mdsize(1)-2)
> BRlon = mm2dlonS(mdsize(0)-2, mdsize(1)-2)
>
> ;generate rectilinear grid for subset region 500m res.
> deg500m = 0.004166667 ;500m expressed in decimal degrees (approx)
> nptsY500 = round(((TLlat - BRlat)/deg500m),3)
> nptsX500 = round(((BRlon - TLlon)/deg500m),3)
> newlat500 = fspan(BRlat, TLlat, nptsY500)
> newlon500 = fspan(TLlon, BRlon, nptsX500)
>
> ;assign attributes to new grid
> newlat500@res = "500m intervals within subset"
> newlat500@units = "degrees north"
> newlon500@res = newlat500@res
> newlon500@units = "degrees east"
>
> ;regrid data to new rectilinear grid coordinates
> newdata = rcm2rgrid_Wrap(mm2dlatS, mm2dlonS, dataS, newlat500, newlon500,1)
> printVarSummary(newdata)
> print("# missing values in regridded data="+num(ismissing(newdata)))
>
> ;---the end
>
> NCL returns:
>
> Variable: data
> Type: float
> Total Size: 12880000 bytes
> 3220000 values
> Number of Dimensions: 2
> Dimensions and sizes: [2300] x [1400]
> Coordinates:
>
> (0) # original data missing values=0
>
> Variable: dataS
> Type: float
> Total Size: 167240 bytes
> 41810 values
> Number of Dimensions: 2
> Dimensions and sizes: [mm2dlatS | 226] x [mm2dlonS | 185]
>
> (0) # spatial subset data missing values=0
>
> Variable: newdata
> Type: float
> Total Size: 235104 bytes
> 58776 values
> Number of Dimensions: 2
> Dimensions and sizes: [MM2DLATs | 237] x [MM2DLONs | 248]
> Coordinates:
> MM2DLATs: [-36.49172..-35.50626]
> MM2DLONs: [147.9893..149.0239]
> Number Of Attributes: 2
> _FillValue : 9.96921e+36
> ncl : rcm2rgrid used for interpolation
>
> (0) # missing values in regridded data=58683
>
>
>
>
> --
> Kathryn Bormann
> PhD Candidate
> Climate Change Research Centre
> University of New South Wales
> Kensington Campus
> Sydney, NSW, 2052
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
>>> --
>>> Kathryn Bormann
>>> PhD Candidate
>>> Climate Change Research Centre
>>> University of New South Wales
>>> Kensington Campus
>>> Sydney, NSW, 2052
> >
>
>
> --
> Kathryn Bormann
> PhD Candidate
> Climate Change Research Centre
> University of New South Wales
> Kensington Campus
> Sydney, NSW, 2052
>
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Wed Jun 15 09:17:48 2011