当前位置: 首页 > news >正文

一阶数字高通滤波器

本文的主要内容包含一阶高通滤波器公式的推导和数字算法的实现以及编程和仿真

1 计算公式推导

1.1.2 算法实现及仿真

利用python实现的代码如下:

import numpy as np
# from scipy.signal import butter, lfilter, freqz
import matplotlib.pyplot as plt
#2pifWc = 2*np.pi*30Tsw = 0.00314           #采样时间
halfdigiW=np.tan(Wc/2*Tsw )b10=1/(halfdigiW+1)
b11=-b10
a10=(1-halfdigiW)/(halfdigiW+1)x=np.linspace(-np.pi,np.pi,2000)    #在[-pi,pi]区间上分割正2000个点  可以理解为信号采样时间为 2*pi/2000s
data=np.zeros_like(x)               #输入信号 保存被干扰的信号
data1=np.zeros_like(x)              #输入信号,保存未被干扰的信号,方便与滤波后的波形进行比较
y1=np.zeros_like(x)                 #一阶滤波输出
y2=np.zeros_like(x)                 #二阶滤波输出 陷波:对希望滤除的特定频率有很好的滤除作用
for i in range(len(x)):data[i] =np.sin( 2 * np.pi * x[i])+0.5*np.sin(30* 2 * np.pi * x[i])  #幅值为1频率为 1Hz的低频信号  + 幅值为0.5 频率为30hz的高频信号data1[i]=0.5*np.sin(130* 2 * np.pi * x[i])# y2[i] = b0*data[i]+b1*data[i-1]+a0*y2[i-1]y1[i] = b10 * data[i] + b11 * data[i - 1] + a10 * y1[i - 1]
y2[0] = 0
y2[1] = 0
y1[0] = 0
y1[1] = 0
#绘原始信号  + 滤波后的信号
plt.subplot(2, 1, 1)
plt.plot(x,data ,label='sig+noise')
plt.plot(x,y1,   'r',label='first order HP')#绘制理想信号  + 滤波后的信号
plt.subplot(2, 1, 2)
plt.plot(x,data1 ,label='sig')
plt.plot(x,y1,   'r',label='first order HP')
plt.grid()
plt.legend()
plt.show()

以下是在python中仿真的波形图:

图1-2 fc=0.5Hz

图1-3 fc=2Hz

图1-4 fc=5Hz

图1-5 fc=10Hz

图1-6 fc=20Hz

图1-7 fc=30Hz

图1-8 fc=40Hz

图1-9 fc=60Hz

图1-10 fc=159Hz

输入的信号是幅值为1 频率为 1Hz的低频信号  加上 幅值为0.5 频率为30Hz的高频信号,采样时间为0.003s,从图1-2到图1-10的仿真波形可以看出,当fc为0.5Hz时,滤后的波形有微小的衰减作用,但几乎和原波形一致,当逐渐增大截止频率fc,对低频的滤除结果越来越强,高频越来越接近高频本身的波形,当截止频率高于高频频率时,高频本身也会被滤掉。当截止频率大于等于1/2采样频率时,输出的是一条直线。

利用C语言实现的代码如下:

