Коррекция геометрических искажений

Идеальные фигуры встречаются только в геометрии.
(Кто-то пошутил)

Не сильно ошибусь, если предположу, что у подавляющего большинства человечества, использующего фото-редакторы, самой востребованной фичей является функция коррекции линейных искажений. Хотя некоторые идут еще дальше: Фотосессия длиной в 4 года. А некоторые, совсем наоборот, утверждают, что желание исправлять геометрию рисунка — это в чистом виде технодрочерство и никакого отношения к фотографии не имеет. Мне же представляется, что все инструменты известные в природе на сегодняшний день для решения этой задачи, не совсем удобные и иногда избыточные с точки зрения возможностей (хотя, все не так, да), если говорить не о технологически точных фото-задачах, а о простом повседневном использовании: ну, там поехал в отпуск, фоткнул красивые руины, дома посмотрел, а там все криво — открыл свою любимую программку (которую сегодня напишем), поправил — мир стал лучше.

Поэтому сегодня пост, о том как самостоятельно изготовить быстрый, легкий и в меру незагруженный тулкит. Будем решать частный, но самый распространенный случай коррекций: исправление искажений проекции: ровно тот случай, что описал Антон, когда плоскость проекции изображения фотографии в момент фиксации не соответствует «правильной» (плоскости объекта съемки и фокальная плоскость не параллельны).

В результате получим полноценный движок, с помощью которого, за считанные секунды, будем править фотографии любой геометрической степени запущенности:

Приложение можно сразу же скачать, текст дальше не читать, а писать багрепорты, например.

Трансформации бывают не только аффинными

В предыдущем посте мы решали задачи кадрирования, поворотов и масштаба проекции изображения на фото-пластину в терминах аффинных преобразований. А так же определили специальный класс IMPQuad, по сути являющийся контейнером произвольных вершин четырех-угольника (пока всегда подразумеваем выпуклый) и солвером для некоторых операций над ними:

public struct IMPQuad {
/// Left bottom point of the quad
public var left_bottom = float2( -1, -1)

/// Left top point of the quad
public var left_top = float2( -1, 1)

/// Right bottom point of the quad
public var right_bottom = float2( 1, -1)

/// Right top point of the quad
public var right_top = float2( 1, 1)
}

Тут есть один странный момент: не смотря на достаточную востребованность на решатели задач преобразования фигур на плоскости (расчет произвольных матриц трансформации и т.п.), ни iOS ни OSX не имеют по умолчанию встроенной библиотеки работы с неаффинными методами. Есть Quartz 2D Reference Collection, в полном объеме реализующий набор линейных трансформаций плоских фигур с сохранением смежности, но в случае задачи коррекции линейных искажений, нам не хватает чего-то подобного для расчета матрицы «деформаций» без сохранения смежности, т.е. когда параллельные прямые трансформируются в пересекающиеся и наоборот. Если говорить о линейной коррекции искажений в плоскости, то задача всегда сводится к проекции аффинных преобразований в 3D на плоскость. Но, во-первых, манипулировать 3D объектом в ручную не всегда удобно, во-вторых не все, что нам нужно деформировать или скорректировать является результатом искажений проекции. Поэтому, обычно коррекцию искажений проекции, как частного случая коррекции геометрических линейных искажений, сводят к задаче деформации фигуры на плоскости. Именно поэтому мы будем работать не с представлением прямоугольника, а произвольным (но выпуклым) четырёх-угольником.

Конечно же, готовых, но сторонних инструментов для решения вот этого всего есть предостаточно: начиная от CGAL и заканчивая openCV, но тащить в один небольшой проект эти могучие либы, мне представляется расточительным. Тем более, что все решение займет не более нескольких строчек кода на базе все тех же классов из акселератора SIMD, который с прошлого года доступен на всех платформах Apple.

Варианты решения

Напомню, что мы работаем с достаточно простым объектом: фото-пластиной, которая в общем виде описывается 4-мя вершинами. В сегодняшнем упражнении нам нужно 4-ре исходных координаты преобразовать в 4-ре координаты целевого объекта. Координаты эти, надо полагать нам известны: получаем мы их, как не трудно предположить, в результате манипуляций пользователя с приложением, к примеру считываем координаты курсора мыши или тачпада (хотя никто не мешает такие данные генерить аналитически и запускать цикл коррекции в потоковом процессе). Поэтому, все что нужно,  нужно решить обратную задачу класса:

