[PROJ] Affine transformation using DERIVEDPROJCRS

[Bart.Duijndam at ziggo.nl]

** RE-SENT AS PLAIN TEXT IN VIEW OF EMAIL SIZE **

Hi Even,

I worked things through, and have included the following three routines in test_crs.cpp, that compile and test okay.
I have a few questions / remarks though:
1. Why is the EPSG definition ID of an affine transform not being exported:  ID["EPSG",9624]],
2. Same for the ID’s of the Affine transformation parameters: 8623, 8624, 8625, 8639, 8640 and 8641.
3. The USAGE description isn’t being included in the generated WKT string and is also missing in the derived transform
   USAGE[
   SCOPE["Engineering survey, topographic mapping."],
   AREA["Netherlands - onshore, including Waddenzee, Dutch Wadden Islands and 12-mile offshore coastal zone."],
   BBOX[50.75,3.2,53.7,7.22]]]  
This is really a pity. Is there an ‘easy’ way to fix this ?
4. Nevertheless, maybe the example below is more ‘meaningful’ than the current one, that is included in the test_crs.file.

Kind regards,

Bart


ProjectedCRSNNPtr createProjected_v2()
{
    auto dbContext = DatabaseContext::create();
    auto factory = AuthorityFactory::create(dbContext, "EPSG");
    auto crs = factory->createCoordinateReferenceSystem("28992");
    auto proj_crs = nn_dynamic_pointer_cast<ProjectedCRS>(crs);
    EXPECT_TRUE(proj_crs != nullptr);

    return NN_CHECK_ASSERT(proj_crs);
}

static DerivedProjectedCRSNNPtr createDerivedProjectedCRS_v2() {

    auto derivingConversion = Conversion::create
    (
        PropertyMap().set(IdentifiedObject::NAME_KEY, "Affine"),
        PropertyMap().set(IdentifiedObject::NAME_KEY, "Affine parametric transformation"),
        std::vector<OperationParameterNNPtr>
        {
            OperationParameter::create(PropertyMap().set(
                IdentifiedObject::NAME_KEY, "A0")),
            OperationParameter::create(PropertyMap().set(
                IdentifiedObject::NAME_KEY, "A1")),
            OperationParameter::create(PropertyMap().set(
                IdentifiedObject::NAME_KEY, "A2")),
            OperationParameter::create(PropertyMap().set(
                IdentifiedObject::NAME_KEY, "B0")),
            OperationParameter::create(PropertyMap().set(
                IdentifiedObject::NAME_KEY, "B1")),
            OperationParameter::create(PropertyMap().set(
                IdentifiedObject::NAME_KEY, "B2")),
        },
        std::vector<ParameterValueNNPtr>
        {
            ParameterValue::create(Length(-257762.015)),
            ParameterValue::create(Scale ( 0.914845326707141)),
            ParameterValue::create(Scale ( 0.403804443019272)),
            ParameterValue::create(Length(-374179.972)),
            ParameterValue::create(Scale (-0.403804443019272)),
            ParameterValue::create(Scale ( 0.914845326707141)),
        }
    );

    // create a DerivedProjectedCRS, but ...
    // first we need to get the BaseCRS from an EPSG code.

    return DerivedProjectedCRS::create(
        PropertyMap().set(IdentifiedObject::NAME_KEY, "Abtswoudsebrug"),
        createProjected_v2(), derivingConversion,
        CartesianCS::createEastingNorthing(UnitOfMeasure::METRE));
}

