LCOV - code coverage report
Current view: top level - home/lbartoletti/qgis_coverage_src/external/poly2tri/sweep - sweep.cc (source / functions) Hit Total Coverage
Test: coverage.info.cleaned Lines: 259 440 58.9 %
Date: 2021-04-10 08:29:14 Functions: 0 0 -
Branches: 0 0 -

           Branch data     Line data    Source code
       1                 :            : /*
       2                 :            :  * Poly2Tri Copyright (c) 2009-2018, Poly2Tri Contributors
       3                 :            :  * https://github.com/jhasse/poly2tri
       4                 :            :  *
       5                 :            :  * All rights reserved.
       6                 :            :  *
       7                 :            :  * Redistribution and use in source and binary forms, with or without modification,
       8                 :            :  * are permitted provided that the following conditions are met:
       9                 :            :  *
      10                 :            :  * * Redistributions of source code must retain the above copyright notice,
      11                 :            :  *   this list of conditions and the following disclaimer.
      12                 :            :  * * Redistributions in binary form must reproduce the above copyright notice,
      13                 :            :  *   this list of conditions and the following disclaimer in the documentation
      14                 :            :  *   and/or other materials provided with the distribution.
      15                 :            :  * * Neither the name of Poly2Tri nor the names of its contributors may be
      16                 :            :  *   used to endorse or promote products derived from this software without specific
      17                 :            :  *   prior written permission.
      18                 :            :  *
      19                 :            :  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
      20                 :            :  * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
      21                 :            :  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
      22                 :            :  * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
      23                 :            :  * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
      24                 :            :  * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
      25                 :            :  * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
      26                 :            :  * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
      27                 :            :  * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
      28                 :            :  * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
      29                 :            :  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
      30                 :            :  */
      31                 :            : #include "sweep.h"
      32                 :            : #include "sweep_context.h"
      33                 :            : #include "advancing_front.h"
      34                 :            : #include "../common/utils.h"
      35                 :            : 
      36                 :            : #include <cassert>
      37                 :            : #include <stdexcept>
      38                 :            : 
      39                 :            : namespace p2t {
      40                 :            : 
      41                 :            : // Triangulate simple polygon with holes
      42                 :         13 : void Sweep::Triangulate(SweepContext& tcx)
      43                 :            : {
      44                 :         13 :   tcx.InitTriangulation();
      45                 :         13 :   tcx.CreateAdvancingFront();
      46                 :            :   // Sweep points; build mesh
      47                 :         13 :   SweepPoints(tcx);
      48                 :            :   // Clean up
      49                 :         13 :   FinalizationPolygon(tcx);
      50                 :         13 : }
      51                 :            : 
      52                 :         13 : void Sweep::SweepPoints(SweepContext& tcx)
      53                 :            : {
      54                 :         91 :   for (size_t i = 1; i < tcx.point_count(); i++) {
      55                 :         78 :     Point& point = *tcx.GetPoint(i);
      56                 :         78 :     Node* node = &PointEvent(tcx, point);
      57                 :        169 :     for (unsigned int i = 0; i < point.edge_list.size(); i++) {
      58                 :         91 :       EdgeEvent(tcx, point.edge_list[i], node);
      59                 :         91 :     }
      60                 :         78 :   }
      61                 :         13 : }
      62                 :            : 
      63                 :         13 : void Sweep::FinalizationPolygon(SweepContext& tcx)
      64                 :            : {
      65                 :            :   // Get an Internal triangle to start with
      66                 :         13 :   Triangle* t = tcx.front()->head()->next->triangle;
      67                 :         13 :   Point* p = tcx.front()->head()->next->point;
      68                 :         13 :   while (!t->GetConstrainedEdgeCW(*p)) {
      69                 :          0 :     t = t->NeighborCCW(*p);
      70                 :            :   }
      71                 :            : 
      72                 :            :   // Collect interior triangles constrained by edges
      73                 :         13 :   tcx.MeshClean(*t);
      74                 :         13 : }
      75                 :            : 
      76                 :         78 : Node& Sweep::PointEvent(SweepContext& tcx, Point& point)
      77                 :            : {
      78                 :         78 :   Node& node = tcx.LocateNode(point);
      79                 :         78 :   Node& new_node = NewFrontTriangle(tcx, point, node);
      80                 :            : 
      81                 :            :   // Only need to check +epsilon since point never have smaller
      82                 :            :   // x value than node due to how we fetch nodes from the front
      83                 :         78 :   if (point.x <= node.point->x + EPSILON) {
      84                 :         39 :     Fill(tcx, node);
      85                 :         39 :   }
      86                 :            : 
      87                 :            :   //tcx.AddNode(new_node);
      88                 :            : 
      89                 :         78 :   FillAdvancingFront(tcx, new_node);
      90                 :         78 :   return new_node;
      91                 :            : }
      92                 :            : 
      93                 :         91 : void Sweep::EdgeEvent(SweepContext& tcx, Edge* edge, Node* node)
      94                 :            : {
      95                 :         91 :   tcx.edge_event.constrained_edge = edge;
      96                 :         91 :   tcx.edge_event.right = (edge->p->x > edge->q->x);
      97                 :            : 
      98                 :         91 :   if (IsEdgeSideOfTriangle(*node->triangle, *edge->p, *edge->q)) {
      99                 :         39 :     return;
     100                 :            :   }
     101                 :            : 
     102                 :            :   // For now we will do all needed filling
     103                 :            :   // TODO: integrate with flip process might give some better performance
     104                 :            :   //       but for now this avoid the issue with cases that needs both flips and fills
     105                 :         52 :   FillEdgeEvent(tcx, edge, node);
     106                 :         52 :   EdgeEvent(tcx, *edge->p, *edge->q, node->triangle, *edge->q);
     107                 :         91 : }
     108                 :            : 
     109                 :        117 : void Sweep::EdgeEvent(SweepContext& tcx, Point& ep, Point& eq, Triangle* triangle, Point& point)
     110                 :            : {
     111                 :        117 :   if (IsEdgeSideOfTriangle(*triangle, ep, eq)) {
     112                 :         52 :     return;
     113                 :            :   }
     114                 :            : 
     115                 :         65 :   Point* p1 = triangle->PointCCW(point);
     116                 :         65 :   Orientation o1 = Orient2d(eq, *p1, ep);
     117                 :         65 :   if (o1 == COLLINEAR) {
     118                 :          0 :     if( triangle->Contains(&eq, p1)) {
     119                 :          0 :       triangle->MarkConstrainedEdge(&eq, p1 );
     120                 :            :       // We are modifying the constraint maybe it would be better to
     121                 :            :       // not change the given constraint and just keep a variable for the new constraint
     122                 :          0 :       tcx.edge_event.constrained_edge->q = p1;
     123                 :          0 :       triangle = &triangle->NeighborAcross(point);
     124                 :          0 :       EdgeEvent( tcx, ep, *p1, triangle, *p1 );
     125                 :          0 :     } else {
     126                 :          0 :       std::runtime_error("EdgeEvent - collinear points not supported");
     127                 :          0 :       assert(0);
     128                 :            :     }
     129                 :          0 :     return;
     130                 :            :   }
     131                 :            : 
     132                 :         65 :   Point* p2 = triangle->PointCW(point);
     133                 :         65 :   Orientation o2 = Orient2d(eq, *p2, ep);
     134                 :         65 :   if (o2 == COLLINEAR) {
     135                 :          0 :     if( triangle->Contains(&eq, p2)) {
     136                 :          0 :       triangle->MarkConstrainedEdge(&eq, p2 );
     137                 :            :       // We are modifying the constraint maybe it would be better to
     138                 :            :       // not change the given constraint and just keep a variable for the new constraint
     139                 :          0 :       tcx.edge_event.constrained_edge->q = p2;
     140                 :          0 :       triangle = &triangle->NeighborAcross(point);
     141                 :          0 :       EdgeEvent( tcx, ep, *p2, triangle, *p2 );
     142                 :          0 :     } else {
     143                 :          0 :       std::runtime_error("EdgeEvent - collinear points not supported");
     144                 :          0 :       assert(0);
     145                 :            :     }
     146                 :          0 :     return;
     147                 :            :   }
     148                 :            : 
     149                 :         65 :   if (o1 == o2) {
     150                 :            :     // Need to decide if we are rotating CW or CCW to get to a triangle
     151                 :            :     // that will cross edge
     152                 :         65 :     if (o1 == CW) {
     153                 :         65 :       triangle = triangle->NeighborCCW(point);
     154                 :         65 :     }       else{
     155                 :          0 :       triangle = triangle->NeighborCW(point);
     156                 :            :     }
     157                 :         65 :     EdgeEvent(tcx, ep, eq, triangle, point);
     158                 :         65 :   } else {
     159                 :            :     // This triangle crosses constraint so lets flippin start!
     160                 :          0 :     FlipEdgeEvent(tcx, ep, eq, triangle, point);
     161                 :            :   }
     162                 :        117 : }
     163                 :            : 
     164                 :        208 : bool Sweep::IsEdgeSideOfTriangle(Triangle& triangle, Point& ep, Point& eq)
     165                 :            : {
     166                 :        208 :   const int index = triangle.EdgeIndex(&ep, &eq);
     167                 :            : 
     168                 :        208 :   if (index != -1) {
     169                 :         91 :     triangle.MarkConstrainedEdge(index);
     170                 :         91 :     Triangle* t = triangle.GetNeighbor(index);
     171                 :         91 :     if (t) {
     172                 :         39 :       t->MarkConstrainedEdge(&ep, &eq);
     173                 :         39 :     }
     174                 :         91 :     return true;
     175                 :            :   }
     176                 :        117 :   return false;
     177                 :        208 : }
     178                 :            : 
     179                 :         78 : Node& Sweep::NewFrontTriangle(SweepContext& tcx, Point& point, Node& node)
     180                 :            : {
     181                 :         78 :   Triangle* triangle = new Triangle(point, *node.point, *node.next->point);
     182                 :            : 
     183                 :         78 :   triangle->MarkNeighbor(*node.triangle);
     184                 :         78 :   tcx.AddToMap(triangle);
     185                 :            : 
     186                 :         78 :   Node* new_node = new Node(point);
     187                 :         78 :   nodes_.push_back(new_node);
     188                 :            : 
     189                 :         78 :   new_node->next = node.next;
     190                 :         78 :   new_node->prev = &node;
     191                 :         78 :   node.next->prev = new_node;
     192                 :         78 :   node.next = new_node;
     193                 :            : 
     194                 :         78 :   if (!Legalize(tcx, *triangle)) {
     195                 :         65 :     tcx.MapTriangleToNodes(*triangle);
     196                 :         65 :   }
     197                 :            : 
     198                 :         78 :   return *new_node;
     199                 :          0 : }
     200                 :            : 
     201                 :         65 : void Sweep::Fill(SweepContext& tcx, Node& node)
     202                 :            : {
     203                 :         65 :   Triangle* triangle = new Triangle(*node.prev->point, *node.point, *node.next->point);
     204                 :            : 
     205                 :            :   // TODO: should copy the constrained_edge value from neighbor triangles
     206                 :            :   //       for now constrained_edge values are copied during the legalize
     207                 :         65 :   triangle->MarkNeighbor(*node.prev->triangle);
     208                 :         65 :   triangle->MarkNeighbor(*node.triangle);
     209                 :            : 
     210                 :         65 :   tcx.AddToMap(triangle);
     211                 :            : 
     212                 :            :   // Update the advancing front
     213                 :         65 :   node.prev->next = node.next;
     214                 :         65 :   node.next->prev = node.prev;
     215                 :            : 
     216                 :            :   // If it was legalized the triangle has already been mapped
     217                 :         65 :   if (!Legalize(tcx, *triangle)) {
     218                 :         65 :     tcx.MapTriangleToNodes(*triangle);
     219                 :         65 :   }
     220                 :            : 
     221                 :         65 : }
     222                 :            : 
     223                 :         78 : void Sweep::FillAdvancingFront(SweepContext& tcx, Node& n)
     224                 :            : {
     225                 :            : 
     226                 :            :   // Fill right holes
     227                 :         78 :   Node* node = n.next;
     228                 :            : 
     229                 :         91 :   while (node->next) {
     230                 :            :     // if HoleAngle exceeds 90 degrees then break.
     231                 :         65 :     if (LargeHole_DontFill(node)) break;
     232                 :         13 :     Fill(tcx, *node);
     233                 :         13 :     node = node->next;
     234                 :            :   }
     235                 :            : 
     236                 :            :   // Fill left holes
     237                 :         78 :   node = n.prev;
     238                 :            : 
     239                 :         78 :   while (node->prev) {
     240                 :            :     // if HoleAngle exceeds 90 degrees then break.
     241                 :         65 :     if (LargeHole_DontFill(node)) break;
     242                 :          0 :     Fill(tcx, *node);
     243                 :          0 :     node = node->prev;
     244                 :            :   }
     245                 :            : 
     246                 :            :   // Fill right basins
     247                 :         78 :   if (n.next && n.next->next) {
     248                 :         52 :     const double angle = BasinAngle(n);
     249                 :         52 :     if (angle < PI_3div4) {
     250                 :         13 :       FillBasin(tcx, n);
     251                 :         13 :     }
     252                 :         52 :   }
     253                 :         78 : }
     254                 :            : 
     255                 :            : // True if HoleAngle exceeds 90 degrees.
     256                 :        130 : bool Sweep::LargeHole_DontFill(const Node* node) const {
     257                 :            : 
     258                 :        130 :   const Node* nextNode = node->next;
     259                 :        130 :   const Node* prevNode = node->prev;
     260                 :        130 :   if (!AngleExceeds90Degrees(node->point, nextNode->point, prevNode->point))
     261                 :         13 :           return false;
     262                 :            : 
     263                 :            :   // Check additional points on front.
     264                 :        117 :   const Node* next2Node = nextNode->next;
     265                 :            :   // "..Plus.." because only want angles on same side as point being added.
     266                 :        117 :   if ((next2Node != NULL) && !AngleExceedsPlus90DegreesOrIsNegative(node->point, next2Node->point, prevNode->point))
     267                 :          0 :           return false;
     268                 :            : 
     269                 :        117 :   const Node* prev2Node = prevNode->prev;
     270                 :            :   // "..Plus.." because only want angles on same side as point being added.
     271                 :        117 :   if ((prev2Node != NULL) && !AngleExceedsPlus90DegreesOrIsNegative(node->point, nextNode->point, prev2Node->point))
     272                 :          0 :           return false;
     273                 :            : 
     274                 :        117 :   return true;
     275                 :        130 : }
     276                 :            : 
     277                 :        130 : bool Sweep::AngleExceeds90Degrees(const Point* origin, const Point* pa, const Point* pb) const {
     278                 :        130 :   const double angle = Angle(origin, pa, pb);
     279                 :        130 :   return ((angle > PI_div2) || (angle < -PI_div2));
     280                 :            : }
     281                 :            : 
     282                 :        169 : bool Sweep::AngleExceedsPlus90DegreesOrIsNegative(const Point* origin, const Point* pa, const Point* pb) const {
     283                 :        169 :   const double angle = Angle(origin, pa, pb);
     284                 :        169 :   return (angle > PI_div2) || (angle < 0);
     285                 :            : }
     286                 :            : 
     287                 :        299 : double Sweep::Angle(const Point* origin, const Point* pa, const Point* pb) const {
     288                 :            :   /* Complex plane
     289                 :            :    * ab = cosA +i*sinA
     290                 :            :    * ab = (ax + ay*i)(bx + by*i) = (ax*bx + ay*by) + i(ax*by-ay*bx)
     291                 :            :    * atan2(y,x) computes the principal value of the argument function
     292                 :            :    * applied to the complex number x+iy
     293                 :            :    * Where x = ax*bx + ay*by
     294                 :            :    *       y = ax*by - ay*bx
     295                 :            :    */
     296                 :        299 :   const double px = origin->x;
     297                 :        299 :   const double py = origin->y;
     298                 :        299 :   const double ax = pa->x- px;
     299                 :        299 :   const double ay = pa->y - py;
     300                 :        299 :   const double bx = pb->x - px;
     301                 :        299 :   const double by = pb->y - py;
     302                 :        299 :   const double x = ax * by - ay * bx;
     303                 :        299 :   const double y = ax * bx + ay * by;
     304                 :        299 :   return atan2(x, y);
     305                 :            : }
     306                 :            : 
     307                 :         52 : double Sweep::BasinAngle(const Node& node) const
     308                 :            : {
     309                 :         52 :   const double ax = node.point->x - node.next->next->point->x;
     310                 :         52 :   const double ay = node.point->y - node.next->next->point->y;
     311                 :         52 :   return atan2(ay, ax);
     312                 :            : }
     313                 :            : 
     314                 :          0 : double Sweep::HoleAngle(const Node& node) const
     315                 :            : {
     316                 :            :   /* Complex plane
     317                 :            :    * ab = cosA +i*sinA
     318                 :            :    * ab = (ax + ay*i)(bx + by*i) = (ax*bx + ay*by) + i(ax*by-ay*bx)
     319                 :            :    * atan2(y,x) computes the principal value of the argument function
     320                 :            :    * applied to the complex number x+iy
     321                 :            :    * Where x = ax*bx + ay*by
     322                 :            :    *       y = ax*by - ay*bx
     323                 :            :    */
     324                 :          0 :   const double ax = node.next->point->x - node.point->x;
     325                 :          0 :   const double ay = node.next->point->y - node.point->y;
     326                 :          0 :   const double bx = node.prev->point->x - node.point->x;
     327                 :          0 :   const double by = node.prev->point->y - node.point->y;
     328                 :          0 :   return atan2(ax * by - ay * bx, ax * bx + ay * by);
     329                 :            : }
     330                 :            : 
     331                 :        169 : bool Sweep::Legalize(SweepContext& tcx, Triangle& t)
     332                 :            : {
     333                 :            :   // To legalize a triangle we start by finding if any of the three edges
     334                 :            :   // violate the Delaunay condition
     335                 :        637 :   for (int i = 0; i < 3; i++) {
     336                 :        481 :     if (t.delaunay_edge[i])
     337                 :         26 :       continue;
     338                 :            : 
     339                 :        455 :     Triangle* ot = t.GetNeighbor(i);
     340                 :            : 
     341                 :        455 :     if (ot) {
     342                 :        221 :       Point* p = t.GetPoint(i);
     343                 :        221 :       Point* op = ot->OppositePoint(t, *p);
     344                 :        221 :       int oi = ot->Index(op);
     345                 :            : 
     346                 :            :       // If this is a Constrained Edge or a Delaunay Edge(only during recursive legalization)
     347                 :            :       // then we should not try to legalize
     348                 :        221 :       if (ot->constrained_edge[oi] || ot->delaunay_edge[oi]) {
     349                 :         52 :         t.constrained_edge[i] = ot->constrained_edge[oi];
     350                 :         52 :         continue;
     351                 :            :       }
     352                 :            : 
     353                 :        169 :       bool inside = Incircle(*p, *t.PointCCW(*p), *t.PointCW(*p), *op);
     354                 :            : 
     355                 :        169 :       if (inside) {
     356                 :            :         // Lets mark this shared edge as Delaunay
     357                 :         13 :         t.delaunay_edge[i] = true;
     358                 :         13 :         ot->delaunay_edge[oi] = true;
     359                 :            : 
     360                 :            :         // Lets rotate shared edge one vertex CW to legalize it
     361                 :         13 :         RotateTrianglePair(t, *p, *ot, *op);
     362                 :            : 
     363                 :            :         // We now got one valid Delaunay Edge shared by two triangles
     364                 :            :         // This gives us 4 new edges to check for Delaunay
     365                 :            : 
     366                 :            :         // Make sure that triangle to node mapping is done only one time for a specific triangle
     367                 :         13 :         bool not_legalized = !Legalize(tcx, t);
     368                 :         13 :         if (not_legalized) {
     369                 :         13 :           tcx.MapTriangleToNodes(t);
     370                 :         13 :         }
     371                 :            : 
     372                 :         13 :         not_legalized = !Legalize(tcx, *ot);
     373                 :         13 :         if (not_legalized)
     374                 :         13 :           tcx.MapTriangleToNodes(*ot);
     375                 :            : 
     376                 :            :         // Reset the Delaunay edges, since they only are valid Delaunay edges
     377                 :            :         // until we add a new triangle or point.
     378                 :            :         // XXX: need to think about this. Can these edges be tried after we
     379                 :            :         //      return to previous recursive level?
     380                 :         13 :         t.delaunay_edge[i] = false;
     381                 :         13 :         ot->delaunay_edge[oi] = false;
     382                 :            : 
     383                 :            :         // If triangle have been legalized no need to check the other edges since
     384                 :            :         // the recursive legalization will handles those so we can end here.
     385                 :         13 :         return true;
     386                 :            :       }
     387                 :        156 :     }
     388                 :        390 :   }
     389                 :        156 :   return false;
     390                 :        169 : }
     391                 :            : 
     392                 :        169 : bool Sweep::Incircle(const Point& pa, const Point& pb, const Point& pc, const Point& pd) const
     393                 :            : {
     394                 :        169 :   const double adx = pa.x - pd.x;
     395                 :        169 :   const double ady = pa.y - pd.y;
     396                 :        169 :   const double bdx = pb.x - pd.x;
     397                 :        169 :   const double bdy = pb.y - pd.y;
     398                 :            : 
     399                 :        169 :   const double adxbdy = adx * bdy;
     400                 :        169 :   const double bdxady = bdx * ady;
     401                 :        169 :   const double oabd = adxbdy - bdxady;
     402                 :            : 
     403                 :        169 :   if (oabd <= 0)
     404                 :         65 :     return false;
     405                 :            : 
     406                 :        104 :   const double cdx = pc.x - pd.x;
     407                 :        104 :   const double cdy = pc.y - pd.y;
     408                 :            : 
     409                 :        104 :   const double cdxady = cdx * ady;
     410                 :        104 :   const double adxcdy = adx * cdy;
     411                 :        104 :   const double ocad = cdxady - adxcdy;
     412                 :            : 
     413                 :        104 :   if (ocad <= 0)
     414                 :         52 :     return false;
     415                 :            : 
     416                 :         52 :   const double bdxcdy = bdx * cdy;
     417                 :         52 :   const double cdxbdy = cdx * bdy;
     418                 :            : 
     419                 :         52 :   const double alift = adx * adx + ady * ady;
     420                 :         52 :   const double blift = bdx * bdx + bdy * bdy;
     421                 :         52 :   const double clift = cdx * cdx + cdy * cdy;
     422                 :            : 
     423                 :         52 :   const double det = alift * (bdxcdy - cdxbdy) + blift * ocad + clift * oabd;
     424                 :            : 
     425                 :         52 :   return det > 0;
     426                 :        169 : }
     427                 :            : 
     428                 :         13 : void Sweep::RotateTrianglePair(Triangle& t, Point& p, Triangle& ot, Point& op) const
     429                 :            : {
     430                 :            :   Triangle* n1, *n2, *n3, *n4;
     431                 :         13 :   n1 = t.NeighborCCW(p);
     432                 :         13 :   n2 = t.NeighborCW(p);
     433                 :         13 :   n3 = ot.NeighborCCW(op);
     434                 :         13 :   n4 = ot.NeighborCW(op);
     435                 :            : 
     436                 :            :   bool ce1, ce2, ce3, ce4;
     437                 :         13 :   ce1 = t.GetConstrainedEdgeCCW(p);
     438                 :         13 :   ce2 = t.GetConstrainedEdgeCW(p);
     439                 :         13 :   ce3 = ot.GetConstrainedEdgeCCW(op);
     440                 :         13 :   ce4 = ot.GetConstrainedEdgeCW(op);
     441                 :            : 
     442                 :            :   bool de1, de2, de3, de4;
     443                 :         13 :   de1 = t.GetDelunayEdgeCCW(p);
     444                 :         13 :   de2 = t.GetDelunayEdgeCW(p);
     445                 :         13 :   de3 = ot.GetDelunayEdgeCCW(op);
     446                 :         13 :   de4 = ot.GetDelunayEdgeCW(op);
     447                 :            : 
     448                 :         13 :   t.Legalize(p, op);
     449                 :         13 :   ot.Legalize(op, p);
     450                 :            : 
     451                 :            :   // Remap delaunay_edge
     452                 :         13 :   ot.SetDelunayEdgeCCW(p, de1);
     453                 :         13 :   t.SetDelunayEdgeCW(p, de2);
     454                 :         13 :   t.SetDelunayEdgeCCW(op, de3);
     455                 :         13 :   ot.SetDelunayEdgeCW(op, de4);
     456                 :            : 
     457                 :            :   // Remap constrained_edge
     458                 :         13 :   ot.SetConstrainedEdgeCCW(p, ce1);
     459                 :         13 :   t.SetConstrainedEdgeCW(p, ce2);
     460                 :         13 :   t.SetConstrainedEdgeCCW(op, ce3);
     461                 :         13 :   ot.SetConstrainedEdgeCW(op, ce4);
     462                 :            : 
     463                 :            :   // Remap neighbors
     464                 :            :   // XXX: might optimize the markNeighbor by keeping track of
     465                 :            :   //      what side should be assigned to what neighbor after the
     466                 :            :   //      rotation. Now mark neighbor does lots of testing to find
     467                 :            :   //      the right side.
     468                 :         13 :   t.ClearNeighbors();
     469                 :         13 :   ot.ClearNeighbors();
     470                 :         13 :   if (n1) ot.MarkNeighbor(*n1);
     471                 :         13 :   if (n2) t.MarkNeighbor(*n2);
     472                 :         13 :   if (n3) t.MarkNeighbor(*n3);
     473                 :         13 :   if (n4) ot.MarkNeighbor(*n4);
     474                 :         13 :   t.MarkNeighbor(ot);
     475                 :         13 : }
     476                 :            : 
     477                 :         13 : void Sweep::FillBasin(SweepContext& tcx, Node& node)
     478                 :            : {
     479                 :         13 :   if (Orient2d(*node.point, *node.next->point, *node.next->next->point) == CCW) {
     480                 :         13 :     tcx.basin.left_node = node.next->next;
     481                 :         13 :   } else {
     482                 :          0 :     tcx.basin.left_node = node.next;
     483                 :            :   }
     484                 :            : 
     485                 :            :   // Find the bottom and right node
     486                 :         13 :   tcx.basin.bottom_node = tcx.basin.left_node;
     487                 :         13 :   while (tcx.basin.bottom_node->next
     488                 :         13 :          && tcx.basin.bottom_node->point->y >= tcx.basin.bottom_node->next->point->y) {
     489                 :          0 :     tcx.basin.bottom_node = tcx.basin.bottom_node->next;
     490                 :            :   }
     491                 :         13 :   if (tcx.basin.bottom_node == tcx.basin.left_node) {
     492                 :            :     // No valid basin
     493                 :         13 :     return;
     494                 :            :   }
     495                 :            : 
     496                 :          0 :   tcx.basin.right_node = tcx.basin.bottom_node;
     497                 :          0 :   while (tcx.basin.right_node->next
     498                 :          0 :          && tcx.basin.right_node->point->y < tcx.basin.right_node->next->point->y) {
     499                 :          0 :     tcx.basin.right_node = tcx.basin.right_node->next;
     500                 :            :   }
     501                 :          0 :   if (tcx.basin.right_node == tcx.basin.bottom_node) {
     502                 :            :     // No valid basins
     503                 :          0 :     return;
     504                 :            :   }
     505                 :            : 
     506                 :          0 :   tcx.basin.width = tcx.basin.right_node->point->x - tcx.basin.left_node->point->x;
     507                 :          0 :   tcx.basin.left_highest = tcx.basin.left_node->point->y > tcx.basin.right_node->point->y;
     508                 :            : 
     509                 :          0 :   FillBasinReq(tcx, tcx.basin.bottom_node);
     510                 :         13 : }
     511                 :            : 
     512                 :          0 : void Sweep::FillBasinReq(SweepContext& tcx, Node* node)
     513                 :            : {
     514                 :            :   // if shallow stop filling
     515                 :          0 :   if (IsShallow(tcx, *node)) {
     516                 :          0 :     return;
     517                 :            :   }
     518                 :            : 
     519                 :          0 :   Fill(tcx, *node);
     520                 :            : 
     521                 :          0 :   if (node->prev == tcx.basin.left_node && node->next == tcx.basin.right_node) {
     522                 :          0 :     return;
     523                 :          0 :   } else if (node->prev == tcx.basin.left_node) {
     524                 :          0 :     Orientation o = Orient2d(*node->point, *node->next->point, *node->next->next->point);
     525                 :          0 :     if (o == CW) {
     526                 :          0 :       return;
     527                 :            :     }
     528                 :          0 :     node = node->next;
     529                 :          0 :   } else if (node->next == tcx.basin.right_node) {
     530                 :          0 :     Orientation o = Orient2d(*node->point, *node->prev->point, *node->prev->prev->point);
     531                 :          0 :     if (o == CCW) {
     532                 :          0 :       return;
     533                 :            :     }
     534                 :          0 :     node = node->prev;
     535                 :          0 :   } else {
     536                 :            :     // Continue with the neighbor node with lowest Y value
     537                 :          0 :     if (node->prev->point->y < node->next->point->y) {
     538                 :          0 :       node = node->prev;
     539                 :          0 :     } else {
     540                 :          0 :       node = node->next;
     541                 :            :     }
     542                 :            :   }
     543                 :            : 
     544                 :          0 :   FillBasinReq(tcx, node);
     545                 :          0 : }
     546                 :            : 
     547                 :          0 : bool Sweep::IsShallow(SweepContext& tcx, Node& node)
     548                 :            : {
     549                 :            :   double height;
     550                 :            : 
     551                 :          0 :   if (tcx.basin.left_highest) {
     552                 :          0 :     height = tcx.basin.left_node->point->y - node.point->y;
     553                 :          0 :   } else {
     554                 :          0 :     height = tcx.basin.right_node->point->y - node.point->y;
     555                 :            :   }
     556                 :            : 
     557                 :            :   // if shallow stop filling
     558                 :          0 :   if (tcx.basin.width > height) {
     559                 :          0 :     return true;
     560                 :            :   }
     561                 :          0 :   return false;
     562                 :          0 : }
     563                 :            : 
     564                 :         52 : void Sweep::FillEdgeEvent(SweepContext& tcx, Edge* edge, Node* node)
     565                 :            : {
     566                 :         52 :   if (tcx.edge_event.right) {
     567                 :          0 :     FillRightAboveEdgeEvent(tcx, edge, node);
     568                 :          0 :   } else {
     569                 :         52 :     FillLeftAboveEdgeEvent(tcx, edge, node);
     570                 :            :   }
     571                 :         52 : }
     572                 :            : 
     573                 :          0 : void Sweep::FillRightAboveEdgeEvent(SweepContext& tcx, Edge* edge, Node* node)
     574                 :            : {
     575                 :          0 :   while (node->next->point->x < edge->p->x) {
     576                 :            :     // Check if next node is below the edge
     577                 :          0 :     if (Orient2d(*edge->q, *node->next->point, *edge->p) == CCW) {
     578                 :          0 :       FillRightBelowEdgeEvent(tcx, edge, *node);
     579                 :          0 :     } else {
     580                 :          0 :       node = node->next;
     581                 :            :     }
     582                 :            :   }
     583                 :          0 : }
     584                 :            : 
     585                 :          0 : void Sweep::FillRightBelowEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
     586                 :            : {
     587                 :          0 :   if (node.point->x < edge->p->x) {
     588                 :          0 :     if (Orient2d(*node.point, *node.next->point, *node.next->next->point) == CCW) {
     589                 :            :       // Concave
     590                 :          0 :       FillRightConcaveEdgeEvent(tcx, edge, node);
     591                 :          0 :     } else{
     592                 :            :       // Convex
     593                 :          0 :       FillRightConvexEdgeEvent(tcx, edge, node);
     594                 :            :       // Retry this one
     595                 :          0 :       FillRightBelowEdgeEvent(tcx, edge, node);
     596                 :            :     }
     597                 :          0 :   }
     598                 :          0 : }
     599                 :            : 
     600                 :          0 : void Sweep::FillRightConcaveEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
     601                 :            : {
     602                 :          0 :   Fill(tcx, *node.next);
     603                 :          0 :   if (node.next->point != edge->p) {
     604                 :            :     // Next above or below edge?
     605                 :          0 :     if (Orient2d(*edge->q, *node.next->point, *edge->p) == CCW) {
     606                 :            :       // Below
     607                 :          0 :       if (Orient2d(*node.point, *node.next->point, *node.next->next->point) == CCW) {
     608                 :            :         // Next is concave
     609                 :          0 :         FillRightConcaveEdgeEvent(tcx, edge, node);
     610                 :          0 :       } else {
     611                 :            :         // Next is convex
     612                 :            :       }
     613                 :          0 :     }
     614                 :          0 :   }
     615                 :            : 
     616                 :          0 : }
     617                 :            : 
     618                 :          0 : void Sweep::FillRightConvexEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
     619                 :            : {
     620                 :            :   // Next concave or convex?
     621                 :          0 :   if (Orient2d(*node.next->point, *node.next->next->point, *node.next->next->next->point) == CCW) {
     622                 :            :     // Concave
     623                 :          0 :     FillRightConcaveEdgeEvent(tcx, edge, *node.next);
     624                 :          0 :   } else{
     625                 :            :     // Convex
     626                 :            :     // Next above or below edge?
     627                 :          0 :     if (Orient2d(*edge->q, *node.next->next->point, *edge->p) == CCW) {
     628                 :            :       // Below
     629                 :          0 :       FillRightConvexEdgeEvent(tcx, edge, *node.next);
     630                 :          0 :     } else{
     631                 :            :       // Above
     632                 :            :     }
     633                 :            :   }
     634                 :          0 : }
     635                 :            : 
     636                 :         52 : void Sweep::FillLeftAboveEdgeEvent(SweepContext& tcx, Edge* edge, Node* node)
     637                 :            : {
     638                 :         65 :   while (node->prev->point->x > edge->p->x) {
     639                 :            :     // Check if next node is below the edge
     640                 :         13 :     if (Orient2d(*edge->q, *node->prev->point, *edge->p) == CW) {
     641                 :         13 :       FillLeftBelowEdgeEvent(tcx, edge, *node);
     642                 :         13 :     } else {
     643                 :          0 :       node = node->prev;
     644                 :            :     }
     645                 :            :   }
     646                 :         52 : }
     647                 :            : 
     648                 :         13 : void Sweep::FillLeftBelowEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
     649                 :            : {
     650                 :         13 :   if (node.point->x > edge->p->x) {
     651                 :         13 :     if (Orient2d(*node.point, *node.prev->point, *node.prev->prev->point) == CW) {
     652                 :            :       // Concave
     653                 :         13 :       FillLeftConcaveEdgeEvent(tcx, edge, node);
     654                 :         13 :     } else {
     655                 :            :       // Convex
     656                 :          0 :       FillLeftConvexEdgeEvent(tcx, edge, node);
     657                 :            :       // Retry this one
     658                 :          0 :       FillLeftBelowEdgeEvent(tcx, edge, node);
     659                 :            :     }
     660                 :         13 :   }
     661                 :         13 : }
     662                 :            : 
     663                 :          0 : void Sweep::FillLeftConvexEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
     664                 :            : {
     665                 :            :   // Next concave or convex?
     666                 :          0 :   if (Orient2d(*node.prev->point, *node.prev->prev->point, *node.prev->prev->prev->point) == CW) {
     667                 :            :     // Concave
     668                 :          0 :     FillLeftConcaveEdgeEvent(tcx, edge, *node.prev);
     669                 :          0 :   } else{
     670                 :            :     // Convex
     671                 :            :     // Next above or below edge?
     672                 :          0 :     if (Orient2d(*edge->q, *node.prev->prev->point, *edge->p) == CW) {
     673                 :            :       // Below
     674                 :          0 :       FillLeftConvexEdgeEvent(tcx, edge, *node.prev);
     675                 :          0 :     } else{
     676                 :            :       // Above
     677                 :            :     }
     678                 :            :   }
     679                 :          0 : }
     680                 :            : 
     681                 :         13 : void Sweep::FillLeftConcaveEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
     682                 :            : {
     683                 :         13 :   Fill(tcx, *node.prev);
     684                 :         13 :   if (node.prev->point != edge->p) {
     685                 :            :     // Next above or below edge?
     686                 :          0 :     if (Orient2d(*edge->q, *node.prev->point, *edge->p) == CW) {
     687                 :            :       // Below
     688                 :          0 :       if (Orient2d(*node.point, *node.prev->point, *node.prev->prev->point) == CW) {
     689                 :            :         // Next is concave
     690                 :          0 :         FillLeftConcaveEdgeEvent(tcx, edge, node);
     691                 :          0 :       } else{
     692                 :            :         // Next is convex
     693                 :            :       }
     694                 :          0 :     }
     695                 :          0 :   }
     696                 :            : 
     697                 :         13 : }
     698                 :            : 
     699                 :          0 : void Sweep::FlipEdgeEvent(SweepContext& tcx, Point& ep, Point& eq, Triangle* t, Point& p)
     700                 :            : {
     701                 :          0 :   Triangle& ot = t->NeighborAcross(p);
     702                 :          0 :   Point& op = *ot.OppositePoint(*t, p);
     703                 :            : 
     704                 :          0 :   if (InScanArea(p, *t->PointCCW(p), *t->PointCW(p), op)) {
     705                 :            :     // Lets rotate shared edge one vertex CW
     706                 :          0 :     RotateTrianglePair(*t, p, ot, op);
     707                 :          0 :     tcx.MapTriangleToNodes(*t);
     708                 :          0 :     tcx.MapTriangleToNodes(ot);
     709                 :            : 
     710                 :          0 :     if (p == eq && op == ep) {
     711                 :          0 :       if (eq == *tcx.edge_event.constrained_edge->q && ep == *tcx.edge_event.constrained_edge->p) {
     712                 :          0 :         t->MarkConstrainedEdge(&ep, &eq);
     713                 :          0 :         ot.MarkConstrainedEdge(&ep, &eq);
     714                 :          0 :         Legalize(tcx, *t);
     715                 :          0 :         Legalize(tcx, ot);
     716                 :          0 :       } else {
     717                 :            :         // XXX: I think one of the triangles should be legalized here?
     718                 :            :       }
     719                 :          0 :     } else {
     720                 :          0 :       Orientation o = Orient2d(eq, op, ep);
     721                 :          0 :       t = &NextFlipTriangle(tcx, (int)o, *t, ot, p, op);
     722                 :          0 :       FlipEdgeEvent(tcx, ep, eq, t, p);
     723                 :            :     }
     724                 :          0 :   } else {
     725                 :          0 :     Point& newP = NextFlipPoint(ep, eq, ot, op);
     726                 :          0 :     FlipScanEdgeEvent(tcx, ep, eq, *t, ot, newP);
     727                 :          0 :     EdgeEvent(tcx, ep, eq, t, p);
     728                 :            :   }
     729                 :          0 : }
     730                 :            : 
     731                 :          0 : Triangle& Sweep::NextFlipTriangle(SweepContext& tcx, int o, Triangle& t, Triangle& ot, Point& p, Point& op)
     732                 :            : {
     733                 :          0 :   if (o == CCW) {
     734                 :            :     // ot is not crossing edge after flip
     735                 :          0 :     int edge_index = ot.EdgeIndex(&p, &op);
     736                 :          0 :     ot.delaunay_edge[edge_index] = true;
     737                 :          0 :     Legalize(tcx, ot);
     738                 :          0 :     ot.ClearDelunayEdges();
     739                 :          0 :     return t;
     740                 :            :   }
     741                 :            : 
     742                 :            :   // t is not crossing edge after flip
     743                 :          0 :   int edge_index = t.EdgeIndex(&p, &op);
     744                 :            : 
     745                 :          0 :   t.delaunay_edge[edge_index] = true;
     746                 :          0 :   Legalize(tcx, t);
     747                 :          0 :   t.ClearDelunayEdges();
     748                 :          0 :   return ot;
     749                 :          0 : }
     750                 :            : 
     751                 :          0 : Point& Sweep::NextFlipPoint(Point& ep, Point& eq, Triangle& ot, Point& op)
     752                 :            : {
     753                 :          0 :   Orientation o2d = Orient2d(eq, op, ep);
     754                 :          0 :   if (o2d == CW) {
     755                 :            :     // Right
     756                 :          0 :     return *ot.PointCCW(op);
     757                 :          0 :   } else if (o2d == CCW) {
     758                 :            :     // Left
     759                 :          0 :     return *ot.PointCW(op);
     760                 :            :   }
     761                 :          0 :   throw std::runtime_error("[Unsupported] Opposing point on constrained edge");
     762                 :          0 : }
     763                 :            : 
     764                 :          0 : void Sweep::FlipScanEdgeEvent(SweepContext& tcx, Point& ep, Point& eq, Triangle& flip_triangle,
     765                 :            :                               Triangle& t, Point& p)
     766                 :            : {
     767                 :          0 :   Triangle& ot = t.NeighborAcross(p);
     768                 :          0 :   Point& op = *ot.OppositePoint(t, p);
     769                 :            : 
     770                 :          0 :   if (InScanArea(eq, *flip_triangle.PointCCW(eq), *flip_triangle.PointCW(eq), op)) {
     771                 :            :     // flip with new edge op->eq
     772                 :          0 :     FlipEdgeEvent(tcx, eq, op, &ot, op);
     773                 :            :     // TODO: Actually I just figured out that it should be possible to
     774                 :            :     //       improve this by getting the next ot and op before the the above
     775                 :            :     //       flip and continue the flipScanEdgeEvent here
     776                 :            :     // set new ot and op here and loop back to inScanArea test
     777                 :            :     // also need to set a new flip_triangle first
     778                 :            :     // Turns out at first glance that this is somewhat complicated
     779                 :            :     // so it will have to wait.
     780                 :          0 :   } else{
     781                 :          0 :     Point& newP = NextFlipPoint(ep, eq, ot, op);
     782                 :          0 :     FlipScanEdgeEvent(tcx, ep, eq, flip_triangle, ot, newP);
     783                 :            :   }
     784                 :          0 : }
     785                 :            : 
     786                 :         13 : Sweep::~Sweep() {
     787                 :            : 
     788                 :            :     // Clean up memory
     789                 :         91 :     for(size_t i = 0; i < nodes_.size(); i++) {
     790                 :         78 :         delete nodes_[i];
     791                 :         78 :     }
     792                 :            : 
     793                 :         13 : }
     794                 :            : 
     795                 :            : }
     796                 :            : 

Generated by: LCOV version 1.14