Кратковременное преобразование Фурье

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

Кратковременное преобразование Фурье ( STFT ) — это преобразование Фурье, используемое для определения синусоидальной частоты и фазового состава локальных участков сигнала по мере его изменения во времени. ^[1] На практике процедура вычисления STFT заключается в разделении более длительного сигнала на более короткие сегменты одинаковой длины, а затем вычислении преобразования Фурье отдельно для каждого более короткого сегмента. Это раскрывает спектр Фурье на каждом более коротком сегменте. Затем обычно строят изменяющиеся спектры как функцию времени, известные как спектрограммы или водопадные графики, например, обычно используемые в дисплеях спектра на основе программно определяемой радиосвязи (SDR). В дисплеях с полной полосой пропускания, охватывающих весь диапазон SDR, на настольных компьютерах обычно используются быстрые преобразования Фурье (БПФ) с 2 ^ 24 точками.

Спектрограмма, визуализирующая результаты STFT слов «девятнадцатый век». Здесь частоты показаны возрастающими по вертикальной оси, а время по горизонтальной оси. Легенда справа показывает, что интенсивность цвета увеличивается с увеличением плотности.

Форвард СТФТ

STFT непрерывного времени

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

${\displaystyle \mathbf {STFT} \{x(t)\}(\tau ,\omega )\equiv X(\tau ,\omega )=\int _{-\infty }^{\infty }x(t)w(t-\tau )e^{-i\omega t}\,dt}$

где – оконная функция , обычно окно Ханна или окно Гаусса с центром вокруг нуля, и – сигнал, который необходимо преобразовать (обратите внимание на разницу между оконной функцией и частотой ). по сути, является преобразованием Фурье , сложной функцией , представляющей фазу и величину сигнала в зависимости от времени и частоты. Часто развертка фазы используется вдоль одной или обеих осей времени и оси частоты, чтобы подавить любой скачок скачка результата фазы STFT. Индекс времени обычно считается « медленным » временем и обычно не выражается в таком высоком разрешении, как время . Учитывая, что STFT по сути представляет собой преобразование Фурье, умноженное на оконную функцию, STFT также называется оконным преобразованием Фурье или зависящим от времени преобразованием Фурье. ${\displaystyle w(\tau )}$ ${\displaystyle x(t)}$ ${\displaystyle w}$ ${\displaystyle \omega }$ ${\displaystyle X(\tau ,\omega )}$ ${\displaystyle x(t)w(t-\tau )}$ ${\displaystyle \tau }$ ${\displaystyle \omega }$ ${\displaystyle \tau }$ ${\displaystyle t}$

STFT дискретного времени

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

${\displaystyle \mathbf {STFT} \{x[n]\}(m,\omega )\equiv X(m,\omega )=\sum _{n=-\infty }^{\infty }x[n]w[n-m]e^{-i\omega n}}$

аналогично с сигналом и окном . В этом случае m дискретна, а ω непрерывна, но в большинстве типичных приложений STFT выполняется на компьютере с использованием быстрого преобразования Фурье , поэтому обе переменные дискретны и квантованы . ${\displaystyle x[n]}$ ${\displaystyle w[n]}$

Квадрат величины STFT дает представление спектрограммы спектральной плотности мощности функции:

${\displaystyle \operatorname {spectrogram} \{x(t)\}(\tau ,\omega )\equiv |X(\tau ,\omega )|^{2}}$

См. также модифицированное дискретное косинусное преобразование (MDCT), которое также является преобразованием Фурье, использующим перекрывающиеся окна.

Скользящая ДПФ

Если требуется только небольшое количество ω или если желательно, чтобы STFT оценивалось для каждого сдвига m окна, тогда STFT может быть более эффективно оценен с использованием алгоритма скользящего ДПФ . ^[2]

Обратный STFT

STFT является обратимым , то есть исходный сигнал может быть восстановлен из преобразования с помощью обратного STFT. Наиболее широко распространенный способ инвертирования STFT — использование метода сложения-перекрытия (OLA) , который также позволяет вносить изменения в комплексный спектр STFT. Это создает универсальный метод обработки сигналов, ^[3] называемый методом перекрытия и добавления с модификациями .

STFT непрерывного времени

Учитывая ширину и определение оконной функции w ( t ), мы изначально требуем, чтобы область оконной функции была масштабирована так, чтобы

${\displaystyle \int _{-\infty }^{\infty }w(\tau )\,d\tau =1.}$

Отсюда легко следует, что

${\displaystyle \int _{-\infty }^{\infty }w(t-\tau )\,d\tau =1\quad \forall \ t}$

