C语言 C编程中的FIR滤波器实现

FIR filter implementation in C programming

提问by D X

Can anyone tell me how to implement an FIR filter using c programming language.

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

回答by Steve314

Designing an FIR filter is NOTa simple topic, but implementing an already-designed filter (assuming you already have the FIR coefficients) isn't too bad. The algorithm is called convolution. Here's a naive implementation...

设计一个 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 [k] * p_in [j];
      }
    }

    p_out [i] = tmp;
  }
}

Basically, the convolution does a moving weighted average of the input signal. The weights are the filter coefficients, which are assumed to sum to 1.0. If the weights sum to something other than 1.0, you get some amplification/attenuation as well as filtering.

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

BTW - it's possible this function has the coefficients array backwards - I haven't double-checked and it's a while since I thought about these things.

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

For how you calculate the FIR coefficients for a particular filter, there's a fair amount of mathematics behind that - you really need a good book on digital signal processing. This oneis available free for a PDF, but I'm not sure how good it is. I have Rorabaughand Orfandis, both published in the mid nineties but these things don't really get obsolete.

对于如何计算特定滤波器的 FIR 系数，背后有大量数学知识 - 您确实需要一本关于数字信号处理的好书。这个可以免费获得 PDF 版本，但我不确定它有多好。我有Rorabaugh和Orfandis，它们都在九十年代中期出版，但这些东西并没有真正过时。

回答by Eric Postpischil

To combine multiple filters:

组合多个过滤器：

Start with an unit impulse (a signal with a 1 in the first position and 0 everywhere else). Apply the first filter. Apply the second filter. Continue until all filters are applied. The result shows how the combined filters convolve the unit impulse (provided the array is long enough that no data was lost), so the values in it are the coefficients for one filter that is the composition of the other filters.

从单位脉冲（第一个位置为 1，其他位置为 0 的信号）开始。应用第一个过滤器。应用第二个过滤器。继续直到应用所有过滤器。结果显示了组合滤波器如何对单位脉冲进行卷积（假设数组足够长，不会丢失数据），因此其中的值是一个滤波器的系数，该滤波器是其他滤波器的组成。

Here is sample code:

这是示例代码：

#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;
}

回答by user3741627

I found this code snippet did not work for me (Visual Studio 2005).

我发现这个代码片段对我不起作用（Visual Studio 2005）。

I did eventually find convolution questions having a great answer:

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

1d linear convolution in ANSI C code?

ANSI C 代码中的一维线性卷积？

For those who don't know - the convolution is exactly the same operation as FIR filtering - the "kernel" is the FIR Filter impulse response and the signal is the input signal.

对于那些不知道的人 - 卷积与 FIR 滤波完全相同的操作 - “内核”是 FIR 滤波器脉冲响应，信号是输入信号。

I hope this helps some poor sap who was looking for FIR code :-)

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