Алгоритмы Винограда малоточечных ДПФ

Содержание

Введение

Представление ДПФ в матричной форме

Алгоритм Рейдера постановка задачи

Расчет ДПФ через циклическую свёртку

Расчет перестановочных матриц

Алгоритм 3-точечного ДПФ

Алгоритм 5-точечного ДПФ

Выводы

Примечания

Список литературы

DSPL-2.0 — свободная библиотека алгоритмов цифровой обработки сигналов

Распространяется под лицензией LGPL v3

Страница проекта на SourceForge

Введение

В предыдущих разделах мы рассмотрели алгоритмы БПФ по основанию два с прореживанием по времени и с прореживанием по частоте. Данные алгоритм очень эффективны, но они имеют существенное ограничение: длина входных данных и выходного вектора должна быть целой степенью двойки N = 2^t .

В данном разделе мы рассмотрим алгоритмы БПФ для произвольной составной длины $N = M \cdot L$ , где и — произвольные положительные целые. Данный подход позволит построить эффективные процедуры дискретного преобразования Фурье практически для всех длин, а не только для N = 2^t . Кроме того, мы покажем, что рассмотренные ранее алгоритмы БПФ по основанию два являются частным случаем алгоритма БПФ составной длины.

Представление ДПФ в матричной форме

Рассмотрим -точечное ДПФ сигнала s(n) , $n = 0\ldots N-1$ :

(1)

где $W_N = \exp(-j 2\pi /N)$ — основание поворотного коэффициента ДПФ.

Выражение (1) можно трактовать как произведение матрицы поворотных коэффициентов на вектор входного сигнала:

(2)

где $\mathbf{S}$ — вектор размерности $[N \times 1]$ спектральных отсчетов S(k)

; $\mathbf{W}$ — матрица размерности $[N \times N]$ поворотных коэффициентов $W_N^{nk}$ ; $\mathbf{s}$ — вектор значений входного сигнала s(n)

; индексы $n,k = 0\ldots N-1$ :

(3)

(4)

Квадратная матрица $\mathbf{W}$ содержит N^2

комплексных экспонент, повернутых на угол $2\pi nk / N$ , и как нетрудно убедится, матрица $\mathbf{W}$ содержит лишь

различных коэффициентов, которые можно отобразить на единичной окружности комплексной плоскости, как это показано на рисунке 1 для $N = 3, \,4, \, 5$ и N = 8

Поворотные коэффициенты ДПФ на комплексной плоскости для и

Рисунок 1. Поворотные коэффициенты ДПФ на комплексной плоскости для $N = 3, \,4, \, 5$ и N = 8

Все элементы матрицы $\mathbf{W}$ -точечного ДПФ, с учетом периодичности комплексной экспоненты могут быть представлены как степени $W_N^{(nk \mod N)}$ по модулю .

Ввиду того, что матрица $\mathbf{W}$ содержит лишь различных комплексных экспонент, она обладает рядом замечательных свойств. Во первых матрица $\mathbf{W}$ является унитарной с весом , т.е. эрмитово сопряжение приводит нас к масштабированной обратной матрице:

(5)

Таким образом обратное преобразование Фурье в матричной форме можно представить как:

(6)

Заметим, что первая строка и первый столбец матрицы $\mathbf{W}$ состоит из единиц, а остальные строки и столбцы имеют целые степени поворотных коэффициентов, которые повторяются ввиду периодичности комплексной экспоненты. Например для N = 5

матрица $\mathbf{W}$ равна:

(7)

Алгоритм Рейдера постановка задачи

Поскольку первая строка матрицы $\mathbf{W}$ содержит только единичные элементы, то первый спектральный отсчет (постоянная составляющая сигнала)

(8)

есть сумма всех входных отсчетов.

С другой стороны, первый столбец матрицы $\mathbf{W}$ также содержит только единичные значения. Это значит, что спектральный отсчет с индексом можно представить как:

(9)

Таким образом, в выражении 9 расчет S(k)

сводится к умножение «подматрицы» $\mathbf{W}$ размерности $[(N-1) \times (N-1)]$ , не содержащей нулевые степени поворотных коэффициентов на «подвектор» входного сигнала, не содержащий нулевой отсчет, как это показано на рисунке 2 при исключении их матричного приведения всех тривиальных множителей первого столбца и первой строки.

Рисунок 2. ДПФ в матричной форме при исключении тривиальных умножений на единицу

Тогда вычисление ДПФ для индексов $k = 1\ldots N-1$ в матричном виде можно представить для N=5 как:

(10)

Можно показать, что при значении длины

равной простому числу (3, 5, 7, 11 и т.д.) в матрице поворотных коэффициентов не останется тривиальных коэффициентов $\pm 1$ и $\pm j$ . При этом размеры ДПФ в виде простых множителей очень важны, потому что эффективный алгоритм для составной длины, рассмотренный в предыдущем разделе требует расчета ДПФ простой длины.

