Дискретное преобразование Фурье (ДПФ)

Содержание
Введение
Дискретное преобразование Фурье (ДПФ) — один из распространенных инструментов спектрального анализа сигналов, широко применяемый в самых разных отраслях науки и техники. При этом разработано множество быстрых алгоритмов для высокой вычислительной эффективности ДПФ.
В данном разделе будет уделено особое внимание переходу от непрерывного интеграла Фурье к дискретно-временному преобразованию Фурье (ДВПФ) и, далее, к дискретному преобразованию Фурье. Понимание данного перехода позволит лучше понять свойства ДПФ и сущность цифрового спектрального анализа в целом.
Пара непрерывного преобразования Фурье (интеграл Фурье) имеет вид:
$$ S(\omega) = \int \limits_{-\infty}^{\infty} s(t) \cdot \exp \left( -j \cdot \omega \cdot t\right) d t; \nonumber $$ $$ s(t) = \frac{1}{2 \pi} \cdot \int \limits_{-\infty}^{\infty} S(\omega) \cdot \exp \left( j \cdot \omega \cdot t\right) d \omega, $$
(1)
где $ S(\omega) $ — спектр сигнала $s(t)$ (в общем случае и сигнал и спектр — комплексные).
Выражения для прямого ДПФ и обратного дискретного преобразования Фурье (ОДПФ) имеют вид:
$$ S_{\textrm{d}}(k) = \sum_{n=0}^{N-1} s(n) \cdot \exp \left( -j \cdot \frac{2\pi}{N} \cdot n \cdot k \right),\text{ } k=0 \ldots N-1; $$ $$ s(n) = \frac{1}{N} \cdot \sum_{k=0}^{N-1} S_{\textrm{d}}(k) \cdot \exp \left( j \cdot \frac{2\pi}{N} \cdot n \cdot k\right), \text{ } n=0 \ldots N-1. $$
(2)
ДПФ ставит в соответствие $N$ отсчетам сигнала $s(n)$, $n = 0\ldots N-1$, $N$ отсчетов комплексного спектра $S_{\textrm{d}}(k)$, $k = 0 \ldots N-1$. Здесь и далее в данном разделе переменная $n$ индексирует временные отсчеты сигнала, а переменная $k$ индексирует спектральные отсчеты ДПФ.
Как в непрерывном, так и в дискретном случаях в выражениях для обратного преобразования имеется нормировочный коэффициент. В случае интеграла Фурье это $\frac{1}{2 \pi}$, в случае ОДПФ – $\frac{1}{N}$.
Нормировочный коэффициент необходим для корректного масштабирования сигнала из частотной области во временную. Нормировочный коэффициент уменьшает амплитуду сигнала на выходе обратного преобразования, для того чтобы она совпадала с амплитудой исходного сигнала. Если последовательно рассчитать прямое преобразование Фурье некоторого сигнала, а после взять обратное преобразование Фурье, то результат обратного преобразования должен полностью совпадать с исходным сигналом.
Дискретизация сигнала по времени. Дискретно-временное преобразование Фурье.
Рассмотрим дискретный сигнал $s_{\text{d}} (t)$ как результат умножения непрерывного сигнала $s(t)$ на решетчатую функцию:
$$ s_\text{d} (t) = s(t) \cdot \sum_{n=0}^{N-1} \delta(t - n \cdot \Delta t) = \sum_{n=0}^{N-1} s(t) \cdot \delta(t - n \cdot \Delta t), $$
(3)
где $\delta(t)$ – дельта-функция:
$$ \delta(t) = \begin{cases} \infty &\text{если $t = 0$;}\\ 0 &\text{если $t \neq 0$,} \end{cases} $$
(4)
$\Delta t$ – интервал дискретизации. Графически процесс дискретизации можно представить, как это показано на рисунке 1.
Рисунок 1. Процесс дискретизации сигнала
Вычислим преобразование Фурье дискретного сигнала $s_{\textrm{d}}(t)$, для этого подставим выражение (3) в выражение для преобразования Фурье (1):
$$ S_{\textrm{d}}(\omega) = \int \limits_{-\infty}^{\infty} s_{\text{d}}(t) \cdot \exp \left( -j \cdot \omega\ \cdot t \right) dt = \int \limits_{-\infty}^{\infty} \sum_{n=0}^{N-1} s(t) \cdot \delta(t - n \cdot \Delta t) \cdot \exp \left( -j \cdot \omega\ \cdot t \right) dt. $$
(5)
Поменяем местами операции суммирования и интегрирования и используем фильтрующее свойство дельта-функции:
$$ \int \limits_{-\infty}^{\infty}s(t) \cdot \delta(t-\tau) dt = s(\tau). $$
(6)
Тогда выражение (5) с учетом (6) принимает вид:
$$ S_{\textrm{d}}( \omega) = \sum_{n=0}^{N-1} {\int \limits_{-\infty}^{\infty} s(t) \cdot \delta(t - n \cdot \Delta t) \cdot \exp \left( -j \cdot \omega\ \cdot t \right) dt} = \sum_{n=0}^{N-1} s(n \cdot \Delta t) \cdot \exp \left( -j \cdot \omega\ \cdot n \cdot \Delta t \right). $$
(7)
Таким образом, мы избавились от интегрирования в бесконечных пределах, заменив его конечным суммированием комплексных экспонент.
Комплексные экспоненты $\exp \left( -j \cdot \omega\ \cdot n \cdot \Delta t \right)$ выражении (7) являются периодическими функциями с периодом:
$$ \Omega(n) = \frac{2 \pi}{n \cdot \Delta t} = 2 \pi \cdot n \cdot F_{\textrm{s}} , \text{ рад/с, } n = 1 \ldots N-1, $$
(8)
где $F_{\rm{s}} = \frac{1}{\Delta t}$ - частота дискретизации сигнала (Гц).
Необходимо отметить, что $n=0$ исключено из выражения (8), так как при $n=0$ комплексная экспонента равна единице. Максимальный период повторения спектра $S_{\textrm{d}}(\omega)$ будет при $n=1$, в этом случае он равен $\Omega_{\text{max}} = \Omega(1) = \frac{2 \pi}{\Delta t} = 2\pi \cdot F_{\rm{s}}$ рад/c.
Таким образом, спектр $S_{\textrm{d}}(\omega)$ дискретного сигнала $s_{\textrm{d}}(t)$, есть $2\pi \cdot F_{\textrm{s}}$ — периодическая функция циклической частоты $\omega$, определенная как (7). Если мы введем нормировку частоты дискретизации $F_{\textrm{s}} = 1$ Гц, то (7) переходит к выражению дискретно-временного преобразования Фурье (ДВПФ):
$$ S(\omega_{\textrm{н}}) = \sum_{n=0}^{N-1} s(n) \cdot \exp \left( -j \cdot \omega_{\textrm{н}} \cdot n \right). $$
(9)
ДВПФ использует только индексы отсчетов входного сигнала $s(n)$ при частоте дискретизации $F_{\textrm{s}} = 1$ Гц. В результате ДВПФ мы получим $2\pi$ периодическую функцию $S_{\textrm{d}}(\omega_{\textrm{н}})$ нормированной циклической частоты $\omega_{\textrm{н}} = \frac{\omega}{F_{\textrm{s}}}$.
Поскольку спектр дискретного сигнала — периодическая функция, то можно рассматривать только один период повторения спектра $S_{\textrm{d}}(\omega)$ при $\omega = \big[0, \: 2\pi \cdot F_{\textrm{s}} \big]$ рад/с или $S_{\textrm{d}}(f)$ при $f = \big[0, \: F_{\rm{s}} \big]$ Гц.

