abstract architecture background contemporary
|

PyQGIS: Get Feature Geometry

ByKonrad Hafen

Geometry is what defines the shape and spatial location of a feature. Features can be points, lines, or polygons, can be anywhere on earth, and can be represented by coordinates in a number of different units (depending on the coordinate reference system, or projection). Points are represented by a single coordinate pair while lines and polygons are represented by multiple coordinate pairs.

Copying, or updating the shape or location of a feature requires you to access the feature’s geometry. This tutorial will show the basics of accessing feature geometry with PyQGIS. This example uses a line feature, but can be easily adapted for lines or polygons.

To start, define the path to the file you want to open, open that file as a QGIS vector layer, then get the features in the layer as a list of objects. 

fn = "path/to/shapefile.shp"
layer = QgsVectorLayer(fn, '', 'ogr')
feats = layer.getFeatures()

You have to get the geometry of each feature individually, so loop through each feature in the list of features. A feature’s geometry can be accessed by calling geometry().

for feat in feats:
    geom = feat.geometry()

If you already know exactly what type your geometry is, you can get the shape’s final geometry with geom.asPolyline() or geom.asMultiline(). However, if you’re unsure of the layer’s geometry type, or if it is single-part or multi-part you’ll need a little more code. It’s always safest to check for the geometry type in the code to avoid any errors.

The geometry’s type can easily be accessed with geom.wkbType(). WKB stands for ‘well-known binary’ format. Geometry type can be either single type — one shape (line, point, polygon) per geometry entry– or multi type — more than one shape per geometry. This can be assessed by passing the geometry type to isSingleType, which returns True or False.

for feat in feats:
    geom = feat.geometry()
    geomSingleType = QgsWkbTypes.isSingleType(geom.wkbType())

Now we’re going to make two checks using if statements. First, we’ll check to be sure the geometry type is a Line Geometry. Second, we’ll check if the geometry is single-part or multi-part using the geomSingleType variable we created above. Note that QgsWkbTypes.LineGeometry can be substituted for points or polygons with PointGeometry or PolygonGeometry.

for feat in feats:
    geom = feat.geometry()
    geomSingleType = QgsWkbTypes.isSingleType(geom.wkbType())
    if geom.type() == QgsWkbTypes.LineGeometry:
        if geomSingleType:

Finally, if the geometry is single-type (single-part) you can access the geometry with asPolyline(), or if it is multi-type (multi-part) with asMultiPolyline(). The print statement prints the resulting geometry.

for feat in feats:
    geom = feat.geometry()
    geomSingleType = QgsWkbTypes.isSingleType(geom.wkbType())
    if geom.type() == QgsWkbTypes.LineGeometry:
        if geomSingleType:
            x = geom.asPolyline()
            print("line:", x)
        else:
            x = geom.asMultiPolyline()
            print("multiline:", x)

Here is the final code block we have created. Watch the video below for a step-by-step demonstration of the process.

fn = "path/to/shapefile.shp"
layer = QgsVectorLayer(fn, '', 'ogr')
feats = layer.getFeatures()

for feat in feats:
    geom = feat.geometry()
    geomSingleType = QgsWkbTypes.isSingleType(geom.wkbType())
    if geom.type() == QgsWkbTypes.LineGeometry:
        if geomSingleType:
            x = geom.asPolyline()
            print("line:", x)
        else:
            x = geom.asMultiPolyline()
            print("multiline:", x)

Similar Posts