Alpha_Q's blog

How do you polar sort?
By Alpha_Q, history, 5 years ago, In English

Say yo want a set of points sorted in some counterclockwise order. What is the best (fast, precise, short) way to do it? The atan2l function is precise and short but way too slow.

Update: One of the previous entries were wrong. I like the following one which hopefully works. 

inline bool up (point p) {
  return p.y > 0 or (p.y == 0 and p.x >= 0);
}

sort(v.begin(), v.end(), [] (point a, point b) {
  return up(a) == up(b) ? a.x * b.y > a.y * b.x : up(a) < up(b);
});
  • Vote: I like it
  • +141
  • Vote: I do not like it

Comments Comments (17)
| Write comment?
»
5 years ago, # |
  Vote: I like it -18 Vote: I do not like it

if you look at each point as a complex number, then imagine you wana sort points by counterclockwise order from P, your sort can be like this:

complex<ld> a[N];
sort(a, a + n, [](complex<ld> x, complex<ld> y) {
    return arg(x - p) < arg(y - p);
});

if you wana know how does it work search about complex in Codeforces might be helpful. I know 2 of them first_one introduce some functions about complex and second_one is about FFT but has useful parts about complex too.

→ Reply
»
5 years ago, # |
  Vote: I like it +26 Vote: I do not like it

This compare function won't change order of the points on X-axis (y=0). Try this: (-1, 0), (1, 0),(-2,0), (2,0), (0,3).

→ Reply
»
5 years ago, # |
  Vote: I like it +18 Vote: I do not like it

(1, 0) == (-1, 0) according to your comparator. But I'm actually interested if there's an elegant way to do this as well.

→ Reply
»
5 years ago, # |
  Vote: I like it +29 Vote: I do not like it

Good points (no pun intended), thanks FrozenBlood and mmaxio. I guess this one works only when no three points are collinear. Nevertheless, I used a different version before.

inline int quad (point p) {
  if (p.x < 0 and p.y < 0) return 0;
  if (p.x >= 0 and p.y < 0) return 1;
  if (p.x >= 0 and p.y >= 0) return 2;
  if (p.x < 0 and p.y >= 0) return 3;
  assert(69 == 420);
}

sort(v.begin(), v.end(), [] (point a, point b) {
  return quad(a) == quad(b) ? a.x * b.y > a.y * b.x : quad(a) < quad(b);
});

This one seems to be correct for all inputs (although a bit longer). So the question remains.

→ Reply
  • »
    »
    5 years ago, # ^ |
      Vote: I like it +54 Vote: I do not like it

    If you want shorter: 

    int ret[2][2] = {{0, 3},{1, 2}}; 
inline int quad(point p) { 
    return ret[p.x>=0][p.y>=0];  
}
    → Reply
    • »
      »
      »
      5 years ago, # ^ |
      ← Rev. 4   Vote: I like it +71 Vote: I do not like it

      I mean, if we're going to golf it ..

      int quad(int x, int y){
   return ((x>=0)^(y>=0))|((y>=0)<<1);
}

      EDIT: What about this cursed version relying on operator precedence rules?

      int quad(int x, int y){
   return x>0^y>0|2*(y>0);
}
      → Reply
»
5 years ago, # |
  Vote: I like it 0 Vote: I do not like it

Auto comment: topic has been updated by Alpha_Q (previous revision, new revision, compare).

→ Reply
»
5 years ago, # |
← Rev. 2   Vote: I like it +43 Vote: I do not like it

I like to use something very similar to your code. First I divide the plane into two half-open $$$180^\circ$$$ angles (by comparing lexicographically with origin), then use the cross product directly. All the relative angles will be smaller than $$$180^\circ$$$, so it'll work.

// in point struct
  long long CrossProd(const point &p) const {
    return x * (long long)p.y - y * (long long)p.x;
  }
  bool operator<(const point &p) const {
    return make_pair(x, y) < make_pair(p.x, p.y);
  }

