/****************************************************************************** * * Project: PROJ * Purpose: ISO19111:2019 implementation * Author: Even Rouault * ****************************************************************************** * Copyright (c) 2018, Even Rouault * * Permission is hereby granted, free of charge, to any person obtaining a * copy of this software and associated documentation files (the "Software"), * to deal in the Software without restriction, including without limitation * the rights to use, copy, modify, merge, publish, distribute, sublicense, * and/or sell copies of the Software, and to permit persons to whom the * Software is furnished to do so, subject to the following conditions: * * The above copyright notice and this permission notice shall be included * in all copies or substantial portions of the Software. * * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER * DEALINGS IN THE SOFTWARE. ****************************************************************************/ #ifndef FROM_PROJ_CPP #define FROM_PROJ_CPP #endif //! @cond Doxygen_Suppress #define DO_NOT_DEFINE_EXTERN_DERIVED_CRS_TEMPLATE //! @endcond #include "proj/crs.hpp" #include "proj/common.hpp" #include "proj/coordinateoperation.hpp" #include "proj/coordinatesystem.hpp" #include "proj/io.hpp" #include "proj/metadata.hpp" #include "proj/util.hpp" #include "proj/internal/coordinatesystem_internal.hpp" #include "proj/internal/internal.hpp" #include "proj/internal/io_internal.hpp" #include "operation/oputils.hpp" #include "proj_constants.h" #include "proj_json_streaming_writer.hpp" #include #include #include #include #include #include #include #include using namespace NS_PROJ::internal; #if 0 namespace dropbox{ namespace oxygen { template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; template<> nn::~nn() = default; }} #endif NS_PROJ_START namespace crs { // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct CRS::Private { BoundCRSPtr canonicalBoundCRS_{}; std::string extensionProj4_{}; bool implicitCS_ = false; bool allowNonConformantWKT1Export_ = false; // for what was initially a COMPD_CS with a VERT_CS with a datum type == // ellipsoidal height / 2002 CompoundCRSPtr originalCompoundCRS_{}; void setImplicitCS(const util::PropertyMap &properties) { const auto pVal = properties.get("IMPLICIT_CS"); if (pVal) { if (const auto genVal = dynamic_cast(pVal->get())) { if (genVal->type() == util::BoxedValue::Type::BOOLEAN && genVal->booleanValue()) { implicitCS_ = true; } } } } }; //! @endcond // --------------------------------------------------------------------------- CRS::CRS() : d(internal::make_unique()) {} // --------------------------------------------------------------------------- CRS::CRS(const CRS &other) : ObjectUsage(other), d(internal::make_unique(*(other.d))) {} // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRS::~CRS() = default; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress /** \brief Return whether the CRS has an implicit coordinate system * (e.g from ESRI WKT) */ bool CRS::hasImplicitCS() const { return d->implicitCS_; } //! @endcond // --------------------------------------------------------------------------- /** \brief Return the BoundCRS potentially attached to this CRS. * * In the case this method is called on a object returned by * BoundCRS::baseCRSWithCanonicalBoundCRS(), this method will return this * BoundCRS * * @return a BoundCRSPtr, that might be null. */ const BoundCRSPtr &CRS::canonicalBoundCRS() PROJ_PURE_DEFN { return d->canonicalBoundCRS_; } // --------------------------------------------------------------------------- /** \brief Return the GeodeticCRS of the CRS. * * Returns the GeodeticCRS contained in a CRS. This works currently with * input parameters of type GeodeticCRS or derived, ProjectedCRS, * CompoundCRS or BoundCRS. * * @return a GeodeticCRSPtr, that might be null. */ GeodeticCRSPtr CRS::extractGeodeticCRS() const { auto raw = extractGeodeticCRSRaw(); if (raw) { return std::dynamic_pointer_cast( raw->shared_from_this().as_nullable()); } return nullptr; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress const GeodeticCRS *CRS::extractGeodeticCRSRaw() const { auto geodCRS = dynamic_cast(this); if (geodCRS) { return geodCRS; } auto projCRS = dynamic_cast(this); if (projCRS) { return projCRS->baseCRS()->extractGeodeticCRSRaw(); } auto compoundCRS = dynamic_cast(this); if (compoundCRS) { for (const auto &subCrs : compoundCRS->componentReferenceSystems()) { auto retGeogCRS = subCrs->extractGeodeticCRSRaw(); if (retGeogCRS) { return retGeogCRS; } } } auto boundCRS = dynamic_cast(this); if (boundCRS) { return boundCRS->baseCRS()->extractGeodeticCRSRaw(); } return nullptr; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress const std::string &CRS::getExtensionProj4() const noexcept { return d->extensionProj4_; } //! @endcond // --------------------------------------------------------------------------- /** \brief Return the GeographicCRS of the CRS. * * Returns the GeographicCRS contained in a CRS. This works currently with * input parameters of type GeographicCRS or derived, ProjectedCRS, * CompoundCRS or BoundCRS. * * @return a GeographicCRSPtr, that might be null. */ GeographicCRSPtr CRS::extractGeographicCRS() const { auto raw = extractGeodeticCRSRaw(); if (raw) { return std::dynamic_pointer_cast( raw->shared_from_this().as_nullable()); } return nullptr; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static util::PropertyMap createPropertyMap(const common::IdentifiedObject *obj) { auto props = util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, obj->nameStr()); if (obj->isDeprecated()) { props.set(common::IdentifiedObject::DEPRECATED_KEY, true); } return props; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRSNNPtr CRS::alterGeodeticCRS(const GeodeticCRSNNPtr &newGeodCRS) const { auto geodCRS = dynamic_cast(this); if (geodCRS) { return newGeodCRS; } auto projCRS = dynamic_cast(this); if (projCRS) { return ProjectedCRS::create(createPropertyMap(this), newGeodCRS, projCRS->derivingConversion(), projCRS->coordinateSystem()); } auto compoundCRS = dynamic_cast(this); if (compoundCRS) { std::vector components; for (const auto &subCrs : compoundCRS->componentReferenceSystems()) { components.emplace_back(subCrs->alterGeodeticCRS(newGeodCRS)); } return CompoundCRS::create(createPropertyMap(this), components); } return NN_NO_CHECK( std::dynamic_pointer_cast(shared_from_this().as_nullable())); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRSNNPtr CRS::alterCSLinearUnit(const common::UnitOfMeasure &unit) const { { auto projCRS = dynamic_cast(this); if (projCRS) { return ProjectedCRS::create( createPropertyMap(this), projCRS->baseCRS(), projCRS->derivingConversion(), projCRS->coordinateSystem()->alterUnit(unit)); } } { auto geodCRS = dynamic_cast(this); if (geodCRS && geodCRS->isGeocentric()) { auto cs = dynamic_cast( geodCRS->coordinateSystem().get()); assert(cs); return GeodeticCRS::create( createPropertyMap(this), geodCRS->datum(), geodCRS->datumEnsemble(), cs->alterUnit(unit)); } } { auto geogCRS = dynamic_cast(this); if (geogCRS && geogCRS->coordinateSystem()->axisList().size() == 3) { return GeographicCRS::create( createPropertyMap(this), geogCRS->datum(), geogCRS->datumEnsemble(), geogCRS->coordinateSystem()->alterLinearUnit(unit)); } } { auto vertCRS = dynamic_cast(this); if (vertCRS) { return VerticalCRS::create( createPropertyMap(this), vertCRS->datum(), vertCRS->datumEnsemble(), vertCRS->coordinateSystem()->alterUnit(unit)); } } { auto engCRS = dynamic_cast(this); if (engCRS) { auto cartCS = util::nn_dynamic_pointer_cast( engCRS->coordinateSystem()); if (cartCS) { return EngineeringCRS::create(createPropertyMap(this), engCRS->datum(), cartCS->alterUnit(unit)); } else { auto vertCS = util::nn_dynamic_pointer_cast( engCRS->coordinateSystem()); if (vertCS) { return EngineeringCRS::create(createPropertyMap(this), engCRS->datum(), vertCS->alterUnit(unit)); } } } } return NN_NO_CHECK( std::dynamic_pointer_cast(shared_from_this().as_nullable())); } //! @endcond // --------------------------------------------------------------------------- /** \brief Return the VerticalCRS of the CRS. * * Returns the VerticalCRS contained in a CRS. This works currently with * input parameters of type VerticalCRS or derived, CompoundCRS or BoundCRS. * * @return a VerticalCRSPtr, that might be null. */ VerticalCRSPtr CRS::extractVerticalCRS() const { auto vertCRS = dynamic_cast(this); if (vertCRS) { return std::dynamic_pointer_cast( shared_from_this().as_nullable()); } auto compoundCRS = dynamic_cast(this); if (compoundCRS) { for (const auto &subCrs : compoundCRS->componentReferenceSystems()) { auto retVertCRS = subCrs->extractVerticalCRS(); if (retVertCRS) { return retVertCRS; } } } auto boundCRS = dynamic_cast(this); if (boundCRS) { return boundCRS->baseCRS()->extractVerticalCRS(); } return nullptr; } // --------------------------------------------------------------------------- /** \brief Returns potentially * a BoundCRS, with a transformation to EPSG:4326, wrapping this CRS * * If no such BoundCRS is possible, the object will be returned. * * The purpose of this method is to be able to format a PROJ.4 string with * a +towgs84 parameter or a WKT1:GDAL string with a TOWGS node. * * This method will fetch the GeographicCRS of this CRS and find a * transformation to EPSG:4326 using the domain of the validity of the main CRS, * and there's only one Helmert transformation. * * @return a CRS. */ CRSNNPtr CRS::createBoundCRSToWGS84IfPossible( const io::DatabaseContextPtr &dbContext, operation::CoordinateOperationContext::IntermediateCRSUse allowIntermediateCRSUse) const { auto thisAsCRS = NN_NO_CHECK( std::static_pointer_cast(shared_from_this().as_nullable())); auto boundCRS = util::nn_dynamic_pointer_cast(thisAsCRS); if (!boundCRS) { boundCRS = canonicalBoundCRS(); } if (boundCRS) { if (boundCRS->hubCRS()->_isEquivalentTo( GeographicCRS::EPSG_4326.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { return NN_NO_CHECK(boundCRS); } } auto compoundCRS = dynamic_cast(this); if (compoundCRS) { const auto &comps = compoundCRS->componentReferenceSystems(); if (comps.size() == 2) { auto horiz = comps[0]->createBoundCRSToWGS84IfPossible( dbContext, allowIntermediateCRSUse); auto vert = comps[1]->createBoundCRSToWGS84IfPossible( dbContext, allowIntermediateCRSUse); if (horiz.get() != comps[0].get() || vert.get() != comps[1].get()) { return CompoundCRS::create(createPropertyMap(this), {horiz, vert}); } } return thisAsCRS; } if (!dbContext) { return thisAsCRS; } const auto &l_domains = domains(); metadata::ExtentPtr extent; if (!l_domains.empty()) { extent = l_domains[0]->domainOfValidity(); } std::string crs_authority; const auto &l_identifiers = identifiers(); // If the object has an authority, restrict the transformations to // come from that codespace too. This avoids for example EPSG:4269 // (NAD83) to use a (dubious) ESRI transformation. if (!l_identifiers.empty()) { crs_authority = *(l_identifiers[0]->codeSpace()); } auto authorities = dbContext->getAllowedAuthorities( crs_authority, metadata::Identifier::EPSG); if (authorities.empty()) { authorities.emplace_back(); } // Vertical CRS ? auto vertCRS = dynamic_cast(this); if (vertCRS) { auto hubCRS = util::nn_static_pointer_cast(GeographicCRS::EPSG_4979); for (const auto &authority : authorities) { try { auto authFactory = io::AuthorityFactory::create( NN_NO_CHECK(dbContext), authority == "any" ? std::string() : authority); auto ctxt = operation::CoordinateOperationContext::create( authFactory, extent, 0.0); ctxt->setAllowUseIntermediateCRS(allowIntermediateCRSUse); // ctxt->setSpatialCriterion( // operation::CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); auto list = operation::CoordinateOperationFactory::create() ->createOperations(hubCRS, thisAsCRS, ctxt); CRSPtr candidateBoundCRS; for (const auto &op : list) { auto transf = util::nn_dynamic_pointer_cast< operation::Transformation>(op); // Only keep transformations that use a known grid if (transf && !transf->hasBallparkTransformation()) { auto gridsNeeded = transf->gridsNeeded(dbContext, true); bool gridsKnown = !gridsNeeded.empty(); for (const auto &gridDesc : gridsNeeded) { if (gridDesc.packageName.empty() && !(!gridDesc.url.empty() && gridDesc.openLicense) && !gridDesc.available) { gridsKnown = false; break; } } if (gridsKnown) { if (candidateBoundCRS) { candidateBoundCRS = nullptr; break; } candidateBoundCRS = BoundCRS::create(thisAsCRS, hubCRS, NN_NO_CHECK(transf)) .as_nullable(); } } } if (candidateBoundCRS) { return NN_NO_CHECK(candidateBoundCRS); } } catch (const std::exception &) { } } return thisAsCRS; } // Geodetic/geographic CRS ? auto geodCRS = util::nn_dynamic_pointer_cast(thisAsCRS); auto geogCRS = extractGeographicCRS(); auto hubCRS = util::nn_static_pointer_cast(GeographicCRS::EPSG_4326); if (geodCRS && !geogCRS) { if (geodCRS->_isEquivalentTo(GeographicCRS::EPSG_4978.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { return thisAsCRS; } hubCRS = util::nn_static_pointer_cast(GeodeticCRS::EPSG_4978); } else if (!geogCRS || geogCRS->_isEquivalentTo( GeographicCRS::EPSG_4326.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { return thisAsCRS; } else { geodCRS = geogCRS; } for (const auto &authority : authorities) { try { auto authFactory = io::AuthorityFactory::create( NN_NO_CHECK(dbContext), authority == "any" ? std::string() : authority); metadata::ExtentPtr extentResolved(extent); if (!extent) { getResolvedCRS(thisAsCRS, authFactory, extentResolved); } auto ctxt = operation::CoordinateOperationContext::create( authFactory, extentResolved, 0.0); ctxt->setAllowUseIntermediateCRS(allowIntermediateCRSUse); // ctxt->setSpatialCriterion( // operation::CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); auto list = operation::CoordinateOperationFactory::create() ->createOperations(NN_NO_CHECK(geodCRS), hubCRS, ctxt); CRSPtr candidateBoundCRS; for (const auto &op : list) { auto transf = util::nn_dynamic_pointer_cast( op); if (transf && !starts_with(transf->nameStr(), "Ballpark geo")) { try { transf->getTOWGS84Parameters(); } catch (const std::exception &) { continue; } if (candidateBoundCRS) { candidateBoundCRS = nullptr; break; } candidateBoundCRS = BoundCRS::create(thisAsCRS, hubCRS, NN_NO_CHECK(transf)) .as_nullable(); } else { auto concatenated = dynamic_cast( op.get()); if (concatenated) { // Case for EPSG:4807 / "NTF (Paris)" that is made of a // longitude rotation followed by a Helmert // The prime meridian shift will be accounted elsewhere const auto &subops = concatenated->operations(); if (subops.size() == 2) { auto firstOpIsTransformation = dynamic_cast( subops[0].get()); auto firstOpIsConversion = dynamic_cast( subops[0].get()); if ((firstOpIsTransformation && firstOpIsTransformation ->isLongitudeRotation()) || (dynamic_cast(thisAsCRS.get()) && firstOpIsConversion)) { transf = util::nn_dynamic_pointer_cast< operation::Transformation>(subops[1]); if (transf && !starts_with(transf->nameStr(), "Ballpark geo")) { try { transf->getTOWGS84Parameters(); } catch (const std::exception &) { continue; } if (candidateBoundCRS) { candidateBoundCRS = nullptr; break; } candidateBoundCRS = BoundCRS::create(thisAsCRS, hubCRS, NN_NO_CHECK(transf)) .as_nullable(); } } } } } } if (candidateBoundCRS) { return NN_NO_CHECK(candidateBoundCRS); } } catch (const std::exception &) { } } return thisAsCRS; } // --------------------------------------------------------------------------- /** \brief Returns a CRS whose coordinate system does not contain a vertical * component * * @return a CRS. */ CRSNNPtr CRS::stripVerticalComponent() const { auto self = NN_NO_CHECK( std::dynamic_pointer_cast(shared_from_this().as_nullable())); auto geogCRS = dynamic_cast(this); if (geogCRS) { const auto &axisList = geogCRS->coordinateSystem()->axisList(); if (axisList.size() == 3) { auto cs = cs::EllipsoidalCS::create(util::PropertyMap(), axisList[0], axisList[1]); return util::nn_static_pointer_cast(GeographicCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, nameStr()), geogCRS->datum(), geogCRS->datumEnsemble(), cs)); } } auto projCRS = dynamic_cast(this); if (projCRS) { const auto &axisList = projCRS->coordinateSystem()->axisList(); if (axisList.size() == 3) { auto cs = cs::CartesianCS::create(util::PropertyMap(), axisList[0], axisList[1]); return util::nn_static_pointer_cast(ProjectedCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, nameStr()), projCRS->baseCRS(), projCRS->derivingConversion(), cs)); } } return self; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress /** \brief Return a shallow clone of this object. */ CRSNNPtr CRS::shallowClone() const { return _shallowClone(); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRSNNPtr CRS::allowNonConformantWKT1Export() const { const auto boundCRS = dynamic_cast(this); if (boundCRS) { return BoundCRS::create( boundCRS->baseCRS()->allowNonConformantWKT1Export(), boundCRS->hubCRS(), boundCRS->transformation()); } auto crs(shallowClone()); crs->d->allowNonConformantWKT1Export_ = true; return crs; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRSNNPtr CRS::attachOriginalCompoundCRS(const CompoundCRSNNPtr &compoundCRS) const { const auto boundCRS = dynamic_cast(this); if (boundCRS) { return BoundCRS::create( boundCRS->baseCRS()->attachOriginalCompoundCRS(compoundCRS), boundCRS->hubCRS(), boundCRS->transformation()); } auto crs(shallowClone()); crs->d->originalCompoundCRS_ = compoundCRS.as_nullable(); return crs; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRSNNPtr CRS::alterName(const std::string &newName) const { auto crs = shallowClone(); auto newNameMod(newName); auto props = util::PropertyMap(); if (ends_with(newNameMod, " (deprecated)")) { newNameMod.resize(newNameMod.size() - strlen(" (deprecated)")); props.set(common::IdentifiedObject::DEPRECATED_KEY, true); } props.set(common::IdentifiedObject::NAME_KEY, newNameMod); crs->setProperties(props); return crs; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRSNNPtr CRS::alterId(const std::string &authName, const std::string &code) const { auto crs = shallowClone(); auto props = util::PropertyMap(); props.set(metadata::Identifier::CODESPACE_KEY, authName) .set(metadata::Identifier::CODE_KEY, code); crs->setProperties(props); return crs; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static bool mustAxisOrderBeSwitchedForVisualizationInternal( const std::vector &axisList) { const auto &dir0 = axisList[0]->direction(); const auto &dir1 = axisList[1]->direction(); if (&dir0 == &cs::AxisDirection::NORTH && &dir1 == &cs::AxisDirection::EAST) { return true; } // Address EPSG:32661 "WGS 84 / UPS North (N,E)" if (&dir0 == &cs::AxisDirection::SOUTH && &dir1 == &cs::AxisDirection::SOUTH) { const auto &meridian0 = axisList[0]->meridian(); const auto &meridian1 = axisList[1]->meridian(); return meridian0 != nullptr && meridian1 != nullptr && std::abs(meridian0->longitude().convertToUnit( common::UnitOfMeasure::DEGREE) - 180.0) < 1e-10 && std::abs(meridian1->longitude().convertToUnit( common::UnitOfMeasure::DEGREE) - 90.0) < 1e-10; } if (&dir0 == &cs::AxisDirection::NORTH && &dir1 == &cs::AxisDirection::NORTH) { const auto &meridian0 = axisList[0]->meridian(); const auto &meridian1 = axisList[1]->meridian(); return meridian0 != nullptr && meridian1 != nullptr && (( // Address EPSG:32761 "WGS 84 / UPS South (N,E)" std::abs(meridian0->longitude().convertToUnit( common::UnitOfMeasure::DEGREE) - 0.0) < 1e-10 && std::abs(meridian1->longitude().convertToUnit( common::UnitOfMeasure::DEGREE) - 90.0) < 1e-10) || // Address EPSG:5482 "RSRGD2000 / RSPS2000" (std::abs(meridian0->longitude().convertToUnit( common::UnitOfMeasure::DEGREE) - 180) < 1e-10 && std::abs(meridian1->longitude().convertToUnit( common::UnitOfMeasure::DEGREE) - -90.0) < 1e-10)); } return false; } // --------------------------------------------------------------------------- bool CRS::mustAxisOrderBeSwitchedForVisualization() const { const CompoundCRS *compoundCRS = dynamic_cast(this); if (compoundCRS) { const auto &comps = compoundCRS->componentReferenceSystems(); if (!comps.empty()) { return comps[0]->mustAxisOrderBeSwitchedForVisualization(); } } const GeographicCRS *geogCRS = dynamic_cast(this); if (geogCRS) { return mustAxisOrderBeSwitchedForVisualizationInternal( geogCRS->coordinateSystem()->axisList()); } const ProjectedCRS *projCRS = dynamic_cast(this); if (projCRS) { return mustAxisOrderBeSwitchedForVisualizationInternal( projCRS->coordinateSystem()->axisList()); } return false; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRSNNPtr CRS::normalizeForVisualization() const { const auto createProperties = [this](const std::string &newName = std::string()) { auto props = util::PropertyMap().set( common::IdentifiedObject::NAME_KEY, !newName.empty() ? newName : nameStr() + " (with axis order normalized for visualization)"); const auto &l_domains = domains(); if (!l_domains.empty()) { auto array = util::ArrayOfBaseObject::create(); for (const auto &domain : l_domains) { array->add(domain); } if (!array->empty()) { props.set(common::ObjectUsage::OBJECT_DOMAIN_KEY, array); } } const auto &l_identifiers = identifiers(); const auto &l_remarks = remarks(); if (l_identifiers.size() == 1) { std::string remarks("Axis order reversed compared to "); remarks += *(l_identifiers[0]->codeSpace()); remarks += ':'; remarks += l_identifiers[0]->code(); if (!l_remarks.empty()) { remarks += ". "; remarks += l_remarks; } props.set(common::IdentifiedObject::REMARKS_KEY, remarks); } else if (!l_remarks.empty()) { props.set(common::IdentifiedObject::REMARKS_KEY, l_remarks); } return props; }; const CompoundCRS *compoundCRS = dynamic_cast(this); if (compoundCRS) { const auto &comps = compoundCRS->componentReferenceSystems(); if (!comps.empty() && comps[0]->mustAxisOrderBeSwitchedForVisualization()) { std::vector newComps; newComps.emplace_back(comps[0]->normalizeForVisualization()); std::string l_name = newComps.back()->nameStr(); for (size_t i = 1; i < comps.size(); i++) { newComps.emplace_back(comps[i]); l_name += " + "; l_name += newComps.back()->nameStr(); } return util::nn_static_pointer_cast( CompoundCRS::create(createProperties(l_name), newComps)); } } const GeographicCRS *geogCRS = dynamic_cast(this); if (geogCRS) { const auto &axisList = geogCRS->coordinateSystem()->axisList(); if (mustAxisOrderBeSwitchedForVisualizationInternal(axisList)) { auto cs = axisList.size() == 2 ? cs::EllipsoidalCS::create(util::PropertyMap(), axisList[1], axisList[0]) : cs::EllipsoidalCS::create(util::PropertyMap(), axisList[1], axisList[0], axisList[2]); return util::nn_static_pointer_cast( GeographicCRS::create(createProperties(), geogCRS->datum(), geogCRS->datumEnsemble(), cs)); } } const ProjectedCRS *projCRS = dynamic_cast(this); if (projCRS) { const auto &axisList = projCRS->coordinateSystem()->axisList(); if (mustAxisOrderBeSwitchedForVisualizationInternal(axisList)) { auto cs = axisList.size() == 2 ? cs::CartesianCS::create(util::PropertyMap(), axisList[1], axisList[0]) : cs::CartesianCS::create(util::PropertyMap(), axisList[1], axisList[0], axisList[2]); return util::nn_static_pointer_cast( ProjectedCRS::create(createProperties(), projCRS->baseCRS(), projCRS->derivingConversion(), cs)); } } return NN_NO_CHECK( std::static_pointer_cast(shared_from_this().as_nullable())); } //! @endcond // --------------------------------------------------------------------------- /** \brief Identify the CRS with reference CRSs. * * The candidate CRSs are either hard-coded, or looked in the database when * authorityFactory is not null. * * Note that the implementation uses a set of heuristics to have a good * compromise of successful identifications over execution time. It might miss * legitimate matches in some circumstances. * * The method returns a list of matching reference CRS, and the percentage * (0-100) of confidence in the match. The list is sorted by decreasing * confidence. *
    *
  • 100% means that the name of the reference entry * perfectly matches the CRS name, and both are equivalent. In which case a * single result is returned. * Note: in the case of a GeographicCRS whose axis * order is implicit in the input definition (for example ESRI WKT), then axis * order is ignored for the purpose of identification. That is the CRS built * from * GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]], * PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]] * will be identified to EPSG:4326, but will not pass a * isEquivalentTo(EPSG_4326, util::IComparable::Criterion::EQUIVALENT) test, * but rather isEquivalentTo(EPSG_4326, * util::IComparable::Criterion::EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS) *
    • *
  • 90% means that CRS are equivalent, but the names are not exactly the * same.
    • *
  • 70% means that CRS are equivalent), but the names do not match at * all.
    • *
  • 25% means that the CRS are not equivalent, but there is some similarity * in * the names.
    • *
