@@ -63,13 +63,21 @@ class RadiallyLessThen {
6363 polarCompare (const Coordinate* o, const Coordinate* p,
6464 const Coordinate* q)
6565 {
66- int orient = Orientation::index (*o, *p, *q);
66+ // To use this comparator in std::sort it must provide a stable ordering, such that
67+ // if cmp(a, b) is true then cmp(b, a) is false. Unfortunately Orientation::index may
68+ // not provide this guarantee when the inputs differ by many orders of magnitude. To
69+ // guard against this, we normalize the order of P and Q before calling OrientationIndex
70+ // and flip the result if the inputs were flipped.
71+ const bool swap = geom::CoordinateLessThan ()(p, q);
72+
73+ const int orient = swap ? Orientation::index (*o, *q, *p) : Orientation::index (*o, *p, *q);
6774
6875 if (orient == Orientation::COUNTERCLOCKWISE) {
69- return 1 ;
76+ return swap ? - 1 : 1 ;
7077 }
78+
7179 if (orient == Orientation::CLOCKWISE) {
72- return -1 ;
80+ return swap ? 1 : -1 ;
7381 }
7482
7583 /* *
@@ -102,7 +110,7 @@ class RadiallyLessThen {
102110 RadiallyLessThen (const Coordinate* c): origin(c) {}
103111
104112 bool
105- operator ()(const Coordinate* p1, const Coordinate* p2)
113+ operator ()(const Coordinate* p1, const Coordinate* p2) const
106114 {
107115 return (polarCompare (origin, p1, p2) == -1 );
108116 }
0 commit comments