Angle Sort
(geometry/angle_sort.hpp)
- View this file on GitHub
- Last update: 2026-10-05 22:23:07+09:00
- Include:
#include "geometry/angle_sort.hpp"
Overview
This header sorts points counterclockwise by their angle around a chosen
origin, without calling atan2.
sort_by_angle(points, origin, start);
auto sorted = angle_sorted(points, origin, start);
sort_by_angle modifies its argument. angle_sorted returns a sorted copy.
Both take $O(N\log N)$ time.
Complexity Notation
-
Nis the number of points being sorted.
sort_by_angle and angle_sorted take $O(N\log N)$ time. The comparator uses
$O(1)$ additional memory per comparison.
Ordering Convention
The default, AngleSortStart::NegativeXAxis, sorts by the usual atan2(y, x)
value from $-\pi$ to $\pi$. Imagine starting just below the negative x-axis
and rotating counterclockwise:
| Order | Direction | Angle |
|---|---|---|
| 1 | down-left | close to $-\pi$ |
| 2 | down | $-\pi/2$ |
| 3 | right | $0$ |
| 4 | up | $\pi/2$ |
| 5 | left | $\pi$ |
In short:
down-left -> down -> right -> up -> up-left -> left
The origin is treated as having angle zero, so it appears with points pointing right. This convention matches Library Checker’s Sort Points by Argument problem.
AngleSortStart::PositiveXAxis uses the more common circular ordering that
starts on the positive x-axis:
right -> up -> left -> down -> back to right
Points on the same ray are ordered by increasing distance from origin.
Equivalent duplicate points may appear in either relative order.
Comparator
AngleLess<T> is the comparator used by the helpers:
m1une::geometry::AngleLess<long long> less(origin);
std::sort(points.begin(), points.end(), less);
For integral coordinates, comparisons use signed 128-bit arithmetic. As with the other exact geometry predicates, coordinate differences and products must fit that type.
Example
#include "geometry/angle_sort.hpp"
#include <iostream>
#include <vector>
int main() {
using Point = m1une::geometry::Point<long long>;
std::vector<Point> points;
points.emplace_back(1, 0);
points.emplace_back(0, 1);
points.emplace_back(-1, 0);
points.emplace_back(0, -1);
m1une::geometry::sort_by_angle(
points,
Point(0, 0),
m1une::geometry::AngleSortStart::PositiveXAxis
);
for (const Point& point : points) {
std::cout << point.x << ' ' << point.y << "\n";
}
}
Depends on
Required by
Verified with
verify/geometry/angle_sort.test.cpp
verify/geometry/centroid.test.cpp
verify/geometry/geometry_algorithms.test.cpp
verify/geometry/rational.test.cpp
Code
#ifndef M1UNE_GEOMETRY_ANGLE_SORT_HPP
#define M1UNE_GEOMETRY_ANGLE_SORT_HPP 1
#include <algorithm>
#include <vector>
#include "point.hpp"
namespace m1une {
namespace geometry {
enum class AngleSortStart {
NegativeXAxis,
PositiveXAxis,
};
template <Coordinate T>
struct AngleLess {
Point<T> origin;
AngleSortStart start;
constexpr explicit AngleLess(
Point<T> origin_value = Point<T>(),
AngleSortStart start_value = AngleSortStart::NegativeXAxis
) : origin(origin_value), start(start_value) {}
constexpr bool operator()(
const Point<T>& first,
const Point<T>& second
) const {
using W = wide_type<T>;
W first_x = W(first.x) - W(origin.x);
W first_y = W(first.y) - W(origin.y);
W second_x = W(second.x) - W(origin.x);
W second_y = W(second.y) - W(origin.y);
W first_distance = first_x * first_x + first_y * first_y;
W second_distance = second_x * second_x + second_y * second_y;
// atan2(0, 0) is treated as angle zero.
if (first_distance == 0) first_x = 1;
if (second_distance == 0) second_x = 1;
auto half = [this](W x, W y) {
if (start == AngleSortStart::PositiveXAxis) {
return y < 0 || (y == 0 && x < 0);
}
return y > 0 || (y == 0 && x < 0);
};
bool first_half = half(first_x, first_y);
bool second_half = half(second_x, second_y);
if (first_half != second_half) return first_half < second_half;
W product = first_x * second_y - first_y * second_x;
if (product != 0) return product > 0;
return first_distance < second_distance;
}
};
// Sorts points counterclockwise by angle around `origin`.
template <Coordinate T>
void sort_by_angle(
std::vector<Point<T>>& points,
Point<T> origin = Point<T>(),
AngleSortStart start = AngleSortStart::NegativeXAxis
) {
std::sort(points.begin(), points.end(), AngleLess<T>(origin, start));
}
// Returns a counterclockwise angle-sorted copy.
template <Coordinate T>
std::vector<Point<T>> angle_sorted(
std::vector<Point<T>> points,
Point<T> origin = Point<T>(),
AngleSortStart start = AngleSortStart::NegativeXAxis
) {
sort_by_angle(points, origin, start);
return points;
}
} // namespace geometry
} // namespace m1une
#endif // M1UNE_GEOMETRY_ANGLE_SORT_HPP#line 1 "geometry/angle_sort.hpp"
#include <algorithm>
#include <vector>
#line 1 "geometry/point.hpp"
#include <cmath>
#include <concepts>
#include <cassert>
#include <type_traits>
#line 1 "geometry/detail/floating_predicate.hpp"
namespace m1une {
namespace geometry {
namespace predicate_detail {
template <typename T>
constexpr T absolute(T value) {
return value < T(0) ? -value : value;
}
template <typename T>
constexpr T max_value(T first, T second) {
return first < second ? second : first;
}
template <typename T>
constexpr T vector_scale(T x, T y) {
return max_value(absolute(x), absolute(y));
}
template <bool Exact, typename T>
constexpr int scaled_sign(T value, T scale, long double eps) {
if constexpr (Exact) {
return (value > T(0)) - (value < T(0));
} else {
const T tolerance = T(eps) * scale;
return (value > tolerance) - (value < -tolerance);
}
}
template <bool Exact, typename T>
constexpr T determinant_scale(T ax, T ay, T bx, T by) {
if constexpr (Exact) {
return T(0);
} else {
return vector_scale(ax, ay) * vector_scale(bx, by);
}
}
template <bool Exact, typename T>
constexpr int determinant_sign(
T ax,
T ay,
T bx,
T by,
long double eps
) {
const T determinant = ax * by - ay * bx;
return scaled_sign<Exact>(
determinant,
determinant_scale<Exact>(ax, ay, bx, by),
eps
);
}
template <bool Exact, typename T>
constexpr int orientation_sign(
T direction_x,
T direction_y,
T offset_x,
T offset_y,
long double eps
) {
const T determinant =
direction_x * offset_y - direction_y * offset_x;
T scale = T(0);
if constexpr (!Exact) {
const T direction_scale =
vector_scale(direction_x, direction_y);
scale = direction_scale * max_value(
direction_scale,
vector_scale(offset_x, offset_y)
);
}
return scaled_sign<Exact>(determinant, scale, eps);
}
template <bool Exact, typename T>
constexpr int dot_sign(
T ax,
T ay,
T bx,
T by,
long double eps
) {
const T value = ax * bx + ay * by;
T scale = T(0);
if constexpr (!Exact) {
scale = vector_scale(ax, ay) * vector_scale(bx, by);
}
return scaled_sign<Exact>(value, scale, eps);
}
} // namespace predicate_detail
} // namespace geometry
} // namespace m1une
#line 10 "geometry/point.hpp"
namespace m1une {
namespace geometry {
template <typename T>
concept Coordinate = !std::same_as<std::remove_cv_t<T>, bool> &&
(std::is_arithmetic_v<T> ||
(std::copyable<T> && std::totally_ordered<T> && requires(T a, T b) {
T(0);
T(1);
static_cast<long double>(a);
{ +a } -> std::same_as<T>;
{ -a } -> std::same_as<T>;
{ a + b } -> std::same_as<T>;
{ a - b } -> std::same_as<T>;
{ a * b } -> std::same_as<T>;
{ a / b } -> std::same_as<T>;
{ a += b } -> std::same_as<T&>;
{ a -= b } -> std::same_as<T&>;
}));
// Custom coordinate types keep their own exact arithmetic.
template <typename T>
concept ExactCoordinate = Coordinate<T> && !std::floating_point<T>;
template <Coordinate T>
using wide_type = std::conditional_t<std::integral<T>, __int128_t,
std::conditional_t<std::floating_point<T>, long double, T>>;
template <Coordinate T>
struct Point {
T x;
T y;
constexpr Point() : x(0), y(0) {}
constexpr Point(T x_value, T y_value) : x(x_value), y(y_value) {}
template <Coordinate U>
explicit constexpr Point(const Point<U>& other)
: x(static_cast<T>(other.x)), y(static_cast<T>(other.y)) {}
constexpr Point& operator+=(const Point& other) {
x += other.x;
y += other.y;
return *this;
}
constexpr Point& operator-=(const Point& other) {
x -= other.x;
y -= other.y;
return *this;
}
constexpr Point operator+() const {
return *this;
}
constexpr Point operator-() const {
return Point(-x, -y);
}
friend constexpr Point operator+(Point left, const Point& right) {
return left += right;
}
friend constexpr Point operator-(Point left, const Point& right) {
return left -= right;
}
friend constexpr bool operator==(const Point&, const Point&) = default;
friend constexpr bool operator<(const Point& left, const Point& right) {
if (left.x != right.x) return left.x < right.x;
return left.y < right.y;
}
};
template <Coordinate T>
constexpr Point<long double> centroid(const Point<T>& point) {
return Point<long double>(point);
}
template <Coordinate T, typename Scalar>
requires (std::is_arithmetic_v<Scalar> || Coordinate<Scalar>)
constexpr auto operator*(const Point<T>& point, Scalar scalar) {
using Result = std::common_type_t<T, Scalar>;
return Point<Result>(
Result(point.x) * Result(scalar),
Result(point.y) * Result(scalar)
);
}
template <typename Scalar, Coordinate T>
requires (std::is_arithmetic_v<Scalar> || Coordinate<Scalar>)
constexpr auto operator*(Scalar scalar, const Point<T>& point) {
return point * scalar;
}
template <Coordinate T, typename Scalar>
requires (std::is_arithmetic_v<Scalar> || Coordinate<Scalar>)
constexpr auto operator/(const Point<T>& point, Scalar scalar) {
using Result = std::common_type_t<T, Scalar>;
return Point<Result>(
Result(point.x) / Result(scalar),
Result(point.y) / Result(scalar)
);
}
template <Coordinate T>
constexpr wide_type<T> dot(const Point<T>& a, const Point<T>& b) {
using W = wide_type<T>;
return W(a.x) * W(b.x) + W(a.y) * W(b.y);
}
template <Coordinate T>
constexpr wide_type<T> cross(const Point<T>& a, const Point<T>& b) {
using W = wide_type<T>;
return W(a.x) * W(b.y) - W(a.y) * W(b.x);
}
template <Coordinate T>
constexpr wide_type<T> cross(
const Point<T>& origin,
const Point<T>& a,
const Point<T>& b
) {
using W = wide_type<T>;
W ax = W(a.x) - W(origin.x);
W ay = W(a.y) - W(origin.y);
W bx = W(b.x) - W(origin.x);
W by = W(b.y) - W(origin.y);
return ax * by - ay * bx;
}
template <Coordinate T>
constexpr wide_type<T> norm2(const Point<T>& point) {
return dot(point, point);
}
template <Coordinate T>
constexpr wide_type<T> distance2(const Point<T>& a, const Point<T>& b) {
using W = wide_type<T>;
W dx = W(a.x) - W(b.x);
W dy = W(a.y) - W(b.y);
return dx * dx + dy * dy;
}
template <Coordinate T>
long double norm(const Point<T>& point) {
return std::hypot(
static_cast<long double>(point.x),
static_cast<long double>(point.y)
);
}
template <Coordinate T>
long double distance(const Point<T>& a, const Point<T>& b) {
return std::hypot(
static_cast<long double>(a.x) - static_cast<long double>(b.x),
static_cast<long double>(a.y) - static_cast<long double>(b.y)
);
}
template <Coordinate T, typename M, typename N>
requires (std::is_arithmetic_v<M> || Coordinate<M>) &&
(std::is_arithmetic_v<N> || Coordinate<N>)
constexpr Point<long double> internal_division_point(
const Point<T>& a,
const Point<T>& b,
M m,
N n
) {
long double first_ratio = static_cast<long double>(m);
long double second_ratio = static_cast<long double>(n);
long double denominator = first_ratio + second_ratio;
assert(denominator != 0);
Point<long double> first(a);
Point<long double> direction = Point<long double>(b) - first;
return first + direction * (first_ratio / denominator);
}
template <Coordinate T, typename M, typename N>
requires (std::is_arithmetic_v<M> || Coordinate<M>) &&
(std::is_arithmetic_v<N> || Coordinate<N>)
constexpr Point<long double> external_division_point(
const Point<T>& a,
const Point<T>& b,
M m,
N n
) {
long double first_ratio = static_cast<long double>(m);
long double second_ratio = static_cast<long double>(n);
long double denominator = first_ratio - second_ratio;
assert(denominator != 0);
Point<long double> first(a);
Point<long double> direction = Point<long double>(b) - first;
return first + direction * (first_ratio / denominator);
}
template <Coordinate T>
constexpr int sign(wide_type<T> value, long double eps = 1e-12L) {
return predicate_detail::scaled_sign<ExactCoordinate<T>>(
value,
wide_type<T>(1),
eps
);
}
template <Coordinate T>
constexpr int orientation(
const Point<T>& a,
const Point<T>& b,
const Point<T>& c,
long double eps = 1e-12L
) {
using W = wide_type<T>;
const W first_x = W(b.x) - W(a.x);
const W first_y = W(b.y) - W(a.y);
const W second_x = W(c.x) - W(a.x);
const W second_y = W(c.y) - W(a.y);
return predicate_detail::orientation_sign<ExactCoordinate<T>>(
first_x,
first_y,
second_x,
second_y,
eps
);
}
template <Coordinate T>
constexpr bool collinear(
const Point<T>& a,
const Point<T>& b,
const Point<T>& c,
long double eps = 1e-12L
) {
return orientation(a, b, c, eps) == 0;
}
template <Coordinate T>
Point<long double> rotate(const Point<T>& point, long double angle) {
long double cosine = std::cos(angle);
long double sine = std::sin(angle);
return Point<long double>(
static_cast<long double>(point.x) * cosine -
static_cast<long double>(point.y) * sine,
static_cast<long double>(point.x) * sine +
static_cast<long double>(point.y) * cosine
);
}
template <Coordinate T>
Point<long double> normalized(const Point<T>& point) {
long double length = norm(point);
assert(length != 0);
return Point<long double>(
static_cast<long double>(point.x) / length,
static_cast<long double>(point.y) / length
);
}
} // namespace geometry
} // namespace m1une
#line 8 "geometry/angle_sort.hpp"
namespace m1une {
namespace geometry {
enum class AngleSortStart {
NegativeXAxis,
PositiveXAxis,
};
template <Coordinate T>
struct AngleLess {
Point<T> origin;
AngleSortStart start;
constexpr explicit AngleLess(
Point<T> origin_value = Point<T>(),
AngleSortStart start_value = AngleSortStart::NegativeXAxis
) : origin(origin_value), start(start_value) {}
constexpr bool operator()(
const Point<T>& first,
const Point<T>& second
) const {
using W = wide_type<T>;
W first_x = W(first.x) - W(origin.x);
W first_y = W(first.y) - W(origin.y);
W second_x = W(second.x) - W(origin.x);
W second_y = W(second.y) - W(origin.y);
W first_distance = first_x * first_x + first_y * first_y;
W second_distance = second_x * second_x + second_y * second_y;
// atan2(0, 0) is treated as angle zero.
if (first_distance == 0) first_x = 1;
if (second_distance == 0) second_x = 1;
auto half = [this](W x, W y) {
if (start == AngleSortStart::PositiveXAxis) {
return y < 0 || (y == 0 && x < 0);
}
return y > 0 || (y == 0 && x < 0);
};
bool first_half = half(first_x, first_y);
bool second_half = half(second_x, second_y);
if (first_half != second_half) return first_half < second_half;
W product = first_x * second_y - first_y * second_x;
if (product != 0) return product > 0;
return first_distance < second_distance;
}
};
// Sorts points counterclockwise by angle around `origin`.
template <Coordinate T>
void sort_by_angle(
std::vector<Point<T>>& points,
Point<T> origin = Point<T>(),
AngleSortStart start = AngleSortStart::NegativeXAxis
) {
std::sort(points.begin(), points.end(), AngleLess<T>(origin, start));
}
// Returns a counterclockwise angle-sorted copy.
template <Coordinate T>
std::vector<Point<T>> angle_sorted(
std::vector<Point<T>> points,
Point<T> origin = Point<T>(),
AngleSortStart start = AngleSortStart::NegativeXAxis
) {
sort_by_angle(points, origin, start);
return points;
}
} // namespace geometry
} // namespace m1une