filter.txt

实现DSP编程中重要的二阶巴特沃思滤波器的实现

#include "stdio.h"
  #include "math.h"
  #define pi 3.1415926 
  
  int ns;
  double fs,nlpass,nlstop,nhpass,nhstop,a[3],b[3],x[130],y[130];

  void biir2lpdes(double fs, double nlpass, double nlstop, double a[], double b[]);

/*
	Fs为采样频率，f1、f2分别对应2路输入信号x1、x2之频率；
	x为2路输入信号x1、x2之和；
	y为滤波器滤波结果；
	biir2lpdes为低通滤波器设计子程序；

*/ 
  void biir2lpdes(double Fs, double nlpass, double nlstop, double a[], double b[])
  
  {
    int i,u,v;
	double wp,omp,gsa,t;
    wp=nlpass*2*pi;
	omp=tan(wp/2.0);
	gsa=omp*omp;
	
	for (i=0; i<=2; i++)
	  {
		u=i%2;
		v=i-1;
		a[i]=gsa*pow(2,u)-sqrt(2.0)*omp*v+pow(-2,u);
	  }
	
	for (i=0; i<=2; i++)
	  {	u=i%2;
		b[i]=gsa*pow(2,u);
	  }
	t=a[0];
	for (i=0; i<=2; i++)
	  {	a[i]=a[i]/t;
		b[i]=b[i]/t;
	  }
  }

  void biir2filter(double x[], double y[], double a[], double b[], int ns)
  {
	int i;
	double c2,c1,c0;
	c2=c1=c0=0.0;
	for (i=0; i<=1; i++)
	  {
	    x[i]=0;
	    y[i]=0;
	  }
    for (i=2; i<=ns-1; i++)
	  {
		c2=x[i]-a[1]*c1-a[2]*c0;
      	y[i]=b[0]*c2+b[1]*c1+b[2]*c0;
      	c0=c1;
      	c1=c2;
      }
  }
 
  
 main()
  {	int i;
	double Fs,f1,f2;
	double x1[130],x2[130];
	
	Fs=1000;
	f1=100;
	f2=400;
	
	ns=130;
	
	for (i=0; i<=1; i++)
	  {	
	  	x[i]=x1[i]=x2[i]=0;
	  }
	for (i=2; i<=ns-1; i++)
	  {
		x1[i]=sin(2*pi*f1*((i-2)/fs));
		x2[i]=sin(2*pi*f2*((i-2)/fs));
				
		x[i]=x1[i]+x2[i];
      }

//	IIR LP
	nlpass=0.1;
	nlstop=0.3;
	biir2lpdes(fs,nlpass,nlstop,a,b);
	
	biir2filter(x,y,a,b,ns);

	for (;;)
	{
	    i=0;  /* 在此处设断点 */	
	}
}
/********************************************************************
** End of File
	本实验为2阶Butterworth型IIR低通滤波器算法程序
	x: 为输入信号，为两正弦波信号之和（频率分别为f1,f2;Fs为采样频率），
	      长度130，32位浮点型；
	y: 为IIR低通滤波器滤波结果，长度130，32位浮点型；
***********************************/