Merge pull request #2068 from rouault/add_proj_get_suggested_operation

Add proj_get_suggested_operation()

4 files changed, 341 insertions, 188 deletions

diff --git a/src/4D_api.cpp b/src/4D_api.cpp
index 52bd7f4f..6e494793 100644
--- a/src/4D_api.cpp
+++ b/src/4D_api.cpp
@@ -179,7 +179,58 @@ double proj_roundtrip (PJ *P, PJ_DIRECTION direction, int n, PJ_COORD *coord) {
     return proj_xyz_dist (org, t);
 }
 
+/**************************************************************************************/
+int pj_get_suggested_operation(PJ_CONTEXT*,
+                               const std::vector<CoordOperation>& opList,
+                               const int iExcluded[2],
+                               PJ_DIRECTION direction,
+                               PJ_COORD coord)
+/**************************************************************************************/
+{
+    // Select the operations that match the area of use
+    // and has the best accuracy.
+    int iBest = -1;
+    double bestAccuracy = std::numeric_limits<double>::max();
+    const int nOperations = static_cast<int>(opList.size());
+    for( int i = 0; i < nOperations; i++ ) {
+        if( i == iExcluded[0] || i == iExcluded[1] ) {
+            continue;
+        }
+        const auto &alt = opList[i];
+        bool spatialCriterionOK = false;
+        if( direction == PJ_FWD ) {
+            if( coord.xyzt.x >= alt.minxSrc &&
+                coord.xyzt.y >= alt.minySrc &&
+                coord.xyzt.x <= alt.maxxSrc &&
+                coord.xyzt.y <= alt.maxySrc) {
+                spatialCriterionOK = true;
+            }
+        } else {
+            if( coord.xyzt.x >= alt.minxDst &&
+                coord.xyzt.y >= alt.minyDst &&
+                coord.xyzt.x <= alt.maxxDst &&
+                coord.xyzt.y <= alt.maxyDst ) {
+                spatialCriterionOK = true;
+            }
+        }
 
+        if( spatialCriterionOK ) {
+            // The offshore test is for the "Test bug 245 (use +datum=carthage)"
+            // of testvarious. The long=10 lat=34 point belongs both to the
+            // onshore and offshore Tunisia area of uses, but is slightly
+            // onshore. So in a general way, prefer a onshore area to a
+            // offshore one.
+            if( iBest < 0 ||
+                (alt.accuracy >= 0 && alt.accuracy < bestAccuracy &&
+                !alt.isOffshore) ) {
+                iBest = i;
+                bestAccuracy = alt.accuracy;
+            }
+        }
+    }
+
+    return iBest;
+}
 
 /**************************************************************************************/
 PJ_COORD proj_trans (PJ *P, PJ_DIRECTION direction, PJ_COORD coord) {
@@ -211,45 +262,11 @@ similarly, but prefers the 2D resp. 3D interfaces if available.
         {
             // Do a first pass and select the operations that match the area of use
             // and has the best accuracy.
-            int iBest = -1;
-            double bestAccuracy = std::numeric_limits<double>::max();
-            for( int i = 0; i < nOperations; i++ ) {
-                if( i == iExcluded[0] || i == iExcluded[1] ) {
-                    continue;
-                }
-                const auto &alt = P->alternativeCoordinateOperations[i];
-                bool spatialCriterionOK = false;
-                if( direction == PJ_FWD ) {
-                    if( coord.xyzt.x >= alt.minxSrc &&
-                        coord.xyzt.y >= alt.minySrc &&
-                        coord.xyzt.x <= alt.maxxSrc &&
-                        coord.xyzt.y <= alt.maxySrc) {
-                        spatialCriterionOK = true;
-                    }
-                } else {
-                    if( coord.xyzt.x >= alt.minxDst &&
-                        coord.xyzt.y >= alt.minyDst &&
-                        coord.xyzt.x <= alt.maxxDst &&
-                        coord.xyzt.y <= alt.maxyDst ) {
-                        spatialCriterionOK = true;
-                    }
-                }
-
-                if( spatialCriterionOK ) {
-                    // The offshore test is for the "Test bug 245 (use +datum=carthage)"
-                    // of testvarious. The long=10 lat=34 point belongs both to the
-                    // onshore and offshore Tunisia area of uses, but is slightly
-                    // onshore. So in a general way, prefer a onshore area to a
-                    // offshore one.
-                    if( iBest < 0 ||
-                        (alt.accuracy >= 0 && alt.accuracy < bestAccuracy &&
-                        !alt.isOffshore) ) {
-                        iBest = i;
-                        bestAccuracy = alt.accuracy;
-                    }
-                }
-            }
-
+            int iBest = pj_get_suggested_operation(P->ctx,
+                                                   P->alternativeCoordinateOperations,
+                                                   iExcluded,
+                                                   direction,
+                                                   coord);
             if( iBest < 0 ) {
                 break;
             }
@@ -968,13 +985,15 @@ static void reproject_bbox(PJ* pjGeogToCrs,
 
 
 /*****************************************************************************/
-static PJ* add_coord_op_to_list(PJ* op,
+static PJ* add_coord_op_to_list(
+                            int idxInOriginalList,
+                            PJ* op,
                             double west_lon, double south_lat,
                             double east_lon, double north_lat,
                             PJ* pjGeogToSrc,
                             PJ* pjGeogToDst,
                             bool isOffshore,
-                            std::vector<PJconsts::CoordOperation>& altCoordOps) {
+                            std::vector<CoordOperation>& altCoordOps) {
 /*****************************************************************************/
 
     double minxSrc;
@@ -997,9 +1016,10 @@ static PJ* add_coord_op_to_list(PJ* op,
         std::string name(c_name ? c_name : "");
 
         const double accuracy = proj_coordoperation_get_accuracy(op->ctx, op);
-        altCoordOps.emplace_back(minxSrc, minySrc, maxxSrc, maxySrc,
-                                    minxDst, minyDst, maxxDst, maxyDst,
-                                    op, name, accuracy, isOffshore);
+        altCoordOps.emplace_back(idxInOriginalList,
+                                 minxSrc, minySrc, maxxSrc, maxySrc,
+                                 minxDst, minyDst, maxxDst, maxyDst,
+                                 op, name, accuracy, isOffshore);
         op = nullptr;
     }
     return op;
@@ -1140,6 +1160,91 @@ PJ  *proj_create_crs_to_crs (PJ_CONTEXT *ctx, const char *source_crs, const char
     return ret;
 }
 
+
+/*****************************************************************************/
+std::vector<CoordOperation> pj_create_prepared_operations(PJ_CONTEXT *ctx,
+                                                     const PJ *source_crs,
+                                                     const PJ *target_crs,
+                                                     PJ_OBJ_LIST* op_list)
+/*****************************************************************************/
+{
+    auto pjGeogToSrc = create_operation_to_geog_crs(ctx, source_crs);
+    if( !pjGeogToSrc )
+    {
+        proj_context_log_debug(ctx,
+            "Cannot create transformation from geographic CRS of source CRS to source CRS");
+        return {};
+    }
+
+    auto pjGeogToDst = create_operation_to_geog_crs(ctx, target_crs);
+    if( !pjGeogToDst )
+    {
+        proj_context_log_debug(ctx,
+            "Cannot create transformation from geographic CRS of target CRS to target CRS");
+        proj_destroy(pjGeogToSrc);
+        return {};
+    }
+
+    try
+    {
+        std::vector<CoordOperation> preparedOpList;
+
+        // Iterate over source->target candidate transformations and reproject
+        // their long-lat bounding box into the source CRS.
+        const auto op_count = proj_list_get_count(op_list);
+        for( int i = 0; i < op_count; i++ )
+        {
+            auto op = proj_list_get(ctx, op_list, i);
+            assert(op);
+            double west_lon = 0.0;
+            double south_lat = 0.0;
+            double east_lon = 0.0;
+            double north_lat = 0.0;
+
+            const char* areaName = nullptr;
+            if( proj_get_area_of_use(ctx, op,
+                        &west_lon, &south_lat, &east_lon, &north_lat, &areaName) )
+            {
+                const bool isOffshore =
+                    areaName && strstr(areaName, "offshore");
+                if( west_lon <= east_lon )
+                {
+                    op = add_coord_op_to_list(i, op,
+                        west_lon, south_lat, east_lon, north_lat,
+                        pjGeogToSrc, pjGeogToDst, isOffshore,
+                        preparedOpList);
+                }
+                else
+                {
+                    auto op_clone = proj_clone(ctx, op);
+
+                    op = add_coord_op_to_list(i, op,
+                        west_lon, south_lat, 180, north_lat,
+                        pjGeogToSrc, pjGeogToDst, isOffshore,
+                        preparedOpList);
+                    op_clone = add_coord_op_to_list(i, op_clone,
+                        -180, south_lat, east_lon, north_lat,
+                        pjGeogToSrc, pjGeogToDst, isOffshore,
+                        preparedOpList);
+                    proj_destroy(op_clone);
+                }
+            }
+
+            proj_destroy(op);
+        }
+
+        proj_destroy(pjGeogToSrc);
+        proj_destroy(pjGeogToDst);
+        return preparedOpList;
+    }
+    catch( const std::exception& )
+    {
+        proj_destroy(pjGeogToSrc);
+        proj_destroy(pjGeogToDst);
+        return {};
+    }
+}
+
 /*****************************************************************************/
 PJ  *proj_create_crs_to_crs_from_pj (PJ_CONTEXT *ctx, const PJ *source_crs, const PJ *target_crs, PJ_AREA *area, const char* const *) {
 /******************************************************************************
@@ -1177,16 +1282,15 @@ PJ  *proj_create_crs_to_crs_from_pj (PJ_CONTEXT *ctx, const PJ *source_crs, cons
             PROJ_GRID_AVAILABILITY_DISCARD_OPERATION_IF_MISSING_GRID);
 
     auto op_list = proj_create_operations(ctx, source_crs, target_crs, operation_ctx);
+    proj_operation_factory_context_destroy(operation_ctx);
 
     if( !op_list ) {
-        proj_operation_factory_context_destroy(operation_ctx);
         return nullptr;
     }
 
     auto op_count = proj_list_get_count(op_list);
     if( op_count == 0 ) {
         proj_list_destroy(op_list);
-        proj_operation_factory_context_destroy(operation_ctx);
 
         proj_context_log_debug(ctx, "No operation found matching criteria");
         return nullptr;
@@ -1199,112 +1303,39 @@ PJ  *proj_create_crs_to_crs_from_pj (PJ_CONTEXT *ctx, const PJ *source_crs, cons
         proj_get_type(source_crs) == PJ_TYPE_GEOCENTRIC_CRS ||
         proj_get_type(target_crs) == PJ_TYPE_GEOCENTRIC_CRS ) {
         proj_list_destroy(op_list);
-        proj_operation_factory_context_destroy(operation_ctx);
         return P;
     }
 
-    auto pjGeogToSrc = create_operation_to_geog_crs(ctx, source_crs);
-    if( !pjGeogToSrc )
+    auto preparedOpList = pj_create_prepared_operations(ctx, source_crs, target_crs,
+                                                   op_list);
+    proj_list_destroy(op_list);
+
+    if( preparedOpList.empty() )
     {
-        proj_list_destroy(op_list);
-        proj_operation_factory_context_destroy(operation_ctx);
-        proj_context_log_debug(ctx,
-            "Cannot create transformation from geographic CRS of source CRS to source CRS");
         proj_destroy(P);
         return nullptr;
     }
 
-    auto pjGeogToDst = create_operation_to_geog_crs(ctx, target_crs);
-    if( !pjGeogToDst )
+    // If there's finally juste a single result, return it directly
+    if( preparedOpList.size() == 1 )
     {
-        proj_list_destroy(op_list);
-        proj_operation_factory_context_destroy(operation_ctx);
-        proj_context_log_debug(ctx,
-            "Cannot create transformation from geographic CRS of target CRS to target CRS");
+        auto retP = preparedOpList[0].pj;
+        preparedOpList[0].pj = nullptr;
         proj_destroy(P);
-        proj_destroy(pjGeogToSrc);
-        return nullptr;
+        return retP;
     }
 
-    try
-    {
-        // Iterate over source->target candidate transformations and reproject
-        // their long-lat bounding box into the source CRS.
-        for( int i = 0; i < op_count; i++ )
-        {
-            auto op = proj_list_get(ctx, op_list, i);
-            assert(op);
-            double west_lon = 0.0;
-            double south_lat = 0.0;
-            double east_lon = 0.0;
-            double north_lat = 0.0;
+    P->alternativeCoordinateOperations = std::move(preparedOpList);
+    // The returned P is rather dummy
+    P->iso_obj = nullptr;
+    P->fwd = nullptr;
+    P->inv = nullptr;
+    P->fwd3d = nullptr;
+    P->inv3d = nullptr;
+    P->fwd4d = nullptr;
+    P->inv4d = nullptr;
 
-            const char* areaName = nullptr;
-            if( proj_get_area_of_use(ctx, op,
-                        &west_lon, &south_lat, &east_lon, &north_lat, &areaName) )
-            {
-                const bool isOffshore =
-                    areaName && strstr(areaName, "offshore");
-                if( west_lon <= east_lon )
-                {
-                    op = add_coord_op_to_list(op,
-                        west_lon, south_lat, east_lon, north_lat,
-                        pjGeogToSrc, pjGeogToDst, isOffshore,
-                        P->alternativeCoordinateOperations);
-                }
-                else
-                {
-                    auto op_clone = proj_clone(ctx, op);
-
-                    op = add_coord_op_to_list(op,
-                        west_lon, south_lat, 180, north_lat,
-                        pjGeogToSrc, pjGeogToDst, isOffshore,
-                        P->alternativeCoordinateOperations);
-                    op_clone = add_coord_op_to_list(op_clone,
-                        -180, south_lat, east_lon, north_lat,
-                        pjGeogToSrc, pjGeogToDst, isOffshore,
-                        P->alternativeCoordinateOperations);
-                    proj_destroy(op_clone);
-                }
-            }
-
-            proj_destroy(op);
-        }
-
-        proj_list_destroy(op_list);
-
-        proj_operation_factory_context_destroy(operation_ctx);
-        proj_destroy(pjGeogToSrc);
-        proj_destroy(pjGeogToDst);
-
-        // If there's finally juste a single result, return it directly
-        if( P->alternativeCoordinateOperations.size() == 1 ) {
-            auto retP = P->alternativeCoordinateOperations[0].pj;
-            P->alternativeCoordinateOperations[0].pj = nullptr;
-            proj_destroy(P);
-            P = retP;
-        } else {
-            // The returned P is rather dummy
-            P->iso_obj = nullptr;
-            P->fwd = nullptr;
-            P->inv = nullptr;
-            P->fwd3d = nullptr;
-            P->inv3d = nullptr;
-            P->fwd4d = nullptr;
-            P->inv4d = nullptr;
-        }
-
-        return P;
-    }
-    catch( const std::exception& )
-    {
-        proj_list_destroy(op_list);
-        proj_operation_factory_context_destroy(operation_ctx);
-        proj_destroy(pjGeogToSrc);
-        proj_destroy(pjGeogToDst);
-        proj_destroy(P);
-        return nullptr;
-    }
+    return P;
 }
 
 PJ *proj_destroy (PJ *P) {
diff --git a/src/iso19111/c_api.cpp b/src/iso19111/c_api.cpp
index 13152662..d7fb35f3 100644
--- a/src/iso19111/c_api.cpp
+++ b/src/iso19111/c_api.cpp
@@ -216,12 +216,17 @@ struct PJ_OBJ_LIST {
 
     explicit PJ_OBJ_LIST(std::vector<IdentifiedObjectNNPtr> &&objectsIn)
         : objects(std::move(objectsIn)) {}
+    virtual ~PJ_OBJ_LIST();
 
     PJ_OBJ_LIST(const PJ_OBJ_LIST &) = delete;
     PJ_OBJ_LIST &operator=(const PJ_OBJ_LIST &) = delete;
     //! @endcond
 };
 
+//! @cond Doxygen_Suppress
+PJ_OBJ_LIST::~PJ_OBJ_LIST() = default;
+//! @endcond
+
 // ---------------------------------------------------------------------------
 
 //! @cond Doxygen_Suppress
@@ -7511,6 +7516,62 @@ void PROJ_DLL proj_operation_factory_context_set_discard_superseded(
 
 // ---------------------------------------------------------------------------
 
+//! @cond Doxygen_Suppress
+/** \brief Opaque object representing a set of operation results. */
+struct PJ_OPERATION_LIST : PJ_OBJ_LIST {
+
+    PJ *source_crs;
+    PJ *target_crs;
+    bool hasPreparedOperation = false;
+    std::vector<CoordOperation> preparedOperations{};
+
+    explicit PJ_OPERATION_LIST(PJ_CONTEXT *ctx, const PJ *source_crsIn,
+                               const PJ *target_crsIn,
+                               std::vector<IdentifiedObjectNNPtr> &&objectsIn);
+    ~PJ_OPERATION_LIST() override;
+
+    PJ_OPERATION_LIST(const PJ_OPERATION_LIST &) = delete;
+    PJ_OPERATION_LIST &operator=(const PJ_OPERATION_LIST &) = delete;
+
+    const std::vector<CoordOperation> &getPreparedOperations(PJ_CONTEXT *ctx);
+};
+
+// ---------------------------------------------------------------------------
+
+PJ_OPERATION_LIST::PJ_OPERATION_LIST(
+    PJ_CONTEXT *ctx, const PJ *source_crsIn, const PJ *target_crsIn,
+    std::vector<IdentifiedObjectNNPtr> &&objectsIn)
+    : PJ_OBJ_LIST(std::move(objectsIn)),
+      source_crs(proj_clone(ctx, source_crsIn)),
+      target_crs(proj_clone(ctx, target_crsIn)) {}
+
+// ---------------------------------------------------------------------------
+
+PJ_OPERATION_LIST::~PJ_OPERATION_LIST() {
+    auto tmpCtxt = proj_context_create();
+    proj_assign_context(source_crs, tmpCtxt);
+    proj_assign_context(target_crs, tmpCtxt);
+    proj_destroy(source_crs);
+    proj_destroy(target_crs);
+    proj_context_destroy(tmpCtxt);
+}
+
+// ---------------------------------------------------------------------------
+
+const std::vector<CoordOperation> &
+PJ_OPERATION_LIST::getPreparedOperations(PJ_CONTEXT *ctx) {
+    if (!hasPreparedOperation) {
+        hasPreparedOperation = true;
+        preparedOperations =
+            pj_create_prepared_operations(ctx, source_crs, target_crs, this);
+    }
+    return preparedOperations;
+}
+
+//! @endcond
+
+// ---------------------------------------------------------------------------
+
 /** \brief Find a list of CoordinateOperation from source_crs to target_crs.
  *
  * The operations are sorted with the most relevant ones first: by
@@ -7562,7 +7623,8 @@ proj_create_operations(PJ_CONTEXT *ctx, const PJ *source_crs,
         for (const auto &op : ops) {
             objects.emplace_back(op);
         }
-        return new PJ_OBJ_LIST(std::move(objects));
+        return new PJ_OPERATION_LIST(ctx, source_crs, target_crs,
+                                     std::move(objects));
     } catch (const std::exception &e) {
         proj_log_error(ctx, __FUNCTION__, e.what());
         return nullptr;
@@ -7571,6 +7633,47 @@ proj_create_operations(PJ_CONTEXT *ctx, const PJ *source_crs,
 
 // ---------------------------------------------------------------------------
 
+/** Return the index of the operation that would be the most appropriate to
+ * transform the specified coordinates.
+ *
+ * This operation may use resources that are not locally available, depending
+ * on the search criteria used by proj_create_operations().
+ *
+ * This could be done by using proj_create_operations() with a punctual bounding
+ * box, but this function is faster when one needs to evaluate on many points
+ * with the same (source_crs, target_crs) tuple.
+ *
+ * @param ctx PROJ context, or NULL for default context
+ * @param operations List of operations returned by proj_create_operations()
+ * @param direction Direction into which to transform the point.
+ * @param coord Coordinate to transform
+ * @return the index in operations that would be used to transform coord. Or -1
+ * in case of error, or no match.
+ *
+ * @since 7.1
+ */
+int proj_get_suggested_operation(PJ_CONTEXT *ctx, PJ_OBJ_LIST *operations,
+                                 PJ_DIRECTION direction, PJ_COORD coord) {
+    SANITIZE_CTX(ctx);
+    auto opList = dynamic_cast<PJ_OPERATION_LIST *>(operations);
+    if (opList == nullptr) {
+        proj_log_error(ctx, __FUNCTION__,
+                       "operations is not a list of operations");
+        return -1;
+    }
+
+    int iExcluded[2] = {-1, -1};
+    const auto &preparedOps = opList->getPreparedOperations(ctx);
+    int idx = pj_get_suggested_operation(ctx, preparedOps, iExcluded, direction,
+                                         coord);
+    if (idx >= 0) {
+        idx = preparedOps[idx].idxInOriginalList;
+    }
+    return idx;
+}
+
+// ---------------------------------------------------------------------------
+
 /** \brief Return the number of objects in the result set
  *
  * @param result Object of type PJ_OBJ_LIST (must not be NULL)
@@ -7892,8 +7995,8 @@ PJ *proj_normalize_for_visualization(PJ_CONTEXT *ctx, const PJ *obj) {
                         }
                     }
                     pjNew->alternativeCoordinateOperations.emplace_back(
-                        minxSrc, minySrc, maxxSrc, maxySrc, minxDst, minyDst,
-                        maxxDst, maxyDst,
+                        alt.idxInOriginalList, minxSrc, minySrc, maxxSrc,
+                        maxySrc, minxDst, minyDst, maxxDst, maxyDst,
                         pj_obj_create(ctx, co->normalizeForVisualization()),
                         co->nameStr(), alt.accuracy, alt.isOffshore);
                 }
diff --git a/src/proj.h b/src/proj.h
index 69aae6d8..4988e1d0 100644
--- a/src/proj.h
+++ b/src/proj.h
@@ -1217,6 +1217,11 @@ PJ PROJ_DLL *proj_list_get(PJ_CONTEXT *ctx,
 
 void PROJ_DLL proj_list_destroy(PJ_OBJ_LIST *result);
 
+int PROJ_DLL proj_get_suggested_operation(PJ_CONTEXT *ctx,
+                                          PJ_OBJ_LIST *operations,
+                                          PJ_DIRECTION direction,
+                                          PJ_COORD coord);
+
 /* ------------------------------------------------------------------------- */
 
 PJ PROJ_DLL *proj_crs_get_geodetic_crs(PJ_CONTEXT *ctx, const PJ *crs);
diff --git a/src/proj_internal.h b/src/proj_internal.h
index 78aff49f..ffc73b06 100644
--- a/src/proj_internal.h
+++ b/src/proj_internal.h
@@ -287,6 +287,55 @@ typedef    PJ_COORD  (* PJ_OPERATOR)    (PJ_COORD, PJ *);
 #define PJD_GRIDSHIFT 3
 #define PJD_WGS84     4   /* WGS84 (or anything considered equivalent) */
 
+struct CoordOperation
+{
+    int idxInOriginalList;
+    double minxSrc = 0.0;
+    double minySrc = 0.0;
+    double maxxSrc = 0.0;
+    double maxySrc = 0.0;
+    double minxDst = 0.0;
+    double minyDst = 0.0;
+    double maxxDst = 0.0;
+    double maxyDst = 0.0;
+    PJ* pj = nullptr;
+    std::string name{};
+    double accuracy = -1.0;
+    bool isOffshore = false;
+
+    CoordOperation(int idxInOriginalListIn,
+                   double minxSrcIn, double minySrcIn, double maxxSrcIn, double maxySrcIn,
+                   double minxDstIn, double minyDstIn, double maxxDstIn, double maxyDstIn,
+                   PJ* pjIn, const std::string& nameIn, double accuracyIn, bool isOffshoreIn):
+        idxInOriginalList(idxInOriginalListIn),
+        minxSrc(minxSrcIn), minySrc(minySrcIn), maxxSrc(maxxSrcIn), maxySrc(maxySrcIn),
+        minxDst(minxDstIn), minyDst(minyDstIn), maxxDst(maxxDstIn), maxyDst(maxyDstIn),
+        pj(pjIn), name(nameIn),
+        accuracy(accuracyIn),
+        isOffshore(isOffshoreIn)
+    {
+    }
+
+    CoordOperation(const CoordOperation&) = delete;
+
+    CoordOperation(CoordOperation&& other):
+        idxInOriginalList(other.idxInOriginalList),
+        minxSrc(other.minxSrc), minySrc(other.minySrc), maxxSrc(other.maxxSrc), maxySrc(other.maxySrc),
+        minxDst(other.minxDst), minyDst(other.minyDst), maxxDst(other.maxxDst), maxyDst(other.maxyDst),
+        name(std::move(other.name)),
+        accuracy(other.accuracy),
+        isOffshore(other.isOffshore) {
+        pj = other.pj;
+        other.pj = nullptr;
+    }
+
+    CoordOperation& operator=(const CoordOperation&) = delete;
+
+    ~CoordOperation()
+    {
+        proj_destroy(pj);
+    }
+};
 
 /* base projection data structure */
 struct PJconsts {
@@ -493,52 +542,6 @@ struct PJconsts {
     /*************************************************************************************
      proj_create_crs_to_crs() alternative coordinate operations
     **************************************************************************************/
-
-    struct CoordOperation
-    {
-        double minxSrc = 0.0;
-        double minySrc = 0.0;
-        double maxxSrc = 0.0;
-        double maxySrc = 0.0;
-        double minxDst = 0.0;
-        double minyDst = 0.0;
-        double maxxDst = 0.0;
-        double maxyDst = 0.0;
-        PJ* pj = nullptr;
-        std::string name{};
-        double accuracy = -1.0;
-        bool isOffshore = false;
-
-        CoordOperation(double minxSrcIn, double minySrcIn, double maxxSrcIn, double maxySrcIn,
-                       double minxDstIn, double minyDstIn, double maxxDstIn, double maxyDstIn,
-                       PJ* pjIn, const std::string& nameIn, double accuracyIn, bool isOffshoreIn):
-            minxSrc(minxSrcIn), minySrc(minySrcIn), maxxSrc(maxxSrcIn), maxySrc(maxySrcIn),
-            minxDst(minxDstIn), minyDst(minyDstIn), maxxDst(maxxDstIn), maxyDst(maxyDstIn),
-            pj(pjIn), name(nameIn),
-            accuracy(accuracyIn),
-            isOffshore(isOffshoreIn)
-        {
-        }
-
-        CoordOperation(const CoordOperation&) = delete;
-
-        CoordOperation(CoordOperation&& other):
-            minxSrc(other.minxSrc), minySrc(other.minySrc), maxxSrc(other.maxxSrc), maxySrc(other.maxySrc),
-            minxDst(other.minxDst), minyDst(other.minyDst), maxxDst(other.maxxDst), maxyDst(other.maxyDst),
-            name(std::move(other.name)),
-            accuracy(other.accuracy),
-            isOffshore(other.isOffshore) {
-            pj = other.pj;
-            other.pj = nullptr;
-        }
-
-        CoordOperation& operator=(const CoordOperation&) = delete;
-
-        ~CoordOperation()
-        {
-            proj_destroy(pj);
-        }
-    };
     std::vector<CoordOperation> alternativeCoordinateOperations{};
     int iCurCoordOp = -1;
 
@@ -873,6 +876,17 @@ std::string PROJ_DLL pj_context_get_user_writable_directory(PJ_CONTEXT *ctx, boo
 void PROJ_DLL pj_context_set_user_writable_directory(PJ_CONTEXT* ctx, const std::string& path);
 std::string PROJ_DLL pj_get_relative_share_proj(PJ_CONTEXT *ctx);
 
+std::vector<CoordOperation> pj_create_prepared_operations(PJ_CONTEXT *ctx,
+                                                     const PJ *source_crs,
+                                                     const PJ *target_crs,
+                                                     PJ_OBJ_LIST* op_list);
+
+int pj_get_suggested_operation(PJ_CONTEXT *ctx,
+                               const std::vector<CoordOperation>& opList,
+                               const int iExcluded[2],
+                               PJ_DIRECTION direction,
+                               PJ_COORD coord);
+
 const PJ_UNITS *pj_list_linear_units();
 
 /* classic public API */