Таким образом необходимо найти эффективный алгоритм для матрично-векторного умножения в выражении (10). Разработка такого алгоритма заняла весьма продолжительное время. В 1968 году Рейдер опубликовал первый важный шаг этого алгоритма: сведение матрично-векторного умножения в выражении (10) к расчету циклической свертки. И спустя 10 лет, в 1978 году, Виноград опубликовал статью, в которой приводил максимально эффективные алгоритмы расчета циклической свертки, с минимальным числом операций умножения.

Расчет ДПФ через циклическую свёртку

При рассмотрении циклической свёртки мы говорили, что ее можно представить в матричном виде как:

(11)

где матрица $\mathbf{B}$ носит название матрицы циркулянта, потому что каждая стока или столбец данной матрицы представляет собой циклическую задержку предыдущей строки или столбца. Из структуры матрицы циркулянта размерности $[M \times M]$ можно заметить, что она содержит только

различных элементов, причем каждая строка и каждый стоблец также содержит

различных элементов.

Обозначим «подматрицу» поворотных коэффициентов в выражении (10) через $\mathbf{Y}$ , а «подвектор» входного сигнала через $\mathbf{x}$ . Тогда произведение матрицы $\mathbf{Y}$ на вектор $\mathbf{x}$ можно представить как:

(12)

Прибавив ко всем элементам вектора $\mathbf{y}$ выражения (12) значение s(0)

получим выходные спектральные отсчеты ДПФ.

Как можно видеть из выражения (12) матрица поворотных коэффициентов $\mathbf{Y}$ не является циркулянтом, но как и циркулянт содержит только N-1 различный элемент, причем как и в циркулянте все строки и столбцы содержат N-1 различный коэффициент $W_N^{nk}$ .

Рейдер в своей работе показал, что если размер ДПФ представляет собой простое число (3, 5, 7, 11 и т.д.), то переставляя строки и столбцы матрицы $\mathbf{Y}$ размерности $[(N-1) \times (N-1)]$ («подматрицы» $\mathbf{W}$ ), мы всегда можем свести ее к матрице циркулянту, а расчет $\mathbf{y = Yx}$ к циклической свертке. Для перестановки столбцов матрицы $\mathbf{Y}$ , умножим обе части уравнения (12) на перестановочную матрицу $\mathbf{P}_0$ слева:

(13)

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

Для перестановки строк матрицы $\mathbf{Y}$ перепишем (13) виде:

(14)

где $\mathbf{P}_1$ — матрица перестановок строк

, а $\mathbf{P}_1^{-1}$ — обратная матрице перестановок $\mathbf{P}_1$ , чтобы не нарушать равенства.

Как мы знаем, матрицы перестановок являются ортогональными, что означает $\mathbf{P}_1^{-1} = \mathbf{P}^T$ и выражение (14) можно записать:

(15)

где $\mathbf{P}_1^{T}\mathbf{x}$ задает перестановку входного вектора $\mathbf{x}$ , а $\mathbf{P}_0 \mathbf{y}$ — перестановку выхода циклической свертки, потому что $\mathbf{P}_0 \mathbf{Y} \mathbf{P}_1$ теперь циркулянт и можно рассчитывать произведение (15) как циклическую свертку.

Например для N = 5 матрицы перестановок равны:

(16)

тогда

(17)

и окончательно алгоритм Рейдера позволяет свести расчет спектральных отсчетов для индексов $k = 1\ldots N-1$ к циклической свертке вида:

(18)

Но перед тем как применить алгоритмы Винограда для коротких циклических сверток, необходимо понять как получить матрицы перестановок $\mathbf{P}_0$ и $\mathbf{P}_1$ . И тут надо отметить, что для ДПФ простой длины N > 3

таких матриц может быть несколько. Например для N = 5

имеется две комбинации перестановочных матриц $\mathbf{P}_0$ и $\mathbf{P}_1$ . Первую мы уже рассмотрели выше, а вторая комбинация матриц имеет вид:

(19)

которая приводит к выражению ДПФ Рейдера вида:

(20)

Для N = 7 комбинаций можно также найти две комбинации перестановочных матриц, а вот для N = 11 и N = 17 таких комбинаций будет уже 4, что приводит к возможности построения различных циркуляционных матриц алгоритма ДПФ. Для простых длин большей длины количество различных комбинаций перестановочных матриц может быть различным, но всегда найдется хотя бы одна такая комбинация перестановочных матриц, которая сведет расчет ДПФ простой длины к циклической свертке согласно алгоритму Рейдера.

Расчет перестановочных матриц

