Workaround for .sel() with twodimensional coordinates

[Caio Stringari]

Hello,

I am working on a dataset that looks like this:

<xarray.Dataset>

Dimensions:     (i: 3, j: 3, u: 259, v: 1257)

Coordinates:

    X          (u, v) float64 459.5 459.5 459.6 459.6 459.7 459.8 459.8 ...

    Y          (u, v) float64 527.5 527.5 527.5 527.5 527.5 527.5 527.5 ...

    time        datetime64[ns] 2016-08-15T02:33:01

  * i           (i) int64 0 1 2

  * j           (j) int64 0 1 2

  * u          (u) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 ...

  * v          (v) int64 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 ...

Data variables:

    red         (u, v) uint8 50 49 48 48 48 49 49 49 49 49 49 49 49 49 48 48 ...

    green     (u, v) uint8 92 91 90 90 90 91 91 91 91 91 91 91 91 91 90 90 ...

    blue       (u, v) uint8 117 116 115 115 115 116 116 116 116 116 116 116 ...

    gray       (u, v) uint8 85 84 83 83 83 84 84 84 84 84 84 84 84 84 83 83 ...

    homography  (i, j) float64 -0.1432 3.836 1.098e+03 -0.5929 3.26 ...

It is a extracted frame from a video recording that had the pixel coordinates (u,v) projected to real-world coordinates (X,Y). While (u,v) are unidimensional, (X,Y) became two dimensional in the process of coordinate transformation due to the image warping. 

Now, I need to extract some pixel intensity values in (very well) known real-world coordinates. I would like to apply something something like ds["gray"].sel(X=x_points,Y=y_points), however, xarray raises ValueError: Coordinate objects must be 1-dimensional.

Have anyone worked on a solution for a similar problem ?

Thanks in advance

Caio E. Stringari

[Stephan Hoyer]

You can efficiently look up nearest-neighbor positions in two dimensions using a KDTree (e.g., scipy.spatial.cKDTree). To use SciPy's KDTree directly with 2D arrays, you'll need to need to flatten them first and lookup positions with np.unravel_index. It might be easier to keep track of indexes if you flatten first, e.g., ds = ds.stack(space=['x', 'y']).

We've contemplated adding a wrapper for a KDTree in xarray to facilitate just this type of indexing, but haven't gotten around to it yet. See https://github.com/pydata/xarray/issues/475 for more discussion.

There is pretty good precedence for the KDTree API with how we are handling pandas.MultiIndex (https://github.com/pydata/xarray/pull/947).

--
You received this message because you are subscribed to the Google Groups "xarray" group.
To unsubscribe from this group and stop receiving emails from it, send an email to xarray+unsubscribe@googlegroups.com.
To post to this group, send email to xar...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/xarray/b0c3c8a3-596e-4456-86b5-eeaf5e142c11%40googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Caio Stringari]

Thanks Stephan !

I had a KDTree search kind of working before posting here but had no knowledge of np.unravel_index, this made my code much simpler. Looking forward for the xarray implementation of kdtree searches.

Caio

To unsubscribe from this group and stop receiving emails from it, send an email to xarray+un...@googlegroups.com.

[Benjamin Root]

Just a heads-up with regards to np.unravel_index(). It currently raises an exception if given an empty array of indexes to unravel (rather than simply returning an empty array). While that wouldn't be an issue in most situations, it is a stupid edge case that has to be accounted for. I have a PR that actually fixes this, but it got hung up over the question of what dtype the (empty!) return array should be -- I kid you not.

Ben Root

To unsubscribe from this group and stop receiving emails from it, send an email to xarray+unsubscribe@googlegroups.com.

To post to this group, send email to xar...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/xarray/32944fe9-b86f-46a9-b2e2-61936822d8b0%40googlegroups.com.