* Other confidence values may be returned by some specialized implementations. * * This is implemented for GeodeticCRS, ProjectedCRS, VerticalCRS and * CompoundCRS. * * @param authorityFactory Authority factory (or null, but degraded * functionality) * @return a list of matching reference CRS, and the percentage (0-100) of * confidence in the match. */ std::list> CRS::identify(const io::AuthorityFactoryPtr &authorityFactory) const { return _identify(authorityFactory); } // --------------------------------------------------------------------------- /** \brief Return CRSs that are non-deprecated substitutes for the current CRS. */ std::list CRS::getNonDeprecated(const io::DatabaseContextNNPtr &dbContext) const { std::list res; const auto &l_identifiers = identifiers(); if (l_identifiers.empty()) { return res; } const char *tableName = nullptr; if (dynamic_cast(this)) { tableName = "geodetic_crs"; } else if (dynamic_cast(this)) { tableName = "projected_crs"; } else if (dynamic_cast(this)) { tableName = "vertical_crs"; } else if (dynamic_cast(this)) { tableName = "compound_crs"; } if (!tableName) { return res; } const auto &id = l_identifiers[0]; auto tmpRes = dbContext->getNonDeprecated(tableName, *(id->codeSpace()), id->code()); for (const auto &pair : tmpRes) { res.emplace_back(io::AuthorityFactory::create(dbContext, pair.first) ->createCoordinateReferenceSystem(pair.second)); } return res; } // --------------------------------------------------------------------------- /** \brief Return a variant of this CRS "promoted" to a 3D one, if not already * the case. * * The new axis will be ellipsoidal height, oriented upwards, and with metre * units. * * @param newName Name of the new CRS. If empty, nameStr() will be used. * @param dbContext Database context to look for potentially already registered * 3D CRS. May be nullptr. * @return a new CRS promoted to 3D, or the current one if already 3D or not * applicable. * @since 6.3 */ CRSNNPtr CRS::promoteTo3D(const std::string &newName, const io::DatabaseContextPtr &dbContext) const { auto upAxis = cs::CoordinateSystemAxis::create( util::PropertyMap().set(IdentifiedObject::NAME_KEY, cs::AxisName::Ellipsoidal_height), cs::AxisAbbreviation::h, cs::AxisDirection::UP, common::UnitOfMeasure::METRE); return promoteTo3D(newName, dbContext, upAxis); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CRSNNPtr CRS::promoteTo3D(const std::string &newName, const io::DatabaseContextPtr &dbContext, const cs::CoordinateSystemAxisNNPtr &verticalAxisIfNotAlreadyPresent) const { const auto createProperties = [this, &newName]() { auto props = util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, !newName.empty() ? newName : nameStr()); const auto &l_domains = domains(); if (!l_domains.empty()) { auto array = util::ArrayOfBaseObject::create(); for (const auto &domain : l_domains) { auto extent = domain->domainOfValidity(); if (extent) { // Propagate only the extent, not the scope, as it might // imply more that we can guarantee with the promotion to // 3D. auto newDomain = common::ObjectDomain::create( util::optional(), extent); array->add(newDomain); } } if (!array->empty()) { props.set(common::ObjectUsage::OBJECT_DOMAIN_KEY, array); } } const auto &l_identifiers = identifiers(); const auto &l_remarks = remarks(); if (l_identifiers.size() == 1) { std::string remarks("Promoted to 3D from "); remarks += *(l_identifiers[0]->codeSpace()); remarks += ':'; remarks += l_identifiers[0]->code(); if (!l_remarks.empty()) { remarks += ". "; remarks += l_remarks; } props.set(common::IdentifiedObject::REMARKS_KEY, remarks); } else if (!l_remarks.empty()) { props.set(common::IdentifiedObject::REMARKS_KEY, l_remarks); } return props; }; const auto geogCRS = dynamic_cast(this); if (geogCRS) { const auto &axisList = geogCRS->coordinateSystem()->axisList(); if (axisList.size() == 2) { const auto &l_identifiers = identifiers(); // First check if there is a Geographic 3D CRS in the database // of the same name. // This is the common practice in the EPSG dataset. if (dbContext && l_identifiers.size() == 1) { auto authFactory = io::AuthorityFactory::create( NN_NO_CHECK(dbContext), *(l_identifiers[0]->codeSpace())); auto res = authFactory->createObjectsFromName( nameStr(), {io::AuthorityFactory::ObjectType::GEOGRAPHIC_3D_CRS}, false); if (!res.empty()) { const auto &firstRes = res.front(); const auto firstResGeog = dynamic_cast(firstRes.get()); const auto &firstResAxisList = firstResGeog->coordinateSystem()->axisList(); if (firstResAxisList[2]->_isEquivalentTo( verticalAxisIfNotAlreadyPresent.get(), util::IComparable::Criterion::EQUIVALENT) && geogCRS->is2DPartOf3D(NN_NO_CHECK(firstResGeog), dbContext)) { return NN_NO_CHECK( util::nn_dynamic_pointer_cast(firstRes)); } } } auto cs = cs::EllipsoidalCS::create( util::PropertyMap(), axisList[0], axisList[1], verticalAxisIfNotAlreadyPresent); return util::nn_static_pointer_cast( GeographicCRS::create(createProperties(), geogCRS->datum(), geogCRS->datumEnsemble(), cs)); } } const auto projCRS = dynamic_cast(this); if (projCRS) { const auto &axisList = projCRS->coordinateSystem()->axisList(); if (axisList.size() == 2) { auto base3DCRS = projCRS->baseCRS()->promoteTo3D(std::string(), dbContext); auto cs = cs::CartesianCS::create(util::PropertyMap(), axisList[0], axisList[1], verticalAxisIfNotAlreadyPresent); return util::nn_static_pointer_cast(ProjectedCRS::create( createProperties(), NN_NO_CHECK( util::nn_dynamic_pointer_cast(base3DCRS)), projCRS->derivingConversion(), cs)); } } const auto boundCRS = dynamic_cast(this); if (boundCRS) { auto base3DCRS = boundCRS->baseCRS()->promoteTo3D( newName, dbContext, verticalAxisIfNotAlreadyPresent); auto transf = boundCRS->transformation(); try { transf->getTOWGS84Parameters(); return BoundCRS::create( base3DCRS, boundCRS->hubCRS()->promoteTo3D(std::string(), dbContext), transf->promoteTo3D(std::string(), dbContext)); } catch (const io::FormattingException &) { return BoundCRS::create(base3DCRS, boundCRS->hubCRS(), transf); } } return NN_NO_CHECK( std::static_pointer_cast(shared_from_this().as_nullable())); } //! @endcond // --------------------------------------------------------------------------- /** \brief Return a variant of this CRS "demoted" to a 2D one, if not already * the case. * * * @param newName Name of the new CRS. If empty, nameStr() will be used. * @param dbContext Database context to look for potentially already registered * 2D CRS. May be nullptr. * @return a new CRS demoted to 2D, or the current one if already 2D or not * applicable. * @since 6.3 */ CRSNNPtr CRS::demoteTo2D(const std::string &newName, const io::DatabaseContextPtr &dbContext) const { const auto geogCRS = dynamic_cast(this); if (geogCRS) { return geogCRS->demoteTo2D(newName, dbContext); } const auto projCRS = dynamic_cast(this); if (projCRS) { return projCRS->demoteTo2D(newName, dbContext); } const auto boundCRS = dynamic_cast(this); if (boundCRS) { auto base2DCRS = boundCRS->baseCRS()->demoteTo2D(newName, dbContext); auto transf = boundCRS->transformation(); try { transf->getTOWGS84Parameters(); return BoundCRS::create( base2DCRS, boundCRS->hubCRS()->demoteTo2D(std::string(), dbContext), transf->demoteTo2D(std::string(), dbContext)); } catch (const io::FormattingException &) { return BoundCRS::create(base2DCRS, boundCRS->hubCRS(), transf); } } const auto compoundCRS = dynamic_cast(this); if (compoundCRS) { const auto &components = compoundCRS->componentReferenceSystems(); if (components.size() >= 2) { return components[0]; } } return NN_NO_CHECK( std::static_pointer_cast(shared_from_this().as_nullable())); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> CRS::_identify(const io::AuthorityFactoryPtr &) const { return {}; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct SingleCRS::Private { datum::DatumPtr datum{}; datum::DatumEnsemblePtr datumEnsemble{}; cs::CoordinateSystemNNPtr coordinateSystem; Private(const datum::DatumPtr &datumIn, const datum::DatumEnsemblePtr &datumEnsembleIn, const cs::CoordinateSystemNNPtr &csIn) : datum(datumIn), datumEnsemble(datumEnsembleIn), coordinateSystem(csIn) { if ((datum ? 1 : 0) + (datumEnsemble ? 1 : 0) != 1) { throw util::Exception("datum or datumEnsemble should be set"); } } }; //! @endcond // --------------------------------------------------------------------------- SingleCRS::SingleCRS(const datum::DatumPtr &datumIn, const datum::DatumEnsemblePtr &datumEnsembleIn, const cs::CoordinateSystemNNPtr &csIn) : d(internal::make_unique(datumIn, datumEnsembleIn, csIn)) {} // --------------------------------------------------------------------------- SingleCRS::SingleCRS(const SingleCRS &other) : CRS(other), d(internal::make_unique(*other.d)) {} // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress SingleCRS::~SingleCRS() = default; //! @endcond // --------------------------------------------------------------------------- /** \brief Return the datum::Datum associated with the CRS. * * This might be null, in which case datumEnsemble() return will not be null. * * @return a Datum that might be null. */ const datum::DatumPtr &SingleCRS::datum() PROJ_PURE_DEFN { return d->datum; } // --------------------------------------------------------------------------- /** \brief Return the datum::DatumEnsemble associated with the CRS. * * This might be null, in which case datum() return will not be null. * * @return a DatumEnsemble that might be null. */ const datum::DatumEnsemblePtr &SingleCRS::datumEnsemble() PROJ_PURE_DEFN { return d->datumEnsemble; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress /** \brief Return the real datum or a synthetized one if a datumEnsemble. */ const datum::DatumNNPtr SingleCRS::datumNonNull(const io::DatabaseContextPtr &dbContext) const { return d->datum ? NN_NO_CHECK(d->datum) : d->datumEnsemble->asDatum(dbContext); } //! @endcond // --------------------------------------------------------------------------- /** \brief Return the cs::CoordinateSystem associated with the CRS. * * @return a CoordinateSystem. */ const cs::CoordinateSystemNNPtr &SingleCRS::coordinateSystem() PROJ_PURE_DEFN { return d->coordinateSystem; } // --------------------------------------------------------------------------- bool SingleCRS::baseIsEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherSingleCRS = dynamic_cast(other); if (otherSingleCRS == nullptr || (criterion == util::IComparable::Criterion::STRICT && !ObjectUsage::_isEquivalentTo(other, criterion, dbContext))) { return false; } if (criterion == util::IComparable::Criterion::STRICT) { const auto &thisDatum = d->datum; const auto &otherDatum = otherSingleCRS->d->datum; if (thisDatum) { if (!thisDatum->_isEquivalentTo(otherDatum.get(), criterion, dbContext)) { return false; } } else { if (otherDatum) { return false; } } const auto &thisDatumEnsemble = d->datumEnsemble; const auto &otherDatumEnsemble = otherSingleCRS->d->datumEnsemble; if (thisDatumEnsemble) { if (!thisDatumEnsemble->_isEquivalentTo(otherDatumEnsemble.get(), criterion, dbContext)) { return false; } } else { if (otherDatumEnsemble) { return false; } } } else { if (!datumNonNull(dbContext)->_isEquivalentTo( otherSingleCRS->datumNonNull(dbContext).get(), criterion, dbContext)) { return false; } } return d->coordinateSystem->_isEquivalentTo( otherSingleCRS->d->coordinateSystem.get(), criterion, dbContext) && getExtensionProj4() == otherSingleCRS->getExtensionProj4(); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void SingleCRS::exportDatumOrDatumEnsembleToWkt( io::WKTFormatter *formatter) const // throw(io::FormattingException) { const auto &l_datum = d->datum; if (l_datum) { l_datum->_exportToWKT(formatter); } else { const auto &l_datumEnsemble = d->datumEnsemble; assert(l_datumEnsemble); l_datumEnsemble->_exportToWKT(formatter); } } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct GeodeticCRS::Private { std::vector velocityModel{}; datum::GeodeticReferenceFramePtr datum_; explicit Private(const datum::GeodeticReferenceFramePtr &datumIn) : datum_(datumIn) {} }; // --------------------------------------------------------------------------- static const datum::DatumEnsemblePtr & checkEnsembleForGeodeticCRS(const datum::GeodeticReferenceFramePtr &datumIn, const datum::DatumEnsemblePtr &ensemble) { const char *msg = "One of Datum or DatumEnsemble should be defined"; if (datumIn) { if (!ensemble) { return ensemble; } msg = "Datum and DatumEnsemble should not be defined"; } else if (ensemble) { const auto &datums = ensemble->datums(); assert(!datums.empty()); auto grfFirst = dynamic_cast(datums[0].get()); if (grfFirst) { return ensemble; } msg = "Ensemble should contain GeodeticReferenceFrame"; } throw util::Exception(msg); } //! @endcond // --------------------------------------------------------------------------- GeodeticCRS::GeodeticCRS(const datum::GeodeticReferenceFramePtr &datumIn, const datum::DatumEnsemblePtr &datumEnsembleIn, const cs::EllipsoidalCSNNPtr &csIn) : SingleCRS(datumIn, checkEnsembleForGeodeticCRS(datumIn, datumEnsembleIn), csIn), d(internal::make_unique(datumIn)) {} // --------------------------------------------------------------------------- GeodeticCRS::GeodeticCRS(const datum::GeodeticReferenceFramePtr &datumIn, const datum::DatumEnsemblePtr &datumEnsembleIn, const cs::SphericalCSNNPtr &csIn) : SingleCRS(datumIn, checkEnsembleForGeodeticCRS(datumIn, datumEnsembleIn), csIn), d(internal::make_unique(datumIn)) {} // --------------------------------------------------------------------------- GeodeticCRS::GeodeticCRS(const datum::GeodeticReferenceFramePtr &datumIn, const datum::DatumEnsemblePtr &datumEnsembleIn, const cs::CartesianCSNNPtr &csIn) : SingleCRS(datumIn, checkEnsembleForGeodeticCRS(datumIn, datumEnsembleIn), csIn), d(internal::make_unique(datumIn)) {} // --------------------------------------------------------------------------- GeodeticCRS::GeodeticCRS(const GeodeticCRS &other) : SingleCRS(other), d(internal::make_unique(*other.d)) {} // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress GeodeticCRS::~GeodeticCRS() = default; //! @endcond // --------------------------------------------------------------------------- CRSNNPtr GeodeticCRS::_shallowClone() const { auto crs(GeodeticCRS::nn_make_shared(*this)); crs->assignSelf(crs); return crs; } // --------------------------------------------------------------------------- /** \brief Return the datum::GeodeticReferenceFrame associated with the CRS. * * @return a GeodeticReferenceFrame or null (in which case datumEnsemble() * should return a non-null pointer.) */ const datum::GeodeticReferenceFramePtr &GeodeticCRS::datum() PROJ_PURE_DEFN { return d->datum_; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress /** \brief Return the real datum or a synthetized one if a datumEnsemble. */ const datum::GeodeticReferenceFrameNNPtr GeodeticCRS::datumNonNull(const io::DatabaseContextPtr &dbContext) const { return NN_NO_CHECK( d->datum_ ? d->datum_ : util::nn_dynamic_pointer_cast( SingleCRS::getPrivate()->datumEnsemble->asDatum(dbContext))); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static datum::GeodeticReferenceFrame *oneDatum(const GeodeticCRS *crs) { const auto &l_datumEnsemble = crs->datumEnsemble(); assert(l_datumEnsemble); const auto &l_datums = l_datumEnsemble->datums(); return static_cast(l_datums[0].get()); } //! @endcond // --------------------------------------------------------------------------- /** \brief Return the PrimeMeridian associated with the GeodeticReferenceFrame * or with one of the GeodeticReferenceFrame of the datumEnsemble(). * * @return the PrimeMeridian. */ const datum::PrimeMeridianNNPtr &GeodeticCRS::primeMeridian() PROJ_PURE_DEFN { if (d->datum_) { return d->datum_->primeMeridian(); } return oneDatum(this)->primeMeridian(); } // --------------------------------------------------------------------------- /** \brief Return the ellipsoid associated with the GeodeticReferenceFrame * or with one of the GeodeticReferenceFrame of the datumEnsemble(). * * @return the PrimeMeridian. */ const datum::EllipsoidNNPtr &GeodeticCRS::ellipsoid() PROJ_PURE_DEFN { if (d->datum_) { return d->datum_->ellipsoid(); } return oneDatum(this)->ellipsoid(); } // --------------------------------------------------------------------------- /** \brief Return the velocity model associated with the CRS. * * @return a velocity model. might be null. */ const std::vector & GeodeticCRS::velocityModel() PROJ_PURE_DEFN { return d->velocityModel; } // --------------------------------------------------------------------------- /** \brief Return whether the CRS is a geocentric one. * * A geocentric CRS is a geodetic CRS that has a Cartesian coordinate system * with three axis, whose direction is respectively * cs::AxisDirection::GEOCENTRIC_X, * cs::AxisDirection::GEOCENTRIC_Y and cs::AxisDirection::GEOCENTRIC_Z. * * @return true if the CRS is a geocentric CRS. */ bool GeodeticCRS::isGeocentric() PROJ_PURE_DEFN { const auto &cs = coordinateSystem(); const auto &axisList = cs->axisList(); return axisList.size() == 3 && dynamic_cast(cs.get()) != nullptr && &axisList[0]->direction() == &cs::AxisDirection::GEOCENTRIC_X && &axisList[1]->direction() == &cs::AxisDirection::GEOCENTRIC_Y && &axisList[2]->direction() == &cs::AxisDirection::GEOCENTRIC_Z; } // --------------------------------------------------------------------------- /** \brief Instantiate a GeodeticCRS from a datum::GeodeticReferenceFrame and a * cs::SphericalCS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datum The datum of the CRS. * @param cs a SphericalCS. * @return new GeodeticCRS. */ GeodeticCRSNNPtr GeodeticCRS::create(const util::PropertyMap &properties, const datum::GeodeticReferenceFrameNNPtr &datum, const cs::SphericalCSNNPtr &cs) { return create(properties, datum.as_nullable(), nullptr, cs); } // --------------------------------------------------------------------------- /** \brief Instantiate a GeodeticCRS from a datum::GeodeticReferenceFrame or * datum::DatumEnsemble and a cs::SphericalCS. * * One and only one of datum or datumEnsemble should be set to a non-null value. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datum The datum of the CRS, or nullptr * @param datumEnsemble The datum ensemble of the CRS, or nullptr. * @param cs a SphericalCS. * @return new GeodeticCRS. */ GeodeticCRSNNPtr GeodeticCRS::create(const util::PropertyMap &properties, const datum::GeodeticReferenceFramePtr &datum, const datum::DatumEnsemblePtr &datumEnsemble, const cs::SphericalCSNNPtr &cs) { auto crs( GeodeticCRS::nn_make_shared(datum, datumEnsemble, cs)); crs->assignSelf(crs); crs->setProperties(properties); properties.getStringValue("EXTENSION_PROJ4", crs->CRS::getPrivate()->extensionProj4_); return crs; } // --------------------------------------------------------------------------- /** \brief Instantiate a GeodeticCRS from a datum::GeodeticReferenceFrame and a * cs::CartesianCS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datum The datum of the CRS. * @param cs a CartesianCS. * @return new GeodeticCRS. */ GeodeticCRSNNPtr GeodeticCRS::create(const util::PropertyMap &properties, const datum::GeodeticReferenceFrameNNPtr &datum, const cs::CartesianCSNNPtr &cs) { return create(properties, datum.as_nullable(), nullptr, cs); } // --------------------------------------------------------------------------- /** \brief Instantiate a GeodeticCRS from a datum::GeodeticReferenceFrame or * datum::DatumEnsemble and a cs::CartesianCS. * * One and only one of datum or datumEnsemble should be set to a non-null value. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datum The datum of the CRS, or nullptr * @param datumEnsemble The datum ensemble of the CRS, or nullptr. * @param cs a CartesianCS * @return new GeodeticCRS. */ GeodeticCRSNNPtr GeodeticCRS::create(const util::PropertyMap &properties, const datum::GeodeticReferenceFramePtr &datum, const datum::DatumEnsemblePtr &datumEnsemble, const cs::CartesianCSNNPtr &cs) { auto crs( GeodeticCRS::nn_make_shared(datum, datumEnsemble, cs)); crs->assignSelf(crs); crs->setProperties(properties); properties.getStringValue("EXTENSION_PROJ4", crs->CRS::getPrivate()->extensionProj4_); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress // Try to format a Geographic/ProjectedCRS 3D CRS as a // GEOGCS[]/PROJCS[],VERTCS[...,DATUM[],...] if we find corresponding objects static bool exportAsESRIWktCompoundCRSWithEllipsoidalHeight( const CRS *self, const GeodeticCRS *geodCRS, io::WKTFormatter *formatter) { const auto &dbContext = formatter->databaseContext(); if (!dbContext) { return false; } const auto l_datum = geodCRS->datumNonNull(formatter->databaseContext()); auto l_alias = dbContext->getAliasFromOfficialName( l_datum->nameStr(), "geodetic_datum", "ESRI"); if (l_alias.empty()) { return false; } auto authFactory = io::AuthorityFactory::create(NN_NO_CHECK(dbContext), std::string()); auto list = authFactory->createObjectsFromName( l_alias, {io::AuthorityFactory::ObjectType::GEODETIC_REFERENCE_FRAME}, false /* approximate=false*/); if (list.empty()) { return false; } auto gdatum = util::nn_dynamic_pointer_cast(list.front()); if (gdatum == nullptr || gdatum->identifiers().empty()) { return false; } const auto &gdatum_ids = gdatum->identifiers(); auto vertCRSList = authFactory->createVerticalCRSFromDatum( "ESRI", "from_geogdatum_" + *gdatum_ids[0]->codeSpace() + '_' + gdatum_ids[0]->code()); if (vertCRSList.size() != 1) { return false; } self->demoteTo2D(std::string(), dbContext)->_exportToWKT(formatter); vertCRSList.front()->_exportToWKT(formatter); return true; } // --------------------------------------------------------------------------- // Try to format a Geographic/ProjectedCRS 3D CRS as a // GEOGCS[]/PROJCS[],VERTCS["Ellipsoid (metre)",DATUM["Ellipsoid",2002],...] static bool exportAsWKT1CompoundCRSWithEllipsoidalHeight( const CRSNNPtr &base2DCRS, const cs::CoordinateSystemAxisNNPtr &verticalAxis, io::WKTFormatter *formatter) { std::string verticalCRSName = "Ellipsoid ("; verticalCRSName += verticalAxis->unit().name(); verticalCRSName += ')'; auto vertDatum = datum::VerticalReferenceFrame::create( util::PropertyMap() .set(common::IdentifiedObject::NAME_KEY, "Ellipsoid") .set("VERT_DATUM_TYPE", "2002")); auto vertCRS = VerticalCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, verticalCRSName), vertDatum.as_nullable(), nullptr, cs::VerticalCS::create(util::PropertyMap(), verticalAxis)); formatter->startNode(io::WKTConstants::COMPD_CS, false); formatter->addQuotedString(base2DCRS->nameStr() + " + " + verticalCRSName); base2DCRS->_exportToWKT(formatter); vertCRS->_exportToWKT(formatter); formatter->endNode(); return true; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void GeodeticCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; const bool isGeographic = dynamic_cast(this) != nullptr; const auto &cs = coordinateSystem(); const auto &axisList = cs->axisList(); const auto oldAxisOutputRule = formatter->outputAxis(); auto l_name = nameStr(); const auto &dbContext = formatter->databaseContext(); if (!isWKT2 && formatter->useESRIDialect() && axisList.size() == 3) { if (!isGeographic) { io::FormattingException::Throw( "Geocentric CRS not supported in WKT1_ESRI"); } // Try to format the Geographic 3D CRS as a GEOGCS[],VERTCS[...,DATUM[]] // if we find corresponding objects if (dbContext) { if (exportAsESRIWktCompoundCRSWithEllipsoidalHeight(this, this, formatter)) { return; } } io::FormattingException::Throw( "Cannot export this Geographic 3D CRS in WKT1_ESRI"); } if (!isWKT2 && formatter->isStrict() && isGeographic && axisList.size() != 2 && oldAxisOutputRule != io::WKTFormatter::OutputAxisRule::NO) { auto geogCRS2D = demoteTo2D(std::string(), dbContext); if (dbContext) { const auto res = geogCRS2D->identify(io::AuthorityFactory::create( NN_NO_CHECK(dbContext), metadata::Identifier::EPSG)); if (res.size() == 1) { const auto &front = res.front(); if (front.second == 100) { geogCRS2D = front.first; } } } if (CRS::getPrivate()->allowNonConformantWKT1Export_) { formatter->startNode(io::WKTConstants::COMPD_CS, false); formatter->addQuotedString(l_name + " + " + l_name); geogCRS2D->_exportToWKT(formatter); const auto oldTOWGSParameters = formatter->getTOWGS84Parameters(); formatter->setTOWGS84Parameters({}); geogCRS2D->_exportToWKT(formatter); formatter->setTOWGS84Parameters(oldTOWGSParameters); formatter->endNode(); return; } auto &originalCompoundCRS = CRS::getPrivate()->originalCompoundCRS_; if (originalCompoundCRS) { originalCompoundCRS->_exportToWKT(formatter); return; } if (formatter->isAllowedEllipsoidalHeightAsVerticalCRS()) { if (exportAsWKT1CompoundCRSWithEllipsoidalHeight( geogCRS2D, axisList[2], formatter)) { 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; if (dbContext) { auto l_alias = dbContext->getAliasFromOfficialName( l_name, "geodetic_crs", "ESRI"); if (!l_alias.empty()) { l_name = l_alias; aliasFound = true; } } if (!aliasFound) { l_name = io::WKTFormatter::morphNameToESRI(l_name); if (!starts_with(l_name, "GCS_")) { l_name = "GCS_" + l_name; } } } } if (!isWKT2 && !formatter->useESRIDialect() && isDeprecated()) { l_name += " (deprecated)"; } formatter->addQuotedString(l_name); const auto &unit = axisList[0]->unit(); formatter->pushAxisAngularUnit(common::UnitOfMeasure::create(unit)); exportDatumOrDatumEnsembleToWkt(formatter); primeMeridian()->_exportToWKT(formatter); formatter->popAxisAngularUnit(); if (!isWKT2) { unit._exportToWKT(formatter); } if (oldAxisOutputRule == io::WKTFormatter::OutputAxisRule::WKT1_GDAL_EPSG_STYLE && isGeocentric()) { formatter->setOutputAxis(io::WKTFormatter::OutputAxisRule::YES); } cs->_exportToWKT(formatter); formatter->setOutputAxis(oldAxisOutputRule); ObjectUsage::baseExportToWKT(formatter); if (!isWKT2 && !formatter->useESRIDialect()) { const auto &extensionProj4 = CRS::getPrivate()->extensionProj4_; if (!extensionProj4.empty()) { formatter->startNode(io::WKTConstants::EXTENSION, false); formatter->addQuotedString("PROJ4"); formatter->addQuotedString(extensionProj4); formatter->endNode(); } } formatter->endNode(); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void GeodeticCRS::addGeocentricUnitConversionIntoPROJString( io::PROJStringFormatter *formatter) const { const auto &axisList = coordinateSystem()->axisList(); const auto &unit = axisList[0]->unit(); if (!unit._isEquivalentTo(common::UnitOfMeasure::METRE, util::IComparable::Criterion::EQUIVALENT)) { if (formatter->getCRSExport()) { io::FormattingException::Throw( "GeodeticCRS::exportToPROJString() only " "supports metre unit"); } formatter->addStep("unitconvert"); formatter->addParam("xy_in", "m"); formatter->addParam("z_in", "m"); { auto projUnit = unit.exportToPROJString(); if (!projUnit.empty()) { formatter->addParam("xy_out", projUnit); formatter->addParam("z_out", projUnit); return; } } const auto &toSI = unit.conversionToSI(); formatter->addParam("xy_out", toSI); formatter->addParam("z_out", toSI); } else if (formatter->getCRSExport()) { formatter->addParam("units", "m"); } } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void GeodeticCRS::_exportToPROJString( io::PROJStringFormatter *formatter) const // throw(io::FormattingException) { const auto &extensionProj4 = CRS::getPrivate()->extensionProj4_; if (!extensionProj4.empty()) { formatter->ingestPROJString( replaceAll(extensionProj4, " +type=crs", "")); formatter->addNoDefs(false); return; } if (!isGeocentric()) { io::FormattingException::Throw( "GeodeticCRS::exportToPROJString() only " "supports geocentric coordinate systems"); } if (!formatter->getCRSExport()) { formatter->addStep("cart"); } else { formatter->addStep("geocent"); } addDatumInfoToPROJString(formatter); addGeocentricUnitConversionIntoPROJString(formatter); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void GeodeticCRS::addDatumInfoToPROJString( io::PROJStringFormatter *formatter) const // throw(io::FormattingException) { const auto &TOWGS84Params = formatter->getTOWGS84Parameters(); bool datumWritten = false; const auto &nadgrids = formatter->getHDatumExtension(); const auto l_datum = datumNonNull(formatter->databaseContext()); if (formatter->getCRSExport() && TOWGS84Params.empty() && nadgrids.empty()) { if (l_datum->_isEquivalentTo( datum::GeodeticReferenceFrame::EPSG_6326.get(), util::IComparable::Criterion::EQUIVALENT)) { datumWritten = true; formatter->addParam("datum", "WGS84"); } else if (l_datum->_isEquivalentTo( datum::GeodeticReferenceFrame::EPSG_6267.get(), util::IComparable::Criterion::EQUIVALENT)) { datumWritten = true; formatter->addParam("datum", "NAD27"); } else if (l_datum->_isEquivalentTo( datum::GeodeticReferenceFrame::EPSG_6269.get(), util::IComparable::Criterion::EQUIVALENT)) { datumWritten = true; if (formatter->getLegacyCRSToCRSContext()) { // We do not want datum=NAD83 to cause a useless towgs84=0,0,0 formatter->addParam("ellps", "GRS80"); } else { formatter->addParam("datum", "NAD83"); } } } if (!datumWritten) { ellipsoid()->_exportToPROJString(formatter); primeMeridian()->_exportToPROJString(formatter); } if (TOWGS84Params.size() == 7) { formatter->addParam("towgs84", TOWGS84Params); } if (!nadgrids.empty()) { formatter->addParam("nadgrids", nadgrids); } } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void GeodeticCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("GeodeticCRS", !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } const auto &l_datum(datum()); if (l_datum) { writer->AddObjKey("datum"); l_datum->_exportToJSON(formatter); } else { writer->AddObjKey("datum_ensemble"); formatter->setOmitTypeInImmediateChild(); datumEnsemble()->_exportToJSON(formatter); } writer->AddObjKey("coordinate_system"); formatter->setOmitTypeInImmediateChild(); coordinateSystem()->_exportToJSON(formatter); ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static util::IComparable::Criterion getStandardCriterion(util::IComparable::Criterion criterion) { return criterion == util::IComparable::Criterion:: EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS ? util::IComparable::Criterion::EQUIVALENT : criterion; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress bool GeodeticCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { if (other == nullptr || !util::isOfExactType(*other)) { return false; } return _isEquivalentToNoTypeCheck(other, criterion, dbContext); } bool GeodeticCRS::_isEquivalentToNoTypeCheck( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { const auto standardCriterion = getStandardCriterion(criterion); // TODO test velocityModel return SingleCRS::baseIsEquivalentTo(other, standardCriterion, dbContext); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static util::PropertyMap createMapNameEPSGCode(const char *name, int code) { return util::PropertyMap() .set(common::IdentifiedObject::NAME_KEY, name) .set(metadata::Identifier::CODESPACE_KEY, metadata::Identifier::EPSG) .set(metadata::Identifier::CODE_KEY, code); } //! @endcond // --------------------------------------------------------------------------- GeodeticCRSNNPtr GeodeticCRS::createEPSG_4978() { return create( createMapNameEPSGCode("WGS 84", 4978), datum::GeodeticReferenceFrame::EPSG_6326, cs::CartesianCS::createGeocentric(common::UnitOfMeasure::METRE)); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static bool hasCodeCompatibleOfAuthorityFactory( const common::IdentifiedObject *obj, const io::AuthorityFactoryPtr &authorityFactory) { const auto &ids = obj->identifiers(); if (!ids.empty() && authorityFactory->getAuthority().empty()) { return true; } for (const auto &id : ids) { if (*(id->codeSpace()) == authorityFactory->getAuthority()) { return true; } } return false; } static bool hasCodeCompatibleOfAuthorityFactory( const metadata::IdentifierNNPtr &id, const io::AuthorityFactoryPtr &authorityFactory) { if (authorityFactory->getAuthority().empty()) { return true; } return *(id->codeSpace()) == authorityFactory->getAuthority(); } //! @endcond // --------------------------------------------------------------------------- /** \brief Identify the CRS with reference CRSs. * * The candidate CRSs are either hard-coded, or looked in the database when * authorityFactory is not null. * * Note that the implementation uses a set of heuristics to have a good * compromise of successful identifications over execution time. It might miss * legitimate matches in some circumstances. * * The method returns a list of matching reference CRS, and the percentage * (0-100) of confidence in the match: *
    *
  • 100% means that the name of the reference entry * perfectly matches the CRS name, and both are equivalent. In which case a * single result is returned. * Note: in the case of a GeographicCRS whose axis * order is implicit in the input definition (for example ESRI WKT), then axis * order is ignored for the purpose of identification. That is the CRS built * from * GEOGCS["GCS_WGS_1984",DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137.0,298.257223563]], * PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]] * will be identified to EPSG:4326, but will not pass a * isEquivalentTo(EPSG_4326, util::IComparable::Criterion::EQUIVALENT) test, * but rather isEquivalentTo(EPSG_4326, * util::IComparable::Criterion::EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS) *
    • *
  • 90% means that CRS are equivalent, but the names are not exactly the * same. *
  • 70% means that CRS are equivalent (equivalent datum and coordinate * system), * but the names are not equivalent.
    • *
  • 60% means that ellipsoid, prime meridian and coordinate systems are * equivalent, but the CRS and datum names do not match.
    • *
  • 25% means that the CRS are not equivalent, but there is some similarity * in * the names.
    • *
