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

Если делать не в режиме Фарроу, а просто как программную реализацию полифазного варианта, то по какому правилу следует выбирать момент и фазу выхода полифазного фильтра чтоб, например, передискретизировать с частоты 4*1/3 в частоту 2 при условии, что полифазный фильтр состоит из 16 ветвей и в каждой ветви 16 коэффициентов.

ничего не понял, особенно что за режим Фарроу (Farrow это структура построения интерполятора, а не какой то алгоритм) %)

Режим Фарроу -эта структура, которая просто была реализована чисто программно.

ну вот и какое отношение может имеет эта структура, к алгоритму расчета дробной ошибки на основе биения частот ?

Да сам Фарроу не имеет, и полифазник не имеет.

X1-X100 отсчеты сигнала
h1-h4 отсчеты ИХ первой фазы
h5-h8 отсчеты ИХ второй фазы
h9-h12 отсчеты ИХ третьей фазы
h13-h16 отсчеты ИХ четвертой фазы

На первую фазу идут Х1 Х5 ...
На вторую фазу идут Х2 Х6 ...
На третью фазу идут Х3 Х7 ...
На четвертую фазу идут Х4 Х8 ...

Если мне надо выход иметь на входной частоте, но сдвинутый по фазе, то тупо выбираю подходящую фазу.
Но меня клинит как организовать съем данных в случае, когда мне надо перейти от входной частоты 4*1/3 к выходной 2.

ну дык так бы и спросили сразу, а не описывая структуры и подчеркивая что у вас другой интерполятор %)

Полифазник как был так и остался, его задача получить сдвиги базового сигнала. А вот из частоты 4*1/3 частоту 2 получить элементарно. собираете NCO на тактовой частоте 4*1/3 и на нем генерируете частоту 2.
получаете моменты отчета в дискретах частоты 4*1/3 и дробные ошибки.

Интервал дискретизации входной - 1/(4*1/3)=0.75
Интервал дискретизации выходной - 1/2=0.5

Педположим необходим 5й отсчет выходного интервала, т.е. в новой шкале времени это момент времени 5*0.5=2.5

Вопрос - при каком номере входного отсчета и с какой фазы брать отсчет?

Когда я работаю интерполятором в виде Фарроу, то нахожу базовый индекс как целая часть отношения 2.5/0.75, определяю дробную часть, рассчитываю коэффициенты и по 4м отсчетам относительно базовой точки рассчитываю выход.
А при полифазнике я должен снимать его выход в момент времени, когда в него зашел соответствующий базовый отсчет?
Когда схема нарисована в железе, там понятно вроде бы. Но я делаю чисто программно, поэтому для меня как программиста фраза типа "собираете NCO на тактовой частоте 4*1/3..." требует осмысление в программной реализации, особенно с учетом того, модель делается в Матлабе, где достаточно применяется подход обработки блоком, а не поотсчетно.

Как программист почитайте
http://www.analog.com/static/imported-file...s/EE183Rev5.pdf

Они там строят передискретизацию с 32-разрядной точностью для акустики ( используя программы MDS для синтеза НЧ-фильтров)

Полифазная реализация — это всего лишь оптимизация прямого алгоритма, которой делает апсэмплинг до частоты 4 (в вашем случае), фильтрует и прореживает. Забудьте про полифазные компоненты, рассмотрите единый фильтр, работающий на частоте 4, и уже через это рассмотрение выведите формулу для "фаз" в полифазной реализации. Она будет примерно такой: phase := (phase + D) % U, где U и D — коэффициенты интерполяции и прореживания.

простой рисунок входной/выходной тактовой, сделанный на коленке вам в этом не поможет ? Все же очевидно и не зависит от реализации интерполятора

UPD. в данном случае вам необходимо интерполировать входной сигнал на интервале между третьим и четвертым отсчетом, посчитать дробную ошибку (в вашем случае, судя по рисунку ее абсолютное значение 0.25, знак зависит от реализации интерполятора), и т.к. интерполяция у вас будет с определенным дискретом (как я понял вы делаете на банке фильтров) выбрать наиболее близкую фазу (в лагранжевых это решается расчетом полинома). как бы все %)

UPD2. Вы интерполируете в 1,5 раза, ИМХО потребуется достаточно хороший интерполятор, что бы не вылезло ничего левого в выходном спектре выше ~70дБ

Пользуясь случаем, хочу спросить - почему в любезно выложенном fontp файле для ресемплинга используется алгоритм "повышение частоты в L раз - фильтрация - фильтрация - понижение частоты в M раз" вместо имхо простого варианта рассчитать на лету и применить полифазный фильтр с нужной фазой? В документе написано

