aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--include/proj/crs.hpp10
-rw-r--r--src/apps/projinfo.cpp4
-rw-r--r--src/iso19111/crs.cpp142
-rw-r--r--src/iso19111/io.cpp17
-rw-r--r--test/unit/test_io.cpp142
5 files changed, 285 insertions, 30 deletions
diff --git a/include/proj/crs.hpp b/include/proj/crs.hpp
index 30134cb4..61844f79 100644
--- a/include/proj/crs.hpp
+++ b/include/proj/crs.hpp
@@ -138,6 +138,8 @@ class PROJ_GCC_DLL CRS : public common::ObjectUsage,
PROJ_INTERNAL CRSNNPtr normalizeForVisualization() const;
+ PROJ_INTERNAL CRSNNPtr allowNonConformantWKT1Export() const;
+
//! @endcond
protected:
@@ -892,6 +894,14 @@ class PROJ_GCC_DLL CompoundCRS final : public CRS,
const std::vector<CRSNNPtr>
&components); // throw InvalidCompoundCRSException
+ //! @cond Doxygen_Suppress
+ PROJ_INTERNAL static CRSNNPtr
+ createLax(const util::PropertyMap &properties,
+ const std::vector<CRSNNPtr> &components,
+ const io::DatabaseContextPtr
+ &dbContext); // throw InvalidCompoundCRSException
+ //! @endcond
+
protected:
// relaxed: standard say SingleCRSNNPtr
PROJ_INTERNAL explicit CompoundCRS(const std::vector<CRSNNPtr> &components);
diff --git a/src/apps/projinfo.cpp b/src/apps/projinfo.cpp
index ed732aaa..cd04433d 100644
--- a/src/apps/projinfo.cpp
+++ b/src/apps/projinfo.cpp
@@ -479,8 +479,8 @@ static void outputObject(
std::cout << "WKT1:GDAL string:" << std::endl;
}
- auto formatter =
- WKTFormatter::create(WKTFormatter::Convention::WKT1_GDAL);
+ auto formatter = WKTFormatter::create(
+ WKTFormatter::Convention::WKT1_GDAL, dbContext);
if (outputOpt.singleLine) {
formatter->setMultiLine(false);
}
diff --git a/src/iso19111/crs.cpp b/src/iso19111/crs.cpp
index b0f5bcf9..77261211 100644
--- a/src/iso19111/crs.cpp
+++ b/src/iso19111/crs.cpp
@@ -93,6 +93,7 @@ struct CRS::Private {
BoundCRSPtr canonicalBoundCRS_{};
std::string extensionProj4_{};
bool implicitCS_ = false;
+ bool allowNonConformantWKT1Export_ = false;
void setImplicitCS(const util::PropertyMap &properties) {
const auto pVal = properties.get("IMPLICIT_CS");
@@ -559,6 +560,18 @@ CRSNNPtr CRS::shallowClone() const { return _shallowClone(); }
//! @cond Doxygen_Suppress
+CRSNNPtr CRS::allowNonConformantWKT1Export() const {
+ auto crs = shallowClone();
+ crs->d->allowNonConformantWKT1Export_ = true;
+ return crs;
+}
+
+//! @endcond
+
+// ---------------------------------------------------------------------------
+
+//! @cond Doxygen_Suppress
+
CRSNNPtr CRS::alterName(const std::string &newName) const {
auto crs = shallowClone();
auto newNameMod(newName);
@@ -1334,30 +1347,49 @@ void GeodeticCRS::_exportToWKT(io::WKTFormatter *formatter) const {
const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2;
const bool isGeographic =
dynamic_cast<const GeographicCRS *>(this) != nullptr;
- formatter->startNode(isWKT2
- ? ((formatter->use2019Keywords() && isGeographic)
- ? io::WKTConstants::GEOGCRS
- : io::WKTConstants::GEODCRS)
- : isGeocentric() ? io::WKTConstants::GEOCCS
- : io::WKTConstants::GEOGCS,
- !identifiers().empty());
- auto l_name = nameStr();
+
const auto &cs = coordinateSystem();
const auto &axisList = cs->axisList();
-
const auto oldAxisOutputRule = formatter->outputAxis();
+ auto l_name = nameStr();
+ const auto &dbContext = formatter->databaseContext();
if (formatter->useESRIDialect()) {
if (axisList.size() != 2) {
io::FormattingException::Throw(
"Only export of Geographic 2D CRS is supported in WKT1_ESRI");
}
+ }
+
+ if (!isWKT2 && formatter->isStrict() && isGeographic &&
+ axisList.size() != 2 &&
+ oldAxisOutputRule != io::WKTFormatter::OutputAxisRule::NO) {
+ if (CRS::getPrivate()->allowNonConformantWKT1Export_) {
+ formatter->startNode(io::WKTConstants::COMPD_CS, false);
+ formatter->addQuotedString(l_name + " + " + l_name);
+ auto geogCRS = demoteTo2D(std::string(), dbContext);
+ geogCRS->_exportToWKT(formatter);
+ geogCRS->_exportToWKT(formatter);
+ formatter->endNode();
+ return;
+ }
+ io::FormattingException::Throw(
+ "WKT1 does not support Geographic 3D CRS.");
+ }
+ formatter->startNode(isWKT2
+ ? ((formatter->use2019Keywords() && isGeographic)
+ ? io::WKTConstants::GEOGCRS
+ : io::WKTConstants::GEODCRS)
+ : isGeocentric() ? io::WKTConstants::GEOCCS
+ : io::WKTConstants::GEOGCS,
+ !identifiers().empty());
+
+ if (formatter->useESRIDialect()) {
if (l_name == "WGS 84") {
l_name = "GCS_WGS_1984";
} else {
bool aliasFound = false;
- const auto &dbContext = formatter->databaseContext();
if (dbContext) {
auto l_alias = dbContext->getAliasFromOfficialName(
l_name, "geodetic_crs", "ESRI");
@@ -1373,11 +1405,6 @@ void GeodeticCRS::_exportToWKT(io::WKTFormatter *formatter) const {
}
}
}
- } else if (!isWKT2 && formatter->isStrict() && isGeographic &&
- axisList.size() != 2 &&
- oldAxisOutputRule != io::WKTFormatter::OutputAxisRule::NO) {
- io::FormattingException::Throw(
- "WKT1 does not support Geographic 3D CRS.");
}
if (!isWKT2 && !formatter->useESRIDialect() && isDeprecated()) {
@@ -3163,6 +3190,36 @@ void ProjectedCRS::_exportToWKT(io::WKTFormatter *formatter) const {
const auto &dbContext = formatter->databaseContext();
auto l_name = nameStr();
+ const auto &l_coordinateSystem = d->coordinateSystem();
+ const auto &axisList = l_coordinateSystem->axisList();
+ if (axisList.size() == 3 && !(isWKT2 && formatter->use2019Keywords())) {
+ if (!formatter->useESRIDialect() &&
+ CRS::getPrivate()->allowNonConformantWKT1Export_) {
+ formatter->startNode(io::WKTConstants::COMPD_CS, false);
+ formatter->addQuotedString(l_name + " + " + baseCRS()->nameStr());
+ auto projCRS2D = demoteTo2D(std::string(), dbContext);
+ if (dbContext) {
+ const auto res =
+ projCRS2D->identify(io::AuthorityFactory::create(
+ NN_NO_CHECK(dbContext), "EPSG"));
+ if (res.size() == 1) {
+ const auto &front = res.front();
+ if (front.second == 100) {
+ projCRS2D = front.first;
+ }
+ }
+ }
+ projCRS2D->_exportToWKT(formatter);
+ baseCRS()
+ ->demoteTo2D(std::string(), dbContext)
+ ->_exportToWKT(formatter);
+ formatter->endNode();
+ return;
+ }
+ io::FormattingException::Throw(
+ "Projected 3D CRS can only be exported since WKT2:2019");
+ }
+
std::string l_alias;
if (formatter->useESRIDialect() && dbContext) {
l_alias = dbContext->getAliasFromOfficialName(l_name, "projected_crs",
@@ -3214,13 +3271,6 @@ void ProjectedCRS::_exportToWKT(io::WKTFormatter *formatter) const {
}
}
- const auto &l_coordinateSystem = d->coordinateSystem();
- const auto &axisList = l_coordinateSystem->axisList();
- if (axisList.size() == 3 && !(isWKT2 && formatter->use2019Keywords())) {
- io::FormattingException::Throw(
- "Projected 3D CRS can only be exported since WKT2:2019");
- }
-
const auto exportAxis = [&l_coordinateSystem, &axisList, &formatter]() {
const auto oldAxisOutputRule = formatter->outputAxis();
if (oldAxisOutputRule ==
@@ -4091,6 +4141,54 @@ CompoundCRSNNPtr CompoundCRS::create(const util::PropertyMap &properties,
// ---------------------------------------------------------------------------
//! @cond Doxygen_Suppress
+
+/** \brief Instantiate a CompoundCRS, a Geographic 3D CRS or a Projected CRS
+ * from a vector of CRS.
+ *
+ * Be a bit "lax", in allowing formulations like EPSG:4326+4326 or
+ * EPSG:32631+4326 to express Geographic 3D CRS / Projected3D CRS.
+ *
+ * @param properties See \ref general_properties.
+ * At minimum the name should be defined.
+ * @param components the component CRS of the CompoundCRS.
+ * @return new CRS.
+ * @throw InvalidCompoundCRSException
+ */
+CRSNNPtr CompoundCRS::createLax(const util::PropertyMap &properties,
+ const std::vector<CRSNNPtr> &components,
+ const io::DatabaseContextPtr &dbContext) {
+
+ if (components.size() == 2) {
+ auto comp0 = components[0].get();
+ auto comp1 = components[1].get();
+ auto comp0Geog = dynamic_cast<const GeographicCRS *>(comp0);
+ auto comp0Proj = dynamic_cast<const ProjectedCRS *>(comp0);
+ auto comp1Geog = dynamic_cast<const GeographicCRS *>(comp1);
+ if ((comp0Geog != nullptr || comp0Proj != nullptr) &&
+ comp1Geog != nullptr) {
+ const auto horizGeog =
+ (comp0Proj != nullptr)
+ ? comp0Proj->baseCRS().as_nullable().get()
+ : comp0Geog;
+ if (horizGeog->_isEquivalentTo(
+ comp1Geog->demoteTo2D(std::string(), dbContext).get())) {
+ return components[0]
+ ->promoteTo3D(std::string(), dbContext)
+ ->allowNonConformantWKT1Export();
+ }
+ throw InvalidCompoundCRSException(
+ "The 'vertical' geographic CRS is not equivalent to the "
+ "geographic CRS of the horizontal part");
+ }
+ }
+
+ return create(properties, components);
+}
+//! @endcond
+
+// ---------------------------------------------------------------------------
+
+//! @cond Doxygen_Suppress
void CompoundCRS::_exportToWKT(io::WKTFormatter *formatter) const {
const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2;
formatter->startNode(isWKT2 ? io::WKTConstants::COMPOUNDCRS
diff --git a/src/iso19111/io.cpp b/src/iso19111/io.cpp
index e516e4b5..59116ad6 100644
--- a/src/iso19111/io.cpp
+++ b/src/iso19111/io.cpp
@@ -1334,7 +1334,7 @@ struct WKTParser::Private {
DerivedVerticalCRSNNPtr buildDerivedVerticalCRS(const WKTNodeNNPtr &node);
- CompoundCRSNNPtr buildCompoundCRS(const WKTNodeNNPtr &node);
+ CRSNNPtr buildCompoundCRS(const WKTNodeNNPtr &node);
BoundCRSNNPtr buildBoundCRS(const WKTNodeNNPtr &node);
@@ -4157,8 +4157,7 @@ WKTParser::Private::buildDerivedVerticalCRS(const WKTNodeNNPtr &node) {
// ---------------------------------------------------------------------------
-CompoundCRSNNPtr
-WKTParser::Private::buildCompoundCRS(const WKTNodeNNPtr &node) {
+CRSNNPtr WKTParser::Private::buildCompoundCRS(const WKTNodeNNPtr &node) {
std::vector<CRSNNPtr> components;
for (const auto &child : node->GP()->children()) {
auto crs = buildCRS(child);
@@ -4166,7 +4165,13 @@ WKTParser::Private::buildCompoundCRS(const WKTNodeNNPtr &node) {
components.push_back(NN_NO_CHECK(crs));
}
}
- return CompoundCRS::create(buildProperties(node), components);
+
+ if (ci_equal(node->GP()->value(), WKTConstants::COMPD_CS)) {
+ return CompoundCRS::createLax(buildProperties(node), components,
+ dbContext_);
+ } else {
+ return CompoundCRS::create(buildProperties(node), components);
+ }
}
// ---------------------------------------------------------------------------
@@ -5851,11 +5856,11 @@ static BaseObjectNNPtr createFromUserInput(const std::string &text,
tokensCode[0], false));
auto crs2(factory->createCoordinateReferenceSystem(
tokensCode[1], false));
- return CompoundCRS::create(
+ return CompoundCRS::createLax(
util::PropertyMap().set(
IdentifiedObject::NAME_KEY,
crs1->nameStr() + " + " + crs2->nameStr()),
- {crs1, crs2});
+ {crs1, crs2}, dbContext);
}
throw;
}
diff --git a/test/unit/test_io.cpp b/test/unit/test_io.cpp
index a227ff31..da26aa2d 100644
--- a/test/unit/test_io.cpp
+++ b/test/unit/test_io.cpp
@@ -2430,6 +2430,53 @@ TEST(wkt_parse, COMPD_CS) {
// ---------------------------------------------------------------------------
+TEST(wkt_parse, COMPD_CS_non_conformant_horizontal_plus_horizontal_as_in_LAS) {
+ auto obj = WKTParser().createFromWKT(
+ "COMPD_CS[\"horizontal + vertical\",\n"
+ " PROJCS[\"WGS 84 / UTM zone 31N\",\n"
+ " GEOGCS[\"WGS 84\",\n"
+ " DATUM[\"World Geodetic System 1984\",\n"
+ " SPHEROID[\"WGS 84\",6378137,298.257223563,\n"
+ " AUTHORITY[\"EPSG\",\"7030\"]],\n"
+ " AUTHORITY[\"EPSG\",\"6326\"]],\n"
+ " PRIMEM[\"Greenwich\",0,\n"
+ " AUTHORITY[\"EPSG\",\"8901\"]],\n"
+ " UNIT[\"degree\",0.0174532925199433,\n"
+ " AUTHORITY[\"EPSG\",\"9122\"]],\n"
+ " AXIS[\"Latitude\",NORTH],\n"
+ " AXIS[\"Longitude\",EAST],\n"
+ " AUTHORITY[\"EPSG\",\"4326\"]],\n"
+ " PROJECTION[\"Transverse_Mercator\"],\n"
+ " PARAMETER[\"latitude_of_origin\",0],\n"
+ " PARAMETER[\"central_meridian\",3],\n"
+ " PARAMETER[\"scale_factor\",0.9996],\n"
+ " PARAMETER[\"false_easting\",500000],\n"
+ " PARAMETER[\"false_northing\",0],\n"
+ " UNIT[\"metre\",1,\n"
+ " AUTHORITY[\"EPSG\",\"9001\"]],\n"
+ " AXIS[\"Easting\",EAST],\n"
+ " AXIS[\"Northing\",NORTH],\n"
+ " AUTHORITY[\"EPSG\",\"32631\"]],\n"
+ " GEOGCS[\"WGS 84\",\n"
+ " DATUM[\"World Geodetic System 1984\",\n"
+ " SPHEROID[\"WGS 84\",6378137,298.257223563,\n"
+ " AUTHORITY[\"EPSG\",\"7030\"]],\n"
+ " AUTHORITY[\"EPSG\",\"6326\"]],\n"
+ " PRIMEM[\"Greenwich\",0,\n"
+ " AUTHORITY[\"EPSG\",\"8901\"]],\n"
+ " UNIT[\"degree\",0.0174532925199433,\n"
+ " AUTHORITY[\"EPSG\",\"9122\"]],\n"
+ " AXIS[\"Latitude\",NORTH],\n"
+ " AXIS[\"Longitude\",EAST],\n"
+ " AUTHORITY[\"EPSG\",\"4326\"]]]");
+ auto crs = nn_dynamic_pointer_cast<ProjectedCRS>(obj);
+ ASSERT_TRUE(crs != nullptr);
+ EXPECT_EQ(crs->nameStr(), "WGS 84 / UTM zone 31N");
+ EXPECT_EQ(crs->coordinateSystem()->axisList().size(), 3U);
+}
+
+// ---------------------------------------------------------------------------
+
TEST(wkt_parse, COORDINATEOPERATION) {
std::string src_wkt;
@@ -9490,6 +9537,101 @@ TEST(io, createFromUserInput) {
EXPECT_EQ(crs->derivingConversion()->getEPSGCode(), 16031);
}
+ // We accept non-conformant EPSG:4326+4326
+ {
+ auto obj = createFromUserInput("EPSG:4326+4326", dbContext);
+ auto crs = nn_dynamic_pointer_cast<GeographicCRS>(obj);
+ ASSERT_TRUE(crs != nullptr);
+ EXPECT_EQ(crs->nameStr(), "WGS 84");
+ EXPECT_EQ(crs->getEPSGCode(), 4979);
+ EXPECT_EQ(crs->coordinateSystem()->axisList().size(), 3U);
+
+ const auto wkt =
+ "COMPD_CS[\"WGS 84 + WGS 84\",\n"
+ " GEOGCS[\"WGS 84\",\n"
+ " DATUM[\"WGS_1984\",\n"
+ " SPHEROID[\"WGS 84\",6378137,298.257223563,\n"
+ " AUTHORITY[\"EPSG\",\"7030\"]],\n"
+ " AUTHORITY[\"EPSG\",\"6326\"]],\n"
+ " PRIMEM[\"Greenwich\",0,\n"
+ " AUTHORITY[\"EPSG\",\"8901\"]],\n"
+ " UNIT[\"degree\",0.0174532925199433,\n"
+ " AUTHORITY[\"EPSG\",\"9122\"]],\n"
+ " AUTHORITY[\"EPSG\",\"4326\"]],\n"
+ " GEOGCS[\"WGS 84\",\n"
+ " DATUM[\"WGS_1984\",\n"
+ " SPHEROID[\"WGS 84\",6378137,298.257223563,\n"
+ " AUTHORITY[\"EPSG\",\"7030\"]],\n"
+ " AUTHORITY[\"EPSG\",\"6326\"]],\n"
+ " PRIMEM[\"Greenwich\",0,\n"
+ " AUTHORITY[\"EPSG\",\"8901\"]],\n"
+ " UNIT[\"degree\",0.0174532925199433,\n"
+ " AUTHORITY[\"EPSG\",\"9122\"]],\n"
+ " AUTHORITY[\"EPSG\",\"4326\"]]]";
+
+ EXPECT_EQ(
+ crs->exportToWKT(WKTFormatter::create(
+ WKTFormatter::Convention::WKT1_GDAL, dbContext)
+ .get()),
+ wkt);
+ }
+
+ // Non consistent
+ EXPECT_THROW(createFromUserInput("EPSG:4326+4258", dbContext),
+ InvalidCompoundCRSException);
+
+ // We accept non-conformant EPSG:32631+4326
+ {
+ auto obj = createFromUserInput("EPSG:32631+4326", dbContext);
+ auto crs = nn_dynamic_pointer_cast<ProjectedCRS>(obj);
+ ASSERT_TRUE(crs != nullptr);
+ EXPECT_EQ(crs->nameStr(), "WGS 84 / UTM zone 31N");
+ EXPECT_EQ(crs->baseCRS()->getEPSGCode(), 4979);
+ EXPECT_EQ(crs->baseCRS()->coordinateSystem()->axisList().size(), 3U);
+ EXPECT_EQ(crs->coordinateSystem()->axisList().size(), 3U);
+
+ const auto wkt =
+ "COMPD_CS[\"WGS 84 / UTM zone 31N + WGS 84\",\n"
+ " PROJCS[\"WGS 84 / UTM zone 31N\",\n"
+ " GEOGCS[\"WGS 84\",\n"
+ " DATUM[\"WGS_1984\",\n"
+ " SPHEROID[\"WGS 84\",6378137,298.257223563,\n"
+ " AUTHORITY[\"EPSG\",\"7030\"]],\n"
+ " AUTHORITY[\"EPSG\",\"6326\"]],\n"
+ " PRIMEM[\"Greenwich\",0,\n"
+ " AUTHORITY[\"EPSG\",\"8901\"]],\n"
+ " UNIT[\"degree\",0.0174532925199433,\n"
+ " AUTHORITY[\"EPSG\",\"9122\"]],\n"
+ " AUTHORITY[\"EPSG\",\"4326\"]],\n"
+ " PROJECTION[\"Transverse_Mercator\"],\n"
+ " PARAMETER[\"latitude_of_origin\",0],\n"
+ " PARAMETER[\"central_meridian\",3],\n"
+ " PARAMETER[\"scale_factor\",0.9996],\n"
+ " PARAMETER[\"false_easting\",500000],\n"
+ " PARAMETER[\"false_northing\",0],\n"
+ " UNIT[\"metre\",1,\n"
+ " AUTHORITY[\"EPSG\",\"9001\"]],\n"
+ " AXIS[\"Easting\",EAST],\n"
+ " AXIS[\"Northing\",NORTH],\n"
+ " AUTHORITY[\"EPSG\",\"32631\"]],\n"
+ " GEOGCS[\"WGS 84\",\n"
+ " DATUM[\"WGS_1984\",\n"
+ " SPHEROID[\"WGS 84\",6378137,298.257223563,\n"
+ " AUTHORITY[\"EPSG\",\"7030\"]],\n"
+ " AUTHORITY[\"EPSG\",\"6326\"]],\n"
+ " PRIMEM[\"Greenwich\",0,\n"
+ " AUTHORITY[\"EPSG\",\"8901\"]],\n"
+ " UNIT[\"degree\",0.0174532925199433,\n"
+ " AUTHORITY[\"EPSG\",\"9122\"]],\n"
+ " AUTHORITY[\"EPSG\",\"4326\"]]]";
+
+ EXPECT_EQ(
+ crs->exportToWKT(WKTFormatter::create(
+ WKTFormatter::Convention::WKT1_GDAL, dbContext)
+ .get()),
+ wkt);
+ }
+
EXPECT_THROW(createFromUserInput(
"urn:ogc:def:crs,crs:EPSG::4979,"
"cs:PROJ::ENh,"