C++ Butterworth N阶滤波器设计

我正在寻找一个像Matlab函数一样计算Butterworth Nth滤波器设计系数的函数:

[bl,al]=butter(but_order,Ws);

和

[bh,ah]=butter(but_order,2*bandwidth(1)/fs,'high');

我发现了许多计算2阶而不是N阶的示例(例如，我使用18阶...)。 -不幸的是，我对DSP一无所知。

您知道任何库或轻松实现此方法的方法吗？当我知道订购时，请切断频率和采样率。我只需要得到B(分子)和A(分母)的 vector 。

还要求该方法在不同的平台(Windows，Linux，...

提前致谢。

最佳答案

可以很容易地找到它(在Debian或Ubuntu中):

$ aptitude search ~dbutterworth | grep lib

这会立即给您答案:

p   librtfilter-dev         - realtime digital filtering library (dev)
p   librtfilter1            - realtime digital filtering library
p   librtfilter1-dbg        - realtime digital filtering library (dbg)

因此，您需要一个名为rtfilter的库。描述:

该库是跨平台的，即可以在Linux，MacOS和Windows下运行。从
official site:

您可以这样安装:

$ sudo aptitude install librtfilter-dev librtfilter1

安装-dev软件包后，您甚至可以在ojit_code上找到一个示例(使用Butterworth过滤器)。这个示例(以及相应的/usr/share/doc/librtfilter1/examples/butterworth.c)也可以找到here。

特别是您对Makefile函数感兴趣。您可以通过以下命令访问此功能的文档:

$ man rtf_create_butterworth

或者您可以阅读here。

您可以指定任何过滤器顺序并将其作为rtf_create_butterworth()参数传递给num_pole函数(据我所知，极数与过滤器顺序相同)。

更新

该库未提供外部API 用于系数计算。它仅提供实际的过滤功能，因此您可以在过滤后使用rtf_create_butterworth()获取样本(数据)。

但是，您可以在库源代码中找到用于系数计算的代码。请参见compute_cheby_iir()函数。此函数是rtf_filter()，因此只能在库本身内部使用。但是，您可以将此功能代码和复制到您的项目中并使用它。另外，不要让此函数的名称使您感到困惑:它是用于Butterworth滤波器和Chebyshev滤波器系数计算的相同算法。

假设您已经为static函数准备了参数:

const double cutoff     = 8.0;          /* cutoff frequency, in Hz */
const double fs         = 512.0;        /* sampling rate, in Hz */

unsigned int nchann     = 1;            /* channels number */
int proctype            = RTF_FLOAT;    /* samples have float type */
double fc               = cutoff / fs;  /* normalized cut-off frequency, Hz */
unsigned int num_pole   = 2;            /* filter order */
int highpass            = 0;            /* lowpass filter */

现在，您要计算过滤器的分子和分母。我已经为您写了 wrapper :

struct coeff {
    double *num;
    double *den;
};

/* TODO: Insert compute_cheby_iir() function here, from library:
 * https://github.com/nbourdau/rtfilter/blob/master/src/common-filters.c#L250
 */

/* Calculate coefficients for Butterworth filter.
 * coeff: contains calculated coefficients
 * Returns 0 on success or negative value on failure.
 */
static int calc_coeff(unsigned int nchann, int proctype, double fc,
                      unsigned int num_pole, int highpass,
                      struct coeff *coeff)
{
    double *num = NULL, *den = NULL;
    double ripple = 0.0;
    int res = 0;

    if (num_pole % 2 != 0)
        return -1;

    num = calloc(num_pole+1, sizeof(*num));
    if (num == NULL)
        return -2;
    den = calloc(num_pole+1, sizeof(*den));
    if (den == NULL) {
        res = -3;
        goto err1;
    }

    /* Prepare the z-transform of the filter */
    if (!compute_cheby_iir(num, den, num_pole, highpass, ripple, fc)) {
        res = -4;
        goto err2;
    }

    coeff->num = num;
    coeff->den = den;

    return 0;

err2:
    free(den);
err1:
    free(num);
    return res;
}

您可以像下面这样使用该包装器:

int main(void)
{
    struct coeff coeff;
    int res;
    int i;

    /* Calculate coefficients */
    res = calc_coeff(nchann, proctype, fc, num_pole, highpass, &coeff);
    if (res != 0) {
        fprintf(stderr, "Error: unable to calculate coefficients: %d\n", res);
        return EXIT_FAILURE;
    }

    /* TODO: Work with calculated coefficients here (coeff.num, coeff.den) */
    for (i = 0; i < num_pole+1; ++i)
        printf("num[%d] = %f\n", i, coeff.num[i]);
    for (i = 0; i < num_pole+1; ++i)
        printf("den[%d] = %f\n", i, coeff.den[i]);

    /* Don't forget to free memory allocated in calc_coeff() */
    free(coeff.num);
    free(coeff.den);

    return EXIT_SUCCESS;
}

如果您对那些系数计算的数学背景感兴趣，请查看DSP Guide, chapter 33。

关于C++ Butterworth N阶滤波器设计，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/29053932/