Или "на лету" это займет больше машинного времени? А там DSP, поэтому это важно? А если на ПЛИС, где можно параллелить процессы, проще будет рассчитать коэффициенты честного полифазного интерполятора без апсемплинга/децимации?

Непонятно, что такое "на лету", "честный полифазный интерполятор без апсэмплинга" и т.д...

Документ вполне справедлив.

Буквально сейчас играюсь - считаю в своем симуляторе полифазные интерполяционные КИХ фильтры из 160 коэффициентов. Рассчитал на середину интервала = 1/2, получился апсемплинг в 2 раза. По тому же правилу рассчитываю коэффициенты для любой другой точки интервала - 1/3, 0.99 и т.п. Может будет понятнее если я употреблю термин "смещенный синк", хотя у меня не совсем синк. Смотрю - величина ошибки интерполяции та же. Вот я и думаю, почему не пошли по пути - получили точку внутри интервала (допустим 0.3766535484757473093947447656), рассчитали для нее нужное количество коэффициентов и сделали свертку (так это называется, кажется)? Как быстро рассчитывать коэффициенты для каждой точки - отдельный вопрос, но у меня есть некоторые имхи. Но в принципе то это должно хорошо работать. А на ПЛИС так ещё и быстро (достаточно для реалтайма).

Для каждой точки считать свои коэффициенты долго, да и не нужно. Именно поэтому и созданы полифазные фильтры: в них всё посчитано заранее.

Видимо я не понимаю что такое полифазные фильтры - я считал что это именно то что я описал, "синк со смещенной фазой". Придется разобраться. Но неужели апсемплировать сигнал в несколько десятков-сотен раз, потом фильтровать его на этой частоте и выкидывать из результата 99% отсчетов будет проще по объему вычислений? В отстаиваемом мной варианте весь вопрос упирается в оптимальную процедуру получения коэффициентов на каждом шаге.

Это одно и то же.

изменение частоты потока данных через интерполяцию/децимацию производиться без ошибок интерполяции, в отличии от использования различных интерполяторов.

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


Причем, на нули, обычно, никто и не умножает (хоть могут и быть варианты). Повышение частоты часто держим в уме для понимания. Поскольку офсеты повторяются периодически и их можно рассчитать заранее. Это же передискретизация "рациональная", а не какая попало. И хотя доказано, что дробью n/m можно приблизить любое число с точностью порядкка 1/(m*m), всё таки n и m используются только в разумных пределах

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

Где-то глубоко внутри я подозреваю что это так и есть, просто наверное сейчас я ещё не дошел до ясного и целостного понимания картины.

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

ЗЫ если я не ошибаюсь, частный случай отстаиваемого мной подхода называется "банк FIR фильтров" - если отношение частот у нас кратное и через N раз отсчеты выходного сигнала будут совпадать с входным, то нам надо хранить только N наборов коэффициентов, а точнее даже N/2 т.к. они симметричные. А если частоты некратные - тогда придется рассчитывать.

fontp, спасибо за некоторые разъяснения. Глупый вопрос:
мы разбавили сигнал нулями раз в 25, потом применяем НЧ фильтр на этой новой частоте. Дошли до определенного отсчета - "слева" у нас все бывшие нули стали нормальными отсчетами (ибо прошли фильтрацию) - и на них уже нужно домножать коэффициенты НЧ фильтра? Мы используем уже рассчитанные новые значения при фильтрации последующих? Если да, то здесь наверное может быть выигрыш по сравнению с "прямым чем-то похожим на синк".


Ну тут ещё вопрос, является ли аудио случаем самых жестких требований к качеству ресемплинга, но речь не об этом. Именно по причине реализации очень хорошего ресемплера (который ресемплит 16 бит даже точнее чем они заданы) интересно вникнуть в суть оптимизации метода. Хотя например при кратной сетке частот 24/25 я не понимаю почему сложнее хранить банк из 24 FIR фильтров. Хотя очевидно, что если исходный сигнал "грязен". т.е. имеет частоты от 20 до 22 кгц, к примеру, то вариант с НЧ фильтром его еще и вычистит параллельно с ресемплингом, а синк оставит эти составляющие.

Нет, входные отсчеты. Используются просто более короткие фильтры чем синк)), но с лучшей прямоугольностью, чем синк усеченный

в вашем оно может быть и так, а вот в спектральном все четко. не придерешься.

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

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

