Создание бесшовного изображения
Введение
В последнее время все большую значимость и масштабы принимает использование геопространственной информации и средств географического ориентирования. В первую очередь, с каждым годом увеличиваются объемы собираемых данных, таких как данные, полученные путем дистанционного зондирования, космические снимки и аэрофотосъемка. Также происходит значительное повышение эффективности и качества в способах получения и обработки информации — применяются более совершенные и точные приборы, улучшаются и разрабатываются новые методы фотограмметрии, т.е. методы для вычисления геометрических свойств объектов по их фотографиям, расширяется использование летательных аппаратов и космических спутников. Кроме того, эти данные получают массовое распространение и становятся доступными большому кругу пользователей.
Область применения полученных данных постоянно расширяется, не ограничиваясь только специфическими задачами государственного уровня, а проникая в тесное взаимодействие практически со всеми сферами общественной жизни. Для управленческой и предпринимательской деятельности существенно необходимо знание о месторасположении для достижения наиболее быстрых и правильных решений. Исследование процессов движения земной коры, рельефа, анализ изменения погоды, прокладывание дорог и планирование строительства, предсказывание природных катастроф, проведение спасательных операций — все это трудно проводить в отсутствие надежных средств ориентирования. Уже невозможно представить жизнь людей без доступности использования подробных карт, панорам различных достопримечательностей и улиц городов, таких средств геопространственного ориентирования, как GPS и ГЛОНАСС. Портативные устройства с использованием карт местности с большим успехом используются для отслеживания положения и вычисления оптимального пути, встраиваются в смартфоны, коммуникаторы и автомобильные бортовые компьютеры.
В частности, в настоящее время осуществляется сбор большого количества данных ДЗЗ (Дистанционное зондирование Земли) с БПЛА (Беспилотный летательный аппарат) и самолетов. В большинстве случаев эти данные представляют собой кадровую съемку, представляющую собой набор отдельных фотографий. Иногда может предоставляться 3D модель рельефа, параметры камеры и глобальные координаты. Таким образом, среди многих направлений в эффективной обработке большого количества геопространственных данных возникает задача формирования мозаики — бесшовного единого изображения, полученного путем сопоставления перекрывающихся изображений, в то же самое время имеющего географическую привязку и максимальное совпадение с рельефом. При этом, встает задача поиска связующих точек изображений на маршруте и межмаршрутных изображений или просто создание панорам.
Существует библиотека GDAL (Geospatial Data Abstraction Library), которая используется многими системами по обработке геопространственных данных. Но сейчас в ней отсутствует алгоритм для полного автоматического создания мозаики, и остро стоит необходимость во включении в базовый функционал библиотеки такого алгоритма с реализацией пользовательского интерфейса в виде плагина для QuantumGIS (свободная кроссплатформенная геоинформационная система).
Постановка задачи
Цель работы: прототип алгоритма для автоматического создания единого бесшовного изображения (мозаики) из отдельных фотографий.
В процессе достижения поставленных целей возникают следующие задачи:
- Быстрый и надежный способ нахождения связующих точек между изображениями.
- Применение метода градиентного спуска как алгоритма для нахождения и оптимизации параметров.
- Реализация прототипа программы для автоматической обработки 2х изображений в случае использования аффинных преобразований.
Вводимые ограничения:
- Отсутствие удаления искажений, вызванных оптикой (дисторсии).
- Применение аффинных преобразований.
- Обработка только 2х изображений.
- Использование метода градиентного спуска.
Проект участвует в программе Google Summer Of Code 2012, планируется дальнейшее развитие и интеграция в GDAL с реализацией пользовательского плагина.
Анализ задачи
Сбор большого количества фотографий, в частности при проведении аэрофотосъемки, приводит к важной задаче по исправлению геометрии фотографий, «сшиванию» в единую мозаику, привязке полученного изображения к глобальным географическим координатам и, что является наиболее важным, к задаче по наибольшей автоматизации этого процесса. В идеале, оператор должен предоставить на вход программе фотографии, при возможности телеметрическую информацию (данные GPS, углы съемки) и получить в итоге готовое изображение. Таким образом, данная автоматизация позволит значительно увеличить производительность труда операторов, взяв на себя большую часть рутинных операций по обработке.
Одним из самых трудоемких этапов является процесс нахождения и сопоставления связующих точек — характерных особенностей на изображениях («feature points»). Этот этап предлагается полностью автоматизировать, исключив участие человека, применив специализированные алгоритмы. Также предлагается оптимизировать параметры для трансформации фотографий в единую мозаику. При этом задача оптимизации в общем случае включает использование телеметрии, то необходимо учитывать также и эту информацию.
Таким образом, главной задачей является создание достаточно гибкой и надежной программы, которая бы требовала минимального участия человека, при этом была бы применима к любым данным и могла выдавать приемлемый результат при отсутствии телеметрических данных, либо какой-то другой известной информации для первоначального ориентирования фотографий.
Анализ предметной области
Для геопространственных систем необходима так называемая ортофотомозаика — единое изображение, имеющее ортогональную проекцию. Для построения ортофото и выпрямления в соответствии с поверхностью Земли используются специальные наземные измерения, применяются алгоритмы триангуляции и стереозрения для восстановления рельефа по перекрывающимся изображениям — происходит создание так называемой DEM (Digital Elevation Model). В том числе применяется уравнивание сети точек по блоку изображений. Также производится учет данных телеметрии для географической привязки фотографий к местности при помощи GCP (Ground Control Points) — это точки, которые имеют координаты как на фотографии, так и соответствующие глобальные координаты на карте.
Анализ существующих решений
Так как задача не является сравнительно новой, есть множество наработок в данной области — концепций, различных подходов и алгоритмов. В специализированных системах используется один из следующих подходов, он будет рассмотрен на примере цифровой фотограмметрической системы PHOTOMOD. Нахождение связующих точек состоит из 5 этапов:
- Сегментация изображений — разбиение изображений на области с примерно одинаковыми фотометрическими характеристиками.
- Полученные области переводятся в аппроксимированные контурные фигуры.
- Для каждой пары перекрывающихся изображений вычисляются коэффициенты корреляции определенным образом выбранного подмножества пар фигур.
- С помощью оптимизированного переборного алгоритма, использующего рассчитанный коэффициент корреляции, а также дополнительные геометрические особенности расположения фигур, производится окончательное сопоставление фигур пары снимков и вычисление приблизительного соответствия точек на паре изображений.
- Результаты этапа 4 поступают на вход иерархическому площадному коррелятору, с помощью которого полностью автоматически вычисляется необходимое число точек. После этого применяются различные критерии отбраковки, тем самым исключаются потенциально ошибочные точки.
Далее происходит расположение фотографий согласно полученным расчетам.
Недостатками такого похода являются такие моменты, как достаточно трудоемкий и аппроксимирующий расчет первоначальных контурных фигур для последующего сопоставления и невозможность найти достаточно большое количество пар этих фигур для повышения точности ориентирования изображений. Более того, для такого коррелятора в ряде случаев затруднен поиск связующих точек, и приходится осуществлять коррекцию вручную, что серьезно замедляет обработку фотографий.
HUGIN Software
Данная программа позволяет производить множество сложных преобразований: создание мозаики из перекрывающихся изображений и особенно панорам, нахождение контрольных точек, HDR (High dynamic range imaging). Проект ориентирован на любительские и профессиональные цели и не предназначен для обработки геопространственной информации.
Ceres-Solver
Применительно к нашей задаче по оптимизации преобразований, существует open-source проект Ceres-Solver. Он представляет собой мощный инструмент для решения нелинейных проблем по минимизации функций. Для использования необходимо задать свою функцию стоимости, а также производные по искомым параметрам, предусмотрена возможность автоматического нахождения частных производных. При этом для сходимости используются алгоритмы, позволяющие достичь решения за минимальное число шагов. В дальнейшем вполне может быть рассмотрена возможность использовать этот проект в качестве алгоритма для оптимизации параметров.
Описание и структура алгоритма
Фотографии, получаемые с БПЛА, обычно снимаются с очень большой частотой, и изображения получаются с сильным перекрытием. Также из-за небольшого изменения угла съемки можно считать, что фотографии представляют собой приближенную ортогональную проекцию. Таким образом, на данном этапе мы можем пренебречь сложными полиномиальными преобразованиями, ограничившись аффинными, которые были выбраны из-за своей простоты и быстроты реализации для демонстрации работы прототипа, и достигнуть вполне приемлемых результатов с относительно небольшой погрешностью в некоторых случаях. При этом обратим внимание, что преобразования не должны быть слишком большими, в противном случае исходные фотографии будут неправдоподобно сильно искажены и не соответствовать реальной поверхности.
Для упрощения алгоритма также принимается во внимание, что искажения, вызванные оптикой фотокамеры (такие как дисторсия объектива), исправляются до работы алгоритма. Это позволяет приблизиться к правильной ортогональной проекции и повысить точность конечной мозаики.
Общая последовательность работы алгоритма предлагается следующей (некоторые пункты в данной работе не рассматриваются):
- Исправление искажений, вызванное оптикой.
- Компрессия изображений в случае сверхвысоких разрешений для ускорения оптимизации.
- Нахождение связующих точек на изображениях.
- Анализ и сопоставление полученных на предыдущем этапе точек, создание пар совпавших точек.
- Использование первоначальной информации о местоположении и углах съемки для предварительной ориентации изображений и определения начальных параметров для проведения более быстрого и надежного способа оптимизации.
- Расчет минимально необходимых параметров трансформации для каждой пары точек для достижения заданной степени погрешности конечной мозаики при помощи итеративного метода градиентного спуска.
- Применение полученных параметров по отношению к точкам к целым изображениям и окончательное формирование единой мозаики.
- Выравнивание яркости изображений и сглаживание швов.
В прототипе реализованы и используются пункты 3, 4, 6, 7, при этом пункты 3-4 - с использованием стандартных средств библиотеки компьютерного зрения OpenCV.
Нахождение связующих (особых) точек («feature points»)
Связующие (или особые, «feature») точки — это точки на изображении, которые отличаются от остальных по какому-либо признаку. Они называются связующими в связи с тем, что используются для связи и дальнейшего сопоставления фотографий друг с другом. Обычно для каждой точки оценивается некоторая уникальным образом характеризующая ее окрестность. Другими словами, такую окрестность какой-либо точки (и соответственно саму такую точку) можно отличить от других окрестностей остальных точек изображения. Существует достаточно много алгоритмов для поиска особых точек, например FAST corner detection, SIFT, SUSAN и т.д. В прототипе используется алгоритм SURF (Speeded up Robust Features).
На этом этапе производится поиск связующих точек на каждом из входных изображений.
Сопоставление и создание пар связующих точек
После нахождения связующих точек (feature points) требуется их сопоставление для получения пар точек. Одна такая пара представляет собой координаты одной и той же точки, но на двух различных перекрывающихся изображениях. Эти пары связующих точек в дальнейшем будут использоваться в методе оптимизации параметров преобразования, являясь как данными для нахождения параметров, так и основой для оценки правильности конечной мозаики. Для сопоставления используются точки, полученные на предыдущем этапе.
Одним из быстрых и надежных алгоритмов является алгоритм «ближайшего соседа» (NNS - Nearest Neighbor Search). Он позволяет в разы увеличить производительность по сравнению с методами простого перебора (brute force). В данной реализации используется специализированная библиотека FLANN (Fast Library for Approximate Nearest Neighbors) для сопоставления точек.
Используемые преобразования
В данной работе в качестве доказательства работы алгоритма используются линейные преобразования изображений (аффинные преобразования) . Преобразования в общем случае для двумерного пространства имеют следующий вид:
(1)
|
---|
где – координаты точки до преобразования,
- преобразованные координаты
- параметры трансформации для 1-го и 2-го изображений соответственно.
Если представить в более удобном матричном виде, получим:
(2)
|
---|
Аффинные преобразования позволяют осуществлять линейные искажения, такие как поворот, масштабирование, растяжение и параллельный перенос.
Для нашего случая наиболее оптимальным и простым с точки зрения наглядности и последующего применения алгоритма является следующий порядок преобразования:
- Осуществление поворота на угол
- Масштабирование и растяжение.
- Параллельный перенос.
Также введем ограничение, что изображение находится в центре его собственной системы координат. Это позволит при расчетах не смещать центр изображения и упростить некоторые вычисления, тем самым позволив алгоритму быстрее находить оптимальные параметры.
Для осуществления преобразований введем вектор параметров :
(3)
|
---|
- угол поворота,
- масштабирование изображения,
- искажение (растяжение и перекос),
- координаты нового центра изображения относительно старой системы координат (параллельный перенос).
В данной работе будем применять такие аффинные преобразования:
- Матрица поворота на угол альфа.
(4)
|
---|
- Матрица для изменения размеров изображения и перекоса.
(5)
|
---|
- Матрица переноса изображения. Координаты новой системы координат в старой системе координат.
(6)
|
---|
Скомбинировав определенные выше матрицы трансформации в этом порядке, получим целостное выражение для трансформации изображений.
(7)
|
---|
где
- пиксельные координаты точки на исходном изображении
- новые, преобразованные координаты точки
- сформированная матрица трансформации, отражающая порядок преобразования, имеет вид:
(8)
|
---|
Применив матрицу (8) для каждой точки фотографии, мы получим новое трансформированное изображение.
Использование данных телеметрии
Анализируя полученные параметры трансформаций, можно с небольшими изменениями использовать данные телеметрии, получаемые при проведении аэрофотосъемки.
- Использование известного угла поворота, подстановка его в неизменном виде в начальные параметры преобразования.
- Использование координат GPS, определяющих местоположение фотографии для расчета относительных параметровво внутренней системе координат для параллельных переносов изображений на свои приближенные позиции.
- Нахождение и сопоставление точек только в тех областях изображений, которые перекрываются, исходя из ориентирования во внутренней системе координат. Данный шаг может значительно увеличить производительность алгоритма.
Также для того, чтобы сохранить географическую привязку, ограничим изменение угла поворота и параллельные переносы:
(9)
|
---|
где - константы для ограничения соответствующих параметров.
Таким образом, телеметрическая информация может способствовать ускорению работы предложенного метода создания мозаики и быть использована практически без дополнительного изменения.
Определение функции стоимости (Cost Function)
В задаче по оптимизации функция стоимости является главной частью алгоритма, определяющей оценку качества составленной мозаики и наискорейшее приближение параметров к оптимальным. Определим функцию как зависимость от квадратов расстояний между соответствующими точками на уже преобразованных изображениях.
Тогда расстояние будет определяться выражением:
(10)
|
---|
Или, если добавить матрицы преобразований (8):
(11)
|
---|
Так как расчет производится для нескольких пар точек, произведем суммирование по m — количеству найденных пар связующих точек.
Функция стоимости для 2х перекрывающихся изображений будет иметь вид:
(12)
|
---|
Применительно к нашей задаче, нахождение минимальной погрешности является частичной оценкой функции стоимости, так как требуется дополнительный учет величины параметров. В случае, когда изображения имеют приблизительную ортопроекцию, нет необходимости в сильном искажении. Обратное будет свидетельствовать о неправильности найденного решения. Таким образом, необходимо ограничить изменение параметров, не позволяя алгоритму изменять их на значительную величину.
Введем дополнительное слагаемое для регулирования параметров, определяемое как сумму квадратов величин значений параметров:
(13)
|
---|
Где - параметр, влияющий на величину изменения параметров. Чем он больше, тем меньшие параметры будет находить алгоритм;
N – количество параметров трансформации (длина вектора).
Таким образом, учитывая все слагаемые, получаем окончательный вид функции стоимости:
(14)
|
---|
Метод градиентного спуска (Gradient descent method)
Этот алгоритм использован в данной работе для нахождения параметров трансформации. Он является простым в реализации, применим к пространствам с неограниченным количеством измерений и позволяет достигать минимума (в некоторых случаях - локального) за некоторое количество шагов (итераций).
Метод предполагает итеративное одновременное обновление (изменение величин) всех параметров в строну наискорейшего уменьшения функции стоимости до тех пор, пока не будет достигнуты заданная точность или схождение алгоритма (параметры практически не будут изменяться за каждую итерацию). Для этого необходимо вычислить частные производные по каждому параметру и произвести изменение в соответствии с формулой:
(15)
|
---|
, где - параметр, влияющий на скорость сходимости алгоритма — чем он больше, тем за меньшее количество итераций метод достигнет оптимума.
Более подробно, применительно к нашим вычислениям, для вектора параметров необходимо вычислить следующие производные:
(16)
|
---|
Необходимо учитывать, что эти частные производные только для изменения параметров одного изображения. В случае двух и более изображений нам необходимо одновременно пересчитывать параметры для каждого изображения, таким образом число производных будет равно , где N — количество параметров, k — число изображений.
Так как в данном случае для двух фотографий число параметров возрастает до 14, на графике работу градиентного спуска сложно отобразить. Упрощенно метод можно проиллюстрировать на примере графика зависимости функции стоимости от вектора в случае 2х параметров:
На вертикальной оси показана величина функции стоимости (, на горизонтальных — величины параметров. Как видно из рисунка, алгоритм начинает свою работу с первоначальных параметров () и соответствующей им «вершины» поверхности (где величина функции стоимости является большой) и за несколько итераций достигает минимума («спускается» к минимуму) в точке со значениями.
Дальнейшие улучшения
Так как представленный проект является прототипом (proof-of-concept), то существуют способы для ускорения и совершенствования работы алгоритма. Основным проблемным местом является метод градиентного спуска, который имеет тенденцию к очень медленному схождению и к нахождению локального минимума. В будущем предполагается реализация быстрого алгоритма Левенберга — Марквардта, направленного на решение задач о наименьших квадратах, а также, что является наиболее важным, переход от обычных аффинных преобразований к перспективным, которые позволят существенно повысить точность конечной мозаики и ее правдоподобность. Также представляется возможным расширение алгоритма для случая с большим количеством изображений.
Средства реализации
Основной алгоритм написан на языке Python в среде Eclipse с использованием плагина PyDev, тестирование проводилось на дистрибутиве Linux Xubuntu. Для использования вспомогательных алгоритмов применялась библиотека компьютерного зрения OpenCV и библиотека NumPy для проведения быстрых вычислений над матрицами.
Заключение
В рамках данной курсовой работы был разработан прототип (proof-of-concept) программы для создания бесшовной мозаики для случая 2х изображений. Были реализованы такие задачи, как:
- Определение функции стоимости для оценки погрешности в сопоставленных изображениях.
- Реализован прототип алгоритма для минимизации параметров трансформации и нахождение оптимальных преобразований методом градиентного спуска в случае 2х изображений.
Прототип успешно выполняет работу по сопоставлению изображений за достаточно небольшое время.
Дальнейшие планируемые изменения и дополнения для интеграции в GDAL в рамках GSoC 2012:
- Реализация алгоритмов поиска и сопоставления точек.
- Расширение алгоритма для большего числа изображений.
- Применение более совершенных методов по минимизации функции стоимости (Conjugate gradient, L-BFGS, алгоритм Левенберга-Марквардта или др.)
- Реализация пользовательского плагина к Quantum GIS.
Список литературы
- Gary Bradski, Adrian Kaehler. Learning OpenCV. Computer Vision with the OpenCV Library- O'Reilly Media - September 2008 – 580 с.
- Robert Laganière. Opencv 2 Computer Vision. Application Programming Cookbook - Packt Publishing – 298c.
- Lengyel E. Mathematics for 3D Game Programming and Computer Graphics, Third Edition - 2011 – 546c.
- Mark Lutz. Programming Python, 4th Edition - December 2010 – 1632c.
- Безменов В.М. Фотограмметрия. Построение и уравнивание аналитической фототриангуляции — КГУ, Казань — 86с.
- Документация по NumPy - (http://numpy.scipy.org/)
- Andrew Ng. Курс онлайн лекций по машинному обучению — (http://ml-class.org)
- А.Ю. Сечин, М.А. Дракин, А.С. Киселева. Беспилотные летательные аппараты: применение в целях аэрофотосъемки для картографирования. «Ракурс», Москва, Россия
- Herbert Bay, Andreas Ess, Tinne Tuytelaars, and Luc Van Gool. Speeded-Up Robust Features (SURF) – 2008 – 14с.
- Marcus ARTHUR, Barbados, Raid AL-TAHIR, and Dexter DAVIS. Rapid Processing of Unmanned Aerial Vehicles Imagery for Disaster Management.
- U. Niethammer, S. Rothmund, U. Schwaderer, J. Zeman, M. Joswig. Open source image-processing tools for low-cost uav-based landslide investigations – 2011