Skip to content

Commit 1212a32

Browse files
author
mhugent
committed
Fix crash with constrained triangulation. Fixes bug #2005, but some numerical problems with constrained delaunay triangulations remain
git-svn-id: http://svn.osgeo.org/qgis/trunk/qgis@11941 c8812cc2-4d05-0410-92ff-de0c093fc19c
1 parent 1f517db commit 1212a32

2 files changed

Lines changed: 197 additions & 4 deletions

File tree

‎src/analysis/interpolation/DualEdgeTriangulation.cc‎

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -959,6 +959,11 @@ QList<int>* DualEdgeTriangulation::getSurroundingTriangles( int pointno )
959959

960960
bool DualEdgeTriangulation::getTriangle( double x, double y, Point3D* p1, int* n1, Point3D* p2, int* n2, Point3D* p3, int* n3 )
961961
{
962+
if ( mPointVector.size() < 3 )
963+
{
964+
return false;
965+
}
966+
962967
if ( p1 && p2 && p3 )
963968
{
964969
Point3D point( x, y, 0 );
@@ -1072,6 +1077,11 @@ bool DualEdgeTriangulation::getTriangle( double x, double y, Point3D* p1, int* n
10721077

10731078
bool DualEdgeTriangulation::getTriangle( double x, double y, Point3D* p1, Point3D* p2, Point3D* p3 )
10741079
{
1080+
if ( mPointVector.size() < 3 )
1081+
{
1082+
return false;
1083+
}
1084+
10751085
if ( p1 && p2 && p3 )
10761086
{
10771087
Point3D point( x, y, 0 );

‎src/analysis/interpolation/qgstininterpolator.cpp‎

Lines changed: 187 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -167,10 +167,12 @@ int QgsTINInterpolator::insertData( QgsFeature* f, bool zCoord, int attr, InputT
167167
}
168168
}
169169

170-
//parse WKB. We cannot use the methods with QgsPoint because they don't contain z-values for 25D types
170+
//parse WKB. It is ugly, but we cannot use the methods with QgsPoint because they don't contain z-values for 25D types
171171
bool hasZValue = false;
172172
double x, y, z;
173173
unsigned char* currentWkbPtr = g->asWkb();
174+
//maybe a structure or break line
175+
Line3D* line = 0;
174176

175177
QGis::WkbType wkbType = g->wkbType();
176178
switch ( wkbType )
@@ -199,12 +201,36 @@ int QgsTINInterpolator::insertData( QgsFeature* f, bool zCoord, int attr, InputT
199201
}
200202
break;
201203
}
204+
case QGis::WKBMultiPoint25D:
205+
hasZValue = true;
206+
case QGis::WKBMultiPoint:
207+
{
208+
currentWkbPtr += ( 1 + sizeof( int ) );
209+
int* npoints = ( int* )currentWkbPtr;
210+
currentWkbPtr += sizeof( int );
211+
for ( int index = 0; index < *npoints; ++index )
212+
{
213+
currentWkbPtr += ( 1 + sizeof( int ) ); //skip endian and point type
214+
x = *(( double* )currentWkbPtr );
215+
currentWkbPtr += sizeof( double );
216+
y = *(( double* )currentWkbPtr );
217+
currentWkbPtr += sizeof( double );
218+
if ( hasZValue ) //skip z-coordinate for 25D geometries
219+
{
220+
z = *(( double* )currentWkbPtr );
221+
currentWkbPtr += sizeof( double );
222+
}
223+
else
224+
{
225+
z = attributeValue;
226+
}
227+
}
228+
break;
229+
}
202230
case QGis::WKBLineString25D:
203231
hasZValue = true;
204232
case QGis::WKBLineString:
205233
{
206-
//maybe a structure or break line
207-
Line3D* line = 0;
208234
if ( type != POINTS )
209235
{
210236
line = new Line3D();
@@ -248,8 +274,165 @@ int QgsTINInterpolator::insertData( QgsFeature* f, bool zCoord, int attr, InputT
248274
}
249275
break;
250276
}
277+
case QGis::WKBMultiLineString25D:
278+
hasZValue = true;
279+
case QGis::WKBMultiLineString:
280+
{
281+
currentWkbPtr += ( 1 + sizeof( int ) );
282+
int* nlines = ( int* )currentWkbPtr;
283+
int* npoints = 0;
284+
currentWkbPtr += sizeof( int );
285+
for ( int index = 0; index < *nlines; ++index )
286+
{
287+
if ( type != POINTS )
288+
{
289+
line = new Line3D();
290+
}
291+
currentWkbPtr += ( sizeof( int ) + 1 );
292+
npoints = ( int* )currentWkbPtr;
293+
currentWkbPtr += sizeof( int );
294+
for ( int index2 = 0; index2 < *npoints; ++index2 )
295+
{
296+
x = *(( double* )currentWkbPtr );
297+
currentWkbPtr += sizeof( double );
298+
y = *(( double* )currentWkbPtr );
299+
currentWkbPtr += sizeof( double );
300+
301+
if ( hasZValue ) //skip z-coordinate for 25D geometries
302+
{
303+
z = *(( double* ) currentWkbPtr );
304+
currentWkbPtr += sizeof( double );
305+
}
306+
else
307+
{
308+
z = attributeValue;
309+
}
310+
311+
if ( type == POINTS )
312+
{
313+
//todo: handle error code -100
314+
mTriangulation->addPoint( new Point3D( x, y, z ) );
315+
}
316+
else
317+
{
318+
line->insertPoint( new Point3D( x, y, z ) );
319+
}
320+
}
321+
if ( type != POINTS )
322+
{
323+
mTriangulation->addLine( line, type == BREAK_LINES );
324+
}
325+
}
326+
break;
327+
}
328+
case QGis::WKBPolygon25D:
329+
hasZValue = true;
330+
case QGis::WKBPolygon:
331+
{
332+
currentWkbPtr += ( 1 + sizeof( int ) );
333+
int* nrings = ( int* )currentWkbPtr;
334+
currentWkbPtr += sizeof( int );
335+
int* npoints;
336+
for ( int index = 0; index < *nrings; ++index )
337+
{
338+
if ( type != POINTS )
339+
{
340+
line = new Line3D();
341+
}
342+
343+
npoints = ( int* )currentWkbPtr;
344+
currentWkbPtr += sizeof( int );
345+
for ( int index2 = 0; index2 < *npoints; ++index2 )
346+
{
347+
x = *(( double* )currentWkbPtr );
348+
currentWkbPtr += sizeof( double );
349+
y = *(( double* )currentWkbPtr );
350+
currentWkbPtr += sizeof( double );
351+
if ( hasZValue ) //skip z-coordinate for 25D geometries
352+
{
353+
z = *(( double* )currentWkbPtr );;
354+
currentWkbPtr += sizeof( double );
355+
}
356+
else
357+
{
358+
z = attributeValue;
359+
}
360+
if ( type == POINTS )
361+
{
362+
//todo: handle error code -100
363+
mTriangulation->addPoint( new Point3D( x, y, z ) );
364+
}
365+
else
366+
{
367+
line->insertPoint( new Point3D( x, y, z ) );
368+
}
369+
}
370+
371+
if ( type != POINTS )
372+
{
373+
mTriangulation->addLine( line, type == BREAK_LINES );
374+
}
375+
}
376+
break;
377+
}
378+
379+
case QGis::WKBMultiPolygon25D:
380+
hasZValue = true;
381+
case QGis::WKBMultiPolygon:
382+
{
383+
currentWkbPtr += ( 1 + sizeof( int ) );
384+
int* npolys = ( int* )currentWkbPtr;
385+
int* nrings;
386+
int* npoints;
387+
currentWkbPtr += sizeof( int );
388+
for ( int index = 0; index < *npolys; ++index )
389+
{
390+
currentWkbPtr += ( 1 + sizeof( int ) ); //skip endian and polygon type
391+
nrings = ( int* )currentWkbPtr;
392+
currentWkbPtr += sizeof( int );
393+
for ( int index2 = 0; index2 < *nrings; ++index2 )
394+
{
395+
if ( type != POINTS )
396+
{
397+
line = new Line3D();
398+
}
399+
npoints = ( int* )currentWkbPtr;
400+
currentWkbPtr += sizeof( int );
401+
for ( int index3 = 0; index3 < *npoints; ++index3 )
402+
{
403+
x = *(( double* )currentWkbPtr );
404+
currentWkbPtr += sizeof( double );
405+
y = *(( double* )currentWkbPtr );
406+
currentWkbPtr += sizeof( double );
407+
if ( hasZValue ) //skip z-coordinate for 25D geometries
408+
{
409+
z = *(( double* )currentWkbPtr );
410+
currentWkbPtr += sizeof( double );
411+
}
412+
else
413+
{
414+
z = attributeValue;
415+
}
416+
if ( type == POINTS )
417+
{
418+
//todo: handle error code -100
419+
mTriangulation->addPoint( new Point3D( x, y, z ) );
420+
}
421+
else
422+
{
423+
line->insertPoint( new Point3D( x, y, z ) );
424+
}
425+
}
426+
if ( type != POINTS )
427+
{
428+
mTriangulation->addLine( line, type == BREAK_LINES );
429+
}
430+
}
431+
}
432+
break;
433+
}
251434
default:
252-
//todo: add the same for multiline, polygon, multipolygon
435+
//should not happen...
253436
break;
254437
}
255438

0 commit comments

Comments
 (0)