В данной статье я рассмотрю два алгоритма, первый — непосредственно кластеризация, второй — построение контура кластера в виде выпуклого многоугольника, прикладная задача для улучшеного восприятия полученного результата.
Кластеризация
Вообще закомство с данной областью довольно поверхностное, поэтому вполне ожидаемо что такой алгоритм давно существует и как-то назван, если кто знает, прошу просветить.
Входные данные
- Координаты точек — просто массив двухмерных координат
- R — Максимальное расстояние между точками в кластере — основной показатель для построения кластера
- k — Коэффициент сетки — степень для числа 2, по которой вычисляется сторона ячейки сети: a = 2^k
Для эффективной работы алгоритма стоит соблюдать следующее условие:
R < 2^k, либо модифицировать алгоритм.
Алгоритм
- Сопоставление координат точек с ячейками регулярной сети.
Xr = x >> k
Yr = y >> k, где k — коэффициент сетки, >> — побитовое смещение. Благодаря отсутствию операции деления на этом шаге получается хороший выигрыш в скорости. Результатом работы на это шаге будет словарь вида:
[координаты ячейки сети]: [ссылки на точки относящиеся к данной ячейке] - Построение кластеров за один проход по точкам:
- Берется очередная точка, если она не привязана к кластеру, создается новый кластер и точка добавляется в него, в самой точке проставляется ссылка на кластер.
- Берется список ячеек в него входят: ячейка с точкой и соседние с ней.
- Перебираются все точки в ячейках из полученного списка.
- Находится расстояние между текущей точкой и проверяемой из списка, если расстояние <= R то возможны два случая
- проверяемая точка пока не входит в кластер => точка вносится в тот же кластер в котором находится текущая точка
- проверяемая точка уже принадлежит какому-то кластеру, отличному от текущего => кластер с большим количеством точек поглощает кластер с меньшим.
Результатом работы будет набор кластеров в каждом из которых расстояние между точками не превышает R, примерно так:
Или более наглядно так:
Точки одного кластера обозначены одним цветом.
Контур
Работает, но результат хочется смотреть в более приглядном виде, самый простой путь, выделить клатер контуром, желательно в виде выпуклого многоугольника.
Главным инструментом служит векторное перемножение по трем точкам:
public static int CrossLinePoint(ClusterPoint l1, ClusterPoint l2, ClusterPoint p)
{
return (p.X - l2.X) * (l2.Y - l1.Y) - (p.Y - l2.Y) * (l2.X - l1.X);
}
Если результат функции меньше 0, то изгиб в средней точке идет по часовой стрелке, если больше 0, то против, если равен 0, изгиба нет. Само значение не важно, важен только знак.
Алгоритм состоит всего из двух шагов.
Выделение описывающего контура
- Находится центр масс кластера, без учета весовых коэффициентов. Т.е. просто складываются все координаты, отдельно по x и y, и делятся на количество точек в кластере, полученные координаты [Xc;Yc] и будут центром масс кластера.
- Ищутся три наиболее удаленных от центра масс точки, которые образуют начальный контур.
- Перебираются оставшиеся точки, и для каждой проверяется вхождение в начальный контур, если точка входит в него, ничего не происходит, если снаружи, становится частью контура.
- При добавлении точки в контур, массив точек контура сортируется, в качестве величины для сортировки берется угол из полярных координат точки, где центром координат выступает центр масс кластера.
- Результатом будет ломанный контур описывающий точки входящие в кластер, причем точки в контуре уже упорядочены по часовой стрелке относительно цетра масс.
Функция для получения угла из полярных координат:
private float GetPolarAngle(ClusterPoint p)
{
float a = 0, x = p.X - center.X, y = p.Y - center.Y;
if (x > 0 && y >= 0)
{
a = (float)Math.Atan(y / x);
}
else if (x > 0 && y < 0)
{
a = (float)(Math.Atan(y / x) + 2 * Math.PI);
}
else if (x < 0)
{
a = (float)(Math.Atan(y / x) + Math.PI);
}
else if (x == 0 && y > 0)
{
a = (float)(Math.PI * .5);
}
else if (x == 0 && y < 0)
{
a = (float)(Math.PI * 1.5);
}
else if (x == 0 && y == 0)
{
a = 0;
}
return (float)(a * 180 / Math.PI);
}
Здесь center — центр масс кластера.
Сортировка по углу нужна для отрисовки точек контура в правильном порядке, и для определения вхождение/невхождения точки в контур, при построении контура. Если оставить точки в случайном порядке, алгоритм работать не будет.
Результат работы:
Зеленая точка, это центр масс.
Все еще не особо приглядно, поэтому переходим к последнему шагу.
Построение выпуклого контура
Самый простой шаг, в котором точки уже упорядочены по часовой стрелке, делается проход по точкам и выкидываются из контура те, в которых изгиб идет против часовой стрелки. Для определения изгиба используется приведенная выше функция CrossLinePoint.
Результат работы:
Тестирование
Графики скорости работы алгоритма.
Поле имеет постоянный размер, плостность точек возрастает:
Поле увеличивается в размере, плотность точек примерно одинакова:
Время обработки сильно зависит от начальных данных. Если задать большое расстояние между точками в кластере на поле с высокой плотностью точек, время вырастает в разы, что логично, т.к. количество сравниваемых точек возрастает. К примеру, на приведенных графиках скорость работы отличается в три раза. На первом видно что для 20.000 точек затрачено 140 секунд, а в случае когда плотность не увеличивается всего 40 секунд для тех же 20.000 точек.
Если вы занимались подобной задачей и у вас есть собственные замеры, буду рад увидеть результаты в комментариях, для сравнения эффективности алгоритмов.
Для тестирования не использовалось построение контура. В виду неоптимизированности кода поиска контура, время на его построение в разы превышает время нахождение самих кластеров.
Что дальше
Алгоритм работает, но только для набора точек имеющих координаты, следующим шагом будет попытка модификации для выделения областей на непрерывных изображениях (рисунки, фотографии и т.п.).
Реализация
Если кого-то заинтересовал алгоритм, есть программа и исходники. Win/C#
Win x86
Sources
Кроме показанных в статье отображений кластеров в виде точек и в виде точек с ограничивающим контуром из линий, в программе возможны еще отображения в виде замкнутых кривых, и закрашенной области из замкнутых кривых.
Используются стандартные средства GDI+.