1 /***************************************************************************
2                          qgstriangle.cpp
3                          -------------------
4     begin                : January 2017
5     copyright            : (C) 2017 by Loïc Bartoletti
6     email                : lbartoletti at tuxfamily dot org
7  ***************************************************************************/
8 
9 /***************************************************************************
10  *                                                                         *
11  *   This program is free software; you can redistribute it and/or modify  *
12  *   it under the terms of the GNU General Public License as published by  *
13  *   the Free Software Foundation; either version 2 of the License, or     *
14  *   (at your option) any later version.                                   *
15  *                                                                         *
16  ***************************************************************************/
17 
18 #include "qgstriangle.h"
19 #include "qgsgeometryutils.h"
20 #include "qgslinestring.h"
21 #include "qgswkbptr.h"
22 
23 #include <memory>
24 
QgsTriangle()25 QgsTriangle::QgsTriangle()
26 {
27   mWkbType = QgsWkbTypes::Triangle;
28 }
29 
QgsTriangle(const QgsPoint & p1,const QgsPoint & p2,const QgsPoint & p3)30 QgsTriangle::QgsTriangle( const QgsPoint &p1, const QgsPoint &p2, const QgsPoint &p3 )
31 {
32   mWkbType = QgsWkbTypes::Triangle;
33 
34   const QVector< double > x { p1.x(), p2.x(), p3.x(), p1.x() };
35   const QVector< double > y { p1.y(), p2.y(), p3.y(), p1.y() };
36   QVector< double > z;
37   if ( p1.is3D() )
38   {
39     z = { p1.z(), p2.z(), p3.z(), p1.z() };
40   }
41   QVector< double > m;
42   if ( p1.isMeasure() )
43   {
44     m = {p1.m(), p2.m(), p3.m(), p1.m() };
45   }
46   setExteriorRing( new QgsLineString( x, y, z, m ) );
47 }
48 
QgsTriangle(const QgsPointXY & p1,const QgsPointXY & p2,const QgsPointXY & p3)49 QgsTriangle::QgsTriangle( const QgsPointXY &p1, const QgsPointXY &p2, const QgsPointXY &p3 )
50 {
51   mWkbType = QgsWkbTypes::Triangle;
52 
53   const QVector< double > x { p1.x(), p2.x(), p3.x(), p1.x() };
54   const QVector< double > y {p1.y(), p2.y(), p3.y(), p1.y() };
55   QgsLineString *ext = new QgsLineString( x, y );
56   setExteriorRing( ext );
57 }
58 
QgsTriangle(const QPointF p1,const QPointF p2,const QPointF p3)59 QgsTriangle::QgsTriangle( const QPointF p1, const QPointF p2, const QPointF p3 )
60 {
61   mWkbType = QgsWkbTypes::Triangle;
62 
63   const QVector< double > x{ p1.x(), p2.x(), p3.x(), p1.x() };
64   const QVector< double > y{ p1.y(), p2.y(), p3.y(), p1.y() };
65   QgsLineString *ext = new QgsLineString( x, y );
66   setExteriorRing( ext );
67 }
68 
operator ==(const QgsTriangle & other) const69 bool QgsTriangle::operator==( const QgsTriangle &other ) const
70 {
71   if ( isEmpty() && other.isEmpty() )
72   {
73     return true;
74   }
75   else if ( isEmpty() || other.isEmpty() )
76   {
77     return false;
78   }
79 
80   return ( ( vertexAt( 0 ) == other.vertexAt( 0 ) ) &&
81            ( vertexAt( 1 ) == other.vertexAt( 1 ) ) &&
82            ( vertexAt( 2 ) == other.vertexAt( 2 ) )
83          );
84 }
85 
operator !=(const QgsTriangle & other) const86 bool QgsTriangle::operator!=( const QgsTriangle &other ) const
87 {
88   return !operator==( other );
89 }
90 
geometryType() const91 QString QgsTriangle::geometryType() const
92 {
93   return QStringLiteral( "Triangle" );
94 }
95 
createEmptyWithSameType() const96 QgsTriangle *QgsTriangle::createEmptyWithSameType() const
97 {
98   auto result = std::make_unique< QgsTriangle >();
99   result->mWkbType = mWkbType;
100   return result.release();
101 }
102 
clear()103 void QgsTriangle::clear()
104 {
105   QgsPolygon::clear();
106   mWkbType = QgsWkbTypes::Triangle;
107 }
108 
clone() const109 QgsTriangle *QgsTriangle::clone() const
110 {
111   return new QgsTriangle( *this );
112 }
113 
fromWkb(QgsConstWkbPtr & wkbPtr)114 bool QgsTriangle::fromWkb( QgsConstWkbPtr &wkbPtr )
115 {
116   clear();
117   if ( !wkbPtr )
118   {
119     return false;
120   }
121 
122   const QgsWkbTypes::Type type = wkbPtr.readHeader();
123   if ( QgsWkbTypes::flatType( type ) != QgsWkbTypes::Triangle )
124   {
125     return false;
126   }
127   mWkbType = type;
128 
129   QgsWkbTypes::Type ringType;
130   switch ( mWkbType )
131   {
132     case QgsWkbTypes::TriangleZ:
133       ringType = QgsWkbTypes::LineStringZ;
134       break;
135     case QgsWkbTypes::TriangleM:
136       ringType = QgsWkbTypes::LineStringM;
137       break;
138     case QgsWkbTypes::TriangleZM:
139       ringType = QgsWkbTypes::LineStringZM;
140       break;
141     default:
142       ringType = QgsWkbTypes::LineString;
143       break;
144   }
145 
146   int nRings;
147   wkbPtr >> nRings;
148   if ( nRings > 1 )
149   {
150     return false;
151   }
152 
153   QgsLineString *line = new QgsLineString();
154   line->fromWkbPoints( ringType, wkbPtr );
155   mExteriorRing.reset( line );
156 
157   return true;
158 }
159 
fromWkt(const QString & wkt)160 bool QgsTriangle::fromWkt( const QString &wkt )
161 {
162   clear();
163 
164   const QPair<QgsWkbTypes::Type, QString> parts = QgsGeometryUtils::wktReadBlock( wkt );
165 
166   if ( QgsWkbTypes::geometryType( parts.first ) != QgsWkbTypes::PolygonGeometry )
167     return false;
168 
169   mWkbType = parts.first;
170 
171   QString secondWithoutParentheses = parts.second;
172   secondWithoutParentheses = secondWithoutParentheses.simplified().remove( ' ' );
173   if ( ( parts.second.compare( QLatin1String( "EMPTY" ), Qt::CaseInsensitive ) == 0 ) ||
174        secondWithoutParentheses.isEmpty() )
175     return true;
176 
177   const QString defaultChildWkbType = QStringLiteral( "LineString%1%2" ).arg( is3D() ? QStringLiteral( "Z" ) : QString(), isMeasure() ? QStringLiteral( "M" ) : QString() );
178 
179   const QStringList blocks = QgsGeometryUtils::wktGetChildBlocks( parts.second, defaultChildWkbType );
180   for ( const QString &childWkt : blocks )
181   {
182     const QPair<QgsWkbTypes::Type, QString> childParts = QgsGeometryUtils::wktReadBlock( childWkt );
183 
184     const QgsWkbTypes::Type flatCurveType = QgsWkbTypes::flatType( childParts.first );
185 
186     if ( flatCurveType == QgsWkbTypes::LineString )
187       mInteriorRings.append( new QgsLineString() );
188     else
189     {
190       clear();
191       return false;
192     }
193     if ( !mInteriorRings.back()->fromWkt( childWkt ) )
194     {
195       clear();
196       return false;
197     }
198   }
199 
200   mExteriorRing.reset( mInteriorRings.takeFirst() );
201   if ( ( mExteriorRing->numPoints() < 3 ) || ( mExteriorRing->numPoints() > 4 ) || ( mExteriorRing->numPoints() == 4 && mExteriorRing->startPoint() != mExteriorRing->endPoint() ) )
202   {
203     clear();
204     return false;
205   }
206 
207   //scan through rings and check if dimensionality of rings is different to CurvePolygon.
208   //if so, update the type dimensionality of the CurvePolygon to match
209   bool hasZ = false;
210   bool hasM = false;
211   if ( mExteriorRing )
212   {
213     hasZ = hasZ || mExteriorRing->is3D();
214     hasM = hasM || mExteriorRing->isMeasure();
215   }
216   if ( hasZ )
217     addZValue( 0 );
218   if ( hasM )
219     addMValue( 0 );
220 
221   return true;
222 }
223 
asGml3(QDomDocument & doc,int precision,const QString & ns,const AxisOrder axisOrder) const224 QDomElement QgsTriangle::asGml3( QDomDocument &doc, int precision, const QString &ns, const AxisOrder axisOrder ) const
225 {
226 
227   QDomElement elemTriangle = doc.createElementNS( ns, QStringLiteral( "Triangle" ) );
228 
229   if ( isEmpty() )
230     return elemTriangle;
231 
232   QDomElement elemExterior = doc.createElementNS( ns, QStringLiteral( "exterior" ) );
233   QDomElement curveElem = exteriorRing()->asGml3( doc, precision, ns, axisOrder );
234   if ( curveElem.tagName() == QLatin1String( "LineString" ) )
235   {
236     curveElem.setTagName( QStringLiteral( "LinearRing" ) );
237   }
238   elemExterior.appendChild( curveElem );
239   elemTriangle.appendChild( elemExterior );
240 
241   return elemTriangle;
242 }
243 
surfaceToPolygon() const244 QgsPolygon *QgsTriangle::surfaceToPolygon() const
245 {
246   return toPolygon();
247 }
248 
toCurveType() const249 QgsCurvePolygon *QgsTriangle::toCurveType() const
250 {
251   std::unique_ptr<QgsCurvePolygon> curvePolygon( new QgsCurvePolygon() );
252   curvePolygon->setExteriorRing( mExteriorRing->clone() );
253 
254   return curvePolygon.release();
255 }
256 
addInteriorRing(QgsCurve * ring)257 void QgsTriangle::addInteriorRing( QgsCurve *ring )
258 {
259   delete ring;
260 }
261 
deleteVertex(QgsVertexId position)262 bool QgsTriangle::deleteVertex( QgsVertexId position )
263 {
264   Q_UNUSED( position )
265   return false;
266 }
267 
insertVertex(QgsVertexId position,const QgsPoint & vertex)268 bool QgsTriangle::insertVertex( QgsVertexId position, const QgsPoint &vertex )
269 {
270   Q_UNUSED( position )
271   Q_UNUSED( vertex )
272   return false;
273 }
274 
moveVertex(QgsVertexId vId,const QgsPoint & newPos)275 bool QgsTriangle::moveVertex( QgsVertexId vId, const QgsPoint &newPos )
276 {
277   if ( isEmpty() )
278     return false;
279 
280   if ( !mExteriorRing || vId.part != 0 || vId.ring != 0 || vId.vertex < 0 || vId.vertex > 4 )
281   {
282     return false;
283   }
284 
285   if ( vId.vertex == 4 )
286   {
287     vId.vertex = 0;
288   }
289 
290   const int n = mExteriorRing->numPoints();
291   const bool success = mExteriorRing->moveVertex( vId, newPos );
292   if ( success )
293   {
294     // If first or last vertex is moved, also move the last/first vertex
295     if ( vId.vertex == 0 )
296       mExteriorRing->moveVertex( QgsVertexId( vId.part, vId.ring, n - 1 ), newPos );
297     clearCache();
298   }
299   return success;
300 }
301 
setExteriorRing(QgsCurve * ring)302 void QgsTriangle::setExteriorRing( QgsCurve *ring )
303 {
304   if ( !ring )
305   {
306     return;
307   }
308 
309   if ( ring->hasCurvedSegments() )
310   {
311     //need to segmentize ring as polygon does not support curves
312     QgsCurve *line = ring->segmentize();
313     delete ring;
314     ring = line;
315   }
316 
317   if ( ( ring->numPoints() > 4 ) || ( ring->numPoints() < 3 ) )
318   {
319     delete ring;
320     return;
321   }
322   else if ( ring->numPoints() == 4 )
323   {
324     if ( !ring->isClosed() )
325     {
326       delete ring;
327       return;
328     }
329   }
330   else if ( ring->numPoints() == 3 )
331   {
332     if ( ring->isClosed() )
333     {
334       delete ring;
335       return;
336     }
337     QgsLineString *lineString = static_cast< QgsLineString *>( ring );
338     if ( !lineString->isClosed() )
339     {
340       lineString->close();
341     }
342     ring = lineString;
343   }
344 
345   mExteriorRing.reset( ring );
346 
347   //set proper wkb type
348   setZMTypeFromSubGeometry( ring, QgsWkbTypes::Triangle );
349 
350   clearCache();
351 }
352 
boundary() const353 QgsCurve *QgsTriangle::boundary() const
354 {
355   if ( !mExteriorRing )
356     return nullptr;
357 
358   return mExteriorRing->clone();
359 }
360 
vertexAt(int atVertex) const361 QgsPoint QgsTriangle::vertexAt( int atVertex ) const
362 {
363   if ( isEmpty() )
364     return QgsPoint();
365 
366   const QgsVertexId id( 0, 0, atVertex );
367   return mExteriorRing->vertexAt( id );
368 }
369 
lengths() const370 QVector<double> QgsTriangle::lengths() const
371 {
372   QVector<double> lengths;
373   if ( isEmpty() )
374     return lengths;
375 
376   lengths.append( vertexAt( 0 ).distance( vertexAt( 1 ) ) ); // c = |AB|
377   lengths.append( vertexAt( 1 ).distance( vertexAt( 2 ) ) ); // a = |BC|
378   lengths.append( vertexAt( 0 ).distance( vertexAt( 2 ) ) ); // b = |AC|
379 
380   return lengths;
381 }
382 
angles() const383 QVector<double> QgsTriangle::angles() const
384 {
385   QVector<double> angles;
386   if ( isEmpty() )
387     return angles;
388 
389   QVector<double> l = lengths();
390 
391   const double a = l[1];
392   const double b = l[2];
393   const double c = l[0];
394 
395   const double a2 = a * a;
396   const double b2 = b * b;
397   const double c2 = c * c;
398 
399   const double alpha = acos( ( b2 + c2 - a2 ) / ( 2 * b * c ) );
400   const double beta = acos( ( a2 + c2 - b2 ) / ( 2 * a * c ) );
401   const double gamma = M_PI - alpha - beta; // acos((a2 + b2 - c2)/(2*a*b)); but ensure that alpha+beta+gamma = 180.0
402 
403   angles.append( alpha );
404   angles.append( beta );
405   angles.append( gamma );
406 
407   return angles;
408 }
409 
isDegenerate()410 bool QgsTriangle::isDegenerate()
411 {
412   if ( isEmpty() )
413     return true;
414 
415   const QgsPoint p1( vertexAt( 0 ) );
416   const QgsPoint p2( vertexAt( 1 ) );
417   const QgsPoint p3( vertexAt( 2 ) );
418   return ( ( ( p1 == p2 ) || ( p1 == p3 ) || ( p2 == p3 ) ) || QgsGeometryUtils::leftOfLine( p1.x(), p1.y(), p2.x(), p2.y(), p3.x(), p3.y() ) == 0 );
419 }
420 
isIsocele(double lengthTolerance) const421 bool QgsTriangle::isIsocele( double lengthTolerance ) const
422 {
423   if ( isEmpty() )
424     return false;
425   const QVector<double> sides = lengths();
426   const bool ab_bc = qgsDoubleNear( sides.at( 0 ), sides.at( 1 ), lengthTolerance );
427   const bool bc_ca = qgsDoubleNear( sides.at( 1 ), sides.at( 2 ), lengthTolerance );
428   const bool ca_ab = qgsDoubleNear( sides.at( 2 ), sides.at( 0 ), lengthTolerance );
429 
430   return ( ab_bc || bc_ca || ca_ab );
431 }
432 
isEquilateral(double lengthTolerance) const433 bool QgsTriangle::isEquilateral( double lengthTolerance ) const
434 {
435   if ( isEmpty() )
436     return false;
437   const QVector<double> sides = lengths();
438   const bool ab_bc = qgsDoubleNear( sides.at( 0 ), sides.at( 1 ), lengthTolerance );
439   const bool bc_ca = qgsDoubleNear( sides.at( 1 ), sides.at( 2 ), lengthTolerance );
440   const bool ca_ab = qgsDoubleNear( sides.at( 2 ), sides.at( 0 ), lengthTolerance );
441 
442   return ( ab_bc && bc_ca && ca_ab );
443 }
444 
isRight(double angleTolerance) const445 bool QgsTriangle::isRight( double angleTolerance ) const
446 {
447   if ( isEmpty() )
448     return false;
449   QVector<double> a = angles();
450   QVector<double>::iterator ita = a.begin();
451   while ( ita != a.end() )
452   {
453     if ( qgsDoubleNear( *ita, M_PI_2, angleTolerance ) )
454       return true;
455     ++ita;
456   }
457   return false;
458 }
459 
isScalene(double lengthTolerance) const460 bool QgsTriangle::isScalene( double lengthTolerance ) const
461 {
462   if ( isEmpty() )
463     return false;
464   return !isIsocele( lengthTolerance );
465 }
466 
altitudes() const467 QVector<QgsLineString> QgsTriangle::altitudes() const
468 {
469   QVector<QgsLineString> alt;
470   if ( isEmpty() )
471     return alt;
472 
473   alt.append( QgsGeometryUtils::perpendicularSegment( vertexAt( 0 ), vertexAt( 2 ), vertexAt( 1 ) ) );
474   alt.append( QgsGeometryUtils::perpendicularSegment( vertexAt( 1 ), vertexAt( 0 ), vertexAt( 2 ) ) );
475   alt.append( QgsGeometryUtils::perpendicularSegment( vertexAt( 2 ), vertexAt( 0 ), vertexAt( 1 ) ) );
476 
477   return alt;
478 }
479 
medians() const480 QVector<QgsLineString> QgsTriangle::medians() const
481 {
482   QVector<QgsLineString> med;
483   if ( isEmpty() )
484     return med;
485 
486   QgsLineString med1;
487   QgsLineString med2;
488   QgsLineString med3;
489   med1.setPoints( QgsPointSequence() << vertexAt( 0 ) << QgsGeometryUtils::midpoint( vertexAt( 1 ), vertexAt( 2 ) ) );
490   med2.setPoints( QgsPointSequence() << vertexAt( 1 ) << QgsGeometryUtils::midpoint( vertexAt( 0 ), vertexAt( 2 ) ) );
491   med3.setPoints( QgsPointSequence() << vertexAt( 2 ) << QgsGeometryUtils::midpoint( vertexAt( 0 ), vertexAt( 1 ) ) );
492   med.append( med1 );
493   med.append( med2 );
494   med.append( med3 );
495 
496   return med;
497 }
498 
bisectors(double lengthTolerance) const499 QVector<QgsLineString> QgsTriangle::bisectors( double lengthTolerance ) const
500 {
501   QVector<QgsLineString> bis;
502   if ( isEmpty() )
503     return bis;
504 
505   QgsLineString bis1;
506   QgsLineString bis2;
507   QgsLineString bis3;
508   const QgsPoint incenter = inscribedCenter();
509   QgsPoint out;
510   bool intersection = false;
511 
512   QgsGeometryUtils::segmentIntersection( vertexAt( 0 ), incenter, vertexAt( 1 ), vertexAt( 2 ), out, intersection, lengthTolerance );
513   bis1.setPoints( QgsPointSequence() <<  vertexAt( 0 ) << out );
514 
515   QgsGeometryUtils::segmentIntersection( vertexAt( 1 ), incenter, vertexAt( 0 ), vertexAt( 2 ), out, intersection, lengthTolerance );
516   bis2.setPoints( QgsPointSequence() <<  vertexAt( 1 ) << out );
517 
518   QgsGeometryUtils::segmentIntersection( vertexAt( 2 ), incenter, vertexAt( 0 ), vertexAt( 1 ), out, intersection, lengthTolerance );
519   bis3.setPoints( QgsPointSequence() <<  vertexAt( 2 ) << out );
520 
521   bis.append( bis1 );
522   bis.append( bis2 );
523   bis.append( bis3 );
524 
525   return bis;
526 }
527 
medial() const528 QgsTriangle QgsTriangle::medial() const
529 {
530   if ( isEmpty() )
531     return QgsTriangle();
532   QgsPoint p1, p2, p3;
533   p1 = QgsGeometryUtils::midpoint( vertexAt( 0 ), vertexAt( 1 ) );
534   p2 = QgsGeometryUtils::midpoint( vertexAt( 1 ), vertexAt( 2 ) );
535   p3 = QgsGeometryUtils::midpoint( vertexAt( 2 ), vertexAt( 0 ) );
536   return QgsTriangle( p1, p2, p3 );
537 }
538 
orthocenter(double lengthTolerance) const539 QgsPoint QgsTriangle::orthocenter( double lengthTolerance ) const
540 {
541   if ( isEmpty() )
542     return QgsPoint();
543   const QVector<QgsLineString> alt = altitudes();
544   QgsPoint ortho;
545   bool intersection;
546   QgsGeometryUtils::segmentIntersection( alt.at( 0 ).pointN( 0 ), alt.at( 0 ).pointN( 1 ), alt.at( 1 ).pointN( 0 ), alt.at( 1 ).pointN( 1 ), ortho, intersection, lengthTolerance );
547 
548   return ortho;
549 }
550 
circumscribedCenter() const551 QgsPoint QgsTriangle::circumscribedCenter() const
552 {
553   if ( isEmpty() )
554     return QgsPoint();
555   double r, x, y;
556   QgsGeometryUtils::circleCenterRadius( vertexAt( 0 ), vertexAt( 1 ), vertexAt( 2 ), r, x, y );
557   return QgsPoint( x, y );
558 }
559 
circumscribedRadius() const560 double QgsTriangle::circumscribedRadius() const
561 {
562   if ( isEmpty() )
563     return 0.0;
564   double r, x, y;
565   QgsGeometryUtils::circleCenterRadius( vertexAt( 0 ), vertexAt( 1 ), vertexAt( 2 ), r, x, y );
566   return r;
567 }
568 
circumscribedCircle() const569 QgsCircle QgsTriangle::circumscribedCircle() const
570 {
571   if ( isEmpty() )
572     return QgsCircle();
573   return QgsCircle( circumscribedCenter(), circumscribedRadius() );
574 }
575 
inscribedCenter() const576 QgsPoint QgsTriangle::inscribedCenter() const
577 {
578   if ( isEmpty() )
579     return QgsPoint();
580 
581   const QVector<double> l = lengths();
582   const double x = ( l.at( 0 ) * vertexAt( 2 ).x() +
583                      l.at( 1 ) * vertexAt( 0 ).x() +
584                      l.at( 2 ) * vertexAt( 1 ).x() ) / perimeter();
585   const double y = ( l.at( 0 ) * vertexAt( 2 ).y() +
586                      l.at( 1 ) * vertexAt( 0 ).y() +
587                      l.at( 2 ) * vertexAt( 1 ).y() ) / perimeter();
588 
589   QgsPoint center( x, y );
590 
591   QgsPointSequence points;
592   points << vertexAt( 0 ) << vertexAt( 1 ) << vertexAt( 2 );
593   QgsGeometryUtils::transferFirstZOrMValueToPoint( points, center );
594 
595   return center;
596 }
597 
inscribedRadius() const598 double QgsTriangle::inscribedRadius() const
599 {
600   if ( isEmpty() )
601     return 0.0;
602   return ( 2.0 * area() / perimeter() );
603 }
604 
inscribedCircle() const605 QgsCircle QgsTriangle::inscribedCircle() const
606 {
607   if ( isEmpty() )
608     return QgsCircle();
609   return QgsCircle( inscribedCenter(), inscribedRadius() );
610 }
611