src/operation/linemerge/LineSequencer.cpp

   1 /**********************************************************************
   2  *
   3  * GEOS - Geometry Engine Open Source
   4  * http://geos.osgeo.org
   5  *
   6  * Copyright (C) 2011 Sandro Santilli <strk@keybit.net>
   7  * Copyright (C) 2006 Refractions Research Inc.
   8  *
   9  * This is free software; you can redistribute and/or modify it under
  10  * the terms of the GNU Lesser General Public Licence as published
  11  * by the Free Software Foundation.
  12  * See the COPYING file for more information.
  13  *
  14  **********************************************************************
  15  *
  16  * Last port: operation/linemerge/LineSequencer.java r378 (JTS-1.12)
  17  *
  18  **********************************************************************/
  19
  20 #include <geos/operation/linemerge/LineSequencer.h>
  21 #include <geos/operation/linemerge/LineMergeEdge.h>
  22 #include <geos/geom/Geometry.h>
  23 #include <geos/geom/GeometryFactory.h>
  24 #include <geos/geom/MultiLineString.h>
  25 #include <geos/geom/Coordinate.h>
  26 #include <geos/geom/CoordinateSequence.h>
  27 #include <geos/geom/LineString.h>
  28 #include <geos/planargraph/Node.h>
  29 #include <geos/planargraph/DirectedEdge.h>
  30 #include <geos/planargraph/Subgraph.h>
  31 #include <geos/planargraph/algorithm/ConnectedSubgraphFinder.h>
  32 #include <geos/util/Assert.h>
  33
  34 #include <cassert>
  35 #include <limits>
  36 #include <vector>
  37
  38 using namespace std;
  39 //using namespace geos::planargraph;
  40 using namespace geos::geom;
  41 //using namespace geos::planargraph::algorithm;
  42
  43 #ifdef _MSC_VER
  44 #pragma warning(disable : 4127)
  45 #endif
  46
  47 namespace geos {
  48 namespace operation { // geos.operation
  49 namespace linemerge { // geos.operation.linemerge
  50
  51 /* static */
  52 bool
  53 LineSequencer::isSequenced(const Geometry* geom)
  54 {
  55         const MultiLineString *mls;
  56
  57         if ( 0 == (mls=dynamic_cast<const MultiLineString *>(geom)) )
  58         {
  59                 return true;
  60         }
  61
  62         // the nodes in all subgraphs which have been completely scanned
  63         Coordinate::ConstSet prevSubgraphNodes;
  64         Coordinate::ConstVect currNodes;
  65
  66         const Coordinate* lastNode = NULL;
  67
  68         for (unsigned int i=0, n=mls->getNumGeometries(); i<n; ++i)
  69         {
  70                 const LineString* lineptr = \
  71                         dynamic_cast<const LineString*>(mls->getGeometryN(i));
  72                 assert(lineptr);
  73                 const LineString& line = *lineptr;
  74
  75
  76                 const Coordinate* startNode = &(line.getCoordinateN(0));
  77                 const Coordinate* endNode = &(line.getCoordinateN(line.getNumPoints() - 1));
  78
  79                 /**
  80                  * If this linestring is connected to a previous subgraph,
  81                  * geom is not sequenced
  82                  */
  83                 if (prevSubgraphNodes.find(startNode) != prevSubgraphNodes.end())
  84                 {
  85                         return false;
  86                 }
  87                 if (prevSubgraphNodes.find(endNode) != prevSubgraphNodes.end())
  88                 {
  89                         return false;
  90                 }
  91
  92                 if (lastNode != NULL)
  93                 {
  94                         if (! startNode->equals2D(*lastNode))
  95                         {
  96                                 // start new connected sequence
  97                                 prevSubgraphNodes.insert(currNodes.begin(),
  98                                                 currNodes.end());
  99                                 currNodes.clear();
 100                         }
 101                 }
 102                 currNodes.push_back(startNode);
 103                 currNodes.push_back(endNode);
 104                 lastNode = endNode;
 105         }
 106         return true;
 107 }
 108
 109 /* private */
 110 bool
 111 LineSequencer::hasSequence(planargraph::Subgraph& graph)
 112 {
 113         int oddDegreeCount = 0;
 114         for (planargraph::NodeMap::container::const_iterator
 115                 it=graph.nodeBegin(), endIt=graph.nodeEnd();
 116                 it!=endIt;
 117                 ++it)
 118         {
 119                 planargraph::Node* node = it->second;
 120                 if (node->getDegree() % 2 == 1)
 121                 oddDegreeCount++;
 122         }
 123         return oddDegreeCount <= 2;
 124 }
 125
 126 void
 127 LineSequencer::delAll( LineSequencer::Sequences& s)
 128 {
 129   for ( Sequences::iterator i=s.begin(), e=s.end(); i!=e; ++i )
 130   {
 131     delete *i;
 132   }
 133 }
 134
 135 /*private*/
 136 LineSequencer::Sequences*
 137 LineSequencer::findSequences()
 138 {
 139         Sequences *sequences = new Sequences();
 140         planargraph::algorithm::ConnectedSubgraphFinder csFinder(graph);
 141         vector<planargraph::Subgraph*> subgraphs;
 142         csFinder.getConnectedSubgraphs(subgraphs);
 143         for (vector<planargraph::Subgraph*>::const_iterator
 144                 it=subgraphs.begin(), endIt=subgraphs.end();
 145                 it!=endIt;
 146                 ++it )
 147         {
 148                 planargraph::Subgraph* subgraph = *it;
 149                 if (hasSequence(*subgraph)) {
 150                         planargraph::DirectedEdge::NonConstList* seq=findSequence(*subgraph);
 151                         sequences->push_back(seq);
 152                 }
 153                 else {
 154                         // if any subgraph cannot be sequenced, abort
 155                         delete subgraph;
 156                         delAll(*sequences);
 157                         delete sequences;
 158                         return NULL;
 159                 }
 160                 delete subgraph;
 161         }
 162         return sequences;
 163 }
 164
 165 /*private*/
 166 void
 167 LineSequencer::addLine(const LineString *lineString)
 168 {
 169         if (factory == NULL) {
 170                 factory = lineString->getFactory();
 171         }
 172         graph.addEdge(lineString);
 173         ++lineCount;
 174 }
 175
 176 /* private */
 177 void
 178 LineSequencer::computeSequence()
 179 {
 180         if (isRun) return;
 181         isRun = true;
 182
 183         Sequences* sequences = findSequences();
 184         if (sequences == NULL) return;
 185
 186         sequencedGeometry = auto_ptr<Geometry>(buildSequencedGeometry(*sequences));
 187         isSequenceableVar = true;
 188
 189         delAll(*sequences);
 190         delete sequences;
 191
 192         // Lines were missing from result
 193         assert(lineCount == sequencedGeometry->getNumGeometries());
 194
 195         // Result is not linear
 196         assert(dynamic_cast<LineString *>(sequencedGeometry.get())
 197                 || dynamic_cast<MultiLineString *>(sequencedGeometry.get()));
 198 }
 199
 200 /*private*/
 201 Geometry*
 202 LineSequencer::buildSequencedGeometry(const Sequences& sequences)
 203 {
 204         auto_ptr<Geometry::NonConstVect> lines(new Geometry::NonConstVect);
 205
 206         for (Sequences::const_iterator
 207                 i1=sequences.begin(), i1End=sequences.end();
 208                 i1 != i1End;
 209                 ++i1)
 210         {
 211                 planargraph::DirectedEdge::NonConstList& seq = *(*i1);
 212                 for(planargraph::DirectedEdge::NonConstList::iterator i2=seq.begin(),
 213                         i2End=seq.end(); i2 != i2End; ++i2)
 214                 {
 215                         const planargraph::DirectedEdge* de = *i2;
 216                         assert(dynamic_cast<LineMergeEdge* >(de->getEdge()));
 217                         LineMergeEdge* e = static_cast<LineMergeEdge* >(de->getEdge());
 218                         const LineString* line = e->getLine();
 219
 220                         // lineToAdd will be a *copy* of input things
 221                         LineString* lineToAdd;
 222
 223                         if ( ! de->getEdgeDirection() && ! line->isClosed() ) {
 224                                 lineToAdd = reverse(line);
 225                         } else {
 226                                 Geometry* lineClone = line->clone();
 227                                 lineToAdd = dynamic_cast<LineString *>(lineClone);
 228                                 assert(lineToAdd);
 229                         }
 230
 231                         lines->push_back(lineToAdd);
 232                 }
 233         }
 234
 235         if ( lines->size() == 0 ) {
 236                 return NULL;
 237         } else {
 238                 Geometry::NonConstVect *l=lines.get();
 239                 lines.release();
 240                 return factory->buildGeometry(l);
 241         }
 242 }
 243
 244 /*static private*/
 245 LineString *
 246 LineSequencer::reverse(const LineString *line)
 247 {
 248         CoordinateSequence* cs=line->getCoordinates();
 249         CoordinateSequence::reverse(cs);
 250         return line->getFactory()->createLineString(cs);
 251 }
 252
 253 /*private static*/
 254 const planargraph::Node*
 255 LineSequencer::findLowestDegreeNode(const planargraph::Subgraph& graph)
 256 {
 257         size_t minDegree = numeric_limits<size_t>::max();
 258         const planargraph::Node* minDegreeNode = NULL;
 259         for (planargraph::NodeMap::container::const_iterator
 260                 it = graph.nodeBegin(), itEnd = graph.nodeEnd();
 261                 it != itEnd;
 262                 ++it )
 263         {
 264                 const planargraph::Node* node = (*it).second;
 265                 if (minDegreeNode == NULL || node->getDegree() < minDegree)
 266                 {
 267                         minDegree = node->getDegree();
 268                         minDegreeNode = node;
 269                 }
 270         }
 271         return minDegreeNode;
 272 }
 273
 274 /*private static*/
 275 const planargraph::DirectedEdge*
 276 LineSequencer::findUnvisitedBestOrientedDE(const planargraph::Node* node)
 277 {
 278         using planargraph::DirectedEdge;
 279         using planargraph::DirectedEdgeStar;
 280
 281         const DirectedEdge* wellOrientedDE = NULL;
 282         const DirectedEdge* unvisitedDE = NULL;
 283         const DirectedEdgeStar* des=node->getOutEdges();
 284         for (DirectedEdge::NonConstVect::const_iterator i=des->begin(),
 285                 e=des->end();
 286                 i!=e;
 287                 ++i)
 288         {
 289                 planargraph::DirectedEdge* de = *i;
 290                 if (! de->getEdge()->isVisited()) {
 291                         unvisitedDE = de;
 292                         if (de->getEdgeDirection()) wellOrientedDE = de;
 293                 }
 294         }
 295         if (wellOrientedDE != NULL)
 296                 return wellOrientedDE;
 297         return unvisitedDE;
 298 }
 299
 300
 301 /*private*/
 302 void
 303 LineSequencer::addReverseSubpath(const planargraph::DirectedEdge *de,
 304                 planargraph::DirectedEdge::NonConstList& deList,
 305                 planargraph::DirectedEdge::NonConstList::iterator lit,
 306                 bool expectedClosed)
 307 {
 308         using planargraph::Node;
 309         using planargraph::DirectedEdge;
 310
 311         // trace an unvisited path *backwards* from this de
 312         Node* endNode = de->getToNode();
 313
 314         Node* fromNode = NULL;
 315         while (true) {
 316                 deList.insert(lit, de->getSym());
 317                 de->getEdge()->setVisited(true);
 318                 fromNode = de->getFromNode();
 319                 const DirectedEdge* unvisitedOutDE = findUnvisitedBestOrientedDE(fromNode);
 320
 321                 // this must terminate, since we are continually marking edges as visited
 322                 if (unvisitedOutDE == NULL) break;
 323                 de = unvisitedOutDE->getSym();
 324         }
 325         if ( expectedClosed ) {
 326                 // the path should end at the toNode of this de,
 327                 // otherwise we have an error
 328                 util::Assert::isTrue(fromNode == endNode, "path not contiguos");
 329                 //assert(fromNode == endNode);
 330         }
 331
 332 }
 333
 334 /*private*/
 335 planargraph::DirectedEdge::NonConstList*
 336 LineSequencer::findSequence(planargraph::Subgraph& graph)
 337 {
 338         using planargraph::DirectedEdge;
 339         using planargraph::Node;
 340         using planargraph::GraphComponent;
 341
 342         GraphComponent::setVisited(graph.edgeBegin(),
 343                         graph.edgeEnd(), false);
 344
 345         const Node* startNode = findLowestDegreeNode(graph);
 346
 347         const DirectedEdge *startDE = *(startNode->getOutEdges()->begin());
 348         const DirectedEdge *startDESym = startDE->getSym();
 349
 350         DirectedEdge::NonConstList *seq = new DirectedEdge::NonConstList();
 351
 352         DirectedEdge::NonConstList::iterator lit=seq->begin();
 353         addReverseSubpath(startDESym, *seq, lit, false);
 354
 355         lit=seq->end();
 356         while (lit != seq->begin()) {
 357                 const DirectedEdge* prev = *(--lit);
 358                 const DirectedEdge* unvisitedOutDE = findUnvisitedBestOrientedDE(prev->getFromNode());
 359                 if (unvisitedOutDE != NULL)
 360                         addReverseSubpath(unvisitedOutDE->getSym(), *seq, lit, true);
 361         }
 362
 363         // At this point, we have a valid sequence of graph DirectedEdges,
 364         // but it is not necessarily appropriately oriented relative to
 365         // the underlying geometry.
 366         DirectedEdge::NonConstList* orientedSeq = orient(seq);
 367
 368         if (orientedSeq != seq) delete seq;
 369
 370         return orientedSeq;
 371 }
 372
 373 /* private */
 374 planargraph::DirectedEdge::NonConstList*
 375 LineSequencer::orient(planargraph::DirectedEdge::NonConstList* seq)
 376 {
 377         using namespace geos::planargraph;
 378
 379         const DirectedEdge* startEdge = seq->front();
 380         const DirectedEdge* endEdge = seq->back();
 381         Node* startNode = startEdge->getFromNode();
 382         Node* endNode = endEdge->getToNode();
 383
 384         bool flipSeq = false;
 385         bool hasDegree1Node = \
 386                 startNode->getDegree() == 1 || endNode->getDegree() == 1;
 387
 388         if (hasDegree1Node)
 389         {
 390                 bool hasObviousStartNode = false;
 391
 392                 // test end edge before start edge, to make result stable
 393                 // (ie. if both are good starts, pick the actual start
 394                 if (endEdge->getToNode()->getDegree() == 1 &&
 395                                 endEdge->getEdgeDirection() == false)
 396                 {
 397                         hasObviousStartNode = true;
 398                         flipSeq = true;
 399                 }
 400                 if (startEdge->getFromNode()->getDegree() == 1 &&
 401                                 startEdge->getEdgeDirection() == true)
 402                 {
 403                         hasObviousStartNode = true;
 404                         flipSeq = false;
 405                 }
 406
 407                 // since there is no obvious start node,
 408                 // use any node of degree 1
 409                 if (! hasObviousStartNode)
 410                 {
 411                         // check if the start node should actually
 412                         // be the end node
 413                         if (startEdge->getFromNode()->getDegree() == 1)
 414                                 flipSeq = true;
 415                         // if the end node is of degree 1, it is
 416                         // properly the end node
 417                 }
 418
 419         }
 420
 421
 422         // if there is no degree 1 node, just use the sequence as is
 423         // (Could insert heuristic of taking direction of majority of
 424         // lines as overall direction)
 425
 426         if (flipSeq)
 427         {
 428                 return reverse(*seq);
 429         }
 430         return seq;
 431 }
 432
 433 /* private */
 434 planargraph::DirectedEdge::NonConstList*
 435 LineSequencer::reverse(planargraph::DirectedEdge::NonConstList& seq)
 436 {
 437         using namespace geos::planargraph;
 438
 439         DirectedEdge::NonConstList* newSeq = new DirectedEdge::NonConstList();
 440         DirectedEdge::NonConstList::iterator it=seq.begin(), itEnd=seq.end();
 441         for (; it!=itEnd; ++it) {
 442                 const DirectedEdge *de = *it;
 443                 newSeq->push_front(de->getSym());
 444         }
 445         return newSeq;
 446 }
 447
 448
 449
 450 } // namespace geos.operation.linemerge
 451 } // namespace geos.operation
 452 } // namespace geos
 453
 454 /**********************************************************************
 455  * $Log$
 456  * Revision 1.10  2006/06/12 10:49:43  strk
 457  * unsigned int => size_t
 458  *
 459  * Revision 1.9  2006/03/24 11:04:44  strk
 460  * Changed assert() with Assert::isTrue in addReverseSubpath
 461  *
 462  * Revision 1.8  2006/03/24 10:44:07  strk
 463  * Fixed to build with -DNDEBUG
 464  *
 465  * Revision 1.7  2006/03/22 10:13:54  strk
 466  * opLinemerge.h split
 467  *
 468  * Revision 1.6  2006/03/21 21:42:54  strk
 469  * planargraph.h header split, planargraph:: classes renamed to match JTS symbols
 470  *
 471  **********************************************************************/