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

Перед видеокамерой установлен вращающийся линейный поляризатор с дистанционным управлением, во время эксперимента получаем по 4 кадра статического изображения при 4^-х угловых положениях поляризатора 0°, 45°, 90° и 135°, а при повороте он вносит смещение в получаемое изображение на несколько точек. На этой установке уже получено много требующего обработки материала, и, чтобы не менять поляризатор на плёночный (тогда бы смещение не мешало), решим задачу по скоррекции смещения. Для этого проведён эксперимент с примерами смещения на изображении чёрно-белой настроечной таблицы:

Рис. 1а — уменьшенный общий вид (затемнение к краю кадра — это границы анализатора),

Рис. 1б — смещение на фрагменте с углом настроечной таблицы

После ряда экспериментов на глаз определяется, что
поворот на 45° даёт смещение ≈3-4 точки вправо и ≈1-2 точки вниз,
поворот на 90° даёт смещение ≈4-5 точек вправо и ≈5-6 точек вниз,
поворот на 135° даёт смещение ≈1,5-2 точки вправо и ≈8-9 точек вниз относительно полученного изображения при положении 0°.
Для дальнейшей обработки поляризационной картины полученные изображения, конечно, надо совместить, с этой целью нужно определить смещение точнее и составить алгоритм для коррекции.

Для изображения фрагмента угла настроечной таблицы достаточно провести одномерное определение смещений по горизонтали и по вертикали.

Рис. 2 — нерезкий край с плавным изменением яркости изображения, показанного на рис. 1б
(экспорт в таблицу Excel проведён с помощью программы Color Table).

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

Так как узлы интерполяции отстоят друг от друга на равное расстояние и упорядочены по величине, то можно использовать интерполяционную формулу Ньютона.
Число точек n + 1. x_n — значение конечной точки, а x₀ — начальной, зададим для удобства отсчёт с x₀ = 0.
Шаг между узлами: x_i+1–x_i = h, который зададим в единицах пикселей h = 1
Прямая интерполяционная формула Ньютона:
P(x) = I₀ + q⋅Δ¹I₀ + q⋅(q – 1)⋅Δ²I₀/2! + ... + (q⋅(q – 1)⋅...⋅(q – n + 1))⋅ΔⁿI₀/n!
Здесь q = (x – x₀)/h, а при взятых h = 1 и x₀ = 0 становится q = x.
I_i = f(x_i) — выборка значений функции яркости, распределённой в пространстве,
а выражения Δ^kI₀ — нулевые конечные разности k порядка, например:
Δ¹I₀ = I₁ – I₀,
Δ²I₀ = Δ¹I₁ – Δ¹I₀ = I₂ – 2⋅I₁ + I₀.

И для нашего случая с 7 точками:
P(x) = I₀ + x⋅Δ¹I₀ + x⋅(x – 1)⋅Δ²I₀/2! + x⋅(x – 1)⋅(x – 2)⋅Δ³I₀/3! + x⋅(x – 1)⋅(x – 2)⋅(x – 3)⋅Δ⁴I₀/4! + x⋅(x – 1)⋅(x – 2)⋅(x – 3)⋅(x – 4)⋅Δ⁵I₀/5! + x⋅(x – 1)⋅(x – 2)⋅(x – 3)⋅(x – 4)⋅(x – 5)⋅Δ⁶I₀/6!

Вычислив по данной формуле интерполяционную функцию, можно использовать как нахождение среднего значения I, так и вычислить первую производную по x — P'_x и вторую — P''_x. Нулевое значение второй производной наряду с максимумом первой производной (по абсолютному значению) в этой окрестности укажут на наиболее резкий перепад яркости. На рис. 3 показан график полученной интерполяции и её двух первых производных:

Рис. 3 —  график исходных значений яркости I(x) (см. рис. 2),
 жирная линия — функция интерполяции P(x) по значениям яркости,
 первая производная этой функции P'_x(x) и
 вторая производная P''_x(x).
В месте наиболее резкого яркостного перепада максимум |P'_x| и P''_x(x)=0, где x≈2,4.

Аналогично было сделано и для остальных 3 кадров. Численно определив пространственные координаты нулевого значения второй производной по x от функции интерполяции P(x) (при максимуме |P'(x)|), запишем результаты координат перепада яркости по горизонтали:
при 0° отсчёт от начала выборки x_0°≈2,4, в котором смещение берётся за 0,
при 45° отсчёт от начала выборки x_45°≈3+2,55 (3 — это смещение в целых пикселях относительно изображения при 0°), и общее смещение 5,55 – x_0° = 3,15 ,
при 90° отсчёт от начала выборки x_90°≈4+2,2, и общее смещение 6,2 – x_0° = 3,8 ,
при 135° отсчёт от начала выборки x_135°≈1+2,9, и общее смещение 3,9 – x_0° = 1,5 ,

А для смещения по вертикали
при 0° отсчёт от начала выборки y_0°≈ – 2,3 (где P''(y≈–2,3)=0), в котором смещение 0 ,
при 45° отсчёт от начала выборки y_45°≈2 – 2,1, и общее смещение вниз 0,1 – y_0° = 2,4 ,
при 90° отсчёт от начала выборки y_90°≈6 – 2,3, и общее смещение вниз  3,7 – y_0° = 6 ,
при 135° отсчёт от начала выборки y_135°≈10 – 3,2, и общее смещение вниз 6,8 – y_0° = 9,1 .

Теперь можно предложить простой линейный алгоритм совмещения сместившихся изображений (хотя в общем случае даже при квадратных фотоэлементах ПЗС матрицы алгоритм учёта смещения для большей точности должен быть более сложным). В одномерном случае при смещении на заданное нецелое число пикселей х_см. новую яркость точки будем вычислять как
I_нов.(x) = (x_см. – int(x_см.))⋅I(х – int(x_см.)) + (int(x_см. + 1) – x_см.)⋅I(х – int(x_см. + 1)),
при x – int(x_см.) ≥ 0,
где int(.) — функция выделения целой части числа.

В двухмерном случае последовательно проводится коррекция смещения по x, а затем по y:
I_нов.(y) = (y_см. – int(y_см.))⋅I(х – int(y_см.)) + (int(y_см. + 1) – y_см.)⋅I(х – int(y_см. + 1)),
при y – int(y_см.) ≥ 0,

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

Результат коррекции смещения оказался не таким успешным, как ожидалось.