matlab filter 函数 C语言实现

C语言实现：

#include <stdio.h>
#include <math.h>
#include <memory.h>
#include <stdlib.h>
#define EPS 0.000001


//filter函数
void filter(const double* x, double* y, int xlen, double* a, double* b, int nfilt, double* zi)//nfilt为系数数组长度
{
    double tmp;
    int i,j;


    //normalization
    if( (*a-1.0>EPS) || (*a-1.0<-EPS) )
{
        tmp=*a;
        for(i=0;i<nfilt;i++)
{
            b[i]/=tmp;
            a[i]/=tmp;
        }
}


    memset(y,0,xlen*sizeof(double));//将y清零，以双浮点为单位


    a[0]=0.0;
    for(i=0;i<xlen;i++)
{
        for(j=0;i>=j&&j<nfilt;j++)
{
            y[i] += (b[j]*x[i-j]-a[j]*y[i-j]);
        }
        if(zi&&i<nfilt-1) y[i] += zi[i];
    }
    a[0]=1.0;

}

对应的是matlab nRet = filter(a,b,x); 

c++ 的 a数组 和 b数组 需要颠倒位置