Vektor réteg vizsgálata

Kezdésként indítsuk el a QGIS Python Konzolját. Írjunk egy rövid kódot, ami betölt egy vektoros réteget egy változóba:

path = "d:/geodata/"    
filename = "place_of_worship"    
layer = QgsVectorLayer(path + filename + ".shp", filename, "ogr")

A ‘path’ változóban adtuk meg a téradatot tartalmazó mappa elérési útját, ez azért szükséges, mert több tértadatot is tárolhatunk itt. A ‘filename’ változóba a fájl nevére hivatkozunk, melyet a QgsVectorLayer második argumentumaként is felhasználunk, így a betöltött réteg is ezt a nevet fogja kapni.

Egyelőre még nem adtuk hozzá a réteget a projektünkhöz, ezért érdemes letesztelni egyáltalán sikeres volt-e a betöltés:

path = "d:/geodata/"    
filename = "place_of_worship"    
layer = QgsVectorLayer(path + filename + ".shp", filename, "ogr")    
if not layer.isValid():    
  print("Layer failed to load!")

Ellenőrizzük, hogy kapunk-e visszajelzést. Ha nem sikerült a betöltés, akkor nem kapunk visszajelzést. Teszteljük a visszajelzést a fájlnév elrontásával. Adjuk hozzá a réteget a palettához:

path = "d:/geodata/"    
filename = "place_of_worship"    
layer = QgsVectorLayer(path + filename + ".shp", filename, "ogr")    
if not layer.isValid():    
  print("Layer failed to load!")  
QgsProject.instance().addMapLayer(layer)

Figyeljük meg, ha újrafuttatjuk a kódot, akkor még egyszer hozzáadja a QGIS a réteget a palettához, kerüljük ezt el egy plusz sorral:

path = "d:/geodata/"    
filename = "place_of_worship"    
layer = QgsVectorLayer(path + filename + ".shp", filename, "ogr")    
if not layer.isValid():    
  print("Layer failed to load!")    
if len(QgsProject.instance().mapLayersByName(filename)) == 0:    
QgsProject.instance().addMapLayer(layer)

Így csak akkor adjuk hozzá a réteget, ha még nincs a palettán a réteg nevével rendelkező réteg.

1. ábra: Magyarország templomai, forrás: OSM
1. ábra: Magyarország templomai, forrás: OSM

Vektor réteg metaadatai

Kérdezzük le az elemek számát egy változóba és irassuk ki:

fcount = layer.featureCount()    
print(fcount)

Nézzük meg a réteg vetületét:

crsname = layer.crs().description()    
print(crsname)

A vetület EPSG azonosítóját:

crscode = layer.crs().authid()    
print(crscode)

Illetve a vetület Proj4 paramétereit:

crsproj4 = layer.crs().toProj4()    
print(crsproj4)

Váltsunk magyarra, s gyűjtsük tovább a metaadatokat, nézzük meg a réteg terjedelmével kapcsolatos tulajdonságokat (attribútumokat):

layerext = layer.extent()    
print(u'"A réteg %s fok széles."' %layerext.width())

A speciális karakterkódolás miatt használtuk ‘u unicode jelzést, ettől függetlenül még mindig hibajelzést kapunk, definiáljuk a kódolást a scriptünk elején a következő sorral:

# -*- coding: utf-8 -*-

A szélesség könnyebben értelmezhető, ha átváltjuk hosszegységbe. Vegyük ki az adatokból levezethető közép-parallelkör szélességét:

print (u'"A templomok alapján a közép parallelkör szélessége %s fok."' %layerext.center()[1])

Egy geodétának a gömbi geometria nem okoz nehézséget. Számítsuk ki az említett parallelkör sugarát:

parradius = numpy.cos(layerext.center()[1]*numpy.pi/180)*6378

Ehhez viszont szükség van a numpy csomagra, ezért szúrjuk be ezt is a scriptünk elejére:

import numpy

Majd számítsuk ki a középponti szöghöz tartozó ívhosszat:

width_in_km = parradius * numpy.pi * 2 * layerext.width() / 360

Azaz, a templomok alapján - amik ránézésre jól lefedik Magyarországot - az ország szélessége közel 500 km. Számítsuk ki az ország magasságát is. Ez egy kicsit egyszerűbb, ugyanis a meridiánok identikus méretűek:

height_in_km = layerext.height() / 360 * 2 * 6378 * numpy.pi

Irassuk ki a befoglaló keret bal felső és jobb alsó sarokkordinátáit:

print(u'"A bal felső sarok koordinátája %s, %s. A jobb alsóé pedig %s, %s."' %(layerext.xMinimum(), layerext.yMaximum(), layerext.xMaximum(), layerext.yMinimum()))

Ennyivel befejeztük a réteg terjedelmével kapcsolatos sorokat. Nézzük meg milyen tulajdonságai vannak még a rétegnek:

dir(layer)

És milyen képességei vannak ennek a rétegnek, azaz mit lehet vele csinálni:

layer.capabilitiesString()

Kérdezzük le a réteg geometriájának típusát:

layer.wkbType()

Egy számot kapunk vissza, melynek a jelentése a következő: 0 = WKBUnknown, 1 = WKBPoint, 2 = WKBLineString, 3 = WKBPolygon, 4 = WKBMultiPoint, 5 = WKBMultiLineString, 6 = WKBMultiPolygon.

© Dr. Wirth Ervin