lower_bound_2d (найти ближайшее)

[Igors]

Добрый день

Есть массив/контейнер эл-тов, найти в нем ближайшее к заданному ключу. Если напр тип эл-та double (или др численный) то ближайший тот эл-т для которого

Код

C++ (Qt)
double Near( const Key & k1, const Key & k2 ) {  return fabs(k1 - k2); } 
 

минимально. Ну конечно возможны вариации: ближайшее большее, меньшее, все и.т.п. - не суть.

В простейшем случае это решается с помощью std::lower_bound, но это не распространяется на 2 и более измерений. Пример: key определен как QPointF. И здесь имеем функтор

Код

C++ (Qt)
double Near( const QPointF & k1, const QPointF  & k2 ) 
{
  QPointF delta = k1 - k2;  
  return delta.x() * delta.x() + delta.y() * delta.y(); 
} 
 

Но как сортировать массив? Первый ключ x(), второй y()? Тогда lower_bound найдет значение далеко от ближайшего, напр ключ QPointF(5, 3)

(4.96, 2)
(4.97, 3) <- этот ближайший
(4.99, 4)
(4.99, 5)
.... 
(5, 100) <- этот будет найден lower_bound

Так как же бум искать?

Спасибо

[Fregloin]

Вам нужно найти ближайшую точку к заданной? Если так, то лучше через вектора. Т.е. найти наименьший вектор от этой точки к сравниваемой.

[Igors]

Вам нужно найти ближайшую точку к заданной? Если так, то лучше через вектора. Т.е. найти наименьший вектор от этой точки к сравниваемой.

Точнее "вектор с наименьшей длиной". Так и записано (см функтор Near выше)

[kamre]

Но как сортировать массив?

Никак не сортировать, нельзя задать линейный порядок таким образом чтобы искать ближайшую точку как для одномерного случая. Обсуждали уже раньше на форуме: http://www.prog.org.ru/topic_24463_0.html

Так как же бум искать?

Если нужно искать быстро, то нужны специальные структуры данных: http://en.wikipedia.org/wiki/Spatial_database#Spatial_index

[Igors]

Если нужно искать быстро, то нужны специальные структуры данных: http://en.wikipedia.org/wiki/Spatial_database#Spatial_index

Т.е. нужно городить деревья, ноды, без этого не обойтись, а просто использовать массив - никак. Правильно ли я Вас понял?

Никак не сортировать, нельзя задать линейный порядок таким образом ..

А нелинейный? Давайте так:

- есть массив эл-ты которого можно переставлять/организовывать как угодно (сортировать или как - Ваше дело). Реализовать быстрый поиск в таком массиве, не создавая доп структур данных (таких как ноды и.т.п)

[kramer]

Самый алгоритмически простой способ - накрыть пространство/плоскость сеткой из прямоугольных бакетов, в каждом из которых хранить индексы входящих в него точек. Тогда ближайшую к данной точку можной искать перебором, начиная с бакета данной точки, и далее по спирали.
Если точки распределены более-менее равномерно, то это довольно эффективный метод.

А так, да, деревья, вот это все.

Для сходной задачи я, было дело, использовал покоординатное индексирование, но там не совсем поиск ближайшей, нужно было все точки на заданном расстоянии от данной найти.

- есть массив эл-ты которого можно переставлять/организовывать как угодно (сортировать или как - Ваше дело). Реализовать быстрый поиск в таком массиве, не создавая доп структур данных (таких как ноды и.т.п)

В общем случае, по-моему, это невозможно. Хотя есть такая штука, неявные К-мерные деревья, но как на них поиск организовать - представления не имею. Еще можно сортировать вдоль заполняющей кривой (Гильбертовой, например) некоторой глубины, но и тут поиск будет отнюдь не тривиальным. Такая сортировка будет, по сути, неявным линейным октодеревом (кваддеревом на плоскости), так что алгоритм поиска будет сходным, но вот индексы соседей замучаетесь вычислять.

[Igors]

Самый алгоритмически простой способ - накрыть пространство/плоскость сеткой из прямоугольных бакетов, в каждом из которых хранить индексы входящих в него точек. Тогда ближайшую к данной точку можной искать перебором, начиная с бакета данной точки, и далее по спирали.
Если точки распределены более-менее равномерно, то это довольно эффективный метод.

Приятно говорить с понимающим человеком. Но это годится для немного но все же др задачи: найти на расстоянии не больше заданного. Вот тогда попав в "пустой" бакет (и проверив соседей) спокойно уходим - требуемых ближайших нет. А насчет "по спирали" - там только проверка соседей геморрой не маленький, напр в 3D их уже 26 