{Ax=B}

В общем виде она сводится к задаче проективного преобразования. Т.е. нам по известным координатам источника и целевой фигуры нужно найти финальную матрицу коэффициентов трансформации. Обычно, для обратимых матриц хорошо работает подход называемый DLT — Direct linear transformation или, к примеру, решение через уравнения вида Ax=B методом la_solve из пакета LAPACK (он входит в состав Accelerate SDK из iOS и OSX). Для решения нам соотвественно нужно подготовить данные для решателя. Получится, что-то вроде:

extension IMPQuad {
func transformTo(destination d:IMPQuad) -> float4x4{
var xy0 = left_bottom
var xy1 = left_top
var xy2 = right_bottom
var xy3 = right_top

var uv0 = d.left_bottom
var uv1 = d.left_top
var uv2 = d.right_bottom
var uv3 = d.right_top

var A:[Float] = [
xy0.x, xy0.y, 1, 0, 0, 0, -uv0.x * xy0.x, -uv0.x * xy0.y,
0, 0, 0, xy0.x, xy0.y, 1, -uv0.y * xy0.x, -uv0.y * xy0.y,
xy1.x, xy1.y, 1, 0, 0, 0, -uv1.x * xy1.x, -uv1.x * xy1.y,
0, 0, 0, xy1.x, xy1.y, 1, -uv1.y * xy1.x, -uv1.y * xy1.y,
xy2.x, xy2.y, 1, 0, 0, 0, -uv2.x * xy2.x, -uv2.x * xy2.y,
0, 0, 0, xy2.x, xy2.y, 1, -uv2.y * xy2.x, -uv2.y * xy2.y,
xy3.x, xy3.y, 1, 0, 0, 0, -uv3.x * xy3.x, -uv3.x * xy3.y,
0, 0, 0, xy3.x, xy3.y, 1, -uv3.y * xy3.x, -uv3.y * xy3.y,
]

var B:[Float] = [
uv0.x,
uv0.y,
uv1.x,
uv1.y,
uv2.x,
uv2.y,
uv3.x,
uv3.y
]

let lA = la_matrix_from_float_buffer( &A, 8, 8, 8, 0, la_attribute_t(LA_ATTRIBUTE_ENABLE_LOGGING) )
let lB = la_matrix_from_float_buffer( &B, 8, 1, 1, 0, la_attribute_t(LA_ATTRIBUTE_ENABLE_LOGGING) )

let lh = la_solve(lA, lB)

var h:[Float] = [Float](count:9, repeatedValue:1)
let status = la_vector_to_float_buffer( &h, 1, lh )

if status >= 0 {

let H:[float3] = [
float3(h[0], h[3], h[6]),
float3(h[1], h[4], h[7]),
float3(h[2], h[5], h[8])
]
let t = float3x3(rows:H)
return float4x4(rows: [
[t[0,0], t[1,0], 0, t[2,0]],
[t[0,1], t[1,1], 0, t[2,1]],
[0 , 0 , 1, 0 ],
[t[0,2], t[1,2], 0, t[2,2]]
])
}
else{
NSLog("IMPWarpSolver: error \(status)")
return float4x4()
}
}
}

Не знаю.., как-то не очень изящно. Мне больше нравится, не столь универсальный, но более короткий способ через присоединенные матрицы. В случае «квадратных» исходных данных, а они у нас «квадрантные», решение поиска матрицы трансформации T сводится к произведению матриц, где матрицы задающие систему линейных уравнений описывающие исходный и целевой четырех-угольники (с учетом требования квадратности дополняем матрицы единицами):

{\begin{bmatrix}x_0 & x_1 & x_2\\y_0 & y_1 & y_2 \\ 1 & 1 & 1\end{bmatrix} \begin{bmatrix}a_0\\ a_1\\ a_2\end{bmatrix}=\begin{bmatrix}x_3\\ y_3\\ 1\end{bmatrix}}

Решение системы, где координаты известны, а не известны коэффициенты:

{A_s = \begin{bmatrix}x_0 & x_1 & x_2 \\y_0 & y_1 & y_2\\1 & 1 & 1 \end{bmatrix}, \\ B_s = \begin{bmatrix}x_3 & y_3 & 1\end{bmatrix}, \\ \begin{bmatrix}a_0\\ a_1\\ a_2\end{bmatrix} = A^{-1}B}