Для того чтобы построить эффективные алгоритмы ДПФ простой длины необходимо научиться рассчитывать перестановочные матрицы $\mathbf{P}_0$ и $\mathbf{P}_1$ для заданной длины . Тогда мы сможем свести расчет ДПФ к циклической свертке, согласно алгоритму Рейдера и использовать эффективные алгоритмы Винограда для расчета циклической свертки.

Для начала вспомним малую теорему Ферма, которая утверждает, что для любого простого и натурального $1 \leq a < N$ выполняется равенство:

(21)

Например для N = 5 имеем:

(22)

Кроме того, из теории чисел известно, что для простого числа существуют примитивные корни , представляющие собой натуральные числа 1 < r < N , такие, что

(23)

Для того, что понять является ли некоторое 1 < r < N примитивным корнем простого числа необходимо рассчитать вектор значений $\mathbf{v} = \left(r^{n} \right) \Mod N$ для $n = 1\ldots N-2$ , и если в данном векторе значений не встречается единица, то значит — примитивный корень .

Рассмотрим пример для N=5 . Возможные значения $r = 2,\,\,3,\,\,4$ , а возможные значения $n = 1,\,\,2,\,\,3,\,\,4$ . Тогда можно составить таблицу, как это показано на рисунке 3

Рисунок 3. Примитивные корни N=5

Из рисунка 3 можно заключить, что r = 2 и r = 3 являются примитивными корнями для простого N = 5 .

Тогда каждый примитивный корень будет порождать матрицу перестановок $\mathbf{P}_1$ . Для того чтобы сгенерировать матрицу $\mathbf{P}_1$ , размерности $[N-1\times N-1]$ необходимо установить в единицу элементы строк с индексами $\big(\left(r^{n} \right) \Mod N \big)-1$ для столбцов $n = 0 \ldots N-2$ (индексация строк и столбцов матрицы начинается с нуля ). Так для примитивного корня r = 2 для N = 5 необходимо установить в единицу следующие элементы матрицы $\mathbf{P}_1$ :

(24)

Таким образом, матрица $\mathbf{P}_1$ равна:

(25)

Матрицу $\mathbf{P}_0$ можно рассчитать при известной матрице $\mathbf{P}_1$ из уравнения:

(26)

Тогда пара перестановочных матриц для примитивного корня r = 2

имеет вид:

(27)

Если взять примитивный корень r = 3

для

, то можно получить еще одну комбинацию перестановочных матриц $\mathbf{P}_0$ и $\mathbf{P}_1$ :

(28)

Заметим, что матрица $\mathbf{Q}$ сохраняет свою структуру и для других размерностей: единичными являются элементы с индексами Q[0,0]

и все элементы на диагонали ниже побочной диагонали[1]:

(29)

Таким образом, мы можем найти различные матрицы перестановок для сведения алгоритма Рейдера к расчету циклической свертки. Применение алгоритма Винограда для коротких циклических сверток позволит получить эффективные алгоритмы малоточечных ДПФ, которые могут быть использованы в алгоритме составной длины.

Алгоритм 3-точечного ДПФ

ДПФ размерности N = 3 точки можно записать в виде алгоритма Рейдера как:

(30)

(31)

Матрица поворотных коэффициентов в алгоритме Рейдера при N=3

представляет собой циркулянт, поэтому никаких перестановок не требуется и мы можем применить алгоритм Винограда для двухточечной циклической свертки:

(32)

который имеет вид:

(33)

Учтем, что $W_3 = \exp(-j2\pi/3)$ , тогда нетрудно показать, что

(34)

(35)

Объединив алгоритм Рейдера (30)–(31) с алгоритмом Винограда циклической свертки (33), с учетом свойств поворотных коэффициентов (34)–(35) получим алгоритм для 3-точечного ДПФ вида:

(36)

Заметим, что слагаемое s(0)

в алгоритме Рейдера (30)–(31) учтено в m_1

, а при расчете S(0)

использовалась ранее посчитанная сумма t_1

. Таким образом 3-точечное ДПФ требует 2 умножения на чисто вещественный множитель $-\frac{1}{2}$ и чисто мнимый $j \sin(\pi/3)$ , а также 6 сложений.

Алгоритм 5-точечного ДПФ

ДПФ размерности N = 5 точки можно записать в виде алгоритма Рейдера как (20), которые требует расчета 4-точечной циклической свертки.

Алгоритм Винограда для 4-точечной циклической свертки представлен в параграфе «Алгоритм Винограда для расчета циклической свертки», там же представлена оптимизированная версия алгоритма с минимальным числом сумматоров [1, стр. 67]:

(37)

Сопоставляя данный алгоритм с циклической сверткой выражения (20) можно сделать заключения:

отсчеты входного сигнала $s(1)\ldots s(4)$ соответствуют значениям циклической свертки $a_0 \ldots a_3$ :

