1 /**********************************************************************
3 * GEOS - Geometry Engine Open Source
4 * http://geos.osgeo.org
6 * Copyright (C) 2010 Sandro Santilli <strk@kbt.io>
8 * This is free software; you can redistribute and/or modify it under
9 * the terms of the GNU Lesser General Public Licence as published
10 * by the Free Software Foundation.
11 * See the COPYING file for more information.
13 **********************************************************************
15 * Last port: original work
17 * Developed by Sandro Santilli (strk@kbt.io)
18 * for Faunalia (http://www.faunalia.it)
19 * with funding from Regione Toscana - Settore SISTEMA INFORMATIVO
20 * TERRITORIALE ED AMBIENTALE - for the project: "Sviluppo strumenti
21 * software per il trattamento di dati geografici basati su QuantumGIS
22 * e Postgis (CIG 0494241492)"
24 **********************************************************************/
26 #include <geos/operation/sharedpaths/SharedPathsOp.h>
27 #include <geos/geom/GeometryFactory.h>
28 #include <geos/geom/Geometry.h>
29 #include <geos/geom/LineString.h>
30 #include <geos/geom/MultiLineString.h>
31 #include <geos/linearref/LinearLocation.h>
32 #include <geos/linearref/LocationIndexOfPoint.h>
33 #include <geos/operation/overlay/OverlayOp.h>
34 #include <geos/util/IllegalArgumentException.h>
36 using namespace geos::geom
;
39 namespace operation
{ // geos.operation
40 namespace sharedpaths
{ // geos.operation.sharedpaths
44 SharedPathsOp::sharedPathsOp(const Geometry
& g1
, const Geometry
& g2
,
45 PathList
& sameDirection
,
46 PathList
& oppositeDirection
)
48 SharedPathsOp
sp(g1
, g2
);
49 sp
.getSharedPaths(sameDirection
, oppositeDirection
);
53 SharedPathsOp::SharedPathsOp(
54 const geom::Geometry
& g1
, const geom::Geometry
& g2
)
60 checkLinealInput(_g1
);
61 checkLinealInput(_g2
);
66 SharedPathsOp::checkLinealInput(const geom::Geometry
& g
)
68 if ( ! dynamic_cast<const LineString
*>(&g
) &&
69 ! dynamic_cast<const MultiLineString
*>(&g
) )
71 throw util::IllegalArgumentException("Geometry is not lineal");
77 SharedPathsOp::getSharedPaths(PathList
& forwDir
, PathList
& backDir
)
80 findLinearIntersections(paths
);
81 for (size_t i
=0, n
=paths
.size(); i
<n
; ++i
)
83 LineString
* path
= paths
[i
];
84 if ( isSameDirection(*path
) ) forwDir
.push_back(path
);
85 else backDir
.push_back(path
);
91 SharedPathsOp::clearEdges(PathList
& edges
)
93 for (PathList::const_iterator
94 i
=edges
.begin(), e
=edges
.end();
104 SharedPathsOp::findLinearIntersections(PathList
& to
)
106 using geos::operation::overlay::OverlayOp
;
108 // TODO: optionally use the tolerance,
109 // snapping _g2 over _g1 ?
111 std::auto_ptr
<Geometry
> full ( OverlayOp::overlayOp(
112 &_g1
, &_g2
, OverlayOp::opINTERSECTION
) );
114 // NOTE: intersection of equal lines yelds splitted lines,
115 // should we sew them back ?
117 for (size_t i
=0, n
=full
->getNumGeometries(); i
<n
; ++i
)
119 const Geometry
* sub
= full
->getGeometryN(i
);
120 const LineString
* path
= dynamic_cast<const LineString
*>(sub
);
122 // NOTE: we're making a copy here, wouldn't be needed
123 // for a simple predicate
124 to
.push_back(_gf
.createLineString(*path
).release());
131 SharedPathsOp::isForward(const geom::LineString
& edge
,
132 const geom::Geometry
& geom
)
134 using namespace geos::linearref
;
138 1. find first point of edge on geom (linearref)
139 2. find second point of edge on geom (linearref)
140 3. if first < second, we're forward
143 1. edge has at least 2 points
144 2. edge first two points are not equal
148 const Coordinate
& pt1
= edge
.getCoordinateN(0);
149 const Coordinate
& pt2
= edge
.getCoordinateN(1);
152 * We move the coordinate somewhat closer, to avoid
153 * vertices of the geometry being checked (geom).
155 * This is mostly only needed when one of the two points
156 * of the edge is an endpoint of a _closed_ geom.
157 * We have an unit test for this...
159 Coordinate pt1i
= LinearLocation::pointAlongSegmentByFraction(pt1
, pt2
, 0.1);
160 Coordinate pt2i
= LinearLocation::pointAlongSegmentByFraction(pt1
, pt2
, 0.9);
162 LinearLocation l1
= LocationIndexOfPoint::indexOf(&geom
, pt1i
);
163 LinearLocation l2
= LocationIndexOfPoint::indexOf(&geom
, pt2i
);
165 return l1
.compareTo(l2
) < 0;
168 } // namespace geos.operation.sharedpaths
169 } // namespace geos::operation