C/C++编程-算法学习-数字滤波器

发布于:2024-07-27 ⋅ 阅读:(30) ⋅ 点赞:(0)

一阶低通滤波器

结论

其基本原理基于以下公式:
o u t p u t [ n ] = α ∗ i n p u t [ n ] + ( 1 − α ) ∗ o u t p u t [ n − 1 ] output[n] = α * input[n] + (1 - α) * output[n - 1] output[n]=αinput[n]+(1α)output[n1]

  • output[n] 是当前的输出值
  • input[n] 是当前的输入值
  • α是滤波系数,取值范围在 0 到 1 之间。α值越大,对输入的响应越迅速,但滤波效果相对较弱;α值越小,滤波效果越强,但对输入的响应越慢
  • output[n - 1] 是上一次的输出值
    例如,如果 alpha = 0.1,输入值在短时间内快速变化,由于 (1 - alpha) 的权重较大,上一次的输出值对当前输出值的影响较大,从而起到平滑和抑制高频变化的作用。

推导1

1. 基本公式推导

对应电路模型一阶RC滤波器
Y ( s ) / X ( s ) = H ( s ) = 1 r c s + 1 = 1 r ⋅ j w c + 1 = 1 τ s + 1 Y(s)/X(s)=H(s) = \frac{1}{rcs +1} = \frac{1}{r·jwc + 1} =\frac{1}{\tau s + 1} Y(s)/X(s)=H(s)=rcs+11=rjwc+11=τs+11
在自控中称为一阶惯性环节
转成时域方程
X ( s ) = Y ( s ) / H ( s ) = Y ( s ) ( τ s + 1 ) X(s) = Y(s)/H(s) = Y(s)(\tau s +1) X(s)=Y(s)/H(s)=Y(s)(τs+1)
x ( t ) = τ y ′ ( t ) + y ( t ) x(t) = \tau y'(t) + y(t) x(t)=τy(t)+y(t)
将导数拆开(使用一阶后向差分法,对上面微分方程进行离散化)
x ( t ) = τ ( y ( t ) − y ( t − T ) T ) + y ( t ) x(t) = \tau (\frac {y(t) - y(t-T)}{T}) + y(t) x(t)=τ(Ty(t)y(tT))+y(t)
整理成可递归迭代函数
y ( t ) = ( 1 − T T + τ ) ⋅ y ( t − T ) + T T + τ x ( t ) y(t) = (1-\frac {T}{T+\tau})·y(t-T) + \frac{T}{T+\tau}x(t) y(t)=(1T+τT)y(tT)+T+τTx(t)
a = T T + τ a = \frac{T}{T+\tau} a=T+τT 则可得一般表达式
y ( t ) = ( 1 − a ) y ( t − T ) + a x ( t ) y(t) = (1-a)y(t-T)+ax(t) y(t)=(1a)y(tT)+ax(t)

2. 截止频率 和 采样频率 推导

实现

#include <stdio.h>
#include <stdlib.h>

// 一阶低通滤波器函数
float lowPassFilter(float input, float prevOutput, float alpha) {
    return alpha * input + (1 - alpha) * prevOutput;
}

int main() {
    float input = 10.0;  // 输入值
    float prevOutput = 5.0;  // 上一次的输出值
    float alpha = 0.2;  // 滤波系数

    for(int i=0; i<20; i++)
    {
        float output = lowPassFilter(input, prevOutput, alpha);

        prevOutput = output;

        printf("filter current result: %f\n", output);
    }

    system("pause");
    return 0;
}

运行结果
在这里插入图片描述

二阶低通滤波器

实现1

#include <stdio.h>
// #include </lib/gcc/x86_64-linux-gnu/9/math.h>
#include <math.h>

// 二阶低通滤波器参数
#define SAMPLING_FREQ 1000  // 采样频率
#define CUTOFF_FREQ 100  // 截止频率

