Skip to content

Commit 4d60d0c

Browse files
committedSep 2, 2016
[FEATURE] Add option to QgsGeometry::smooth to not smooth
segments shorter than a certain threshold or sharp corners with an angle exceeding a threshold Expose the angle threshold to processing smooth algorithm Also: - optimise QgsGeometry::smooth for new geometry classes - Fix smooth does not work with geometries containing Z/M
1 parent 01da222 commit 4d60d0c

File tree

15 files changed

+429
-86
lines changed

15 files changed

+429
-86
lines changed
 

‎doc/api_break.dox

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -581,6 +581,7 @@ method to determine if a geometry is valid.</li>
581581
<li>static bool compare( const QgsPolyline& p1, const QgsPolyline& p2, double epsilon ) has been renamed to comparePolylines</li>
582582
<li>static bool compare( const QgsPolygon& p1, const QgsPolygon& p2, double epsilon ) has been renamed to comparePolygons</li>
583583
<li>static bool compare( const QgsMultiPolygon& p1, const QgsMultiPolygon& p2, double epsilon ) has been renamed to compareMultiPolygons</li>
584+
<li>smoothLine and smoothPolygon are no longer public API (use smooth() instead)</li>
584585
</ul>
585586

586587
\subsection qgis_api_break_3_0_QgsGeometryAnalyzer QgsGeometryAnalyzer

‎python/core/geometry/qgsgeometry.sip

Lines changed: 4 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -849,14 +849,12 @@ class QgsGeometry
849849
* @param offset fraction of line to create new vertices along, between 0 and 1.0
850850
* eg the default value of 0.25 will create new vertices 25% and 75% along each line segment
851851
* of the geometry for each iteration. Smaller values result in "tighter" smoothing.
852+
* @param minimumDistance minimum segment length to apply smoothing to
853+
* @param maxAngle maximum angle at node (0-180) at which smoothing will be applied
852854
* @note added in 2.9
853855
*/
854-
QgsGeometry smooth( const unsigned int iterations, const double offset ) const;
855-
856-
/** Smooths a polygon using the Chaikin algorithm*/
857-
QgsPolygon smoothPolygon( const QgsPolygon &polygon, const unsigned int iterations = 1, const double offset = 0.25 ) const;
858-
/** Smooths a polyline using the Chaikin algorithm*/
859-
QgsPolyline smoothLine( const QgsPolyline &polyline, const unsigned int iterations = 1, const double offset = 0.25 ) const;
856+
QgsGeometry smooth( const unsigned int iterations = 1, const double offset = 0.25,
857+
double minimumDistance = -1.0, double maxAngle = 180.0 ) const;
860858

861859
/** Creates and returns a new geometry engine
862860
*/

‎python/plugins/processing/algs/help/qgis.yaml

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -444,6 +444,15 @@ qgis:singlesidedbuffer: >
444444

445445
The mitre limit parameter is only applicable for mitre join styles, and controls the maximum distance from the buffer to use when creating a mitred join.
446446

447+
qgis:smoothgeometry: >
448+
This algorithm smooths the geometries in a line or polygon layer. It creates a new layer with the same features as the ones in the input layer, but with geometries containing a higher number of vertices and corners in the geometries smoothed out.
449+
450+
The iterations parameter dictates how many smoothing iterations will be applied to each geometry. A higher number of iterations results in smoother geometries with the cost of greater number of nodes in the geometries.
451+
452+
The offset parameter controls how "tightly" the smoothed geometries follow the original geometries. Smaller values results in a tighter fit, and larger values will create a looser fit.
453+
454+
The maximum angle parameter can be used to prevent smoothing of nodes with large angles. Any node where the angle of the segments to either side is larger then this will not be smoothed. Eg setting the maximum angle to 90 degrees or lower would preserve right angles in the geometry.
455+
447456
qgis:snappointstogrid: >
448457
This algorithm modifies the position of points in a vector layer, so they fall in the coordinates of a grid.
449458

‎python/plugins/processing/algs/qgis/Smooth.py

Lines changed: 13 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -38,6 +38,7 @@ class Smooth(GeoAlgorithm):
3838
INPUT_LAYER = 'INPUT_LAYER'
3939
OUTPUT_LAYER = 'OUTPUT_LAYER'
4040
ITERATIONS = 'ITERATIONS'
41+
MAX_ANGLE = 'MAX_ANGLE'
4142
OFFSET = 'OFFSET'
4243

4344
def defineCharacteristics(self):
@@ -50,37 +51,37 @@ def defineCharacteristics(self):
5051
self.tr('Iterations'), default=1, minValue=1, maxValue=10))
5152
self.addParameter(ParameterNumber(self.OFFSET,
5253
self.tr('Offset'), default=0.25, minValue=0.0, maxValue=0.5))
54+
self.addParameter(ParameterNumber(self.MAX_ANGLE,
55+
self.tr('Maximum node angle to smooth'), default=180.0, minValue=0.0, maxValue=180.0))
5356
self.addOutput(OutputVector(self.OUTPUT_LAYER, self.tr('Smoothed')))
5457

