Importante
unireLa traduzione è uno sforzo comunitario you can join. Questa pagina è attualmente tradotta al 100.00%.
7. Gestione della Geometria
Suggerimento
I frammenti di codice di questa pagina necessitano delle seguenti importazioni se sei è al di fuori della console pyqgis:
1from qgis.core import (
2 QgsGeometry,
3 QgsGeometryCollection,
4 QgsPoint,
5 QgsPointXY,
6 QgsWkbTypes,
7 QgsProject,
8 QgsFeatureRequest,
9 QgsVectorLayer,
10 QgsDistanceArea,
11 QgsUnitTypes,
12 QgsCoordinateTransform,
13 QgsCoordinateReferenceSystem
14)
Punti, linee e poligoni che rappresentano un elemento spaziale sono comunemente indicati come geometrie. In QGIS sono rappresentati con la QgsGeometry
class.
Alcune volte una geometria é effettivamente una collezione di geometrie (parti singole) piú semplici. Se contiene un tipo di geometria semplice, la chiameremo punti multipli, string multi linea o poligoni multipli. Ad esempio, un Paese formato da piú isole puó essere rappresentato come un poligono multiplo.
Le coordinate delle geometrie possono essere in qualsiasi sistema di riferimento delle coordinate (CRS). Quando si estraggono delle caratteristiche da un vettore, le geometrie associate avranno le coordinate nel CRS del vettore.
La descrizione e le specifiche di tutte le possibili geometrie e relazioni sono disponibili per dettagli avanzati nel documento OGC Simple Feature Access Standards .
7.1. Costruzione della Geometria
PyQGIS offre diverse opzioni per la creazione di una geometria:
dalle coordinate
1gPnt = QgsGeometry.fromPointXY(QgsPointXY(1,1)) 2print(gPnt) 3gLine = QgsGeometry.fromPolyline([QgsPoint(1, 1), QgsPoint(2, 2)]) 4print(gLine) 5gPolygon = QgsGeometry.fromPolygonXY([[QgsPointXY(1, 1), 6 QgsPointXY(2, 2), QgsPointXY(2, 1)]]) 7print(gPolygon)
Le coordinate vengono fornite utilizzando la classe
QgsPoint
o la classeQgsPointXY
. La differenza tra queste classi è cheQgsPoint
supporta le dimensioni M e Z.Una Polilinea (Linestring) è rappresentata da un elenco di punti.
Un Poligono è rappresentato da un elenco di anelli lineari (cioè di linee chiuse). Il primo anello è l’anello esterno (confine), mentre eventuali anelli successivi sono buchi nel poligono. Da notare che, a differenza di alcuni programmi, QGIS chiude l’anello per voi, quindi non è necessario duplicare il primo punto come ultimo.
Le geometrie a parti multiple vanno ad un livello successivo: punti multipli é una lista di punti, una stringa multi linea é una linea di linee ed un poligono multiplo é una lista di poligoni.
da well-known text (WKT)
geom = QgsGeometry.fromWkt("POINT(3 4)") print(geom)
da well-known binary (WKB)
1g = QgsGeometry() 2wkb = bytes.fromhex("010100000000000000000045400000000000001440") 3g.fromWkb(wkb) 4 5# print WKT representation of the geometry 6print(g.asWkt())
7.2. Accedere alla Geometria
Per prima cosa, devi trovare il tipo di geometria. Il metodo wkbType()
è quello da utilizzare. Restituisce un valore dell’enumerazione QgsWkbTypes.Type
.
1print(gPnt.wkbType())
2# output: 'WkbType.Point'
3print(gLine.wkbType())
4# output: 'WkbType.LineString'
5print(gPolygon.wkbType())
6# output: 'WkbType.Polygon'
In alternativa, si può usare il metodo type()
che restituisce un valore dell’enumerazione QgsWkbTypes.GeometryType
.
print(gLine.type())
# output: 'GeometryType.Line'
Puoi usare la funzione displayString()
per ottenere un tipo di geometria leggibile in chiaro.
1print(QgsWkbTypes.displayString(gPnt.wkbType()))
2# output: 'Point'
3print(QgsWkbTypes.displayString(gLine.wkbType()))
4# output: 'LineString'
5print(QgsWkbTypes.displayString(gPolygon.wkbType()))
6# output: 'Polygon'
Esiste anche una funzione di aiuto isMultipart()
per scoprire se una geometria è multipart o meno.
Per estrarre informazioni da una geometria esistono funzioni di selezione per ogni tipologia di vettore. Ecco un esempio di utilizzo di queste funzioni di selezione:
1print(gPnt.asPoint())
2# output: <QgsPointXY: POINT(1 1)>
3print(gLine.asPolyline())
4# output: [<QgsPointXY: POINT(1 1)>, <QgsPointXY: POINT(2 2)>]
5print(gPolygon.asPolygon())
6# output: [[<QgsPointXY: POINT(1 1)>, <QgsPointXY: POINT(2 2)>, <QgsPointXY: POINT(2 1)>, <QgsPointXY: POINT(1 1)>]]
Nota
Le tuple (x,y) non sono tuple reali, ma oggetti QgsPoint
, i cui valori sono accessibili con i metodi x()
e y()
.
Per le geometrie multiparte esistono funzioni di accesso simili: asMultiPoint()
, asMultiPolyline()
e asMultiPolygon()
.
È possibile iterare su tutte le parti di una geometria, indipendentemente dal tipo di geometria. Ad esempio
geom = QgsGeometry.fromWkt( 'MultiPoint( 0 0, 1 1, 2 2)' )
for part in geom.parts():
print(part.asWkt())
Point (0 0)
Point (1 1)
Point (2 2)
geom = QgsGeometry.fromWkt( 'LineString( 0 0, 10 10 )' )
for part in geom.parts():
print(part.asWkt())
LineString (0 0, 10 10)
gc = QgsGeometryCollection()
gc.fromWkt('GeometryCollection( Point(1 2), Point(11 12), LineString(33 34, 44 45))')
print(gc[1].asWkt())
Point (11 12)
È anche possibile modificare ciascuna parte della geometria usando il metodo QgsGeometry.parts()
.
1geom = QgsGeometry.fromWkt( 'MultiPoint( 0 0, 1 1, 2 2)' )
2for part in geom.parts():
3 part.transform(QgsCoordinateTransform(
4 QgsCoordinateReferenceSystem("EPSG:4326"),
5 QgsCoordinateReferenceSystem("EPSG:3111"),
6 QgsProject.instance())
7 )
8
9print(geom.asWkt())
MultiPoint ((-10334728.12541878595948219 -5360106.25905461423099041),(-10462135.16126426123082638 -5217485.4735023295506835),(-10589399.84444035589694977 -5072021.45942386891692877))
7.3. Predicati ed Operazioni delle Geometrie
QGIS utilizza la libreria GEOS per operazioni geometriche avanzate come i predicati geometrici (contains()
, intersects()
, …) e operazioni di insieme (combine()
, difference()
, …). Può anche calcolare le proprietà geometriche delle geometrie, come l’area (nel caso dei poligoni) o le lunghezze (per poligoni e linee).
Vediamo un esempio che combina l’iterazione degli elementi in un dato layer e l’esecuzione di alcuni calcoli geometrici basati sulle loro geometrie. Il codice seguente calcolerà e stamperà l’area e il perimetro di ogni countrie nel layer ``countries”” del nostro progetto QGIS.
Il codice seguente presuppone che layer
sia un oggetto QgsVectorLayer
che ha una tipologia Poligono.
1# let's access the 'countries' layer
2layer = QgsProject.instance().mapLayersByName('countries')[0]
3
4# let's filter for countries that begin with Z, then get their features
5query = '"name" LIKE \'Z%\''
6features = layer.getFeatures(QgsFeatureRequest().setFilterExpression(query))
7
8# now loop through the features, perform geometry computation and print the results
9for f in features:
10 geom = f.geometry()
11 name = f.attribute('NAME')
12 print(name)
13 print('Area: ', geom.area())
14 print('Perimeter: ', geom.length())
1Zambia
2Area: 62.82279065343119
3Perimeter: 50.65232014052552
4Zimbabwe
5Area: 33.41113559136517
6Perimeter: 26.608288555013935
Ora hai calcolato e stampato le aree e i perimetri delle geometrie. Tuttavia, puoi subito notare che i valori sono strani. Questo perché le aree e i perimetri non tengono conto dei SR quando vengono calcolati con i metodi area()
e length()
della classe QgsGeometry
. Per un calcolo più preciso di aree e distanze, occorre utilizzare la classe QgsDistanceArea
, che può eseguire calcoli basati su ellissoidi:
Il codice seguente presuppone che layer
sia un oggetto QgsVectorLayer
che ha una tipologia Poligono.
1d = QgsDistanceArea()
2d.setEllipsoid('WGS84')
3
4layer = QgsProject.instance().mapLayersByName('countries')[0]
5
6# let's filter for countries that begin with Z, then get their features
7query = '"name" LIKE \'Z%\''
8features = layer.getFeatures(QgsFeatureRequest().setFilterExpression(query))
9
10for f in features:
11 geom = f.geometry()
12 name = f.attribute('NAME')
13 print(name)
14 print("Perimeter (m):", d.measurePerimeter(geom))
15 print("Area (m2):", d.measureArea(geom))
16
17 # let's calculate and print the area again, but this time in square kilometers
18 print("Area (km2):", d.convertAreaMeasurement(d.measureArea(geom), QgsUnitTypes.AreaSquareKilometers))
1Zambia
2Perimeter (m): 5539361.250294601
3Area (m2): 751989035032.9031
4Area (km2): 751989.0350329031
5Zimbabwe
6Perimeter (m): 2865021.3325076113
7Area (m2): 389267821381.6008
8Area (km2): 389267.8213816008
In alternativa, potresti voler conoscere la distanza tra due punti.
1d = QgsDistanceArea()
2d.setEllipsoid('WGS84')
3
4# Let's create two points.
5# Santa claus is a workaholic and needs a summer break,
6# lets see how far is Tenerife from his home
7santa = QgsPointXY(25.847899, 66.543456)
8tenerife = QgsPointXY(-16.5735, 28.0443)
9
10print("Distance in meters: ", d.measureLine(santa, tenerife))
Puoi trovare molti esempi di algoritmi che sono inclusi in QGIS ed utilizzare questi metodi per analizzare e trasformare i dati vettoriali. Di seguito i link al codice di alcuni di questi.
Distanza e area utilizzando la classe
QgsDistanceArea
: Distance matrix algorithm.