Pull to refresh

Comments 21

Было бы интересно узнать результаты прослушивания на слышимость артефактов. Например, пре-эхо.

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

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

Можно сравнить визуально. Возьмём высокодобротный фильтр низких частот,
АЧХ

Тогда его фазолинейная импульсная характеристика будет выглядеть так

а минимально-фазовая — так:


Видно, что в первом случае и фронт более резкий, и пульсации затухают сильнее.

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

Не там искали.

КИХ фильтры в обработке звука обычно не используют, так как ухо нечувствительно к фазе сигнала, только в редких случаях прибегают к Allpass фильтрам для глубокого баса, но это редкость и мало кто это способен хоть как-то услышать.

Спасибо за интересный и полезный материал!

Позвольте пару вопросов от неспециалиста:

  1. Как эффективно применять полученный КИХ-фильтр?
    Простейшая реализация - прямая свёртка входного сигнала с фильтром. Это даёт почти нулевую задержку, но само вычисление занимает слишком много времени. Ускорить вычисление можно с помощью БПФ (умножая Фурье-образ фрагмента входного сигнала на Фурье-образ фильтра). Но алгоритм БПФ предполагает, что входной сигнал зациклен, а в нашем случае это не так, т.к. мы просто выдёргиваем произвольный фрагмент из временного ряда. В результате такой фильтрации конец и начало фрагмента будут "смешаны" друг с другом. Как этого избежать? Использовать оконную функцию (например, Ганна)? А как потом восстановить фильтрованный сигнал из фрагментов?

  2. Пробовали ли Вы в качестве функции пропускания (для АЧХ) использовать интеграл атомарной функции? Она обладает двумя замечательными свойствами: её производная представляет собой сумму двух таких же функций (с точностью до масштаба и смещения); в точках "стыковки" все её производные равны нулю, т.е. стык получается идеально гладким.

1.а. Чтобы сделать свёртку через FFT, не нужны оконные функции, а нужен метод Overlap–add. Его идея в следующем: допустим, у нас FIR из 3-х отсчётов, а использовать будем FFT на 8 отсчётов. Тогда мы этот FIR дополняем 5-ю нулями до 8-и, а входной сигнал делим на порции из 5-и отсчётов, каждый из которых дополняем 3-мя нулями (тоже до 8). Спектры перемножаем, и результаты складываем с нахлёстом в 3 отсчёта.

1.б. Для свёртки без задержки есть алгоритм zero-delay convolution, идея которого тоже проста. Допустим, у нас FIR в 16 отсчётов. Тогда мы его разбиваем на 4,4 и 8 отсчётов. Первые 4 вычисляем прямой свёрткой, остальные — через FFT (что даст задержку соответственно в 4 и 8 отсчётов). Результаты складываем.

2. Нет, не слышал о таких функциях.

Спасибо за интересную и доступную подачу материала. Раз вы упомянули про QMF, может знаете что кроме научных статей можно почитать про расчет FIR фильтра для полифазного фильтрбанка?

Почитать не могу ничего посоветовать, кроме того что искать надо на англоязычных источниках. А у вас есть какая-то конкретная задача или просто интересно?

Ух. Есть у меня pet project - открытая реализация atrac кодека. Для atrac3plus нужно реализовать анализирующий фильтрбанк для вот такой реализации синтезирующего https://github.com/FFmpeg/FFmpeg/blob/master/libavcodec/atrac3plusdsp.c#L502-L645

Попытка "с наскоку" взять и адаптировать код pqf из другого кодека для моего случая (M=16, 384 fir коэффициента) не очень сработала - хотя для DC или синуса результат очень похож на правду, и после обратного синтеза получилось что то похожее на near perfect reconstruction, но для дельта функции или просто прямоугольного фронта после синтеза получается уже совсем не похожий на входной сигнал.

Что хочется для себя прояснить.

  1. Как по прототипу понять является ли фильтрбанк near perfect reconstruction или perfect reconstruction?

  2. Из разных статей понял, что есть случаи когда анализирующий и синтезирующий фильтрбанк используют разные прототипы. Есть подозрения что это именно мой случай. Можно ли это как то узнать из анализа прототипа? И как получить из одного другой?

