信号处理.FIR滤波器在C程序中的实现



有人能告诉我如何使用c编程语言实现FIR滤波器吗。

设计FIR滤波器不是一个简单的话题,但实现一个已经设计好的滤波器(假设你已经有了FIR系数)也不算太糟。该算法称为卷积。这是一个天真的实现。。。

void convolve (double *p_coeffs, int p_coeffs_n,
               double *p_in, double *p_out, int n)
{
  int i, j, k;
  double tmp;
  for (k = 0; k < n; k++)  //  position in output
  {
    tmp = 0;
    for (i = 0; i < p_coeffs_n; i++)  //  position in coefficients array
    {
      j = k - i;  //  position in input
      if (j >= 0)  //  bounds check for input buffer
      {
        tmp += p_coeffs [i] * p_in [j];
      }
    }
    p_out [k] = tmp;
  }
}

基本上,卷积对输入信号进行移动加权平均。权重是滤波器系数,假设其总和为1.0。如果权重之和不是1.0,则会得到一些放大/衰减以及滤波。

顺便说一句,这个函数的系数数组可能是向后的——我还没有仔细检查过,我已经有一段时间没有考虑过这些事情了。

对于如何计算特定滤波器的FIR系数,这背后有相当多的数学知识——你真的需要一本关于数字信号处理的好书。这本书是免费的PDF版,但我不确定它有多好。我有Rorabaugh和Orfandis,都是在90年代中期出版的,但这些东西并没有真正过时。

组合多个过滤器:

从单位脉冲开始(第一个位置为1,其他位置为0的信号)。应用第一个过滤器。应用第二个过滤器。继续,直到应用所有筛选器。结果显示了组合滤波器如何卷积单位脉冲(前提是阵列足够长,没有数据丢失),因此其中的值是一个滤波器的系数,该滤波器是其他滤波器的组成。

这是示例代码:

#include <stdio.h>
#include <string.h>

#define NumberOf(a) (sizeof (a) / sizeof *(a))

/*  Convolve Signal with Filter.
    Signal must contain OutputLength + FilterLength - 1 elements.  Conversely,
    if there are N elements in Signal, OutputLength may be at most
    N+1-FilterLength.
*/
static void convolve(
    float *Signal,
    float *Filter, size_t FilterLength,
    float *Output, size_t OutputLength)
{
    for (size_t i = 0; i < OutputLength; ++i)
    {
        double sum = 0;
        for (size_t j = 0; j < FilterLength; ++j)
            sum += Signal[i+j] * Filter[FilterLength - 1 - j];
        Output[i] = sum;
    }
}

int main(void)
{
    //  Define a length for buffers that is long enough for this demonstration.
    #define LongEnough  128

    //  Define some sample filters.
    float Filter0[] = { 1, 2, -1 };
    float Filter1[] = { 1, 5, 7, 5, 1 };
    size_t Filter0Length = NumberOf(Filter0);
    size_t Filter1Length = NumberOf(Filter1);

    //  Define a unit impulse positioned so it captures all of the filters.
    size_t UnitImpulsePosition = Filter0Length - 1 + Filter1Length - 1;
    float UnitImpulse[LongEnough];
    memset(UnitImpulse, 0, sizeof UnitImpulse);
    UnitImpulse[UnitImpulsePosition] = 1;

    //  Calculate a filter that is Filter0 and Filter1 combined.
    float CombinedFilter[LongEnough];
    //  Set N to number of inputs that must be used.
    size_t N = UnitImpulsePosition + 1 + Filter0Length - 1 + Filter1Length - 1;
    //  Subtract to find number of outputs of first convolution, then convolve.
    N -= Filter0Length - 1;
    convolve(UnitImpulse,    Filter0, Filter0Length, CombinedFilter, N);
    //  Subtract to find number of outputs of second convolution, then convolve.
    N -= Filter1Length - 1;
    convolve(CombinedFilter, Filter1, Filter1Length, CombinedFilter, N);
    //  Remember size of resulting filter.
    size_t CombinedFilterLength = N;
    //  Display filter.
    for (size_t i = 0; i < CombinedFilterLength; ++i)
        printf("CombinedFilter[%zu] = %g.n", i, CombinedFilter[i]);

    //  Define two identical signals.
    float Buffer0[LongEnough];
    float Buffer1[LongEnough];
    for (size_t i = 0; i < LongEnough; ++i)
    {
        Buffer0[i] = i;
        Buffer1[i] = i;
    }

    //  Convolve Buffer0 by using the two filters separately.
    //  Start with buffer length.
    N = LongEnough;
    //  Subtract to find number of outputs of first convolution, then convolve.
    N -= Filter0Length - 1;
    convolve(Buffer0, Filter0, Filter0Length, Buffer0, N);
    //  Subtract to find number of outputs of second convolution, then convolve.
    N -= Filter1Length - 1;
    convolve(Buffer0, Filter1, Filter1Length, Buffer0, N);
    //  Remember the length of the result.
    size_t ResultLength = N;

    //  Convolve Buffer1 with the combined filter.
    convolve(Buffer1, CombinedFilter, CombinedFilterLength, Buffer1, ResultLength);

    //  Show the contents of Buffer0 and Buffer1, and their differences.
    for (size_t i = 0; i < ResultLength; ++i)
    {
        printf("Buffer0[%zu] = %g.  Buffer1[%zu] = %g.  Difference = %g.n",
            i, Buffer0[i], i, Buffer1[i], Buffer0[i] - Buffer1[i]);
    }
    return 0;
}

我发现这个代码片段不适合我(Visual Studio 2005)。

我最终发现卷积问题有一个很好的答案:

ANSI C代码中的1d线性卷积?

对于那些不知道的人来说——卷积与FIR滤波完全相同——"内核"是FIR滤波器的脉冲响应,信号是输入信号。

我希望这能帮助一些可怜的sap,他们正在寻找FIR代码:-)