Pull to refresh

Comments 63

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

Кодек MP2/3 смотрит с недоумением и с перепуга пытается начать крутить песни по кругу в попытке сделать их периодическими...

Сложно себе представить работу такого фильтра в реальном масштабе времени. Ведь надо много вычислений O(n^2)

На этом месте MP3 заплакал и ушёл, прихватив с собой половину современной электроники, а среди математиков и программистов началась повальное пьянство.

Надо заранее знать самую низкую частоту. Получившаяся формула ряде Фурье будет отрисовывать сигнал с периодом равным размеру выборки умноженную на период дискретизации.

Интересно у с тут получается. Что значит "самую низкую"? Если у вас периодический сигнал, то волшебным образом "самая низкая частота" равна его периоду. А если непериодический, то разложение Фурье совершенно точно по вашим же словам к нему не применимо, и в рамках предлагаемой логики вопрос не имеет смысла. Да и вообще, нижняя граница частоты непериодического сигнала равна его длительности, что судя по всему делает спектральную обработку таких сигналов если не невозможной, то невероятно сложной.

Самая низкая - это разве не нулевая (постоянная составляющая)?

 Что значит "самую низкую"? Если у вас периодический сигнал, то волшебным образом "самая низкая частота" равна его периоду. 

Ситуация такая:

Есть таблица чисел от измерения сигнала. Про этот сигнал мало что известно (ровным счетом ничего).

Не известна даже частота его самой низкочастотной гармоники. Не ясно сколько отсчетов брать для преобразования DFT.

Как определить период самой длинной гармоники, если, например, записана выборка всего длительностью t=0.75*T от этого самого длинного периода T?

Тут преобразование Фурье не поможет.

Тут перепутана лошадь и повозка, сдаётся мне. Нужно как раз наоборот - получить спектр и выбрать необходимую частоту, ниже которой сделать сигнала уже не интересна для ожидаемого результата обработки. Очевидно, что даже не нужно строить полный спектр, достаточно построить только интересующий низ и с интересной дискретностью. Короче, такой проблемы нет, а если она по какой-то причине есть, то её решение известно и относительно несложно.

Ряд Фурье для периодических сигналов. Интеграл Фурье - для непериодических. Периодический сигнал имеет дискретный спектр, непериодический сигнал имеет непрерывный (сплошной) спектр.

FIR фильтры, DFT в цифровой обработке сигналов? Назад в прошлое?

Используйте БИХ-фильтры и fft, их уже минимум 80 лет как раз для этого изобрели (или 220 лет назад, смотря как считать), и работают они намного быстрее.

По крайней мере БИХ-фильтром можно реализовать примитивную низкочастотную фильтрацию за три арифметических действия очень легко:

OUT(t)=(1-a)*OUT(t-1)+a*IN(t)

a \in [0,1]

Чем меньше a, тем у'же полоса НЧ фильтрации. Если нужен более резкий срез в частотной области, добавляют элементов в правой части из OUT и IN, для деталей почитайте про проектирование БИХ-фильтров, есть и готовые программки для расчета коэффициентов для заданных характеристик фильтров.

Ведь 10 операций сложение/умножение для непосредственного получения фильтрованного сигнала при БИХ-фильтрации выполняются значительно быстрее, чем 10^5 .. 10^6 операций при DFT или 10^3..10^4 при FIR.

Вообще в описанной задаче больше подошло бы двойное или тройное экспоненциальное сглаживание.

Кстати, то, что написали Вы, иначе называется просто экспоненциальным сглаживанием.

Все три варианта отлично подходят для потоковой обработки данных.

Есть и более сложные (вычислительно) фильтры. Например, Ходрика-Прескотта

Вообще, тема анализа временных рядов (в т.ч. и с учетом сезонности) исследована достаточно глубоко.

