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,进行观察:
文章来源:https://www.toymoban.com/news/detail-819493.html
可以看到,频率为4Hz的信号幅度衰减超过一半,而1Hz的衰减非常少。说明效果还是不错的。文章来源地址https://www.toymoban.com/news/detail-819493.html
到了这里,关于一阶低通滤波的C语言实现(简单易移植)的文章就介绍完了。如果您还想了解更多内容,请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章,希望大家以后多多支持TOY模板网!