АНАЛИЗ ТОЧНОСТИ РЕАЛИЗАЦИИ АЛГОРИТМОВ БПФ

8. АНАЛИЗ ТОЧНОСТИ РЕАЛИЗАЦИИ АЛГОРИТМОВ БПФ

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

Методику оценки вычислительных ошибок БПФ рассмотрим на примере реализации БПФ по основанию 2 и с прореживанием по времени. Рассматриваемая методика может быть применена и для анализа ошибок других алгоритмов БПФ. Будем предполагать, что: обрабатываемые числа представляются с помощью b₁+1 разрядов, а коэффициенты - с помощью b₂+1 разрядов с учетом знака; для аппроксимации произведений используется операция округления; масштабирование промежуточных результатов производится на входе каждой операции «бабочка» путем сдвига чисел на один разряд вправо (деление на два); входные данные пронормированы таким образом, что, и подчиняются равномерному закону распределения т. е. имеют математическое ожидание равное нулю, и дисперсию  равную 1/3. Следовательно, среднеквадратическое значение (СК3) входной последовательности равно также 1/3. В соответствии с теоремой Парсеваля

выходная последовательность X(k) БПФ будет иметь СК3 N/3, где N-размер преобразования. С целью уяснения методики анализа ошибок получим алгоритм БПФ в аналитическом виде. Для этого в выражение для N=2^v -точечного ДПФ (5.1)

(8.1)

необходимо подставить двоичное разложение коэффициентов п и k:

(8.2)

в результате алгоритм БПФ можно представить, как ранее убедились, в виде v+1-ступенчатого процесса.

На ступени т=0 производится двоично-инверсная перестановка входной последовательности

На каждой из остальных v ступеней (т=1, 2,….,v) производится преобразование типа «бабочка» выходной последовательности предыдущей ступени.

Так, на ступени т = 1 производится преобразование последовательности

 

X₀(n₁,…,п_v):

На ступени m=2 -преобразование последовательности Х₁ (п₁,..., n_v-1, k₁)

На m-й ступени

 (8.3)

Так постепенно в двоичном представлении индекса последовательности Х_m с увеличением т происходит замена коэффициентов n_i на k_j

Наконец при т = v

Выходная последовательность последней ступени является искомой:

 

X(k)=X(k_v2^v-1 +... +k¹)=X_v(k_v 2^v-1 +k_v-12^v-2+... +k₁).

Представим индекс элемента т-й ступени в виде

 

(n₁,..., n_v-m,k_m,,..., k₁)=2^v-1n₁+2^v-2n₂+… +2^mn_v-m+2^m-1k_m+…+k₁=2^ml+q (8.4)

 

Тогда число Xm(2ml+q) можно рассматривать как q-й элемент l-го блока т-й ступени.

Пример б. Рассмотрим описанную процедуру синтеза алгоритма БПФ с прореживанием по времени на примере 16-точечного ДПФ. В этом случае v=4. Индексы n и k представим следующим образом: n=n₄2³+n₃2²+n₂2+n₁, k=k₄2³ +k₃2² +k₂2+k_1.

Подставляя n и k в выражение для 16-точечного ДПФ, получаем

 

X(k₄k₃k₂k₁)=    (n₄n_Зn₂n₁)

Теперь распишем алгоритм по ступеням:

т = 0 - инверсия входной последовательности:

 