5558
def processAlgorithm(self, progress):
5659
layer = dataobjects.getObjectFromUri(
5760
self.getParameterValue(self.INPUT_LAYER))
5861
iterations = self.getParameterValue(self.ITERATIONS)
5962
offset = self.getParameterValue(self.OFFSET)
63+
max_angle = self.getParameterValue(self.MAX_ANGLE)
6064

6165
writer = self.getOutputFromName(
6266
self.OUTPUT_LAYER).getVectorWriter(
6367
layer.fields().toList(),
6468
layer.wkbType(),
6569
layer.crs())
6670

67-
outFeat = QgsFeature()
68-
6971
features = vector.features(layer)
7072
total = 100.0 / len(features)
7173

72-
for current, inFeat in enumerate(features):
73-
inGeom = inFeat.geometry()
74-
attrs = inFeat.attributes()
74+
for current, input_feature in enumerate(features):
75+
output_feature = input_feature
76+
if input_feature.geometry():
77+
output_geometry = input_feature.geometry().smooth(iterations, offset, -1, max_angle)
78+
if not output_geometry:
79+
raise GeoAlgorithmExecutionException(
80+
self.tr('Error smoothing geometry'))
7581

76-
outGeom = inGeom.smooth(iterations, offset)
77-
if not outGeom:
78-
raise GeoAlgorithmExecutionException(
79-
self.tr('Error smoothing geometry'))
82+
output_feature.setGeometry(output_geometry)
8083

81-
outFeat.setGeometry(outGeom)
82-
outFeat.setAttributes(attrs)
83-
writer.addFeature(outFeat)
84+
writer.addFeature(output_feature)
8485
progress.setPercentage(int(current * total))
8586

