Coverage Report

Created: 2025-11-24 06:58

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/PROJ/src/pipeline.cpp
Line
Count
Source
1
/*******************************************************************************
2
3
                       Transformation pipeline manager
4
5
                    Thomas Knudsen, 2016-05-20/2016-11-20
6
7
********************************************************************************
8
9
    Geodetic transformations are typically organized in a number of
10
    steps. For example, a datum shift could be carried out through
11
    these steps:
12
13
    1. Convert (latitude, longitude, ellipsoidal height) to
14
       3D geocentric cartesian coordinates (X, Y, Z)
15
    2. Transform the (X, Y, Z) coordinates to the new datum, using a
16
       7 parameter Helmert transformation.
17
    3. Convert (X, Y, Z) back to (latitude, longitude, ellipsoidal height)
18
19
    If the height system used is orthometric, rather than ellipsoidal,
20
    another step is needed at each end of the process:
21
22
    1. Add the local geoid undulation (N) to the orthometric height
23
       to obtain the ellipsoidal (i.e. geometric) height.
24
    2. Convert (latitude, longitude, ellipsoidal height) to
25
       3D geocentric cartesian coordinates (X, Y, Z)
26
    3. Transform the (X, Y, Z) coordinates to the new datum, using a
27
       7 parameter Helmert transformation.
28
    4. Convert (X, Y, Z) back to (latitude, longitude, ellipsoidal height)
29
    5. Subtract the local geoid undulation (N) from the ellipsoidal height
30
       to obtain the orthometric height.
31
32
    Additional steps can be added for e.g. change of vertical datum, so the
33
    list can grow fairly long. None of the steps are, however, particularly
34
    complex, and data flow is strictly from top to bottom.
35
36
    Hence, in principle, the first example above could be implemented using
37
    Unix pipelines:
38
39
    cat my_coordinates | geographic_to_xyz | helmert | xyz_to_geographic >
40
my_transformed_coordinates
41
42
    in the grand tradition of Software Tools [1].
43
44
    The proj pipeline driver implements a similar concept: Stringing together
45
    a number of steps, feeding the output of one step to the input of the next.
46
47
    It is a very powerful concept, that increases the range of relevance of the
48
    proj.4 system substantially. It is, however, not a particularly intrusive
49
    addition to the PROJ.4 code base: The implementation is by and large
50
completed by adding an extra projection called "pipeline" (i.e. this file),
51
which handles all business, and a small amount of added functionality in the
52
    pj_init code, implementing support for multilevel, embedded pipelines.
53
54
    Syntactically, the pipeline system introduces the "+step" keyword (which
55
    indicates the start of each transformation step), and reintroduces the +inv
56
    keyword (indicating that a given transformation step should run in reverse,
57
i.e. forward, when the pipeline is executed in inverse direction, and vice
58
versa).
59
60
    Hence, the first transformation example above, can be implemented as:
61
62
    +proj=pipeline +step proj=cart +step proj=helmert <ARGS> +step proj=cart
63
+inv
64
65
    Where <ARGS> indicate the Helmert arguments: 3 translations (+x=..., +y=...,
66
    +z=...), 3 rotations (+rx=..., +ry=..., +rz=...) and a scale factor
67
(+s=...). Following geodetic conventions, the rotations are given in arcseconds,
68
    and the scale factor is given as parts-per-million.
69
70
    [1] B. W. Kernighan & P. J. Plauger: Software tools.
71
        Reading, Massachusetts, Addison-Wesley, 1976, 338 pp.
72
73
********************************************************************************
74
75
Thomas Knudsen, thokn@sdfe.dk, 2016-05-20
76
77
********************************************************************************
78
* Copyright (c) 2016, 2017, 2018 Thomas Knudsen / SDFE
79
*
80
* Permission is hereby granted, free of charge, to any person obtaining a
81
* copy of this software and associated documentation files (the "Software"),
82
* to deal in the Software without restriction, including without limitation
83
* the rights to use, copy, modify, merge, publish, distribute, sublicense,
84
* and/or sell copies of the Software, and to permit persons to whom the
85
* Software is furnished to do so, subject to the following conditions:
86
*
87
* The above copyright notice and this permission notice shall be included
88
* in all copies or substantial portions of the Software.
89
*
90
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
91
* OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
92
* FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
93
* THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
94
* LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
95
* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
96
* DEALINGS IN THE SOFTWARE.
97
*
98
********************************************************************************/
99
100
#include <math.h>
101
#include <stack>
102
#include <stddef.h>
103
#include <string.h>
104
#include <vector>
105
106
#include "geodesic.h"
107
#include "proj.h"
108
#include "proj_internal.h"
109
110
PROJ_HEAD(pipeline, "Transformation pipeline manager");
111
PROJ_HEAD(pop, "Retrieve coordinate value from pipeline stack");
112
PROJ_HEAD(push, "Save coordinate value on pipeline stack");
113
114
/* Projection specific elements for the PJ object */
115
namespace { // anonymous namespace
116
117
struct Step {
118
    PJ *pj = nullptr;
119
    bool omit_fwd = false;
120
    bool omit_inv = false;
121
122
    Step(PJ *pjIn, bool omitFwdIn, bool omitInvIn)
123
78.5k
        : pj(pjIn), omit_fwd(omitFwdIn), omit_inv(omitInvIn) {}
124
    Step(Step &&other)
125
80.5k
        : pj(std::move(other.pj)), omit_fwd(other.omit_fwd),
126
80.5k
          omit_inv(other.omit_inv) {
127
80.5k
        other.pj = nullptr;
128
80.5k
    }
129
    Step(const Step &) = delete;
130
    Step &operator=(const Step &) = delete;
131
132
159k
    ~Step() { proj_destroy(pj); }
133
};
134
135
struct Pipeline {
136
    char **argv = nullptr;
137
    char **current_argv = nullptr;
138
    std::vector<Step> steps{};
139
    std::stack<double> stack[4];
140
};
141
142
struct PushPop {
143
    bool v1;
144
    bool v2;
145
    bool v3;
146
    bool v4;
147
};
148
} // anonymous namespace
149
150
static void pipeline_forward_4d(PJ_COORD &point, PJ *P);
151
static void pipeline_reverse_4d(PJ_COORD &point, PJ *P);
152
static PJ_XYZ pipeline_forward_3d(PJ_LPZ lpz, PJ *P);
153
static PJ_LPZ pipeline_reverse_3d(PJ_XYZ xyz, PJ *P);
154
static PJ_XY pipeline_forward(PJ_LP lp, PJ *P);
155
static PJ_LP pipeline_reverse(PJ_XY xy, PJ *P);
156
157
0
static void pipeline_reassign_context(PJ *P, PJ_CONTEXT *ctx) {
158
0
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
159
0
    for (auto &step : pipeline->steps)
160
0
        proj_assign_context(step.pj, ctx);
161
0
}
162
163
970k
static void pipeline_forward_4d(PJ_COORD &point, PJ *P) {
164
970k
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
165
2.57M
    for (auto &step : pipeline->steps) {
166
2.57M
        if (!step.omit_fwd) {
167
2.57M
            if (!step.pj->inverted)
168
2.33M
                pj_fwd4d(point, step.pj);
169
239k
            else
170
239k
                pj_inv4d(point, step.pj);
171
2.57M
            if (point.xyzt.x == HUGE_VAL) {
172
5.09k
                break;
173
5.09k
            }
174
2.57M
        }
175
2.57M
    }
176
970k
}
177
178
0
static void pipeline_reverse_4d(PJ_COORD &point, PJ *P) {
179
0
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
180
0
    for (auto iterStep = pipeline->steps.rbegin();
181
0
         iterStep != pipeline->steps.rend(); ++iterStep) {
182
0
        const auto &step = *iterStep;
183
0
        if (!step.omit_inv) {
184
0
            if (step.pj->inverted)
185
0
                pj_fwd4d(point, step.pj);
186
0
            else
187
0
                pj_inv4d(point, step.pj);
188
0
            if (point.xyzt.x == HUGE_VAL) {
189
0
                break;
190
0
            }
191
0
        }
192
0
    }
193
0
}
194
195
0
static PJ_XYZ pipeline_forward_3d(PJ_LPZ lpz, PJ *P) {
196
0
    PJ_COORD point = {{0, 0, 0, 0}};
197
0
    point.lpz = lpz;
198
0
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
199
0
    for (auto &step : pipeline->steps) {
200
0
        if (!step.omit_fwd) {
201
0
            point = pj_approx_3D_trans(step.pj, PJ_FWD, point);
202
0
            if (point.xyzt.x == HUGE_VAL) {
203
0
                break;
204
0
            }
205
0
        }
206
0
    }
207
208
0
    return point.xyz;
209
0
}
210
211
0
static PJ_LPZ pipeline_reverse_3d(PJ_XYZ xyz, PJ *P) {
212
0
    PJ_COORD point = {{0, 0, 0, 0}};
213
0
    point.xyz = xyz;
214
0
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
215
0
    for (auto iterStep = pipeline->steps.rbegin();
216
0
         iterStep != pipeline->steps.rend(); ++iterStep) {
217
0
        const auto &step = *iterStep;
218
0
        if (!step.omit_inv) {
219
0
            point = proj_trans(step.pj, PJ_INV, point);
220
0
            if (point.xyzt.x == HUGE_VAL) {
221
0
                break;
222
0
            }
223
0
        }
224
0
    }
225
226
0
    return point.lpz;
227
0
}
228
229
0
static PJ_XY pipeline_forward(PJ_LP lp, PJ *P) {
230
0
    PJ_COORD point = {{0, 0, 0, 0}};
231
0
    point.lp = lp;
232
0
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
233
0
    for (auto &step : pipeline->steps) {
234
0
        if (!step.omit_fwd) {
235
0
            point = pj_approx_2D_trans(step.pj, PJ_FWD, point);
236
0
            if (point.xyzt.x == HUGE_VAL) {
237
0
                break;
238
0
            }
239
0
        }
240
0
    }
241
242
0
    return point.xy;
243
0
}
244
245
0
static PJ_LP pipeline_reverse(PJ_XY xy, PJ *P) {
246
0
    PJ_COORD point = {{0, 0, 0, 0}};
247
0
    point.xy = xy;
248
0
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
249
0
    for (auto iterStep = pipeline->steps.rbegin();
250
0
         iterStep != pipeline->steps.rend(); ++iterStep) {
251
0
        const auto &step = *iterStep;
252
0
        if (!step.omit_inv) {
253
0
            point = pj_approx_2D_trans(step.pj, PJ_INV, point);
254
0
            if (point.xyzt.x == HUGE_VAL) {
255
0
                break;
256
0
            }
257
0
        }
258
0
    }
259
260
0
    return point.lp;
261
0
}
262
263
18.9k
static PJ *destructor(PJ *P, int errlev) {
264
18.9k
    if (nullptr == P)
265
0
        return nullptr;
266
267
18.9k
    if (nullptr == P->opaque)
268
0
        return pj_default_destructor(P, errlev);
269
270
18.9k
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
271
272
18.9k
    free(pipeline->argv);
273
18.9k
    free(pipeline->current_argv);
274
275
18.9k
    delete pipeline;
276
18.9k
    P->opaque = nullptr;
277
278
18.9k
    return pj_default_destructor(P, errlev);
279
18.9k
}
280
281
/* count the number of args in pipeline definition, and mark all args as used */
282
18.9k
static size_t argc_params(paralist *params) {
283
18.9k
    size_t argc = 0;
284
934k
    for (; params != nullptr; params = params->next) {
285
915k
        argc++;
286
915k
        params->used = 1;
287
915k
    }
288
18.9k
    return ++argc; /* one extra for the sentinel */
289
18.9k
}
290
291
/* Sentinel for argument list */
292
static const char *argv_sentinel = "step";
293
294
/* turn paralist into argc/argv style argument list */
295
18.9k
static char **argv_params(paralist *params, size_t argc) {
296
18.9k
    char **argv;
297
18.9k
    size_t i = 0;
298
18.9k
    argv = static_cast<char **>(calloc(argc, sizeof(char *)));
299
18.9k
    if (nullptr == argv)
300
0
        return nullptr;
301
934k
    for (; params != nullptr; params = params->next)
302
915k
        argv[i++] = params->param;
303
18.9k
    argv[i++] = const_cast<char *>(argv_sentinel);
304
18.9k
    return argv;
305
18.9k
}
306
307
/* Being the special operator that the pipeline is, we have to handle the    */
308
/* ellipsoid differently than usual. In general, the pipeline operation does */
309
/* not need an ellipsoid, but in some cases it is beneficial nonetheless.    */
310
/* Unfortunately we can't use the normal ellipsoid setter in pj_init, since  */
311
/* it adds a +ellps parameter to the global args if nothing else is specified*/
312
/* This is problematic since that ellipsoid spec is then passed on to the    */
313
/* pipeline children. This is rarely what we want, so here we implement our  */
314
/* own logic instead. If an ellipsoid is set in the global args, it is used  */
315
/* as the pipeline ellipsoid. Otherwise we use GRS80 parameters as default.  */
316
/* At last we calculate the rest of the ellipsoid parameters and             */
317
/* re-initialize P->geod.                                                    */
318
17.7k
static void set_ellipsoid(PJ *P) {
319
17.7k
    paralist *cur, *attachment;
320
17.7k
    int err = proj_errno_reset(P);
321
322
    /* Break the linked list after the global args */
323
17.7k
    attachment = nullptr;
324
68.6k
    for (cur = P->params; cur != nullptr; cur = cur->next)
325
        /* cur->next will always be non 0 given argv_sentinel presence, */
326
        /* but this is far from being obvious for a static analyzer */
327
68.6k
        if (cur->next != nullptr &&
328
68.6k
            strcmp(argv_sentinel, cur->next->param) == 0) {
329
17.7k
            attachment = cur->next;
330
17.7k
            cur->next = nullptr;
331
17.7k
            break;
332
17.7k
        }
333
334
    /* Check if there's any ellipsoid specification in the global params. */
335
    /* If not, use GRS80 as default                                       */
336
17.7k
    if (0 != pj_ellipsoid(P)) {
337
110
        P->a = 6378137.0;
338
110
        P->f = 1.0 / 298.257222101;
339
110
        P->es = 2 * P->f - P->f * P->f;
340
341
        /* reset an "unerror": In this special use case, the errno is    */
342
        /* not an error signal, but just a reply from pj_ellipsoid,      */
343
        /* telling us that "No - there was no ellipsoid definition in    */
344
        /* the PJ you provided".                                         */
345
110
        proj_errno_reset(P);
346
110
    }
347
17.7k
    P->a_orig = P->a;
348
17.7k
    P->es_orig = P->es;
349
350
17.7k
    if (pj_calc_ellipsoid_params(P, P->a, P->es) == 0)
351
17.7k
        geod_init(P->geod, P->a, P->f);
352
353
    /* Re-attach the dangling list */
354
    /* Note: cur will always be non 0 given argv_sentinel presence, */
355
    /* but this is far from being obvious for a static analyzer */
356
17.7k
    if (cur != nullptr)
357
17.7k
        cur->next = attachment;
358
17.7k
    proj_errno_restore(P, err);
359
17.7k
}
360
361
18.9k
PJ *OPERATION(pipeline, 0) {
362
18.9k
    int i, nsteps = 0, argc;
363
18.9k
    int i_pipeline = -1, i_first_step = -1, i_current_step;
364
18.9k
    char **argv, **current_argv;
365
366
18.9k
    if (P->ctx->pipelineInitRecursiongCounter == 5) {
367
        // Can happen for a string like:
368
        // proj=pipeline step "x="""," u=" proj=pipeline step ste=""[" u="
369
        // proj=pipeline step ste="[" u=" proj=pipeline step ste="[" u="
370
        // proj=pipeline step ste="[" u=" proj=pipeline step ste="[" u="
371
        // proj=pipeline step ste="[" u=" proj=pipeline step ste="[" u="
372
        // proj=pipeline step ste="[" u=" proj=pipeline p step ste="[" u="
373
        // proj=pipeline step ste="[" u=" proj=pipeline step ste="[" u="
374
        // proj=pipeline step ste="[" u=" proj=pipeline step ""x="""""""""""
375
        // Probably an issue with the quoting handling code
376
        // But doesn't hurt to add an extra safety check
377
0
        proj_log_error(P, _("Pipeline: too deep recursion"));
378
0
        return destructor(
379
0
            P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX); /* ERROR: nested pipelines */
380
0
    }
381
382
18.9k
    P->fwd4d = pipeline_forward_4d;
383
18.9k
    P->inv4d = pipeline_reverse_4d;
384
18.9k
    P->fwd3d = pipeline_forward_3d;
385
18.9k
    P->inv3d = pipeline_reverse_3d;
386
18.9k
    P->fwd = pipeline_forward;
387
18.9k
    P->inv = pipeline_reverse;
388
18.9k
    P->destructor = destructor;
389
18.9k
    P->reassign_context = pipeline_reassign_context;
390
391
    /* Currently, the pipeline driver is a raw bit mover, enabling other
392
     * operations */
393
    /* to collaborate efficiently. All prep/fin stuff is done at the step
394
     * levels.
395
     */
396
18.9k
    P->skip_fwd_prepare = 1;
397
18.9k
    P->skip_fwd_finalize = 1;
398
18.9k
    P->skip_inv_prepare = 1;
399
18.9k
    P->skip_inv_finalize = 1;
400
401
18.9k
    P->opaque = new (std::nothrow) Pipeline();
402
18.9k
    if (nullptr == P->opaque)
403
0
        return destructor(P, PROJ_ERR_INVALID_OP /* ENOMEM */);
404
405
18.9k
    argc = (int)argc_params(P->params);
406
18.9k
    auto pipeline = static_cast<struct Pipeline *>(P->opaque);
407
18.9k
    pipeline->argv = argv = argv_params(P->params, argc);
408
18.9k
    if (nullptr == argv)
409
0
        return destructor(P, PROJ_ERR_INVALID_OP /* ENOMEM */);
410
411
18.9k
    pipeline->current_argv = current_argv =
412
18.9k
        static_cast<char **>(calloc(argc, sizeof(char *)));
413
18.9k
    if (nullptr == current_argv)
414
0
        return destructor(P, PROJ_ERR_OTHER /*ENOMEM*/);
415
416
    /* Do some syntactical sanity checking */
417
877k
    for (i = 0; i < argc && argv[i] != nullptr; i++) {
418
859k
        if (0 == strcmp(argv_sentinel, argv[i])) {
419
166k
            if (-1 == i_pipeline) {
420
0
                proj_log_error(P, _("Pipeline: +step before +proj=pipeline"));
421
0
                return destructor(P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX);
422
0
            }
423
166k
            if (0 == nsteps)
424
17.8k
                i_first_step = i;
425
166k
            nsteps++;
426
166k
            continue;
427
166k
        }
428
429
692k
        if (0 == strcmp("proj=pipeline", argv[i])) {
430
18.9k
            if (-1 != i_pipeline) {
431
0
                proj_log_error(P, _("Pipeline: Nesting only allowed when child "
432
0
                                    "pipelines are wrapped in '+init's"));
433
0
                return destructor(
434
0
                    P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX); /* ERROR: nested
435
                                                             pipelines */
436
0
            }
437
18.9k
            i_pipeline = i;
438
673k
        } else if (0 == nsteps && 0 == strncmp(argv[i], "proj=", 5)) {
439
            // Non-sensical to have proj= in the general pipeline parameters.
440
            // Would not be a big issue in itself, but this makes bad
441
            // performance in parsing hostile pipelines more likely, such as the
442
            // one of
443
            // https://bugs.chromium.org/p/oss-fuzz/issues/detail?id=41290
444
1.13k
            proj_log_error(
445
1.13k
                P, _("Pipeline: proj= operator before first step not allowed"));
446
1.13k
            return destructor(P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX);
447
672k
        } else if (0 == nsteps && 0 == strncmp(argv[i], "o_proj=", 7)) {
448
            // Same as above.
449
11
            proj_log_error(
450
11
                P,
451
11
                _("Pipeline: o_proj= operator before first step not allowed"));
452
11
            return destructor(P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX);
453
11
        }
454
692k
    }
455
17.8k
    nsteps--; /* Last instance of +step is just a sentinel */
456
457
17.8k
    if (-1 == i_pipeline)
458
0
        return destructor(
459
0
            P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX); /* ERROR: no pipeline def */
460
461
17.8k
    if (0 == nsteps)
462
54
        return destructor(
463
54
            P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX); /* ERROR: no pipeline def */
464
465
17.7k
    set_ellipsoid(P);
466
467
    /* Now loop over all steps, building a new set of arguments for each init */
468
17.7k
    i_current_step = i_first_step;
469
96.2k
    for (i = 0; i < nsteps; i++) {
470
88.3k
        int j;
471
88.3k
        int current_argc = 0;
472
88.3k
        int err;
473
88.3k
        PJ *next_step = nullptr;
474
475
        /* Build a set of setup args for the current step */
476
88.3k
        proj_log_trace(P, "Pipeline: Building arg list for step no. %d", i);
477
478
        /* First add the step specific args */
479
454k
        for (j = i_current_step + 1; 0 != strcmp("step", argv[j]); j++)
480
366k
            current_argv[current_argc++] = argv[j];
481
482
88.3k
        i_current_step = j;
483
484
        /* Then add the global args */
485
582k
        for (j = i_pipeline + 1; 0 != strcmp("step", argv[j]); j++)
486
494k
            current_argv[current_argc++] = argv[j];
487
488
88.3k
        proj_log_trace(P, "Pipeline: init - %s, %d", current_argv[0],
489
88.3k
                       current_argc);
490
860k
        for (j = 1; j < current_argc; j++)
491
772k
            proj_log_trace(P, "    %s", current_argv[j]);
492
493
88.3k
        err = proj_errno_reset(P);
494
495
88.3k
        P->ctx->pipelineInitRecursiongCounter++;
496
88.3k
        next_step = pj_create_argv_internal(P->ctx, current_argc, current_argv);
497
88.3k
        P->ctx->pipelineInitRecursiongCounter--;
498
88.3k
        proj_log_trace(P, "Pipeline: Step %d (%s) at %p", i, current_argv[0],
499
88.3k
                       next_step);
500
501
88.3k
        if (nullptr == next_step) {
502
            /* The step init failed, but possibly without setting errno. If so,
503
             * we say "malformed" */
504
9.83k
            int err_to_report = proj_errno(P);
505
9.83k
            if (0 == err_to_report)
506
0
                err_to_report = PROJ_ERR_INVALID_OP_WRONG_SYNTAX;
507
9.83k
            proj_log_error(P, _("Pipeline: Bad step definition: %s (%s)"),
508
9.83k
                           current_argv[0],
509
9.83k
                           proj_context_errno_string(P->ctx, err_to_report));
510
9.83k
            return destructor(P, err_to_report); /* ERROR: bad pipeline def */
511
9.83k
        }
512
78.5k
        next_step->parent = P;
513
514
78.5k
        proj_errno_restore(P, err);
515
516
        /* Is this step inverted? */
517
872k
        for (j = 0; j < current_argc; j++) {
518
794k
            if (0 == strcmp("inv", current_argv[j])) {
519
                /* if +inv exists in both global and local args the forward
520
                 * operation should be used */
521
25.6k
                next_step->inverted = next_step->inverted == 0 ? 1 : 0;
522
25.6k
            }
523
794k
        }
524
525
78.5k
        bool omit_fwd = pj_param(P->ctx, next_step->params, "bomit_fwd").i != 0;
526
78.5k
        bool omit_inv = pj_param(P->ctx, next_step->params, "bomit_inv").i != 0;
527
78.5k
        pipeline->steps.emplace_back(next_step, omit_fwd, omit_inv);
528
529
78.5k
        proj_log_trace(P, "Pipeline at [%p]:    step at [%p] (%s) done", P,
530
78.5k
                       next_step, current_argv[0]);
531
78.5k
    }
532
533
    /* Require a forward path through the pipeline */
534
41.7k
    for (auto &step : pipeline->steps) {
535
41.7k
        PJ *Q = step.pj;
536
41.7k
        if (step.omit_fwd) {
537
69
            continue;
538
69
        }
539
41.6k
        if (Q->inverted) {
540
12.5k
            if (Q->inv || Q->inv3d || Q->inv4d) {
541
12.4k
                continue;
542
12.4k
            }
543
127
            proj_log_error(
544
127
                P, _("Pipeline: Inverse operation for %s is not available"),
545
127
                Q->short_name);
546
127
            return destructor(P, PROJ_ERR_OTHER_NO_INVERSE_OP);
547
29.1k
        } else {
548
29.1k
            if (Q->fwd || Q->fwd3d || Q->fwd4d) {
549
29.1k
                continue;
550
29.1k
            }
551
0
            proj_log_error(
552
0
                P, _("Pipeline: Forward operation for %s is not available"),
553
0
                Q->short_name);
554
0
            return destructor(P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX);
555
29.1k
        }
556
41.6k
    }
557
558
    /* determine if an inverse operation is possible */
559
41.4k
    for (auto &step : pipeline->steps) {
560
41.4k
        PJ *Q = step.pj;
561
41.4k
        if (step.omit_inv || pj_has_inverse(Q)) {
562
41.0k
            continue;
563
41.0k
        } else {
564
348
            P->inv = nullptr;
565
348
            P->inv3d = nullptr;
566
348
            P->inv4d = nullptr;
567
348
            break;
568
348
        }
569
41.4k
    }
570
571
    /* Replace PJ_IO_UNITS_WHATEVER with input/output units of neighbouring
572
     * steps where */
573
    /* it make sense. It does in most cases but not always, for instance */
574
    /*      proj=pipeline step proj=unitconvert xy_in=deg xy_out=rad step ... */
575
    /* where the left-hand side units of the first step shouldn't be changed to
576
     * RADIANS */
577
    /* as it will result in deg->rad conversions in cs2cs and other
578
     * applications.
579
     */
580
581
41.5k
    for (i = nsteps - 2; i >= 0; --i) {
582
33.7k
        auto pj = pipeline->steps[i].pj;
583
33.7k
        if (pj_left(pj) == PJ_IO_UNITS_WHATEVER &&
584
10.7k
            pj_right(pj) == PJ_IO_UNITS_WHATEVER) {
585
9.30k
            const auto right_pj = pipeline->steps[i + 1].pj;
586
9.30k
            const auto right_pj_left = pj_left(right_pj);
587
9.30k
            const auto right_pj_right = pj_right(right_pj);
588
9.30k
            if (right_pj_left != right_pj_right ||
589
8.15k
                right_pj_left != PJ_IO_UNITS_WHATEVER) {
590
8.15k
                pj->left = right_pj_left;
591
8.15k
                pj->right = right_pj_left;
592
8.15k
            }
593
9.30k
        }
594
33.7k
    }
595
596
41.5k
    for (i = 1; i < nsteps; i++) {
597
33.7k
        auto pj = pipeline->steps[i].pj;
598
33.7k
        if (pj_left(pj) == PJ_IO_UNITS_WHATEVER &&
599
3.98k
            pj_right(pj) == PJ_IO_UNITS_WHATEVER) {
600
2.50k
            const auto left_pj = pipeline->steps[i - 1].pj;
601
2.50k
            const auto left_pj_left = pj_left(left_pj);
602
2.50k
            const auto left_pj_right = pj_right(left_pj);
603
2.50k
            if (left_pj_left != left_pj_right ||
604
1.64k
                left_pj_right != PJ_IO_UNITS_WHATEVER) {
605
1.64k
                pj->left = left_pj_right;
606
1.64k
                pj->right = left_pj_right;
607
1.64k
            }
608
2.50k
        }
609
33.7k
    }
610
611
    /* Check that units between each steps match each other, fail if they don't
612
     */
613
39.2k
    for (i = 0; i + 1 < nsteps; i++) {
614
31.7k
        enum pj_io_units curr_step_output = pj_right(pipeline->steps[i].pj);
615
31.7k
        enum pj_io_units next_step_input = pj_left(pipeline->steps[i + 1].pj);
616
617
31.7k
        if (curr_step_output == PJ_IO_UNITS_WHATEVER ||
618
28.7k
            next_step_input == PJ_IO_UNITS_WHATEVER)
619
3.95k
            continue;
620
621
27.8k
        if (curr_step_output != next_step_input) {
622
311
            proj_log_error(
623
311
                P, _("Pipeline: Mismatched units between step %d and %d"),
624
311
                i + 1, i + 2);
625
311
            return destructor(P, PROJ_ERR_INVALID_OP_WRONG_SYNTAX);
626
311
        }
627
27.8k
    }
628
629
7.49k
    proj_log_trace(
630
7.49k
        P, "Pipeline: %d steps built. Determining i/o characteristics", nsteps);
631
632
    /* Determine forward input (= reverse output) data type */
633
7.49k
    P->left = pj_left(pipeline->steps.front().pj);
634
635
    /* Now, correspondingly determine forward output (= reverse input) data type
636
     */
637
7.49k
    P->right = pj_right(pipeline->steps.back().pj);
638
7.49k
    return P;
639
7.80k
}
640
641
25.8k
static void push(PJ_COORD &point, PJ *P) {
642
25.8k
    if (P->parent == nullptr)
643
0
        return;
644
645
25.8k
    struct Pipeline *pipeline =
646
25.8k
        static_cast<struct Pipeline *>(P->parent->opaque);
647
25.8k
    struct PushPop *pushpop = static_cast<struct PushPop *>(P->opaque);
648
649
25.8k
    if (pushpop->v1)
650
8.40k
        pipeline->stack[0].push(point.v[0]);
651
25.8k
    if (pushpop->v2)
652
4.70k
        pipeline->stack[1].push(point.v[1]);
653
25.8k
    if (pushpop->v3)
654
11.4k
        pipeline->stack[2].push(point.v[2]);
655
25.8k
    if (pushpop->v4)
656
0
        pipeline->stack[3].push(point.v[3]);
657
25.8k
}
658
659
102k
static void pop(PJ_COORD &point, PJ *P) {
660
102k
    if (P->parent == nullptr)
661
0
        return;
662
663
102k
    struct Pipeline *pipeline =
664
102k
        static_cast<struct Pipeline *>(P->parent->opaque);
665
102k
    struct PushPop *pushpop = static_cast<struct PushPop *>(P->opaque);
666
667
102k
    if (pushpop->v1 && !pipeline->stack[0].empty()) {
668
6.37k
        point.v[0] = pipeline->stack[0].top();
669
6.37k
        pipeline->stack[0].pop();
670
6.37k
    }
671
672
102k
    if (pushpop->v2 && !pipeline->stack[1].empty()) {
673
2.34k
        point.v[1] = pipeline->stack[1].top();
674
2.34k
        pipeline->stack[1].pop();
675
2.34k
    }
676
677
102k
    if (pushpop->v3 && !pipeline->stack[2].empty()) {
678
9.73k
        point.v[2] = pipeline->stack[2].top();
679
9.73k
        pipeline->stack[2].pop();
680
9.73k
    }
681
682
102k
    if (pushpop->v4 && !pipeline->stack[3].empty()) {
683
0
        point.v[3] = pipeline->stack[3].top();
684
0
        pipeline->stack[3].pop();
685
0
    }
686
102k
}
687
688
17.6k
static PJ *setup_pushpop(PJ *P) {
689
17.6k
    auto pushpop =
690
17.6k
        static_cast<struct PushPop *>(calloc(1, sizeof(struct PushPop)));
691
17.6k
    P->opaque = pushpop;
692
17.6k
    if (nullptr == P->opaque)
693
0
        return destructor(P, PROJ_ERR_OTHER /*ENOMEM*/);
694
695
17.6k
    if (pj_param_exists(P->params, "v_1"))
696
1.09k
        pushpop->v1 = true;
697
698
17.6k
    if (pj_param_exists(P->params, "v_2"))
699
571
        pushpop->v2 = true;
700
701
17.6k
    if (pj_param_exists(P->params, "v_3"))
702
14.0k
        pushpop->v3 = true;
703
704
17.6k
    if (pj_param_exists(P->params, "v_4"))
705
14
        pushpop->v4 = true;
706
707
17.6k
    P->left = PJ_IO_UNITS_WHATEVER;
708
17.6k
    P->right = PJ_IO_UNITS_WHATEVER;
709
710
17.6k
    return P;
711
17.6k
}
712
713
8.28k
PJ *OPERATION(push, 0) {
714
8.28k
    P->fwd4d = push;
715
8.28k
    P->inv4d = pop;
716
717
8.28k
    return setup_pushpop(P);
718
8.28k
}
719
720
9.40k
PJ *OPERATION(pop, 0) {
721
9.40k
    P->inv4d = push;
722
9.40k
    P->fwd4d = pop;
723
724
9.40k
    return setup_pushpop(P);
725
9.40k
}