Skip to content
Permalink
Browse files

Port Points Displacement algorithm to new API

Rework algorithm to use same approach as points displacement
renderer. Also maintain Z/M values, and add tests and docs.
  • Loading branch information
nyalldawson committed Aug 8, 2017
1 parent 69c991e commit c8ac7841f0436692cb955790548a3b1ee75e6902
@@ -360,7 +360,7 @@ qgis:pointsalonglines: >
An optional start and end offset can be specified, which controls how far from the start and end of the geometry the points should be created.

qgis:pointsdisplacement:

Offsets nearby point features by moving nearby points by a preset amount to minimize overlapping features.

qgis:pointslayerfromtable: >
This algorithm generates a points layer based on the values from an input table.
@@ -27,27 +27,26 @@
__revision__ = '$Format:%H$'

import math
from qgis.core import (QgsApplication,
QgsFeatureRequest,
QgsFeature,
QgsFeatureSink,
from qgis.core import (QgsFeatureSink,
QgsGeometry,
QgsPointXY,
QgsProcessingUtils)
from processing.tools import dataobjects
QgsSpatialIndex,
QgsRectangle,
QgsProcessing,
QgsProcessingParameterFeatureSource,
QgsProcessingParameterNumber,
QgsProcessingParameterBoolean,
QgsProcessingParameterFeatureSink)
from processing.algs.qgis.QgisAlgorithm import QgisAlgorithm
from processing.core.parameters import ParameterVector
from processing.core.parameters import ParameterNumber
from processing.core.parameters import ParameterBoolean
from processing.core.outputs import OutputVector


class PointsDisplacement(QgisAlgorithm):

INPUT_LAYER = 'INPUT_LAYER'
INPUT = 'INPUT'
DISTANCE = 'DISTANCE'
PROXIMITY = 'PROXIMITY'
HORIZONTAL = 'HORIZONTAL'
OUTPUT_LAYER = 'OUTPUT_LAYER'
OUTPUT = 'OUTPUT'

def group(self):
return self.tr('Vector geometry tools')
@@ -56,14 +55,17 @@ def __init__(self):
super().__init__()

def initAlgorithm(self, config=None):
self.addParameter(ParameterVector(self.INPUT_LAYER,
self.tr('Input layer'), [dataobjects.TYPE_VECTOR_POINT]))
self.addParameter(ParameterNumber(self.DISTANCE,
self.tr('Displacement distance'),
0.00001, 999999999.999990, 0.00015))
self.addParameter(ParameterBoolean(self.HORIZONTAL,
self.tr('Horizontal distribution for two point case')))
self.addOutput(OutputVector(self.OUTPUT_LAYER, self.tr('Displaced'), datatype=[dataobjects.TYPE_VECTOR_POINT]))
self.addParameter(QgsProcessingParameterFeatureSource(self.INPUT,
self.tr('Input layer'), [QgsProcessing.TypeVectorPoint]))
self.addParameter(QgsProcessingParameterNumber(self.PROXIMITY,
self.tr('Minimum distance to other points'), type=QgsProcessingParameterNumber.Double,
minValue=0.00001, defaultValue=0.00015))
self.addParameter(QgsProcessingParameterNumber(self.DISTANCE,
self.tr('Displacement distance'), type=QgsProcessingParameterNumber.Double,
minValue=0.00001, defaultValue=0.00015))
self.addParameter(QgsProcessingParameterBoolean(self.HORIZONTAL,
self.tr('Horizontal distribution for two point case')))
self.addParameter(QgsProcessingParameterFeatureSink(self.OUTPUT, self.tr('Displaced'), QgsProcessing.TypeVectorPoint))

def name(self):
return 'pointsdisplacement'
@@ -72,64 +74,110 @@ def displayName(self):
return self.tr('Points displacement')

def processAlgorithm(self, parameters, context, feedback):
radius = self.getParameterValue(self.DISTANCE)
horizontal = self.getParameterValue(self.HORIZONTAL)
output = self.getOutputFromName(self.OUTPUT_LAYER)
source = self.parameterAsSource(parameters, self.INPUT, context)
proximity = self.parameterAsDouble(parameters, self.PROXIMITY, context)
radius = self.parameterAsDouble(parameters, self.DISTANCE, context)
horizontal = self.parameterAsBool(parameters, self.HORIZONTAL, context)