и

${\displaystyle x(t)=x(t)\int _{-\infty }^{\infty }w(t-\tau )\,d\tau =\int _{-\infty }^{\infty }x(t)w(t-\tau )\,d\tau .}$

Непрерывное преобразование Фурье

${\displaystyle X(\omega )=\int _{-\infty }^{\infty }x(t)e^{-i\omega t}\,dt.}$

Подставляя x ( t ) сверху:

${\displaystyle X(\omega )=\int _{-\infty }^{\infty }\left[\int _{-\infty }^{\infty }x(t)w(t-\tau )\,d\tau \right]\,e^{-i\omega t}\,dt}$

${\displaystyle =\int _{-\infty }^{\infty }\int _{-\infty }^{\infty }x(t)w(t-\tau )\,e^{-i\omega t}\,d\tau \,dt.}$

Изменение порядка интеграции:

${\displaystyle X(\omega )=\int _{-\infty }^{\infty }\int _{-\infty }^{\infty }x(t)w(t-\tau )\,e^{-i\omega t}\,dt\,d\tau }$

${\displaystyle =\int _{-\infty }^{\infty }\left[\int _{-\infty }^{\infty }x(t)w(t-\tau )\,e^{-i\omega t}\,dt\right]\,d\tau }$

${\displaystyle =\int _{-\infty }^{\infty }X(\tau ,\omega )\,d\tau .}$

Таким образом, преобразование Фурье можно рассматривать как своего рода фазово-когерентную сумму всех STFT x ( t ). Поскольку обратное преобразование Фурье

${\displaystyle x(t)={\frac {1}{2\pi }}\int _{-\infty }^{\infty }X(\omega )e^{+i\omega t}\,d\omega ,}$

тогда x ( t ) можно восстановить из X (τ,ω) как

${\displaystyle x(t)={\frac {1}{2\pi }}\int _{-\infty }^{\infty }\int _{-\infty }^{\infty }X(\tau ,\omega )e^{+i\omega t}\,d\tau \,d\omega .}$

или

${\displaystyle x(t)=\int _{-\infty }^{\infty }\left[{\frac {1}{2\pi }}\int _{-\infty }^{\infty }X(\tau ,\omega )e^{+i\omega t}\,d\omega \right]\,d\tau .}$

По сравнению с приведенным выше можно увидеть, что оконная «зернистость» или «вейвлет» x ( t )

${\displaystyle x(t)w(t-\tau )={\frac {1}{2\pi }}\int _{-\infty }^{\infty }X(\tau ,\omega )e^{+i\omega t}\,d\omega .}$

обратное преобразование Фурье X (τ,ω) при фиксированном τ.

Альтернативное определение обратного преобразования, действительное только вблизи τ:

${\displaystyle x(t)={\frac {1}{w(t-\tau )}}{\frac {1}{2\pi }}\int _{-\infty }^{\infty }X(\tau ,\omega )e^{+i\omega t}\,d\omega .}$

В общем случае оконная функция имеет следующие свойства: ${\displaystyle w(t)}$

а) четная симметрия: ; ${\displaystyle w(t)=w(-t)}$

(б) невозрастающая (для положительного времени): если ; ${\displaystyle w(t)\geq w(s)}$ ${\displaystyle |t|\leq |s|}$

(в) компактный носитель: равен нулю, когда |t| большой. ${\displaystyle w(t)}$

Проблемы с разрешением

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

Сравнение разрешения STFT. Левый имеет лучшее временное разрешение, а правый — лучшее частотное разрешение.

Это одна из причин создания вейвлет-преобразования и многоразрешительного анализа , которые могут дать хорошее временное разрешение для высокочастотных событий и хорошее частотное разрешение для низкочастотных событий. Комбинация лучше всего подходит для многих реальных сигналов.

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

Можно рассматривать STFT для изменения размера окна как двумерную область (время и частота), как показано в примере ниже, которую можно рассчитать, изменяя размер окна. Однако это уже не строго частотно-временное представление — ядро ​​не является постоянным на протяжении всего сигнала.

Примеры

Когда исходная функция:

${\displaystyle X(t,f)=\int _{-\infty }^{\infty }w(t-\tau )x(\tau )e^{-j2\pi f\tau }d\tau }$

У нас может быть простой пример:

w(t) = 1 для |t| меньше или равно B

w(t) = 0 в противном случае

Б = окно

Теперь исходную функцию кратковременного преобразования Фурье можно изменить следующим образом:

${\displaystyle X(t,f)=\int _{t-B}^{t+B}x(\tau )e^{-j2\pi f\tau }d\tau }$

