Merge pull request #2364 from rouault/fix_issue_NAD83_NAD83_2011

createOperations(): several fixes affecting NAD83 -> NAD83(2011)

1 files changed, 56 insertions, 31 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));
                             }
                         }
                     }