layer = QgsProcessingUtils.mapLayerFromString(self.getParameterValue(self.INPUT_LAYER), context)
(sink, dest_id) = self.parameterAsSink(parameters, self.OUTPUT, context,
source.fields(), source.wkbType(), source.sourceCrs())

writer = output.getVectorWriter(layer.fields(), layer.wkbType(), layer.crs(), context)
features = source.getFeatures()

features = QgsProcessingUtils.getFeatures(layer, context)
total = 100.0 / source.featureCount() if source.featureCount() else 0

total = 100.0 / layer.featureCount() if layer.featureCount() else 0
def searchRect(p):
return QgsRectangle(p.x() - proximity, p.y() - proximity, p.x() + proximity, p.y() + proximity)

duplicates = dict()
index = QgsSpatialIndex()

# NOTE: this is a Python port of QgsPointDistanceRenderer::renderFeature. If refining this algorithm,
# please port the changes to QgsPointDistanceRenderer::renderFeature also!

clustered_groups = []
group_index = {}
group_locations = {}
for current, f in enumerate(features):
wkt = f.geometry().exportToWkt()
if wkt not in duplicates:
duplicates[wkt] = [f.id()]
if feedback.isCanceled():
break

if not f.hasGeometry():
continue

point = f.geometry().asPoint()

other_features_within_radius = index.intersects(searchRect(point))
if not other_features_within_radius:
index.insertFeature(f)
group = [f]
clustered_groups.append(group)
group_index[f.id()] = len(clustered_groups) - 1
group_locations[f.id()] = point
else:
duplicates[wkt].extend([f.id()])
# find group with closest location to this point (may be more than one within search tolerance)
min_dist_feature_id = other_features_within_radius[0]
min_dist = group_locations[min_dist_feature_id].distance(point)
for i in range(1, len(other_features_within_radius)):
candidate_id = other_features_within_radius[i]
new_dist = group_locations[candidate_id].distance(point)
if new_dist < min_dist:
min_dist = new_dist
min_dist_feature_id = candidate_id

group_index_pos = group_index[min_dist_feature_id]
group = clustered_groups[group_index_pos]

# calculate new centroid of group
old_center = group_locations[min_dist_feature_id]
group_locations[min_dist_feature_id] = QgsPointXY((old_center.x() * len(group) + point.x()) / (len(group) + 1.0),
(old_center.y() * len(group) + point.y()) / (len(group) + 1.0))
# add to a group
clustered_groups[group_index_pos].append(f)
group_index[f.id()] = group_index_pos

feedback.setProgress(int(current * total))

current = 0
total = 100.0 / len(duplicates) if duplicates else 1
total = 100.0 / len(clustered_groups) if clustered_groups else 1
feedback.setProgress(0)

fullPerimeter = 2 * math.pi

for (geom, fids) in list(duplicates.items()):
count = len(fids)
for group in clustered_groups:
if feedback.isCanceled():
break

count = len(group)
if count == 1:
f = next(layer.getFeatures(QgsFeatureRequest().setFilterFid(fids[0])))
writer.addFeature(f, QgsFeatureSink.FastInsert)
sink.addFeature(group[0], QgsFeatureSink.FastInsert)
else:
angleStep = fullPerimeter / count
if count == 2 and horizontal:
currentAngle = math.pi / 2
else:
currentAngle = 0

old_point = QgsGeometry.fromWkt(geom).asPoint()
old_point = group_locations[group[0].id()]

for f in group:
if feedback.isCanceled():
break

request = QgsFeatureRequest().setFilterFids(fids).setFlags(QgsFeatureRequest.NoGeometry)
for f in layer.getFeatures(request):
sinusCurrentAngle = math.sin(currentAngle)
cosinusCurrentAngle = math.cos(currentAngle)
dx = radius * sinusCurrentAngle
dy = radius * cosinusCurrentAngle

