Geoprocesamiento con Python usando la librería OGR

En el día de hoy me decidí a buscar la manera de incorporar en mis scripts de python las librerías GDAL/OGR y encontré un curso muy interesante de la Universidad de Utah sobre geoprocesamiento con Python y SIG Open Source. El curso está orientado a Windows y Mac pero como yo tengo instaladas las librerías de GDAL/OGR en mi Linux Debian no tuve problemas en probar algunas cosas con Python. Las librerías que pueden ser usadas en los scripts son 5 (gdal, ogr, osr, gdal_array y gdalconst) pero acá sólo me limitaré a un ejemplo sencillo de uso con la librería ogr (Simple Feature Library) y archivos vectoriales.

Para ello, voy a determinar el tipo, la extensión (donde corresponda) que abarca un shapefile (formato de ESRI) y las coordenadas de punto sin necesidad de usar software SIG (ArcGis, QGIS o GRASS); sólo un script muy sencillo en texto plano cargando la librería OGR. Para realizar una tarea como ésta, además de cargar la librería (from osgeo import ogr), el script debe posicionarse en el directorio donde están los shapefiles e invocar los métodos para identificar el driver que corresponde al formato vectorial usado (‘ESRI Shapefile’), abrir el vectorial escogido, obtener la capa correspondiente y, finalmente, invocar los métodos para determinar la extensión y obtener las coordenadas. El script es el siguiente:

#!/usr/bin/env python
# -*- coding: utf-8

from os import chdir,system
from osgeo import ogr

system("clear")

#Se cambia al directorio de los shapefiles
chdir('/home/zeito/Desktop/PRUEBAS_QGIS')

#Establece el driver a usar
driver = ogr.GetDriverByName('ESRI Shapefile')

#Abre archivo tipo polígono
dataSource1 = driver.Open('cuenca.shp', 0)

#Abre archivo tipo punto
dataSource2 = driver.Open('centroide.shp', 0)

#Obtiene la fuente de datos
layer1 = dataSource1.GetLayer()
layer2 = dataSource2.GetLayer()

#Dos modelos de despliegue de la extensión
extent = layer1.GetExtent()

print 'Dos modelos de despliegue de la extensión (Polígono cuenca)\nExtent:', extent
print 'UL:', extent[0], extent[3]
print 'LR:', extent[1], extent[2]

#Obtiene el rasgo
feature1 = layer1.GetFeature(0)
feature2 = layer2.GetFeature(0)

#Obtiene la geometría
geometry1 = feature1.GetGeometryRef()
geometry2 = feature2.GetGeometryRef()

#Imprime geometría y coordenadas de punto como tupla
print "\nGeometría y vértices de polígono\n", geometry1
print "\nGeometría y coordenadas de centroide\n", geometry2

#Individualiza coordenadas
x = geometry2.GetX()
y = geometry2.GetY()
z = geometry2.GetZ()

print "X = ", x, "Y = ", y, "Z = ", z

cuya ejecución produce esta salida:

No obstante, todavía no conozco los métodos para acceder a las coordenadas de los vértices del polígono que despliega geometry1 (para los vectoriales tipo punto es directo y está en el script); lo verdaderamente importante.

Esta entrada fue publicada en Código Python, GDAL/OGR, SIG, Software Libre. Guarda el enlace permanente.

2 respuestas a Geoprocesamiento con Python usando la librería OGR

  1. Ignacio Torres dijo:

    Hola. No se si te llegará este mail directamente Necesito escribir una base de datos, en una planilla excel, para luego guardarla como texte delimitado por coma para poder cargarla en Quantum Gis El problema es que no sé como o de donde sacar las cordenadas de las calles que estoy escribiendo en excel para que Qgis pueda referenciarlas. Es decir yo escribo la dirección con la altura pero la coordenada geográfica como la pongo? Gracias, espero me puedas asesorar

    El 10 de octubre de 2012 09:01, WordPress.com

  2. Pingback: Cómo extraer las coordenadas x, y de un polígono mediante el módulo OGR de Python |

Responder

Por favor, inicia sesión con uno de estos métodos para publicar tu comentario:

Logo de WordPress.com

Estás comentando usando tu cuenta de WordPress.com. Cerrar sesión / Cambiar )

Imagen de Twitter

Estás comentando usando tu cuenta de Twitter. Cerrar sesión / Cambiar )

Foto de Facebook

Estás comentando usando tu cuenta de Facebook. Cerrar sesión / Cambiar )

Google+ photo

Estás comentando usando tu cuenta de Google+. Cerrar sesión / Cambiar )

Conectando a %s