aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--src/iso19111/coordinateoperation.cpp87
-rw-r--r--test/unit/test_operation.cpp95
2 files changed, 149 insertions, 33 deletions
diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp
index ef262cc8..6593d3e9 100644
--- a/src/iso19111/coordinateoperation.cpp
+++ b/src/iso19111/coordinateoperation.cpp
@@ -608,6 +608,7 @@ struct CoordinateOperation::Private {
interpolationCRS_(other.interpolationCRS_),
sourceCoordinateEpoch_(other.sourceCoordinateEpoch_),
targetCoordinateEpoch_(other.targetCoordinateEpoch_),
+ hasBallparkTransformation_(other.hasBallparkTransformation_),
strongRef_(other.strongRef_ ? internal::make_unique<CRSStrongRef>(
*(other.strongRef_))
: nullptr) {}
@@ -11038,7 +11039,7 @@ struct CoordinateOperationFactory::Private {
const crs::CRSNNPtr &targetCRS, Context &context);
private:
- static constexpr bool allowEmptyIntersection = true;
+ static constexpr bool disallowEmptyIntersection = true;
static void
buildCRSIds(const crs::CRSNNPtr &crs, Private::Context &context,
@@ -12553,7 +12554,8 @@ createGeodToGeodPROJBased(const crs::CRSNNPtr &geodSrc,
static CoordinateOperationNNPtr createHorizVerticalPROJBased(
const crs::CRSNNPtr &sourceCRS, const crs::CRSNNPtr &targetCRS,
const operation::CoordinateOperationNNPtr &horizTransform,
- const operation::CoordinateOperationNNPtr &verticalTransform) {
+ const operation::CoordinateOperationNNPtr &verticalTransform,
+ bool checkExtent) {
auto geogDst = util::nn_dynamic_pointer_cast<crs::GeographicCRS>(targetCRS);
assert(geogDst);
@@ -12579,10 +12581,16 @@ static CoordinateOperationNNPtr createHorizVerticalPROJBased(
verticalTransform->coordinateOperationAccuracies(),
verticalTransform->hasBallparkTransformation());
} else {
- bool dummy = false;
+ bool emptyIntersection = false;
auto ops = std::vector<CoordinateOperationNNPtr>{horizTransform,
verticalTransform};
- auto extent = getExtent(ops, true, dummy);
+ auto extent = getExtent(ops, true, emptyIntersection);
+ if (checkExtent && emptyIntersection) {
+ std::string msg(
+ "empty intersection of area of validity of concatenated "
+ "operations");
+ throw InvalidOperationEmptyIntersection(msg);
+ }
auto properties = util::PropertyMap();
properties.set(common::IdentifiedObject::NAME_KEY,
computeConcatenatedName(ops));
@@ -12844,7 +12852,7 @@ CoordinateOperationFactory::Private::createOperationsGeogToGeog(
}
auto op = ConcatenatedOperation::createComputeMetadata(
- steps, !allowEmptyIntersection);
+ steps, disallowEmptyIntersection);
op->setHasBallparkTransformation(!sameDatum);
res.emplace_back(op);
return res;
@@ -13113,8 +13121,11 @@ void CoordinateOperationFactory::Private::createOperationsWithDatumPivot(
}
logTrace("transformation " + debugStr);
#endif
- res.emplace_back(ConcatenatedOperation::createComputeMetadata(
- subOps, !allowEmptyIntersection));
+ try {
+ res.emplace_back(ConcatenatedOperation::createComputeMetadata(
+ subOps, disallowEmptyIntersection));
+ } catch (const InvalidOperationEmptyIntersection &) {
+ }
}
};
@@ -13700,7 +13711,7 @@ CoordinateOperationFactory::Private::createOperationsGeogToVertFromGeoid(
context, targetOp, vertDst, tmp);
assert(!tmp.empty());
auto ret = ConcatenatedOperation::createComputeMetadata(
- {op, tmp.front()}, !allowEmptyIntersection);
+ {op, tmp.front()}, disallowEmptyIntersection);
return ret;
};
@@ -13858,7 +13869,7 @@ std::vector<CoordinateOperationNNPtr> CoordinateOperationFactory::Private::
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
{opFirst, opsSecond.front()},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
}
}
}
@@ -14043,7 +14054,7 @@ void CoordinateOperationFactory::Private::createOperationsGeodToGeod(
auto opSecond =
Conversion::createGeographicGeocentric(interm_crs, targetCRS);
res.emplace_back(ConcatenatedOperation::createComputeMetadata(
- {opFirst, opSecond}, !allowEmptyIntersection));
+ {opFirst, opSecond}, disallowEmptyIntersection));
} else {
// Apply previous case in reverse way
std::vector<CoordinateOperationNNPtr> resTmp;
@@ -14109,7 +14120,7 @@ void CoordinateOperationFactory::Private::createOperationsDerivedTo(
for (const auto &opSecond : opsSecond) {
try {
res.emplace_back(ConcatenatedOperation::createComputeMetadata(
- {opFirst, opSecond}, !allowEmptyIntersection));
+ {opFirst, opSecond}, disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
}
@@ -14168,7 +14179,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
ConcatenatedOperation::createComputeMetadata(
{NN_NO_CHECK(opIntermediate),
boundSrc->transformation()},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
} else {
@@ -14191,7 +14202,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
subops.emplace_back(boundSrc->transformation());
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
- subops, !allowEmptyIntersection));
+ subops, disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
}
@@ -14231,7 +14242,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
subops.emplace_back(opLast);
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
- subops, !allowEmptyIntersection));
+ subops, disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
}
@@ -14282,7 +14293,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
try {
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
- {opFirst, transf}, !allowEmptyIntersection));
+ {opFirst, transf}, disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
}
@@ -14308,7 +14319,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
{opFirst, boundSrc->transformation()},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
}
@@ -14329,15 +14340,19 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
for (const auto &opFirst : opsFirst) {
for (const auto &opLast : opsLast) {
// Exclude artificial transformations from the hub
- // to the target CRS
- if (!opLast->hasBallparkTransformation()) {
+ // to the target CRS, if it is the only one.
+ if (opsLast.size() > 1 ||
+ !opLast->hasBallparkTransformation()) {
try {
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
{opFirst, opLast},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
+ } else {
+ // std::cerr << "excluded " << opLast->nameStr() <<
+ // std::endl;
}
}
}
@@ -14395,7 +14410,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
try {
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
- {op, conv}, !allowEmptyIntersection));
+ {op, conv}, disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
}
@@ -14417,10 +14432,13 @@ void CoordinateOperationFactory::Private::createOperationsBoundToGeog(
ConcatenatedOperation::
createComputeMetadata(
{opFirst, opLast},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
} catch (
const InvalidOperationEmptyIntersection &) {
}
+ } else {
+ // std::cerr << "excluded " << opLast->nameStr() <<
+ // std::endl;
}
}
}
@@ -14625,7 +14643,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToBound(
}
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
- ops, !allowEmptyIntersection));
+ ops, disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
}
@@ -14672,7 +14690,7 @@ void CoordinateOperationFactory::Private::createOperationsBoundToBound(
try {
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
- {opFirst, opLast}, !allowEmptyIntersection));
+ {opFirst, opLast}, disallowEmptyIntersection));
} catch (const InvalidOperationEmptyIntersection &) {
}
}
@@ -14844,7 +14862,7 @@ void CoordinateOperationFactory::Private::createOperationsCompoundToGeog(
res.emplace_back(
ConcatenatedOperation::createComputeMetadata(
{opsFirst.front(), opLast},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
} catch (const std::exception &) {
}
}
@@ -15082,7 +15100,7 @@ void CoordinateOperationFactory::Private::createOperationsCompoundToGeog(
try {
auto op = createHorizVerticalPROJBased(
sourceCRS, targetCRS, horizTransform,
- verticalTransform);
+ verticalTransform, disallowEmptyIntersection);
res.emplace_back(op);
} catch (const std::exception &) {
}
@@ -15117,9 +15135,13 @@ void CoordinateOperationFactory::Private::createOperationsToGeod(
for (const auto &op : sourceToGeog3DOps) {
auto newOp = op->shallowClone();
setCRSs(newOp.get(), sourceCRS, intermGeog3DCRS);
- res.emplace_back(ConcatenatedOperation::createComputeMetadata(
- {newOp, geog3DToTargetOps.front()},
- !allowEmptyIntersection));
+ try {
+ res.emplace_back(
+ ConcatenatedOperation::createComputeMetadata(
+ {newOp, geog3DToTargetOps.front()},
+ disallowEmptyIntersection));
+ } catch (const InvalidOperationEmptyIntersection &) {
+ }
}
}
}
@@ -15186,11 +15208,14 @@ void CoordinateOperationFactory::Private::createOperationsCompoundToCompound(
createOperations(intermGeogSrc, intermGeogDst, context);
for (const auto &op1 : opsSrcToGeog) {
if (op1->hasBallparkTransformation()) {
+ // std::cerr << "excluded " << op1->nameStr() << std::endl;
continue;
}
for (const auto &op2 : opsGeogSrcToGeogDst) {
for (const auto &op3 : opsGeogToTarget) {
if (op3->hasBallparkTransformation()) {
+ // std::cerr << "excluded " << op3->nameStr() <<
+ // std::endl;
continue;
}
try {
@@ -15204,7 +15229,7 @@ void CoordinateOperationFactory::Private::createOperationsCompoundToCompound(
CoordinateOperationNNPtr>{op1,
op2,
op3},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
} catch (const std::exception &) {
}
}
@@ -15340,14 +15365,14 @@ void CoordinateOperationFactory::Private::createOperationsBoundToCompound(
createComputeMetadata(
{opSrc, intermOps.front(),
opDst},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
}
} else {
res.emplace_back(
ConcatenatedOperation::
createComputeMetadata(
{opSrc, opDst},
- !allowEmptyIntersection));
+ disallowEmptyIntersection));
}
}
}
diff --git a/test/unit/test_operation.cpp b/test/unit/test_operation.cpp
index 23f3489b..c7c0a262 100644
--- a/test/unit/test_operation.cpp
+++ b/test/unit/test_operation.cpp
@@ -6801,7 +6801,7 @@ TEST(operation, nadgrids_with_pm) {
dst = authFactory->createCoordinateReferenceSystem("4258");
list = CoordinateOperationFactory::create()->createOperations(
NN_NO_CHECK(src), dst, ctxt);
- ASSERT_EQ(list.size(), 1U);
+ ASSERT_GE(list.size(), 1U);
EXPECT_EQ(list[0]->exportToPROJString(PROJStringFormatter::create().get()),
"+proj=pipeline "
"+step +inv +proj=tmerc +lat_0=39.6666666666667 +lon_0=1 "
@@ -6819,7 +6819,7 @@ TEST(operation, nadgrids_with_pm) {
ASSERT_TRUE(crsFromWkt);
list = CoordinateOperationFactory::create()->createOperations(
NN_NO_CHECK(crsFromWkt), dst, ctxt);
- ASSERT_EQ(list.size(), 1U);
+ ASSERT_GE(list.size(), 1U);
EXPECT_EQ(list[0]->exportToPROJString(PROJStringFormatter::create().get()),
"+proj=pipeline "
"+step +inv +proj=tmerc +lat_0=39.6666666666667 +lon_0=1 "
@@ -7621,6 +7621,97 @@ TEST(
// ---------------------------------------------------------------------------
+TEST(operation,
+ compoundCRS_with_boundGeogCRS_and_geoid_to_geodCRS_NAD2011_ctxt) {
+
+ auto dbContext = DatabaseContext::create();
+
+ const char *wktSrc =
+ "COMPD_CS[\"NAD83 / California zone 5 (ftUS) + "
+ "NAVD88 height - Geoid12B (ftUS)\","
+ " PROJCS[\"NAD83 / California zone 5 (ftUS)\","
+ " GEOGCS[\"NAD83\","
+ " DATUM[\"North_American_Datum_1983\","
+ " SPHEROID[\"GRS 1980\",6378137,298.257222101,"
+ " AUTHORITY[\"EPSG\",\"7019\"]],"
+ " TOWGS84[0,0,0,0,0,0,0],"
+ " AUTHORITY[\"EPSG\",\"6269\"]],"
+ " PRIMEM[\"Greenwich\",0,"
+ " AUTHORITY[\"EPSG\",\"8901\"]],"
+ " UNIT[\"degree\",0.0174532925199433,"
+ " AUTHORITY[\"EPSG\",\"9122\"]],"
+ " AUTHORITY[\"EPSG\",\"4269\"]],"
+ " PROJECTION[\"Lambert_Conformal_Conic_2SP\"],"
+ " PARAMETER[\"standard_parallel_1\",35.46666666666667],"
+ " PARAMETER[\"standard_parallel_2\",34.03333333333333],"
+ " PARAMETER[\"latitude_of_origin\",33.5],"
+ " PARAMETER[\"central_meridian\",-118],"
+ " PARAMETER[\"false_easting\",6561666.667],"
+ " PARAMETER[\"false_northing\",1640416.667],"
+ " UNIT[\"US survey foot\",0.3048006096012192,"
+ " AUTHORITY[\"EPSG\",\"9003\"]],"
+ " AXIS[\"X\",EAST],"
+ " AXIS[\"Y\",NORTH],"
+ " AUTHORITY[\"EPSG\",\"2229\"]],"
+ "VERT_CS[\"NAVD88 height - Geoid12B (ftUS)\","
+ " VERT_DATUM[\"North American Vertical Datum 1988\",2005,"
+ " AUTHORITY[\"EPSG\",\"5103\"]],"
+ " UNIT[\"US survey foot\",0.3048006096012192,"
+ " AUTHORITY[\"EPSG\",\"9003\"]],"
+ " AXIS[\"Gravity-related height\",UP],"
+ " AUTHORITY[\"EPSG\",\"6360\"]]]";
+ auto objSrc =
+ WKTParser().attachDatabaseContext(dbContext).createFromWKT(wktSrc);
+ auto srcCRS = nn_dynamic_pointer_cast<CompoundCRS>(objSrc);
+ ASSERT_TRUE(srcCRS != nullptr);
+
+ auto authFactoryEPSG = AuthorityFactory::create(dbContext, "EPSG");
+ // NAD83(2011) geocentric
+ auto dstCRS = authFactoryEPSG->createCoordinateReferenceSystem("6317");
+
+ auto authFactory = AuthorityFactory::create(dbContext, std::string());
+ auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0);
+ ctxt->setGridAvailabilityUse(
+ CoordinateOperationContext::GridAvailabilityUse::
+ IGNORE_GRID_AVAILABILITY);
+ ctxt->setSpatialCriterion(
+ CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION);
+ auto list = CoordinateOperationFactory::create()->createOperations(
+ NN_NO_CHECK(srcCRS), dstCRS, ctxt);
+ bool found = false;
+ for (const auto &op : list) {
+ if (op->nameStr() ==
+ "Inverse of unnamed + "
+ "Transformation from NAD83 to WGS84 + "
+ "Ballpark geographic offset from WGS 84 to NAD83(2011) + "
+ "Transformation from NAVD88 height (ftUS) to NAVD88 height + "
+ "Inverse of NAD83(2011) to NAVD88 height (1) + "
+ "Conversion from NAD83(2011) (geog3D) to NAD83(2011) "
+ "(geocentric)") {
+ found = true;
+ EXPECT_EQ(
+ op->exportToPROJString(PROJStringFormatter::create().get()),
+ "+proj=pipeline "
+ "+step +proj=unitconvert +xy_in=us-ft +xy_out=m "
+ "+step +inv +proj=lcc +lat_0=33.5 +lon_0=-118 "
+ "+lat_1=35.4666666666667 +lat_2=34.0333333333333 "
+ "+x_0=2000000.0001016 +y_0=500000.0001016 +ellps=GRS80 "
+ "+step +proj=unitconvert +z_in=us-ft +z_out=m "
+ "+step +proj=vgridshift +grids=us_noaa_g2012bu0.tif "
+ "+multiplier=1 "
+ "+step +proj=cart +ellps=GRS80");
+ }
+ }
+ EXPECT_TRUE(found);
+ if (!found) {
+ for (const auto &op : list) {
+ std::cerr << op->nameStr() << std::endl;
+ }
+ }
+}
+
+// ---------------------------------------------------------------------------
+
TEST(operation, geocent_to_compoundCRS) {
auto objSrc = PROJStringParser().createFromPROJString(
"+proj=geocent +datum=WGS84 +units=m +type=crs");