Коллеги,есть несколько вопросов по фильтрам описанным в книге Лайонса в главе 13.24.1 Нерекурсивные АГФ.
Там предлагается строить фильтр по следующей структуре(см.рисунок)
Но если построить АЧХ фильтра первого каскада для М = 5 то получим следующую картинку. Имхо,фильтр с такой АЧХ вообще нельзя использовать в качестве децимирующего
Подскажите кто в курсе в чём здесь секрет.

Так у всех ИГФ такая АЧХ. Там же Лайонс предлагает дополнять корректирующим фильтром, в конце раздела. Я в Матлабе тоже потыкался, ничего хорошего не получил.
Если вы децимируете, значит, вам нужна ограниченная полоса. Там, где спад еще небольшой.

Хорошо,давайте посмотрим что получается. На рисунке по оси частот число 0,5 это есть частота Найквиста после децимации. Т.е. если на входе фильтра будет присутствовать сигнал на частоте 0,6 то на выходе фильтра получим образ на частоте 0,4 (или 0,8 Найквиста после децимации) причём с подавлением всего -20дБ! Это очень плохо.
А при рекурсивной реализации ведь никто не заставляет использовать первый порядок,а например 5-й, тогда образ будет задавлен на 60дБ,но коэффициент усиления будет офигенный. А у нерекурсивного на мой взгляд можно нормировать сигнал между каскадами,тогда не будет большого коэффициента усиления.
Вот я и ищу разумный компромисс, пока не очень получается.

Матлаб показывает диапазон только до половины частоты дискретизации.

Я не вижу разницы между рекурсивным и нерекурсивным фильтром. Я и не помню нерекурсивного в книге. "Цифровая обработка сигналов"?

Значит Вы ничего не поняли про -20дБ.

Ясно.

По-моему, мы уже говорили. Сначала задавите полосу, потом децимируйте.

Наш разговор немного напоминает разговор слепого с глухим. Если не можете ответить на поставленный вопрос то я немного озвучу задачу. Мне нужно продецимировать сигнал полосой 50МГц в 16 раз из 1200 MSPS в 75 MSPS. АЧХ CIC фильтра пятого порядка на картинке. Видно офигенное подавление боковых лепестков что хорошо, но так же видно и коэффициент передачи 150дБ! Неслабенький такой. Второй минус - неравномерность АЧХ в полосе около 50дБ. Только не говорите мне что эту неравномерность можно легко скомпенсировать дополнительным фильтром. Думаю сначала децимировать на 8 (тогда неравномерность будет всего 10 дБ,решаемо), а потом поставить ещё один дециматор на 2 на half-band фильтре. Но и в этом случае коэффициент усиления остаётся очень большим. Вот и ищу компромиссы...

Я много моделировал CIC в симулинке, правда децимация степень двойки, так же как и другие параметры и это усиление компенсируется банальным отбрасыванием младших разрядов. Главная проблема действительно в создании компенсирующего фильтра с обратной sinc характеристикой, поэтому на практике (в литературе встречал много раз) этим фильтром децимируют до новой частоты дискретизации большей не менее 8 раз полосы сигнала. Остальное режут обычным ФНЧ. Задача CIC не полностью отфильтровать сигнал, а свернуть исходную полосу таким образом, чтобы полоса сигнала попала в нули функции, не исказив исходный сигнал. Ну а вообще этот процесс лучше всего описан не у Лайонса, а на dsplib.

1200 это 4 фазы по 300МГц?

Всё правильно, но вот насчёт не уверен. Если звено состоит из интегратор-ГФ-дециматор то на выходе отбросить разряды не проблема. А если интегратор-дециматор-ГФ то не уверен. Помнится,где-то читал что так делать нельзя.

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

2. В вашем случае нужен cic с прореживанием на 8 и ких с прореживанием на 2.

3. Картинки ваши - лажа. Ибо соответствуют фильтру аж 32-го (???) порядка.

+1. Где то ошибка в порядке. Если у вас 4 фазы по 300МГц, то :
1. Делаете фильтр скользящего среднего '{1,1,1,1} - простые сумматоры.
2. Ставите 6 звеньев подряд. В итоге имеете полосу по -3дб ~57.587МГц, затухание на границе полосы задерживания зеркала ~-101 дб.
3. Децимируете на 4. Переходите на 300МГц.
4. Фильтруете обычным FIR/bireciprocal IIR/half-band FIR
5. Децимируете на 4.

