C++で実装するAngular Sweep（角度スイープ）アルゴリズムの解説

本記事では、与えられた半径の円に含めることができる点の最大数を求めるアルゴリズムについて解説します。具体的には、半径 r の円と2次元平面上の点の集合が与えられたとき、その円の内部（円周上は除く）に収まる点の数の最大値を見つけるという問題です。

この問題を効率的に解くための最も有効な手法が「Angular Sweep（角度スイープ）アルゴリズム」です。

アルゴリズムの基本的な考え方

各点 P を基準としたとき、半径 r の円が別の点 Q を覆える角度の範囲は「区間」として表現できます。点 Q が P から距離 2r 以内にあれば、P の周りで円を回転させる際に Q が円内に入る角度範囲 [α, β] が存在します。この区間の開始・終了をイベントとしてソートし、走査（スイープ）することで、同時に覆盖できる点の最大数を求められます。計算量は O(n² log n) です。

アルゴリズムの手順

• 与えられた n 個の点から作られる ⁿC₂ 通りのペアすべてについて、点間距離を事前に計算しておきます。

• 任意の点 P を基準に選び、その点を中心に半径 r の円を回転させたときに、円内に含まれる点の最大数を角度スイープによって求めます。

• すべての基準点について得られた結果のうち、最大の値を問題の答えとして返します。

C++による実装例

#include <bits/stdc++.h>
using namespace std;
#define MAX_POINTS 500
typedef complex<double> Point;
Point arr[MAX_POINTS];
double dis[MAX_POINTS][MAX_POINTS];
int getPointsInside(int i, double r, int n) {
    vector <pair<double, bool> > angles;
    for (int j = 0; j < n; j++) {
        if (i != j && dis[i][j] <= 2*r) {
            double B = acos(dis[i][j]/(2*r));
            double A = arg(arr[j]-arr[i]);
            double alpha = A-B;
            double beta = A+B;
            angles.push_back(make_pair(alpha, true));
            angles.push_back(make_pair(beta, false));
        }
    }
    sort(angles.begin(), angles.end());
    int count = 1, res = 1;
    vector <pair<double, bool>>::iterator it;
    for (it=angles.begin(); it!=angles.end(); ++it) {
        if ((*it).second)
            count++;
        else
            count--;
        if (count > res)
            res = count;
    }
    return res;
}
int maxPoints(Point arr[], int n, int r) {
    for (int i = 0; i < n-1; i++)
        for (int j=i+1; j < n; j++)
            dis[i][j] = dis[j][i] = abs(arr[i]-arr[j]);
    int ans = 0;
    for (int i = 0; i < n; i++)
        ans = max(ans, getPointsInside(i, r, n));
    return ans;
}
int main() {
    Point arr[] = {Point(6.47634, 7.69628), Point(5.16828, 4.79915), Point(6.69533, 6.20378)};
    int r = 1;
    int n = sizeof(arr)/sizeof(arr[0]);
    cout << "The maximum number of points are: " << maxPoints(arr, n, r);
    return 0;
}

実行結果

The maximum number of points are: 2

コードのポイント

• acos(dis[i][j]/(2*r)) により、点 Q が円に出入りする境界角度 B を算出しています。

• arg(arr[j]-arr[i]) は複素数を利用して、ベクトル PQ の偏角 A を求めています。

• 区間の開始（true）でカウントを増やし、終了（false）で減らすことで、スイープ中の最大覆盖数を記録します。