Găsirea Poligoanelor Învecinate Dintr-un Strat

Notă

Dacă utilizați QGIS3, există o modalitate mai simplă și mai bună de a face această analiză folosind funcția de agregare. Parcurgeți, în acest sens, următoarea postare Identificarea poligoanelor învecinate folosind funcția de agregare din QGIS

Există unele cazuri de utilizare, când se dorește găsirea tuturor poligoanelor învecinate, pentru fiecare dintre poligoanele dintr-un strat. Cu un mic script în Python, putem realiza acest lucru și multe altele, în QGIS. Aici este un exemplu de script, pe care îl puteți folosi pentru a găsi toate poligoanele care au granița comună, pentru fiecare dintre poligoanele stratului și, de asemenea, să le adăugați numele în tabela de atribute. Ca bonus, script-ul însumează un atribut, la alegere, din toate poligoanele învecinate.

Privire de ansamblu asupra activității

Pentru a demonstra modul în care funcționează script-ul, vom folosi un strat poligonal de țări și vom găsi țările care împart aceeași frontieră. De asemenea, vom calcula totalul populației din țările învecinate.

Obținerea datelor

Vom folosi setul de date Admin 0 - Țări de la Natural Earth.

Vom folosi setul de date Admin 0 - fișierul shape al țărilor.

Sursa de date [NATURALEARTH]

Obținerea script-ului

Descărcați script-ul neighbors.py și salvați-l pe discul dumneavoastră.

Procedura

  1. Încărcați stratul ne_10m_admin_0_countries mergând la Layer ‣ Add Vector Layer.

../_images/1178.png
  1. Script-ul folosește 2 câmpuri pentru a efectua acțiunea. Un câmp pentru nume, și un câmp pe care doriți să-l însumați. Utilizați instrumentul Identify pentru a face clic pe orice element, pentru a-i examina atributele. În acest caz, câmpul pentru nume este NAME, și dorim să însumăm estimările de populație din câmpul POP_EST.

../_images/2147.png
  1. Mergeți la Plugins ‣ Python Console.

../_images/388.png
  1. În fereastra Python Console, faceți clic pe butonul Show Editor.

../_images/454.png
  1. În panoul Editor, faceți clic pe butonul Open file, navigați la script-ul descărcat neighbors.py și apăsați Open.

../_images/547.png
  1. O dată ce script-ul este încărcat, este posibil să doriți să schimbați valorile _NAME_FIELD și _SUM_FIELD pentru a se potrivi cu atributele din propriul strat. Dacă lucrați cu stratul ne_10m_admin_0_countries, le puteți lăsa așa cum sunt. Faceți clic pe butonul Save din panoul Editor, dacă ați făcut vreo schimbare. Apoi apăsați butonul Run script pentru a rula script-ul.

../_images/645.png
  1. O dată ce script-ul se încheie, faceți clic-dreapta pe stratul ne_10m_admin_0_countries, apoi selectați Open Attribute Table.

../_images/745.png
  1. Veți observa 2 noi atribute, numite NEIGHBORS și SUM. Acestea au fost adăugate de către script.

../_images/844.png

Pentru referință, mai jos este script-ul complet. Îl puteți modifica pentru a se potrivi cerințelor dumneavoastră.

################################################################################
# Copyright 2014 Ujaval Gandhi
#
#This program is free software; you can redistribute it and/or
#modify it under the terms of the GNU General Public License
#as published by the Free Software Foundation; either version 2
#of the License, or (at your option) any later version.
#
#This program is distributed in the hope that it will be useful,
#but WITHOUT ANY WARRANTY; without even the implied warranty of
#MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#GNU General Public License for more details.
#
#You should have received a copy of the GNU General Public License
#along with this program; if not, write to the Free Software
#Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
#
################################################################################
from qgis.utils import iface
from PyQt4.QtCore import QVariant

# Replace the values below with values from your layer.
# For example, if your identifier field is called 'XYZ', then change the line
# below to _NAME_FIELD = 'XYZ'
_NAME_FIELD = 'NAME'
# Replace the value below with the field name that you want to sum up.
# For example, if the # field that you want to sum up is called 'VALUES', then
# change the line below to _SUM_FIELD = 'VALUES'
_SUM_FIELD = 'POP_EST'

# Names of the new fields to be added to the layer
_NEW_NEIGHBORS_FIELD = 'NEIGHBORS'
_NEW_SUM_FIELD = 'SUM'

layer = iface.activeLayer()

# Create 2 new fields in the layer that will hold the list of neighbors and sum
# of the chosen field.
layer.startEditing()
layer.dataProvider().addAttributes(
        [QgsField(_NEW_NEIGHBORS_FIELD, QVariant.String),
         QgsField(_NEW_SUM_FIELD, QVariant.Int)])
layer.updateFields()
# Create a dictionary of all features
feature_dict = {f.id(): f for f in layer.getFeatures()}

# Build a spatial index
index = QgsSpatialIndex()
for f in feature_dict.values():
    index.insertFeature(f)

# Loop through all features and find features that touch each feature
for f in feature_dict.values():
    print 'Working on %s' % f[_NAME_FIELD]
    geom = f.geometry()
    # Find all features that intersect the bounding box of the current feature.
    # We use spatial index to find the features intersecting the bounding box
    # of the current feature. This will narrow down the features that we need
    # to check neighboring features.
    intersecting_ids = index.intersects(geom.boundingBox())
    # Initalize neighbors list and sum
    neighbors = []
    neighbors_sum = 0
    for intersecting_id in intersecting_ids:
        # Look up the feature from the dictionary
        intersecting_f = feature_dict[intersecting_id]

        # For our purpose we consider a feature as 'neighbor' if it touches or
        # intersects a feature. We use the 'disjoint' predicate to satisfy
        # these conditions. So if a feature is not disjoint, it is a neighbor.
        if (f != intersecting_f and
            not intersecting_f.geometry().disjoint(geom)):
            neighbors.append(intersecting_f[_NAME_FIELD])
            neighbors_sum += intersecting_f[_SUM_FIELD]
    f[_NEW_NEIGHBORS_FIELD] = ','.join(neighbors)
    f[_NEW_SUM_FIELD] = neighbors_sum
    # Update the layer with new attribute values.
    layer.updateFeature(f)

layer.commitChanges()
print 'Processing complete.'

If you want to give feedback or share your experience with this tutorial, please comment below. (requires GitHub account)