Comments 94
В статье сначала речь идет о классическом Фурье преобразовании бесконечного сигнала, но в реализации внезапно говорится про прямоугольное окно, то есть это уже другое преобразование (оконное). Это разные преобразования с различными свойствами и применимостью, а при неверном выборе окна можно совершенно фантастические эффекты получить. Очевидно, если размер окна меньше периода исследуемого сигнала, результат не имеет смысла. Кроме того, существуют такие известные проблемы, как явление (эффект) Гиббса, делающий невозможным восстановление исходного сигнала даже в теории, если только сигнал не удовлетворяет определенным ограничениям. В статье же не приведено никакой оценки самого сигнала и валидности выбранных параметров и результатов, так что и полученные результаты могут быть совершенно неверны и построенные на их основании выводы тоже. Если автор хотел показать практическое применение Фурье преобразования, то нужно это делать корректно, а не методом тыка, авось красивый график получится и все равно, имеет ли он хоть какой-то смысл.
Спасибо за интерес к материалам!
Цель этой статьи - поделиться мощным рабочим рабочим инструментом, который был выработан в результате многолетнего труда. Исследовать же его свойства на различных размерах окон и типах сигналов можно самостоятельно, скачав приложение Solfeggio (для Windows) или воспользовавшись исходными кодами.
Это не только красивые графики - это надёжный рабочий метод интерполяции, который по ряду показателей обходит другие подходы и алгоритмы для высокоточного анализа сигналов.
Чтобы частотные пики гармоник сигнала, которым не повезло угодить точно в сетку частот БПФ, не "растекались" по соседним линиям спектра, обычно применяют различные оконные функции. Их придумали в ассортименте: Хемминга, Ханнига, Блекмана, и т.д.
Какую именно функцию лучше использовать - зависит от того, что у вас за сигнал и что вы хотите из него добиться. Одно дело - нарисовать красивую картинку для качественной оценки. Совсем другое - выполнить какие-то измерения и получить по возможности точные результаты, имея оцифрованный сигнал в качестве исходных данных.
Исследование близко расположенных частотных составляющих в сигнале. То, что в статье называется "резонанс". Я бы это назвал не резонансом, а биениями, но не суть. Их можно исследовать, просто увеличивая длительность оцифрованного сигнала. Чем длиннее выборка, тем лучше разрешение в частотной области у преобразования Фурье. Хотите исследовать гармоники сигнала с точностью до 0.001 Гц? Оцифровывайте его в течение 1000 сек, и будет вам счастье. Правда, тут все хорошо в меру. Если вы увеличите длину выборки до миллионов и миллиардов сэмплов, то БПФ уже даст неприятные погрешности.
Чтобы этого не случилось, можно применить такой вычислительный трюк, как "растяжка спектра" (поиск в гугле по словам zoom FFT). Продолжительность вашей выборки (в секундах) он уменьшить не позволит, (в примере выше - записывать сигнал на протяжении 1000 сек все равно придется), но размеры выборки (в сэмплах) сократятся до вполне вменяемых величин. Рекомендую попробовать.
Из трюков и фокусов, связанных с исследовпнием биений, можно упомянуть еще преобразование Гильберта. (При определенных условиях и предположениях его можно использовать для разделения модулированного по амплитуде сигнала, на несущую и огибающую). А биения - это, по сути, и есть амплитудная модуляция. И преобразование Гильберта в данных условиях неплохо работает.
Благодарю за интерес к статье и рекомендации!
ФАИ обеспечивает более высокую точность при прямоугольном окне, чем применение каких-либо других оконных функций.
Качественные результаты получаются не только на сгенерированных идеальных сигналах, но и вполне реальных, например, аудиосигнале с микрофона. Как минимум, по частоте обеспечивается точность гитарных тюнеров при сопоставимых размерах окна.
Автор постоянно упоминает БПФ, но это лишь способ быстрого вычисления коэффициентов дискретного преобразования Фурье. Ну и любые попытки восстановить частоты и/или значения сигнала в промежуточных точках - это попытка высосать информацию из пальца.
Спасибо за интерес к статье!
Возможно, вы пока что не вполне уловили её суть... Рекомендую скачать приложение Solfeggio и сравнить результаты на идеальных и реальных сигналах (например, с микрофона) в результате чистого БПФ и связки БПФ+ФАИ.
Возможно, вы пока что не вполне уловили её суть...
С каким из двух моих утверждений вы не согласны?
Со вторым
Ну и любые попытки восстановить частоты и/или значения сигнала в промежуточных точках - это попытка высосать информацию из пальца.
Это утверждение отчасти может быть верно для белого шума, но не для разрежённых сигналов с чёткими или растёкшимися в результате выполнения БПФ амплитудными пиками.
Это утверждение отчасти может быть верно для белого шума
Т.е. вы утверждаете, что иногда мое утверждение неверно для белого шума? Поскольку у вашей статьи есть тег "Математика", я бы с удовольствием увидел математически корректное доказательство вашего утверждения ;)
но не для разрежённых сигналов с чёткими или растёкшимися в результате выполнения БПФ амплитудными пиками.
Вы понимаете, что чисто гармонический сигнал (Фурье интеграл - дельта-функция) - это когда синусоида с постоянной амплитудой в произвольный момент времени от минус до плюс бесконечность?
Не берусь ничего утверждать, но, возможно, ваше утверждение категорически верно для идеального белого шума, но в остальных случаях лишь отчасти, с определённой долей вероятности и погрешности.
У меня прикладной взгляд в данном вопросе - этот алгоритм с высокой точностью работает на реальных сигналах, даже если кто-то говорит, что это невозможно ;)
Всем скептикам рекомендую скачать приложение и воочию протестировать ФАИ на реальных сигналах.
Исходные коды полностью открыты, можно удостовериться, что больше никакой дополнительной обработки к сигналам не применяется.
Всем скептикам рекомендую скачать приложение и воочию протестировать ФАИ на реальных сигналах.
Исходные коды полностью открыты, можно удостовериться, что больше никакой дополнительной обработки к сигналам не применяется.
Вы извините, но восстанавливать математический алгоритм по коду я не буду. И использовать код, который непонятно как работает (или врет) тоже не буду. Еще раз напомню, что это вы поставили тег "математика". Где она? А на нет и суда нет. Ну если хотя бы то была библиотека на C++ или пакет на Питоне, чтобы можно было нормально протестировать так чтобы глаза не вытекли....
Ну и любые попытки восстановить частоты и/или значения сигнала в промежуточных точках — это попытка высосать информацию из пальца.Это не так, теорема Котельникова как раз и рассматривает однозначное восстановление непрерывного сигнала из набора дискретных отсчётов.
Тут конечное число отсчетов и еще нет ограниченного спектра.
Первое - лишь приблизительно, Второе - где цикличность у затухающего сигнала?
Первое — с чего это вдруг?
Я с удовольствием увижу физически реализуемый фильтр, убивающий 100% частот, выше заданной.
Второе — цикличность не у самого исходного сигнала, а у мат. аппарата, его рассматривающего — точно там же, где цикличность суммы конечного числа синусоид с кратными частотами. Достаточно вычислить обратное ДПФ аналитически, а не численно, чтобы убедиться в этом.
Так, давайте на конкретном примере. Пусть мы измерили сигнал 4 раза:
t=0, x=0
t=1, x=1
t=2, x=2
t=3 x=3
Чему будет равен сигнал при t=4? А при t=1.5? А ведь это еще без шума тракта и шума дискретизации.
Я с удовольствием увижу физически реализуемый фильтр, убивающий 100% частот, выше заданной.А 100% и не надо, достаточно ниже шумов квантования и прочих, потому что при всём желании вы не можете численно работать с неограниченной точностью. Ну и при любом раскладе в физическом мире верхняя частота любого сигнала ограничена скоростью света.
Так, давайте на конкретном примере. Чему будет равен сигнал при t=4? А при t=1.5?Вы в этом примере очень хитро подменили реальный сигнал и его мат.модель — предсказывать будущее ещё никто не научился.
Ну а в мат.модели — если вы уверены, что все прочие х равны нулю, а сигнал удовлетворяет требованиям Котельникова — то:
Ну или вы знаете, что все прочие x равны -1:
Ну или вы знаете, что сигнал периодичен и точно знаете период — то например:
Подставляйте интересующее вас x в функцию, получите конкретное значение.
А 100% и не надо, достаточно ниже шумов квантования и прочих, потому что при всём желании вы не можете численно работать с неограниченной точностью. Ну и при любом раскладе в физическом мире верхняя частота любого сигнала ограничена скоростью света.
Вы представляете сложность даже такого фильтра? Ну а если я хочу смотреть на нефильтрованный сигнал? Есть померянные с некоторой точностью отсчеты, а между ними - неизвестность. Хотите точнее - собирайте больше данных. Чаще, дольше. Та же ситуация, что и в мат.статистике.
Вы в этом примере очень хитро подменили реальный сигнал и его мат.модель — предсказывать будущее ещё никто не научился.
Я там и про t=1.5 спросил. И по Вашим графикам (Mathematica?) видно, что при разных предположениях о сигнале, даже при ограниченном спектре, в промежуточных точках получаются разные предсказания. Также видно, что идеальная фильтрация сильно искажает форму сигнала x(t)=t.
что при разных предположениях о сигналеНу это логично, а как иначе-то? Если брать сигнал с магнитофона после кнопки «Вкл.» и до кнопки «Выкл.», то вполне логично дополнять его нулями (то бишь тишиной) с обоих сторон. Если у микрофона мембрана физически не в состоянии реагировать на частоты выше определённой — то вполне логично считать сигнал с него спектрально ограниченным.
Также видно, что идеальная фильтрация сильно искажает форму сигнала x(t)=tОчевидно что сигнал с бесконечным спектром под критерии Котельникова не подходит, и для работы с таким сигналами используются другие мат.модели и другие инструменты. Равно как и в реальном физическом мире вы такой сигнал ничем воспроизвести не сможете — из-за физических ограничений элементной базы.
Ну это логично, а как иначе-то?
Я же не спорю, это автор пытается добыть из данных дополнительную информацию, изначально там отсутствующую.
Ничего автор не добывает.
Алгоритм в результате наоборот убирает избыточную и малозначимую информацию, с довольно высокой точностью реконструируя идеальный сигнал.
На сильно зашумленных сигналах он, конечно, чудес не сотворит, но с реальными (музыка, речь) справляется достойно.
У вас при обратном Фурье-преобразовании значения в исходных точках восстанавливаются точно?
В заключительной части статьи добавлен раздел "песочница", который содержит ссылку на реализацию алгоритма в онлайн-компиляторе. Там можно самому поиграть с данными и получить непосредственное представление, какая информация теряется, а какая восстанавливается, и насколько точно.
Вы представляете сложность даже такого фильтра?
Это уже давно не проблема для АЦП, работающих в звуковом диапазоне. АЦП с оверсемплингом, после АЦП - цифровой фильтр. Перед АЦП - обычный аналоговый фильтр.
Это уже давно не проблема для АЦП, работающих в звуковом диапазоне. АЦП с оверсемплингом, после АЦП - цифровой фильтр. Перед АЦП - обычный аналоговый фильтр.
На самом деле @Refridgerator, как я только сейчас понял, меня обманул. Если поставить физически реализуемый фильтр перед АЦП, после АЦП, то да, финальный сигнал можно полностью восстановить. Но это не изначальный сигнал!
А что мешает выполнить преобразование Фурье не для дискретного набора частот (как в БПФ) а для нужного набора частот (скажем от 0 до fs/2 с шагом df)?
Вычислительная сложность. БПФ значительно сокращает объёмы вычислений и расходы памяти (то есть обладает наименьшей избыточностью по данным показателям, чем другие методы вычисления преобразования Фурье), при этом обеспечивая обратимость преобразования (возможность однозначного восстановления исходного сигнала до точности округлений).
В теории можно рассчитать коэффициенты преобразования Фурье с огромной точностью, но для этого потребуется много вычислительных ресурсов, что накладывает серьёзные ограничения на их вычисление в режиме реального времени.
Тут основное ограничение - на длину выборки сигнала.
Дискретное преобразование Фурье устроено так, что если мы хотим исследовать в подробностях составляющую сигнала с частотой f, то в нашей выборке обязано находиться целое количество периодов такого сигнала, который имеет частоту f.
Второе ограничение связано с БПФ.
БПФ - это, по сути, хитрый вычислительный трюк, который позволяет очень быстро выполнить дискретное преобразование Фурье, при выполнении очень мааааленького дополнительного ограничения: длина выборки сигнала обязана быть степенью числа 2. Но зато преобразование выполнится ну ооочень быстро.
То есть если выборка у нас получилась длиной 1024 смп, то все отлично, вычисляем спектр за сущие микросекунды. А если (не дай бог) для того, чтобы исследовать в подробностях интересующую вас в сигнале частоту f , приходится брать выборку в 1025 или в 1023 смп, то происходит страшное... Бпф не работает, и для вычисления такого же спектра требуется в десятки раз больше времени.
Значит, способов решения проблемы тут два.
подобрать частоту квантования АЦП так, чтобы нужная частота сигнала уложилась в подходящее для бпф количество сэмплов. Это не всегда удобно и не всегда возможно технически.
Использовать БПФ только для получения общей картины спектра сигнала. А для подробного и детального исследования его составляющих лучше применять ДПФ. Да, это конечно, будет медленно. Но вам ведь не нужно вычислять все коэфф Фурье, весь их полный набор для вашей выборки сигнала! Достаточно вычислить только один единственный коэффициэнт, который соответствует интересующей вас сейчас частоте f. Или набору частот (их все равно будет сильно меньше, чем длина выборки) Для современных компьютеров - вполне посильная задача. Тогда и не придется умничать с интерполяциями. С точки зрения математики, интерполяция, экстраполяция - это, по сути, игра в угадайку. Иногда это может сработать и вы угадаете подходящее значение (если окажутся верны исходные предположения о сигнале). А может так случиться, что и не угадаете.
Если изначально выборка удовлетворяет теореме Найквиста-Шеннона-Котельникова, то её длину можно искусственно варьировать, интерполируя сигнал во временной области и выполняя передискретизацию с большей частотой, что будет приводить к возрастанию спектрального разрешения за счёт увеличения объёмов вычислений и потребления памяти.
Преимущество ФАИ в том, что алгоритм позволяет весьма точно «угадать» характеристики сразу нескольких достаточно отдалённых гармоник сигнала при малом объёме дополнительных вычислений.
Другими словами, ФАИ - это метод высокоточного сжатия разрежённых сигналов с минимальными потерями полезной информации.
передискретизацию с большей частотой, что будет приводить к возрастанию спектрального разрешения
Ну, вообще, для увеличения спектрального разрешения надо увеличивать длительность выборки. Увеличением частоты дискретизации можно расширить спектр в область высоких частот.
А так да, можно вырезать из выборки кусок, чтобы туда укладывалось целое число периодов сигнала, а потом передискретизировать\интерполировать число точек к ближайшей степени двойки.
Спасибо, верно подметили!
Не совсем корректно выразил мысль: скажем, у нас есть выборка из 1024 отсчётов, передискретизируем её с частотой в два раза выше, чтобы получить уже 2048 отсчётов. При этом размер выборки и спектральное разрешение удвоится, но временная длительность выборки останется прежней.
Ну вот если под спектральным разрешением понимать величину шага между гармониками, то она не увеличится. Если длительность выборки 1 сек, то шаг между гармониками 1 Гц. Если у нас при этом 1024 отсчета в выборке, то в спектре мы сможем получить 512 гармоник с шагом 1 Гц. Если мы увеличим частоту дискретизации в 2 раза, так что на выборке длит. 1 сек мы получим 2048 отсчетов, то в спектре сможем получить 1024 гармоники, но с тем же с шагом 1 Гц. То есть увеличение числа отсчетов (частоты дискретизации) при сохранении длительности выборки просто расширит спектр в область верхних частот. Но шаг между гармониками останется прежним, то есть спектральное разрешение не улучшится.
Чтобы улучшить спектральное разрешение надо увеличивать длительность выборки. Если длительность будет 10 сек, то шаг между гармониками будет 0.1 Гц.
Ну это понятно, за период первой гармоники в спектре берется длительность выборки. То есть на длительности выборки первая гармоника уложится 1 раз, вторая - 2 раза и т.д.
Ещё можно отметить, поскольку ФАИ позволяет довольно точно восстановить гармоническую векторную составляющую разрежённого сигнала (получить исходные коэффициенты косинусоид близкие к аналоговым), на этой основе легко выполнить передискретизацию данного сигнала в произвольное число раз (даже нецелое), а также при необходимости частотную фильтрацию или даже точечную модификацию конкретных гармонических составляющих (приглушить их или усилить, сдвинуть по частоте и даже фазово).
То есть открывается большая свобода действий по обработке сигнала...
БПФ - это, по сути, хитрый вычислительный трюк, который позволяет очень длина выборки сигнала обязана быть степенью числа 2. Но зато преобразование выполнится ну ооочень быстро.
Не-не. Главное, чтобы не простое. Но, да, чем лучше факторизуется, тем быстрее.
Вдогонку, можно даже и простое: https://en.wikipedia.org/wiki/Fast_Fourier_transform
А если (не дай бог) для того, чтобы исследовать в подробностях интересующую вас в сигнале частоту f, приходится брать выборку в 1025 или в 1023 смп, то происходит страшное…Ничего страшного не происходит, для этого есть prime-factor и Bluestein алгоритмы. Они просто сложнее в реализации, Bluestein дополнительно требует увеличения точности для больших n (>50000).
Частоту сигнала намного точнее можно вычислять при помощи автокорреляции (тоже через БПФ). Пик автокорреляции даст период с точностью до нескольких отсчётов, а уж период в частоту перевести легко. Результат получается на порядки точнее, чем по спектру, особенно для низких частот. Собственно, чем ниже частота, тем выше точность по автокорреляции и ниже по спектру. Если, к примеру, писать программу-тюнер для гитары, то я бы советовал использовать автокорреляцию.
Есть и подводные камни:
Нужно использовать дополнительные эвристики для отсеивания пиков на гармониках.
Для смешанного сигнала (несколько частот) понадобится сначала отделить нужную частоту при помощи фильтра частот.
(Без претензий к ТС, просто вдруг кому-то пригодится.)
Спасибо за отклик и интерес к статье!
Даже на низких частотах ФАИ обеспечивает точность сравнимую с автокорреляционными музыкальными тюнерами при эквивалентных длинах выборки (сам проверял). Возможно, какая-то разница будет наблюдаться на сверхнизких частотах, но это требует дальнейшей проверки и исследования.
может быть, существует какой-то метод интерполяции в частотной? Хотя на момент написания статьи математического обоснования этой гипотезе нетСуществует уже почти 100 лет — читайте про sinc и periodic-sinc функции, алгоритмы передискретизации, ну и вообще любой доступный учебник по обработке сигналов.
Интересное наблюдение: резонанс очень напоминает амплитудную модуляцию низкочастотного сигнала высокочастотным, например, резонанс 440 и 444 герц похож на 4 герцевый сигнал (или частоты такого порядка), модулированной частотой примерно в 442 герца. Данное предположение требует исследования, уточнений и дальнейшей проверки.Данное предположение описывается и доказывается в школьном курсе тригонометрии (разложение произведения синусов на сумму и наоборот).
Спасибо за интерес к статье!
Не знаю, в какой области работают sinc и periodic-sinc функции, но алгоритмы передискретизации точно работают во временной (применяются к исходной выборке), речь же идёт про алгоритм интерполяции в частотной области (применяется к спектру после БПФ).
Я почему хотел увидеть алгоритм автора статьи, чтобы хотя бы убедиться, что при обратном преобразовании сигнал остался неизменным или хотя бы вещественным.
Ну а поскольку восстановление исходного сигнала в задачи автора не входило, то это у него и не получится.
Ну а поскольку восстановление исходного сигнала в задачи автора не входило, то это у него и не получится.
Сигналы с разрежёнными гармониками вроде звуковых восстанавливаются с довольно высокой точностью, при этом происходит очистка от шумов (малозначимых гармоник).
Да, это неоднозначное восстановление с потерями, но и преимущество для определённых задач.
Так автор только сейчас (см.ниже) свою весьма специфическую цель обозначил. А, то, скачайте, попробуйте...
Алгоритм представлен в статье. С разрежёнными гармониками в идеальных и реальных сигналах он справляется очень даже неплохо.
Вы бы его еще на brainfack'е написали :)
Но... Чем вам C# не угодил? Это современнный (не побоюсь этого слова, передовой) язык программирования широкого предназначения. Синтаксически близок к Java, во многих аспектах опережает её, также есть немало общих черт с C++.
Да и алгоритм не настолько сложен (по сути, около десятка формул), чтобы даже интуитивно при внимательном рассмотрении с ним не разобраться. Если же есть какие-то затруднения или вопросы, то можно спросить, уточнить...
:)
Я хотел бы только заметить, что тот мой комментарий, где я в первый раз отказался восстанавливать алгоритм из кода, понравился еще двум читателям. Т.е. или нам троим заниматься reverse engineering'ом, или вам эти формулы опубликовать. Или у вас их нет, а есть только код?
Код содержит все ключевые формулы.
:)
Добавил ссылку в разделе «песочница» на онлайн-компилятор с реализацией алгоритма, где можно его исследовать на различных сигналах.
Лучший и худший - понятия довольно относительные, которые зависят от критериев оценки.
Первоочередная цель при создании приложения Solfeggio - точное распознавание нот из звуковых сигналов в режиме реального времени. С этой задачей ФАИ справляется на уровне автокорреляционных тюнеров для музыкальных инструментов, при этом в вычислительном плане метод не требователен к ресурсам.
Возможно, существуют ещё более совершенные подходы для решения подобного рода задач, но мне они не известны.
Первоочередная цель при создании приложения Solfeggio — точное распознавание нот из звуковых сигналов в режиме реального времени. С этой задачей ФАИ справляется на уровне автокорреляционных тюнеров для музыкальных инструментовВот именно, что «распознавание нот» и «тюнер для музыкальных инструментов» — это совершенно разные задачи, и в обоих случаях есть более совершенные в плане точности алгоритмы, в том числе и без FFT вообще. FFT даёт линейную шкалу частот, в то время как распределение частот в нотах логарифмическое. Для распознавания нот банк узкополосных фильтров с выделением огибающей через преобразование Гильберта будет точнее.
Для распознавания нот банк узкополосных фильтров с выделением огибающей через преобразование Гильберта будет точнее.
Большой вопрос к вычислительной сложножности данного подхода... А точность можно существенно повысить и для БПФ, опять же за счёт повышения требований к вычислительным ресурсам.
«распознавание нот» и «тюнер для музыкальных инструментов» — это совершенно разные задачи
Тюнер - это точное распознование одной наиболее значимой ноты в сигнале, то есть частный случай более общей задачи по точному распознованию множества нот.
А в распознавании нот частотная точность вообще не важна, потому что частотные модуляции возникают как при игре на музыкальных инструментах, так и при различных интерференционных эффектах. А что действительно важно, так это:
1) последовательность нот;
2) отсутствие ложных срабатываний.
Причём на входе мы имеем полифонический сигнал с кучей шума, где одна мелодия на одном инструменте звучит поверх другой на другом инструменте поверх гармонической аккордовой последовательности на третьем. Вот, например, реальный спектр реального музыкального сигнала длиной 8192 отсчётов:
Какие ноты тут звучат, в какой последовательности, и сколько их? Интерполяция здесь вообще ничего не решает.
На всякий случай уточню, изначальная задача для приложения Solfeggio имеет наиболее общий характер: максимально точное распознавание частоты всех нот в аудиосигнале при коротких временных размерах окон (чтобы разобрать порядок звучания в реальном времени).
Например, если музыкант играет на клавишном инструменте, то программа должна показывать с достаточно высокой точностью, на какие клавиши он сейчас нажимает.
Для отдельно звучащих инструментов алгоритм справляется весьма неплохо. Помимо прочего программу можно использовать в качестве продвинутого музыкального тюнера, который позволяет производить настройку инструмента не только по одному обертону, но и сразу по нескольким.
Большой вопрос к вычислительной сложности данного подхода...Вычислительная сложность не имеет никакого значения, если вы не знаете, как именно это сделать. Здесь сложность прежде всего алгоритмическая, которая на несколько порядков выше, чем написать FFT.Transform(data).
Добавлен раздел «песочница» со ссылкой на онлайн-компилятор, где можно самому пощупать алгоритм и поиграть с входными данными сигнала.
Мне будет очень интересно познакомиться с альтернативными подходами, которые обеспечат примерно такой же качественный уровень распознавания при сопоставимой алгоритмической сложности...
Не-не, имею в виду готовую реализацию на каком-либо языке программирования, чтоб можно было запустить и сравнить.
Так и сам могу кого угодно отослать по умным запросам, мол, изучайте, реализуйте и будет лучше, чем у вас. Поверьте, за годы исследований довелось познакомиться с немалым количеством материалов по теме обработки сигналов, но во многих случаях это теории или частные случаи. Так, чтобы рабочий код, это редкость.
Да, в открытом доступе есть реализации гитарных тюнеров на основе автокорреляции, но и они не лишены недостатков.
В статье представлен рабочий код, так что от «конкурентов» и альтернатив ожидается того же ;)
Благодарю за ссылку. Для определённых задач вейвлеты очень даже годный инструмент, но вычислительно он на порядки дороже БПФ (из моих изысканий), что ограничивает их область применения, в особенности при обработке сигналов в реальном времени.
Не оспариваю тот факт, что различными подходами можно увеличивать точность распознавания гармоник, но опять же это требует дополнительных ресурсов.
Из моего опыта, связка БПФ и ФАИ - это в практическом смысле наиболее оптимальный путь достижения довольно высокой точности распознавания при минимальных затратах вычислительных ресурсов.
Возможно, если эту связку использовать для подготовительного анализа сигнала, а потом уже уточнять свойства гармоник другими методами, то можно будет достичь и ещё более впечатляющих результатов...
Но для моих задач точности хватало, поэтому такими экспериментами не занимался. Проект же развивается на энтузиазме, мне просто интересно разбираться в этой теме, поэтому делюсь результатами в открытом доступе :)
Используемая вами реализация FFT, кстати, на больших размерах будет врать из-за накопления ошибок в роторе.
Не исключаю, что если задействовать векторные инструкции процессора (как в FFTW), то возможно достичь большей скорости вычислений на тех же объёмах.
Плюс в теории можно выполнить предрасчёт ротора, чтобы закэшировать значения (правда, в данном случае нужно проверять на конкретном железе, поскольку перерасчёт значений в процессоре может оказаться быстрее, чем обращения за ними к оперативной памяти).
магическим математическим смыслом
Меня очень цепляет этот загадочный смысл. Не могу строго объяснить эти во многом интуитивные формулы и почему алгоритм работает, но он работает и весьма точно!)
---
Используемая вами реализация FFT, кстати, на больших размерах будет врать из-за накопления ошибок в роторе.
Вместо умножений там можно возводить в степень, будет чуть медленнее, но, предположительно, точнее
//rotor *= rotorBase;
rotor = Complex.Pow(rotorBase, i + 1);Не исключаю, что если задействовать векторные инструкции процессора (как в FFTW), то возможно достичь большей скорости вычислений на тех же объёмах.Не только поэтому — c# очень медленно работает с типом Complex, потому что это отдельный класс, а не встроенный тип данных. Ну и алгоритмическая оптимизация — Radix-4, Split-Radix, etc.
Вместо умножений там можно возводить в степень
Меня очень цепляет этот загадочный смысл. Не могу строго объяснить эти во многом интуитивные формулы и почему алгоритм работаетЗагадочный смысл где-то возле линейной интерполяции и точки пересечения двух прямых.
Загадочный смысл где-то возле линейной интерполяции и точки пересечения двух прямых.
Линейной, да не совсем. Если бы всё было настолько линейно и очевидно, то, скорее всего, способ был бы давно замечен и широко известен в открытом доступе.
Ведь существует немало альтернативных публичных подходов к интерполяции: амплитудная параболическая и по Гауссу, фазовый сдвиг при смещении кадров, различные классы оконных функций...
Но для моих задач соотношение точности и производительности для ФАИ наиболее оптимальное.
С некоторыми вариантами интерполяции можно ознакомиться по ссылкам
Improving FFT Frequency Measurement Resolution by Parabolic and Gaussian Interpolation (cern.ch)
Guitar Pitch Shifter - Algorithm
Pitch Shifting Using The Fourier Transform | Stephan Bernsee's Blog (zynaptiq.com)
Линейность - многозначительный термин: если имеется в виду линейная функция от многих переменных, то не вполне уверен, что все формулы удовлетворяют определению таких функций. Но это терминолгические тонкости, сути дела они не меняют :)
В математике много недоказанных гипотез разной сложности. Некоторые проверены вычислительно до огромных чисел, но всё ещё остаются гипотезами.
Причём, из теоремы Гёделя о неполноте следует, что вполне могут существовать и неразрешимые утверждения: их нельзя ни доказать ни опровергнуть, ни даже узнать наверняка, разрешимые они или нет...
Слабое место математики: можно ли доказать всё, что истинно? [Veritasium] - YouTube
Что же касается нашего случая, то обоснование для него может поискать любой желающий, возможно, оно существует ;) отчасти для этого и создана данная статья. Кроме того, у меня есть подозрения, что некоторые формулы можно даже слегка улучшить для большей точности, тоже вопрос открытый...
Ой, зря вы так, некоторые ролики дают визуализацию и представление о сути таких вещей, которую по формулам уловить чрезвычайно трудно. Это как лекции с очень качественным видеорядом. Есть ряд математических каналов, при просмотре которых мой интерес к красоте этой области знаний очень усилился...
А теорема Гёделя имеет фундаментальное отношение ко всем разделам математики.
Теорема Гёделя возникла из чисто философского вопроса об истинности математики в целом. Инструменты мат.анализа из неё никак не выводятся.
Меня до сих пор вдохновляют некоторые ролики о комплексных числах... Чего только стоит серия "everything is rotation"
и другие по запросу
Во многом благодаря таким роликам в программе Solfeggio появился модуль по геометрии звука (2D цветок и 3D лента)
В своё время курс по геометрическому введению в теорию групп и другие ролики Алексея Савватеева очень поспособствовали развитию моего математического мышления или представления, кто как назовёт :-)
Геометрия и группы. Алексей Савватеев. Лекция 1.1. Мотивирующий пример - YouTube
В общем-то алгоритм рисования геометрий в Solfeggio работает на том же принципе, что и в ролике, только он переводит звук в изображение. Единственное, на данный момент не пробовал реализовывать обратное преобразование, чтобы перевести визуальное изображение в звук.
Если рассмотреть в 3D, то окружность - это ортографическая проекция спирали (пружинки) с торца, а сбоку получится обычная косинусоида/синусоида. То есть полокружность - это половина периода косинусоиды.
Мы привыкли, что комплексные числа - это что-то лишь про комплекную плоскость, но на самом деле более общими являются эйлеровы числа (тригонометрическая или показательная формы комплексных чисел), которые за счёт их периодичности позволяют задавать координаты на трёхмерных объектах вроде бесконечно протяжённых спиралей.
Классическое комплексное число в алгебраической форме - это ортографическая проекция на комплексную плоскость множества эйлеровых чисел с одинаковым периодом.
По-моему, интересная тема для исследований... ))
В общем-то алгоритм рисования геометрий в Solfeggio работает на том же принципе, что и в ролике, только он переводит звук в изображение.Не-не, вы мне зубы не заговаривайте, тут задача прямо противоположная. Чтобы из звука один канал направить на координату X, а другой на Y, много ума не надо.
То есть полокружность — это половина периода косинусоидыЭто дуга, в полуокружности ещё замыкающая линия с другой стороны.
Существует метод рисования, когда один канал направляют на ось X, а другой на Y.
Но тут принципиально другой подход и работа с моносигналом. Каждая гармоника - это три параметра: амплитуда, фаза и частота, то есть набор векторов, которые вращаются с разной угловой скоростью при их временном сдвиге по оси Z. При суммировании и продвижении вглубь, они образуют различные трёхмерные спиралевидные фигуры необычных очертаний (некоторые их торцевые ортографические и перспективные проекции представлены на картинках).
Для единичной синусоиды получается трёхмерная спираль-пружинка, растяжение и число оборотов которой зависит от соотношения длины временного интервала и частоты, радиус от амплитуды, а начальный поворот от фазы гармоники.
В вашем случае получилось довольно сложное выражение, поскольку полуокружность замкнута линией (чтобы описать такое нужны высокочастотные компоненты), но если её убрать, то полуокружность или дуга любой длины опишется единичной гармоникой (при верном подборе временного интервала к частоте).
код метода для рисования геометрий
public static IEnumerable<Point> DrawGeometry(Bin[] peaks, int sampleSize, double sampleRate,
double centerX = 0d, double centerY = 0d)
{
var radiusMax = Math.Max(Math.Min(centerX, centerY), 1d);
for (var i = 0; i < 4 * sampleSize; i++)
{
var a = 0d;
var b = 0d;
for (var j = 0; j < peaks.Length; j++)
{
var peak = peaks[j];
var w = Pi.Double * peak.Frequency;
var t = 0.5 * i / sampleRate;
var phase = peak.Phase + w * t;
var radius = peak.Magnitude * radiusMax;
a += radius * Math.Cos(phase);
b += radius * Math.Sin(phase);
}
yield return new(a + centerX, b + centerY);
}
}
если её убратьСтало быть, ничего сложнее дуги окружности вы смоделировать в состоянии. Что, собственно, и подтверждает мой исходный тезис.
В вашем случае получилось довольно сложное выражениеНикакое оно не сложное — там просто 40 коэффициентов, которые получены из обобщённой формулы, на получение которой ушло пара минут. Так-то понятно, что ряд получится бесконечным.
Стало быть, ничего сложнее дуги окружности вы смоделировать в состоянии. Что, собственно, и подтверждает мой исходный тезис.
Вы сами себе формулируете узкие тезисы, ищете им подтверждения и этим ограничиваетесь. Просили дугу - её и получили, каких-то других запросов не поступало.
Кроме того, на скриншотах есть и более сложные варианты фигур, чем дуга, но их почему-то вы не берёте в расчёт.
Благодарю за видео, красивое! Правда, насколько понял, метод рисования там двухканальный, а не вращательный, как в серии роликов “everything is rotation” и программе Solfeggio.
Не умаляю ваших школьных достижений, но, во-первых, задача состояла в математически обоснованном преобразовании звука в геометрические объекты, во-вторых, это далеко не метод тыка, хоть и матчасть не слишком-то и сложная.
Ещё хочется отметить, что статья создана без цели с кем-то посоревноваться в понимании матчасти или что-то важное доказать другим) главные её задачи - привлечь в внимание к теме и поделиться текущими наработками в виде доступном для читателей разного уровня: от школьника до специалиста.
Также мне было интересно подвергнуть материалы внешней критике. Например, благодаря вашему участию были замечены небольшие огрехи в коде. Из дискуссий тоже сформировал для себя некоторые выводы. Так что большое спасибо вам лично и другим участникам за уделённое время и внимание! Как минимум для меня это ценно.
Применять ли информацию из статьи дальше в практической деятельности, думаю, читатель уже решит сам.
В заключительной части статьи добавлен раздел ПЕСОЧНИЦА, который содержит ссылку на реализацию алгоритма ФАИ в онлайн-компиляторе, где можно можно поиграть с ним на различных сгенерированных сигналах
var harmonics = new Harmonic[]
{
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 9.00 },
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 3900.00 },
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 4500.00 },
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 7000.00 },
};F = 15.62 Hz, M = 0.91, P = -1.43
F = 3900.01 Hz, M = 1.00, P = -0.98Если слегка уплотнить спектр, то тоже начинает заметно ошибаться:
var harmonics = new Harmonic[]
{
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 900.00 },
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 950.00 },
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 1000.00 },
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 1050.00 },
new() { Magnitude = 1.0, BasisFunc = Cos, Phase = 0.0, Frequency = 1100.00 },
};F = 899.17 Hz, M = 1.02, P = 0.51
F = 953.12 Hz, M = 1.01, P = -0.25
F = 1000.00 Hz, M = 1.12, P = 0.39
F = 1046.88 Hz, M = 1.00, P = 1.03
F = 1100.97 Hz, M = 1.02, P = 0.24С плотными спектрами действительно начинают возникать погрешности, поскольку спектральные утечки начинают заметно накладываться друг на друга, постепенно переходя в резонансные биения. Но даже при таких условиях результаты зачастую точнее, чем при обычном БПФ.
По высокочастотным компонентам уточню чуть позже, возможно, есть ограничения или что-то не учтено.
Спасибо, что глянули! Там в методе EnumeratePeaks небольшая ошибочка закралась.
Вместо
var count = spectrum.Count / 2 - 3;нужно
var count = spectrum.Count - 3;тогда высокочастотные компоненты не обрезаются. Подправлю это!
Для работы с плотными спектрами рекомендуется уплотнить шаг частотной решётки либо путём понижения частоты дискретизации, либо увеличением размера выборки. Например, для частоты 16000 Гц и длины выборки в 1024 отчсёта шаг решётки составляет 15.625 Гц, а длительность 64 мс, для комбинаций 8000/1024 и 16000/2048 уже 7.8125 Гц и 128 мс.
Тогда на ваших данных погрешности по частоте не превышают 0.2 Гц.
спойлер
F = 899.91 Hz, M = 0.99, P = 0.73
F = 949.90 Hz, M = 1.01, P = 0.79
F = 1000.00 Hz, M = 1.00, P = 0.78
F = 1050.13 Hz, M = 1.01, P = 0.77
F = 1100.10 Hz, M = 0.99, P = -0.34Там в методе EnumeratePeaks небольшая ошибочка закралась.
Вот поэтому я и не стал восстанавливать математический алгоритм по коду.
var count = spectrum.Count / 2 — 3;Не нашёл такого в песочнице. Ну ладно, пытаем дальше, уменьшаем размер блока:
var sampleRate = 48000d;
var frameSize = 32;
var harmonics = new Harmonic[]
{
new() { Magnitude = 1.00, BasisFunc = Cos, Phase = 1.00, Frequency = 2444.00 }
};↓
F = 2325.56 Hz, M = 0.96, P = -1.26Погрешность больше 100 Гц.
Судя по тому, что тут график, а не точки, это интеграл, а не ряд?
Погрешность больше 100 Гц.
Бывает, что онлайн компилятор глючит при обновлении кода, но для ваших входных данных у меня получился немного другой результат
F = 2423.47 Hz, M = 0.71, P = 0.19То есть погрешность составила около 22 Гц.
Это моё субъективное мнение, но для шага дискретной решётки в
48000 / 32 = 1500 Гц (!)Погрешность в 100 Гц ещё терпима, а 22 Гц - это очень неплохой результат.
Можно понизить часоту дискретизации, скажем, в 6 раз при прежнем размере блока
8000 / 32 = 250 Гц
F = 2447.52 Hz, M = 1.00, P = 1.56Тогда для шага решётки в 250 Гц погрешность составит порядка 22 / 6 ~ 3.67 Гц
Загадки быстрого преобразования Фурье