+proj=push/pop: add a +bank=<name> parameter to interleave push/pop o…

…perations

The motivation is to be able to do transformation that first do a geoid
transformation followed by a Helmert one.

The following example shows a ETRS89 to S-JTSK/05 + Baltic 1957 height
transformation that chains a geoid model registered against ETRS89 with a
Helmert transformation.

```
       +proj=pipeline \
            +step +proj=push +v_3 +omit_inv \                         # Save ETRS89 ellipsoidal height
            +step +inv +proj=vgridshift +grids=CR2005.tif +multiplier=1 \  # Apply geoid
            +step +proj=push +v_3 +bank=baltic_height \               # Save Baltic 1957 height
            +step +proj=pop +v_3 +omit_inv \                          # Restore ETRS89 ellipsoidal height
            +step +proj=cart +ellps=GRS80 \                           # Helmert transformation
            +step +inv +proj=helmert +x=572.213 +y=85.334 +z=461.94 \
                +rx=-4.9732 +ry=-1.529 +rz=-5.2484 +s=3.5378 +convention=coordinate_frame \
            +step +inv +proj=cart +ellps=bessel \
            +step +proj=pop +v_3 +bank=baltic_height                  # Restore Baltic 1957 height
```

Without that, the current alternative is a ugly workaround involving
using the v_4 component and axisswap
```
       +proj=pipeline \
            +step +proj=push +v_4 \                                   # Save v_4 as we are going to use it as a temp variable
            +step +proj=push +v_3 +omit_inv \                         # Save ETRS89 ellipsoidal height
            +step +inv +proj=vgridshift +grids=CR2005.tif +multiplier=1 \  # Apply geoid
            +step +proj=push +v_3 +omit_fwd \                         # On reverse path, restore initial Baltic height
            +step +proj=axisswap +order=1,2,4,3 +omit_inv \           # On forward path, save Baltic height in v_4 component...
            +step +proj=pop +v_3 +omit_inv \                          # On forward parth, restore initial ellipsoidal height
            +step +proj=cart +ellps=GRS80 \                           # Helmert transformation
            +step +inv +proj=helmert +x=572.213 +y=85.334 +z=461.94 \
                +rx=-4.9732 +ry=-1.529 +rz=-5.2484 +s=3.5378 +convention=coordinate_frame \
            +step +inv +proj=cart +ellps=bessel \
            +step +proj=axisswap +order=1,2,4,3 +omit_inv \           # On forward path, restore Baltic height from v_4 component...
            +step +proj=pop +v_3 +omit_fwd \                          # On reverse path, save initial Baltic height
            +step +proj=pop +v_4
```

docs/source/operations/conversions/pop.rst

@@ -86,6 +86,13 @@ Parameters
 
    Retrieves the fourth coordinate component from the pipeline stack
 
+.. option:: +bank=<bank_name>
+
+   .. versionadded:: 9.4.0
+
+   Restore the above components from a named "bank" instead of the pipeline
+   stack. This must be paired with a corresponding :ref:`push <push>` with the same bank name.
+
 
 Further reading
 ################################################################################

docs/source/operations/conversions/push.rst

@@ -85,6 +85,32 @@ Parameters
 
    Stores the fourth coordinate component on the pipeline stack
 
+.. option:: +bank=<bank_name>
+
+   .. versionadded:: 9.4.0
+
+   Store the above saved components into a named "bank" instead of the pipeline
+   stack. They will be restored when using :ref:`pop <pop>` with the same bank name.
+   Push/pop operations with bank names can be interleaved with other push/pop
+   operations with or without bank names.
+
+   The following example shows a ETRS89 to S-JTSK/05 + Baltic 1957 height
+   transformation that chains a geoid model registered against ETRS89 with a
+   Helmert transformation.
+
+   ::
+
+       +proj=pipeline \
+            +step +proj=push +v_3 +omit_inv \                               # Save ETRS89 ellipsoidal height
+            +step +inv +proj=vgridshift +grids=CR2005.tif +multiplier=1 \   # Apply geoid
+            +step +proj=push +v_3 +bank=baltic_height \                     # Save Baltic 1957 height
+            +step +proj=pop +v_3 +omit_inv \                                # Restore ETRS89 ellipsoidal height
+            +step +proj=cart +ellps=GRS80 \                                 # Helmert transformation
+            +step +inv +proj=helmert +x=572.213 +y=85.334 +z=461.94 \
+                +rx=-4.9732 +ry=-1.529 +rz=-5.2484 +s=3.5378 +convention=coordinate_frame \
+            +step +inv +proj=cart +ellps=bessel \
+            +step +proj=pop +v_3 +bank=baltic_height                        # Restore Baltic 1957 height
+
 
 Further reading
 ################################################################################

src/pipeline.cpp

@@ -97,6 +97,8 @@ Thomas Knudsen, thokn@sdfe.dk, 2016-05-20
 *
 ********************************************************************************/
 