new_point = QgsPointXY(old_point.x() + dx, old_point.y() + dy)
out_feature = QgsFeature()
out_feature.setGeometry(QgsGeometry.fromPoint(new_point))
out_feature.setAttributes(f.attributes())
# we want to keep any existing m/z values
point = f.geometry().geometry().clone()
point.setX(old_point.x() + dx)
point.setY(old_point.y() + dy)
f.setGeometry(QgsGeometry(point))

writer.addFeature(out_feature, QgsFeatureSink.FastInsert)
sink.addFeature(f, QgsFeatureSink.FastInsert)
currentAngle += angleStep

current += 1
feedback.setProgress(int(current * total))

del writer
return {self.OUTPUT: dest_id}
@@ -94,6 +94,7 @@
from .PointDistance import PointDistance
from .PointOnSurface import PointOnSurface
from .PointsAlongGeometry import PointsAlongGeometry
from .PointsDisplacement import PointsDisplacement
from .PointsInPolygon import PointsInPolygon
from .PointsLayerFromTable import PointsLayerFromTable
from .PointsToPaths import PointsToPaths
@@ -153,7 +154,6 @@
# from .GeometryConvert import GeometryConvert
# from .FieldsCalculator import FieldsCalculator
# from .FieldPyculator import FieldsPyculator
# from .PointsDisplacement import PointsDisplacement
# from .PointsFromPolygons import PointsFromPolygons
# from .PointsFromLines import PointsFromLines
# from .SetVectorStyle import SetVectorStyle
@@ -192,7 +192,7 @@ def getAlgs(self):
# GeometryConvert(), FieldsCalculator(),
# FieldsPyculator(),
#
# RasterLayerStatistics(), PointsDisplacement(),
# RasterLayerStatistics(),
# PointsFromPolygons(),
# PointsFromLines(),
# SetVectorStyle(), SetRasterStyle(),
@@ -261,6 +261,7 @@ def getAlgs(self):
PointDistance(),
PointOnSurface(),
PointsAlongGeometry(),
PointsDisplacement(),
PointsInPolygon(),
PointsLayerFromTable(),
PointsToPaths(),
@@ -0,0 +1,26 @@
<GMLFeatureClassList>
<GMLFeatureClass>
<Name>displace_points</Name>
<ElementPath>displace_points</ElementPath>
<!--POINT-->
<GeometryType>1</GeometryType>
<SRSName>EPSG:4326</SRSName>
<DatasetSpecificInfo>
<FeatureCount>4</FeatureCount>
<ExtentXMin>1.00000</ExtentXMin>
<ExtentXMax>4.00000</ExtentXMax>
<ExtentYMin>1.00000</ExtentYMin>
<ExtentYMax>1.00000</ExtentYMax>
</DatasetSpecificInfo>
<PropertyDefn>
<Name>id</Name>
<ElementPath>id</ElementPath>
<Type>Integer</Type>
</PropertyDefn>
<PropertyDefn>
<Name>id2</Name>
<ElementPath>id2</ElementPath>
<Type>Integer</Type>
</PropertyDefn>
</GMLFeatureClass>
</GMLFeatureClassList>
@@ -0,0 +1,41 @@
<?xml version="1.0" encoding="utf-8" ?>
<ogr:FeatureCollection
xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
xsi:schemaLocation="http://ogr.maptools.org/ displace_points.xsd"
xmlns:ogr="http://ogr.maptools.org/"
xmlns:gml="http://www.opengis.net/gml">
<gml:boundedBy>
<gml:Box>
<gml:coord><gml:X>1</gml:X><gml:Y>1</gml:Y></gml:coord>
<gml:coord><gml:X>4</gml:X><gml:Y>1</gml:Y></gml:coord>
</gml:Box>
</gml:boundedBy>

