陷波器的离散化及仿真验证
一、陷波器在连续域的传递函数
1、最基本的陷波器传函
(1)
其中,wo是所谓“中心频率”,也就是你想要“陷掉”的频率。而 ζ 则是“陷阱”的宽度。
根据公式可以发现,当输入信号频率很小(s=0)或者很大( s=+∞)的时候,上面式子的值是1;当输入信号频率刚好等于 s=jωo的时候,分子是0,所以增益变成0,那这个频率的信号当然就全都被衰减掉了。
由上图可见,ζ越大,则弦波带宽越宽,但弦波频率处的衰减越小。
2、三参数陷波器传函
(2)
其中,ωo是陷波频率(即凹陷的中心频率),ζ1和ζ2是陷波系数
陷波滤波器重点关注的参数一般有三个:
(1)陷波频率(ωo rad/s可转换Hz)
(2)陷波深度(depth为衰减倍数)
例如对于100Hz频率处的衰减深度是100,那经过该滤波器后,幅值衰减100倍。
(3)陷波宽度(△f单位Hz)
即中心频率两侧,幅值衰减-3dB时,对应的两个频率的差值。
ζ1和ζ2与陷波深度depth和陷波宽度△f(Hz)的关系表示如下:
二、陷波器传函在MATLAB中的表达及其离散化
1、以最基本的陷波器传函为例
a、MATLAB中编写如下m文件:
syms w0 s Ts z zeta % 定义符号变量
G1 =(s^2+w0^2)/(s^2+2*w0*zeta*s+w0^2) %传递函数
sys_s2c = 2*(z-1)/Ts/(z+1);
G2 = subs(G1,s,sys_s2c) %离散化 tustin变换
G3 = collect(G2,z) % 将表达式G2中的以z为变量合并相同次幂;
得到连续域和Z域的传递函数如下:
G1 =(s^2 + w0^2)/(s^2 + 2*zeta*s*w0 + w0^2)
G2 =(w0^2 + (2*z - 2)^2/(Ts^2*(z + 1)^2))/(w0^2 + (2*z - 2)^2/(Ts^2*(z + 1)^2) + (2*w0*zeta*(2*z - 2))/(Ts*(z + 1)))
G3 =((Ts^2*w0^2 + 4)*z^2 + (2*Ts^2*w0^2 - 8)*z + Ts^2*w0^2 + 4)/((Ts^2*w0^2 + 4*zeta*Ts*w0 + 4)*z^2 + (2*Ts^2*w0^2 - 8)*z + Ts^2*w0^2 - 4*zeta*Ts*w0 + 4)
根据离散化的方法:
分子分母同时除以A0,得到差分方程系数
差分方程为:
y(k) = b0*x(k) + b1*x(k-1) +b2*x(k-2) - a1*y(k-1) - a2*y(k-2);
则
B0 = Ts^2*w0^2 + 4;
B1 = 2*Ts^2*w0^2 - 8;
B2 = B0; (3)
A0 = Ts^2*w0^2 + 4*zeta*Ts*w0 + 4;
A1 = B1;
A2 = Ts^2*w0^2 - 4*zeta*Ts*w0 + 4;
可得差分方程系数
a0 = 1
b0 = B0/A0
b1 = B1/A0 (4)
b2 = B2/A0
a1 = A1/A0
a2 = A2/A0
b、代入实际参数,求得差分方程系数
f0 = 100; %目标频率
w0 = 2*pi*f0;
Ts = 10e-6; %离散化周期
zeta = 0.5; %带宽% Control object funciton
gxnum2 = [1 0 w0^2]; % 传函分子
gxden3 = [1 2*w0*zeta w0^2]; % 传函分母
sys2 = tf(gxnum2, gxden3) % 传函
zpk(sys2, 's'); % 将传函用零极点格式表示
dsys2 = c2d(sys2, Ts, 'tustin') % s到z,即,连续到离散域的变换
zpk(dsys2, 'z'); % 将传函用零极点格式表示
[num2, den2] = tfdata(dsys2, 'v'); % 提取传递函数分子、分母的z次幂的系数,并保存到数组
得到传递函数如下:
需要注意的是,在上图Z域传递函数dsys2的各项系数是经四舍五入的,如果直接用这些系数到Simulink或编写程序进行仿真,得到的结果是错误的!!!需要使用更加精确的系数。这些系数可以将式3和式4代入MATLAB中直接求解得到。
b0 = 0.996868276853708 b1 = -1.993697199313698
b2 = 0.996868276853708 a1 = -1.993697199313698
a2 = 0.993736553707416
2、以三参数陷波器为例
三、仿真验证
1、使用经四舍五入的系数仿真
结果错误
2、使用正确的系数
结果正确:
四、C程序验证
1、方式一
直接将差分方程用代码呈现:
double ADValue;//Notching Filter Coefficient
#define Notching_filter_100Hz_a0 1
#define Notching_filter_100Hz_a1 -1.998704711158930
#define Notching_filter_100Hz_a2 0.998744164397945
#define Notching_filter_100Hz_b0 0.999372082198973
#define Notching_filter_100Hz_b1 -1.998704711158930
#define Notching_filter_100Hz_b2 0.999372082198973// Control law 2p2z data define
typedef struct IIR_2OR_DATA_TAG{double coeff_a0;double coeff_a1;double coeff_a2;double coeff_b0;double coeff_b1;double coeff_b2;double filter_out;double filter_y1;double filter_y2;double filter_u1;double filter_u2;
} IIR_2OR_DATA_DEF;IIR_2OR_DATA_DEF notching_data = {Notching_filter_100Hz_a0,Notching_filter_100Hz_a1, Notching_filter_100Hz_a2, Notching_filter_100Hz_b0, Notching_filter_100Hz_b1, Notching_filter_100Hz_b2, 0, 0, 0, 0 ,0
};double iir_2or_func(IIR_2OR_DATA_DEF *filter_date, double target)
{////y(out) = b0*x(k) + b1*x(k-1) +b2*x(k-2) - a1*y(k-1) - a2*y(k-2);//filter_date->filter_out = (filter_date->coeff_b0 * (target)) \+ (filter_date->coeff_b1 * filter_date->filter_u1) \+ (filter_date->coeff_b2 * filter_date->filter_u2) \- (filter_date->coeff_a1 * filter_date->filter_y1) \- (filter_date->coeff_a2 * filter_date->filter_y2);//// Update last data//filter_date->filter_y2 = filter_date->filter_y1;filter_date->filter_y1 = filter_date->filter_out;filter_date->filter_u2 = filter_date->filter_u1;filter_date->filter_u1 = target;//// Return Value//return(filter_date->filter_out);
}
得到结果正确:
2、方式二
参考文献3中的方法,增加一个中间变量w,来实现。
double ADValue;//Notching Filter Coefficient #define Notching_filter_100Hz_GAIN 1
#define Notching_filter_100Hz_A1 -1.998704711158930
#define Notching_filter_100Hz_A2 0.998744164397945
#define Notching_filter_100Hz_B0 0.999372082198973
#define Notching_filter_100Hz_B1 -1.998704711158930
#define Notching_filter_100Hz_B2 0.999372082198973// Control law 2p2z data define
typedef struct IIR_2OR_DATA_TAG{double coeff_GAIN;double coeff_B0;double coeff_B1;double coeff_B2;double coeff_A1;double coeff_A2;double filter_out;double filter_W0;double filter_W1;double filter_W2;
} IIR_2OR_DATA_DEF;IIR_2OR_DATA_DEF notching_data = {Notching_filter_100Hz_GAIN,Notching_filter_100Hz_B0, Notching_filter_100Hz_B1, Notching_filter_100Hz_B2, Notching_filter_100Hz_A1, Notching_filter_100Hz_A2, 0, 0, 0, 0
};double iir_2or_func(IIR_2OR_DATA_DEF *filter_date, double target)
{//// w0 = x(0) - A1 * W1 - A2 * W2//filter_date->filter_W0 = (target) - filter_date->coeff_A1 * filter_date->filter_W1 \- filter_date->coeff_A2 * filter_date->filter_W2;//// Y(0) = Gain * (B0 * W0 + B1 * W1 + B2 * W2)//filter_date->filter_out = filter_date->coeff_GAIN * ( filter_date->coeff_B0 * filter_date->filter_W0 \+ filter_date->coeff_B1 * filter_date->filter_W1 \+ filter_date->coeff_B2 * filter_date->filter_W2 );//// Update last data//filter_date->filter_W2 = filter_date->filter_W1;filter_date->filter_W1 = filter_date->filter_W0;//// Return Value//return(filter_date->filter_out);
}
仿真结果:
与方法一完全一样。
参考文献:
1、Simulink 窄带陷波滤波器(Notch filter)仿真到代码生成
2、温故知新(五)——三参数陷波滤波器离散化推导及MATLAB实现
3、陷波器及其算法(基于C语言)
PLECS仿真:notch_filter.plecs
陷波器的离散化及仿真验证相关推荐
- 50hz 60hz 级联 陷波器,心电信号50Hz陷波器的FPGA实现
收稿日期: 2009 - 03 - 23 心电信号 50 Hz 陷波器的 FPGA 实现 林 霖 , 张志德 (南方医科大学 生物医学工程学院 , 广东广州 510515) [中图分类号]TH772 ...
- 有源滤波器: 基于UAF42的50Hz陷波器仿真
上一小节,我们设计出一个基于UAF42的50Hz陷波器.在本 小节,我们将使用免费的仿真软件TINA对这个电路进行仿真分析.具体原理图如下所示. 其中3.1831M欧的电阻用两个E96标准的电阻串联组 ...
- 基础补充——什么是陷波器?作用是什么?
在伺服控制系统中,常因为机械安装等因素导致系统出现振动,甚至有时候在电气控制频率与机械谐振频率相等时,会出现发散性震荡,损坏设备,造成严重的经济损失.目前常用的解决措施是在系统中加入一个与谐振频率对应 ...
- 滤波器基础03——Sallen-Key滤波器、多反馈滤波器与Bainter陷波器
滤波器基础系列博客,传送门: 滤波器基础01--滤波器的种类与特性 滤波器基础02--滤波器的传递函数与性能参数 滤波器基础03--Sallen-Key滤波器.多反馈滤波器与Bainter陷波器 滤波 ...
- 心电图ECG常用滤波器之陷波器
1.介绍 什么是陷波器? 其实就是一直特殊的带通滤波器,只不过一般它的频带极窄,可以理解为就是滤除某种特定频率的滤波器. 为什么心电图ECG需要陷波器? 由于我们所处的环境,到处都有50/60Hz的交 ...
- 50hz 60hz 级联 陷波器_自适应陷波器级联神经网络抗干扰算法
随着卫星导航系统的快速发展,诸多军事及民用领域越来越依赖卫星导航系统.但是由于全球定位系统(GlobalPositioningSystem,GPS)信号频率和调制特征公开而且信噪比很低,GPS卫星信号 ...
- 锁相环(单相+陷波器)入门理解
新手数控工程师一枚,大概学习了一下锁相环,以下是我的简单总结. 锁相环工作原理 和电网相连接的一些设备,往往都需要对电网相角的精准估计才能实现更好的控制效果.比如双闭环的PFC,电压环能输出Iref, ...
- c语言实现陷波器算法,50Hz数字陷波器的设计.doc
50Hz数字陷波器的设计.doc 四川理工学院毕业设计(论文) 50Hz数字陷波器的设计 学 生:孙全成 学 号:04021030312 专 业:通信工程 班 级:2004.3 指导教师:徐永俊 四川 ...
- 使用数字陷波器滤除工频信号
在实际测量时经常会受到工频信号(交流50Hz)的干扰,有时干扰还很大,有用信号完全被淹没了.可以应用数字陷波器来消除工频信号的干扰. 数字陷波器函数如下: 函数:iirnotch 功能:数字陷波器设计 ...
最新文章
- 元气骑士机器人修好后怎么用_《元气骑士》五大“难度”挑战,从手速到恶搞很嗨,还能解锁皮肤...
- 位置式PID与增量式PID的介绍和代码实现
- android简单的自定义按钮,Android 自定义button简单示例
- 「C++」C++ Primer Plus 笔记:第三章 处理数据
- 如何在PFSense中设置故障转移和负载平衡
- mysql开发中遇到的坑_mysql优化过程中遇见的坑(mysql优化问题特别注意)
- Edge浏览器怎么打开开发者模式
- 使用SQL SSIS和BIML自动化Salesforce数据复制
- Spring Boot项目中使用RestTemplate时出现乱码时的解决方案
- java cc_6.JavaCC官方入门指南-例1
- python写水仙花数
- oracle中job的retry次数,Oracle 19c注意事项: DBMS_JOB 行为变化
- 2015Esri全球用户大会top10的QA
- 线性表表长是否要算入头结点
- manul_css.css:1 Failed to load resource: the server responded with a status of 404 ()
- 无效的月份oracle,Oracle插入失败:无效的月份
- 对于Transformer 模型----可以从哪些地方进行创新和改进
- 关于登录时验证码无法显示
- 手机怎样压缩图片大小?手机照片内存怎么缩小?
- ssssssssss
热门文章
- 计算机操作系统——LINUX的C语言编程与shell编程
- 笨方法学Python3复习
- UX术语详解:任务流,用户流,流程图以及其它全新术语
- 学区摇号软件设计_重学 Java 设计模式:实战观察者模式「模拟类似小客车指标摇号过程,监听消息通知用户中签场景」...
- Mysql 常用 时间函数
- Android:利用sdk中的build-tools对包进行签名
- 【Linux】usermod 命令的使用
- 内核spinlock raw_spin_lock spin_lock_bh
- word文档纯字数统计_如何在您的Word文档中插入字数统计
- 详解二叉树的前序遍历