本文整理汇总了Python中qgis.core.QgsProcessingUtils.createSpatialIndex方法的典型用法代码示例。如果您正苦于以下问题:Python QgsProcessingUtils.createSpatialIndex方法的具体用法?Python QgsProcessingUtils.createSpatialIndex怎么用?Python QgsProcessingUtils.createSpatialIndex使用的例子?那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类qgis.core.QgsProcessingUtils
的用法示例。
在下文中一共展示了QgsProcessingUtils.createSpatialIndex方法的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的Python代码示例。
示例1: regularMatrix
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def regularMatrix(self, context, inLayer, inField, targetLayer, targetField,
nPoints, feedback):
index = QgsProcessingUtils.createSpatialIndex(targetLayer, context)
inIdx = inLayer.fields().lookupField(inField)
distArea = QgsDistanceArea()
first = True
features = QgsProcessingUtils.getFeatures(inLayer, context)
total = 100.0 / QgsProcessingUtils.featureCount(inLayer, context)
for current, inFeat in enumerate(features):
inGeom = inFeat.geometry()
inID = str(inFeat.attributes()[inIdx])
featList = index.nearestNeighbor(inGeom.asPoint(), nPoints)
if first:
first = False
data = ['ID']
for i in range(len(featList)):
data.append('DIST_{0}'.format(i + 1))
self.writer.addRecord(data)
data = [inID]
for i in featList:
request = QgsFeatureRequest().setFilterFid(i)
outFeat = next(targetLayer.getFeatures(request))
outGeom = outFeat.geometry()
dist = distArea.measureLine(inGeom.asPoint(),
outGeom.asPoint())
data.append(str(float(dist)))
self.writer.addRecord(data)
feedback.setProgress(int(current * total))
示例2: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
layerA = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(Difference.INPUT), context)
layerB = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(Difference.OVERLAY), context)
geomType = QgsWkbTypes.multiType(layerA.wkbType())
writer = self.getOutputFromName(
Difference.OUTPUT).getVectorWriter(layerA.fields(), geomType, layerA.crs(), context)
outFeat = QgsFeature()
index = QgsProcessingUtils.createSpatialIndex(layerB, context)
selectionA = QgsProcessingUtils.getFeatures(layerA, context)
total = 100.0 / layerA.featureCount() if layerA.featureCount() else 0
for current, inFeatA in enumerate(selectionA):
geom = inFeatA.geometry()
diff_geom = QgsGeometry(geom)
attrs = inFeatA.attributes()
intersections = index.intersects(geom.boundingBox())
request = QgsFeatureRequest().setFilterFids(intersections).setSubsetOfAttributes([])
for inFeatB in layerB.getFeatures(request):
tmpGeom = inFeatB.geometry()
if diff_geom.intersects(tmpGeom):
diff_geom = QgsGeometry(diff_geom.difference(tmpGeom))
try:
outFeat.setGeometry(diff_geom)
outFeat.setAttributes(attrs)
writer.addFeature(outFeat, QgsFeatureSink.FastInsert)
except:
QgsMessageLog.logMessage(self.tr('Feature geometry error: One or more output features ignored due to invalid geometry.'), self.tr('Processing'), QgsMessageLog.WARNING)
continue
feedback.setProgress(int(current * total))
del writer
示例3: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
layerPoints = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.POINTS), context)
layerHubs = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.HUBS), context)
fieldName = self.getParameterValue(self.FIELD)
units = self.UNITS[self.getParameterValue(self.UNIT)]
if layerPoints.source() == layerHubs.source():
raise GeoAlgorithmExecutionException(
self.tr('Same layer given for both hubs and spokes'))
fields = layerPoints.fields()
fields.append(QgsField('HubName', QVariant.String))
fields.append(QgsField('HubDist', QVariant.Double))
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fields, QgsWkbTypes.LineString, layerPoints.crs(),
context)
index = QgsProcessingUtils.createSpatialIndex(layerHubs, context)
distance = QgsDistanceArea()
distance.setSourceCrs(layerPoints.crs())
distance.setEllipsoid(QgsProject.instance().ellipsoid())
# Scan source points, find nearest hub, and write to output file
features = QgsProcessingUtils.getFeatures(layerPoints, context)
total = 100.0 / layerPoints.featureCount() if layerPoints.featureCount() else 0
for current, f in enumerate(features):
src = f.geometry().boundingBox().center()
neighbors = index.nearestNeighbor(src, 1)
ft = next(layerHubs.getFeatures(QgsFeatureRequest().setFilterFid(neighbors[0]).setSubsetOfAttributes([fieldName], layerHubs.fields())))
closest = ft.geometry().boundingBox().center()
hubDist = distance.measureLine(src, closest)
attributes = f.attributes()
attributes.append(ft[fieldName])
if units == 'Feet':
attributes.append(hubDist * 3.2808399)
elif units == 'Miles':
attributes.append(hubDist * 0.000621371192)
elif units == 'Kilometers':
attributes.append(hubDist / 1000.0)
elif units != 'Meters':
attributes.append(sqrt(
pow(src.x() - closest.x(), 2.0) +
pow(src.y() - closest.y(), 2.0)))
else:
attributes.append(hubDist)
feat = QgsFeature()
feat.setAttributes(attributes)
feat.setGeometry(QgsGeometry.fromPolyline([src, closest]))
writer.addFeature(feat, QgsFeatureSink.FastInsert)
feedback.setProgress(int(current * total))
del writer
示例4: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
filename = self.getParameterValue(self.INPUT)
inputLayer = QgsProcessingUtils.mapLayerFromString(filename, context)
method = self.getParameterValue(self.METHOD)
filename2 = self.getParameterValue(self.INTERSECT)
selectLayer = QgsProcessingUtils.mapLayerFromString(filename2, context)
predicates = self.getParameterValue(self.PREDICATE)
precision = self.getParameterValue(self.PRECISION)
oldSelection = set(inputLayer.selectedFeatureIds())
inputLayer.removeSelection()
index = QgsProcessingUtils.createSpatialIndex(inputLayer, context)
if 'disjoint' in predicates:
disjoinSet = []
for feat in QgsProcessingUtils.getFeatures(inputLayer, context):
disjoinSet.append(feat.id())
geom = QgsGeometry()
selectedSet = []
features = QgsProcessingUtils.getFeatures(selectLayer, context)
total = 100.0 / QgsProcessingUtils.featureCount(selectLayer, context)
for current, f in enumerate(features):
geom = vector.snapToPrecision(f.geometry(), precision)
bbox = geom.boundingBox()
bbox.grow(0.51 * precision)
intersects = index.intersects(bbox)
request = QgsFeatureRequest().setFilterFids(intersects).setSubsetOfAttributes([])
for feat in inputLayer.getFeatures(request):
tmpGeom = vector.snapToPrecision(feat.geometry(), precision)
res = False
for predicate in predicates:
if predicate == 'disjoint':
if tmpGeom.intersects(geom):
try:
disjoinSet.remove(feat.id())
except:
pass # already removed
else:
res = getattr(tmpGeom, predicate)(geom)
if res:
selectedSet.append(feat.id())
break
feedback.setProgress(int(current * total))
if 'disjoint' in predicates:
selectedSet = selectedSet + disjoinSet
if method == 1:
selectedSet = list(oldSelection.union(selectedSet))
elif method == 2:
selectedSet = list(oldSelection.difference(selectedSet))
inputLayer.selectByIds(selectedSet)
self.setOutputValue(self.OUTPUT, filename)
示例5: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
polyLayer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.POLYGONS), context)
pointLayer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.POINTS), context)
fieldName = self.getParameterValue(self.FIELD)
fieldIdx = pointLayer.fields().lookupField(self.getParameterValue(self.WEIGHT))
fields = polyLayer.fields()
fields.append(QgsField(fieldName, QVariant.Int))
(idxCount, fieldList) = vector.findOrCreateField(polyLayer,
polyLayer.fields(), fieldName)
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fields, polyLayer.wkbType(),
polyLayer.crs(), context)
spatialIndex = QgsProcessingUtils.createSpatialIndex(pointLayer, context)
ftPoint = QgsFeature()
outFeat = QgsFeature()
geom = QgsGeometry()
features = QgsProcessingUtils.getFeatures(polyLayer, context)
total = 100.0 / QgsProcessingUtils.featureCount(polyLayer, context)
for current, ftPoly in enumerate(features):
geom = ftPoly.geometry()
engine = QgsGeometry.createGeometryEngine(geom.geometry())
engine.prepareGeometry()
attrs = ftPoly.attributes()
count = 0
points = spatialIndex.intersects(geom.boundingBox())
if len(points) > 0:
feedback.setProgressText(str(len(points)))
request = QgsFeatureRequest().setFilterFids(points).setSubsetOfAttributes([fieldIdx])
fit = pointLayer.getFeatures(request)
ftPoint = QgsFeature()
while fit.nextFeature(ftPoint):
tmpGeom = QgsGeometry(ftPoint.geometry())
if engine.contains(tmpGeom.geometry()):
weight = str(ftPoint.attributes()[fieldIdx])
try:
count += float(weight)
except:
# Ignore fields with non-numeric values
pass
outFeat.setGeometry(geom)
if idxCount == len(attrs):
attrs.append(count)
else:
attrs[idxCount] = count
outFeat.setAttributes(attrs)
writer.addFeature(outFeat)
feedback.setProgress(int(current * total))
del writer
示例6: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
polyLayer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.POLYGONS), context)
pointLayer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.POINTS), context)
fieldName = self.getParameterValue(self.FIELD)
classFieldName = self.getParameterValue(self.CLASSFIELD)
fields = polyLayer.fields()
fields.append(QgsField(fieldName, QVariant.Int))
classFieldIndex = pointLayer.fields().lookupField(classFieldName)
(idxCount, fieldList) = vector.findOrCreateField(polyLayer,
polyLayer.fields(), fieldName)
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fields, polyLayer.wkbType(),
polyLayer.crs(), context)
spatialIndex = QgsProcessingUtils.createSpatialIndex(pointLayer, context)
ftPoint = QgsFeature()
outFeat = QgsFeature()
geom = QgsGeometry()
features = QgsProcessingUtils.getFeatures(polyLayer, context)
total = 100.0 / polyLayer.featureCount() if polyLayer.featureCount() else 0
for current, ftPoly in enumerate(features):
geom = ftPoly.geometry()
engine = QgsGeometry.createGeometryEngine(geom.geometry())
engine.prepareGeometry()
attrs = ftPoly.attributes()
classes = set()
points = spatialIndex.intersects(geom.boundingBox())
if len(points) > 0:
request = QgsFeatureRequest().setFilterFids(points).setSubsetOfAttributes([classFieldIndex])
fit = pointLayer.getFeatures(request)
ftPoint = QgsFeature()
while fit.nextFeature(ftPoint):
tmpGeom = QgsGeometry(ftPoint.geometry())
if engine.contains(tmpGeom.geometry()):
clazz = ftPoint.attributes()[classFieldIndex]
if clazz not in classes:
classes.add(clazz)
outFeat.setGeometry(geom)
if idxCount == len(attrs):
attrs.append(len(classes))
else:
attrs[idxCount] = len(classes)
outFeat.setAttributes(attrs)
writer.addFeature(outFeat, QgsFeatureSink.FastInsert)
feedback.setProgress(int(current * total))
del writer
示例7: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
vlayerA = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.INPUT), context)
vlayerB = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.INPUT2), context)
geomType = QgsWkbTypes.multiType(vlayerA.wkbType())
fields = vector.combineVectorFields(vlayerA, vlayerB)
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fields, geomType, vlayerA.crs(), context)
outFeat = QgsFeature()
index = QgsProcessingUtils.createSpatialIndex(vlayerB, context)
selectionA = QgsProcessingUtils.getFeatures(vlayerA, context)
total = 100.0 / vlayerA.featureCount() if vlayerA.featureCount() else 0
for current, inFeatA in enumerate(selectionA):
feedback.setProgress(int(current * total))
geom = inFeatA.geometry()
atMapA = inFeatA.attributes()
intersects = index.intersects(geom.boundingBox())
request = QgsFeatureRequest().setFilterFids(intersects)
engine = None
if len(intersects) > 0:
# use prepared geometries for faster intersection tests
engine = QgsGeometry.createGeometryEngine(geom.geometry())
engine.prepareGeometry()
for inFeatB in vlayerB.getFeatures(request):
tmpGeom = inFeatB.geometry()
if engine.intersects(tmpGeom.geometry()):
atMapB = inFeatB.attributes()
int_geom = QgsGeometry(geom.intersection(tmpGeom))
if int_geom.wkbType() == QgsWkbTypes.Unknown or QgsWkbTypes.flatType(int_geom.geometry().wkbType()) == QgsWkbTypes.GeometryCollection:
int_com = geom.combine(tmpGeom)
int_geom = QgsGeometry()
if int_com:
int_sym = geom.symDifference(tmpGeom)
int_geom = QgsGeometry(int_com.difference(int_sym))
if int_geom.isEmpty() or not int_geom.isGeosValid():
raise GeoAlgorithmExecutionException(
self.tr('GEOS geoprocessing error: One or '
'more input features have invalid '
'geometry.'))
try:
if int_geom.wkbType() in wkbTypeGroups[wkbTypeGroups[int_geom.wkbType()]]:
outFeat.setGeometry(int_geom)
attrs = []
attrs.extend(atMapA)
attrs.extend(atMapB)
outFeat.setAttributes(attrs)
writer.addFeature(outFeat, QgsFeatureSink.FastInsert)
except:
raise GeoAlgorithmExecutionException(
self.tr('Feature geometry error: One or more '
'output features ignored due to invalid '
'geometry.'))
del writer
示例8: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
filename = self.getParameterValue(self.INPUT)
layer = QgsProcessingUtils.mapLayerFromString(filename, context)
filename = self.getParameterValue(self.INTERSECT)
selectLayer = QgsProcessingUtils.mapLayerFromString(filename, context)
predicates = self.getParameterValue(self.PREDICATE)
precision = self.getParameterValue(self.PRECISION)
index = QgsProcessingUtils.createSpatialIndex(layer, context)
output = self.getOutputFromName(self.OUTPUT)
writer = output.getVectorWriter(layer.fields(), layer.wkbType(), layer.crs(), context)
if 'disjoint' in predicates:
disjoinSet = []
for feat in QgsProcessingUtils.getFeatures(layer, context):
disjoinSet.append(feat.id())
selectedSet = []
features = QgsProcessingUtils.getFeatures(selectLayer, context)
total = 100.0 / selectLayer.featureCount() if selectLayer.featureCount() else 0
for current, f in enumerate(features):
geom = vector.snapToPrecision(f.geometry(), precision)
bbox = geom.boundingBox()
bbox.grow(0.51 * precision)
intersects = index.intersects(bbox)
request = QgsFeatureRequest().setFilterFids(intersects).setSubsetOfAttributes([])
for feat in layer.getFeatures(request):
tmpGeom = vector.snapToPrecision(feat.geometry(), precision)
res = False
for predicate in predicates:
if predicate == 'disjoint':
if tmpGeom.intersects(geom):
try:
disjoinSet.remove(feat.id())
except:
pass # already removed
else:
res = getattr(tmpGeom, predicate)(geom)
if res:
selectedSet.append(feat.id())
break
feedback.setProgress(int(current * total))
if 'disjoint' in predicates:
selectedSet = selectedSet + disjoinSet
features = QgsProcessingUtils.getFeatures(layer, context)
total = 100.0 / layer.featureCount() if layer.featureCount() else 0
for current, f in enumerate(features):
if f.id() in selectedSet:
writer.addFeature(f, QgsFeatureSink.FastInsert)
feedback.setProgress(int(current * total))
del writer
示例9: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
layer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.VECTOR), context)
pointCount = int(self.getParameterValue(self.POINT_NUMBER))
minDistance = float(self.getParameterValue(self.MIN_DISTANCE))
bbox = layer.extent()
idxLayer = QgsProcessingUtils.createSpatialIndex(layer, context)
fields = QgsFields()
fields.append(QgsField('id', QVariant.Int, '', 10, 0))
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fields, QgsWkbTypes.Point, layer.crs(), context)
nPoints = 0
nIterations = 0
maxIterations = pointCount * 200
total = 100.0 / pointCount
index = QgsSpatialIndex()
points = dict()
random.seed()
while nIterations < maxIterations and nPoints < pointCount:
rx = bbox.xMinimum() + bbox.width() * random.random()
ry = bbox.yMinimum() + bbox.height() * random.random()
pnt = QgsPointXY(rx, ry)
geom = QgsGeometry.fromPoint(pnt)
ids = idxLayer.intersects(geom.buffer(5, 5).boundingBox())
if len(ids) > 0 and \
vector.checkMinDistance(pnt, index, minDistance, points):
request = QgsFeatureRequest().setFilterFids(ids).setSubsetOfAttributes([])
for f in layer.getFeatures(request):
tmpGeom = f.geometry()
if geom.within(tmpGeom):
f = QgsFeature(nPoints)
f.initAttributes(1)
f.setFields(fields)
f.setAttribute('id', nPoints)
f.setGeometry(geom)
writer.addFeature(f)
index.insertFeature(f)
points[nPoints] = pnt
nPoints += 1
feedback.setProgress(int(nPoints * total))
nIterations += 1
if nPoints < pointCount:
QgsMessageLog.logMessage(self.tr('Can not generate requested number of random points. '
'Maximum number of attempts exceeded.'), self.tr('Processing'), QgsMessageLog.INFO)
del writer
示例10: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
layer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.POINTS), context)
output = self.getOutputValue(self.OUTPUT)
spatialIndex = QgsProcessingUtils.createSpatialIndex(layer, context)
neighbour = QgsFeature()
distance = QgsDistanceArea()
sumDist = 0.00
A = layer.extent()
A = float(A.width() * A.height())
features = QgsProcessingUtils.getFeatures(layer, context)
count = QgsProcessingUtils.featureCount(layer, context)
total = 100.0 / count
for current, feat in enumerate(features):
neighbourID = spatialIndex.nearestNeighbor(
feat.geometry().asPoint(), 2)[1]
request = QgsFeatureRequest().setFilterFid(neighbourID).setSubsetOfAttributes([])
neighbour = next(layer.getFeatures(request))
sumDist += distance.measureLine(neighbour.geometry().asPoint(),
feat.geometry().asPoint())
feedback.setProgress(int(current * total))
do = float(sumDist) / count
de = float(0.5 / math.sqrt(count / A))
d = float(do / de)
SE = float(0.26136 / math.sqrt(count ** 2 / A))
zscore = float((do - de) / SE)
data = []
data.append('Observed mean distance: ' + str(do))
data.append('Expected mean distance: ' + str(de))
data.append('Nearest neighbour index: ' + str(d))
data.append('Number of points: ' + str(count))
data.append('Z-Score: ' + str(zscore))
self.createHTML(output, data)
self.setOutputValue(self.OBSERVED_MD, float(data[0].split(': ')[1]))
self.setOutputValue(self.EXPECTED_MD, float(data[1].split(': ')[1]))
self.setOutputValue(self.NN_INDEX, float(data[2].split(': ')[1]))
self.setOutputValue(self.POINT_COUNT, float(data[3].split(': ')[1]))
self.setOutputValue(self.Z_SCORE, float(data[4].split(': ')[1]))
示例11: linearMatrix
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def linearMatrix(self, context, inLayer, inField, targetLayer, targetField,
matType, nPoints, feedback):
if matType == 0:
self.writer.addRecord(['InputID', 'TargetID', 'Distance'])
else:
self.writer.addRecord(['InputID', 'MEAN', 'STDDEV', 'MIN', 'MAX'])
index = QgsProcessingUtils.createSpatialIndex(targetLayer, context)
inIdx = inLayer.fields().lookupField(inField)
outIdx = targetLayer.fields().lookupField(targetField)
distArea = QgsDistanceArea()
features = QgsProcessingUtils.getFeatures(inLayer, context)
total = 100.0 / QgsProcessingUtils.featureCount(inLayer, context)
for current, inFeat in enumerate(features):
inGeom = inFeat.geometry()
inID = str(inFeat.attributes()[inIdx])
featList = index.nearestNeighbor(inGeom.asPoint(), nPoints)
distList = []
vari = 0.0
request = QgsFeatureRequest().setFilterFids(featList).setSubsetOfAttributes([outIdx])
for outFeat in targetLayer.getFeatures(request):
outID = outFeat.attributes()[outIdx]
outGeom = outFeat.geometry()
dist = distArea.measureLine(inGeom.asPoint(),
outGeom.asPoint())
if matType == 0:
self.writer.addRecord([inID, str(outID), str(dist)])
else:
distList.append(float(dist))
if matType != 0:
mean = sum(distList) / len(distList)
for i in distList:
vari += (i - mean) * (i - mean)
vari = math.sqrt(vari / len(distList))
self.writer.addRecord([inID, str(mean),
str(vari), str(min(distList)),
str(max(distList))])
feedback.setProgress(int(current * total))
示例12: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
layerA = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.INPUT_A), context)
layerB = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.INPUT_B), context)
sameLayer = self.getParameterValue(self.INPUT_A) == self.getParameterValue(self.INPUT_B)
fieldList = layerA.fields()
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fieldList, QgsWkbTypes.LineString, layerA.crs(),
context)
spatialIndex = QgsProcessingUtils.createSpatialIndex(layerB, context)
outFeat = QgsFeature()
features = QgsProcessingUtils.getFeatures(layerA, context)
total = 100.0 / QgsProcessingUtils.featureCount(layerA, context)
for current, inFeatA in enumerate(features):
inGeom = inFeatA.geometry()
attrsA = inFeatA.attributes()
outFeat.setAttributes(attrsA)
inLines = [inGeom]
lines = spatialIndex.intersects(inGeom.boundingBox())
engine = None
if len(lines) > 0:
# use prepared geometries for faster intersection tests
engine = QgsGeometry.createGeometryEngine(inGeom.geometry())
engine.prepareGeometry()
if len(lines) > 0: # hasIntersections
splittingLines = []
request = QgsFeatureRequest().setFilterFids(lines).setSubsetOfAttributes([])
for inFeatB in layerB.getFeatures(request):
# check if trying to self-intersect
if sameLayer:
if inFeatA.id() == inFeatB.id():
continue
splitGeom = inFeatB.geometry()
if engine.intersects(splitGeom.geometry()):
splittingLines.append(splitGeom)
if len(splittingLines) > 0:
for splitGeom in splittingLines:
splitterPList = vector.extractPoints(splitGeom)
outLines = []
split_geom_engine = QgsGeometry.createGeometryEngine(splitGeom.geometry())
split_geom_engine.prepareGeometry()
while len(inLines) > 0:
inGeom = inLines.pop()
inPoints = vector.extractPoints(inGeom)
if split_geom_engine.intersects(inGeom.geometry()):
try:
result, newGeometries, topoTestPoints = inGeom.splitGeometry(splitterPList, False)
except:
QgsMessageLog.logMessage(self.tr('Geometry exception while splitting'),
self.tr('Processing'), QgsMessageLog.WARNING)
result = 1
# splitGeometry: If there are several intersections
# between geometry and splitLine, only the first one is considered.
if result == 0: # split occurred
if inPoints == vector.extractPoints(inGeom):
# bug in splitGeometry: sometimes it returns 0 but
# the geometry is unchanged
outLines.append(inGeom)
else:
inLines.append(inGeom)
for aNewGeom in newGeometries:
inLines.append(aNewGeom)
else:
outLines.append(inGeom)
else:
outLines.append(inGeom)
inLines = outLines
for aLine in inLines:
if len(aLine.asPolyline()) > 2 or \
(len(aLine.asPolyline()) == 2 and
aLine.asPolyline()[0] != aLine.asPolyline()[1]):
# sometimes splitting results in lines of zero length
outFeat.setGeometry(aLine)
writer.addFeature(outFeat)
feedback.setProgress(int(current * total))
del writer
示例13: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
layerA = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.INPUT), context)
layerB = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.OVERLAY), context)
geomType = QgsWkbTypes.multiType(layerA.wkbType())
fields = vector.combineVectorFields(layerA, layerB)
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fields, geomType, layerA.crs(), context)
featB = QgsFeature()
outFeat = QgsFeature()
indexA = QgsProcessingUtils.createSpatialIndex(layerB, context)
indexB = QgsProcessingUtils.createSpatialIndex(layerA, context)
featuresA = QgsProcessingUtils.getFeatures(layerA, context)
featuresB = QgsProcessingUtils.getFeatures(layerB, context)
total = 100.0 / (QgsProcessingUtils.featureCount(layerA, context) * QgsProcessingUtils.featureCount(layerB, context))
count = 0
for featA in featuresA:
geom = featA.geometry()
diffGeom = QgsGeometry(geom)
attrs = featA.attributes()
intersects = indexA.intersects(geom.boundingBox())
request = QgsFeatureRequest().setFilterFids(intersects).setSubsetOfAttributes([])
for featB in layerB.getFeatures(request):
tmpGeom = featB.geometry()
if diffGeom.intersects(tmpGeom):
diffGeom = QgsGeometry(diffGeom.difference(tmpGeom))
try:
outFeat.setGeometry(diffGeom)
outFeat.setAttributes(attrs)
writer.addFeature(outFeat)
except:
QgsMessageLog.logMessage(self.tr('Feature geometry error: One or more output features ignored due to invalid geometry.'),
self.tr('Processing'), QgsMessageLog.WARNING)
continue
count += 1
feedback.setProgress(int(count * total))
length = len(layerA.fields())
for featA in featuresB:
geom = featA.geometry()
diffGeom = QgsGeometry(geom)
attrs = featA.attributes()
attrs = [NULL] * length + attrs
intersects = indexB.intersects(geom.boundingBox())
request = QgsFeatureRequest().setFilterFids(intersects).setSubsetOfAttributes([])
for featB in layerA.getFeatures(request):
tmpGeom = featB.geometry()
if diffGeom.intersects(tmpGeom):
diffGeom = QgsGeometry(diffGeom.difference(tmpGeom))
try:
outFeat.setGeometry(diffGeom)
outFeat.setAttributes(attrs)
writer.addFeature(outFeat)
except:
QgsMessageLog.logMessage(self.tr('Feature geometry error: One or more output features ignored due to invalid geometry.'),
self.tr('Processing'), QgsMessageLog.WARNING)
continue
count += 1
feedback.setProgress(int(count * total))
del writer
示例14: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
target = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.TARGET), context)
join = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.JOIN), context)
predicates = self.getParameterValue(self.PREDICATE)
precision = self.getParameterValue(self.PRECISION)
summary = self.getParameterValue(self.SUMMARY) == 1
keep = self.getParameterValue(self.KEEP) == 1
sumList = self.getParameterValue(self.STATS).lower().split(',')
targetFields = target.fields()
joinFields = join.fields()
fieldList = QgsFields()
if not summary:
joinFields = vector.testForUniqueness(targetFields, joinFields)
seq = list(range(len(targetFields) + len(joinFields)))
targetFields.extend(joinFields)
targetFields = dict(list(zip(seq, targetFields)))
else:
numFields = {}
for j in range(len(joinFields)):
if joinFields[j].type() in [QVariant.Int, QVariant.Double, QVariant.LongLong, QVariant.UInt, QVariant.ULongLong]:
numFields[j] = []
for i in sumList:
field = QgsField(i + str(joinFields[j].name()), QVariant.Double, '', 24, 16)
fieldList.append(field)
field = QgsField('count', QVariant.Double, '', 24, 16)
fieldList.append(field)
joinFields = vector.testForUniqueness(targetFields, fieldList)
targetFields.extend(fieldList)
seq = list(range(len(targetFields)))
targetFields = dict(list(zip(seq, targetFields)))
fields = QgsFields()
for f in list(targetFields.values()):
fields.append(f)
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fields, target.wkbType(), target.crs(), context)
outFeat = QgsFeature()
inFeatB = QgsFeature()
inGeom = QgsGeometry()
index = QgsProcessingUtils.createSpatialIndex(join, context)
mapP2 = dict()
features = QgsProcessingUtils.getFeatures(join, context)
for f in features:
mapP2[f.id()] = QgsFeature(f)
features = QgsProcessingUtils.getFeatures(target, context)
total = 100.0 / target.featureCount() if target.featureCount() else 0
for c, f in enumerate(features):
atMap1 = f.attributes()
outFeat.setGeometry(f.geometry())
inGeom = vector.snapToPrecision(f.geometry(), precision)
none = True
joinList = []
if inGeom.type() == QgsWkbTypes.PointGeometry:
bbox = inGeom.buffer(10, 2).boundingBox()
else:
bbox = inGeom.boundingBox()
bbox.grow(0.51 * precision)
joinList = index.intersects(bbox)
if len(joinList) > 0:
count = 0
for i in joinList:
inFeatB = mapP2[i]
inGeomB = vector.snapToPrecision(inFeatB.geometry(), precision)
res = False
for predicate in predicates:
res = getattr(inGeom, predicate)(inGeomB)
if res:
break
if res:
count = count + 1
none = False
atMap2 = inFeatB.attributes()
if not summary:
atMap = atMap1
atMap2 = atMap2
atMap.extend(atMap2)
atMap = dict(list(zip(seq, atMap)))
break
else:
for j in list(numFields.keys()):
numFields[j].append(atMap2[j])
if summary and not none:
atMap = atMap1
for j in list(numFields.keys()):
for k in sumList:
if k == 'sum':
atMap.append(sum(self._filterNull(numFields[j])))
elif k == 'mean':
#.........这里部分代码省略.........
示例15: processAlgorithm
# 需要导入模块: from qgis.core import QgsProcessingUtils [as 别名]
# 或者: from qgis.core.QgsProcessingUtils import createSpatialIndex [as 别名]
def processAlgorithm(self, parameters, context, feedback):
lineLayer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.LINES), context)
polyLayer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.POLYGONS), context)
lengthFieldName = self.getParameterValue(self.LEN_FIELD)
countFieldName = self.getParameterValue(self.COUNT_FIELD)
(idxLength, fieldList) = vector.findOrCreateField(polyLayer,
polyLayer.fields(), lengthFieldName)
(idxCount, fieldList) = vector.findOrCreateField(polyLayer, fieldList,
countFieldName)
writer = self.getOutputFromName(self.OUTPUT).getVectorWriter(fieldList, polyLayer.wkbType(),
polyLayer.crs(), context)
spatialIndex = QgsProcessingUtils.createSpatialIndex(lineLayer, context)
ftLine = QgsFeature()
ftPoly = QgsFeature()
outFeat = QgsFeature()
inGeom = QgsGeometry()
outGeom = QgsGeometry()
distArea = QgsDistanceArea()
features = QgsProcessingUtils.getFeatures(polyLayer, context)
total = 100.0 / QgsProcessingUtils.featureCount(polyLayer, context)
hasIntersections = False
for current, ftPoly in enumerate(features):
inGeom = ftPoly.geometry()
attrs = ftPoly.attributes()
count = 0
length = 0
hasIntersections = False
lines = spatialIndex.intersects(inGeom.boundingBox())
engine = None
if len(lines) > 0:
hasIntersections = True
# use prepared geometries for faster intersection tests
engine = QgsGeometry.createGeometryEngine(inGeom.geometry())
engine.prepareGeometry()
if hasIntersections:
request = QgsFeatureRequest().setFilterFids(lines).setSubsetOfAttributes([])
for ftLine in lineLayer.getFeatures(request):
tmpGeom = ftLine.geometry()
if engine.intersects(tmpGeom.geometry()):
outGeom = inGeom.intersection(tmpGeom)
length += distArea.measureLength(outGeom)
count += 1
outFeat.setGeometry(inGeom)
if idxLength == len(attrs):
attrs.append(length)
else:
attrs[idxLength] = length
if idxCount == len(attrs):
attrs.append(count)
else:
attrs[idxCount] = count
outFeat.setAttributes(attrs)
writer.addFeature(outFeat)
feedback.setProgress(int(current * total))
del writer