2

In order to apply shapely.affinity.scale (see documentation) I have had to transform a series of geodetic coordinates to projected coordinates, see python - UserWarning: Geometry is in a geographic CRS. Results from 'buffer' are likely incorrect - Stack Overflow and UserWarning: Geometry is in a geographic CRS. Results from 'buffer' are likely incorrect for context

All my coordinates are within a bounding box defined as 35.5, -11.0, 71.5, 31.5 (South, West, North, East) covering most of Europe.

Transformation of e.g. 54.1767° latitude and 7.8911° longitude from EPSG:4326 (WGS 84) to EPSG:4936 (ETRS89) results in X: 3705858.6428144597 and Y: 513643.6629113175. But the reverse transformation results in 0° latitude and 7.8911° longitude! So something is wrong.

Transformation from EPSG:4326 to EPSG:3034 (ETRS89-extended / LCC Europe) and back to EPSG:4326 works just fine.

Both coordinate systems (EPSG:4936 and EPSG:3034) cover the same area (WGS84 bounds: -16.1 32.88, 40.18 84.73) and the coordinates used fall within these bounds.

What limitation of coordinates transformations am I missing?

UPDATE: a Minimal, Reproducible Example in Python

In [1]: import geopandas as gpd 
        import pandas as pd 
        import shapely 

In [2]: point = shapely.Point([7.8911,54.1767]) 

In [3]: gdf = gpd.GeoDataFrame(index=[0], crs="epsg:4326", 
                               geometry=[point]) 

In [4]: # transform point from WGS 84 to ETRS89 
        gdf["geometry"].to_crs(4936)
Out[4]: 0    POINT (3705858.64281 513643.66291)
        Name: geometry, dtype: geometry

In [5]: # transform point from WGS 84 to ETRS89 and back to WGS 84 
        gdf["geometry"].to_crs(4936).to_crs(4326) 
Out[5]: 0    POINT (7.89110 0.00000)
        Name: geometry, dtype: geometry

In [6]: # transform point from WGS 84 to ETRS89-extended 
        gdf["geometry"].to_crs(3034)
Out[6]: 0    POINT (3866814.669 3036011.045)
        Name: geometry, dtype: geometry

In [7]: # transform point from WGS 84 to ETRS89-extended and back to WGS 84 
        gdf["geometry"].to_crs(3034).to_crs(4326) 
Out[7]: 0    POINT (7.89110 54.17670)
        Name: geometry, dtype: geometry
Improve this question
5
  • simply all kinds of transformation come with some loss on accuracy, in particular when changing the so-called geodetic datum. Transforming twice will therefor double the error.
    – MakePeaceGreatAgain
    Mar 7 at 21:24
  • @MakePeaceGreatAgain I understand that, but 54.34707° to 0° latitude is not a rounding error.
    – marianoju
    Mar 7 at 21:33
  • 1
    On the inverse the values need to be swapped. What you have labeled as X is the northing values and Y is the easting value.
    – mkennedy
    Mar 19 at 23:30
  • @mkennedy I do not understand your comment. But I have simplified the MRE to make my original question easier to read. Which values exactly do you think need to be swapped?
    – marianoju
    Mar 20 at 22:06
  • 1
    OH, I should have checked the WKIDs. 4936 is the XYZ (3D geocentric/ECEF) version of ETRS89 so you should be printing 3 values for the output coordinates. I think dropping the 3rd value is what causes the very different results in the inverse calculation. I'd assumed 4936 was a UTM or GK projected CRS.
    – mkennedy
    Mar 21 at 15:15

1 Answer 1

Reset to default
1

Both EPSG:4326 (WGS 84) and EPSG:3034 (ETRS89-extended / LCC Europe) are 2D cartesian coordinate systems. EPSG:4936 (ETRS89) is a (geocentric) 3D cartesian coordinate system:

Cartesian 3D CS (geocentric). Axes: geocentric X Y Z. Orientations: X and Y in equatorial plane, X positive through intersection with prime meridian, Y through 0°N 90°E. Z axis parallel to mean earth rotation axis and positive towards North Pole. UoM: m.

The transformation of coordinates into a 3D coordinate system requires three parameters. A missing parameter leads to an error during the transformation.

Thanks to mkennedy for pointing me in the right direction!

Improve this answer

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.