EPSG Helmert and Molodensky-Badekas methods in the Geographic 2D domain: use the push/pop v_3 operator to preserve the Z component

1 files changed, 95 insertions, 0 deletions

diff --git a/src/iso19111/io.cpp b/src/iso19111/io.cpp
index 0e9eac02..6edce579 100644
--- a/src/iso19111/io.cpp
+++ b/src/iso19111/io.cpp
@@ -4837,6 +4837,12 @@ const std::string &PROJStringFormatter::toString() const {
                    step.paramValues[1].keyEquals("xy_out") &&
                    step.paramValues[0].value == step.paramValues[1].value) {
             iter = d->steps_.erase(iter);
+        } else if (step.name == "push" && step.inverted) {
+            step.name = "pop";
+            step.inverted = false;
+        } else if (step.name == "pop" && step.inverted) {
+            step.name = "push";
+            step.inverted = false;
         } else {
             ++iter;
         }
@@ -4912,6 +4918,30 @@ const std::string &PROJStringFormatter::toString() const {
                 break;
             }
 
+            // push v_x followed by pop v_x is a no-op.
+            if (curStep.name == "pop" && prevStep.name == "push" &&
+                !curStep.inverted && !prevStep.inverted &&
+                curStepParamCount == 1 && prevStepParamCount == 1 &&
+                curStep.paramValues[0].key == prevStep.paramValues[0].key) {
+                ++iterCur;
+                d->steps_.erase(iterPrev, iterCur);
+                changeDone = true;
+                break;
+            }
+
+            // pop v_x followed by push v_x is, almost, a no-op. For our
+            // purposes,
+            // we consider it as a no-op for better pipeline optimizations.
+            if (curStep.name == "push" && prevStep.name == "pop" &&
+                !curStep.inverted && !prevStep.inverted &&
+                curStepParamCount == 1 && prevStepParamCount == 1 &&
+                curStep.paramValues[0].key == prevStep.paramValues[0].key) {
+                ++iterCur;
+                d->steps_.erase(iterPrev, iterCur);
+                changeDone = true;
+                break;
+            }
+
             // unitconvert (xy) followed by its inverse is a no-op
             if (curStep.name == "unitconvert" &&
                 prevStep.name == "unitconvert" && !curStep.inverted &&
@@ -5072,6 +5102,71 @@ const std::string &PROJStringFormatter::toString() const {
                 }
             }
 
+            // axisswap order=2,1, pop/push v_3, axisswap order=2,1 -> can
+            // suppress axisswap
+            if (i + 1 < d->steps_.size() && prevStep.name == "axisswap" &&
+                (curStep.name == "push" || curStep.name == "pop") &&
+                prevStepParamCount == 1 &&
+                prevStep.paramValues[0].equals("order", "2,1") &&
+                curStepParamCount == 1 && curStep.paramValues[0].key == "v_3") {
+                auto iterNext = iterCur;
+                ++iterNext;
+                auto &nextStep = *iterNext;
+                if (nextStep.name == "axisswap" &&
+                    nextStep.paramValues.size() == 1 &&
+                    nextStep.paramValues[0].equals("order", "2,1")) {
+                    d->steps_.erase(iterPrev);
+                    d->steps_.erase(iterNext);
+                    changeDone = true;
+                    break;
+                }
+            }
+
+            // push v_3, axisswap order=2,1, pop v_3 -> can suppress push/pop
+            if (i + 1 < d->steps_.size() && prevStep.name == "push" &&
+                prevStepParamCount == 1 &&
+                prevStep.paramValues[0].key == "v_3" &&
+                curStep.name == "axisswap" && curStepParamCount == 1 &&
+                curStep.paramValues[0].equals("order", "2,1")) {
+                auto iterNext = iterCur;
+                ++iterNext;
+                auto &nextStep = *iterNext;
+                if (nextStep.name == "pop" &&
+                    nextStep.paramValues.size() == 1 &&
+                    nextStep.paramValues[0].key == "v_3") {
+                    d->steps_.erase(iterPrev);
+                    d->steps_.erase(iterNext);
+                    changeDone = true;
+                    break;
+                }
+            }
+
+            // unitconvert xy_in=A xy_out=B, pop/push v_3, unitconvert xy_in=B
+            // xy_out=A -> can suppress unitconvert
+            if (i + 1 < d->steps_.size() && prevStep.name == "unitconvert" &&
+                (curStep.name == "push" || curStep.name == "pop") &&
+                prevStepParamCount == 2 &&
+                prevStep.paramValues[0].key == "xy_in" &&
+                prevStep.paramValues[1].key == "xy_out" &&
+                curStepParamCount == 1 && curStep.paramValues[0].key == "v_3") {
+                auto iterNext = iterCur;
+                ++iterNext;
+                auto &nextStep = *iterNext;
+                if (nextStep.name == "unitconvert" &&
+                    nextStep.paramValues.size() == 2 &&
+                    nextStep.paramValues[0].key == "xy_in" &&
+                    nextStep.paramValues[1].key == "xy_out" &&
+                    nextStep.paramValues[0].value ==
+                        prevStep.paramValues[1].value &&
+                    nextStep.paramValues[1].value ==
+                        prevStep.paramValues[0].value) {
+                    d->steps_.erase(iterPrev);
+                    d->steps_.erase(iterNext);
+                    changeDone = true;
+                    break;
+                }
+            }
+
             // for practical purposes WGS84 and GRS80 ellipsoids are
             // equivalents (cartesian transform between both lead to differences
             // of the order of 1e-14 deg..).