Resampling 44100 -> 12800, нужна помощь по реализации. Опции

[evil_laugh]

Приветствую уважаемых гуру цифровой обработки

Прошу помочь мне с алгоритмом, который бы осуществлял изменение частоты дискретизации сигнала с 44100 до 12800 Гц. Сигнал - 16битные отсчёты с микрофона, обычный PCM.

Требования к передискретизатору следующие:

а) пропускание частоты сигнала в полосе от 0 до 6 кГц с неравномерностью АЧХ не более 3 дБ
б) подавление сигнала в полосе 6200 Гц - 6400 Гц не меньше, чем на 60 дБ
в) и самое главное - крайне жёсткие требования к объёму используемой оперативной памяти и процессорного времени. Процессор, для которого пишется алгоритм, помимо передискретизации загружен ещё тонной задач. Оперативной памяти для алгоритма может быть выделено не более 2 кБ. Время - как всегда, чем меньше, тем лучше.

Расскажу о том, какие исследования я провёл в этой области.

1. Естественно, что чтобы сделать изменение Fд 44100->12800, нужно сперва повысить частоту в 128 раз, потом понизить в 441 раз. Фильтр, который нужен для такой передискретизации в один каскад, имеет больше 7000 коэффициентов (рассматриваем КИХ-фильтрацию), поэтому вариант с таким алгоритмом отпал сразу.
   
2. Следующим моим шагом было создание многокаскадной системы: сперва повышение частоты в 2*4*16=128 раз, потом её уменьшение в 21*7*3=441 раз. Естественно, использовалось полифазное разбиение фильтров. На входе алгоритма был 441 отсчёт, на выходе - 128. Всё работало и обеспечивало требуемые параметры АЧХ, однако требовалось много места (нужно же где-то хранить 441*2*8 промежуточных отсчётов сигнала на высокой частоте) и много процессорного времени. Поэтому от этого алгоритма тоже пришлось отказаться. Было поставлено требование: только однокаскадная передискретизация.
   
3. Попутно, кстати, была прочитана горка литературы по фильтрам Фарроу, Ланцоша и прочим полиномиальным фильтрам, но ввиду почти полного отсутствия вменяемых реализаций с графиками и прочими данными - я понял только общий принцип и то, что обмануть никого не удастся и всё равно придётся реализовывать огроменный фильтр.
   
4. Дальше попробовал сделать так: разработать фильтр с меньшим количеством отсчётов, а при фильтрации проводить интерполяцию полученных значений. При этом было решено поступиться требованиями к переходной полосе фильтра и сделать так, чтобы на 6200 Гц обеспечивалось подавление не ниже 15 дБ, но все дальнейшие частоты выше какой-то величины (у меня получилось 7 кГц) всё равно должны были подавляться не меньше, чем на 60 дБ. Алгоритм был написан на Си и показал свою работоспособность. Казалось бы, всё в порядке и можно писать программу для DSP, но меня не покидало ощущение того, что можно попробовать что-то ещё. Поэтому плавно переезжаем к...
   
5. ...фильтрации в частотной области, которую я решил рассмотреть как вариант реализации передискретизатора. Что нужно? Правильно, взять ДПФ исходных 441 точек сигнала, занулить ненужные гармоники, выкинуть ненужные 441-128 отсчётов из середины спектра и золотой ключик у меня в кармане. Однако снова возник ряд проблем. Первая - это "как быстро взять ДПФ 441 точек?". Дополнение нулями до 512 точек, судя по всему, не прокатывает - для передискретизации в частотной области придётся выкинуть нецелое количество точек. Ладно, для этого вопроса решение было найдено путём штудирования книжки Блейхута "Быстрые алгоритмы цифровой обработки", откуда я выудил алгоритмы Кули-Таки и алгоритм Гуда-Томаса (который прекрасно работает для числа 441, ибо это 3*3*7*7 - произведение простых чисел). Решение спорное, впрочем, ибо архитектура DSP, ИМХО, ну ни разу не приспособлена для подобного рода выкрутасов со всеми этими адресными тасованиями, тем более в моём дилетантском исполнении. Впрочем, возникает вторая проблема, которая затмевает собой первую: как избежать эффекта Гиббса в данной ситуации? Напрашивается всем известный overlap-save, однако не получится ли так, что время, затрачиваемое на вычисление "лишних" точек, которые потом выкидываются, сведёт весь выигрыш от используемого алгоритма на нет?