8687
del writer
Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,15 @@
1+
<GMLFeatureClassList>
2+
<GMLFeatureClass>
3+
<Name>smoothed_lines</Name>
4+
<ElementPath>smoothed_lines</ElementPath>
5+
<GeometryType>2</GeometryType>
6+
<SRSName>EPSG:4326</SRSName>
7+
<DatasetSpecificInfo>
8+
<FeatureCount>7</FeatureCount>
9+
<ExtentXMin>-1.00000</ExtentXMin>
10+
<ExtentXMax>11.00000</ExtentXMax>
11+
<ExtentYMin>-3.00000</ExtentYMin>
12+
<ExtentYMax>5.00000</ExtentYMax>
13+
</DatasetSpecificInfo>
14+
</GMLFeatureClass>
15+
</GMLFeatureClassList>
Lines changed: 48 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,48 @@
1+
<?xml version="1.0" encoding="utf-8" ?>
2+
<ogr:FeatureCollection
3+
xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
4+
xsi:schemaLocation=""
5+
xmlns:ogr="http://ogr.maptools.org/"
6+
xmlns:gml="http://www.opengis.net/gml">
7+
<gml:boundedBy>
8+
<gml:Box>
9+
<gml:coord><gml:X>-1</gml:X><gml:Y>-3</gml:Y></gml:coord>
10+
<gml:coord><gml:X>11</gml:X><gml:Y>5</gml:Y></gml:coord>
11+
</gml:Box>
12+
</gml:boundedBy>
13+
14+
<gml:featureMember>
15+
<ogr:smoothed_lines fid="lines.0">
16+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>6,2 8.25,2.0 9.0,2.25 9.0,2.75 9.5,3.5 11,5</gml:coordinates></gml:LineString></ogr:geometryProperty>
17+
</ogr:smoothed_lines>
18+
</gml:featureMember>
19+
<gml:featureMember>
20+
<ogr:smoothed_lines fid="lines.1">
21+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>-1,-1 1,-1</gml:coordinates></gml:LineString></ogr:geometryProperty>
22+
</ogr:smoothed_lines>
23+
</gml:featureMember>
24+
<gml:featureMember>
25+
<ogr:smoothed_lines fid="lines.2">
26+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>2,0 2.0,1.5 2.25,2.0 2.75,2.0 3.0,2.25 3,3</gml:coordinates></gml:LineString></ogr:geometryProperty>
27+
</ogr:smoothed_lines>
28+
</gml:featureMember>
29+
<gml:featureMember>
30+
<ogr:smoothed_lines fid="lines.3">
31+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>3,1 5,1</gml:coordinates></gml:LineString></ogr:geometryProperty>
32+
</ogr:smoothed_lines>
33+
</gml:featureMember>
34+
<gml:featureMember>
35+
<ogr:smoothed_lines fid="lines.4">
36+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>7,-3 10,-3</gml:coordinates></gml:LineString></ogr:geometryProperty>
37+
</ogr:smoothed_lines>
38+
</gml:featureMember>
39+
<gml:featureMember>
40+
<ogr:smoothed_lines fid="lines.5">
41+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>6,-3 10,1</gml:coordinates></gml:LineString></ogr:geometryProperty>
42+
</ogr:smoothed_lines>
43+
</gml:featureMember>
44+
<gml:featureMember>
45+
<ogr:smoothed_lines fid="lines.6">
46+
</ogr:smoothed_lines>
47+
</gml:featureMember>
48+
</ogr:FeatureCollection>
Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,15 @@
1+
<GMLFeatureClassList>
2+
<GMLFeatureClass>
3+
<Name>smoothed_lines_max_angle</Name>
4+
<ElementPath>smoothed_lines_max_angle</ElementPath>
5+
<GeometryType>2</GeometryType>
6+
<SRSName>EPSG:4326</SRSName>
7+
<DatasetSpecificInfo>
8+
<FeatureCount>7</FeatureCount>
9+
<ExtentXMin>-1.00000</ExtentXMin>
10+
<ExtentXMax>11.00000</ExtentXMax>
11+
<ExtentYMin>-3.00000</ExtentYMin>
12+
<ExtentYMax>5.00000</ExtentYMax>
13+
</DatasetSpecificInfo>
14+
</GMLFeatureClass>
15+
</GMLFeatureClassList>
Lines changed: 48 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,48 @@
1+
<?xml version="1.0" encoding="utf-8" ?>
2+
<ogr:FeatureCollection
3+
xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
4+
xsi:schemaLocation=""
5+
xmlns:ogr="http://ogr.maptools.org/"
6+
xmlns:gml="http://www.opengis.net/gml">
7+
<gml:boundedBy>
8+
<gml:Box>
9+
<gml:coord><gml:X>-1</gml:X><gml:Y>-3</gml:Y></gml:coord>
10+
<gml:coord><gml:X>11</gml:X><gml:Y>5</gml:Y></gml:coord>
11+
</gml:Box>
12+
</gml:boundedBy>
13+
14+
<gml:featureMember>
15+
<ogr:smoothed_lines_max_angle fid="lines.0">
16+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>6,2 9,2 9.0,2.75 9.5,3.5 11,5</gml:coordinates></gml:LineString></ogr:geometryProperty>
17+
</ogr:smoothed_lines_max_angle>
18+
</gml:featureMember>
19+
<gml:featureMember>
20+
<ogr:smoothed_lines_max_angle fid="lines.1">
21+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>-1,-1 1,-1</gml:coordinates></gml:LineString></ogr:geometryProperty>
22+
</ogr:smoothed_lines_max_angle>
23+
</gml:featureMember>
24+
<gml:featureMember>
25+
<ogr:smoothed_lines_max_angle fid="lines.2">
26+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>2,0 2,2 3,2 3,3</gml:coordinates></gml:LineString></ogr:geometryProperty>
27+
</ogr:smoothed_lines_max_angle>
28+
</gml:featureMember>
29+
<gml:featureMember>
30+
<ogr:smoothed_lines_max_angle fid="lines.3">
31+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>3,1 5,1</gml:coordinates></gml:LineString></ogr:geometryProperty>
32+
</ogr:smoothed_lines_max_angle>
33+
</gml:featureMember>
34+
<gml:featureMember>
35+
<ogr:smoothed_lines_max_angle fid="lines.4">
36+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>7,-3 10,-3</gml:coordinates></gml:LineString></ogr:geometryProperty>
37+
</ogr:smoothed_lines_max_angle>
38+
</gml:featureMember>
39+
<gml:featureMember>
40+
<ogr:smoothed_lines_max_angle fid="lines.5">
41+
<ogr:geometryProperty><gml:LineString srsName="EPSG:4326"><gml:coordinates>6,-3 10,1</gml:coordinates></gml:LineString></ogr:geometryProperty>
42+
</ogr:smoothed_lines_max_angle>
43+
</gml:featureMember>
44+
<gml:featureMember>
45+
<ogr:smoothed_lines_max_angle fid="lines.6">
46+
</ogr:smoothed_lines_max_angle>
47+
</gml:featureMember>
48+
</ogr:FeatureCollection>

‎python/plugins/processing/tests/testdata/qgis_algorithm_tests.yaml

Lines changed: 28 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -965,3 +965,31 @@ tests:
965965
OUTPUT:
966966
name: expected/simplify_grid_lines.gml
967967
type: vector
968+
969+
- algorithm: qgis:smoothgeometry
970+
name: Smooth (lines)
971+
params:
972+
INPUT_LAYER:
973+
name: lines.gml
974+
type: vector
975+
ITERATIONS: 1
976+
MAX_ANGLE: 180.0
977+
OFFSET: 0.25
978+
results:
979+
OUTPUT_LAYER:
980+
name: expected/smoothed_lines.gml
981+
type: vector
982+
983+
- algorithm: qgis:smoothgeometry
984+
name: Smooth (lines, with max angle)
985+
params:
986+
INPUT_LAYER:
987+
name: lines.gml
988+
type: vector
989+
ITERATIONS: 1
990+
MAX_ANGLE: 60.0
991+
OFFSET: 0.25
992+
results:
993+
OUTPUT_LAYER:
994+
name: expected/smoothed_lines_max_angle.gml
995+
type: vector

‎src/core/geometry/qgsgeometry.cpp

Lines changed: 124 additions & 61 deletions
Original file line numberDiff line numberDiff line change
@@ -2138,55 +2138,56 @@ bool QgsGeometry::compare( const QgsMultiPolygon &p1, const QgsMultiPolygon &p2,
21382138
return true;
21392139
}
21402140