Другой пример:

Используя следующий пример сигнала , который состоит из набора из четырех синусоидальных сигналов, соединенных последовательно. Каждая форма волны состоит только из одной из четырех частот (10, 25, 50, 100 Гц ). Определение : ${\displaystyle x(t)}$ ${\displaystyle x(t)}$

${\displaystyle x(t)={\begin{cases}\cos(2\pi 10t)&0\,\mathrm {s} \leq t<5\,\mathrm {s} \\\cos(2\pi 25t)&5\,\mathrm {s} \leq t<10\,\mathrm {s} \\\cos(2\pi 50t)&10\,\mathrm {s} \leq t<15\,\mathrm {s} \\\cos(2\pi 100t)&15\,\mathrm {s} \leq t<20\,\mathrm {s} \\\end{cases}}}$

Затем он дискретизируется с частотой 400 Гц. Были получены следующие спектрограммы:

Окно в 25 мс позволяет нам определить точное время изменения сигналов, но точные частоты определить сложно. На другом конце шкалы окно в 1000 мс позволяет точно видеть частоты, но время между изменениями частот размыто.

Другие примеры:

${\displaystyle w(t)=exp(\sigma -t^{2})}$

Обычно мы вызываем функцию Гаусса или функцию Габора. Когда мы его используем, кратковременное преобразование Фурье называется «преобразованием Габора». ${\displaystyle exp(\sigma -t^{2})}$

Объяснение

Это также можно объяснить с помощью дискретизации и частоты Найквиста .

Возьмите окно N выборок из произвольного сигнала с действительным знаком с частотой дискретизации f _s . Преобразование Фурье дает N комплексных коэффициентов. Из этих коэффициентов полезна только половина (последний N/2 является комплексно-сопряженным первым N/2 в обратном порядке, поскольку это сигнал с действительным значением).

Эти коэффициенты N/2 представляют частоты от 0 до f _s /2 (Найквист), а два последовательных коэффициента разнесены на f _s / N Гц.

Для увеличения частотного разрешения окна необходимо уменьшить частотный интервал коэффициентов. Есть только две переменные, но уменьшение f ₍ и сохранение N постоянным) приведет к увеличению размера окна, поскольку теперь в единицу времени поступает меньше выборок. Другая альтернатива — увеличить N , но это снова приведет к увеличению размера окна. Таким образом, любая попытка увеличить разрешение по частоте приводит к увеличению размера окна и, следовательно, к уменьшению разрешения по времени — и наоборот.

Частота Рэлея

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

Частота Рэлея — это минимальная частота, которую можно разрешить с помощью временного окна конечной продолжительности. ^{[4] [5]}

Учитывая временное окно длиной Т секунд, минимальная частота, которую можно разрешить, составляет 1/Т Гц.

Частота Рэлея является важным фактором при применении кратковременного преобразования Фурье (STFT), а также любого другого метода гармонического анализа сигнала конечной длины записи. ^{[6] [7]}

Приложение

STFT используется для анализа аудиосигнала во времени.

STFT, а также стандартные преобразования Фурье и другие инструменты часто используются для анализа музыки. Спектрограмма может, например, отображать частоту на горизонтальной оси: самые низкие частоты слева , а самые высокие справа. Высота каждой полосы (дополненной цветом) представляет амплитуду частот в этом диапазоне. Измерение глубины представляет собой время, когда каждый новый столбец представлял собой отдельное преобразование. Аудиоинженеры используют этот вид визуальных эффектов для получения информации об аудиообразце, например, для определения частот определенных шумов (особенно при использовании с более высоким частотным разрешением) или для поиска частот, которые могут быть более или менее резонансными в пространстве, где находится сигнал был записан. Эту информацию можно использовать для эквалайзера или настройки других звуковых эффектов.

Выполнение

Оригинальная функция

${\displaystyle X(t,f)=\int _{-\infty }^{\infty }w(t-\tau )x(\tau )e^{-j2\pi f\tau }d\tau }$

Преобразование в дискретную форму:

${\displaystyle t=n\Delta _{t},f=m\Delta _{f},\tau =p\Delta _{t}}$

${\displaystyle X(n\Delta _{t},m\Delta _{f})=\sum _{-\infty }^{\infty }w((n-p)\Delta _{t})x(p\Delta _{t})e^{-j2\pi pm\Delta _{t}\Delta _{f}}\Delta _{t}}$

Предположим, что

${\displaystyle w(t)\cong 0{\text{ for }}|t|>B,{\frac {B}{\Delta _{t}}}=Q}$

Тогда мы можем записать исходную функцию в

