diff --git a/src/pipeline.cpp b/src/pipeline.cpp index c4a9d2c048..07aa90d323 100644 --- a/src/pipeline.cpp +++ b/src/pipeline.cpp @@ -100,6 +100,7 @@ Thomas Knudsen, thokn@sdfe.dk, 2016-05-20 #include #include #include +#include #include "geodesic.h" #include "proj.h" @@ -115,6 +116,8 @@ struct pj_opaque { char **argv; char **current_argv; PJ **pipeline; + bool *push[4]; + bool *pop[4]; }; } // anonymous namespace @@ -132,26 +135,59 @@ static PJ_LP pipeline_reverse (PJ_XY xy, PJ *P); static PJ_COORD pipeline_forward_4d (PJ_COORD point, PJ *P) { int i, first_step, last_step; + pj_opaque *opaque = static_cast(P->opaque); + std::stack stack[4]; first_step = 1; - last_step = static_cast(P->opaque)->steps + 1; + last_step = opaque->steps + 1; - for (i = first_step; i != last_step; i++) - point = proj_trans (static_cast(P->opaque)->pipeline[i], PJ_FWD, point); + for (i = first_step; i != last_step; i++) { + for (int j=0; j<4; j++) { + if (opaque->push[j][i-1]) { + stack[j].push(point.v[j]); + } + } + + point = proj_trans (opaque->pipeline[i], PJ_FWD, point); + + for (int j=0; j<4; j++) { + if (opaque->pop[j][i-1] && stack[j].size() > 0) { + point.v[j] = stack[j].top(); + stack[j].pop(); + } + } + } return point; } static PJ_COORD pipeline_reverse_4d (PJ_COORD point, PJ *P) { int i, first_step, last_step; + pj_opaque *opaque = static_cast(P->opaque); + std::stack stack[4]; first_step = static_cast(P->opaque)->steps; last_step = 0; - for (i = first_step; i != last_step; i--) + for (i = first_step; i != last_step; i--) { + /* pop first, which in reverse is push */ + for (int j=0; j<4; j++) { + if (opaque->pop[j][i-1]) { + stack[j].push(point.v[j]); + } + } + point = proj_trans (static_cast(P->opaque)->pipeline[i], PJ_INV, point); + /* then push, which in reverse is pop */ + for (int j=0; j<4; j++) { + if (opaque->push[j][i-1] && stack[j].size() > 0) { + point.v[j] = stack[j].top(); + stack[j].pop(); + } + } + } return point; } @@ -231,13 +267,22 @@ static PJ *destructor (PJ *P, int errlev) { static PJ *pj_create_pipeline (PJ *P, size_t steps) { - + pj_opaque *opaque = static_cast(P->opaque); /* Room for the pipeline: An array of PJ * with room for sentinels at both ends */ - static_cast(P->opaque)->pipeline = static_cast(pj_calloc (steps + 2, sizeof(PJ *))); - if (nullptr==static_cast(P->opaque)->pipeline) + opaque->pipeline = static_cast(pj_calloc (steps + 2, sizeof(PJ *))); + if (nullptr==opaque->pipeline) return nullptr; - static_cast(P->opaque)->steps = (int)steps; + opaque->steps = (int)steps; + + for (int i=0; i<4; i++) { + opaque->push[i] = static_cast(pj_calloc(steps, sizeof(bool))); + opaque->pop[i] = static_cast(pj_calloc(steps, sizeof(bool))); + for (int j=0; j<(int)steps; j++) { + opaque->push[i][j] = false; + opaque->pop[i][j] = false; + } + } return P; } @@ -427,8 +472,10 @@ PJ *OPERATION(pipeline,0) { if (nullptr==pj_create_pipeline (P, nsteps)) return destructor (P, ENOMEM); + pj_opaque *opaque = static_cast(P->opaque); set_ellipsoid(P); + /* Now loop over all steps, building a new set of arguments for each init */ i_current_step = i_first_step; for (i = 0; i < nsteps; i++) { @@ -470,12 +517,23 @@ PJ *OPERATION(pipeline,0) { proj_errno_restore (P, err); - /* Is this step inverted? */ - for (j = 0; j < current_argc; j++) + /* take care of pipeline specific parameters */ + for (j = 0; j < current_argc; j++) { + /* Is this step inverted? */ if (0==strcmp("inv", current_argv[j])) { /* if +inv exists in both global and local args the forward operation should be used */ next_step->inverted = next_step->inverted == 0 ? 1 : 0; } + /* do we need to push or pop something on the coordinate stacks? */ + if (strcmp("push=1", current_argv[j])==0) opaque->push[0][i] = true; + if (strcmp("push=2", current_argv[j])==0) opaque->push[1][i] = true; + if (strcmp("push=3", current_argv[j])==0) opaque->push[2][i] = true; + if (strcmp("push=4", current_argv[j])==0) opaque->push[3][i] = true; + if (strcmp("pop=1", current_argv[j])==0) opaque->pop[0][i] = true; + if (strcmp("pop=2", current_argv[j])==0) opaque->pop[1][i] = true; + if (strcmp("pop=3", current_argv[j])==0) opaque->pop[2][i] = true; + if (strcmp("pop=4", current_argv[j])==0) opaque->pop[3][i] = true; + } static_cast(P->opaque)->pipeline[i+1] = next_step; diff --git a/test/gie/4D-API_cs2cs-style.gie b/test/gie/4D-API_cs2cs-style.gie index 215971a0af..bb98b25eca 100644 --- a/test/gie/4D-API_cs2cs-style.gie +++ b/test/gie/4D-API_cs2cs-style.gie @@ -294,6 +294,35 @@ accept 12 56 expect 1335.8339 7522.963 ------------------------------------------------------------------------------- +------------------------------------------------------------------------------- +Test Pipeline Coordinate Stack +------------------------------------------------------------------------------- +operation +proj=pipeline + +step +proj=utm +zone=32 +push=1 + +step +proj=utm +zone=33 +pop=1 +inv + +accept 12 56 0 2020 +expect 12 56 0 2020 +roundtrip 10 + +operation +proj=pipeline + +step +proj=latlon # dummy step + +step +proj=utm +zone=32 +push=1 + +step +proj=utm +zone=33 +pop=1 +inv + +step +proj=affine # dummy step + +accept 12 56 0 2020 +expect 12 56 0 2020 +roundtrip 10 + +operation +proj=pipeline + +step +inv +proj=eqearth ++push=2 + +step +proj=laea +pop=2 + +accept 900000 6000000 0 2020 +expect 896633.0226 6000000 0 2020 + +------------------------------------------------------------------------------- ------------------------------------------------------------------------------- Test bugfix of https://github.com/OSGeo/proj.4/issues/1002