(38)

Поворотные коэффициенты $W_5\ldots W_5^4$ соответствуют значениям циклической свертки $b_0 \ldots b_3$ :

(39)

Коэффициенты $r_1 \ldots r_7$ зависят от поворотных коэффициентов и могут быть просчитаны однократно до использования алгоритма.
Выходные значения свертки $c_0\ldots c_3$ используются для расчета спектральных отсчетов согласно алгоритму Рейдера как:

(40)

Учтем также, что

(41)

Переобозначим промежуточные суммы $s_1 \ldots s_7$ при расчете циклической свертки как $t_1 \ldots t_7$ , чтобы не путать промежуточные суммы и сигнал.

Обобщив все замечания, алгоритм 5-точечного ДПФ можно записать как

(42)

Из выражения (42) можно выразить коэффициенты r_1

, а также разности поворотных коэффициентов W_5 - W_5^4

, как это показано на рисунке 4.

Рисунок 4. Пересчет поворотных коэффициентов

Можно видеть, что

(43)

а также

(44)

Тогда

(45)

Заметим, что поворотные множители r_3

являются вещественными, а $r_5 \ldots r_7$ — чисто мнимыми. Это существенно упрощает расчет произведений $m_0 \ldots m_4$ .

Также можно слагаемое s(0) при расчете спектральных отсчетов учесть в произведении m_0 , которое используется в q_0 и q_1 . Окончательно, алгоритм Винограда 5-точечного ДПФ может быть представлен как:

(46)

Полученный алгоритм требует 17 сложений и 5 умножений на предварительно рассчитанные множители $r_3 \ldots r_7$ .

Выводы

В данном разделе мы рассмотрели алгоритм Винограда расчета коротких ДПФ по простому основанию. Ядра коротких ДПФ могут быть использованы для построения быстрых алгоритмов составной длины, отличной от степени двойки. Это позволяет расширить применение алгоритмов БПФ и построить эффективные процедуры практически для любой длины исходной выборки.

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

В данном разделе приведен алгоритм расчета перестановочных матриц для приведения ДПФ к расчету циклической свертки, а также приведены алгоритмы для длин N = 3 и N = 5 .

Смотри также

Алгоритм БПФ составной длины
Алгоритм Винограда для расчета циклической свертки
Алгоритм БПФ с прореживанием по времени
Алгоритм БПФ с прореживанием по частоте

Примечания

[1] Замечание автора. Структура матрицы $\mathbf{Q}$ найдена мной, как автором, эмпирическим путем, который заключался в том, что анализировались свойства перестановочных матриц $\mathbf{P}_0$ и $\mathbf{P}_1$ и выявлении свойства $\mathbf{Q}=\mathbf{P}_0\mathbf{P}_1$ . В результате структура $\mathbf{Q}$ сохранялась для всех перестановочных матриц, которые позволяли свести ДПФ к циклической свертке. Мне неизвестно доказательство того, что данный алгоритм работоспособен для всех простых длин ДПФ. Однако я проверил данный алгоритм численно для всех простых длин до N=59 включительно и для всех проверенных длин, алгоритм расчета перестановочных матриц давал правильный результат. Для длин больше N = 59 я не смог проверить, ввиду численной неустойчивости алгоритма расчета примитивных корней. Если кто-то из читателей найдет аналитическое доказательство корректности подхода, равно как найдет опровержение, то прошу незамедлительно сообщить.

Список литературы

[1] Нуссбаумер Г. Быстрое преобразование Фурье и алгоритмы вычисления сверток. Москва, Радио и связь, 1985.

[2] James W. Cooley and John W. Tukey An Algorithm for the Machine Calculation of Complex Fourier Series. Mathematics of Computation, Vol. 19, no. 2, April 1965, pp. 297-301.

[3] Оппенгейм А., Шаффер Р. Цифровая обработка сигналов. Москва, Техносфера, 2012. 1048 с. ISBN 978-5-94836-329-5

[4] Bracewell R. The Fourier Transform and Its Applications McGraw-Hills, 1986, 474 c. ISBN 0-07-007-015-6

[5] Гольденберг Л.М., Матюшкин Б.Д., Поляк М.Н. Цифровая обработка сигналов. Москва, Радио и связь, 1985.

[6] Блейхут Р. Быстрые алгоритмы цифровой обработки сигналов. Москва, Мир, 1989, 448 c.

[7] Сергиенко А.Б. Цифровая обработка сигналов СПб, Питер, 2002.

[8] Рабинер, Л., Гоулд, Б. Теория и применение цифровой обработки сигналов. Москва, Мир, 1978. 848 с.

Последнее изменение страницы: 20.01.2024 (19:21:33)

Страница создана Latex to HTML translator ver. 5.20.11.14