Насчет экспоненциального сглаживания согласен. А если попробовать вместо a=-0.5, а вместо (1-a)=0.5 то получите фильтр высоких частот. В общем, БИХ - очень широкий и эффективный класс фильтров, где можно сравнительно небольшим числом операций решать достаточно сложные задачи.

БИХ фильтр искажает ФЧХ, что в ряде задач недопустимо. Выходом может быть блочная обработка и filtfilt, но при этом невозможно обрабатывать сигнал в реальном времени.

Никогда особно не вникал в эту проблему (чаще всего работал с реальными данными, а у них фаза всегда 0, поэтому 'фазовых искажений' нет в принципе), но по построению для комплексной области это зависит от фазы коэффициентов. Если исходный сигнал задан с фазой (в комплексном виде или в квадратурных компонентах ), то скорее всего можно сделать фазовый доворот комплекснозначных коэффициентов a_i, b_i для необходимой компенсации фазы или ФЧХ. При необходимости это можно делать в динамике, 10 комплексных умножений много времени не сожрут.

Фаза у любых реальных сигналов, у которых может гулять частота, будет смещаться в зависимости от ФЧХ БИХ фильтра. Тут компенсировать фазу может оказаться сложновато. Обычно БИХ применяются там, где важна амплитуда, а не фаза; либо там, где частота сигнала находится вдали от сильных искажений.

Не всегда БИХ/КИХ фильтры подходят под задачу. БПФ, например, неплохо подходит там, где нужно отфильтровать сигнал и получить фазу сигнала одновременно. Приходится подбирать инструменты под задачу и все таки это не "назад в прошлое".

В принципе, да, в каком-то смысле согласен. Но есть одно но - терминологическая путаница. Фазы у реальнозначных сигналов нет, потому она и смещаться никуда не может. Эта терминологическая ошибка - исключительно проблема их неудачного представления в виде ряда Фурье.

Для реальных сигналов представление в виде ряда Фурье используется традиционно, но оно не совсем корректно. Отсюда и возникающие проблемы с кажущейся фазой (которой реально нет, ведь для реальнозначных сигналов есть лишь задержки). Если использовать более корректное для этого случая преобразование Хартли (где в качестве базовой функции выступает cas(x)=cos(x)+sin(x), которое еще и считается в 2 раза быстрее, чем БПФ), проблем с фазой (которой нет) не возникнет, как того и ожидается для реальнозначных сигналов, у которых фаза нулевая всегда по определению. Просто преобразование Хартли корректный, но не очень общепринятый метод анализа.

Ну а в электрических схемах нужно использовать квадратурные компоненты сигнала (или пару ток-напряжение), для них все честно и корректно, и фазовая информация не теряется.

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

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

Есть достаточно простой приём для коррекции ФЧХ любых IIR и FIR - после фильтрации сделать ещё одну обратив во времени полученный сигнал. ФЧХ скомпенсируется, а частотные свойства фильтра улучшатся.

Сезонные данные так не анализируют- потому что, к примеру, солнце всходит и заходит в разное время каждый день. Отфильтровать можно с помощью оконного фильтра Гаусса (как любой оконный фильтр, легко применим и к потоковым данным с некоторым аккумулятором). А для полноценного разложения на константу, тренд и сезонную компоненту (или компоненты) есть специальные методы, см. https://www.statsmodels.org/dev/examples/notebooks/generated/stl_decomposition.html

P.S. Давно уже придумано оконное преобразование Фурье, или вейвлеты, которые намного лучше применимы в мире не-бесконечных сигналов.

