resp.c

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

/*计算数字滤波器的单位冲激响应和单位阶跃响应
若输入序列为delta(n),则为单位冲激响应；
若输入序列为u(n),则为单位阶跃响应*/
/*x,y,lx,ly:输入序列参数；b,a,m,n:滤波器参数；*/
void resp(x,y,lx,ly,b,a,m,n)
/*x:双精度实型一维数组，长度为lx,存放滤波器输入序列；
y:双精度实型一维数组，长度为ly,存放滤波器输出序列；
lx:输入序列长度；
ly:输出序列长度；
b:双精度实型一维数组，长度为（m+1）,存放滤波器分子多项式系数b(i);
a:双精度实型一维数组，长度为（n+1）,存放滤波器分母多项式系数a(i);
m:滤波器分子多项式阶数；
n:滤波器分母多项式阶数；*/
int lx,ly,m,n;
double x[],y[],b[],a[];
{	int k,i,i1;
	double sum;
	for(k=0;k<ly;k++)
	{	sum=0.0;
		for(i=0;i<=m;i++)
		{	if((k-i)>=0)
			{	i1=	(k-i)<lx?(k-i):(lx-1); /*输入数据x若为delta(n)，则x[lx-1]=0.0;  若为u(n),则x[lx-1]=1.0;即x[n]=x[lx-1]其中(n>=lx);*/
			     sum=sum+b[i]*x[i1];  /*sum_{k=0}^{m} b(i)*x(k-i);*/
			}
		}
		for(i=1;i<=n;i++)
		{	if((k-i)>=0)
				sum=sum-a[i]*y[k-i];  /*sum-=sum_{k=0}^{m} a(i)*y(k-i);*/
		}
		y[k]=sum;
	}
}