参考文章: https://www.sohu.com/a/233360876_466960
参考文章:https://blog.youkuaiyun.com/qq_34040067/article/details/114048308
我是摘抄自网络, 对其进行了优化排版, 使其易于理解.
4.接下来是低通滤波公式推导(高通滤波公式可以用相同的思路推导出来)
PS:非专业分析,仅供参考,如有错误请指正。
先定义符号意义:
Q -电荷量(库伦),
C -电容(F),
R -电阻(欧姆),
U -电压(V),
fs -采样率(sampling rate, 次/秒)
fc -截止频率(Hz),
PI -圆周率(3.14…),
t -时间(秒)
这里额外插一句, 先定义符号的意义是非常有助于理解接下来的数学公式的. 很多的文章写的公式一大堆, 符号意义都不解释,
这样的文章就是垃圾.建议以后各位一定要在文章中写公式中符号的意义和单位.
RC高/低通滤波电路
RC带通滤波电路
根据RC低通电路图可知,输出电压等于电容电压
如果采样电压与电容电压不一致,电容就会充/放电,
假设电容充电后电压增加Ur(放电则Ur为负数):
符号:Uo-输出电压,U-电容电压, Ui-输入电压
式1
Uo=U+Ur
因为 C=Q/U, 所以U=Q/C
根据 U=Q/C, Q=I*t,I=U/R, t = 1/fs 得:
式2
Ur=(Ui-U) / (R*C*fs)
代入1得:
式3
Uo= U+ (Ui-U) / (R*C*fs)
由于我们是软件滤波,不需要真正的电容电阻, 所以令k=1/(RCfs),代入得:
Uo=U+(Ui-U)*k
简化:
式4 Uo=k*Ui+(1-k)*U
到这里已经得出滤波公式,如何通过截止频率得出k值呢
根据 截止频率公式
fc=1/(2*Pi*R*C)
得:
RC=1/(2*Pi*fc)
代入k=1 / (R * C * fs)得:
式5 k=(2*Pi*fc)/fs
C#版完整的,低通,高通滤波代码
using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
namespace Filter
{
/// <summary>
/// 模仿RC电路滤波的数字滤波算法
/// 参考文章 https://blog.youkuaiyun.com/qq_41848097/article/details/106124194
/// </summary>
///
public class RCFilter
{
float lastValue;
float fc;
float fs;
float k;
public RCFilter(float fc=200, float fs=100000 )
{
this.fc = fc;
this.fs = fs;
this.k = (float)(fc * 2 * Math.PI / fs);
}
/// <summary>
/// RC低通滤波
/// Y(n)= a* X(n) + (1-a) * Y(n-1);
/// a = (fc * 2π) / fl
/// fc :截止频率
/// fl :采样频率
/// </summary>
/// <param name="data">数据源</param>
/// <param name="fc">截止频率</param>
/// <param name="fs">采样频率</param>
public float 低通滤波(ref float data )
{
//float k = (float)(fc * 2 * Math.PI / fs); //滤波系数
lastValue = this.k * data + (1 - this.k) * lastValue;
return lastValue;
}
/// <summary>
/// RC高通滤波
/// Y(n)= a* X(n) + (1-a) * Y(n-1);
/// a = (fc * 2π) / fl
/// fc :截止频率
/// fl :采样频率
/// </summary>
/// <param name="data">数据源</param>
/// <param name="fc">截止频率</param>
/// <param name="fs">采样频率</param>
public float 高通滤波(ref float data)
{
//float k = (float)(fc * 2 * Math.PI / fs); //滤波系数
lastValue = this.k * data + (1 - this.k) * lastValue;
//return data - lastValue;//如果直接返回滤波结果,滤波后图像是倒转的,在心电图等一些场合,需要将图像再镜像过来
return -( data - lastValue);//滤波结果
}
/// <summary>
/// RC低通滤波
/// Y(n)= a* X(n) + (1-a) * Y(n-1);
/// a = (fc * 2π) / fl
/// fc :截止频率
/// fl :采样频率
/// </summary>
/// <param name="data">数据源</param>
/// <param name="fc">截止频率</param>
/// <param name="fs">采样频率</param>
public static double[] 低通滤波(ref double[] data, double fc, double fs)
{
//————————————————
//版权声明:本文为优快云博主「C# 学习者」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
//原文链接:https://blog.youkuaiyun.com/qq_41848097/article/details/106124194
double k = fc * 2 * Math.PI / fs; //滤波系数
double[] outdata = new double[data.Length];
outdata[0] = data[0];
double lastValue = data[0];//必须取第一个, 取特殊值,0或1都不对.
for (int i = 1; i < data.Length; i++)
{
lastValue = k * data[i] + (1 - k) * lastValue;
outdata[i] = lastValue;
}
return data;
}
/// <summary>
/// RC低通滤波
/// Y(n)= a* X(n) + (1-a) * Y(n-1);
/// a = (fc * 2π) / fl
/// fc :截止频率
/// fl :采样频率
/// </summary>
/// <param name="data">数据源</param>
/// <param name="fc">截止频率</param>
/// <param name="fs">采样频率</param>
public static float[] 低通滤波 (ref float[] data, float fc, float fs)
{
//————————————————
//版权声明:本文为优快云博主「C# 学习者」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
//原文链接:https://blog.youkuaiyun.com/qq_41848097/article/details/106124194
if (data == null || data.Length == 0)
{
return data;
}
float k = (float)( fc * 2 * Math.PI / fs); //滤波系数
float[] outdata = new float[data.Length];
outdata[0] = data[0];
float lastValue = data[0];//必须取第一个, 取特殊值,0或1都不对.
for (int i = 1; i < data.Length; i++)
{
lastValue = k * data[i] + (1 - k) * lastValue;
outdata[i] = lastValue;
}
return data;
}
/// <summary>
/// RC高通滤波
/// Y(n)= a* X(n) + (1-a) * Y(n-1);
/// a = (fc * 2π) / fl
/// fc :截止频率
/// fl :采样频率
/// </summary>
/// <param name="data">数据源</param>
/// <param name="fc">截止频率</param>
/// <param name="fs">采样频率</param>
public static float[] 高通滤波(ref float[] data, float fc, float fs)
{
//————————————————
//版权声明:本文为优快云博主「C# 学习者」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
//原文链接:https://blog.youkuaiyun.com/qq_41848097/article/details/106124194
if (data == null || data.Length == 0)
{
return data;
}
float k = (float)(fc * 2 * Math.PI / fs); //滤波系数
float[] outdata = new float[data.Length];
outdata[0] = data[0];
float lastValue = data[0];//必须取第一个, 取特殊值,0或1都不对.
for (int i = 1; i < data.Length; i++)
{
lastValue = k * data[i] + (1 - k) * lastValue;
outdata[i] = data[i] - lastValue; // 高通这里是原始值减去上次的值, 低通是直接输出lastvalue
//outdata[i] = data[i] - lastValue;//如果直接返回滤波结果,滤波后图像是倒转的,在心电图等一些场合,需要将图像再镜像过来
outdata[i] = 0-( data[i] - lastValue); //滤波结果
}
return data;
}
/// <summary>
/// RC高通滤波
/// Y(n)= a* X(n) + (1-a) * Y(n-1);
/// a = (fc * 2π) / fl
/// fc :截止频率
/// fl :采样频率
/// </summary>
/// <param name="data">数据源</param>
/// <param name="fc">截止频率</param>
/// <param name="fs">采样频率</param>
public static double[] 高通滤波(ref double[] data, double fc, double fs)
{
//————————————————
//版权声明:本文为优快云博主「C# 学习者」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
//原文链接:https://blog.youkuaiyun.com/qq_41848097/article/details/106124194
if (data == null || data.Length == 0)
{
return data;
}
double k = (double)(fc * 2 * Math.PI / fs); //滤波系数
double[] outdata = new double[data.Length];
outdata[0] = data[0];
double lastValue = data[0];//必须取第一个, 取特殊值,0或1都不对.
for (int i = 1; i < data.Length; i++)
{
lastValue = k * data[i] + (1 - k) * lastValue;
outdata[i] = data[i] - lastValue; // 高通这里是原始值减去上次的值, 低通是直接输出lastvalue
//outdata[i] = data[i] - lastValue;//如果直接返回滤波结果,滤波后图像是倒转的,在心电图等一些场合,需要将图像再镜像过来
outdata[i] = 0 - (data[i] - lastValue); //滤波结果
}
return data;
}
}
}
C#版有两种用法, 一种是直接调用静态方法处理整个数组, 另外一种用法是new 一个新对象,
下面的用法是new 新对象的用法.
static RCFilter rcFilterA1;
static RCFilter rcFilterB1;
static RCFilter rcFilterA2;
static RCFilter rcFilterB2;
public static void 初始化滤波器(ref float fs)
{
rcFilterA1 = new RCFilter(GD.低通滤波截止频率, fs);//一阶滤波
rcFilterA2 = new RCFilter(GD.低通滤波截止频率, fs);//二阶滤波
rcFilterB1 = new RCFilter(GD.低通滤波截止频率, fs);
rcFilterB2 = new RCFilter(GD.低通滤波截止频率, fs);
}
public static void main(){
GD.CurrentCycle_A[i] = rcFilterA1.低通滤波(ref GD.CurrentCycle_A[i]);
GD.CurrentCycle_A[i] = rcFilterA2.低通滤波(ref GD.CurrentCycle_A[i]);
GD.CurrentCycle_B[i] = rcFilterB1.低通滤波(ref GD.CurrentCycle_B[i]);
GD.CurrentCycle_B[i] = rcFilterB2.低通滤波(ref GD.CurrentCycle_B[i]);
}
静态变量一次性使用的方法
float[] data = new float[]{1,2,3,3,43,43,.....};
RCFilter.低通滤波(ref data,100,50);//将data中高于50hz的信号都过滤掉
RCFilter.低通滤波(ref data,100,60);//将data中高于60hz的信号都过滤掉
RCFilter.高通滤波(ref data,100,10);//将data中低于10hz的信号都过滤掉
RCFilter.高通滤波(ref data,100,20);//将data中低于20hz的信号都过滤掉
如果即执行了高通又执行了低通, 就变成带通了.
注意RC滤波会引起相位后移, 如果需要计算相位的,请注意相位位移问题. 我这个里面执行的次数越多相位后移的也越多.