Заметающая прямая
Также: Sweep line · plane sweep · scanline · заметание плоскостью · Бентли — Оттман · Bentley–Ottmann · Шеймос — Хоуи · Shamos–Hoey · площадь объединения прямоугольников · задача Кли
Вертикальная прямая идёт по плоскости слева направо и останавливается только в событиях — концах отрезков и точках их встречи. В каждый момент она помнит, что пересекает, в порядке снизу вверх, и сравнивает только соседей. Так за O((n + k) log n) находятся все k пересечений n отрезков, а с деревом отрезков — площадь объединения прямоугольников. На береговой линии мира — места, где упрощение заставило берег пересечь самого себя.
Береговые линии мира в масштабе 1:50 миллионов — это 1422 линии и 60 409 точек. Чтобы карта рисовалась быстрее, её упрощают: алгоритм Рамера — Дугласа — Пекера выбрасывает точки, которые отстоят от прямой меньше чем на допуск. При допуске 0,2° — это около 22 километров — от 58 987 отрезков остаётся 6306, а 942 острова исчезают совсем. Но берег от этого портится: в 72 местах отрезки, которые раньше были далеко друг от друга, начинают пересекаться или касаться. Пар отрезков почти двадцать миллионов. Прямая, которая проходит карту с запада на восток, находит все 72 места, проверив 10 554 пары, — и при этом в каждый момент помнит не больше 69 отрезков.
Вертикальная линия идёт по карте, подсвечены отрезки, которые она сейчас пересекает, красные кружки — найденные места. Под картой — каждое место крупно: тонкая линия — исходный берег, толстая — упрощённый, красные отрезки — те, что пересеклись. Исходные данные чисты: без упрощения прямая проверяет 109 785 пар соседей и не находит ни одного пересечения. Уже при 0,02° — двух километрах — появляются три. При 1° пересечений меньше, чем при 0,2°, потому что исчезает 1345 островов из 1422: пересекаться почти нечему.
Идея
Проверить все пары n отрезков — n(n − 1)/2 проверок. Но два отрезка, которые лежат в разных концах карты, пересечься не могут, и хочется сравнивать только тех, кто рядом. «Рядом» можно понять так: проведём вертикальную прямую; отрезки, которые она пересекает, упорядочены по высоте точки пересечения. Два отрезка, между которыми на этой прямой есть третий, в этом месте не встречаются.
Будем двигать вертикальную прямую слева направо и хранить состояние — отрезки, которые она пересекает, снизу вверх. Порядок меняется только в событиях: в концах отрезков (отрезок входит в состояние или уходит из него) и в точках, где два отрезка встречаются (они меняются местами). Между событиями сравнивать нечего. В событии проверяются только пары, которые стали соседями; если соседи встречаются правее, точка встречи добавляется в очередь событий.
Состояние — упорядоченный список, но ключ у отрезка — высота на текущей вертикали — меняется по мере движения. Это не мешает: пока прямая между событиями, порядок не меняется, и двоичный поиск по нему работает, если в каждом сравнении вычислять высоту заново. Очередь событий — куча точек в порядке «слева направо, при равных x — снизу вверх».
Пусть q — общая точка двух отрезков. Тогда прямая остановится в q, и в этом событии найдутся все отрезки, проходящие через q.
Если q — конец какого-нибудь отрезка, это событие с самого начала. Пусть q лежит внутри всех отрезков, которые через неё проходят, и пусть p — последнее событие перед q. Между p и q событий нет, поэтому состояние на этом промежутке одно и то же. Возьмём вертикаль чуть левее q. Отрезки через q на ней идут подряд: если бы между двумя из них лежал отрезок r, не проходящий через q, то при приближении вертикали к q его высота оказалась бы зажата между высотами двух отрезков, которые сходятся к q, — и r прошёл бы через q. Значит, какие-то два отрезка через q — соседи в состоянии после p.
Соседями они стали в каком-то событии не позже p: при вставке одного из них, при удалении отрезка между ними или когда в точке встречи поменялся порядок. В каждом из этих случаев алгоритм проверяет пары, ставшие соседями, — значит, проверил и эту и положил q в очередь. Наконец, в событии q все отрезки, которые через неё проходят, — это ровно те отрезки состояния, у которых высота на вертикали через q равна высоте q; они идут подряд, и двоичный поиск находит их все.
Из теоремы видно, почему хватает двух проверок на событие. Отрезки, которые проходят через точку события, уже встретились в ней; второй раз два отрезка встречаются, только если лежат на одной прямой, а такое наложение находят события концов — конец одного лежит на другом. Значит, новые соседи — это только две пары на краях: нижний из прошедших через точку со своим соседом снизу и верхний — со своим соседом сверху. Если через точку ничего не проходит дальше, проверяются два отрезка, между которыми стало пусто.
Это пролив между двумя островами при допуске 0,2°. В самом узком месте он шириной около 2,4 километра, а допуск — около 22: от одного острова остался четырёхугольник, от конца другого — три отрезка, и они пересеклись дважды. На семи отрезках 10 событий: 8 концов и две точки встречи. Прямая проверяет 9 пар соседей из 21 возможной; на целой карте разница — в тысячи раз.
Что происходит в событии
Пусть прямая дошла до точки p. Двоичный поиск по состоянию находит отрезки, которые проходят через p: положение отрезка относительно точки — знак векторного произведения, без вычисления высоты. Среди них — те, что в p кончаются, и те, что через p идут дальше. Если вместе с отрезками, которые в p начинаются, их больше одного — найдена точка встречи. Дальше кончающиеся уходят, а проходящие и начинающиеся встают на освободившееся место в порядке наклона: правее p выше тот, у кого наклон круче. Отрезки, которые пересеклись в p, при этом сами собой меняются местами — отдельной операции «поменять» не нужно.
Два вырожденных случая решает порядок точек. Вертикальный отрезок не пересекает вертикаль в одной точке, но если сравнивать точки сначала по x, а потом по y, его нижний конец — левый, а сам он считается «круче всех»: в своей нижней точке он встаёт выше остальных. Это то же самое, что повернуть карту на бесконечно малый угол. А несколько отрезков, сходящихся в одной точке, обрабатываются одним событием — все сразу.
Точки встречи — дроби: у отрезков с целыми концами пересечение имеет вид (X/W, Y/W). Сравнивать их в числах с плавающей точкой опасно: две почти совпадающие точки могут встать не в том порядке, и событие обработается дважды или пропадёт. В листинге на Python это Fraction, в Go — целые числитель и знаменатель и сравнение произведений через 128-битное умножение: в тысячных долях градуса произведения могут доходить до 1027, далеко за пределы 64 бит.
Код
Пары соседних отрезков одной линии сходятся в общей вершине — это не ошибка, если они не легли друг на друга. При 0,2° из 72 пар 40 — линия сама с собой, 32 — две разные линии, как два острова в проливе выше. Проверочный подсчёт перебором пар по сетке даёт те же числа при каждом допуске.
Кстати, упрощение начинается раньше, чем кажется: если просто округлить координаты этой карты до сотых долей градуса — около километра, — появятся 24 пары касающихся или пересёкшихся отрезков, почти все на крошечных островах. Поэтому данные здесь хранятся в тысячных.
Площадь объединения прямоугольников
Прямая годится не только для пересечений. Пусть на карте лежат прямоугольники со сторонами вдоль осей — кадры спутниковой съёмки, — и нужна площадь, покрытая хотя бы одним. Событиями будут левые и правые края. Между двумя соседними событиями прямая накрыта одним и тем же множеством отрезков, и площадь полосы — накрытая длина, умноженная на ширину полосы. Остаётся быстро знать накрытую длину, когда на прямую добавляется или с неё убирается отрезок [y1, y2].
Для этого — дерево отрезков над сжатыми координатами y: лист — промежуток между соседними различными y. Отрезок раскладывается на O(log n) узлов дерева, как любой запрос. В узле два числа: cnt — сколько добавленных отрезков легли на этот узел целиком, и cov — накрытая длина внутри узла. Если cnt > 0, cov — вся длина узла, иначе — сумма cov детей (у листа — 0). Ответ «сколько накрыто» лежит в корне.
Если каждый отрезок удаляется только после того, как был добавлен, то cnt никогда не становится отрицательным, и cov в корне равна длине объединения отрезков, которые сейчас на прямой.
Разложение промежутка [a, b) на узлы дерева зависит только от a и b. Поэтому удаление вычитает единицу ровно из тех узлов, куда добавление её прибавило, и cnt[v] всегда равен числу текущих отрезков, в разложение которых входит v, — это неотрицательное число.
Докажем индукцией снизу вверх, что cov[v] — длина той части узла v, которую накрывают узлы разложений, лежащие в поддереве v (включая сам v). Для листа это 0 или его длина в зависимости от cnt. Для внутреннего узла: если на него лёг отрезок целиком, накрыт весь узел; иначе накрытое внутри — объединение накрытого в двух детях, которые не пересекаются, и длины складываются. Пересчёт идёт по всем узлам, у которых изменился cnt, и по всем их предкам; у остальных узлов поддерево не изменилось. У корня предков нет, а каждый отрезок на прямой раскладывается на узлы внутри дерева — значит, cov в корне и есть длина объединения.
Выходит O(n log n) на n прямоугольников. Отложенных операций, как в обычном дереве с прибавлением на отрезке, не нужно: счётчик никогда не спускается к детям. Этим решением Джон Бентли в 1977 году ответил на двумерный вариант вопроса, который задал Виктор Кли. Задача «Покрытие снимками» ниже — ровно оно.
Сколько это стоит
Событий не больше 2n + k: концы отрезков и k точек встречи. Каждое — извлечение из кучи, двоичный поиск по состоянию из не более n отрезков и вставка, то есть O(log n) сравнений, плюс работа, пропорциональная числу отрезков через точку. Всего O((n + k) log n). Вставка в массив в листингах сдвигает хвост и формально стоит O(n), но состояние на карте маленькое — меридиан пересекает самое большее 85 отрезков исходного берега из 58 987. Для гарантированной оценки состояние держат в сбалансированном дереве поиска.
Перебор пар растёт как n², прямая — почти линейно. В Python перебор проигрывает с самого начала: на 2259 отрезках — в семьдесят с лишним раз. В Go перебор простых умножений быстрее прямой на сотнях отрезков, примерно до полутора тысяч, а на всей карте прямая выигрывает в три десятка раз. Если нужно только узнать, есть ли хоть одно пересечение, точки встречи не нужны вовсе: алгоритм Шеймоса — Хоуи останавливается на первой найденной паре, событиями у него только концы, и он работает за O(n log n) — это задача «Граница без петель».
Множитель k убрать нельзя: пересечения нужно хотя бы перечислить, а их бывает до n²/2. Но O((n + k) log n) — не предел: Бернард Шазель и Герберт Эдельсбруннер в 1992 году получили O(n log n + k). Когда пересечений очень много, проще бывает сетка: в среднем она быстра, но гарантий не даёт.
Где легко ошибиться
- Числа с плавающей точкой в событиях. Две близкие точки встречи можно поставить в очередь в неверном порядке или посчитать разными, когда это одна точка. Тогда состояние перестаёт быть упорядоченным, и двоичный поиск начинает промахиваться — ошибка всплывает далеко от места, где случилась. Для проверки корректности геометрии нужна точная арифметика.
- Проверка только при вставке. Когда отрезок уходит, его соседи снизу и сверху становятся соседями друг друга, и эту пару тоже нужно проверить. Без этого теряются пересечения, которые «загораживал» короткий отрезок.
- Повторные события. Одну и ту же точку встречи находят несколько раз — каждый раз, когда её пара снова становится соседями. В очередь она должна попасть один раз: помогает множество уже поставленных точек или проверка «правее текущей».
- Вертикальные отрезки и несколько отрезков в одной точке. Если сравнивать точки по x, а потом по y, и обрабатывать все отрезки через точку одним событием, отдельные случаи не нужны. Пропустить это — значит получить алгоритм, который работает на случайных данных и ломается на настоящих: в береговых линиях есть и вертикальные отрезки, и вершины, лежащие на чужих отрезках.
- Соседние отрезки одной линии. Они всегда встречаются в общей вершине; считать это пересечением нельзя, а вот наложение — когда линия идёт назад по себе — нужно.
Где работает
Картография. Упрощение границ — обычная операция перед публикацией карты, и программы после него ищут появившиеся самопересечения: в mapshaper починка пересечений после упрощения включена по умолчанию, а найденные пересечения можно подсветить на карте. Проверка «граница не пересекает сама себя» входит в определение правильного многоугольника в геоинформационных системах.
Вычислительная геометрия. В библиотеке CGAL пересечения кривых на плоскости считаются заметающей прямой за O((n + k) log n), с наложениями, вертикальными отрезками и касаниями. На том же приёме — алгоритм Стивена Форчуна для диаграммы Вороного: прямая идёт через точки, а состояние — «береговая линия» из дуг парабол.
Графика. Заливка многоугольника на растре идёт строка за строкой: для каждой строки известен список рёбер, которые она пересекает, отсортированный по x, и закрашиваются промежутки между парами. Это та же прямая, только горизонтальная и с шагом в пиксель.
Задачи на отрезки и интервалы. Многие задачи о прямоугольниках и интервалах — сколько точек в прямоугольнике, какая наибольшая глубина перекрытия, какой отрезок первым накроет точку — решаются обходом событий в отсортированном порядке со структурой данных, которая помнит текущее состояние. Пример — подсчёт населения прямоугольников в статье о сжатии координат.
История
Майкл Шеймос и Дэн Хоуи в 1976 году показали, что узнать, пересекаются ли хоть какие-то два из n отрезков, можно за O(n log n), и сделали это прямой, идущей по плоскости. В 1977 году Виктор Кли спросил, можно ли длину объединения n отрезков на прямой найти быстрее O(n log n), и Джон Бентли в неопубликованных заметках решил двумерный вариант — площадь объединения прямоугольников — заметающей прямой с деревом отрезков. В 1979 году Бентли и Томас Оттман опубликовали алгоритм, который находит все k пересечений за O((n + k) log n). Оптимальную оценку O(n log n + k) получили Шазель и Эдельсбруннер в 1992 году, а в 1995 году Иван Балабан — с памятью O(n). Стивен Форчун в 1986 году построил заметающей прямой диаграмму Вороного.
Задачи
Проверка, что упрощённая граница не пересекает сама себя, — Шеймос и Хоуи; площадь, покрытая кадрами, — прямая с деревом отрезков; подсчёт мостов и тоннелей на улицах — все пересечения внутренностей отрезков.
Граница без петель Покрытие снимками Мосты и тоннелиДанные: Береговые линии: Natural Earth, общественное достояние