Determinar si shapefile y de trama se superponen en python

Question

Estoy construyendo un script en python usando osgeo/gdal.

Tengo un conjunto de shapefiles y un conjunto de GeoTiff archivos ráster.

Me gustaría tener mi script ignorar los shapefiles si no se intersecan con la trama de la zona.

El shapefile no es un rectángulo, así que no puedo comparar simplemente el xmin/xmax,ymin/ymax valores devueltos por capa.GetExtent(). Necesito el real polígono que representa la forma general y, a continuación, algunos de forma de determinar si ese polígono se cruza con la trama de la plaza.

Yo estaba pensando que de alguna manera podría combinar todos los polígonos en el shapefile en una función y, a continuación, lea la geometría en esa función, y luego comparar esa información a la trama medida. Sin embargo, estoy seguro de concreto para la ejecución de este.

Mi pregunta se reduce a 2 puntos:

1. Cómo extraer la frontera del polígono información de shapefile?
2. Cómo determinar si el polígono se cruza una determinada zona de la plaza?

Cualquier visión se agradece.

Answer 1

La siguiente secuencia de comandos determina el cuadro delimitador de trama y crea basándose en el cuadro delimitador de una geometría.

import ogr, gdal

raster = gdal.Open('sample.tif')
vector = ogr.Open('sample.shp')

# Get raster geometry
transform = raster.GetGeoTransform()
pixelWidth = transform[1]
pixelHeight = transform[5]
cols = raster.RasterXSize
rows = raster.RasterYSize

xLeft = transform[0]
yTop = transform[3]
xRight = xLeft+cols*pixelWidth
yBottom = yTop-rows*pixelHeight

ring = ogr.Geometry(ogr.wkbLinearRing)
ring.AddPoint(xLeft, yTop)
ring.AddPoint(xLeft, yBottom)
ring.AddPoint(xRight, yTop)
ring.AddPoint(xRight, yBottom)
ring.AddPoint(xLeft, yTop)
rasterGeometry = ogr.Geometry(ogr.wkbPolygon)
rasterGeometry.AddGeometry(ring)

Siguiente, la geometría de las vectorial de polígono se determina. Esto responde a su primera pregunta.

# Get vector geometry
layer = vector.GetLayer()
feature = layer.GetFeature(0)
vectorGeometry = feature.GetGeometryRef()

Por último, la geometría del vector y raster son probados para la intersección (devuelve True o False). Esto responde a tu segunda pregunta.

print rasterGeometry.Intersect(vectorGeometry)