Х₀(n₁ 2³+n₂ 2²+n₃ 2+n₄)=x(n₄ 2³+n_3'2²+n₂2+n₁);

 

т= l - преобразование «бабочка» последовательности Х₀:

 

т = 2 - преобразование «бабочка» последовательности X₁:

m=3 - преобразование «бабочка» последовательности Х₂:

m=4 - преобразование «бабочка» последовательности Х_З:

Искомая выходная последовательность

Рисунок 7

Направленный граф полученного в примере 16-точечного БПФ с указанием источников элементарных ошибок округления произведений и масштабирования, а также путей их распространения изображен на рис.7. Ошибки округления появляются в местах умножения комплексных чисел на нетривиальные поворачивающие множители.

Ошибки масштабирования появляются на обоих входах каждой операции «бабочка» из-за сдвига входных чисел на один разряд вправо (деление на два).

Элементарные ошибки округления и масштабирования, возникающие на различных ступенях алгоритма, приводят к тому, что отсчеты ДПФ X(k) на выходе вычисляются неточно.

Обозначим ошибку вычисления k-ro отсчета ДПФ через e(k)=X'(k) -X(k), где Х'(k)-вычисленное значение отсчета.

Анализ траекторий распространения элементарных ошибок, возникающих на различных ступенях алгоритма, позволяет сделать следующие выводы:

1. Элементарная ошибка с дисперсией , возникающая на т-й ступени алгоритма. вызывает ошибку с дисперсией  в 2^v-m выходных точках БПФ.

2. Дисперсия ошибки каждого выходного отсчета алгоритма, обусловленной округлением произведений, равна сумме ошибок, 2^v-m из которых вызывается ошибками т-й ступени.

Перейдем к анализу арифметических ошибок. В силу того, что математическое ожидание всех элементарных ошибок равно нулю, математическое ожидание результирующей ошибки E(e(k))=0 для всех k. Тогда СКЗ ошибок ДПФ будет определяться только дисперсией элементарных ошибок.

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

Вычислим СКЗ ошибок e(k), k=0,...,N-1. Из анализа алгоритма БПФ (8.3) и соответствующего направленного графа следует, что нетривиальные умножения, связанные с вычислением отсчета Х (k), появляются, начиная со ступени

(8.5)

где.  Это объясняется тем, что первое появление

единицы в двоичном представлении , , соответствует умножению на коэффициент -1 на ступени т, тогда на следующей т+ l-й ступени наличие коэффициента  = 1 приведет к умножению на j, а уже на т+2 и далее ступенях все умножения будут нетривиальными. Это хорошо проследить, анализируя выражение (8.3), описывающее т-ю ступень алгоритма. Можно показать путем вычислений s (k), что для определения отсчетов с номерами k=0, N/4; N/2; 3N/4 не требуется выполнение нетривиальных умножений, все умножения здесь на ; . Ошибка округления этих отсчетов равна нулю. А для вычисления СКЗ ошибок округления отсчетов с другими k воспользуемся вторым выводом. В результате получим

(8.6)

Пример 7. Рассмотрим вычисление ошибок округления отсчетов 16-точечного БПФ. Для этого выпишем для каждого номера отсчета его двоичное представление . Затем пользуясь выражением (8.6), вычисляем ошибки округления. Результаты расчетов приводятся в табл.1. Заметим, что отсчеты с номерами k=4, 8, 12 вычисляются точно так, как s(k)>4, где 4-максимальное число возможных ступеней алгоритма, т. е. в формировании отсчетов с номерами 0, 4, 8, 12 участвуют умножения только на тривиальные множители  (см. также граф на рис.5).

Вычислим СКЗ ошибки, усредненное по всем k. На первых двух ступенях алгоритма (т=1,2) все умножения являются точными, так как множители тривиальны ; на выходах т-й ступени (т = 3,…, у) появляется N произведений в  блоках, причем четыре произведения в каждом блоке являются точными в результате умножения также на тривиальные множители.

Таблица 1

Тогда, пользуясь выводом 1, получаем

(8.7)

В случае 16-точечного БПФ

Определим СКЗ выходной ошибки алгоритма, обусловленной масштабированием промежуточных результатов. Ошибка масштабирования комплексного числа на m-й ступени алгоритма имеет нулевое среднее и дисперсию . Множитель  появляется из-за масштабирования на m предыдущих ступенях.

В соответствии с выводом 2 СКЗ индивидуальной ошибки равно

(8.8)

в случае 16-точечного БПФ

Усредненное по всем выходам k СК3 также равно

(8.9)

Среднеквадратическое значение суммарной ошибки, обусловленной округлением и масштабированием, вычисляется как сумма отдельных СК3:

(8.10)

СК3 суммарной ошибки, усредненное по всем выходам алгоритма, равно

 (8.11)

В случае 1б-точечного БПФ

Результаты вычисления СК3 суммарных ошибок также приведены в табл.1. Точность алгоритма принято оценивать по отношению «СК3 ошибки / СК3 выходного сигнала». В данном случае оно составляет

____«СКЗ ошибки»______ = (8.12)

«СКЗ выходного сигнала»

В случае 1б-точечного БПФ

___«СК3 ошибки»_______ =

«СК3 выходного сигнала»

Как видно из полученных выражений, точность алгоритма зависит от двух параметров: длины преобразования N и разрядности представления чисел b_1. Если известны требования по точности вычисления и размер преобразования, разрядность процессора БПФ можно определить из следующего выражения, полученного логарифмированием выражения (8.12):

 (8.13)

Теперь перейдем к анализу ошибок БПФ, вызванных неточным представлением значений поворачивающих множителей  Пусть

 - точные значения коэффициентов Фурье, а '(k) - неточные, полученные при условии представления коэффициентов конечным числом разрядов

.

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

(8.14)

где

(8.15)

Дисперсия элементарных ошибок округления  комплексных коэффициентов равна

Индивидуальная ошибка БПФ, обусловленная округлением коэффициентов, равна

(8.16)

Вычитая (8.15) из (8.14), получаем

 члены высшего порядка.

Пренебрегая членами высшего порядка и подставляя последнее выражение в (8.16), получаем следующее СК3 ошибки, вызванной квантованием коэффициентов:

(8.17)

По теореме Парсеваля СК3 выходного сигнала равно

(8.18)

Отсюда «СК3 ошибки / СК3 выходного сигнала» = (v/б) 2^-2b₂. В случае 16-точечного БПФ «СК3 ошибки/СК3 выходного сигнала» = (2/3) 2^-2b₂.

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

ЗАКЛЮЧЕНИЕ

До середины 1960-х для представления спектрального разложения использовались точные формулы для нахождения параметров синусов и косинусов. Соответствующие вычисления требовали как минимум N**2 (комплексных) умножений. Таким образом, даже сегодня высокоскоростному компьютеру потребовалось бы очень много времени для анализа даже небольшого временного ряда (для 8,000 наблюдений потребовалось бы, по меньшей мере 64 миллиона умножений).

Ситуация кардинально изменилась с открытием так называемого алгоритма быстрого преобразования Фурье, или БПФ для краткости. Достаточно сказать, что при применении алгоритма БПФ время выполнения спектрального анализа ряда длины N стало пропорционально N*log2(N) что конечно является огромным прогрессом.

Однако недостаток стандартного алгоритма БПФ состоит в том, что число данных ряда должно быть равным степени 2 (т.е. 16, 64, 128, 256,...). Обычно это приводит к необходимости добавлять нули во временной ряд, который, как описано выше, в большинстве случаев не меняет характерные пики периодограммы или оценки спектральной плотности. Тем не менее, в некоторых случаях, когда единица времени значительна, добавление констант во временной ряд может сделать результаты более громоздкими.

СПИСОК ИСПОЛЬЗОВАННОЙ ЛИТЕРАТУРЫ

1.       Цифровая обработка сигналов: Учебн. Пособие для вузов/Л. М. Гольденберг, Б.Д. Матюшкин, М. Н. Поляк. – 2изд., перераб. и доп. – М.: Радио и связь, 1990. – 256 с.: ил.

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

3.       Макклеплан Д.Х., Рейдер Ч. М. Применение теории чисел в цифровой обработке сигналов: Пер. с англ. – М.: Радио и связь, 1983.-264с.

4.       Гольденберг Л.М., Матюшкин Б.Д., Поляк М. Н. Цифровая обработка сигналов: Справочник.- М.: Радио и связь, 1985. – 312с., ил.