%%%%%%%%%%%%%%%%%%心电信号降噪
目前创新互联已为上千的企业提供了网站建设、域名、虚拟主机、网站改版维护、企业网站设计、新罗网站维护等服务,公司将坚持客户导向、应用为本的策略,正道将秉承"和谐、参与、激情"的文化,与客户和合作伙伴齐心协力一起成长,共同发展。
%%%%%%%%%%%%%%%Birge-Massart策略阈值降噪
%基于小波变换的心电信号的降噪
ecg=fopen('100.dat','r');% 调用心电数据库 r为只读,ecg是打开文件的识别符
N=1201;%常数赋值,要读数据个数
data=fread(ecg,N,'int16'); %从一个流中读N个数据,数据格式是int16,16进制整数
data=data/10000;%数据缩小10000倍
fclose(ecg);%关闭打开的文件
x=data;%把数据转赋给x变量
wavename='db5'; %db5是小波名
level=4;%4级分解
[c,l]=wavedec(x,level,wavename); %4级小波分解,c保存各级分解系数,l是薄记矩阵,保存各级的系数的个数
alpha=1.5; %1.5用于信号压缩,3用于降噪
sorh='h'; %为硬阈值
[thr,nkeep]=wdcbm(c,l,alpha);%使用Birgé-Massart策略计算一维小波分解或压缩的阈值thr和各级的系数个数nkeep
[xc,cxc,lxc,perf0,perfl2]=wdencmp('lvd',c,l,wavename,level,thr,sorh); %小波压缩重构后的图像
t1=0:0.004:(length(x)-1)*0.004;%一行数据
figure(4);%打开一个图形窗口
subplot(211); %子图1
plot(t1,x);%画图形
title('从人体采集的原始的ECG信号');%加上子图名称
subplot(212);%子图2
plot(t1,xc);%画图形
title('Birge-Massart策略阈值降噪后的ECG信号(wname=db5 level=4)');%加上子图名称
load noisdopp %读入信号
x=noisdopp
subplot(211);plot(x);title('原始信号');
s1=wden(x,'minimaxi','s','one',5,'db4'); %用db4小波进行软阈值去噪
subplot(212);plot(s1);title('使用sym4小波最小极大阈值软阈值降噪后的信号');
s %自己定义
%画出原始信号
subplot(221);
plot(s);
title('原始信号');
ylabel('幅值A');
%用db3小波对信号进行3层分解并提取系数
[c,l]=wavedec(s,3,'db3');
a3=appcoef(c,l,'db3',3);
d3=detcoef(c,l,3);
d2=detcoef(c,l,2);
d1=detcoef(c,l,1);
%强制消噪处理
dd3=zeros(1,length(d3));
dd2=zeros(1,length(d2));
dd1=zeros(1,length(d1));
c1=[a3 dd3 dd2 dd1];
s1=waverec(c1,l,'db3');
subplot(222);
plot(s1);
title('强制消噪信号');
%默认阈值进行消噪
%用ddencmp函数获得信号的默认阈值
[thr,sorh,keepapp]=ddencmp('den','wv',s);
s2=wdencmp('gbl',c,l,'db3',3,thr,sorh,keepapp);
subplot(223);
plot(s2);
title('默认阈值消噪');
xlabel('样本序号n');
ylabel('幅值A');
%用给定的软阈值进行消噪
softd1=wthresh(d1,'s',1.465);
softd2=wthresh(d2,'s',1.823);
softd3=wthresh(d3,'s',2.768);
c2=[a3 softd3 softd2 softd1];
s3=waverec(c2,l,'db3');
subplot(224);
plot(s3);
title('给定软阈值消噪');
哈哈!傻眼了吧,理论和实际不可兼得啊!对于你这个问题有如下建议:
小波去噪的试验,十个有九个都喜欢用正、余弦函数,但由于小波函数的几何特征,其实不易得到满意效果,你只能选择线性强的小波基,即对称性强并且光滑的小波,嘿嘿,貌似能做DWT的所有小波基中只有sym8的对称性最强(参看;msgtype=2),这样你将会得到其降噪效果好于其它小波基的谬论,实际应用中如果原始信号本身特征就不规律对称则有可能其它小波基的处理效果会好于sym8,这需要试错的方法试验,是你这个原始信号对称的前提假设造成的,但不可认为sym8适合你这个特殊例子就得出结论,你这试验只适合你的假设,其实并不适合所有情形,个人认为其实不具说服力和可比性。如同我家钥匙开不了你家门,是不能得出我家钥匙不如你家钥匙的结论滴。
当你的噪声特征与信号的特征的频率相近时,小波也无能为力,它不是神,也一样分不出啥是噪声,所以其一是可增加分解层数,你这个信号只有100个数,5层已经很高了,再增大也没啥用了,可能会过多显示小波基的特征,造成扭曲失真(如果用SWT会好很多,但需要自己编制函数);其二是参数SCAL可以改为伸缩的sln,而不是固定的one,这样分解层数和SCAL都将起作用,你可以试着改改玩,效果还行。
分解5层
分解8层
你可以试试只分解一层的状况,频率最低的几处噪声会保留下来哦!