Быстрое преобразование Хафа (англ. Fast Hough Transform, сокр. БПХ) — модификация преобразования Хафа, позволяющая параметрически идентифицировать прямые (а также с дополнительной модификацией отрезки и четырехугольники) с меньшей вычислительной сложностью за счет использования факта самопересечения рассматриваемых дискретных прямых.

Алгоритм впервые был предложен М. Л. Брейди в 1992 году,^[1] в дальнейшем переизобретался несколько раз различными авторами.^[2][3]

Пусть дано изображение размера ${\displaystyle m\times n}$. Рассмотрим ${\displaystyle m^{2}}$ диадических прямых (прямых специального вида), состоящих из ${\displaystyle n}$ пикселей на изображении каждая (один пиксель в столбце).

Пусть ${\displaystyle I^{(s,t)}(i)}$ — интенсивность ${\displaystyle i}$ -го пикселя, принадлежащего диадической прямой, заданной параметрами ${\displaystyle s,t}$; ${\displaystyle S^{(s,t)}}$ — Хаф-образ этой диадической прямой.

Образ дискретного преобразование Хафа определяется следующей формулой:

${\displaystyle S^{(s,t)}=\sum _{i=0}^{n-1}(I^{(s,t)}(i))}$

Прямое вычисление всех значений ${\displaystyle S^{(s,t)}}$ требует ${\displaystyle O}$${\displaystyle \left(m^{2}n\right)}$ операций: перебор по ${\displaystyle m}$ различных значений параметров ${\displaystyle s}$, ${\displaystyle t}$, перебор для каждой пары ${\displaystyle n}$ значений ${\displaystyle i}$.

В свою очередь, алгоритм БПХ, основанный на учёте пересечений отрезков друг с другом, требует ${\displaystyle \Theta }$${\displaystyle \left(mn\log m\right)}$ действий, на одну прямую необходимо ${\displaystyle \Theta }$${\displaystyle \left({\frac {n\log m}{m}}\right)}$ операций (для квадратных изображений ${\displaystyle \Theta }$${\displaystyle \left(\log n\right)}$). По теореме, сформулированной Т. М. Ханиповым,^[4] невозможно сложить диадические прямые с асимптотически меньшей вычислительной сложностью.

Алгоритм основан на принципе «разделяй и властвуй». Задача состоит в нахождении сумм значений пикселей вдоль отрезков, соединяющих «левый» и «правый» края изображения. Изображение делится пополам, в каждой части задача решается независимо. Для получения итоговой суммы на каждом из отрезков складываются ответы на «левой» и «правой» его половинах.

В алгоритме БПХ пиксели прямых произвольного вида дискретно приближаются диадическими прямыми. Пиксели диадического приближения прямой на изображении размера ${\displaystyle n\times n}$ удалены от исходной прямой не более, чем на ${\displaystyle O}$${\displaystyle (\log n)}$ пикселей.^[5]

Отрезки параметризуются центрами соединяемых пикселей. Поэтому разбиение отрезка на подотрезки лишь приближенно составляет исходный отрезок. Ошибка аппроксимации с шагами рекурсии накапливается, однако не более чем на ${\displaystyle O}$${\displaystyle (\log n)}$ пикселей.^[5] Построенная таким образом дискретизация отрезка до пикселей называется диадической аппроксимацией.

Далее паттерном называется множество пикселей ${\displaystyle p=\{(x_{i},i)\}_{i=0}^{2^{k}-1}}$, содержащее по элементу в каждой вертикали изображения. Отклонением паттерна будет называться значение ${\displaystyle t(p)=x_{2^{k}-1}-x_{0}}$, а координатой — значение ${\displaystyle s(p)=x_{0}}$. Сдвигом паттерна будет называться множество $${\displaystyle p\nearrow (a,b)=\lbrace \left(x_{i}+a,y_{i}+b\right)\ |\ (x_{i},y_{i})\in p\rbrace }$$Порождающие диадические паттерны ширины ${\displaystyle 2^{k}}$ и наклона ${\displaystyle t<2^{k}}$ определяются рекурсивно. При ${\displaystyle k=t=0}$ паттерн состоит из одного пикселя, а при ${\displaystyle k>0}$ выражается через ${\displaystyle k-1}$.

