From ea73297ea426eac1dcc0133c4cd730ff029e26a8 Mon Sep 17 00:00:00 2001
From: Even Rouault <even.rouault@spatialys.com>
Date: Mon, 30 Dec 2019 00:50:08 +0100
Subject: deformation: add support for +grids= for GeoTIFF grids

This option is to load grid(s) that have both the horizontal and
vertical velocities in the same file.

Can be tested with the following grid
https://github.com/rouault/sample_proj_gtiff_grids/blob/master/nkgrf03vel_realigned.tif
converted from the original .ct2 and .gtx with
```
gdal_translate nkgrf03vel_realigned.vrt nkgrf03vel_realigned.tif -co COMPRESS=DEFLATE -co PREDICTOR=3 -co BLOCKYSIZE=241 -co INTERLEAVE=BAND
```

where nkgrf03vel_realigned.vrt is
```
<VRTDataset rasterXSize="223" rasterYSize="241">
  <SRS dataAxisToSRSAxisMapping="2,1">GEOGCRS["Unknown based on GRS80",
    DATUM["Unknown based on GRS80",
        ELLIPSOID["GRS 1980",6378137,298.257222101,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["geodetic latitude (Lat)",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic longitude (Lon)",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]]]
</SRS>
  <GeoTransform>  2.9166666666666670e+00,  1.6666666666666666e-01,  0.0000000000000000e+00,  7.3041666666666686e+01,  0.0000000000000000e+00, -8.3333333333333329e-02</GeoTransform>
  <Metadata>
    <MDI key="DESCRIPTION"></MDI>
    <MDI key="area_of_use">Nordic and Baltic countries</MDI>
    <MDI key="AREA_OR_POINT">Point</MDI>
    <MDI key="TIFFTAG_COPYRIGHT">The Nordic Geodetic Commission. Creative Commons Attribution 4.0 https://creativecommons.org/licenses/by/4.0/</MDI>
    <MDI key="TIFFTAG_DATETIME">2019:12:30 00:00:00</MDI>
    <MDI key="TIFFTAG_IMAGEDESCRIPTION">Deformation model covering the Nordic and Baltic countries. Used in transformations between global reference frames and the local realisations of ETRS89 in the Nordic and Baltic countries</MDI>
    <MDI key="TYPE">VELOCITY</MDI>
  </Metadata>
  <VRTRasterBand dataType="Float32" band="1">
    <Description>east_velocity</Description>
    <UnitType>mm/year</UnitType>
    <SimpleSource>
      <SourceFilename relativeToVRT="0">/home/even/proj/proj-datumgrid/europe/nkgrf03vel_realigned_xy.ct2</SourceFilename>
      <SourceBand>2</SourceBand> <!-- GDAL exposes the first physical component of the file (longitude offset normally, here east velocity) as the second band -->
      <SrcRect xOff="0" yOff="0" xSize="223" ySize="241" />
      <DstRect xOff="0" yOff="0" xSize="223" ySize="241" />
    </SimpleSource>
  </VRTRasterBand>
  <VRTRasterBand dataType="Float32" band="2">
    <Description>north_velocity</Description>
    <UnitType>mm/year</UnitType>
    <SimpleSource>
      <SourceFilename relativeToVRT="0">/home/even/proj/proj-datumgrid/europe/nkgrf03vel_realigned_xy.ct2</SourceFilename>
      <SourceBand>1</SourceBand> <!-- and the second physical component (latitude offset normally, here north velocity) as the first band -->
      <SrcRect xOff="0" yOff="0" xSize="223" ySize="241" />
      <DstRect xOff="0" yOff="0" xSize="223" ySize="241" />
    </SimpleSource>
  </VRTRasterBand>
  <VRTRasterBand dataType="Float32" band="3">
    <Description>up_velocity</Description>
    <UnitType>mm/year</UnitType>
    <SimpleSource>
      <SourceFilename relativeToVRT="0">/home/even/proj/proj-datumgrid/europe/nkgrf03vel_realigned_z.gtx</SourceFilename>
      <SourceBand>1</SourceBand>
      <SrcRect xOff="0" yOff="0" xSize="223" ySize="241" />
      <DstRect xOff="0" yOff="0" xSize="223" ySize="241" />
    </SimpleSource>
  </VRTRasterBand>
</VRTDataset>
```
---
 test/gie/deformation.gie | 29 +++++++++++++++++++++++++++++
 1 file changed, 29 insertions(+)

(limited to 'test')

diff --git a/test/gie/deformation.gie b/test/gie/deformation.gie
index 37f62e6e..848b9e89 100644
--- a/test/gie/deformation.gie
+++ b/test/gie/deformation.gie
@@ -11,6 +11,35 @@ The input coordinate is located at lon=60, lam=-160 - somewhere in Alaska.
 
 <gie>
 
+-------------------------------------------------------------------------------
+Test with an extract of nkgrf03vel_realigned with ctable2+gtx
+-------------------------------------------------------------------------------
+operation    +proj=pipeline
+             +step +proj=cart +ellps=GRS80
+             +step +proj=deformation
+             +xy_grids=tests/nkgrf03vel_realigned_xy_extract.ct2
+             +z_grids=tests/nkgrf03vel_realigned_z_extract.gtx +ellps=GRS80 +dt=1
+             +step +proj=cart +ellps=GRS80 +inv
+-------------------------------------------------------------------------------
+tolerance   0.05 mm
+accept      21.5            63                   0
+expect      21.5000000049   62.9999999937        0.0083
+roundtrip   5
+
+-------------------------------------------------------------------------------
+Test with an extract of nkgrf03vel_realigned with GeoTIFF
+-------------------------------------------------------------------------------
+operation    +proj=pipeline
+             +step +proj=cart +ellps=GRS80
+             +step +proj=deformation
+             +grids=tests/nkgrf03vel_realigned_extract.tif +ellps=GRS80 +dt=1
+             +step +proj=cart +ellps=GRS80 +inv
+-------------------------------------------------------------------------------
+tolerance   0.05 mm
+accept      21.5            63                   0
+expect      21.5000000049   62.9999999937        0.0083
+roundtrip   5
+
 -------------------------------------------------------------------------------
 Test the +dt parameter
 -------------------------------------------------------------------------------
-- 
cgit v1.2.3