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

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

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

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

Вперед 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 может быть более эффективно оценен с использованием алгоритма скользящего DFT . ^[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|}$

(c) компактный носитель: равен нулю, когда |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 _s (и сохранение 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}}$

Преобразование Chirp 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}}}$

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

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

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

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

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

Ссылки

1. Sejdić E.; Djurović I.; Jiang J. (2009). «Представление частотно-временных характеристик с использованием концентрации энергии: обзор последних достижений». Цифровая обработка сигналов . 19 (1): 153–183. Bibcode : 2009DSP....19..153S. doi : 10.1016/j.dsp.2007.12.004.
2. Э. Якобсен и Р. Лайонс, Скользящее ДПФ, Журнал обработки сигналов , т. 20, выпуск 2, стр. 74–80 (март 2003 г.).
3. Джонт Б. Аллен (июнь 1977 г.). «Кратковременный спектральный анализ, синтез и модификация с помощью дискретного преобразования Фурье». Труды IEEE по акустике, речи и обработке сигналов . ASSP-25 (3): 235–238. doi :10.1109/TASSP.1977.1162950.
4. Кляйнфельд, Дэвид; Митра, Парта П. (март 2014 г.). «Спектральные методы функциональной визуализации мозга». Cold Spring Harbor Protocols . 2014 (3): pdb.top081075. doi :10.1101/pdb.top081075. PMID  24591695.
5. «Что означает «недостаточное заполнение для требуемого разрешения по частоте»? – FieldTrip toolbox».
6. Zeitler M, Fries P, Gielen S (2008). «Предвзятая конкуренция через вариации амплитуды гамма-колебаний». J Comput Neurosci . 25 (1): 89–107. doi :10.1007/s10827-007-0066-2. PMC 2441488 . PMID  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
• Растянутое во времени кратковременное преобразование Фурье для частотно-временного анализа сверхширокополосных сигналов
• Лицензированный BSD класс Matlab для выполнения STFT и обратного STFT
• LTFAT – бесплатный (GPL) набор инструментов Matlab/Octave для работы с кратковременными преобразованиями Фурье и частотно-временным анализом
• Сонограмма видимой речи – Бесплатное (GPL) программное обеспечение для кратковременных преобразований Фурье и частотно-временного анализа
• Национальный университет Тайваня, частотно-временной анализ и вейвлет-преобразование 2021 г., профессор кафедры электротехники Цзянь-Цзюнь Дин