sharpetronics/QGIS
1
0
Fork 0
mirror of https://github.com/qgis/QGIS.git synced 2025-02-25 00:58:06 -05:00
Code Issues Packages Projects Releases Wiki Activity

[FEATURE] A reshape tool to apply to line/polygon geometries. The part of a geometry between the first and last intersection of the reshape line will be replaced

Browse Source 
git-svn-id: http://svn.osgeo.org/qgis/trunk@11500 c8812cc2-4d05-0410-92ff-de0c093fc19c
This commit is contained in:
mhugent 2009-08-25 19:54:28 +00:00
parent 61604a3260
commit 20591bac19
9 changed files with 719 additions and 24 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

BIN
images/themes/default/mActionReshape.png Normal file
View File

Binary file not shown.
Side by Side

After

Width:  |  Height:  |  Size: 1.8 KiB

5
python/core/qgsgeometry.sip
View File

@ -204,6 +204,11 @@ not disjoint with existing polygons of the feature*/
@return 0 in case of success, 1 if geometry has not been split, error else*/
int splitGeometry(const QList<QgsPoint>& splitLine, QList<QgsGeometry*>& newGeometries, bool topological, QList<QgsPoint>& topologyTestPoints);
/**Replaces a part of this geometry with another line
@return 0 in case of success
@note: this function was added in version 1.3*/
int reshapeGeometry( const QList<QgsPoint>& reshapeWithLine );
/**Changes this geometry such that it does not intersect the other geometry
@param other geometry that should not be intersect
@return 0 in case of success*/

1
src/app/CMakeLists.txt
View File

