Por Rober J
Importación de módulos
Al comienzo de un script se debe siempre importar los módulos que vayan a usarse:
## Módulos imprescindibles con las funciones básicas de QGIS from qgis.core import * from qgis.utils import *from PyQt5.QtCore import * from PyQt5.QtGui import *
from qgis.PyQt.QtCore import QSettings, QTranslator, QcoreApplication from qgis.PyQt.QtGui import Qicon from qgis.PyQt.QtWidgets import QAction, QfileDialog## recomendable
import os import shutil import processing from qgis.analysis import * ## Herramientas de analisis raster
Cargar capas ráster
Cargar una capa ráster implica crear un objeto QgsRasterLayer que contendrá los datos del archivo ráster. A dicho objeto se le podrán aplicar los distintos métodos de su clase para operar con sus datos.
Hay dos maneras de cargar las capas ráster usando PyQGIS dependiendo de si queremos visualizarlas o no. Si estamos usando Python dentro de QGIS probablemente queramos añadir las capas al panel de capas, pero en caso de usar un IDE externo (VSCode, PyCharm…) esto no será posible y solo podremos acceder a los datos de las capas sin poder visualizarlas
A diferencia de las capas vectoriales, no será necesario especificar el proveedor puesto que para ráster el único disponible es GDAL (Geospatial Data Abstraction Library) y es el equivalente a OGR para datos vectoriales, es decir, admite multitud de formatos ráster como TIFF, IMG…
Con el siguiente código, además de cargar en una variable llamada capa_raster los datos de una capa ráster, dicha capa se añadirá al panel de capas de QGIS
capa_raster = iface.addRasterLayer('ruta_capa','nombre_capa')
Esta función no funcionará si trabajamos con PyQGIS en una IDE externa.
Con el siguiente código se cargan solo los datos de la capa ráster en una variable nueva llamada capa_raster, de modo que se podrá usar para llevar a cabo geoprocesos pero no se visualizarán.
capa_raster = QgsRasterLayer('ruta_capa','nombre_capa')
Se puede usar tanto en una IDE como dentro de QGIS si no queremos añadir la capa a la visualización.
Al igual que con los datos vectoriales, debe usarse el método .isValid() sobre el objeto que guarda la capa. Combinado con una sentencia if podremos generar un mensaje en caso de error:
if not capa_raster.isValid():
print("Error al cargar la capa")
Extraer información de la capa
Una vez cargada la capa ráster en un objeto QgsRasterLayer podremos obtener información de ella a través de distintos métodos:
capa_raster.extent()
capa_raster.rasterType()
Podremos obtener 3 resultados:
- 0 = Escala de grises o indefinido (monobanda)
- 1 = Monobanda
- 2 = Multibanda
capa_raster.width() ## columnas capa_raster.height() ## filas
capa_raster.bandCount()
capa_raster.metadata()
El siguiente código usa el método .identify() sobre dataProvider() indicando un objeto QgsPoint. Si es correcto (las coordenadas del punto se encuentran dentro de la extensión espacial del ráster) imprime el valor del píxel que se encuentra en el objeto de tipo QgsIdentifyResult().
coordx = 5000 coordy = -3500 punto = QgsPointXY(coordx, coordy)valor_pixel = rlayer.dataProvider().identify(punto, QgsRaster.IdentifyFormatValue)
if valor_pixel.isValid(): print(valor_pixel.results())
También se puede extraer el valor de un píxel en forma de tupla con el método .sample() sobre dataProvider():
coordx = 5000 coordy = -3500 n_banda = 1 punto = QgsPointXY(coordx, coordy)valor_pixel, res = capa_raster.dataProvider().sample(punto, n_banda) print(valor_pixel)
Calculadora ráster
Las operaciones matemáticas sobre los píxeles de las capas ráster en QGIS se llevan a cabo en la calculadora ráster. En PyQGIS, la función detrás de esta herramienta es QgsRasterCalculator(). Sus parámetros son, por orden:
QgsRasterCalculator(expresion, capa salida, formato, extension, ancho, alto, capas de entrada)
| Parámetro | Descripción | Tipo de objeto |
|---|---|---|
| expresión | Expresión lógica con las operaciones que queremos realizar | Cadena de texto |
| capa de salida | Ruta con el nombre y la extensión de la capa en la que guardar los resultados | Cadena de texto |
| formato | Formato de la capa de salida | Cadena de texto |
| extensión | Coordenadas de la capa de salida | QgsRectangle |
| ancho (columnas) | Número de columnas de la capa de salida | Entero |
| alto (filas) | Número de filas de la capa de salida | Entero |
| capas de entrada | Lista de objetos que almacenan las capas que se usarán en la calculadora | QgsRasterCalculatorEntry |
El proceso para crear la lista con los objetos QgsRasterCalculatorEntry junto con el resto de parámetros y la ejecución de la calculadora es el siguiente:
## Cargar capas mi_raster = QgsRasterLayer(ruta, "nombre") mi_raster2 = QgsRasterLayer(ruta, "nombre2")tipo = 'GTiff' extension = mi_raster.extent() filas = mi_raster.height() columnas = mi_raster.width() expresion = 'capa1@1 * capa2@1' ## Multiplicar los valores de las dos capas de entrada capa_salida = 'C:\ruta\' ## Ruta con el nombre de la capa de salida y la extension capas_entrada = [] ## Lista vacia que almacenara los objetos de entrada
capa1 = QgsRasterCalculatorEntry() ## Crear el objeto de entrada capa1.ref = 'capa1@1' ## Crear referencia para la expresion capa1.raster = mi_raster ## Seleccion de la capa
capa1.bandNumber = 1 ## Definir el nº de banda a usar capas_entrada.append(capa1) ## Añadir a la lista vaciacapa2 = QgsRasterCalculatorEntry() capa2.ref = 'capa2@1' capa2.raster = mi_raster2 capa2.bandNumber = 3 capas_entrada.append(capa2)
calculos = QgsRasterCalculator(expresion, capa_salida, tipo, extension, columnas, filas, capas_entrada)
calculos.processCalculation()
Si una de las capas de entrada fuese una capa de máscara, es decir, con valores 1 y 0, al multiplicarlas recortaríamos la otra capa, pues los píxeles que se multipliquen por 0 perderían su valor.
El ejemplo anterior toma dos capas de entrada para sumar sus valores, pero podrían usarse más o usar tan solo una para reclasificar sus valores y convertir variables continuas en discretas. En tal caso habrá que adecuar la expresión:
## Los pixeles cuyo valor sea mayor a 1000 pasaran a valer 1 y el resto 0 expresion = 'capa@1 > 1000'expresion = 'capa@1 * 0.7'
expresion = '( capa1@1 < 100 ) * 1 + ( capa1@1 >= 100 AND capa1@1 < 500 ) * 2 + ( capa1@1 >= 500 AND capa1@1 < 1000 ) * 3 + ( capa1@1 >= 1000 ) * 4'
He dejado un ejemplo completo de reclasificación en GitHub en el que se lleva a cabo la reclasificación por intervalos para todos los TIF de una carpeta.
Cargar capas vectoriales
Principalmente existen dos funciones para cargar las capas vectoriales usando PyQGIS dependiendo de si queremos visualizarlas o no. En cualquier caso, lo que se hace es crear un objeto de tipo QgsVectorLayer, es decir, un objeto o variable que almacena una capa.
Crear esta clase de objetos es fundamental para poder acceder después a sus entidades, y a partir de ellas a sus atributos y su geometría.
Con la función iface.addVectorLayer() además de cargar en una variable nueva los datos de una capa vectorial, dicha capa se añadirá al panel de capas de QGIS.
Esta función no funcionará si trabajamos con PyQGIS en una IDE externa (VSCode, PyCharm…) .
capa_vectorial = iface.addVectorLayer('ruta_capa','nombre_capa','proveedor)
Con la función QgsVectorLayer() se cargan solo los datos de la capa vectorial en una variable nueva, de modo que se podrá usar para llevar a cabo geoprocesos pero no se visualizarán.
capa_vectorial = QgsVectorLayer('ruta_capa','nombre_capa','proveedor')
Se puede usar tanto en una IDE como dentro de QGIS si no queremos añadir la capa a la visualización.
Proveedores
El procedimiento para cargar capas vectoriales varía según el proveedor, puesto que el origen de datos varía. En cualquier caso, funcionan igual tanto con QgsVectorLayer() como con iface.addVectorLayer():
Es el proveedor para las capas vectoriales de los tipos de archivo más comunes (.shp, .kml o .geojson) y suele usarse para utilizar archivos que se encuentran en el sistema de archivos de nuestro ordenador.
capa_vectorial = QgsVectorLayer('ruta_capa','nombre_capa','ogr')
Sirve para cargar información vectorial almacenada en bases de datos de PostgreSQL o PostGIS.
# Conectar con la base de datos
uri = QgsDataSourceUri()
uri.setConnection("host", "puerto", "nombre_database", "usuario", "contraseña")
uri.setDataSource("esquema", "tabla", "campo_geometría")
capa = QgsVectorLayer(uri.uri(), "Nombre_capa", "postgres")
Para conectar con una base de datos de tipo MySQL se debe usar otro proveedor como OGR:
# Conectar con la base de datos uri = "MySQL:nombre_database,host=nombre_host,port=numero_puerto,user=nombre_usuario,password=contraseña|layername=nombre_tabla"capa = QgsVectorLayer(uri,"Nombre_capa","ogr" ),"Nombre_capa","ogr"
Para cargar servicios WFS (Web Feature Services) habrá que indicar una URL como origen de los datos
uri = "URL del servicio WFS" capa = QgsVectorLayer(uri, "Nombre_capa", "WFS")
Estos servicios también pueden cargarse usando el módulo de Python urllib:
import urllib
params = {
'service': 'WFS',
'version': '1.0.0', 'request': 'petición',
'typename': 'nombre_capa', 'srsname': "EPSG:XXXX"
}
uri =
'ruta_servicio?'+urllib.parse.unquote(urllib.parse.ur
lencode(params))
capa = QgsVectorLayer(uri, "Nombre_capa", "WFS")
Los servicios WFS están estandarizados por la OGC (Open Geospatial Consortium) y todos poseen los siguientes tipos de peticiones:
- GetCapabilities – obtener metadatos
- GetFeature – obtener entidades
- DescribeFeatureType – obtener esquema XML de los servicios del servidor del WFS
⚠ En QGIS 2 el proveedor que debe especificarse al cargar capas WFS debe ser OGR
Se puede también cargar archivos CSV o de texto delimitado indicando el caracter separador (espacios, comas, puntos, punto y comas…) y las columnas de las coordenadas X e Y
uri =
"ruta_archivo.extensión?delimiter={}&xField={}
&yField={}".format("delimitador", "nombre_campoX", "nombre_campoY")
capa = QgsVectorLayer(uri, "layer name you like", "delimitedtext")
Las capas ‘memory’ son capas temporales y vacías que solo existen mientras se ejecuta un script. Son útiles porque nos evita el andar creando archivos innecesarios y actúan como un lienzo en blanco en el que almacenar nuevas geometrías.
En el ejemplo se especifica que la capa temporal sea de tipo multipolígono, su SRC 25830, la creación de un campo id de tipo entero y la creación del índice espacial:
uri = "MultiPolygon?crs=epsg:25830"+"&field=id:integer""&index=yes" capa = iface.addVectorLayer(uri, "capa_temporal", "memory")
Al haber usado iface.addVectorLayer(), la capa temporal se cargará al panel de capas de QGIS y duraría hasta que cerráramos la sesión. Para guardar la capa temporal como capa nueva hay que usar la función QgsVectorFileWriter.writeAsVectorFormat()
Validez de la capa
Para saber si una capa espacial se ha cargado correctamente se debe usar el método .isValid() sobre dicha capa. Lo que hace es devolver un booleano en el que True significa que la capa es correcta y False si no lo es.
Lo que suele hacerse es combinarlo con estructuras de control if – else para obtener un mensaje según si se ha cargado bien o no:
capa = QgsVectorLayer('ruta_capa','nombre_capa','proveedor')
if not capa.isValid():
print('Error al cargar la capa')
else:
print('La capa se cargó correctamente')
Capa activa
Consiste en seleccionar una capa que pasa a ser la ‘capa activa’ para que sea el objetivo de alguna herramienta de QGIS.
La instrucción siguiente crea una variable que contendrá los datos de la capa que se encuentra activa en QGIS. Esto significa que:
- solo podremos utilizar la capa activa si utilizamos Python dentro de QGIS (no funciona en IDEs externos)
- los datos de la variable creada cambiarán según la capa que establezcamos como activa
capa = iface.activelayer()
Vale tanto para capas vectoriales como ráster y solo puede haber una capa activa a la vez.
Copiar capas
Para copiar capas tendremos que cargar la capa que queramos copiar y usarla como valor de entrada en la función de escritura de capas vectoriales:
capa = QgsVectorLayer('ruta_capa','nombre_capa','proveedor')
CRS = capa.crs()
ruta_salida = 'C:\\ruta\\nombre.extension'
QgsVectorFileWriter.writeAsVectorFormat(capa, ruta_salida, CRS, driverName = 'ESRI Shapefile')
En la documentación oficial hay más información acerca de la escritura de archivos vectoriales.
Las capas espaciales vectoriales se componen de entidades: objetos geométricos que representan fenómenos territoriales. Cada uno de estos objetos ocupa un registro en la tabla de atributos de dicha capa, y entre estos atributos se encuentra la geometría del objeto.
En PyQGIS se distinguen dos clases de variables para estos elementos:
- Variables de tipo QgsFeature que almacenan las entidades de una capa
- Variables de tipo QgsGeometry que almacenan la geometría de las entidades
QgsFeature
Para trabajar con la información vectorial es necesario crear un objeto que indexe las entidades de una capa cargada previamente. Dicho objeto será de tipo QgsFeature, y para obtenerlo se debe usar el método .getFeatures() sobre un objeto QgsVectorLayer.
Sin embargo, para poder acceder a cada una de esas entidades indexadas será necesario iterar sobre ellas:
capa_vectorial = QgsVectorLayer('capa.shp','nombre_capa','ogr') ## Objeto Layer
features = capa_vectorial.getFeatures() ## Objeto Feature
Con esto ya podríamos iterar sobre las entidades de la capa usando estructuras de control para obtener
for entidad in features:
print(entidad.id()) ## ver el identificador
print(entidad.attributes()) ## ver los valores de todos los campos
print(entidad[0]) ## ver el valor de la primera columna usando el índice de posición
print(entidad['nombre']) ## ver el valor de la columna usando su nombre
print(entidad.geometry()) ## ver la geometría
QgsGeometry
Una vez accedemos a las entidades podemos extraer su geometría para obtener información en base a ella como puede ser el área, la longitud o el perímetro de nuestras entidades.
Para ello es necesario crear una nueva variable de tipo QgsGeometry que almacene la geometría recogida por el método .geometry() sobre cada entidad y aplicarle un método compatible con el tipo de geometría:
capa_vectorial = QgsVectorLayer('capa.shp','nombre_capa','ogr')
features = capa_vectorial.getFeatures()
for entidad in features:
geom = entidad.geometry() ## Creación del objeto con la geometría
print(geom.area()) ## ver el área de una geometría de tipo polígono
print(geom.length()) ## ver la longitud de una geometría de tipo línea o el perímetro si es de tipo polígono
print(geom.asPoint()) ## ver las coordenadas de una geometría de tipo punto
Además, la geometría es también el producto de los geoprocesos, que nos permitirán generar geometrías nuevas a partir de las existentes.
Por atributos
La selección de entidades según sus atributos o consultas temáticas consiste en seleccionar aquellas entidades cuyos atributos coincidan con nuestra consulta.
Para ello se crean expresiones SQL con la función QgsExpression() para usarlas en una petición mediante la función QgsFeatureRequest(). Dicha petición se emplea para acceder a la geometría de aquellas entidades cuyos atributos casen con la consulta.
En el siguiente ejemplo se seleccionan de una capa con ríos aquellos de longitud superior a 100 km para luego iterar sobre ellos:
capa_vectorial = QgsVectorLayer('rios.shp','red_fluvial','ogr')
expresion = QgsExpression('longitud_km > 100')
peticion = QgsFeatureRequest(expresion)
for entidad in capa_vectorial.getFeatures(peticion):
# operaciones a llevar a cabo con cada río
Otra forma es seleccionar entidades según su id usando el método .selectByIds(). Esto resulta útil pues
capa_rios = QgsVectorLayer('rios.shp','red_fluvial','ogr')
expresion = QgsExpression('longitud_km > 100')
peticion = QgsFeatureRequest(expresion)
entidades = capa_rios.getFeatures(peticion)
lista_id = [rio.id() for rio in entidades]
capa_rios.selectByIds(lista_id)
Por localización
La selección por localización o consultas espaciales consisten en seleccionar aquellas entidades que cumplen con algún tipo de relación topológica, es decir, necesitaremos hacer uso de la geometría de las entidades para llevarlas a cabo.
En PyQGIS podemos hacerlo de dos maneras:
Esto es, utilizar una AOI (Area de Interés) para seleccionar todo aquello que se encuentre dentro de ella. La AOI podemos obtenerla de varias maneras, como tomando el ‘extent’ o recuadro delimitador de una capa (sus coordenadas X e Y máximas y mínimas) o dibujándola nosotros mismos.
En el siguiente ejemplo se utiliza como AOI el ‘extent’ de una capa para que se seleccionen las entidades de otra capa:
# Cargar la capa de la que seleccionar las entidades
capa = QgsVectorLayer('capa1.shp','Mi capa','ogr')
AOI = QgsVectorLayer('capa2.shp','Area de interes','ogr').extent()
peticion = QgsFeatureRequest()
peticion.setFilterRect(AOI)
for entidad in capa.getFeatures(peticion):
print(feature.id())
En este ejemplo la petición estaba vacía, pero se puede combinar filtros espaciales con consultas de atributos para obtener solo algunas entidades de entre todas las que se encuentren en la AOI. En el siguiente ejemplo se obtiene solo entidades de uso de suelo urbano dentro del filtro espacial:
# Cargar la capa de la que seleccionar las entidades
capa = QgsVectorLayer('capa1.shp','Mi capa','ogr')
AOI = QgsVectorLayer('capa2.shp','Area de interes','ogr').extent()
expresion = QgsExpression('USO ILIKE '%urbano%'')
usos_urbanos = QgsFeatureRequest(expresion)
usos_urbanos.setFilterRect(AOI)
for entidad in capa.getFeatures(usos_urbanos):
print(feature.id())
Consiste en quedarnos con aquellas entidades de una capa que cumplen algún tipo de relación espacial con las entidades de otra capa: se encuentran dentro, fuera, tocando, ocupando el mismo espacio…
Para ello tenemos que acceder a la geometría de las entidades de las capas y usar métodos de la clase geometry sobre dichas geometrías para comprobar las relaciones espaciales:
# Cargar la capa de la que seleccionar las entidades
capa1 = QgsVectorLayer("capa1.shp",'Mi capa','ogr')
capa2 = QgsVectorLayer("capa2.shp",'Mi capa 2','ogr')
for f1 in capa1.getFeatures(): # Bucle que itera sobre la geometría de las entidades de la capa 1
for f2 in capa2.getFeatures(): # Bucle para iterar sobre la geometría decada entidad de la capa 2 por cada una de la capa 1
if f1.geometry().intersects(f2.geometry()): # Estructura if para comprobar si hay intersección entre las entidades
print(f1.id(),'intersects',f2.id()) # Acción a realizar si se cumple dicha relación
En el ejemplo anterior se ha usado el método para averiguar si existe intersección, pero hay muchos más predicados espaciales como pueden ser:
| geometría1.intersects(geometría2) | averiguar si existe intersección entre las geometrías |
| geometría1.disjoints(geometría2) | averiguar si no existe intersección |
| geometría1.contains(geometría2) | averiguar si la geometría 1 contiene a la geometría 2 |
| geometría1.crosses(geometría2) | averiguar si la geometría 1 cruza en algún punto con la geometría 2 |
Recomiendo el post sobre relaciones espaciales para conocerlas más a fondo.
Guardar en capa nueva
Para crear capas vectoriales nuevas hay que usar la función QgsVectorFileWriter.writeAsVectorFormat() indicando dentro de los paréntesis los siguientes parámetros, por orden:
| Datos | Una variable u objeto de tipo QgsVectorLayer, es decir, la variable que contiene los datos que queremos guardar en la capa nueva |
| Ruta + nombre + extensión | La ruta junto con el nombre del archivo a guardar y su extensión |
| Codificación | La codificación de caracteres para esta capa: utf-8, latin1, iso… |
| SRC | Sistema de coordenadas asociado a la capa nueva |
| Formato | Tipo de capa o archivo que estamos creando: shapefile, geojson, kml… |
| Modo | Especificar si en el parámetro datos estamos usando entidades seleccionadas (True/False) |
Por ejemplo, con el siguiente código se guarda en una capa nueva una selección de usos del suelo que se encuentran dentro del ‘extent’ de otra capa.
from qgis.core import * from qgis.utils import *ruta = 'C:\ruta\'
capa1 = QgsVectorLayer(ruta+'usos.shp','Mi capa','ogr') capa2 = QgsVectorLayer(ruta+'aoi.shp','Area de interes','ogr').extent() expresion = QgsExpression('USO ILIKE '%urbano%') usos_urbanos = QgsFeatureRequest(expresion) usos_urbanos.setFilterRect(capa2) seleccion = capa.getFeatures(usos_urbanos)
ids = [i.id() for i in seleccion] capa1.selectByIds(ids)
QgsVectorFileWriter.writeAsVectorFormat(capa1, ruta+capa.name()+'_nuevo.shp', 'UTF-8', capa1.crs(), driverName="ESRI Shapefile", onlySelected = True)
Añadir campos
Los atributos de las capas vectoriales se organizan en forma de tabla, en la que cada columna es un campo y cada fila una entidad. Por tanto, los valores de los campos para cada fila serán los atributos de dicha fila.
Se pueden crear campos al crear capas temporales. El siguiente ejemplo crea un campo de cada tipo (text, integer, double y date) en una nueva capa temporal (Fuente: Geoinnova):
URI = "MultiPolygon?field=id:integer&index=yes&field=decimal:double&field=tex to:string&field=fecha:date" capa_temporal = iface.addVectorLayer(URI,"temp", "memory")
También pueden crearse campos en capas ya existentes activando la edición de dicha capa y encadenando los siguientes métodos de la clase QgsVectorLayer:
capa.dataProvider().addAttributes([QgsField])
- método .dataProvider() ## para modificar tablas de atributos
- método .addAttributes([lista_objetos_QgsField]) ## para añadir nuevos campos
Los objetos de tipo campo o QgsField se crean con la siguiente sintaxis:
QgsField("nombre_campo", QVariant.Text) ## cadenas de texto
QgsField("nombre_campo", QVariant.Integer) ## numeros enteros
QgsField("nombre_campo", QVariant.Double) ## numeros decimales
QgsField("nombre_campo", QVariant.Date) ## valores de tipo fecha
Al completo, el procedimiento sería el siguiente:
## Cargar la capa
capa = iface.addVectorLayer('ruta_capa','nombre_capa','proveedor)
capa.startEditing()
capa.dataProvider().addAttributes([QgsField("texto", QVariant.String),
QgsField("entero", QVariant.Int), QgsField("decimal",
QVariant.Double),QgsField("fecha", QVariant.Date)])
capa.updateFields()
capa.commitChanges()
Para ello habrá que obtener un objeto de la clase QgsField con los campos de la capa de la que queremos copiar dichos campos. Después solo habrá que crear una capa nueva, activar su edición y copiarle dichos campos:
## Cargar capa
capa = iface.addVectorLayer('ruta_capa','nombre_capa','proveedor)
crs = capa.crs().postgisSrid()
campos = capa.dataProvider().fields() ## Obtención de campos
URI = "MultiPolygon?crs=epsg:"+str(crs)
nueva_capa = iface.addVectorLayer('URI','nueva','memory')
nueva_capa.startEditing()
nueva_capa.dataProvider().addAttributes(campos.toList())
nueva_capa.commitChanges()
Una de las operaciones más comunes en cuanto a la creación de nuevos campos es crear IDs que identifiquen de forma única a cada entidad:
## Cargar la capa y empezar la edición capa = QgsVectorLayer(output,nombre,'ogr') capa.startEditing()capa.dataProvider().addAttributes([QgsField('ID', QVariant.Int)]) capa.updateFields()
entidades = capa.getFeatures()
contador = 0
for entidad in entidades:
## Establecer el valor del atributo ID para la entidad a partir del valor del contador entidad['ID'] = contador ## Actualizar la entidad capa.updateFeature(entidad) ## Aumentar el contador contador += 1capa.commitChanges()
Modificar atributos
Con el método .setAttributes([]) podemos actuar sobre un objeto de tipo entidad en el momento de crearla, es decir, cuando generamos nuevos objetos de la clase QgsFeature:
## cargar capa
capa = iface.addVectorLayer('ruta_capa','nombre_capa','proveedor)
lista_valores = [30, 40, 60]
for f in range(1,10):
## nueva entidad
entidad = QgsFeature()
## asignar los valores a las columnas
entidad.setAttributes(lista_valores)
## añadir la nueva entidad a la capa
capa.addFeatures([entidad])
capa.commitChanges()
La lista de valores para la nueva entidad sigue el orden de posición de los campos. Si se quisiera cambiar, por ejemplo, la tercera columna, habría que dejar vacías las posiciones 0 y 1 de la lista.
Mediante .changeAttributeValue() actuamos sobre un objeto de la clase QgsVectorLayer, es decir, sobre una capa indicando el id de la fila y la posición de la columna en la que actuar:
capa = iface.addVectorLayer('ruta_capa','nombre_capa','proveedor)
capa.startEditing()
nuevo_valor = 50
id_entidad = 3
id_campo = 0
capa.changeAttributeValue(id_entidad, id_campo, nuevo_valor)
capa.commitChanges()
Editar geometrías
La edición de geometrías implica operar sobre los objetos QgsGeometry, es decir, la geometría de las entidades, para borrarlas, modificarlas o añadir nuevas geometrías que se deriven de algún geoproceso sobre ellas.
Por ejemplo, al hacer un buffer tomamos las geometrías de una capa y se les aplica un geoproceso que genera nuevos polígonos en base a un valor de distancia que indicamos en la función. Esos nuevos polígonos deben almacenarse en otra capa que se encuentre en edición.
Para comenzar a editar las geometrías es necesario activar la edición de la capa de entrada, y al finalizar la edición deben guardarse los cambios:
capa = QgsVectorLayer('capa1.shp','Mi capa','ogr') # Objeto tipo layer con los datos de la capa
capa.startEditing() # Activación de la edición de la capa
## Geoprocesos ##
capa.commitChanges() # Guardar los cambios
Los principales métodos para editar la geometría de una capa son:
objeto_qgsfeature.setGeometry(nueva_geometría) ## Establecer la geometría de una entidad objeto_qgsvectorlayer.addFeatures([objeto_qgsfeature]) ## Añadir la entidad a una capa
A diferencia de las selecciones, los resultados de los geoprocesos deben almacenarse en capas nuevas que se encuentren en edición, ya que estamos hablando de geometrías nuevas que deben almacenarse en algún lado.
A continuación tenéis ejemplos con distintos geoprocesos básicos:
El buffer o área de influencia se calcula tomando una geometría y aplicándole el método .buffer() junto a una serie de parámetros (ver la documentación oficial).
El procedimiento para llevarlo a cabo consiste en iterar sobre las geometrías de las entidades de la capa, en este caso para que se aplique el buffer a cada una de ellas y la geometría resultante se almacene en un objeto QgsFeature que pase a formar parte de una capa temporal. Además, se crea un contador para darle un ID a cada buffer que se cree:
## Cargar capa vectorial de entrada y obtener su sistema de coordenadas ruta = 'C:\\ruta\\' capa_entrada = QgsVectorLayer(ruta, "entidades_de_entrada", "ogr") CRS = capa_entrada.crs().postgisSrid()uri = "MultiPolygon?crs=epsg:"+ str(CRS) + "&field=id:integer""&index=yes" capa_buffers = QgsVectorLayer(uri, "buffers", "memory")
capa_buffers.startEditing() ## Comenzar la edición de la capa temporal id = 0 ## Valor id que se le asignará a cada nueva entidad
for f in capa_entrada.getFeatures(): ## Bucle para obtener las entidades de la capa de entrada geom = f.geometry() ## Acceso a la geometría de cada entidad buffer = geom.buffer(200,-1) ## Parámetros del buffer: 200 metros, nº segmentos automático entidad = QgsFeature() ## Creación de objeto de tipo Feature vacío entidad.setGeometry(buffer) ## Añadirle al objeto la geometría del buffer entidad.setAttributes([id]) ## Añadirle al objeto el valor de la variable id capa_buffers.addFeatures([entidad]) ## Añadirle a la capa temporal la entidad id += 1 ## Aumentar el id al acabar cada paso del bucle
capa_buffers.commitChanges() ## Guardar los cambios de la capa temporal
En el caso de la intersección, el método a aplicar sobre los objetos QgsGeometry es .intersection() (ver la documentación oficial)
No debe confundirse con el método .intersect(), porque a diferencia de éste ahora estaremos generando geometrías nuevas, es decir, objetos de la clase QgsGeometry.
Con el siguiente código se crea una capa nueva con la geometría de las intersecciones entre las entidades de dos capas de entrada, es decir, con la superficie que ambas capas tienen en común (por tanto deberán ser dos capas de tipo poligonal):
## Módulos from qgis.core import * from qgis.utils import *ruta1 = 'C:\ruta' ruta2 = 'C:\ruta2' capa_entrada1 = QgsVectorLayer(ruta1, "entidades_de_entrada", "ogr") capa_entrada2 = QgsVectorLayer(ruta2, "entidades_de_entrada2", "ogr") CRS = capa_entrada1.crs().postgisSrid()
uri = "MultiPolygon?crs=epsg:"+ str(CRS) + "&field=id:integer""&index=yes" capa_intersecciones = iface.addVectorLayer(uri, "intersecciones", "memory")
capa_intersecciones.startEditing() ## Comenzar la edición de la capa temporal id = 0 ## Valor id que se le asignará a cada nueva entidad
for f1 in capa_entrada1.getFeatures(): ## Bucle para obtener las entidades de la capa de entrada 1 for f2 in capa_entrada2.getFeatures(): ## Bucle para obtener las entidades de la capa de entrada 2 if f1.geometry().intersects(f2.geometry()): interseccion_geom = f1.geometry().intersection(f2.geometry()) ## Creación de la geometría de la intersección entidad = QgsFeature() ## Creación de objeto de tipo Feature vacío entidad.setGeometry(interseccion_geom) ## Añadirle al objeto la geometría de la intersección entidad.setAttributes([id]) ## Añadirle al objeto el valor de la variable id capa_intersecciones.addFeatures([entidad]) ## Añadirle a la capa temporal la entidad id += 1 ## Aumentar el id al acabar cada paso del bucle
capa_intersecciones.commitChanges() ## Guardar los cambios de la capa temporal
Con la diferencia obtendríamos la superficie de cada entidad de una primera capa que no coincide espacialmente con ninguna entidad de la segunda capa. El procedimiento es casi igual que con la intersección (necesario que las capas sean polígonos):
from qgis.core import * from qgis.utils import *ruta1 = 'C:\ruta1' ruta2 = 'C:\ruta2' capa_entrada1 = QgsVectorLayer(ruta1, "entidades_de_entrada", "ogr") capa_entrada2 = QgsVectorLayer(ruta2, "entidades_de_entrada2", "ogr") CRS = capa_entrada1.crs().postgisSrid()
uri = "MultiPolygon?crs=epsg:"+ str(CRS) + "&field=id:integer""&index=yes" capa_diferencia = iface.addVectorLayer(uri, "diferencias", "memory")
capa_diferencia.startEditing()
id = 0for f1 in capa_entrada1.getFeatures(): diff_geom = f1.geometry() for f2 in capa_entrada2.getFeatures(): if diff_geom.intersects(f2.geometry()): diff_geom = diff_geom.difference(f2.geometry()) entidad = QgsFeature() entidad.setGeometry(diff_geom) entidad.setAttributes([id]) capa_diferencia.addFeatures([entidad]) id += 1 capa_diferencia.commitChanges()
En este caso, lo que se hace es crear un objeto con la geometría de cada entidad de la primera capa y modificarlo o ‘recortarlo’ con las geometrías de las entidades de la segunda capa. Después, la geometría resultante es almacenada en la capa temporal junto a su id.