aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--data/sql/grid_alternatives.sql16
-rw-r--r--scripts/reference_exported_symbols.txt1
-rw-r--r--src/apply_gridshift.cpp6
-rw-r--r--src/apps/proj.cpp8
-rw-r--r--src/gridcatalog.cpp4
-rw-r--r--src/iso19111/c_api.cpp28
-rw-r--r--src/iso19111/crs.cpp49
-rw-r--r--src/iso19111/io.cpp10
-rw-r--r--src/nad_cvt.cpp43
-rw-r--r--src/proj.h3
-rw-r--r--src/proj_internal.h4
-rw-r--r--src/transform.cpp16
-rw-r--r--test/cli/CMakeLists.txt1
-rw-r--r--test/cli/Makefile.am7
-rw-r--r--test/cli/ntv2_out.dist3
-rwxr-xr-xtest/cli/testntv28
-rwxr-xr-xtest/cli/testproj55
-rw-r--r--test/cli/testproj_out.dist1
-rw-r--r--test/unit/pj_transform_test.cpp34
-rw-r--r--test/unit/test_c_api.cpp30
-rw-r--r--test/unit/test_crs.cpp40
-rw-r--r--test/unit/test_io.cpp35
22 files changed, 376 insertions, 26 deletions
diff --git a/data/sql/grid_alternatives.sql b/data/sql/grid_alternatives.sql
index 6a893c5e..070a0a03 100644
--- a/data/sql/grid_alternatives.sql
+++ b/data/sql/grid_alternatives.sql
@@ -1106,6 +1106,22 @@ INSERT INTO grid_alternatives(original_grid_name,
'proj-datumgrid-oceania',
NULL, NULL, NULL, NULL);
+INSERT INTO grid_alternatives(original_grid_name,
+ proj_grid_name,
+ proj_grid_format,
+ proj_method,
+ inverse_direction,
+ package_name,
+ url, direct_download, open_license, directory)
+ SELECT grid_name,
+ 'AUSGeoid98.gtx',
+ 'GTX',
+ 'vgridshift',
+ 1,
+ 'proj-datumgrid-oceania',
+ NULL, NULL, NULL, NULL FROM grid_transformation WHERE
+ grid_name LIKE '%DAT.htm' AND name LIKE 'GDA94 to AHD height%';
+
-- Netherlands / RDNAP (non-free grids)
INSERT INTO grid_alternatives(original_grid_name,
diff --git a/scripts/reference_exported_symbols.txt b/scripts/reference_exported_symbols.txt
index e20f29b3..3f4cd51b 100644
--- a/scripts/reference_exported_symbols.txt
+++ b/scripts/reference_exported_symbols.txt
@@ -801,6 +801,7 @@ proj_context_use_proj4_init_rules
proj_convert_conversion_to_other_method
proj_coord
proj_coord_error()
+proj_coordoperation_create_inverse
proj_coordoperation_get_accuracy
proj_coordoperation_get_grid_used
proj_coordoperation_get_grid_used_count
diff --git a/src/apply_gridshift.cpp b/src/apply_gridshift.cpp
index e6c7b2b1..39e7c3b5 100644
--- a/src/apply_gridshift.cpp
+++ b/src/apply_gridshift.cpp
@@ -115,7 +115,7 @@ int pj_apply_gridshift_2( PJ *defn, int inverse,
/* Determine which grid is the correct given an input coordinate. */
/************************************************************************/
-static struct CTABLE* find_ctable(projCtx ctx, PJ_LP input, int grid_count, PJ_GRIDINFO **tables) {
+struct CTABLE* find_ctable(projCtx ctx, PJ_LP input, int grid_count, PJ_GRIDINFO **tables) {
int itable;
/* keep trying till we find a table that works */
@@ -210,7 +210,7 @@ int pj_apply_gridshift_3( projCtx ctx, PJ_GRIDINFO **gridlist, int gridlist_coun
ct = find_ctable(ctx, input, gridlist_count, gridlist);
if( ct != nullptr )
{
- output = nad_cvt( input, inverse, ct );
+ output = nad_cvt( ctx, input, inverse, ct, gridlist_count, gridlist);
if ( output.lam != HUGE_VAL && debug_count++ < 20 )
pj_log( ctx, PJ_LOG_DEBUG_MINOR, "pj_apply_gridshift(): used %s", ct->id );
@@ -356,7 +356,7 @@ PJ_LP proj_hgrid_apply(PJ *P, PJ_LP lp, PJ_DIRECTION direction) {
}
inverse = direction == PJ_FWD ? 0 : 1;
- out = nad_cvt(lp, inverse, ct);
+ out = nad_cvt(P->ctx, lp, inverse, ct, P->gridlist_count, P->gridlist);
if (out.lam == HUGE_VAL || out.phi == HUGE_VAL)
pj_ctx_set_errno(P->ctx, PJD_ERR_GRID_AREA);
diff --git a/src/apps/proj.cpp b/src/apps/proj.cpp
index 7fe08023..09c8a81d 100644
--- a/src/apps/proj.cpp
+++ b/src/apps/proj.cpp
@@ -490,6 +490,14 @@ int main(int argc, char **argv) {
exit(0);
}
+ // Ugly hack. See https://github.com/OSGeo/PROJ/issues/1782
+ if( Proj->right == PJ_IO_UNITS_WHATEVER && Proj->descr &&
+ strncmp(Proj->descr, "General Oblique Transformation",
+ strlen("General Oblique Transformation")) == 0 )
+ {
+ Proj->right = PJ_IO_UNITS_PROJECTED;
+ }
+
if (inverse) {
if (!Proj->inv)
emess(3,"inverse projection not available");
diff --git a/src/gridcatalog.cpp b/src/gridcatalog.cpp
index 15d81dd7..9b94fef8 100644
--- a/src/gridcatalog.cpp
+++ b/src/gridcatalog.cpp
@@ -164,7 +164,7 @@ int pj_gc_apply_gridshift( PJ *defn, int inverse,
return PJD_ERR_FAILED_TO_LOAD_GRID;
}
- output_after = nad_cvt( input, inverse, gi->ct );
+ output_after = nad_cvt( defn->ctx, input, inverse, gi->ct, 0, nullptr );
if( output_after.lam == HUGE_VAL )
{
if( defn->ctx->debug_level >= PJ_LOG_DEBUG_MAJOR )
@@ -213,7 +213,7 @@ int pj_gc_apply_gridshift( PJ *defn, int inverse,
return PJD_ERR_FAILED_TO_LOAD_GRID;
}
- output_before = nad_cvt( input, inverse, gi->ct );
+ output_before = nad_cvt( defn->ctx, input, inverse, gi->ct, 0, nullptr );
if( output_before.lam == HUGE_VAL )
{
if( defn->ctx->debug_level >= PJ_LOG_DEBUG_MAJOR )
diff --git a/src/iso19111/c_api.cpp b/src/iso19111/c_api.cpp
index 9db9e5b3..5e2ac522 100644
--- a/src/iso19111/c_api.cpp
+++ b/src/iso19111/c_api.cpp
@@ -7769,6 +7769,34 @@ PJ *proj_normalize_for_visualization(PJ_CONTEXT *ctx, const PJ *obj) {
// ---------------------------------------------------------------------------
+/** \brief Returns a PJ* coordinate operation object which represents the
+ * inverse operation of the specified coordinate operation.
+ *
+ * @param ctx PROJ context, or NULL for default context
+ * @param obj Object of type CoordinateOperation (must not be NULL)
+ * @return a new PJ* object to free with proj_destroy() in case of success, or
+ * nullptr in case of error
+ * @since 6.3
+ */
+PJ *proj_coordoperation_create_inverse(PJ_CONTEXT *ctx, const PJ *obj) {
+
+ SANITIZE_CTX(ctx);
+ auto co = dynamic_cast<const CoordinateOperation *>(obj->iso_obj.get());
+ if (!co) {
+ proj_log_error(ctx, __FUNCTION__,
+ "Object is not a CoordinateOperation");
+ return nullptr;
+ }
+ try {
+ return pj_obj_create(ctx, co->inverse());
+ } catch (const std::exception &e) {
+ proj_log_debug(ctx, __FUNCTION__, e.what());
+ return nullptr;
+ }
+}
+
+// ---------------------------------------------------------------------------
+
/** \brief Returns the number of steps of a concatenated operation.
*
* The input object must be a concatenated operation.
diff --git a/src/iso19111/crs.cpp b/src/iso19111/crs.cpp
index 49cc050f..46aeee20 100644
--- a/src/iso19111/crs.cpp
+++ b/src/iso19111/crs.cpp
@@ -49,6 +49,7 @@
#include <algorithm>
#include <cassert>
+#include <cmath>
#include <cstring>
#include <iostream>
#include <memory>
@@ -595,12 +596,44 @@ CRSNNPtr CRS::alterId(const std::string &authName,
//! @cond Doxygen_Suppress
-static bool isAxisListNorthEast(
+static bool mustAxisOrderBeSwitchedForVisualizationInternal(
const std::vector<cs::CoordinateSystemAxisNNPtr> &axisList) {
const auto &dir0 = axisList[0]->direction();
const auto &dir1 = axisList[1]->direction();
- return (&dir0 == &cs::AxisDirection::NORTH &&
- &dir1 == &cs::AxisDirection::EAST);
+ if (&dir0 == &cs::AxisDirection::NORTH &&
+ &dir1 == &cs::AxisDirection::EAST) {
+ return true;
+ }
+
+ // Address EPSG:32661 "WGS 84 / UPS North (N,E)"
+ if (&dir0 == &cs::AxisDirection::SOUTH &&
+ &dir1 == &cs::AxisDirection::SOUTH) {
+ const auto &meridian0 = axisList[0]->meridian();
+ const auto &meridian1 = axisList[1]->meridian();
+ return meridian0 != nullptr && meridian1 != nullptr &&
+ std::abs(meridian0->longitude().convertToUnit(
+ common::UnitOfMeasure::DEGREE) -
+ 180.0) < 1e-10 &&
+ std::abs(meridian1->longitude().convertToUnit(
+ common::UnitOfMeasure::DEGREE) -
+ 90.0) < 1e-10;
+ }
+
+ // Address EPSG:32761 "WGS 84 / UPS South (N,E)"
+ if (&dir0 == &cs::AxisDirection::NORTH &&
+ &dir1 == &cs::AxisDirection::NORTH) {
+ const auto &meridian0 = axisList[0]->meridian();
+ const auto &meridian1 = axisList[1]->meridian();
+ return meridian0 != nullptr && meridian1 != nullptr &&
+ std::abs(meridian0->longitude().convertToUnit(
+ common::UnitOfMeasure::DEGREE) -
+ 0.0) < 1e-10 &&
+ std::abs(meridian1->longitude().convertToUnit(
+ common::UnitOfMeasure::DEGREE) -
+ 90.0) < 1e-10;
+ }
+
+ return false;
}
// ---------------------------------------------------------------------------
@@ -616,12 +649,14 @@ bool CRS::mustAxisOrderBeSwitchedForVisualization() const {
const GeographicCRS *geogCRS = dynamic_cast<const GeographicCRS *>(this);
if (geogCRS) {
- return isAxisListNorthEast(geogCRS->coordinateSystem()->axisList());
+ return mustAxisOrderBeSwitchedForVisualizationInternal(
+ geogCRS->coordinateSystem()->axisList());
}
const ProjectedCRS *projCRS = dynamic_cast<const ProjectedCRS *>(this);
if (projCRS) {
- return isAxisListNorthEast(projCRS->coordinateSystem()->axisList());
+ return mustAxisOrderBeSwitchedForVisualizationInternal(
+ projCRS->coordinateSystem()->axisList());
}
return false;
@@ -655,7 +690,7 @@ CRSNNPtr CRS::normalizeForVisualization() const {
const GeographicCRS *geogCRS = dynamic_cast<const GeographicCRS *>(this);
if (geogCRS) {
const auto &axisList = geogCRS->coordinateSystem()->axisList();
- if (isAxisListNorthEast(axisList)) {
+ if (mustAxisOrderBeSwitchedForVisualizationInternal(axisList)) {
auto cs = axisList.size() == 2
? cs::EllipsoidalCS::create(util::PropertyMap(),
axisList[1], axisList[0])
@@ -670,7 +705,7 @@ CRSNNPtr CRS::normalizeForVisualization() const {
const ProjectedCRS *projCRS = dynamic_cast<const ProjectedCRS *>(this);
if (projCRS) {
const auto &axisList = projCRS->coordinateSystem()->axisList();
- if (isAxisListNorthEast(axisList)) {
+ if (mustAxisOrderBeSwitchedForVisualizationInternal(axisList)) {
auto cs =
axisList.size() == 2
? cs::CartesianCS::create(util::PropertyMap(), axisList[1],
diff --git a/src/iso19111/io.cpp b/src/iso19111/io.cpp
index e0e6152a..c704e1c1 100644
--- a/src/iso19111/io.cpp
+++ b/src/iso19111/io.cpp
@@ -2683,7 +2683,9 @@ WKTParser::Private::buildGeodeticCRS(const WKTNodeNNPtr &node) {
auto cs = buildCS(csNode, node, angularUnit);
auto ellipsoidalCS = nn_dynamic_pointer_cast<EllipsoidalCS>(cs);
if (ellipsoidalCS) {
- assert(!ci_equal(nodeName, WKTConstants::GEOCCS));
+ if (ci_equal(nodeName, WKTConstants::GEOCCS)) {
+ throw ParsingException("ellipsoidal CS not expected in GEOCCS");
+ }
try {
auto crs = GeographicCRS::create(props, datum, datumEnsemble,
NN_NO_CHECK(ellipsoidalCS));
@@ -7154,6 +7156,12 @@ void PROJStringFormatter::stopInversion() {
// the current end of steps
for (auto iter = startIter; iter != d->steps_.end(); ++iter) {
iter->inverted = !iter->inverted;
+ for (auto &paramValue : iter->paramValues) {
+ if (paramValue.key == "omit_fwd")
+ paramValue.key = "omit_inv";
+ else if (paramValue.key == "omit_inv")
+ paramValue.key = "omit_fwd";
+ }
}
// And reverse the order of steps in that range as well.
std::reverse(startIter, d->steps_.end());
diff --git a/src/nad_cvt.cpp b/src/nad_cvt.cpp
index 79441d0a..e8b8e9b7 100644
--- a/src/nad_cvt.cpp
+++ b/src/nad_cvt.cpp
@@ -7,10 +7,12 @@
#include "proj_internal.h"
#include <math.h>
+#include <limits>
+
#define MAX_ITERATIONS 10
#define TOL 1e-12
-PJ_LP nad_cvt(PJ_LP in, int inverse, struct CTABLE *ct) {
+PJ_LP nad_cvt(projCtx ctx, PJ_LP in, int inverse, struct CTABLE *ct, int grid_count, PJ_GRIDINFO **tables) {
PJ_LP t, tb,del, dif;
int i = MAX_ITERATIONS;
const double toltol = TOL*TOL;
@@ -40,18 +42,37 @@ PJ_LP nad_cvt(PJ_LP in, int inverse, struct CTABLE *ct) {
do {
del = nad_intr(t, ct);
- /* This case used to return failure, but I have
- changed it to return the first order approximation
- of the inverse shift. This avoids cases where the
- grid shift *into* this grid came from another grid.
- While we aren't returning optimally correct results
- I feel a close result in this case is better than
- no result. NFW
- To demonstrate use -112.5839956 49.4914451 against
- the NTv2 grid shift file from Canada. */
- if (del.lam == HUGE_VAL)
+ /* In case we are (or have switched) on the null grid, exit now */
+ if( del.lam == 0 && del.phi == 0 )
break;
+ /* We can possibly go outside of the initial guessed grid, so try */
+ /* to fetch a new grid into which iterate... */
+ if (del.lam == HUGE_VAL)
+ {
+ if( grid_count == 0 )
+ break;
+ PJ_LP lp;
+ lp.lam = t.lam + ct->ll.lam;
+ lp.phi = t.phi + ct->ll.phi;
+ auto newCt = find_ctable(ctx, lp, grid_count, tables);
+ if( newCt == nullptr || newCt == ct )
+ break;
+ pj_log(ctx, PJ_LOG_DEBUG_MINOR, "Switching from grid %s to grid %s",
+ ct->id,
+ newCt->id);
+ ct = newCt;
+ t.lam = lp.lam - ct->ll.lam;
+ t.phi = lp.phi - ct->ll.phi;
+ tb = in;
+ tb.lam -= ct->ll.lam;
+ tb.phi -= ct->ll.phi;
+ tb.lam = adjlon (tb.lam - M_PI) + M_PI;
+ dif.lam = std::numeric_limits<double>::max();
+ dif.phi = std::numeric_limits<double>::max();
+ continue;
+ }
+
dif.lam = t.lam - del.lam - tb.lam;
dif.phi = t.phi + del.phi - tb.phi;
t.lam -= dif.lam;
diff --git a/src/proj.h b/src/proj.h
index fc309542..b8b90277 100644
--- a/src/proj.h
+++ b/src/proj.h
@@ -1091,6 +1091,9 @@ int PROJ_DLL proj_coordoperation_get_towgs84_values(PJ_CONTEXT *ctx,
int value_count,
int emit_error_if_incompatible);
+PJ PROJ_DLL *proj_coordoperation_create_inverse(PJ_CONTEXT *ctx, const PJ *obj);
+
+
int PROJ_DLL proj_concatoperation_get_step_count(PJ_CONTEXT *ctx,
const PJ *concatoperation);
diff --git a/src/proj_internal.h b/src/proj_internal.h
index 3ca927a3..3e219682 100644
--- a/src/proj_internal.h
+++ b/src/proj_internal.h
@@ -863,8 +863,10 @@ int pj_deriv(PJ_LP, double, const PJ *, struct DERIVS *);
int pj_factors(PJ_LP, const PJ *, double, struct FACTORS *);
/* nadcon related protos */
+struct CTABLE* find_ctable(projCtx_t *ctx, PJ_LP input, int grid_count, PJ_GRIDINFO **tables);
+
PJ_LP nad_intr(PJ_LP, struct CTABLE *);
-PJ_LP nad_cvt(PJ_LP, int, struct CTABLE *);
+PJ_LP nad_cvt(projCtx_t *ctx, PJ_LP in, int inverse, struct CTABLE *ct, int grid_count, PJ_GRIDINFO **tables);
struct CTABLE *nad_init(projCtx_t *ctx, char *);
struct CTABLE *nad_ctable_init( projCtx_t *ctx, struct projFileAPI_t* fid );
int nad_ctable_load( projCtx_t *ctx, struct CTABLE *, struct projFileAPI_t* fid );
diff --git a/src/transform.cpp b/src/transform.cpp
index d111d835..781c0061 100644
--- a/src/transform.cpp
+++ b/src/transform.cpp
@@ -231,6 +231,14 @@ static int geographic_to_projected (PJ *P, long n, int dist, double *x, double *
return 0;
}
+ // Ugly hack. See https://github.com/OSGeo/PROJ/issues/1782
+ if( P->right == PJ_IO_UNITS_WHATEVER && P->descr &&
+ strncmp(P->descr, "General Oblique Transformation",
+ strlen("General Oblique Transformation")) == 0 )
+ {
+ P->right = PJ_IO_UNITS_PROJECTED;
+ }
+
for( i = 0; i <n; i++ )
{
PJ_XY projected_loc;
@@ -342,6 +350,14 @@ static int projected_to_geographic (PJ *P, long n, int dist, double *x, double *
return 0;
}
+ // Ugly hack. See https://github.com/OSGeo/PROJ/issues/1782
+ if( P->right == PJ_IO_UNITS_WHATEVER && P->descr &&
+ strncmp(P->descr, "General Oblique Transformation",
+ strlen("General Oblique Transformation")) == 0 )
+ {
+ P->right = PJ_IO_UNITS_PROJECTED;
+ }
+
/* Fallback to the original PROJ.4 API 2d inversion - inv */
for( i = 0; i < n; i++ ) {
PJ_XY projected_loc;
diff --git a/test/cli/CMakeLists.txt b/test/cli/CMakeLists.txt
index 4e1ab75a..feef5bf0 100644
--- a/test/cli/CMakeLists.txt
+++ b/test/cli/CMakeLists.txt
@@ -7,6 +7,7 @@ set(PROJINFO_BIN "projinfo")
set(CCT_BIN "cct")
proj_add_test_script_sh("test27" PROJ_BIN)
proj_add_test_script_sh("test83" PROJ_BIN)
+proj_add_test_script_sh("testproj" PROJ_BIN)
proj_add_test_script_sh("testvarious" CS2CS_BIN)
proj_add_test_script_sh("testdatumfile" CS2CS_BIN "connu")
proj_add_test_script_sh("testIGNF" CS2CS_BIN "ntf_r93.gsb")
diff --git a/test/cli/Makefile.am b/test/cli/Makefile.am
index 758352c6..fb7169e4 100644
--- a/test/cli/Makefile.am
+++ b/test/cli/Makefile.am
@@ -10,6 +10,7 @@ CCTEXE = $(EXEPATH)/cct
# PROJ.4 test scripts
TEST27 = $(THIS_DIR)/test27
TEST83 = $(THIS_DIR)/test83
+TESTPROJ = $(THIS_DIR)/testproj
TESTNTV2 = $(THIS_DIR)/testntv2
TESTVARIOUS = $(THIS_DIR)/testvarious
TESTFLAKY = $(THIS_DIR)/testflaky
@@ -24,6 +25,7 @@ EXTRA_DIST = pj_out27.dist pj_out83.dist td_out.dist \
testIGNF proj_outIGNF.dist \
testprojinfo testprojinfo_out.dist \
testcct testcct_out.dist \
+ testproj testproj_out.dist \
CMakeLists.txt
testprojinfo-check:
@@ -35,6 +37,9 @@ test27-check:
test83-check:
$(TEST83) $(PROJEXE)
+testproj-check:
+ $(TESTPROJ) $(PROJEXE)
+
testvarious-check:
PROJ_LIB=$(PROJ_LIB) $(TESTVARIOUS) $(CS2CSEXE)
@@ -56,4 +61,4 @@ testntv2-check:
testcct-check:
PROJ_LIB=$(PROJ_LIB) $(TESTCCT) $(CCTEXE)
-check-local: testprojinfo-check test27-check test83-check testvarious-check testdatumfile-check testign-check testntv2-check testcct-check
+check-local: testprojinfo-check test27-check test83-check testproj-check testvarious-check testdatumfile-check testign-check testntv2-check testcct-check
diff --git a/test/cli/ntv2_out.dist b/test/cli/ntv2_out.dist
index 9a97f9cf..650a69d8 100644
--- a/test/cli/ntv2_out.dist
+++ b/test/cli/ntv2_out.dist
@@ -9,3 +9,6 @@ Try with NTv2 and NTv1 together ... falls back to NTv1
99d00'00.000"W 65d00'00.000"N 0.0 99d0'1.5885"W 65d0'1.3482"N 0.000
111d00'00.000"W 46d00'00.000"N 0.0 111d0'3.1897"W 45d59'59.7489"N 0.000
111d00'00.000"W 47d30'00.000"N 0.0 111d0'2.7989"W 47d29'59.9896"N 0.000
+##############################################################
+Switching between NTv2 subgrids
+-112.5839956 49.4914451 0 -112.58307487 49.49145197 0.00000000
diff --git a/test/cli/testntv2 b/test/cli/testntv2
index 73371dbe..2a31304e 100755
--- a/test/cli/testntv2
+++ b/test/cli/testntv2
@@ -52,6 +52,14 @@ $EXE +proj=latlong +ellps=clrk66 +nadgrids=ntv2_0.gsb,ntv1_can.dat,conus \
111d00'00.000"W 46d00'00.000"N 0.0
111d00'00.000"W 47d30'00.000"N 0.0
EOF
+
+echo "##############################################################" >> ${OUT}
+echo Switching between NTv2 subgrids >> ${OUT}
+# Initial guess is in ALraymnd, going to parent CAwest afterwards
+$EXE +proj=latlong +datum=NAD83 +to +proj=latlong +ellps=clrk66 +nadgrids=ntv2_0.gsb -E -d 8 >>${OUT} <<EOF
+-112.5839956 49.4914451 0
+EOF
+
#
##############################################################################
# Done!
diff --git a/test/cli/testproj b/test/cli/testproj
new file mode 100755
index 00000000..8686224e
--- /dev/null
+++ b/test/cli/testproj
@@ -0,0 +1,55 @@
+:
+# Script to test proj exe
+#
+TEST_CLI_DIR=`dirname $0`
+EXE=$1
+
+usage()
+{
+ echo "Usage: ${0} <path to 'proj' program>"
+ echo
+ exit 1
+}
+
+if test -z "${EXE}"; then
+ EXE=../../src/cs2cs
+fi
+
+if test ! -x ${EXE}; then
+ echo "*** ERROR: Can not find '${EXE}' program!"
+ exit 1
+fi
+
+if test -z "${PROJ_LIB}"; then
+ export PROJ_LIB="`dirname $0`/../../data"
+fi
+
+echo "============================================"
+echo "Running ${0} using ${EXE}:"
+echo "============================================"
+
+OUT=testproj_out
+#
+echo "doing tests into file ${OUT}, please wait"
+#
+$EXE +ellps=WGS84 +proj=ob_tran +o_proj=latlon +o_lon_p=0.0 +o_lat_p=90.0 +lon_0=360.0 +to_meter=0.0174532925199433 +no_defs -E -f '%.3f' >${OUT} <<EOF
+2 49
+EOF
+
+#
+# do 'diff' with distribution results
+echo "diff ${OUT} with testproj_out.dist"
+diff -u -b ${OUT} ${TEST_CLI_DIR}/testproj_out.dist
+if [ $? -ne 0 ] ; then
+ echo ""
+ echo "PROBLEMS HAVE OCCURRED"
+ echo "test file ${OUT} saved"
+ echo
+ exit 100
+else
+ echo "TEST OK"
+ echo "test file ${OUT} removed"
+ echo
+ /bin/rm -f ${OUT}
+ exit 0
+fi
diff --git a/test/cli/testproj_out.dist b/test/cli/testproj_out.dist
new file mode 100644
index 00000000..fad0a014
--- /dev/null
+++ b/test/cli/testproj_out.dist
@@ -0,0 +1 @@
+2 49 2.000 49.000
diff --git a/test/unit/pj_transform_test.cpp b/test/unit/pj_transform_test.cpp
index b3a061b4..1f4473c1 100644
--- a/test/unit/pj_transform_test.cpp
+++ b/test/unit/pj_transform_test.cpp
@@ -614,4 +614,38 @@ TEST(proj_api_h, pj_set_finder) {
pj_set_finder(nullptr);
}
+// ---------------------------------------------------------------------------
+
+TEST(pj_transform_test, ob_tran_to_meter_as_dest) {
+ auto src = pj_init_plus(
+ "+ellps=WGS84 +a=57.29577951308232 +proj=eqc +lon_0=0.0 +no_defs");
+ auto dst = pj_init_plus("+ellps=WGS84 +proj=ob_tran +o_proj=latlon "
+ "+o_lon_p=0.0 +o_lat_p=90.0 +lon_0=360.0 "
+ "+to_meter=0.0174532925199433 +no_defs");
+ double x = 2 * DEG_TO_RAD;
+ double y = 49 * DEG_TO_RAD;
+ EXPECT_EQ(pj_transform(src, dst, 1, 0, &x, &y, nullptr), 0);
+ EXPECT_NEAR(x, 2 * DEG_TO_RAD, 1e-12) << x / DEG_TO_RAD;
+ EXPECT_NEAR(y, 49 * DEG_TO_RAD, 1e-12) << y / DEG_TO_RAD;
+ pj_free(src);
+ pj_free(dst);
+}
+
+// ---------------------------------------------------------------------------
+
+TEST(pj_transform_test, ob_tran_to_meter_as_srouce) {
+ auto src = pj_init_plus("+ellps=WGS84 +proj=ob_tran +o_proj=latlon "
+ "+o_lon_p=0.0 +o_lat_p=90.0 +lon_0=360.0 "
+ "+to_meter=0.0174532925199433 +no_defs");
+ auto dst = pj_init_plus(
+ "+ellps=WGS84 +a=57.29577951308232 +proj=eqc +lon_0=0.0 +no_defs");
+ double x = 2 * DEG_TO_RAD;
+ double y = 49 * DEG_TO_RAD;
+ EXPECT_EQ(pj_transform(src, dst, 1, 0, &x, &y, nullptr), 0);
+ EXPECT_NEAR(x, 2 * DEG_TO_RAD, 1e-12) << x / DEG_TO_RAD;
+ EXPECT_NEAR(y, 49 * DEG_TO_RAD, 1e-12) << y / DEG_TO_RAD;
+ pj_free(src);
+ pj_free(dst);
+}
+
} // namespace
diff --git a/test/unit/test_c_api.cpp b/test/unit/test_c_api.cpp
index f87f6589..22e2ac11 100644
--- a/test/unit/test_c_api.cpp
+++ b/test/unit/test_c_api.cpp
@@ -3486,6 +3486,36 @@ TEST_F(CApi, proj_normalize_for_visualization_on_crs) {
// ---------------------------------------------------------------------------
+TEST_F(CApi, proj_coordoperation_create_inverse) {
+
+ auto P = proj_create(
+ m_ctxt, "+proj=pipeline +step +proj=axisswap +order=2,1 +step "
+ "+proj=unitconvert +xy_in=deg +xy_out=rad +step +proj=push "
+ "+v_3 +step +proj=cart +ellps=evrst30 +step +proj=helmert "
+ "+x=293 +y=836 +z=318 +rx=0.5 +ry=1.6 +rz=-2.8 +s=2.1 "
+ "+convention=position_vector +step +inv +proj=cart "
+ "+ellps=WGS84 +step +proj=pop +v_3 +step +proj=unitconvert "
+ "+xy_in=rad +xy_out=deg +step +proj=axisswap +order=2,1");
+ ObjectKeeper keeper_P(P);
+ ASSERT_NE(P, nullptr);
+ auto Pinversed = proj_coordoperation_create_inverse(m_ctxt, P);
+ ObjectKeeper keeper_Pinversed(Pinversed);
+ ASSERT_NE(Pinversed, nullptr);
+
+ auto projstr = proj_as_proj_string(m_ctxt, Pinversed, PJ_PROJ_5, nullptr);
+ ASSERT_NE(projstr, nullptr);
+ EXPECT_EQ(std::string(projstr),
+ "+proj=pipeline +step +proj=axisswap +order=2,1 +step "
+ "+proj=unitconvert +xy_in=deg +xy_out=rad +step +proj=push +v_3 "
+ "+step +proj=cart +ellps=WGS84 +step +inv +proj=helmert +x=293 "
+ "+y=836 +z=318 +rx=0.5 +ry=1.6 +rz=-2.8 +s=2.1 "
+ "+convention=position_vector +step +inv +proj=cart "
+ "+ellps=evrst30 +step +proj=pop +v_3 +step +proj=unitconvert "
+ "+xy_in=rad +xy_out=deg +step +proj=axisswap +order=2,1");
+}
+
+// ---------------------------------------------------------------------------
+
TEST_F(CApi, proj_get_remarks) {
auto co = proj_create_from_database(m_ctxt, "EPSG", "8048",
PJ_CATEGORY_COORDINATE_OPERATION, false,
diff --git a/test/unit/test_crs.cpp b/test/unit/test_crs.cpp
index e896853f..33d67e0a 100644
--- a/test/unit/test_crs.cpp
+++ b/test/unit/test_crs.cpp
@@ -943,6 +943,16 @@ TEST(crs, EPSG_32661_projected_north_pole_north_east) {
ASSERT_TRUE(op != nullptr);
EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()),
proj_string);
+
+ auto opNormalized = op->normalizeForVisualization();
+ auto proj_string_normalized =
+ "+proj=pipeline +step +proj=unitconvert +xy_in=deg +xy_out=rad "
+ "+step "
+ "+proj=stere +lat_0=90 +lon_0=0 +k=0.994 +x_0=2000000 +y_0=2000000 "
+ "+ellps=WGS84";
+ EXPECT_EQ(
+ opNormalized->exportToPROJString(PROJStringFormatter::create().get()),
+ proj_string_normalized);
}
// ---------------------------------------------------------------------------
@@ -964,6 +974,16 @@ TEST(crs, EPSG_5041_projected_north_pole_east_north) {
ASSERT_TRUE(op != nullptr);
EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()),
proj_string);
+
+ auto opNormalized = op->normalizeForVisualization();
+ auto proj_string_normalized =
+ "+proj=pipeline +step +proj=unitconvert +xy_in=deg +xy_out=rad "
+ "+step "
+ "+proj=stere +lat_0=90 +lon_0=0 +k=0.994 +x_0=2000000 +y_0=2000000 "
+ "+ellps=WGS84";
+ EXPECT_EQ(
+ opNormalized->exportToPROJString(PROJStringFormatter::create().get()),
+ proj_string_normalized);
}
// ---------------------------------------------------------------------------
@@ -985,6 +1005,16 @@ TEST(crs, EPSG_32761_projected_south_pole_north_east) {
ASSERT_TRUE(op != nullptr);
EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()),
proj_string);
+
+ auto opNormalized = op->normalizeForVisualization();
+ auto proj_string_normalized =
+ "+proj=pipeline +step +proj=unitconvert +xy_in=deg +xy_out=rad "
+ "+step "
+ "+proj=stere +lat_0=-90 +lon_0=0 +k=0.994 +x_0=2000000 +y_0=2000000 "
+ "+ellps=WGS84";
+ EXPECT_EQ(
+ opNormalized->exportToPROJString(PROJStringFormatter::create().get()),
+ proj_string_normalized);
}
// ---------------------------------------------------------------------------
@@ -1006,6 +1036,16 @@ TEST(crs, EPSG_5042_projected_south_pole_east_north) {
ASSERT_TRUE(op != nullptr);
EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()),
proj_string);
+
+ auto opNormalized = op->normalizeForVisualization();
+ auto proj_string_normalized =
+ "+proj=pipeline +step +proj=unitconvert +xy_in=deg +xy_out=rad "
+ "+step "
+ "+proj=stere +lat_0=-90 +lon_0=0 +k=0.994 +x_0=2000000 +y_0=2000000 "
+ "+ellps=WGS84";
+ EXPECT_EQ(
+ opNormalized->exportToPROJString(PROJStringFormatter::create().get()),
+ proj_string_normalized);
}
// ---------------------------------------------------------------------------
diff --git a/test/unit/test_io.cpp b/test/unit/test_io.cpp
index 38dfc2b4..15ab8706 100644
--- a/test/unit/test_io.cpp
+++ b/test/unit/test_io.cpp
@@ -5577,6 +5577,17 @@ TEST(wkt_parse, invalid_GEOCCS) {
"NORTH],AXIS[\"longitude\",EAST]]"),
ParsingException);
+ // ellipsoidal CS is invalid in a GEOCCS
+ EXPECT_THROW(WKTParser().createFromWKT(
+ "GEOCCS[\"WGS 84\",DATUM[\"World Geodetic System 1984\","
+ "ELLIPSOID[\"WGS 84\",6378274,298.257223564,"
+ "LENGTHUNIT[\"metre\",1]]],"
+ "CS[ellipsoidal,2],AXIS[\"geodetic latitude (Lat)\",north,"
+ "ANGLEUNIT[\"degree\",0.0174532925199433]],"
+ "AXIS[\"geodetic longitude (Lon)\",east,"
+ "ANGLEUNIT[\"degree\",0.0174532925199433]]]"),
+ ParsingException);
+
// 3 axis required
EXPECT_THROW(WKTParser().createFromWKT(
"GEOCCS[\"x\",DATUM[\"x\",SPHEROID[\"x\",1,0.5]],PRIMEM["
@@ -6856,6 +6867,30 @@ TEST(io, projstringformatter_optim_hgridshift_vgridshift_hgridshift_inv) {
"+step +proj=pop +v_1 +v_2");
}
+ // Test omit_fwd->omit_inv when inversing the pipeline
+ {
+ auto fmt = PROJStringFormatter::create();
+ fmt->startInversion();
+ fmt->ingestPROJString("+proj=hgridshift +grids=foo +omit_fwd");
+ fmt->stopInversion();
+
+ EXPECT_EQ(fmt->toString(),
+ "+proj=pipeline "
+ "+step +inv +proj=hgridshift +grids=foo +omit_inv");
+ }
+
+ // Test omit_inv->omit_fwd when inversing the pipeline
+ {
+ auto fmt = PROJStringFormatter::create();
+ fmt->startInversion();
+ fmt->ingestPROJString("+proj=hgridshift +grids=foo +omit_inv");
+ fmt->stopInversion();
+
+ EXPECT_EQ(fmt->toString(),
+ "+proj=pipeline "
+ "+step +inv +proj=hgridshift +grids=foo +omit_fwd");
+ }
+
// Variant with first hgridshift inverted, and second forward
{
auto fmt = PROJStringFormatter::create();