<gml:featureMember>
<ogr:displace_points fid="points.4">
<ogr:geometryProperty><gml:Point srsName="EPSG:4326"><gml:coordinates>4,1</gml:coordinates></gml:Point></ogr:geometryProperty>
<ogr:id>5</ogr:id>
<ogr:id2>1</ogr:id2>
</ogr:displace_points>
</gml:featureMember>
<gml:featureMember>
<ogr:displace_points fid="points.3">
<ogr:geometryProperty><gml:Point srsName="EPSG:4326"><gml:coordinates>4,1</gml:coordinates></gml:Point></ogr:geometryProperty>
<ogr:id>4</ogr:id>
<ogr:id2>2</ogr:id2>
</ogr:displace_points>
</gml:featureMember>
<gml:featureMember>
<ogr:displace_points fid="points.0">
<ogr:geometryProperty><gml:Point srsName="EPSG:4326"><gml:coordinates>1,1</gml:coordinates></gml:Point></ogr:geometryProperty>
<ogr:id>1</ogr:id>
<ogr:id2>2</ogr:id2>
</ogr:displace_points>
</gml:featureMember>
<gml:featureMember>
<ogr:displace_points fid="points.1">
<ogr:id>2</ogr:id>
<ogr:id2>1</ogr:id2>
</ogr:displace_points>
</gml:featureMember>
</ogr:FeatureCollection>
@@ -0,0 +1,26 @@
<GMLFeatureClassList>
<GMLFeatureClass>
<Name>displaced_points</Name>
<ElementPath>displaced_points</ElementPath>
<!--POINT-->
<GeometryType>1</GeometryType>
<SRSName>EPSG:4326</SRSName>
<DatasetSpecificInfo>
<FeatureCount>3</FeatureCount>
<ExtentXMin>1.00000</ExtentXMin>
<ExtentXMax>4.86603</ExtentXMax>
<ExtentYMin>0.50000</ExtentYMin>
<ExtentYMax>2.00000</ExtentYMax>
</DatasetSpecificInfo>
<PropertyDefn>
<Name>id</Name>
<ElementPath>id</ElementPath>
<Type>Integer</Type>
</PropertyDefn>
<PropertyDefn>
<Name>id2</Name>
<ElementPath>id2</ElementPath>
<Type>Integer</Type>
</PropertyDefn>
</GMLFeatureClass>
</GMLFeatureClassList>
@@ -0,0 +1,35 @@
<?xml version="1.0" encoding="utf-8" ?>
<ogr:FeatureCollection
xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
xsi:schemaLocation=""
xmlns:ogr="http://ogr.maptools.org/"
xmlns:gml="http://www.opengis.net/gml">
<gml:boundedBy>
<gml:Box>
<gml:coord><gml:X>1</gml:X><gml:Y>0</gml:Y></gml:coord>
<gml:coord><gml:X>4</gml:X><gml:Y>2</gml:Y></gml:coord>
</gml:Box>
</gml:boundedBy>

<gml:featureMember>
<ogr:displaced_points fid="points.4">
<ogr:geometryProperty><gml:Point srsName="EPSG:4326"><gml:coordinates>4,2</gml:coordinates></gml:Point></ogr:geometryProperty>
<ogr:id>5</ogr:id>
<ogr:id2>1</ogr:id2>
</ogr:displaced_points>
</gml:featureMember>
<gml:featureMember>
<ogr:displaced_points fid="points.3">
<ogr:geometryProperty><gml:Point srsName="EPSG:4326"><gml:coordinates>4,0</gml:coordinates></gml:Point></ogr:geometryProperty>
<ogr:id>4</ogr:id>
<ogr:id2>2</ogr:id2>
</ogr:displaced_points>
</gml:featureMember>
<gml:featureMember>
<ogr:displaced_points fid="points.0">
<ogr:geometryProperty><gml:Point srsName="EPSG:4326"><gml:coordinates>1,1</gml:coordinates></gml:Point></ogr:geometryProperty>
<ogr:id>1</ogr:id>
<ogr:id2>2</ogr:id2>
</ogr:displaced_points>
</gml:featureMember>
</ogr:FeatureCollection>
@@ -2905,3 +2905,17 @@ tests:
- 'Mean value: 147.17197994967066'
- 'Standard deviation: 43.9618116337985'
- 'Sum of the squares: 252304334.52061242'

- algorithm: qgis:pointsdisplacement
name: Point displacement
params:
DISTANCE: 1.0
HORIZONTAL: false
INPUT:
name: custom/displace_points.gml
type: vector
PROXIMITY: 2.0
results:
OUTPUT:
name: expected/displaced_points.gml
type: vector

0 comments on commit c8ac784

Please sign in to comment.
You can’t perform that action at this time.