Тогда присоединенная матрица будет произведением исходной и произведением диагональной матрицы составленной из вектора решения:

{D_s = \begin{bmatrix}a_0 & 0 & 0 \\0 & a_1 & 0\\0 & 0 & a_2 \end{bmatrix}, \\C_s = A_sD_s}

Точно также решаем систему для результирующего четырех-угольника и находим искомую матрицу трансформации:

 {T = C_dC_s^{-1}} , где  {C_d} — присоединенная матрица целевой фигуры,  {C_s} — присоединенная матрица исходного фигуры четырех-угольника

Теперь запишем вот все вот это в виде кода на Swift с использованием simd-классов для представления матриц:

extension IMPQuad {

/// Присоединенная матрица линейной системы уровнений описывающая четырех-угольник
public var adjugate:float3x3 {
get{
let A = float3x3(rows:[
[left_bottom.x,left_top.x,right_bottom.x],
[left_bottom.y,left_top.y,right_bottom.y],
[ 1, 1, 1]
])
let B = float3(right_top.x,right_top.y,1)
let X = A.inverse * B
// C = (Ai)*B
return A * float3x3(diagonal: X)
}
}

///
/// Вычисление матрицы трансформации
///
public func transformTo(destination d:IMPQuad) -> float4x4 {
let t = d.adjugate * adjugate.inverse

/// адаптируем 3x3 под наши практические нужды - переводим в эквивалетную запись
/// для системы 4x4
return float4x4(rows: [
[t[0,0], t[1,0], 0, t[2,0]],
[t[0,1], t[1,1], 0, t[2,1]],
[0 , 0 , 1, 0 ],
[t[0,2], t[1,2], 0, t[2,2]]
])
}

}

Графический шейдер трансформации

Все, что нам осталось реализовать в контексте реализации фильтра коррекции — это написать шейдер и обертку к нему:

vertex IMPVertexOut vertex_warpTransformation(
const device IMPVertex* vertex_array [[ buffer(0) ]],
const device float4x4 &homography_model [[ buffer(1) ]],
unsigned int vid [[ vertex_id ]])
{
IMPVertex in = vertex_array[vid];
float3 position = float3(in.position);

IMPVertexOut out;
//
// единственная операция преобразования - перемножение матрицы трансформации
// и координат текущей позиции
//
out.position = homography_model * float4(position,1);

out.texcoord = float2(float3(in.texcoord).xy);

return out;
}

fragment float4 fragment_passthrough(
IMPVertexOut in [[stage_in]],
texture2d<float, access::sample> texture [[ texture(0) ]]
)
{
// линейный семплер текстуры, хотя может быть и любой самописный -
// для повышения качества интерполяции
constexpr sampler s(address::clamp_to_edge, filter::linear, coord::normalized);
return texture.sample(s, in.texcoord.xy);
}

Обертка для исполнения рендеринга реализована в классе IMPWarpFilter. Код UI-приложения для интерактивной коррекции изображений можно скачать из проекта ImageMetalling, под номером 12. В коде контроллера основного окна можно почитать особенности имплементации реализованных нами фильтров.

PS

Соглашусь, пожалуй, что истории типа Фотосессия длиной в 4 года — годный и неприкрытый снобизмом фан, но практичнее иногда сделать вот так:

dolboed-loop

Буквально за несколько секунд.

 


Авторы блога не преследуют целей быть предельно корректным, но если заметили явную ашипку, если написали явную глупость, если что-то не понятно, есть идеи улучшения: комментируйте или пишите на: imagemetalling [*] gmail.com.

Реклама

Коррекция геометрических искажений: Один комментарий

Добавить комментарий

Заполните поля или щелкните по значку, чтобы оставить свой комментарий:

Логотип WordPress.com

Для комментария используется ваша учётная запись WordPress.com. Выход / Изменить )

Фотография Twitter

Для комментария используется ваша учётная запись Twitter. Выход / Изменить )

Фотография Facebook

Для комментария используется ваша учётная запись Facebook. Выход / Изменить )

Google+ photo

Для комментария используется ваша учётная запись Google+. Выход / Изменить )

Connecting to %s