Оглавление.

Компьютерное распознавание и порождение речи

Глава 2. Основы цифровой обработки звуковых сигналов.

2.3. Преобразование Фурье.

Теперь настало время более подробно поговорить о разложении Фурье, дискретном преобразовании Фурье и алгоритмах его быстрого вычисления.

В 1807 году французский математик и физик Жан-Батист Жозеф Фурье выдвинул идею о том, что любой непрерывный периодический сигнал может быть представлен как сумма должным образом выбранных синусоидальных волн. И хотя он оказался прав не на все сто процентов, ибо суммированием синусоид нельзя получить сигнал с уголковыми вершинами, все же подобным разложением можно получить очень хорошо приближенные значения даже для таких сигналов, а уже тем паче для других.

Нас будет интересовать частный случай подобного разложения, а именно дискретное преобразование Фурье (ДПФ или в английском варианте FFT – Fast Fourier Transformation), потому что мы имеем дело не с непрерывным (аналоговым) сигналом, а с дискретным сигналом, представляющим собой набор цифровых отсчетов, следующих с определенной частотой дискретизации.

На вход ДПФ подается сигнал, представляющий собой последовательность амплитуд в N следующих друг за другом с одинаковыми достаточно короткими промежутками моментах времени. На выходе ДПФ получается последовательность N/2+1 амплитуд следующих друг за другом частот синуса и последовательность N/2+1 амплитуд следующих друг за другом частот косинуса. Принято говорить, что входной сигнал находится в домене времени, а оба выходных сигнала – в домене частоты (или частотном домене).

Сигналы домена времени обозначают строчными буквами : x[], y[] и т.д. Сигналы же частотного домена обозначают большими буквами соответствующими сигналу домена времени. То есть сигналу домена времени x[] соответствует сигнал частотного домена X[], а сигналу y[] – сигнал Y[]. Сигнал частотного домена - например X[]- состоит из двух выборок длиной N/2+1. Они называются вещественной частью из X[] , которая обозначается как ReX[], и мнимой частью из X[] , которая обозначается как ImX[]. Значения ReX[]- амплитуды волн косинуса, а значения ImX[] – амплитуды волн синуса. Если сигналы домена времени нумеруются от x[0] до x[N-1], то частотного домена – от Re[0] до Re[N/2] и от Im[0] до Im[N/2].

Используемые в ДПФ волны синуса и косинуса называются базисными функциями ДПФ. Базисные функции – набор волн синуса и косинуса с амплитудой равной единице. На выходе ДПФ мы получаем набор амплитуд, которые мы должны приложить к определенной базисной функции - волне косинуса или синуса, - и получить в результате масштабируемые волны косинуса и синуса, которые при сложении дают сигнал времени.

Уравнения для базисных функций таковы:

ck[i]=cos(2πki/N) волна косинуса,

sk[i]=sin(2πki/N) волна синуса,

для точек i от 0 до N-1 домена времени, для точек k от 0 до N/2 домена частоты. То есть точке k в частотном домене соответствует частота 2πki/N. Таким образом дискретное преобразование Фурье состоит в нахождении таких точек X[i], что

	N/2		N/2
x[i]=	∑	ReX[i]cos(2πki/N) +	∑	ImX[i]sin(2πki/N)
	k=0		k=0

Эта формула называется обратным преобразованием Фурье. Для нахождения же самих коэффициентов X[i] применяется формула прямого преобразования Фурье:

	N-1
ReX [k]=	∑	x[i] cos(2πki/N)
	i=0

	N-1
ImX [k]=	∑	x[i] sin(2πki/N)
	i=0

Или в общем виде:

	N-1		N-1
X[k]=	∑	x[i]cos(2πki/N) +j	∑	x[i]sin(2πki/N)
	i=0		i=0

Здесь j- мнимая единица

По формуле Эйлера e^jx=cosx+jsinx.

Отсюда мы получаем комплексное представление прямого ДПФ:

	N-1
X[k]=	∑	x[i]e-j2π ki/N		(1)
	i=0

Обратное ДПФ будет вычисляться по формуле

	N-1
x[i]=	∑	X[k]e-j2π ki/N		(2)
	i=0

Величина

W_N^k=exp{-j2πk/N} называется поворачивающим множителем.

Тогда прямое преобразование Фурье:

	N-1
X[k]=	∑	x[i]WNk		(3)
	i=0

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

Для любознательных рассмотрим подробно один из таких алгоритмов. Остальные же могут спокойно сразу переходить к следующему параграфу и использовать при программировании готовый алгоритм, не вдаваясь в то, как он работает. Готовый алгоритм можно скачать в Интернете, например п адресу http:\\psi-logic.narod.ru.

Из приведенной ранее формулы Эйлера e^jx=cosx+jsinx легко вывести равенство

e^j2πN = 1 (4) для любого целого N.

Пользуясь этой формулой можно вывести следующие свойства поворачивающего множителя:

1) W_N^k периодичен по k и n с периодом N. (5)

2) W_N^k =- W_N^k-N/2 (6)

Далее определим две последовательности: {x_чет} и {x_неч} через последовательность {x}следующим образом:

x_чет[i] = x[2i], x_неч[i] = x[2i+1], где i=(1,N/2-1)

Пусть к этим последовательностям применено дискретное преобразование Фурье, обозначим его результаты как последовательности {X_чет} и {X_неч} по N/2 элементов каждая.

Пользуясь свойствами поворачивающего множителя можно доказать, что

X[k]=X_чет[k]+ W_N^kX_неч[k], при k=(0,N/2-1) (8)

X[k]=X_чет[k-N/2]+ W_N^k-N/2X_неч[k-N/2], при k=(N/2,N-1)