* * @param authorityFactory Authority factory (or null, but degraded * functionality) * @return a list of matching reference CRS, and the percentage (0-100) of * confidence in the match. */ std::list> GeodeticCRS::identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; const auto &thisName(nameStr()); io::DatabaseContextPtr dbContext = authorityFactory ? authorityFactory->databaseContext().as_nullable() : nullptr; const bool l_implicitCS = hasImplicitCS(); const auto crsCriterion = l_implicitCS ? util::IComparable::Criterion::EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS : util::IComparable::Criterion::EQUIVALENT; if (authorityFactory == nullptr || authorityFactory->getAuthority().empty() || authorityFactory->getAuthority() == metadata::Identifier::EPSG) { const GeographicCRSNNPtr candidatesCRS[] = {GeographicCRS::EPSG_4326, GeographicCRS::EPSG_4267, GeographicCRS::EPSG_4269}; for (const auto &crs : candidatesCRS) { const bool nameEquivalent = metadata::Identifier::isEquivalentName( thisName.c_str(), crs->nameStr().c_str()); const bool nameEqual = thisName == crs->nameStr(); const bool isEq = _isEquivalentTo(crs.get(), crsCriterion, dbContext); if (nameEquivalent && isEq && (!authorityFactory || nameEqual)) { res.emplace_back(util::nn_static_pointer_cast(crs), nameEqual ? 100 : 90); return res; } else if (nameEqual && !isEq && !authorityFactory) { res.emplace_back(util::nn_static_pointer_cast(crs), 25); return res; } else if (isEq && !authorityFactory) { res.emplace_back(util::nn_static_pointer_cast(crs), 70); return res; } } } std::string geodetic_crs_type; if (isGeocentric()) { geodetic_crs_type = "geocentric"; } else { auto geogCRS = dynamic_cast(this); if (geogCRS) { if (coordinateSystem()->axisList().size() == 2) { geodetic_crs_type = "geographic 2D"; } else { geodetic_crs_type = "geographic 3D"; } } } if (authorityFactory) { const auto thisDatum(datumNonNull(dbContext)); auto searchByDatumCode = [this, &authorityFactory, &res, &geodetic_crs_type, crsCriterion, &dbContext](const common::IdentifiedObjectNNPtr &l_datum) { for (const auto &id : l_datum->identifiers()) { try { auto tempRes = authorityFactory->createGeodeticCRSFromDatum( *id->codeSpace(), id->code(), geodetic_crs_type); for (const auto &crs : tempRes) { if (_isEquivalentTo(crs.get(), crsCriterion, dbContext)) { res.emplace_back(crs, 70); } } } catch (const std::exception &) { } } }; auto searchByEllipsoid = [this, &authorityFactory, &res, &thisDatum, &geodetic_crs_type, l_implicitCS, &dbContext]() { const auto &thisEllipsoid = thisDatum->ellipsoid(); const auto ellipsoids = thisEllipsoid->identifiers().empty() ? authorityFactory->createEllipsoidFromExisting( thisEllipsoid) : std::list{thisEllipsoid}; for (const auto &ellps : ellipsoids) { for (const auto &id : ellps->identifiers()) { try { auto tempRes = authorityFactory->createGeodeticCRSFromEllipsoid( *id->codeSpace(), id->code(), geodetic_crs_type); for (const auto &crs : tempRes) { const auto crsDatum(crs->datumNonNull(dbContext)); if (crsDatum->ellipsoid()->_isEquivalentTo( ellps.get(), util::IComparable::Criterion::EQUIVALENT, dbContext) && crsDatum->primeMeridian()->_isEquivalentTo( thisDatum->primeMeridian().get(), util::IComparable::Criterion::EQUIVALENT, dbContext) && (!l_implicitCS || coordinateSystem()->_isEquivalentTo( crs->coordinateSystem().get(), util::IComparable::Criterion::EQUIVALENT, dbContext))) { res.emplace_back(crs, 60); } } } catch (const std::exception &) { } } } }; const auto searchByDatumOrEllipsoid = [&authorityFactory, &res, &thisDatum, searchByDatumCode, searchByEllipsoid]() { if (!thisDatum->identifiers().empty()) { searchByDatumCode(thisDatum); } else { auto candidateDatums = authorityFactory->createObjectsFromName( thisDatum->nameStr(), {io::AuthorityFactory::ObjectType:: GEODETIC_REFERENCE_FRAME}, false); const size_t sizeBefore = res.size(); for (const auto &candidateDatum : candidateDatums) { searchByDatumCode(candidateDatum); } if (sizeBefore == res.size()) { searchByEllipsoid(); } } }; const bool insignificantName = thisName.empty() || ci_equal(thisName, "unknown") || ci_equal(thisName, "unnamed"); if (insignificantName) { searchByDatumOrEllipsoid(); } else if (hasCodeCompatibleOfAuthorityFactory(this, authorityFactory)) { // If the CRS has already an id, check in the database for the // official object, and verify that they are equivalent. for (const auto &id : identifiers()) { if (hasCodeCompatibleOfAuthorityFactory(id, authorityFactory)) { try { auto crs = io::AuthorityFactory::create( authorityFactory->databaseContext(), *id->codeSpace()) ->createGeodeticCRS(id->code()); bool match = _isEquivalentTo(crs.get(), crsCriterion, dbContext); res.emplace_back(crs, match ? 100 : 25); return res; } catch (const std::exception &) { } } } } else { bool gotAbove25Pct = false; for (int ipass = 0; ipass < 2; ipass++) { const bool approximateMatch = ipass == 1; auto objects = authorityFactory->createObjectsFromName( thisName, {io::AuthorityFactory::ObjectType::GEODETIC_CRS}, approximateMatch); for (const auto &obj : objects) { auto crs = util::nn_dynamic_pointer_cast(obj); assert(crs); auto crsNN = NN_NO_CHECK(crs); if (_isEquivalentTo(crs.get(), crsCriterion, dbContext)) { if (crs->nameStr() == thisName) { res.clear(); res.emplace_back(crsNN, 100); return res; } const bool eqName = metadata::Identifier::isEquivalentName( thisName.c_str(), crs->nameStr().c_str()); res.emplace_back(crsNN, eqName ? 90 : 70); gotAbove25Pct = true; } else { res.emplace_back(crsNN, 25); } } if (!res.empty()) { break; } } if (!gotAbove25Pct) { searchByDatumOrEllipsoid(); } } const auto &thisCS(coordinateSystem()); // Sort results res.sort([&thisName, &thisDatum, &thisCS, &dbContext](const Pair &a, const Pair &b) { // First consider confidence if (a.second > b.second) { return true; } if (a.second < b.second) { return false; } // Then consider exact name matching const auto &aName(a.first->nameStr()); const auto &bName(b.first->nameStr()); if (aName == thisName && bName != thisName) { return true; } if (bName == thisName && aName != thisName) { return false; } // Then datum matching const auto aDatum(a.first->datumNonNull(dbContext)); const auto bDatum(b.first->datumNonNull(dbContext)); const auto thisEquivADatum(thisDatum->_isEquivalentTo( aDatum.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)); const auto thisEquivBDatum(thisDatum->_isEquivalentTo( bDatum.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)); if (thisEquivADatum && !thisEquivBDatum) { return true; } if (!thisEquivADatum && thisEquivBDatum) { return false; } // Then coordinate system matching const auto &aCS(a.first->coordinateSystem()); const auto &bCS(b.first->coordinateSystem()); const auto thisEquivACs(thisCS->_isEquivalentTo( aCS.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)); const auto thisEquivBCs(thisCS->_isEquivalentTo( bCS.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)); if (thisEquivACs && !thisEquivBCs) { return true; } if (!thisEquivACs && thisEquivBCs) { return false; } // Then dimension of the coordinate system matching const auto thisCSAxisListSize = thisCS->axisList().size(); const auto aCSAxistListSize = aCS->axisList().size(); const auto bCSAxistListSize = bCS->axisList().size(); if (thisCSAxisListSize == aCSAxistListSize && thisCSAxisListSize != bCSAxistListSize) { return true; } if (thisCSAxisListSize != aCSAxistListSize && thisCSAxisListSize == bCSAxistListSize) { return false; } // Favor the CRS whole ellipsoid names matches the ellipsoid // name (WGS84...) const bool aEllpsNameEqCRSName = metadata::Identifier::isEquivalentName( aDatum->ellipsoid()->nameStr().c_str(), a.first->nameStr().c_str()); const bool bEllpsNameEqCRSName = metadata::Identifier::isEquivalentName( bDatum->ellipsoid()->nameStr().c_str(), b.first->nameStr().c_str()); if (aEllpsNameEqCRSName && !bEllpsNameEqCRSName) { return true; } if (bEllpsNameEqCRSName && !aEllpsNameEqCRSName) { return false; } // Arbitrary final sorting criterion return aName < bName; }); // If there are results with 90% confidence, only keep those if (res.size() >= 2 && res.front().second == 90) { std::list newRes; for (const auto &pair : res) { if (pair.second == 90) { newRes.push_back(pair); } else { break; } } return newRes; } } return res; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> GeodeticCRS::_identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; auto resTemp = identify(authorityFactory); for (const auto &pair : resTemp) { res.emplace_back(pair.first, pair.second); } return res; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct GeographicCRS::Private { cs::EllipsoidalCSNNPtr coordinateSystem_; explicit Private(const cs::EllipsoidalCSNNPtr &csIn) : coordinateSystem_(csIn) {} }; //! @endcond // --------------------------------------------------------------------------- GeographicCRS::GeographicCRS(const datum::GeodeticReferenceFramePtr &datumIn, const datum::DatumEnsemblePtr &datumEnsembleIn, const cs::EllipsoidalCSNNPtr &csIn) : SingleCRS(datumIn, datumEnsembleIn, csIn), GeodeticCRS(datumIn, checkEnsembleForGeodeticCRS(datumIn, datumEnsembleIn), csIn), d(internal::make_unique(csIn)) {} // --------------------------------------------------------------------------- GeographicCRS::GeographicCRS(const GeographicCRS &other) : SingleCRS(other), GeodeticCRS(other), d(internal::make_unique(*other.d)) {} // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress GeographicCRS::~GeographicCRS() = default; //! @endcond // --------------------------------------------------------------------------- CRSNNPtr GeographicCRS::_shallowClone() const { auto crs(GeographicCRS::nn_make_shared(*this)); crs->assignSelf(crs); return crs; } // --------------------------------------------------------------------------- /** \brief Return the cs::EllipsoidalCS associated with the CRS. * * @return a EllipsoidalCS. */ const cs::EllipsoidalCSNNPtr &GeographicCRS::coordinateSystem() PROJ_PURE_DEFN { return d->coordinateSystem_; } // --------------------------------------------------------------------------- /** \brief Instantiate a GeographicCRS from a datum::GeodeticReferenceFrameNNPtr * and a * cs::EllipsoidalCS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datum The datum of the CRS. * @param cs a EllipsoidalCS. * @return new GeographicCRS. */ GeographicCRSNNPtr GeographicCRS::create(const util::PropertyMap &properties, const datum::GeodeticReferenceFrameNNPtr &datum, const cs::EllipsoidalCSNNPtr &cs) { return create(properties, datum.as_nullable(), nullptr, cs); } // --------------------------------------------------------------------------- /** \brief Instantiate a GeographicCRS from a datum::GeodeticReferenceFramePtr * or * datum::DatumEnsemble and a * cs::EllipsoidalCS. * * One and only one of datum or datumEnsemble should be set to a non-null value. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datum The datum of the CRS, or nullptr * @param datumEnsemble The datum ensemble of the CRS, or nullptr. * @param cs a EllipsoidalCS. * @return new GeographicCRS. */ GeographicCRSNNPtr GeographicCRS::create(const util::PropertyMap &properties, const datum::GeodeticReferenceFramePtr &datum, const datum::DatumEnsemblePtr &datumEnsemble, const cs::EllipsoidalCSNNPtr &cs) { GeographicCRSNNPtr crs( GeographicCRS::nn_make_shared(datum, datumEnsemble, cs)); crs->assignSelf(crs); crs->setProperties(properties); properties.getStringValue("EXTENSION_PROJ4", crs->CRS::getPrivate()->extensionProj4_); crs->CRS::getPrivate()->setImplicitCS(properties); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress /** \brief Return whether the current GeographicCRS is the 2D part of the * other 3D GeographicCRS. */ bool GeographicCRS::is2DPartOf3D(util::nn other, const io::DatabaseContextPtr &dbContext) PROJ_PURE_DEFN { const auto &axis = d->coordinateSystem_->axisList(); const auto &otherAxis = other->d->coordinateSystem_->axisList(); if (!(axis.size() == 2 && otherAxis.size() == 3)) { return false; } const auto &firstAxis = axis[0]; const auto &secondAxis = axis[1]; const auto &otherFirstAxis = otherAxis[0]; const auto &otherSecondAxis = otherAxis[1]; if (!(firstAxis->_isEquivalentTo( otherFirstAxis.get(), util::IComparable::Criterion::EQUIVALENT) && secondAxis->_isEquivalentTo( otherSecondAxis.get(), util::IComparable::Criterion::EQUIVALENT))) { return false; } try { const auto thisDatum = datumNonNull(dbContext); const auto otherDatum = other->datumNonNull(dbContext); return thisDatum->_isEquivalentTo( otherDatum.get(), util::IComparable::Criterion::EQUIVALENT); } catch (const util::InvalidValueTypeException &) { // should not happen really, but potentially thrown by // Identifier::Private::setProperties() assert(false); return false; } } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress bool GeographicCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { if (other == nullptr || !util::isOfExactType(*other)) { return false; } const auto standardCriterion = getStandardCriterion(criterion); if (GeodeticCRS::_isEquivalentToNoTypeCheck(other, standardCriterion, dbContext)) { return true; } if (criterion != util::IComparable::Criterion::EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS) { return false; } const auto axisOrder = coordinateSystem()->axisOrder(); if (axisOrder == cs::EllipsoidalCS::AxisOrder::LONG_EAST_LAT_NORTH || axisOrder == cs::EllipsoidalCS::AxisOrder::LAT_NORTH_LONG_EAST) { const auto &unit = coordinateSystem()->axisList()[0]->unit(); return GeographicCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, nameStr()), datum(), datumEnsemble(), axisOrder == cs::EllipsoidalCS::AxisOrder::LONG_EAST_LAT_NORTH ? cs::EllipsoidalCS::createLatitudeLongitude(unit) : cs::EllipsoidalCS::createLongitudeLatitude(unit)) ->GeodeticCRS::_isEquivalentToNoTypeCheck(other, standardCriterion, dbContext); } if (axisOrder == cs::EllipsoidalCS::AxisOrder::LONG_EAST_LAT_NORTH_HEIGHT_UP || axisOrder == cs::EllipsoidalCS::AxisOrder::LAT_NORTH_LONG_EAST_HEIGHT_UP) { const auto &angularUnit = coordinateSystem()->axisList()[0]->unit(); const auto &linearUnit = coordinateSystem()->axisList()[2]->unit(); return GeographicCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, nameStr()), datum(), datumEnsemble(), axisOrder == cs::EllipsoidalCS::AxisOrder:: LONG_EAST_LAT_NORTH_HEIGHT_UP ? cs::EllipsoidalCS:: createLatitudeLongitudeEllipsoidalHeight( angularUnit, linearUnit) : cs::EllipsoidalCS:: createLongitudeLatitudeEllipsoidalHeight( angularUnit, linearUnit)) ->GeodeticCRS::_isEquivalentToNoTypeCheck(other, standardCriterion, dbContext); } return false; } //! @endcond // --------------------------------------------------------------------------- GeographicCRSNNPtr GeographicCRS::createEPSG_4267() { return create(createMapNameEPSGCode("NAD27", 4267), datum::GeodeticReferenceFrame::EPSG_6267, cs::EllipsoidalCS::createLatitudeLongitude( common::UnitOfMeasure::DEGREE)); } // --------------------------------------------------------------------------- GeographicCRSNNPtr GeographicCRS::createEPSG_4269() { return create(createMapNameEPSGCode("NAD83", 4269), datum::GeodeticReferenceFrame::EPSG_6269, cs::EllipsoidalCS::createLatitudeLongitude( common::UnitOfMeasure::DEGREE)); } // --------------------------------------------------------------------------- GeographicCRSNNPtr GeographicCRS::createEPSG_4326() { return create(createMapNameEPSGCode("WGS 84", 4326), datum::GeodeticReferenceFrame::EPSG_6326, cs::EllipsoidalCS::createLatitudeLongitude( common::UnitOfMeasure::DEGREE)); } // --------------------------------------------------------------------------- GeographicCRSNNPtr GeographicCRS::createOGC_CRS84() { util::PropertyMap propertiesCRS; propertiesCRS .set(metadata::Identifier::CODESPACE_KEY, metadata::Identifier::OGC) .set(metadata::Identifier::CODE_KEY, "CRS84") .set(common::IdentifiedObject::NAME_KEY, "WGS 84 (CRS84)"); return create(propertiesCRS, datum::GeodeticReferenceFrame::EPSG_6326, cs::EllipsoidalCS::createLongitudeLatitude( // Long Lat ! common::UnitOfMeasure::DEGREE)); } // --------------------------------------------------------------------------- GeographicCRSNNPtr GeographicCRS::createEPSG_4979() { return create( createMapNameEPSGCode("WGS 84", 4979), datum::GeodeticReferenceFrame::EPSG_6326, cs::EllipsoidalCS::createLatitudeLongitudeEllipsoidalHeight( common::UnitOfMeasure::DEGREE, common::UnitOfMeasure::METRE)); } // --------------------------------------------------------------------------- GeographicCRSNNPtr GeographicCRS::createEPSG_4807() { auto ellps(datum::Ellipsoid::createFlattenedSphere( createMapNameEPSGCode("Clarke 1880 (IGN)", 7011), common::Length(6378249.2), common::Scale(293.4660212936269))); auto cs(cs::EllipsoidalCS::createLatitudeLongitude( common::UnitOfMeasure::GRAD)); auto datum(datum::GeodeticReferenceFrame::create( createMapNameEPSGCode("Nouvelle Triangulation Francaise (Paris)", 6807), ellps, util::optional(), datum::PrimeMeridian::PARIS)); return create(createMapNameEPSGCode("NTF (Paris)", 4807), datum, cs); } // --------------------------------------------------------------------------- /** \brief Return a variant of this CRS "demoted" to a 2D one, if not already * the case. * * * @param newName Name of the new CRS. If empty, nameStr() will be used. * @param dbContext Database context to look for potentially already registered * 2D CRS. May be nullptr. * @return a new CRS demoted to 2D, or the current one if already 2D or not * applicable. * @since 6.3 */ GeographicCRSNNPtr GeographicCRS::demoteTo2D(const std::string &newName, const io::DatabaseContextPtr &dbContext) const { const auto &axisList = coordinateSystem()->axisList(); if (axisList.size() == 3) { const auto &l_identifiers = identifiers(); // First check if there is a Geographic 2D CRS in the database // of the same name. // This is the common practice in the EPSG dataset. if (dbContext && l_identifiers.size() == 1) { auto authFactory = io::AuthorityFactory::create( NN_NO_CHECK(dbContext), *(l_identifiers[0]->codeSpace())); auto res = authFactory->createObjectsFromName( nameStr(), {io::AuthorityFactory::ObjectType::GEOGRAPHIC_2D_CRS}, false); if (!res.empty()) { const auto &firstRes = res.front(); auto firstResAsGeogCRS = util::nn_dynamic_pointer_cast(firstRes); if (firstResAsGeogCRS && firstResAsGeogCRS->is2DPartOf3D( NN_NO_CHECK(this), dbContext)) { return NN_NO_CHECK(firstResAsGeogCRS); } } } auto cs = cs::EllipsoidalCS::create(util::PropertyMap(), axisList[0], axisList[1]); return GeographicCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, !newName.empty() ? newName : nameStr()), datum(), datumEnsemble(), cs); } return NN_NO_CHECK(std::dynamic_pointer_cast( shared_from_this().as_nullable())); } // --------------------------------------------------------------------------- //! @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) { const auto &extensionProj4 = CRS::getPrivate()->extensionProj4_; if (!extensionProj4.empty()) { formatter->ingestPROJString( replaceAll(extensionProj4, " +type=crs", "")); formatter->addNoDefs(false); return; } if (!formatter->omitProjLongLatIfPossible() || primeMeridian()->longitude().getSIValue() != 0.0 || !formatter->getTOWGS84Parameters().empty() || !formatter->getHDatumExtension().empty()) { formatter->addStep("longlat"); bool done = false; if (formatter->getLegacyCRSToCRSContext() && formatter->getHDatumExtension().empty() && formatter->getTOWGS84Parameters().empty()) { const auto l_datum = datumNonNull(formatter->databaseContext()); if (l_datum->_isEquivalentTo( datum::GeodeticReferenceFrame::EPSG_6326.get(), util::IComparable::Criterion::EQUIVALENT)) { done = true; formatter->addParam("ellps", "WGS84"); } else if (l_datum->_isEquivalentTo( datum::GeodeticReferenceFrame::EPSG_6269.get(), util::IComparable::Criterion::EQUIVALENT)) { done = true; // We do not want datum=NAD83 to cause a useless towgs84=0,0,0 formatter->addParam("ellps", "GRS80"); } } if (!done) { addDatumInfoToPROJString(formatter); } } if (!formatter->getCRSExport()) { addAngularUnitConvertAndAxisSwap(formatter); } } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void GeographicCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("GeographicCRS", !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } const auto &l_datum(datum()); if (l_datum) { writer->AddObjKey("datum"); l_datum->_exportToJSON(formatter); } else { writer->AddObjKey("datum_ensemble"); formatter->setOmitTypeInImmediateChild(); datumEnsemble()->_exportToJSON(formatter); } writer->AddObjKey("coordinate_system"); formatter->setOmitTypeInImmediateChild(); coordinateSystem()->_exportToJSON(formatter); ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct VerticalCRS::Private { std::vector geoidModel{}; std::vector velocityModel{}; }; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static const datum::DatumEnsemblePtr & checkEnsembleForVerticalCRS(const datum::VerticalReferenceFramePtr &datumIn, const datum::DatumEnsemblePtr &ensemble) { const char *msg = "One of Datum or DatumEnsemble should be defined"; if (datumIn) { if (!ensemble) { return ensemble; } msg = "Datum and DatumEnsemble should not be defined"; } else if (ensemble) { const auto &datums = ensemble->datums(); assert(!datums.empty()); auto grfFirst = dynamic_cast(datums[0].get()); if (grfFirst) { return ensemble; } msg = "Ensemble should contain VerticalReferenceFrame"; } throw util::Exception(msg); } //! @endcond // --------------------------------------------------------------------------- VerticalCRS::VerticalCRS(const datum::VerticalReferenceFramePtr &datumIn, const datum::DatumEnsemblePtr &datumEnsembleIn, const cs::VerticalCSNNPtr &csIn) : SingleCRS(datumIn, checkEnsembleForVerticalCRS(datumIn, datumEnsembleIn), csIn), d(internal::make_unique()) {} // --------------------------------------------------------------------------- VerticalCRS::VerticalCRS(const VerticalCRS &other) : SingleCRS(other), d(internal::make_unique(*other.d)) {} // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress VerticalCRS::~VerticalCRS() = default; //! @endcond // --------------------------------------------------------------------------- CRSNNPtr VerticalCRS::_shallowClone() const { auto crs(VerticalCRS::nn_make_shared(*this)); crs->assignSelf(crs); return crs; } // --------------------------------------------------------------------------- /** \brief Return the datum::VerticalReferenceFrame associated with the CRS. * * @return a VerticalReferenceFrame. */ const datum::VerticalReferenceFramePtr VerticalCRS::datum() const { return std::static_pointer_cast( SingleCRS::getPrivate()->datum); } // --------------------------------------------------------------------------- /** \brief Return the geoid model associated with the CRS. * * Geoid height model or height correction model linked to a geoid-based * vertical CRS. * * @return a geoid model. might be null */ const std::vector & VerticalCRS::geoidModel() PROJ_PURE_DEFN { return d->geoidModel; } // --------------------------------------------------------------------------- /** \brief Return the velocity model associated with the CRS. * * @return a velocity model. might be null. */ const std::vector & VerticalCRS::velocityModel() PROJ_PURE_DEFN { return d->velocityModel; } // --------------------------------------------------------------------------- /** \brief Return the cs::VerticalCS associated with the CRS. * * @return a VerticalCS. */ const cs::VerticalCSNNPtr VerticalCRS::coordinateSystem() const { return util::nn_static_pointer_cast( SingleCRS::getPrivate()->coordinateSystem); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress /** \brief Return the real datum or a synthetized one if a datumEnsemble. */ const datum::VerticalReferenceFrameNNPtr VerticalCRS::datumNonNull(const io::DatabaseContextPtr &dbContext) const { return NN_NO_CHECK( util::nn_dynamic_pointer_cast( SingleCRS::datumNonNull(dbContext))); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void VerticalCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; formatter->startNode(isWKT2 ? io::WKTConstants::VERTCRS : formatter->useESRIDialect() ? io::WKTConstants::VERTCS : io::WKTConstants::VERT_CS, !identifiers().empty()); auto l_name = nameStr(); const auto &dbContext = formatter->databaseContext(); if (formatter->useESRIDialect()) { bool aliasFound = false; if (dbContext) { auto l_alias = dbContext->getAliasFromOfficialName( l_name, "vertical_crs", "ESRI"); if (!l_alias.empty()) { l_name = l_alias; aliasFound = true; } } if (!aliasFound) { l_name = io::WKTFormatter::morphNameToESRI(l_name); } } formatter->addQuotedString(l_name); const auto l_datum = datum(); if (formatter->useESRIDialect() && l_datum && l_datum->getWKT1DatumType() == "2002") { bool foundMatch = false; if (dbContext) { auto authFactory = io::AuthorityFactory::create( NN_NO_CHECK(dbContext), std::string()); auto list = authFactory->createObjectsFromName( l_datum->nameStr(), {io::AuthorityFactory::ObjectType::GEODETIC_REFERENCE_FRAME}, false /* approximate=false*/); if (!list.empty()) { auto gdatum = util::nn_dynamic_pointer_cast(list.front()); if (gdatum) { gdatum->_exportToWKT(formatter); foundMatch = true; } } } if (!foundMatch) { // We should export a geodetic datum, but we cannot really do better l_datum->_exportToWKT(formatter); } } else { exportDatumOrDatumEnsembleToWkt(formatter); } const auto &cs = SingleCRS::getPrivate()->coordinateSystem; const auto &axisList = cs->axisList(); if (formatter->useESRIDialect()) { // Seems to be a constant value... formatter->startNode(io::WKTConstants::PARAMETER, false); formatter->addQuotedString("Vertical_Shift"); formatter->add(0.0); formatter->endNode(); formatter->startNode(io::WKTConstants::PARAMETER, false); formatter->addQuotedString("Direction"); formatter->add( axisList[0]->direction() == cs::AxisDirection::UP ? 1.0 : -1.0); formatter->endNode(); } if (!isWKT2) { axisList[0]->unit()._exportToWKT(formatter); } const auto oldAxisOutputRule = formatter->outputAxis(); if (oldAxisOutputRule == io::WKTFormatter::OutputAxisRule::WKT1_GDAL_EPSG_STYLE) { formatter->setOutputAxis(io::WKTFormatter::OutputAxisRule::YES); } cs->_exportToWKT(formatter); formatter->setOutputAxis(oldAxisOutputRule); if (isWKT2 && formatter->use2019Keywords() && !d->geoidModel.empty()) { const auto &model = d->geoidModel[0]; formatter->startNode(io::WKTConstants::GEOIDMODEL, false); formatter->addQuotedString(model->nameStr()); model->formatID(formatter); formatter->endNode(); } ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void VerticalCRS::_exportToPROJString( io::PROJStringFormatter *formatter) const // throw(io::FormattingException) { auto geoidgrids = formatter->getVDatumExtension(); if (!geoidgrids.empty()) { formatter->addParam("geoidgrids", geoidgrids); } auto &axisList = coordinateSystem()->axisList(); if (!axisList.empty()) { auto projUnit = axisList[0]->unit().exportToPROJString(); if (projUnit.empty()) { formatter->addParam("vto_meter", axisList[0]->unit().conversionToSI()); } else { formatter->addParam("vunits", projUnit); } } } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void VerticalCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("VerticalCRS", !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } const auto &l_datum(datum()); if (l_datum) { writer->AddObjKey("datum"); l_datum->_exportToJSON(formatter); } else { writer->AddObjKey("datum_ensemble"); formatter->setOmitTypeInImmediateChild(); datumEnsemble()->_exportToJSON(formatter); } writer->AddObjKey("coordinate_system"); formatter->setOmitTypeInImmediateChild(); coordinateSystem()->_exportToJSON(formatter); if (!d->geoidModel.empty()) { const auto &model = d->geoidModel[0]; writer->AddObjKey("geoid_model"); auto objectContext2(formatter->MakeObjectContext(nullptr, false)); writer->AddObjKey("name"); writer->Add(model->nameStr()); if (model->identifiers().empty()) { const auto &interpCRS = model->interpolationCRS(); if (interpCRS) { writer->AddObjKey("interpolation_crs"); interpCRS->_exportToJSON(formatter); } } model->formatID(formatter); } ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void VerticalCRS::addLinearUnitConvert( io::PROJStringFormatter *formatter) const { auto &axisList = coordinateSystem()->axisList(); if (!axisList.empty()) { if (axisList[0]->unit().conversionToSI() != 1.0) { formatter->addStep("unitconvert"); formatter->addParam("z_in", "m"); auto projVUnit = axisList[0]->unit().exportToPROJString(); if (projVUnit.empty()) { formatter->addParam("z_out", axisList[0]->unit().conversionToSI()); } else { formatter->addParam("z_out", projVUnit); } } } } //! @endcond // --------------------------------------------------------------------------- /** \brief Instantiate a VerticalCRS from a datum::VerticalReferenceFrame and a * cs::VerticalCS. * * @param properties See \ref general_properties. * At minimum the name should be defined. The GEOID_MODEL property can be set * to a TransformationNNPtr object. * @param datumIn The datum of the CRS. * @param csIn a VerticalCS. * @return new VerticalCRS. */ VerticalCRSNNPtr VerticalCRS::create(const util::PropertyMap &properties, const datum::VerticalReferenceFrameNNPtr &datumIn, const cs::VerticalCSNNPtr &csIn) { return create(properties, datumIn.as_nullable(), nullptr, csIn); } // --------------------------------------------------------------------------- /** \brief Instantiate a VerticalCRS from a datum::VerticalReferenceFrame or * datum::DatumEnsemble and a cs::VerticalCS. * * One and only one of datum or datumEnsemble should be set to a non-null value. * * @param properties See \ref general_properties. * At minimum the name should be defined. The GEOID_MODEL property can be set * to a TransformationNNPtr object. * @param datumIn The datum of the CRS, or nullptr * @param datumEnsembleIn The datum ensemble of the CRS, or nullptr. * @param csIn a VerticalCS. * @return new VerticalCRS. */ VerticalCRSNNPtr VerticalCRS::create(const util::PropertyMap &properties, const datum::VerticalReferenceFramePtr &datumIn, const datum::DatumEnsemblePtr &datumEnsembleIn, const cs::VerticalCSNNPtr &csIn) { auto crs(VerticalCRS::nn_make_shared(datumIn, datumEnsembleIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); const auto geoidModelPtr = properties.get("GEOID_MODEL"); if (geoidModelPtr) { auto transf = util::nn_dynamic_pointer_cast( *geoidModelPtr); if (transf) { crs->d->geoidModel.emplace_back(NN_NO_CHECK(transf)); } } return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress bool VerticalCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherVertCRS = dynamic_cast(other); // TODO test geoidModel and velocityModel return otherVertCRS != nullptr && SingleCRS::baseIsEquivalentTo(other, criterion, dbContext); } //! @endcond // --------------------------------------------------------------------------- /** \brief Identify the CRS with reference CRSs. * * The candidate CRSs are looked in the database when * authorityFactory is not null. * * Note that the implementation uses a set of heuristics to have a good * compromise of successful identifications over execution time. It might miss * legitimate matches in some circumstances. * * The method returns a list of matching reference CRS, and the percentage * (0-100) of confidence in the match. * 100% means that the name of the reference entry * perfectly matches the CRS name, and both are equivalent. In which case a * single result is returned. * 90% means that CRS are equivalent, but the names are not exactly the same. * 70% means that CRS are equivalent (equivalent datum and coordinate system), * but the names are not equivalent. * 25% means that the CRS are not equivalent, but there is some similarity in * the names. * * @param authorityFactory Authority factory (if null, will return an empty * list) * @return a list of matching reference CRS, and the percentage (0-100) of * confidence in the match. */ std::list> VerticalCRS::identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; const auto &thisName(nameStr()); if (authorityFactory) { const io::DatabaseContextNNPtr &dbContext = authorityFactory->databaseContext(); const bool insignificantName = thisName.empty() || ci_equal(thisName, "unknown") || ci_equal(thisName, "unnamed"); if (hasCodeCompatibleOfAuthorityFactory(this, authorityFactory)) { // If the CRS has already an id, check in the database for the // official object, and verify that they are equivalent. for (const auto &id : identifiers()) { if (hasCodeCompatibleOfAuthorityFactory(id, authorityFactory)) { try { auto crs = io::AuthorityFactory::create( dbContext, *id->codeSpace()) ->createVerticalCRS(id->code()); bool match = _isEquivalentTo( crs.get(), util::IComparable::Criterion::EQUIVALENT, dbContext); res.emplace_back(crs, match ? 100 : 25); return res; } catch (const std::exception &) { } } } } else if (!insignificantName) { for (int ipass = 0; ipass < 2; ipass++) { const bool approximateMatch = ipass == 1; auto objects = authorityFactory->createObjectsFromName( thisName, {io::AuthorityFactory::ObjectType::VERTICAL_CRS}, approximateMatch); for (const auto &obj : objects) { auto crs = util::nn_dynamic_pointer_cast(obj); assert(crs); auto crsNN = NN_NO_CHECK(crs); if (_isEquivalentTo( crs.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { if (crs->nameStr() == thisName) { res.clear(); res.emplace_back(crsNN, 100); return res; } res.emplace_back(crsNN, 90); } else { res.emplace_back(crsNN, 25); } } if (!res.empty()) { break; } } } // Sort results res.sort([&thisName](const Pair &a, const Pair &b) { // First consider confidence if (a.second > b.second) { return true; } if (a.second < b.second) { return false; } // Then consider exact name matching const auto &aName(a.first->nameStr()); const auto &bName(b.first->nameStr()); if (aName == thisName && bName != thisName) { return true; } if (bName == thisName && aName != thisName) { return false; } // Arbitrary final sorting criterion return aName < bName; }); // Keep only results of the highest confidence if (res.size() >= 2) { const auto highestConfidence = res.front().second; std::list newRes; for (const auto &pair : res) { if (pair.second == highestConfidence) { newRes.push_back(pair); } else { break; } } return newRes; } } return res; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> VerticalCRS::_identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; auto resTemp = identify(authorityFactory); for (const auto &pair : resTemp) { res.emplace_back(pair.first, pair.second); } return res; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct DerivedCRS::Private { SingleCRSNNPtr baseCRS_; operation::ConversionNNPtr derivingConversion_; Private(const SingleCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn) : baseCRS_(baseCRSIn), derivingConversion_(derivingConversionIn) {} // For the conversion make a _shallowClone(), so that we can later set // its targetCRS to this. Private(const Private &other) : baseCRS_(other.baseCRS_), derivingConversion_(other.derivingConversion_->shallowClone()) {} }; //! @endcond // --------------------------------------------------------------------------- // DerivedCRS is an abstract class, that virtually inherits from SingleCRS // Consequently the base constructor in SingleCRS will never be called by // that constructor. clang -Wabstract-vbase-init and VC++ underline this, but // other // compilers will complain if we don't call the base constructor. DerivedCRS::DerivedCRS(const SingleCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::CoordinateSystemNNPtr & #if !defined(COMPILER_WARNS_ABOUT_ABSTRACT_VBASE_INIT) cs #endif ) : #if !defined(COMPILER_WARNS_ABOUT_ABSTRACT_VBASE_INIT) SingleCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), cs), #endif d(internal::make_unique(baseCRSIn, derivingConversionIn)) { } // --------------------------------------------------------------------------- DerivedCRS::DerivedCRS(const DerivedCRS &other) : #if !defined(COMPILER_WARNS_ABOUT_ABSTRACT_VBASE_INIT) SingleCRS(other), #endif d(internal::make_unique(*other.d)) { } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress DerivedCRS::~DerivedCRS() = default; //! @endcond // --------------------------------------------------------------------------- /** \brief Return the base CRS of a DerivedCRS. * * @return the base CRS. */ const SingleCRSNNPtr &DerivedCRS::baseCRS() PROJ_PURE_DEFN { return d->baseCRS_; } // --------------------------------------------------------------------------- /** \brief Return the deriving conversion from the base CRS to this CRS. * * @return the deriving conversion. */ const operation::ConversionNNPtr DerivedCRS::derivingConversion() const { return d->derivingConversion_->shallowClone(); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress const operation::ConversionNNPtr & DerivedCRS::derivingConversionRef() PROJ_PURE_DEFN { return d->derivingConversion_; } //! @endcond // --------------------------------------------------------------------------- bool DerivedCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherDerivedCRS = dynamic_cast(other); const auto standardCriterion = getStandardCriterion(criterion); if (otherDerivedCRS == nullptr || !SingleCRS::baseIsEquivalentTo(other, standardCriterion, dbContext)) { return false; } return d->baseCRS_->_isEquivalentTo(otherDerivedCRS->d->baseCRS_.get(), criterion, dbContext) && d->derivingConversion_->_isEquivalentTo( otherDerivedCRS->d->derivingConversion_.get(), standardCriterion, dbContext); } // --------------------------------------------------------------------------- void DerivedCRS::setDerivingConversionCRS() { derivingConversionRef()->setWeakSourceTargetCRS( baseCRS().as_nullable(), std::static_pointer_cast(shared_from_this().as_nullable())); } // --------------------------------------------------------------------------- void DerivedCRS::baseExportToWKT(io::WKTFormatter *formatter, const std::string &keyword, const std::string &baseKeyword) const { formatter->startNode(keyword, !identifiers().empty()); formatter->addQuotedString(nameStr()); const auto &l_baseCRS = d->baseCRS_; formatter->startNode(baseKeyword, formatter->use2019Keywords() && !l_baseCRS->identifiers().empty()); formatter->addQuotedString(l_baseCRS->nameStr()); l_baseCRS->exportDatumOrDatumEnsembleToWkt(formatter); if (formatter->use2019Keywords() && !(formatter->idOnTopLevelOnly() && formatter->topLevelHasId())) { l_baseCRS->formatID(formatter); } formatter->endNode(); formatter->setUseDerivingConversion(true); derivingConversionRef()->_exportToWKT(formatter); formatter->setUseDerivingConversion(false); coordinateSystem()->_exportToWKT(formatter); ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void DerivedCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext(className(), !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } writer->AddObjKey("base_crs"); baseCRS()->_exportToJSON(formatter); writer->AddObjKey("conversion"); formatter->setOmitTypeInImmediateChild(); derivingConversionRef()->_exportToJSON(formatter); writer->AddObjKey("coordinate_system"); formatter->setOmitTypeInImmediateChild(); coordinateSystem()->_exportToJSON(formatter); ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct ProjectedCRS::Private { GeodeticCRSNNPtr baseCRS_; cs::CartesianCSNNPtr cs_; Private(const GeodeticCRSNNPtr &baseCRSIn, const cs::CartesianCSNNPtr &csIn) : baseCRS_(baseCRSIn), cs_(csIn) {} inline const GeodeticCRSNNPtr &baseCRS() const { return baseCRS_; } inline const cs::CartesianCSNNPtr &coordinateSystem() const { return cs_; } }; //! @endcond // --------------------------------------------------------------------------- ProjectedCRS::ProjectedCRS( const GeodeticCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::CartesianCSNNPtr &csIn) : SingleCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), DerivedCRS(baseCRSIn, derivingConversionIn, csIn), d(internal::make_unique(baseCRSIn, csIn)) {} // --------------------------------------------------------------------------- ProjectedCRS::ProjectedCRS(const ProjectedCRS &other) : SingleCRS(other), DerivedCRS(other), d(internal::make_unique(other.baseCRS(), other.coordinateSystem())) {} // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress ProjectedCRS::~ProjectedCRS() = default; //! @endcond // --------------------------------------------------------------------------- CRSNNPtr ProjectedCRS::_shallowClone() const { auto crs(ProjectedCRS::nn_make_shared(*this)); crs->assignSelf(crs); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- /** \brief Return the base CRS (a GeodeticCRS, which is generally a * GeographicCRS) of the ProjectedCRS. * * @return the base CRS. */ const GeodeticCRSNNPtr &ProjectedCRS::baseCRS() PROJ_PURE_DEFN { return d->baseCRS(); } // --------------------------------------------------------------------------- /** \brief Return the cs::CartesianCS associated with the CRS. * * @return a CartesianCS */ const cs::CartesianCSNNPtr &ProjectedCRS::coordinateSystem() PROJ_PURE_DEFN { return d->coordinateSystem(); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void ProjectedCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; const auto &l_identifiers = identifiers(); // Try to perfectly round-trip ESRI projectedCRS if the current object // perfectly matches the database definition 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())) { auto projCRS2D = demoteTo2D(std::string(), dbContext); if (dbContext) { const auto res = projCRS2D->identify(io::AuthorityFactory::create( NN_NO_CHECK(dbContext), metadata::Identifier::EPSG)); if (res.size() == 1) { const auto &front = res.front(); if (front.second == 100) { projCRS2D = front.first; } } } if (formatter->useESRIDialect() && dbContext) { // Try to format the ProjecteD 3D CRS as a // PROJCS[],VERTCS[...,DATUM[]] // if we find corresponding objects if (exportAsESRIWktCompoundCRSWithEllipsoidalHeight( this, baseCRS().as_nullable().get(), formatter)) { return; } } if (!formatter->useESRIDialect() && CRS::getPrivate()->allowNonConformantWKT1Export_) { formatter->startNode(io::WKTConstants::COMPD_CS, false); formatter->addQuotedString(l_name + " + " + baseCRS()->nameStr()); projCRS2D->_exportToWKT(formatter); baseCRS() ->demoteTo2D(std::string(), dbContext) ->_exportToWKT(formatter); formatter->endNode(); return; } auto &originalCompoundCRS = CRS::getPrivate()->originalCompoundCRS_; if (!formatter->useESRIDialect() && originalCompoundCRS) { originalCompoundCRS->_exportToWKT(formatter); return; } if (!formatter->useESRIDialect() && formatter->isAllowedEllipsoidalHeightAsVerticalCRS()) { if (exportAsWKT1CompoundCRSWithEllipsoidalHeight( projCRS2D, axisList[2], formatter)) { 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", "ESRI"); } if (!isWKT2 && formatter->useESRIDialect() && !l_identifiers.empty() && *(l_identifiers[0]->codeSpace()) == "ESRI" && dbContext) { try { const auto definition = dbContext->getTextDefinition( "projected_crs", "ESRI", l_identifiers[0]->code()); if (starts_with(definition, "PROJCS")) { auto crsFromFromDef = io::WKTParser() .attachDatabaseContext(dbContext) .createFromWKT(definition); if (_isEquivalentTo( dynamic_cast(crsFromFromDef.get()), util::IComparable::Criterion::EQUIVALENT)) { formatter->ingestWKTNode( io::WKTNode::createFrom(definition)); return; } } } catch (const std::exception &) { } } else if (!isWKT2 && formatter->useESRIDialect() && !l_alias.empty()) { try { auto res = io::AuthorityFactory::create(NN_NO_CHECK(dbContext), "ESRI") ->createObjectsFromName( l_alias, {io::AuthorityFactory::ObjectType::PROJECTED_CRS}, false); if (res.size() == 1) { const auto definition = dbContext->getTextDefinition( "projected_crs", "ESRI", res.front()->identifiers()[0]->code()); if (starts_with(definition, "PROJCS")) { if (_isEquivalentTo( dynamic_cast(res.front().get()), util::IComparable::Criterion::EQUIVALENT)) { formatter->ingestWKTNode( io::WKTNode::createFrom(definition)); return; } } } } catch (const std::exception &) { } } const auto exportAxis = [&l_coordinateSystem, &axisList, &formatter]() { const auto oldAxisOutputRule = formatter->outputAxis(); if (oldAxisOutputRule == io::WKTFormatter::OutputAxisRule::WKT1_GDAL_EPSG_STYLE) { if (&axisList[0]->direction() == &cs::AxisDirection::EAST && &axisList[1]->direction() == &cs::AxisDirection::NORTH) { formatter->setOutputAxis(io::WKTFormatter::OutputAxisRule::YES); } } l_coordinateSystem->_exportToWKT(formatter); formatter->setOutputAxis(oldAxisOutputRule); }; if (!isWKT2 && !formatter->useESRIDialect() && starts_with(nameStr(), "Popular Visualisation CRS / Mercator")) { formatter->startNode(io::WKTConstants::PROJCS, !l_identifiers.empty()); formatter->addQuotedString(nameStr()); formatter->setTOWGS84Parameters({0, 0, 0, 0, 0, 0, 0}); baseCRS()->_exportToWKT(formatter); formatter->setTOWGS84Parameters({}); formatter->startNode(io::WKTConstants::PROJECTION, false); formatter->addQuotedString("Mercator_1SP"); formatter->endNode(); formatter->startNode(io::WKTConstants::PARAMETER, false); formatter->addQuotedString("central_meridian"); formatter->add(0.0); formatter->endNode(); formatter->startNode(io::WKTConstants::PARAMETER, false); formatter->addQuotedString("scale_factor"); formatter->add(1.0); formatter->endNode(); formatter->startNode(io::WKTConstants::PARAMETER, false); formatter->addQuotedString("false_easting"); formatter->add(0.0); formatter->endNode(); formatter->startNode(io::WKTConstants::PARAMETER, false); formatter->addQuotedString("false_northing"); formatter->add(0.0); formatter->endNode(); axisList[0]->unit()._exportToWKT(formatter); exportAxis(); derivingConversionRef()->addWKTExtensionNode(formatter); ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); return; } formatter->startNode(isWKT2 ? io::WKTConstants::PROJCRS : io::WKTConstants::PROJCS, !l_identifiers.empty()); if (formatter->useESRIDialect()) { if (l_alias.empty()) { l_name = io::WKTFormatter::morphNameToESRI(l_name); } else { l_name = l_alias; } } if (!isWKT2 && !formatter->useESRIDialect() && isDeprecated()) { l_name += " (deprecated)"; } formatter->addQuotedString(l_name); const auto &l_baseCRS = d->baseCRS(); const auto &geodeticCRSAxisList = l_baseCRS->coordinateSystem()->axisList(); if (isWKT2) { formatter->startNode( (formatter->use2019Keywords() && dynamic_cast(l_baseCRS.get())) ? io::WKTConstants::BASEGEOGCRS : io::WKTConstants::BASEGEODCRS, formatter->use2019Keywords() && !l_baseCRS->identifiers().empty()); formatter->addQuotedString(l_baseCRS->nameStr()); l_baseCRS->exportDatumOrDatumEnsembleToWkt(formatter); // insert ellipsoidal cs unit when the units of the map // projection angular parameters are not explicitly given within those // parameters. See // http://docs.opengeospatial.org/is/12-063r5/12-063r5.html#61 if (formatter->primeMeridianOrParameterUnitOmittedIfSameAsAxis()) { geodeticCRSAxisList[0]->unit()._exportToWKT(formatter); } l_baseCRS->primeMeridian()->_exportToWKT(formatter); if (formatter->use2019Keywords() && !(formatter->idOnTopLevelOnly() && formatter->topLevelHasId())) { l_baseCRS->formatID(formatter); } formatter->endNode(); } else { const auto oldAxisOutputRule = formatter->outputAxis(); formatter->setOutputAxis(io::WKTFormatter::OutputAxisRule::NO); l_baseCRS->_exportToWKT(formatter); formatter->setOutputAxis(oldAxisOutputRule); } formatter->pushAxisLinearUnit( common::UnitOfMeasure::create(axisList[0]->unit())); formatter->pushAxisAngularUnit( common::UnitOfMeasure::create(geodeticCRSAxisList[0]->unit())); derivingConversionRef()->_exportToWKT(formatter); formatter->popAxisAngularUnit(); formatter->popAxisLinearUnit(); if (!isWKT2) { axisList[0]->unit()._exportToWKT(formatter); } exportAxis(); if (!isWKT2 && !formatter->useESRIDialect()) { const auto &extensionProj4 = CRS::getPrivate()->extensionProj4_; if (!extensionProj4.empty()) { formatter->startNode(io::WKTConstants::EXTENSION, false); formatter->addQuotedString("PROJ4"); formatter->addQuotedString(extensionProj4); formatter->endNode(); } else { derivingConversionRef()->addWKTExtensionNode(formatter); } } ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); return; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void ProjectedCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("ProjectedCRS", !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } writer->AddObjKey("base_crs"); formatter->setAllowIDInImmediateChild(); formatter->setOmitTypeInImmediateChild(); baseCRS()->_exportToJSON(formatter); writer->AddObjKey("conversion"); formatter->setOmitTypeInImmediateChild(); derivingConversionRef()->_exportToJSON(formatter); writer->AddObjKey("coordinate_system"); formatter->setOmitTypeInImmediateChild(); coordinateSystem()->_exportToJSON(formatter); ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- void ProjectedCRS::_exportToPROJString( io::PROJStringFormatter *formatter) const // throw(io::FormattingException) { const auto &extensionProj4 = CRS::getPrivate()->extensionProj4_; if (!extensionProj4.empty()) { formatter->ingestPROJString( replaceAll(extensionProj4, " +type=crs", "")); formatter->addNoDefs(false); return; } derivingConversionRef()->_exportToPROJString(formatter); } // --------------------------------------------------------------------------- /** \brief Instantiate a ProjectedCRS from a base CRS, a deriving * operation::Conversion * and a coordinate system. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param baseCRSIn The base CRS, a GeodeticCRS that is generally a * GeographicCRS. * @param derivingConversionIn The deriving operation::Conversion (typically * using a map * projection method) * @param csIn The coordniate system. * @return new ProjectedCRS. */ ProjectedCRSNNPtr ProjectedCRS::create(const util::PropertyMap &properties, const GeodeticCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::CartesianCSNNPtr &csIn) { auto crs = ProjectedCRS::nn_make_shared( baseCRSIn, derivingConversionIn, csIn); crs->assignSelf(crs); crs->setProperties(properties); crs->setDerivingConversionCRS(); properties.getStringValue("EXTENSION_PROJ4", crs->CRS::getPrivate()->extensionProj4_); crs->CRS::getPrivate()->setImplicitCS(properties); return crs; } // --------------------------------------------------------------------------- bool ProjectedCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { return other != nullptr && util::isOfExactType(*other) && DerivedCRS::_isEquivalentTo(other, criterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress ProjectedCRSNNPtr ProjectedCRS::alterParametersLinearUnit(const common::UnitOfMeasure &unit, bool convertToNewUnit) const { return create( createPropertyMap(this), baseCRS(), derivingConversion()->alterParametersLinearUnit(unit, convertToNewUnit), coordinateSystem()); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void ProjectedCRS::addUnitConvertAndAxisSwap(io::PROJStringFormatter *formatter, bool axisSpecFound) const { const auto &axisList = d->coordinateSystem()->axisList(); const auto &unit = axisList[0]->unit(); const auto *zUnit = axisList.size() == 3 ? &(axisList[2]->unit()) : nullptr; if (!unit._isEquivalentTo(common::UnitOfMeasure::METRE, util::IComparable::Criterion::EQUIVALENT) || (zUnit && !zUnit->_isEquivalentTo(common::UnitOfMeasure::METRE, util::IComparable::Criterion::EQUIVALENT))) { auto projUnit = unit.exportToPROJString(); const double toSI = unit.conversionToSI(); if (!formatter->getCRSExport()) { formatter->addStep("unitconvert"); formatter->addParam("xy_in", "m"); if (zUnit) formatter->addParam("z_in", "m"); if (projUnit.empty()) { formatter->addParam("xy_out", toSI); } else { formatter->addParam("xy_out", projUnit); } if (zUnit) { auto projZUnit = zUnit->exportToPROJString(); const double zToSI = zUnit->conversionToSI(); if (projZUnit.empty()) { formatter->addParam("z_out", zToSI); } else { formatter->addParam("z_out", projZUnit); } } } else { if (projUnit.empty()) { formatter->addParam("to_meter", toSI); } else { formatter->addParam("units", projUnit); } } } else if (formatter->getCRSExport() && !formatter->getLegacyCRSToCRSContext()) { formatter->addParam("units", "m"); } if (!axisSpecFound && !formatter->getCRSExport()) { const auto &dir0 = axisList[0]->direction(); const auto &dir1 = axisList[1]->direction(); if (!(&dir0 == &cs::AxisDirection::EAST && &dir1 == &cs::AxisDirection::NORTH) && // For polar projections, that have south+south direction, // we don't want to mess with axes. dir0 != dir1) { const char *order[2] = {nullptr, nullptr}; 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] = "1"; else if (&dir == &cs::AxisDirection::SOUTH) order[i] = "-2"; else if (&dir == &cs::AxisDirection::NORTH) order[i] = "2"; } if (order[0] && order[1]) { formatter->addStep("axisswap"); char orderStr[10]; sprintf(orderStr, "%.2s,%.2s", order[0], order[1]); formatter->addParam("order", orderStr); } } else { const auto &name0 = axisList[0]->nameStr(); const auto &name1 = axisList[1]->nameStr(); const bool northingEasting = ci_starts_with(name0, "northing") && ci_starts_with(name1, "easting"); // case of EPSG:32661 ["WGS 84 / UPS North (N,E)]" // case of EPSG:32761 ["WGS 84 / UPS South (N,E)]" if (((&dir0 == &cs::AxisDirection::SOUTH && &dir1 == &cs::AxisDirection::SOUTH) || (&dir0 == &cs::AxisDirection::NORTH && &dir1 == &cs::AxisDirection::NORTH)) && northingEasting) { formatter->addStep("axisswap"); formatter->addParam("order", "2,1"); } } } } //! @endcond // --------------------------------------------------------------------------- /** \brief Identify the CRS with reference CRSs. * * The candidate CRSs are either hard-coded, or looked in the database when * authorityFactory is not null. * * Note that the implementation uses a set of heuristics to have a good * compromise of successful identifications over execution time. It might miss * legitimate matches in some circumstances. * * The method returns a list of matching reference CRS, and the percentage * (0-100) of confidence in the match. The list is sorted by decreasing * confidence. * * 100% means that the name of the reference entry * perfectly matches the CRS name, and both are equivalent. In which case a * single result is returned. * 90% means that CRS are equivalent, but the names are not exactly the same. * 70% means that CRS are equivalent (equivalent base CRS, conversion and * coordinate system), but the names are not equivalent. * 60% means that CRS have strong similarity (equivalent base datum, conversion * and coordinate system), but the names are not equivalent. * 50% means that CRS have similarity (equivalent base ellipsoid and * conversion), * but the coordinate system do not match (e.g. different axis ordering or * axis unit). * 25% means that the CRS are not equivalent, but there is some similarity in * the names. * * For the purpose of this function, equivalence is tested with the * util::IComparable::Criterion::EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS, that is * to say that the axis order of the base GeographicCRS is ignored. * * @param authorityFactory Authority factory (or null, but degraded * functionality) * @return a list of matching reference CRS, and the percentage (0-100) of * confidence in the match. */ std::list> ProjectedCRS::identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; const auto &thisName(nameStr()); io::DatabaseContextPtr dbContext = authorityFactory ? authorityFactory->databaseContext().as_nullable() : nullptr; std::list> baseRes; const auto &l_baseCRS(baseCRS()); const auto l_datum = l_baseCRS->datumNonNull(dbContext); const bool significantNameForDatum = !ci_starts_with(l_datum->nameStr(), "unknown") && l_datum->nameStr() != "unnamed"; const auto &ellipsoid = l_baseCRS->ellipsoid(); auto geogCRS = dynamic_cast(l_baseCRS.get()); if (geogCRS && geogCRS->coordinateSystem()->axisOrder() == cs::EllipsoidalCS::AxisOrder::LONG_EAST_LAT_NORTH) { baseRes = GeographicCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, geogCRS->nameStr()), geogCRS->datum(), geogCRS->datumEnsemble(), cs::EllipsoidalCS::createLatitudeLongitude( geogCRS->coordinateSystem()->axisList()[0]->unit())) ->identify(authorityFactory); } else { baseRes = l_baseCRS->identify(authorityFactory); } int zone = 0; bool north = false; auto computeConfidence = [&thisName](const std::string &crsName) { return crsName == thisName ? 100 : metadata::Identifier::isEquivalentName( crsName.c_str(), thisName.c_str()) ? 90 : 70; }; const auto &conv = derivingConversionRef(); const auto &cs = coordinateSystem(); if (baseRes.size() == 1 && baseRes.front().second >= 70 && (authorityFactory == nullptr || authorityFactory->getAuthority().empty() || authorityFactory->getAuthority() == metadata::Identifier::EPSG) && conv->isUTM(zone, north) && cs->_isEquivalentTo( cs::CartesianCS::createEastingNorthing(common::UnitOfMeasure::METRE) .get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { auto computeUTMCRSName = [](const char *base, int l_zone, bool l_north) { return base + toString(l_zone) + (l_north ? "N" : "S"); }; if (baseRes.front().first->_isEquivalentTo( GeographicCRS::EPSG_4326.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { std::string crsName( computeUTMCRSName("WGS 84 / UTM zone ", zone, north)); res.emplace_back( ProjectedCRS::create( createMapNameEPSGCode(crsName.c_str(), (north ? 32600 : 32700) + zone), GeographicCRS::EPSG_4326, conv->identify(), cs), computeConfidence(crsName)); return res; } else if (((zone >= 1 && zone <= 22) || zone == 59 || zone == 60) && north && baseRes.front().first->_isEquivalentTo( GeographicCRS::EPSG_4267.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { std::string crsName( computeUTMCRSName("NAD27 / UTM zone ", zone, north)); res.emplace_back( ProjectedCRS::create( createMapNameEPSGCode(crsName.c_str(), (zone >= 59) ? 3370 + zone - 59 : 26700 + zone), GeographicCRS::EPSG_4267, conv->identify(), cs), computeConfidence(crsName)); return res; } else if (((zone >= 1 && zone <= 23) || zone == 59 || zone == 60) && north && baseRes.front().first->_isEquivalentTo( GeographicCRS::EPSG_4269.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { std::string crsName( computeUTMCRSName("NAD83 / UTM zone ", zone, north)); res.emplace_back( ProjectedCRS::create( createMapNameEPSGCode(crsName.c_str(), (zone >= 59) ? 3372 + zone - 59 : 26900 + zone), GeographicCRS::EPSG_4269, conv->identify(), cs), computeConfidence(crsName)); return res; } } const bool l_implicitCS = hasImplicitCS(); const auto addCRS = [&](const ProjectedCRSNNPtr &crs, const bool eqName) { const auto &l_unit = cs->axisList()[0]->unit(); if (_isEquivalentTo(crs.get(), util::IComparable::Criterion:: EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS, dbContext) || (l_implicitCS && l_unit._isEquivalentTo( crs->coordinateSystem()->axisList()[0]->unit(), util::IComparable::Criterion::EQUIVALENT) && l_baseCRS->_isEquivalentTo( crs->baseCRS().get(), util::IComparable::Criterion:: EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS, dbContext) && derivingConversionRef()->_isEquivalentTo( crs->derivingConversionRef().get(), util::IComparable::Criterion::EQUIVALENT, dbContext))) { if (crs->nameStr() == thisName) { res.clear(); res.emplace_back(crs, 100); } else { res.emplace_back(crs, eqName ? 90 : 70); } } else if (ellipsoid->_isEquivalentTo( crs->baseCRS()->ellipsoid().get(), util::IComparable::Criterion::EQUIVALENT, dbContext) && derivingConversionRef()->_isEquivalentTo( crs->derivingConversionRef().get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { if ((l_implicitCS && l_unit._isEquivalentTo( crs->coordinateSystem()->axisList()[0]->unit(), util::IComparable::Criterion::EQUIVALENT)) || cs->_isEquivalentTo(crs->coordinateSystem().get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { if (!significantNameForDatum || l_datum->_isEquivalentTo( crs->baseCRS()->datumNonNull(dbContext).get(), util::IComparable::Criterion::EQUIVALENT)) { res.emplace_back(crs, 70); } else { res.emplace_back(crs, 60); } } else { res.emplace_back(crs, 50); } } else { res.emplace_back(crs, 25); } }; if (authorityFactory) { const bool insignificantName = thisName.empty() || ci_equal(thisName, "unknown") || ci_equal(thisName, "unnamed"); bool foundEquivalentName = false; if (hasCodeCompatibleOfAuthorityFactory(this, authorityFactory)) { // If the CRS has already an id, check in the database for the // official object, and verify that they are equivalent. for (const auto &id : identifiers()) { if (hasCodeCompatibleOfAuthorityFactory(id, authorityFactory)) { try { auto crs = io::AuthorityFactory::create( authorityFactory->databaseContext(), *id->codeSpace()) ->createProjectedCRS(id->code()); bool match = _isEquivalentTo( crs.get(), util::IComparable::Criterion:: EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS, dbContext); res.emplace_back(crs, match ? 100 : 25); return res; } catch (const std::exception &) { } } } } else if (!insignificantName) { for (int ipass = 0; ipass < 2; ipass++) { const bool approximateMatch = ipass == 1; auto objects = authorityFactory->createObjectsFromNameEx( thisName, {io::AuthorityFactory::ObjectType::PROJECTED_CRS}, approximateMatch); for (const auto &pairObjName : objects) { auto crs = util::nn_dynamic_pointer_cast( pairObjName.first); assert(crs); auto crsNN = NN_NO_CHECK(crs); const bool eqName = metadata::Identifier::isEquivalentName( thisName.c_str(), pairObjName.second.c_str()); foundEquivalentName |= eqName; addCRS(crsNN, eqName); if (res.back().second == 100) { return res; } } if (!res.empty()) { break; } } } const auto lambdaSort = [&thisName](const Pair &a, const Pair &b) { // First consider confidence if (a.second > b.second) { return true; } if (a.second < b.second) { return false; } // Then consider exact name matching const auto &aName(a.first->nameStr()); const auto &bName(b.first->nameStr()); if (aName == thisName && bName != thisName) { return true; } if (bName == thisName && aName != thisName) { return false; } // Arbitrary final sorting criterion return aName < bName; }; // Sort results res.sort(lambdaSort); if (!hasCodeCompatibleOfAuthorityFactory(this, authorityFactory) && !foundEquivalentName && (res.empty() || res.front().second < 50)) { std::set> alreadyKnown; for (const auto &pair : res) { const auto &ids = pair.first->identifiers(); assert(!ids.empty()); alreadyKnown.insert(std::pair( *(ids[0]->codeSpace()), ids[0]->code())); } auto self = NN_NO_CHECK(std::dynamic_pointer_cast( shared_from_this().as_nullable())); auto candidates = authorityFactory->createProjectedCRSFromExisting(self); for (const auto &crs : candidates) { const auto &ids = crs->identifiers(); assert(!ids.empty()); if (alreadyKnown.find(std::pair( *(ids[0]->codeSpace()), ids[0]->code())) != alreadyKnown.end()) { continue; } addCRS(crs, insignificantName); } res.sort(lambdaSort); } // Keep only results of the highest confidence if (res.size() >= 2) { const auto highestConfidence = res.front().second; std::list newRes; for (const auto &pair : res) { if (pair.second == highestConfidence) { newRes.push_back(pair); } else { break; } } return newRes; } } return res; } // --------------------------------------------------------------------------- /** \brief Return a variant of this CRS "demoted" to a 2D one, if not already * the case. * * * @param newName Name of the new CRS. If empty, nameStr() will be used. * @param dbContext Database context to look for potentially already registered * 2D CRS. May be nullptr. * @return a new CRS demoted to 2D, or the current one if already 2D or not * applicable. * @since 6.3 */ ProjectedCRSNNPtr ProjectedCRS::demoteTo2D(const std::string &newName, const io::DatabaseContextPtr &dbContext) const { const auto &axisList = coordinateSystem()->axisList(); if (axisList.size() == 3) { auto cs = cs::CartesianCS::create(util::PropertyMap(), axisList[0], axisList[1]); const auto &l_baseCRS = baseCRS(); const auto geogCRS = dynamic_cast(l_baseCRS.get()); const auto newBaseCRS = geogCRS ? util::nn_static_pointer_cast( geogCRS->demoteTo2D(std::string(), dbContext)) : l_baseCRS; return ProjectedCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, !newName.empty() ? newName : nameStr()), newBaseCRS, derivingConversion(), cs); } return NN_NO_CHECK(std::dynamic_pointer_cast( shared_from_this().as_nullable())); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> ProjectedCRS::_identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; auto resTemp = identify(authorityFactory); for (const auto &pair : resTemp) { res.emplace_back(pair.first, pair.second); } return res; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress InvalidCompoundCRSException::InvalidCompoundCRSException(const char *message) : Exception(message) {} // --------------------------------------------------------------------------- InvalidCompoundCRSException::InvalidCompoundCRSException( const std::string &message) : Exception(message) {} // --------------------------------------------------------------------------- InvalidCompoundCRSException::~InvalidCompoundCRSException() = default; // --------------------------------------------------------------------------- InvalidCompoundCRSException::InvalidCompoundCRSException( const InvalidCompoundCRSException &) = default; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct CompoundCRS::Private { std::vector components_{}; }; //! @endcond // --------------------------------------------------------------------------- CompoundCRS::CompoundCRS(const std::vector &components) : CRS(), d(internal::make_unique()) { d->components_ = components; } // --------------------------------------------------------------------------- CompoundCRS::CompoundCRS(const CompoundCRS &other) : CRS(other), d(internal::make_unique(*other.d)) {} // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress CompoundCRS::~CompoundCRS() = default; //! @endcond // --------------------------------------------------------------------------- CRSNNPtr CompoundCRS::_shallowClone() const { auto crs(CompoundCRS::nn_make_shared(*this)); crs->assignSelf(crs); return crs; } // --------------------------------------------------------------------------- /** \brief Return the components of a CompoundCRS. * * @return the components. */ const std::vector & CompoundCRS::componentReferenceSystems() PROJ_PURE_DEFN { return d->components_; } // --------------------------------------------------------------------------- /** \brief Instantiate a CompoundCRS from a vector of CRS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param components the component CRS of the CompoundCRS. * @return new CompoundCRS. * @throw InvalidCompoundCRSException */ CompoundCRSNNPtr CompoundCRS::create(const util::PropertyMap &properties, const std::vector &components) { if (components.size() < 2) { throw InvalidCompoundCRSException( "compound CRS should have at least 2 components"); } auto comp0 = components[0].get(); auto comp0Bound = dynamic_cast(comp0); if (comp0Bound) { comp0 = comp0Bound->baseCRS().get(); } auto comp0Geog = dynamic_cast(comp0); auto comp0Proj = dynamic_cast(comp0); auto comp0Eng = dynamic_cast(comp0); auto comp1 = components[1].get(); auto comp1Bound = dynamic_cast(comp1); if (comp1Bound) { comp1 = comp1Bound->baseCRS().get(); } auto comp1Vert = dynamic_cast(comp1); auto comp1Eng = dynamic_cast(comp1); // Loose validation based on // http://docs.opengeospatial.org/as/18-005r4/18-005r4.html#34 bool ok = false; if ((comp0Geog && comp0Geog->coordinateSystem()->axisList().size() == 2 && (comp1Vert || (comp1Eng && comp1Eng->coordinateSystem()->axisList().size() == 1))) || (comp0Proj && comp0Proj->coordinateSystem()->axisList().size() == 2 && (comp1Vert || (comp1Eng && comp1Eng->coordinateSystem()->axisList().size() == 1))) || (comp0Eng && comp0Eng->coordinateSystem()->axisList().size() <= 2 && comp1Vert)) { // Spatial compound coordinate reference system ok = true; } else { bool isComp0Spatial = comp0Geog || comp0Proj || comp0Eng || dynamic_cast(comp0) || dynamic_cast(comp0); if (isComp0Spatial && dynamic_cast(comp1)) { // Spatio-temporal compound coordinate reference system ok = true; } else if (isComp0Spatial && dynamic_cast(comp1)) { // Spatio-parametric compound coordinate reference system ok = true; } } if (!ok) { throw InvalidCompoundCRSException( "components of the compound CRS do not belong to one of the " "allowed combinations of " "http://docs.opengeospatial.org/as/18-005r4/18-005r4.html#34"); } auto compoundCRS(CompoundCRS::nn_make_shared(components)); compoundCRS->assignSelf(compoundCRS); compoundCRS->setProperties(properties); if (!properties.get(common::IdentifiedObject::NAME_KEY)) { std::string name; for (const auto &crs : components) { if (!name.empty()) { name += " + "; } const auto &l_name = crs->nameStr(); if (!l_name.empty()) { name += l_name; } else { name += "unnamed"; } } util::PropertyMap propertyName; propertyName.set(common::IdentifiedObject::NAME_KEY, name); compoundCRS->setProperties(propertyName); } return compoundCRS; } // --------------------------------------------------------------------------- //! @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 &components, const io::DatabaseContextPtr &dbContext) { if (components.size() == 2) { auto comp0 = components[0].get(); auto comp1 = components[1].get(); auto comp0Geog = dynamic_cast(comp0); auto comp0Proj = dynamic_cast(comp0); auto comp0Bound = dynamic_cast(comp0); if (comp0Geog == nullptr && comp0Proj == nullptr) { if (comp0Bound) { const auto *baseCRS = comp0Bound->baseCRS().get(); comp0Geog = dynamic_cast(baseCRS); comp0Proj = dynamic_cast(baseCRS); } } auto comp1Geog = dynamic_cast(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"); } // Detect a COMPD_CS whose VERT_CS is for ellipoidal heights auto comp1Vert = util::nn_dynamic_pointer_cast(components[1]); if (comp1Vert != nullptr && comp1Vert->datum() && comp1Vert->datum()->getWKT1DatumType() == "2002") { const auto &axis = comp1Vert->coordinateSystem()->axisList()[0]; std::string name(components[0]->nameStr()); if (!(axis->unit()._isEquivalentTo( common::UnitOfMeasure::METRE, util::IComparable::Criterion::EQUIVALENT) && &(axis->direction()) == &(cs::AxisDirection::UP))) { name += " (" + comp1Vert->nameStr() + ')'; } auto newVertAxis = cs::CoordinateSystemAxis::create( util::PropertyMap().set(IdentifiedObject::NAME_KEY, cs::AxisName::Ellipsoidal_height), cs::AxisAbbreviation::h, axis->direction(), axis->unit()); return components[0] ->promoteTo3D(name, dbContext, newVertAxis) ->attachOriginalCompoundCRS(create( properties, comp0Bound ? std::vector{comp0Bound->baseCRS(), components[1]} : components)); } } return create(properties, components); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void CompoundCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; const auto &l_components = componentReferenceSystems(); if (!isWKT2 && formatter->useESRIDialect() && l_components.size() == 2) { l_components[0]->_exportToWKT(formatter); l_components[1]->_exportToWKT(formatter); } else { formatter->startNode(isWKT2 ? io::WKTConstants::COMPOUNDCRS : io::WKTConstants::COMPD_CS, !identifiers().empty()); formatter->addQuotedString(nameStr()); for (const auto &crs : l_components) { crs->_exportToWKT(formatter); } ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void CompoundCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("CompoundCRS", !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } writer->AddObjKey("components"); { auto componentsContext(writer->MakeArrayContext(false)); for (const auto &crs : componentReferenceSystems()) { crs->_exportToJSON(formatter); } } ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- void CompoundCRS::_exportToPROJString( io::PROJStringFormatter *formatter) const // throw(io::FormattingException) { for (const auto &crs : componentReferenceSystems()) { auto crs_exportable = dynamic_cast(crs.get()); if (crs_exportable) { crs_exportable->_exportToPROJString(formatter); } } } // --------------------------------------------------------------------------- bool CompoundCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherCompoundCRS = dynamic_cast(other); if (otherCompoundCRS == nullptr || (criterion == util::IComparable::Criterion::STRICT && !ObjectUsage::_isEquivalentTo(other, criterion, dbContext))) { return false; } const auto &components = componentReferenceSystems(); const auto &otherComponents = otherCompoundCRS->componentReferenceSystems(); if (components.size() != otherComponents.size()) { return false; } for (size_t i = 0; i < components.size(); i++) { if (!components[i]->_isEquivalentTo(otherComponents[i].get(), criterion, dbContext)) { return false; } } return true; } // --------------------------------------------------------------------------- /** \brief Identify the CRS with reference CRSs. * * The candidate CRSs are looked in the database when * authorityFactory is not null. * * Note that the implementation uses a set of heuristics to have a good * compromise of successful identifications over execution time. It might miss * legitimate matches in some circumstances. * * The method returns a list of matching reference CRS, and the percentage * (0-100) of confidence in the match. The list is sorted by decreasing * confidence. * * 100% means that the name of the reference entry * perfectly matches the CRS name, and both are equivalent. In which case a * single result is returned. * 90% means that CRS are equivalent, but the names are not exactly the same. * 70% means that CRS are equivalent (equivalent horizontal and vertical CRS), * but the names are not equivalent. * 25% means that the CRS are not equivalent, but there is some similarity in * the names. * * @param authorityFactory Authority factory (if null, will return an empty * list) * @return a list of matching reference CRS, and the percentage (0-100) of * confidence in the match. */ std::list> CompoundCRS::identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; const auto &thisName(nameStr()); const auto &components = componentReferenceSystems(); bool l_implicitCS = components[0]->hasImplicitCS(); if (!l_implicitCS) { const auto projCRS = dynamic_cast(components[0].get()); if (projCRS) { l_implicitCS = projCRS->baseCRS()->hasImplicitCS(); } } const auto crsCriterion = l_implicitCS ? util::IComparable::Criterion::EQUIVALENT_EXCEPT_AXIS_ORDER_GEOGCRS : util::IComparable::Criterion::EQUIVALENT; if (authorityFactory) { const io::DatabaseContextNNPtr &dbContext = authorityFactory->databaseContext(); const bool insignificantName = thisName.empty() || ci_equal(thisName, "unknown") || ci_equal(thisName, "unnamed"); bool foundEquivalentName = false; if (hasCodeCompatibleOfAuthorityFactory(this, authorityFactory)) { // If the CRS has already an id, check in the database for the // official object, and verify that they are equivalent. for (const auto &id : identifiers()) { if (hasCodeCompatibleOfAuthorityFactory(id, authorityFactory)) { try { auto crs = io::AuthorityFactory::create( dbContext, *id->codeSpace()) ->createCompoundCRS(id->code()); bool match = _isEquivalentTo(crs.get(), crsCriterion, dbContext); res.emplace_back(crs, match ? 100 : 25); return res; } catch (const std::exception &) { } } } } else if (!insignificantName) { for (int ipass = 0; ipass < 2; ipass++) { const bool approximateMatch = ipass == 1; auto objects = authorityFactory->createObjectsFromName( thisName, {io::AuthorityFactory::ObjectType::COMPOUND_CRS}, approximateMatch); for (const auto &obj : objects) { auto crs = util::nn_dynamic_pointer_cast(obj); assert(crs); auto crsNN = NN_NO_CHECK(crs); const bool eqName = metadata::Identifier::isEquivalentName( thisName.c_str(), crs->nameStr().c_str()); foundEquivalentName |= eqName; if (_isEquivalentTo(crs.get(), crsCriterion, dbContext)) { if (crs->nameStr() == thisName) { res.clear(); res.emplace_back(crsNN, 100); return res; } res.emplace_back(crsNN, eqName ? 90 : 70); } else { res.emplace_back(crsNN, 25); } } if (!res.empty()) { break; } } } const auto lambdaSort = [&thisName](const Pair &a, const Pair &b) { // First consider confidence if (a.second > b.second) { return true; } if (a.second < b.second) { return false; } // Then consider exact name matching const auto &aName(a.first->nameStr()); const auto &bName(b.first->nameStr()); if (aName == thisName && bName != thisName) { return true; } if (bName == thisName && aName != thisName) { return false; } // Arbitrary final sorting criterion return aName < bName; }; // Sort results res.sort(lambdaSort); if (identifiers().empty() && !foundEquivalentName && (res.empty() || res.front().second < 50)) { std::set> alreadyKnown; for (const auto &pair : res) { const auto &ids = pair.first->identifiers(); assert(!ids.empty()); alreadyKnown.insert(std::pair( *(ids[0]->codeSpace()), ids[0]->code())); } auto self = NN_NO_CHECK(std::dynamic_pointer_cast( shared_from_this().as_nullable())); auto candidates = authorityFactory->createCompoundCRSFromExisting(self); for (const auto &crs : candidates) { const auto &ids = crs->identifiers(); assert(!ids.empty()); if (alreadyKnown.find(std::pair( *(ids[0]->codeSpace()), ids[0]->code())) != alreadyKnown.end()) { continue; } if (_isEquivalentTo(crs.get(), crsCriterion, dbContext)) { res.emplace_back(crs, insignificantName ? 90 : 70); } else { res.emplace_back(crs, 25); } } res.sort(lambdaSort); // If there's a single candidate at 90% confidence with same name, // then promote it to 100% if (res.size() == 1 && res.front().second == 90 && thisName == res.front().first->nameStr()) { res.front().second = 100; } } // If we didn't find a match for the CompoundCRS, check if the // horizontal and vertical parts are not themselves well known. if (identifiers().empty() && res.empty() && components.size() == 2) { auto candidatesHorizCRS = components[0]->identify(authorityFactory); auto candidatesVertCRS = components[1]->identify(authorityFactory); if (candidatesHorizCRS.size() == 1 && candidatesVertCRS.size() == 1 && candidatesHorizCRS.front().second >= 70 && candidatesVertCRS.front().second >= 70) { auto newCRS = CompoundCRS::create( util::PropertyMap().set( common::IdentifiedObject::NAME_KEY, candidatesHorizCRS.front().first->nameStr() + " + " + candidatesVertCRS.front().first->nameStr()), {candidatesHorizCRS.front().first, candidatesVertCRS.front().first}); const bool eqName = metadata::Identifier::isEquivalentName( thisName.c_str(), newCRS->nameStr().c_str()); res.emplace_back( newCRS, std::min(thisName == newCRS->nameStr() ? 100 : eqName ? 90 : 70, std::min(candidatesHorizCRS.front().second, candidatesVertCRS.front().second))); } } // Keep only results of the highest confidence if (res.size() >= 2) { const auto highestConfidence = res.front().second; std::list newRes; for (const auto &pair : res) { if (pair.second == highestConfidence) { newRes.push_back(pair); } else { break; } } return newRes; } } return res; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> CompoundCRS::_identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; auto resTemp = identify(authorityFactory); for (const auto &pair : resTemp) { res.emplace_back(pair.first, pair.second); } return res; } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct PROJ_INTERNAL BoundCRS::Private { CRSNNPtr baseCRS_; CRSNNPtr hubCRS_; operation::TransformationNNPtr transformation_; Private(const CRSNNPtr &baseCRSIn, const CRSNNPtr &hubCRSIn, const operation::TransformationNNPtr &transformationIn); inline const CRSNNPtr &baseCRS() const { return baseCRS_; } inline const CRSNNPtr &hubCRS() const { return hubCRS_; } inline const operation::TransformationNNPtr &transformation() const { return transformation_; } }; BoundCRS::Private::Private( const CRSNNPtr &baseCRSIn, const CRSNNPtr &hubCRSIn, const operation::TransformationNNPtr &transformationIn) : baseCRS_(baseCRSIn), hubCRS_(hubCRSIn), transformation_(transformationIn) {} //! @endcond // --------------------------------------------------------------------------- BoundCRS::BoundCRS(const CRSNNPtr &baseCRSIn, const CRSNNPtr &hubCRSIn, const operation::TransformationNNPtr &transformationIn) : d(internal::make_unique(baseCRSIn, hubCRSIn, transformationIn)) { } // --------------------------------------------------------------------------- BoundCRS::BoundCRS(const BoundCRS &other) : CRS(other), d(internal::make_unique(other.d->baseCRS(), other.d->hubCRS(), other.d->transformation())) {} // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress BoundCRS::~BoundCRS() = default; //! @endcond // --------------------------------------------------------------------------- BoundCRSNNPtr BoundCRS::shallowCloneAsBoundCRS() const { auto crs(BoundCRS::nn_make_shared(*this)); crs->assignSelf(crs); return crs; } // --------------------------------------------------------------------------- CRSNNPtr BoundCRS::_shallowClone() const { return shallowCloneAsBoundCRS(); } // --------------------------------------------------------------------------- /** \brief Return the base CRS. * * This is the CRS into which coordinates of the BoundCRS are expressed. * * @return the base CRS. */ const CRSNNPtr &BoundCRS::baseCRS() PROJ_PURE_DEFN { return d->baseCRS_; } // --------------------------------------------------------------------------- // The only legit caller is BoundCRS::baseCRSWithCanonicalBoundCRS() void CRS::setCanonicalBoundCRS(const BoundCRSNNPtr &boundCRS) { d->canonicalBoundCRS_ = boundCRS; } // --------------------------------------------------------------------------- /** \brief Return a shallow clone of the base CRS that points to a * shallow clone of this BoundCRS. * * The base CRS is the CRS into which coordinates of the BoundCRS are expressed. * * The returned CRS will actually be a shallow clone of the actual base CRS, * with the extra property that CRS::canonicalBoundCRS() will point to a * shallow clone of this BoundCRS. Use this only if you want to work with * the base CRS object rather than the BoundCRS, but wanting to be able to * retrieve the BoundCRS later. * * @return the base CRS. */ CRSNNPtr BoundCRS::baseCRSWithCanonicalBoundCRS() const { auto baseCRSClone = baseCRS()->_shallowClone(); baseCRSClone->setCanonicalBoundCRS(shallowCloneAsBoundCRS()); return baseCRSClone; } // --------------------------------------------------------------------------- /** \brief Return the target / hub CRS. * * @return the hub CRS. */ const CRSNNPtr &BoundCRS::hubCRS() PROJ_PURE_DEFN { return d->hubCRS_; } // --------------------------------------------------------------------------- /** \brief Return the transformation to the hub RS. * * @return transformation. */ const operation::TransformationNNPtr & BoundCRS::transformation() PROJ_PURE_DEFN { return d->transformation_; } // --------------------------------------------------------------------------- /** \brief Instantiate a BoundCRS from a base CRS, a hub CRS and a * transformation. * * @param baseCRSIn base CRS. * @param hubCRSIn hub CRS. * @param transformationIn transformation from base CRS to hub CRS. * @return new BoundCRS. */ BoundCRSNNPtr BoundCRS::create(const CRSNNPtr &baseCRSIn, const CRSNNPtr &hubCRSIn, const operation::TransformationNNPtr &transformationIn) { auto crs = BoundCRS::nn_make_shared(baseCRSIn, hubCRSIn, transformationIn); crs->assignSelf(crs); const auto &l_name = baseCRSIn->nameStr(); if (!l_name.empty()) { crs->setProperties(util::PropertyMap().set( common::IdentifiedObject::NAME_KEY, l_name)); } return crs; } // --------------------------------------------------------------------------- /** \brief Instantiate a BoundCRS from a base CRS and TOWGS84 parameters * * @param baseCRSIn base CRS. * @param TOWGS84Parameters a vector of 3 or 7 double values representing WKT1 * TOWGS84 parameter. * @return new BoundCRS. */ BoundCRSNNPtr BoundCRS::createFromTOWGS84(const CRSNNPtr &baseCRSIn, const std::vector &TOWGS84Parameters) { auto transf = operation::Transformation::createTOWGS84(baseCRSIn, TOWGS84Parameters); return create(baseCRSIn, transf->targetCRS(), transf); } // --------------------------------------------------------------------------- /** \brief Instantiate a BoundCRS from a base CRS and nadgrids parameters * * @param baseCRSIn base CRS. * @param filename Horizontal grid filename * @return new BoundCRS. */ BoundCRSNNPtr BoundCRS::createFromNadgrids(const CRSNNPtr &baseCRSIn, const std::string &filename) { const auto sourceGeographicCRS = baseCRSIn->extractGeographicCRS(); auto transformationSourceCRS = sourceGeographicCRS ? NN_NO_CHECK(std::static_pointer_cast(sourceGeographicCRS)) : baseCRSIn; if (sourceGeographicCRS != nullptr && sourceGeographicCRS->primeMeridian()->longitude().getSIValue() != 0.0) { transformationSourceCRS = GeographicCRS::create( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, sourceGeographicCRS->nameStr() + " (with Greenwich prime meridian)"), datum::GeodeticReferenceFrame::create( util::PropertyMap().set( common::IdentifiedObject::NAME_KEY, sourceGeographicCRS->datumNonNull(nullptr)->nameStr() + " (with Greenwich prime meridian)"), sourceGeographicCRS->datumNonNull(nullptr)->ellipsoid(), util::optional(), datum::PrimeMeridian::GREENWICH), cs::EllipsoidalCS::createLatitudeLongitude( common::UnitOfMeasure::DEGREE)); } std::string transformationName = transformationSourceCRS->nameStr(); transformationName += " to WGS84"; return create( baseCRSIn, GeographicCRS::EPSG_4326, operation::Transformation::createNTv2( util::PropertyMap().set(common::IdentifiedObject::NAME_KEY, transformationName), transformationSourceCRS, GeographicCRS::EPSG_4326, filename, std::vector())); } // --------------------------------------------------------------------------- bool BoundCRS::isTOWGS84Compatible() const { return dynamic_cast(d->hubCRS().get()) != nullptr && ci_equal(d->hubCRS()->nameStr(), "WGS 84"); } // --------------------------------------------------------------------------- std::string BoundCRS::getHDatumPROJ4GRIDS() const { if (ci_equal(d->hubCRS()->nameStr(), "WGS 84")) { return d->transformation()->getNTv2Filename(); } return std::string(); } // --------------------------------------------------------------------------- std::string BoundCRS::getVDatumPROJ4GRIDS() const { if (dynamic_cast(d->baseCRS().get()) && ci_equal(d->hubCRS()->nameStr(), "WGS 84")) { return d->transformation()->getHeightToGeographic3DFilename(); } return std::string(); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void BoundCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; if (isWKT2) { formatter->startNode(io::WKTConstants::BOUNDCRS, false); formatter->startNode(io::WKTConstants::SOURCECRS, false); d->baseCRS()->_exportToWKT(formatter); formatter->endNode(); formatter->startNode(io::WKTConstants::TARGETCRS, false); d->hubCRS()->_exportToWKT(formatter); formatter->endNode(); formatter->setAbridgedTransformation(true); d->transformation()->_exportToWKT(formatter); formatter->setAbridgedTransformation(false); formatter->endNode(); } else { auto vdatumProj4GridName = getVDatumPROJ4GRIDS(); if (!vdatumProj4GridName.empty()) { formatter->setVDatumExtension(vdatumProj4GridName); d->baseCRS()->_exportToWKT(formatter); formatter->setVDatumExtension(std::string()); return; } auto hdatumProj4GridName = getHDatumPROJ4GRIDS(); if (!hdatumProj4GridName.empty()) { formatter->setHDatumExtension(hdatumProj4GridName); d->baseCRS()->_exportToWKT(formatter); formatter->setHDatumExtension(std::string()); return; } if (!isTOWGS84Compatible()) { io::FormattingException::Throw( "Cannot export BoundCRS with non-WGS 84 hub CRS in WKT1"); } auto params = d->transformation()->getTOWGS84Parameters(); if (!formatter->useESRIDialect()) { formatter->setTOWGS84Parameters(params); } d->baseCRS()->_exportToWKT(formatter); formatter->setTOWGS84Parameters(std::vector()); } } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void BoundCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("BoundCRS", !identifiers().empty())); writer->AddObjKey("source_crs"); d->baseCRS()->_exportToJSON(formatter); writer->AddObjKey("target_crs"); d->hubCRS()->_exportToJSON(formatter); writer->AddObjKey("transformation"); formatter->setOmitTypeInImmediateChild(); formatter->setAbridgedTransformation(true); d->transformation()->_exportToJSON(formatter); formatter->setAbridgedTransformation(false); } //! @endcond // --------------------------------------------------------------------------- void BoundCRS::_exportToPROJString( io::PROJStringFormatter *formatter) const // throw(io::FormattingException) { auto crs_exportable = dynamic_cast(d->baseCRS_.get()); if (!crs_exportable) { io::FormattingException::Throw( "baseCRS of BoundCRS cannot be exported as a PROJ string"); } auto vdatumProj4GridName = getVDatumPROJ4GRIDS(); if (!vdatumProj4GridName.empty()) { formatter->setVDatumExtension(vdatumProj4GridName); crs_exportable->_exportToPROJString(formatter); formatter->setVDatumExtension(std::string()); } else { auto hdatumProj4GridName = getHDatumPROJ4GRIDS(); if (!hdatumProj4GridName.empty()) { formatter->setHDatumExtension(hdatumProj4GridName); crs_exportable->_exportToPROJString(formatter); formatter->setHDatumExtension(std::string()); } else { if (isTOWGS84Compatible()) { auto params = transformation()->getTOWGS84Parameters(); formatter->setTOWGS84Parameters(params); } crs_exportable->_exportToPROJString(formatter); formatter->setTOWGS84Parameters(std::vector()); } } } // --------------------------------------------------------------------------- bool BoundCRS::_isEquivalentTo(const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherBoundCRS = dynamic_cast(other); if (otherBoundCRS == nullptr || (criterion == util::IComparable::Criterion::STRICT && !ObjectUsage::_isEquivalentTo(other, criterion, dbContext))) { return false; } const auto standardCriterion = getStandardCriterion(criterion); return d->baseCRS_->_isEquivalentTo(otherBoundCRS->d->baseCRS_.get(), criterion, dbContext) && d->hubCRS_->_isEquivalentTo(otherBoundCRS->d->hubCRS_.get(), criterion, dbContext) && d->transformation_->_isEquivalentTo( otherBoundCRS->d->transformation_.get(), standardCriterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> BoundCRS::_identify(const io::AuthorityFactoryPtr &authorityFactory) const { typedef std::pair Pair; std::list res; if (!authorityFactory) return res; std::list resMatchOfTransfToWGS84; const io::DatabaseContextNNPtr &dbContext = authorityFactory->databaseContext(); if (d->hubCRS_->_isEquivalentTo(GeographicCRS::EPSG_4326.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { auto resTemp = d->baseCRS_->identify(authorityFactory); std::string refTransfPROJString; bool refTransfPROJStringValid = false; auto refTransf = d->transformation_->normalizeForVisualization(); try { refTransfPROJString = refTransf->exportToPROJString( io::PROJStringFormatter::create().get()); refTransfPROJString = replaceAll( refTransfPROJString, " +rx=0 +ry=0 +rz=0 +s=0 +convention=position_vector", ""); refTransfPROJStringValid = true; } catch (const std::exception &) { } bool refIsNullTransform = false; if (isTOWGS84Compatible()) { auto params = transformation()->getTOWGS84Parameters(); if (params == std::vector{0, 0, 0, 0, 0, 0, 0}) { refIsNullTransform = true; } } for (const auto &pair : resTemp) { const auto &candidateBaseCRS = pair.first; auto projCRS = dynamic_cast(candidateBaseCRS.get()); auto geodCRS = projCRS ? projCRS->baseCRS().as_nullable() : util::nn_dynamic_pointer_cast( candidateBaseCRS); if (geodCRS) { auto context = operation::CoordinateOperationContext::create( authorityFactory, nullptr, 0.0); context->setSpatialCriterion( operation::CoordinateOperationContext::SpatialCriterion:: PARTIAL_INTERSECTION); auto ops = operation::CoordinateOperationFactory::create() ->createOperations(NN_NO_CHECK(geodCRS), GeographicCRS::EPSG_4326, context); bool foundOp = false; for (const auto &op : ops) { auto opNormalized = op->normalizeForVisualization(); std::string opTransfPROJString; bool opTransfPROJStringValid = false; const auto &opName = op->nameStr(); if (starts_with( opName, operation::BALLPARK_GEOCENTRIC_TRANSLATION) || starts_with(opName, operation::NULL_GEOGRAPHIC_OFFSET)) { if (refIsNullTransform) { res.emplace_back(create(candidateBaseCRS, d->hubCRS_, transformation()), pair.second); foundOp = true; break; } continue; } try { opTransfPROJString = opNormalized->exportToPROJString( io::PROJStringFormatter::create().get()); opTransfPROJStringValid = true; opTransfPROJString = replaceAll(opTransfPROJString, " +rx=0 +ry=0 +rz=0 +s=0 " "+convention=position_vector", ""); } catch (const std::exception &) { } if ((refTransfPROJStringValid && opTransfPROJStringValid && refTransfPROJString == opTransfPROJString) || opNormalized->_isEquivalentTo( refTransf.get(), util::IComparable::Criterion::EQUIVALENT, dbContext)) { resMatchOfTransfToWGS84.emplace_back( create(candidateBaseCRS, d->hubCRS_, NN_NO_CHECK(util::nn_dynamic_pointer_cast< operation::Transformation>(op))), pair.second); foundOp = true; break; } } if (!foundOp) { res.emplace_back( create(candidateBaseCRS, d->hubCRS_, transformation()), std::min(70, pair.second)); } } } } return !resMatchOfTransfToWGS84.empty() ? resMatchOfTransfToWGS84 : res; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct DerivedGeodeticCRS::Private {}; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress DerivedGeodeticCRS::~DerivedGeodeticCRS() = default; //! @endcond // --------------------------------------------------------------------------- DerivedGeodeticCRS::DerivedGeodeticCRS( const GeodeticCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::CartesianCSNNPtr &csIn) : SingleCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), GeodeticCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), DerivedCRS(baseCRSIn, derivingConversionIn, csIn), d(nullptr) {} // --------------------------------------------------------------------------- DerivedGeodeticCRS::DerivedGeodeticCRS( const GeodeticCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::SphericalCSNNPtr &csIn) : SingleCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), GeodeticCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), DerivedCRS(baseCRSIn, derivingConversionIn, csIn), d(nullptr) {} // --------------------------------------------------------------------------- DerivedGeodeticCRS::DerivedGeodeticCRS(const DerivedGeodeticCRS &other) : SingleCRS(other), GeodeticCRS(other), DerivedCRS(other), d(nullptr) {} // --------------------------------------------------------------------------- CRSNNPtr DerivedGeodeticCRS::_shallowClone() const { auto crs(DerivedGeodeticCRS::nn_make_shared(*this)); crs->assignSelf(crs); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- /** \brief Return the base CRS (a GeodeticCRS) of a DerivedGeodeticCRS. * * @return the base CRS. */ const GeodeticCRSNNPtr DerivedGeodeticCRS::baseCRS() const { return NN_NO_CHECK(util::nn_dynamic_pointer_cast( DerivedCRS::getPrivate()->baseCRS_)); } // --------------------------------------------------------------------------- /** \brief Instantiate a DerivedGeodeticCRS from a base CRS, a deriving * conversion and a cs::CartesianCS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param baseCRSIn base CRS. * @param derivingConversionIn the deriving conversion from the base CRS to this * CRS. * @param csIn the coordinate system. * @return new DerivedGeodeticCRS. */ DerivedGeodeticCRSNNPtr DerivedGeodeticCRS::create( const util::PropertyMap &properties, const GeodeticCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::CartesianCSNNPtr &csIn) { auto crs(DerivedGeodeticCRS::nn_make_shared( baseCRSIn, derivingConversionIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- /** \brief Instantiate a DerivedGeodeticCRS from a base CRS, a deriving * conversion and a cs::SphericalCS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param baseCRSIn base CRS. * @param derivingConversionIn the deriving conversion from the base CRS to this * CRS. * @param csIn the coordinate system. * @return new DerivedGeodeticCRS. */ DerivedGeodeticCRSNNPtr DerivedGeodeticCRS::create( const util::PropertyMap &properties, const GeodeticCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::SphericalCSNNPtr &csIn) { auto crs(DerivedGeodeticCRS::nn_make_shared( baseCRSIn, derivingConversionIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void DerivedGeodeticCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; if (!isWKT2) { io::FormattingException::Throw( "DerivedGeodeticCRS can only be exported to WKT2"); } formatter->startNode(io::WKTConstants::GEODCRS, !identifiers().empty()); formatter->addQuotedString(nameStr()); auto l_baseCRS = baseCRS(); formatter->startNode((formatter->use2019Keywords() && dynamic_cast(l_baseCRS.get())) ? io::WKTConstants::BASEGEOGCRS : io::WKTConstants::BASEGEODCRS, !baseCRS()->identifiers().empty()); formatter->addQuotedString(l_baseCRS->nameStr()); auto l_datum = l_baseCRS->datum(); if (l_datum) { l_datum->_exportToWKT(formatter); } else { auto l_datumEnsemble = datumEnsemble(); assert(l_datumEnsemble); l_datumEnsemble->_exportToWKT(formatter); } l_baseCRS->primeMeridian()->_exportToWKT(formatter); formatter->endNode(); formatter->setUseDerivingConversion(true); derivingConversionRef()->_exportToWKT(formatter); formatter->setUseDerivingConversion(false); coordinateSystem()->_exportToWKT(formatter); ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } //! @endcond // --------------------------------------------------------------------------- void DerivedGeodeticCRS::_exportToPROJString( io::PROJStringFormatter *) const // throw(io::FormattingException) { throw io::FormattingException( "DerivedGeodeticCRS cannot be exported to PROJ string"); } // --------------------------------------------------------------------------- bool DerivedGeodeticCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherDerivedCRS = dynamic_cast(other); return otherDerivedCRS != nullptr && DerivedCRS::_isEquivalentTo(other, criterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> DerivedGeodeticCRS::_identify(const io::AuthorityFactoryPtr &factory) const { return CRS::_identify(factory); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct DerivedGeographicCRS::Private {}; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress DerivedGeographicCRS::~DerivedGeographicCRS() = default; //! @endcond // --------------------------------------------------------------------------- DerivedGeographicCRS::DerivedGeographicCRS( const GeodeticCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::EllipsoidalCSNNPtr &csIn) : SingleCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), GeographicCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), DerivedCRS(baseCRSIn, derivingConversionIn, csIn), d(nullptr) {} // --------------------------------------------------------------------------- DerivedGeographicCRS::DerivedGeographicCRS(const DerivedGeographicCRS &other) : SingleCRS(other), GeographicCRS(other), DerivedCRS(other), d(nullptr) {} // --------------------------------------------------------------------------- CRSNNPtr DerivedGeographicCRS::_shallowClone() const { auto crs(DerivedGeographicCRS::nn_make_shared(*this)); crs->assignSelf(crs); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- /** \brief Return the base CRS (a GeodeticCRS) of a DerivedGeographicCRS. * * @return the base CRS. */ const GeodeticCRSNNPtr DerivedGeographicCRS::baseCRS() const { return NN_NO_CHECK(util::nn_dynamic_pointer_cast( DerivedCRS::getPrivate()->baseCRS_)); } // --------------------------------------------------------------------------- /** \brief Instantiate a DerivedGeographicCRS from a base CRS, a deriving * conversion and a cs::EllipsoidalCS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param baseCRSIn base CRS. * @param derivingConversionIn the deriving conversion from the base CRS to this * CRS. * @param csIn the coordinate system. * @return new DerivedGeographicCRS. */ DerivedGeographicCRSNNPtr DerivedGeographicCRS::create( const util::PropertyMap &properties, const GeodeticCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::EllipsoidalCSNNPtr &csIn) { auto crs(DerivedGeographicCRS::nn_make_shared( baseCRSIn, derivingConversionIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void DerivedGeographicCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; if (!isWKT2) { io::FormattingException::Throw( "DerivedGeographicCRS can only be exported to WKT2"); } formatter->startNode(formatter->use2019Keywords() ? io::WKTConstants::GEOGCRS : io::WKTConstants::GEODCRS, !identifiers().empty()); formatter->addQuotedString(nameStr()); auto l_baseCRS = baseCRS(); formatter->startNode((formatter->use2019Keywords() && dynamic_cast(l_baseCRS.get())) ? io::WKTConstants::BASEGEOGCRS : io::WKTConstants::BASEGEODCRS, !l_baseCRS->identifiers().empty()); formatter->addQuotedString(l_baseCRS->nameStr()); l_baseCRS->exportDatumOrDatumEnsembleToWkt(formatter); l_baseCRS->primeMeridian()->_exportToWKT(formatter); formatter->endNode(); formatter->setUseDerivingConversion(true); derivingConversionRef()->_exportToWKT(formatter); formatter->setUseDerivingConversion(false); coordinateSystem()->_exportToWKT(formatter); ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } //! @endcond // --------------------------------------------------------------------------- void DerivedGeographicCRS::_exportToPROJString( io::PROJStringFormatter *formatter) const // throw(io::FormattingException) { const auto &l_conv = derivingConversionRef(); const auto &methodName = l_conv->method()->nameStr(); for (const char *substr : {"PROJ ob_tran o_proj=longlat", "PROJ ob_tran o_proj=lonlat", "PROJ ob_tran o_proj=latlon", "PROJ ob_tran o_proj=latlong"}) { if (starts_with(methodName, substr)) { l_conv->_exportToPROJString(formatter); return; } } if (ci_equal(methodName, PROJ_WKT2_NAME_METHOD_POLE_ROTATION_GRIB_CONVENTION)) { l_conv->_exportToPROJString(formatter); return; } throw io::FormattingException( "DerivedGeographicCRS cannot be exported to PROJ string"); } // --------------------------------------------------------------------------- bool DerivedGeographicCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherDerivedCRS = dynamic_cast(other); return otherDerivedCRS != nullptr && DerivedCRS::_isEquivalentTo(other, criterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> DerivedGeographicCRS::_identify(const io::AuthorityFactoryPtr &factory) const { return CRS::_identify(factory); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct DerivedProjectedCRS::Private {}; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress DerivedProjectedCRS::~DerivedProjectedCRS() = default; //! @endcond // --------------------------------------------------------------------------- DerivedProjectedCRS::DerivedProjectedCRS( const ProjectedCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::CoordinateSystemNNPtr &csIn) : SingleCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), DerivedCRS(baseCRSIn, derivingConversionIn, csIn), d(nullptr) {} // --------------------------------------------------------------------------- DerivedProjectedCRS::DerivedProjectedCRS(const DerivedProjectedCRS &other) : SingleCRS(other), DerivedCRS(other), d(nullptr) {} // --------------------------------------------------------------------------- CRSNNPtr DerivedProjectedCRS::_shallowClone() const { auto crs(DerivedProjectedCRS::nn_make_shared(*this)); crs->assignSelf(crs); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- /** \brief Return the base CRS (a ProjectedCRS) of a DerivedProjectedCRS. * * @return the base CRS. */ const ProjectedCRSNNPtr DerivedProjectedCRS::baseCRS() const { return NN_NO_CHECK(util::nn_dynamic_pointer_cast( DerivedCRS::getPrivate()->baseCRS_)); } // --------------------------------------------------------------------------- /** \brief Instantiate a DerivedProjectedCRS from a base CRS, a deriving * conversion and a cs::CS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param baseCRSIn base CRS. * @param derivingConversionIn the deriving conversion from the base CRS to this * CRS. * @param csIn the coordinate system. * @return new DerivedProjectedCRS. */ DerivedProjectedCRSNNPtr DerivedProjectedCRS::create( const util::PropertyMap &properties, const ProjectedCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::CoordinateSystemNNPtr &csIn) { auto crs(DerivedProjectedCRS::nn_make_shared( baseCRSIn, derivingConversionIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void DerivedProjectedCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; if (!isWKT2 || !formatter->use2019Keywords()) { io::FormattingException::Throw( "DerivedProjectedCRS can only be exported to WKT2:2019"); } formatter->startNode(io::WKTConstants::DERIVEDPROJCRS, !identifiers().empty()); formatter->addQuotedString(nameStr()); { auto l_baseProjCRS = baseCRS(); formatter->startNode(io::WKTConstants::BASEPROJCRS, !l_baseProjCRS->identifiers().empty()); formatter->addQuotedString(l_baseProjCRS->nameStr()); auto l_baseGeodCRS = l_baseProjCRS->baseCRS(); auto &geodeticCRSAxisList = l_baseGeodCRS->coordinateSystem()->axisList(); formatter->startNode( dynamic_cast(l_baseGeodCRS.get()) ? io::WKTConstants::BASEGEOGCRS : io::WKTConstants::BASEGEODCRS, !l_baseGeodCRS->identifiers().empty()); formatter->addQuotedString(l_baseGeodCRS->nameStr()); l_baseGeodCRS->exportDatumOrDatumEnsembleToWkt(formatter); // insert ellipsoidal cs unit when the units of the map // projection angular parameters are not explicitly given within those // parameters. See // http://docs.opengeospatial.org/is/12-063r5/12-063r5.html#61 if (formatter->primeMeridianOrParameterUnitOmittedIfSameAsAxis() && !geodeticCRSAxisList.empty()) { geodeticCRSAxisList[0]->unit()._exportToWKT(formatter); } l_baseGeodCRS->primeMeridian()->_exportToWKT(formatter); formatter->endNode(); l_baseProjCRS->derivingConversionRef()->_exportToWKT(formatter); formatter->endNode(); } formatter->setUseDerivingConversion(true); derivingConversionRef()->_exportToWKT(formatter); formatter->setUseDerivingConversion(false); coordinateSystem()->_exportToWKT(formatter); ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } //! @endcond // --------------------------------------------------------------------------- bool DerivedProjectedCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherDerivedCRS = dynamic_cast(other); return otherDerivedCRS != nullptr && DerivedCRS::_isEquivalentTo(other, criterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct TemporalCRS::Private {}; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress TemporalCRS::~TemporalCRS() = default; //! @endcond // --------------------------------------------------------------------------- TemporalCRS::TemporalCRS(const datum::TemporalDatumNNPtr &datumIn, const cs::TemporalCSNNPtr &csIn) : SingleCRS(datumIn.as_nullable(), nullptr, csIn), d(nullptr) {} // --------------------------------------------------------------------------- TemporalCRS::TemporalCRS(const TemporalCRS &other) : SingleCRS(other), d(nullptr) {} // --------------------------------------------------------------------------- CRSNNPtr TemporalCRS::_shallowClone() const { auto crs(TemporalCRS::nn_make_shared(*this)); crs->assignSelf(crs); return crs; } // --------------------------------------------------------------------------- /** \brief Return the datum::TemporalDatum associated with the CRS. * * @return a TemporalDatum */ const datum::TemporalDatumNNPtr TemporalCRS::datum() const { return NN_NO_CHECK(std::static_pointer_cast( SingleCRS::getPrivate()->datum)); } // --------------------------------------------------------------------------- /** \brief Return the cs::TemporalCS associated with the CRS. * * @return a TemporalCS */ const cs::TemporalCSNNPtr TemporalCRS::coordinateSystem() const { return util::nn_static_pointer_cast( SingleCRS::getPrivate()->coordinateSystem); } // --------------------------------------------------------------------------- /** \brief Instantiate a TemporalCRS from a datum and a coordinate system. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datumIn the datum. * @param csIn the coordinate system. * @return new TemporalCRS. */ TemporalCRSNNPtr TemporalCRS::create(const util::PropertyMap &properties, const datum::TemporalDatumNNPtr &datumIn, const cs::TemporalCSNNPtr &csIn) { auto crs(TemporalCRS::nn_make_shared(datumIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void TemporalCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; if (!isWKT2) { io::FormattingException::Throw( "TemporalCRS can only be exported to WKT2"); } formatter->startNode(io::WKTConstants::TIMECRS, !identifiers().empty()); formatter->addQuotedString(nameStr()); datum()->_exportToWKT(formatter); coordinateSystem()->_exportToWKT(formatter); ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void TemporalCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("TemporalCRS", !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } writer->AddObjKey("datum"); formatter->setOmitTypeInImmediateChild(); datum()->_exportToJSON(formatter); writer->AddObjKey("coordinate_system"); formatter->setOmitTypeInImmediateChild(); coordinateSystem()->_exportToJSON(formatter); ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- bool TemporalCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherTemporalCRS = dynamic_cast(other); return otherTemporalCRS != nullptr && SingleCRS::baseIsEquivalentTo(other, criterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct EngineeringCRS::Private {}; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress EngineeringCRS::~EngineeringCRS() = default; //! @endcond // --------------------------------------------------------------------------- EngineeringCRS::EngineeringCRS(const datum::EngineeringDatumNNPtr &datumIn, const cs::CoordinateSystemNNPtr &csIn) : SingleCRS(datumIn.as_nullable(), nullptr, csIn), d(internal::make_unique()) {} // --------------------------------------------------------------------------- EngineeringCRS::EngineeringCRS(const EngineeringCRS &other) : SingleCRS(other), d(internal::make_unique(*(other.d))) {} // --------------------------------------------------------------------------- CRSNNPtr EngineeringCRS::_shallowClone() const { auto crs(EngineeringCRS::nn_make_shared(*this)); crs->assignSelf(crs); return crs; } // --------------------------------------------------------------------------- /** \brief Return the datum::EngineeringDatum associated with the CRS. * * @return a EngineeringDatum */ const datum::EngineeringDatumNNPtr EngineeringCRS::datum() const { return NN_NO_CHECK(std::static_pointer_cast( SingleCRS::getPrivate()->datum)); } // --------------------------------------------------------------------------- /** \brief Instantiate a EngineeringCRS from a datum and a coordinate system. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datumIn the datum. * @param csIn the coordinate system. * @return new EngineeringCRS. */ EngineeringCRSNNPtr EngineeringCRS::create(const util::PropertyMap &properties, const datum::EngineeringDatumNNPtr &datumIn, const cs::CoordinateSystemNNPtr &csIn) { auto crs(EngineeringCRS::nn_make_shared(datumIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void EngineeringCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; formatter->startNode(isWKT2 ? io::WKTConstants::ENGCRS : io::WKTConstants::LOCAL_CS, !identifiers().empty()); formatter->addQuotedString(nameStr()); if (isWKT2 || !datum()->nameStr().empty()) { datum()->_exportToWKT(formatter); } if (!isWKT2) { coordinateSystem()->axisList()[0]->unit()._exportToWKT(formatter); } const auto oldAxisOutputRule = formatter->outputAxis(); formatter->setOutputAxis(io::WKTFormatter::OutputAxisRule::YES); coordinateSystem()->_exportToWKT(formatter); formatter->setOutputAxis(oldAxisOutputRule); ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void EngineeringCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("EngineeringCRS", !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } writer->AddObjKey("datum"); formatter->setOmitTypeInImmediateChild(); datum()->_exportToJSON(formatter); writer->AddObjKey("coordinate_system"); formatter->setOmitTypeInImmediateChild(); coordinateSystem()->_exportToJSON(formatter); ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- bool EngineeringCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherEngineeringCRS = dynamic_cast(other); return otherEngineeringCRS != nullptr && SingleCRS::baseIsEquivalentTo(other, criterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct ParametricCRS::Private {}; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress ParametricCRS::~ParametricCRS() = default; //! @endcond // --------------------------------------------------------------------------- ParametricCRS::ParametricCRS(const datum::ParametricDatumNNPtr &datumIn, const cs::ParametricCSNNPtr &csIn) : SingleCRS(datumIn.as_nullable(), nullptr, csIn), d(nullptr) {} // --------------------------------------------------------------------------- ParametricCRS::ParametricCRS(const ParametricCRS &other) : SingleCRS(other), d(nullptr) {} // --------------------------------------------------------------------------- CRSNNPtr ParametricCRS::_shallowClone() const { auto crs(ParametricCRS::nn_make_shared(*this)); crs->assignSelf(crs); return crs; } // --------------------------------------------------------------------------- /** \brief Return the datum::ParametricDatum associated with the CRS. * * @return a ParametricDatum */ const datum::ParametricDatumNNPtr ParametricCRS::datum() const { return NN_NO_CHECK(std::static_pointer_cast( SingleCRS::getPrivate()->datum)); } // --------------------------------------------------------------------------- /** \brief Return the cs::TemporalCS associated with the CRS. * * @return a TemporalCS */ const cs::ParametricCSNNPtr ParametricCRS::coordinateSystem() const { return util::nn_static_pointer_cast( SingleCRS::getPrivate()->coordinateSystem); } // --------------------------------------------------------------------------- /** \brief Instantiate a ParametricCRS from a datum and a coordinate system. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param datumIn the datum. * @param csIn the coordinate system. * @return new ParametricCRS. */ ParametricCRSNNPtr ParametricCRS::create(const util::PropertyMap &properties, const datum::ParametricDatumNNPtr &datumIn, const cs::ParametricCSNNPtr &csIn) { auto crs(ParametricCRS::nn_make_shared(datumIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void ParametricCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; if (!isWKT2) { io::FormattingException::Throw( "ParametricCRS can only be exported to WKT2"); } formatter->startNode(io::WKTConstants::PARAMETRICCRS, !identifiers().empty()); formatter->addQuotedString(nameStr()); datum()->_exportToWKT(formatter); coordinateSystem()->_exportToWKT(formatter); ObjectUsage::baseExportToWKT(formatter); formatter->endNode(); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void ParametricCRS::_exportToJSON( io::JSONFormatter *formatter) const // throw(io::FormattingException) { auto writer = formatter->writer(); auto objectContext( formatter->MakeObjectContext("ParametricCRS", !identifiers().empty())); writer->AddObjKey("name"); auto l_name = nameStr(); if (l_name.empty()) { writer->Add("unnamed"); } else { writer->Add(l_name); } writer->AddObjKey("datum"); formatter->setOmitTypeInImmediateChild(); datum()->_exportToJSON(formatter); writer->AddObjKey("coordinate_system"); formatter->setOmitTypeInImmediateChild(); coordinateSystem()->_exportToJSON(formatter); ObjectUsage::baseExportToJSON(formatter); } //! @endcond // --------------------------------------------------------------------------- bool ParametricCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherParametricCRS = dynamic_cast(other); return otherParametricCRS != nullptr && SingleCRS::baseIsEquivalentTo(other, criterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress struct DerivedVerticalCRS::Private {}; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress DerivedVerticalCRS::~DerivedVerticalCRS() = default; //! @endcond // --------------------------------------------------------------------------- DerivedVerticalCRS::DerivedVerticalCRS( const VerticalCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::VerticalCSNNPtr &csIn) : SingleCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), VerticalCRS(baseCRSIn->datum(), baseCRSIn->datumEnsemble(), csIn), DerivedCRS(baseCRSIn, derivingConversionIn, csIn), d(nullptr) {} // --------------------------------------------------------------------------- DerivedVerticalCRS::DerivedVerticalCRS(const DerivedVerticalCRS &other) : SingleCRS(other), VerticalCRS(other), DerivedCRS(other), d(nullptr) {} // --------------------------------------------------------------------------- CRSNNPtr DerivedVerticalCRS::_shallowClone() const { auto crs(DerivedVerticalCRS::nn_make_shared(*this)); crs->assignSelf(crs); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- /** \brief Return the base CRS (a VerticalCRS) of a DerivedVerticalCRS. * * @return the base CRS. */ const VerticalCRSNNPtr DerivedVerticalCRS::baseCRS() const { return NN_NO_CHECK(util::nn_dynamic_pointer_cast( DerivedCRS::getPrivate()->baseCRS_)); } // --------------------------------------------------------------------------- /** \brief Instantiate a DerivedVerticalCRS from a base CRS, a deriving * conversion and a cs::VerticalCS. * * @param properties See \ref general_properties. * At minimum the name should be defined. * @param baseCRSIn base CRS. * @param derivingConversionIn the deriving conversion from the base CRS to this * CRS. * @param csIn the coordinate system. * @return new DerivedVerticalCRS. */ DerivedVerticalCRSNNPtr DerivedVerticalCRS::create( const util::PropertyMap &properties, const VerticalCRSNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const cs::VerticalCSNNPtr &csIn) { auto crs(DerivedVerticalCRS::nn_make_shared( baseCRSIn, derivingConversionIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress void DerivedVerticalCRS::_exportToWKT(io::WKTFormatter *formatter) const { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; if (!isWKT2) { bool useBaseMethod = true; const DerivedVerticalCRS *dvcrs = this; while (true) { // If the derived vertical CRS is obtained through simple conversion // methods that just do unit change or height/depth reversal, export // it as a regular VerticalCRS const int methodCode = dvcrs->derivingConversionRef()->method()->getEPSGCode(); if (methodCode == EPSG_CODE_METHOD_CHANGE_VERTICAL_UNIT || methodCode == EPSG_CODE_METHOD_CHANGE_VERTICAL_UNIT_NO_CONV_FACTOR || methodCode == EPSG_CODE_METHOD_HEIGHT_DEPTH_REVERSAL) { dvcrs = dynamic_cast(baseCRS().get()); if (dvcrs == nullptr) { break; } } else { useBaseMethod = false; break; } } if (useBaseMethod) { VerticalCRS::_exportToWKT(formatter); return; } io::FormattingException::Throw( "DerivedVerticalCRS can only be exported to WKT2"); } baseExportToWKT(formatter, io::WKTConstants::VERTCRS, io::WKTConstants::BASEVERTCRS); } //! @endcond // --------------------------------------------------------------------------- void DerivedVerticalCRS::_exportToPROJString( io::PROJStringFormatter *) const // throw(io::FormattingException) { throw io::FormattingException( "DerivedVerticalCRS cannot be exported to PROJ string"); } // --------------------------------------------------------------------------- bool DerivedVerticalCRS::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherDerivedCRS = dynamic_cast(other); return otherDerivedCRS != nullptr && DerivedCRS::_isEquivalentTo(other, criterion, dbContext); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress std::list> DerivedVerticalCRS::_identify(const io::AuthorityFactoryPtr &factory) const { return CRS::_identify(factory); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress template struct DerivedCRSTemplate::Private {}; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress template DerivedCRSTemplate::~DerivedCRSTemplate() = default; //! @endcond // --------------------------------------------------------------------------- template DerivedCRSTemplate::DerivedCRSTemplate( const BaseNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const CSNNPtr &csIn) : SingleCRS(baseCRSIn->datum().as_nullable(), nullptr, csIn), BaseType(baseCRSIn->datum(), csIn), DerivedCRS(baseCRSIn, derivingConversionIn, csIn), d(nullptr) {} // --------------------------------------------------------------------------- template DerivedCRSTemplate::DerivedCRSTemplate( const DerivedCRSTemplate &other) : SingleCRS(other), BaseType(other), DerivedCRS(other), d(nullptr) {} // --------------------------------------------------------------------------- template const typename DerivedCRSTemplate::BaseNNPtr DerivedCRSTemplate::baseCRS() const { auto l_baseCRS = DerivedCRS::getPrivate()->baseCRS_; return NN_NO_CHECK(util::nn_dynamic_pointer_cast(l_baseCRS)); } // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress template CRSNNPtr DerivedCRSTemplate::_shallowClone() const { auto crs(DerivedCRSTemplate::nn_make_shared(*this)); crs->assignSelf(crs); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- template typename DerivedCRSTemplate::NNPtr DerivedCRSTemplate::create( const util::PropertyMap &properties, const BaseNNPtr &baseCRSIn, const operation::ConversionNNPtr &derivingConversionIn, const CSNNPtr &csIn) { auto crs(DerivedCRSTemplate::nn_make_shared( baseCRSIn, derivingConversionIn, csIn)); crs->assignSelf(crs); crs->setProperties(properties); crs->setDerivingConversionCRS(); return crs; } // --------------------------------------------------------------------------- template const char *DerivedCRSTemplate::className() const { return DerivedCRSTraits::CRSName().c_str(); } // --------------------------------------------------------------------------- static void DerivedCRSTemplateCheckExportToWKT(io::WKTFormatter *formatter, const std::string &crsName, bool wkt2_2019_only) { const bool isWKT2 = formatter->version() == io::WKTFormatter::Version::WKT2; if (!isWKT2 || (wkt2_2019_only && !formatter->use2019Keywords())) { io::FormattingException::Throw(crsName + " can only be exported to WKT2" + (wkt2_2019_only ? ":2019" : "")); } } // --------------------------------------------------------------------------- template void DerivedCRSTemplate::_exportToWKT( io::WKTFormatter *formatter) const { DerivedCRSTemplateCheckExportToWKT(formatter, DerivedCRSTraits::CRSName(), DerivedCRSTraits::wkt2_2019_only); baseExportToWKT(formatter, DerivedCRSTraits::WKTKeyword(), DerivedCRSTraits::WKTBaseKeyword()); } // --------------------------------------------------------------------------- template bool DerivedCRSTemplate::_isEquivalentTo( const util::IComparable *other, util::IComparable::Criterion criterion, const io::DatabaseContextPtr &dbContext) const { auto otherDerivedCRS = dynamic_cast(other); return otherDerivedCRS != nullptr && DerivedCRS::_isEquivalentTo(other, criterion, dbContext); } //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static const std::string STRING_DerivedEngineeringCRS("DerivedEngineeringCRS"); const std::string &DerivedEngineeringCRSTraits::CRSName() { return STRING_DerivedEngineeringCRS; } const std::string &DerivedEngineeringCRSTraits::WKTKeyword() { return io::WKTConstants::ENGCRS; } const std::string &DerivedEngineeringCRSTraits::WKTBaseKeyword() { return io::WKTConstants::BASEENGCRS; } template class DerivedCRSTemplate; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static const std::string STRING_DerivedParametricCRS("DerivedParametricCRS"); const std::string &DerivedParametricCRSTraits::CRSName() { return STRING_DerivedParametricCRS; } const std::string &DerivedParametricCRSTraits::WKTKeyword() { return io::WKTConstants::PARAMETRICCRS; } const std::string &DerivedParametricCRSTraits::WKTBaseKeyword() { return io::WKTConstants::BASEPARAMCRS; } template class DerivedCRSTemplate; //! @endcond // --------------------------------------------------------------------------- //! @cond Doxygen_Suppress static const std::string STRING_DerivedTemporalCRS("DerivedTemporalCRS"); const std::string &DerivedTemporalCRSTraits::CRSName() { return STRING_DerivedTemporalCRS; } const std::string &DerivedTemporalCRSTraits::WKTKeyword() { return io::WKTConstants::TIMECRS; } const std::string &DerivedTemporalCRSTraits::WKTBaseKeyword() { return io::WKTConstants::BASETIMECRS; } template class DerivedCRSTemplate; //! @endcond // --------------------------------------------------------------------------- } // namespace crs NS_PROJ_END