Если берете FIR, то операцию 4, можно делать с учетом операции 5. Будет небольшая экономия.

ЗЫ. Количество звеньев можно уменьшить

UPD. Ошибся чуток -101 будет для децимации на 2. Если считать на полосу 50МГц и децимацию на 4, то будет так :
Затухание на полосе 50МГц = -2.287дб, наименьшее затухание в полосе задерживания - 67.83дб. Но это лечится либо добиванием еще одного звена, либо заменой одного звена(пары звеньев) на фильтр {1,1,1} (тогда результат будет -2.1дб и -73.2дб соответственно)

ЗЗЫ. полоса 50МГц на 75MSPS? Она комплексная? Тогда верхняя частота в каждом канале будет 25, а не 50. Еще меньше влияние фильтра.

Наконец сеть появилась.
Денис, спасибо! Это решение лучшее из всех рассмотренных мной. А как посмотреть АЧХ 6 звеньев подряд?

я делаю так (это маткад, он нагляднее).

Спорить не буду ибо считаю Вас одним из самых компетентных посетителей конференции.На картинке изображено три секции CIC - значит фильтр третьего порядка. Я в своём фильтре задаю 5 секций, поэтому и думал что работаю с фильтром пятого порядка. Как из параметров R,M,N получить номер порядка?

Order = R*M


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

Ой, с маткадом уже лет 8 не работаю... А в матлабе есть такая возможность?


Ни фига не пойму. А как же Лайонс , врёт что ли? Только от М... А М у него - количество секций.

Спасибо!

Порядок звена да, а порядок фильтра нет. Тем и плохи CIC фильтры, что без децимации внутри не работают. Точнее вырождаются в фильтры скользящего среднего (допускающие усекание разрядов кстати). В приложении теория "на пальцах"

Спасибо!

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

Это Лайонс и есть, 1 в 1. Можно то же прочитать на русском в книге.

Спасибо за лестный отзыв о моей компетентности. Что касается порядка cic:

порядок определяется количеством нулей передаточной функции. 3 нуля - 3-й порядок (N=4 коэффициента эквивалентного однородного), 32 - 32-го (N=33 к-та).
числа 3, 5 и т д означают количество каскадов или кратность нулей передаточной функции cic.

Смысл использования cic при децимации/интерполяции заключается в том, что такой фильтр имеет нули точно на частотах Fd/N*k k = 1 ... N-1 и при децимации в N раз накладываются на частоты 0 и Fd. Избирательные свойства однородного фильтра весьма убоги, поэтому используют последовательное соединение нескольких таких фильтров, что в свою очередь приводит к значительной неравномерности в полосе полезного сигнала.

В вашем случае нужна децимация в 16 раз (1200/75) при этом вы получаете полезный сигнал с 1.5-оверсэмплингом т е нужен фильтр с 16-ю к-тами. Результат можно скорректировать (если есть смысл) ких-фильтром без децимации.

Что получается?

p=ones(1,16);% однородный эквивалент cic1
h=conv(conv(conv(conv(p,p),p),p),p);%эквивалент cic5

plot(0:1200/10000:1200-1200/10000,20*log10(abs(freqz(h,1,10000,'whole'))))

Видно, что на частотах 0... 25мгц подавление наложенных частот будет не менее 40 дб. Но и неравномерность в полосе будет 8 дб. Если такое подавление недостаточно - делаем cic c децимацией на 8 и 2 с ких. итд

Теперь про разрядность:

К-т усиления одного звена будет N, к звеньев будут иметь к-т передачи N^k. Поэтому разрядность аккумуляторов интеграторов выбирают ceil(k*log2(N))+ разрядность исходных данных. После интеграторов результат прореживается в N раз и обрабатывается каскадом гребенок которая вырождается в фильтр y(n)=x(n)-x(n-1). Затем результат приводится к нужной разрядности через правый сдвиг на М разрядов. Обычно, М=ceil(k*log2(N)), т е к исходной разрядности. К-т передачи такого фильтра будет (N^k) / (2^ceil(k*log2(N)))

