QGIS API Documentation 3.37.0-Master (fdefdf9c27f)
NormVecDecorator.cpp
Go to the documentation of this file.
1/***************************************************************************
2 NormVecDecorator.cpp
3 --------------------
4 copyright : (C) 2004 by Marco Hugentobler
6 ***************************************************************************/
7
8/***************************************************************************
9 * *
10 * This program is free software; you can redistribute it and/or modify *
11 * it under the terms of the GNU General Public License as published by *
12 * the Free Software Foundation; either version 2 of the License, or *
13 * (at your option) any later version. *
14 * *
15 ***************************************************************************/
16
17#include <QApplication>
18
19#include "qgsfeedback.h"
20#include "qgslogger.h"
21#include "qgsfields.h"
22#include "qgspoint.h"
23
24#include "NormVecDecorator.h"
25
26
28{
29 //remove all the normals
30 if ( !mNormVec->isEmpty() )
31 {
32 for ( int i = 0; i < mNormVec->count(); i++ )
33 {
34 delete ( *mNormVec )[i];
35 }
36 }
37
38 delete mNormVec;
39 delete mPointState;
40 delete mTIN;
41}
42
44{
45 if ( mTIN )
46 {
47 int pointno;
48 pointno = mTIN->addPoint( p );
49
50 if ( pointno == -100 )//a numerical error occurred
51 {
52// QgsDebugError("warning, numerical error");
53 return -100;
54 }
55
56 //recalculate the necessary normals
57 if ( alreadyestimated )
58 {
59 estimateFirstDerivative( pointno );
60 //update also the neighbours of the new point
61 const QList<int> list = mTIN->surroundingTriangles( pointno );
62 auto it = list.constBegin();//iterate through the list and analyze it
63 while ( it != list.constEnd() )
64 {
65 int point;
66 point = *it;
67 if ( point != -1 )
68 {
70 }
71 ++it;
72 ++it;
73 ++it;
74 ++it;
75 }
76 }
77 return pointno;
78 }
79
80 return -1;
81}
82
83bool NormVecDecorator::calcNormal( double x, double y, QgsPoint &result )
84{
85 if ( !alreadyestimated )
86 {
88 alreadyestimated = true;
89 }
90
91 if ( mInterpolator )
92 {
93 const bool b = mInterpolator->calcNormVec( x, y, result );
94 return b;
95 }
96 else
97 {
98 QgsDebugError( QStringLiteral( "warning, null pointer" ) );
99 return false;
100 }
101}
102
103bool NormVecDecorator::calcNormalForPoint( double x, double y, int pointIndex, Vector3D *result )
104{
105 if ( !alreadyestimated )
106 {
108 alreadyestimated = true;
109 }
110
111 if ( result )
112 {
113 int numberofbreaks = 0;//number of breaklines around the point
114 int ffirstbp = -1000; //numbers of the points related to the first breakline
115 int lfirstbp = -1000;
116 bool pointfound = false;//is set to true, if the triangle with the point in is found
117 int numberofruns = 0;//number of runs of the loop. This integer can be used to prevent endless loops
118 const int limit = 100000;//ater this number of iterations, the method is terminated
119
120 result->setX( 0 );
121 result->setY( 0 );
122 result->setZ( 0 );
123
124 const QList<int> vlist = surroundingTriangles( pointIndex );//get the value list
125 if ( vlist.empty() )//an error occurred in 'getSurroundingTriangles'
126 {
127 return false;
128 }
129
130 if ( ( ( vlist.count() ) % 4 ) != 0 ) //number of items in vlist has to be a multiple of 4
131 {
132 QgsDebugError( QStringLiteral( "warning, wrong number of items in vlist" ) );
133 return false;
134 }
135
136 auto it = vlist.constBegin();
137
138 bool firstrun;
139
140 while ( true )//endless loop to analyze vlist
141 {
142 numberofruns++;
143 if ( numberofruns > limit )
144 {
145 QgsDebugError( QStringLiteral( "warning, a probable endless loop is detected" ) );
146 return false;
147 }
148
149 int p1, p2, p3, line;
150 firstrun = false;//flag which tells, if it is the first run with a breakline
151 p1 = ( *it );
152 ++it;
153 p2 = ( *it );
154 ++it;
155 p3 = ( *it );
156 ++it;
157 line = ( *it );
158
159
160 if ( numberofbreaks > 0 )
161 {
162
163 if ( p1 != -1 && p2 != -1 && p3 != -1 )
164 {
165 if ( MathUtils::pointInsideTriangle( x, y, point( p1 ), point( p2 ), point( p3 ) ) )
166 {
167 pointfound = true;
168 }
169
170 Vector3D addvec( 0, 0, 0 );
171 MathUtils::normalFromPoints( point( p1 ), point( p2 ), point( p3 ), &addvec );
172 result->setX( result->getX() + addvec.getX() );
173 result->setY( result->getY() + addvec.getY() );
174 result->setZ( result->getZ() + addvec.getZ() );
175 }
176 }
177
178 if ( line == -10 )//we found a breakline
179 {
180
181 if ( numberofbreaks == 0 )//it is the first breakline
182 {
183 firstrun = true;
184 ffirstbp = p2;//set the marks to recognize the breakline later
185 lfirstbp = p3;
186 }
187
188 if ( p2 == ffirstbp && p3 == lfirstbp && !firstrun )//we are back at the first breakline
189 {
190 if ( !pointfound )//the point with coordinates x, y was in no triangle
191 {
192 QgsDebugError( QStringLiteral( "warning: point (x,y) was in no triangle" ) );
193 return false;
194 }
195 result->standardise();
196 break;
197 }
198
199 if ( numberofbreaks > 0 && pointfound )//we found the second break line and the point is between the first and the second
200 {
201 result->standardise();
202 numberofbreaks++;//to make the distinction between endpoints and points on a breakline easier
203 break;
204 }
205
206 result->setX( 0 );
207 result->setY( 0 );
208 result->setZ( 0 );
209 numberofbreaks++;
210 }
211
212 ++it;
213 if ( it == vlist.constEnd() )//restart at the beginning of the loop
214 {
215 it = vlist.constBegin();
216 }
217 }
218 return true;
219 }
220 else
221 {
222 QgsDebugError( QStringLiteral( "warning, null pointer" ) );
223 return false;
224 }
225
226}
227
228bool NormVecDecorator::calcPoint( double x, double y, QgsPoint &result )
229{
230
231 if ( !alreadyestimated )
232 {
234 alreadyestimated = true;
235 }
236
237 if ( mInterpolator )
238 {
239 const bool b = mInterpolator->calcPoint( x, y, result );
240 return b;
241 }
242 else
243 {
244 QgsDebugError( QStringLiteral( "warning, null pointer" ) );
245 return false;
246 }
247}
248
249bool NormVecDecorator::getTriangle( double x, double y, QgsPoint &p1, Vector3D *v1, QgsPoint &p2, Vector3D *v2, QgsPoint &p3, Vector3D *v3 )
250{
251 if ( v1 && v2 && v3 )
252 {
253 int nr1 = 0;
254 int nr2 = 0;
255 int nr3 = 0;
256
257 if ( TriDecorator::triangleVertices( x, y, p1, nr1, p2, nr2, p3, nr3 ) )//everything alright
258 {
259 if ( ( *mNormVec )[ nr1 ] && ( *mNormVec )[ nr2 ] && ( *mNormVec )[ nr3 ] )
260 {
261 v1->setX( ( *mNormVec )[ nr1 ]->getX() );
262 v1->setY( ( *mNormVec )[nr1 ]->getY() );
263 v1->setZ( ( *mNormVec )[nr1 ]->getZ() );
264
265 v2->setX( ( *mNormVec )[nr2 ]->getX() );
266 v2->setY( ( *mNormVec )[nr2 ]->getY() );
267 v2->setZ( ( *mNormVec )[nr2 ]->getZ() );
268
269 v3->setX( ( *mNormVec )[nr3 ]->getX() );
270 v3->setY( ( *mNormVec )[nr3 ]->getY() );
271 v3->setZ( ( *mNormVec )[nr3 ]->getZ() );
272 }
273 else
274 {
275 QgsDebugError( QStringLiteral( "warning, null pointer" ) );
276 return false;
277 }
278 return true;
279 }
280 else
281 {
282 return false;
283 }
284 }
285
286 else
287 {
288 QgsDebugError( QStringLiteral( "warning, null pointer" ) );
289 return false;
290 }
291}
292
294{
295 if ( pointno >= 0 )
296 {
297 return mPointState->at( pointno );
298 }
299 else
300 {
301 QgsDebugError( QStringLiteral( "warning, number below 0" ) );
302 return mPointState->at( 0 );//just to avoid a compiler warning
303 }
304}
305
306
307bool NormVecDecorator::getTriangle( double x, double y, QgsPoint &p1, int &ptn1, Vector3D *v1, PointState *state1, QgsPoint &p2, int &ptn2, Vector3D *v2, PointState *state2, QgsPoint &p3, int &ptn3, Vector3D *v3, PointState *state3 )
308{
309 if ( v1 && v2 && v3 && state1 && state2 && state3 )
310 {
311 if ( TriDecorator::triangleVertices( x, y, p1, ptn1, p2, ptn2, p3, ptn3 ) )//everything alright
312 {
313 v1->setX( ( *mNormVec )[( ptn1 )]->getX() );
314 v1->setY( ( *mNormVec )[( ptn1 )]->getY() );
315 v1->setZ( ( *mNormVec )[( ptn1 )]->getZ() );
316
317 ( *state1 ) = ( *mPointState )[( ptn1 )];
318
319 v2->setX( ( *mNormVec )[( ptn2 )]->getX() );
320 v2->setY( ( *mNormVec )[( ptn2 )]->getY() );
321 v2->setZ( ( *mNormVec )[( ptn2 )]->getZ() );
322
323 ( *state2 ) = ( *mPointState )[( ptn2 )];
324
325 v3->setX( ( *mNormVec )[( ptn3 )]->getX() );
326 v3->setY( ( *mNormVec )[( ptn3 )]->getY() );
327 v3->setZ( ( *mNormVec )[( ptn3 )]->getZ() );
328
329 ( *state3 ) = ( *mPointState )[( ptn3 )];
330
331 return true;
332 }
333 else
334 {
335 return false;
336 }
337
338 }
339 else
340 {
341 QgsDebugError( QStringLiteral( "warning, null pointer" ) );
342 return false;
343 }
344}
345
347{
348 if ( pointno == -1 )
349 {
350 return false;
351 }
352
353 Vector3D part;
354 Vector3D total;
355 total.setX( 0 );
356 total.setY( 0 );
357 total.setZ( 0 );
358 int numberofbreaks = 0;//number of counted breaklines
359 double weights = 0;//sum of the weights
360 double currentweight = 0;//current weight
361 PointState status;
362
363 const QList<int> vlist = surroundingTriangles( pointno );//get the value list
364
365 if ( vlist.empty() )
366 {
367 //something went wrong in getSurroundingTriangles, set the normal to (0,0,0)
368 if ( mNormVec->size() <= mNormVec->count() )//allocate more memory if necessary
369 {
370 QgsDebugMsgLevel( QStringLiteral( "resizing mNormVec from %1 to %2" ).arg( mNormVec->size() ).arg( mNormVec->size() + 1 ), 2 );
371 mNormVec->resize( mNormVec->size() + 1 );
372 }
373
374 //todo:resize mNormVec if necessary
375
376 if ( !( ( *mNormVec )[pointno] ) ) //insert a pointer to a Vector3D, if there is none at this position
377 {
378 Vector3D *vec = new Vector3D( total.getX(), total.getY(), total.getZ() );
379 mNormVec->insert( pointno, vec );
380 }
381 else
382 {
383 ( *mNormVec )[pointno]->setX( 0 );
384 ( *mNormVec )[pointno]->setY( 0 );
385 ( *mNormVec )[pointno]->setZ( 0 );
386 }
387 return false;
388 }
389
390 if ( ( vlist.count() % 4 ) != 0 ) //number of items in vlist has to be a multiple of 4
391 {
392 QgsDebugError( QStringLiteral( "warning, wrong number of items in vlist" ) );
393 return false;
394 }
395
396 auto it = vlist.constBegin();//iterate through the list and analyze it
397 while ( it != vlist.constEnd() )
398 {
399 int p1, p2, p3, flag;
400 part.setX( 0 );
401 part.setY( 0 );
402 part.setZ( 0 );
403
404 currentweight = 0;
405
406 p1 = ( *it );
407 ++it;
408 p2 = ( *it );
409 ++it;
410 p3 = ( *it );
411 ++it;
412 flag = ( *it );
413
414 if ( flag == -10 )//we found a breakline.
415 {
416 numberofbreaks++;
417 }
418
419 if ( p1 != -1 && p2 != -1 && p3 != -1 )//don't calculate normal, if a point is a virtual point
420 {
421 MathUtils::normalFromPoints( point( p1 ), point( p2 ), point( p3 ), &part );
422 const double dist1 = point( p3 )->distance3D( *point( p1 ) );
423 const double dist2 = point( p3 )->distance3D( *point( p2 ) );
424 //don't add the normal if the triangle is horizontal
425 if ( ( point( p1 )->z() != point( p2 )->z() ) || ( point( p1 )->z() != point( p3 )->z() ) )
426 {
427 currentweight = 1 / ( dist1 * dist1 * dist2 * dist2 );
428 total.setX( total.getX() + part.getX()*currentweight );
429 total.setY( total.getY() + part.getY()*currentweight );
430 total.setZ( total.getZ() + part.getZ()*currentweight );
431 weights += currentweight;
432 }
433 }
434 ++it;
435 }
436
437 if ( total.getX() == 0 && total.getY() == 0 && total.getZ() == 0 )//we have a point surrounded by horizontal triangles
438 {
439 total.setZ( 1 );
440 }
441 else
442 {
443 total.setX( total.getX() / weights );
444 total.setY( total.getY() / weights );
445 total.setZ( total.getZ() / weights );
446 total.standardise();
447 }
448
449
450 if ( numberofbreaks == 0 )
451 {
452 status = Normal;
453 }
454 else if ( numberofbreaks == 1 )
455 {
456 status = EndPoint;
457 }
458 else
459 {
460 status = BreakLine;
461 }
462
463 //insert the new calculated vector
464 if ( mNormVec->size() <= mNormVec->count() )//allocate more memory if necessary
465 {
466 mNormVec->resize( mNormVec->size() + 1 );
467 }
468
469 if ( !( ( *mNormVec )[pointno] ) ) //insert a pointer to a Vector3D, if there is none at this position
470 {
471 Vector3D *vec = new Vector3D( total.getX(), total.getY(), total.getZ() );
472 mNormVec->insert( pointno, vec );
473 }
474 else
475 {
476 ( *mNormVec )[pointno]->setX( total.getX() );
477 ( *mNormVec )[pointno]->setY( total.getY() );
478 ( *mNormVec )[pointno]->setZ( total.getZ() );
479 }
480
481 //insert the new status
482
483 if ( pointno >= mPointState->size() )
484 {
485 mPointState->resize( mPointState->size() + 1 );
486 }
487
488 ( *mPointState )[pointno] = status;
489
490 return true;
491}
492
493//weighted method of little
495{
496 const int numberPoints = pointsCount();
497 for ( int i = 0; i < numberPoints; i++ )
498 {
499 if ( feedback )
500 {
501 feedback->setProgress( 100.0 * static_cast< double >( i ) / numberPoints );
502 }
504 }
505
506 return true;
507}
508
510{
511 if ( mTIN )
512 {
513 if ( alreadyestimated )
514 {
517 }
518 else
519 {
521 }
522 }
523 else
524 {
525 QgsDebugError( QStringLiteral( "warning, null pointer" ) );
526 }
527}
528
530{
531 if ( pointno >= 0 )
532 {
533 ( *mPointState )[pointno] = s;
534 }
535 else
536 {
537 QgsDebugError( QStringLiteral( "warning, pointno>0" ) );
538 }
539}
540
541bool NormVecDecorator::swapEdge( double x, double y )
542{
543 if ( mTIN )
544 {
545 bool b = false;
546 if ( alreadyestimated )
547 {
548 const QList<int> list = pointsAroundEdge( x, y );
549 if ( !list.empty() )
550 {
551 b = mTIN->swapEdge( x, y );
552 for ( const int i : list )
553 {
555 }
556 }
557 }
558 else
559 {
560 b = mTIN->swapEdge( x, y );
561 }
562 return b;
563 }
564 else
565 {
566 QgsDebugError( QStringLiteral( "warning, null pointer" ) );
567 return false;
568 }
569}
570
572{
573 if ( !mTIN )
574 {
575 return false;
576 }
577 return mTIN->saveTriangulation( sink, feedback );
578}
579
581{
582 if ( !mTIN )
583 {
584 return QgsMesh();
585 }
586 return mTIN->triangulationToMesh( feedback );
587}
588
bool getTriangle(double x, double y, QgsPoint &p1, Vector3D *v1, QgsPoint &p2, Vector3D *v2, QgsPoint &p3, Vector3D *v3)
Finds out, in which triangle a point with coordinates x and y is and assigns the triangle points to p...
QVector< PointState > * mPointState
Vector who stores, it a point is not on a breakline, if it is a normal point of the breakline or if i...
bool alreadyestimated
Is true, if the normals already have been estimated.
bool calcNormal(double x, double y, QgsPoint &result) override
Calculates the normal at a point on the surface and assigns it to 'result'. Returns true in case of s...
bool saveTriangulation(QgsFeatureSink *sink, QgsFeedback *feedback=nullptr) const override
Saves the triangulation features to a feature sink.
~NormVecDecorator() override
QgsMesh triangulationToMesh(QgsFeedback *feedback=nullptr) const override
Returns a QgsMesh corresponding to the triangulation.
void eliminateHorizontalTriangles() override
Eliminates the horizontal triangles by swapping or by insertion of new points. If alreadyestimated is...
bool calcPoint(double x, double y, QgsPoint &result) override
Calculates x-, y and z-value of the point on the surface and assigns it to 'result'.
TriangleInterpolator * mInterpolator
Association with an interpolator object.
void setState(int pointno, PointState s)
Sets the state (BreakLine, Normal, EndPoint) of a point.
QVector< Vector3D * > * mNormVec
Vector that stores the normals for the points. If 'estimateFirstDerivatives()' was called and there i...
bool swapEdge(double x, double y) override
Swaps the edge which is closest to the point with x and y coordinates (if this is possible) and force...
PointState getState(int pointno) const
Returns the state of the point with the number 'pointno'.
bool calcNormalForPoint(double x, double y, int pointIndex, Vector3D *result)
Calculates the normal of a triangle-point for the point with coordinates x and y. This is needed,...
bool estimateFirstDerivatives(QgsFeedback *feedback=nullptr)
This method adds the functionality of estimating normals at the data points. Return true in the case ...
bool estimateFirstDerivative(int pointno)
Estimates the first derivative a point. Return true in case of success and false otherwise.
PointState
Enumeration for the state of a point. Normal means, that the point is not on a BreakLine,...
int addPoint(const QgsPoint &p) override
Adds a point to the triangulation.
An interface for objects which accept features via addFeature(s) methods.
Base class for feedback objects to be used for cancellation of something running in a worker thread.
Definition: qgsfeedback.h:44
void setProgress(double progress)
Sets the current progress for the feedback object.
Definition: qgsfeedback.h:61
Point geometry type, with support for z-dimension and m-values.
Definition: qgspoint.h:49
double distance3D(double x, double y, double z) const
Returns the Cartesian 3D distance between this point and a specified x, y, z coordinate.
Definition: qgspoint.h:437
virtual QgsMesh triangulationToMesh(QgsFeedback *feedback=nullptr) const =0
Returns a QgsMesh corresponding to the triangulation.
virtual void eliminateHorizontalTriangles()=0
Eliminates the horizontal triangles by swapping.
virtual bool saveTriangulation(QgsFeatureSink *sink, QgsFeedback *feedback=nullptr) const =0
Saves the triangulation features to a feature sink.
virtual int addPoint(const QgsPoint &point)=0
Adds a point to the triangulation.
virtual bool swapEdge(double x, double y)=0
Reads the content of a taff-file.
virtual QList< int > surroundingTriangles(int pointno)=0
Returns a value list with the information of the triangles surrounding (counterclockwise) a point.
QgsPoint * point(int i) const override
Returns a pointer to the point with number i.
QgsTriangulation * mTIN
Association with a Triangulation object.
Definition: TriDecorator.h:65
QList< int > pointsAroundEdge(double x, double y) override
Returns a value list with the numbers of the four points, which would be affected by an edge swap.
QList< int > surroundingTriangles(int pointno) override
Returns a value list with the information of the triangles surrounding (counterclockwise) a point.
int pointsCount() const override
Returns the number of points.
bool triangleVertices(double x, double y, QgsPoint &p1, int &n1, QgsPoint &p2, int &n2, QgsPoint &p3, int &n3) override
Finds out in which triangle the point with coordinates x and y is and assigns the numbers of the vert...
virtual bool calcNormVec(double x, double y, QgsPoint &result)=0
Calculates the normal vector and assigns it to vec.
virtual bool calcPoint(double x, double y, QgsPoint &result)=0
Performs a linear interpolation in a triangle and assigns the x-,y- and z-coordinates to point.
Class Vector3D represents a 3D-Vector, capable to store x-,y- and z-coordinates in double values.
Definition: Vector3D.h:36
void standardise()
Standardises the vector.
Definition: Vector3D.cpp:24
void setX(double x)
Sets the x-component of the vector.
Definition: Vector3D.h:106
double getY() const
Returns the y-component of the vector.
Definition: Vector3D.h:96
double getX() const
Returns the x-component of the vector.
Definition: Vector3D.h:91
void setY(double y)
Sets the y-component of the vector.
Definition: Vector3D.h:111
double getZ() const
Returns the z-component of the vector.
Definition: Vector3D.h:101
void setZ(double z)
Sets the z-component of the vector.
Definition: Vector3D.h:116
bool ANALYSIS_EXPORT pointInsideTriangle(double x, double y, QgsPoint *p1, QgsPoint *p2, QgsPoint *p3)
Returns true, if the point with coordinates x and y is inside (or at the edge) of the triangle p1,...
Definition: MathUtils.cpp:610
void ANALYSIS_EXPORT normalFromPoints(QgsPoint *p1, QgsPoint *p2, QgsPoint *p3, Vector3D *vec)
Calculates the normal vector of the plane through the points p1, p2 and p3 and assigns the result to ...
Definition: MathUtils.cpp:562
#define QgsDebugMsgLevel(str, level)
Definition: qgslogger.h:39
#define QgsDebugError(str)
Definition: qgslogger.h:38
Mesh - vertices, edges and faces.