Generalize generalize_proj_crs_create_bound_vertical_crs_to_WGS84()

In recent commits, we added a generalize_proj_crs_create_bound_vertical_crs_to_WGS84()
function, but there are situations where more accurate results can be obtained, if
instead of specifying WGS84 as the hub CRS, the user can specify the exact hub CRS.
For example the GEOID2018 grid is against NAD83(2011).
So replace this function with proj_crs_create_bound_vertical_crs()

2 files changed, 26 insertions, 14 deletions

diff --git a/src/iso19111/c_api.cpp b/src/iso19111/c_api.cpp
index 4fedbe05..f5f7ba55 100644
--- a/src/iso19111/c_api.cpp
+++ b/src/iso19111/c_api.cpp
@@ -1810,7 +1810,8 @@ PJ *proj_crs_create_bound_crs_to_WGS84(PJ_CONTEXT *ctx, const PJ *crs,
 
 // ---------------------------------------------------------------------------
 
-/** \brief Returns a BoundCRS, with a transformation to EPSG:4979 using a grid.
+/** \brief Returns a BoundCRS, with a transformation to a hub geographic 3D crs
+ * (use EPSG:4979 for WGS84 for example), using a grid.
  *
  * The returned object must be unreferenced with proj_destroy() after
  * use.
@@ -1818,33 +1819,43 @@ PJ *proj_crs_create_bound_crs_to_WGS84(PJ_CONTEXT *ctx, const PJ *crs,
  *
  * @param ctx PROJ context, or NULL for default context
  * @param vert_crs Object of type VerticalCRS (must not be NULL)
+ * @param hub_geographic_3D_crs Object of type Geographic 3D CRS (must not be
+ * NULL)
  * @param grid_name Grid name (typically a .gtx file)
  * @return Object that must be unreferenced with proj_destroy(), or NULL
  * in case of error.
  * @since 7.0
  */
-PJ *proj_crs_create_bound_vertical_crs_to_WGS84(PJ_CONTEXT *ctx,
-                                                const PJ *vert_crs,
-                                                const char *grid_name) {
+PJ *proj_crs_create_bound_vertical_crs(PJ_CONTEXT *ctx, const PJ *vert_crs,
+                                       const PJ *hub_geographic_3D_crs,
+                                       const char *grid_name) {
     SANITIZE_CTX(ctx);
     assert(vert_crs);
+    assert(hub_geographic_3D_crs);
     assert(grid_name);
     auto l_crs = std::dynamic_pointer_cast<VerticalCRS>(vert_crs->iso_obj);
     if (!l_crs) {
-        proj_log_error(ctx, __FUNCTION__, "Object is not a VerticalCRS");
+        proj_log_error(ctx, __FUNCTION__, "vert_crs is not a VerticalCRS");
+        return nullptr;
+    }
+    auto hub_crs =
+        std::dynamic_pointer_cast<CRS>(hub_geographic_3D_crs->iso_obj);
+    if (!hub_crs) {
+        proj_log_error(ctx, __FUNCTION__, "hub_geographic_3D_crs is not a CRS");
         return nullptr;
     }
     try {
         auto nnCRS = NN_NO_CHECK(l_crs);
+        auto nnHubCRS = NN_NO_CHECK(hub_crs);
         auto transformation =
             Transformation::createGravityRelatedHeightToGeographic3D(
                 PropertyMap().set(IdentifiedObject::NAME_KEY,
-                                  "unknown to WGS84 ellipsoidal height"),
-                nnCRS, GeographicCRS::EPSG_4979, nullptr,
-                std::string(grid_name), std::vector<PositionalAccuracyNNPtr>());
-        return pj_obj_create(
-            ctx,
-            BoundCRS::create(nnCRS, GeographicCRS::EPSG_4979, transformation));
+                                  "unknown to " + hub_crs->nameStr() +
+                                      " ellipsoidal height"),
+                nnCRS, nnHubCRS, nullptr, std::string(grid_name),
+                std::vector<PositionalAccuracyNNPtr>());
+        return pj_obj_create(ctx,
+                             BoundCRS::create(nnCRS, nnHubCRS, transformation));
     } catch (const std::exception &e) {
         proj_log_error(ctx, __FUNCTION__, e.what());
         if (ctx->cpp_context) {
diff --git a/src/proj_experimental.h b/src/proj_experimental.h
index 63b858d5..c6d5bc45 100644
--- a/src/proj_experimental.h
+++ b/src/proj_experimental.h
@@ -317,9 +317,10 @@ PJ PROJ_DLL *proj_crs_create_bound_crs_to_WGS84(PJ_CONTEXT *ctx,
                                                         const PJ *crs,
                                                         const char *const *options);
 
-PJ PROJ_DLL *proj_crs_create_bound_vertical_crs_to_WGS84(PJ_CONTEXT *ctx,
-                                                         const PJ* vert_crs,
-                                                         const char* grid_name);
+PJ PROJ_DLL *proj_crs_create_bound_vertical_crs(PJ_CONTEXT *ctx,
+                                                const PJ* vert_crs,
+                                                const PJ* hub_geographic_3D_crs,
+                                                const char* grid_name);
 
 /* BEGIN: Generated by scripts/create_c_api_projections.py*/
 PJ PROJ_DLL *proj_create_conversion_utm(