В софте типа гитарных процесоров фильтр - Ких - импульсная реакция колонки или "гитарноного кабинета" как вычисляется ких фильтром на скользящем БПФ в РЕАЛЬНОМ врмени, причем импульсная реакция колонки или длина КИХ может быть скажем 65к точек и все работает в реальном времени на умеренно мощных процессорах. И это уже примерно как 7..10 лет существует. Подробнее вот ссылка https://guitarclan.com/free-cabinet-impulse-responses/ Добавлю именно в реальном времени и именно на бесконечном потоке данных и о ужас даже с задержкой близкой к нескольким миллисекундам, теория математики этого отлично разработана и опубликована в открытых источниках, довольно сложна, но это так на самом деле КИХ на БПФ организованный определеным хитрым способом почти не имеет задержки из-за какзалось обязательного буферизирования данных для блока БПФ - а вот и нет, математики это решили! (хитрое каскадиование и переменных длины буфер от 1 мс до нескольких секунд). Гитаристы играют в живую через плгуины VST с фильтрами КИХ на БПФ с длиной импульса скажем в 0.5 сек или 22000 отсчетов с задержкой 5..10 мс! (или 256 отсчетов) Да задержка меньше длины КИХ фильтра в примерно 100 раз. Казалось бы невозможно - но это факт и алгоритмы для этого легко гугляться. (фильтр на скользащем БПФ)

Добавлю именно в реальном времени и именно на бесконечном потоке данных и о ужас даже с задержкой близкой к нескольким миллисекундам, теория математики этого отлично разработана и опубликована в открытых источниках, довольно сложна, но это так на самом деле КИХ на БПФ организованный определеным хитрым способом почти не имеет задержки из-за какзалось обязательного буферизирования данных для блока БПФ - а вот и нет, математики это решили! (хитрое каскадиование и переменных длины буфер от 1 мс до нескольких секунд).

Бесконечного потока данных не бывает, в лучшем случае, есть непрерывный. А вот на границах разрывов сигнала что происходит? (ничего хорошего). Далее, сначала вы клянетесь про «задержкой близкой к нескольким миллисекундам», а потом «буфер от 1 мс до нескольких секунд» - тут разница на три порядка. После начала регистрации сигнала или его резком изменении сколько времени заполнения буфера ждать? Более того, и в микросекундном интервале бывает достаточно отсчетов для обработки сигналов, и наоборот (вы годовые изменения с миллисекундным буфером изучать собрались?), так что считать в единицах времени, а не в количестве отсчетов - вообще бессмысленно. Не надо ссылаться на «теорию математики», надо с фактами поаккуратнее обращаться.

длиной импульса скажем в 0.5 сек или 22000 отсчетов с задержкой 5..10 мс! 

Задержка имеет место быть в начале, сначала должны пройти эти 22000 отсчетов для вычисления отклика, а потом отфильтрованный сигнал будет обновляться с приходом каждого нового отсчета. Ну или как там будет успевать обсчитывать процессор.
Для простого скользящего среднего сначала приходят эти 22000 отсчетов, вычисляется среднее, а потом с приходом нового отсчета к полученному прибавляется этот отсчет
Xn+1 /22000 и вычитается первый отсчет X0/22000 .

близко! вы поняли но не до конца... ких фильтр длиной скажем 128000 можно разбить на 1000 фильтров по 128 отсчетов (фактически все 1000 фильтров по 128000 но не равных нулю коэф. ких будет по 128 штук), 1-й ких имеет импульсную реакцию отличную от 0 на дискретах с 0 до 127, 2-й ких фильтр имеет импульсную реакцию отличную от 0 на отсчетах от 128 до 255, и тд.

1-й фильтр считается, да, с задержкой для заполнения буффера, 128 дискретов, второй и все последующие начинают считаться в паралель и в момент начала расчета 1-го фильтра - с учётом того что ких там где 0 НЕ считается и ресурсы не тратятся, результаты суммируются с учетом задержке при разбиении фильтра на 1000 фильтров. И так все время. При такой схеме хадержка будет не 128000 а только 128 дискретов. Есть и более сложные схемы где задержка будет только 1 дискрет. Попробуйте придумать сами!

Неважно, как вы считаете, окончательный результат будет после того, как придут все 128000 отсчетов. Если результат готов через 128 дискрет, то остальные 128000-128 отсчетов никак не используются, просто потому, что их еще не существует в природе. Это если вы хотите посчитать среднее за 128000 отсчетов.

