Преобразования ФНЧ-ПФ

Все что касается фильтрации
Аватара пользователя
DimKaKiber
Сообщения: 29
Зарегистрирован: 26 мар 2020, 15:00
Откуда: Томск
Контактная информация:

Преобразования ФНЧ-ПФ

Сообщение DimKaKiber »

Добрый день! Понадобыилось реализовать библиотеку для расчеты БИХ фильтров в режиме реального времени. Будут перестраиваться в зависимости от входного ряда.

Пишу на С++ с вкраплениями С (если это важно). Все действия реализовал строго по Вашим статьям, но споткнулся о преобразования ФНЧ-ПФ. Здесь на выходе получается в два раза больше коэффициентов полинома. Это по статьям. Я не могу добиться такого эффекта ((
На вход дробно-рациональной подставновки подаю вектора коэффициентов аналогового прототипа единой длины.
Вектора a и b подаю размером 3:

a = [1,0,Wc2]
b = [0,B,0]

В результате получаю только первую половину коэффициентов. Остальное нули. Подскажите, пожалуйста, все ли правильно я делаю?

Аватара пользователя
DimKaKiber
Сообщения: 29
Зарегистрирован: 26 мар 2020, 15:00
Откуда: Томск
Контактная информация:

Re: Коэффициенты БИХ фильтров

Сообщение DimKaKiber »

Здравствуйте! Работаю над библиотекой для расчета БИХ фильтров и споткнулся о преобразование ФНЧ-ПФ.
Следую указаниям из Ваших статей. Не могу понять какие вектора подавать на вход дробно-рациональной подстановки. Она работает верно - Ваш пример обсчитывается правильно + преобразования ФНЧ-ФНЧ и ФНЧ-ФВЧ также выдают корректные результаты.

С преобразованием же ФНЧ-ПФ ступор - получаю неверные значения коэффициентов числителя и ровно половину коэффициентов знаменателя. Аналоговый прототип эллиптического фильтра (по Вашей статье) рассчитывается корректно.

Реализация подстановки (матрицы - мой класс. Здесь используется только для удобства создания/удаления):

Код: Выделить всё

void IIRcreator::rationalSubstitution(std::vector<double> &numerator,std::vector<double> &denumerator, double *a,int  aSize,double *b,int bSize)
{
 if (Grade == 0) return;
 if (a == nullptr || b == nullptr) return;
 if (numerator.size() == 0 || denumerator.size() == 0) return;

 int n = Grade+1,
     m = Grade+1; //полиномы для подстановок у меня все в первой степени изначально, поэтому формулу N*P-1 вполне справедливо модернизировать так как я ее записал.

 if (Type == FilterType::BPF || Type == FilterType::BSF) //Р = 2. Это значит, что длина строки становится равной 2*N+1.
  {
   m = 2*Grade+1;
  }

 Matrix NUM(n,m),
        DEN(n,m);
 NUM.data[0][0] = 1.0;
 DEN.data[0][0] = 1.0;
 int num_size = 1,
     denum_size = 1,
     buffSize = 2*m;
 double *buff = (double*)calloc(buffSize,sizeof(double));
 for (int i = 1; i < n; ++i)
  {
   //линейная свертка для строки  NUM
   linearConv(buff,NUM.data[i-1],b,buffSize,num_size,bSize,num_size);
   memcpy(NUM.data[i],buff,sizeof(double)*num_size);
   //линейная свертка для строки DEN
   linearConv(buff,DEN.data[i-1],a,buffSize,denum_size,aSize,denum_size);
   memcpy(DEN.data[i],buff,sizeof(double)*denum_size);
  }

 Matrix ND(n,m);
 int nd_size = 0;
     denum_size = Grade+1;
     num_size = 1;
 for (int i = 0; i < n; ++i,--denum_size,++num_size)
  {
   linearConv(buff,DEN.data[Grade-i],NUM.data[i],buffSize,denum_size,num_size,nd_size);
   memcpy(ND.data[i],buff,sizeof(double)*nd_size);
  }

 ND.debugMatrix("ND");
 //построчное умножение на вектора коэффициентов.
 for (int i = 0; i < n; ++i)
  {
   for (int j = 0; j < m; ++j)
    {
     NUM.data[i][j] = ND.data[i][j]*numerator[i];
     DEN.data[i][j] = ND.data[i][j]*denumerator[i];
    }
  }

 NUM.debugMatrix("Num");
 DEN.debugMatrix("DEN");

 //теперь суммирую значения по столбцам и получаю новые коэффициенты.
 memset(numerator.data(),0,sizeof(double)*numerator.size());
 memset(denumerator.data(),0,sizeof(double)*denumerator.size());
 for (int i = 0; i < n; ++i)
  {
   for (int j = 0; j < n; ++j)
    {
     numerator[i] += NUM.data[j][i];
     denumerator[i] += DEN.data[j][i];
    }
  }

 free(buff); buff = nullptr;
}

Для подстановки на вход подаю вектора:
a = [1,0,Wp^2]
b = [0,B,0]

Направьте, пожалуйста, в нужную сторону? Уже третий день не могу понять как все правильно сделать.....

Аватара пользователя
Бахурин Сергей
Администратор
Сообщения: 1114
Зарегистрирован: 05 окт 2010, 19:55
Контактная информация:

Re: Преобразования ФНЧ-ПФ

Сообщение Бахурин Сергей »

1. Да коэффициентов должно быть в два раза больше, потому что подставляются полиномы второй степени и степени перемножаются в результате.
2. Проверяйте функцию подстановки если она не перемножает полиномы.
3. Можете взять мою си библиотеку, в которой все это уже есть (хотя возможно не до конца задокументировано), но функции подстановки и частотных преобразований документированы. Там же может посмотреть исходники.

Аватара пользователя
DimKaKiber
Сообщения: 29
Зарегистрирован: 26 мар 2020, 15:00
Откуда: Томск
Контактная информация:

Re: Преобразования ФНЧ-ПФ

Сообщение DimKaKiber »

Большое спасибо за отклик!

Как оказалось - у меня была ошибка внутри функции дробно-рациональной подстановки. Неверно формировал матрицу ND. Не заметил, что размер строк начинается с половины результирующей длины. Для квадратной матрицы все работало верно. Для большего количества столбцов (чем строк) полезла неприятная ошибка.

Спасибо за предложение попользоваться Вашей библиотекой, но на Вашем ресурсе итак все разжевано и привык знать как мой код работает. Удобнее потом ошибки исправлять и дорабатывать. Да и эллиптический фильтр крови попил пока реализовывал - после такого свой код не брошу :))

