Improve API documentation of GridMapping

[heptaflar]

Describe the bug
The API documentation of GridMapping is not clear in the context of a Coords2DGridMapping: It is not clear what x and y is. It is not clear what i and j is. And it is not clear how i and j are related to x and y. It is not clear what is_j_axis_up means because it is not clear what up means.

To Reproduce
Imagine a satellite in equatorial orbit instead of a polar orbit.

Additional context
Test case was provided to Konstantin.

[konstntokas]

Thank you for pointing this out. I agree that the naming schema in the GridMapping objects does look a bit odd and could benefit from a review.

@forman, do you already foresee any backward compatibility issues for downstream applications if we were to adjust the naming schema? It would be good to understand the potential impact before making changes.

[heptaflar]

Thank you Konstantin. Maybe terms like along track, across track, northward, eastward are quite useful in this context. Also note that the code computing resolution in x and y directions is wrong for equatorial satellite orbits. Anyway, I preferred computation of resolution along these directions being optional. In my use case I know them and there may be other use cases where the resolution is known.

[konstntokas]

Also note that the code computing resolution in x and y directions is wrong for equatorial satellite orbits.

Note that for 2d coordinate grid-mappings, the code just estimates the resolution. For this, it takes the 4 corner points and the number of pixels along each edge and estimates the resolution with:

    x_res = abs(np.mean([x00, x10]) - np.mean([x11, x01])) / (width - 1)
    y_res = abs(np.mean([y00, y01]) - np.mean([y11, y10])) / (height - 1)

@heptaflar If you have a better estimation in mind, please share. However, I want an estimation, which works on only a few pixels so that the grid-mapping generation stays fast.

[konstntokas]

I would then propose the following naming changes in the objects field of GridMapping:

old name	new name
x_coords	eastward_coords
y_coords	northward_coords
xy_coords	spatial_coords
xy_coords_chunks	spatial_coords_chunks
xy_var_names	spatial_var_names
xy_dim_names	spatial_dim_names
xy_bbox	bbox
x_min	eastward_min
y_min	northward_min
x_max	eastward_max
y_max	northward_max
xy_res	spatial_res
x_res	easterward_res
y_res	northward_res
is_j_axis_up	is_northward_axis_up
ij_to_xy_transform	?
xy_to_ij_transform	?
ij_bbox	pixel_bbox
ij_bboxes	pixel_bboxes
xy_bboxes	bboxes

@heptaflar and @forman what are your thoughts on this?

[forman]

I would apply any renaming at all. A gridmapping should remain an abstract concept of mapping a 2D raster to x,y, coordinates, see CF conventions where it has its roots. A gridmapping describes both rectified/regular coordinates and unrectified coordines such as often seen in L1/L2 satellite imagery.

Instead I believe it is better to clearly document the x,y and i,j coordinates systems and specify all used terms, especially what up means. Also always refer to CF.

[konstntokas]

I enhanced the GridMapping API and linked to the CF convention: https://xcube-dev.github.io/xcube-resampling/api/#xcube_resampling.gridmapping.GridMapping

[heptaflar]

Also note that the code computing resolution in x and y directions is wrong for equatorial satellite orbits.

Note that for 2d coordinate grid-mappings, the code just estimates the resolution. For this, it takes the 4 corner points and the number of pixels along each edge and estimates the resolution with:

    x_res = abs(np.mean([x00, x10]) - np.mean([x11, x01])) / (width - 1)
    y_res = abs(np.mean([y00, y01]) - np.mean([y11, y10])) / (height - 1)

@heptaflar If you have a better estimation in mind, please share. However, I want an estimation, which works on only a few pixels so that the grid-mapping generation stays fast.

Dear @konstntokas ,

I do not have a better estimation in mind. Your quick estimation is OK for points 0_ and 1_ making a step along the y axis and points _0 and _1 making a step along the x axis. This is the usual case for polar orbit imagery.

But for equatorial imagery the y axis is along the equator (longitudinal) and the x axis is across the equator (latitudinal). Then the code above does not work anymore. Instead, you need to compute:

x_res = abs(np.mean([x00, x01]) - np.mean([x11, x10])) / (width - 1)
y_res = abs(np.mean([y00, y10]) - np.mean([y11, y01])) / (height - 1)

where width is the latitudinal dimension size (i.e., across track as in the polar orbit case) and height is the longitudinal dimension size (i.e., along track as in the polar orbit case). I did not check the code, just guessed it.

Use case: I rotate a polar orbit into an equatorial orbit centered on the prime meridian by 3D Euler rotations to minimize distortions by the equirectangular projection of Earth's spherical geometry and to avoid problems with poles and anti-meridian in the image.

[heptaflar]

I would then propose the following naming changes in the objects field of GridMapping:

old name new name
x_coords eastward_coords
y_coords northward_coords
xy_coords spatial_coords
xy_coords_chunks spatial_coords_chunks
xy_var_names spatial_var_names
xy_dim_names spatial_dim_names
xy_bbox bbox
x_min eastward_min
y_min northward_min
x_max eastward_max
y_max northward_max
xy_res spatial_res
x_res easterward_res
y_res northward_res
is_j_axis_up is_northward_axis_up
ij_to_xy_transform ?
xy_to_ij_transform ?
ij_bbox pixel_bbox
ij_bboxes pixel_bboxes
xy_bboxes bboxes
@heptaflar and @forman what are your thoughts on this?

I am happy with @forman 's advice. My concern is that ij and xy do not necessarily have the same (or a similar) orientation. They can be perpendicular or reflected, in principle. Especially in this case terms like north, east, south, west and up are quite confusing.