Далее из уравнений (8) можно легко получить:

X[k]=X_чет[k]+ W_N^kX_неч[k], при k=(0,N/2-1) (9)

X[N/2+k]=X_чет[k]- W_N^kX_неч[k], при k=(0,N/2-1)

Формула (9) позволяет сократить число умножений вдвое, если вычислять X[k] не последовательно, а попарно: X[0] и Х[N/2], X[1] и Х[N/2+1] и так далее. Также отпадает необходимость хранить вычисленные значения X_чет[k] и X_неч[k] после использования при вычислении очередной пары и одно вычисление W_N^k можно использовать для расчета двух элементов последовательности {X}. В результате этих нехитрых манипуляций скорость время вычисления преобразования Фурье существенно сокращается.

Рассмотрим БПФ для различных N. Введем дополнительные обозначения.

X_{R}[k]– это k-й элемент последовательности из R элементов,

X_{R/2}чет - k-й элемент подпоследовательности элементов с четными номерами из последовательности R элементов,

X_{R/2}неч- k-й элемент подпоследовательности элементов с нечетными номерами из последовательности R элементов.

Для N=1:

	N-1
X{1}[k] =	∑	x{1}[i]WNk = x{1}[0]e0 = x{1}[0],

то есть ДПФ над одним числом дает само это число.

Для N=2, пользуясь (9) получим:

X_{2}[0] = X_{1}чет[0]+X_{1}неч[0]W₄⁰ = x[0]+x[1],

X{2}[1] = X_{1}чет[0]-X_{1}неч[0]W₄⁰ = x[0]-x[1]

Для N=4:

X{4}[0] = X_{2}чет[0]+X_{2}неч[0]W₄⁰ = x[0]+x[2]+ x[1]+x[3],

X{4}[2] = X_{2}чет[0]-X_{2}неч[0]W₄⁰ = x[0]+x[2]- x[1]-x[3],

X{4}[1] = X_{2}чет[1]+X_{2}неч[1]W₂¹ = x[0]-x[2]-j(x[1]-x[3]),

X{4}[3] = X_{2}чет[1]-X_{2}неч[1]W₂¹ = x[0]-x[2]+j(x[1]-x[3])

И так далее. При вычислении следует учесть, что

W_N=exp{-j2π/N}= exp{-j2πL/NL}= W_NL^L,
то есть, например вычислив один раз W256 , мы можем получить:

W₁₂₈= W₂₅₆²
W₆₄= W₂₅₆⁴ = W₁₂₈²
W₃₂= W₂₅₆⁸ = W₆₄²
W₁₆= W₂₅₆¹⁶ = W₃₂²
W₈= W₂₅₆³² = W₁₆²
W₄= W₂₅₆⁶⁴ = W₈²
W₂= W₂₅₆¹²⁸ = W₄² ,

то есть мы избавлены от необходимости полного вычисления поворотного множителя на каждом уровне.

Теперь имея полный математический аппарат перейдем к описанию самого алгоритма БПФ.

Рассмотрим случай для N=2^T, где T – целое.

Сначала разделим входящую последовательность x на две подпоследовательности – с четными и нечетными номерами. Одну помещаем слева, другую – справа. Затем перенумеровываем каждую из подпоследовательностей и снова разделяем каждую из них на две подпоследовательности: с четными и нечетными номерами и помещаем одну подпоследовательность слева, а другую справа ( в рамках «материнской» подпоследовательности). Эту своеобразную сортировку продолжаем до тех пор пока размер каждой подпоследовательности не достигнет 2 элементов. Это достигается ровно в T этапов. В результате( что можно обосновать математически) номер каждого исходного элемента в двоичной системе переворачивается. То есть, например для однобайтных номеров двоичный номер 00100001 станет номером 10000100, номер 10100100 – номером 00100101 и так далее. На этом и основан применяемый здесь алгоритм сортировки: мы имеем заранее подготовленный массив из 256 элементов, в котором для номера i содержится элемент со значением равным перевернутому в двоичной системе i.. Для, например, четырехбайтных чисел мы заносим каждый байт в однобайтную переменную, выбираем для них из массива четыре однобайтных значения перевернутых номеров и формируем из них новый четырехбайтный номер: перевернутый байт номер 0 становится байтом номер 3, байт номер 1 – байтом номер 2, байт номер 2 – байтом номер 1, байт номер 3 – байтом номер 0.

Описанная выше сортировка подготавливает нас к основному алгоритму БПФ.

На рисунке (ниже) изображен порядок использования отсортированных элементов на втором этапе для случая N=8. Для других N порядок строится по аналогии. Алгоритм будет состоять из трех циклов, вложенных друг в друга. Внешний цикл состоит из T итераций (если N=2T) , соответствующих количеству уровней обработки. Средний цикл содержит количество итераций (переменная k), равное количеству получающихся на данном уровне подпоследовательностей. Здесь же вычисляется WNk. Количество итераций самого внутреннего цикла соответствует количеству единичных элементов, входящих в результирующую подпоследовательность. Здесь же элементы результирующей подпоследовательности вычисляются по формулам (9):

X[k]=X_чет[k]+ W_N^kX_неч[k]

X[N/2+k]=X_чет[k]- W_N^kX_неч[k]

При работе с комплексными числами учитываем, что если x₁ и x₂ – комплексные числа, такие что x₁=a₁+jb₁, x2=a2+jb2, S-действительное число, то

x₁+x₂=(a₁+a₂)+j(b₁+b₂),

x₁x₂=(a₁a₂-b₁b₂)+j(a₁b₂+b₁a₂)

exp{-jS}=cosS-jsinS.

Распознавание речи. Обработка текстов на естественном языке.