// 计算滤波器系数
void calculateFilterCoefficients(double *a, double *b) {
    double omega = 2 * M_PI * CUTOFF_FREQ / SAMPLING_FREQ;
    double alpha = sin(omega) / (2 * 0.707);
    double beta = cos(omega);

    double a0 = 1 + alpha;
    double a1 = -2 * beta;
    double a2 = 1 - alpha;

    double b0 = (1 - beta) / 2;
    double b1 = 1 - beta;
    double b2 = (1 - beta) / 2;

    *a = a0;
    *(a + 1) = a1;
    *(a + 2) = a2;

    *b = b0;
    *(b + 1) = b1;
    *(b + 2) = b2;
}

// 二阶低通滤波函数
double lowPassFilter(double input, double *prevInputs, double *prevOutputs, double *a, double *b) {
    double output = *b * input + *b * prevInputs[0] + *b * prevInputs[1] - *a * prevOutputs[0] - *a * prevOutputs[1];

    prevInputs[1] = prevInputs[0];
    prevInputs[0] = input;

    prevOutputs[1] = prevOutputs[0];
    prevOutputs[0] = output;

    return output;
}

int main() {
    double a[3], b[3];
    calculateFilterCoefficients(a, b);

    double prevInputs[2] = {0};
    double prevOutputs[2] = {0};

    double input = 10;  // 输入值,可根据实际情况修改
    double filteredOutput = lowPassFilter(input, prevInputs, prevOutputs, a, b);

    printf("滤波后的输出: %f\n", filteredOutput);

    return 0;
}



/**
如果你在使用gcc编译含数学函数的 C 程序时,出现undefined reference to 'sin'、undefined reference to 'cos'等错误,一般是由于缺少库造成的。因为在 Ubuntu 系统中,gcc的数学函数(如sin、cos等)是定义在libm.so里面的,而数学库不在默认路径下。通过添加-lm选项,就可以告诉编译器到正确的库中查找这些函数。

注意:在使用cmake进行编译时,需要添加命令target_link_libraries(your_target_name m)来链接数学库,其中your_target_name是你的目标名称。
*/

经我实际验证ubuntu20,的math库在如下路径
dpkg -l | grep math
在这里插入图片描述

在这里插入图片描述
在这里插入图片描述

实现2

#include <stdio.h>
#include <math.h>

// 二阶低通滤波器系数
typedef struct {
    double a0, a1, a2, b1, b2;
} FilterCoefficients;

// 计算二阶低通滤波器系数
void calculateFilterCoefficients(double cutoffFrequency, double samplingFrequency, FilterCoefficients *coefficients) {
    double omega = 2.0 * 3.14159 * cutoffFrequency / samplingFrequency;
    double cosOmega = cos(omega);
    double sinOmega = sin(omega);
    double alpha = sinOmega / (2.0 * 0.707);

    double a0 =  1 + alpha;
    double a1 = -2 * cosOmega;
    double a2 =  1 - alpha;
    double b1 = -2 * cosOmega;
    double b2 =  1 - alpha;

    coefficients->a0 = 1.0 / a0;
    coefficients->a1 = a1 / a0;
    coefficients->a2 = a2 / a0;
    coefficients->b1 = b1 / a0;
    coefficients->b2 = b2 / a0;
}

// 二阶低通滤波器函数
void secondOrderLowPassFilter(double input[], double output[], int length, FilterCoefficients coefficients) {
    output[0] = input[0];
    output[1] = coefficients.a0 * input[1] + coefficients.a1 * input[0] + coefficients.b1 * output[0];

    for (int i = 2; i < length; i++) {
        output[i] = coefficients.a0 * input[i] + coefficients.a1 * input[i - 1] + coefficients.a2 * input[i - 2]
                    - coefficients.b1 * output[i - 1] - coefficients.b2 * output[i - 2];
    }
}

int main() {
    double cutoffFrequency = 10.0;  // 截止频率
    double samplingFrequency = 50.0;  // 采样频率
    FilterCoefficients coefficients;

    calculateFilterCoefficients(cutoffFrequency, samplingFrequency, &coefficients);

    double input[] = {1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0};
    double output[10];
    int length = 10;

    secondOrderLowPassFilter(input, output, length, coefficients);

    for (int i = 0; i < length; i++) {
        printf("Output[%d] = %f\n", i, output[i]);
    }

    return 0;
}

实验结果:

在这里插入图片描述


网站公告

今日签到

点亮在社区的每一天
去签到