irpas技术客

一阶低通滤波的C语言实现(简单易移植)_月落三千雪_c语言实现低通滤波器

未知 4282

一阶低通滤波的C语言实现 0 引言1 一阶低通滤波器模型2 matlab 实现2.1 matlab 代码2.2 总结 3 c语言实现4 matlab 查看波形频率(快速傅里叶变换,FFT)

0 引言

一阶低通滤波器(Low Pass Filter,LPF),核心参数为截止频率fc,该算法可以保留截止频率以内的信号,而衰减截止频率之外的信号。主要用于去除高频噪声。

1 一阶低通滤波器模型

一阶低通滤波公式如下:

也可以写作:

其中:

参数说明:y(n)为本次滤波输出值,y(n-1)为上次滤波输出值,x(n)为本次采样值。Ts为采样周期,fc为截止频率。α范围为[0,1]

2 matlab 实现

我们假设,现在有一个信号,它包含了频率为1Hz(幅值为3)和4Hz(幅值为1)的两个正弦波。

如图所示:黄色为叠加后的信号。 我们假定采样周期为0.02s,且需要保留1Hz的信号(那么截止频率应该在1Hz~4Hz之间,先取1Hz)。

2.1 matlab 代码 t=(0:0.02:5); x=3*sin(2*pi*1*t)+1*sin(2*pi*4*t);%采集的混合信号 plot(t,x); hold on fc=2;%截止频率 Ts=0.02;%采样周期 b=2*pi*fc*Ts; alpha=b/(b+1);%alpha y=zeros(1,251); y(1)=x(1) for n= 2:251 y(n)=y(n-1)+alpha*(x(n)-y(n-1));%一阶滤波 end plot(t,y);

我们将滤波后的波形和理论的1Hz的波形对比:

你会发现:

确实有效果波形幅度变小了波形滞后了

因此,一阶低通滤波,存在滞后的现象。那么我们可以适当修改一下α,比如,将截止频率设置为2Hz,再看看效果。

你会发现:

幅值变化不大相位滞后不大 2.2 总结

对于一阶低通滤波器,截止频率应在有效信号频率和杂波频率之间。 越接近有效信号频率,则滤波效果越强,但存在相位更加滞后,幅度衰减越厉害的情况。 因此,应该多次设值,从而寻找一个能够兼顾响应速度与滤波效果的值。

3 c语言实现

虽然在网上看到别人写的c语言实现,但基本是原理说明,搬过来没办法直接使用,还要修修改改的,就很难受。 我这里模拟采集信号的过程,写一个c语言的实现。主要添加了一个模拟输入的效果,数据其实就是从matlab里面复制出来的:x=3*sin(2*pi*1*t)+1*sin(2*pi*4*t)。

