Что получается в результате преобразования фурье
Практическое применение преобразования Фурье для обработки сигналов
Книги и публикации по цифровой обработке сигналов пишут авторы зачастую не догадывающиеся и не понимающие задач, стоящих перед разработчиками. Особенно это касается систем, работающих в реальном времени. Эти авторы отводят себе скромную роль бога, существующего вне времени и пространства, что вызывает некоторое недоумение у читателей подобной литературы. Данная публикация имеет целью развеять недоумения, возникающие у большинства разработчиков, и помочь им преодолеть «порог вхождения», для этих целей в тексте сознательно используется аналогии и терминология сферы программирования.
Данный опус не претендует на полноту и связность изложения.
Добавлено после прочтения комментариев.
Публикаций о том как делать БПФ немеряно, а о том как сделать БПФ, преобразовать спектр, и собрать сигнал заново, да еще и в реальном времени, явно не хватает. Автор пытается восполнить этот пробел.
Часть первая, обзорная
Существуют два основных способа построения дискретных линейных динамических систем. В литературе, такие системы принято называть цифровыми фильтрами, которые подразделяются на два основных типа: фильтры с конечной импульсной характеристикой (КИХ) и фильтры с бесконечной импульсной характеристикой (БИХ).
Алгоритмическая сущность фильтра с КИХ заключается в дискретном вычислении интеграла свертки:
Где x(t) – входной сигнал
y(t) – выходной сигнал
h(t) – импульсная характеристика фильтра или реакция фильтра на дельта функцию. Импульсная характеристика является обратным преобразованием Фурье комплексной частотной характеристики фильтра K(f).
Для формирования ясной картины у читателя, приведем пример дискретного вычисления интеграла свертки на языке С в реальном времени.
Вызывая данную функцию через определенные интервалы времени T и передавая ей в качестве аргумента входной сигнал, на выходе мы получим выходной сигнал, соответствующий реакции фильтра с импульсной характеристикой вида:
h(t)=1 при 0 >alfa);, но в этом случае происходит потеря alfa значащих разрядов. Рекуррентное выражение фильтра, из примера кода, построено таким образом, чтобы избежать потери значащих разрядов. Именно конечная точность вычислений может испортить всю прелесть цифрового фильтра с бесконечной импульсной характеристикой. Особенно это заметно на фильтрах высоких порядков, отличающихся высокой добротностью. В реальных динамических системах такая проблема не возникает, наша Матрица производит вычисления с невероятной для нас точностью.
Синтезу подобных фильтров посвящена масса литературы, также имеются готовые программные продукты (см. выше).
Часть вторая. Фурье – фильтр
Из вузовских курсов (у вашего покорного слуги это был курс ОТЭЦ) многие собравшие помнят два основных подхода к анализу линейных динамических систем: анализ во временной области и анализ в частотной области. Анализ во временной области — это решение дифференциальных уравнений, интегралы свертки и Дюамеля. Эти методы анализа дискретно воплотились в цифровых фильтрах БИХ и КИХ.
Но существует частотный подход к анализу линейных динамических систем. Иногда его называют операторным. В качестве операторов используются преобразование Фурье, Лапласа и т.п. Далее мы будем говорить только о преобразовании Фурье.
Данный метод анализа не получил широкого распространения при построении цифровых фильтров. Автору не удалось найти вменяемых практических рекомендаций по построению подобных фильтров на русском языке. Единственное краткое упоминание такого фильтра в практической литературе [Рабинер Л., Гоулд Б., Теория и применение цифровой обработки сигналов 1978], но в данной книге рассмотрение подобного фильтра очень поверхностно. В указанной книге данная схема построения фильтра называется: «свертка в реальном времени методом БПФ», что, по моему скромному мнению, совершенно не отражает сути, название должно быть коротким, иначе времени на отдых не останется.
Реакция линейной динамической системы есть обратное преобразование Фурье от произведения изображения по Фурье входного сигнала x(t) на комплексный коэффициент передачи K(f):
В практическом плане, данное аналитическое выражение предполагает следующий порядок действий: берем преобразование Фурье от входного сигнала, умножаем результат на комплексный коэффициент передачи, выполняем обратное преобразование Фурье, результатом которого является выходной сигнал. В реальном дискретном времени такой порядок действий выполнить невозможно. Как брать интеграл по времени от минус до плюс бесконечности?! Его можно взять только находясь вне времени…
В дискретном мире для выполнения преобразования Фурье существует инструмент — алгоритм быстрого преобразования Фурье (БПФ). Именно его мы и будем использовать при реализации нашего Фурье-фильтра. Аргументом функции БПФ является массив временных отсчетов из 2^n элементов, результатом два массива длинной 2^n элементов соответствующие действительной и мнимой части преобразования Фурье. Дискретной особенностью алгоритма БПФ является то, что входной сигнал считается периодичным с интервалом 2^n. Это накладывает некоторые ограничения на алгоритм Фурье-фильтра. Если взять последовательность выборок входного сигнала, провести от них БПФ, умножить результат БПФ на комплексный коэффициент передачи фильтра и выполнить обратное преобразование …ничего получится! Выходной сигнал будет иметь огромные нелинейные искажения в окрестности стыков выборок.
Для решения этой проблемы необходимо применить два приема:
Такие функции широко применяются в технике цифровой обработки сигналов, и называть их принято — окнами. По скромному мнению автора лучшим, с практической точки зрения, является окно имени Хана:
На рисунке приведены графики иллюстрирующие свойства окна Хана длинной 2^n=256. Экземпляры окна построены с половинным перекрытием k=128. Как видно все оговоренные выше свойства имеются в наличии.
По просьбам трудящихся, на следующем рисунке приведена схема вычислений Фурье-фильтра, при длине выборки 2^n=8, количество выборок 3. На подобных рисунках очень сложно отобразить процесс вычислений, особенно тяжело показать его цикличность, поэтому мы и ограничились количеством выборок равным трем.
Входной сигнал разбивается на блоки длинной 2^n=8 с перекрытием 50%, от каждого блока берется БПФ, результаты БПФ подвергаются нужной трансформации, берется обратное БПФ, результат обратного БПФ скалярно умножается на окно, после умножения блоки складываются с перекрытием.
При выполнение трансформаций спектра, не стоит забывать о главном свойстве массива БПФ действительных сигналов, первая половина массива БПФ комплексно сопряжена со второй половиной, т.е Re[i]=Re[(1
Фурье-вычисления для сравнения изображений
Традиционная техника “начального уровня”, сравнения текущего изображения с эталоном основывается на рассмотрении изображений как двумерных функций яркости (дискретных двумерных матриц интенсивности). При этом измеряется либо расстояние между изображениями, либо мера их близости.
Как правило, для вычисления расстояний между изображениями используется формула, являющаяся суммой модулей или квадратов разностей интенсивности:
Если помимо простого сравнения двух изображений требуется решить задачу обнаружения позиции фрагмента одного изображения в другом, то классический метод “начального уровня”, заключающийся в переборе всех координат и вычисления расстояния по указанной формуле, как правило, терпит неудачу практического использования из-за требуемого большого количества вычислений.
Одним из методов, позволяющих значительно сократить количество вычислений, является применение Фурье преобразований и дискретных Фурье преобразований для расчёта меры совпадения двух изображений при различных смещениях их между собой. Вычисления при этом происходят одновременно для различных комбинаций сдвигов изображений относительно друг друга.
Наличие большого числа библиотек, реализующих Фурье преобразований (во всевозможных вариантах быстрых версий), делает реализацию алгоритмов сравнения изображений не очень сложной задачей для программирования.
Постановка задачи
образец
в изображении
Корреляция как мера между изображениями
m(X,Y) = SUM ( X[i,j] * Y[i,j] ) / ( SQRT ( SUM X[i,j] ^2 ) * SQRT ( SUM Y[i,j] ^2 ) )
Данная величина получена из операции скалярного произведения векторов (рассматривая изображения как векторы в многомерном пространстве). И даже более — эта же формула представляет собой и стандартную статистическую формулу критерия для гипотезы о совпадении двух вероятностных распределений.
Примечание:
При вычислении корреляции между фрагментами изображений, если одно изображение меньше другого, будем делить только на значение норм у пересекающийся частей.
Свёртка двух функций
Преобразование Фурье
Преобразование Фурье (ℱ) — операция, сопоставляющая одной функции вещественной переменной другую функцию, также вещественной переменной. Эта новая функция описывает коэффициенты («амплитуды») при разложении исходной функции на элементарные составляющие — гармонические колебания с разными частотами.
Преобразование Фурье функции f вещественной переменной является интегральным и задаётся следующей формулой:
Разные источники могут давать определения, отличающиеся от приведённого выше выбором коэффициента перед интегралом, а также знака «−» в показателе экспоненты. Но все свойства будут те же, хотя вид некоторых формул может измениться.
Кроме того, существуют разнообразные обобщения данного понятия.
Многомерное преобразование Фурье
Преобразование Фурье функций, заданных на пространстве ℝ^n, определяется формулой:
Обратное преобразование в этом случае задается формулой:
Как и прежде, в разных источниках определения многомерного преобразования Фурье могут отличаться выбором константы перед интегралом.
Дискретное преобразование Фурье
Дискретное преобразование Фурье (в англоязычной литературе DFT, Discrete Fourier Transform) — это одно из преобразований Фурье, широко применяемых в алгоритмах цифровой обработки сигналов (его модификации применяются в сжатии звука в MP3, сжатии изображений в JPEG и др.), а также в других областях, связанных с анализом частот в дискретном (к примеру, оцифрованном аналоговом) сигнале. Дискретное преобразование Фурье требует в качестве входа дискретную функцию. Такие функции часто создаются путём дискретизации (выборки значений из непрерывных функций). Дискретные преобразования Фурье помогают решать дифференциальные уравнения в частных производных и выполнять такие операции, как свёртки. Дискретные преобразования Фурье также активно используются в статистике, при анализе временных рядов. Существуют многомерные дискретные преобразования Фурье.
Формулы дискретных преобразований
Дискретное преобразование Фурье является линейным преобразованием, которое переводит вектор временных отсчётов в вектор спектральных отсчётов той же длины. Таким образом преобразование может быть реализовано как умножение симметричной квадратной матрицы на вектор:
Фурье-преобразования для вычисления свёртки
Одним из замечательных свойств преобразований Фурье является возможность быстрого вычисления корреляции двух функций определённых, либо на действительном аргументе (при использовании классической формулы), либо на конечном кольце (при использовании дискретных преобразований).
И хотя подобные свойства присущи многим линейным преобразованиям, для практического применения, для вычисления операции свёртки, согласно данному нами определению, используется формула
Фурье-преобразования для вычисления корреляции
Пусть (t) равна корреляции получаемой в результате сдвига одного вектора, относительно другого на шаг t
Тогда, как уже показано ранее, выполняется
Фурье-преобразования для решения задачи
Упрощение формул для решения поставленной задачи
Алгоритм поиска фрагмента в полном изображении
Примеры реализации
Реализованные алгоритмы являются частью библиотеки с открытым исходным кодом FFTTools. Интернет-адрес: github.com/dprotopopov/FFTTools
Преобразование Фурье в действии: точное определение частоты сигнала и выделение нот
Начнём с пианино. Очень упрощёно этот музыкальный инструмент представляет собой набор белых и чёрных клавиш, при нажатии на каждую из которых извлекается определённый звук заранее заданной частоты от низкого до высокого. Конечно, каждый клавишный инструмент имеет свою уникальную тембральную окраску звучания, благодаря которой мы можем отличить, например, аккордеон от фортепиано, но если грубо обобщить, то каждая клавиша представляет собой просто генератор синусоидальных акустических волн определённой частоты.
Когда музыкант играет композицию, то он поочерёдно или одновременно зажимает и отпускает клавиши, в результате чего несколько синусоидальных сигналов накладываются друг на друга образуя рисунок. Именно этот рисунок воспринимается нами как мелодия, благодаря чему мы без труда узнаём одно произведение, исполняемое на различных инструментах в разных жанрах или даже непрофессионально напеваемое человеком.
Наглядная иллюстрация нотного рисунка
Определение частоты (режим гитарного тюнера)
Обратная задача состоит в том, чтобы разобрать звучащую музыкальную композицию на ноты. То есть разложить суммарный акустический сигнал, улавливаемый ухом, на исходные синусоиды. По сути, этот процесс и представляет собой прямое преобразование Фурье. А нажатие на клавиши и извлечение звука есть процесс обратного преобразования Фурье.
Математически в первом случае происходит разложение сложной периодической (на некотором временном интервале) функции в ряд более элементарных ортогональных функций (синусоид и косинусоид). А во втором их обратное суммирование, то есть синтез сложного сигнала.
Ортогональность, в некотором роде, обозначает несмешиваемость функций. Например, если мы возьмём несколько кусочков цветного пластилина и склеим их, то потом всё же сможем разобрать, какие цвета были изначально, но если хорошенько перемешаем несколько баночек гуашевых красок, то точно восстановить исходные цвета без дополнительной информации уже будет невозможно.
(!) Важно понимать, когда мы берёмся анализировать реальный сигнал с помощью преобразования Фурье, мы идеализируем ситуацию и исходим из предположения, что он периодический на текущем временном интервале и состоит из элементарных синусоид. Зачастую это именно так, поскольку акустические сигналы, как правило, имеют гармоническую природу, но вообще возможны и более сложные случаи. Любые наши допущения о природе сигнала обычно ведут к частичным искажениям и погрешностям, но без этого выделить полезную информацию из него крайне сложно.
Теперь опишем весь процесс анализа более подробно:
1. Всё начинается с того, что звуковые волны колеблют мембрану микрофона, который преобразует их в аналоговые колебания электрического тока.
2. Затем происходит дискретизация аналогового электрического сигнала в цифровую форму. На этом моменте стоит остановиться подробно.
Поскольку аналоговый сигнал математически состоит из бесконечного непрерывного во времени множества точек-значений амплитуды, в процессе измерения мы можем выделить из него лишь конечный ряд значений в дискретные моменты времени, то есть, по сути, выполнить квантование по времени…
Как правило, значения-отсчёты берутся через небольшие равные временные промежутки, то есть с определённой частотой, например, 16000 или 22000 Гц. Однако в общем случае дискретные отсчёты могут идти и неравномерно, но это усложняет математический аппарат анализа, поэтому на практике обычно не применяется.
Существует важная теорема Котельникова-Найквиста-Шеннона, которая гласит, что аналоговый периодический сигнал, имеющий конечный (ограниченный по ширине) спектр, может быть однозначно восстановлен без искажений и потерь по своим отсчётам, взятым с частотой, большей или равной удвоенной верхней частоте спектра (называемой частотой дискретизации или Найквиста).
Для этого восстановления необходимо применить специальные интерполирующие функции, но проблема в том, что при использовании данных функций вычисления нужно выполнять на бесконечном временном интервале, что на практике невозможно. Поэтому в реальной жизни нельзя сколь угодно повысить частоту дискретизации искусственным образом без искажений даже если изначально она удовлетворяет теореме Котельникова-Найквиста-Шеннона. Для этой операции применяются фильтры Фарроу.
Также дискретизация происходит не только по времени, но и по уровню значений амплитуды, поскольку компьютер способен манипулировать лишь ограниченным множеством чисел. Это также вносит небольшие погрешности.
3. На следующем этапе происходит само дискретное прямое преобразование Фурье.
Мы выделяем короткий кадр (интервал) композиции, состоящий из дискретных отсчётов, который условно считаем периодическим и применяем к нему преобразование Фурье. В результате преобразования получаем массив комплексных чисел, содержащий информацию об амплитудном и фазовом спектрах анализируемого кадра. Причём спектры также являются дискретными с шагом равным (частота дискретизации)/(количество отсчётов). То есть чем больше мы берём отсчётов, тем более точное разрешение получаем по частоте. Однако при постоянной частоте дискретизации увеличивая число отсчётов, мы увеличиваем анализируемый временной интервал, а поскольку в реальных музыкальных произведениях ноты имеют различную длительность звучания и могут быстро сменять друг друга, происходит их наложение, поэтому амплитуда длительных нот «затмевает» собой амплитуду коротких. С другой стороны для гитарных тюнеров такой способ увеличения разрешения по частоте подходит хорошо, поскольку нота, как правило, звучит долго и одна.
Существует также довольно простой трюк для увеличения разрешения по частоте — нужно исходный дискретный сигнал заполнить нулями между отсчётами. Однако в результате такого заполнения сильно искажается фазовый спектр, но зато увеличивается разрешение амплитудного. Также возможно применение фильтров Фарроу и искусственное увеличение частоты дискретизации, однако и оно вносит искажения в спектры.
Длительность кадра обычно составляет приблизительно от 30 мс до 1 с. Чем он короче, тем лучшее разрешение мы получаем по времени, но худшее по частоте, чем сэмпл длиннее, тем лучшее по частоте, но худшее по времени. Это очень напоминает принцип неопределённости Гейзенберга из квантовой механики..и не с проста, как гласит Википедия, соотношение неопределенностей в квантовой механике в математическом смысле есть прямое следствие свойств преобразования Фурье…
Интересно и то, что в результате анализа сэмпла одиночного синусоидального сигнала амплитудный спектр очень напоминает дифракционную картинку…
Синусоидальный сигнал, ограниченный прямоугольным окном, и его «дифракция»
Дифракция световых волн 
На практике это нежелательный эффект, затрудняющий анализ сигналов, поэтому его стараются понизить путём применения оконных функций. Таких функций придумано немало, ниже представлены реализации некоторых из них, а также сравнительное влияние на спектр одиночного синусоидального сигнала.
Применяется оконная функция ко входному кадру очень просто:
Что касается компьютеров, в своё время был разработан алгоритм быстрого преобразования Фурье, который минимизирует число математических операций, необходимых для его вычисления. Единственное требование алгоритма состоит в том, чтобы число отсчётов было кратно степени двойки (256, 512, 1024 и так далее).
Ниже его классическая рекурсивная реализация на языке C#.
Существует две разновидности алгоритма БПФ — с прореживанием по времени и по частоте, но оба дают идентичный результат. Функции принимают массив комплексных чисел, заполненный реальными значениями амплитуд сигнала во временной области, а после своего выполнения возвращают массив комплексных чисел, содержащий информацию об амплитудном и фазовом спектрах. Стоит помнить, что реальная и мнимая части комплексного числа — это далеко не то же самое, что его амплитуда и фаза!
magnitude = Math.Sqrt(x.Real*x.Real + x.Imaginary*x.Imaginary)
phase = Math.Atan2(x.Imaginary, x.Real)
Результирующий массив комплексных чисел заполнен полезной информацией ровно на половину, другая половина является лишь зеркальным отражением первой и спокойно может быть исключена из рассмотрения. Если вдуматься, то этот момент хорошо иллюстрирует теорему Котельникова-Найквиста-Шеннона, о том, что частота дискретизации должна быть не меньше максимальной удвоенной частоты сигнала…
Также существует разновидность алгоритма БПФ без рекурсии по Кули-Тьюки, которая часто применяется на практике, но она чуть более сложна для восприятия.
Сразу после вычисления преобразования Фурье удобно нормализовать амплитудный спектр:
Это приведёт к тому, что величина значений амплитуды получится одного порядка не зависимо от размеров сэмпла.
Вычислив амплитудный и частотный спектры, легко производить обработку сигнала, например, применять частотную фильтрацию или производить сжатие. По сути, таким образом можно сделать эквалайзер: выполнив прямое преобразование Фурье, легко увеличить или уменьшить амплитуду определённой области частот, после чего выполнить обратное преобразование Фурье (хотя работа настоящих эквалайзеров обычно основана на другом принципе — фазовом сдвиге сигнала). Да и сжать сигнал очень просто — нужно всего лишь сделать словарь, где ключом является частота, а значением соответствующее комплексное число. В словарь нужно занести лишь те частоты, амплитуда сигнала на которых превышает какой-то минимальный порог. Информация о «тихих» частотах, не слышимых ухом, будет потеряна, но получится ощутимое сжатие при сохранении приемлемого качества звучания. Отчасти этот принцип лежит в основе многих кодеков.
4. Точное определение частоты
Дискретное преобразование Фурье даёт нам дискретный спектр, где каждое значение амплитуды отстоит от соседних на равные промежутки по частоте. И если частота в сигнале кратна шагу равному (частота дискретизации)/(количество отсчётов), то мы получим выраженный остроконечный пик, но если частота сигнала лежит где-то между границами шага ближе к середине у нас выйдет пик со «срезанной» вершиной и нам будет затруднительно сказать, что же там за частота. Очень может быть что в сигнале присутствуют две частоты лежащие рядом друг с другом. В этом и заключается ограничение разрешения по частоте. Так же как на фотоснимке с низким разрешением мелкие предметы склеиваются и становятся неразличимы, так же и тонкие детали спектра могут теряться.
Но частоты музыкальных нот лежат далеко не на сетке шагов преобразования Фурье, а для повседневных задач настройки музыкальных инструментов и распознавания нот необходимо знать именно точную частоту. Более того, на низких октавах при разрешении от 1024 отсчётов и ниже сетка частот Фурье становится настолько редкой, что попросту на одном шаге начинают умещаться несколько нот и определить какая же на самом деле из них играет становится фактически невозможно.
Чтобы как-то обойти это ограничение иногда применяют аппроксимирующие функции, например, параболические.
www.ingelec.uns.edu.ar/pds2803/Materiales/Articulos/AnalisisFrecuencial/04205098.pdf
mgasior.web.cern.ch/mgasior/pap/biw2004_poster.pdf
Но всё это искусственные меры, которые улучшая одни показатели могут давать искажения в других.
Существует ли более естественный путь для точного определения частоты?
Да, и скрыт он как раз-таки в использовании фазового спектра сигнала, которым часто пренебрегают.
Данный метод уточнения частоты сигнала, основан на вычислении задержки фаз у спектров двух кадров, наложенных друг на друга, но немного сдвинутых во времени.
На C# реализация метода выглядит довольно просто:
Применение также несложное:
Обычно исходные кадры сдвинуты на 1/16 или 1/32 своей длины, то есть ShiftsPerFrame равно 16 или 32.
В результате мы получим словарь частота-амплитуда, где значения частот будут довольно близки к реальным. Однако «срезанные пики» всё ещё будут наблюдаться, хоть и менее выражено. Чтобы устранить этот недостаток, можно просто «дорисовать» их.
Нотный анализ музыкальных произведений открывает ряд интересных возможностей. Ведь имея в наличии готовый нотный рисунок, можно осуществлять поиск других музыкальных композиций со схожим рисунком.
Например, одно и то же произведение может быть исполнено на другом инструменте, в различной манере, с другим тембром, либо транспонировано по октавам, однако нотный рисунок останется похожим, что позволит найти различные варианты исполнения одного и того же произведения. Это очень напоминает игру «угадай мелодию».
В некоторых случаях подобный анализ поможет выявить плагиат в музыкальных произведениях. Также по нотному рисунку, теоретически, можно искать произведения определённого настроения или жанра, что поднимает поиск на новый уровень.
В этой статье изложены основные принципы точного определения частот акустических сигналов и выделения нот. А также показана некоторая тонкая интуитивная связь дискретного преобразования Фурье с квантовой физикой, что подталкивает на размышления о единой картине мира.