Еще, как оказалось (это уже у Вас в коде подсмотрел) вектор a{} у меня зеркальным должен был быть: a = [Wc2,0,1]. Только, если судить по формуле в статье изначально все верно подавал.

Хотел бы еще вопросы задать:
1. На форуме вычитал что при проектировании ФВЧ и РФ предварительно надо сделать преобразование ФНЧ-ФНЧ - это обязательный процесс и после него не надо осуществлять нормировки никакие? Сразу гоню результат на дальнейшее преобразование в ФВЧ или РФ?
2. Предискажение частот обязательно при проектировании каждого из фильтров?

Еще раз большое спасибо за ответ!

Аватара пользователя
Бахурин Сергей
Администратор
Сообщения: 1114
Зарегистрирован: 05 окт 2010, 19:55
Контактная информация:

Re: Преобразования ФНЧ-ПФ

Сообщение Бахурин Сергей »

порядок коэффициентов вектора выбирается в соответствии с реализацией полинома. У меня в реализации a[0] это коэффициент при нулевой степени a[1] при первой степени и т.д. Другие пакеты (матлаб октава) наоборот a[0] это при старшей степени.

Да при расчете режекторного фильтра я делаю преобразование для того чтобы в заданной полосе [w0 w1] гарантировать подавление Rs дБ. Т.о. существуют формулы пересчета частоты на которой фильтра имеет Rs и я делаю перенос в этой частоты на 1 рад/c так после преобразования режекторный фильтр обеспечивает подавление. Формулы я выводил, нигде их не встречал, если надо поделюсь.

Что имеется ввиду предыскажение частот?

Аватара пользователя
DimKaKiber
Сообщения: 29
Зарегистрирован: 26 мар 2020, 15:00
Откуда: Томск
Контактная информация:

Re: Преобразования ФНЧ-ПФ

Сообщение DimKaKiber »

Спасибо. Пока что режекторный фильтр не применял к своим задачам - мне нужно выделать низкие частоты в диапазоне 20-40Гц, поэтому пользуюсь ПФ в основном. Для общего развития было бы интересно.

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

Предискажение отсюда взял - http://www.dsplib.ru/content/filters/bu ... terex.html, формула (1). Но вопрос, наверное, глупый, т.к. тут же есть и объяснение.

Аватара пользователя
Бахурин Сергей
Администратор
Сообщения: 1114
Зарегистрирован: 05 окт 2010, 19:55
Контактная информация:

