How do I convert WKT geometry from 3D to 2D

Question

I'm working with kml files using Python gdal.ogr module. When I extract feature's geometry using feature.GetGeometryRef(), I get:

POLYGON ((44.184069865521963 41.971859752507207 0,44.355248797385123 41.976261187945859 0, ...

However, I do not need this Z value - I need just a plain 2D WKT. After hours of searching I found FlattenTo2D() method, but feature.GetGeometryRef().FlattenTo2D() returns None

So how do I convert 3D WKT to 2D?

Answer 1

Using of FlattenTo2D() method is correct way, but the method does not returns new object (like you expected) - it changes current object. So you should use the method like this:

g = ogr.CreateGeometryFromWkt("POINT(1 1 1)")
g.FlattenTo2D()
print(g.ExportToWkt())
>> 'POINT (1 1)'

Answer 2

WKT is only a text markup language for representing vector geometry, therefore you cannot convert 3D WKT to 2D (= text to text), you need to change the geometry

The problem with ogr is that it seems that all the geometries are 3D by defaut

1) Creation of a 2D point with ogr

point = ogr.Geometry(ogr.wkbPoint)
# 2D point
point.AddPoint(1198054.34, 648493.09)
print point.ExportToWkt()
POINT (1198054.34 648493.09 0)
point.GetGeometryType()==ogr.wkbPoint
False
point.GetGeometryType()==ogr.wkbPoint25D
True

2) But, with the creation of a 2D point with Shapely

from shapely.geometry import Point
pt = Point(1198054.34, 648493.09)
# conversion to an ogr geometry
point = ogr.CreateGeometryFromWkb(pt.wkb)
print(point.ExportToWkt())
POINT (1198054.34 648493.09)
point.GetGeometryType()==ogr.wkbPoint
True
point.GetGeometryType()==ogr.wkbPoint25D
False

3) Therefore it is possible to modify the geometry and an easy solution is to use the GeoJson format (feature.ExportToJson())

# The resulting WKT 
feature.ExportToWkt()
'Polygon ((250325.46051841106964275 142166.43077902274671942 0, 250775.29772301684715785 142019.67537893861299381 0, 250906.10144917882280424 141828.25529187230858952 0, 250555.1646228906174656 141550.69616562619921751 0, 250309.5088444888824597 141238.04335675123729743 0, 249933.04933992517180741 141777.20993532129796222 0, 250325.46051841106964275 142166.43077902274671942 0))'

json_geom = feature.ExportToJson()
print json_geom
'{ "type": "Polygon", "coordinates": [ [ [ 250325.460518411069643, 142166.430779022746719, 0.0 ], [ 250775.297723016847158, 142019.675378938612994, 0.0 ], [ 250906.101449178822804, 141828.25529187230859, 0.0 ], [ 250555.164622890617466, 141550.696165626199218, 0.0 ], [ 250309.50884448888246, 141238.043356751237297, 0.0 ], [ 249933.049339925171807, 141777.209935321297962, 0.0 ], [ 250325.460518411069643, 142166.430779022746719, 0.0 ] ] ] }'
# convert json format to a Python dictionary
import json
geom = json.loads(json_geom)
print geom['coordinates']
[[[250325.4605184111, 142166.4307790227, 100.0], [250775.2977230168, 142019.6753789386, 150.0], [250906.1014491788, 141828.2552918723, 200.0], [250555.1646228906, 141550.6961656262, 230.0], [250309.5088444889, 141238.0433567512, 340.0], [249933.0493399252, 141777.2099353213, 220.0], [250325.4605184111, 142166.4307790227, 100.0]]]
# slice the coordinates to eliminate 3D
new_coords =  [[i[:2] for i in geom['coordinates'][0]]]
print new_coords
[[[250325.4605184111, 142166.4307790227], [250775.2977230168, 142019.6753789386], [250906.1014491788, 141828.2552918723], [250555.1646228906, 141550.6961656262], [250309.5088444889, 141238.0433567512], [249933.0493399252, 141777.2099353213], [250325.4605184111, 142166.4307790227]]]
new = { "type": "Polygon",'coordinates': new_coords}
print new
{'type': 'Polygon', 'coordinates': [[[250325.4605184111, 142166.4307790227], [250775.2977230168, 142019.6753789386], [250906.1014491788, 141828.2552918723], [250555.1646228906, 141550.6961656262], [250309.5088444889, 141238.0433567512], [249933.0493399252, 141777.2099353213], [250325.4605184111, 142166.4307790227]]]}
# recreate the ogr geometry
poly = ogr.CreateGeometryFromJson(json.dumps(new))
poly.ExportToWkt() 
'POLYGON ((250325.460518411098747 142166.430779022688512,250775.29772301678895 142019.675378938612994,250906.1014491787937 141828.25529187230859,250555.164622890588362 141550.696165626199218,250309.508844488911564 141238.043356751208194,249933.049339925200911 141777.209935321297962,250325.460518411098747 142166.430779022688512))'

You can also use geomet to transform the correct GeoJson to WKT

from geomet import wkt
wkt.dumps(new, decimals=4)
'POLYGON ((250325.4605 142166.4308, 250775.2977 142019.6754, 250906.1014 141828.2553, 250555.1646 141550.6962, 250309.5088 141238.0434, 249933.0493 141777.2099, 250325.4605 142166.4308))'

Answer 3

Though quite an old post but I was having the same issue and this helped me:

ring = ogr.Geometry(ogr.wkbLinearRing)

ring.AddPoint(coordinates['xmin'], coordinates['ymin'])
ring.AddPoint(coordinates['xmin'], coordinates['ymax'])
ring.AddPoint(coordinates['xmax'], coordinates['ymax'])
ring.AddPoint(coordinates['xmax'], coordinates['ymin'])
ring.CloseRings()

geom = ogr.Geometry(ogr.wkbPolygon)
geom.AddGeometry(ring)

geom.FlattenTo2D()

Thanks to programcreek