// when sorting
sort(v.begin(), v.end(), [] (point a, point b) {
  const point origin{0, 0};
  bool ba = a < origin, bb = b < origin;
  if (ba != bb) { return ba < bb; }
  return a.CrossProd(b) > 0;
}

Looks quite long (but nice), but we can make it shorter and messier, too:

sort(v.begin(), v.end(), [] (point a, point b) {
  bool ba = make_pair(a.x,a.y) < pair<int,int>(), bb = make_pair(b.x,b.y) < pair<int,int>();
  if (ba != bb) { return ba < bb; }
  return a.x * b.y > a.y * b.x;
}

You can also use atan2 without worrying about precision issues with a simple hack. If the angles differ by more than, say, $$$0.1$$$, you compare the angles directly, otherwise you directly use the cross product.

BTW why try to make the code as short as possible? If you're not code golfing, you should rather aim for readability. It pays off when debugging (at least it works for me).

→ Reply
  • »
    »
    5 years ago, # ^ |
    ← Rev. 2   Vote: I like it 0 Vote: I do not like it

    Thanks for sharing. You're right about readability. I guess I just wanted to look at some elegant ways to do it from different people.

    Also since I use std::pair directly as my point structure, this is less messy for me.

    sort(v.begin(), v.end(), [] (point a, point b) {
  bool x = a < point(), y = b < point();
  return x == y ? a.x * b.y > a.y * b.x : x < y;
});
    → Reply
»
5 years ago, # |
← Rev. 5   Vote: I like it 0 Vote: I do not like it

(EDIT: I think I was mistaken in the first part, sorry.)

What I do is very similar to what mnbvmar described. I use function

bool IsInUpper(Point p) {
  return p.y > 0 || (p.y == 0 && p.x > 0);
}

and then sort points as you did replacing your quad function with my IsInUpper function.

→ Reply
  • »
    »
    5 years ago, # ^ |
    ← Rev. 2   Vote: I like it 0 Vote: I do not like it

    This is exactly what I was looking for, thank you.

    I actually tried the first snippet for the first time today, and getting AC misled me to believe it's correct for all inputs. It seems like I only needed to differentiate between the two sides of x-axis to get it right.

    I'll update the blog with the following snippet only. 

    inline bool up (point p) {
  return p.y > 0 or (p.y == 0 and p.x >= 0);
}

sort(v.begin(), v.end(), [] (point a, point b) {
  return up(a) == up(b) ? a.x * b.y > a.y * b.x : up(a) < up(b);
});
    → Reply
»
5 years ago, # |
  Vote: I like it 0 Vote: I do not like it 
bool half(pt p, pt v){
    if(v.x==0 && v.y==0) return (p.y>0 || (p.y==0 && p.x<0));
    else return (cross_prod(v, p)<0 || (cross_prod(v, p)==0 && dot_prod(v, p)<0));
}
void polar_sort(vector<pt> &vec, pt sort_around, pt point_first){
    sort(vec.begin(), vec.end(), [sort_around, point_first](pt v, pt w){
        return make_tuple(half(v-sort_around, point_first), 0)<
            make_tuple(half(w-sort_around, point_first), cross_prod(v-sort_around, w-sort_around));
    });
}

this code is adapted from the one given here

→ Reply
»
5 years ago, # |
  Vote: I like it -7 Vote: I do not like it

I'd add that a.x * b.y > a.y * b.x works fine as long as all your points are in the same strict half plane whose border contains $$$(0,0)$$$.

→ Reply
»
5 years ago, # |
  Vote: I like it 0 Vote: I do not like it

In the end, in Goodbye E, I ended up fixing one point P, splitting the other points into those on the left side of the line origin-P and those on the right side, and sorting those on each side by their cosine with respect to the ray origin->P. Basically, I removed the circularity that complicates polar sorting and the rest is just 2 pointers.

→ Reply
»
4 years ago, # |
← Rev. 3   Vote: I like it -59 Vote: I do not like it

sorry I mistook something...just forget what I said

→ Reply