TEST(crs, derivedProjectedCRS_WKT2_2019_v2) {

    auto expected_v2 =
        "DERIVEDPROJCRS[\"Abtswoudsebrug\",\n"
        "    BASEPROJCRS[\"Amersfoort / RD New\",\n"
        "        BASEGEOGCRS[\"Amersfoort\",\n"
        "            DATUM[\"Amersfoort\",\n"
        "                ELLIPSOID[\"Bessel 1841\",6377397.155,299.1528128,\n"
        "                    LENGTHUNIT[\"metre\",1]]],\n"
        "            PRIMEM[\"Greenwich\",0,\n"
        "                ANGLEUNIT[\"degree\",0.0174532925199433]]],\n"
        "        CONVERSION[\"RD New\",\n"
        "            METHOD[\"Oblique Stereographic\",\n"
        "                ID[\"EPSG\",9809]],\n"
        "            PARAMETER[\"Latitude of natural origin\",52.1561605555556,\n"
        "                ANGLEUNIT[\"degree\",0.0174532925199433],\n"
        "                ID[\"EPSG\",8801]],\n"
        "            PARAMETER[\"Longitude of natural origin\",5.38763888888889,\n"
        "                ANGLEUNIT[\"degree\",0.0174532925199433],\n"
        "                ID[\"EPSG\",8802]],\n"
        "            PARAMETER[\"Scale factor at natural origin\",0.9999079,\n"
        "                SCALEUNIT[\"unity\",1],\n"
        "                ID[\"EPSG\",8805]],\n"
        "            PARAMETER[\"False easting\",155000,\n"
        "                LENGTHUNIT[\"metre\",1],\n"
        "                ID[\"EPSG\",8806]],\n"
        "            PARAMETER[\"False northing\",463000,\n"
        "                LENGTHUNIT[\"metre\",1],\n"
        "                ID[\"EPSG\",8807]]]],\n"
        "    DERIVINGCONVERSION[\"Affine\",\n"
        "        METHOD[\"Affine parametric transformation\"],\n"
        "        PARAMETER[\"A0\",-257762.015,\n"
        "            LENGTHUNIT[\"metre\",1,\n"
        "                ID[\"EPSG\",9001]]],\n"
        "        PARAMETER[\"A1\",0.914845326707141,\n"
        "            SCALEUNIT[\"unity\",1,\n"
        "                ID[\"EPSG\",9201]]],\n"
        "        PARAMETER[\"A2\",0.403804443019272,\n"
        "            SCALEUNIT[\"unity\",1,\n"
        "                ID[\"EPSG\",9201]]],\n"
        "        PARAMETER[\"B0\",-374179.972,\n"
        "            LENGTHUNIT[\"metre\",1,\n"
        "                ID[\"EPSG\",9001]]],\n"
        "        PARAMETER[\"B1\",-0.403804443019272,\n"
        "            SCALEUNIT[\"unity\",1,\n"
        "                ID[\"EPSG\",9201]]],\n"
        "        PARAMETER[\"B2\",0.914845326707141,\n"
        "            SCALEUNIT[\"unity\",1,\n"
        "                ID[\"EPSG\",9201]]]],\n"
        "    CS[Cartesian,2],\n"
        "        AXIS[\"(E)\",east,\n"
        "            ORDER[1],\n"
        "            LENGTHUNIT[\"metre\",1,\n"
        "                ID[\"EPSG\",9001]]],\n"
        "        AXIS[\"(N)\",north,\n"
        "            ORDER[2],\n"
        "            LENGTHUNIT[\"metre\",1,\n"
        "                ID[\"EPSG\",9001]]]]";
//        "    USAGE[\",\n"
//        "        SCOPE[\"Engineering survey, topographic mapping.\"],\n"
//        "        AREA[\"Netherlands - onshore, including Waddenzee, Dutch Wadden Islands and 12-mile offshore coastal zone.\"],\n"
//        "        BBOX[50.75, 3.2, 53.7, 7.22]]]";

    auto crs = createDerivedProjectedCRS_v2();
    EXPECT_EQ(
        crs->exportToWKT(
            WKTFormatter::create(WKTFormatter::Convention::WKT2_2019).get()), expected_v2);

    EXPECT_TRUE(crs->isEquivalentTo(crs.get()));
    EXPECT_TRUE(crs->shallowClone()->isEquivalentTo(crs.get()));
    EXPECT_FALSE(crs->isEquivalentTo(createUnrelatedObject().get()));

    auto geodCRS = crs->extractGeodeticCRS();
    EXPECT_TRUE(geodCRS != nullptr);
}