Make proj_create_crs_to_crs() use proj_obj_create_operations() and use area of use argument, and make createFromUserInput() recognize init=epsg: / init=IGNF: in legacy mode, that is when proj_context_get_use_proj4_init_rules() is used

10 files changed, 471 insertions, 118 deletions

diff --git a/include/proj/io.hpp b/include/proj/io.hpp
index c649fa9f..3b5019c1 100644
--- a/include/proj/io.hpp
+++ b/include/proj/io.hpp
@@ -597,7 +597,8 @@ class PROJ_GCC_DLL WKTNode {
 
 PROJ_DLL util::BaseObjectNNPtr
 createFromUserInput(const std::string &text,
-                    const DatabaseContextPtr &dbContext);
+                    const DatabaseContextPtr &dbContext,
+                    bool usePROJ4InitRules = false);
 
 // ---------------------------------------------------------------------------
 
@@ -654,6 +655,8 @@ class PROJ_GCC_DLL PROJStringParser {
     PROJ_DLL PROJStringParser &
     attachDatabaseContext(const DatabaseContextPtr &dbContext);
 
+    PROJ_DLL PROJStringParser &setUsePROJ4InitRules(bool enable);
+
     PROJ_DLL std::vector<std::string> warningList() const;
 
     PROJ_DLL util::BaseObjectNNPtr createFromPROJString(
diff --git a/src/c_api.cpp b/src/c_api.cpp
index fbba0f51..ba1b9534 100644
--- a/src/c_api.cpp
+++ b/src/c_api.cpp
@@ -276,6 +276,18 @@ PJ_GUESSED_WKT_DIALECT proj_context_guess_wkt_dialect(PJ_CONTEXT *ctx,
 
 // ---------------------------------------------------------------------------
 
+//! @cond Doxygen_Suppress
+static const char *getOptionValue(const char *option,
+                                  const char *keyWithEqual) noexcept {
+    if (ci_starts_with(option, keyWithEqual)) {
+        return option + strlen(keyWithEqual);
+    }
+    return nullptr;
+}
+//! @endcond
+
+// ---------------------------------------------------------------------------
+
 /** \brief Instanciate an object from a WKT string, PROJ string or object code
  * (like "EPSG:4326", "urn:ogc:def:crs:EPSG::4326",
  * "urn:ogc:def:coordinateOperation:EPSG::1671").
@@ -287,7 +299,17 @@ PJ_GUESSED_WKT_DIALECT proj_context_guess_wkt_dialect(PJ_CONTEXT *ctx,
  *
  * @param ctx PROJ context, or NULL for default context
  * @param text String (must not be NULL)
- * @param options should be set to NULL for now
+ * @param options null-terminated list of options, or NULL. Currently
+ * supported options are:
+ * <ul>
+ * <li>USE_PROJ4_INIT_RULES=YES/NO. Defaults to NO. When set to YES,
+ * init=epsg:XXXX syntax will be allowed and will be interpreted according to
+ * PROJ.4 and PROJ.5 rules, that is geodeticCRS will have longitude, latitude
+ * order and will expect/output coordinates in radians. ProjectedCRS will have
+ * easting, northing axis order (except the ones with Transverse Mercator South
+ * Orientated projection). In that mode, the epsg:XXXX syntax will be also
+ * interprated the same way.</li>
+ * </ul>
  * @return Object that must be unreferenced with proj_obj_unref(), or NULL in
  * case of error.
  */
@@ -298,8 +320,20 @@ PJ_OBJ *proj_obj_create_from_user_input(PJ_CONTEXT *ctx, const char *text,
     (void)options;
     auto dbContext = getDBcontextNoException(ctx, __FUNCTION__);
     try {
+        bool usePROJ4InitRules = false;
+        for (auto iter = options; iter && iter[0]; ++iter) {
+            const char *value;
+            if ((value = getOptionValue(*iter, "USE_PROJ4_INIT_RULES="))) {
+                usePROJ4InitRules = ci_equal(value, "YES");
+            } else {
+                std::string msg("Unknown option :");
+                msg += *iter;
+                proj_log_error(ctx, __FUNCTION__, msg.c_str());
+                return nullptr;
+            }
+        }
         auto identifiedObject = nn_dynamic_pointer_cast<IdentifiedObject>(
-            createFromUserInput(text, dbContext));
+            createFromUserInput(text, dbContext, usePROJ4InitRules));
         if (identifiedObject) {
             return PJ_OBJ::create(ctx, NN_NO_CHECK(identifiedObject));
         }
@@ -809,18 +843,6 @@ const char *proj_obj_get_id_code(PJ_OBJ *obj, int index) {
 
 // ---------------------------------------------------------------------------
 
-//! @cond Doxygen_Suppress
-static const char *getOptionValue(const char *option,
-                                  const char *keyWithEqual) noexcept {
-    if (ci_starts_with(option, keyWithEqual)) {
-        return option + strlen(keyWithEqual);
-    }
-    return nullptr;
-}
-//! @endcond
-
-// ---------------------------------------------------------------------------
-
 /** \brief Get a WKT representation of an object.
  *
  * The returned string is valid while the input obj parameter is valid,
@@ -986,26 +1008,28 @@ const char *proj_obj_as_proj_string(PJ_OBJ *obj, PJ_PROJ_STRING_TYPE type,
 /** \brief Return the area of use of an object.
  *
  * @param obj Object (must not be NULL)
- * @param p_west_lon Pointer to a double to receive the west longitude (in
- * degrees). Or NULL. If the returned value is -1000, the bounding box is
+ * @param p_west_lon_degree Pointer to a double to receive the west longitude
+ * (in degrees). Or NULL. If the returned value is -1000, the bounding box is
  * unknown.
- * @param p_south_lat Pointer to a double to receive the south latitude (in
- * degrees). Or NULL. If the returned value is -1000, the bounding box is
+ * @param p_south_lat_degree Pointer to a double to receive the south latitude
+ * (in degrees). Or NULL. If the returned value is -1000, the bounding box is
  * unknown.
-  * @param p_east_lon Pointer to a double to receive the east longitude (in
- * degrees). Or NULL. If the returned value is -1000, the bounding box is
+ * @param p_east_lon_degree Pointer to a double to receive the east longitude
+ * (in degrees). Or NULL. If the returned value is -1000, the bounding box is
  * unknown.
- * @param p_north_lat Pointer to a double to receive the north latitude (in
- * degrees). Or NULL. If the returned value is -1000, the bounding box is
+ * @param p_north_lat_degree Pointer to a double to receive the north latitude
+ * (in degrees). Or NULL. If the returned value is -1000, the bounding box is
  * unknown.
  * @param p_area_name Pointer to a string to receive the name of the area of
  * use. Or NULL. *p_area_name is valid while obj is valid itself.
  * @return TRUE in case of success, FALSE in case of error or if the area
  * of use is unknown.
  */
-int proj_obj_get_area_of_use(PJ_OBJ *obj, double *p_west_lon,
-                             double *p_south_lat, double *p_east_lon,
-                             double *p_north_lat, const char **p_area_name) {
+int proj_obj_get_area_of_use(PJ_OBJ *obj, double *p_west_lon_degree,
+                             double *p_south_lat_degree,
+                             double *p_east_lon_degree,
+                             double *p_north_lat_degree,
+                             const char **p_area_name) {
     if (p_area_name) {
         *p_area_name = nullptr;
     }
@@ -1031,32 +1055,32 @@ int proj_obj_get_area_of_use(PJ_OBJ *obj, double *p_west_lon,
         auto bbox =
             dynamic_cast<const GeographicBoundingBox *>(geogElements[0].get());
         if (bbox) {
-            if (p_west_lon) {
-                *p_west_lon = bbox->westBoundLongitude();
+            if (p_west_lon_degree) {
+                *p_west_lon_degree = bbox->westBoundLongitude();
             }
-            if (p_south_lat) {
-                *p_south_lat = bbox->southBoundLatitude();
+            if (p_south_lat_degree) {
+                *p_south_lat_degree = bbox->southBoundLatitude();
             }
-            if (p_east_lon) {
-                *p_east_lon = bbox->eastBoundLongitude();
+            if (p_east_lon_degree) {
+                *p_east_lon_degree = bbox->eastBoundLongitude();
             }
-            if (p_north_lat) {
-                *p_north_lat = bbox->northBoundLatitude();
+            if (p_north_lat_degree) {
+                *p_north_lat_degree = bbox->northBoundLatitude();
             }
             return true;
         }
     }
-    if (p_west_lon) {
-        *p_west_lon = -1000;
+    if (p_west_lon_degree) {
+        *p_west_lon_degree = -1000;
     }
-    if (p_south_lat) {
-        *p_south_lat = -1000;
+    if (p_south_lat_degree) {
+        *p_south_lat_degree = -1000;
     }
-    if (p_east_lon) {
-        *p_east_lon = -1000;
+    if (p_east_lon_degree) {
+        *p_east_lon_degree = -1000;
     }
-    if (p_north_lat) {
-        *p_north_lat = -1000;
+    if (p_north_lat_degree) {
+        *p_north_lat_degree = -1000;
     }
     return true;
 }
@@ -3743,21 +3767,21 @@ void proj_operation_factory_context_set_desired_accuracy(
 /** \brief Set the desired area of interest for the resulting coordinate
  * transformations.
  *
- * For an area of interest crossing the anti-meridian, west_lon will be
- * greater than east_lon.
+ * For an area of interest crossing the anti-meridian, west_lon_degree will be
+ * greater than east_lon_degree.
  *
  * @param ctxt Operation factory context. must not be NULL
- * @param west_lon West longitude (in degrees).
- * @param south_lat South latitude (in degrees).
- * @param east_lon East longitude (in degrees).
- * @param north_lat North latitude (in degrees).
+ * @param west_lon_degree West longitude (in degrees).
+ * @param south_lat_degree South latitude (in degrees).
+ * @param east_lon_degree East longitude (in degrees).
+ * @param north_lat_degree North latitude (in degrees).
  */
 void proj_operation_factory_context_set_area_of_interest(
-    PJ_OPERATION_FACTORY_CONTEXT *ctxt, double west_lon, double south_lat,
-    double east_lon, double north_lat) {
+    PJ_OPERATION_FACTORY_CONTEXT *ctxt, double west_lon_degree,
+    double south_lat_degree, double east_lon_degree, double north_lat_degree) {
     assert(ctxt);
-    ctxt->operationContext->setAreaOfInterest(
-        Extent::createFromBBOX(west_lon, south_lat, east_lon, north_lat));
+    ctxt->operationContext->setAreaOfInterest(Extent::createFromBBOX(
+        west_lon_degree, south_lat_degree, east_lon_degree, north_lat_degree));
 }
 
 // ---------------------------------------------------------------------------
diff --git a/src/io.cpp b/src/io.cpp
index c2ca484d..0d220e13 100644
--- a/src/io.cpp
+++ b/src/io.cpp
@@ -4066,10 +4066,22 @@ BaseObjectNNPtr WKTParser::Private::build(const WKTNodeNNPtr &node) {
  *     determine the appropriate best match.</li>
  * </ul>
  *
+ * @param text One of the above mentionned text format
+ * @param dbContext Database context, or nullptr (in which case database
+ * lookups will not work)
+ * @param usePROJ4InitRules When set to true,
+ * init=epsg:XXXX syntax will be allowed and will be interpreted according to
+ * PROJ.4 and PROJ.5 rules, that is geodeticCRS will have longitude, latitude
+ * order and will expect/output coordinates in radians. ProjectedCRS will have
+ * easting, northing axis order (except the ones with Transverse Mercator South
+ * Orientated projection). In that mode, the epsg:XXXX syntax will be also
+ * interprated the same way.
  * @throw ParsingException
  */
 BaseObjectNNPtr createFromUserInput(const std::string &text,
-                                    const DatabaseContextPtr &dbContext) {
+                                    const DatabaseContextPtr &dbContext,
+                                    bool usePROJ4InitRules) {
+
     for (const auto &wktConstants : WKTConstants::constants()) {
         if (ci_starts_with(text, wktConstants)) {
             return WKTParser().attachDatabaseContext(dbContext).createFromWKT(
@@ -4085,6 +4097,7 @@ BaseObjectNNPtr createFromUserInput(const std::string &text,
         strncmp(textWithoutPlusPrefix, "title=", strlen("title=")) == 0) {
         return PROJStringParser()
             .attachDatabaseContext(dbContext)
+            .setUsePROJ4InitRules(usePROJ4InitRules)
             .createFromPROJString(text);
     }
 
@@ -5293,6 +5306,7 @@ const DatabaseContextPtr &PROJStringFormatter::databaseContext() const {
 
 struct PROJStringParser::Private {
     DatabaseContextPtr dbContext_{};
+    bool usePROJ4InitRules_ = false;
     std::vector<std::string> warningList_{};
 
     std::string projString_{};
@@ -5400,6 +5414,22 @@ PROJStringParser::attachDatabaseContext(const DatabaseContextPtr &dbContext) {
 
 // ---------------------------------------------------------------------------
 
+/** \brief Set how init=epsg:XXXX syntax should be interpreted.
+ *
+ * @param enable When set to true,
+ * init=epsg:XXXX syntax will be allowed and will be interpreted according to
+ * PROJ.4 and PROJ.5 rules, that is geodeticCRS will have longitude, latitude
+ * order and will expect/output coordinates in radians. ProjectedCRS will have
+ * easting, northing axis order (except the ones with Transverse Mercator South
+ * Orientated projection).
+ */
+PROJStringParser &PROJStringParser::setUsePROJ4InitRules(bool enable) {
+    d->usePROJ4InitRules_ = enable;
+    return *this;
+}
+
+// ---------------------------------------------------------------------------
+
 /** \brief Return the list of warnings found during parsing.
  */
 std::vector<std::string> PROJStringParser::warningList() const {
@@ -6817,6 +6847,9 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
     std::string vunits;
     std::string vto_meter;
 
+    d->steps_.clear();
+    d->title_.clear();
+    d->globalParamValues_.clear();
     d->projString_ = projString;
     PROJStringSyntaxParser(projString, d->steps_, d->globalParamValues_,
                            d->title_, vunits, vto_meter);
@@ -6853,6 +6886,99 @@ PROJStringParser::createFromPROJString(const std::string &projString) {
                                             : -1));
     }
 
+    // +init=xxxx:yyyy syntax
+    if (d->steps_.size() == 1 && d->steps_[0].isInit &&
+        d->steps_[0].paramValues.size() == 0) {
+
+        // Those used to come from a text init file
+        // We only support them in compatibility mode
+        if (ci_starts_with(d->steps_[0].name, "epsg:") ||
+            ci_starts_with(d->steps_[0].name, "IGNF:")) {
+            bool usePROJ4InitRules = d->usePROJ4InitRules_;
+            if (!usePROJ4InitRules) {
+                PJ_CONTEXT *ctx = proj_context_create();
+                if (ctx) {
+                    usePROJ4InitRules =
+                        proj_context_get_use_proj4_init_rules(ctx) == TRUE;
+                    proj_context_destroy(ctx);
+                }
+            }
+            if (!usePROJ4InitRules) {
+                throw ParsingException("init=epsg:/init=IGNF: syntax not "
+                                       "supported in non-PROJ4 emulation mode");
+            }
+            auto obj =
+                createFromUserInput(d->steps_[0].name, d->dbContext_, true);
+            auto geogCRS = dynamic_cast<GeographicCRS *>(obj.get());
+            if (geogCRS) {
+                // Override with longitude latitude in radian
+                return GeographicCRS::create(
+                    PropertyMap().set(IdentifiedObject::NAME_KEY,
+                                      geogCRS->nameStr()),
+                    geogCRS->datum(), geogCRS->datumEnsemble(),
+                    EllipsoidalCS::createLongitudeLatitude(
+                        UnitOfMeasure::RADIAN));
+            }
+            auto projCRS = dynamic_cast<ProjectedCRS *>(obj.get());
+            if (projCRS) {
+                // Override with easting northing orer
+                auto conv = projCRS->derivingConversionRef();
+                if (conv->method()->getEPSGCode() !=
+                    EPSG_CODE_METHOD_TRANSVERSE_MERCATOR_SOUTH_ORIENTATED) {
+                    return ProjectedCRS::create(
+                        PropertyMap().set(IdentifiedObject::NAME_KEY,
+                                          projCRS->nameStr()),
+                        projCRS->baseCRS(), conv,
+                        CartesianCS::createEastingNorthing(
+                            projCRS->coordinateSystem()
+                                ->axisList()[0]
+                                ->unit()));
+                }
+            }
+            return obj;
+        }
+
+        paralist *init = pj_mkparam(("init=" + d->steps_[0].name).c_str());
+        if (!init) {
+            throw ParsingException("out of memory");
+        }
+        PJ_CONTEXT *ctx = proj_context_create();
+        if (!ctx) {
+            pj_dealloc(init);
+            throw ParsingException("out of memory");
+        }
+        paralist *list = pj_expand_init(ctx, init);
+        proj_context_destroy(ctx);
+        if (!list) {
+            pj_dealloc(init);
+            throw ParsingException("cannot expand " + projString);
+        }
+        std::string expanded;
+        bool first = true;
+        bool has_init_term = false;
+        for (auto t = list; t;) {
+            if (!expanded.empty()) {
+                expanded += ' ';
+            }
+            if (first) {
+                // first parameter is the init= itself
+                first = false;
+            } else if (starts_with(t->param, "init=")) {
+                has_init_term = true;
+            } else {
+                expanded += t->param;
+            }
+
+            auto n = t->next;
+            pj_dealloc(t);
+            t = n;
+        }
+
+        if (!has_init_term) {
+            return createFromPROJString(expanded);
+        }
+    }
+
     int iFirstGeogStep = -1;
     int iSecondGeogStep = -1;
     int iProjStep = -1;
diff --git a/src/pj_ctx.c b/src/pj_ctx.c
index 6626d5ee..54e2cfb7 100644
--- a/src/pj_ctx.c
+++ b/src/pj_ctx.c
@@ -85,6 +85,7 @@ projCtx pj_get_default_ctx()
         default_context.app_data = NULL;
         default_context.fileapi = pj_get_default_fileapi();
         default_context.cpp_context = NULL;
+        default_context.use_proj4_init_rules = FALSE;
 
         if( getenv("PROJ_DEBUG") != NULL )
         {
@@ -114,6 +115,7 @@ projCtx pj_ctx_alloc()
     memcpy( ctx, pj_get_default_ctx(), sizeof(projCtx_t) );
     ctx->last_errno = 0;
     ctx->cpp_context = NULL;
+    ctx->use_proj4_init_rules = FALSE;
 
     return ctx;
 }
diff --git a/src/proj.h b/src/proj.h
index 39d16b30..4b599eba 100644
--- a/src/proj.h
+++ b/src/proj.h
@@ -335,6 +335,8 @@ typedef struct projCtx_t PJ_CONTEXT;
 PJ_CONTEXT PROJ_DLL *proj_context_create (void);
 PJ_CONTEXT PROJ_DLL *proj_context_destroy (PJ_CONTEXT *ctx);
 
+void PROJ_DLL proj_context_use_proj4_init_rules(PJ_CONTEXT *ctx, int enable);
+int PROJ_DLL proj_context_get_use_proj4_init_rules(PJ_CONTEXT *ctx);
 
 /* Manage the transformation definition object PJ */
 PJ PROJ_DLL *proj_create (PJ_CONTEXT *ctx, const char *definition);
@@ -342,11 +344,14 @@ PJ PROJ_DLL *proj_create_argv (PJ_CONTEXT *ctx, int argc, char **argv);
 PJ PROJ_DLL *proj_create_crs_to_crs(PJ_CONTEXT *ctx, const char *srid_from, const char *srid_to, PJ_AREA *area);
 PJ PROJ_DLL *proj_destroy (PJ *P);
 
-/* Setter-functions for the opaque PJ_AREA struct */
-/* Uncomment these when implementing support for area-based transformations.
-void proj_area_bbox(PJ_AREA *area, LP ll, LP ur);
-void proj_area_description(PJ_AREA *area, const char *descr);
-*/
+
+PJ_AREA PROJ_DLL *proj_area_create(void);
+void PROJ_DLL proj_area_set_bbox(PJ_AREA *area,
+                                 double west_lon_degree,
+                                 double south_lat_degree,
+                                 double east_lon_degree,
+                                 double north_lat_degree);
+void PROJ_DLL proj_area_destroy(PJ_AREA* area);
 
 /* Apply transformation to observation - in forward or inverse direction */
 enum PJ_DIRECTION {
@@ -1194,10 +1199,10 @@ const char PROJ_DLL* proj_obj_get_id_auth_name(PJ_OBJ *obj, int index);
 const char PROJ_DLL* proj_obj_get_id_code(PJ_OBJ *obj, int index);
 
 int PROJ_DLL proj_obj_get_area_of_use(PJ_OBJ *obj,
-                                      double* p_west_lon,
-                                      double* p_south_lat,
-                                      double* p_east_lon,
-                                      double* p_north_lat,
+                                      double* p_west_lon_degree,
+                                      double* p_south_lat_degree,
+                                      double* p_east_lon_degree,
+                                      double* p_north_lat_degree,
                                       const char **p_area_name);
 
 /** \brief WKT version. */
@@ -1278,10 +1283,10 @@ void PROJ_DLL proj_operation_factory_context_set_desired_accuracy(
 
 void PROJ_DLL proj_operation_factory_context_set_area_of_interest(
                                             PJ_OPERATION_FACTORY_CONTEXT *ctxt,
-                                            double west_lon,
-                                            double south_lat,
-                                            double east_lon,
-                                            double north_lat);
+                                            double west_lon_degree,
+                                            double south_lat_degree,
+                                            double east_lon_degree,
+                                            double north_lat_degree);
 
 /** Specify how source and target CRS extent should be used to restrict
   * candidate operations (only taken into account if no explicit area of
diff --git a/src/proj_4D_api.c b/src/proj_4D_api.c
index 56694aae..72e1a2d6 100644
--- a/src/proj_4D_api.c
+++ b/src/proj_4D_api.c
@@ -32,6 +32,9 @@
 #include <stdio.h>
 #include <stdlib.h>
 #include <string.h>
+#ifndef _MSC_VER
+#include <strings.h>
+#endif
 
 #include "proj.h"
 #include "proj_internal.h"
@@ -634,6 +637,75 @@ indicator, as in {"+proj=utm", "+zone=32"}, or leave it out, as in {"proj=utm",
     return P;
 }
 
+/** Create an area of use */
+PJ_AREA * proj_area_create(void) {
+    return pj_calloc(1, sizeof(PJ_AREA));
+}
+
+/** Assign a bounding box to an area of use. */
+void proj_area_set_bbox(PJ_AREA *area,
+                                 double west_lon_degree,
+                                 double south_lat_degree,
+                                 double east_lon_degree,
+                                 double north_lat_degree) {
+    area->bbox_set = TRUE;
+    area->west_lon_degree = west_lon_degree;
+    area->south_lat_degree = south_lat_degree;
+    area->east_lon_degree = east_lon_degree;
+    area->north_lat_degree = north_lat_degree;
+}
+
+/** Free an area of use */
+void proj_area_destroy(PJ_AREA* area) {
+    pj_dealloc(area);
+}
+
+/************************************************************************/
+/*                  proj_context_use_proj4_init_rules()                 */
+/************************************************************************/
+
+void proj_context_use_proj4_init_rules(PJ_CONTEXT *ctx, int enable) {
+    if( ctx == NULL ) {
+        ctx = pj_get_default_ctx();
+    }
+    ctx->use_proj4_init_rules = enable;
+}
+
+/************************************************************************/
+/*                              EQUAL()                                 */
+/************************************************************************/
+
+static int EQUAL(const char* a, const char* b) {
+#ifdef _MSC_VER
+    return _stricmp(a, b) == 0;
+#else
+    return strcasecmp(a, b) == 0;
+#endif
+}
+
+/************************************************************************/
+/*                  proj_context_get_use_proj4_init_rules()             */
+/************************************************************************/
+
+int proj_context_get_use_proj4_init_rules(PJ_CONTEXT *ctx) {
+    const char* val = getenv("PROJ_USE_PROJ4_INIT_RULES");
+
+    if( ctx == NULL ) {
+        ctx = pj_get_default_ctx();
+    }
+
+    if( val ) {
+        if( EQUAL(val, "yes") || EQUAL(val, "on") || EQUAL(val, "true") ) {
+            return TRUE;
+        }
+        if( EQUAL(val, "no") || EQUAL(val, "off") || EQUAL(val, "false") ) {
+            return FALSE;
+        }
+        pj_log(ctx, PJ_LOG_ERROR, "Invalid value for PROJ_USE_PROJ4_INIT_RULES");
+    }
+
+    return ctx->use_proj4_init_rules;
+}
 
 
 /*****************************************************************************/
@@ -643,42 +715,88 @@ PJ  *proj_create_crs_to_crs (PJ_CONTEXT *ctx, const char *srid_from, const char
     systems.
 
     srid_from and srid_to should be the value part of a +init=... parameter
-    set, i.e. "epsg:25833" or "IGNF:AMST63". Any projection definition that
+    set, i.e. "EPSG:25833" or "IGNF:AMST63". Any projection definition that
     can be found in a init-file in PROJ_LIB is a valid input to this function.
 
-    For now the function mimics the cs2cs app: An input and an output CRS is
-    given and coordinates are transformed via a hub datum (WGS84). This
-    transformation strategy is referred to as "early-binding" by the EPSG. The
-    function can be extended to support "late-binding" transformations in the
-    future without affecting users of the function.
-
-    An "area of use" can be specified in area. In the current version of this
-    function is has no function, but is added in anticipation of a
-    "late-binding" implementation in the future. The idea being, that if a user
-    supplies an area of use, the more accurate transformation between two given
-    systems can be chosen.
+    An "area of use" can be specified in area. When it is supplied, the more
+    accurate transformation between two given systems can be chosen.
 
     Example call:
 
-        PJ *P = proj_create_crs_to_crs(0, "epsg:25832", "epsg:25833", NULL);
+        PJ *P = proj_create_crs_to_crs(0, "EPSG:25832", "EPSG:25833", NULL);
 
 ******************************************************************************/
     PJ *P;
-    char buffer[512];
-    size_t len;
+    PJ_OBJ* src;
+    PJ_OBJ* dst;
+    PJ_OPERATION_FACTORY_CONTEXT* operation_ctx;
+    PJ_OBJ_LIST* op_list;
+    PJ_OBJ* op;
+    const char* proj_string;
+    const char* const optionsProj4Mode[] = { "USE_PROJ4_INIT_RULES=YES", NULL };
+    const char* const* optionsImportCRS =
+        proj_context_get_use_proj4_init_rules(ctx) ? optionsProj4Mode : NULL;
+
+    src = proj_obj_create_from_user_input(ctx, srid_from, optionsImportCRS);
+    if( !src ) {
+        return NULL;
+    }
+
+    dst = proj_obj_create_from_user_input(ctx, srid_to, optionsImportCRS);
+    if( !dst ) {
+        proj_obj_unref(src);
+        return NULL;
+    }
+
+    operation_ctx = proj_create_operation_factory_context(ctx, NULL);
+    if( !operation_ctx ) {
+        proj_obj_unref(src);
+        proj_obj_unref(dst);
+        return NULL;
+    }
+
+    if( area && area->bbox_set ) {
+        proj_operation_factory_context_set_area_of_interest(
+                                            operation_ctx,
+                                            area->west_lon_degree,
+                                            area->south_lat_degree,
+                                            area->east_lon_degree,
+                                            area->north_lat_degree);
+    }
+
+    proj_operation_factory_context_set_grid_availability_use(
+        operation_ctx, PROJ_GRID_AVAILABILITY_DISCARD_OPERATION_IF_MISSING_GRID);
+
+    op_list = proj_obj_create_operations(src, dst, operation_ctx);
+
+    proj_operation_factory_context_unref(operation_ctx);
+    proj_obj_unref(src);
+    proj_obj_unref(dst);
+
+    if( !op_list ) {
+        return NULL;
+    }
 
-    /* area not in use yet, suppressing warning */
-    (void)area;
+    if( proj_obj_list_get_count(op_list) == 0 ) {
+        proj_obj_list_unref(op_list);
+        return NULL;
+    }
+
+    op = proj_obj_list_get(op_list, 0);
+    proj_obj_list_unref(op_list);
+    if( !op ) {
+        return NULL;
+    }
+
+    proj_string = proj_obj_as_proj_string(op, PJ_PROJ_5, NULL);
+    if( !proj_string) {
+        proj_obj_unref(op);
+        return NULL;
+    }
 
-    strcpy(buffer, "+proj=pipeline +step +init=");
-    len = strlen(buffer);
-    strncat(buffer + len, srid_from, sizeof(buffer)-1-len);
-    len += strlen(buffer + len);
-    strncat(buffer + len, " +inv +step +init=", sizeof(buffer)-1-len);
-    len += strlen(buffer + len);
-    strncat(buffer + len, srid_to, sizeof(buffer)-1-len);
+    P = proj_create(ctx, proj_string);
 
-    P = proj_create(ctx, buffer);
+    proj_obj_unref(op);
 
     return P;
 }
diff --git a/src/proj_symbol_rename.h b/src/proj_symbol_rename.h
index d0dd09eb..353473b5 100644
--- a/src/proj_symbol_rename.h
+++ b/src/proj_symbol_rename.h
@@ -101,14 +101,19 @@
 #define pj_transform internal_pj_transform
 #define proj_angular_input internal_proj_angular_input
 #define proj_angular_output internal_proj_angular_output
+#define proj_area_create internal_proj_area_create
+#define proj_area_destroy internal_proj_area_destroy
+#define proj_area_set_bbox internal_proj_area_set_bbox
 #define proj_context_create internal_proj_context_create
 #define proj_context_delete_cpp_context internal_proj_context_delete_cpp_context
 #define proj_context_destroy internal_proj_context_destroy
 #define proj_context_errno internal_proj_context_errno
 #define proj_context_get_database_path internal_proj_context_get_database_path
+#define proj_context_get_use_proj4_init_rules internal_proj_context_get_use_proj4_init_rules
 #define proj_context_guess_wkt_dialect internal_proj_context_guess_wkt_dialect
 #define proj_context_set internal_proj_context_set
 #define proj_context_set_database_path internal_proj_context_set_database_path
+#define proj_context_use_proj4_init_rules internal_proj_context_use_proj4_init_rules
 #define proj_coord internal_proj_coord
 #define proj_coord_error internal_proj_coord_error
 #define proj_coordoperation_get_accuracy internal_proj_coordoperation_get_accuracy
@@ -132,7 +137,6 @@
 #define proj_factors internal_proj_factors
 #define proj_free_int_list internal_proj_free_int_list
 #define proj_free_string_list internal_proj_free_string_list
-#define proj_geocentric_latitude internal_proj_geocentric_latitude
 #define proj_geod internal_proj_geod
 #define proj_get_authorities_from_database internal_proj_get_authorities_from_database
 #define proj_get_codes_from_database internal_proj_get_codes_from_database
diff --git a/src/projects.h b/src/projects.h
index e34fc9e0..8ab6f478 100644
--- a/src/projects.h
+++ b/src/projects.h
@@ -224,10 +224,11 @@ struct PJ_REGION_S {
 };
 
 struct PJ_AREA {
-    int     id;         /* Area ID in the EPSG database */
-    LP      ll;         /* Lower left corner of bounding box */
-    LP      ur;         /* Upper right corner of bounding box */
-    char    descr[64];  /* text representation of area */
+    int     bbox_set;
+    double  west_lon_degree;
+    double  south_lat_degree;
+    double  east_lon_degree;
+    double  north_lat_degree;
 };
 
 struct projCtx_t;
@@ -596,6 +597,7 @@ struct projCtx_t {
     void    *app_data;
     struct projFileAPI_t *fileapi;
     struct projCppContext* cpp_context; /* internal context for C++ code */
+    int     use_proj4_init_rules;
 };
 
 /* classic public API */
diff --git a/test/unit/gie_self_tests.cpp b/test/unit/gie_self_tests.cpp
index 67a44fe5..8b90a22d 100644
--- a/test/unit/gie_self_tests.cpp
+++ b/test/unit/gie_self_tests.cpp
@@ -49,24 +49,11 @@ TEST( gie, cart_selftest ) {
     PJ_COORD a, b, obs[2];
     PJ_COORD coord[2];
 
-    PJ_INFO info;
-    PJ_PROJ_INFO pj_info;
-    PJ_GRID_INFO grid_info;
-    PJ_INIT_INFO init_info;
-
-    PJ_FACTORS factors;
-
-    const PJ_OPERATIONS *oper_list;
-    const PJ_ELLPS *ellps_list;
-    const PJ_UNITS *unit_list;
-    const PJ_PRIME_MERIDIANS *pm_list;
-
     int err;
     size_t n, sz;
     double dist, h, t;
     const char * const args[3] = {"proj=utm", "zone=32", "ellps=GRS80"};
     char arg[50] = {"+proj=utm; +zone=32; +ellps=GRS80"};
-    char buf[40];
 
     /* An utm projection on the GRS80 ellipsoid */
     P = proj_create (PJ_DEFAULT_CTX, arg);
@@ -246,10 +233,32 @@ TEST( gie, cart_selftest ) {
 
     /* Clean up  after proj_trans_* tests */
     proj_destroy (P);
+}
+
+// ---------------------------------------------------------------------------
+
+class gieTest : public ::testing::Test {
+
+    static void DummyLogFunction(void *, int, const char *) {}
+
+  protected:
+    void SetUp() override {
+        m_ctxt = proj_context_create();
+        proj_log_func(m_ctxt, nullptr, DummyLogFunction);
+    }
 
+    void TearDown() override { proj_context_destroy(m_ctxt); }
+
+    PJ_CONTEXT *m_ctxt = nullptr;
+};
+
+// ---------------------------------------------------------------------------
+
+TEST_F( gieTest, proj_create_crs_to_crs ) {
     /* test proj_create_crs_to_crs() */
-    P = proj_create_crs_to_crs(PJ_DEFAULT_CTX, "epsg:25832", "epsg:25833", NULL);
+    auto P = proj_create_crs_to_crs(PJ_DEFAULT_CTX, "epsg:25832", "epsg:25833", NULL);
     ASSERT_TRUE( P != nullptr );
+    PJ_COORD a, b;
 
     a.xy.x =  700000.0;
     a.xy.y = 6000000.0;
@@ -257,14 +266,39 @@ TEST( gie, cart_selftest ) {
     b.xy.y = 5999669.3036037628;
 
     a = proj_trans(P, PJ_FWD, a);
-    ASSERT_LE(dist, 1e-7);
+    EXPECT_NEAR( a.xy.x, b.xy.x, 1e-9);
+    EXPECT_NEAR( a.xy.y, b.xy.y, 1e-9);
     proj_destroy(P);
 
-    /* let's make sure that only entries in init-files results in a usable PJ */
+    /* we can also allow PROJ strings as a usable PJ */
     P = proj_create_crs_to_crs(PJ_DEFAULT_CTX, "proj=utm +zone=32 +datum=WGS84", "proj=utm +zone=33 +datum=WGS84", NULL);
-    ASSERT_TRUE( P == nullptr );
+    ASSERT_TRUE( P != nullptr );
     proj_destroy(P);
 
+    EXPECT_TRUE( proj_create_crs_to_crs(m_ctxt, "invalid", "EPSG:25833", NULL) == nullptr );
+    EXPECT_TRUE( proj_create_crs_to_crs(m_ctxt, "EPSG:25832", "invalid", NULL) == nullptr );
+}
+
+// ---------------------------------------------------------------------------
+
+TEST( gie, info_functions ) {
+    PJ_INFO info;
+    PJ_PROJ_INFO pj_info;
+    PJ_GRID_INFO grid_info;
+    PJ_INIT_INFO init_info;
+
+    PJ_FACTORS factors;
+
+    const PJ_OPERATIONS *oper_list;
+    const PJ_ELLPS *ellps_list;
+    const PJ_UNITS *unit_list;
+    const PJ_PRIME_MERIDIANS *pm_list;
+
+    char buf[40];
+    PJ *P;
+    char arg[50] = {"+proj=utm; +zone=32; +ellps=GRS80"};
+    PJ_COORD a;
+
     /* ********************************************************************** */
     /*                          Test info functions                           */
     /* ********************************************************************** */
@@ -345,7 +379,7 @@ TEST( gie, cart_selftest ) {
     proj_destroy(P);
 
     /* Check that proj_list_* functions work by looping through them */
-    n = 0;
+    size_t n = 0;
     for (oper_list = proj_list_operations(); oper_list->id; ++oper_list) n++;
     ASSERT_NE(n, 0U);
 
@@ -361,11 +395,15 @@ TEST( gie, cart_selftest ) {
     for (pm_list = proj_list_prime_meridians(); pm_list->id; ++pm_list) n++;
     ASSERT_NE(n, 0U);
 
+}
+
+// ---------------------------------------------------------------------------
 
+TEST( gie, io_predicates ) {
     /* check io-predicates */
 
     /* angular in on fwd, linear out */
-    P = proj_create (PJ_DEFAULT_CTX, "+proj=cart +ellps=GRS80");
+    auto P = proj_create (PJ_DEFAULT_CTX, "+proj=cart +ellps=GRS80");
     ASSERT_TRUE( P != nullptr );
     ASSERT_TRUE(proj_angular_input (P, PJ_FWD));
     ASSERT_FALSE(proj_angular_input (P, PJ_INV));
@@ -425,7 +463,7 @@ TEST( gie, cart_selftest ) {
     /* Test that pj_fwd* and pj_inv* returns NaNs when receiving NaN input */
     P = proj_create(PJ_DEFAULT_CTX, "+proj=merc");
     ASSERT_TRUE( P != nullptr );
-    a = proj_coord(NAN, NAN, NAN, NAN);
+    auto a = proj_coord(NAN, NAN, NAN, NAN);
     a = proj_trans(P, PJ_FWD, a);
     ASSERT_TRUE ( ( std::isnan(a.v[0]) && std::isnan(a.v[1]) && std::isnan(a.v[2]) && std::isnan(a.v[3]) ) );
 
diff --git a/test/unit/test_io.cpp b/test/unit/test_io.cpp
index 6ccd5578..99f58739 100644
--- a/test/unit/test_io.cpp
+++ b/test/unit/test_io.cpp
@@ -7448,12 +7448,43 @@ TEST(io, projparse_projected_title) {
 
 TEST(io, projparse_init) {
 
+    // Not allowed in non-compatibillity mode
+    EXPECT_THROW(PROJStringParser().createFromPROJString("init=epsg:4326"),
+                 ParsingException);
+
+    {
+        // EPSG:4326 is normally latitude-longitude order with degree,
+        // but in compatibillity mode it will be long-lat radian
+        auto dbContext = DatabaseContext::create();
+        auto obj = createFromUserInput("init=epsg:4326", dbContext, true);
+        auto crs = nn_dynamic_pointer_cast<GeographicCRS>(obj);
+        ASSERT_TRUE(crs != nullptr);
+        EXPECT_TRUE(crs->coordinateSystem()->isEquivalentTo(
+            EllipsoidalCS::createLongitudeLatitude(UnitOfMeasure::RADIAN)
+                .get()));
+    }
+
+    {
+        // EPSG:3040 is normally northing-easting order, but in compatibillity
+        // mode it will be easting-northing
+        auto dbContext = DatabaseContext::create();
+        auto obj = createFromUserInput("init=epsg:3040", dbContext, true);
+        auto crs = nn_dynamic_pointer_cast<ProjectedCRS>(obj);
+        ASSERT_TRUE(crs != nullptr);
+        EXPECT_TRUE(crs->coordinateSystem()->isEquivalentTo(
+            CartesianCS::createEastingNorthing(UnitOfMeasure::METRE).get()));
+    }
+
     {
-        auto obj = PROJStringParser().createFromPROJString("init=epsg:4326");
+        auto obj =
+            PROJStringParser().createFromPROJString("init=ITRF2000:ITRF2005");
         auto co = nn_dynamic_pointer_cast<CoordinateOperation>(obj);
         ASSERT_TRUE(co != nullptr);
         EXPECT_EQ(co->exportToPROJString(PROJStringFormatter::create().get()),
-                  "+init=epsg:4326");
+                  "+proj=helmert +x=-0.0001 +y=0.0008 +z=0.0058 +rx=0 +ry=0 "
+                  "+rz=0 +s=-0.0004 +dx=0.0002 +dy=-0.0001 +dz=0.0018 +drx=0 "
+                  "+dry=0 +drz=0 +ds=-8e-06 +t_epoch=2000 "
+                  "+convention=position_vector");
     }
 
     {