C# transform coordinates using specific transform code using GDAL

Question

I have an application that can import / export geospatial data to / from an ESRI SDE Database (SQLServer). The data is transformed to WGS84 before it is inserted into the DB and again from WGS84 to the projected Coordinate System the data was originally in when exported.

I use EPSG codes for the coordinate system and the transformation but when I come to transform the data I can't seem to specify the transformation using the ESPG code.

I use GDAL v3.4.1 packages (GDAL & GDAL.Native) from NUGET

Here is the code I use for creating the CoordinateTransformation -

public static CoordinateTransformation getTransformationParameters(int TransformationSR, int FromSRCode, int ToSRCode)
    {
        // Create SpatialReference from ESPG codes one will always be 4326 depending on direction
        SpatialReference FromSR = new SpatialReference(null);
        FromSR.ImportFromEPSG(FromSRCode);
        SpatialReference ToSR = new SpatialReference(null);
        ToSR.ImportFromEPSG(ToSRCode);

        // Create TranformationOptions 
        CoordinateTransformationOptions TransformationOptions = new CoordinateTransformationOptions();
        TransformationOptions.SetOperation($@"urn:ogc:def:coordinateOperation:EPSG::{TransformationSR}");
        CoordinateTransformation Transformation = new CoordinateTransformation(FromSR, ToSR, TransformationOptions);

        return Transformation;
    }

If I remove the TransformationOptions and let GDAL use the normal computed pipeline it transforms the coordinates using whatever it thinks the best transformation is, but I have to be able to transform the coordinate using very specific transformations and I just can't get this working.

The ESPG codes that I am using are as follows -

FromSRCode - 4326 (WGS84)
ToSRCode - 23090 (ED50 TM0)
TransformationSR - 1311 (ED50 to WGS 84 (18))

I'd greatly appricate any help.

Answer 1

It has to go through a GEOGCS conversion first, in my example from ED50 TM0 (EPSG:23090) > ED50 (EPSG:4230) then run the URN pipeline transformation from ED50 (EPSG:4230) > WGS84 (4326).

public static CoordinateTransformation getTransformationParameters(SpatialReference FromSRCode, SpatialReference ToSRCode, int TransformationSR = 0)
    {
        if (TransformationSR == 0)
            return new CoordinateTransformation(FromSRCode, ToSRCode);

        CoordinateTransformationOptions TransformationOptions = new CoordinateTransformationOptions();

        if (TransformationOptions.SetOperation($"urn:ogc:def:coordinateOperation:EPSG::{TransformationSR}"))
            return new CoordinateTransformation(FromSRCode, ToSRCode, TransformationOptions);
        else
            throw new Exception("Exception creating CoordinateTranformation Pipeline");
    }

To create the correct transformation between ED50TM0 (EPSG:23030) and WGS84 (EPSG:4326) -

    private static Geometry transform(CoordinateTransformation coordinateTransformation, Geometry sqlGeom)//, out List<double> z, out List<double> m, out double[] translationArrayX, out double[] translationArrayY)
    {
        if (sqlGeom.Transform(coordinateTransformation) == Ogr.OGRERR_NONE)
            return sqlGeom;
        throw new Exception("ERROR when tranforming");
    }

    public static Geometry transformToWGS84(SpatialReference ProjectSR, int TransformationSR, wkbGeometryType geomType, SqlGeometry sqlGeom)
    {
        Geometry GeogTransform = transform(getTransformationParameters(ProjectSR, ProjectSR.CloneGeogCS()), sqlGeom.toOSRGeometry(ProjectSR));
        Geometry ProjTransform = transform(getTransformationParameters(ProjectSR.CloneGeogCS(), WGS84CRS, TransformationSR), GeogTransform);
        return ProjTransform;
    }

One last thing, when creating any SpatialReference use Traditional GIS Axis order for the coordinates -

public static SpatialReference getSpatialReference(int projectSR)
    {
        SpatialReference spatialReference = new SpatialReference("");
        if (projectSR != 0)
        {
            if (spatialReference.SetFromUserInput($"urn:ogc:def:crs:EPSG::{projectSR}") == Ogr.OGRERR_NONE)
            {
                spatialReference.SetAxisMappingStrategy(AxisMappingStrategy.OAMS_TRADITIONAL_GIS_ORDER);
                return spatialReference;
            }
            throw new Exception(string.Format("ERROR getting spatial reference {0}", projectSR));
        }
        return spatialReference;
    }