clear
Fs=360;N=2048;
load('9.mat');
x1=y;
load('11.mat');
x2=y;
figure
subplot(2,1,1);plot(x1);grid;title('原始心电信号x1.9.mat');
subplot(2,1,2);plot(x2);title('原始心电信号x2.9.mat');grid
%(1)谱分析
mf1=fft(x1,N);%进行频谱变换(傅里叶变换)
mf2=fft(x2,N);
mag1=abs(mf1);
mag2=abs(mf2);
f1=(0:length(mf1)-1)*Fs/length(mf1); %进行频率变换
f2=(0:length(mf2)-1)*Fs/length(mf2);
figure
plot(f1,mag1);axis([0,370,0,190]);grid; %画出频谱图
xlabel('频率(HZ)');ylabel('幅值');title('心电信号x1幅度谱');
figure
plot(f2,mag2);axis([0,370,0,190]);grid; %画出频谱图
xlabel('频率(HZ)');ylabel('幅值');title('心电信号x2幅度图');
figure
plot(f1,angle(mf1)/pi);axis([0 360 -1.3 1.3]);grid;
xlabel('频率(HZ)');ylabel('相位');title('心电信号x1相位谱');
figure
plot(f1,angle(mf2)/pi);axis([0 360 -1.3 1.3]);grid;
xlabel('频率(HZ)');ylabel('相位');title('心电信号x2相位谱');
wn1=x1;
wn2=x2;
P1=10*log10(abs(fft(wn1).^2)/N);
P2=10*log10(abs(fft(wn2).^2)/N);
f1=(0:length(P1)-1)/length(P1);
f2=(0:length(P2)-1)/length(P2);
figure
plot(f1,P1);grid
xlabel('归一化频率');ylabel('功率(dB)');title('心电信号x1的功率谱');
figure
plot(f2,P2);grid
xlabel('归一化频率');ylabel('功率(dB)');title('心电信号x2的功率谱');
%(2)相关分析
avr1=mean(x1);avr2=mean(x2);
fprintf('信号x1的均值= %f\n',avr1);fprintf('信号x2的均值= %f\n',avr2);
var1=var(x1);var2=var(x2);
fprintf('信号x1的方差= %f\n',var1);fprintf('信号x2的方差= %f\n',var2);
rx1=xcorr(x1,'biased');
rx2=xcorr(x2,'biased');
rx1x2=xcorr(x1,x2);
figure
plot(rx1);grid;title('心电信号x1的自相关函数');
figure
plot(rx2);grid;title('心电信号x2的自相关函数');
figure
plot(rx1x2);grid;title('心电信号x1,x2的互相关函数');
%(3)数字滤波器设计
SNR=10*log(100/8); % 2%是能量比
x11=awgn(x1,SNR);
figure
subplot(2,1,1), plot(x1);title('原信号x1');% 加入噪声后有毛刺,但2%的噪声有点小,毛刺不明显。
subplot(2,1,2), plot(x11);title('加高斯白噪信号');
dt=1/1023.5; % 取樣時距(或週期),sec
t=(0:dt:2)'; % 建立一個0-2秒的時間向量
y_high=sin(2*pi*1000*t)/10; % 100 Hz的高頻訊號,振福1/10
y_out=x1+y_high; % 基頻載上一組高頻的輸出。
figure
plot(t,y_out);grid;title('加入高频噪声1000hz');
%——————IIR零相移数字滤波器纠正基线漂移——————-
Wp=1.4*2/Fs; %通带截止频率
Ws=0.6*2/Fs; %阻带截止频率
devel=0.005; %通带纹波
Rp=20*log10((1+devel)/(1-devel)); %通带纹波系数
Rs=20; %阻带衰减
[N, Wn]=ellipord(Wp,Ws,Rp,Rs,'s'); %求椭圆滤波器的阶次
[b, a]=ellip(N,Rp,Rs,Wn,'high'); %求椭圆滤波器的系数
[hw,w]=freqz(b,a,512);
result =filter(b,a,x1);
figure
freqz(b,a);title('IIR数字滤波器幅频曲线');
figure;subplot(2,1,1); plot(x1);
xlabel('t(s)');ylabel('幅值');title('原始信号');grid
subplot(2,1,2); plot(result);
xlabel('t(s)');ylabel('幅值');title('线性滤波后信号');grid
figure
N=512;
subplot(2,1,1);plot(abs(fft(x1))*2/N);
xlabel('频率(Hz)');ylabel('幅值');title('原始信号x1频谱');grid;
subplot(2,1,2);plot(abs(fft(result))*2/N);
xlabel('频率(Hz)');ylabel('幅值');title('线性滤波去掉基线漂移频谱');grid;
%(5)维纳滤波器去除工频干扰:
dt=1/1023.5; % 取樣時距(或週期),sec
t=(0:dt:2)'; % 建立一個0-2秒的時間向量
y=sin(2*pi*50*t)/10; % 100 Hz的高頻訊號,振福1/10
x_noise=x1+y; % 基頻載上一組高頻的輸出。
figure
plot(t,x_noise);grid;title('加入噪声50hz');
%维纳滤波
Mlag=100;
N=100;%维纳滤波器长度
Rxn=xcorr(x_noise,Mlag,'biased');
Rxnx=xcorr(x1,x_noise,Mlag,'biased');%产生输入信号与原始信号的互相关函数
rxnx=zeros(N,1);
rxnx(:)=Rxnx(101:101+N-1);
Rxx=zeros(N,N);%产生输入信号自相关矩阵
Rxx=diag(Rxn(101)*ones(1,N));
for i=2:N
c=Rxn(101+i)*ones(1,N+1-i);
Rxx=Rxx+diag(c,i-1)+diag(c,-i+1);
end
Rxx;
h=zeros(N,1);
h=inv(Rxx)*rxnx; %计算维纳滤波器的h(n)
yn=filter(h,1,x_noise); %将输入信号通过维纳滤波器
figure
plot(yn);title('经过维纳滤波器后信号');
ynmean=mean(yn) %计算经过维纳滤波器后信号均值
ynms=mean(yn.^2) %计算经过维纳滤波器后信号均方值
ynvar=var(yn,1) %计算经过维纳滤波器后信号方差
Ryn=xcorr(yn,Mlag,'biased'); % 计算经过维纳滤波器后信号自相关函数
Y=fft(yn); %计算经过维纳滤波器后信号序列的快速离散傅里叶变换
Y1=fft(x_noise);
Py=Y.*conj(Y)/600; % 计算信号频谱
Py1=Y.*conj(Y)/600;
figure
subplot(211)
semilogy(t,Py) %绘制在半对数坐标系下频谱图像
title('经过维纳滤波器后信号在半对数坐标系下频谱图像')
xlabel('频率(HZ)','color','b');ylabel('幅度','color','b')
subplot(212)
semilogy(t,Py1) %绘制在半对数坐标系下频谱图像
title('噪声信号在半对数坐标系下频谱图像')
xlabel('频率(HZ)','color','b');ylabel('幅度','color','b')
pyn=periodogram(yn);%计算经过维纳滤波器后信号的功率谱密度
pyn1=periodogram(x_noise);
figure
subplot(211)
semilogy(pyn)%绘制在半对数坐标系下功率谱密度图像
title('经过维纳滤波器后信号在半对数坐标系下功率谱密度图像');
xlabel('频率','color','b');ylabel('幅度','color','b')
subplot(212)
semilogy(pyn);title('噪声信号在半对数坐标系下功率谱密度图像');
xlabel('频率','color','b');ylabel('幅度','color','b')
没有合适的资源?快使用搜索试试~ 我知道了~
心电信号基于matlab心电信号特征提取+分析处理【含Matlab源码 289期】.zip
共25个文件
jpg:19个
mat:4个
doc:1个
1.该资源内容由用户上传,如若侵权请联系客服进行举报
2.虚拟产品一经售出概不退款(资源遇到问题,请及时私信上传者)
2.虚拟产品一经售出概不退款(资源遇到问题,请及时私信上传者)
版权申诉
5星 · 超过95%的资源 13 下载量 165 浏览量
2021-12-31
12:15:18
上传
评论 30
收藏 1.24MB ZIP 举报
温馨提示
CSDN海神之光上传的代码均可运行,亲测可用,直接替换数据即可,适合小白; 1、代码压缩包内容 主函数:main.m; 调用函数:其他m文件;无需运行 运行结果效果图; 2、代码运行版本 Matlab 2019b;若运行有误,根据提示修改;若不会,私信博主; 3、运行操作步骤 步骤一:将所有文件放到Matlab的当前文件夹中; 步骤二:双击打开main.m文件; 步骤三:点击运行,等程序运行完得到结果; 4、仿真咨询 如需其他服务,可私信博主或扫描博客文章底部QQ名片; 4.1 博客或资源的完整代码提供 4.2 期刊或参考文献复现 4.3 Matlab程序定制 4.4 科研合作 功率谱估计: 故障诊断分析: 雷达通信:雷达LFM、MIMO、成像、定位、干扰、检测、信号分析、脉冲压缩 滤波估计:SOC估计 目标定位:WSN定位、滤波跟踪、目标定位 生物电信号:肌电信号EMG、脑电信号EEG、心电信号ECG 通信系统:DOA估计、编码译码、变分模态分解、管道泄漏、滤波器、数字信号处理+传输+分析+去噪、数字信号调制、误码率、信号估计、DTMF、信号检测识别融合
资源推荐
资源详情
资源评论
收起资源包目录
【心电信号】基于matlab心电信号特征提取+分析处理【含Matlab源码 289期】.zip (25个子文件)
【心电信号】基于matlab心电信号特征提取+分析处理【含Matlab源码 289期】
运行结果13.jpg 37KB
Untitled.m 5KB
100-2-3.mat 3KB
运行结果6.jpg 43KB
运行结果4.jpg 51KB
运行结果7.jpg 43KB
运行结果3.jpg 36KB
运行结果19.jpg 0B
105-2-3.mat 3KB
运行结果14.jpg 38KB
运行结果1.jpg 39KB
运行结果12.jpg 41KB
运行结果11.jpg 29KB
11.mat 3KB
运行结果5.jpg 51KB
运行结果15.jpg 32KB
运行结果16.jpg 41KB
9.mat 4KB
生物医学信号检测与处理 期末综合测评.doc 747KB
运行结果8.jpg 36KB
运行结果9.jpg 42KB
运行结果10.jpg 42KB
运行结果17.jpg 25KB
运行结果2.jpg 37KB
运行结果18.jpg 41KB
共 25 条
- 1
海神之光
- 粉丝: 3w+
- 资源: 2085
上传资源 快速赚钱
- 我的内容管理 展开
- 我的资源 快来上传第一个资源
- 我的收益 登录查看自己的收益
- 我的积分 登录查看自己的积分
- 我的C币 登录后查看C币余额
- 我的收藏
- 我的下载
- 下载帮助
安全验证
文档复制为VIP权益,开通VIP直接复制
信息提交成功
- 1
- 2
- 3
前往页