最近在做直流电机的毕设中，由于需要采集转速，电流，电压，温度等参数，常规的采集容易受到干扰，所以特意复习了一下关于数字滤波有关的知识，并作出相应的整理。本文对数字滤波进行简单介绍，讲解七种常用的滤波算法并用C语言实现，并比较其优缺点。由于篇幅原因未能在嵌入式系统中实验，读者可以自行验证。本篇篇幅较长，建议收藏。
所用工具：
1、测试软件：Visual C++6.0
知识概括：
通过本篇文章您将学到：
1、数字滤波算法基础知识
很多嵌入式系统都需要通过AD转换方式采集模拟信号，当干扰作用于模拟信号之后，AD转换结果就会偏离真实值。如果仅采样一次，是无法确定该结果是否可信的，必须多次采样，得到一个AD转换的数据序列，通过某种处理后，才能得到一个可信度较高的结果。这种用软件算法从采样数据序列中提取逼近真值数据的方法，通常称为数字滤波。数字滤波有硬件滤波的功效，却节省了硬件投资。实现数字滤波功能的软件算法称为数字滤波算法，由于数字滤波算法的灵活性，其效果往往是硬件滤波电路达不到的。
数字滤波的不足之处就是需要消耗一定的CPU机时，在进行实时信号处理时必须充分考虑到这一点。当有用信号为高频信号，干扰信号为低频信号时，需要采用高通滤波，而实时高通滤波算法收到CPU速度的限制，往往力不从心，而硬件高通滤波器却很容易处理这种实时信号。当有用信号为低频信号（如温度，湿度，重量，水位等），其频率通常很低（甚至接近0.001Hz），而干扰信号频率先对较高（如各种电磁干扰），需要采用低通滤波。由于硬件低通滤波体积往往较大，与当前电子产品微型化趋势不相容。而当前CPU的速度在执行实时低通滤波算法时完全可以胜任，故数字滤波主要应用领域为实时低通滤波。在某些情况下，数字滤波也用来进行带通滤波。
对于非实时信息处理，由于不受处理速度的制约，数字滤波算法广泛应用于各种音频数据，图片数据和说数据的加工处理。
如果用硬件手段（FPGA/CPLD）来实现数字滤波，其处理速度将有质的飞越，滤波算法采用数学工具软件MATLAB来设计，滤波功能也将更加丰富多彩，应用频率范围也大大扩展。以下介绍在嵌入式系统中广泛使用的几种数字滤波算法，这些算法都是用软件进行实时低频信号处理，达到提高系统抗干扰能力的目的。
生活经验告诉我们，很多物理量的变化是需要一定时间的，相邻两次采样值之间的变化有一个合理的限度（即有用信号高频分量的最高频率有限）。例如在热处理车间的大型电炉里，工作的温度变化不可能在短时间（采样间隔）内发生剧烈变化，如果相邻两次采样值之间的变化未超过预定范围，说明该采样值未受到明显干扰，可以采用。如果相邻两次采样值之间的变化超过预定范围，说明该采样值受到了明显干扰（毛刺型突发干扰），不能采用。对于不可信得采样数据，必须输出一个合理的替代数据，以保证采样数据序列的连续性和完整性。在要求不高的场合，可以用前一个采样数据来替代这个受干扰的采样数据。
限幅滤波C语言代码如下：
#include "stdafx.h"
int adcData[10]={10,12,15,13,14,25,16,16,20,9}; //手动定义ADC采集的值
#define D 10 //相邻两次采样值之间最大允许变化值
int history; //上次采样值
int main(int argc, char* argv[])
{
int neo; //定义本次采样值
history=adcData[0]; //定义上次采样值初值
for(int i=0;i<10;i++) //循环访问数组adcData
{
neo=adcData[i]; //将数组里的元素作为本次采样值
if(i>0) //i=0时，无上次采样值，故判断
{
if((neo-history>D)||(history-neo>D)) //限幅判断
neo=history; //如果超过最大变化值，则新值变为上次采样值
history=neo; //本次采样值变为上次采样值，开始下一轮循环
}
printf("%d\n",neo); //输出滤波后的值
}
return 0;
}
//循环neo输出结果为：10,12,15,13,14,14,16,16,20,20
//可见第六位和第10位由于相差超过10，故被前一个值替换
当被检测物理量处于明显变化阶段时，相邻两个采样值之间的差别也比较明显，这时用前一个采样值来替代当前受干扰的采样值存在较大的误差。因此，可以利用被检测物理量的变化趋势来进行预测，选择一个更加合理的数据来替代当前受干扰的采样值，使得替代误差大大减小。
为了掌握被检测物理量当前变化趋势，至少需要两个以上的采样值历史记录。我们采用最简单的线性预测算法，只有两个历史采样值就可以了。线性预测算法的含义为“当前采样值的变化趋势与前一次采样值的变化趋势相同”，即采样值在短期内保持相同的上升趋势和下降趋势。由于代码类似，只是加一些限制条件，故这里省略。
限幅滤波的关键是合理确定采样周期和相邻两次采样值的最大允许变化值。采样周期可以通过检测精度（硬件是AD模块分辨率，软件是采样频率）和物理量最大变化速率（比如温度每秒钟对多上升多少摄氏度）来计算。另外，采样周期不能太短，必须大于干扰的持续时间，一面一次干扰造成两次以上的采样值不准确。
限幅滤波本质上是低通滤波，由于毛刺型突发干扰为高频干扰，故可以被很好地滤除。另一方面，限幅滤波采用限幅手段来达到滤波目的，对于幅度比较小的噪声干扰无能为力，即不能滤除随机干扰，因而只能用于对精度要求不是很高的场合。
限幅滤波有一个潜在的隐患，在连续两个以上采样周期受到强干扰后，系统可能会不稳定。因此，限幅滤波只能在基本没有干扰的场合下采样，他只能滤除极个别偶发的毛刺型干扰。在环境恶化，干扰频繁的情况下，可以采用中值滤波算法来处理。该算法描述是：连续进行奇数次采样，然后将采样得到的数据样本进行排序，取中间数据样本作为有效采样值。例如连续采样5次，得到5个采样数据样本（例如17、17、29、27、16），然后进行排序（16、17、17、27、29），最后取中间的（第三个，即17）采样数据样本作为有效采样值输出。由于受到干扰的采样值偏离有效采样值，排序后必然处于两端的位置，只要受到干扰的采样数据样本个数小于总采样数据样本数目的一半，就可以确保中值采样数据样本的有效性。
在采样中值滤波算法时，需要确定两个采样周期：一个是主采样周期（T0）和子采样周期（T1），每个主采样周期输出一个有效采样值，具体的时间长度由系统精度和有效信号的高频分量决定。每个子采样周期进行一次采样操作，它是由干扰频繁程度和连续采样次数N决定。
采用中值滤波算法必须满足以下条件：
被采样的物理量本身在连续N次子采样周期期间是基本稳定的，其变化量小于系统精度要求，可以忽略不计。
每次干扰的最长持续时间已知，相邻两次干扰的间隔时间虽然不定，但远大于一次干扰的持续时间。这样，可以使得连续采样的数据样本中被干扰的样本数目不会超过总样本数目N的一半，确保中值样本的有效性。
中值滤波C语言代码如下：
#include "stdafx.h"
int adcData[20]={10,12,15,13,14,25,16,16,20,9, //手动定义ADC采集的值
13,14,25,16,16,20,9,10,12,15};
#define N 5 //定义N个数取中值
int value_buff[N]; //定义中值计算缓冲区
int main(int argc, char* argv[])
{
int neo; //定义滤波计算完输出的值
char temp; //定义冒泡排序法临时变量
for(int major_col=0;major_col<20;major_col+=N) //定义主采样周期，每个包括N个子采样周期
{
for(int minor_col=0;minor_col<N;minor_col++) //采集每个子采样周期的值
{
value_buff[minor_col]=adcData[minor_col+major_col];
}
//冒泡排序法
for(int i=0;i<N-1;i++)
{
for(int j=0;j<N-i-1;j++)
{
if(value_buff[j]>value_buff[j+1])
{
temp=value_buff[j];
value_buff[j]=value_buff[j+1];
value_buff[j+1]=temp;
}
}
}
neo=value_buff[(N-1)/2]; //输出的值取缓冲区内的中值
printf("%d\n",neo);
}
return 0;
}
//循环neo输出结果为：13,16,16,12
//可见四个主采样周期的值都是子采样周期数据的中值
中值滤波本质上也是低通滤波，由于毛刺型突发干扰为高频干扰，故可以被很好地滤除。由于中值滤波进行N次采样才输出一次有效值（即采样输出比N：1），抗突发干扰能力比限幅滤波要提高很多，能够在干扰比较频繁的场合正常工作。另外，在N个数据样本中，中值样本通常受到的噪声干扰也比较小，故中值滤波算法也具备一定的抗随机干扰能力（其能力与N正相关）。
中值滤波的效果比限幅滤波提高很多，但必须满足主采样周期远大于子采样周期的条件。而子采样周期受到毛刺型突发干扰持续时间的限制不能太短，故主采样周期必然较长，即被检测物理量的基波频率一定比较低（慢变信号），其应用场合受到一定的限制。
由于随机干扰（噪声型干扰）随着数据样本的增加其统计平均值区域零，故对被检测物理量进行连续多次采样，然后求其算数平均值作为有效采样值，就可以达到抑制随机干扰的效果。连续采样次数越多，抑制随机干扰的效果越好。这种滤波算法就是算术平均滤波算法。
算术平均滤波C语言代码如下：
#include "stdafx.h"
float adcData[20]={10,12,15,13,14,25,16,16,20,9, //手动定义ADC采集的值
13,14,25,16,16,20,9,10,12,15};
#define N 5 //定义N个数取算术平均
float value_buff[N]; //定义算术平均计算缓冲区
int main(int argc, char* argv[])
{
float neo=0.0; //定义滤波计算完输出的值
for(int major_col=0;major_col<20;major_col+=N) //定义主采样周期，每个包括N个子采样周期
{
for(int minor_col=0;minor_col<N;minor_col++) //采集每个子采样周期的值
{
value_buff[minor_col]=adcData[minor_col+major_col];
neo+=value_buff[minor_col];
}
neo=neo/N; //输出算数平均
printf("%0.2f\n",neo);
}
return 0;
}
//循环neo输出结果为：12.80 19.76 20.75 17.35
//可见四个主采样周期的值都是子采样周期数据的算数平均
算术平均滤波对随机干扰的抑制效果与连续采样次数N密切相关，当采样次数N增加到一定程度后，被测物理量本身产生的变化量就会超过允许范围，再也不能忽略不计了，因此连续采样次数N是不能任意增加的。其应用场合与中值滤波的应用场合类似，必须满足主采样周期远大于子采样周期的条件。而子采样周期受到毛刺型突发干扰持续时间的限制不能太短，故主采样周期必然较长，即被检测物理量的基波频率一定比较低（慢变信号），其应用场合受到一定限制。
虽然算术平均滤波对随机干扰的抑制效果比较好，但却不能消除毛刺型突发干扰。只要有一个采样数据样本受到毛刺型突发干扰，算术平均值将明显偏离真实值，所以一般需要结合其他限制条件去对算术平均滤波算法进行优化。
算术平均滤波不能消除毛刺型突发干扰，只是通过平均操作将其影响削弱。因毛刺型突发干扰使采样值远离真实值，针对毛刺型突发干扰这一特点，可以比较容易地将其剔除，不参加平均值计算，从而使平均滤波的输出值更接近真实值。算法原理如下：连续采集N次，将其累加求和，同时找出其中的最大值和最小值，再从累加和中减去最大值和最小值，将剩余的N-2个采样值求平均，即得有效采样值。这种滤波算法就是去极值平均滤波算法。
算术平均滤波C语言代码如下：
#include "stdafx.h"
float adcData[20]={10,12,15,13,14,25,16,16,20,9, //手动定义ADC采集的值
13,14,25,16,16,20,9,10,12,15};
#define N 5 //定义N个数取算术平均
float value_buff[N]; //定义算术平均计算缓冲区
float nvalue_buff[N-2]; //定义去极值之后的缓冲区
int main(int argc, char* argv[])
{
float neo=0.0; //定义滤波计算完输出的值
float temp;
for(int major_col=0;major_col<20;major_col+=N) //定义主采样周期，每个包括N个子采样周期
{
for(int minor_col=0;minor_col<N;minor_col++) //采集每个子采样周期的值
{
value_buff[minor_col]=adcData[minor_col+major_col];
}
//冒泡排序法
for(int i=0;i<N-1;i++)
{
for(int j=0;j<N-i-1;j++)
{
if(value_buff[j]>value_buff[j+1])
{
temp=value_buff[j];
value_buff[j]=value_buff[j+1];
value_buff[j+1]=temp;
}
}
}
for(int k=0;k<N-2;k++) //去极值并且放入缓冲区
{
nvalue_buff[k]=value_buff[k+1];
neo+=nvalue_buff[k];
}
neo=neo/(N-2); //输出值为去极值之后的算数平均值
printf("%0.2f\n",neo);
}
return 0;
}
//循环neo输出结果为：13.00 21.67 22.56 19.85
//可见四个主采样周期的值都是子采样周期数据去极值后的算数平均
去极值平均滤波同时具有消除毛刺型突发干扰和噪声型随机干扰的能力，在低频信号采集系统中被广泛运用。比如在各种文艺表演大赛的评奖计分中经常听到主持人所说的“去掉一个最高分，去掉一个最低分，最终的得分是xxx”，这实际上就是采用了去极值平均滤波算法。
如果连续N次采样中没有受到毛刺型干扰的数据样本，被剔除的将是两个随机误差最大的数据样本，滤波效果仍然很好。去极值平均滤波算法也有效果不佳的时候，如果连续N次采样中有两个以上的正毛刺数据样本或者两个以上的负毛刺数据样本，将至少有一个毛刺数据样本不能剔除，最终平均值必然与真实值存在明显偏差。
去极值平均滤波的应用条件与算术平均滤波的应用条件相同，主要用来对低频信号进行软件滤波。
前面介绍的两种平均滤波算法有一个共同点，即每个主采样周期取得的一个有效采样值是由连续若干个子采样值计算出来的，为此，这些算法的主采样周期都是比较长的。当被检测物理量变化快时，较长的主采样周期不能及时采集到被检测物理量中的高频信息，致使系统反应迟钝，实时性难以保证。要保证系统的实时性，必须要按技术指标的要求相应缩短主采样周期，而子采样周期不能随意缩短（受干扰持续时间的制约），只好减少连续采样次数，连续采样次数的减少又直接降低了滤波效果。既然如此，干脆直接按技术指标的要求确定一个合理的采样周期，每个采样周期只进行一次采样操作（即采样输出比1:1），然后将当前的采样值和之前连续的若干次采样值（最近历史采样值）一起进行平均，将得到的平均值作为当前有效采样值投入使用。由于参与平均运算的历史采样值个数固定且内容不断更新，相当于一个滑动的时间窗口，故将这种平均滤波方式称为“滑动平均滤波”，窗口的大小用采样数据样本个数m来表示。按这种方式进行采样就不再存在主采样周期和子采样周期的问题。
滑动平均滤波C语言代码如下：
#include "stdafx.h"
float adcData[20]={10,12,15,13,14,25,16,16,20,9, //手动定义ADC采集的值
13,14,25,16,16,20,9,10,12,15};
#define N 5 //定义N个数取算术平均
float value_buff[N]; //定义算术平均计算缓冲区
int main(int argc, char* argv[])
{
float neo=0.0; //定义滤波计算完输出的值
for(int step=0;step<20-N;step++) //缓冲区的循环
{
for(int count=0;count<N;count++) //缓冲区内N个数的赋值
{
value_buff[count]=adcData[count+step];
neo+=value_buff[count];
}
neo=neo/N; //输出滑动平均后的值，一共20-N个
printf("%0.2f\n",neo);
}
return 0;
}
//循环neo输出结果为：12.80 18.36 20.27 20.85 22.37 21.67 19.13 18.23
// 19.85 19.37 20.67 22.33 21.67 18.53 17.11
//可见添加滑动窗之后会更加平滑
//理应输入输出是1:1，这里为了简便代码，故这样输出，只要知道思想就可以了，具体算法具体使用
在滑动平均滤波算法中，窗口中的m个采样数据样本以平等身份参与运算，这对抑制随机干扰比较有利，但将有超时严重的滞后效果，降低系统反应速度。为了提高系统的反应速度，滤波算法的输出应该及时反映当前采样值中包含的有效信息，即增加即时数据样本和近期数据样本的权重，降低早期数据样本的权重。为此，为滑动窗口中的各个数据样本分配不同的“加权系数”，进行加权平均运算。这种滤波算法称为“滑动加权滤波”。
滑动加权滤波C语言代码如下：
#include "stdafx.h"
float adcData[20]={10,12,15,13,14,25,16,16,20,9, //手动定义ADC采集的值
13,14,25,16,16,20,9,10,12,15};
#define N 5 //定义N个数取算术平均
float value_buff[N]; //定义算术平均计算缓冲区
float perc_buff[N]={0.1,0.2,0.1,0.3,0.3}; //定义缓冲区权重
int main(int argc, char* argv[])
{
for(int step=0;step<20-N;step++) //缓冲区的循环
{
float neo=0.0; //定义滤波计算完输出的值
for(int count=0;count<N;count++) //缓冲区内N个数的赋值
{
value_buff[count]=adcData[count+step];
value_buff[count]*=perc_buff[count]; //每个值乘对应的权重
neo+=value_buff[count]; //累加
}
printf("%0.2f\n",neo);
}
return 0;
}
//循环neo输出结果为：13.00 17.20 17.80 16.20 18.80 16.00 13.40 14.60
// 16.80 17.20 16.20 18.80 16.00 12.50 13.10
//可见添加滑动窗之后会更加平滑
将一阶滞后滤波器的微分方程用差分方程来表示，便可以用软件算法来模拟硬件一阶滞后滤波器（RC低通滤波器）的功能。
由上式可以看出，本次滤波的输出值主要取决于上次滤波的输出值，本次采样值对本次滤波输出的贡献是比较小的，但多少有些修正作用。这种算法便模拟了具有较大惯性的一阶滞后滤波器的功能。
当目标参数为变化很慢的物理量时（如大型储水池的水位信号），采用一阶滞后滤波算法是很有效的。
另一方面，他不能滤除高于二分之一采样频率的干扰信号。一阶滞后滤波算法程序流程与加权平均滤波相似，而加权系数只有两个：a和（1-a）。
一阶滞后滤波C语言代码如下：
#include "stdafx.h"
float adcData[10]={10,12,15,13,14,25,16,16,20,9}; //手动定义ADC采集的值
#define D 10 //相邻两次采样值之间最大允许变化值
#define A 0.85 //滤波系数
float history=0.0; //上次采样值
int main(int argc, char* argv[])
{
float neo=0.0; //定义本次采样值
history=adcData[0]; //定义上次采样值初值
for(int i=0;i<10;i++) //循环访问数组adcData
{
neo=adcData[i]; //将数组里的元素作为本次采样值
neo=A*neo+(1-A)*history; //滤波加权
printf("%0.2f\n",neo); //输出滤波后的值
}
return 0;
}
//循环neo输出结果为：10.00 11.70 14.25 12.55 13.40 22.75 15.10 15.10 18.50 9.15
//可见滤波系数的取值对滤波效果尤为重要
在设置滤波系数时需要同时兼顾系统的平稳性和灵敏性，而两者在一阶滞后滤波算法中属于不能同事达到最优，只能取得相对最优解，所以找到一个平衡点很重要，这也需要具体工程具体应用。当数据变化较为快速时需要以灵敏性优先考虑，当数据趋于稳定时需要以平稳性优先考虑。
下面将从采样输出比，康毛刺干扰能力和抗噪声干扰能力三个方面进行对比，如下表所示：
每种数字滤波算法都有自己的特点和应用条件，应该根据实际应用场合来合理选择数字滤波算法，并根据技术参数要求来确定采样周期。以上这些也只算是比较基础的滤波方式，可以相互结合使用，但也要考虑采集数据的类型。读者掌握这些之后可以在学习一些更为升入的滤波方式，比如卡尔曼，维纳等，技多不压身。
PS：最最最重要的，码字不易，求各位大佬看完点个赞，感谢！！！