2141-
QgsGeometry QgsGeometry::smooth( const unsigned int iterations, const double offset ) const
2141+
QgsGeometry QgsGeometry::smooth( const unsigned int iterations, const double offset, double minimumDistance, double maxAngle ) const
21422142
{
2143-
switch ( wkbType() )
2143+
if ( d->geometry->isEmpty() )
2144+
return QgsGeometry();
2145+
2146+
QgsGeometry geom = *this;
2147+
if ( QgsWkbTypes::isCurvedType( wkbType() ) )
2148+
geom = QgsGeometry( d->geometry->segmentize() );
2149+
2150+
switch ( QgsWkbTypes::flatType( geom.wkbType() ) )
21442151
{
21452152
case QgsWkbTypes::Point:
2146-
case QgsWkbTypes::Point25D:
21472153
case QgsWkbTypes::MultiPoint:
2148-
case QgsWkbTypes::MultiPoint25D:
21492154
//can't smooth a point based geometry
2150-
return QgsGeometry( *this );
2155+
return geom;
21512156

21522157
case QgsWkbTypes::LineString:
2153-
case QgsWkbTypes::LineString25D:
21542158
{
2155-
QgsPolyline line = asPolyline();
2156-
return QgsGeometry::fromPolyline( smoothLine( line, iterations, offset ) );
2159+
QgsLineString* lineString = static_cast< QgsLineString* >( d->geometry );
2160+
return QgsGeometry( smoothLine( *lineString, iterations, offset, minimumDistance, maxAngle ) );
21572161
}
21582162

21592163
case QgsWkbTypes::MultiLineString:
2160-
case QgsWkbTypes::MultiLineString25D:
21612164
{
2162-
QgsMultiPolyline multiline = asMultiPolyline();
2163-
QgsMultiPolyline resultMultiline;
2164-
QgsMultiPolyline::const_iterator lineIt = multiline.constBegin();
2165-
for ( ; lineIt != multiline.constEnd(); ++lineIt )
2165+
QgsMultiLineString* multiLine = static_cast< QgsMultiLineString* >( d->geometry );
2166+
2167+
QgsMultiLineString* resultMultiline = new QgsMultiLineString();
2168+
for ( int i = 0; i < multiLine->numGeometries(); ++i )
21662169
{
2167-
resultMultiline << smoothLine( *lineIt, iterations, offset );
2170+
resultMultiline->addGeometry( smoothLine( *( static_cast< QgsLineString* >( multiLine->geometryN( i ) ) ), iterations, offset, minimumDistance, maxAngle ) );
21682171
}
2169-
return QgsGeometry::fromMultiPolyline( resultMultiline );
2172+
return QgsGeometry( resultMultiline );
21702173
}
21712174

21722175
case QgsWkbTypes::Polygon:
2173-
case QgsWkbTypes::Polygon25D:
21742176
{
2175-
QgsPolygon poly = asPolygon();
2176-
return QgsGeometry::fromPolygon( smoothPolygon( poly, iterations, offset ) );
2177+
QgsPolygonV2* poly = static_cast< QgsPolygonV2* >( d->geometry );
2178+
return QgsGeometry( smoothPolygon( *poly, iterations, offset, minimumDistance, maxAngle ) );
21772179
}
21782180

21792181
case QgsWkbTypes::MultiPolygon:
2180-
case QgsWkbTypes::MultiPolygon25D:
21812182
{
2182-
QgsMultiPolygon multipoly = asMultiPolygon();
2183-
QgsMultiPolygon resultMultipoly;
2184-
QgsMultiPolygon::const_iterator polyIt = multipoly.constBegin();
2185-
for ( ; polyIt != multipoly.constEnd(); ++polyIt )
2183+
QgsMultiPolygonV2* multiPoly = static_cast< QgsMultiPolygonV2* >( d->geometry );
2184+
2185+
QgsMultiPolygonV2* resultMultiPoly = new QgsMultiPolygonV2();
2186+
for ( int i = 0; i < multiPoly->numGeometries(); ++i )
21862187
{
2187-
resultMultipoly << smoothPolygon( *polyIt, iterations, offset );
2188+
resultMultiPoly->addGeometry( smoothPolygon( *( static_cast< QgsPolygonV2* >( multiPoly->geometryN( i ) ) ), iterations, offset, minimumDistance, maxAngle ) );
21882189
}
2189-
return QgsGeometry::fromMultiPolygon( resultMultipoly );
2190+
return QgsGeometry( resultMultiPoly );
21902191
}
21912192

21922193
case QgsWkbTypes::Unknown:
@@ -2195,58 +2196,120 @@ QgsGeometry QgsGeometry::smooth( const unsigned int iterations, const double off
21952196
}
21962197
}
21972198

2198-
inline QgsPoint interpolatePointOnLine( const QgsPoint& p1, const QgsPoint& p2, const double offset )
2199+
inline QgsPointV2 interpolatePointOnLine( const QgsPointV2& p1, const QgsPointV2& p2, const double offset )
21992200
{
22002201
double deltaX = p2.x() - p1.x();
22012202
double deltaY = p2.y() - p1.y();
2202-
return QgsPoint( p1.x() + deltaX * offset, p1.y() + deltaY * offset );
2203+
return QgsPointV2( p1.x() + deltaX * offset, p1.y() + deltaY * offset );
22032204
}
22042205

2205-
QgsPolyline QgsGeometry::smoothLine( const QgsPolyline& polyline, const unsigned int iterations, const double offset ) const
2206+
QgsLineString* smoothCurve( const QgsLineString& line, const unsigned int iterations,
2207+
const double offset, double squareDistThreshold, double maxAngleRads,
2208+
bool isRing )
22062209
{
2207-
QgsPolyline result = polyline;
2210+
QScopedPointer< QgsLineString > result( new QgsLineString( line ) );
22082211
for ( unsigned int iteration = 0; iteration < iterations; ++iteration )
22092212
{
2210-
QgsPolyline outputLine = QgsPolyline();
2211-
outputLine.reserve( 2 * ( result.count() - 1 ) );
2212-
for ( int i = 0; i < result.count() - 1; i++ )
2213+
QgsPointSequence outputLine;
2214+
outputLine.reserve( 2 * ( result->numPoints() - 1 ) );
2215+
bool skipFirst = false;
2216+
bool skipLast = false;
2217+
if ( isRing )
22132218
{
2214-
const QgsPoint& p1 = result.at( i );
2215-
const QgsPoint& p2 = result.at( i + 1 );
2216-
outputLine << ( i == 0 ? result.at( i ) : interpolatePointOnLine( p1, p2, offset ) );
2217-
outputLine << ( i == result.count() - 2 ? result.at( i + 1 ) : interpolatePointOnLine( p1, p2, 1.0 - offset ) );
2219+
QgsPointV2 p1 = result->pointN( result->numPoints() - 2 );
2220+
QgsPointV2 p2 = result->pointN( 0 );
2221+
QgsPointV2 p3 = result->pointN( 1 );
2222+
double angle = QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(),
2223+
p3.x(), p3.y() );
2224+
angle = qAbs( M_PI - angle );
2225+
skipFirst = angle > maxAngleRads;
22182226
}
2219-
result = outputLine;
2220-
}
2221-
return result;
2222-
}
2223-
2224-
QgsPolygon QgsGeometry::smoothPolygon( const QgsPolygon& polygon, const unsigned int iterations, const double offset ) const
2225-
{
2226-
QgsPolygon resultPoly;
2227-
QgsPolygon::const_iterator ringIt = polygon.constBegin();
2228-
for ( ; ringIt != polygon.constEnd(); ++ringIt )
2229-
{
2230-
QgsPolyline resultRing = *ringIt;
2231-
for ( unsigned int iteration = 0; iteration < iterations; ++iteration )
2227+
for ( int i = 0; i < result->numPoints() - 1; i++ )
22322228
{
2233-
QgsPolyline outputRing = QgsPolyline();
2234-
outputRing.reserve( 2 * ( resultRing.count() - 1 ) + 1 );
2235-
for ( int i = 0; i < resultRing.count() - 1; ++i )
2229+
QgsPointV2 p1 = result->pointN( i );
2230+
QgsPointV2 p2 = result->pointN( i + 1 );
2231+
2232+
double angle = M_PI;
2233+
if ( i == 0 && isRing )
2234+
{
2235+
QgsPointV2 p3 = result->pointN( result->numPoints() - 2 );
2236+
angle = QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(),
2237+
p3.x(), p3.y() );
2238+
}
2239+
else if ( i < result->numPoints() - 2 )
2240+
{
2241+
QgsPointV2 p3 = result->pointN( i + 2 );
2242+
angle = QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(),
2243+
p3.x(), p3.y() );
2244+
}
2245+
else if ( i == result->numPoints() - 2 && isRing )
22362246
{
2237-
const QgsPoint& p1 = resultRing.at( i );
2238-
const QgsPoint& p2 = resultRing.at( i + 1 );
2239-
outputRing << interpolatePointOnLine( p1, p2, offset );
2240-
outputRing << interpolatePointOnLine( p1, p2, 1.0 - offset );
2247+
QgsPointV2 p3 = result->pointN( 1 );
2248+
angle = QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(),
2249+
p3.x(), p3.y() );
22412250
}
2242-
//close polygon
2243-
outputRing << outputRing.at( 0 );
22442251

