Реализация нерекурсивного фильтра CIC Опции

[_Anatoliy]

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

Эскизы прикрепленных изображений

 

[ViKo]

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

[_Anatoliy]

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

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

[ViKo]

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

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

[_Anatoliy]

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

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

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

Ясно.

[ViKo]

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

[_Anatoliy]

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

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

Эскизы прикрепленных изображений

 

[bogaev_roman]

Видно офигенное подавление боковых лепестков, но так же видно и коэффициент передачи 150дБ!

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

[des00]

Мне нужно продецимировать сигнал полосой 50МГц в 16 раз из 1200 MSPS в 75 MSPS.

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

[_Anatoliy]

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

Всё правильно, но вот насчёт

Код

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

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

[thermit]

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

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

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

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

[des00]

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МГц в 16 раз из 1200 MSPS в 75 MSPS

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

[_Anatoliy]

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

[des00]

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

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

Эскизы прикрепленных изображений

 

[_Anatoliy]

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

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

Эскизы прикрепленных изображений

 

[des00]

Как из параметров R,M,N получить номер порядка?

Order = R*M

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

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

[_Anatoliy]

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

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

Order = R*M

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

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

Спасибо!

Эскизы прикрепленных изображений

 

[des00]

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

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

Прикрепленные файлы

 cic_filter_aplication.7z ( 922.98 килобайт ) Кол-во скачиваний: 47

 

[_Anatoliy]

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

Спасибо!

[des00]

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

[ViKo]

В приложении теория "на пальцах"

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

[thermit]

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

Спасибо за лестный отзыв о моей компетентности. Что касается порядка 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)))

[_Anatoliy]

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

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

[thermit]

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

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

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

[_Anatoliy]

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

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

[_Anatoliy]

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

[ViKo]

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

[_Anatoliy]

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

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

[ViKo]

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

[_Anatoliy]

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

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

[soldat_shveyk]

Скажите, а как правильно посмотреть суммарную АЧХ такой цепочки : fir1-децимация-fir2?

Анатолий, здесь все очень просто.
АЧХ децимируюещй системы (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);

[_Anatoliy]

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

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

Код

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

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

[thermit]

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

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

дальше freqz

[_Anatoliy]

h=conv(h1,8*upsample(h2,8));
дальше freqz

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

[soldat_shveyk]

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

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

А я пытался наоборот - децимировать h1 и естественно ничего не получалось.

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

[_Anatoliy]

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

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

[_Anatoliy]

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

Эскизы прикрепленных изображений

 

[Кнкн]

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

дальше freqz

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

[thermit]

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

Можно.

[honeycomb0]

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

Мне нужно продецимировать сигнал полосой 50МГц в 16 раз из 1200 MSPS в 75 MSPS.

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

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

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

[honeycomb0]

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

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