${\displaystyle X(n\Delta _{t},m\Delta _{f})=\sum _{p=n-Q}^{n+Q}w((n-p)\Delta _{t})x(p\Delta _{t})e^{-j2\pi pm\Delta _{t}\Delta _{f}}\Delta _{t}}$

Прямая реализация

Ограничения

а. Критерий Найквиста (устранение эффекта наложения):

${\displaystyle \Delta _{t}<{\frac {1}{2\Omega }}}$, где полоса пропускания ${\displaystyle \Omega }$ ${\displaystyle x(\tau )w(t-\tau )}$

Метод на основе БПФ

Ограничение

а. , где целое число ${\displaystyle \Delta _{t}\Delta _{f}={\tfrac {1}{N}}}$ ${\displaystyle N}$

б. ${\displaystyle N\geq 2Q+1}$

в. Критерий Найквиста (устранение эффекта наложения):

${\displaystyle \Delta _{t}<{\frac {1}{2\Omega }}}$, это полоса пропускания ${\displaystyle \Omega }$ ${\displaystyle x(\tau )w(t-\tau )}$

${\displaystyle X(n\Delta _{t},m\Delta _{f})=\sum _{p=n-Q}^{n+Q}w((n-p)\Delta _{t})x(p\Delta _{t})e^{-{\frac {2\pi jpm}{N}}}\Delta _{t}}$

${\displaystyle {\text{if we have }}q=p-(n-Q),{\text{ then }}p=(n-Q)+q}$

${\displaystyle X(n\Delta _{t},m\Delta _{f})=\Delta _{t}e^{\frac {2\pi j(Q-n)m}{N}}\sum _{q=0}^{N-1}x_{1}(q)e^{-{\frac {2\pi jqm}{N}}}}$

${\displaystyle {\text{where }}x_{1}(q)={\begin{cases}w((Q-q)\Delta _{t})x((n-Q+q)\Delta _{t})&0\leq q\leq 2Q\\0&2Q<q<N\end{cases}}}$

Рекурсивный метод

Ограничение

а. , где целое число ${\displaystyle \Delta _{t}\Delta _{f}={\tfrac {1}{N}}}$ ${\displaystyle N}$

б. ${\displaystyle N\geq 2Q+1}$

в. Критерий Найквиста (устранение эффекта наложения):

${\displaystyle \Delta _{t}<{\frac {1}{2\Omega }}}$, это полоса пропускания ${\displaystyle \Omega }$ ${\displaystyle x(\tau )w(t-\tau )}$

д. Только для реализации прямоугольного STFT

Прямоугольное окно накладывает ограничение

${\displaystyle w((n-p)\Delta _{t})=1}$

Замена дает:

${\displaystyle {\begin{aligned}X(n\Delta _{t},m\Delta _{f})&=\sum _{p=n-Q}^{n+Q}w((n-p)\Delta _{t})&x(p\Delta _{t})e^{-{\frac {j2\pi pm}{N}}}\Delta _{t}\\&=\sum _{p=n-Q}^{n+Q}&x(p\Delta _{t})e^{-{\frac {j2\pi pm}{N}}}\Delta _{t}\\\end{aligned}}}$

Замена переменной n -1 на n :

${\displaystyle X((n-1)\Delta _{t},m\Delta _{f})=\sum _{p=n-1-Q}^{n-1+Q}x(p\Delta _{t})e^{-{\frac {j2\pi pm}{N}}}\Delta _{t}}$

Рассчитаем с помощью N -точечного БПФ: ${\displaystyle X(\min {n}\Delta _{t},m\Delta _{f})}$

${\displaystyle X(n_{0}\Delta _{t},m\Delta _{f})=\Delta _{t}e^{\frac {j2\pi (Q-n_{0})m}{N}}\sum _{q=0}^{N-1}x_{1}(q)e^{-j{\frac {2\pi qm}{N}}},\qquad n_{0}=\min {(n)}}$

где

${\displaystyle x_{1}(q)={\begin{cases}x((n-Q+q)\Delta _{t})&q\leq 2Q\\0&q>2Q\end{cases}}}$

Применение рекурсивной формулы для расчета ${\displaystyle X(n\Delta _{t},m\Delta _{f})}$

${\displaystyle X(n\Delta _{t},m\Delta _{f})=X((n-1)\Delta _{t},m\Delta _{f})-x((n-Q-1)\Delta _{t})e^{-{\frac {j2\pi (n-Q-1)m}{N}}}\Delta _{t}+x((n+Q)\Delta _{t})e^{-{\frac {j2\pi (n+Q)m}{N}}}\Delta _{t}}$