#include <stdio.h> #include <string.h> /*************************** 模拟输入数据 ********************************/ float input_data[251] = { 0.00, 0.86, 1.59, 2.10, 2.35, 2.35, 2.18, 1.94, 1.76, 1.73, 1.90, 2.26, 2.75, 3.24, 3.63, 3.80, 3.70, 3.30, 2.68, 1.93, 1.18, 0.54, 0.11, -0.10, -0.11, 0.00, 0.11, 0.10, -0.11, -0.54, -1.18, -1.93, -2.68, -3.30, -3.70, -3.80, -3.63, -3.24, -2.75, -2.26, -1.90, -1.73, -1.76, -1.94, -2.18, -2.35, -2.35, -2.10, -1.59, -0.86, 0.00, 0.86, 1.59, 2.10, 2.35, 2.35, 2.18, 1.94, 1.76, 1.73, 1.90, 2.26, 2.75, 3.24, 3.63, 3.80, 3.70, 3.30, 2.68, 1.93, 1.18, 0.54, 0.11, -0.10, -0.11, 0.00, 0.11, 0.10, -0.11, -0.54, -1.18, -1.93, -2.68, -3.30, -3.70, -3.80, -3.63, -3.24, -2.75, -2.26, -1.90, -1.73, -1.76, -1.94, -2.18, -2.35, -2.35, -2.10, -1.59, -0.86, 0.00, 0.86, 1.59, 2.10, 2.35, 2.35, 2.18, 1.94, 1.76, 1.73, 1.90, 2.26, 2.75, 3.24, 3.63, 3.80, 3.70, 3.30, 2.68, 1.93, 1.18, 0.54, 0.11, -0.10, -0.11, 0.00, 0.11, 0.10, -0.11, -0.54, -1.18, -1.93, -2.68, -3.30, -3.70, -3.80, -3.63, -3.24, -2.75, -2.26, -1.90, -1.73, -1.76, -1.94, -2.18, -2.35, -2.35, -2.10, -1.59, -0.86, 0.00, 0.86, 1.59, 2.10, 2.35, 2.35, 2.18, 1.94, 1.76, 1.73, 1.90, 2.26, 2.75, 3.24, 3.63, 3.80, 3.70, 3.30, 2.68, 1.93, 1.18, 0.54, 0.11, -0.10, -0.11, 0.00, 0.11, 0.10, -0.11, -0.54, -1.18, -1.93, -2.68, -3.30, -3.70, -3.80, -3.63, -3.24, -2.75, -2.26, -1.90, -1.73, -1.76, -1.94, -2.18, -2.35, -2.35, -2.10, -1.59, -0.86, 0.00, 0.86, 1.59, 2.10, 2.35, 2.35, 2.18, 1.94, 1.76, 1.73, 1.90, 2.26, 2.75, 3.24, 3.63, 3.80, 3.70, 3.30, 2.68, 1.93, 1.18, 0.54, 0.11, -0.10, -0.11, 0.00, 0.11, 0.10, -0.11, -0.54, -1.18, -1.93, -2.68, -3.30, -3.70, -3.80, -3.63, -3.24, -2.75, -2.26, -1.90, -1.73, -1.76, -1.94, -2.18, -2.35, -2.35, -2.10, -1.59, -0.86, 0.00}; float get_data(void) { static int i = 0; return (input_data[i++]); if (i == 251) //轮回 i = 0; } float fc = 2.0f; //截止频率 float Ts = 0.02f; //采样周期 float pi = 3.14159f; //π float alpha = 0; //滤波系数 /************************ 滤波器初始化 alpha *****************************/ void low_pass_filter_init(void) { float b = 2.0 * pi * fc * Ts; alpha = b / (b + 1); } float low_pass_filter(float value) { static float out_last = 0; //上一次滤波值 float out; /***************** 如果第一次进入,则给 out_last 赋值 ******************/ static char fisrt_flag = 1; if (fisrt_flag == 1) { fisrt_flag = 0; out_last = value; } /*************************** 一阶滤波 *********************************/ out = out_last + alpha * (value - out_last); out_last = out; return out; } void main(void) { float result[251]; low_pass_filter_init(); for (int i = 0; i < 251; i++) { result[i] = low_pass_filter(get_data()); printf("%f,", result[i]); } }

我这里使用的是VSCode,将终端输出的数据直接复制到matlab加点代码就可以观察了。(excel也行,但插入数据后会坨在一格里,要点击分列)。 效果如下:

4 matlab 查看波形频率(快速傅里叶变换,FFT)

有时候,我们想先分析我们的数据的频率组成,或者观察滤波效果,我们可以观察他们的频率。这里用的是快速傅里叶变换。 这里分析X=3*sin(2*pi*1*t)+1*sin(2*pi*4*t),代码如下:

t=(0:0.02:5);% 时间向量 X=3*sin(2*pi*1*t)+1*sin(2*pi*4*t);%采集的混合信号 Ts = 0.02; % 采样周期 Fs = 1/Ts; % 采样频率 L = length(t); % 获取信号长度 n = 2^nextpow2(L);%首先从原始信号长度确定下一个 2 次幂的新输入长度。这将用尾随零填充信号 X 以改善 fft 的性能。 %%%%%%%%%%%%%%%%%%%%%%% 快速傅里叶变换 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Y = fft(X,n);%信号的傅里叶变换 P2 = abs(Y/n);%计算双侧频谱 P2 P1 = P2(1:n/2+1);%基于 P2 和偶数信号长度 L 计算单侧频谱 P1。 P1(2:end-1) = 2*P1(2:end-1);%双侧边变成单侧侧边,幅度叠加 f = Fs*(0:(n/2))/n;%定义频域 f plot(f,P1) title('Single-Sided Amplitude Spectrum of X(t)') xlabel('f (Hz)') ylabel('|P1(f)|')

我们将之前代码的滤波输出(fc=2)直接赋值给X,进行观察:

可以看到,频率为4Hz的信号幅度衰减超过一半,而1Hz的衰减非常少。说明效果还是不错的。


1.本站遵循行业规范,任何转载的稿件都会明确标注作者和来源;2.本站的原创文章,会注明原创字样,如未注明都非原创,如有侵权请联系删除!;3.作者投稿可能会经我们编辑修改或补充;4.本站不提供任何储存功能只提供收集或者投稿人的网盘链接。

标签: #c语言实现低通滤波器