Skip to content

Commit

Permalink
createOperations(): support one of source/target to be a CoordinateMe…
Browse files Browse the repository at this point in the history 
…tadata

but not both.

Note also that the core of createOperations() does not take into account
the source/target coordinate epoch in its logic to infer pipelines. It
just keep those values in its metadata, so that PJ* objects constructed
from them can later override the PJ_COORD.xyzt.t component.

So basically this is currently useful only for transformations between a
dynamic CRS and a static CRS (or the reverse) with time dependent Helmert
transformations.
  • Loading branch information
@rouault
rouault committed Jan 14, 2023
1 parent f6ec2b2 commit e70473e
Show file tree
Hide file tree
Showing 11 changed files with 385 additions and 27 deletions.
10 changes: 8 additions & 2 deletions docs/source/development/reference/functions.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -139,6 +139,9 @@ paragraph for more details.
- more generally any string accepted by :c:func:`proj_create` representing
a CRS
Starting with PROJ 9.2, source_crs or target_crs can be a CoordinateMetadata
with an associated coordinate epoch (but only one of them, not both).
An "area of use" can be specified in area. When it is supplied, the more
accurate transformation between two given systems can be chosen.
Expand All @@ -164,9 +167,9 @@ paragraph for more details.
:param ctx: Threading context.
:type ctx: :c:type:`PJ_CONTEXT` *
:param `source_crs`: Source CRS.
:param `source_crs`: Source CRS or CoordinateMetadata.
:type `source_crs`: `const char*`
:param `target_crs`: Destination SRS.
:param `target_crs`: Destination SRS or CoordinateMetadata
:type `target_crs`: `const char*`
:param `area`: Descriptor of the desired area for the transformation.
:type `area`: :c:type:`PJ_AREA` *
Expand All @@ -182,6 +185,9 @@ paragraph for more details.
This is the same as :c:func:`proj_create_crs_to_crs` except that the source and
target CRS are passed as PJ* objects which must be of the CRS variety.
Starting with PROJ 9.2, source_crs or target_crs can be a CoordinateMetadata
with an associated coordinate epoch (but only one of them, not both).
:param `options`: a list of NUL terminated options, or NULL.
The list of supported options is:
Expand Down
37 changes: 37 additions & 0 deletions include/proj/coordinateoperation.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,6 +49,12 @@ class DerivedCRS;
class ProjectedCRS;
} // namespace crs

namespace coordinates {
class CoordinateMetadata;
using CoordinateMetadataPtr = std::shared_ptr<CoordinateMetadata>;
using CoordinateMetadataNNPtr = util::nn<CoordinateMetadataPtr>;
} // namespace coordinates