${\displaystyle P_{t}^{k}={\begin{cases}\{(0,0)\}&,k=0\\\ P_{t/2}^{k-1}\ \ \ \cup \ P_{t/2}^{k-1}\ \ \ \ \nearrow (t/2,2^{k})&,k>0,\ \ t\equiv 0\ (mod2)\\P_{(t-1)/2}^{k-1}\cup \ P_{(t-1)/2}^{k-1}\nearrow ((t+1)/2,2^{k})&,k>0,\ \ t\equiv 1\ (mod2)\end{cases}}}$

Преимущественно горизонтальные, направленные вверх, диадические прямые получаются из всех сдвинутых по вертикали порождающих паттернов ${\displaystyle P_{s,t}^{k}=P_{t}^{k}\nearrow (s,0)}$, построенных со всеми возможными координатами на изображении ${\displaystyle s\in \{0,...,2^{k}-1\},t\in \{0,...,2^{k}-1\}}$.

Для приближенного вычисления преобразования Хафа требуется найти суммы ${\displaystyle S_{s,t}=\sum _{(x,y)\in P_{s,t}^{k}}I_{x,y}}$ по всем диадическим прямым на изображении. В этой сумме ${\displaystyle 2^{k}\cdot 2^{k}}$ прямых по ${\displaystyle 2^{k}}$ точек каждая. За счет рекурсивного перехода это суммирование сводится к подсчету отдельно левых, отдельно правых половин, что позволяет свести вычисление к подсчету ${\displaystyle 2\cdot 2^{k}\cdot 2^{k-1}}$ сумм по ${\displaystyle 2^{k-1}}$ точек каждая.

Рассмотрим бинарные слова, состоящие из чисел 0 и 1. Множество диадических слов определяется рекурсивно. ${\displaystyle w}$ будет называться диадическим словом, если ${\displaystyle w}$ имеет вид ${\displaystyle a0a}$ или ${\displaystyle a1a}$, где ${\displaystyle a}$ — пустое или диадическое слово. Все диадические слова длины не больше трех: 0, 1, 000, 010, 101, 111.

Для каждого диадического слова ${\displaystyle {\overline {b_{1}b_{2}...b_{n}}}}$рассматривается кумулятивная сумма ${\displaystyle \{a_{i}\}_{0}^{2^{k}-1},a_{i}=\sum _{0}^{i}b_{i}}$. Будем говорить, что последовательность пикселей ${\displaystyle \{(a_{i},i)\}_{0}^{2^{k}-1}}$ — диадическая прямая, соединяющая центры пикселей ${\displaystyle (0,0)}$и ${\displaystyle (a_{2^{k}-1},2^{k}-1)}$.

Существует ровно ${\displaystyle 2^{k}}$диадических прямых длины ${\displaystyle 2^{k}}$. По одной на каждую пару ${\displaystyle (0,0)}$и ${\displaystyle (h,2^{k}-1),0\leqslant h<2^{k}}$.

Алгоритм БПХ устроен следующим образом:^[6]

Исходное состояние матрицы — исходное изображение размера ${\displaystyle 2^{p}\times 2^{q}}$. Затем происходит вычисление на ${\displaystyle q-1}$-их уровнях по очереди, начиная с первого: на ${\displaystyle i}$-том уровне матрица в текущем состоянии делится на группы по принципу равенства целой части координаты второй оси после деления на ${\displaystyle 2^{i}}$; получаются подматрицы ${\displaystyle 2^{p}\times 2^{i}}$; объединим соседние в пары (без пересечений, это возможно, так как размер матрицы — степень двойки) и назовем в этой паре первой ту подматрицу, что расположена на меньших координатах по второй координате в матрице, а другую — второй; вместо первых в каждой паре записывается ее сумма с соответствующей ей второй, а вместо второй — сумму первой и второй с циклическим сдвигом на один влево. Таким образом считается Хаф-образ таких прямых, что для любой пары точек с координатами ${\displaystyle (x_{1},y_{1}),(x_{2},y_{2})}$ с этой прямой выполняется ${\displaystyle (x_{1}\neq x_{2})\rightarrow ((y_{1}\neq y_{2})\wedge (x_{1}<x_{2}\Leftrightarrow y_{1}<y_{2}))}$, с помощью аппроксимации диадическими прямыми. Для того, чтобы посчитать образ для остальных прямых, достаточно поворачивать изображение и проводить ту же самую операцию, а полученные результаты сложить.

Полученные таким образом на каждом уровне матрицы являются элементами БПХ-пирамиды. Формальное описание БПХ-пирамиды: Нулевым уровнем БПХ-пирамиды является исходное изображение (размером ${\displaystyle m\times n}$, а последним — Хаф-образ, содержащий суммы по диадическим прямым длины ${\displaystyle n=2^{p}}$. Для описания ${\displaystyle q}$-ого уровня пирамиды исходное изображение разбивается на горизонтальные полосы ${\displaystyle [0,n-1]\times [k\cdot 2^{q},(k+1)\cdot 2^{q}-1]}$, где ${\displaystyle k}$ — номер полосы, ${\displaystyle k\in [0,2^{p-q}-1]}$. Для каждой полосы в ${\displaystyle q}$-ом уровне БПХ-пирамиды хранятся суммы по всевозможным паттернам полосы, имеющим длину ${\displaystyle 2^{q}}$ и параметры ${\displaystyle s\in [0,n-1],t\in [0,2^{q}-1]}$. Количество таких паттернов для полосы — ${\displaystyle n\cdot 2^{q}}$, поэтому ${\displaystyle q}$-й уровень пирамиды занимает столько же памяти, сколько исходное изображение.

Инвариантность количества затрачиваемой памяти и возможность хранить каждый уровень в матрице того же размера, что оригинальное изображение, без потери интерпретируемости дает следующее свойство: можно хранить БПХ-пирамиду в форме матрицы с размерностью, на одну большую, чем размерность исходного изображения (по одной оси — количество уровней, ${\displaystyle \log _{2}n}$), по остальным — размер изображения).^[7]

Пример реализации на языке python:

import numpy as np
W = 2**5
H = 2**5
img = np.random.random([H, W])