Понимаете, обсуждать то, что у Вас что-то похожее на синк, лучше чем что-то непохожее это бессмысленно. Это означает только что Вы умеете синтезировать более лучшие фильтры НЧ. Это уже схоластика

fontp, давайте попробуем так: предложите условия задачи: исходная частота, новая частота, исходная частота Найквиста, необходимая точность при ресемплинге. Желательно из аудио-примеров. Для начала с рациональным отношением частот и не слишком большими n и m. Я рассчитаю банк FIR, определю количество их коэффициентов. Вы посмеетесь и скажете что полифазник гораздо меньше по вычислениям

А чем Ваш "банк FIR" не будет тем же полифазником?

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

Можно каскадно, а можно и "в лоб". Ничего, что фильтр длинный. Стандартные частоты в аудио требуют ↑ в пределах пары сотен раз (напр., 44100 → 96000). Умножаем это на типичную длину полифазной компоненты (20-30) и получаем всего десяток тысяч отсчётов. Умножаем синк на окно — фильтр готов.

Мне кажется имелось в виду что для "очень хорошего фильтра" недостаточно "типичной длины 20-30" а нужно побольше. Насколько - зависит от степени "хорошести" фильтра. И все это потом ещё умножить на количество различных фаз. Имхо вопрос из серии: насколько экономнее по вычислениям прогнать сигнал через 2 последовательных относительно коротких фильтра чем через один длинный, обеспечивающий такую же прямоугольность.

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

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

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

Детали тоже интересны, особенно если остальное более-менее понятно. Например, объединение в один полифазный фильтр двух с разным количеством фаз (в документе кратко описано, надо разобраться) и тот самый алгоритм трехкаскадного сабсемплинга от MDS (явно с документе не описан, только сказано что при количестве каскадов 4 и более ну нереально сложно сделать оптимальный алгоритм).


Это если исходный сигнал дискретизирован "впритык", то есть там есть спектр почти до Найквиста. А если дискретизирован с большим запасом, то имхо прямоугольность можно не обеспечивать, только большое подавление в полосе задержания, а переходную полосу можно растянуть насколько запас позволит. Здесь я чуть ли не в первый раз согласен с des00

Там на каждой ступени фильтр оптимален (по Ремезу). А в целом возможен только перебор - перераспределяя количество тапов по каскадам ( с учетом заданной разрядности поскольку общая х-ка как раз страдает).
Если Вам такие детали интересны, я мог бы забросить программу для ознакомления на опыты))
У них, кстати, есть тоже неплохой AN но там более широко чем передискретизация, там фильтры multirate
http://mds.com/system/resources/BAhbBlsHOg...ate_article.pdf

Мне интересны такие детали. Только например если Вы забросите код на чем-нибудь типа VisualDSP++ или как её там или чем-то подобном, на хитром специальном скрипте, то вряд ли я смогу его посмотреть.
У меня много вопросов возникает по сабжу вообще и этому документу в частности (с учетом того что я только несколько часов назад узнал про метод расчета фильтров через синк с окном). Например - если это DSP, то есть готовая прошитая микросхема, которая больше ничего не делает (или там есть возможность дописать ещё исполняемый код?) то зачем оптимизация в каждом случае? Имхо достаточно оптимизировать самый затратный случай а остальное делать чтобы было не дольше.

Это очень хорошо! Может в документе было упоминание об этом, сейчас ещё просмотрю. Но этот AN весьма интересен.
UPD спасибо, скачал ПДФ с сайта MDS. Что мне мешало раньше зайти на их сайт (по ссылке из документа)?....
Есть несколько вопросов по этой теме, но думаю до изучения AN они преждевременны. Например, насчет полиномиальной или сплайновой интерполяции после двух каскадов сабсемплигна - не будет ли быстрее при той же точности? Хотя надо для начала узнать все ли отсчеты сабсемплируются.

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

Речь шла о коммерческих программах синтеза фильтров для Виндуз, упомянутых в статье. Это экзешники а не скрипты))
Предполагается, что фильтры конструируются однажды заранее, а не налету.
Программы довольно старые и когда-то были доступны в Сети, они должны быть в закромах.
Оптимизация там делается в обоих возможных смыслах.
Если это DSP то нужно минимизировать число операций. Процессор один и он должен последовательно выполнить все фильтры разных каскадов, естественно с разной частотой.
Если это FPGA то делается попытка минимизировать число арифметических устройств в конвейере каскадов(т.е. в каком-то смысле вентилей при некоторых предположениях о сложности умножителей)

