[ncl-talk] test if lat-lon point is within gridded model domain

Tue Nov 21 16:58:32 MST 2017

Hi Bill,

Thanks for clarifying some of that, and alerting me to other use cases that
leverage the current functionality. Perhaps an alternative solution to the
one you propose would be to allow an attribute to be attached to opt that
would trigger this automatic bounds checking, and replace any (i,j) outside
the bounds of the domain with missing values. (And of course, whatever
action is taken, some clarification of the behavior of this function in the
online documentation would be helpful.)

Jared

On Tue, Nov 21, 2017 at 4:42 PM, Bill Ladwig <ladwig at ucar.edu> wrote:

> Wrf_user_ll_to_ij does not do any bounds checking to ensure points are in
> the WRF domain because it can also be generically used to determine grid
> point indexes based only on a map projection and geographic points.  The
> underlying function is wrf_ll_to_ij (https://www.ncl.ucar.edu/Docu
> ment/Functions/Built-in/wrf_ll_to_ij.shtml), where the wrf_user_ll_to_ij
> version does the extraction of the map projection parameters for you.  So,
> it's possible to use the underlying function without any actual WRF data
> files.  I'd be leery of taking away the ability of the function to answer
> questions like "how big of a grid do I need to include this lat,lon
> point?".  So, we can probably create a new function, say
> "in_wrf_domain(...)", that applies the same bounds checking that you are
> doing on the result.  I'll make a ticket for it.
>
> Regards,
>
> Bill
>
> On Tue, Nov 21, 2017 at 2:24 PM, Jared Lee <jaredlee at ucar.edu> wrote:
>
>> Also, as one more follow-up, I think I found an acceptable work-around to
>> the issue of identifying which observations are within the WRF domain. I
>> ran all the observation (lon,lat) pairs through wrf_user_ll_to_ij, and then
>> found the indices of the (i,j) pairs that were actually within the WRF
>> domain (that's the workaround for the bug in wrf_user_ll_to_ij). I've
>> checked the remaining non-Alaska observations that make it through that
>> process, and those are indeed within the domain. So I think that works for
>> my needs.
>>
>> Jared
>>
>> On Tue, Nov 21, 2017 at 11:34 AM, Jared Lee <jaredlee at ucar.edu> wrote:
>>
>>> Hi Mary,
>>>
>>> I'm attaching both the mapgrid.png and another map of the domain I made
>>> earlier from my own script. (The landmask plot was made from a wrfout file
>>> that contained all the static fields, such as LANDMASK and HGT; the wrfout
>>> files I'm reading in are subsetted with only a small number of variables to
>>> minimize storage requirements.)
>>>
>>> I have also written a demo script to show what I believe are bugs both
>>> in the wrf_user_ll_to_ij and rcm2points_Wrap functions. That script points
>>> to a sample wrfout and madis file on glade. I've attached that script here
>>> as well. The wrf_user_ll_to_ij function returns nonsense (i,j) values for
>>> (lon,lat) points outside the WRF domain. rcm2points also interpolates WRF
>>> data to many of these points outside the WRF domain. I set up some print
>>> statements in stages of thinning the data to progressively highlight the
>>> problems with both functions.
>>>
>>> Preferred behavior (from my perspective) would be to return missing
>>> values for (lat,lon) points outside the WRF domain, both for
>>> wrf_user_ll_to_ij and rcm2points. At the very least, wrf_user_ll_to_ij
>>> should not give non-existing (i,j) coordinates. And rcm2points should only
>>> be an interpolation function, not an extrapolation function.
>>>
>>> Jared
>>>
>>>
>>>
>>> On Mon, Nov 20, 2017 at 9:59 PM, Mary Haley <haley at ucar.edu> wrote:
>>>
>>>> Hi Jared,
>>>>
>>>> It would help if I could see a picture of your WRF grid. You said that
>>>> it is polar stereographic, yet it crosses the dateline.
>>>>
>>>> I think that gc_inout might still work for you, but I agree there could
>>>> be issues if it crosses the dateline. However, it's possible we could break
>>>> this up into two grids, but I first need to see what it looks like.
>>>>
>>>> Can you take the attached script, plug in your WRF file, run it, and
>>>> send me the mapgrid.png file?
>>>>
>>>> Or, even better, if have the WRF and MADIS files on the NCAR supers,
>>>> can you point me to them?
>>>>
>>>> Thanks,
>>>>
>>>> --Mary
>>>>
>>>>
>>>>
>>>>
>>>> On Mon, Nov 20, 2017 at 6:17 PM, Jared Lee <jaredlee at ucar.edu> wrote:
>>>>
>>>>> Hi Dave,
>>>>>
>>>>> Thanks for the suggestion, but I'm not optimistic that that will work
>>>>> for my problem. I think finding the boundaries of an arbitrary WRF grid as
>>>>> a spherical polygon would be a non-trivial task; I don't know how I would
>>>>> even go about doing that.
>>>>>
>>>>> Thinking about it, I can probably approximately solve this problem in
>>>>> a few steps with a combination of wrf_user_ll_to_ij and then gc_latlon (by
>>>>> calculating the distance between obs and "nearest grid points" along the
>>>>> edge of the WRF domain and then applying a distance threshold). I'll try
>>>>> that tomorrow. I think I've found a bug with wrf_user_ll_to_ij, though (it
>>>>> gives me non-existent (i,j) values for (lon,lat) points outside the
>>>>> domain), and I'll submit a separate ticket on that.
>>>>>
>>>>> Jared
>>>>>
>>>>> On Mon, Nov 20, 2017 at 5:39 PM, Dave Allured - NOAA Affiliate <
>>>>> dave.allured at noaa.gov> wrote:
>>>>>
>>>>>> Jared,
>>>>>>
>>>>>> See whether the function gc_inout will work in your case.  You would
>>>>>> need to get the boundaries of the WRF grid as a spherical polygon.  There
>>>>>> may already be published WRF polygons as shape files; I am not familiar
>>>>>> with that.  Then use the function gc_inout to determine inside or
>>>>>> outside for each obs location.
>>>>>>
>>>>>> In essence you would be applying two mask tests to each obs point,
>>>>>> because you also say you want to select within a lat/lon box.  Use gc_inout
>>>>>> to test against the larger WRF grid, and use simple arithmetic comparisons
>>>>>> for the lat/lon box.
>>>>>>
>>>>>> Caution, there is a little trap here.  Do not use gc_inout for
>>>>>> testing simple four-sided lat/lon bounding boxes.  You would probably get
>>>>>> something unexpected.  gc_inout uses spherical polygons on the earth's
>>>>>> surface.  Such polygons have curved edges on the simple lat/lon plane.
>>>>>>
>>>>>> --Dave
>>>>>>
>>>>>>
>>>>>> On Mon, Nov 20, 2017 at 4:26 PM, Jared Lee <jaredlee at ucar.edu> wrote:
>>>>>>
>>>>>>> Hi, I have a WRF domain on a polar stereographic grid (which also
>>>>>>> happens to span the dateline), and I want to interpolate model data to
>>>>>>> MADIS observation locations that are within my domain. Is there a
>>>>>>> straightforward way to evaluate which observation lat/lon points are
>>>>>>> outside the WRF grid?
>>>>>>>
>>>>>>> Because it's a polar stereographic grid, doing a first pass to
>>>>>>> eliminate MADIS stations by comparing to the min/max latitude and longitude
>>>>>>> of the WRF grid still leaves a ton of geographic area that would be inside
>>>>>>> that lat/lon box but outside the WRF grid (and that's complicated further
>>>>>>> by the WRF domain spanning the dateline).
>>>>>>>
>>>>>>> I tried using rcm2points to do horizontal interpolation, but that
>>>>>>> function is still interpolating to numerous grid points that are thousands
>>>>>>> of km outside my WRF domain. I also tried feeding wrf_user_ll_to_ij some
>>>>>>> lat/lon values for points well outside my domain, and it gives me nonsense
>>>>>>> (but non-missing) values for nearest grid points.
>>>>>>>
>>>>>>> Anyone have any ideas? Or is there a function that already exists to
>>>>>>> do this that I'm apparently not seeing? It seems like a function to return
>>>>>>> the (i,j) values of the four surrounding grid points would be an ideal way
>>>>>>> to accomplish this and be useful in additional contexts and applications.
>>>>>>>
>>>>>>> Jared
>>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> ===============================
>>>>> Jared A. Lee, Ph.D.
>>>>> Project Scientist I
>>>>> Research Applications Laboratory
>>>>> National Center for Atmospheric Research
>>>>> Boulder, Colorado, USA
>>>>>
>>>>> Member, AMS Planning Commission
>>>>>
>>>>> Email: jaredlee at ucar.edu (w)
>>>>> Phone: 303.497.8485 <(303)%20497-8485> (w)
>>>>> Web: https://staff.ucar.edu/users/jaredlee
>>>>> ===============================
>>>>>
>>>>> _______________________________________________
>>>>> ncl-talk mailing list
>>>>> ncl-talk at ucar.edu
>>>>> List instructions, subscriber options, unsubscribe:
>>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>>
>>>>>
>>>>
>>>
>>>
>>> --
>>> ===============================
>>> Jared A. Lee, Ph.D.
>>> Project Scientist I
>>> Research Applications Laboratory
>>> National Center for Atmospheric Research
>>> Boulder, Colorado, USA
>>>
>>> Member, AMS Planning Commission
>>>
>>> Email: jaredlee at ucar.edu (w)
>>> Phone: 303.497.8485 <(303)%20497-8485> (w)
>>> Web: https://staff.ucar.edu/users/jaredlee
>>> ===============================
>>>
>>
>>
>>
>> --
>> ===============================
>> Jared A. Lee, Ph.D.
>> Project Scientist I
>> Research Applications Laboratory
>> National Center for Atmospheric Research
>> Boulder, Colorado, USA
>>
>> Member, AMS Planning Commission
>>
>> Email: jaredlee at ucar.edu (w)
>> Phone: 303.497.8485 <(303)%20497-8485> (w)
>> Web: https://staff.ucar.edu/users/jaredlee
>> ===============================
>>
>> _______________________________________________
>> ncl-talk mailing list
>> ncl-talk at ucar.edu
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>>
>

-- 
===============================
Jared A. Lee, Ph.D.
Project Scientist I
Research Applications Laboratory
National Center for Atmospheric Research
Boulder, Colorado, USA

Member, AMS Planning Commission

Email: jaredlee at ucar.edu (w)
Phone: 303.497.8485 (w)
Web: https://staff.ucar.edu/users/jaredlee
===============================
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171121/00378ebc/attachment.html>