/** osgeo.proj.operation namespace
\brief Coordinate operations (relationship between any two coordinate
Expand Down Expand Up @@ -193,6 +199,11 @@ class PROJ_GCC_DLL CoordinateOperation : public common::ObjectUsage,
const std::vector<metadata::PositionalAccuracyNNPtr> &accuracies);
PROJ_INTERNAL void setHasBallparkTransformation(bool b);

PROJ_INTERNAL void
setSourceCoordinateEpoch(const util::optional<common::DataEpoch> &epoch);
PROJ_INTERNAL void
setTargetCoordinateEpoch(const util::optional<common::DataEpoch> &epoch);

PROJ_INTERNAL void
setProperties(const util::PropertyMap
&properties); // throw(InvalidValueTypeException)
Expand Down Expand Up @@ -1899,12 +1910,28 @@ class PROJ_GCC_DLL CoordinateOperationContext {
PROJ_DLL const std::vector<std::pair<std::string, std::string>> &
getIntermediateCRS() const;

PROJ_DLL void
setSourceCoordinateEpoch(const util::optional<common::DataEpoch> &epoch);

PROJ_DLL const util::optional<common::DataEpoch> &
getSourceCoordinateEpoch() const;

PROJ_DLL void
setTargetCoordinateEpoch(const util::optional<common::DataEpoch> &epoch);

PROJ_DLL const util::optional<common::DataEpoch> &
getTargetCoordinateEpoch() const;

PROJ_DLL static CoordinateOperationContextNNPtr
create(const io::AuthorityFactoryPtr &authorityFactory,
const metadata::ExtentPtr &extent, double accuracy);

PROJ_DLL CoordinateOperationContextNNPtr clone() const;

protected:
PROJ_INTERNAL CoordinateOperationContext();
PROJ_INTERNAL
CoordinateOperationContext(const CoordinateOperationContext &);
INLINED_MAKE_UNIQUE

private:
Expand Down Expand Up @@ -1942,6 +1969,16 @@ class PROJ_GCC_DLL CoordinateOperationFactory {
const crs::CRSNNPtr &targetCRS,
const CoordinateOperationContextNNPtr &context) const;

PROJ_DLL std::vector<CoordinateOperationNNPtr> createOperations(
const coordinates::CoordinateMetadataNNPtr &sourceCoordinateMetadata,
const crs::CRSNNPtr &targetCRS,
const CoordinateOperationContextNNPtr &context) const;

PROJ_DLL std::vector<CoordinateOperationNNPtr> createOperations(
const crs::CRSNNPtr &sourceCRS,
const coordinates::CoordinateMetadataNNPtr &targetCoordinateMetadata,
const CoordinateOperationContextNNPtr &context) const;

PROJ_DLL static CoordinateOperationFactoryNNPtr create();

protected:
Expand Down
7 changes: 7 additions & 0 deletions scripts/reference_exported_symbols.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -614,6 +614,7 @@ osgeo::proj::operation::Conversion::createWagnerV(osgeo::proj::util::PropertyMap
osgeo::proj::operation::Conversion::identify() const
osgeo::proj::operation::Conversion::inverse() const
osgeo::proj::operation::Conversion::isUTM(int&, bool&) const
osgeo::proj::operation::CoordinateOperationContext::clone() const
osgeo::proj::operation::CoordinateOperationContext::~CoordinateOperationContext()
osgeo::proj::operation::CoordinateOperationContext::create(std::shared_ptr<osgeo::proj::io::AuthorityFactory> const&, std::shared_ptr<osgeo::proj::metadata::Extent> const&, double)
osgeo::proj::operation::CoordinateOperationContext::getAllowBallparkTransformations() const
Expand All @@ -625,7 +626,9 @@ osgeo::proj::operation::CoordinateOperationContext::getDiscardSuperseded() const
osgeo::proj::operation::CoordinateOperationContext::getGridAvailabilityUse() const
osgeo::proj::operation::CoordinateOperationContext::getIntermediateCRS() const
osgeo::proj::operation::CoordinateOperationContext::getSourceAndTargetCRSExtentUse() const
osgeo::proj::operation::CoordinateOperationContext::getSourceCoordinateEpoch() const
osgeo::proj::operation::CoordinateOperationContext::getSpatialCriterion() const
osgeo::proj::operation::CoordinateOperationContext::getTargetCoordinateEpoch() const
osgeo::proj::operation::CoordinateOperationContext::getUsePROJAlternativeGridNames() const
osgeo::proj::operation::CoordinateOperationContext::setAllowBallparkTransformations(bool)
osgeo::proj::operation::CoordinateOperationContext::setAllowUseIntermediateCRS(osgeo::proj::operation::CoordinateOperationContext::IntermediateCRSUse)
Expand All @@ -635,13 +638,17 @@ osgeo::proj::operation::CoordinateOperationContext::setDiscardSuperseded(bool)
osgeo::proj::operation::CoordinateOperationContext::setGridAvailabilityUse(osgeo::proj::operation::CoordinateOperationContext::GridAvailabilityUse)
osgeo::proj::operation::CoordinateOperationContext::setIntermediateCRS(std::vector<std::pair<std::string, std::string>, std::allocator<std::pair<std::string, std::string> > > const&)
osgeo::proj::operation::CoordinateOperationContext::setSourceAndTargetCRSExtentUse(osgeo::proj::operation::CoordinateOperationContext::SourceTargetCRSExtentUse)
osgeo::proj::operation::CoordinateOperationContext::setSourceCoordinateEpoch(osgeo::proj::util::optional<osgeo::proj::common::DataEpoch> const&)
osgeo::proj::operation::CoordinateOperationContext::setSpatialCriterion(osgeo::proj::operation::CoordinateOperationContext::SpatialCriterion)
osgeo::proj::operation::CoordinateOperationContext::setTargetCoordinateEpoch(osgeo::proj::util::optional<osgeo::proj::common::DataEpoch> const&)
osgeo::proj::operation::CoordinateOperationContext::setUsePROJAlternativeGridNames(bool)
osgeo::proj::operation::CoordinateOperation::~CoordinateOperation()
osgeo::proj::operation::CoordinateOperation::coordinateOperationAccuracies() const
osgeo::proj::operation::CoordinateOperationFactory::~CoordinateOperationFactory()
osgeo::proj::operation::CoordinateOperationFactory::create()
osgeo::proj::operation::CoordinateOperationFactory::createOperation(dropbox::oxygen::nn<std::shared_ptr<osgeo::proj::crs::CRS> > const&, dropbox::oxygen::nn<std::shared_ptr<osgeo::proj::crs::CRS> > const&) const
osgeo::proj::operation::CoordinateOperationFactory::createOperations(dropbox::oxygen::nn<std::shared_ptr<osgeo::proj::coordinates::CoordinateMetadata> > const&, dropbox::oxygen::nn<std::shared_ptr<osgeo::proj::crs::CRS> > const&, dropbox::oxygen::nn<std::unique_ptr<osgeo::proj::operation::CoordinateOperationContext, std::default_delete<osgeo::proj::operation::CoordinateOperationContext> > > const&) const
osgeo::proj::operation::CoordinateOperationFactory::createOperations(dropbox::oxygen::nn<std::shared_ptr<osgeo::proj::crs::CRS> > const&, dropbox::oxygen::nn<std::shared_ptr<osgeo::proj::coordinates::CoordinateMetadata> > const&, dropbox::oxygen::nn<std::unique_ptr<osgeo::proj::operation::CoordinateOperationContext, std::default_delete<osgeo::proj::operation::CoordinateOperationContext> > > const&) const
osgeo::proj::operation::CoordinateOperationFactory::createOperations(dropbox::oxygen::nn<std::shared_ptr<osgeo::proj::crs::CRS> > const&, dropbox::oxygen::nn<std::shared_ptr<osgeo::proj::crs::CRS> > const&, dropbox::oxygen::nn<std::unique_ptr<osgeo::proj::operation::CoordinateOperationContext, std::default_delete<osgeo::proj::operation::CoordinateOperationContext> > > const&) const
osgeo::proj::operation::CoordinateOperation::hasBallparkTransformation() const
osgeo::proj::operation::CoordinateOperation::interpolationCRS() const
Expand Down
24 changes: 24 additions & 0 deletions src/4D_api.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -374,6 +374,8 @@ similarly, but prefers the 2D resp. 3D interfaces if available.
P->iCurCoordOp = iBest;
}
PJ_COORD res = coord;
if( alt.pj->hasCoordinateEpoch )
coord.xyzt.t = alt.pj->coordinateEpoch;
if( direction == PJ_FWD )
pj_fwd4d( res, alt.pj );
else
Expand Down Expand Up @@ -439,6 +441,8 @@ similarly, but prefers the 2D resp. 3D interfaces if available.
}

P->iCurCoordOp = 0; // dummy value, to be used by proj_trans_get_last_used_operation()
if( P->hasCoordinateEpoch )
coord.xyzt.t = P->coordinateEpoch;
if (direction == PJ_FWD)
pj_fwd4d (coord, P);
else
Expand Down Expand Up @@ -1691,9 +1695,29 @@ static PJ* add_coord_op_to_list(
return op;
}

namespace {
struct ObjectKeeper {
PJ *m_obj = nullptr;
explicit ObjectKeeper(PJ *obj) : m_obj(obj) {}
~ObjectKeeper() { proj_destroy(m_obj); }
ObjectKeeper(const ObjectKeeper &) = delete;
ObjectKeeper& operator=(const ObjectKeeper &) = delete;
};
} // namespace

/*****************************************************************************/
static PJ* create_operation_to_geog_crs(PJ_CONTEXT* ctx, const PJ* crs) {
/*****************************************************************************/

std::unique_ptr<ObjectKeeper> keeper;
if( proj_get_type(crs) == PJ_TYPE_COORDINATE_METADATA ) {
auto tmp = proj_get_source_crs(ctx, crs);
assert(tmp);
keeper.reset(new ObjectKeeper(tmp));
crs = tmp;
}
(void)keeper;

// Create a geographic 2D long-lat degrees CRS that is related to the
// CRS
auto geodetic_crs = proj_crs_get_geodetic_crs(ctx, crs);
Expand Down
66 changes: 55 additions & 11 deletions src/apps/projinfo.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,7 @@

#include <proj/common.hpp>
#include <proj/coordinateoperation.hpp>
#include <proj/coordinates.hpp>
#include <proj/crs.hpp>
#include <proj/io.hpp>
#include <proj/metadata.hpp>
Expand All @@ -49,6 +50,7 @@
#include "proj/internal/internal.hpp" // for split

using namespace NS_PROJ::common;
using namespace NS_PROJ::coordinates;
using namespace NS_PROJ::crs;
using namespace NS_PROJ::io;
using namespace NS_PROJ::metadata;
Expand Down Expand Up @@ -835,22 +837,53 @@ static void outputOperations(
CoordinateOperationContext::IntermediateCRSUse::NEVER,
promoteTo3D, normalizeAxisOrder, outputOpt.quiet);
auto sourceCRS = nn_dynamic_pointer_cast<CRS>(sourceObj);
CoordinateMetadataPtr sourceCoordinateMetadata;
if (!sourceCRS) {
std::cerr << "source CRS string is not a CRS" << std::endl;
std::exit(1);
sourceCoordinateMetadata =
nn_dynamic_pointer_cast<CoordinateMetadata>(sourceObj);
if (!sourceCoordinateMetadata) {
std::cerr
<< "source CRS string is not a CRS or a CoordinateMetadata"
<< std::endl;
std::exit(1);
}
if (!sourceCoordinateMetadata->coordinateEpoch().has_value()) {
sourceCRS = sourceCoordinateMetadata->crs().as_nullable();
sourceCoordinateMetadata.reset();
}
}

auto targetObj =
buildObject(dbContext, targetCRSStr, "crs", "target CRS", false,
CoordinateOperationContext::IntermediateCRSUse::NEVER,
promoteTo3D, normalizeAxisOrder, outputOpt.quiet);
auto targetCRS = nn_dynamic_pointer_cast<CRS>(targetObj);
CoordinateMetadataPtr targetCoordinateMetadata;
if (!targetCRS) {
std::cerr << "target CRS string is not a CRS" << std::endl;
targetCoordinateMetadata =
nn_dynamic_pointer_cast<CoordinateMetadata>(targetObj);
if (!targetCoordinateMetadata) {
std::cerr
<< "target CRS string is not a CRS or a CoordinateMetadata"
<< std::endl;
std::exit(1);
}
if (!targetCoordinateMetadata->coordinateEpoch().has_value()) {
targetCRS = targetCoordinateMetadata->crs().as_nullable();
targetCoordinateMetadata.reset();
}
}

if (sourceCoordinateMetadata != nullptr &&
targetCoordinateMetadata != nullptr) {
std::cerr << "CoordinateMetadata with epoch to CoordinateMetadata "
"with epoch not supported currently."
<< std::endl;
std::exit(1);
}

if (dbContext && !promoteTo3D) {
// TODO: handle promotion of CoordinateMetadata
if (sourceCRS && targetCRS && dbContext && !promoteTo3D) {
// Auto-promote source/target CRS if it is specified by its name,
// if it has a known 3D version of it and that the other CRS is 3D.
// e.g projinfo -s "WGS 84 + EGM96 height" -t "WGS 84"
Expand All @@ -875,8 +908,6 @@ static void outputOperations(
}
}

auto nnSourceCRS = NN_NO_CHECK(sourceCRS);
auto nnTargetCRS = NN_NO_CHECK(targetCRS);
std::vector<CoordinateOperationNNPtr> list;
size_t spatialCriterionPartialIntersectionResultCount = 0;
bool spatialCriterionPartialIntersectionMoreRelevant = false;
Expand All @@ -888,6 +919,22 @@ static void outputOperations(
: nullptr;
auto ctxt =
CoordinateOperationContext::create(authFactory, bboxFilter, 0);

const auto createOperations = [&]() {
if (sourceCoordinateMetadata) {
return CoordinateOperationFactory::create()->createOperations(
NN_NO_CHECK(sourceCoordinateMetadata),
NN_NO_CHECK(targetCRS), ctxt);
} else if (targetCoordinateMetadata) {
return CoordinateOperationFactory::create()->createOperations(
NN_NO_CHECK(sourceCRS),
NN_NO_CHECK(targetCoordinateMetadata), ctxt);
} else {
return CoordinateOperationFactory::create()->createOperations(
NN_NO_CHECK(sourceCRS), NN_NO_CHECK(targetCRS), ctxt);
}
};

ctxt->setSpatialCriterion(spatialCriterion);
ctxt->setSourceAndTargetCRSExtentUse(crsExtentUse);
ctxt->setGridAvailabilityUse(gridAvailabilityUse);
Expand All @@ -899,18 +946,15 @@ static void outputOperations(
if (minimumAccuracy >= 0) {
ctxt->setDesiredAccuracy(minimumAccuracy);
}
list = CoordinateOperationFactory::create()->createOperations(
nnSourceCRS, nnTargetCRS, ctxt);
list = createOperations();
if (!spatialCriterionExplicitlySpecified &&
spatialCriterion == CoordinateOperationContext::SpatialCriterion::
STRICT_CONTAINMENT) {
try {
ctxt->setSpatialCriterion(
CoordinateOperationContext::SpatialCriterion::
PARTIAL_INTERSECTION);
auto list2 =
CoordinateOperationFactory::create()->createOperations(
nnSourceCRS, nnTargetCRS, ctxt);
auto list2 = createOperations();
spatialCriterionPartialIntersectionResultCount = list2.size();
if (spatialCriterionPartialIntersectionResultCount == 1 &&
list.size() == 1 &&
Expand Down
Loading

0 comments on commit e70473e

Please sign in to comment.