filterc example.c

自己编的数值计算方法的c语言实现源程序。

#include"stdio.h"
#include"filterc.c"
/*求H(z)的单位冲激响应；
H(z)=0.001836(1+z^-1)^4/(1-1.4996z^-1+0.8482z^-2)(1-1.5548z^-1+0.6493z^-2)*/
void main()
{	int i,j,n,ns,len,nblk;
	double b[2][3]={{0.001836,0.003672,0.001836},{1.0,2.0,1.0}};
	double a[2][3]={{1.0,-1.4996,0.8482},{1.0,-1.5548,0.6493}};
	double px[2][3],py[2][3];
	double x[25],data[100];
	FILE *fp;
	n=2;	/*每节的阶数*/
	ns=2;   /*节数L*/
	len=25;	/*输入和输出序列长度*/
	nblk=4;	/*块*/
	/*输入长度为len的单位冲激响应*/
	x[0]=1.0;
	for(i=1;i<len;i++)
	{	x[i]=0.0;	}
	/*初始化px,py;*/
	for(i=0;i<ns;i++)
	for(j=0;j<=n;j++)
	{	px[i][j]=0.0;
		py[i][j]=0.0;
	}
	/*计算级联型IIR数字滤波；*/
	for(j=0;j<nblk;j++)
	{	filterc(b,a,n,ns,x,len,px,py);
		for(i=0;i<len;i++)
		{	data[j*len+i]=x[i];
			x[i]=0.0;
		}
	}
	/*输出结果*/
	printf("Unit Impulse Response.\n");
	if(( fp=fopen("filterc.txt","w"))==NULL)
	{	printf("cannot open file 'filterc.txt'!\n");
		exit(0);
	}
	for(i=0;i<100;i++)
	{	fprintf(fp,"%lf ",data[i]); }
	fclose(fp);
}