А так, да, деревья, вот это все.

А мне кажется можно. Пусть есть 100 эл-тов (0..99), напр тех же QPointF. Схема поиска:

- полагаем начальное расстояние = infinity

- берем средний эл-т A[50]. Массив (A) упорядочен так что эл-ты [0..49] имеют x "не больше" чем A[50]. а эл-ты [51..99] "не меньше". Вычисляем расстояние от запросной точки до A[50]. Ну пока мы ничего не можем сказать, просто улучшили расстояние

- повторяем поиск рекурсивно в 2 половинках чередуя ось. Напр в первой половинке [0..49] эл-т A[25] делит по Y. Теперь если расстояние по Y до A[25] оказалось больше имеющегося - эл-ты [0..25] пропускаем,

Что Вы об этом думаете?

[kramer]

Приятно говорить с понимающим человеком.

Спасибо.

Но это годится для немного но все же др задачи: найти на расстоянии не больше заданного. Вот тогда попав в "пустой" бакет (и проверив соседей) спокойно уходим - требуемых ближайших нет. А насчет "по спирали" - там только проверка соседей геморрой не маленький, напр в 3D их уже 26

Да, вообще, пространственный поиск - это непростая задача. Однако алгоритм "поиск ближайшей точки" практически всегда можно свести к алгоритму "поиск точек в радиусе", главное суметь выбрать начальный радиус. Вообще, с этой сеткой алгоритм тоже не прост - есть масса нюансов. Но с точки зрения алгоритмизации, на мой взгляд, он один из самых легких для осознания.

А мне кажется можно. Пусть есть 100 эл-тов (0..99), напр тех же QPointF. Схема поиска:

У вас получилось как раз медианное К-мерное дерево - т.е. intrusive algorithm именно таким образом сортирует исходный массив. Но там все же есть само дерево, которое хранит отношения родитель-потомок. Я что-то не могу сходу сообразить, нафига оно нужно, но раз есть, значит нужно, наверное. Загляните на вики, почитайте про kd-trees, там и алгоритм поиска описан. В частности, отсечения по расстояниям до соответствующей оси там описаны именно так, как у вас. Я этот алгоритм реализовывал, он в принципе несложный, 100 строк кода всего. Однако хочу сразу предупредить насчет граблей - классическое медианное К-дерево хорошо работает только для однородного распределения точек. Т.е. если у вас облако точек, вытянутое вдоль одной из осей, или, не дай бог, два локализованных облака точек на некотором расстоянии друг от друга, эффективность падает в разы, если не в сотни раз.

[kramer]

По поводу покоординатного индексирования - вот статейка, которой я руководствовался, http://www1.cs.columbia.edu/CAVE/publications/pdfs/Nene_TR95.pdf

Если отсутствие накладных расходов по памяти у вас не принципиальное требование, а просто нежелание городить огород с деревьями, и не слишком высокие требования к эффективности по построению - я думаю, вам вполне подойдет. Вообще, он разрабатывался для пространств высокой размерности, поэтому на плоскости будет не слишком эффективен, но зато бесплатный бонус - для одиннадцатимерного пространства будет работать всего в 11 раз медленнее

