【肌电信号】基于matlab带通滤波肌电信号处理【含Matlab源码 965期】

一、获取代码方式

获取代码方式1：
通过订阅紫极神光博客付费专栏，凭支付凭证，私信博主，可获得此代码。

获取代码方式2：
完整代码已上传我的资源：【肌电信号】基于matlab带通滤波肌电信号处理【含Matlab源码 965期】

备注：
订阅紫极神光博客付费专栏，可免费获得1份代码（有效期为订阅日起，三天内有效）；

二、带通滤波器简介

1 滤波器
滤波器按照频带划分可以分为：低通滤波器(LPF)、高通滤波器(HPF)、带通滤波器(BPF)、带阻滤波器(BEF)。其中射频天线领域主要采用带通滤波器(BPF)。

2 带通滤波器(BPF)

①：RF BPF：从天线中取出期望频带的有效信号，滤除不需要频带的电磁波信号和噪声。
②：RF BPF：用下级混频器(降频混频器)进行频率转换时，防止在IF频带中引入干扰信号，用BPF将不期望频率成分滤除。
③：IF BPF：从混频器进行频率转换所产生的许多频率成分中，用BPF只取出期望的IF信号，也具有限制接收信号频带的作用。
④：IF BPF：用混频器(升频转换器)进行频率转换时，用BPF防止在RF频带中引入干扰信号，也具有限制发射信号频带的作用。
⑤：RF BPF：从混频器进行频率转换所产生的许多频率成分中，用BPF只取出期望的RF信号。
⑥：RF BPF：用BPF除去功率放大器放大时失真所产生的寄生部分，以及大功率信号的输入在天线开关时所产生的寄生部分，以免寄生成分通过天线发射出去。

3 带通滤波器(BPF)选择
(1)对于低频电路，选择LC滤波器。LC滤波器通常用于基频领域。
(2)对于100MHz~1000MH频率范围，选择SAW滤波器。SAW滤波器通常用于IF BPF领域。
(3)对于2GHz频率电路，选择介质滤波器。介质滤波器通常用于RF BPF领域。

4 带通滤波器(BPF)详解
核心思想：BPF基本上是谐振电路。谐振电路分为串联谐振电路和并联谐振电路，对于串联谐振电路，若L和C无损耗，则在谐振频率时阻抗变为零；对于并联谐振电路，谐振频率时阻抗变为无穷大。对于串联谐振电路可知，L越大，衰减特性越陡；对于并联谐振电路可知，L越小C越大，衰减特性越陡。
谐振电路的谐振频率：

串联谐振与并联谐振的组合电路：

5 介质滤波器

6 微带线制作滤波器

三、部分源代码

