From 805b187edd4c9246629e9aeb062118f8c2de2dfe Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 29 Oct 2019 11:48:17 +0100 Subject: Little tweaks in implicit 2D/3D geog conversions in compoundCRS to geogCRS transformations --- test/unit/test_operation.cpp | 36 +++++++++++++++++++++++++----------- 1 file changed, 25 insertions(+), 11 deletions(-) (limited to 'test/unit/test_operation.cpp') diff --git a/test/unit/test_operation.cpp b/test/unit/test_operation.cpp index 739d8ec3..79541d88 100644 --- a/test/unit/test_operation.cpp +++ b/test/unit/test_operation.cpp @@ -7129,6 +7129,11 @@ TEST(operation, compoundCRS_to_geogCRS_2D_promote_to_3D_context) { auto authFactory = AuthorityFactory::create(DatabaseContext::create(), "EPSG"); auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0); + ctxt->setSpatialCriterion( + CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); + ctxt->setGridAvailabilityUse( + CoordinateOperationContext::GridAvailabilityUse:: + IGNORE_GRID_AVAILABILITY); // NAD83 + NAVD88 height auto srcObj = createFromUserInput("EPSG:4269+5703", authFactory->databaseContext(), false); @@ -7136,17 +7141,26 @@ TEST(operation, compoundCRS_to_geogCRS_2D_promote_to_3D_context) { ASSERT_TRUE(src != nullptr); auto nnSrc = NN_NO_CHECK(src); auto dst = authFactory->createCoordinateReferenceSystem("4269"); // NAD83 - dst = dst->promoteTo3D(std::string(), authFactory->databaseContext()); - { - auto list = CoordinateOperationFactory::create()->createOperations( - nnSrc, dst, ctxt); - ASSERT_GE(list.size(), 8U); - } - { - auto list = CoordinateOperationFactory::create()->createOperations( - dst, nnSrc, ctxt); - ASSERT_GE(list.size(), 8U); - } + + auto listCompoundToGeog2D = + CoordinateOperationFactory::create()->createOperations(nnSrc, dst, + ctxt); + // The checked value is not that important, but in case this changes, + // likely due to a EPSG upgrade, worth checking + ASSERT_EQ(listCompoundToGeog2D.size(), 469U); + + auto listGeog2DToCompound = + CoordinateOperationFactory::create()->createOperations(dst, nnSrc, + ctxt); + ASSERT_EQ(listGeog2DToCompound.size(), listCompoundToGeog2D.size()); + + auto listCompoundToGeog3D = + CoordinateOperationFactory::create()->createOperations( + nnSrc, + dst->promoteTo3D(std::string(), authFactory->databaseContext()), + ctxt); + // It includes a trailing ballpart transformation + ASSERT_EQ(listCompoundToGeog3D.size(), listCompoundToGeog2D.size() + 1); } // --------------------------------------------------------------------------- -- cgit v1.2.3 From ffc865a41aa540673eaedb2552565cf9f8d18679 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Tue, 29 Oct 2019 22:20:24 +0100 Subject: Vertical transformations: improve situations similar to transforming from 'NAVD88 (ftUS)' to X, where we now consider the available transformations from 'NAVD88' to X that might exist in the database --- test/unit/test_operation.cpp | 121 +++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 121 insertions(+) (limited to 'test/unit/test_operation.cpp') diff --git a/test/unit/test_operation.cpp b/test/unit/test_operation.cpp index 79541d88..7fbbf70e 100644 --- a/test/unit/test_operation.cpp +++ b/test/unit/test_operation.cpp @@ -7165,6 +7165,127 @@ TEST(operation, compoundCRS_to_geogCRS_2D_promote_to_3D_context) { // --------------------------------------------------------------------------- +TEST(operation, compoundCRS_to_geogCRS_with_vertical_unit_change) { + auto authFactory = + AuthorityFactory::create(DatabaseContext::create(), "EPSG"); + auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0); + ctxt->setSpatialCriterion( + CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); + ctxt->setGridAvailabilityUse( + CoordinateOperationContext::GridAvailabilityUse:: + IGNORE_GRID_AVAILABILITY); + // NAD83(2011) + NAVD88 height (ftUS) + auto srcObj = createFromUserInput("EPSG:6318+6360", + authFactory->databaseContext(), false); + auto src = nn_dynamic_pointer_cast(srcObj); + ASSERT_TRUE(src != nullptr); + auto nnSrc = NN_NO_CHECK(src); + auto dst = + authFactory->createCoordinateReferenceSystem("6319"); // NAD83(2011) 3D + + auto listCompoundToGeog = + CoordinateOperationFactory::create()->createOperations(nnSrc, dst, + ctxt); + ASSERT_TRUE(!listCompoundToGeog.empty()); + + // NAD83(2011) + NAVD88 height + auto srcObjCompoundVMetre = createFromUserInput( + "EPSG:6318+5703", authFactory->databaseContext(), false); + auto srcCompoundVMetre = nn_dynamic_pointer_cast(srcObjCompoundVMetre); + ASSERT_TRUE(srcCompoundVMetre != nullptr); + auto listCompoundMetreToGeog = + CoordinateOperationFactory::create()->createOperations( + NN_NO_CHECK(srcCompoundVMetre), dst, ctxt); + + // Check that we get the same and similar results whether we start from + // regular NAVD88 height or its ftUs variant + ASSERT_EQ(listCompoundToGeog.size(), listCompoundMetreToGeog.size()); + + EXPECT_EQ(listCompoundToGeog[0]->nameStr(), + "Transformation from NAVD88 height (ftUS) to NAVD88 height + " + + listCompoundMetreToGeog[0]->nameStr()); + EXPECT_EQ( + listCompoundToGeog[0]->exportToPROJString( + PROJStringFormatter::create(PROJStringFormatter::Convention::PROJ_5, + authFactory->databaseContext()) + .get()), + replaceAll(listCompoundMetreToGeog[0]->exportToPROJString( + PROJStringFormatter::create( + PROJStringFormatter::Convention::PROJ_5, + authFactory->databaseContext()) + .get()), + "+step +proj=unitconvert +xy_in=deg +xy_out=rad", + "+step +proj=unitconvert +xy_in=deg +z_in=us-ft +xy_out=rad " + "+z_out=m")); + + // Check reverse path + auto listGeogToCompound = + CoordinateOperationFactory::create()->createOperations(dst, nnSrc, + ctxt); + EXPECT_EQ(listGeogToCompound.size(), listCompoundToGeog.size()); +} + +// --------------------------------------------------------------------------- + +TEST( + operation, + compoundCRS_to_geogCRS_with_vertical_unit_change_and_complex_horizontal_change) { + auto authFactory = + AuthorityFactory::create(DatabaseContext::create(), "EPSG"); + auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0); + ctxt->setSpatialCriterion( + CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); + ctxt->setGridAvailabilityUse( + CoordinateOperationContext::GridAvailabilityUse:: + IGNORE_GRID_AVAILABILITY); + // NAD83(2011) + NAVD88 height (ftUS) + auto srcObj = createFromUserInput("EPSG:6318+6360", + authFactory->databaseContext(), false); + auto src = nn_dynamic_pointer_cast(srcObj); + ASSERT_TRUE(src != nullptr); + auto nnSrc = NN_NO_CHECK(src); + auto dst = + authFactory->createCoordinateReferenceSystem("7665"); // WGS84(G1762) 3D + + auto listCompoundToGeog = + CoordinateOperationFactory::create()->createOperations(nnSrc, dst, + ctxt); + + // NAD83(2011) + NAVD88 height + auto srcObjCompoundVMetre = createFromUserInput( + "EPSG:6318+5703", authFactory->databaseContext(), false); + auto srcCompoundVMetre = nn_dynamic_pointer_cast(srcObjCompoundVMetre); + ASSERT_TRUE(srcCompoundVMetre != nullptr); + auto listCompoundMetreToGeog = + CoordinateOperationFactory::create()->createOperations( + NN_NO_CHECK(srcCompoundVMetre), dst, ctxt); + + // Check that we get the same and similar results whether we start from + // regular NAVD88 height or its ftUs variant + ASSERT_EQ(listCompoundToGeog.size(), listCompoundMetreToGeog.size()); + + ASSERT_GE(listCompoundToGeog.size(), 1U); + + EXPECT_EQ(listCompoundToGeog[0]->nameStr(), + "Transformation from NAVD88 height (ftUS) to NAVD88 height + " + + listCompoundMetreToGeog[0]->nameStr()); + EXPECT_EQ( + listCompoundToGeog[0]->exportToPROJString( + PROJStringFormatter::create(PROJStringFormatter::Convention::PROJ_5, + authFactory->databaseContext()) + .get()), + replaceAll(listCompoundMetreToGeog[0]->exportToPROJString( + PROJStringFormatter::create( + PROJStringFormatter::Convention::PROJ_5, + authFactory->databaseContext()) + .get()), + "+step +proj=unitconvert +xy_in=deg +xy_out=rad", + "+step +proj=unitconvert +xy_in=deg +z_in=us-ft +xy_out=rad " + "+z_out=m")); +} + +// --------------------------------------------------------------------------- + TEST(operation, compoundCRS_from_WKT2_to_geogCRS_3D_context) { auto authFactory = AuthorityFactory::create(DatabaseContext::create(), "EPSG"); -- cgit v1.2.3 From 05a7bcfa56a03437b2ba73616a6bc21c9347d2a7 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Wed, 30 Oct 2019 11:28:43 +0100 Subject: Rework importing of Vertical unit change from EPSG db, add support for Height Depth Reversal and use it in createOperations() --- test/unit/test_operation.cpp | 169 ++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 167 insertions(+), 2 deletions(-) (limited to 'test/unit/test_operation.cpp') diff --git a/test/unit/test_operation.cpp b/test/unit/test_operation.cpp index 7fbbf70e..f7a86eb4 100644 --- a/test/unit/test_operation.cpp +++ b/test/unit/test_operation.cpp @@ -6872,6 +6872,43 @@ TEST(operation, vertCRS_to_vertCRS) { EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()), "+proj=affine +s33=0.999998"); } + + auto vertCRSMetreUp = + nn_dynamic_pointer_cast(WKTParser().createFromWKT( + "VERTCRS[\"my height\",VDATUM[\"my datum\"],CS[vertical,1]," + "AXIS[\"gravity-related height (H)\",up," + "LENGTHUNIT[\"metre\",1]]]")); + ASSERT_TRUE(vertCRSMetreUp != nullptr); + + auto vertCRSMetreDown = + nn_dynamic_pointer_cast(WKTParser().createFromWKT( + "VERTCRS[\"my depth\",VDATUM[\"my datum\"],CS[vertical,1]," + "AXIS[\"depth (D)\",down,LENGTHUNIT[\"metre\",1]]]")); + ASSERT_TRUE(vertCRSMetreDown != nullptr); + + auto vertCRSMetreDownFtUS = + nn_dynamic_pointer_cast(WKTParser().createFromWKT( + "VERTCRS[\"my depth (ftUS)\",VDATUM[\"my datum\"],CS[vertical,1]," + "AXIS[\"depth (D)\",down,LENGTHUNIT[\"US survey " + "foot\",0.304800609601219]]]")); + ASSERT_TRUE(vertCRSMetreDownFtUS != nullptr); + + { + auto op = CoordinateOperationFactory::create()->createOperation( + NN_CHECK_ASSERT(vertCRSMetreUp), NN_CHECK_ASSERT(vertCRSMetreDown)); + ASSERT_TRUE(op != nullptr); + EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()), + "+proj=axisswap +order=1,2,-3"); + } + + { + auto op = CoordinateOperationFactory::create()->createOperation( + NN_CHECK_ASSERT(vertCRSMetreUp), + NN_CHECK_ASSERT(vertCRSMetreDownFtUS)); + ASSERT_TRUE(op != nullptr); + EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()), + "+proj=affine +s33=-3.28083333333333"); + } } // --------------------------------------------------------------------------- @@ -7202,7 +7239,7 @@ TEST(operation, compoundCRS_to_geogCRS_with_vertical_unit_change) { ASSERT_EQ(listCompoundToGeog.size(), listCompoundMetreToGeog.size()); EXPECT_EQ(listCompoundToGeog[0]->nameStr(), - "Transformation from NAVD88 height (ftUS) to NAVD88 height + " + + "Inverse of NAVD88 height to NAVD88 height (ftUS) + " + listCompoundMetreToGeog[0]->nameStr()); EXPECT_EQ( listCompoundToGeog[0]->exportToPROJString( @@ -7267,7 +7304,7 @@ TEST( ASSERT_GE(listCompoundToGeog.size(), 1U); EXPECT_EQ(listCompoundToGeog[0]->nameStr(), - "Transformation from NAVD88 height (ftUS) to NAVD88 height + " + + "Inverse of NAVD88 height to NAVD88 height (ftUS) + " + listCompoundMetreToGeog[0]->nameStr()); EXPECT_EQ( listCompoundToGeog[0]->exportToPROJString( @@ -7286,6 +7323,134 @@ TEST( // --------------------------------------------------------------------------- +TEST(operation, compoundCRS_to_geogCRS_with_height_depth_reversal) { + auto authFactory = + AuthorityFactory::create(DatabaseContext::create(), "EPSG"); + auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0); + ctxt->setSpatialCriterion( + CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); + ctxt->setGridAvailabilityUse( + CoordinateOperationContext::GridAvailabilityUse:: + IGNORE_GRID_AVAILABILITY); + // NAD83(2011) + NAVD88 depth + auto srcObj = createFromUserInput("EPSG:6318+6357", + authFactory->databaseContext(), false); + auto src = nn_dynamic_pointer_cast(srcObj); + ASSERT_TRUE(src != nullptr); + auto nnSrc = NN_NO_CHECK(src); + auto dst = + authFactory->createCoordinateReferenceSystem("6319"); // NAD83(2011) 3D + + auto listCompoundToGeog = + CoordinateOperationFactory::create()->createOperations(nnSrc, dst, + ctxt); + ASSERT_TRUE(!listCompoundToGeog.empty()); + + // NAD83(2011) + NAVD88 height + auto srcObjCompoundVMetre = createFromUserInput( + "EPSG:6318+5703", authFactory->databaseContext(), false); + auto srcCompoundVMetre = nn_dynamic_pointer_cast(srcObjCompoundVMetre); + ASSERT_TRUE(srcCompoundVMetre != nullptr); + auto listCompoundMetreToGeog = + CoordinateOperationFactory::create()->createOperations( + NN_NO_CHECK(srcCompoundVMetre), dst, ctxt); + + // Check that we get the same and similar results whether we start from + // regular NAVD88 height or its depth variant + ASSERT_EQ(listCompoundToGeog.size(), listCompoundMetreToGeog.size()); + + EXPECT_EQ(listCompoundToGeog[0]->nameStr(), + "Inverse of NAVD88 height to NAVD88 depth + " + + listCompoundMetreToGeog[0]->nameStr()); + EXPECT_EQ( + listCompoundToGeog[0]->exportToPROJString( + PROJStringFormatter::create(PROJStringFormatter::Convention::PROJ_5, + authFactory->databaseContext()) + .get()), + replaceAll(listCompoundMetreToGeog[0]->exportToPROJString( + PROJStringFormatter::create( + PROJStringFormatter::Convention::PROJ_5, + authFactory->databaseContext()) + .get()), + "+step +proj=unitconvert +xy_in=deg +xy_out=rad", + "+step +proj=unitconvert +xy_in=deg +xy_out=rad " + "+step +proj=axisswap +order=1,2,-3")); + + // Check reverse path + auto listGeogToCompound = + CoordinateOperationFactory::create()->createOperations(dst, nnSrc, + ctxt); + EXPECT_EQ(listGeogToCompound.size(), listCompoundToGeog.size()); +} + +// --------------------------------------------------------------------------- + +TEST( + operation, + compoundCRS_to_geogCRS_with_vertical_unit_change_and_height_depth_reversal) { + auto authFactory = + AuthorityFactory::create(DatabaseContext::create(), "EPSG"); + auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0); + ctxt->setSpatialCriterion( + CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); + ctxt->setGridAvailabilityUse( + CoordinateOperationContext::GridAvailabilityUse:: + IGNORE_GRID_AVAILABILITY); + // NAD83(2011) + NAVD88 depth (ftUS) + auto srcObj = createFromUserInput("EPSG:6318+6358", + authFactory->databaseContext(), false); + auto src = nn_dynamic_pointer_cast(srcObj); + ASSERT_TRUE(src != nullptr); + auto nnSrc = NN_NO_CHECK(src); + auto dst = + authFactory->createCoordinateReferenceSystem("6319"); // NAD83(2011) 3D + + auto listCompoundToGeog = + CoordinateOperationFactory::create()->createOperations(nnSrc, dst, + ctxt); + ASSERT_TRUE(!listCompoundToGeog.empty()); + + // NAD83(2011) + NAVD88 height + auto srcObjCompoundVMetre = createFromUserInput( + "EPSG:6318+5703", authFactory->databaseContext(), false); + auto srcCompoundVMetre = nn_dynamic_pointer_cast(srcObjCompoundVMetre); + ASSERT_TRUE(srcCompoundVMetre != nullptr); + auto listCompoundMetreToGeog = + CoordinateOperationFactory::create()->createOperations( + NN_NO_CHECK(srcCompoundVMetre), dst, ctxt); + + // Check that we get the same and similar results whether we start from + // regular NAVD88 height or its depth (ftUS) variant + ASSERT_EQ(listCompoundToGeog.size(), listCompoundMetreToGeog.size()); + + EXPECT_EQ(listCompoundToGeog[0]->nameStr(), + "Inverse of NAVD88 height (ftUS) to NAVD88 depth (ftUS) + " + "Inverse of NAVD88 height to NAVD88 height (ftUS) + " + + listCompoundMetreToGeog[0]->nameStr()); + EXPECT_EQ( + listCompoundToGeog[0]->exportToPROJString( + PROJStringFormatter::create(PROJStringFormatter::Convention::PROJ_5, + authFactory->databaseContext()) + .get()), + replaceAll(listCompoundMetreToGeog[0]->exportToPROJString( + PROJStringFormatter::create( + PROJStringFormatter::Convention::PROJ_5, + authFactory->databaseContext()) + .get()), + "+step +proj=unitconvert +xy_in=deg +xy_out=rad", + "+step +proj=unitconvert +xy_in=deg +xy_out=rad " + "+step +proj=axisswap +order=1,2,-3 " + "+step +proj=unitconvert +z_in=us-ft +z_out=m")); + + // Check reverse path + auto listGeogToCompound = + CoordinateOperationFactory::create()->createOperations(dst, nnSrc, + ctxt); + EXPECT_EQ(listGeogToCompound.size(), listCompoundToGeog.size()); +} + +// --------------------------------------------------------------------------- + TEST(operation, compoundCRS_from_WKT2_to_geogCRS_3D_context) { auto authFactory = AuthorityFactory::create(DatabaseContext::create(), "EPSG"); -- cgit v1.2.3 From 81e06f42c7552494bcb3f466b0b1317341187679 Mon Sep 17 00:00:00 2001 From: Even Rouault Date: Wed, 30 Oct 2019 18:20:17 +0100 Subject: Add a test to check we can use a VerticalCRS from its name only without the EPSG code --- test/unit/test_operation.cpp | 53 ++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 53 insertions(+) (limited to 'test/unit/test_operation.cpp') diff --git a/test/unit/test_operation.cpp b/test/unit/test_operation.cpp index f7a86eb4..d61c74c0 100644 --- a/test/unit/test_operation.cpp +++ b/test/unit/test_operation.cpp @@ -7202,6 +7202,59 @@ TEST(operation, compoundCRS_to_geogCRS_2D_promote_to_3D_context) { // --------------------------------------------------------------------------- +TEST(operation, compoundCRS_from_wkt_without_id_to_geogCRS) { + auto authFactory = + AuthorityFactory::create(DatabaseContext::create(), "EPSG"); + auto ctxt = CoordinateOperationContext::create(authFactory, nullptr, 0.0); + ctxt->setSpatialCriterion( + CoordinateOperationContext::SpatialCriterion::PARTIAL_INTERSECTION); + ctxt->setGridAvailabilityUse( + CoordinateOperationContext::GridAvailabilityUse:: + IGNORE_GRID_AVAILABILITY); + auto wkt = + "COMPOUNDCRS[\"NAD83(2011) + NAVD88 height\",\n" + " GEOGCRS[\"NAD83(2011)\",\n" + " DATUM[\"NAD83 (National Spatial Reference System 2011)\",\n" + " ELLIPSOID[\"GRS 1980\",6378137,298.257222101,\n" + " LENGTHUNIT[\"metre\",1]]],\n" + " PRIMEM[\"Greenwich\",0,\n" + " ANGLEUNIT[\"degree\",0.0174532925199433]],\n" + " CS[ellipsoidal,2],\n" + " AXIS[\"geodetic latitude (Lat)\",north,\n" + " ORDER[1],\n" + " ANGLEUNIT[\"degree\",0.0174532925199433]],\n" + " AXIS[\"geodetic longitude (Lon)\",east,\n" + " ORDER[2],\n" + " ANGLEUNIT[\"degree\",0.0174532925199433]]],\n" + " VERTCRS[\"NAVD88 height\",\n" + " VDATUM[\"North American Vertical Datum 1988\"],\n" + " CS[vertical,1],\n" + " AXIS[\"gravity-related height (H)\",up,\n" + " LENGTHUNIT[\"metre\",1]]]]"; + auto srcObj = + createFromUserInput(wkt, authFactory->databaseContext(), false); + auto src = nn_dynamic_pointer_cast(srcObj); + ASSERT_TRUE(src != nullptr); + auto dst = + authFactory->createCoordinateReferenceSystem("6319"); // NAD83(2011) + + auto list = CoordinateOperationFactory::create()->createOperations( + NN_NO_CHECK(src), dst, ctxt); + // NAD83(2011) + NAVD88 height + auto srcRefObj = createFromUserInput("EPSG:6318+5703", + authFactory->databaseContext(), false); + auto srcRef = nn_dynamic_pointer_cast(srcRefObj); + ASSERT_TRUE(srcRef != nullptr); + ASSERT_TRUE( + src->isEquivalentTo(srcRef.get(), IComparable::Criterion::EQUIVALENT)); + auto listRef = CoordinateOperationFactory::create()->createOperations( + NN_NO_CHECK(srcRef), dst, ctxt); + + EXPECT_EQ(list.size(), listRef.size()); +} + +// --------------------------------------------------------------------------- + TEST(operation, compoundCRS_to_geogCRS_with_vertical_unit_change) { auto authFactory = AuthorityFactory::create(DatabaseContext::create(), "EPSG"); -- cgit v1.2.3