-
-
Notifications
You must be signed in to change notification settings - Fork 3k
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Ported 'Join by location' from fTools to Processing
- Loading branch information
1 parent
f9c2b91
commit 8daf8c3
Showing
3 changed files
with
225 additions
and
1 deletion.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,209 @@ | ||
# -*- coding: utf-8 -*- | ||
|
||
""" | ||
*************************************************************************** | ||
SpatialJoin.py | ||
--------------------- | ||
Date : October 2013 | ||
Copyright : (C) 2013 by Joshua Arnott | ||
Email : josh at snorfalorpagus dot net | ||
*************************************************************************** | ||
* * | ||
* 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. * | ||
* * | ||
*************************************************************************** | ||
""" | ||
from processing.core.GeoAlgorithmExecutionException import GeoAlgorithmExecutionException | ||
__author__ = 'Joshua Arnott' | ||
__date__ = 'October 2013' | ||
__copyright__ = '(C) 2013, Joshua Arnott' | ||
# This will get replaced with a git SHA1 when you do a git archive | ||
__revision__ = '$Format:%H$' | ||
|
||
from processing.core.GeoAlgorithm import GeoAlgorithm | ||
from PyQt4.QtCore import * | ||
from PyQt4.QtGui import * | ||
from qgis.core import * | ||
from processing.parameters.ParameterVector import ParameterVector | ||
from processing.parameters.ParameterSelection import ParameterSelection | ||
from processing.parameters.ParameterString import ParameterString | ||
from processing.tools import dataobjects | ||
from processing.outputs.OutputVector import OutputVector | ||
from processing.tools import vector | ||
from processing.core.ProcessingLog import ProcessingLog | ||
import os | ||
|
||
def myself(L): | ||
#median computation | ||
nVal = len(L) | ||
if nVal == 1: | ||
return L[0] | ||
L.sort() | ||
#test for list length | ||
medianVal = 0 | ||
if nVal > 1: | ||
if ( nVal % 2 ) == 0: | ||
#index begin at 0 | ||
#remove 1 to index in standard median computation | ||
medianVal = 0.5 * ( (L[ (nVal) / 2 - 1]) + (L[ (nVal) / 2 ] )) | ||
else: | ||
medianVal = L[ (nVal + 1) / 2 - 1] | ||
return medianVal | ||
|
||
class SpatialJoin(GeoAlgorithm): | ||
''' | ||
Join attributes by location | ||
Port of the spatial join algorithm from fTools to the Processing Toolbox. | ||
''' | ||
|
||
INPUT1 = "INPUT1" | ||
INPUT2 = "INPUT2" | ||
SUMMARY = "SUMMARY" | ||
STATS = "STATS" | ||
KEEP = "KEEP" | ||
OUTPUT = "OUTPUT" | ||
|
||
SUMMARYS = [ | ||
'Take attributes of the first located feature', | ||
'Take summary of intersecting features' | ||
] | ||
|
||
KEEPS = [ | ||
'Only keep matching records', | ||
'Keep all records (including non-matching target records)' | ||
] | ||
|
||
#=========================================================================== | ||
# def getIcon(self): | ||
# return QIcon(os.path.dirname(__file__) + "/icons/join_location.png") | ||
#=========================================================================== | ||
|
||
def defineCharacteristics(self): | ||
self.name = "Join attributes by location" | ||
self.group = "Vector general tools" | ||
self.addParameter(ParameterVector(SpatialJoin.INPUT1, "Target vector layer", [ParameterVector.VECTOR_TYPE_ANY])) | ||
self.addParameter(ParameterVector(SpatialJoin.INPUT2, "Join vector layer", [ParameterVector.VECTOR_TYPE_ANY])) | ||
self.addParameter(ParameterSelection(self.SUMMARY, "Attribute summary", self.SUMMARYS, 0)) | ||
self.addParameter(ParameterString(self.STATS, "Statistics for summary (comma separated)", "sum,mean,min,max,med")) | ||
self.addParameter(ParameterSelection(self.KEEP, "Output table", self.KEEPS, 0)) | ||
self.addOutput(OutputVector(SpatialJoin.OUTPUT, "Output layer")) | ||
|
||
def processAlgorithm(self, progress): | ||
progress.setText('Analysing inputs...') | ||
|
||
summary = self.getParameterValue(self.SUMMARY) == 1 | ||
sumList = self.getParameterValue(self.STATS).upper().replace(' ','').split(',') | ||
keep = self.getParameterValue(self.KEEP) == 1 | ||
|
||
input1 = self.getParameterValue(self.INPUT1) | ||
layer1 = dataobjects.getObjectFromUri(input1) | ||
provider1 = layer1.dataProvider() | ||
fieldList1 = provider1.fields() | ||
|
||
input2 = self.getParameterValue(self.INPUT2) | ||
layer2 = dataobjects.getObjectFromUri(input2) | ||
provider2 = layer2.dataProvider() | ||
fieldList2 = provider2.fields() | ||
|
||
fieldList = QgsFields() | ||
if not summary: | ||
fieldList2 = vector.testForUniqueness(fieldList1, fieldList2) | ||
seq = range(0, len(fieldList1) + len(fieldList2)) | ||
fieldList1.extend(fieldList2) | ||
fieldList1 = dict(zip(seq, fieldList1)) | ||
else: | ||
numFields = {} | ||
for j in xrange(len(fieldList2)): | ||
if fieldList2[j].type() == QVariant.Int or fieldList2[j].type() == QVariant.Double: | ||
numFields[j] = [] | ||
for i in sumList: | ||
field = QgsField(i + unicode(fieldList2[j].name()), QVariant.Double, "real", 24, 16, "Summary field" ) | ||
fieldList.append(field) | ||
field = QgsField("COUNT", QVariant.Double, "real", 24, 0, "Summary field" ) | ||
fieldList.append(field) | ||
fieldList2 = vector.testForUniqueness(fieldList1, fieldList) | ||
fieldList1.extend(fieldList) | ||
seq = range(0, len(fieldList1)) | ||
fieldList1 = dict(zip(seq, fieldList1)) | ||
|
||
sRs = provider1.crs() | ||
progress.setPercentage(13) | ||
fields = QgsFields() | ||
for f in fieldList1.values(): | ||
fields.append(f) | ||
output = self.getOutputFromName(self.OUTPUT) | ||
writer = output.getVectorWriter(fields, provider1.geometryType(), sRs) | ||
|
||
inFeat = QgsFeature() | ||
outFeat = QgsFeature() | ||
inFeatB = QgsFeature() | ||
inGeom = QgsGeometry() | ||
|
||
progress.setPercentage(15) | ||
start = 15.00 | ||
add = 85.00 / provider1.featureCount() | ||
|
||
progress.setText('Creating spatial index...') | ||
index = vector.spatialindex(layer2) | ||
progress.setText('Processing spatial join...') | ||
fit1 = provider1.getFeatures() | ||
while fit1.nextFeature(inFeat): | ||
inGeom = inFeat.geometry() | ||
atMap1 = inFeat.attributes() | ||
outFeat.setGeometry(inGeom) | ||
none = True | ||
joinList = [] | ||
if inGeom.type() == QGis.Point: | ||
joinList = index.intersects( inGeom.buffer(10,2).boundingBox() ) | ||
if len(joinList) > 0: check = 0 | ||
else: check = 1 | ||
else: | ||
joinList = index.intersects( inGeom.boundingBox() ) | ||
if len(joinList) > 0: check = 0 | ||
else: check = 1 | ||
if check == 0: | ||
count = 0 | ||
for i in joinList: | ||
provider2.getFeatures( QgsFeatureRequest().setFilterFid( int(i) ) ).nextFeature( inFeatB ) | ||
if inGeom.intersects(inFeatB.geometry()): | ||
count = count + 1 | ||
none = False | ||
atMap2 = inFeatB.attributes() | ||
if not summary: | ||
atMap = atMap1 | ||
atMap2 = atMap2 | ||
atMap.extend(atMap2) | ||
atMap = dict(zip(seq, atMap)) | ||
break | ||
else: | ||
for j in numFields.keys(): | ||
numFields[j].append(atMap2[j]) | ||
if summary and not none: | ||
atMap = atMap1 | ||
for j in numFields.keys(): | ||
for k in sumList: | ||
if k == "SUM": atMap.append(sum(numFields[j])) | ||
elif k == "MEAN": atMap.append(sum(numFields[j]) / count) | ||
elif k == "MIN": atMap.append(min(numFields[j])) | ||
elif k == "MED": atMap.append(myself(numFields[j])) | ||
else: atMap.append(max(numFields[j])) | ||
numFields[j] = [] | ||
atMap.append(count) | ||
atMap = dict(zip(seq, atMap)) | ||
if none: | ||
outFeat.setAttributes(atMap1) | ||
else: | ||
outFeat.setAttributes(atMap.values()) | ||
if keep: # keep all records | ||
writer.addFeature(outFeat) | ||
else: # keep only matching records | ||
if not none: | ||
writer.addFeature(outFeat) | ||
start = start + add | ||
progress.setPercentage(start) | ||
|
||
del writer |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters