Дискретное преобразование Хартли

Математическое преобразование, связанное с Фурье

Дискретное преобразование Хартли (DHT) — это преобразование Фурье дискретных периодических данных, похожее на дискретное преобразование Фурье (DFT), с аналогичными применениями в обработке сигналов и смежных областях. Его главное отличие от DFT заключается в том, что оно преобразует действительные входы в действительные выходы без внутреннего вовлечения комплексных чисел . Так же, как DFT является дискретным аналогом непрерывного преобразования Фурье (FT), DHT является дискретным аналогом непрерывного преобразования Хартли (HT), введенного Ральфом В. Л. Хартли в 1942 году. ^[1]

Поскольку существуют быстрые алгоритмы для DHT, аналогичные быстрому преобразованию Фурье (FFT), DHT был первоначально предложен Рональдом Н. Брейсвеллом в 1983 году ^[2] как более эффективный вычислительный инструмент в общем случае, когда данные являются чисто реальными. Однако впоследствии было высказано мнение, что специализированные алгоритмы FFT для реальных входов или выходов обычно можно найти с немного меньшим количеством операций, чем любой соответствующий алгоритм для DHT.

Определение

Формально дискретное преобразование Хартли представляет собой линейную обратимую функцию H : R ⁿ → R ⁿ (где R обозначает множество действительных чисел ). N действительных чисел x ₀ , ..., x _{N −1} преобразуются в N действительных чисел H ₀ , ..., H _{N −1} по формуле

${\displaystyle H_{k}=\sum _{n=0}^{N-1}x_{n}\operatorname {cas} \left({\frac {2\pi }{N}}nk\right)=\sum _{n=0}^{N-1}x_{n}\left[\cos \left({\frac {2\pi }{N}}nk\right)+\sin \left({\frac {2\pi }{N}}nk\right)\right]\quad \quad k=0,\dots ,N-1.}$

Эту комбинацию иногда обозначают как cas( z ) , и ее не следует путать с cis( z ) = e ^iz = cos( z ) + i  sin( z ) , или e ^{− iz} = cis(− z ), которые фигурируют в определении DFT (где i — мнимая единица ). ${\displaystyle \cos(z)+\sin(z)}$ ${\displaystyle ={\sqrt {2}}\cos \left(z-{\frac {\pi }{4}}\right)}$

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

Характеристики

Преобразование можно интерпретировать как умножение вектора ( x ₀ , ...., x _{N −1} ) на матрицу размером N на N ; следовательно, дискретное преобразование Хартли является линейным оператором . Матрица обратима; обратное преобразование, которое позволяет восстановить x _n из H _k , является просто DHT H _{k ,} умноженным на 1/ N . То есть DHT является своим собственным обратным ( инволютивным ) с точностью до общего масштабного коэффициента.

DHT можно использовать для вычисления DFT, и наоборот. Для действительных входных данных x _n выходное DFT X _k имеет действительную часть ( H _k + H _N−k )/2 и мнимую часть ( H _N−k − H _k )/2. Наоборот, DHT эквивалентно вычислению DFT x _{n ,} умноженного на 1 +  i , а затем взятию действительной части результата.

Как и в случае с DFT, циклическая свертка z = x ∗ y двух векторов x = ( x _n ) и y = ( y _n ) для получения вектора z = ( z _n ), все длиной N , становится простой операцией после DHT. В частности, предположим, что векторы X , Y , и Z обозначают DHT для x , y , и z соответственно. Тогда элементы Z задаются как:

${\displaystyle {\begin{matrix}Z_{k}&=&\left[X_{k}\left(Y_{k}+Y_{N-k}\right)+X_{N-k}\left(Y_{k}-Y_{N-k}\right)\right]/2\\Z_{N-k}&=&\left[X_{N-k}\left(Y_{k}+Y_{N-k}\right)-X_{k}\left(Y_{k}-Y_{N-k}\right)\right]/2\end{matrix}}}$

где мы берем все векторы периодическими по N ( X _N = X ₀ и т. д.). Таким образом, так же, как DFT преобразует свертку в поточечное умножение комплексных чисел ( пары действительных и мнимых частей), DHT преобразует свертку в простую комбинацию пар действительных частотных компонентов. Обратное DHT затем дает желаемый вектор z . Таким образом, быстрый алгоритм для DHT (см. ниже) дает быстрый алгоритм для свертки. (Это немного дороже, чем соответствующая процедура для DFT, не включая стоимость преобразований ниже, потому что попарная операция выше требует 8 действительных арифметических операций по сравнению с 6 для комплексного умножения. Это количество не включает деление на 2, которое может быть поглощено, например, в 1/ N нормализацию обратного DHT.)

Быстрые алгоритмы

Как и для DFT, оценка определения DHT напрямую потребует O( N ² ) арифметических операций (см. нотацию Big O ). Однако существуют быстрые алгоритмы, похожие на FFT, которые вычисляют тот же результат всего за O( N log N ) операций. Почти каждый алгоритм FFT, от Cooley–Tukey до prime-factor , Winograd (1985) ^[3] и Bruun (1993), ^[4], имеет прямой аналог для дискретного преобразования Хартли. (Однако несколько более экзотических алгоритмов FFT, таких как QFT, еще не были исследованы в контексте DHT.)

В частности, аналог DHT алгоритма Кули–Тьюки широко известен как алгоритм быстрого преобразования Хартли (FHT), и был впервые описан Брейсвеллом в 1984 году. ^[5] Этот алгоритм FHT, по крайней мере, применительно к размерам N , равным степени двойки , является предметом патента США номер 4,646,256, выданного в 1987 году Стэнфордскому университету . Стэнфорд передал этот патент в общественное достояние в 1994 году (Брейсвелл, 1995). ^[6]

Как упоминалось выше, алгоритмы DHT обычно немного менее эффективны (с точки зрения количества операций с плавающей точкой ), чем соответствующий алгоритм DFT (FFT), специализированный для действительных входов (или выходов). Это впервые было высказано Соренсеном и др. (1987) ^[7] и Дюамелем и Веттерли (1987). ^[8] Последние авторы получили то, что, по-видимому, является наименьшим опубликованным числом операций для DHT размеров степени двойки, используя алгоритм с разделенным основанием (похожий на алгоритм FFT с разделенным основанием ), который разбивает DHT длины N на DHT длины N /2 и два действительных входных DFT ( не DHT) длины N /4. Таким образом, они утверждали, что DHT длины степени двойки можно вычислить, в лучшем случае, с помощью 2 сложений больше, чем соответствующее число арифметических операций для действительного входного DFT.

На современных компьютерах производительность определяется в большей степени соображениями кэша и конвейера ЦП , чем строгими подсчетами операций, и небольшая разница в арифметической стоимости вряд ли будет значительной. Поскольку алгоритмы FHT и FFT с реальным входом имеют схожие вычислительные структуры, ни один из них, по-видимому, не имеет существенного априорного преимущества в скорости (Popović  [sr] и Šević, 1994). ^[9] С практической точки зрения, высокооптимизированные библиотеки FFT с реальным входом доступны из многих источников (например, от поставщиков ЦП, таких как Intel ), тогда как высокооптимизированные библиотеки DHT менее распространены.

С другой стороны, избыточные вычисления в БПФ из-за реальных входов сложнее устранить для больших простых N , несмотря на существование O( N log N ) алгоритмов комплексных данных для таких случаев, поскольку избыточности скрыты за сложными перестановками и/или фазовыми вращениями в этих алгоритмах. Напротив, стандартный алгоритм БПФ простого размера, алгоритм Рейдера , может быть напрямую применен к DHT реальных данных примерно в два раза меньше, чем эквивалентный комплексный БПФ (Frigo и Johnson, 2005). ^[10] С другой стороны, также возможна не основанная на DHT адаптация алгоритма Рейдера для ДПФ с реальными входами (Chu и Burrus , 1982). ^[11]

Многомерное дискретное преобразование Хартли (MD-DHT)

rD-DHT (MD-DHT с размерами «r») определяется как

${\displaystyle X(k_{1},k_{2},...,k_{r})=\sum _{n_{1}=0}^{N_{1}-1}\sum _{n_{2}=0}^{N_{2}-1}\dots \sum _{n_{r}=0}^{N_{r}-1}x(n_{1},n_{2},...,n_{r}){\rm {cas}}({\frac {2\pi n_{1}k_{1}}{N_{1}}}+\dots +{\frac {2\pi n_{r}k_{r}}{N_{r}}}),}$

с и где ${\displaystyle k_{i}=0,1,\ldots ,N_{i}-1}$ ${\displaystyle {\rm {cas}}(x)=\cos(x)+\sin(x).}$

Подобно случаю 1-D, как действительное и симметричное преобразование, MD-DHT проще, чем MD-DFT. Во-первых, обратное DHT идентично прямому преобразованию с добавлением масштабного коэффициента;

и во-вторых, поскольку ядро ​​является действительным, оно избегает вычислительной сложности комплексных чисел . Кроме того, DFT напрямую получается из DHT с помощью простой аддитивной операции (Брейсвелл, 1983). ^[2]

MD-DHT широко используется в таких областях, как обработка изображений и оптических сигналов. Конкретные приложения включают компьютерное зрение, телевидение высокой четкости и телеконференции, области, которые обрабатывают или анализируют движущиеся изображения (Zeng, 2000). ^[12]

Быстрые алгоритмы для MD-DHT

По мере того, как скорость вычислений продолжает расти, большие многомерные проблемы становятся вычислительно осуществимыми, требуя необходимости быстрых многомерных алгоритмов. Ниже приведены три таких алгоритма.

В целях обеспечения разделимости для повышения эффективности мы рассматриваем следующее преобразование (Брейсвелл, 1983), ^[2]

${\displaystyle {\hat {X}}(k_{1},k_{2},...,k_{r})=\sum _{n_{1}=0}^{N_{1}-1}\sum _{n_{2}=0}^{N_{2}-1}\dots \sum _{n_{r}=0}^{N_{r}-1}x(n_{1},n_{2},...,n_{r}){\rm {cas}}({\frac {2\pi n_{1}k_{1}}{N_{1}}})\dots {\rm {cas}}({\frac {2\pi n_{r}k_{r}}{N_{r}}}).}$

^{В работе Бортфельда (1995) [13]} было показано, что эти два явления могут быть связаны несколькими дополнениями. Например, в 3-D,

${\displaystyle X(k_{1},k_{2},k_{3})={\frac {1}{2}}[{\hat {X}}(k_{1},k_{2},-k_{3})+{\hat {X}}(k_{1},-k_{2},k_{3})+{\hat {X}}(-k_{1},k_{2},k_{3})-{\hat {X}}(-k_{1},-k_{2},-k_{3})].}$

Для , тогда могут быть реализованы алгоритмы строка-столбец. Этот метод широко используется из-за простоты таких алгоритмов RC, но они не оптимизированы для общих пространств MD. ${\displaystyle {\hat {X}}}$

Были разработаны и другие быстрые алгоритмы, такие как radix-2, radix-4 и split radix. Например, Boussakta (2000) ^[14] разработал 3-D векторный radix,

${\displaystyle X(k_{1},k_{2},...,k_{r})=\sum _{n_{1}=0}^{N-1}\sum _{n_{2}=0}^{N-1}\sum _{n_{r}=0}^{N-1}x(n_{1},n_{2},n_{3}){\rm {cas}}({\frac {2\pi }{N}}(n_{1}k_{1}+n_{2}k_{2}+n_{3}k_{3}))}$

${\displaystyle =\sum _{n_{1}:even}\sum _{n_{2}:even}\sum _{n_{3}:even}+\sum _{n_{1}:even}\sum _{n_{2}:even}\sum _{n_{3}:odd}+\sum _{n_{1}:even}\sum _{n_{2}:odd}\sum _{n_{3}:even}}$

${\displaystyle +\sum _{n_{1}:even}\sum _{n_{2}:odd}\sum _{n_{3}:odd}+\sum _{n_{1}:odd}\sum _{n_{2}:even}\sum _{n_{3}:even}+\sum _{n_{1}:odd}\sum _{n_{2}:even}\sum _{n_{3}:odd}}$

${\displaystyle +\sum _{n_{1}:odd}\sum _{n_{2}:odd}\sum _{n_{3}:even}+\sum _{n_{1}:odd}\sum _{n_{2}:odd}\sum _{n_{3}:odd}.}$

^{Также в работе Boussakta (2000) [14]} было показано , что этот алгоритм 3D-векторного радикса использует умножения и сложения по сравнению с умножениями и сложениями из подхода «строка-столбец». Недостатком является то, что реализацию этих алгоритмов радиксного типа трудно обобщить для сигналов произвольных измерений. ${\displaystyle ({\frac {7}{4}})N^{3}\log _{2}N}$ ${\displaystyle ({\frac {31}{8}})N^{3}\log _{2}N}$ ${\displaystyle 3N^{3}\log _{2}N}$ ${\displaystyle ({\frac {9}{2}})N^{3}\log _{2}N+3N^{2}}$

Теоретико-числовые преобразования также использовались для решения MD-DHT, поскольку они выполняют чрезвычайно быстрые свертки. В Boussakta (1988), ^[15] было показано, как разложить преобразование MD-DHT в форму, состоящую из сверток:

Для случая 2-D (случай 3-D также рассматривается в указанной ссылке),

${\displaystyle X(k,l)=\sum _{n=0}^{N-1}\sum _{m=0}^{M-1}x(n,m){\rm {cas}}({\frac {2\pi nk}{N}}+{\frac {2\pi ml}{M}}),\;}$ ${\displaystyle k=0,1,\ldots ,N-1}$, ${\displaystyle l=0,1,\ldots ,M-1}$

можно разложить на одномерные и двумерные круговые свертки следующим образом:

${\displaystyle X(k,l)={\begin{cases}X_{1}(k,0)\\X_{2}(0,l)\\X_{3}(k,l)\end{cases}}}$

где

${\displaystyle X_{1}(k,0)=\sum _{n=0}^{N-1}(\sum _{m=0}^{M-1}x(n,m)){\rm {cas}}({\frac {2\pi nk}{N}}),\;}$ ${\displaystyle k=0,1,\ldots ,N-1}$

${\displaystyle X_{2}(0,l)=\sum _{m=0}^{M-1}(\sum _{n=0}^{N-1}x(n,m)){\rm {cas}}({\frac {2\pi ml}{M}}),\;}$ ${\displaystyle l=1,2,\dots ,M-1}$

${\displaystyle X_{3}(k,l)=\sum _{n=0}^{N-1}\sum _{m=0}^{M-1}x(n,m){\rm {cas}}({\frac {2\pi nk}{N}}+{\frac {2\pi ml}{M}})\;,}$

${\displaystyle k=1,2,\ldots ,N-1}$

${\displaystyle l=1,2,\ldots ,M-1.}$

Развиваясь дальше, ${\displaystyle X_{3}}$

${\displaystyle X_{3}(k,l)=\sum _{n=0}^{N-1}x(n,0){\rm {cas}}({\frac {2\pi nk}{N}})+\sum _{m=1}^{M-1}x(0,m){\rm {cas}}({\frac {2\pi ml}{M}})}$

${\displaystyle +\sum _{n=1}^{N-1}\sum _{m=1}^{M-1}x(n,m){\rm {cas}}({\frac {2\pi nk}{N}}+{\frac {2\pi ml}{M}}).}$

На этом этапе мы представляем преобразование чисел Ферма (ПФЧ). t- ^е число Ферма задается как , причем . Хорошо известные числа Ферма для ( является простым для ), (Буссакта, 1988). ^[15] Преобразование чисел Ферма задается как ${\displaystyle F_{t}=2^{b}+1}$ ${\displaystyle b=2^{t}}$ ${\displaystyle t=0,1,2,3,4,5,6}$ ${\displaystyle F_{t}}$ ${\displaystyle 0\leq t\leq 4}$

${\displaystyle X(k,l)=\sum _{n=0}^{N-1}\sum _{m=0}^{M-1}x(n,m)\alpha _{1}^{nk}\alpha _{2}^{ml}\mod F_{t}}$

причем . и являются корнями единицы порядка и соответственно . ${\displaystyle k=0,\ldots ,N-1,l=0,\ldots ,M-1}$ ${\displaystyle \alpha _{1}}$ ${\displaystyle \alpha _{2}}$ ${\displaystyle N}$ ${\displaystyle M}$ ${\displaystyle (\alpha _{1}^{N}=\alpha _{2}^{M}=1\mod F_{t})}$

Возвращаясь к разложению, последний член для обозначим как , тогда ${\displaystyle X_{3}(k,l)}$ ${\displaystyle X_{4}(k,l)}$

${\displaystyle X_{4}(k,l)=\sum _{n=1}^{N-1}\sum _{m=1}^{M-1}x(n,m){\rm {cas}}({\frac {2\pi nk}{N}}+{\frac {2\pi ml}{M}}),}$

${\displaystyle k=1,2,\ldots ,N-1}$

${\displaystyle l=1,2,\ldots ,M-1.}$

Если и являются примитивными корнями чисел и (которые гарантированно существуют, если и являются простыми числами ), то и отображаются в Итак, отображая и в и , получаем следующее: ${\displaystyle g_{1}}$ ${\displaystyle g_{2}}$ ${\displaystyle N}$ ${\displaystyle M}$ ${\displaystyle M}$ ${\displaystyle N}$ ${\displaystyle g_{1}}$ ${\displaystyle g_{2}}$ ${\displaystyle (n,m)}$ ${\displaystyle (g_{1}^{n}\mod N,g_{2}^{m}\mod M).}$ ${\displaystyle n,m,k}$ ${\displaystyle l}$ ${\displaystyle g_{1}^{-n},g_{2}^{-m},g_{1}^{k}}$ ${\displaystyle g_{2}^{l}}$

${\displaystyle X_{4}(g_{1}^{k},g_{2}^{l})=\sum _{n=0}^{N-2}\sum _{m=0}^{M-2}x(g_{1}^{-n},g_{2}^{-m}){\rm {cas}}({\frac {2\pi g_{1}^{(-n+k)}}{N}}+{\frac {2\pi g_{2}^{(-m+l)}}{M}}),}$

${\displaystyle k=0,1,\ldots ,N-2}$

${\displaystyle l=0,1,\ldots ,M-2}$.

Что теперь является круговой сверткой . С , , и , имеем ${\displaystyle Y(k,l)=X_{4}(g_{1}^{k},g_{2}^{l})}$ ${\displaystyle y(n,m)=x(g_{1}^{-n},g_{2}^{-m})}$ ${\displaystyle h(n,m)={\rm {cas}}({\frac {2\pi g_{1}^{n}}{N}}+{\frac {2\pi g_{2}^{m}}{M}})}$

${\displaystyle Y(k,l)=\sum _{n=0}^{N-2}\sum _{m=0}^{M-2}y(n,m)h(<k-n>_{N},<l-m>_{M})}$

${\displaystyle Y(k,l)=FNT^{-1}\{FNT[y(n,m)]\otimes FNT[h(n,m)]}$

где обозначает почленное умножение. В (Boussakta, 1988) ^[15] также утверждалось , что этот алгоритм уменьшает количество умножений в 8–20 раз по сравнению с другими алгоритмами DHT за счет небольшого увеличения количества операций сдвига и сложения, которые, как предполагается, проще умножений. Недостатком этого алгоритма является ограничение, что каждое измерение преобразования имеет примитивный корень . ${\displaystyle \otimes }$

Ссылки

1. Хартли, Ральф ВЛ (март 1942 г.). «Более симметричный анализ Фурье, применяемый к проблемам передачи». Труды IRE . 30 (3): 144–150. doi :10.1109/JRPROC.1942.234333. S2CID  51644127.
2. Брейсвелл, Рональд Н. (1983). «Дискретное преобразование Хартли». Журнал оптического общества Америки . 73 (12): 1832–1835. doi :10.1364/josa.73.001832. S2CID  120611904.
3. Соренсен, Хенрик В.; Джонс, Дуглас Л .; Беррус, Чарльз Сидней ; Хайдеман, Майкл Т. (1985). «О вычислении дискретного преобразования Хартли». Труды IEEE по акустике, речи и обработке сигналов . ASSP-33 (4): 1231–1238.
4. Бини, Дарио Андреа; Боццо, Энрико (1993). «Быстрое дискретное преобразование с помощью собственных многочленов». Компьютеры и математика с приложениями . 26 (9): 35–52. doi : 10.1016/0898-1221(93)90004-f .
5. Брейсвелл, Рональд Н. (1984). «Быстрое преобразование Хартли». Труды IEEE . 72 (8): 1010–1018. doi :10.1109/proc.1984.12968. S2CID  21988816.
6. Брейсвелл, Рональд Н. (1995). «Вычисления с преобразованием Хартли». Компьютеры в физике . 9 (4): 373–379. Bibcode :1995ComPh...9..373B. doi : 10.1063/1.168534 .
7. Соренсен, Хенрик В.; Джонс, Дуглас Л .; Хайдеман, Майкл Т.; Беррус, Чарльз Сидни (1987). «Алгоритмы быстрого преобразования Фурье с действительными значениями». Труды IEEE по акустике, речи и обработке сигналов . ASSP-35 (6): 849–863.
8. Дюамель, Пьер; Веттерли, Мартин (1987). «Улучшенные алгоритмы преобразования Фурье и Хартли: применение к циклической свертке реальных данных». Труды IEEE по акустике, речи и обработке сигналов . ASSP-35: 818–824.
9. Поповић [Попович], Миодраг [Миодраг] [на сербском языке] ; Шевич, Драгутин (1994). «Новый взгляд на сравнение быстрых преобразований Хартли и Фурье». Транзакции IEEE по обработке сигналов . 42 (8): 2178–2182. Бибкод : 1994ITSP...42.2178P. дои : 10.1109/78.301854.
10. Фриго, Маттео; Джонсон, Стивен Г. (2005). «Проектирование и реализация FFTW3» (PDF) . Труды IEEE . 93 (2): 216–231. CiteSeerX 10.1.1.66.3097 . doi :10.1109/jproc.2004.840301. S2CID  6644892. }
11. Чу, Шуни; Беррус, Чарльз Сидни (1982). «Алгоритм FTT [ sic ] на основе простого множителя с использованием распределенной арифметики». Труды IEEE по акустике, речи и обработке сигналов . 30 (2): 217–227. doi :10.1109/tassp.1982.1163875.
12. Цзэн, Юнхан; Би, Гоан; Лейман, Абдул Рахим (2000). «Алгоритмы полиномиального преобразования для многомерного дискретного преобразования Хартли». Международный симпозиум IEEE по схемам и системам (V): 517–520.
13. Бортфельд, Томас; Динтер, Вольфганг (1995). «Вычисление многомерных преобразований Хартли с использованием одномерных преобразований Фурье». Труды IEEE по обработке сигналов . 43 (5): 1306–1310. Bibcode : 1995ITSP...43.1306B. doi : 10.1109/78.382424.
14. Boussakta, Said; Alshibami, Osama (2000). «Быстрый алгоритм для трехмерного дискретного преобразования Хартли». Международная конференция по акустике, речи и обработке сигналов '00 (4): 2302–2305.
15. Boussakta, Said; Holt, Alan GJ (1988). "Быстрое многомерное дискретное преобразование Хартли с использованием преобразования чисел Ферма". Труды IEE G - Электронные схемы и системы . 135 (6): 235–237. doi :10.1049/ip-g-1.1988.0036.

Дальнейшее чтение

• Брейсвелл, Рональд Н. (1986). Преобразование Хартли (1-е изд.). Oxford University Press . ISBN 978-0-19503969-6.
• Буссакта, Саид; Холт, Алан Г. Дж. (1988). «Быстрое многомерное дискретное преобразование Хартли с использованием преобразования чисел Ферма». Труды IEE G — Электронные схемы и системы . 135 (6): 235–237. doi :10.1049/ip-g-1.1988.0036.
• Хонг, Джонатан; Веттерли, Мартин ; Дюамель, Пьер (1994). «Преобразования базового поля со свойством свертки» (PDF) . Труды IEEE . 82 (3): 400–412. doi :10.1109/5.272145.
• О'Нил, Марк А. (1988). «Быстрее, чем быстрый Фурье». BYTE . 13 (4): 293–300.
• Olnejniczak, Kraig J.; Heydt, Gerald T. (март 1994). «Сканирование специального раздела преобразования Хартли». Труды IEEE . 82 : 372–380.(Примечание. Содержит обширную библиографию.)