Custom datum transformation in Python

Question

I need to transform a shapefile from Deir ez Zor / Levant Stereographic (EPSG:22780) to UTM 36N using 7-parameters Helmet transformation on Python script.

In ArcGIS, I can achieve the following using this code:

GEOGTRAN["from levant to wgs",
    GEOGCS["GCS_levant",
        DATUM["D_levant",
            SPHEROID["Clarke_1880_IGN",6378249.2,293.46602]],
        PRIMEM["Greenwich",0.0],
        UNIT["Degree",0.0174532925199433]],
    GEOGCS["GCS_WGS_1984",
        DATUM["D_WGS_1984",
            SPHEROID["WGS_1984",6378137.0,298.257223563]],
        PRIMEM["Greenwich",0.0],
        UNIT["Degree",0.0174532925199433]],
    METHOD["Position_Vector"],
    PARAMETER["X_Axis_Translation",xxxxx],
    PARAMETER["Y_Axis_Translation",xxxxx],
    PARAMETER["Z_Axis_Translation",xxxxx],
    PARAMETER["X_Axis_Rotation",xxxxx],
    PARAMETER["Y_Axis_Rotation",xxxxx],
    PARAMETER["Z_Axis_Rotation",xxxxx],
    PARAMETER["Scale_Difference",xxxxx]]

I experimented with several commands such as shpproj with no luck.

Answer 1

If all else fails - read in the file using Shapely and reproject each vertex using the pyproj library and write back to new shapefile. Example use of pyproj, where you can edit the custom transformation:

from pyproj import Transformer
from pyproj import CRS
dlon = 33.0
dlat = 40.0
projWGS84=CRS.from_epsg(4326)
projUTM36N=CRS.from_epsg(32636)
proj22780="+proj=sterea +lat_0=34.2 +lon_0=39.15 +k=0.9995341 +x_0=0 +y_0=0 +a=6378249.2 +b=6356515 +towgs84=-190.421,8.532,238.69,0,0,0,0 +units=m +no_defs"
#sourced from https://epsg.io/22780 
#alternatively use:
#proj22780=CRS.from_epsg(22780)
print ("WGS84",dlon,dlat)
transformer = Transformer.from_crs(projWGS84, projUTM36N)
x,y=transformer.transform(dlon,dlat)
print ("UTM36N",x,y)
transformer = Transformer.from_crs(projUTM36N, proj22780)
x,y=transformer.transform(x,y)
print ("EPSG22780",x,y)

Has output:

WGS84 33.0 40.0
UTM36N 1154576.4163343282 3673130.830795592
EPSG22780 79280.46018819655 -132668.34403660995