Coordinate contains check

[groutr]

Implements a useful utility function for checking of a coordinate is contained within the bounding box of a raster. The use case that prompted this was testing to see which rasters covered a specific location(s). Testing this before reading seemed a much simpler route than checking to see if the read pixel was a nodata value.

It's a small utility function designed to make life a little more simple and the coord module seemed like the perfect place for this to live.

[sgillies]

@groutr how about a different numerical approach? If we apply the inverse of the raster's transform to the x, y coords we get column and row indices, right? Then we can compare these to the raster's width and height. It would account for rotation, too, unless I'm mistaken.

[groutr]

@sgillies Using the transform is definitely a more precise method. In my use case, I wasn't dealing with rotated rasters and an approximate check was more than good enough. I was attempting to stay consistent with the other functions in coords.py by not requiring a transform.

I guess the question to decide is do we want an approximate check using just the bounds, or a precise check using the transform. I don't have a strong opinion either way. The goal is to avoid unnecessary reads.

[sgillies]

@groutr I'm seeing more and more rotated rasters at my new day job. Umbra's open data, for example, is generally rotated.

# rio info "https://umbra-open-data-catalog.s3.amazonaws.com/sar-data/tasks/Cal%20Poly,%20San%20Luis%20Obispo,%20CA/7ab8228e-42b0-47b9-ab45-c7d6fe96c200/2023-12-01-19-20-51_UMBRA-08/2023-12-01-19-20-51_UMBRA-08_GEC.tif" | jq '.transform'
[
  0.1443650591189317,
  0.6522549621710828,
  710240.0900495084,
  0.6522549621710314,
  -0.14436505911896663,
  3907811.078473875,
  0.0,
  0.0,
  1.0
]

I think it would be good if new API features were rotation correct.

[groutr]

I can make the check use transforms. Should it live in the transform.py module or stay in coords.py?

Also, I'm trying to dig into why rowcol transforms, especially for scalar coordinates seem slow and if there is a way to speed it up.

On my system

# t = identity transform
%timeit t * (-77, 84)
459 ns ± 4.65 ns per loop (mean ± std. dev. of 7 runs, 1,000,000 loops each)

%timeit riotransform.rowcol(t, -77, 84)
23.1 µs ± 295 ns per loop (mean ± std. dev. of 7 runs, 10,000 loops each)

A = riotransform.AffineTransform(t)
%timeit A.rowcol(-77, 84)
20.7 µs ± 219 ns per loop (mean ± std. dev. of 7 runs, 10,000 loops each)

%timeit A._transform([-77], [84], 0, transform_direction=riotransform.TransformDirection.reverse)
2.52 µs ± 13.7 ns per loop (mean ± std. dev. of 7 runs, 100,000 loops each)

EDIT: I'm aware that for a contains check, I can just use the transform directly to avoid all the overhead of the rowcol function. I'm curious where the overhead is coming from though.

[groutr]

In looking at this some more, I realized that TransformMethodsMixin has an .index() method which will return the the row, col coordinates of a point. Why this isn't called rowcol to be consistent with the rest of transform module is a discussion for later...

I think a good solution here might be an additional method, TransformMethodsMixin.contains which calls .index(), does the bounds checking, and returns a boolean. This approach should be transparent to all supported transforms.

[sgillies]

@groutr no, let's not extend TransformMethodsMixin.