Вкратце, на каждое измерение вам нужно будет завести по два массива индексов, прямой и обратный, затем отсортировать их по соответствующим координатам. Поиск работает так - для заданной точки вырезаете кусок массива прямых индексов по первой оси, в который входят точки на некотором расстоянии dist от нее, потом через обратный массив для второй оси выкидываете из них те, что не входят в такой же кусок для второй оси. Получаются точки, лежащие внутри квадрата (x-dist, y-dist)(x+dist, y+dist). Для них уже тупо перебор. Понятно, это тоже поиск внутри заданного радиуса, но начальный радиус для поиска можно найти так dist=max(min(dx), min(dy), где dx, dy - расстояния между соседними точками в отсортированном массиве.

UPD: Наврал в конце, все проще - для данной точки (x,y): dist = max(abs(lower_bound_x - x), abs(upper_bound_x - x)), а еще лучше минимум из таких максимумов по всем осям. Получается один лишний вызов lower_bound на ось.

[Igors]

Однако хочу сразу предупредить насчет граблей - классическое медианное К-дерево хорошо работает только для однородного распределения точек. Т.е. если у вас облако точек, вытянутое вдоль одной из осей, или, не дай бог, два локализованных облака точек на некотором расстоянии друг от друга, эффективность падает в разы, если не в сотни раз.

Давайте подробнее. Вот ситуация

здесь много точек массива (куча/сгусток)
**********
**********
**********


                                  х(запросная точка)

                                                     
                                                                     *(еще 1 точка массива)

Ну вот оказался один урод далеко, фишка так легла. Дерево делится "по числу эл-тов" поэтому ось/плоскость деления все время будет находиться в "куче", мин расстояние до запросной точки остается большим, размер по оси ее превысить не может. В результате дерево переберет все (или почти) точки. Думаю это то самое "в сотни раз" о котором Вы говорите.

Но все это только если не задано исходное (достаточно малое) мин расстояние. Иначе все норм, имеем пресловутый log N - ведь при каждом делении одна из половинок кучи имеет расстояние по оси больше заданного.

Др ситуации
Два облака: ну если они соразмеримы по числу точек - все норм

Вытянуто по оси - скорость поиска вдвое упадет, т.к. деление по одной из осей ничего не дает. Но для этого есть решение: делить не абы как, а по большей из осей. Но тогда индекс деления надо где-то хранить (напр в самой точке) - а это доп данные

Оказывается "идеально сбалансированное дерево" - необязательно хорошо/здорово 

[kramer]

Давайте подробнее. Вот ситуация

Насколько смогу. Я этим занимался два года назад, тонкостей, естественно, уже не помню.

Ну вот оказался один урод далеко, фишка так легла.

Да-да, вот это - один из самых мерзких случаев. При локальных облаках, даже соизмеримых по числу точек, производительность падает несколько меньше, но тоже порядочно. По той же причине - если первоначальный кандидат лежит в одном облаке, а реально ближайшая точка - в другом - придется перебрать все дерево. Бороться с этим можно кластерингом - т.е. вычислять обходом вширь идентификатор принадлежности к такому локальному облаку для каждой точки, и для каждого облака строить отдельное дерево. Но обход вширь дорогого стоит, мне просто эта информация нужна была для других целей, поэтому я ей воспользовался, так сказать, попутно. Но и в этом случае тонкостей тоже масса.

Что до разбиения по максимальной оси - опытным путем у меня не получилось получить прироста, а при некоторых конфигурациях производительность даже падала. Реализуется это элементарно, одним ifdef'ом, разница в 3 строчки. Я долго пытался понять, почему это происходит, но эффект проявляется только на больших облаках, а средств визуализировать дерево у меня тогда еще не было. В конце концов плюнул, оставил (i%3) - работает, и ладно.

Сбалансированное дерево гарантирует лучшую "среднюю" производительность, однако в крайних случаях проседает сильнее, чем несбалансированное.

А скажите, какой у вас порядок количества точек, хотя бы вилку?

[Igors]

Если отсутствие накладных расходов по памяти у вас не принципиальное требование, а просто нежелание городить огород с деревьями, и не слишком высокие требования к эффективности по построению - я думаю, вам вполне подойдет. Вообще, он разрабатывался для пространств высокой размерности, поэтому на плоскости будет не слишком эффективен, но зато бесплатный бонус - для одиннадцатимерного пространства будет работать всего в 11 раз медленнее

Вкратце, на каждое измерение вам нужно будет завести по два массива индексов, прямой и обратный, затем отсортировать их по соответствующим координатам. Поиск работает так - для заданной точки вырезаете кусок массива прямых индексов по первой оси, в который входят точки на некотором расстоянии dist от нее, потом через обратный массив для второй оси выкидываете из них те, что не входят в такой же кусок для второй оси. Получаются точки, лежащие внутри квадрата (x-dist, y-dist)(x+dist, y+dist). Для них уже тупо перебор. Понятно, это тоже поиск внутри заданного радиуса, но начальный радиус для поиска можно найти так dist=max(min(dx), min(dy), где dx, dy - расстояния между соседними точками в отсортированном массиве.

UPD: Наврал в конце, все проще - для данной точки (x,y): dist = max(abs(lower_bound_x - x), abs(upper_bound_x - x)), а еще лучше минимум из таких максимумов по всем осям. Получается один лишний вызов lower_bound на ось.

Можно еще проще

- сортируем по одной из осей, напр по X
- с помощью lower_bound получаем ближайшую по X
- дальше просматриваем массив вперед и назад

Код

C++ (Qt)
int index = std::distance(pt.begin(), std::lower_bound(pt.begin(), pt.end(), CompareX));
int minIndex = index;
double minD = Near(cntr, pt[index]);
for (int i = index + 1; i < pt.size(); ++i) {
 double dx = pt[i].x() - cntr.x();
 if (dx * dx > minD) break;
 double d = Near(cntr, pt[i]);
 if (d > minD) continue;
 minD = d;
 minIndex = i;
}

Просмотр назад аналогичен, опускаем

Насколько я понял, в статье предлагается улучшить этот алгоритм за счет замены Near на более быструю проверку по индексу. Согласен, какое-то ускорение есть, но, во-первых, оно не так уж велико (расчетов в Near мало, можно еще уменьшить). Во-вторых в коде выше нам может повезти если minD быстро уменьшится, а по методе статьи - все равно тратиться на поиск по Y.  В любом случае "ну очень быстрого поиска" ожидать не приходится - число проверяемых может быть просто велико. Но зато простота реализации.

Тут такой интересный вопрос: а как удачнее выбрать ось для первого поиска? (от этого производительность сильно зависит). А еще лучше "кластеризоваться"...

[kramer]

- сортируем по одной из осей, напр по X
- с помощью lower_bound получаем ближайшую по X
- дальше просматриваем массив вперед и назад

Тут плохая асимптотика, линейная в худшем случае, а он вполне вероятен.

Еще короче суть можно изложить так - для каждой оси получаем фрагмент прямого (отсортированного) массива лежащий между плоскостями (x-dist, x+dist), а потом с помощью обратного массива ищем пересечение этих фрагментов, т.е. точки, которые входят во все такие фрагменты по всем осям. Обратный массив нужен для получения из индекса точки в отсортированном массиве ее индекса в изначальном. Асимптотика - логарифм + перебор в полученном массиве.

По опыту - работает для сравнительно малых радиусов примерно также, как KD-дерево. Проверял, правда, только на малом числе точек, до тысячи примерно, такова специфика задачи. Зато пишется все за полчаса, все на стандартных контейнерах и алгоритмах из STL.

Тут такой интересный вопрос: а как удачнее выбрать ось для первого поиска? (от этого производительность сильно зависит). А еще лучше "кластеризоваться"...

На первый потолочный взгляд - выбирать наибольшую по линейным размерам, при равномерном распределении на ней будет меньше всего точек в интересующем нас регионе.

[Igors]

А скажите, какой у вас порядок количества точек, хотя бы вилку?

"Общий случай" т.к. размер сцены не ограничен. На практике мульены - дело рядовое.

Тут плохая асимптотика, линейная в худшем случае, а он вполне вероятен.

асимптотика - это "O(..)"?  То для бездельников-теоретиков  

Еще короче суть можно изложить так - для каждой оси получаем фрагмент прямого (отсортированного) массива лежащий между плоскостями (x-dist, x+dist), а потом с помощью обратного массива ищем пересечение этих фрагментов, т.е. точки, которые входят во все такие фрагменты по всем осям. Обратный массив нужен для получения из индекса точки в отсортированном массиве ее индекса в изначальном. Асимптотика - логарифм + перебор в полученном массиве.

Вовсе нет. Только логарифм (точнее 2 * log) только по первой оси, по остальным плюс перебор того что попало в диапазон первой. Напр по X имеем 10K точек. Значит в цикле до 10K будем проверять попадает ли индекс(ы). Пусть проверки быстрые (2 обращения к массиву) - но 10K есть 10K, сожрет прилично

И все это только если задано dist - но его запросто может не быть. В общем - слабовато

По опыту - работает для сравнительно малых радиусов примерно также, как KD-дерево. Проверял, правда, только на малом числе точек, до тысячи примерно, такова специфика задачи.

Ну при малых радиусах и поиск по одной оси работает прилично  

[kramer]

асимптотика - это "O(..)"?  То для бездельников-теоретиков

Это всего лишь способ оценки роста сложности. Если точек может быть миллион - это становится важным.

Вовсе нет. Только логарифм (точнее 2 * log) ...
И все это только если задано dist - но его запросто может не быть. В общем - слабовато

Как посчитать радиус поиска для заданной точки - я там выше писал. Обойдется в лишний вызов lower_bound на ось минус один, так как один можно переиспользовать. В асимптотике константные множители по-моему не учитываются, на то она и асимптотика. Она характеризует рост сложности при росте количества элементов, а не саму сложность. Т.е. O(N) > O(K*log(N)) для любого конечного K.

Ну при малых радиусах и поиск по одной оси работает прилично

Сортировка по одной оси, на мой взгляд, вообще ничего не дает, поскольку самая дальняя вдоль оси точка может вполне оказаться ближайшей по расстоянию, и наооборот, так что вам все равно придется перебрать весь массив.