def calc_sums(img, xmin, xmax):
    res = np.zeros([W, xmax-xmin])
    if xmax - xmin == 1:
        res[:, 0] = img[:, xmin]
    else:
        mid = (xmin + xmax) // 2
        ans1 = calc_sums(img, xmin, mid)
        ans2 = calc_sums(img, mid, xmax)
        for x in range(W):
            for shift in range(xmax-xmin):
                res[x, shift] = ans1[x, shift//2] + ans2[(x + shift//2 + shift%2) % W, shift//2]
    return res
res = calc_sums(img, 0, W)

Алгоритм реализован в библиотеке opencv^[8] и может быть использован, например, для быстрого нахождения точки схода.^[9]

Решение этой задачи подразумевает использование алгоритма для двумерного случая.

Хаф-образ плоскостей будет так же трехмерным (плоскость задается через три координаты вектора, перпендикулярного ей). Пусть он будет задан параметризацией ${\displaystyle s,t,t_{2}}$, где ${\displaystyle s}$ — координата ${\displaystyle x}$ пересечения плоскости с границей изображения на луче ${\displaystyle Ox}$, ${\displaystyle s+t}$ — координата ${\displaystyle x}$ точки пересечения с границей изображения, параллельной лучу ${\displaystyle Ox}$ в плоскости ${\displaystyle xOz}$, ${\displaystyle t_{2}}$ — разница координат ${\displaystyle x}$ второй и первой точки пересечения плоскости с границами изображения. Первая точка находится в пересечении содержащих границу изображения плоскостей ${\displaystyle xOz}$ и плоскости, параллельной ${\displaystyle xOy}$. Вторая точка — в пересечении содержащих границу изображения плоскостей, параллельных ${\displaystyle xOy}$ и ${\displaystyle xOz}$.

Будем называть плоскость преимущественно-перпендикулярной к оси координат, если нормаль к ней образует с данной осью меньший угол, чем с двумя другими. Будем рассматривать только плоскости, преимущественно-перпендикулярные к оси y. Они разделяются на 4 типа наклонов, как это показано на рисунке 4. Без ограничения общности будем считать, что рассматриваемые плоскости имеют тип I.

Построение Хаф-образа методом перебора плоскостей имеет асимптотическую сложность ${\displaystyle O}$${\displaystyle \left(mnk\cdot (mn+nk+mk)\right)}$ (количество плоскостей, умноженное на количество операций для суммирования одной плоскости), где m, n, k — размеры изображения в каждой размерности.

Быстрым преобразованием Хафа в данном случае будет следующий алгоритм:

1. Для каждой плоскости перпендикулярной оси ${\displaystyle Oz}$ с координатой ${\displaystyle z}$ по этой оси считается быстрое преобразование Хафа, а результат кладется в трехмерное пространство по координатам ${\displaystyle s,t,z}$.
2. Для каждой плоскости в получившемся трехмерном пространстве перпендикулярной оси ${\displaystyle Ot}$ с координатой ${\displaystyle t}$ по этой оси считается быстрое преобразование Хафа, а результат кладется в куб по координатам ${\displaystyle s,t,t_{2}}$.

Сложность такого преобразования складывается из сложности первого шага (${\displaystyle \Theta }$${\displaystyle \left(k\cdot mn\log m\right)}$) и сложности второго шага (${\displaystyle \Theta }$${\displaystyle \left(n\cdot mk\log m\right)}$), которые считаются, как произведение количества рассматриваемых плоскостей и количества операций на одну плоскость. Итого, ${\displaystyle \Theta }$${\displaystyle \left(mnk\log m\right)}$, в пересчете на одну плоскость ${\displaystyle \Theta }$${\displaystyle (\log m)}$.

Хаф-образ трехмерной прямой будет четырехмерным (по два параметра на каждую из двух точек на прямой). Пусть он будет задан параметризацией ${\displaystyle s,t,s_{2},t_{2}}$. ${\displaystyle s,s_{2}}$ — координаты x, y точки на плоскости ${\displaystyle xOy}$, ${\displaystyle s+t,s_{2}+t_{2}}$ — координаты x, y точки пересечения прямой с границей изображения, параллельной плоскости ${\displaystyle xOy}$.

Построение Хаф-образа методом перебора трехмерных прямых имеет асимптотическую сложность ${\displaystyle O}$${\displaystyle \left((mn)^{2}\cdot k\right)}$ (количество прямых умноженное на количество операций для суммирования одной прямой), где m, n, k — размеры изображения в каждой размерности.

Быстрое преобразование Хафа для такого случая формулируется подобно определению для двумерного случая. В двумерном случае возможность сдвига была только вдоль одной оси, теперь же сдвиг будет вдоль одной оси, вдоль второй оси и вдоль двух осей одновременно.

Подсчет паттернов длины два требует ${\displaystyle {\frac {k}{2}}}$ (количество групп суммируемых плоскостей), умноженное на ${\displaystyle 4mn}$ (сложность проводимых сложений для каждой группы), операций. Подсчет паттернов длины 4 требует ${\displaystyle {\frac {k}{4}}\cdot 16mn}$ операций. Длины паттернов ${\displaystyle d}$ — ${\displaystyle {\frac {k}{d}}\cdot d^{2}mn=mnkd=mnk\cdot 2^{i}}$, где ${\displaystyle i}$ определяется, как ${\displaystyle \log _{2}(d)}$, то есть номер рассматриваемой длины паттерна. Суммируя ${\displaystyle \log _{2}(k)}$ слагаемых (количество возможных для рассматриваемого изображения длин паттернов) с использованием формулы суммы геометрической прогрессии, получается сложность алгоритма ${\displaystyle \Theta }$${\displaystyle \left(\sum _{i=1}^{\log _{2}k}mnk\cdot 2^{i}\right)\equiv }$${\displaystyle \Theta }$${\displaystyle \left({\frac {2mnk\cdot (2^{\log _{2}k}-1)}{2-1}}\right)\equiv }$${\displaystyle \Theta }$${\displaystyle \left(mnk^{2}\right)}$, а сложность на одну прямую ${\displaystyle \Theta }$${\displaystyle \left({\frac {k^{2}}{mn}}\right)}$. При ${\displaystyle m=n=k}$ сложность будет константной.

Несмотря на то, что количество операций на одну прямую константное при равном размере изображения в каждой размерности, для всех прямых надо потратить ${\displaystyle n^{4}}$. Но если не нужны все прямые, а нужна только часть их, то можно совершить предпосчет первых ${\displaystyle x}$ шагов^[10] , сохранить их в память, а затем посчитать суммы только для тех прямых, которые нужны.

Этот концепт был закреплен в методе четырех русских. Метод назван в честь первоткрывателей В. Арлазаров, М. Кронрод, Е. Диниц, И. Фараджев.

В оригинальном алгоритме БПХ для трехмерных прямых совершается на каждом из ${\displaystyle \log _{2}k}$ уровней подсчет для прямых определенной длины. С другой стороны можно сделать предпосчет только для первых ${\displaystyle x}$ шагов, а затем производить подсчет для нужных прямых.

Для того, чтобы определить оптимальное количество шагов предпосчета, необходимо решить такое уравнение (${\displaystyle a}$ — количество прямых, которое необходимо найти алгоритму):

${\displaystyle \Theta \left(mnk2^{x}\right)=\Theta \left(a2^{\log _{2}k-x}\right)}$

Слева — количество операций для выполнения предпосчета. Справа — количество операций на нахождение сумм вдоль запрошенных прямых. Пусть необходимо найти все прямые, тогда ${\displaystyle a=mnk}$, тогда решением уравнения будет ${\displaystyle x={\frac {1}{2}}\log _{2}k=\log _{2}{\sqrt {k}}}$, а левая и правая части равны ${\displaystyle \Theta \left(mnk{\sqrt {k}}\right)}$, что в ${\displaystyle {\sqrt {k}}}$ меньше, чем без предпосчета. Тем не менее за уменьшение количества операций необходимо заплатить памятью в том же размере, что занимает Хаф-образ (свойство инвариантности занимаемой памяти на каждом уровне подсчета алгоритмом БПХ).

Принцип вычисления основан на использовании значений не только последнего уровня БПХ-пирамиды (то есть непосредственно самого Хаф-образа), но и других уровней БПХ-пирамиды.

Задача делится на две подзадачи:

1. Найти диадическую прямую проходящую через два заданных пикселя
2. Из суммы значений вдоль этой прямой выделить ту часть суммы, которая относится к паттерну между заданными пикселями

Без ограничения общности считаем, что ${\displaystyle y_{2}\geq y_{1}}$. Здесь будут рассмотрены только преимущественно вертикальные паттерны с наклоном вправо, то есть ${\displaystyle x_{2}\geq x_{1}}$ и ${\displaystyle y_{2}-y_{1}\geq x_{2}-x_{1}}$. Также используется ${\displaystyle (s,t)}$-параметризация и значение ${\displaystyle p}$ равное ${\displaystyle \log _{2}n}$, где ${\displaystyle n}$ — размер изображения по оси ${\displaystyle y}$.

Пусть двоичное разложение параметра ${\displaystyle t}$ диадической прямой выглядит как ${\displaystyle t=\sum _{i=0}^{p-1}t_{i}2^{i}}$ Тогда паттерн может быть записан так (${\displaystyle \left[\cdot \right]}$ — округление до ближайшего целого.):

${\displaystyle x=s+\sum _{r=0}^{p-1}t_{r}\left[{\frac {s^{r}y}{2^{p}-1}}\right]}$

${\displaystyle \Delta x=x_{2}-x_{1}}$ считается из данных задачи. ${\displaystyle \delta _{r}^{y_{1},y_{2}}=\left(D_{r}(y_{2})-D_{r}(y_{1})\right)}$ — число сдвигов рассматриваемого паттерна в полосе ${\displaystyle (y_{1},y_{2})}$, которое также известно. Таким образом необходимо только восстановить биты ${\displaystyle t_{0},t_{1},\ldots ,t_{p-1}}$.

Для восстановления используется «жадный» алгоритм: Все биты ${\displaystyle t_{r}}$ сначала равны нулю. Так как ${\displaystyle i>j\rightarrow \delta _{i}^{y_{1},y_{2}}\geq \delta _{j}^{y_{1},y_{2}}}$, поэтому перебор ведется от большего числа сдвигов к меньшему, от уровня ${\displaystyle p-1}$ до уровня 0. Если ${\displaystyle \Delta y\geq \delta _{r}^{y_{1},y_{2}}}$, то соответствующий этому уровню бит устанавливается в 1, а ${\displaystyle \Delta y}$ уменьшается на ${\displaystyle \delta _{r}^{y_{1},y_{2}}}$. Шаг повторяется, пока ${\displaystyle \Delta y}$ не обратится в 0.

По ${\displaystyle t_{i}}$ вычисляется значение параметра ${\displaystyle t}$. Через этот параметр вычисляется параметр ${\displaystyle s}$ по следующей формуле:

${\displaystyle s=D(y_{1},t)-\sum _{r=0}^{p-1}t_{r}D_{r}(y_{1})}$

${\displaystyle \#\lbrace \delta _{r}^{y_{1},y_{2}}\rbrace =p=\log _{2}n}$, поэтому сложность алгоритма ${\displaystyle \Theta }$${\displaystyle (\log n)}$.^[7]

Обращаясь к рисунку, можно увидеть, что произвольный отрезок на прямой длины ${\displaystyle n}$ вычисляется, нахождением минимального количества диадических паттернов, в которых содержатся части от начала прямой до конца заданного отрезка включительно, и минимальное количество паттернов, которые содержат отрезок от начала прямой до начала заданного отрезка, исключая первый пиксель исходного отрезка. То есть нужно найти суммы для двух отрезков с началом на границе изображения и разными координатами конца.

Для подсчета суммы по такого типа отрезка длины ${\displaystyle l}$ (его двоичное разложение ${\displaystyle \sum _{j=0}^{p}\left(l_{j}\cdot 2^{j}\right)}$) ${\displaystyle \sum _{q=0}^{p}\left[l_{q}\cdot H^{s,t}\left(q,\sum _{j=q+1}^{p}l_{j}\cdot 2^{j-q}\right)\right]}$, где ${\displaystyle H^{s,t}\left(q,k\right)}$ — сумма по паттерну ${\displaystyle 2^{q}}$ в ${\displaystyle k}$-ой полосе ${\displaystyle q}$-ого уровня БПХ=пирамиды для прямой с параметрами ${\displaystyle (s,t)}$.

Внутренняя сумма не требует полного вычисления на каждом шаге, поскольку получается из предыдущей за константное время. Таким образом, трудоемкость алгоритма пропорциональна числу слагаемых во внешней сумме, то есть составляет ${\displaystyle \Theta }$${\displaystyle (\log n)}$. Так как результат вычисляется для двух отрезков этого типа, то результирующая сложность алгоритма так же ${\displaystyle \Theta }$${\displaystyle (\log n)}$. Причем стоит заметить, что пиксель может быть многоканальным.^[7]

Отрезок может быть составлен из минимального количества паттернов внутри отрезка. Для поиска таких паттернов нужно смотреть уровни БПХ-пирамиды, начиная с последних и заканчивая первыми. Сразу можно отсеять те паттерны, которые не входят в отрезок. Если попадается паттерн полностью лежащий внутри отрезка, то его сумма включается в искомую сумму, а его разбиения на следующих уровнях не рассматриваются. Этот способ более вычислительно сложный, чем первый, так как требует перебора всех паттернов, которых более ${\displaystyle n}$.

Аналогично вычислению суммы по отрезку для вычисления суммы по четырехугольнику из промежуточных вычислений Хаф-образа для плоскостей, иначе говоря БПХ-пирамидой для случая плоскостей.

Предполагая, что параметры ${\displaystyle (s,t,t_{2})}$ плоскости на которой находится заданный четырехугольник известны, искомая сумма вычисляется по формуле включения-исключения беря сумму по четырем прямоугольникам, одна вершина которых — угловая вершина диадической плоскости (обозначим ее буквой ${\displaystyle E}$, а сегменты с этой вершиной угловыми сегментами). Обозначим координаты точек ближней и дальней от ${\displaystyle E}$ вершин заданного четырехугольника через ${\displaystyle (i_{1},j_{1})}$ и ${\displaystyle (i_{2},j_{2})}$ соответственно. Суммы обозначенных угловых сегментов с вершинами в ${\displaystyle (i_{1}-1,j_{1}-1)}$ и ${\displaystyle (i_{2},j_{2})}$ берутся со знаком плюс, а суммы тех, что с вершинами в ${\displaystyle (i_{1}-1,j_{2})}$ и ${\displaystyle (i_{2},j_{1}-1)}$, берутся со знаком минус.

Чтобы найти сумму по произвольному угловому сегменту, необходимо разбить его на сегменты, присутствующие в БПХ-пирамиде. Необходимо рассмотреть двоичные разложения ширины и высоты сегмента. Аналогично одномерному случаю сегмент разбивается по горизонтали на ${\displaystyle \log _{2}n}$ вертикальных полос и по вертикали не более чем на ${\displaystyle \log _{2}n}$ горизонтальных полос. Их пересечение даст не более чем ${\displaystyle \log _{2}n}$ сегментов, присутствующих в трехмерной БПХ-пирамиде. Таким образом, сложность вычисления суммы по произвольному сегменту составляет ${\displaystyle \Theta }$${\displaystyle (\log n)}$ операций.^[7]

Хотя у апроксимации прямой диадическим паттерном есть некоторая ошибка, впрочем эксперименты показывают, что это ошибка достаточно незначительна для того, чтобы в решениях задач можно было заменить традиционный алгоритм преобразования Хафа на алгоритм БПХ.^[11]

Применяя М-оценки к задаче линейной регрессии, можно получить радиально-базисные функции. Они составляют «непрерывное» изображение, которое в свою очередь дискретизуется в двумерную гистограмму.

Далее выполняеncя свёртка изображения с дискретизованным ядром ${\displaystyle K}$, соответствующим выбранной М-оценке. На основе полученного вычисляется Хаф-образ с помощью БПХ. Координата максимума в пространстве параметров — и будет являться искомой М-оценкой.

Задача формулируется следующим образом: необходимо найти гиперплоскость, которая разделяет изображение на 2 класса. Изображение представляется в виде нормированной гистограммы изображения ${\displaystyle P}$.

${\displaystyle L}$ — искомая гиперплоскость, делящая изображения на два класса ${\displaystyle A_{L},B_{L}}$, ${\displaystyle \Omega }$ — класс всех элементов гистограммы.

Используемые аддитивные статистики ${\displaystyle X\subset \Omega }$ (${\displaystyle x_{i}}$ — ${\displaystyle i}$-ая координата ${\displaystyle x}$):

• ${\displaystyle D^{0}(X)=\sum _{x\in X}P(x)}$
• ${\displaystyle D_{i}^{1}(X)=\sum _{x\in X}x_{i}P(x)}$
• ${\displaystyle D_{ij}^{2}(X)=\sum _{x\in X}x_{i}x_{j}P(x)}$

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

Наивный алгоритм: Существует ${\displaystyle O(n^{2})}$ дискретных прямых, пересекающих гистограмму с линейным размером ${\displaystyle n}$. Для каждой из них требуется выполнить ${\displaystyle O(n^{2})}$ операций по вычислению весов, матриц ковариации и, наконец, значений критерия. Таким образом, вычислительная сложность наивного алгоритма составляет ${\displaystyle O(n^{4})}$ операций. Аналогично можно показать, что для трехмерного случая вычислительная сложность алгоритма будет составлять ${\displaystyle O(n^{6})}$.

На этом этапе применяется кумулятивное суммирование: в элемент слоя с номером ${\displaystyle i}$ выходного изображения записывается сумма соответствующих элементов всех слоев входного изображения с индексом, не превосходящим ${\displaystyle i}$.

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

Алгоритм:

1. вычислим набор изображений, содержащих значения необходимых аддитивных статистик для каждого из элементов входной гистограммы (${\displaystyle D^{0},D_{i}^{1},D_{ij}^{2}}$) (6 в двумерном случае, 10 в трехмерном)
2. вычисление кумулятивной суммы вдоль каждой из осей, получим кортеж изображений. Для любого изображения этого кортежа, относящегося к размерности ${\displaystyle d}$, сумма по любой гиперплоскости, преимущественно перпендикулярной оси с индексом ${\displaystyle d}$, равна соответствующей аддитивной статистике, вычисленной по полупространству, включающему начало координат и ограниченному выбранной гиперплоскостью. Зная значение аддитивной статистики для одного полупространства, легко получить значение той же статистики и для второго — вычитанием из статистики, посчитанной по всему изображению.
3. теперь, вычислив аддитивные статистики по всем разделениям гиперплоскостей, можно считать значения критерия для выбора оптимальной кластеризации.

Если просто суммировать по всем гиперплоскостям то в двумерном случае сложность ${\displaystyle O(mn^{3})}$, в трехмерном ${\displaystyle O(mn^{5})}$. (В ${\displaystyle k}$-мерном ${\displaystyle O(mn^{2k-1})}$)

Суммирование по гиперплоскостям (прямым в двумерном случае, плоскостям в трехмерном…) можно сделать с помощью БПХ. Это помогает снизить сложность с ${\displaystyle O(n^{2k-1})}$ до ${\displaystyle O(n^{k}\log {\left(n\right)})}$ для каждого из ${\displaystyle m}$ изображений. То есть теперь сложность в двумерном случае ${\displaystyle O(mn^{2}\log {\left(n\right)})}$, в трехмерном ${\displaystyle O(mn^{3}\log {\left(n\right)})}$.

Таким образом, финальный алгоритм:

1. Кумулятивная суммация
2. Подсчет аддитивных статистик
3. БПХ
4. Поиск максимума в пространстве Хафа

Сложность: время ${\displaystyle O(mn^{k}\log {\left(n\right)})}$, память ${\displaystyle O(mn^{k})}$.

В двумерном случае более подробно:

1. Кумулятивная суммация: ${\displaystyle O(n^{2})}$
2. Подготовка к вычислению аддитивных статистик: ${\displaystyle O(n^{2})}$
3. БПХ: ${\displaystyle O\left(n^{2}\log {\left(n\right)}\right)}$
4. Поиск максимума в пространстве Хафа: ${\displaystyle O(n^{2})}$

Итоговая сложность: ${\displaystyle O(n^{2}\log {\left(n\right)})}$

В трехмерном случае более подробно:

1. Кумулятивная суммация: ${\displaystyle O(n^{3})}$
2. Подготовка к вычислению аддитивных статистик: ${\displaystyle O(n^{3})}$
3. БПХ: ${\displaystyle O\left(n^{3}\log {\left(n\right)}\right)}$
4. Поиск максимума в пространстве Хафа: ${\displaystyle O(n^{3})}$

Итоговая сложность: ${\displaystyle O(n^{3}\log {\left(n\right)})}$

Далее представлены только некоторые задачи, которые могут быть решены с использованием преобразования Хафа.

• Отслеживание равномерно движущихся объектов, используя покадровую разность изображений. Эти объекты на своих треках оставляют выраженные прямые линии.^[12][13]
• Детекция точки схода на изображении. Точка схода — это точка на плоскости изображения в которой пересекаются проекции параллельных в трехмерной сцене прямых.^[14][15]
• Томографическое восстановление. Процедуру формирования проекций изображения анализируемого объекта с использованием рентгеновских лучей принято моделировать преобразованиями Хафа и Радона, а получение трехмерной структуры исследуемого объекта часто сводится к решению обратного преобразования Хафа или Радона.^[16]
• Анализ медицинских изображений.^[17]
• В реализации алгоритмов слепой калибровки радиальной дисторсии при условии нахождения прямолинейных объектов на сцене. Через оптимизацию нового функционала Хаф-образа подбираются параметры компенсации радиальной дисторсии.^[18]
• Определение степени сбития камеры. Основано на подсчете БПХ от эпиполярной картины и поиске прямой, на которой лежат точки интересующих прямых на эпиполярной картине.
• Распознавание рукописного текста.^[19]
• Определение наклона шрифта. Основано на том, что в шрифте есть символы состоящие из прямых отрезков, находящихся под одним углом, вдоль такого угла Хаф-образ будет большее значение.^[20]
• Распознавание штрих-кодов.^[21][22]
• Определение степени сходства форм.^[23]
• Векторизации трехмерных изображения.^[24]
• Детекция треков спутников по изображениям с длинной выдержкой.^[25]
• Обнаружение радиолокационных целей.^[26][27]
• Анализа деформаций подземного профиля.^[28]
• Анализ структуры топологии микросхем по фотографиям.^[29]
• Подсчет количества осей транспортного средства по трекам детектора колес изображений, полученных с камеры, снятой сбоку.^[30]
• Трехмерная реконструкции плоских граней прозрачных минералов по набору изображений.^[31]
• Анализ SAR снимков.^[32]