@@ -168,7 +168,61 @@ static QVector3D _calculateNormal( const QgsCurve *curve, double originX, double
168
168
return normal;
169
169
}
170
170
171
- void QgsTessellator::addPolygon ( const QgsPolygonV2 &polygon, float extrusionHeight )
171
+
172
+ static void _normalVectorToXYVectors ( const QVector3D &pNormal, QVector3D &pXVector, QVector3D &pYVector )
173
+ {
174
+ // Here we define the two perpendicular vectors that define the local
175
+ // 2D space on the plane. They will act as axis for which we will
176
+ // calculate the projection coordinates of a 3D point to the plane.
177
+ if ( pNormal.z () > 0.001 || pNormal.z () < -0.001 )
178
+ {
179
+ pXVector = QVector3D ( 1 , 0 , -pNormal.x () / pNormal.z () );
180
+ }
181
+ else if ( pNormal.y () > 0.001 || pNormal.y () < -0.001 )
182
+ {
183
+ pXVector = QVector3D ( 1 , -pNormal.x () / pNormal.y (), 0 );
184
+ }
185
+ else
186
+ {
187
+ pXVector = QVector3D ( -pNormal.y () / pNormal.x (), 1 , 0 );
188
+ }
189
+ pXVector.normalize ();
190
+ pYVector = QVector3D::normal ( pNormal, pXVector );
191
+ }
192
+
193
+
194
+ static void _ringToPoly2tri ( const QgsCurve *ring, const QgsPoint &ptFirst, const QVector3D &pXVector, const QVector3D &pYVector, std::vector<p2t::Point *> &polyline, QHash<p2t::Point *, float > &zHash )
195
+ {
196
+ QgsVertexId::VertexType vt;
197
+ QgsPoint pt;
198
+
199
+ const int pCount = ring->numPoints ();
200
+ double x0 = ptFirst.x (), y0 = ptFirst.y (), z0 = ( std::isnan ( ptFirst.z () ) ? 0 : ptFirst.z () );
201
+
202
+ polyline.reserve ( pCount );
203
+
204
+ for ( int i = 0 ; i < pCount - 1 ; ++i )
205
+ {
206
+ ring->pointAt ( i, pt, vt );
207
+ const float z = std::isnan ( pt.z () ) ? 0 : pt.z ();
208
+ QVector3D tempPt ( pt.x () - x0, pt.y () - y0, z - z0 );
209
+ const float x = QVector3D::dotProduct ( tempPt, pXVector );
210
+ const float y = QVector3D::dotProduct ( tempPt, pYVector );
211
+
212
+ const bool found = std::find_if ( polyline.begin (), polyline.end (), [x, y]( p2t::Point *&p ) { return *p == p2t::Point ( x, y ); } ) != polyline.end ();
213
+
214
+ if ( found )
215
+ {
216
+ continue ;
217
+ }
218
+
219
+ p2t::Point *pt2 = new p2t::Point ( x, y );
220
+ polyline.push_back ( pt2 );
221
+ zHash[pt2] = z;
222
+ }
223
+ }
224
+
225
+ static bool _check_intersecting_rings ( const QgsPolygonV2 &polygon )
172
226
{
173
227
// At this point we assume that input polygons are valid according to the OGC definition.
174
228
// This means e.g. no duplicate points, polygons are simple (no butterfly shaped polygon with self-intersection),
@@ -198,24 +252,28 @@ void QgsTessellator::addPolygon( const QgsPolygonV2 &polygon, float extrusionHei
198
252
for ( int j = i + 1 ; j < geomRings.count (); ++j )
199
253
{
200
254
if ( geomRings[i].intersects ( geomRings[j] ) )
201
- {
202
- // skip the polygon - it would cause a crash inside poly2tri library
203
- qDebug () << " polygon rings intersect each other - skipping" ;
204
- return ;
205
- }
255
+ return false ;
206
256
}
207
257
}
258
+ return true ;
259
+ }
260
+
261
+
262
+ void QgsTessellator::addPolygon ( const QgsPolygonV2 &polygon, float extrusionHeight )
263
+ {
264
+ if ( !_check_intersecting_rings ( polygon ) )
265
+ {
266
+ // skip the polygon - it would cause a crash inside poly2tri library
267
+ qDebug () << " polygon rings intersect each other - skipping" ;
268
+ return ;
269
+ }
208
270
209
271
const QgsCurve *exterior = polygon.exteriorRing ();
210
272
211
273
QList< std::vector<p2t::Point *> > polylinesToDelete;
212
274
QHash<p2t::Point *, float > z;
213
275
214
276
std::vector<p2t::Point *> polyline;
215
- polyline.reserve ( exterior->numPoints () );
216
-
217
- QgsVertexId::VertexType vt;
218
- QgsPoint pt;
219
277
220
278
const QVector3D pNormal = _calculateNormal ( exterior, mOriginX , mOriginY );
221
279
const int pCount = exterior->numPoints ();
@@ -224,6 +282,7 @@ void QgsTessellator::addPolygon( const QgsPolygonV2 &polygon, float extrusionHei
224
282
if ( pCount == 4 )
225
283
{
226
284
QgsPoint pt;
285
+ QgsVertexId::VertexType vt;
227
286
for ( int i = 0 ; i < 3 ; i++ )
228
287
{
229
288
exterior->pointAt ( i, pt, vt );
@@ -234,65 +293,28 @@ void QgsTessellator::addPolygon( const QgsPolygonV2 &polygon, float extrusionHei
234
293
}
235
294
else
236
295
{
237
-
238
296
if ( !qgsDoubleNear ( pNormal.length (), 1 , 0.001 ) )
239
- {
240
- return ;
241
- }
242
-
243
- const QgsPoint ptFirst ( exterior->startPoint () );
244
- QVector3D pOrigin ( ptFirst.x (), ptFirst.y (), std::isnan ( ptFirst.z () ) ? 0 : ptFirst.z () );
245
- QVector3D pXVector;
246
- // Here we define the two perpendicular vectors that define the local
247
- // 2D space on the plane. They will act as axis for which we will
248
- // calculate the projection coordinates of a 3D point to the plane.
249
- if ( pNormal.z () > 0.001 || pNormal.z () < -0.001 )
250
- {
251
- pXVector = QVector3D ( 1 , 0 , -pNormal.x () / pNormal.z () );
252
- }
253
- else if ( pNormal.y () > 0.001 || pNormal.y () < -0.001 )
254
- {
255
- pXVector = QVector3D ( 1 , -pNormal.x () / pNormal.y (), 0 );
256
- }
257
- else
258
- {
259
- pXVector = QVector3D ( -pNormal.y () / pNormal.x (), 1 , 0 );
260
- }
261
- pXVector.normalize ();
262
- const QVector3D pYVector = QVector3D::normal ( pNormal, pXVector );
297
+ return ; // this should not happen - pNormal should be normalized to unit length
263
298
264
- for ( int i = 0 ; i < pCount - 1 ; ++i )
265
- {
266
- exterior->pointAt ( i, pt, vt );
267
- QVector3D tempPt ( pt.x (), pt.y (), ( std::isnan ( pt.z () ) ? 0 : pt.z () ) );
268
- const float x = QVector3D::dotProduct ( tempPt - pOrigin, pXVector );
269
- const float y = QVector3D::dotProduct ( tempPt - pOrigin, pYVector );
270
-
271
- const bool found = std::find_if ( polyline.begin (), polyline.end (), [x, y]( p2t::Point *&p ) { return *p == p2t::Point ( x, y ); } ) != polyline.end ();
299
+ QVector3D pXVector, pYVector;
300
+ _normalVectorToXYVectors ( pNormal, pXVector, pYVector );
272
301
273
- if ( found )
274
- {
275
- continue ;
276
- }
277
-
278
- p2t::Point *pt2 = new p2t::Point ( x, y );
279
- polyline.push_back ( pt2 );
280
-
281
- z[pt2] = std::isnan ( pt.z () ) ? 0 : pt.z ();
282
- }
302
+ const QgsPoint ptFirst ( exterior->startPoint () );
303
+ _ringToPoly2tri ( exterior, ptFirst, pXVector, pYVector, polyline, z );
283
304
polylinesToDelete << polyline;
284
305
285
306
// TODO: robustness (no nearly duplicate points, invalid geometries ...)
286
307
308
+ double x0 = ptFirst.x (), y0 = ptFirst.y ();
287
309
if ( polyline.size () == 3 && polygon.numInteriorRings () == 0 )
288
310
{
289
311
for ( std::vector<p2t::Point *>::iterator it = polyline.begin (); it != polyline.end (); it++ )
290
312
{
291
313
p2t::Point *p = *it;
292
314
const double zPt = z[p];
293
- QVector3D nPoint = pOrigin + pXVector * p->x + pYVector * p->y ;
294
- const double fx = nPoint.x () - mOriginX ;
295
- const double fy = nPoint.y () - mOriginY ;
315
+ QVector3D nPoint = pXVector * p->x + pYVector * p->y ;
316
+ const double fx = nPoint.x () - mOriginX + x0 ;
317
+ const double fy = nPoint.y () - mOriginY + y0 ;
296
318
const double fz = extrusionHeight + ( std::isnan ( zPt ) ? 0 : zPt );
297
319
mData << fx << fz << -fy;
298
320
if ( mAddNormals )
@@ -307,31 +329,14 @@ void QgsTessellator::addPolygon( const QgsPolygonV2 &polygon, float extrusionHei
307
329
for ( int i = 0 ; i < polygon.numInteriorRings (); ++i )
308
330
{
309
331
std::vector<p2t::Point *> holePolyline;
310
- holePolyline.reserve ( exterior->numPoints () );
311
332
const QgsCurve *hole = polygon.interiorRing ( i );
312
- for ( int j = 0 ; j < hole->numPoints () - 1 ; ++j )
313
- {
314
- hole->pointAt ( j, pt, vt );
315
- QVector3D tempPt ( pt.x (), pt.y (), ( std::isnan ( pt.z () ) ? 0 : pt.z () ) );
316
-
317
- const float x = QVector3D::dotProduct ( tempPt - pOrigin, pXVector );
318
- const float y = QVector3D::dotProduct ( tempPt - pOrigin, pYVector );
319
-
320
- const bool found = std::find_if ( polyline.begin (), polyline.end (), [x, y]( p2t::Point *&p ) { return *p == p2t::Point ( x, y ); } ) != polyline.end ();
321
-
322
- if ( found )
323
- {
324
- continue ;
325
- }
326
333
327
- p2t::Point *pt2 = new p2t::Point ( x, y );
328
- holePolyline.push_back ( pt2 );
334
+ _ringToPoly2tri ( hole, ptFirst, pXVector, pYVector, holePolyline, z );
329
335
330
- z[pt2] = std::isnan ( pt.z () ) ? 0 : pt.z ();
331
- }
332
336
cdt->AddHole ( holePolyline );
333
337
polylinesToDelete << holePolyline;
334
338
}
339
+
335
340
try
336
341
{
337
342
cdt->Triangulate ();
@@ -344,11 +349,11 @@ void QgsTessellator::addPolygon( const QgsPolygonV2 &polygon, float extrusionHei
344
349
for ( int j = 0 ; j < 3 ; ++j )
345
350
{
346
351
p2t::Point *p = t->GetPoint ( j );
347
- float zPt = z[p];
348
- QVector3D nPoint = pOrigin + pXVector * p->x + pYVector * p->y ;
349
- float fx = nPoint.x () - mOriginX ;
350
- float fy = nPoint.y () - mOriginY ;
351
- float fz = extrusionHeight + ( std::isnan ( zPt ) ? 0 : zPt );
352
+ const double zPt = z[p];
353
+ QVector3D nPoint = pXVector * p->x + pYVector * p->y ;
354
+ const double fx = nPoint.x () - mOriginX + x0 ;
355
+ const double fy = nPoint.y () - mOriginY + y0 ;
356
+ const double fz = extrusionHeight + ( std::isnan ( zPt ) ? 0 : zPt );
352
357
mData << fx << fz << -fy;
353
358
if ( mAddNormals )
354
359
mData << pNormal.x () << pNormal.z () << - pNormal.y ();
0 commit comments