Ну и в целом тема интересная, есть много вопросов, но к сожалению математику происходящего понимаю плохо (

Тема действительно интересная, но я с ней ещё не сталкивался, к сожалению. Решал похожую задачу, но во временной области — для нелинейной обработки сигнала через FFT. Тут, как я понимаю, аналогичная ситуация, только восстанавливать в единицу нужен. Ещё, поскольку кодек со сжатием, часть высокочастотного сигнала безвозвратно теряется и полного восстановления сигнала получить всё равно не получится. Поэтому для тестирования нужно брать сигнал с идеально ровным спектром и после анализа-синтеза смотреть его на предмет появления пульсаций в рабочей полосе спектра.

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


Выглядит так, как будто авторы специально усложнили алгоритм для реверс-инжиниринга.

Анализирующий pqf это первое преобразование при кодировании, и соответственно синтезирующий pqf - последнее при декодировании. Все потери происходят сильно позже, уже после перехода в частотную область. Таким образом потерь я пока просто не допускаю, просто подаю результат работы своего анализирующего фильтрбанка на вход синтезирующего из ffmpeg. И ожидаю как минимум "near perfect reconstruction".

Коэффициенты часто записываются не последовательно для более простого кода их применения, так почти во всех кодеках, в том числе и opensource. Там же хитрое прореживание совмещенное со сверткой с частью прототипа. В dcaenc боле понятно https://github.com/FFmpeg/FFmpeg/blob/master/libavcodec/dcaenc.c#L319-L367 И вот форма этого прототипа сильно влияет на свойства фильтрбанка. В том же dca есть 2 фильтрбанка, для PR и NPR.

NPR обеспечивает подавление алиасинга только в соседних полосах, часто этого достаточно, так как за пределами соседних полос ослабление уже обычно более 90db, и этим можно пренебречь для задач аудио. Например вот тут https://www.researchgate.net/publication/221380083_A_method_to_convert_near-perfect_into_perfect_reconstruction_FIR_prototype_filters_for_modulated_filter_banks
описано получение из near perfect reconstruction, perfect reconstruction. Забавно, что perfect reconstruction уже не выглядит как sinc(x)/x а имеет некоторые "резкие изломы".

Почитаю еще раз, вашу предыдущую статью может какие идеи возникнут...

А предыстория у вашего кода есть? Нельзя же просто так с нуля взять и аудио-кодек написать, ещё и совместимый с проприетарным.

Предыстория? Я как то порывался статью написать, но все не мог себя заставить)) Аудио и видео кодеки всегда были мне интересны, еще до того как программировать стал. Потом как то сидел в отпуске и захотелось попробовать себя в чем то новом (по работе я не связан с DSP никак), стал читать код ffmpeg, увидел что там есть декодер для ATRAC, а когда то я был фанатом минидиска, ну тут все и совпало - решено написать енкодер по коду декодера из ffmpeg. ATRAC1 мне показался достаточно простым (стек из 2х QMF + MDCT + квантование в частотной области, ну еще опционально переключение размера окна), ну и стал постепенно изучать область и реализовывать. Особых сложностей с ATRAC1 в общем то не было, и достаточно быстро получился рабочий код, дающий в целом не плохой результат. Ну так и втянулся.

Потом взялся за ATRAC3 - он не сильно сложнее, другой фильтрбанк, но все еще имеющий в основе стек из QMF + MDCT, есть кодирование хаффмана. Из интересного, нет переключения размера окна, но есть gain modulation - возможность изменить уровень сигнала в пределах окна перед mdct в кодировщике и соответственно обратным образом изменить его после imdct в декодере, эта фича, кстати, до конца мной не была реализована. Тем не менее для LP2 режима (это что то в районе 132Kbit/s) качество получилось сравнимым с mp3 - 128.

Тем временем код на github заметили другие любители минидиска, и с удивлением для себя я обнаружил, что оказывается проект используют в том числе для заливки треков на минидиск теми кто по каким то причинам до сих пор ими пользуются)) А кто то собрал его под WebAssembly и предлагал как web сервис))

Тем не менее мне интереснее разбираться дальше, и ATRAC3Plus тут как раз интересен тем что использует Generalized Harmonic Analysis (GHA) для извлечения информации о тональных составляющих - в других кодеках такого подхода я не встречал, и было бы интересно поиграться с этим. Но ATRAC3Plus использует pqf для разбивки входного сигнала на 16 суб-полос перед mdct, и подозреваю что разные прототипы для синтеза и анализа, и тут я застреваю, поскольку дальше без понимания математики происходящего в pqf никак.

