aboutsummaryrefslogtreecommitdiff
path: root/src/iso19111/crs.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'src/iso19111/crs.cpp')
-rw-r--r--src/iso19111/crs.cpp176
1 files changed, 110 insertions, 66 deletions
diff --git a/src/iso19111/crs.cpp b/src/iso19111/crs.cpp
index 7c8fcd81..be593878 100644
--- a/src/iso19111/crs.cpp
+++ b/src/iso19111/crs.cpp
@@ -1608,7 +1608,7 @@ GeodeticCRS::velocityModel() PROJ_PURE_DEFN {
// ---------------------------------------------------------------------------
-/** \brief Return whether the CRS is a geocentric one.
+/** \brief Return whether the CRS is a Cartesian geocentric one.
*
* A geocentric CRS is a geodetic CRS that has a Cartesian coordinate system
* with three axis, whose direction is respectively
@@ -1629,6 +1629,31 @@ bool GeodeticCRS::isGeocentric() PROJ_PURE_DEFN {
// ---------------------------------------------------------------------------
+/** \brief Return whether the CRS is a Spherical planetocentric one.
+ *
+ * A Spherical planetocentric CRS is a geodetic CRS that has a spherical
+ * (angular) coordinate system with 2 axis, which represent geocentric latitude/
+ * longitude or longitude/geocentric latitude.
+ *
+ * Such CRS are typically used in use case that apply to non-Earth bodies.
+ *
+ * @return true if the CRS is a Spherical planetocentric CRS.
+ *
+ * @since 8.2
+ */
+bool GeodeticCRS::isSphericalPlanetocentric() PROJ_PURE_DEFN {
+ const auto &cs = coordinateSystem();
+ const auto &axisList = cs->axisList();
+ return axisList.size() == 2 &&
+ dynamic_cast<cs::SphericalCS *>(cs.get()) != nullptr &&
+ ((ci_equal(axisList[0]->nameStr(), "planetocentric latitude") &&
+ ci_equal(axisList[1]->nameStr(), "planetocentric longitude")) ||
+ (ci_equal(axisList[0]->nameStr(), "planetocentric longitude") &&
+ ci_equal(axisList[1]->nameStr(), "planetocentric latitude")));
+}
+
+// ---------------------------------------------------------------------------
+
/** \brief Instantiate a GeodeticCRS from a datum::GeodeticReferenceFrame and a
* cs::SphericalCS.
*
@@ -1971,6 +1996,61 @@ void GeodeticCRS::addGeocentricUnitConversionIntoPROJString(
// ---------------------------------------------------------------------------
//! @cond Doxygen_Suppress
+void GeodeticCRS::addAngularUnitConvertAndAxisSwap(
+ io::PROJStringFormatter *formatter) const {
+ const auto &axisList = coordinateSystem()->axisList();
+
+ formatter->addStep("unitconvert");
+ formatter->addParam("xy_in", "rad");
+ if (axisList.size() == 3 && !formatter->omitZUnitConversion()) {
+ formatter->addParam("z_in", "m");
+ }
+ {
+ const auto &unitHoriz = axisList[0]->unit();
+ const auto projUnit = unitHoriz.exportToPROJString();
+ if (projUnit.empty()) {
+ formatter->addParam("xy_out", unitHoriz.conversionToSI());
+ } else {
+ formatter->addParam("xy_out", projUnit);
+ }
+ }
+ if (axisList.size() == 3 && !formatter->omitZUnitConversion()) {
+ const auto &unitZ = axisList[2]->unit();
+ auto projVUnit = unitZ.exportToPROJString();
+ if (projVUnit.empty()) {
+ formatter->addParam("z_out", unitZ.conversionToSI());
+ } else {
+ formatter->addParam("z_out", projVUnit);
+ }
+ }
+
+ const char *order[2] = {nullptr, nullptr};
+ const char *one = "1";
+ const char *two = "2";
+ for (int i = 0; i < 2; i++) {
+ const auto &dir = axisList[i]->direction();
+ if (&dir == &cs::AxisDirection::WEST) {
+ order[i] = "-1";
+ } else if (&dir == &cs::AxisDirection::EAST) {
+ order[i] = one;
+ } else if (&dir == &cs::AxisDirection::SOUTH) {
+ order[i] = "-2";
+ } else if (&dir == &cs::AxisDirection::NORTH) {
+ order[i] = two;
+ }
+ }
+ if (order[0] && order[1] && (order[0] != one || order[1] != two)) {
+ formatter->addStep("axisswap");
+ char orderStr[10];
+ sprintf(orderStr, "%.2s,%.2s", order[0], order[1]);
+ formatter->addParam("order", orderStr);
+ }
+}
+//! @endcond
+
+// ---------------------------------------------------------------------------
+
+//! @cond Doxygen_Suppress
void GeodeticCRS::_exportToPROJString(
io::PROJStringFormatter *formatter) const // throw(io::FormattingException)
{
@@ -1982,19 +2062,38 @@ void GeodeticCRS::_exportToPROJString(
return;
}
- if (!isGeocentric()) {
- io::FormattingException::Throw(
- "GeodeticCRS::exportToPROJString() only "
- "supports geocentric coordinate systems");
- }
+ if (isGeocentric()) {
+ if (!formatter->getCRSExport()) {
+ formatter->addStep("cart");
+ } else {
+ formatter->addStep("geocent");
+ }
- if (!formatter->getCRSExport()) {
- formatter->addStep("cart");
+ addDatumInfoToPROJString(formatter);
+ addGeocentricUnitConversionIntoPROJString(formatter);
+ } else if (isSphericalPlanetocentric()) {
+ if (!formatter->getCRSExport()) {
+ formatter->addStep("geoc");
+ addDatumInfoToPROJString(formatter);
+ addAngularUnitConvertAndAxisSwap(formatter);
+ } else {
+ io::FormattingException::Throw(
+ "GeodeticCRS::exportToPROJString() not supported on spherical "
+ "planetocentric coordinate systems");
+ // The below code now works as input to PROJ, but I'm not sure we
+ // want to propagate this, given that we got cs2cs doing conversion
+ // in the wrong direction in past versions.
+ /*formatter->addStep("longlat");
+ formatter->addParam("geoc");
+
+ addDatumInfoToPROJString(formatter);*/
+ }
} else {
- formatter->addStep("geocent");
+ io::FormattingException::Throw(
+ "GeodeticCRS::exportToPROJString() only "
+ "supports geocentric or spherical planetocentric "
+ "coordinate systems");
}
- addDatumInfoToPROJString(formatter);
- addGeocentricUnitConversionIntoPROJString(formatter);
}
//! @endcond
@@ -2859,61 +2958,6 @@ GeographicCRS::demoteTo2D(const std::string &newName,
// ---------------------------------------------------------------------------
//! @cond Doxygen_Suppress
-void GeographicCRS::addAngularUnitConvertAndAxisSwap(
- io::PROJStringFormatter *formatter) const {
- const auto &axisList = coordinateSystem()->axisList();
-
- formatter->addStep("unitconvert");
- formatter->addParam("xy_in", "rad");
- if (axisList.size() == 3 && !formatter->omitZUnitConversion()) {
- formatter->addParam("z_in", "m");
- }
- {
- const auto &unitHoriz = axisList[0]->unit();
- const auto projUnit = unitHoriz.exportToPROJString();
- if (projUnit.empty()) {
- formatter->addParam("xy_out", unitHoriz.conversionToSI());
- } else {
- formatter->addParam("xy_out", projUnit);
- }
- }
- if (axisList.size() == 3 && !formatter->omitZUnitConversion()) {
- const auto &unitZ = axisList[2]->unit();
- auto projVUnit = unitZ.exportToPROJString();
- if (projVUnit.empty()) {
- formatter->addParam("z_out", unitZ.conversionToSI());
- } else {
- formatter->addParam("z_out", projVUnit);
- }
- }
-
- const char *order[2] = {nullptr, nullptr};
- const char *one = "1";
- const char *two = "2";
- for (int i = 0; i < 2; i++) {
- const auto &dir = axisList[i]->direction();
- if (&dir == &cs::AxisDirection::WEST) {
- order[i] = "-1";
- } else if (&dir == &cs::AxisDirection::EAST) {
- order[i] = one;
- } else if (&dir == &cs::AxisDirection::SOUTH) {
- order[i] = "-2";
- } else if (&dir == &cs::AxisDirection::NORTH) {
- order[i] = two;
- }
- }
- if (order[0] && order[1] && (order[0] != one || order[1] != two)) {
- formatter->addStep("axisswap");
- char orderStr[10];
- sprintf(orderStr, "%.2s,%.2s", order[0], order[1]);
- formatter->addParam("order", orderStr);
- }
-}
-//! @endcond
-
-// ---------------------------------------------------------------------------
-
-//! @cond Doxygen_Suppress
void GeographicCRS::_exportToPROJString(
io::PROJStringFormatter *formatter) const // throw(io::FormattingException)
{