2245-
resultRing = outputRing;
2252+
skipLast = angle < M_PI - maxAngleRads || angle > M_PI + maxAngleRads;
2253+
2254+
// don't apply distance threshold to first or last segment
2255+
if ( i == 0 || i >= result->numPoints() - 2
2256+
|| QgsGeometryUtils::sqrDistance2D( p1, p2 ) > squareDistThreshold )
2257+
{
2258+
if ( !isRing )
2259+
{
2260+
if ( !skipFirst )
2261+
outputLine << ( i == 0 ? result->pointN( i ) : interpolatePointOnLine( p1, p2, offset ) );
2262+
if ( !skipLast )
2263+
outputLine << ( i == result->numPoints() - 2 ? result->pointN( i + 1 ) : interpolatePointOnLine( p1, p2, 1.0 - offset ) );
2264+
else
2265+
outputLine << p2;
2266+
}
2267+
else
2268+
{
2269+
// ring
2270+
if ( !skipFirst )
2271+
outputLine << interpolatePointOnLine( p1, p2, offset );
2272+
else if ( i == 0 )
2273+
outputLine << p1;
2274+
if ( !skipLast )
2275+
outputLine << interpolatePointOnLine( p1, p2, 1.0 - offset );
2276+
else
2277+
outputLine << p2;
2278+
}
2279+
}
2280+
skipFirst = skipLast;
22462281
}
2247-
resultPoly << resultRing;
2282+
2283+
if ( isRing && outputLine.at( 0 ) != outputLine.at( outputLine.count() - 1 ) )
2284+
outputLine << outputLine.at( 0 );
2285+
2286+
result->setPoints( outputLine );
2287+
}
2288+
return result.take();
2289+
}
2290+
2291+
QgsLineString* QgsGeometry::smoothLine( const QgsLineString& line, const unsigned int iterations, const double offset, double minimumDistance, double maxAngle ) const
2292+
{
2293+
double maxAngleRads = maxAngle * M_PI / 180.0;
2294+
double squareDistThreshold = minimumDistance > 0 ? minimumDistance * minimumDistance : -1;
2295+
return smoothCurve( line, iterations, offset, squareDistThreshold, maxAngleRads, false );
2296+
}
2297+
2298+
QgsPolygonV2* QgsGeometry::smoothPolygon( const QgsPolygonV2& polygon, const unsigned int iterations, const double offset, double minimumDistance, double maxAngle ) const
2299+
{
2300+
double maxAngleRads = maxAngle * M_PI / 180.0;
2301+
double squareDistThreshold = minimumDistance > 0 ? minimumDistance * minimumDistance : -1;
2302+
QScopedPointer< QgsPolygonV2 > resultPoly( new QgsPolygonV2 );
2303+
2304+
resultPoly->setExteriorRing( smoothCurve( *( static_cast< const QgsLineString*>( polygon.exteriorRing() ) ), iterations, offset,
2305+
squareDistThreshold, maxAngleRads, true ) );
2306+
2307+
for ( int i = 0; i < polygon.numInteriorRings(); ++i )
2308+
{
2309+
resultPoly->addInteriorRing( smoothCurve( *( static_cast< const QgsLineString*>( polygon.interiorRing( i ) ) ), iterations, offset,
2310+
squareDistThreshold, maxAngleRads, true ) );
22482311
}
2249-
return resultPoly;
2312+
return resultPoly.take();
22502313
}
22512314

