CC BY
146
СЛЕПАЯ КОМПЕНСАЦИЯ РАДИАЛЬНОЙ ДИСТОРСИИ НА ОДИНОЧНОМ ИЗОБРАЖЕНИИ С ИСПОЛЬЗОВАНИЕМ БЫСТРОГО ПРЕОБРАЗОВАНИЯ ХАФА
И.А. Кунина, С.А. Гладилин, Д.П. Николаев Институт проблем передачи информации РАН, Москва, Россия
Аннотация
В работе представлен способ автоматической компенсации радиальной дисторсии - искажения, характерного для широкоугольных объективов. Предложенный алгоритм оценивает параметры искажения, используя только единичное изображение, полученное из неизвестного источника. При этом не используются какие-либо калибровочные объекты, но предполагается, что исходная сцена содержит прямые линии. Суть метода заключается в поиске радиальной дисторсии с таким набором параметров, что ее исправление даст изображение, на котором «в целом» линии выглядят более прямыми. Для оценки общей искривленности линий в данной работе предлагается использовать быстрое преобразование Хафа, при этом собственно линии на изображении не выделяются. Предложенный алгоритм был проверен на реальных данных, полученных камерами с откалиброванными объективами, имеющими различную степень радиальной дисторсии. Для формальной оценки работы алгоритма был предложен функционал качества компенсации геометрического искажения, нечувствительный к возможной плохой обусловленности задачи определения коэффициентов модели дисторсии.
Ключевые слова: цифровая обработка изображений, анализ изображений, конструкция оптической системы, радиальная дисторсия, автоматическая калибровка, быстрое преобразование Хафа.
Цитирование: Кунина, И.А. Слепая компенсация радиальной дисторсии на одиночном изображении с использованием быстрого преобразования Хафа / И.А. Кунина, С.А. Гладилин, Д.П. Николаев // Компьютерная оптика. - 2016. - Т. 40, № 3. - С. 395-403. - Б01: 10.18287/24126179-2016-40-3-395-403.
Введение
В настоящее время во многих отраслях экономики все большее распространение получают интеллектуальные системы обработки изображений и видеопотоков. При этом для практических приложений оптимальным является использование фото- и видеокамер с оптикой невысокого качества, характеризующейся в том числе значительными дисторсиями, компенсация которых осуществляется в процессе последующей обработки изображений. Для выполнения такой компенсации необходима информация о параметрах искажения. В отдельных случаях она может содержаться в паспортных данных камеры, но чаще всего для определения этих параметров необходимо выполнять калибровку оптической системы.
В данной работе рассматривается задача компенсации радиальной дисторсии - несовершенства оптической системы, вызванного сферической формой линз объективов. В результате действия радиальной дисторсии нарушается геометрическое подобие между объектом и его изображением, прямые линии на изображении становятся кривыми (за исключением прямых, проходящих через оптический центр кадра), при этом степень искривления увеличивается от центра к периферии. Изображения, полученные с широкоугольных камер, не предназначенных для измерительных целей, как правило, подвержены в первую очередь именно радиальной дисторсии.
Наиболее распространенным способом компенсации дисторсий камеры является предварительная калибровка с использованием специального калибровочного объекта, помещаемого в поле зрения камеры
[1, 2]. В качестве калибровочного объекта может выступать как периодическая структура [1], так и случайная текстура с определенными статистическими свойствами [2]. Характеристики такого объекта должны быть известны заранее.
Также находят свое применение и методы, не требующие специального калибровочного объекта, но использующие несколько зарегистрированных изображений одной и той же сцены. Такие методы опираются на априорную информацию о движении камеры или геометрии сцены и решают задачу калибровки с учетом ограничений эпиполярной геометрии [3, 4].
Однако на практике часто возникает ситуация, когда описанные методы неприменимы, т.к. отсутствует возможность получить с камеры требуемые изображения. Так бывает, если доступно только изображение сцены, а информация об использованной камере отсутствует. В этом случае параметры дисторсии могут быть оценены на основе анализа самого изображения. Такой аналитический способ исправления геометрических искажений получил название автоматической (или слепой) калибровки. К таким способам относят анализ изображения в частотной области [5], исследование структуры изображения [6] и другие. Однако наибольшее распространение получили методы, основанные на предположении, что типичная сцена содержит большое количество прямых линий, которые остаются прямыми при центральном проецировании, но приобретают кривизну вследствие радиальной дисторсии. Такое предположение характерно для многих практических задач, например, оптической навигации летательных аппаратов [7-9]. В работах [10-12] предлагается делать пробные ком-
пенсации радиальной дисторсии в различных предположениях о ее параметрах и оценивать кривизну линий на восстановленном изображении через вычисление его Хаф-образа. В работе [13] предлагается похожий подход, но вместо Хаф-образа вычисляется гистограмма ориентированных градиентов. В работе [14] авторы не делают пробных восстановлений, а детектируют линии и аппроксимируют извлеченные кривые дугами окружностей с последующим вычислением параметров радиальной дисторсии.
В данной работе представлено развитие алгоритма [12], основанного на быстром преобразовании Хафа (БПХ) [15, 16].
1. Алгоритм компенсации радиальной дисторсии
Предлагаемый метод, как и его предшественник, выполняет пробные компенсации радиальной дистор-сии в различных предположениях о ее параметрах и затем с помощью аппарата БПХ оценивает качество исправленного изображения. Учет изменения элемента длины кривой при исправлении радиальной дисторсии позволил повысить точность алгоритма и перейти к трехпараметрической модели радиальной дисторсии.
1.1. Модель дисторсии
Будем считать, что оптическая ось линзовой системы пересекает изображение в его геометрическом центре. Расположим начало координат в центре, а оси х и у направим вправо и вниз соответственно. В качестве математической модели дисторсии воспользуемся классической моделью Брауна [17]. В ней координаты на радиально искаженном изображении выражаются через известные координаты на оригинальном (не содержащем искажения) изображении следующим образом:
Г = / (г) = (1 + !>/2'' )г
1=1
•ха = х • О /г , (1)
Уа = У • га / г
где г = ^х2 + у2, к = 1...и - параметры дисторсии,
(х,у) - оригинальное положение точки, (ха,уа) - положение точки в результате действия радиальной дисторсии.
Если га< г, то такой эффект называют отрицательной («бочкообразной») дисторсией, в противном случае (при га> г) - положительной («подушкообразной»). При этом отрицательная дисторсия характерна для широкоугольных объективов. Положительная дистор-сия, характерная для телеобъективов, является менее распространённой в непрофессиональной среде (где возможна утеря информации о параметрах объектива), и в настоящей работе рассматриваться не будет.
Значения коэффициентов к зависят от фокусного расстояния камеры, но эта зависимость не будет обсуждаться в данной работе, поскольку изменение фокусного расстояния изменяет не саму радиальную дисторсию, а только ее параметризацию.
С ростом номера члена ряда / коэффициенты разложения к у реальных оптических систем быстро спадают, поэтому можно рассматривать только несколько первых членов разложения. Мы ограничим себя рассмотрением параметров радиальной дистор-сии до п = 3.
1.2. Компенсация дисторсии при известных параметрах
Пусть дано изображение, искаженное радиальной дисторсией с известными параметрами. Рассмотрим задачу формирования изображения с компенсированной радиальной дисторсией. Для каждого пикселя с координатами (х,у) формируемого изображения I осуществим преобразование координат по формуле (1). Значение яркости в пикселе с координатами (ха, уа) на искаженном изображении 1а будет искомым значением для пикселя на формируемом изображении. При этом, поскольку изображение 1а задано только в узлах целочисленной сетки координат, для нецелых значений ха и уа целочисленные значения должны быть получены путем интерполяции функции 1а.
Обратим внимание, что для обратного к радиальной дисторсии преобразования изображения оказалось достаточным использование только прямого преобразования координат (1) без необходимости его обращения.
1.3. Метод оценки параметров радиальной дисторсии
Как было сказано выше, алгоритм восстановления будет опираться на оценку качества исправления дис-торсии при всевозможных значениях ее параметров на некоторой сетке. Опишем теперь собственно метод оценки качества при известных пробных параметрах дисторсии. Для этого нам понадобится изложить некоторые детали использования аппарата БПХ.
1.3.1. Идея метода
Преобразование Хафа заключается в суммировании значений пикселей изображения вдоль всевозможных пересекающих его прямых. Вычисление Хаф-образа неоднократно предлагалось для оценки прямолинейности линий на изображении [10, 11], однако не нашло развития из-за кажущейся вычислительной сложности наиболее распространенного алгоритма [18, 19]. В настоящей работе применяется этот же подход, но предлагается использовать быстрое преобразование Хафа.
Пусть исходное изображение имеет размер п*ш. При быстром преобразовании Хафа преимущественно вертикальные (с тангенсом угла отклонения от вертикали, ограниченным единицей) прямые исходного изображения параметризованы точкой пересечения с верхним краем изображения 5 и точкой пересечения с нижним краем изображения яИ. Величину 1, характеризующую угол наклона прямой, будем называть тангентой, 1е [-п, п]. В Хаф-пространстве каждой преимущественно вертикальной прямой соответствует точка с координатами (я, 1). Определим матрицу И" размером ш*(2п-1), значение каждого элемента
которой равно сумме значений пикселей исходного изображения вдоль соответствующей прямой. Для преимущественно горизонтальных прямых используется аналогичная параметризация с той лишь разницей, что используются левый и правый края изображения. Полученная матрица сумм Ик°Г имеет размер и*(2от-1). Пара матриц (Иуег, ИЬог) образует полный Хаф-образ И.
Предположим, что на исходном изображении имеется (возможно, прерывистый) яркий длинный отрезок, лежащий на прямой, параметризованной как ^о, /0), и что в его окрестности нет других столь же ярких и длинных отрезков. В этом случае значение суммы яркостей вдоль прямой /0) будет велико, а суммы вдоль окрестных прямых - малы. Таким образом, на Хаф-образе будет наблюдаться резкий пик в точке s0 на прямой / = /0. Предположим теперь, что на исходном изображении имеется яркий криволинейный отрезок. Так как существует много прямых с близкими параметрами, касательных к нему, то суммы вдоль них дадут близкие значения и на Хаф-образе будет наблюдаться пологий «холм».
Поскольку на изображениях реальных сцен, как правило, присутствуют прямые линии, а шанс, что при восстановлении радиальной дисторсии с произвольными параметрами какая-то кривая случайно восстановится до прямой, невелик, то лучшими параметрами радиальной дисторсии на изображении можно считать такие, при которых на восстановленном изображении выше прямолинейность линий, то есть его Хаф-образ содержит резкие пики, а не размытые «холмы».
Опишем теперь конкретный алгоритм вычисления функционала оценки качества исправления дисторсии, основанный на быстром преобразовании Хафа.
1.3.2. Исправление дисторсии с сохранением интеграла интенсивности
Пусть имеется изображение сцены, содержащей прямые линии, подверженное радиальной дисторсии, и значения параметров дисторсии. Поставим задачу оценки, искажено ли изображение радиальной дис-торсией именно с этими параметрами. Для этого выделим на изображении границы путем вычисления модуля градиента яркости в каждой точке, а затем выполним пробное восстановление изображения, как описано в подпараграфе 1.2. Теперь выполним быстрое преобразование Хафа восстановленного изображения. Как было указано выше, в случае, если данные параметры радиальной дисторсии - правильные, то яркие длинные отрезки прямых на восстановленном изображении дадут резкие пики на Хаф-образе, причем чем ярче и длиннее отрезки прямых, тем заметнее будет этот пик на фоне пологих «холмов», соответствующих криволинейным границам. Однако обратим внимание на то, что преобразование радиальной дисторсии не сохраняет длины и площади, поэтому интеграл интенсивности отрезка при преобразовании может произвольно измениться в зависимости от положения на изображении. Более того, если исправле-
ние радиальной дисторсии превращает криволинейный участок в прямолинейный при сохранении расстояния между концами отрезка, то его длина уменьшается, что приводит к уменьшению интегральной интенсивности и падению амплитуды соответствующего ему пика в Хаф-пространстве, а также приводит к опосредованному «штрафу» правильных параметров. Таким образом, для увеличения точности оценки прямолинейности линий по соответствующим им пикам в Хаф-пространстве актуальным является разработка такого способа реконструкции радиально искаженного изображения, которое при восстановлении геометрии одновременно сохраняет интегральную интенсивность границ на изображении.
Пусть имеется искаженное изображение /л. Создадим пустое (все значения которого равны нулю) изображение /'. Пусть имеется обратное к (1) преобразование координат, позволяющее для каждого пикселя (Хс1,ул) искаженного изображения найти координаты (х,у) на изображении /'. Тогда для каждого пикселя искаженного изображения /л выполним операцию
IХх у) = I'(х у) + /л(ха, уа).
Предложенное преобразование изображений удовлетворяет сформулированному требованию сохранения интегральной интенсивности вдоль кривых, хотя и непригодно для визуализации, поскольку не сохраняет значение в точке. Поэтому в данной работе оно используется для пробного восстановления изображений, а стандартное, описанное в подпарагра-фе 1.2, используется для построения выходного изображения по найденным параметрам дисторсии.
Как было указано выше, пробное восстановление, сохраняющее интегральную интенсивность вдоль кривых, требует вычисления преобразования координат, обратного к (1). Однако в случае немонотонности (1) обратного преобразования не существует. Поэтому разумно наложить ограничения на возможные значения параметров дисторсии Ь таким образом, чтобы преобразование (1) было монотонным на участке [0, г*], где г* - такой радиус, что / (г) оказывается не меньше некоторого заранее заданного гс„/. Окружность с этим радиусом будем называть критической.
Если в качестве критической брать окружность с радиусом Тех:, описанную вокруг изображения, то ограничение на параметры дисторсии оказывается слишком жестким в том смысле, что изображения с сильной дисторсией из тестового набора не описываются с разумной точностью полиномиальной моделью с малым числом членов. Если же радиус критической окружности слишком мал (например, равен радиусу гм вписанной в изображение окружности), то найденная оценка параметров исправления дисторсии оказывается малоосмысленной, поскольку на значимой периферийной области изображения искажения могут оказаться сколь угодно велики. В нашей работе радиус критической окружности был выбран равным Тоги = 0,7-Тех/. При этом отношение сторон у изображений составляло 3:4, что дает гш = 0,6-Тех/. Для других
соотношений сторон разумным выбором для радиуса критической окружности, по-видимому, будет значение, делящее отрезок [гы, гх] в близких пропорциях.
Обратим внимание на то, что параметры дистор-сионного искажения влияют на глобальный масштаб изображения, поэтому при восстановлении изображения возможен выход части сцены за его границы, что нарушит интегральную интенсивность по кадру. Введем понятие неподвижной окружности. Назовем так невырожденную окружность с радиусом, не изменяющимся при пробном исправлении радиальной дисторсии. При некоторых параметрах дисторсии такой окружности может не существовать. Теперь при неизменных параметрах дисторсии введем дополнительное масштабирование с коэффициентом к0. Для любой окружности не больше критической для данных параметров дисторсии можно подобрать ко так, что выбранная окружность будет неподвижной (для преобразования с масштабированием). Далее будем выбирать к0 так, чтобы неподвижной была критическая окружность.
Заметим теперь, что монотонность преобразования /(кот) вплоть до Гати не гарантирует сохранность интегральной интенсивности на изображении внутри критической окружности, если последняя выходит за пределы изображения. Для сохранения интеграла внутри усеченной окружности дополнительно потребуем /(ког)<ког на отрезке [ты, гаги]. Теперь для гарантии сохранения интегральной интенсивности на всём изображении достаточно обнулить значения модуля градиента вне критической окружности.
1.3.3. Вычисление обратного преобразования координат при радиальной дисторсии
Предложенная процедура для пробного восстановления радиально искаженного изображения требует вычисления преобразования координат, обратного к (1). Пусть заданы координаты пикселя на искаженном изображении (ха, у) и параметры радиальной дисторсии.
Положим га = х2 + у^ и перепишем уравнение (1) в следующем виде:
rd = f '(r ) = (1 + Jk'f 2 )r
(2)
Найдём масштабирующее преобразование с коэффициентом ко, сохраняющее радиус критической окружности:
k0f (rcrit) rcrit.
(3)
C учетом масштабирования с коэффициентом k0 и замены k,- = koki для i > 0 преобразование координат приобретает вид:
rd = f (r) = (Ykr2i )r .
(4)
Уравнение (4) не имеет аналитического решения относительно г при п > 1, поэтому будем искать решение численно. А поскольку обратное преобразование с фиксированными параметрами дисторсии требуется вычислять многократно (трудоемкость пропорциональна
площади изображения), разумно единожды построить приближенную модель обратного преобразования и пользоваться ей вместо численного решения уравнения (2) в каждой точке. В связи с тем, что преобразование радиальной дисторсии для рассматриваемых объективов оказалось достаточно гладким, а точки относительно центра радиального искажения смещены симметрично, значение функции г =/-1(га) численно оценивалось для т равноотстоящих друг от друга значений га на участке [0, га„,] и использовалась линейная интерполяция для вычисления г для всех точек на входном изображении. В проведенных ниже численных экспериментах т бралось равным 300.
1.3.4. Угловой дескриптор изображения Как уже было сказано выше, прямые линии на изображении порождают на строках Хаф-образа резкие пики, а кривые - более размытые. Чтобы избавиться от размытых пиков и сохранить более резкие, вычтем из Хаф-образа его же, но сглаженный Гауссовым фильтром вдоль оси t, оставив только неотрицательные значения:
Р"ег = тах(0, И™ - О, (о) * И™), Рш = тах(0, Иш - О, (о)* Ик"г).
(5)
где О,(о)- ядро Гауссова фильтра со стандартным отклонением о. В приведенных ниже экспериментах о бралось равным 5,0 при ширине изображения 360 пикселей.
Теперь определим векторы ^ уег длиной 2п-1 и Е1юг длиной 2т-1 так, чтобы в каждом их элементе содержалась оценка прямолинейности для множества краев определенной тангенты ,. В качестве такой оценки возьмем дисперсию значений в строке , отфильтрованных Хаф-образов Р™г и Р["г.
Выбор дисперсии в качестве оценки прямолинейности объясняется тем, что размытые «холмы» обладают меньшей дисперсией значений, чем острые пики. Также дисперсия растет при увеличении числа пиков.
Учтем теперь тот факт, что прямые, проходящие через центр изображения, не искривляются при любой радиальной дисторсии вовсе, а проходящие близко к центру искривляются меньше, чем далекие. Поскольку за счет дискретизации и других возможных искажений амплитуда пиков на Хаф-образе оказывается зашумленной, разумно приписать меньший удельный вес пикам, соответствующим прямым с меньшим расстоянием до центра изображения.
С учетом сказанного оценка прямолинейности для пучка параллельных прямых каждой тангенты , выглядит следующим образом:
РГ = о2( Р^ ),
yver \ ' s,t
Fhor=s2(Ph°rwh°r),
Yflws,t (Ypw
s2(P,W ) = -
Ж
ж,
(6)
(7)
где №'ег и Wh"г - расстояние от прямой (. ,) до оптического центра изображения. Такую редукцию рас-
i=1
,=0
2
s
s
s
сматриваемых Хаф-образов И"ег и Иког к вектору р = р'ог^рког назовем угловым дескриптором.
1.3.5. Оценка исправления дисторсии по угловому дескриптору
На угловом дескрипторе мы ожидаем увидеть пики от пучков квазипараллельных краев, причем от их прямолинейности зависит как амплитуда пиков, так и скорость падения. В качестве меры выраженности пиков, учитывающей оба эти эффекта, воспользуемся величиной энтропии сигнала.
Пусть р имеет к уровней значений, а частота появления уровня к равна Р(Рк). Тогда энтропия сигнала р определяется следующим образом:
К
Е(Р ) = -£Р( Рк )1о§ Р(Рк). (8)
к=1
Это выражение использовалось в работе как оценка качества исправления дисторсии - чем меньше энтропия Е, тем сильнее выражены пики в угловом дескрипторе, что, как ожидается, связано с тем, что прямые линии на изображении менее искажены.
1.4. Общее описание алгоритма слепой компенсации Представленный алгоритм использует схему исчерпывающего перебора в ш-мерном пространстве параметров радиальной дисторсии (эксперименты проводились при п = 3) в заданных границах. Он состоит из следующих шагов:
1. Сформировать контурное изображение, выделив
границы на входном изображении путем вычисления модуля градиента яркости в каждой точке изображения.
2. Обнулить пиксели контурного изображения вне
критической окружности (см. п. 1.3.2).
3. Выполнить перебор в пространстве параметров в за-
данных границах. Для каждого набора параметров:
3.1. Проверить налагаемые в п. 1.3.2 ограничения на дисторсионное преобразование при данном наборе параметров; в случае, если ограничения не выполняются, перейти к следующему набору параметров.
3.2. Вычислить обратное преобразование координат (п. 1.3.3).
3.3. Выполнить обратное дисторсионное преобразование контурного изображения с сохранением интеграла интенсивности (п. 1.3.2).
3.4. Сформировать угловой дескриптор (п. 1.3.4).
3.5. Оценить исправление дисторсии по угловому дескриптору (п. 1.3.5).
4. Найти набор параметров, минимизирующий значе-
ние целевой функции (8).
5. Получить выходное изображение путем компенса-
ции радиальной дисторсии входного изображения, используя найденные параметры (п. 1.2).
2. Экспериментальная проверка качества работы предложенного алгоритма
Тестирование алгоритма проводилось на изображениях документов и трехмерных сцен внутри и вне помещений, содержащих прямые линии (рис. 2а—г).
Изображения были сняты камерами с 6 различными объективами (табл. 1). Также для контроля алгоритм применялся к изображениям, приведенным в [12]. Для каждого набора данных, исключая взятые из [12], с помощью программного модуля MATLAB Single Camera Calibration, использующего алгоритм Жанга [1], были получены референтные параметры радиальной дисторсии. В качестве калибровочного объекта использовался жесткий транспарант со стороной 1 метр с шахматным паттерном шагом 10 см.
Цель эксперимента состояла в количественной оценке точности определения параметров дисторсии для различных модификаций предложенного алгоритма на множестве тестовых изображений. Следует заметить, что простые метрики в пространстве параметров модели дисторсии - не лучший выбор для количественной оценки. Во-первых, нет очевидного способа сопоставить цену равных отклонений для разных коэффициентов модели. Во-вторых, что более важно, два полинома с существенно разными наборами коэффициентов могут мало различаться на интересующем нас отрезке. С практической точки зрения цена ошибки алгоритма компенсации любой дис-торсии напрямую связана с отклонениями координат точек после исправления относительно идеальных координат в области изображения. Далее будет описан количественный способ оценки выраженности остаточной дисторсии, а также приведены результаты работы алгоритма согласно этой методике.
2.1. Способ оценки выраженности дисторсии и точности ее компенсации
Пусть нам известно исходное положение точки r = (х, y) до искажения и наблюдаемое положение точки rd = (xd, yd) в результате действия некоторого искажения. При этом координаты восстановленного положения точки равны (х', y ) = f l(rd), где f-1 -оцениваемая функция компенсации дисторсии. Будем рассматривать евклидову норму отклонения координат при восстановлении в координатах идеального изображения:
Df (х,y) = || г - f -1(rd)||2 (9)
Введем на изображении регулярную прямоугольную сетку N*M (рис. 1а) и определим оценку качества компенсации дисторсии как среднюю норму отклонения в узлах сетки:
1 N,M
d f =-УД f (х, y,). (10)
f NM ~ix f '
Как было уже упомянуто выше, для искаженных изображений понятие масштаба плохо определимо. Кроме того, образ, восстановленный из искаженного прямоугольного изображения (рис. 1б), не будет прямоугольным (рис. 1г). На практике это приводит к необходимости выбора масштаба (и, в общем случае, сдвига) прямоугольника выходного изображения, оптимального с точки зрения конкретной прикладной
задачи. С этой точки зрения общее растяжение или сжатие изображения при восстановлении не должно влиять на оценку качества восстановления. Поэтому модифицируем оценку качества следующим образом:
А/,р(х,у)=|| г -Ьр(/-\га))||2
• 1 , \ , (11) ¿/=тп ш5А / ■ р(х,у)
где Ьр - преобразование масштабирования и сдвига (в нашем случае - только масштабирования) с параметрами р.
Аналогично определим степень выраженности исходной дисторсии изображения ¿0 как оценку качества восстановления в отсутствие всякой компенсации (рис. 1 в):
Л0 = ¿/ |/(Г)=г . (12)
а)
б)
в)
г)
Рис. 1. Оценка алгоритма компенсации радиальной дисторсии: а) используемая сетка, б) результат действия радиальной дисторсии, в) преобразование с такими же
параметрами дисторсии, но с оптимальным масштабированием, используемым для вычисления ¿0,
г) компенсация дисторсии с оптимальным масштабированием, используемым для вычисления ¿/
Теперь мы можем ввести относительную оценку качества компенсации дисторсии:
= 10 -[1 - (¿/ /(¿0 +£»], (13)
где 10 - коэффициент, используемый для приведения оценки к 10-балльной шкале, а £ - параметр, регулирующий толерантность оценки восстановления малых дисторсий. Его значение может меняться в зависимости от требований к точности восстановления. Введенная таким образом относительная оценка позволяет получать сравнимые результаты на наборах изображений с существенно различными уровнями искажений.
В данной работе при вычислении (13) использовалась равномерная сетка 36*48 узлов и £ = 1,0 при размере изображения 360*480. 2.2. Обсуждение экспериментальных результатов
Проведена серия экспериментов, в которых сравнивались результаты работы разных вариантов предложенного алгоритма. В качестве контрольных использовались воспроизведенные результаты работы алгоритма, представленного в [12]. Каждая последующая модификация алгоритма кумулятивно применя-
лась к предыдущей версии. В первом эксперименте предложенный в [12] алгоритм был модернизирован путем замены используемой в нем модели радиальной дисторсии на двупараметрическую, а метода обращения радиальной дисторсии при пробном восстановлении на описанный в п. 1.3.2 за исключением масштабирования, сохраняющего радиус критической окружности. Такое масштабирование было добавлено во втором эксперименте. В третьем эксперименте используемый в [12] функционал оценки исправления дисторсии был заменен на предложенный в п. 1.3.5 за исключением приписывания различных весов пикам на Хаф-образе в зависимости от расстояния прямой до центра изображения (п. 1.3.4). Такое приписывание весов было добавлено в четвертом эксперименте. В пятом эксперименте модель радиальной дисторсии была заменена на трехпараметрическую.
Характеристики оптических систем, использованных для получения тестовых данных, представлены в табл. 1. Результаты экспериментов представлены в табл. 2 и на рис. 2д-з.
Табл. 1. Характеристики используемого набора данных
Номер объектива Количество снимков Выраженность дисторсии ¿0
1 10 1,25
2 10 7,12
3 10 3,61
4 10 6,99
5 10 5,15
6 15 6,47
Заключение
В настоящей работе предложен алгоритм восстановления изображения, подверженного радиальной дисторсии, на основе анализа только самого этого изображения без требований априорной информации о параметрах камеры. Алгоритм был проверен на реальных изображениях, полученных с объективов с разной выраженностью радиальной дисторсии, и показал практическую применимость. Точность исправления радиальной дисторсии значительно превысила точность предшественника, предложенного в [12], в том числе за счет использования модели дисторсии с большим числом параметров.
Предложенный алгоритм основан на минимизации функционала оценки прямолинейности краев, построенного над результатом быстрого преобразования Хафа изображения пробного исправления дисторсии для модуля градиента входного изображения. Как было продемонстрировано в [12], функционалы такого вида могут оказаться невыпуклыми, поэтому для решения оптимизационной задачи применен метод полного перебора в пространстве параметров радиальной дисторсии. В связи с тем, что исправление радиальной дисторсии по единичному изображению из неизвестного источника не представляется задачей, решаемой для больших потоков данных, к нему не предъявляется требование высокого быстродействия, поэтому полный перебор вполне обоснован.
Табл. 2. Результаты компенсации радиальной дисторсии различными модификациями алгоритма
(по 10-балльной шкале, см. п. 3.1)
Среднее качество Среднее
№ Версия алгоритма по используемому набору данных качество
1 2 3 4 5 6 по версии
0 Контроль: алгоритм, предложенный в [12] 8,4 2,6 8,6 0,8 4,6 2,0 4,50
1 Способ пробного восстановления заменен на описанный в п. 1.3.2 без масштабирования; модель дисторсии расширена до двухпараметрической 8,0 5,4 7,9 2,3 1,6 1,9 4,52
2 Добавлено масштабирование 3,0 6,0 3,0 6,3 2,7 5,75 4,46
3 Функционал оценки исправления дисторсии заменен на предложенный в п. 1.3.5 без взвешивания пиков на Хаф-образе 6,3 6,2 7,8 6,5 7,1 4,0 6,32
4 Добавлено взвешивание 9,0 7,8 9,0 7,8 7,9 7,75 8,21
5 Модель дисторсии расширена до трехпараметрической 9,4 8,2 9,6 8,5 8,5 6,5 8,45
ш
Рис. 2. Пример работы алгоритма: а—г) входные изображения, д-з) изображения с исправленной радиальной дисторсией,
полученные финальной версией разрабатываемого алгоритма
Также в работе подробно рассмотрены основные этапы алгоритма и приведены оценки качества работы различных версий алгоритма согласно предложенному критерию. Для наглядности были приведены количественные и визуальные данные, показывающие, что предложенный метод успешно справляется как с вариативностью сцен, так и с различными условиями съемки и степенью искажений.
В дальнейшем планируется ускорить процесс поиска параметров радиального искажения и адаптировать предложенный алгоритм для работы с видеофрагментами, используя информацию с отличных друг от друга кадров для одновременной оптимизации параметров радиального искажения.
Благодарности
Исследование выполнено за счет гранта Российского научного фонда (проект № 14-50-00150).
Литература
1. Zhang, Z. A flexible new technique for camera calibration / Z. Zhang // IEEE Transactions on Pattern Analysis and Machine Intelligence. - 2000. - Vol. 22(11). - P. 1330-1334. -DOI: 10.1109/34.888718.
2. Миронова, Т.В. Анализ деформаций, оптических неод-нородностей и дисторсионных искажений с помощью искусственных спеклов в цифровой фотографии: дисс. к.ф.-м.н., Физический институт им. П.Н. Лебедева РАН, Москва, 2005.
3. Hartley, R. Self-calibration of stationary cameras / R. Hartley // International Journal of Computer Vision. -1997. - Vol. 22(1). - P. 5-23. - DOI: 10.1023/A:1007957826135.
4. Stein, G. Accurate internal camera calibration using rotation, with analysis of sources of error / G. Stein // Computer Vision: Proceedings, Fifth International Conference on. - 1995. - P. 230-236. - DOI: 10.1109/ICCV.1995.466781.
5. Farid, H. Blind removal of lens distortion / H. Farid, A.C. Popescu // Journal of the Optical Society of America A. - 2001. - Vol. 18(9). - P. 2072-2078. - DOI: 10.1364/J0SAA.18.002072.
6. Zhang, Z. Camera calibration with lens distortion from low-rank textures / Z. Zhang, Y. Matsushita, Y. Ma // Computer Vision and Pattern Recognition (CVPR), 2011 IEEE Conference on. - 2011. - P. 2321-2328. - DOI: 10.1109/CVPR.2011.5995548.
7. Karpenko, S. UAV Control on the Basis of 3D Landmark Bearing-Only Observations / S. Karpenko, I. Konovalenko, A. Miller, B. Miller, D. Nikolaev // Sensors. - 2015. -Vol. 15(12). - P. 29802-29820. - DOI: 10.3390/s151229768.
8. Karpenko, S. Visual navigation of the UAVs on the basis of 3D natural landmarks / S. Karpenko, I. Konovalenko, A. Miller, B. Miller, D. Nikolaev // Proceedings of SPIE: Eighth International Conference on Machine Vision (ICMV 2015). - 2015. - Vol. 9875. - 98751I (10 p.).
9. Konovalenko, I. UAV navigation on the basis of the feature points detection on underlying surface / I. Konovalenko, A. Miller, B. Miller, D. Nikolaev // Proceedings of the 29th European Conference on Modeling and Simulation (ECMS 2015), Albena (Varna), Bulgaria. - 2015. - P. 26-29. - DOI: 10.7148/2015-0499.
10. Rosten, E. Camera distortion self-calibration using the plumb-line constraint and minimal Hough entropy / E. Rosten, R. Loveland // Machine Vision and Applications. - 2011. - Vol. 22(1). - P. 77-85. - DOI: 10.1007/s00138-009-0196-9.
11. Alemán-Flores, M. Automatic Lens Distortion Correction Using One-Parameter Division Models / M. Alemán-Flores, L. Alvarez, L. Gomez, D. Santana-Cedres // Image Processing On Line. - 2014. - Vol. 4. - P. 327-343. - DOI: 10.1007/s10851-012-0342-2.
12. Карпенко, С.М. Метод восстановления изображений, подверженных радиальной дисторсии / С.М. Карпенко, С.А. Гладилин, Д.П. Николаев // Сборник трудов конференции «Информационные технологии и системы (ИТиС'08)». - 2008. - С. 502-505.
13. Kanuki, Y. Automatic compensation of radial distortion by minimizing entropy of histogram of oriented gradients / Y. Kanuki, N. Ohta, A. Nagai // Pattern Recognition (ACPR), 2013 2nd IAPR Asian Conference on. - 2013. -P. 912-916. - DOI: 10.1109/ACPR.2013.167.
14. Wang, A. A simple method of radial distortion correction with centre of distortion estimation / A. Wang, T. Qiu,
L. Shao // Journal of Mathematical Imaging and Vision. -2009. - Vol. 35(3). - P. 165-172. - DOI: 10.1007/s10851-009-0162-1.
15. Nikolaev, D. Hough transform: underestimated tool in the computer vision field / D. Nikolaev, S. Karpenko, I. Nikolaev, P. Nikolayev // Proceedings of the 22th European Conference on Modelling and Simulation. - 2008. - P. 238-246.
16. Brady, M.L. A fast discrete approximation algorithm for the Radon transform / M.L. Brady // SIAM Journal on Computing. - 1998. - Vol. 27(1). - P. 107-119. - DOI: 10.1137/S0097539793256673.
17. Brown, D.C. Close-range camera calibration / D.C. Brown // Photogrammetric Engineering. - 1971. - Vol. 37(8). - P. 855866.
18. Duda, R.O. Use of the Hough transformation to detect lines and curves in pictures / R.O. Duda, P.E. Hart // Communications of the ACM. - 1972. - Vol. 15(1). - P. 11-15. -DOI: 10.1145/361237.361242.
19. U.S. Patent 3,069,654 G06K 9/46, G01T 5/02, G01T 5/00, 382/281, 382/202, 342/176, 342/190. Method and means for recognizing complex patterns / P. Hough, filed of March 25, 1960, published of December 18, 1962.
Сведения об авторах
Кунина Ирина Андреевна, младший научный сотрудник ИППИ РАН. Окончила в 2014 г. НИТУ МИСиС. Область научных интересов: обработка изображений, компьютерное зрение, калибровка зрительных систем. Email: kunina@iitp.ru .
Гладилин Сергей Александрович, к. ф.-м. н., научный сотрудник ИППИ РАН. Окончил в 2002 г. МГУ им. М.В. Ломоносова. Область научных интересов: зрительные системы, обработка изображений, распознавание образов. E-mail: gladilin@iitp.ru .
Николаев Дмитрий Петрович, к. ф.-м. н., заведующий лабораторией в ИППИ РАН. Окончил в 2000 г. МГУ им. М.В. Ломоносова. Область научных интересов: машинное зрение, быстрые алгоритмы обработки изображений, распознавание образов. E-mail: dimonstr@iitp.ru .
ГРНТИ: 28.23.15.
Поступила в редакцию 10 мая 2016 г. Окончательный вариант - 4 июня 2016 г.
BLIND RADIAL DISTORTION COMPENSATION IN A SINGLE IMAGE USING FAST HOUGH TRANSFORM
I. A. Kunina, S. A. Gladilin, D. P. Nikolaev Institute for Information Transmission Problems RAS, Moscow, Russia
Abstract
In this paper, we present an automatic technique for compensation of radial distortion, which is characteristic of wide-angle lenses. The proposed method estimates distortion parameters using a single image from unknown source. No calibration objects are required, but it is assumed that the original scene contains straight lines. The method is based on finding such radial distortion parameters, that maximize total length of linear segments. We employ a fast Hough transform to estimate the overall curvature of lines without selecting any. The proposed algorithm is tested on real images obtained using calibrated camera lenses with different radial distortion.
For the formal evaluation of the algorithm, we propose a quality measure for geometric distortion compensation, which works correctly even in the case when the problem of determining the coefficients is ill-conditioned.
Keywords: digital image processing, image analysis, lens system design, radial distortion, automatic calibration, fast Hough transform.
Citation: Kunina IA, Gladilin SA, Nikolaev DP. Blind compensation of radial distortion in a single image using fast Hough transform. Computer Optics 2016; 40(3): 395-403. - DOI: 10.18287/2412-6179-2016-40-3-395-403.
Acknowledgements: This research was conducted at the IITP RAS and supported by Russian Science Foundation Grant #14-50-00150.
References
[1] Zhang Z. A flexible new technique for camera calibration. Pattern Analysis and Machine Intelligence, IEEE Transactions on 2000; 22(11): 1330-1334. DOI: 10.1109/34.888718.
[2] Mironova TV. Analysis of deformations, optical inhomogeneities and distortion in digital photography using artificial speckles [In Russian]: dis. ... cand. sci. (phys.-math.). Moscow: Lebedev Physical Institute; 2005.
[3] Hartley RI. Self-calibration of stationary cameras. International Journal of Computer Vision 1997; 22(1): 5-23. DOI: 10.1023/A:1007957826135.
[4] Stein GP. Accurate internal camera calibration using rotation, with analysis of sources of error. Computer Vision, Proceedings, Fifth International Conference on; 1995. DOI: 10.1109/ICCV.1995.466781.
[5] Farid H, Popescu AC. Blind removal of lens distortion. JOSA A 2001; 18(9): 2072-2078. DOI: 10.1364/J0SAA.18.002072.
[6] Zhang Z, Matsushita Y, Ma Y. Camera calibration with lens distortion from low-rank textures. Computer Vision and Pattern Recognition (CVPR), 2011 IEEE Conference on; 2011. DOI: 10.1109/CVPR.2011.5995548.
[7] Karpenko S, Konovalenko I, Miller A, Miller B, Nikolaev D. UAV Control on the Basis of 3D Landmark Bearing-Only Observations. Sensors 2015; 15(12): 29802-29820. DOI: 10.3390/s151229768.
[8] Karpenko S, Konovalenko I, Miller A, Miller B, Nikolaev D. Visual navigation of the UAVs on the basis of 3D natural landmarks. Eighth International Conference on Machine Vision, International Society for Optics and Photonics 2015; 9875: 9875I. DOI: 10.1117/12.2228793.
[9] Konovalenko I, Miller A, Miller B, Nikolaev D. UAV navigation on the basis of the feature points detection on underlying surface. Proceedings of the 29th European Conference on Modeling and Simulation (ECMS 2015), Albena (Varna), Bulgaria, 2015: 26-29. DOI: 10.7148/2015-0499.
[10] Rosten E, Loveland R. Camera distortion self-calibration using the plumb-line constraint and minimal Hough entropy. Machine Vision and Applications 2011; 22(1): 77-85. DOI: 10.1007/s00138-009-0196-9.
[11] Alemán-Flores, M, Alvarez L, Gomez L, Santana-Cedres D. Automatic Lens Distortion Correction Using One-Parameter Division Models. Image Processing On Line 2014; 4: 327-343. DOI: 10.1007/s10851-012-0342-2.
[12] Karpenko, SM, Gladilin DP, Nikolaev DP. Method for recovery of radial distorted images [In Russian]. Proceedings of the conference Information Technologies and Systems ITaS, 2008: 502-505.
[13] Kanuki Y, Ohta N, Nagai A. Automatic compensation of radial distortion by minimizing entropy of histogram of oriented gradients. Pattern Recognition (ACPR), 2013 2nd IAPR Asian Conference on, 2013: 912-916. DOI: 10.1109/ACPR.2013.167.
[14] Wang A, Qiu T, Shao L. A simple method of radial distortion correction with centre of distortion estimation. Journal of Mathematical Imaging and Vision 2009; 35(3): 165-172. DOI: 10.1007/s10851-009-0162-1.
[15] Nikolaev D, Karpenko S, Nikolaev I, Nikolayev P. Hough transform: underestimated tool in the computer vision field. Proceedings of the 22th European Conference on Modelling and Simulation 2008; 238-246.
[16] Brady ML. A fast discrete approximation algorithm for the Radon transform. SIAM Journal on Computing 1998; 27(1): 107119. DOI: 10.1137/S0097539793256673.
[17] Brown DC. Close-range camera calibration. Photogrammetric engineering 1971; 37(8): 855-866.
[18] Duda RO, Hart PE. Use of the Hough transformation to detect lines and curves in pictures. Communications of the ACM 1972; 15(1): 11-15. DOI: 10.1145/361237.361242.
[19] Hough PVC. Method and means for recognizing complex patterns. US Patent 3069654 of Dec 18, 1962.
Authors' information
Irina Andreevna Kunina, a research assistant at the IITP RAS. Graduated from NUST MISiS in 2014. Research interests are image processing, computer vision, calibration of visual systems. E-mail: kunina@iitp.ru .
Sergey Alexandrovich Gladilin, Ph. D. in Physics and Mathematics, a researcher at the IITP RAS. Graduated from Lomonosov MSU in 2002. Research interests are visual systems, image processing and pattern recognition. E-mail: gladilin@iitp. ru .
Dmitry Petrovich Nikolaev, Ph. D. in Physics and Mathematics, a head of the laboratory at the IITP RAS. Graduated from Lomonosov M.V. in 2000. Research interests are machine vision, algorithms for fast image processing, pattern recognition. E-mail: dimonstr@iitp.ru .
Received May 10, 2016. The final version - June 4, 2016.
