Main question
Using OGR in Python, how can I extract the integer-based geometry type code for a specific feature?
I know how to do this when dealing with a whole layer. For example, reading a layer from a GDB file would go like this:
from osgeo import ogr
gdb_path = '/path/to/my/file.gdb'
layer_name = 'my_input_layer'
driver_name = "OpenFileGDB"
driver = ogr.GetDriverByName(driver_name)
data_source = driver.Open(gdb_path, 0)
layer = data_source.GetLayerByName(layer_name)
layer_geom_type = layer.GetGeomType()
print(layer_geom_type)
# 2005
# This integer indicates that the layer I was reading contained MultiLineString geometries.
Getting the text-based output
However, I can't seem to find the analogous operation for dealing with the actual OGR geometries/features. Since OGR geometries don't have the GetGeomType
method, I can only get the text-based descriptions using the GetGeometryName
method. But can't get the integer-based codes. Here's my best approximation:
mls_m_wkt = "MULTILINESTRING M((0 0 100, 1 0 200, 2 0 300),(3 0 400, 4 0 500))"
mls_m_ogr = ogr.CreateGeometryFromWkt(mls_m_wkt)
geom_type = mls_m_ogr.GetGeometryName()
print(geom_type)
# MULTILINESTRING
# Notice that the output is a string with the text 'MULTILINESTRING' instead of
# the integer `2005`.
Is it even possible to extract the integer-based code directly from an OGR geometry object?
Getting the raw WKB
I know that the codes I'm interested in extracting are related to the WKB (well-known binary) codes for each of the geometries. On Wikipedia, you can find a table that tells you which integer codes are associated with which geometry types (see a snapshot from that table copied below). And I can even extract the raw WKB data from the geometry with a simple mls_m_ogr.ExportToWkb()
. But I'm not sure where to go from there and how to infer those geometry type codes from the raw WKB.
2 Answers 2
After poking around a little longer, I found out that while OGR geometries don't have access to the GetGeomType
method, they DO have access to the GetGeometryType
method, which outputs exactly what I'm looking for. I'm not sure why there is this discrepancy in the code base, but I'm happy that there's a way to solve my problem!
mls_m_wkt = "MULTILINESTRING M((0 0 100, 1 0 200, 2 0 300),(3 0 400, 4 0 500))"
mls_m_ogr = ogr.CreateGeometryFromWkt(mls_m_wkt)
geom_type = mls_m_ogr.GetGeometryType()
print(geom_type)
# 2005
# Now we see the integer code, as desired.
I'll mark my own question as answered. Who knows, maybe someone else will come up against this in the future and find this useful!
If you want solutions with the ogr raw WKB:
from osgeo import ogr
mls_m_wkt = "MULTILINESTRING M((0 0 100, 1 0 200, 2 0 300),(3 0 400, 4 0 500))"
mls_m_ogr = ogr.CreateGeometryFromWkt(mls_m_wkt)
mls_m_wkb = mls_m_ogr.ExportToWkb()
wk = mls_m_wkb.hex() # transform bytearray to hexstring
a) with the struct module (used in parse_wkb for example)
From Wikipedia:
The first byte indicates the byte order for the data:
mls_m_wkb[0]
0 # 00 : big endian
The next 4 bytes are a 32-bit unsigned integer for the geometry type
from struct import unpack
geomType = unpack('>i', mls_m_wkb[1:][:4])[0] # > = big endian, i = integer)
print(geomType)
2005
b) with the binascii module (used in Decode SRID information from a PostGIS WKB string in Python using the standard library adapted from the QuickWKT QGIS plugin)
import binascii
def decodeBinary(wkb):
"""Decode the binary wkb and return as a hex string"""
value = binascii.a2b_hex(wkb)
value = value[::-1]
value = binascii.b2a_hex(value)
return value.decode("UTF-8")
wkb= decodeBinary(wk[2:10])
geomType = int("0x" + decodeBinary(wkb), 0)
print(geomType)
2005
-
Hi @gene, thanks so much for the answer! I'm not sure why, but when I run these methods, I'm not getting the expected
2005
result. For method (a) I get-720961536
and for method (b) I get3574005760
. And I'm running this starting from a fresh new python kernel, so nothing wacky is happening with variables containing data stored from prior calculations and being accidentally pulled in. Super strange... Any idea of what's up? I'm running withpython 3.12.2
, andgdal 3.8.4
.Felipe D.– Felipe D.2024年04月03日 16:26:05 +00:00Commented Apr 3, 2024 at 16:26
Explore related questions
See similar questions with these tags.