22522315
QgsGeometry QgsGeometry::convertToPoint( bool destMultipart ) const

‎src/core/geometry/qgsgeometry.h

Lines changed: 33 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -43,6 +43,7 @@ class QgsVectorLayer;
4343
class QgsMapToPixel;
4444
class QPainter;
4545
class QgsPolygonV2;
46+
class QgsLineString;
4647

4748
/** Polyline is represented as a vector of points */
4849
typedef QVector<QgsPoint> QgsPolyline;
@@ -888,14 +889,12 @@ class CORE_EXPORT QgsGeometry
888889
* @param offset fraction of line to create new vertices along, between 0 and 1.0
889890
* eg the default value of 0.25 will create new vertices 25% and 75% along each line segment
890891
* of the geometry for each iteration. Smaller values result in "tighter" smoothing.
892+
* @param minimumDistance minimum segment length to apply smoothing to
893+
* @param maxAngle maximum angle at node (0-180) at which smoothing will be applied
891894
* @note added in 2.9
892895
*/
893-
QgsGeometry smooth( const unsigned int iterations = 1, const double offset = 0.25 ) const;
894-
895-
/** Smooths a polygon using the Chaikin algorithm*/
896-
QgsPolygon smoothPolygon( const QgsPolygon &polygon, const unsigned int iterations = 1, const double offset = 0.25 ) const;
897-
/** Smooths a polyline using the Chaikin algorithm*/
898-
QgsPolyline smoothLine( const QgsPolyline &polyline, const unsigned int iterations = 1, const double offset = 0.25 ) const;
896+
QgsGeometry smooth( const unsigned int iterations = 1, const double offset = 0.25,
897+
double minimumDistance = -1.0, double maxAngle = 180.0 ) const;
899898

