aboutsummaryrefslogtreecommitdiff
path: root/src
diff options
context:
space:
mode:
Diffstat (limited to 'src')
-rw-r--r--src/iso19111/coordinateoperation.cpp44
-rw-r--r--src/iso19111/io.cpp440
-rw-r--r--src/pipeline.cpp111
-rw-r--r--src/pj_list.h2
-rw-r--r--src/proj_internal.h1
5 files changed, 254 insertions, 344 deletions
diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp
index d36d901f..d964cdc1 100644
--- a/src/iso19111/coordinateoperation.cpp
+++ b/src/iso19111/coordinateoperation.cpp
@@ -5489,7 +5489,7 @@ void Conversion::_exportToPROJString(
if (isUTM(zone, north)) {
bConversionDone = true;
formatter->addStep("utm");
- if( useApprox ) {
+ if (useApprox) {
formatter->addParam("approx");
}
formatter->addParam("zone", zone);
@@ -8138,6 +8138,18 @@ void Transformation::_exportToPROJString(
double z =
parameterValueNumericAsSI(EPSG_CODE_PARAMETER_Z_AXIS_TRANSLATION);
+ if (methodEPSGCode == EPSG_CODE_METHOD_COORDINATE_FRAME_GEOGRAPHIC_2D ||
+ methodEPSGCode ==
+ EPSG_CODE_METHOD_TIME_DEPENDENT_COORDINATE_FRAME_GEOGRAPHIC_2D ||
+ methodEPSGCode == EPSG_CODE_METHOD_POSITION_VECTOR_GEOGRAPHIC_2D ||
+ methodEPSGCode ==
+ EPSG_CODE_METHOD_TIME_DEPENDENT_POSITION_VECTOR_GEOGRAPHIC_2D ||
+ methodEPSGCode ==
+ EPSG_CODE_METHOD_GEOCENTRIC_TRANSLATION_GEOGRAPHIC_2D) {
+ formatter->addStep("push");
+ formatter->addParam("v_3");
+ }
+
setupPROJGeodeticSourceCRS(formatter, sourceCRS(), "Helmert");
formatter->addStep("helmert");
@@ -8204,6 +8216,18 @@ void Transformation::_exportToPROJString(
setupPROJGeodeticTargetCRS(formatter, targetCRS(), "Helmert");
+ if (methodEPSGCode == EPSG_CODE_METHOD_COORDINATE_FRAME_GEOGRAPHIC_2D ||
+ methodEPSGCode ==
+ EPSG_CODE_METHOD_TIME_DEPENDENT_COORDINATE_FRAME_GEOGRAPHIC_2D ||
+ methodEPSGCode == EPSG_CODE_METHOD_POSITION_VECTOR_GEOGRAPHIC_2D ||
+ methodEPSGCode ==
+ EPSG_CODE_METHOD_TIME_DEPENDENT_POSITION_VECTOR_GEOGRAPHIC_2D ||
+ methodEPSGCode ==
+ EPSG_CODE_METHOD_GEOCENTRIC_TRANSLATION_GEOGRAPHIC_2D) {
+ formatter->addStep("pop");
+ formatter->addParam("v_3");
+ }
+
return;
}
@@ -8250,6 +8274,14 @@ void Transformation::_exportToPROJString(
double pz = parameterValueNumericAsSI(
EPSG_CODE_PARAMETER_ORDINATE_3_EVAL_POINT);
+ if (methodEPSGCode ==
+ EPSG_CODE_METHOD_MOLODENSKY_BADEKAS_PV_GEOGRAPHIC_2D ||
+ methodEPSGCode ==
+ EPSG_CODE_METHOD_MOLODENSKY_BADEKAS_CF_GEOGRAPHIC_2D) {
+ formatter->addStep("push");
+ formatter->addParam("v_3");
+ }
+
setupPROJGeodeticSourceCRS(formatter, sourceCRS(),
"Molodensky-Badekas");
@@ -8273,6 +8305,14 @@ void Transformation::_exportToPROJString(
setupPROJGeodeticTargetCRS(formatter, targetCRS(),
"Molodensky-Badekas");
+ if (methodEPSGCode ==
+ EPSG_CODE_METHOD_MOLODENSKY_BADEKAS_PV_GEOGRAPHIC_2D ||
+ methodEPSGCode ==
+ EPSG_CODE_METHOD_MOLODENSKY_BADEKAS_CF_GEOGRAPHIC_2D) {
+ formatter->addStep("pop");
+ formatter->addParam("v_3");
+ }
+
return;
}
@@ -9978,7 +10018,7 @@ struct FilterResults {
!gridDesc.available) {
gridsAvailable = false;
}
- if (gridDesc.packageName.empty()) {
+ if (gridDesc.packageName.empty() && !gridDesc.available) {
gridsKnown = false;
}
}
diff --git a/src/iso19111/io.cpp b/src/iso19111/io.cpp
index 60c28201..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..).
@@ -5844,14 +5939,6 @@ struct PROJStringParser::Private {
CRSNNPtr buildBoundOrCompoundCRSIfNeeded(int iStep, CRSNNPtr crs);
UnitOfMeasure buildUnit(Step &step, const std::string &unitsParamName,
const std::string &toMeterParamName);
- CoordinateOperationNNPtr buildHelmertTransformation(
- int iStep, int iFirstAxisSwap = -1, int iFirstUnitConvert = -1,
- int iFirstGeogStep = -1, int iSecondGeogStep = -1,
- int iSecondAxisSwap = -1, int iSecondUnitConvert = -1);
- CoordinateOperationNNPtr buildMolodenskyTransformation(
- int iStep, int iFirstAxisSwap = -1, int iFirstUnitConvert = -1,
- int iFirstGeogStep = -1, int iSecondGeogStep = -1,
- int iSecondAxisSwap = -1, int iSecondUnitConvert = -1);
enum class AxisType { REGULAR, NORTH_POLE, SOUTH_POLE };
@@ -7170,261 +7257,6 @@ CRSNNPtr PROJStringParser::Private::buildProjectedCRS(
return crs;
}
-// ---------------------------------------------------------------------------
-
-static bool isDatumDefiningParam(const std::string &param) {
- return (param == "datum" || param == "ellps" || param == "a" ||
- param == "b" || param == "rf" || param == "f" || param == "R");
-}
-
-// ---------------------------------------------------------------------------
-
-CoordinateOperationNNPtr PROJStringParser::Private::buildHelmertTransformation(
- int iStep, int iFirstAxisSwap, int iFirstUnitConvert, int iFirstGeogStep,
- int iSecondGeogStep, int iSecondAxisSwap, int iSecondUnitConvert) {
- auto &step = steps_[iStep];
- auto datum = buildDatum(step, std::string());
- auto cs = CartesianCS::createGeocentric(UnitOfMeasure::METRE);
-
- auto mapWithUnknownName = createMapWithUnknownName();
-
- auto sourceCRS =
- iFirstGeogStep >= 0
- ? util::nn_static_pointer_cast<crs::CRS>(
- buildGeographicCRS(iFirstGeogStep, iFirstUnitConvert,
- iFirstAxisSwap, true, false))
- : util::nn_static_pointer_cast<crs::CRS>(
- GeodeticCRS::create(mapWithUnknownName, datum, cs));
- auto targetCRS =
- iSecondGeogStep >= 0
- ? util::nn_static_pointer_cast<crs::CRS>(
- buildGeographicCRS(iSecondGeogStep, iSecondUnitConvert,
- iSecondAxisSwap, true, false))
- : util::nn_static_pointer_cast<crs::CRS>(
- GeodeticCRS::create(mapWithUnknownName, datum, cs));
-
- double x = 0;
- double y = 0;
- double z = 0;
- double rx = 0;
- double ry = 0;
- double rz = 0;
- double s = 0;
- double dx = 0;
- double dy = 0;
- double dz = 0;
- double drx = 0;
- double dry = 0;
- double drz = 0;
- double ds = 0;
- double t_epoch = 0;
- bool rotationTerms = false;
- bool timeDependent = false;
- bool conventionFound = false;
- bool positionVectorConvention = false;
-
- struct Params {
- double *pValue;
- const char *name;
- bool *pPresent;
- };
- const Params knownParams[] = {
- {&x, "x", nullptr},
- {&y, "y", nullptr},
- {&z, "z", nullptr},
- {&rx, "rx", &rotationTerms},
- {&ry, "ry", &rotationTerms},
- {&rz, "rz", &rotationTerms},
- {&s, "s", &rotationTerms},
- {&dx, "dx", &timeDependent},
- {&dy, "dy", &timeDependent},
- {&dz, "dz", &timeDependent},
- {&drx, "drx", &timeDependent},
- {&dry, "dry", &timeDependent},
- {&drz, "drz", &timeDependent},
- {&ds, "ds", &timeDependent},
- {&t_epoch, "t_epoch", &timeDependent},
- {nullptr, "exact", nullptr},
- };
-
- for (const auto &param : step.paramValues) {
- if (isDatumDefiningParam(param.key)) {
- continue;
- }
- if (param.key == "convention") {
- if (param.value == "position_vector") {
- positionVectorConvention = true;
- conventionFound = true;
- } else if (param.value == "coordinate_frame") {
- positionVectorConvention = false;
- conventionFound = true;
- } else {
- throw ParsingException("unsupported convention");
- }
- } else {
- bool found = false;
- for (auto &&knownParam : knownParams) {
- if (param.key == knownParam.name) {
- found = true;
- if (knownParam.pValue)
- *(knownParam.pValue) = getNumericValue(param.value);
- if (knownParam.pPresent)
- *(knownParam.pPresent) = true;
- break;
- }
- }
- if (!found) {
- throw ParsingException("unsupported keyword for Helmert: " +
- param.key);
- }
- }
- }
-
- rotationTerms |= timeDependent;
- if (rotationTerms && !conventionFound) {
- throw ParsingException("missing convention");
- }
-
- std::vector<metadata::PositionalAccuracyNNPtr> emptyAccuracies;
-
- auto transf = ([&]() {
- if (!rotationTerms) {
- return Transformation::createGeocentricTranslations(
- mapWithUnknownName, sourceCRS, targetCRS, x, y, z,
- emptyAccuracies);
- } else if (positionVectorConvention) {
- if (timeDependent) {
- return Transformation::createTimeDependentPositionVector(
- mapWithUnknownName, sourceCRS, targetCRS, x, y, z, rx, ry,
- rz, s, dx, dy, dz, drx, dry, drz, ds, t_epoch,
- emptyAccuracies);
- } else {
- return Transformation::createPositionVector(
- mapWithUnknownName, sourceCRS, targetCRS, x, y, z, rx, ry,
- rz, s, emptyAccuracies);
- }
- } else {
- if (timeDependent) {
- return Transformation::
- createTimeDependentCoordinateFrameRotation(
- mapWithUnknownName, sourceCRS, targetCRS, x, y, z, rx,
- ry, rz, s, dx, dy, dz, drx, dry, drz, ds, t_epoch,
- emptyAccuracies);
- } else {
- return Transformation::createCoordinateFrameRotation(
- mapWithUnknownName, sourceCRS, targetCRS, x, y, z, rx, ry,
- rz, s, emptyAccuracies);
- }
- }
- })();
-
- if (step.inverted) {
- return util::nn_static_pointer_cast<CoordinateOperation>(
- transf->inverse());
- } else {
- return util::nn_static_pointer_cast<CoordinateOperation>(transf);
- }
-}
-
-// ---------------------------------------------------------------------------
-
-CoordinateOperationNNPtr
-PROJStringParser::Private::buildMolodenskyTransformation(
- int iStep, int iFirstAxisSwap, int iFirstUnitConvert, int iFirstGeogStep,
- int iSecondGeogStep, int iSecondAxisSwap, int iSecondUnitConvert) {
- auto &step = steps_[iStep];
-
- double dx = 0;
- double dy = 0;
- double dz = 0;
- double da = 0;
- double df = 0;
-
- struct Params {
- double *pValue;
- const char *name;
- };
- const Params knownParams[] = {
- {&dx, "dx"}, {&dy, "dy"}, {&dz, "dz"}, {&da, "da"}, {&df, "df"},
- };
- bool abridged = false;
-
- for (const auto &param : step.paramValues) {
- if (isDatumDefiningParam(param.key)) {
- continue;
- } else if (param.key == "abridged") {
- abridged = true;
- } else {
- bool found = false;
- for (auto &&knownParam : knownParams) {
- if (param.key == knownParam.name) {
- found = true;
- if (knownParam.pValue)
- *(knownParam.pValue) = getNumericValue(param.value);
- break;
- }
- }
- if (!found) {
- throw ParsingException("unsupported keyword for Molodensky: " +
- param.key);
- }
- }
- }
-
- auto datum = buildDatum(step, std::string());
- auto sourceCRS = iFirstGeogStep >= 0
- ? buildGeographicCRS(iFirstGeogStep, iFirstUnitConvert,
- iFirstAxisSwap, true, false)
- : buildGeographicCRS(iStep, -1, -1, true, false);
-
- const auto &ellps = sourceCRS->ellipsoid();
- const double a = ellps->semiMajorAxis().getSIValue();
- const double rf = ellps->computedInverseFlattening();
- const double target_a = a + da;
- const double target_rf = 1.0 / (1.0 / rf + df);
-
- auto mapWithUnknownName = createMapWithUnknownName();
-
- auto target_ellipsoid =
- Ellipsoid::createFlattenedSphere(mapWithUnknownName, Length(target_a),
- Scale(target_rf))
- ->identify();
- auto target_datum = GeodeticReferenceFrame::create(
- mapWithUnknownName, target_ellipsoid, util::optional<std::string>(),
- PrimeMeridian::GREENWICH);
-
- auto targetCRS = util::nn_static_pointer_cast<crs::CRS>(
- iSecondGeogStep >= 0
- ? buildGeographicCRS(iSecondGeogStep, iSecondUnitConvert,
- iSecondAxisSwap, true, false)
- : GeographicCRS::create(mapWithUnknownName, target_datum,
- EllipsoidalCS::createLongitudeLatitude(
- UnitOfMeasure::DEGREE)));
-
- auto sourceCRS_as_CRS = util::nn_static_pointer_cast<crs::CRS>(sourceCRS);
-
- std::vector<metadata::PositionalAccuracyNNPtr> emptyAccuracies;
-
- auto transf = ([&]() {
- if (abridged) {
- return Transformation::createAbridgedMolodensky(
- mapWithUnknownName, sourceCRS_as_CRS, targetCRS, dx, dy, dz, da,
- df, emptyAccuracies);
- } else {
- return Transformation::createMolodensky(
- mapWithUnknownName, sourceCRS_as_CRS, targetCRS, dx, dy, dz, da,
- df, emptyAccuracies);
- }
- })();
-
- if (step.inverted) {
- return util::nn_static_pointer_cast<CoordinateOperation>(
- transf->inverse());
- } else {
- return util::nn_static_pointer_cast<CoordinateOperation>(transf);
- }
-}
-
//! @endcond
// ---------------------------------------------------------------------------
@@ -7655,10 +7487,6 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
int iSecondUnitConvert = -1;
int iFirstAxisSwap = -1;
int iSecondAxisSwap = -1;
- int iHelmert = -1;
- int iFirstCart = -1;
- int iSecondCart = -1;
- int iMolodensky = -1;
bool unexpectedStructure = d->steps_.empty();
for (int i = 0; i < static_cast<int>(d->steps_.size()); i++) {
const auto &stepName = d->steps_[i].name;
@@ -7689,27 +7517,6 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
unexpectedStructure = true;
break;
}
- } else if (stepName == "helmert") {
- if (iHelmert >= 0) {
- unexpectedStructure = true;
- break;
- }
- iHelmert = i;
- } else if (stepName == "cart") {
- if (iFirstCart < 0) {
- iFirstCart = i;
- } else if (iSecondCart < 0) {
- iSecondCart = i;
- } else {
- unexpectedStructure = true;
- break;
- }
- } else if (stepName == "molodensky") {
- if (iMolodensky >= 0) {
- unexpectedStructure = true;
- break;
- }
- iMolodensky = i;
} else if (isProjectedStep(stepName)) {
if (iProjStep >= 0) {
unexpectedStructure = true;
@@ -7722,19 +7529,13 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
}
}
- if (iHelmert < 0 && iMolodensky < 0 && !d->steps_.empty()) {
+ if (!d->steps_.empty()) {
// CRS candidate
if ((d->steps_.size() == 1 &&
d->getParamValue(d->steps_[0], "type") != "crs") ||
(d->steps_.size() > 1 && d->getGlobalParamValue("type") != "crs")) {
unexpectedStructure = true;
}
- } else if (iHelmert >= 0 &&
- (d->hasParamValue(d->steps_[iHelmert], "theta") ||
- d->hasParamValue(d->steps_[iHelmert], "exact") ||
- d->hasParamValue(d->steps_[iHelmert], "transpose") ||
- d->hasParamValue(d->steps_[iHelmert], "towgs84"))) {
- unexpectedStructure = true;
}
struct Logger {
@@ -7771,6 +7572,7 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
if (d->ctx_) {
d->ctx_->curStringInCreateFromPROJString = projString;
}
+ auto log_level = proj_log_level(pj_context, PJ_LOG_NONE);
auto pj = pj_create_internal(
pj_context, (projString.find("type=crs") != std::string::npos
? projString + " +disable_grid_presence_check"
@@ -7780,6 +7582,7 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
d->ctx_->curStringInCreateFromPROJString.clear();
}
valid = pj != nullptr;
+ proj_log_level(pj_context, log_level);
// Remove parameters not understood by PROJ.
if (valid && d->steps_.size() == 1) {
@@ -7875,8 +7678,7 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
if (!unexpectedStructure) {
if (iFirstGeogStep == 0 && !d->steps_[iFirstGeogStep].inverted &&
- iSecondGeogStep < 0 && iProjStep < 0 && iHelmert < 0 &&
- iFirstCart < 0 && iMolodensky < 0 &&
+ iSecondGeogStep < 0 && iProjStep < 0 &&
(iFirstUnitConvert < 0 || iSecondUnitConvert < 0) &&
(iFirstAxisSwap < 0 || iSecondAxisSwap < 0)) {
const bool ignoreVUnits = false;
@@ -7893,8 +7695,7 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
}
if (iProjStep >= 0 && !d->steps_[iProjStep].inverted &&
(iFirstGeogStep < 0 || iFirstGeogStep + 1 == iProjStep) &&
- iMolodensky < 0 && iSecondGeogStep < 0 && iFirstCart < 0 &&
- iHelmert < 0) {
+ iSecondGeogStep < 0) {
if (iFirstGeogStep < 0)
iFirstGeogStep = iProjStep;
const bool ignoreVUnits = true;
@@ -7920,47 +7721,6 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
}
}
}
-
- if (iProjStep < 0 && iHelmert > 0 && iMolodensky < 0 &&
- (iFirstGeogStep < 0 || iFirstGeogStep == iFirstCart - 1 ||
- (iFirstGeogStep == iSecondCart + 1 && iSecondGeogStep < 0)) &&
- iFirstCart == iHelmert - 1 && iSecondCart == iHelmert + 1 &&
- (iSecondGeogStep < 0 || iSecondGeogStep == iSecondCart + 1) &&
- !d->steps_[iFirstCart].inverted &&
- d->steps_[iSecondCart].inverted && iFirstAxisSwap < iHelmert &&
- iFirstUnitConvert < iHelmert &&
- (iSecondAxisSwap < 0 || iSecondAxisSwap > iHelmert) &&
- (iSecondUnitConvert < 0 || iSecondUnitConvert > iHelmert)) {
- return d->buildHelmertTransformation(
- iHelmert, iFirstAxisSwap, iFirstUnitConvert,
- iFirstGeogStep >= 0 && iFirstGeogStep == iFirstCart - 1
- ? iFirstGeogStep
- : iFirstCart,
- iFirstGeogStep == iSecondCart + 1
- ? iFirstGeogStep
- : iSecondGeogStep == iSecondCart + 1 ? iSecondGeogStep
- : iSecondCart,
- iSecondAxisSwap, iSecondUnitConvert);
- }
-
- if (iProjStep < 0 && iHelmert < 0 && iMolodensky > 0 &&
- (iFirstGeogStep < 0 || iFirstGeogStep == iMolodensky - 1 ||
- (iFirstGeogStep == iMolodensky + 1 && iSecondGeogStep < 0)) &&
- (iSecondGeogStep < 0 || iSecondGeogStep == iMolodensky + 1) &&
- iFirstAxisSwap < iMolodensky && iFirstUnitConvert < iMolodensky &&
- (iSecondAxisSwap < 0 || iSecondAxisSwap > iMolodensky) &&
- (iSecondUnitConvert < 0 || iSecondUnitConvert > iMolodensky)) {
- return d->buildMolodenskyTransformation(
- iMolodensky, iFirstAxisSwap, iFirstUnitConvert,
- iFirstGeogStep >= 0 && iFirstGeogStep == iMolodensky - 1
- ? iFirstGeogStep
- : iMolodensky,
- iFirstGeogStep == iMolodensky + 1
- ? iFirstGeogStep
- : iSecondGeogStep == iMolodensky + 1 ? iSecondGeogStep
- : iMolodensky,
- iSecondAxisSwap, iSecondUnitConvert);
- }
}
auto props = PropertyMap();
diff --git a/src/pipeline.cpp b/src/pipeline.cpp
index c4a9d2c0..39563c65 100644
--- a/src/pipeline.cpp
+++ b/src/pipeline.cpp
@@ -100,6 +100,7 @@ Thomas Knudsen, thokn@sdfe.dk, 2016-05-20
#include <math.h>
#include <stddef.h>
#include <string.h>
+#include <stack>
#include "geodesic.h"
#include "proj.h"
@@ -107,6 +108,8 @@ Thomas Knudsen, thokn@sdfe.dk, 2016-05-20
#include "proj_internal.h"
PROJ_HEAD(pipeline, "Transformation pipeline manager");
+PROJ_HEAD(pop, "Retrieve coordinate value from pipeline stack");
+PROJ_HEAD(push, "Save coordinate value on pipeline stack");
/* Projection specific elements for the PJ object */
namespace { // anonymous namespace
@@ -115,9 +118,16 @@ struct pj_opaque {
char **argv;
char **current_argv;
PJ **pipeline;
+ std::stack<double> *stack[4];
};
-} // anonymous namespace
+struct pj_opaque_pushpop {
+ bool v1;
+ bool v2;
+ bool v3;
+ bool v4;
+};
+} // anonymous namespace
static PJ_COORD pipeline_forward_4d (PJ_COORD point, PJ *P);
@@ -217,7 +227,7 @@ static PJ *destructor (PJ *P, int errlev) {
if (nullptr==P->opaque)
return pj_default_destructor (P, errlev);
- /* Deallocate each pipeine step, then pipeline array */
+ /* Deallocate each pipeline step, then pipeline array */
if (nullptr!=static_cast<struct pj_opaque*>(P->opaque)->pipeline)
for (i = 0; i < static_cast<struct pj_opaque*>(P->opaque)->steps; i++)
proj_destroy (static_cast<struct pj_opaque*>(P->opaque)->pipeline[i+1]);
@@ -226,6 +236,9 @@ static PJ *destructor (PJ *P, int errlev) {
pj_dealloc (static_cast<struct pj_opaque*>(P->opaque)->argv);
pj_dealloc (static_cast<struct pj_opaque*>(P->opaque)->current_argv);
+ for (i=0; i<4; i++)
+ delete static_cast<struct pj_opaque*>(P->opaque)->stack[i];
+
return pj_default_destructor(P, errlev);
}
@@ -384,6 +397,10 @@ PJ *OPERATION(pipeline,0) {
if (nullptr==P->opaque)
return destructor(P, ENOMEM);
+ /* initialize stack */
+ for (i=0; i<4; i++)
+ static_cast<struct pj_opaque*>(P->opaque)->stack[i] = new std::stack<double>;
+
argc = (int)argc_params (P->params);
static_cast<struct pj_opaque*>(P->opaque)->argv = argv = argv_params (P->params, argc);
if (nullptr==argv)
@@ -467,6 +484,7 @@ PJ *OPERATION(pipeline,0) {
proj_log_error (P, "Pipeline: Bad step definition: %s (%s)", current_argv[0], pj_strerrno (err_to_report));
return destructor (P, err_to_report); /* ERROR: bad pipeline def */
}
+ next_step->parent = P;
proj_errno_restore (P, err);
@@ -551,3 +569,92 @@ PJ *OPERATION(pipeline,0) {
P->right = pj_right (static_cast<struct pj_opaque*>(P->opaque)->pipeline[nsteps]);
return P;
}
+
+static PJ_COORD push(PJ_COORD point, PJ *P) {
+ if (P->parent == nullptr)
+ return point;
+
+ struct pj_opaque *pipeline = static_cast<struct pj_opaque*>(P->parent->opaque);
+ struct pj_opaque_pushpop *opaque = static_cast<struct pj_opaque_pushpop*>(P->opaque);
+
+ if (opaque->v1)
+ pipeline->stack[0]->push(point.v[0]);
+ if (opaque->v2)
+ pipeline->stack[1]->push(point.v[1]);
+ if (opaque->v3)
+ pipeline->stack[2]->push(point.v[2]);
+ if (opaque->v4)
+ pipeline->stack[3]->push(point.v[3]);
+
+ return point;
+}
+
+static PJ_COORD pop(PJ_COORD point, PJ *P) {
+ if (P->parent == nullptr)
+ return point;
+
+ struct pj_opaque *pipeline = static_cast<struct pj_opaque*>(P->parent->opaque);
+ struct pj_opaque_pushpop *opaque = static_cast<struct pj_opaque_pushpop*>(P->opaque);
+
+ if (opaque->v1 && !pipeline->stack[0]->empty()) {
+ point.v[0] = pipeline->stack[0]->top();
+ pipeline->stack[0]->pop();
+ }
+
+ if (opaque->v2 && !pipeline->stack[1]->empty()) {
+ point.v[1] = pipeline->stack[1]->top();
+ pipeline->stack[1]->pop();
+ }
+
+ if (opaque->v3 && !pipeline->stack[2]->empty()) {
+ point.v[2] = pipeline->stack[2]->top();
+ pipeline->stack[2]->pop();
+ }
+
+ if (opaque->v4 && !pipeline->stack[3]->empty()) {
+ point.v[3] = pipeline->stack[3]->top();
+ pipeline->stack[3]->pop();
+ }
+
+ return point;
+}
+
+
+
+static PJ *setup_pushpop(PJ *P) {
+ P->opaque = static_cast<struct pj_opaque_pushpop*>(pj_calloc (1, sizeof(struct pj_opaque_pushpop)));
+ if (nullptr==P->opaque)
+ return destructor(P, ENOMEM);
+
+ if (pj_param_exists(P->params, "v_1"))
+ static_cast<struct pj_opaque_pushpop*>(P->opaque)->v1 = true;
+
+ if (pj_param_exists(P->params, "v_2"))
+ static_cast<struct pj_opaque_pushpop*>(P->opaque)->v2 = true;
+
+ if (pj_param_exists(P->params, "v_3"))
+ static_cast<struct pj_opaque_pushpop*>(P->opaque)->v3 = true;
+
+ if (pj_param_exists(P->params, "v_4"))
+ static_cast<struct pj_opaque_pushpop*>(P->opaque)->v4 = true;
+
+ P->left = PJ_IO_UNITS_WHATEVER;
+ P->right = PJ_IO_UNITS_WHATEVER;
+
+ return P;
+}
+
+
+PJ *OPERATION(push, 0) {
+ P->fwd4d = push;
+ P->inv4d = pop;
+
+ return setup_pushpop(P);
+}
+
+PJ *OPERATION(pop, 0) {
+ P->inv4d = push;
+ P->fwd4d = pop;
+
+ return setup_pushpop(P);
+}
diff --git a/src/pj_list.h b/src/pj_list.h
index 3592dcc4..8ab4cdc0 100644
--- a/src/pj_list.h
+++ b/src/pj_list.h
@@ -116,6 +116,8 @@ PROJ_HEAD(pconic, "Perspective Conic")
PROJ_HEAD(patterson, "Patterson Cylindrical")
PROJ_HEAD(pipeline, "Transformation pipeline manager")
PROJ_HEAD(poly, "Polyconic (American)")
+PROJ_HEAD(pop, "Retrieve coordinate value from pipeline stack")
+PROJ_HEAD(push, "Save coordinate value on pipeline stack")
PROJ_HEAD(putp1, "Putnins P1")
PROJ_HEAD(putp2, "Putnins P2")
PROJ_HEAD(putp3, "Putnins P3")
diff --git a/src/proj_internal.h b/src/proj_internal.h
index b7e3de9e..14b69492 100644
--- a/src/proj_internal.h
+++ b/src/proj_internal.h
@@ -321,6 +321,7 @@ struct PJconsts {
const char *descr = nullptr; /* From pj_list.h or individual PJ_*.c file */
paralist *params = nullptr; /* Parameter list */
char *def_full = nullptr; /* Full textual definition (usually 0 - set by proj_pj_info) */
+ PJconsts *parent = nullptr; /* Parent PJ of pipeline steps - nullptr if not a pipeline step */
/* For debugging / logging purposes */
char *def_size = nullptr; /* Shape and size parameters extracted from params */