Кто-то из вас вероятно путает задержку при реализации свёртки и задержку как одну из характеристик непосредственно самой импульсной характеристики.

"Задержка при реализации" - это что имеется в виду? Которая зависит от того, как быстро считает компьютер?

Имеется в виду:

а) обычная реализация, так называемая "линейная" или "прямая" свёртка - два цикла один внутри другого с перемножением и сложением. Сложность n^2 , задержки нет.

б) "Быстрая" свёртка через FFT - прямое преобразование, комплексное перемножение спектров, обратное преобразование. Сложность n*log(n) , задержка равна размеру FFT.

в) "Overlap-Add" алгоритм - вариация предыдущего пункта для случая, когда один сигнал длинный, а второй короткий (ядро свёртки). Позволяет уменьшить размер FFT. Сложность и задержки те же.

г) "Zero-delay fast convolution", комбинация а) и в). Сложность чуть выше предыдущего, задержки нет. Ядро свёртки разбивается по степеням двойки, например 128=16+16+32+64. Первые 16 считаются прямой свёрткой, остальные через FFT.

Я реализовывал все эти варианты.

Такого рода "задержка" будет и на аналоговых фильтрах. Сюрприз.

Вовсе не сюрприз. Это физика.

--Сложно себе представить работу такого фильтра в реальном масштабе времени. Ведь надо много вычислений O(n^2)

Вам ещё предстоит открыть для себя быстрое преобразование Фурье (FFT: fast Fourier transform).

--Преобразование Фурье не подходит для непериодических сигналов.

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

--Если вычислить DFT, обнулить высокие частоты и вычислить обратное DFT,
то не меняется фаза тех низкочастотных гармоник, которые мы оставили.

Когда частотное представление получено, его можно умножить на любую функцию частоты, не обязательно ступеньку. В частности, подобрав должным образом коэффициенты, можно получить тот самый FIR.