#ifndef _MHPF1W_F_H_
#define _MHPF1W_F_H_#include <stdint.h>
struct MHpf1W_F
{/*初始化*/struct{void (*Set)(struct MHpf1W_F *self, float cutFreq, int samFreq);  //设置截止和采样频率void (*VaryCutFreq)(struct MHpf1W_F *self, float cutFreq);       //改变截止频率float cutFreq;                                                   //截止频率float samFreq;                                                   //采样频率} Init;/*采样计算*/struct{int (*In)(struct MHpf1W_F *self, int Xn);int out_y;                                                        //输出值} Prd;/*变量 中间变量  系数等,由初始参数 初始化计算得出*/struct{float Ts;                                                         //采样周期int a0, b0, b1;                                                   //差分系数int Xn_1, Yn_1;} pri;
};
void MHpf1W_F_Create(struct MHpf1W_F *self);
#endif
//创建方式
// struct MHpf1W_F mlp;
// MHpf1W_F_Create(&mlp);
// mlp.Init.Set(&mlp,2, 1000);  fc=20Hz  fs=1000Hz#include "MHpf1W_F.h"
#include <string.h>
#include "math.h"static const float PI = 3.1415926535897932384626f;
#define MID(a,min,max) (a= (a<min)?min:(a<max)?a:max)
#define Q15M(a,b)  ((a*b)>>15)
/******************************************************************************** 函 数 名         : _Update* 函数功能         : 各系数计算,参数更新* 输入参数         : 滤波器对象 struct MHpf1W_F *self* 返 回 值         : 无*******************************************************************************/
static void _Update(struct MHpf1W_F *self)
{float halfdigiW,tgAnaWT ;halfdigiW = PI *  self->Init.cutFreq * self->pri.Ts;tgAnaWT = tan(halfdigiW);                       //ignore the 1/Tself->pri.b0 = 1/(tgAnaWT+1)*32768;             //转成Q15格式self->pri.b1 = -self->pri.b0;                   //转成Q15格式self->pri.a0 =((1-tgAnaWT)/(tgAnaWT+1))*32768;  //转成Q15格式self->pri.Xn_1 = 0;self->pri.Yn_1 = 0;
}
/******************************************************************************** 函 数 名         : InitSet* 函数功能         : 初始化* 输入参数         : cutFreq----截至频率*                   samFreq----计算机采样频率* 返 回 值         : 无*******************************************************************************/
static void InitSet(struct MHpf1W_F *self, float cutFreq, int samFreq)
{self->Init.cutFreq = MID(cutFreq  , 0.0f , samFreq*0.5f);;                       //截止频率self->Init.samFreq = samFreq;                       //采样频率self->pri.Ts = 1.0f / self->Init.samFreq;           //采样周期    1/Ts_Update(self);
}
/******************************************************************************** 函 数 名         : InitVaryCutF* 函数功能         : 改变截止频率* 输入参数         : cutFreq----截至频率* 返 回 值         : 无*******************************************************************************/
static void InitVaryCutF(struct MHpf1W_F *self, float cutFreq)
{self->Init.cutFreq = cutFreq;_Update(self);
}
/******************************************************************************** 函 数 名         : PrdIn* 函数功能         : 本次输出结果计算* 输入参数         : Xn----本次输入值* 返 回 值         :   本次滤波后的值* 计算公式         :Y(n)=b0*X(n)+b1*X(n-1)+a0*Y(n-1)*******************************************************************************/
static int PrdIn(struct MHpf1W_F *self, int Xn)
{/*Y(n)=b0*X(n)+b1*X(n-1)+b2*X(n-2)+a0*Y(n-1)+a2*Y(n-2)*/self->Prd.out_y = Q15M(self->pri.b0 , Xn ) + \Q15M(self->pri.b1 , self->pri.Xn_1 ) + \Q15M(self->pri.a0 , self->pri.Yn_1 ) ;self->pri.Yn_1 = self->Prd.out_y;self->pri.Xn_1 = Xn;return self->Prd.out_y;
}
/******************************************************************************** 函 数 名         : MHpf1W_F_Create* 函数功能         : 创建对象 初始化* 输入参数         : self对象* 返 回 值         : 无*******************************************************************************/
void MHpf1W_F_Create(struct MHpf1W_F *self)
{memset(self, 0, sizeof(struct MHpf1W_F));self->Init.Set = InitSet;self->Init.VaryCutFreq = InitVaryCutF;self->Prd.In = PrdIn;
}

单片机+匿名科创地面站的软件输出波形如下:

图1-11 fc=0.5Hz

图4-12 fc=2Hz

图1-13 fc=5Hz

图1-14 fc=10Hz

图1-15 fc=20Hz

图1-16 fc=30Hz

图1-17 fc=40Hz

图1-18 fc=60Hz

图1-18 fc=100Hz

图1-18 fc=200Hz

                                  图1-18 fc=500Hz

单片机模拟输入的信号是幅值为1000 频率为 1Hz的低频信号加上幅值为500 频率为30Hz的高频信号,采样时间为0.001s,从图1-11到图1-18的波形可以看出,输出变化的规律现象和python仿真的规律和现象是一致的,同样当截止频率大于等于1/2采样频率时,输出的是一条直线。

http://www.lryc.cn/news/351353.html

相关文章:

  • Linux多线程系列2: 模拟封装简易语言级线程库,线程互斥和锁,线程同步和条件变量,线程其他知识点
  • VUE3-form表单保存附件与基本信息
  • 无线网络安全技术基础
  • sheng的学习笔记-docker部署Greenplum
  • 【投稿资讯】区块链会议CCF A -- SP 2025 截止6.6、11.14 附录用率
  • C++哪些函数不能被声明为虚函数
  • vue中数据已经改变了,但是table里面内容没更新渲染!
  • 头歌实践教学平台:Junit实训入门篇
  • matlab使用教程(80)—修改图形对象的透明度
  • mysql bin 日志转成sql
  • 河南道路与桥梁乙级资质申请:注册证书与职称证书准备
  • 3D工业视觉
  • 使用auth_basic模块进行基础认证
  • 深度解析物联网平台:优化数据点位管理的实战策略
  • Spring常见问题
  • MiniMax Golang2轮面试,期望薪资25K
  • MyBatis系统学习篇 - MyBatis的缓存
  • K-means聚类模型
  • 免费分享一套微信小程序旅游推荐(智慧旅游)系统(SpringBoot后端+Vue管理端)【论文+源码+SQL脚本】,帅呆了~~
  • Matlab 2023b学习笔记1——界面认识
  • C++ sort排序的总和应用题
  • [力扣]——231.2的幂
  • 【css】引入背景图时候,路径写入@会报错
  • 【有手就行】使用你自己的声音做语音合成,CPU都能跑,亲测有效
  • 《ESP8266通信指南》番外-(附完整代码)ESP8266获取DHT11接入(基于Lua)
  • [IMX6ULL驱动开发]-Linux对中断的处理(一)
  • PHP基础学习笔记(面向对象OOP)
  • Mysql超详细安装配置教程(保姆级图文)
  • HR招聘测评,如何判断候选人的团队协作能力?
  • [STM32-HAL库]Flash库-HAL库-复杂数据读写-STM32CUBEMX开发-HAL库开发系列-主控STM32F103C6T6