Перейти к содержанию
Поиск
    

Определение параметров фильтра

Ответить

strap89    0

strap89
Опубликовано · Жалоба

Здравствуйте!

Помогите разобраться пожалуйста. Есть цифровой БИХ-фильтр 2-го порядка, режекторный.

Есть алгоритм расчета его коэффициентов. И есть пара коэффициентов, определяющих параметры фильтра (это мои предположения).

Не могли бы Вы подсказать мне, какой тип фильтра тут использован и что определяют его параметры. Сам алгоритм:

void makefiltr2coeff(filtr2& flt,double Tsam,double f0,double f1,double ksi0,double ksi1)

{

if(!f0||!f1||!Tsam||!ksi0||!ksi1)

return;

Tsam*=scaleTsam;

double om0=2.*M_PI*f0;

double om1=2.*M_PI*f1;

double tau0=.5/tan(Tsam*om0/2.);

double tau1=.5/tan(Tsam*om1/2.);

double d=1.+(4.*tau1*(tau1+ksi1));if(!d) return;

double b0=(1.+(4.*tau0*(tau0+ksi0)))/d;

double b1=2.*(1.-4.*tau0*tau0)/d;

double b2=(1.+(4.*tau0*(tau0-ksi0)))/d;

double a1=2.*(1.-4.*tau1*tau1)/d;

double a2=(1.+(4.*tau1*(tau1-ksi1)))/d;

a0=1;

Сама реализация фильтра есть, фильтр работает, но формулу расчета коэффициентов мне найти не удалось.

Как соответственно и информацию, что означают его параметры.

Поделиться сообщением

Ссылка на сообщение
Поделиться на другие сайты

Самурай    12

Самурай
Опубликовано · Жалоба
Здравствуйте!

Помогите разобраться пожалуйста. Есть цифровой БИХ-фильтр 2-го порядка, режекторный.

Есть алгоритм расчета его коэффициентов. И есть пара коэффициентов, определяющих параметры фильтра (это мои предположения).

Не могли бы Вы подсказать мне, какой тип фильтра тут использован и что определяют его параметры. Сам алгоритм:

void makefiltr2coeff(filtr2& flt,double Tsam,double f0,double f1,double ksi0,double ksi1)

{

if(!f0||!f1||!Tsam||!ksi0||!ksi1)

return;

Tsam*=scaleTsam;

double om0=2.*M_PI*f0;

double om1=2.*M_PI*f1;

double tau0=.5/tan(Tsam*om0/2.);

double tau1=.5/tan(Tsam*om1/2.);

double d=1.+(4.*tau1*(tau1+ksi1));if(!d) return;

double b0=(1.+(4.*tau0*(tau0+ksi0)))/d;

double b1=2.*(1.-4.*tau0*tau0)/d;

double b2=(1.+(4.*tau0*(tau0-ksi0)))/d;

double a1=2.*(1.-4.*tau1*tau1)/d;

double a2=(1.+(4.*tau1*(tau1-ksi1)))/d;

a0=1;

Сама реализация фильтра есть, фильтр работает, но формулу расчета коэффициентов мне найти не удалось.

Как соответственно и информацию, что означают его параметры.

 

 

Ваш IIR 2-ого порядка получен через билинейное преобразование аналогового фильтра-прототипа общего вида:

 

?H(s)=\frac{s^2+2\xi_0\omega_0 s+\omega_0^2}{s^2+2\xi_1\omega_1 s+\omega_1^2}

 

т.е. через замену ?s=\frac{2}{T_{sam}}(\frac{1-z^{-1}}{1+z^{-1}}) и ?\omega_{0,1}=\frac{2}{T_{sam}}tg(\frac{\omega_{0,1}T_{sam}}{2})

 

Выводить не буду, там все тривиально, но больно много писанины...

 

Для режекторного фильтра необходимо выполнение условий: f0 = f1 = частоте подавления, ksi0 < ksi1, оба эти коэффициента влияют на уровень затухания и ширину полосы подавления, если ksi0=0, подавление максимальное.

 

 

Поделиться сообщением

Ссылка на сообщение
Поделиться на другие сайты

strap89    0

strap89
Опубликовано · Жалоба
Ваш IIR 2-ого порядка получен через билинейное преобразование аналогового фильтра-прототипа общего вида:

 

?H(s)=\frac{s^2+2\xi_0\omega_0 s+\omega_0^2}{s^2+2\xi_1\omega_1 s+\omega_1^2}

 

т.е. через замену ?s=\frac{2}{T_{sam}}(\frac{1-z^{-1}}{1+z^{-1}}) и ?\omega_{0,1}=\frac{2}{T_{sam}}tg(\frac{\omega_{0,1}T_{sam}}{2})

 

Выводить не буду, там все тривиально, но больно много писанины...

 

Для режекторного фильтра необходимо выполнение условий: f0 = f1 = частоте подавления, ksi0 < ksi1, оба эти коэффициента влияют на уровень затухания и ширину полосы подавления, если ksi0=0, подавление максимальное.

Спасибо за объяснение.

Поделиться сообщением

Ссылка на сообщение
Поделиться на другие сайты

Ruslan1    16

Ruslan1
Опубликовано · Жалоба

Вопрос про усиление: у меня по этим формулам получается та же форма АЧХ, но только с ослаблением сигнала, от 0 dB и вниз, а если через Матлабовскую функцию- то все как у оригинала H(s), где-то на +30dB.

Какую константу я забыл?

 

Имею очень похожую формулу:

H(S) =    

  s^2 + 31.1 s + 799.4
  ---------------------
  s^2 + 7.109 s + 25.27

 

если делаю в Матлабе: 

b =     1   31.1  799.4
a =   1    7.109   25.27
fs = 1000
[bd,ad] = bilinear(b,a,fs)

то имею коэффициенты фильтра

Quote

bd =     1.0121   -1.9925    0.9812

ad =    1.0000   -1.9929    0.9929

И по ним fvtool строит именно то, что у аналогового прототипа (та же форма и то же усиление).

 

Но вот если иду через вышеприведенные формулы (чтобы совсем не трогать код, "вытащил" из моих констант эти кси и омега): 

%  s^2 + 31.1 s + 799.4
%  ---------------------
%  s^2 + 7.109 s + 25.27

clear all
M_PI = pi

om1 = sqrt(25.27);
om0= sqrt(799.4)
ksi0 = 31.1 / (2*om0);
ksi1 = 7.109 / (2*om1);
fs = 1000;   % sample rate = 1000 SPS
Tsam=1/fs;

tau0=.5/tan(Tsam*om0/2.);
tau1=.5/tan(Tsam*om1/2.);
d=1.+(4.*tau1*(tau1+ksi1));
b0=(1.+(4.*tau0*(tau0+ksi0)))/d;
b1=2.*(1.-4.*tau0*tau0)/d;
b2=(1.+(4.*tau0*(tau0-ksi0)))/d;
a1=2.*(1.-4.*tau1*tau1)/d;
a2=(1.+(4.*tau1*(tau1-ksi1)))/d;
a0=1;

то имею коэффициенты:
 

Quote

 

bd =    0.0320   -0.0630    0.0310

ad =     1.0000   -1.9929    0.9929

вижу [bd], уменьшенный примерно в 31 раз  если сравнить с результатом функции матлаба.

Поделиться сообщением

Ссылка на сообщение
Поделиться на другие сайты

Самурай    12

Самурай
Опубликовано (изменено) · Жалоба

Все дело в том, что в формулах ТС есть некоторые упрощения, которые не влияют на результат в частном случае, когда w0 = w1. Для общего случая, формулы нужно слегка поправить :

 

double om0=2.*M_PI*f0;
double om1=2.*M_PI*f1;

 

double tau0=.5/tan(Tsam*om0/2.);
double tau1=.5/tan(Tsam*om1/2.);

 

double d=1.+(4.*tau1*(tau1+ksi1));
double d0=1.+(4.*tau0*(tau0+ksi0));

 

double k=(d0/d)*(tau1/tau0)^2;

 

double b0=k*(1.+(4.*tau0*(tau0+ksi0)))/d0;
double b1=k*2.*(1.-4.*tau0*tau0)/d0;
double b2=k*(1.+(4.*tau0*(tau0-ksi0)))/d0;

 

double a1=2.*(1.-4.*tau1*tau1)/d;
double a2=(1.+(4.*tau1*(tau1-ksi1)))/d;
a0=1;

Изменено пользователем Самурай

Поделиться сообщением

Ссылка на сообщение
Поделиться на другие сайты

Ruslan1    16

Ruslan1
Опубликовано · Жалоба
22 hours ago, Самурай said:

Все дело в том, что в формулах ТС есть некоторые упрощения, которые не влияют на результат в частном случае, когда w0 = w1. Для общего случая, формулы нужно слегка поправить :

...

 

Уважаемый Самурай, все получилось!

Спасибо огромное!!

 

С помощью Вашей подсказки и Википедии получилось сделать универсальный вариант, с любыми нумераторами-деноминаторами 2-й степени (в стартовом коде было a0=b0=1).

Оставлю тут , может пригодится кому-нибудь: 

% second-order bilinear transform
% input:
% Ha(s) = ( b0*s^2 + b1*s + b2 ) / (a0*s^2 + a1*s + a2)
%
% output: (digital biquad filter, a0 = 1)
% Hd(z) = (b0 + b1*z^-1 + b2*z^-2) / (a0 + a1*z^-1 + a2* z^-2)
%
% reference: 
% #1 (General second-order biquad transformation):  https://en.wikipedia.org/wiki/Bilinear_transform
% #2 (correction for case w0 != w1):  https://electronix.ru/forum/index.php?app=forums&module=forums&controller=topic&id=148759
clear all
%% input data
b_in = [23 31.1 799.4];
a_in = [19 7.109 25.27];
fs = 1000;

%% read data from input arrays
b0 = b_in(1);
b1 = b_in(2);
b2 = b_in(3);
a0 = a_in(1);
a1 = a_in(2);
a2 = a_in(3);

%% calculation
w0 = sqrt(b2);
w1= sqrt(a2)
T = 1/fs;
K0 = w0/tan(w0*T/2);  
K1 = w1/tan(w1*T/2);
%not normalled coefficients
bn0 = b0*K0^2 + b1*K0+b2;
bn1 = 2*b2-2*b0*K0^2;
bn2 = b0*K0^2 - b1*K0 + b2;
an0 = a0*K1^2 + a1*K1 + a2;
an1 = 2*a2 - 2*a0*K1^2;
an2 = a0*K1^2 - a1*K1 + a2;
% correction for denominator
k = (K1/K0)^2;
%k=1
% Normally the constant term in the denominator must be normalized to 1
bd0 = k* bn0 / an0;
bd1 = k* bn1 / an0;
bd2 = k* bn2 / an0;
ad0 = an0 / an0;
ad1 = an1 / an0;
ad2 = an2 / an0;
%% output data
bd = [bd0 bd1 bd2]
ad = [ad0 ad1 ad2]

%% test
[bd_matlab,ad_matlab] = bilinear(b_in,a_in,fs);
ad 
ad_matlab 
bd 
bd_matlab
fvtool (bd_matlab,ad_matlab);
fvtool (bd,ad)

 

Поделиться сообщением

Ссылка на сообщение
Поделиться на другие сайты

Присоединяйтесь к обсуждению

Вы можете написать сейчас и зарегистрироваться позже. Если у вас есть аккаунт, авторизуйтесь, чтобы опубликовать от имени своего аккаунта.

Гость
Ответить в этой теме...

×   Вставлено с форматированием.   Вставить как обычный текст

  Разрешено использовать не более 75 эмодзи.

×   Ваша ссылка была автоматически встроена.   Отображать как обычную ссылку

×   Ваш предыдущий контент был восстановлен.   Очистить редактор

×   Вы не можете вставлять изображения напрямую. Загружайте или вставляйте изображения по ссылке.

×
×
×
  • Создать...