А вообще, сигналы можно фильтровать с помощью (почти) любого разложения по полной ортогональной системе функций ([напр.](https://core.ac.uk/download/36719692.pdf)). И даже необязательно ортогоноальной (см. вейвлеты).

Вам ещё предстоит открыть для себя быстрое преобразование Фурье (FFT: fast Fourier transform).

FFT работает только с сигналами размер выборки которых равен степеням двойки. А тут размер выборки непредсказуемый. Поэтому FFT не подходит.

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

И кстати выборку всегда можно передискретизировать к степени двойки.

Если вычислить DFT, обнулить высокие частоты и вычислить обратное DFT, то не меняется фаза тех низкочастотных гармоник, которые мы оставили.

Умножение в частотной области (на ступеньку с занулением определённых частот) == свёртка во временной, ровно то что КИХ фильтры и делают без преобразования туда обратно.

Фаза и так не портится, ФЧХ абсолютна линейна и соответствует просто задержке на половину длины фильра.

Ваше "скользящее среднее" просто надо считать не от X[i]..X[i+N] а от X[i-N/2] до X[i+N/2], и внезапно ничего никуда смещаться не будет.

Дааа, вот он хабр где комментарии можно как вторую статью читать :)

При этом интерес представляет самая низкая частотная составляющая этого
сигнала. Та которая с периодом 24 часа. Надо определить амплитуду и фазу
самой первой гармоники.

Понятно, что надо как-то обобщить этот сигнал. Пропустить его через фильтр нижних частот.

Если мы пропустим сигнал через цифровой фильтр скользящего среднего (FIR
фильтр), то выходной сигнал сместится по фазе вправо и реальная фаза
окажется искажена. Нужен другой способ.

Почему для выделения первой гармоники Вы используете ФНЧ, а не узкополосный полосовой фильтр? Простейший аналог из схемотехники - фильтр на основе одиночного колебательного контура. На частоте резонанса он не дает фазового сдвига.

Почему для выделения первой гармоники Вы используете ФНЧ, а не узкополосный полосовой фильтр?

Потому что частота первой гармоники в реальности не известна.

а для этого есть автокорреляция,

и затем гражданин Герцель для DFT с одним умножением на отсчёт, без тригонометрии.

А какую тогда Вы хотели установить частоту среза ФНЧ?

А какую тогда Вы хотели установить частоту среза ФНЧ?

Очень низкую.
Fc = 1.0/(30.0*3600.0 [s] )= 9.259259259259259e-06 [Hz]

Исходя из чего принято такое значение? Ведь тогда первая гармоника из Вашего же примера попадает за границу полосы пропускания. Это далеко от принципов оптимальной фильтрации сигналов, используемых для максимизации соотношения сигнал/шум. Ну а про фазовые сдвиги сигнала за пределами полосы пропускания и говорить нечего - чем круче фильтр, тем хуже будет ситуация с фазой. А насколько можно понять из Вашей же постановки задачи - фаза для Вас важна.

Исходя из чего принято такое значение?

Я хочу экспериментально увидеть в преобразовании Фурье, что период первой гармоники графика естественной освещенности на самом деле не точно равен 24 часа, а немного меняется в течении года, так как Земля движется вокруг Солнца по эллипсу и угол доворота на Солнце каждый день разный.

То есть Вам нужна только временная характеристика сигнала - частота или фаза? А в статье Вы писали:

Надо определить амплитуду и фазу самой первой гармоники.

Я думаю, что такой "разброд" мнений комментаторов вызван именно нечеткой формулировкой постановки задачи.

Лучше это сначала смоделировать, как там меняются углы освещенности при движении по эллипсу, сгенерировать отсчеты и подсунуть эту выборку Фурье.
Если результат совпадет с ожидаемым, то хорошо.

задача сводится к изменению первой гармоники с точностью в районе сотых долей процента. При этом вариация частоты полностью компенсируется всего на 365 циклах. Выглядит сложно. Мне почему-то кажется, что придется отслеживать период одного колебания точнее. Например щелевой колпачок на фотодатчик. Тогда разность фаз между сигналом и им же, сдвинутым на сколько-то дней, будет явно видна.

Безумная идея, не судите строго:

Фотодатчик, естественно неподвижный относительно земли. На него надеваем щелевой колпачок. нацелим его скажем на юг, чтобы в моменты срабатывания сигнал был повеселее. внутри вся камера с фотодатчиком должна быть зачернена. Замер проводим строго в одно и то же время суток по таймеру. Если щель достаточно узка, а датчик от неё далёк, то он по идее должен выдавать пики на частоте средних солнечных суток (по таймеру же), амплитудно промодулированных солнцем (оно будет то входить в щель, то недоходить/переходить. Частота модуляции наверное будет искомыми отклонениями. очевидный минус - погодой эта частота тоже будет промодулирована.

А если попробовать вычесть друг из друга сигналы? За второй сигнал взять первый, но сдвинутый на полгода/год. Будет похоже на биения? И взять Фурье уже от этих биений - может выйдет точнее. Очередная безумная мысль, прошу прощения

Фурье тут скорее всего плохо проходит - низкая точность. Мы например такие вещи параметрическими спектральными оценками считаем.

А еще, обнуление бинов тождественно накладыванию прямоугольного окна, что во временной области даст бесконечный КИХ.

Если мы пропустим сигнал через цифровой фильтр скользящего среднего (FIR фильтр), то выходной сигнал сместится по фазе вправо и реальная фаза окажется искажена.

Это если вы пытаетесь получать результат в реальном времени. Потому что если пытаться усреднять сигнал на каком-то отрезке времени, то первые значения усредненного сигнала будут готовы после окончания этого отрезка. Если вы проводите постобработку, то всегда можно отнести значения отфильтрованного сигнала к "правильным" моментам времени.
То есть сместить их назад на половину периода усреднения.

Надо определить амплитуду и фазу самой первой гармоники.

И что означает "фаза" в данном контексте? Вот у вас есть синусоидальный сигнал, вы в произвольный момент времени начали измерения.. Так фаза в этом случае будет зависеть от того в какую точку синусоиды попало начало измерения. То есть понятие фазы в данном случае бессмысленно. Фаза имеет смысл когда есть два сигнала одинаковой (ну или близкой) частоты - опорный, по которому стартуют измерения и исследуемый сигнал, фаза которого вычисляется относительно опорного.
Ну или вы стартуете измерения в 0:00 часов по местному времени и измеряете освещенность от восхода солнца. Тогда фаза сигнала освещенности будет отсчитываться относительно этих 0:00 часов.
Тут и изображать ФНЧ не требуется. Вы производите измерения в течение 24 часов, получаете выборку, затем через Фурье вычисляете амплитуду и фазу первой гармоники.

И что означает "фаза" в данном контексте? 


Проще говоря, найти момент времени, когда происходит максимум усреднённого сигнала.

Ну или вы стартуете измерения в 0:00 часов по местному времени и измеряете освещенность от восхода солнца. Тогда фаза сигнала освещенности будет отсчитываться относительно этих 0:00 часов.

Нет. Я измеряю освещенность не по местному времени, а по UTC+0. Поэтому фаза как раз меняется при перемещении по параллелям.

Ну это без разницы, всегда можно внести поправку на часовой пояс.
Сутки длятся 24 часа, за это время происходят все события - восход Солнца, макс. высота, заход Солнца. Потом все повторяется. То есть период равен 24 часа.
Но время захода и восхода немного сдвигается каждый день. Максимальная высота Солнца всегда в полдень.
Тогда измерения ведутся 24 часа (длительность выборки 24 часа). Начало отсчета выбирается произвольно, какое требуется. Потом делается Фурье (FFT например, да любое), вычисляется амплитуда и фаза первой гармоники. Потом измерения повторяются, фаза при этом будет наверное ползти, потому что время восхода\захода поменяется.
Освещенность прямо пропорциональна высоте Солнца над горизонтом? Так это табличные значения вроде?
Можно для тестирования модели погонять эти табличные значения, подсунуть их в Фурье анализ.

Но время захода и восхода немного сдвигается каждый день. 

Для неподвижного наблюдателя эти времена сдвигаются в разные стороны. Либо от либо к друг другу.

Потом измерения повторяются, фаза при этом будет наверное ползти, потому что время восхода\захода поменяется.

Вот именно что фаза для неподвижного наблюдателя никуда ползти не должна.

Освещенность прямо пропорциональна высоте Солнца над горизонтом? Так это табличные значения вроде?

Это вообще по астрономическим формулам чисто аналитически вычисляется.

сутки длятся не 24 часа. в этом и суть

А сколько длятся сутки?

... Вследствие того, что орбита, по которой Земля движется вокруг Солнца, не является окружностью, а ось суточного вращения Земли имеет наклон к плоскости орбиты (из-за чего на Земле происходит смена времён года), истинное солнечное время неравномерно. Максимальная разница в продолжительности истинных солнечных суток в течение года составляет примерно 50 с.

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

Верно. Но с другой стороны, облачность может быть и одинаковой в обеих половинах суток и ничего не сместить.

Абсолютно все утверждения из статьи ошибочны и не соответствуют действительности. Автор не виноват, в информационном пространстве вокруг преобразования Фурье магии и городских легенд намного больше, чем вменяемой информации, в том числе и в ВУЗ-ах, преподаватели в которых ни разу не использовали теорию на практике. И теория эта совсем не простая для настоящего понимания её сути, а не просто фантазий на тему.


В частности, преобразования Фурье есть 4 вида:


а) непрерывный сигнал, непрерывный спектр;
б) непрерывный сигнал, дискретный спектр;
в) дискретный сигнал, непрерывный спектр;
г) дискретный сигнал, дискретный спектр;