Re: Преобразования ФНЧ-ПФ

Сообщение Бахурин Сергей »

что значит много новых синусоид появляется при амплитудной модуляции?

Можно еще SSB перенос сделать. Останется верхняя боковая или нижняя боковая полоса.

Аватара пользователя
DimKaKiber
Сообщения: 29
Зарегистрирован: 26 мар 2020, 15:00
Откуда: Томск
Контактная информация:

Re: Преобразования ФНЧ-ПФ

Сообщение DimKaKiber »

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

Поэтому постараюсь проиллюстрировать. Работаю с фонокардиограммами плода, у которых соотношение сигнал шум минимальное до обработки. После некоторой достаточно муторной разноплановой магии (в т.ч. и фильтрации) в случае наличия полезного сигнала (тон или тоны сердца плода в зависимости) из шума все-таки достаю и хочу дать послушать полученный сигнал акушерам или пациентам, которые пользуются прибором и приложением со специфическими алгоритмами. Сигнал получается низкоамплитудным (да и Бог тут с ним...не проблема совсем) и в собственном диапазоне частот (20-40 Гц). Этот сигнал, при прослушке только "чувствуется" в относительно нормальных наушниках на уровне ощущений - колебания их мембраны на перепонке. При этом его приходится достаточно хорошо усиливать. Что хорошо, но только в offline. При последовательном приеме и воспроизведении пакетов такая роскошь мне не светит - буду искажать сигнал. Соответственно задача - сделать его нормально слышимым, но с таким же звучанием что я его слышу сейчас. Поэкспериментировал с линейным переносом спектра - становится лучше слышно, но все искажается хрюками (до 150 Гц несущая) и дальше превращается в звуки сонара.

А много синусоид - это вот такая картина как на приложенных изображениях (частота дискретизации 44,1 кГц. Привел сигнал к ней из-за того что Adnroid гарантированно воспроизводит такие сигналы):
1. Исходное после Фурье фильтра с 50% перекрытием пакетов в потоке (оставил нужный частотный диапазон, затем по результату обратного БПФ прошелся окном Бартлета-Хана).
Вид после фильтрации
Вид после фильтрации


2. Тоже самое только перенес все на несущую (с оглядкой на частотное разрешение окна БПФ) 150 Гц. Т.е. условно частоты своего диапазона перенес на +150Гц.
После линейного переноса несущей частоты на 150Гц
После линейного переноса несущей частоты на 150Гц
3. Тоже самое что и (2), только несущая - 1,5 кГц.
После линейного переноса несущей частоты на 1,5 кГц
После линейного переноса несущей частоты на 1,5 кГц
По теории прекрасно понимаю почему все так происходит, но все-таки верю что должна быть какая то математическая магия которая даст мне возможность насладиться неизменным звуком на более высоких частотах.

Аватара пользователя
Бахурин Сергей
Администратор
Сообщения: 1114
Зарегистрирован: 05 окт 2010, 19:55
Контактная информация:

Re: Преобразования ФНЧ-ПФ

Сообщение Бахурин Сергей »

Честно говоря я не совсем понимаю что значит неизменный звук на более высокой частоте? Если я перенесу ноту ДО на частоту ноты ЛЯ, это будет низменная нота ДО или нет? Сможете ли вы узнать ноту ДО на частоте ноты ЛЯ? То что вы показали это всего лишь амплитудная модуляция вашего НЧ сигнала с переносом на 150 и 1500 Гц. Очевидно звук будет абсолютно разный. Мне еще не совсем понятно в чем заключается информационная составляющая? например я понимаю, что ЭКГ сердца анализируют по форме пиков и выделяют на графике болезни если они есть. Также применяют различные методы спектрального аналаза ЭКГ сигналов.

Аватара пользователя
DimKaKiber
Сообщения: 29
Зарегистрирован: 26 мар 2020, 15:00
Откуда: Томск
Контактная информация:

Re: Преобразования ФНЧ-ПФ

Сообщение DimKaKiber »

Я это понимаю как звук в области высоких частот, имеющий амплитудное представление как тот же звук, но в области низких (рисунок 1). Не огибающая его, а именно само амплитудное представление. Как при передискретизации. Точек больше, а вид сигнала неизменен.
Эта операция мне нужна только лишь затем, чтобы врач надев наушники услышал, что он установил датчик в область сердца, а не кишечника, например. От графического представления врачи отказываются. Все операции с сигналом прекрасно проходят и в его "родной" области частот :)

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

Ответить