900899
/** Creates and returns a new geometry engine
901900
*/
@@ -941,6 +940,34 @@ class CORE_EXPORT QgsGeometry
941940
QgsGeometry convertToLine( bool destMultipart ) const;
942941
/** Try to convert the geometry to a polygon */
943942
QgsGeometry convertToPolygon( bool destMultipart ) const;
943+
944+
/** Smooths a polyline using the Chaikin algorithm
945+
* @param line line to smooth
946+
* @param iterations number of smoothing iterations to run. More iterations results
947+
* in a smoother geometry
948+
* @param offset fraction of line to create new vertices along, between 0 and 1.0
949+
* eg the default value of 0.25 will create new vertices 25% and 75% along each line segment
950+
* of the geometry for each iteration. Smaller values result in "tighter" smoothing.
951+
* @param minimumDistance minimum segment length to apply smoothing to
952+
* @param maxAngle maximum angle at node (0-180) at which smoothing will be applied
953+
*/
954+
QgsLineString* smoothLine( const QgsLineString & line, const unsigned int iterations = 1, const double offset = 0.25,
955+
double minimumDistance = -1, double maxAngle = 180.0 ) const;
956+
957+
/** Smooths a polygon using the Chaikin algorithm
958+
* @param polygon polygon to smooth
959+
* @param iterations number of smoothing iterations to run. More iterations results
960+
* in a smoother geometry
961+
* @param offset fraction of segment to create new vertices along, between 0 and 1.0
962+
* eg the default value of 0.25 will create new vertices 25% and 75% along each line segment
963+
* of the geometry for each iteration. Smaller values result in "tighter" smoothing.
964+
* @param minimumDistance minimum segment length to apply smoothing to
965+
* @param maxAngle maximum angle at node (0-180) at which smoothing will be applied
966+
*/
967+
QgsPolygonV2* smoothPolygon( const QgsPolygonV2 &polygon, const unsigned int iterations = 1, const double offset = 0.25,
968+
double minimumDistance = -1, double maxAngle = 180.0 ) const;
969+
970+
944971
}; // class QgsGeometry
945972

946973
Q_DECLARE_METATYPE( QgsGeometry )

‎src/core/geometry/qgsgeometryutils.cpp

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -824,6 +824,13 @@ double QgsGeometryUtils::lineAngle( double x1, double y1, double x2, double y2 )
824824
return normalizedAngle( a );
825825
}
826826