Теория обычно даётся для а), но численно доступно только г). Отсутствие понимания того, как все эти типы соотносятся, и приводят к большинству заблуждений.


Далее нужно разобраться, что такое свёртка, как она соотносится с корреляцией, с которой многие уже знакомы через мат.статистику, в чём отличие линейной свёртки от циклической, и в чём смысл мантры "свёртка во временной области равно умножению в частотной и наоборот". Да, низкочастотная фильтрация осуществляется умножением высоких частот на ноль или около того.


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


Ну и напоследок — скользящее среднее это не FIR, а его частный случай. И его, между прочим, можно вычислить за O(1).

Абсолютно все утверждения из статьи ошибочны и не соответствуют действительности.

Тем не менее я смог вычислить то, что хотел: нашел фазу усредненного сигнала от освещенности. Значит достаточно понимания, чтобы решать реальные задачи.

Выше уже намекали на то, что подобная формулировка некорректна. Если вы решили какую-то свою задачу - это замечательно, но не отменяет полной оторванности содержания вашей статьи от реальности.

Вот простой пример, как максимум первой гармоники может отличаться от максимума их суммы:

А ведь здесь их только три и даже нет сдвигов по фазе. Поэтому, если ваша задача состояла в поиске фильтрованного максимума — то найденное вами решение, скажем так, слегка неточное.

