Coordinate Transformations

[gschulze]

This PR implements all coordinate transformation methods used in GeoTIFF:

• Tie Point and Pixel Scale
• Affine Transform
• Multiple Tie Points (Ground Control Points)

In addition, the implementation considers the raster type (RasterPixelIsArea or RasterPixelIsPoint) to determine the origin of the raster space coordinate system.

The Tie Point and Pixel Scale and Affine Transform methods are relatively straightforward. The Multiple Tie Points transformation roughly works as follows:

• Create a Delaunay triangulation of the tie points in raster space
• Apply the resulting triangulation to the tie points in model space
• When doing coordinate transformations, use the containing triangle of a coordinate to determine the three points to be used for the transformation
• Determine the relative position of the coordinate within the triangle in the model space coordinate system
• Use the relative position to linearly interpolate the coordinate within the corresponding triangle in the raster space coordinate system

@weiji14 please have a look.

[weiji14]

Just a couple of quick comments on the organization of the code. I'll try to find some time to do a deeper dive over the weekend!

[kylebarron]

I'm not familiar with how GCPs work, so if @weiji14 has time to review I'll defer to him

[weiji14]

Did an initial look through the changes and looks good so far! Unit tests are solid (you're lucky I did a exchange in Salzvurg, Austria many years ago, so familiar with the geography 😃), and I do appreciate the changes to get_pixel_value to return Some/None.

Leaving a few small comments or now, I'd be keen on a bit more documentation around the public enums/structs if possible (a link to the relevant subsection in the GeoTIFF spec will do). Will squeeze out some time next week to do a deeper review on the affine and Delaunay triangulation math.

[weiji14]

Ok, I've checked the affine transform code and it looks sound. The Delaunay triangulation logic is a bit over my head, but I trust that it is ok (if not, someone will hopefully notice and file a bug report).

Don't want to hold this back any longer because I'm getting into a busy period, so I'll pre-approve this. I have some questions below around the get_value_at method, but hopefully those are non-issues.

[weiji14]

Just want to double-check the logic here because. If I'm not mistaken, the code is:

1. Transforming model coordinates to raster coordinates via transform_to_raster
2. Applying an offset of -0.5 (in the case of PixelisPoint) to the raster coordinates?

Is it correct to apply an offset of -0.5 (which should be a model offset) to raster coordinates?

[gschulze]

@weiji14 The 0.5 offset refers to half a pixel, which only makes sense in raster space. See https://docs.ogc.org/is/19-008r4/19-008r4.html#_raster_space for reference. This issue also has some nice illustrations: opengeospatial/ogcapi-coverages#92

I added a reference to the corresponding section in the standard in the code.

[gschulze]

If I remember correctly, I also cross-checked with QGIS regarding handling of the raster type.

[weiji14]

Ah ok, I got mixed up on raster space (which should be pixels) and model space (which are real world coordinates in metres or other units). I kept thinking raster offset should be in metres and model offset is in pixels, but it should be the other way around. Thanks for providing the links for clarification, and for adding an inline comment to document this in the code!

[gschulze]

Thanks @weiji14 for your comments. I've been busy the last two weeks, but hopefully, I'll find time to address the remaining issues next week.

[gschulze]

@weiji14 Let me know if you want any further changes; otherwise, I will merge this PR.

[weiji14]

Let's go ahead and merge, thanks for adding the extra docs and bug fixes!