¿Crear un nuevo shapefile de puntos cuando las características de la cuadrícula vectorial tienen un atributo común?

Question

Tengo un shapefile muy grande que se compone de muchas celdas. Si hay 4 celdas juntas que tienen el mismo valor de atributo, debo localizar una geometría de puntos, como la siguiente:

Cada lado de la celda tiene 625 metros de longitud. Cualquier idea es bienvenida (QGIS, ArcGIS, o código en PyQGIS...).

Answer 1

Se podría convertir a raster con un tamaño de celda de 1250 (2 veces el tamaño de celda de la cuadrícula vectorial). A continuación, se volvería a convertir a vector, se seleccionarían las celdas que tienen el atributo deseado y se tomarían sus centroides. A continuación, realizaría una unión espacial entre los centroides y la capa original, tomando los atributos de todas las celdas originales que tocan los centroides para eliminar cualquier centroide erróneo que se haya creado debido a una mayoría de atributos deseables en lugar de un bloque completo de cuatro celdas como se muestra a continuación.

Answer 2

Propongo una solución utilizando QGIS (PyQGIS). Partiendo de esta cuadrícula vectorial, que ha sido categorizada utilizando el mismo color para los valores de los atributos comunes del "MY_FIELD" campo):

puede ejecutar el siguiente código desde el Consola Python :

from qgis.core import *
layer = QgsVectorLayer('C:/Users/path_to_file/filename.shp', 'test', 'ogr') # Edit the filepath
comparing_attr = 'MY_FIELD' # Edit the field name

# Create the output layer
crs = layer.crs().toWkt()
outLayer = QgsVectorLayer('Point?crs='+ crs, 'common' , 'memory')
prov = outLayer.dataProvider()
fields = layer.pendingFields()
for field in fields:
    if field.name() == comparing_attr:
        field_type = field.typeName()
        field_length = field.length()
        new_field = [field]
prov.addAttributes(new_field)
outLayer.updateFields()

attrs_dict = {}
for feat in layer.getFeatures():
    tmp_coords = feat.geometry().asPolygon()
    points = []
    for i in tmp_coords:
        points.extend(i)
    coords=list(set(points))

    if feat[comparing_attr] not in attrs_dict.keys():
        attrs_dict[feat[comparing_attr]] = coords
    else:
        attrs_dict[feat[comparing_attr]].extend(coords)

for key, vals in attrs_dict.iteritems():
    common=[x for x in vals if vals.count(x) >= 4]
    already_processed = []
    for item in common:
        if item not in already_processed:
            geom=QgsGeometry.fromPoint(item)
            outGeom = QgsFeature()
            outGeom.setGeometry(geom)
            outGeom.setAttributes([key])
            prov.addFeatures([outGeom])
            already_processed.append(item)

# Add the layer to the Layers panel
QgsMapLayerRegistry.instance().addMapLayer(outLayer)

y devolverá una capa de memoria puntual que almacena el valor del atributo común:

Lo único que tienes que hacer es adaptar adecuadamente el código a tu caso concreto, es decir, tienes que especificar la ruta del archivo de la capa y el nombre del atributo de interés.