В эксперименте, что в тексте выделяется фаза только первой гармоники. Она представляет ключевой интерес.

Остальные гармоники не имеют существенного значения. Всё что после k=1 это просто мусор.

Это не мусор и не может быть мусором. Спектральное разложение — это математический инструмент, а не объективная реальность. Объективную реальность дифференциальные уравнения описывают, а здесь у нас есть множество неочевидных нюансов, которые могут входить в конфликт с интуитивным пониманием вопроса. В частности белый шум имеет равномерный спектр — что значит, что его вклад во все частоты одинаков и его нельзя отфильтровать просто НЧ фильтрацией. НЧ-составляющая шума от этого никуда не денется.


Другой момент — сигналы формируются не только аддитивно, но и мультипликативно. Если взять синусоиду одной частоты и промодулировать её амплитуду синусоидой другой, более низкой частоты — Фурье-анализ не покажет эти частоты. Он покажет их сумму и разность, в полном соответствии с математикой и тригонометрией, и сдвиги фаз их относительно друг друга тоже изменятся.


Поэтому в подобного рода задачах используются понятия "мгновенная амплитуда" и "мгновенная частота" (она же фаза). Сначала сигнал фильтруется полосовым (а не нч!) фильтром. Затем он преобразуется в аналитический (комплексный) сигнал через преобразование Гильберта, из которого получаются необходимые мгновенные значения и которые уже можно анализировать привычными способами.

Если у вас отсчёты сигнала идут равномерно и известен период, то проблем-то особых нет. Берёте ДПФ от отрезка сигнала за один период, затем усредняете квадраты амплитуд гармоник по многим следующим отрезкам. Далее измеряете амплитуду первой гармоники. Параллельно усредняете фазу. Чем больше отрезков, тем точнее результат.

Существует ли простой ASIC чип, который с одной стороны получает I2S/PDM интерфейс и аппаратно вычисляет FFT после каждого семпла на линии LR?

То есть вычислять FFT непрерывно в потоковом режиме по скользящему окну.

Чтобы можно было по SPI вычитывать результат работы FFT (массив мощности/ фазы для каждой частоты в спектре входного аудио сигнала).

Скорее всего это какой-то аппаратный аудиокодек с функцией аппаратного FFT. Очень не хотелось бы делать это с нуля на FPGA.

Кто-нибудь сталкивался с подобными микросхемами?

Хочу демодулировать акустический GFSK сигнал. И передавать бинарные данные звуком. Посылать конфиги. Обновлять прошивку через микрофон.

Sign up to leave a comment.

Articles