Большое спасибо за развёрнутый ответ! Красиво Вы свёртки посчитали,а я каждую отдельно считал...

Кстати, а почему у Вас первый ноль на частоте 75МГц? Я думал нужно чтобы первый ноль был на частоте 75/2 == 37,5МГц.

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

А,так вот почему у меня получилось 50 дБ неравномерности...Спасибо!

Потихоньку двигаюсь дальше. Скажите, а как правильно посмотреть суммарную АЧХ такой цепочки : fir1-децимация-fir2? Допустим в качестве fir1 выступает фильтр с h1=conv(conv(conv(conv(p,p),p),p),p);
Сигнал с его выхода децимирую на 8 и подаю на второй фильтр с ИХ h2.

Так, БПФ от результирующей импульсной характеристики и даст АЧХ.

Даст, и что? Вы вопрос читали внимательно? Как учесть децимацию? Чтобы работать со вторым фильтром нужно перенести АЧХ первого фильтра в новый частотный домен (Fs2 = Fs1 / 8). И вопрос как это сделать через импульсную характеристику(во временнОй области)?

Я вижу, мои ответы вас раздражают. Побережем нервы.

Абсолютно нет. Я спокоен как селёдка.

Анатолий, здесь все очень просто.
АЧХ децимируюещй системы (DDC в наших случаях) мы смотрим на частоте дискретизации АЦП.
С первым фильтром fir1 все понятно, он работает до децимации на 2.
Чтобы привести частоте дискретизации АЦП АЧХ второго фильтра fir2, который стоит после децимации на 2, мы должны его импульсную характристику h2 интероплировать на 2.
После этого берем свертку от h1 и интерполированной h2, и затем АЧХ от результата свертки.

Общее правило такое, если перед каким-то фильтром есть децимация на N, то для приведения его к частоте дискретизации АЦП, надо его импульсную характеристику интерполировать на N.
И так с каждым каскадом DDC. Потом свертка всех интероплированных результатов.


h1 = [-1 4 10 4 -1] ./ 16;
h2 = [0.00142192172987358 -0.00222268508122243 -0.0137433922997260 -0.0322604508637608 -0.0458124362883413 -0.0346645776223707 0.0163109682969986 0.102976693990851 0.197188326792324 0.258848499145519 0.258848499145519 0.197188326792324 0.102976693990851 0.0163109682969986 -0.0346645776223707 -0.0458124362883413 -0.0322604508637608 -0.0137433922997260 -0.00222268508122243 0.00142192172987358];

h_stage1 = conv(h1, h1); % Каскад 1 (без децимации

h_stage2 = h2;
h_stage2 = upsample(h_stage2, 2); % С учетом децимации на 2.

% Результирующая импульсная характеристика.
h_ddc = conv(h_stage1, h_stage2);
fvtool(h_ddc);

Большое спасибо за квалифицированный ответ! А я пытался наоборот - децимировать h1 и естественно ничего не получалось. Поехал дальше.

p.s.
А вот эта строчка зачем?

Зачем сворачивать самого с собой?

h=conv(h1,8*upsample(h2,8));

дальше freqz

Вот это ответ! Проще не бывает. Спасибо!

В моем примере стояло два фильтра h1 перед децимацией на 2, забыл это упомянуть.


Децимировать h надо если перед фильтром есть интеропляция, например когда делаете DUC.

Спасибо! Я понял свою ошибку.

Вот такая сквозная АЧХ моего DDC получилась. Спасибо всем участникам обсуждения!

Скажите,пожалуйста, можно ли это сделать в частотной области?

Можно.

Вот что-то читаю и никак не могу понять...


Полосу сигнала в 50МГц при частоте дискретизации 75МSPS??? Где-то ошибка?


Тут имелось ввиду "задержка" а не порядок? Порядок фильтра обеспечивает подавление лепестков, насколько я это понимаю (http://www.dsplib.ru/content/cic/cic.html). Т.е. кол-во нулей при разных порядках фильтра остается тем же, отличается только амплитуды "лепестков" в АЧХ. Разве не так?

Первый вопрос снят... Наверняка имелось ввиду полоса с центром на 0Гц.