Цифровой фильтр

From

Stas Polozov (2:5000/166.2)

To

Dmitry Luzin

Date

2002-10-03T15:16:24Z

Area

RU.ALGORITHMS

Ничего, Dmitry, сиди, где сидишь. Я не могу тебя видеть, но я слышу,
как дрожат твои пальцы.

24 Сен 02 12:24, Daniel Kamperov -> Dmitry Luzin:

 DK> Ну, вот, простой (но медленный) алгоритм.
 DK>
 DK> Путём некоторого усложнения передаточной функции можно получить
 DK> фильтры более высокого порядка.
  Некоторого усложнения, хех :)
  Короче, пора мне вмешаться :)
  Значит, делаем БИХ фильтр.
  Если конкретно, то фильтр Баттерворта, т.е. фильтр с оптимальными затуханиями
своих звеньев, т.е. его АЧХ не имеет выбросов.
  *Основное положение*: при последовательном соединении 2-х фильтров порядка n
и k получаем фильтр порядка n+k. Что реально означает "последовательное
соединение"? Это значит, ты сначала фильтруешь фильтром n-го порядка, а потом
то, что получилось - фильтром k-го порядка (вобще, по причине того, что фильтр
- линейная система, очередность фильтрации не важна).
  Порядок выбираешь в зависимости от нужной тебе крутизны среза АЧХ.
Возрастание порядка фильтра на 1 соответствует увеличению крутизны среза на
-6Дб/октава или -20 Дб/декада.
  Надо помнить о недостатке БИХ фильтра - наличии переходного процесса в начале
фильтрации.
  Вот исходник фильтра восьмого порядка (-160Дб/декада). Фильтрация сделана
двумя последовательными фильтрами четвертого порядка. Обрати внимание: для
начала фильтрации требуется знать 8 предидущих отфильтрованных значений. То,
что ты забиваешь их нулями, и дает этот неприятный переходный процесс.
  *Ахтунг*: это ФНЧ, чтобы получить ФВЧ, ты должен фильтровать каждым из
фильтров четвертого порядка четвертую производную, умноженную на период
дискретизации в 4 степени.
  Если ничего не понятно, могу залить главу из своего диплома.

-[***]- Начало myfilter_x8.cpp -[***]-
#include <iostream.h>
#include <iomanip.h>
#include <math.h>
#include <stdlib.h>

#define d1 0.19509032201612826784828486847702 //постоянные фильтра
#define d2 0.55557023301960222474283081394853
#define d3 0.83146961230254523707878837761791
#define d4 0.98078528040323044912618223613424

#define fd 93.005952380952380952380952380952 //частота дискретизации
#define K 30000 //число отсчетов

void myfilter (double *data, double *data_filtered1, double *data_filtered2,
double *data_filtered_last1, double *data_filtered_last2, double f0, int N)
{double a_1, b_1, c_1, d_1, a_2, b_2, c_2, d_2, *data_filtered_tmp1,
*data_filtered_tmp2;

a_1 = 2*exp(-d1*2*M_PI*f0/fd)*cos(2*M_PI*f0/fd*sqrt(1-d1*d1))
 +2*exp(-d2*2*M_PI*f0/fd)*cos(2*M_PI*f0/fd*sqrt(1-d2*d2));

b_1 = exp(-2*d1*2*M_PI*f0/fd)+exp(-2*d2*2*M_PI*f0/fd)+
4*exp(-1*(d1+d2)*2*M_PI*f0/fd)*cos(2*M_PI*f0/fd*sqrt(1-d1*d1))*
cos(2*M_PI*f0/fd*sqrt(1-d2*d2));

c_1 = 2*exp(-2*M_PI*f0/fd*(2*d1+d2))*cos(2*M_PI*f0/fd*sqrt(1-d2*d2))+
2*exp(-2*M_PI*f0/fd*(d1+2*d2))*cos(2*M_PI*f0/fd*sqrt(1-d1*d1));

d_1 = exp(-2*2*M_PI*f0/fd*(d1+d2));
/*-------------------------------------------------------*/
a_2 = 2*exp(-d3*2*M_PI*f0/fd)*cos(2*M_PI*f0/fd*sqrt(1-d3*d3))
 +2*exp(-d4*2*M_PI*f0/fd)*cos(2*M_PI*f0/fd*sqrt(1-d4*d4));

b_2 = exp(-2*d3*2*M_PI*f0/fd)+exp(-2*d4*2*M_PI*f0/fd)+
4*exp(-1*(d3+d4)*2*M_PI*f0/fd)*cos(2*M_PI*f0/fd*sqrt(1-d3*d3))*
cos(2*M_PI*f0/fd*sqrt(1-d4*d4));

c_2 = 2*exp(-2*M_PI*f0/fd*(2*d3+d4))*cos(2*M_PI*f0/fd*sqrt(1-d4*d4))+
2*exp(-2*M_PI*f0/fd*(d3+2*d4))*cos(2*M_PI*f0/fd*sqrt(1-d3*d3));

d_2 = exp(-2*2*M_PI*f0/fd*(d3+d4));

data_filtered_tmp1=new double [N+4];
data_filtered_tmp2=new double [N+4];

for (int i=0;i<4;i++)
{data_filtered_tmp1[i]=data_filtered_last1[i];
data_filtered_tmp2[i]=data_filtered_last2[i];
//cout<<data_filtered_tmp[i];
}
//cout<<d_1<<" "<<d_2<<"\n";

for (int i=4;i<N+4;i++)
{data_filtered_tmp1[i]=data_filtered_tmp1[i-1]*a_1-data_filtered_tmp1[i-2]*b_1+
data_filtered_tmp1[i-3]
*c_1-data_filtered_tmp1[i-4]*d_1+(1-a_1+b_1-c_1+d_1)*data[i-4];

data_filtered_tmp2[i]=data_filtered_tmp2[i-1]*a_2-data_filtered_tmp2[i-2]*b_2+d
ata_filtered_tmp2[i-3]
*c_2-data_filtered_tmp2[i-4]*d_2+(1-a_2+b_2-c_2+d_2)*data_filtered_tmp1[i];
}

for(int i=4;i<N+4;i++)
{
data_filtered1[i-4]=data_filtered_tmp1[i];
data_filtered2[i-4]=data_filtered_tmp2[i];
}

delete data_filtered_tmp1;
delete data_filtered_tmp2;
}

void main()
{double in[K], out1[K], out2[K], out_last1[4], out_last2[4];

for (int i=0;i<K;i++)
in[i]=10;

for (int i=0;i<4;i++)
{
out_last1[i]=0;
out_last2[i]=0;
}

myfilter(in,out1,out2,out_last1,out_last2,2,K);

for(int i=0;i<K;i++)
cout<<i<<" "<<setprecision(15)<<out1[i]<<" "<<out2[i]<<"\n";

}
-[***]- Конец myfilter_x8.cpp  -[***]-

... np: Marillion - Three Boats Down From The ...
--- GoldED+/W32 1.1.4.7
 * Origin:  Hollow ambitions in a hollow mind (2:5000/166.2)