Собственно, у меня вопрос к уважаемому сообществу: что можно ещё придумать, чтобы выполнить задачу? Приветствуются пинки ногами, крики "пшёл вон, нуб", а также ссылки на авторитетные источники, личный опыт и другую информацию. Напомню задачу: 44100 -> 12800. пропускание 0-6 кГц с неравномерностью 3 дБ, подавление 6200-6400 минимум 60 дБ, минимум процессорного времени, 2 кБ памяти.

Спасибо


Сообщение отредактировал evil_laugh - Dec 18 2012, 16:03

[evil_laugh]

Serg76, похоже, Вы познали Дао. У меня только один вопрос: КАК? Можно, пожалуйста, подробнее про алгоритм?

[Serg76]

Serg76, похоже, Вы познали Дао. У меня только один вопрос: КАК? Можно, пожалуйста, подробнее про алгоритм?

да как и описывал выше на основе полифазника, здесь ничего сверхестественного, правильно подберите нужное окно. я взял распространенное Наттала, правда пришлось все-таки увеличить порядок фильтра до 300, вроде как в -70 дБ влезли (пардон, -60 дБ) А вообще, -60 дБ - это еще по-Божески, где-то попадалось на просторах Инета ресамплер как раз звукового диапазона с подавлением чуть ли не -120 дБ.

Сообщение отредактировал Serg76 - Dec 20 2012, 19:02

[evil_laugh]

да как и описывал выше на основе полифазника, здесь ничего сверхестественного, правильно подберите нужное окно. я взял распространенное Наттала, правда пришлось все-таки увеличить порядок фильтра до 300, вроде как в -70 дБ влезли (пардон, -60 дБ)

То есть Вы делаете простую однокаскадную передискретизацию, но с таким вот хитрым фильтром? А на какую частоту дискретизации фильтр считался? Я так понимаю, Вы не делали "честное" повышение частоты дискретизации до 44100*128, а сделали как-то по-другому? А как?

[Serg76]

То есть Вы делаете простую однокаскадную передискретизацию, но с таким вот хитрым фильтром? А на какую частоту дискретизации фильтр считался? Я так понимаю, Вы не делали "честное" повышение частоты дискретизации до 44100*128, а сделали как-то по-другому? А как?

Да, однокаскадная децимация. В полифазной фильтрации нет ничего "хитрого", она известна давно, есть даже статья в википедии http://ru.wikipedia.org/wiki/%D0%9F%D0%B5%...%86%D0%B8%D1%8F
довольно понятно все описано. никакого повышения/понижения делать не надо. преимущество в том, что все вычисления производятся на выходной частоте дискретизации.

[evil_laugh]

Да, однокаскадная децимация. В полифазной фильтрации нет ничего "хитрого", она известна давно, есть даже статья в википедии http://ru.wikipedia.org/wiki/%D0%9F%D0%B5%...%86%D0%B8%D1%8F
довольно понятно все описано. никакого повышения/понижения делать не надо. преимущество в том, что все вычисления производятся на выходной частоте дискретизации.

Это мне понятно. Но фильтр же рассчитывается на промежуточную частоту дискретизации? Полифазная фильтрация - это метод разбиения одного большого фильтра на несколько маленьких с низкой частотой дискретизации. Но исходный-то фильтр (который разбиваем) работает на большой частоте?

[Serg76]

Это мне понятно. Но фильтр же рассчитывается на промежуточную частоту дискретизации? Полифазная фильтрация - это метод разбиения одного большого фильтра на несколько маленьких с низкой частотой дискретизации. Но исходный-то фильтр (который разбиваем) работает на большой частоте?

нет, все фильтры работают на низкой частоте. изначально характеристика фильтра рассчитывается на нормированную частоту 12800/44100, а затем для каждого элементарного фильтра пересчитываем эту характеристику с учетом смещения. в процессе фильтрации для каждого выходного отсчета вычисляем отклики двух соседних фильтров и проводим линейную интерполяцию между ними. поэтому выходной спектр и "обогащается" гармониками, которых не было во входном сигнале.

[alex_os]

нет, все фильтры работают на низкой частоте. изначально характеристика фильтра рассчитывается на нормированную частоту 12800/44100, а затем для каждого элементарного фильтра пересчитываем эту характеристику с учетом смещения. в процессе фильтрации для каждого выходного отсчета вычисляем отклики двух соседних фильтров и проводим линейную интерполяцию между ними. поэтому выходной спектр и "обогащается" гармониками, которых не было во входном сигнале.

Т.е метод эквивалентный предложенному ув. termit ? А сколько фаз и какая длина каждой фазы ?