Thanks Dennis, the Natgrid routines appear to do what
I want them to, and your example was very helpful.
I've also learned that Natgrid seems to be incorporated
in the standard NCL package (at least in recent versions),
so I didn't need to do any sort of installation.
Don
On Thu, 18 Nov 2004 15:39:54 -0700 (MST), "Dennis Shea"
<shea@cgd.ucar.edu> said:
> >
> >Hello, a couple of weeks ago I sought help on extracting
> >point data from a grid, given a lat/lon that didn't necessarily
> >correspond to the grid points. I was led to linint2_points() and
> >have made "some" progress, but of course, I had to complicate
> >the initial requirements a little bit :)
> >
> >Given a regular 3x3 grid varying from 45N to 47N and 112W to 114W,
> >and a temperature field, I can easily extract/interpolate point
> >data as follows:
> >
> >;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> >begin
> >
> > temperature = (/ (/100, 100, 50/), \
> > (/100, 100, 100/), \
> > (/100, 100, 100/) /)
> >
> > lat = (/45, 46, 47/)
> > lon = (/-114, -113, -112/)
> >
> > ; Here's the lat/lon I want station data for
> > y0 = 45.5
> > x0 = -112.5
> >
> > station_temperature = linint2_points(lon, lat, temperature, \
> > False, x0, y0, 0)
> >
> > print(station_temperature)
> >
> >end
> >;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> >
> >However, what if my grid is NOT so "regular?" Suppose, as is the
> >case with some of my WRF output, a "row" of points might start at
> >45N and end at 45.2N? Does linint2_points() allow for this? In
> >the following, I've naively tried by placing all of my lat/lon grid
> >points in arrays and passing those to linint2_points(), as follows
> >
> >;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> >begin
> >
> > temperature = (/ (/100, 100, 100/), \
> > (/100, 50, 100/), \
> > (/100, 100, 100/) /)
> >
> > lat_points = (/ (/45, 45.1, 45.2 /), \
> > (/46, 46.1, 46.2 /), \
> > (/47, 47.1, 47.2 /) /)
> >
> > lon_points = (/ (/-114, -113, -112/), \
> > (/-114.1, -113.1, -112.1/), \
> > (/-114.2, -113.2, -112.2/) /)
> >
> > ; Here's the lat/lon of the location I want station data for
> > y0 = 45.5
> > x0 = -112.5
> >
> > station_temperature = linint2_points(lon_points, lat_points, \
> > temperature, False, x0, y0, 0)
> >
> > print(station_temperature)
> >
> >end
> >;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> >
> >
> >but I get an error message saying:
> >
> >fatal:linint2_points: If xi is not one-dimensional, then it must have
> >one less dimension than fi
> >
> >(Note that above, "xi" corresponds to "lon_points" and "fi" corresponds
> >to "temperature").
> >
> >Indeed, the documentation confirms that xi must have one less dimension
> >than fi, but
> >I don't understand why this would be the case (which suggests I don't
> >really
> >understand the purpose of multi-dimension arguments for xi and yi). Is
> >there any place
> >I can go for a little more explanation? A Google search on
> >"linint2_points" doesn't
> >turn up much :)
>
> ============
> Hello
>
> linint2_points was developed to interpolate from
> a "regular" grid to arbitrary locations within the regular grid.
> By "regular", I mean a grid where any point on the grid
> can be accessed via one-dimensional (1d) coordinate variables.
> To be considered a "coordinate variabe" it must be a 1D
> array of monotonically {in/de}creasing values. They need not
> be equally spaced.
>
> The WRF grid is curvilinear. It requires 2D lat/lon arrays to
> describe the grid. It is not so easy to develop a fast, general
> purpose code to interpolate from any arbitrary grid to a point.
> Suffice it to say, some grids are *very* unusual.
>
> I asked the MMM/WRF 'user support' person if they had
> an interpolation function. Here is the response:
>
> > We probably do not have the exact thing you
> > ask for. We have some routines that converts
> > x, y to lat, long, and vice versa. But the x,y
> > here refers to the coordinate on projected maps.
> > If that would be useful to you, I will pass them to you.
>
> ---
> NCL's Natgrid interpolation library may work.
> I have never used Natgrid. Also, I have heard it can
> be slow. The following is untested
>
> http://ngwww.ucar.edu/ngdoc/ng/ref/ncl/functions/nnptinit.html
> http://ngwww.ucar.edu/ngdoc/ng/ref/ncl/functions/nnpnt.html
>
> Let's say lat2d(lat,lon) and lon2d(lat,lon) are the grid coordinates.
> for the variable x(ntim,lat,lon). Let latSta(N), lonSta(N),
> xSta(ntim,N) be the station observations
>
> lat1d = ndtooned(lat2d)
> lon1d = ndtooned(lon2d)
>
> do nt=0,ntim-1
> nnpntinit(lon1d, lat1d, ndtooned(x(nt,:,:)) )
> xSta = nnpnts(lonSta,latSta) ; xSta(N)
> :
> end do
>
> Good Luck
> D
>
--
** NOT ON CAMPUS UNTIL SPRING 2005 **
**************************************************************
Don Morton http://www.cs.umt.edu/~morton/
Department of Computer Science The University of Montana
Missoula, MT 59812 | Voice (406) 243-4975 | Fax (406) 243-5139
May your trails be crooked, winding, lonesome, dangerous, leading to the
most amazing view. May your mountains rise into and above the clouds.
-Edward Abbey, naturalist and author (1927-1989)
_______________________________________________
ncl-talk mailing list
ncl-talk@ucar.edu
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
This archive was generated by hypermail 2b29 : Wed Dec 01 2004 - 10:32:22 MST