Чирп-Z-преобразование

Ограничение

${\displaystyle \exp {(-j2\pi pm\Delta _{t}\Delta _{f})}=\exp {(-j\pi p^{2}\Delta _{t}\Delta _{f})}\cdot \exp {(j\pi (p-m)^{2}\Delta _{t}\Delta _{f})}\cdot \exp {(-j\pi m^{2}\Delta _{t}\Delta _{f})}}$

так

${\displaystyle X(n\Delta _{t},m\Delta _{f})=\Delta _{t}\sum _{p=n-Q}^{n+Q}w((n-p)\Delta _{t})x(p\Delta _{t})e^{-j2\pi pm\Delta _{t}\Delta _{f}}}$

${\displaystyle X(n\Delta _{t},m\Delta _{f})=\Delta _{t}e^{-j2\pi m^{2}\Delta _{t}\Delta _{f}}\sum _{p=n-Q}^{n+Q}w((n-p)\Delta _{t})x(p\Delta _{t})e^{-j\pi p^{2}\Delta _{t}\Delta _{f}}e^{j\pi (p-m)^{2}\Delta _{t}\Delta _{f}}}$

Сравнение реализации

Смотрите также

• Спектральный анализ методом наименьших квадратов
• Оценка спектральной плотности
• Частотно-временной анализ
• Частотно-временное представление
• Метод переназначения

Другие частотно-временные преобразования:

• Функция распределения конической формы
• Преобразование с постоянной Q
• Дробное преобразование Фурье
• Преобразование Габора
• Преобразование Ньюленда
• S-преобразование
• Вейвлет-преобразование
• Преобразование чирплета

Рекомендации

1. Сейдич Э.; Джурович И.; Цзян Дж. (2009). «Представление частотно-временных характеристик с использованием концентрации энергии: обзор последних достижений». Цифровая обработка сигналов . 19 (1): 153–183. дои : 10.1016/j.dsp.2007.12.004.
2. Э. Якобсен и Р. Лайонс, Скользящая ДПФ, журнал Signal Processing Magazine, том. 20, выпуск 2, стр. 74–80 (март 2003 г.).
3. Джонт Б. Аллен (июнь 1977 г.). «Кратковременный спектральный анализ, синтез и модификация с помощью дискретного преобразования Фурье». Транзакции IEEE по акустике, речи и обработке сигналов . АССП-25 (3): 235–238. дои :10.1109/ТАССП.1977.1162950.
4. https://physical.ucsd.edu/neurophysical/publications/Cold%20Spring%20Harb%20Protoc-2014-Kleinfeld-pdb.top081075.pdf ^{[ пустой URL-адрес PDF ]}
5. «Что означает «заполнение недостаточно для требуемого частотного разрешения»? - Набор инструментов FieldTrip» .
6. Цайтлер М., Фрис П., Гилен С. (2008). «Предвзятая конкуренция за счет изменения амплитуды гамма-колебаний». J Вычислительная Неврология . 25 (1): 89–107. дои : 10.1007/s10827-007-0066-2. ПМЦ 2441488 . ПМИД  18293071. 
7. Вингерден, Марин ван; Винк, Мартин; Ланкельма, Ян; Пеннарц, Сириэль М.А. (19 мая 2010 г.). «Фазовая блокировка тета-диапазона орбитофронтальных нейронов во время ожидания вознаграждения». Журнал неврологии . 30 (20): 7078–7087. doi : 10.1523/JNEUROSCI.3860-09.2010. ISSN  0270-6474. ПМЦ 6632657 . ПМИД  20484650. 

Внешние ссылки

• DiscreteTFDs - программное обеспечение для расчета кратковременного преобразования Фурье и других частотно-временных распределений.
• Singular Spectral Analysis – MultiTaper Method Toolkit – бесплатная программа для анализа коротких, зашумленных временных рядов.
• Набор инструментов kSpectra для Mac OS X от SpectraWorks
• Растянутое по времени кратковременное преобразование Фурье для частотно-временного анализа сверхширокополосных сигналов
• Класс Matlab с лицензией BSD для выполнения STFT и обратного STFT.
• LTFAT — бесплатный (GPL) набор инструментов Matlab/Octave для работы с кратковременными преобразованиями Фурье и частотно-временным анализом.
• Видимая речь сонограммы - бесплатное (GPL) бесплатное программное обеспечение для кратковременных преобразований Фурье и частотно-временного анализа.
• Национальный Тайваньский университет, частотно-временной анализ и вейвлет-преобразование 2021, профессор Цзянь-Цзюнь Дин, факультет электротехники