Python OGR: convertir una geometría de puntos en un polígono

Question

Estoy escribiendo una rutina para crear la envolvente del casco convexo alrededor de shapefiles de islas. Podría leer los datos, extraer los vértices de los polígonos en una geometría y luego calcular el casco convexo.

Ahora tengo que guardarlo en un shapefile de polígonos.

¿Hay alguna manera de transformar la geometría (llamada convexHull en mi código) devuelto por el método convexHull en una geometría de polígono; o estoy obligado a analizar todos los vértices de convexHull y crear un anillo de un polígono?

def doConvexHull(infile, outfile):
    inH = ogr.Open(infile, 0)
    if inH is None:
        usage("Could not open file {0}. Exit.".format(infile))
    layer = inH.GetLayer()

    # get all polygons
    thisGeometry = ogr.Geometry(ogr.wkbPoint)
    for index in xrange(layer.GetFeatureCount()):
        feature = layer.GetFeature(index)
        geometry = feature.GetGeometryRef()
        ring = geometry.GetGeometryRef(0)
        points = ring.GetPointCount()
        for p in xrange(points):
            lon, lat, z = ring.GetPoint(p)
            thisGeometry.AddPoint(lon, lat)

    convexHull = thisGeometry.ConvexHull()

    drv = ogr.GetDriverByName( "ESRI Shapefile" )
    ds = drv.CreateDataSource( outfile )
    if ds is None:
        usage("Could not create file {0}".format( outfile) )

    lyrname = "convexHull_${0}".format( layer.GetName() )
    lyr = ds.CreateLayer( lyrname, layer.GetSpatialRef(), ogr.wkbPolygon )
    thisFeature = ogr.Feature( layer.GetLayerDefn() )
    thisFeature.SetGeometry( convexHull )
    lyr.CreateFeature( thisFeature ): error from geometry

Answer 1

En este caso lo que ocurre es que estás creando una geometría de tipo Punto. Las geometrías de tipo punto sólo pueden contener un punto real.

thisGeometry.GetGeometryCount() #1

ConvexHull() crea una geometría del mínimo número de dimensiones que contiene toda la geometría del objeto padre. Así, si hay un punto, el casco convexo será un punto, para dos puntos será una línea (asumiendo que no hay puntos coincidentes), y para tres o más será un polígono (asumiendo que los puntos no están en una línea y no son coincidentes).

Para comprobarlo, lo más sencillo es utilizar la función ogr.GeometryTypeToName método:

ogr.GeometryTypeToName(convexHull.GetGeometryType())

La forma más fácil de obtener el casco convexo de un grupo de islas sería crear una GeometryCollection y obtener el casco convexo a partir de ella, así:

geom_collection = ogr.Geometry(ogr.wkbGeometryCollection)
for feature in layer: #layers (and datasets) have a list interface
    geometry = feature.GetGeometryRef()
    geom_collection.AddGeometry(geometry)

convex_hull = geom_collection.ConvexHull()

Y este casco convexo debería ser un polígono (ya que se crea a partir de una colección de polígonos), lo que significa que debería poder utilizarlo como entrada para su conjunto de datos de salida.

Si aún así quisieras obtener el casco convexo de múltiples puntos tendrías que crear una geometría de tipo MultiPoint y luego añadir los puntos a ese objeto MultiPoint individualmente.

thisGeometry = ogr.Geometry(ogr.wkbMultiPoint)
for feature in layer:
    geometry = feature.GetGeometryRef()
    for ring in geometry: #This is assuming that the geometry is a polygon
        for point in ring:
            thisGeometry.AddGeometry(point)

convexHull = thisGeometry.ConvexHull()

Answer 2

El método OGR convexHull devuelve un polígono, así que la respuesta a tu pregunta es "no". No es necesario analizar los vértices ni convertir la geometría.

Me doy cuenta de que está definiendo thisGeometry sea del tipo ogr.wkbPoint . la geometría wkbPoint es para almacenar puntos individuales y no múltiples puntos o polígonos. Tampoco es necesario analizar cada punto individual, sino que puede añadir geometrías enteras a su colección al por mayor. Para tu caso de uso te sugiero que pruebes:

thisGeometry = ogr.Geometry(ogr.wkbGeometryCollection)
for index in range(layer.GetFeatureCount()):
    feature = layer.GetFeature(index)
    geometry = feature.GetGeometryRef()
    thisGeometry.AddGeometry(geometry)

convexHull = thisGeometry.ConvexHull()

El resultado de esto debería ser un polígono normal de casco convexo.

Answer 3

Si puedes limitar la capacitancia de carga del cable y del amplificador, puedes utilizar el modelo de tensión. Para ello es necesario que no haya ningún cable, ya que los cables suelen añadir 20~30 pF/pie.

En caso contrario, hay que utilizar el modelo de carga y se necesita un amplificador de carga. Así es como funcionan los acelerómetros con cables del sensor al amplificador.