Tags: Регистрация облаков точек SVD Евклидово расстояние
Введение
Евклидово расстояние
Определение минимального расстояния между парами из множества точек
Регистрация облаков точек. Постановка задачи
Оценка соответствия особых точек
Регистрация 3D облаков точек, используя SVD
Полезные ссылки
Введение
Регистрация облаков точек – процесс совмещения нескольких облаков точек одного объекта в единую систему координат. Цель регистрации — найти преобразование, которое оптимально позиционирует два облака (см. Алгоритмы регистрации облаков точек). При этом необходимо оценить насколько точно определяются параметры преобразования.
Сходство или различие между объектами классификации устанавливается в зависимости от выбранного метрического расстояния между ними. Если каждый объект описывается i-свойствами (признаками), то он может быть представлен как точка в i-мерном пространстве, и сходство с другими объектами будет определяться как соответствующее расстояние. При классификации используются различные меры расстояния между объектами.
Евклидово расстояние
Евклидово расстояние — наиболее часто используемая мера расстояния. Она является геометрическим расстоянием в многомерном пространстве и вычисляется следующим образом:
где:
- L – расстояние между объектами A и B;
- Ai – значение i-свойства объекта A;
- Bi – значение i-свойства объекта B.
Применение Евклидова расстояния оправдано в случае, если свойства (признаки) объекта измерены в одних единицах, при этом одинаково важны для классификации и однородны по физическому смыслу. Пример — расстояния между двумя точками в 2D и 3D пространствах:
Формула Евклидова расстояния для Rn является обобщением теоремы Пифагора.
Взвешенное евклидово расстояние применяется в тех случаях, когда каждому i-свойству удается приписать некоторый «вес» wi, пропорционально степени важности признака в задаче классификации:
Нормализация признака. Естественная, с геометрической точки зрения, евклидова мера расстояния может оказаться бессмысленной, если признаки измерены в разных единицах. Чтобы исправить положение, прибегают к нормированию каждого признака.
Процент несогласия — эта мера расстояния используется в тех случаях, когда свойства (признаки) объекта являются категориальными. Например, первый признак объекта – пол, второй – возраст, третий – место работы. Представим значения свойств (признаков) объекта в виде вектора значений. Первый вектор – (муж, 20 лет, учитель), второй вектор – (муж, 28 лет, менеджер). Процент несогласия равен 2/3. Эти вектора различаются на 66.6%.
Определение минимального расстояния между парами из множества точек
Постановка задачи. Пусть задано в 2D пространстве множество точек. Необходимо найти пару точек, Евклидово расстояние между которыми будет минимальным.
Решение задачи. Определяем все возможные комбинации пар точек. Для каждой из них находим Евклидово расстояние. Перебором определяем пару точек с минимальным Евклидовым расстоянием.
Эта простая задача лежит в основе решения более сложной задаче, рассмотренной ниже.
Регистрация облаков точек. Постановка задачи
При создании панорамного фото или в фотограмметрии часто требуется решить задачу соединения нескольких облаков точек в одно. Эта задача обычно решается в 3D пространстве. Проиллюстрируем подход к ее решению на простейшем примере соединения 2-х облаков точек в 2D пространстве.
Постановка задачи
Задано 2 облака точек. Координаты одного облака определены в глобальной системе координат (XY), координаты второго облака – в локальной системе координат (xy).
Известно 2 общие точки (точка 1 и точка 2), по которым облака могут быть совмещены. Эти точки определены как в локальной так и в глобальной системах координат.
Необходимо определить зависимости для переопределения координат точек второго облака из локальной в глобальную систему координат.
По 2-м точкам всегда можно определить взаимное положение систем координат на плоскости. По 3-м точкам можно определить взаимное положение систем координат в в 3D пространстве. Подробнее см. Выделение параметров формы и положения.
Решение задачи
Положение локальной системы координат (ЛСК) относительно глобальной системы координат (ГСК) определяются последовательностью из 3-х преобразований, записанных через матрицы преобразований в однородных координатах:
Из 3-х преобразований можно организовать различную комбинацию последовательностей преобразований. Но не каждая из них приведет к желаемому результату (см., например, Задачу 3).
Уравнения для перехода из системы координат xy в систему координат XY запишутся в следующем виде:
Находим параметры (α, m, n) через решение систем уравнений с заданными координатами точек облаков 1 и 2
Три из 4-х уравнений используем для нахождения параметров (α, m, n). После определения этих параметров можно использовать уравнения (выделены красным) для переопределения координат точек второго облака из локальной в глобальную систему координат.
Эта задача не определяет соответствие между точками изображений в панораме или фотограмметрии. Она лишь демонстрирует подход к решению подобных задач. Более детально ознакомиться с геометрическими преобразованиями в 3D пространстве, имеющими отношение к созданию панорамного фото или фотограмметрии, можно в статьях:
Оценка соответствия особых точек
При регистрации облаков точек по особым точкам определяются матрицы преобразования. И от того, насколько корректно выбраны особые точки, зависит точность регистрации облаков. Ниже рассмотрен алгоритм, который позволяет оценить степень соответствия особых точек.
Выше была рассмотрена задача, в которой 2 общие точки облаков (точки 1 и 2) считались однозначно заданными. Усложним задачу. Представим, что у нас не 2, а 4 общие точки облаков. При этом, некоторые из точек могут быть определены ошибочно, как общие точки облаков. Необходимо отфильтровать ошибочные точки, а корректные точки использовать для соединения 2-х облаков.
Алгоритм решения задачи.
Суть алгоритма. Если параметры (α, m, n) определились одинаковыми по крайней мере через 2 пары точек, то точки в этих парах определены корректно.
Подробное описание алгоритма:
- Из 4-х точек выделяем все возможные комбинации пар точек (6 комбинаций)
1-2 1-3 2-3 1-4 2-4 3-4
- Для каждой пары находим параметры (α, m, n) через решение систем уравнений с заданными координатами точек этих (см. Регистрация облаков точек. Постановка задачи).
- Определяем меру соответствия каждой пары точек с другой парой через Евклидово расстояние.
Каждая пара определяется 3-мя параметрами (α, m, n) в 3-мерном пространстве. Параметры (m, n) и (α) измерены в разных единицах — в единицах длины и градусах (или радианах). Поэтому определяем два Евклидовых расстояния:
- Если расстояния L1 и L2 находятся в пределах допуска, то точки в обеих сравниваемых парах соответствуют, а параметры, которые по ним определены, можно использовать для переопределения координат точек из локальной системы координат в глобальную.
Регистрация 3D облаков точек, используя SVD
Перевод с технической обработкой статьи Finding the optimal/best rotation and translation between two sets of corresponding 3D point data
Матрица поворота рассчитывается с помощью сингулярного разложения (SVD-(Singular Value Decomposition). Сингулярным разложением матрицы H размера m×n является:
H = U S Vt, где
S – матрица размера m×n с неотрицательными элементами (элементы, лежащие на главной диагонали, – сингулярные числа, а элементы, не лежащие на главной диагонали, – нулевые);
Матрицы U (размера m) и V (размера n) – унитарные матрицы, состоящие из левых и правых сингулярных векторов соответственно (Vt – транспонированная матрица к V). Так как H – квадратная матрица, то матрицы U, S и V такого же размера, что и H.
Для того, чтобы найти матрицу поворота, необходимо перенести оба набора точек в центр координат, что позволит избавиться от смещения.
Допустим, есть матрицы PA1, PA2, PA3 (точки 1-го облака) и матрицы PB1, PB,2, PB3 (точки 2-го облака), которые имеют вид:
Тогда центроиды (средние точки) могут быть получены следующим образом:
Затем подсчитываем матрицу следующего вида:
Применив к H сингулярное разложение, определяем U и V, после чего можем рассчитать матрицу поворота:
Смещение t находим по формуле:
Представленное решение будет работать как с данными без шума, так и с шумом. Для совмещения облаков точек A и B требуется как минимум N=3 пары соответствующих точек, чтобы определить матрицы поворота R и перемещения t.
Если данные зашумлены, то точек используется N>3. для получения средневзвешенных значений искомых параметров и определения среднеквадратичной ошибки (RMSE):
Метод наименьших квадратов (МНК) может использоваться для «решения» переопределенных систем уравнений (когда количество уравнений превышает количество неизвестных).
Программу для тестирования алгоритма можно скачать по ссылке. Основные этапы ее работы описывается ниже:
- Создаем матрицы вращения и перемещения. Параметры матриц выбираем наугад (Random).
# Random rotation and translation
R = np.random.rand(3,3)
t = np.random.rand(3,1)
print("Random rotation")
print(R)
print("Random translation")
print(t)
Матрица поворота, созданная наугад, не есть корректной с точки зрения жесткости, поскольку только 3 параметра (угла) из 9 параметров независимые.
Преобразование называют жестким (еще движения или Евклидовым), если оно сохраняет форму и размер. В этом отличие от аффинного преобразования, которое включает еще масштабирование и сдвиг, Евклидово преобразование есть частным случаем аффинного преобразования. В свою очередь, преобразование поворота есть частным случаем Евклидового преобразования (см. подробнее Матрицы преобразований).
- Делаем R правильной матрицей вращения — усиливаем жестокость преобразования, используя SVD
# make R a proper rotation matrix, force orthonormal
U, S, Vt = np.linalg.svd(R) # Singular Value Decomposition
R = U@Vt # @ (at) operator for matrix multiplication
print("U")
print(U)
print("S")
print(S)
print("Vt")
print(Vt)
print("Proper rotation")
print(R)
- Создаем матрицу из 3D координат n=5 точек облаков A и B. Параметры матрицы A выбираем наугад (Random). Параметры матрицы B определяем через поворот R и перемещение t облака A
n = 5 # number of points
A = np.random.rand(3, n) # 3 rows and 5 columns
B = R@A + t
print("A и B")
n = 5 # number of points
A = np.random.rand(3, n) # 3 rows and 5 columns
B = R@A + t
print("A и B")
print(A)
print(B)
- Восстанавливаем R и t по множествам точек облаков A и B через вызов функции
ret_R, ret_t = rigid_transform_3D(A, B)
- Через эту функцию сначала определяем центроиды облаков A и B, как среднеарифметические значения координат
centroid_A = np.mean(A, axis=1) # The arithmetic mean is the sum of the elements
centroid_B = np.mean(B, axis=1) # along the axis divided by the number of elements.
print("centroid_A и centroid_B")
print(centroid_A)
print(centroid_B)
- Меняем форму массива
centroid_A = centroid_A.reshape(-1, 1) # The reshape() function is used to give a new shape
centroid_B = centroid_B.reshape(-1, 1) # to an array without changing its data.
print("centroid_A и centroid_B after reshape()" )
print(centroid_A)
print(centroid_B)
- Определяем матрицы преобразования с центроидами в начале координат
# subtract mean
Am = A - centroid_A
Bm = B - centroid_B
Bmt = np.transpose(Bm) # Reverse or permute the axes of an array
print("Am, Bm и Bmt")
print(Am)
print(Bm)
print(Bmt)
print("matrix multiplication H = Am @ Bmt")
H = Am @ Bmt
print(H)
# find rotation
U, S, Vt = np.linalg.svd(H) # Singular Value Decomposition
print("U")
print(U)
print("S")
print(S)
print("Vt")
print(Vt)
VtT = Vt.T
print("VtT")
print(VtT)
UT = U.T
print("UT")
print(UT)
R = VtT @ UT
print("R = VtT @ UT")
print(R)
t = -R @ centroid_A + centroid_B
print("t = -R @ centroid_A + centroid_B")
print(t)
- Возврат в основную программу
# Compare the recovered R and t with the original B2 = (ret_R@A) + ret_t print(B2)
# Find the root mean squared error err = B2 - B err = err * err err = np.sum(err) rmse = np.sqrt(err/n)
Полезные ссылки:
- Мера расстояния
- Finding the optimal/best rotation and translation between two sets of corresponding 3D point data
- Singular value decomposition
- Среднеквадратичная ошибка (RMSE)
- Разработка программной системы для совмещения облака точек
- Распознавание образов и анализ сцен
- 3D reconstruction
Автор: Николай Свирневский