827+
double QgsGeometryUtils::angleBetweenThreePoints( double x1, double y1, double x2, double y2, double x3, double y3 )
828+
{
829+
double angle1 = atan2( y1 - y2, x1 - x2 );
830+
double angle2 = atan2( y3 - y2, x3 - x2 );
831+
return normalizedAngle( angle1 - angle2 );
832+
}
833+
827834
double QgsGeometryUtils::linePerpendicularAngle( double x1, double y1, double x2, double y2 )
828835
{
829836
double a = lineAngle( x1, y1, x2, y2 );

‎src/core/geometry/qgsgeometryutils.h

Lines changed: 13 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -198,6 +198,19 @@ class CORE_EXPORT QgsGeometryUtils
198198
*/
199199
static double lineAngle( double x1, double y1, double x2, double y2 );
200200

201+
/** Calculates the angle between the lines AB and BC, where AB and BC described
202+
* by points a, b and b, c.
203+
* @param x1 x-coordinate of point a
204+
* @param y1 y-coordinate of point a
205+
* @param x2 x-coordinate of point b
206+
* @param y2 y-coordinate of point b
207+
* @param x3 x-coordinate of point c
208+
* @param y3 y-coordinate of point c
209+
* @returns angle between lines in radians. Returned value is undefined if two or more points are equal.
210+
*/
211+
static double angleBetweenThreePoints( double x1, double y1, double x2, double y2,
212+
double x3, double y3 );
213+
201214
/** Calculates the perpendicular angle to a line joining two points. Returned angle is in radians,
202215
* clockwise from the north direction.
203216
* @param x1 x-coordinate of line start

‎tests/src/core/testqgsgeometry.cpp

Lines changed: 50 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -3539,6 +3539,45 @@ void TestQgsGeometry::smoothCheck()
35393539
<< QgsPoint( 10.0, 7.5 ) << QgsPoint( 12.5, 10.0 ) << QgsPoint( 20.0, 10.0 );
35403540
QVERIFY( QgsGeometry::compare( line, expectedLine ) );
35413541

3542+
//linestring, with min distance
3543+
wkt = "LineString(0 0, 10 0, 10 10, 15 10, 15 20)";
3544+
geom = QgsGeometry::fromWkt( wkt );
3545+
result = geom.smooth( 1, 0.25, 6 );
3546+
line = result.asPolyline();
3547+
expectedLine.clear();
3548+
expectedLine << QgsPoint( 0, 0 ) << QgsPoint( 7.5, 0 ) << QgsPoint( 10.0, 2.5 )
3549+
<< QgsPoint( 10.0, 7.5 ) << QgsPoint( 15, 12.5 ) << QgsPoint( 15.0, 20.0 );
3550+
QVERIFY( QgsGeometry::compare( line, expectedLine ) );
3551+
3552+
//linestring, with max angle
3553+
wkt = "LineString(0 0, 10 0, 15 5, 25 -5, 30 -5 )";
3554+
geom = QgsGeometry::fromWkt( wkt );
3555+
result = geom.smooth( 1, 0.25, 0, 50 );
3556+
line = result.asPolyline();
3557+
expectedLine.clear();
3558+
expectedLine << QgsPoint( 0, 0 ) << QgsPoint( 7.5, 0 ) << QgsPoint( 11.25, 1.25 )
3559+
<< QgsPoint( 15.0, 5.0 ) << QgsPoint( 22.5, -2.5 ) << QgsPoint( 26.25, -5 ) << QgsPoint( 30, -5 );
3560+
QVERIFY( QgsGeometry::compare( line, expectedLine ) );
3561+
3562+
//linestring, with max angle, other direction
3563+
wkt = "LineString( 30 -5, 25 -5, 15 5, 10 0, 0 0 )";
3564+
geom = QgsGeometry::fromWkt( wkt );
3565+
result = geom.smooth( 1, 0.25, 0, 50 );
3566+
line = result.asPolyline();
3567+
expectedLine.clear();
3568+
expectedLine << QgsPoint( 30, -5 ) << QgsPoint( 26.25, -5 ) << QgsPoint( 22.5, -2.5 )
3569+
<< QgsPoint( 15.0, 5.0 ) << QgsPoint( 11.25, 1.25 ) << QgsPoint( 7.5, 0 ) << QgsPoint( 0, 0 );
3570+
QVERIFY( QgsGeometry::compare( line, expectedLine ) );
3571+
3572+
//linestring, max angle, first corner sharp
3573+
wkt = "LineString(0 0, 10 0, 10 10 )";
3574+
geom = QgsGeometry::fromWkt( wkt );
3575+
result = geom.smooth( 1, 0.25, 0, 50 );
3576+
line = result.asPolyline();
3577+
expectedLine.clear();
3578+
expectedLine << QgsPoint( 0, 0 ) << QgsPoint( 10, 0 ) << QgsPoint( 10, 10 );
3579+
QVERIFY( QgsGeometry::compare( line, expectedLine ) );
3580+
35423581
wkt = "MultiLineString ((0 0, 10 0, 10 10, 20 10),(30 30, 40 30, 40 40, 50 40))";
35433582
geom = QgsGeometry::fromWkt( wkt );
35443583
result = geom.smooth( 1, 0.25 );
@@ -3564,7 +3603,17 @@ void TestQgsGeometry::smoothCheck()
35643603
<< QgsPoint( 2.0, 3.5 ) << QgsPoint( 2.0, 2.5 ) << QgsPoint( 2.5, 2.0 ) );
35653604
QVERIFY( QgsGeometry::compare( poly, expectedPolygon ) );
35663605

3567-
//multipolygon
3606+
//polygon with max angle - should be unchanged
3607+
wkt = "Polygon ((0 0, 10 0, 10 10, 0 10, 0 0))";
3608+
geom = QgsGeometry::fromWkt( wkt );
3609+
result = geom.smooth( 1, 0.25, -1, 50 );
3610+
poly = result.asPolygon();
3611+
expectedPolygon.clear();
3612+
expectedPolygon << ( QgsPolyline() << QgsPoint( 0, 0 ) << QgsPoint( 10, 0 ) << QgsPoint( 10.0, 10 )
3613+
<< QgsPoint( 0, 10 ) << QgsPoint( 0, 0 ) );
3614+
QVERIFY( QgsGeometry::compare( poly, expectedPolygon ) );
3615+
3616+
//multipolygon)
35683617
wkt = "MultiPolygon (((0 0, 10 0, 10 10, 0 10, 0 0 )),((2 2, 4 2, 4 4, 2 4, 2 2)))";
35693618
geom = QgsGeometry::fromWkt( wkt );
35703619
result = geom.smooth( 1, 0.1 );

‎tests/src/core/testqgsgeometryutils.cpp

Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -50,6 +50,7 @@ class TestQgsGeometryUtils: public QObject
5050
void testCircleCenterRadius_data();
5151
void testCircleCenterRadius();
5252
void testSqrDistToLine();
53+
void testAngleThreePoints();
5354
};
5455

5556

@@ -518,6 +519,26 @@ void TestQgsGeometryUtils::testSqrDistToLine()
518519
QGSCOMPARENEAR( sqrDist, 11.83, 0.01 );
519520
}
520521

522+
void TestQgsGeometryUtils::testAngleThreePoints()
523+
{
524+
QgsPoint p1( 0, 0 );
525+
QgsPoint p2( 1, 0 );
526+
QgsPoint p3( 1, 1 );
527+
QGSCOMPARENEAR( QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(), p3.x(), p3.y() ), M_PI / 2.0, 0.00000001 );
528+
p3 = QgsPoint( 1, -1 );
529+
QGSCOMPARENEAR( QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(), p3.x(), p3.y() ), 3 * M_PI / 2.0, 0.00000001 );
530+
p3 = QgsPoint( 2, 0 );
531+
QGSCOMPARENEAR( QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(), p3.x(), p3.y() ), M_PI, 0.00000001 );
532+
p3 = QgsPoint( 0, 0 );
533+
QGSCOMPARENEAR( QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(), p3.x(), p3.y() ), 0.0, 0.00000001 );
534+
p3 = QgsPoint( 1, 0 );
535+
//undefined, but want no crash
536+
( void )QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(), p3.x(), p3.y() );
537+
p2 = QgsPoint( 0, 0 );
538+
( void )QgsGeometryUtils::angleBetweenThreePoints( p1.x(), p1.y(), p2.x(), p2.y(), p3.x(), p3.y() );
539+
}
540+
541+
521542

522543
QTEST_MAIN( TestQgsGeometryUtils )
523544
#include "testqgsgeometryutils.moc"

0 commit comments

Comments
 (0)
Please sign in to comment.