¿Determinar si un polígono se encuentra completamente dentro de otro polígono o sólo en parte utilizando QGIS?

Question

Yo utilizo principalmente QGIS.

Tengo dos capas que contienen polígonos: una capa (A) contiene sólo un polígono y la otra contiene varios polígonos (capa (B)).

Puedo ejecutar la herramienta de intersección y QGIS me devuelve los polígonos de la capa (B) que se encuentran total o parcialmente dentro del polígono de la capa (A).

Pero, ¿es posible que QGIS devuelva una clasificación de los polígonos de intersección de la capa (B) por "completamente" dentro del polígono de la capa (A) y "parcialmente" dentro del polígono de la capa (A).

Tal vez esto podría hacerse calculando y comparando el área de los polígonos de la capa(B) antes y después de la intersección, pero me gustaría evitar eso por razones.

¿Y cómo puedo automatizar esto con un script?

Answer 1

En primer lugar, puede seleccionar los polígonos de B que se encuentran completamente dentro de la capa de polígonos A utilizando la función st_contiene función. Clasifíquelos como completamente dentro del polígono A.

A continuación, puede seleccionar todos los polígonos de B que intersecan a A utilizando la función st_intersects pero esto también seleccionará los polígonos que están completamente dentro de A, por lo que tendrá que eliminar los resultados de la primera consulta/clasificación.

En lo que respecta a QGIS, puede utilizar el plugin Spatial Query para hacerlo,

Answer 2

Para automatizar con un script puede utilizar el siguiente código:

mapcanvas = iface.mapCanvas()

layers = mapcanvas.layers()

feats_A = [ feat for feat in layers[0].getFeatures() ]

feats_B = layers[1].getFeatures().next()

for i, feat in enumerate(feats_A):
    if feat.geometry().intersects(feats_B.geometry()):
        if feats_B.geometry().contains(feat.geometry()):
            print "total overlapping, ", "feature: ", i
        else:
            areaA = feat.geometry().intersection(feats_B.geometry()).area()
            areaB = feats.geometry().area()
            print "partial overlapping, feature: %d, perc: %.2f" % (i, (areaA/areaB)*100)

Lo pruebo con dos shapefiles para una condición similar a la de tu pregunta:

Después de ejecutar el código en la consola Python de QGIS, se puede observar en la imagen de arriba, que se imprimieron sólo los tres casos posibles y sus respectivos ids y porcentajes (excepto la superposición total: 100 %)

Answer 3

Si no, puedes hacerlo con la calculadora de campo de la capa B con la expresión:

within( $geometry,geometry( get_feature('LayerA, 'id', 1))) 

Siempre que tu capa A se llame 'CapaA' y que tenga un campo 'id' con valor '1'. Por supuesto, puede modificar todos estos parámetros.

La estructura de get_feature es la siguiente: get_feature(capa, atributo, valor).

La expresión devolverá 1 para los polígonos completamente dentro y 0 para los demás.