Простите, нет опыта правильной контекстной интерпретации некоторых слов. Если Вам не трудно, поделитесь пожалуйста экзешниками. Доступа в "Свои" нет, к FTP соответственно тоже. Есть мыло Ivana_2004 на мэйлру

Скачал демки с сайта MDS, вникаю. И продолжаю задавать глупые вопросы: допустим мы делаем полифазный ресемплинг из 44100 в 48000. Фильтр каждой фазы считаем с оптимальнейшим соотношением прямоугольность/подавление чтобы и полезный сигнал не испортить и отражения подавить. Но в итоге у нас в выходном сигнале полоса "лишнего" будет максимум 20000-24000. Так почему бы не сделать так: считать относительно короткие полифазные фильтры с достаточным подавлением и плохой прямоугольностью (широкой переходной полосой) а потом "вытянуть" получившийся плавный спад результирующей АЧХ ВЧ фильтром выходного сигнала?

Так и делают во всех аппаратных ddc. Одна-две ступени скользящего среднего и на выходе fir.Последний с хорошей прямоугольностью, учитывая искажения cic-фильтров первых ступеней. И оптимизация перебором по параметрам - вариантам распределения коэффициентов сабсэмплинга и возможной длины фильтров по каскадам. Analog дает программы такого расчета для своих ddc. Только до -100 db они не тянут, а для других приложений , радиочастотных например,получаются вполне адекватные фильтры

Спасибо, буду пытаться осознать.
Продолжаю, в попытках прояснить картину. Ответ на этот вопрос я наверное могу получить самостоятельно, но нужен соответствующий софт (Матлаб и т.п.) и умение с ним работать, к тому же мне хочется прояснить теоретические аспекты: есть сигнал с равномерным (или почти) спектром, представляющий из себя шумоподобный или шум, дискретизирован с 44100. Значимый спектр тянется вплоть до 22050. Производится его ресемплинг в 48000 с помощью некоего ресемплера - пусть полифазного фильтра. Полоса пропускания фильтра берется 0-20000. Вопрос: какой (спасенный) спектр будет у ресемплированного сигнала (значимый, не считая тихих шумов)? Варианты ответа:
1) 0-24000 - врядли, частоты 22050-24000 не должны никак присутствовать при адекватном ресемплинге
2) 0-22050 - очень вероятно, если НЕ умножаем на бывшие нули при полифазной фильтрации
3) 0-20000 - очень вероятно, если умножаем на бывшие нули при полифазной фильтрации
4) свой вариант ответа
В данный момент я предполагаю, что в зависимости от реализации ресемплера будет или 2) или 3). Причем, если мы сделаем ресемплер по принципу минимального количества операций (не умножая на бывшие нули), то будет 2).

UPD простите за суету, постепенно приходит понимание Просто к теме фильтрации пришел из аппроксимации, поэтому некоторые сложности.
С "умножаем на бывшие нули" - это я сглупил, перепутал с БИХ-фильтром (кстати, почему бы их не применить?)
А с полосой - все верно, если мы рассчитываем фильтр 0-20000 то у нас эта полоса и останется при любой децимации до частоты дискретизации не ниже 40000. Просто можно мыслить в категориях что мы считаем НЧ фильтр 0-20000 на высокой частоте а потом разбивать его на разные фазы соответствующие ненулевым отсчетам входного сигнала на этой высокой частоте. А если оставаться на исходной частоте (44100), то можно:
1) рассчитывать интерполяционный фильтр, т.е. фильтр 0-22050 - но для любой фазы (что я и делал) - и тогда это аналогично расчету 0-22050 на высокой частоте и по идее должно приводить к тем же самым коэффициентам полифазного фильтра
2) рассчитывать ШП НЧ фильтр 0-20000 - но тут пока непонятно откуда брать разные его фазы для последующего получения 48000 (если это вообще возможно при такой постановке задачи).

Читаю Multirate Filter Design и думаю: какая красивая идея - делать узкополосный НЧ фильтр не в лоб, а сначала отрезать побольше с широкой переходной полосой, заменить сколько можно отсчетов нулями и финально фильтровать на той же частоте экономя операции на этих нулях!...

Если у фильтра передискретизации полоса пропускания 0–20000, а переходная полоса 20000–22050, то сигнал от 0 до 20000 будет передан практически без искажений, сигнал от 20000 до 22050 можно при желании восстановить корректирующим фильтром, а частот выше 22050 в исходном сигнале и не было — непонятно, что там "спасать".