A number of improvements to the grid handling

[busstoptaktik]

Two interpolation tests in ntv2/mod.rs (succesfully validated against ESRI's ntv2_cvt program)

Improve check for NTv2 file header sanity

Improve readability and correct vector capacity in ntv2 parse_subgrid_grid()

Reduce visibility of SubGridHeader

Clean up the use block in ntv2 subgrid.rs

In the inverse case, don't iterate in reverse order over the grid collection in deformation & gridshift operators (we still want to hit the same grid as in the forward case)

[busstoptaktik]

@Rennzie does this look OK to you?

[Rennzie]

🤩 Looks great @busstoptaktik !

[Rennzie]

@busstoptaktik I've integrated this Ntv2 support into geodesy-wasm now. All is working smoothly 🎉 .

One thing I've noticed in my tests is the diff between RG and PROJ when converting an EPSG:27700 coordinate to EPSG:3857. See below:

   const definition = `
  | tmerc inv lat_0=49 lon_0=-2 k_0=0.9996012717 x_0=400000 y_0=-100000 ellps=airy
  | gridshift grids=OSTN15_NTv2_OSGBtoETRS.gsb
  | webmerc lat_0=0 lon_0=0 x_0=0 y_0=0 ellps=WGS84
  `;

  const centralCambridgeCoords = [
    [544748.5367636156, 258372.49178149243, 9.61],
    [544750.8907704452, 258365.94195330486, 9.61],
  ];

 // Results from: echo 544748.5367636156 258372.49178149243 9.61 | cct +proj=pipeline  +step +inv +proj=tmerc +lat_0=49 +lon_0=-2 +k_0=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +step +proj=hgridshift +grids=OSTN15_NTv2_OSGBtoETRS.gsb +step +proj=webmerc +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84
  const expected = [
    [13004.3086, 6837202.7637, 9.61],
    [13007.8289, 6837191.9623, 9.61],
  ];

// Forward:
13004.3041 6837202.7506 9.6100
13007.8245 6837191.9494 9.6100

// Diff Expected:
-0.004527604605755187 -0.013138627633452415 0
-0.0044116277422290295 -0.01289624534547329 0

I'm curious to know if the precision is about what you'd expect given the off by half differences you'd noticed between the PROJ interpolation that of RG?

Incase you're curious you can find OSTN15_NTv2_OSGBtoETRS.gsb used here

[busstoptaktik]

I'm curious to know if the precision is about what you'd expect given the off by half differences you'd noticed between the PROJ interpolation that of RG?

@Rennzie no, these differences are probably too large, and I'm staring myself to blindness today, checking interpolation code between ESRI's ntv2_cvt, PROJ's hgridshift and Rust Geodesy's gridshift.

PROJ and ESRI agree on the tenth to hundredth of a millimeter, while RG is sightly off, only agreeing with the other two on the millimeter level.

I have even reimplemented ESRI's code in RG, with next to no difference, so I assume the problem is in the calculation of the cell relative coordinates.

I'm looking into it, but it is a nasty problem, not the least because the grid is only single precision, so where and when you convert to f64 in the process will also make a difference.

My comparison routine is this:

$ echo 2.1686 41.3874 0 0 | cct -d 5 proj=pipeline step proj=hgridshift grids=./path/to/100800401.gsb step proj=utm zone=31 --

> 430391.39938  4581973.06191

$ ntv2_cvt -f .\path\to\100800401.gsb  41.3874 2.1686 | proj -f %9.5f +proj=utm +zone=31

> 430391.39938    4581973.06190

$ echo  41.3874 2.1686 | cargo run -- "geo:in | gridshift grids=100800401.gsb | utm zone=31"

> 430391.39952 4581973.06152

In summary:

PROJ: 430391.39938  4581973.06191
ESRI: 430391.39938  4581973.06190
RG:   430391.39952  4581973.06152

[Rennzie]

Thanks for looking into it. Please let me know if I can help in any way.