Повторение сигнала во времени. Дискретное преобразование Фурье
Для программной реализации алгоритмов цифровой обработки требуются как дискретные отсчеты сигнала, так и дискретные отсчеты спектра. Известно что дискретным, или, как еще говорят, линейчатым спектром, обладают периодические сигналы. При этом дискретный спектр получается путем разложения в ряд Фурье периодического сигнала. Значит, чтобы получить дискретный спектр, надо сделать исходный дискретный сигнал периодическим путем повторения данного сигнала во времени бесконечное количество раз с некоторым периодом $T$. Тогда спектр периодического сигнала будет содержать дискретные гармоники, кратные $\Delta \omega = \frac{2\pi}{T}$ рад/c.
Графически процесс повторения сигнала во времени представлен на рисунке 2.
Рисунок 2. Повторение сигнала во времени
Черным показан исходный сигнал, красным — его повторения через некоторый период $T$.
Повторять сигнал можно с различным периодом $T$, однако необходимо, чтобы период повторения был больше или равен длительности сигнала $T \ge N \cdot \Delta t$, чтобы сигнал и его периодические повторения не перекрывались во времени. При этом минимальный период повторения сигнала $T_{\text{min}}$, при котором сигнал и его повторения не накладываются друг на друга, равен
$$ T_{\text{min}} = N \cdot \Delta t, {\textrm{ сек}}. $$
(10)
Повторение сигнала с минимальным периодом $T_{\text{min}} = N \cdot \Delta t$ показано на рисунке 3.
Рисунок 3. Повторение сигнала с минимальным периодом
При повторении сигнала с минимальным периодом $T_{\text{min}}$ получим линейчатый спектр сигнала, состоящий из гармоник, кратных:
$$ \Delta \omega = \frac{2\pi}{T_{\text{min}}} = \frac{2\pi}{N \cdot \Delta t} \text{ (рад/с)}. $$
(11)
Таким образом, мы можем продискретизировать спектр $S_{\textrm{d}}(\omega)$ дискретного сигнала $s(n)$ на одном периоде повторения $2\pi \cdot F_{\textrm{s}}$ с шагом $\Delta \omega = \frac{2\pi}{N \cdot \Delta t}$ и получим
$$ \frac{ 2\pi \cdot F_{\textrm{s}} }{\Delta \omega} = \frac{ \frac{ 2 \pi }{ \Delta t } } { \frac{ 2\pi }{ N \cdot \Delta t} } = N $$
(12)
отсчетов спектра.
Учтем вышесказанное в выражении (7):
$$ S_{\textrm{d}}( k \cdot \Delta \omega) = \sum_{n=0}^{N-1} s(n \cdot \Delta t) \exp \left( -j \cdot n \cdot \Delta t \cdot k \cdot \Delta \omega \right) = $$ $$ = \sum_{n=0}^{N-1} s(n \cdot \Delta t) \exp \left( -j \cdot n \cdot \Delta t \cdot k \cdot \frac{2\pi}{N \cdot \Delta t} \right) = $$ $$ = \sum_{n=0}^{N-1} s(n \cdot \Delta t) \exp \left( -j \cdot \frac{2\pi}{N} \cdot n \cdot k \right), \text{ } k = 0 \ldots N-1. $$
(13)
Если опустить в выражении (13) шаг дискретизации по времени $\Delta t$ и по частоте $\Delta \omega$, то получим окончательное выражение для ДПФ:
$$ S_{\textrm{d}}(k)= \sum_{n=0}^{N-1} s(n) \cdot \exp \left( -j \cdot \frac{2\pi}{N} \cdot n \cdot k \right), \text{ } k = 0 \ldots N-1. $$
(14)
ДПФ ставит в соответствие $N$ отсчетам дискретного сигнала $s(n)$, $N$ отсчетов дискретного спектра $S_{\textrm{d}}(k)$, при этом предполагается, что и сигнал и спектр являются периодическими и анализируются на одном периоде повторения.
Детально свойства ДПФ будут рассмотрены в следующих разделах.
Обратное дискретное преобразование Фурье
Аналогично (3) можно записать выражение для дискретного спектра через решетчатую функцию:
$$ S_{\textrm{D}}(\omega) = \sum_{k=0}^{N-1}S_{\textrm{d}}(\omega) \cdot \delta(\omega - k \cdot \Delta \omega), $$
(15)
где $S_{\textrm{D}}(\omega)$ - дискретные отсчеты спектра $S_{\textrm{d}}(\omega)$ на одном периоде повторения $2 \pi \cdot F_{\textrm{s}}$ рад/с.
Подставим (15) в выражение для обратного преобразования Фурье (1):
$$ s(t) = C \cdot \int \limits_{-\infty}^{\infty}S_{\textrm{D}}(\omega) \cdot \exp(j \cdot \omega \cdot t) d \omega = C \cdot \int \limits_{-\infty}^{\infty}\sum_{k=0}^{N-1}S_{\textrm{d}}(\omega) \cdot \delta(\omega - k \cdot \Delta \omega) \cdot \exp(j \cdot \omega \cdot t) d \omega , $$
(16)
где $C$ — коэффициент пропорциональности, задача которого — обеспечить равенство по амплитуде исходного дискретного сигнала и результата ОДПФ. Коэффициент пропорциональности $C$ учитывает коэффициент $\frac{1}{2\pi}$ обратного преобразования Фурье.
Поменяем местами операции суммирования и интегрирования и учтем фильтрующее свойство дельта-функции:
$$ s(t)= C \cdot \sum_{k=0}^{N-1}\int \limits_{-\infty}^{\infty}S_{\textrm{d}}(\omega) \cdot \delta(\omega - k\cdot \Delta \omega) \cdot \exp(j \cdot \omega \cdot t) d \omega = C \cdot \sum_{k=0}^{N-1}S_{\textrm{d}}(k\cdot \Delta\omega) \cdot \exp(j \cdot k \cdot \Delta \omega \cdot t). $$
(17)
Возьмем дискретные отсчеты $s(t)$ через интервал $\Delta t = \frac{1}{F_{\textrm{s}}}$ сек, тогда (17) можно записать:
$$ s(n \cdot \Delta t) = C \cdot \sum_{k=0}^{N-1}S_{\textrm{d}}(k\cdot \Delta\omega) \cdot \exp(j \cdot k \cdot \Delta \omega \cdot n \cdot \Delta t), \text{ } n = 0 \ldots N-1. $$
(18)
Учтем (11) и получим:
$$ s(n \cdot \Delta t) = C \cdot \sum_{k=0}^{N-1}S_{\textrm{d}}(k\cdot \Delta\omega) \cdot \exp \left(j \cdot \frac{2\pi}{N} \cdot n \cdot k \right), \text{ } n = 0\ldots N-1. $$
(19)
Опустив в (19) интервалы дискретизации по частоте и по времени, получим выражение для ОДПФ, в котором остался неизвестным коэффициент пропорциональности $C$:
$$ s(n) = C \cdot \sum_{k=0}^{N-1}S_{\textrm{d}}(k) \cdot \exp \left(j \cdot \frac{2\pi}{N} \cdot n \cdot k \right), \text{ } n = 0 \ldots N-1. $$
(20)
Для того чтобы рассчитать коэффициент $C$, необходимо в выражение для ОДПФ (20) подставить выражение для ДПФ (14):
$$ s(n) = C \cdot \sum_{k=0}^{N-1} \sum_{m=0}^{N-1} s(m) \cdot \exp \left(-j \cdot \frac{2\pi}{N} \cdot m \cdot k \right) \cdot \exp \left(j \cdot \frac{2\pi}{N} \cdot n \cdot k \right), \text{ } n = 0 \ldots N-1. $$
(21)
Поменяем местами в (21) порядок суммирования и объединим экспоненты:
$$ s(n) = C \cdot \sum_{m=0}^{N-1} \cdot s(m) \cdot \sum_{k=0}^{N-1} \exp \left(j \cdot \frac{2\pi}{N} \cdot (n-m) \cdot k \right), \text{ } n = 0 \ldots N-1. $$
(22)
Рассмотрим подробнее сумму комплексных экспонент входящую в (22).
При $m=n$ имеем:
$$ \sum_{k=0}^{N-1} \exp \left(j \cdot \frac{2\pi}{N} \cdot (n-n) \cdot k \right) = \sum_{k=0}^{N-1} 1 = N. $$
(23)
Теперь рассмотрим эту же сумму при $m \neq n$.
Пусть $L=n-m$, тогда:
$$ \sum_{k=0}^{N-1} \exp \left(j \cdot \frac{2\pi}{N} \cdot (n-m) \cdot k \right) = \sum_{k=0}^{N-1} \exp \left(j \cdot \frac{2\pi}{N} \cdot L \cdot k \right). $$
(24)
Каждая комплексная экспонента, входящая в сумму (24), есть вектор на комплексной плоскости единичной длины, повернутый на угол:
$$ \phi(k) = \frac{2\pi \cdot L}{N} \cdot k, \text{ } k = 0 \ldots N-1. $$
(25)
Таких векторов будет $N$, и они повернуты относительно друг друга на одинаковые углы $\Delta\phi$.
Поскольку все векторы выходят из одной точки (из 0 комплексной плоскости) и повернуты относительно друг друга на одинаковые углы $\Delta\phi$, то их сумма равна нулю. Это иллюстрируется рисунком 4 для $N=4$.
Рисунок 4. Сумма комплексных экспонент
Из рисунка 4 можно заключить, что для всех $L=n-m$, $n=0\ldots N-1$, $m=0\ldots N-1$, $m \neq n$ суммы комплексных экспонент (24) равны нулю.
Тогда в выражении (22) от суммы по $m=0 \ldots N-1$ останется только слагаемое при $n=m$, и (22) можно представить в виде:
$$ s(n) = C \cdot \sum_{m=0}^{N-1} \cdot s(m) \cdot \sum_{k=0}^{N-1} \exp \left(j \cdot \frac{2\pi}{N} \cdot (n-m) \cdot k \right) = C \cdot s(n) \cdot N, \text{ } n = 0 \ldots N-1. $$
(26)
Из (26) можно заключить, что $C = \frac{1}{N}$.
Таким образом, окончательное выражение для ОДПФ имеет вид:
$$ s(n) = \frac{1}{N} \cdot \sum_{k=0}^{N-1} S_{\textrm{d}}(k) \cdot \exp \left( j \cdot \frac{2 \pi}{N} \cdot n \cdot k\right), \text{ } n = 0 \ldots N-1. $$
(27)
Выводы
В данном разделе мы осуществили переход от интеграла Фурье к прямому и обратному дискретному преобразованию Фурье путем дискретизации сигнала как по времени, так и по частоте.
Было показано, что спектр $S_{\textrm{d}}(\omega)$ дискретного сигнала - периодическая функция с периодом $2 \pi \cdot F_{\textrm{s}}$ рад/с или $F_{\textrm{s}}$ Гц.
Также было введено понятие дискретно-временного преобразования Фурье $S(\omega_ {\textrm{н}})$ как $2\pi$ периодической функции, нормированной к частоте дискретизации циклической частоты $\omega_{\textrm{н}}$.
ДПФ получено путем дискретизации непрерывной функции $S_{\textrm{d}}(\omega)$ на одном периоде повторения, что эквивалентно периодическому повторению дискретного сигнала во времени с периодом $\frac{N}{F_{\textrm{s}}}$ сек.

Ваши комментарии, вопросы и пожелания вы можете оставить на форуме или в гостевой книге.
Список литературы
[1] Сергиенко А.Б. Цифровая обработка сигналов Питер, 2002.

[2] Оппенгейм А., Шафер Р. Цифровая обработка сигналов Техносфера, 2006.

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

[4] Bracewell R.N. The Fourier Transform and its Applications. McGraw Hill, Singapor, 2000.

[5] Oppenheim Alan V. and Schafer Ronald W. Discrete-Time Signal Processing Second Edition. Prentice-Hall, New Jersey, 1999.

[6] Robert J. Marks II The Joy of Fourier: Analysis, Sampling Theory, Systems, Multidimensions, Stochastic Processes, Random Variables, Signal Recovery, Pocs, Time Scales, & Applications. Baylor University, 2006.

[7] Nussbaumer Henri J. Fast Fourier Transform and Convolution Algorithms. Second Corrected and Updated Edition. Springer-Verlag, 1982.

Oбнаружили ошибку в тексте? Выделите ее мышкой и нажмите