Вот такая история)

но для дельта функции или просто прямоугольного фронта после синтеза получается уже совсем не похожий на входной сигнал.

А должен? Вы же понимаете, что такой сигнал невозможен, если мы говорим о результате работы АЦП. Т.е. он в принципе не в домене звукового кодека лежит.

Эксперименты с фильтрбанками для которых я могу найти заведомо корректную реализацию говорят что да. Для фильтрбанка DTS кодека (реализация есть в ffmpeg), например, ошибка где то в 5м знаке после запятой была.

В статье повстречались мелкие ошибки по синтаксису и орфографии:

1) -"Переход между этими областями осуществляется через преобразование Фурье, которая между функциями во временной и в частотной области задаёт однозначное соответствие." Здесь определение "которая" женского рода, а должно быть среднего, так как относится к слову "преобразование", которое среднего рода. Тут бы не помешало буковку "а" заменить на "о", в окончании слова "которая".

2) -"По сути, представляют собой всё те аналоговые фильтры, но «моделируемые» в цифре..." Здесь недостаёт частички "же" после "всё те", для более удобоваримого прочтения.

3) -"Среди множества разработанных окон одним из наиболее оптимальных и удобным для аналитических вычислений является окно Нуттала, определяемая как..." Тут "определяемая" - по смыслу соотносится к существительному среднего рода "окно", а раз так, то и окончание должно быть так же среднего рода: "определяемое", а не женского - "определяемая". Ну или как вариант, можно добавить слово "функция".

4) -"Используя полином 5-ой степени, с двумя нулевыми производными на границах сопряжения." Здесь недосказанный непорядок в первом слове. "Используя", подразумевает некое дальнейшее действие. Но его нет до самой точки в окончании предложения. А раз нет дальнейшего действия, то целесообразнее изменить окончание у первого слова, например на вариант "используем", тогда получится предложение с констатацией факта. А иначе смысл как бы повисает в воздухе, а логичного завершения в данном предложении не следует.

5) "-Используя рациональный полином — даёт более гладкую характеристику и более простую обратную функцию, чем у кубической:" Здесь первому слову "используя" (что делая здесь и сейчас), глагол "даёт" не соответствует. Просится изменение первого слова на "использование", с коррекцией окончаний в последующих двух словах: -Использование рационального полинома — даёт более гладкую характеристику...

6) -"Для реализации таких фильтров в «мультиампинг»-системах можно использовать готовые решения, позволяющих загружать заранее посчитанные импульсные файлы." Здесь следует изменить окончание в слове "позволяющих", на "позволяющие", так как оно относится к существительному женского рода во множественном числе: "решения". Готовые решения позволяющие загружать... (так получится правильно). А пока читается несколько коряво: -"Готовые решения позволяющих загружать..."

Большое спасибо за замечания! Все эти глупые ошибки — следствие злоупотребления автодополнением в современных IDE при программировании.

Выберем размер массива...

Массива чего?

Каждому элементу в нём будет соответствовать частота от 0 (постоянная составляющая) до samplerate/2 (частота Найквиста), равномерно распределённых в линейном масштабе.

Переведите это на русский, пожалуйста...

Здесь каждую нечётную частоту мы повернули на 180° ...

Частоты можно вращать? (Да вращал я ваши частоты знаете ли...)

...для того чтобы после БПФ максимум импульса был расположен по центру.

Центру чего, простите?

Массива чего?
Ответ дан в следующей Вашей же цитате:
Каждому элементу в нём будет соответствовать частота
Алгоритмы FFT обычно работают inplace, то есть результаты преобразования сохраняются в массив с исходными данными. Поэтому вопрос о том, находятся там отсчёты во временной или частотной области — это вопрос интерпретации.

Частоты можно вращать?
Мне это казалось допустимой формулировкой, потому что речь идёт о комплексных числах. Исправил на фазу.

Центру чего, простите?
Центру массива.

Сижу на паре по цифровой обработке сигналов и благодарю вас за эту статью, потому что теперь наконец лучше понял, что же это за звери такие, КИХ-фильтры)

Sign up to leave a comment.