@ -39,6 +39,7 @@ SET(QGIS_APP_SRCS
qgsmaptoolmovefeature.cpp
qgsmaptoolmovevertex.cpp
qgsmaptoolnodetool.cpp
qgsmaptoolreshape.cpp
qgsmaptoolselect.cpp
qgsmaptoolsimplify.cpp
qgsmaptoolsplitfeatures.cpp

25
src/app/qgisapp.cpp
View File

@ -171,6 +171,7 @@
#include "qgsmaptoolnodetool.h"
#include "qgsmaptoolpan.h"
#include "qgsmaptoolselect.h"
#include "qgsmaptoolreshape.h"
#include "qgsmaptoolsplitfeatures.h"
#include "qgsmaptoolvertexedit.h"
#include "qgsmaptoolzoom.h"
@ -473,6 +474,7 @@ QgisApp::~QgisApp()
delete mMapTools.mCaptureLine;
delete mMapTools.mCapturePolygon;
delete mMapTools.mMoveFeature;
delete mMapTools.mReshapeFeatures;
delete mMapTools.mSplitFeatures;
delete mMapTools.mSelect;
delete mMapTools.mVertexAdd;
@ -660,6 +662,12 @@ void QgisApp::createActions()
connect( mActionMoveFeature, SIGNAL( triggered() ), this, SLOT( moveFeature() ) );
mActionMoveFeature->setEnabled( false );
mActionReshapeFeatures = new QAction( getThemeIcon( "mActionReshape.png" ), tr( "Reshape Features" ), this );
shortcuts->registerAction( mActionReshapeFeatures );
mActionReshapeFeatures->setStatusTip( tr( "Reshape Features" ) );
connect( mActionReshapeFeatures, SIGNAL( triggered() ), this, SLOT( reshapeFeatures() ) );
mActionReshapeFeatures->setEnabled( false );
mActionSplitFeatures = new QAction( getThemeIcon( "mActionSplitFeatures.png" ), tr( "Split Features" ), this );
shortcuts->registerAction( mActionSplitFeatures );
mActionSplitFeatures->setStatusTip( tr( "Split Features" ) );
@ -1062,6 +1070,8 @@ void QgisApp::createActionGroups()
mMapToolGroup->addAction( mActionCapturePolygon );
mActionMoveFeature->setCheckable( true );
mMapToolGroup->addAction( mActionMoveFeature );
mActionReshapeFeatures->setCheckable( true );
mMapToolGroup->addAction( mActionReshapeFeatures );
mActionSplitFeatures->setCheckable( true );
mMapToolGroup->addAction( mActionSplitFeatures );
mMapToolGroup->addAction( mActionDeleteSelected );
@ -1175,6 +1185,7 @@ void QgisApp::createMenus()
mEditMenu->addAction( mActionAddIsland );
mEditMenu->addAction( mActionDeleteRing );
mEditMenu->addAction( mActionDeletePart );
mEditMenu->addAction( mActionReshapeFeatures );
mEditMenu->addAction( mActionSplitFeatures );
mEditMenu->addAction( mActionMergeFeatures );
mEditMenu->addAction( mActionNodeTool );
@ -1383,6 +1394,7 @@ void QgisApp::createToolBars()
mAdvancedDigitizeToolBar->addAction( mActionAddIsland );
mAdvancedDigitizeToolBar->addAction( mActionDeleteRing );
mAdvancedDigitizeToolBar->addAction( mActionDeletePart );
mAdvancedDigitizeToolBar->addAction( mActionReshapeFeatures );
mAdvancedDigitizeToolBar->addAction( mActionSplitFeatures );
mAdvancedDigitizeToolBar->addAction( mActionMergeFeatures );
mAdvancedDigitizeToolBar->addAction( mActionNodeTool );
@ -1605,6 +1617,7 @@ void QgisApp::setTheme( QString theThemeName )
mActionCaptureLine->setIcon( getThemeIcon( "/mActionCaptureLine.png" ) );
mActionCapturePolygon->setIcon( getThemeIcon( "/mActionCapturePolygon.png" ) );
mActionMoveFeature->setIcon( getThemeIcon( "/mActionMoveFeature.png" ) );
mActionReshapeFeatures->setIcon( getThemeIcon( "/mActionReshape.png" ) );
mActionSplitFeatures->setIcon( getThemeIcon( "/mActionSplitFeatures.png" ) );
mActionDeleteSelected->setIcon( getThemeIcon( "/mActionDeleteSelected.png" ) );
mActionAddVertex->setIcon( getThemeIcon( "/mActionAddVertex.png" ) );
@ -1735,6 +1748,8 @@ void QgisApp::createCanvas()
mActionCapturePolygon->setVisible( false );
mMapTools.mMoveFeature = new QgsMapToolMoveFeature( mMapCanvas );
mMapTools.mMoveFeature->setAction( mActionMoveFeature );
mMapTools.mReshapeFeatures = new QgsMapToolReshape( mMapCanvas );
mMapTools.mReshapeFeatures->setAction( mActionReshapeFeatures );
mMapTools.mSplitFeatures = new QgsMapToolSplitFeatures( mMapCanvas );
mMapTools.mSplitFeatures->setAction( mActionSplitFeatures );
mMapTools.mSelect = new QgsMapToolSelect( mMapCanvas );
@ -4337,6 +4352,11 @@ void QgisApp::splitFeatures()
mMapCanvas->setMapTool( mMapTools.mSplitFeatures );
}
void QgisApp::reshapeFeatures()
{
mMapCanvas->setMapTool( mMapTools.mReshapeFeatures );
}
void QgisApp::capturePoint()
{
if ( mMapCanvas && mMapCanvas->isDrawing() )
@ -5634,6 +5654,7 @@ void QgisApp::activateDeactivateLayerRelatedActions( QgsMapLayer* layer )
mActionMoveVertex->setEnabled( false );
mActionAddRing->setEnabled( false );
mActionAddIsland->setEnabled( false );
mActionReshapeFeatures->setEnabled( false );
mActionSplitFeatures->setEnabled( false );
mActionSimplifyFeature->setEnabled( false );
mActionDeleteRing->setEnabled( false );
@ -5650,6 +5671,7 @@ void QgisApp::activateDeactivateLayerRelatedActions( QgsMapLayer* layer )
{
mActionCaptureLine->setEnabled( true );
mActionCaptureLine->setVisible( true );
mActionReshapeFeatures->setEnabled( true );
mActionSplitFeatures->setEnabled( true );
mActionSimplifyFeature->setEnabled( true );
mActionDeletePart->setEnabled( true );
@ -5659,6 +5681,7 @@ void QgisApp::activateDeactivateLayerRelatedActions( QgsMapLayer* layer )
{
mActionCaptureLine->setEnabled( false );
mActionCaptureLine->setVisible( false );
mActionReshapeFeatures->setEnabled( false );
mActionSplitFeatures->setEnabled( false );
mActionSimplifyFeature->setEnabled( false );
mActionDeletePart->setEnabled( false );
@ -5679,6 +5702,7 @@ void QgisApp::activateDeactivateLayerRelatedActions( QgsMapLayer* layer )
mActionCapturePolygon->setVisible( true );
mActionAddRing->setEnabled( true );
mActionAddIsland->setEnabled( true );
mActionReshapeFeatures->setEnabled( true );
mActionSplitFeatures->setEnabled( true );
mActionSimplifyFeature->setEnabled( true );
mActionDeleteRing->setEnabled( true );
@ -5690,6 +5714,7 @@ void QgisApp::activateDeactivateLayerRelatedActions( QgsMapLayer* layer )
mActionCapturePolygon->setVisible( false );
mActionAddRing->setEnabled( false );
mActionAddIsland->setEnabled( false );
mActionReshapeFeatures->setEnabled( false );
mActionSplitFeatures->setEnabled( false );
mActionSimplifyFeature->setEnabled( false );
mActionDeleteRing->setEnabled( false );

4
src/app/qgisapp.h
View File

@ -497,6 +497,8 @@ class QgisApp : public QMainWindow
void deleteSelected();
//! activates the move feature tool
void moveFeature();
//! activates the reshape features tool
void reshapeFeatures();
//! activates the split features tool
void splitFeatures();
//! activates the add vertex tool
@ -725,6 +727,7 @@ class QgisApp : public QMainWindow
QAction *mActionCapturePolygon;
QAction *mActionDeleteSelected;
QAction *mActionMoveFeature;
QAction *mActionReshapeFeatures;
QAction *mActionSplitFeatures;
QAction *mActionAddVertex;
QAction *mActionDeleteVertex;
@ -847,6 +850,7 @@ class QgisApp : public QMainWindow
QgsMapTool* mCaptureLine;
QgsMapTool* mCapturePolygon;
QgsMapTool* mMoveFeature;
QgsMapTool* mReshapeFeatures;
QgsMapTool* mSplitFeatures;
QgsMapTool* mSelect;
QgsMapTool* mVertexAdd;

128
src/app/qgsmaptoolreshape.cpp Normal file
View File

@ -0,0 +1,128 @@
/***************************************************************************
qgsmaptoolreshape.cpp
---------------------------
begin : Juli 2009
copyright : (C) 2009 by Marco Hugentobler
email : marco dot hugentobler at karto dot baug dot ethz dot ch
***************************************************************************
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
***************************************************************************/
/* $Id$ */
#include "qgsmaptoolreshape.h"
#include "qgsgeometry.h"
#include "qgsmapcanvas.h"
#include "qgsrubberband.h"
#include "qgsvectorlayer.h"
#include <QMessageBox>
#include <QMouseEvent>
QgsMapToolReshape::QgsMapToolReshape( QgsMapCanvas* canvas ): QgsMapToolCapture( canvas, QgsMapToolCapture::CaptureLine )
{
}
QgsMapToolReshape::~QgsMapToolReshape()
{
}
void QgsMapToolReshape::canvasReleaseEvent( QMouseEvent * e )
{
//check if we operate on a vector layer //todo: move this to a function in parent class to avoid duplication
QgsVectorLayer *vlayer = dynamic_cast <QgsVectorLayer*>( mCanvas->currentLayer() );
if ( !vlayer )
{
QMessageBox::information( 0, tr( "Not a vector layer" ),
tr( "The current layer is not a vector layer" ) );
return;
}
if ( !vlayer->isEditable() )
{
QMessageBox::information( 0, tr( "Layer not editable" ),
tr( "Cannot edit the vector layer. To make it editable, go to the file item "
"of the layer, right click and check 'Allow Editing'." ) );
return;
}
//add point to list and to rubber band
int error = addVertex( e->pos() );
if ( error == 1 )
{
//current layer is not a vector layer
return;
}
else if ( error == 2 )
{
//problem with coordinate transformation
QMessageBox::information( 0, tr( "Coordinate transform error" ),
tr( "Cannot transform the point to the layers coordinate system" ) );
return;
}
if ( e->button() == Qt::LeftButton )
{
mCapturing = TRUE;
}
else if ( e->button() == Qt::RightButton )
{
mCapturing = FALSE;
delete mRubberBand;
mRubberBand = 0;
//find out bounding box of mCaptureList
if(mCaptureList.size() < 1)
{
return;
}
QgsPoint firstPoint = mCaptureList.at(0);
QgsRectangle bbox(firstPoint.x(), firstPoint.y(), firstPoint.x(), firstPoint.y());
for(int i = 1; i < mCaptureList.size(); ++i)
{
bbox.combineExtentWith(mCaptureList.at(i).x(), mCaptureList.at(i).y());
}
//query all the features that intersect bounding box of capture line
vlayer->select(QgsAttributeList(), bbox, true, false);
QgsFeature f;
int reshapeReturn;
bool reshapeDone = false;
vlayer->beginEditCommand( tr( "Reshape" ) );
while(vlayer->nextFeature(f))
{
//query geometry
//call geometry->reshape(mCaptureList)
//register changed geometry in vector layer
QgsGeometry* geom = f.geometry();
if(geom)
{
reshapeReturn = geom->reshapeGeometry(mCaptureList);
if(reshapeReturn == 0)
{
vlayer->changeGeometry(f.id(), geom);
reshapeDone = true;
}
}
}
if(reshapeDone)
{
vlayer->endEditCommand();
}
else
{
vlayer->destroyEditCommand();
}
mCaptureList.clear();
mCanvas->refresh();
}
}

31
src/app/qgsmaptoolreshape.h Normal file
View File

@ -0,0 +1,31 @@
/***************************************************************************
qgsmaptoolreshape.h
---------------------
begin : Juli 2009
copyright : (C) 2009 by Marco Hugentobler
email : marco.hugentobler at karto dot baug dot ethz dot ch
***************************************************************************
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
***************************************************************************/
/* $Id$ */
#ifndef QGSMAPTOOLRESHAPE_H
#define QGSMAPTOOLRESHAPE_H
#include "qgsmaptoolcapture.h"
/**A map tool that draws a line and splits the features cut by the line*/
class QgsMapToolReshape: public QgsMapToolCapture
{
public:
QgsMapToolReshape( QgsMapCanvas* canvas );
virtual ~QgsMapToolReshape();
void canvasReleaseEvent( QMouseEvent * e );
};
#endif

518
src/core/qgsgeometry.cpp
View File

@ -17,6 +17,7 @@ email : morb at ozemail dot com dot au
#include <limits>
#include <cstdarg>
#include <cstdio>
#include <cmath>
#include "qgis.h"
#include "qgsgeometry.h"
@ -3219,6 +3220,113 @@ int QgsGeometry::splitGeometry( const QList<QgsPoint>& splitLine, QList<QgsGeome
return returnCode;
}
/**Replaces a part of this geometry with another line*/
int QgsGeometry::reshapeGeometry( const QList<QgsPoint>& reshapeWithLine )
{
int returnCode = 0;
if ( type() == QGis::Point )
{
return 1; //cannot reshape points
}
GEOSGeometry* reshapeLineGeos = createGeosLineString( reshapeWithLine.toVector() );
//make sure this geos geometry is up-to-date
if ( !mGeos || mDirtyGeos )
{
exportWkbToGeos();
}
//single or multi?
int numGeoms = GEOSGetNumGeometries(mGeos);
if(numGeoms == -1)
{
return 1;
}
bool isMultiGeom = (numGeoms > 1);
bool isLine = ( type() == QGis::Line);
//polygon or multipolygon?
if ( !isMultiGeom )
{
GEOSGeometry* reshapedGeometry;
if(isLine)
{
reshapedGeometry = reshapeLine( mGeos, reshapeLineGeos );
}
else
{
reshapedGeometry = reshapePolygon( mGeos, reshapeLineGeos );
}
GEOSGeom_destroy(reshapeLineGeos);
if(reshapedGeometry)
{
GEOSGeom_destroy( mGeos );
mGeos = reshapedGeometry;
mDirtyWkb = true;
return 0;
}
else
{
return 1;
}
}
else
{
//call reshape for each geometry part and replace mGeos with new geometry if reshape took place
bool reshapeTookPlace = false;
GEOSGeometry* currentReshapeGeometry = 0;
GEOSGeometry** newGeoms = new GEOSGeometry*[numGeoms];
for(int i = 0; i < numGeoms; ++i)
{
if(isLine)
{
currentReshapeGeometry = reshapeLine( GEOSGetGeometryN(mGeos, i), reshapeLineGeos);
}
else
{
currentReshapeGeometry = reshapePolygon( GEOSGetGeometryN(mGeos, i), reshapeLineGeos);
}
if(currentReshapeGeometry)
{
newGeoms[i] = currentReshapeGeometry;
reshapeTookPlace = true;
}
else
{
newGeoms[i] = GEOSGeom_clone( GEOSGetGeometryN(mGeos, i) );
}
}
GEOSGeom_destroy(reshapeLineGeos);
GEOSGeometry* newMultiGeom = GEOSGeom_createCollection(GEOS_MULTIPOLYGON, newGeoms, numGeoms);
delete[] newGeoms;
if( ! newMultiGeom )
{
return 3;
}
if(reshapeTookPlace)
{
GEOSGeom_destroy( mGeos );
mGeos = newMultiGeom;
mDirtyWkb = true;
return 0;
}
else
{
GEOSGeom_destroy(newMultiGeom);
return 1;
}
}
return 1;
}
int QgsGeometry::makeDifference( QgsGeometry* other )
{
//make sure geos geometry is up to date
@ -4801,28 +4909,13 @@ int QgsGeometry::splitLinearGeometry( GEOSGeometry *splitLine, QList<QgsGeometry
QVector<GEOSGeometry*> testedGeometries;
GEOSGeometry* intersectGeom = 0;
//Create a small buffer around the original geometry
//and intersect candidate line segments with the buffer.
//Then we use the ratio intersection length / segment length to
//decide if the line segment belongs to the original geometry or
//if it is part of the splitting line
double bufferDistance = 0.0000001;
for ( int i = 0; i < GEOSGetNumGeometries( mergedLines ); i++ )
{
const GEOSGeometry *testing = GEOSGetGeometryN( mergedLines, i );
intersectGeom = GEOSIntersection( mGeos, GEOSBuffer( testing, bufferDistance, DEFAULT_QUADRANT_SEGMENTS ) );
double len;
GEOSLength( intersectGeom, &len );
double testingLen;
GEOSLength( testing, &testingLen );
double ratio = len / testingLen;
//the ratios for geometries that belong to the original line are usually close to 1
if ( ratio >= 0.5 && ratio <= 1.5 )
if ( lineContainedInLine( testing, mGeos ) == 1 )
{
testedGeometries << GEOSGeom_clone( testing );
}
GEOSGeom_destroy( intersectGeom );
}
mergeGeometriesMultiTypeSplit( testedGeometries );
@ -4962,6 +5055,331 @@ int QgsGeometry::splitPolygonGeometry( GEOSGeometry* splitLine, QList<QgsGeometr
return 0;
}
GEOSGeometry* QgsGeometry::reshapePolygon( const GEOSGeometry* polygon, const GEOSGeometry* reshapeLineGeos )
{
//go through outer shell and all inner rings and check if there is exactly one intersection of a ring and the reshape line
int nIntersections = 0;
int lastIntersectingRing = -2;
const GEOSGeometry* lastIntersectingGeom = 0;
int nRings = GEOSGetNumInteriorRings(polygon);
if(nRings < 0)
{
return 0;
}
//does outer ring intersect?
const GEOSGeometry* outerRing = GEOSGetExteriorRing(polygon);
if( GEOSIntersects( outerRing, reshapeLineGeos ) == 1)
{
++nIntersections;
lastIntersectingRing = -1;
lastIntersectingGeom = outerRing;
}
//do inner rings intersect?
const GEOSGeometry* innerRings[nRings];
for(int i = 0; i < nRings; ++i)
{
innerRings[i] = GEOSGetInteriorRingN(polygon, i);
if(GEOSIntersects( innerRings[i], reshapeLineGeos) == 1)
{
++nIntersections;
lastIntersectingRing = i;
lastIntersectingGeom = innerRings[i];
}
}
if(nIntersections != 1) //reshape line is only allowed to intersect one ring
{
return 0;
}
//we have one intersecting ring, let's try to reshape it
GEOSGeometry* reshapeResult = reshapeLine( lastIntersectingGeom, reshapeLineGeos);
if(!reshapeResult)
{
return 0;
}
//if reshaping took place, we need to reassemble the polygon and its rings
GEOSGeometry* newRing = 0;
const GEOSCoordSequence* reshapeSequence = GEOSGeom_getCoordSeq(reshapeResult);
GEOSCoordSequence* newCoordSequence = GEOSCoordSeq_clone(reshapeSequence);
GEOSGeom_destroy(reshapeResult);
newRing = GEOSGeom_createLinearRing(newCoordSequence);
if(!newRing)
{
return 0;
}
GEOSGeometry* newOuterRing = 0;
if(lastIntersectingRing == -1)
{
newOuterRing = newRing;
}
else
{
newOuterRing = GEOSGeom_clone(outerRing);
}
GEOSGeometry** newInnerRings = new GEOSGeometry*[nRings];
for(int i = 0; i < nRings; ++i)
{
if( lastIntersectingRing == i)
{
newInnerRings[i] = newRing;
}
else
{
newInnerRings[i] = GEOSGeom_clone(innerRings[i]);
}
}
GEOSGeometry* reshapedPolygon = GEOSGeom_createPolygon( newOuterRing, newInnerRings, nRings);
delete[] newInnerRings;
if(!reshapedPolygon)
{
return 0;
}
return reshapedPolygon;
}
GEOSGeometry* QgsGeometry::reshapeLine( const GEOSGeometry* line, const GEOSGeometry* reshapeLineGeos )
{
if ( !line || !reshapeLineGeos )
{
return 0;
}
//make sure there are at least two instersction between line and reshape geometry
GEOSGeometry* intersectGeom = GEOSIntersection( line, reshapeLineGeos);
bool atLeastTwoIntersections = (GEOSGeomTypeId(intersectGeom) == GEOS_MULTIPOINT && GEOSGetNumGeometries(intersectGeom) > 1);
GEOSGeom_destroy(intersectGeom);
if(!atLeastTwoIntersections)
{
return 0;
}
bool isRing = false;
if(GEOSGeomTypeId(line) == GEOS_LINEARRING)
{
isRing = true;
}
//begin and end point of original line
const GEOSCoordSequence* lineCoordSeq = GEOSGeom_getCoordSeq( line );
if ( !lineCoordSeq )
{
return 0;
}
unsigned int lineCoordSeqSize;
if ( GEOS_DLL GEOSCoordSeq_getSize( lineCoordSeq, &lineCoordSeqSize ) == 0 )
{
return 0;
}
if ( lineCoordSeqSize < 2 )
{
return 0;
}
//first and last vertex of line
double x1, y1, x2, y2;
GEOS_DLL GEOSCoordSeq_getX( lineCoordSeq, 0, &x1 );
GEOS_DLL GEOSCoordSeq_getY( lineCoordSeq, 0, &y1 );
GEOS_DLL GEOSCoordSeq_getX( lineCoordSeq, lineCoordSeqSize - 1, &x2 );
GEOS_DLL GEOSCoordSeq_getY( lineCoordSeq, lineCoordSeqSize - 1, &y2 );
GEOSGeometry* beginLineVertex = createGeosPoint( QgsPoint( x1, y1 ) );
GEOSGeometry* endLineVertex = createGeosPoint( QgsPoint( x2, y2 ) );
//node line and reshape line
GEOSGeometry* nodedGeometry = nodeGeometries( reshapeLineGeos, line );
if ( !nodedGeometry )
{
GEOSGeom_destroy( beginLineVertex );
GEOSGeom_destroy( endLineVertex );
return 0;
}
//and merge them together
GEOSGeometry *mergedLines = GEOSLineMerge( nodedGeometry );
GEOSGeom_destroy( nodedGeometry );
if ( !mergedLines )
{
GEOSGeom_destroy( beginLineVertex );
GEOSGeom_destroy( endLineVertex );
return 0;
}
int numMergedLines = GEOSGetNumGeometries( mergedLines );
if ( numMergedLines < 2 ) //some special cases. Normally it is >2
{
GEOSGeom_destroy( beginLineVertex );
GEOSGeom_destroy( endLineVertex );
if ( numMergedLines == 1 ) //reshape line is from begin to endpoint. So we keep the reshapeline
{
return GEOSGeom_clone( reshapeLineGeos );
}
else
{
return 0;
}
}
QList<GEOSGeometry*> resultLineParts; //collection with the line segments that will be contained in result
QList<GEOSGeometry*> probableParts; //parts where we can decide on inclusion only after going through all the candidates
for ( int i = 0; i < numMergedLines; ++i )
{
const GEOSGeometry* currentGeom;
currentGeom = GEOSGetGeometryN( mergedLines, i );
const GEOSCoordSequence* currentCoordSeq = GEOSGeom_getCoordSeq( currentGeom );
unsigned int currentCoordSeqSize;
GEOS_DLL GEOSCoordSeq_getSize( currentCoordSeq, &currentCoordSeqSize );
if ( currentCoordSeqSize < 2 )
{
continue;
}
//get the two endpoints of the current line merge result
double xBegin, xEnd, yBegin, yEnd;
GEOS_DLL GEOSCoordSeq_getX( currentCoordSeq, 0, &xBegin );
GEOS_DLL GEOSCoordSeq_getY( currentCoordSeq, 0, &yBegin );
GEOS_DLL GEOSCoordSeq_getX( currentCoordSeq, currentCoordSeqSize - 1, &xEnd );
GEOS_DLL GEOSCoordSeq_getY( currentCoordSeq, currentCoordSeqSize - 1, &yEnd );
GEOSGeometry* beginCurrentGeomVertex = createGeosPoint( QgsPoint( xBegin, yBegin ) );
GEOSGeometry* endCurrentGeomVertex = createGeosPoint( QgsPoint( xEnd, yEnd ) );
//check how many endpoints of the line merge result are on the (original) line
int nEndpointsOnOriginalLine = 0;
if ( pointContainedInLine( beginCurrentGeomVertex, line ) == 1 )
{
nEndpointsOnOriginalLine += 1;
}
if ( pointContainedInLine( endCurrentGeomVertex, line ) == 1 )
{
nEndpointsOnOriginalLine += 1;
}
//check how many endpoints equal the endpoints of the original line
int nEndpointsSameAsOriginalLine = 0;
if ( GEOSEquals( beginCurrentGeomVertex, beginLineVertex ) == 1 || GEOSEquals( beginCurrentGeomVertex, endLineVertex ) == 1 )
{
nEndpointsSameAsOriginalLine += 1;
}
if ( GEOSEquals( endCurrentGeomVertex, beginLineVertex ) == 1 || GEOSEquals( endCurrentGeomVertex, endLineVertex ) == 1 )
{
nEndpointsSameAsOriginalLine += 1;
}
//check if the current geometry overlaps the original geometry (GEOSOverlap does not seem to work with linestrings)
bool currentGeomOverlapsOriginalGeom = false;
bool currentGeomOverlapsReshapeLine = false;
if ( QgsGeometry::lineContainedInLine( currentGeom, line ) == 1 )
{
currentGeomOverlapsOriginalGeom = true;
}
if ( QgsGeometry::lineContainedInLine( currentGeom, reshapeLineGeos ) == 1 )
{
currentGeomOverlapsReshapeLine = true;
}
//logic to decide if this part belongs to the result
if ( nEndpointsSameAsOriginalLine == 1 && nEndpointsOnOriginalLine == 2 && currentGeomOverlapsOriginalGeom )
{
resultLineParts.push_back( GEOSGeom_clone( currentGeom ) );
}
//for closed rings, we take one segment from the candidate list
else if(isRing && nEndpointsOnOriginalLine == 2 && currentGeomOverlapsOriginalGeom)
{
probableParts.push_back(GEOSGeom_clone(currentGeom));
}
else if ( nEndpointsOnOriginalLine == 2 && !currentGeomOverlapsOriginalGeom )
{
resultLineParts.push_back( GEOSGeom_clone( currentGeom ) );
}
else if ( nEndpointsSameAsOriginalLine == 2 && !currentGeomOverlapsOriginalGeom )
{
resultLineParts.push_back( GEOSGeom_clone( currentGeom ) );
}
else if ( currentGeomOverlapsOriginalGeom && currentGeomOverlapsReshapeLine )
{
resultLineParts.push_back( GEOSGeom_clone( currentGeom ) );
}
GEOSGeom_destroy( beginCurrentGeomVertex );
GEOSGeom_destroy( endCurrentGeomVertex );
}
//add the longest segment from the probable list for rings (only used for polygon rings)
if(isRing)
{
GEOSGeometry* maxGeom = 0; //the longest geometry in the probabla list
GEOSGeometry* currentGeom = 0;
double maxLength = -DBL_MAX;
double currentLength = 0;
for(int i = 0; i < probableParts.length(); ++i)
{
currentGeom = probableParts.at(i);
GEOSLength( currentGeom, &currentLength);
if(currentLength > maxLength)
{
maxLength = currentLength;
GEOSGeom_destroy(maxGeom);
maxGeom = currentGeom;
}
else
{
GEOSGeom_destroy(currentGeom);
}
}
resultLineParts.push_back(maxGeom);
}
GEOSGeom_destroy( beginLineVertex );
GEOSGeom_destroy( endLineVertex );
GEOSGeom_destroy( mergedLines );
GEOSGeometry* result = 0;
if ( resultLineParts.size() < 1 )
{
return 0;
}
if ( resultLineParts.size() == 1 ) //the whole result was reshaped
{
result = resultLineParts[0];
}
else //>1
{
GEOSGeometry* lineArray[resultLineParts.size()];
for ( int i = 0; i < resultLineParts.size(); ++i )
{
lineArray[i] = resultLineParts[i];
}
//create multiline from resultLineParts
GEOSGeometry* multiLineGeom = GEOSGeom_createCollection( GEOS_MULTILINESTRING, lineArray, resultLineParts.size() );
//then do a linemerge with the newly combined partstrings
result = GEOSLineMerge( multiLineGeom );
GEOSGeom_destroy( multiLineGeom );
}
//now test if the result is a linestring. Otherwise something went wrong
if( GEOSGeomTypeId(result) != GEOS_LINESTRING)
{
GEOSGeom_destroy(result);
return 0;
}
return result;
}
int QgsGeometry::topologicalTestPointsSplit( const GEOSGeometry* splitLine, QList<QgsPoint>& testPoints ) const
{
//Find out the intersection points between splitLineGeos and this geometry.
@ -5020,7 +5438,7 @@ int QgsGeometry::topologicalTestPointsSplit( const GEOSGeometry* splitLine, QLis
return 0;
}
GEOSGeometry *QgsGeometry::nodeGeometries( const GEOSGeometry *splitLine, GEOSGeometry *geom ) const
GEOSGeometry *QgsGeometry::nodeGeometries( const GEOSGeometry *splitLine, const GEOSGeometry *geom )
{
if ( !splitLine || !geom )
{
@ -5028,27 +5446,81 @@ GEOSGeometry *QgsGeometry::nodeGeometries( const GEOSGeometry *splitLine, GEOSGe
}
GEOSGeometry *geometryBoundary = 0;
bool deleteBoundary = false;
if ( GEOSGeomTypeId( geom ) == GEOS_POLYGON || GEOSGeomTypeId( geom ) == GEOS_MULTIPOLYGON )
{
geometryBoundary = GEOSBoundary( geom );
deleteBoundary = true;
}
else
{
geometryBoundary = geom;
geometryBoundary = GEOSGeom_clone( geom );
}
GEOSGeometry *splitLineClone = GEOSGeom_clone( splitLine );
GEOSGeometry *unionGeometry = GEOSUnion( splitLineClone, geometryBoundary );
GEOSGeom_destroy( splitLineClone );
if ( deleteBoundary )
GEOSGeom_destroy( geometryBoundary );
GEOSGeom_destroy( geometryBoundary );
return unionGeometry;
}
int QgsGeometry::lineContainedInLine( const GEOSGeometry* line1, const GEOSGeometry* line2 )
{
if ( !line1 || !line2 )
{
return -1;
}
double bufferDistance = 0.00001;
GEOSGeometry* bufferGeom = GEOSBuffer( line2, bufferDistance, DEFAULT_QUADRANT_SEGMENTS );
if ( !bufferGeom )
{
return -2;
}
GEOSGeometry* intersectionGeom = GEOSIntersection( bufferGeom, line1 );
//compare ratio between line1Length and intersectGeomLength (usually close to 1 if line1 is contained in line2)
double intersectGeomLength;
double line1Length;
GEOSLength( intersectionGeom, &intersectGeomLength );
GEOSLength( line1, &line1Length );
GEOSGeom_destroy( bufferGeom );
GEOSGeom_destroy( intersectionGeom );
double intersectRatio = line1Length / intersectGeomLength;
if ( intersectRatio > 0.9 && intersectRatio < 1.1 )
{
return 1;
}
return 0;
}
int QgsGeometry::pointContainedInLine( const GEOSGeometry* point, const GEOSGeometry* line )
{
if ( !point || !line )
{
return -1;
}
double bufferDistance = 0.0000001;
GEOSGeometry* lineBuffer = GEOSBuffer( line, bufferDistance, 8 );
if ( !lineBuffer )
{
return -2;
}
bool contained = false;
if ( GEOSContains( lineBuffer, point ) == 1 )
{
contained = true;
}
GEOSGeom_destroy( lineBuffer );
return contained;
}
int QgsGeometry::numberOfGeometries( GEOSGeometry* g ) const
{
if ( !g )

31
src/core/qgsgeometry.h
View File

@ -252,6 +252,11 @@ class CORE_EXPORT QgsGeometry
bool topological,
QList<QgsPoint>& topologyTestPoints );
/**Replaces a part of this geometry with another line
@return 0 in case of success
@note: this function was added in version 1.3*/
int reshapeGeometry( const QList<QgsPoint>& reshapeWithLine );
/**Changes this geometry such that it does not intersect the other geometry
@param other geometry that should not be intersect
@return 0 in case of success*/
@ -439,9 +444,33 @@ class CORE_EXPORT QgsGeometry
@return 0 in case of success*/
int topologicalTestPointsSplit( const GEOSGeometry* splitLine, QList<QgsPoint>& testPoints ) const;
/**Creates a new line from an original line and a reshape line. The part of the input line from the first to the last intersection with the \
reshape line will be replaced. The calling function takes ownership of the result.
@param origLine the original line
@param reshapeLineGeos the reshape line
@return the reshaped line or 0 in case of error*/
static GEOSGeometry* reshapeLine( const GEOSGeometry* origLine, const GEOSGeometry* reshapeLineGeos );
/**Creates a new polygon replacing the part from the first to the second intersection with the reshape line. As a polygon ring is always closed,
the method keeps the longer part of the existing boundary
@param polygon geometry to reshape
@param reshapeLineGeos the reshape line
@return the reshaped polygon or 0 in case of error*/
static GEOSGeometry* reshapePolygon( const GEOSGeometry* polygon, const GEOSGeometry* reshapeLineGeos );
/**Nodes together a split line and a (multi-) polygon geometry in a multilinestring
@return the noded multiline geometry or 0 in case of error. The calling function takes ownership of the node geometry*/
GEOSGeometry* nodeGeometries( const GEOSGeometry *splitLine, GEOSGeometry *poly ) const;
static GEOSGeometry* nodeGeometries( const GEOSGeometry *splitLine, const GEOSGeometry *poly );
/**Tests if line1 is completely contained in line2. This method works with a very small buffer around line2 and therefore is less prone
to numerical errors as the corresponding geos method*/
static int lineContainedInLine( const GEOSGeometry* line1, const GEOSGeometry* line2 );
/**Tests if a point is on a line. This method works with a very small buffer and is thus less prone to numerical problems as the direct geos functions.
@param point the point to test
@param line line to test
@return 0 not contained, 1 if contained, <0 in case of error*/
static int pointContainedInLine( const GEOSGeometry* point, const GEOSGeometry* line );
/**Returns number of single geometry in a geos geometry. Is save for geos 2 and 3*/
int numberOfGeometries( GEOSGeometry* g ) const;