clear all;
file_name = 'Data1.txt';    % 文件名
fid =fopen(file_name,'r');  %打开文件
h=fgets(fid);               %去掉第一行；
C=textscan(fid,'%*f%f%f%f%f');%去掉第一列；并将其他行列写入C
fclose(fid); %关闭文件
Data=[C{1}';C{2}';C{3}';C{4}'];%将C中的数据取出放入数组Data1
fid =fopen(['new_' file_name],'w'); %打开新文件夹
fprintf(fid,'%-.10e    %-.10e   %-.10e   %-.10e\r\n  ',Data); %将取出的数据放入新文件夹，纯数据
fclose(fid);%关闭新文件夹
%load 'new_Data1.txt';
t=['new_' file_name];
[newData]=textread(t,'','delimiter',' ');x1=newData(:,1);
x2=newData(:,2);
x3=newData(:,3);
x4=newData(:,4);figure(1);
subplot(4,1,1)
plot(x1);
axis([0 512 0.5 1.7]);grid on;
% figure;
subplot(4,1,2)
plot(x2);
axis([0 512 0.5 1.7]);grid on;
% figure;
subplot(4,1,3)
plot(x3);
axis([0 512 0.5 1.7]);grid on;
% figure;
subplot(4,1,4)
plot(x4);
axis([0 512 0.5 1.7]);grid on;
N=512;
M=1; %调节因子fs=1000;
%数据点数
%时间序列
n=0:N-1;
f=n*fs/N;    %频率序列
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%信号
figure(2);y1=fft(x1,N);    %对信号进行快速Fourier变换
mag1=abs(y1);     %求得Fourier变换后的振幅
subplot(4,1,1) ,plot(f(1:N/2),mag1(1:N/2));
axis([0 500 0 20]);grid on;y2=fft(x2,N);    %对信号进行快速Fourier变换
mag2=abs(y2);    %求得Fourier变换后的振幅
subplot(4,1,2), plot(f(1:N/2),mag2(1:N/2));
axis([0 500 0 20]);grid on;y3=fft(x3,N);    %对信号进行快速Fourier变换
mag3=abs(y3);     %求得Fourier变换后的振幅
subplot(4,1,3) ,plot(f(1:N/2),mag3(1:N/2));
axis([0 500 0 20]);grid on;y4=fft(x4,N);    %对信号进行快速Fourier变换
mag4=abs(y4);     %求得Fourier变换后的振幅
subplot(4,1,4), plot(f(1:N/2),mag4(1:N/2));
axis([0 500 0 20]);grid on;
%=====%part3%===========%带通滤波器%========%
Wp = [50 150]/1000; Ws = [10 500]/1000;
Rp = 1; Rs =50;
[n,Wn] = buttord(Wp,Ws,Rp,Rs)
[b1,a1] = butter(n,Wn);
figure(3)
freqz(b1,a1)
title('n=4 Butterworth Bandpass Filter');
%=====%part3%===========%带通滤波器%========%
%=====%part4%===========%50hz陷波器%=======%
Wp = [10 500]/1000; Ws = [49 51]/1000;
Rp = 1; Rs = 50;
[n,Wn] = buttord(Wp,Ws,Rp,Rs)
[b2,a2] = butter(n,Wn,'stop');
figure(4)
freqz(b2,a2)
title('n=2 Butterworth Bandstop Filter');
%=====%part4%===========%50hz陷波器%=======%
%=====%part5%=============%数字滤波%========%
x1=filter(b1,a1,x1);
x1=filter(b2,a2,x1);
x2=filter(b1,a1,x2);
x2=filter(b2,a2,x2);
x3=filter(b1,a1,x3);
x3=filter(b2,a2,x3);
x4=filter(b1,a1,x4);
x4=filter(b2,a2,x4);figure(5);
subplot(4,1,1)
plot(x1);grid on;
axis([0 512 -0.5 0.5])
subplot(4,1,2)
plot(x2);grid on;
axis([0 512 -0.5 0.5])
subplot(4,1,3)
plot(x3);grid on;
axis([0 512 -0.5 0.5])
subplot(4,1,4)
plot(x4);grid on;
axis([0 512 -0.5 0.5])%=====%part5%=============%数字滤波%========%
figure(6);
y1=fft(x1,N);    %对信号进行快速Fourier变换
mag1=abs(y1);     %求得Fourier变换后的振幅
subplot(4,1,1) ,plot(f(1:N/2),mag1(1:N/2));
axis([0 500 0 11]);
grid on;y2=fft(x2,N);    %对信号进行快速Fourier变换
mag2=abs(y2);    %求得Fourier变换后的振幅
subplot(4,1,2), plot(f(1:N/2),mag2(1:N/2));
axis([0 500 0 10]);
grid on;y3=fft(x3,N);    %对信号进行快速Fourier变换
mag3=abs(y3);     %求得Fourier变换后的振幅
subplot(4,1,3) ,plot(f(1:N/2),mag3(1:N/2));
axis([0 500 0 10]);
grid on;y4=fft(x4,N);    %对信号进行快速Fourier变换
mag4=abs(y4);     %求得Fourier变换后的振幅
subplot(4,1,4), plot(f(1:N/2),mag4(1:N/2));
axis([0 500 0 10]);
grid on;
%=====%part6%=============%验证效果%========%
iemg1=sum(abs(x1))/length(x1);
iemg2=sum(abs(x2))/length(x2);
iemg3=sum(abs(x3))/length(x3);
iemg4=sum(abs(x4))/length(x4);

四、运行结果

五、matlab版本及参考文献

1 matlab版本
2014a

2 参考文献
[1] 沈再阳.精通MATLAB信号处理[M].清华大学出版社，2015.
[2]高宝建,彭进业,王琳,潘建寿.信号与系统——使用MATLAB分析与实现[M].清华大学出版社，2020.
[3]王文光,魏少明,任欣.信号处理与系统分析的MATLAB实现[M].电子工业出版社，2018.