+#include <array>
+#include <map>
 #include <math.h>
 #include <stack>
 #include <stddef.h>
@@ -137,13 +139,17 @@ struct Pipeline {
     char **current_argv = nullptr;
     std::vector<Step> steps{};
     std::stack<double> stack[4];
+    std::map<std::string, int> bank_names_to_idx{};
+    std::vector<std::array<double, 4>> bank{};
 };
 
 struct PushPop {
-    bool v1;
-    bool v2;
-    bool v3;
-    bool v4;
+    std::string bank_name{};
+    int bank_index = -1;
+    bool v1 = false;
+    bool v2 = false;
+    bool v3 = false;
+    bool v4 = false;
 };
 } // anonymous namespace
 
@@ -642,55 +648,112 @@ static void push(PJ_COORD &point, PJ *P) {
     if (P->parent == nullptr)
         return;
 
-    struct Pipeline *pipeline =
-        static_cast<struct Pipeline *>(P->parent->opaque);
-    struct PushPop *pushpop = static_cast<struct PushPop *>(P->opaque);
-
-    if (pushpop->v1)
-        pipeline->stack[0].push(point.v[0]);
-    if (pushpop->v2)
-        pipeline->stack[1].push(point.v[1]);
-    if (pushpop->v3)
-        pipeline->stack[2].push(point.v[2]);
-    if (pushpop->v4)
-        pipeline->stack[3].push(point.v[3]);
+    Pipeline *pipeline = static_cast<struct Pipeline *>(P->parent->opaque);
+    PushPop *pushpop = static_cast<struct PushPop *>(P->opaque);
+
+    if (!pushpop->bank_name.empty() && pushpop->bank_index < 0) {
+        // Check if we already registered this bank name
+        const auto iter = pipeline->bank_names_to_idx.find(pushpop->bank_name);
+        if (iter != pipeline->bank_names_to_idx.end()) {
+            // Yes, then store its known index, so we don't have to fetch
+            // it anymore
+            pushpop->bank_index = iter->second;
+        } else {
+            // Unkown bank name. Add a new one, and register it in the pipeline.
+            const int bank_index = static_cast<int>(pipeline->bank.size());
+            pipeline->bank.resize(pipeline->bank.size() + 1);
+            pushpop->bank_index = bank_index;
+            pipeline->bank_names_to_idx[pushpop->bank_name] = bank_index;
+        }
+    }
+    if (pushpop->bank_index >= 0) {
+        if (pushpop->v1)
+            pipeline->bank[pushpop->bank_index][0] = point.v[0];
+        if (pushpop->v2)
+            pipeline->bank[pushpop->bank_index][1] = point.v[1];
+        if (pushpop->v3)
+            pipeline->bank[pushpop->bank_index][2] = point.v[2];
+        if (pushpop->v4)
+            pipeline->bank[pushpop->bank_index][3] = point.v[3];
+    } else {
+        if (pushpop->v1)
+            pipeline->stack[0].push(point.v[0]);
+        if (pushpop->v2)
+            pipeline->stack[1].push(point.v[1]);
+        if (pushpop->v3)
+            pipeline->stack[2].push(point.v[2]);
+        if (pushpop->v4)
+            pipeline->stack[3].push(point.v[3]);
+    }
 }
 
 static void pop(PJ_COORD &point, PJ *P) {
     if (P->parent == nullptr)
         return;
 
-    struct Pipeline *pipeline =
-        static_cast<struct Pipeline *>(P->parent->opaque);
-    struct PushPop *pushpop = static_cast<struct PushPop *>(P->opaque);
+    Pipeline *pipeline = static_cast<struct Pipeline *>(P->parent->opaque);
+    PushPop *pushpop = static_cast<struct PushPop *>(P->opaque);
 
-    if (pushpop->v1 && !pipeline->stack[0].empty()) {
-        point.v[0] = pipeline->stack[0].top();
-        pipeline->stack[0].pop();
+    if (!pushpop->bank_name.empty() && pushpop->bank_index < 0) {
+        const auto iter = pipeline->bank_names_to_idx.find(pushpop->bank_name);
+        if (iter != pipeline->bank_names_to_idx.end()) {
+            pushpop->bank_index = iter->second;
+        } else {
+            proj_context_errno_set(P->ctx,
+                                   PROJ_ERR_INVALID_OP_ILLEGAL_ARG_VALUE);
+            proj_log_error(P, _("Pipeline: pop: Bank %s is unknown"),
+                           pushpop->bank_name.c_str());
+            point = proj_coord_error();
+            return;
+        }
     }
 
-    if (pushpop->v2 && !pipeline->stack[1].empty()) {
-        point.v[1] = pipeline->stack[1].top();
-        pipeline->stack[1].pop();
-    }
+    if (pushpop->bank_index >= 0) {
+        if (pushpop->v1)
+            point.v[0] = pipeline->bank[pushpop->bank_index][0];
+        if (pushpop->v2)
+            point.v[1] = pipeline->bank[pushpop->bank_index][1];
+        if (pushpop->v3)
+            point.v[2] = pipeline->bank[pushpop->bank_index][2];
+        if (pushpop->v4)
+            point.v[3] = pipeline->bank[pushpop->bank_index][3];
+    } else {
+        if (pushpop->v1 && !pipeline->stack[0].empty()) {
+            point.v[0] = pipeline->stack[0].top();
+            pipeline->stack[0].pop();
+        }
 
-    if (pushpop->v3 && !pipeline->stack[2].empty()) {
-        point.v[2] = pipeline->stack[2].top();
-        pipeline->stack[2].pop();
-    }
+        if (pushpop->v2 && !pipeline->stack[1].empty()) {
+            point.v[1] = pipeline->stack[1].top();
+            pipeline->stack[1].pop();
+        }
+
+        if (pushpop->v3 && !pipeline->stack[2].empty()) {
+            point.v[2] = pipeline->stack[2].top();
+            pipeline->stack[2].pop();
+        }
 
-    if (pushpop->v4 && !pipeline->stack[3].empty()) {
-        point.v[3] = pipeline->stack[3].top();
-        pipeline->stack[3].pop();
+        if (pushpop->v4 && !pipeline->stack[3].empty()) {
+            point.v[3] = pipeline->stack[3].top();
+            pipeline->stack[3].pop();
+        }
     }
 }
 
+static PJ *pj_pushpop_destructor(PJ *P, int errlev) {
+    if (nullptr == P)
+        return nullptr;
+
+    delete static_cast<PushPop *>(P->opaque);
+    P->opaque = nullptr;
+
+    return pj_default_destructor(P, errlev);
+}
+
 static PJ *setup_pushpop(PJ *P) {
-    auto pushpop =
-        static_cast<struct PushPop *>(calloc(1, sizeof(struct PushPop)));
+    auto pushpop = new PushPop;
     P->opaque = pushpop;
-    if (nullptr == P->opaque)
-        return destructor(P, PROJ_ERR_OTHER /*ENOMEM*/);
+    P->destructor = pj_pushpop_destructor;
 
     if (pj_param_exists(P->params, "v_1"))
         pushpop->v1 = true;
@@ -704,6 +767,10 @@ static PJ *setup_pushpop(PJ *P) {
     if (pj_param_exists(P->params, "v_4"))
         pushpop->v4 = true;
 
+    if (pj_param(P->ctx, P->params, "tbank").i) {
+        pushpop->bank_name = pj_param(P->ctx, P->params, "sbank").s;
+    }
+
     P->left = PJ_IO_UNITS_WHATEVER;
     P->right = PJ_IO_UNITS_WHATEVER;
 

test/gie/4D-API_cs2cs-style.gie

@@ -393,6 +393,46 @@ operation   +proj=pop +v_3
 accept      12 56 0 0
 expect      12 56 0 0
 
+-------------------------------------------------------------------------------
+# Test push/pop named banks
+-------------------------------------------------------------------------------
+
+# Simulate ETRS89 to S-JTSK/05 + Baltic 1957 height
+operation  +proj=pipeline \
+            +step +proj=push +v_3 +omit_inv \
+            +step +proj=affine +zoff=-50 \  # Simulate geoid grid
+            +step +proj=push +v_3 +bank=my_bank \
+            +step +proj=pop +v_3 +omit_inv \
+            +step +proj=cart +ellps=GRS80 \
+            +step +inv +proj=helmert +x=572.213 +y=85.334 +z=461.94 \
+                +rx=-4.9732 +ry=-1.529 +rz=-5.2484 +s=3.5378 +convention=coordinate_frame \
+            +step +inv +proj=cart +ellps=bessel \
+            +step +proj=pop +v_3 +bank=my_bank
+tolerance   1 mm
+
+accept      15              50                 1000
+expect      15.0011679017   50.0007533757      950.0000
+roundtrip   1
+
+# Several banks
+operation  +proj=pipeline \
+            +step +proj=push +v_1 +v_3 +bank=bank1 \
+            +step +proj=set +v_3=1000 \
+            +step +proj=push +v_1 +v_3 +bank=bank2 \
+            +step +proj=set +v_1=10 \
+            +step +proj=pop +v_1 +v_3 +bank=bank1 \
+            +step +proj=pop +v_3 +bank=bank2
+accept      1 2 3
+expect      1 2 1000
+
+
+# Error case: wrong bank name in pop
+operation  +proj=pipeline \
+            +step +proj=push +v_3 +bank=some_bank \
+            +step +proj=pop +v_3 +bank=unknown_bank
+accept      0 0 0
+expect      failure errno invalid_op_illegal_arg_value
+
 -------------------------------------------------------------------------------
 # Test Pipeline +omit_inv
 -------------------------------------------------------------------------------