19
19
#include " qgsgeometry.h"
20
20
#include " qgspoint.h"
21
21
#include " qgspolygon.h"
22
+ #include " qgis_sip.h"
22
23
23
24
#include " poly2tri/poly2tri.h"
24
25
@@ -110,6 +111,28 @@ static void _makeWalls( const QgsCurve &ring, bool ccw, float extrusionHeight, Q
110
111
}
111
112
}
112
113
114
+ static QVector3D _calculateNormal ( const QgsCurve *curve )
115
+ {
116
+ QgsVertexId::VertexType vt;
117
+ QgsPoint pt1, pt2;
118
+ QVector3D normal ( 0 , 0 , 0 );
119
+
120
+ for ( int i = 0 ; i < curve->numPoints () - 1 ; i++ )
121
+ {
122
+ QgsPoint pt1, pt2;
123
+ curve->pointAt ( i, pt1, vt );
124
+ curve->pointAt ( i + 1 , pt2, vt );
125
+
126
+ normal .setX ( normal .x () + ( pt1.y () - pt2.y () ) * ( pt1.z () + pt2.z () ) );
127
+ normal .setY ( normal .y () + ( pt1.z () - pt2.z () ) * ( pt1.x () + pt2.x () ) );
128
+ normal .setZ ( normal .z () + ( pt1.x () - pt2.x () ) * ( pt1.y () + pt2.y () ) );
129
+ }
130
+
131
+ normal .normalize ();
132
+
133
+ return normal ;
134
+ }
135
+
113
136
void QgsTessellator::addPolygon ( const QgsPolygonV2 &polygon, float extrusionHeight )
114
137
{
115
138
// At this point we assume that input polygons are valid according to the OGC definition.
@@ -159,20 +182,8 @@ void QgsTessellator::addPolygon( const QgsPolygonV2 &polygon, float extrusionHei
159
182
QgsVertexId::VertexType vt;
160
183
QgsPoint pt;
161
184
162
- QVector3D pNormal ( 0 , 0 , 0 );
185
+ QVector3D pNormal = _calculateNormal ( exterior );
163
186
const int pCount = exterior->numPoints ();
164
- for ( int i = 0 ; i < pCount - 1 ; i++ )
165
- {
166
- QgsPoint pt1, pt2;
167
- exterior->pointAt ( i, pt1, vt );
168
- exterior->pointAt ( i + 1 , pt2, vt );
169
-
170
- pNormal.setX ( pNormal.x () + ( pt1.y () - pt2.y () ) * ( pt1.z () + pt2.z () ) );
171
- pNormal.setY ( pNormal.y () + ( pt1.z () - pt2.z () ) * ( pt1.x () + pt2.x () ) );
172
- pNormal.setZ ( pNormal.z () + ( pt1.x () - pt2.x () ) * ( pt1.y () + pt2.y () ) );
173
- }
174
-
175
- pNormal.normalize ();
176
187
177
188
// Polygon is a triangle
178
189
if ( pCount == 4 )
@@ -185,76 +196,43 @@ void QgsTessellator::addPolygon( const QgsPolygonV2 &polygon, float extrusionHei
185
196
if ( mAddNormals )
186
197
mData << pNormal.x () << pNormal.z () << - pNormal.y ();
187
198
}
188
-
189
- return ;
190
- }
191
-
192
- if ( pNormal.length () < 0.999 || pNormal.length () > 1.001 )
193
- {
194
- return ;
195
- }
196
-
197
- const QgsPoint ptFirst ( exterior->startPoint () );
198
- QVector3D pOrigin ( ptFirst.x (), ptFirst.y (), ptFirst.z () ), pXVector;
199
- // Here we define the two perpendicular vectors that define the local
200
- // 2D space on the plane. They will act as axis for which we will
201
- // calculate the projection coordinates of a 3D point to the plane.
202
- if ( pNormal.z () > 0.001 || pNormal.z () < -0.001 )
203
- {
204
- pXVector = QVector3D ( 1 , 0 , -pNormal.x () / pNormal.z () );
205
- }
206
- else if ( pNormal.y () > 0.001 || pNormal.y () < -0.001 )
207
- {
208
- pXVector = QVector3D ( 1 , -pNormal.x () / pNormal.y (), 0 );
209
199
}
210
200
else
211
201
{
212
- pXVector = QVector3D ( -pNormal.y () / pNormal.x (), 1 , 0 );
213
- }
214
- pXVector.normalize ();
215
- QVector3D pYVector = QVector3D::normal ( pNormal, pXVector );
216
-
217
- for ( int i = 0 ; i < pCount - 1 ; ++i )
218
- {
219
- exterior->pointAt ( i, pt, vt );
220
- QVector3D tempPt ( pt.x (), pt.y (), ( qIsNaN ( pt.z () ) ? 0 : pt.z () ) );
221
- const float x = QVector3D::dotProduct ( tempPt - pOrigin, pXVector );
222
- const float y = QVector3D::dotProduct ( tempPt - pOrigin, pYVector );
223
-
224
- const bool found = std::find_if ( polyline.begin (), polyline.end (), [x, y]( p2t::Point *&p ) { return *p == *new p2t::Point ( x, y ); } ) != polyline.end ();
225
202
226
- if ( found )
203
+ if ( ! qgsDoubleNear ( pNormal. length (), 1 , 0.001 ) )
227
204
{
228
- continue ;
205
+ return ;
229
206
}
230
207
231
- p2t::Point *pt2 = new p2t::Point ( x, y );
232
- polyline.push_back ( pt2 );
233
-
234
- z[pt2] = qIsNaN ( pt.z () ) ? 0 : pt.z ();
235
- }
236
- polylinesToDelete << polyline;
237
-
238
- if ( polyline.size () < 3 )
239
- return ;
240
-
241
- p2t::CDT *cdt = new p2t::CDT ( polyline );
208
+ const QgsPoint ptFirst ( exterior->startPoint () );
209
+ QVector3D pOrigin ( ptFirst.x (), ptFirst.y (), ptFirst.z () ), pXVector;
210
+ // Here we define the two perpendicular vectors that define the local
211
+ // 2D space on the plane. They will act as axis for which we will
212
+ // calculate the projection coordinates of a 3D point to the plane.
213
+ if ( pNormal.z () > 0.001 || pNormal.z () < -0.001 )
214
+ {
215
+ pXVector = QVector3D ( 1 , 0 , -pNormal.x () / pNormal.z () );
216
+ }
217
+ else if ( pNormal.y () > 0.001 || pNormal.y () < -0.001 )
218
+ {
219
+ pXVector = QVector3D ( 1 , -pNormal.x () / pNormal.y (), 0 );
220
+ }
221
+ else
222
+ {
223
+ pXVector = QVector3D ( -pNormal.y () / pNormal.x (), 1 , 0 );
224
+ }
225
+ pXVector.normalize ();
226
+ QVector3D pYVector = QVector3D::normal ( pNormal, pXVector );
242
227
243
- // polygon holes
244
- for ( int i = 0 ; i < polygon.numInteriorRings (); ++i )
245
- {
246
- std::vector<p2t::Point *> holePolyline;
247
- holePolyline.reserve ( exterior->numPoints () );
248
- const QgsCurve *hole = polygon.interiorRing ( i );
249
- for ( int j = 0 ; j < hole->numPoints () - 1 ; ++j )
228
+ for ( int i = 0 ; i < pCount - 1 ; ++i )
250
229
{
251
- hole ->pointAt ( j , pt, vt );
230
+ exterior ->pointAt ( i , pt, vt );
252
231
QVector3D tempPt ( pt.x (), pt.y (), ( qIsNaN ( pt.z () ) ? 0 : pt.z () ) );
253
-
254
232
const float x = QVector3D::dotProduct ( tempPt - pOrigin, pXVector );
255
233
const float y = QVector3D::dotProduct ( tempPt - pOrigin, pYVector );
256
234
257
- const bool found = std::find_if ( polyline.begin (), polyline.end (), [x, y]( p2t::Point *&p ) { return *p == * new p2t::Point ( x, y ); } ) != polyline.end ();
235
+ const bool found = std::find_if ( polyline.begin (), polyline.end (), [x, y]( p2t::Point *&p ) { return *p == p2t::Point ( x, y ); } ) != polyline.end ();
258
236
259
237
if ( found )
260
238
{
@@ -266,62 +244,96 @@ void QgsTessellator::addPolygon( const QgsPolygonV2 &polygon, float extrusionHei
266
244
267
245
z[pt2] = qIsNaN ( pt.z () ) ? 0 : pt.z ();
268
246
}
269
- cdt->AddHole ( holePolyline );
270
- polylinesToDelete << holePolyline;
271
- }
247
+ polylinesToDelete << polyline;
272
248
273
- // TODO: robustness (no nearly duplicate points, invalid geometries ...)
249
+ if ( polyline.size () < 3 )
250
+ return ;
274
251
275
- if ( polyline.size () == 3 )
276
- {
277
- for ( std::vector<p2t::Point *>::iterator it = polyline.begin (); it != polyline.end (); it++ )
252
+ p2t::CDT *cdt = new p2t::CDT ( polyline );
253
+
254
+ // polygon holes
255
+ for ( int i = 0 ; i < polygon.numInteriorRings (); ++i )
278
256
{
279
- p2t::Point *p = *it;
280
- const double zPt = z[p];
281
- QVector3D nPoint = pOrigin + pXVector * p->x + pYVector * p->y ;
282
- const double fx = nPoint.x () - mOriginX ;
283
- const double fy = nPoint.y () - mOriginY ;
284
- const double fz = extrusionHeight + ( qIsNaN ( zPt ) ? 0 : zPt );
285
- mData << fx << fz << -fy;
286
- if ( mAddNormals )
287
- mData << pNormal.x () << pNormal.z () << - pNormal.y ();
288
- }
257
+ std::vector<p2t::Point *> holePolyline;
258
+ holePolyline.reserve ( exterior->numPoints () );
259
+ const QgsCurve *hole = polygon.interiorRing ( i );
260
+ for ( int j = 0 ; j < hole->numPoints () - 1 ; ++j )
261
+ {
262
+ hole->pointAt ( j, pt, vt );
263
+ QVector3D tempPt ( pt.x (), pt.y (), ( qIsNaN ( pt.z () ) ? 0 : pt.z () ) );
289
264
290
- return ;
291
- }
265
+ const float x = QVector3D::dotProduct ( tempPt - pOrigin, pXVector ) ;
266
+ const float y = QVector3D::dotProduct ( tempPt - pOrigin, pYVector );
292
267
293
- try
294
- {
295
- cdt->Triangulate ();
296
- }
297
- catch ( ... )
298
- {
299
- qDebug () << " Triangulation failed. Skipping polygon..." ;
300
- return ;
301
- }
268
+ const bool found = std::find_if ( polyline.begin (), polyline.end (), [x, y]( p2t::Point *&p ) { return *p == p2t::Point ( x, y ); } ) != polyline.end ();
269
+
270
+ if ( found )
271
+ {
272
+ continue ;
273
+ }
302
274
303
- std::vector<p2t::Triangle *> triangles = cdt->GetTriangles ();
275
+ p2t::Point *pt2 = new p2t::Point ( x, y );
276
+ polyline.push_back ( pt2 );
304
277
305
- for ( size_t i = 0 ; i < triangles.size (); ++i )
306
- {
307
- p2t::Triangle *t = triangles[i];
308
- for ( int j = 0 ; j < 3 ; ++j )
278
+ z[pt2] = qIsNaN ( pt.z () ) ? 0 : pt.z ();
279
+ }
280
+ cdt->AddHole ( holePolyline );
281
+ polylinesToDelete << holePolyline;
282
+ }
283
+
284
+ // TODO: robustness (no nearly duplicate points, invalid geometries ...)
285
+
286
+ if ( polyline.size () == 3 && polygon.numInteriorRings () == 0 )
309
287
{
310
- p2t::Point *p = t->GetPoint ( j );
311
- float zPt = z[p];
312
- QVector3D nPoint = pOrigin + pXVector * p->x + pYVector * p->y ;
313
- float fx = nPoint.x () - mOriginX ;
314
- float fy = nPoint.y () - mOriginY ;
315
- float fz = extrusionHeight + ( qIsNaN ( zPt ) ? 0 : zPt );
316
- mData << fx << fz << -fy;
317
- if ( mAddNormals )
318
- mData << pNormal.x () << pNormal.z () << - pNormal.y ();
288
+ for ( std::vector<p2t::Point *>::iterator it = polyline.begin (); it != polyline.end (); it++ )
289
+ {
290
+ p2t::Point *p = *it;
291
+ const double zPt = z[p];
292
+ QVector3D nPoint = pOrigin + pXVector * p->x + pYVector * p->y ;
293
+ const double fx = nPoint.x () - mOriginX ;
294
+ const double fy = nPoint.y () - mOriginY ;
295
+ const double fz = extrusionHeight + ( qIsNaN ( zPt ) ? 0 : zPt );
296
+ mData << fx << fz << -fy;
297
+ if ( mAddNormals )
298
+ mData << pNormal.x () << pNormal.z () << - pNormal.y ();
299
+ }
300
+ }
301
+ else
302
+ {
303
+ try
304
+ {
305
+ cdt->Triangulate ();
306
+ }
307
+ catch ( ... )
308
+ {
309
+ qDebug () << " Triangulation failed. Skipping polygon..." ;
310
+ return ;
311
+ }
312
+
313
+ std::vector<p2t::Triangle *> triangles = cdt->GetTriangles ();
314
+
315
+ for ( size_t i = 0 ; i < triangles.size (); ++i )
316
+ {
317
+ p2t::Triangle *t = triangles[i];
318
+ for ( int j = 0 ; j < 3 ; ++j )
319
+ {
320
+ p2t::Point *p = t->GetPoint ( j );
321
+ float zPt = z[p];
322
+ QVector3D nPoint = pOrigin + pXVector * p->x + pYVector * p->y ;
323
+ float fx = nPoint.x () - mOriginX ;
324
+ float fy = nPoint.y () - mOriginY ;
325
+ float fz = extrusionHeight + ( qIsNaN ( zPt ) ? 0 : zPt );
326
+ mData << fx << fz << -fy;
327
+ if ( mAddNormals )
328
+ mData << pNormal.x () << pNormal.z () << - pNormal.y ();
329
+ }
330
+ }
319
331
}
320
- }
321
332
322
- delete cdt;
323
- for ( int i = 0 ; i < polylinesToDelete.count (); ++i )
324
- qDeleteAll ( polylinesToDelete[i] );
333
+ delete cdt;
334
+ for ( int i = 0 ; i < polylinesToDelete.count (); ++i )
335
+ qDeleteAll ( polylinesToDelete[i] );
336
+ }
325
337
326
338
// add walls if extrusion is enabled
327
339
if ( extrusionHeight != 0 )
0 commit comments