小波变换语音消噪改进阈值.docx

1、小波变换语音消噪改进阈值改进阈值函数进行语音信号消噪，但是在程序运行过程中频频报错。本人经验不足调试不出，希望求得各位指导。改进函数表达式附图clear all; clc; close all;fs=8000; %语音信号采样频率为8000xx=wavread(lw1.wav);x1=xx(:,1);%取单声道t=(0:length(x1)-1)/8000;y1=fft(x1,2048); %对信号做2048点FFT变换f=fs*(0:1023)/2048;figure(1)plot(t,x1) %做原始语音信号的时域图形y=awgn(x1,10,measured); %加10db的高斯白噪声

2、snr,mse=snrmse(x1,y)%求得信噪比 均方误差figure(2)plot(t,y) %做加噪语音信号的时域图形c,l=wavedec(y,3,db1);%多尺度一维分解%用db1小波对信号进行3层分解并提取系数a3=appcoef(c,l,db1,3);%a2=appcoef(c,l,db1,2);%a1=appcoef(c,l,db1,1);d3=detcoef(c,l,3);d2=detcoef(c,l,2);d1=detcoef(c,l,1);thr1=thselect(d1,rigrsure);%阈值获取，使用Stein的无偏风险估计原理thr2=thselect(d2

3、,rigrsure);thr3=thselect(d3,rigrsure);%利用改进阈值函数进行去噪处理gd1=Garrote_gg(d1,thr1);gd2=Garrote_gg(d2,thr2);gd3=Garrote_gg(d3,thr3);c1=a3 gd3 gd2 gd1;y1=waverec(c2,l,db1);%多尺度重构snr,mse=snrmse(x1,y1)%求得信噪比 均方误差figure(3);plot(t,y1);function gd=Garrote_gg(a,b)%a为信号分解后的小波系数，b为获得的阈值m=0.2*(a*a)-(b*b);if (abs(a)=

4、b) gd=sign(a)*(abs(a)-b/exp(m);else (abs(a)=b) gd=sign(a)*(abs(a)-b/exp(m);else gd=zeros(size(a);endy1=waverec(c1,l,db1);%多尺度重构snr,mse=snrmsel(x1,y1) %求得信噪比 均方误差figure(3);plot(t,y1);小波去噪软阈值和硬阈值的matlab仿真程序硬阈值、软阈值这里有一段不知道有用没%设置信噪比和随机种子值snr=4;init=2055615866;%产生原始信号sref和高斯白噪声污染的信号ssref,s=wnoise(1,11,sn

5、r,init);%用db1小波对原始信号进行3层分解并提取系数c,l=wavedec(s,3,db1);a3=appcoef(c,l,db1,3);d3=detcoef(c,l,3);d2=detcoef(c,l,2);d1=detcoef(c,l,1);thr=1;%进行硬阈值处理ythard1=wthresh(d1,h,thr);ythard2=wthresh(d2,h,thr);ythard3=wthresh(d3,h,thr);c2=a3 ythard3 ythard2 ythard1;s3=waverec(c2,l,db1);%进行软阈值处理ytsoftd1=wthresh(d1,s

6、,thr);ytsoftd2=wthresh(d2,s,thr);ytsoftd3=wthresh(d3,s,thr);c3=a3 ytsoftd3 ytsoftd2 ytsoftd1;s4=waverec(c3,l,db1);%对上述信号进行图示subplot(5,1,1);plot(sref);title(参考信号);subplot(5,1,2);plot(s);title(染噪信号);subplot(5,1,3);plot(s3);title(硬阈值处理);subplot(5,1,4);plot(s4);title(软阈值处理);matlab小波除噪，为何硬阈值和软阈值除躁信噪比一样了？

7、matlab小波除噪，为何硬阈值和软阈值除躁信噪比一样了？load leleccum;index=1:1024;f1=leleccum(index); % 产生含噪信号init=2055615866;randn(seed,init);f2=f1+18*randn(size(x);snr=SNR_singlech(f1,f2) %信噪比subplot(2,2,1);plot(f1);title(含噪信号); %axis(1,1024,-1,1);subplot(2,2,2);plot(f2);title(含噪信号); %axis(1,1024,-1,1);%用db5小波对原始信号进行3层分解并提

8、取系数c,l=wavedec(f2,3,db6);a3=appcoef(c,l,db6,3);d3=detcoef(c,l,3);d2=detcoef(c,l,2);d1=detcoef(c,l,1);sigma=wnoisest(c,l,1);thr=wbmpen(c,l,sigma,2);%进行硬阈值处理ythard1=wthresh(d1,h,thr);ythard2=wthresh(d2,h,thr);ythard3=wthresh(d3,h,thr);c2=a3 ythard3 ythard2 ythard1;f3=waverec(c2,l,db6);%进行软阈值处理ytsoftd1

9、=wthresh(d1,s,thr);ytsoftd2=wthresh(d2,s,thr);ytsoftd3=wthresh(d3,s,thr);c3=a3 ytsoftd3 ytsoftd2 ytsoftd1;f4=waverec(c3,l,db6);%对上述信号进行图示subplot(2,2,3);plot(f3);title(硬阈值处理);%axis(1,1024,-1,1);subplot(2,2,4);plot(f4);title(软阈值处理);%axis(1,1024,-1,1);snr=SNR_singlech(f1,f3)snr=SNR_singlech(f1,f4)信噪比函数

10、SNR_singlech(I,In)function snr=SNR_singlech(I,In)% 计算信噪比函数% I：riginal signal% In:noisy signal(ie. original signal + noise signal)snr=0;Ps=sum(sum(I-mean(mean(I).2);%signal powerPn=sum(sum(I-In).2); %noise powersnr=10*log10(Ps/Pn);小波去噪程序Matlab小波去噪（默认，强制，给定三种情况）% 利用小波分析对监测采集的信号进行去噪处理，恢复原始信号%小波分析进行去噪有3

11、中方法：%1、默认阈值去噪处理。该方法利用函数ddencmp( )生成信号的默认阈值，然后利用函数wdencmp( )进行去噪处理；%2、给定阈值去噪处理。在实际的去噪处理过程中，阈值往往可通过经验公式获得，且这种阈值比默认阈值的可信度高。在进行阈值量化处理时可利用函数wthresh( );%3、强制去噪处理。该方法是将小波分解结构中的高频系数全部置0，即滤掉所有高频部分，然后对信号进行小波重构。这种方法比较简单，且去噪后的信号比较平滑，但是容易丢失信号中的有用成分。% 利用小波分析对监测采集的水轮机信号进行去噪处理，恢复原始信号%Program Start% 载入监测所得信号load def

12、ault.txt;%装载采集的信号x= default;lx=length(x);t=0:1:length(x)-1;% 绘制监测所得信号subplot(2,2,1);plot(t,x);title(原始信号);grid onset(gcf,color,w)set(gca,fontname,times New Roman)set(gca,fontsize,14.0)% 用db1小波对原始信号进行3层分解并提取小波系数c,l=wavedec(x,3,db1);%sym8ca3=appcoef(c,l,db1,3);%低频部分cd3=detcoef(c,l,3);%高频部分cd2=detcoef(

13、c,l,2);%高频部分cd1=detcoef(c,l,1);%高频部分% 对信号进行强制去噪处理并图示cdd3=zeros(1,length(cd3);cdd2=zeros(1,length(cd2);cdd1=zeros(1,length(cd1);c1=ca3,cdd3,cdd2,cdd1;x1=waverec(c1,1,db1);subplot(2,2,2);plot(x1);title(强制去噪后信号);grid onset(gcf,color,w)set(gca,fontname,times New Roman)set(gca,fontsize,14.0)% 默认阈值对信号去噪并图

14、示%用ddencmp( )函数获得信号的默认阈值，使用wdencmp( )函数实现去噪过程thr,sorh,keepapp=ddencmp(den,wv,x);x2=wdencmp(gbl,c,l,db1,3,thr,sorh,keepapp);subplot(2,2,3);plot(x2);title(默认阈值去噪后信号);grid onset(gcf,color,w)set(gca,fontname,times New Roman)set(gca,fontsize,14.0)% 给定的软阈值进行去噪处理并图示cd1soft=wthresh(cd1,x,1.465);%经验给出软阈值数cd2

15、soft=wthresh(cd2,x,1.823); %经验给出软阈值数cd3soft=wthresh(cd3,x,2.768); %经验给出软阈值数c2=ca3,cd3soft,cd2soft,cd1soft;x3=waverec(c2,1,db1);subplot(2,2,4);plot(x3);title(给定软阈值去噪后信号);grid onset(gcf,color,w)set(gca,fontname,times New Roman)set(gca,fontsize,14.0)以上就是三种小波去噪的原程序。但红色标注的地方，我却运行不过去。提示分别问：1，? Error using

16、 = horzcatCAT arguments dimensions are not consistent.Error in = tt at 33c1=ca3,cdd3,cdd2,cdd1;2，? Error using = wthresh at 30Invalid argument value.Error in = tt at 54cd1soft=wthresh(cd1,x,1.465);%经验给出软阈值数请问为什么？关于小波新阈值去噪改进方法里，matlab程序编写的讨论新阈值能有效地克服软阈值去噪方法中由于估计值与真实值之间的恒定偏差而带来的去噪误差，也能有效地抑制硬阈值去噪方法中易产生

17、的信号振荡现象。但不懂它程序的编写。 请教高手帮忙解决一下！谢谢了！function y = wthresh(x,sorh,t)%WTHRESH Perform soft or hard thresholding.% Y = WTHRESH(X,SORH,T) returns soft (if SORH = s)% or hard (if SORH = h) T-thresholdingof the input% vector or matrix X. T is the threshold value.% Y = WTHRESH(X,s,T) returns Y = SIGN(X).(|X|-

18、T)+, soft% thresholding is shrinkage.% Y = WTHRESH(X,h,T) returns Y = X.1_(|X|T), hard% thresholding is cruder.% See also WDEN, WDENCMP, WPDENCMP.% M. Misiti, Y. Misiti, G. Oppenheim, J.M. Poggi 12-Mar-96.% Last Revision: 24-Jul-2007.% Copyright 1995-2007 The MathWorks, Inc.% $Revision: 1.11.4.3 $sw

19、itch sorhcase s tmp = (abs(x)-t); tmp = (tmp+abs(tmp)/2; y = sign(x).*tmp;case h y = x.*(abs(x)t);otherwise error(Wavelet:FunctionArgVal:Invalid_ArgVal, . Invalid argument value.)end这是传统的阈值去噪法，只要你有新的算法，直接改就行了，但是我对该法有些疑惑，wthresh（x,sorh,thr）,若我用赵瑞珍法算出了每个尺度（假设J=3）的阈值，那么还用wthresh（x,sorh,thr）去噪时，thr是不是应该

20、在之前写成thr=thr1,thr2,thr3呢？求助高手：Matlab小波去噪中遇见的问题小波程序如下：function y=f1(x);X=xlsread(shoulian.xls); %读入表格数据x=X(:,1); %将表格数据的第一列赋值给x轴y=X(:,2); %将表格数据的第二列赋值给y轴subplot(2,2,1);plot(x,y,-m*) %画出图形hold on;title(原始图形);grid;%画出网格%用db4小波对原始信号进行3层分解并提取系数还没有写好:要增加高频系数图c,l=wavedec(y,3,db4); %对y进行3层分解ca3=appcoef(c,l,

21、db4,3); %提取系数cd3=detcoef(c,l,3);cd2=detcoef(c,l,2);cd1=detcoef(c,l,1);%对信号进行强制性去噪处理并图示结果cdd3=zeros(1,length(cd3);cdd2=zeros(1,length(cd2);cdd1=zeros(1,length(cd1);c1=ca3 cdd3 cdd2 cdd1;s1=waverec(c1,l,db4);subplot(2,2,2);plot(x,s1,-m*);hold on;title(强制去噪后的信号);grid;%用默认阈值对信号进行去噪处理并图示结果%用ddencmp()函数获得

22、信号的默认阈值,使用wdencmp()命令函数实现去噪过程thr,sorh,keepapp=ddencmp(den,wv,y);s2=wdencmp(gbl,c,l,db4,3,thr,sorh,keepapp);subplot(2,2,3);plot(x,s2,-m*);hold on;title(默认阈值去噪后的信号);grid;%用给定的软阈值进行去噪处理cdd1soft=wthresh(cd1,s,1.465);cdd2soft=wthresh(cd2,s,1.823);cdd3soft=wthresh(cd3,s,2.768);c2=ca3 cdd3soft cdd2soft cdd

23、1soft;s3=waverec(c2,l,db4);subplot(2,2,4);plot(x,s3,-m*);title(给定软阈值去噪后的信号);gridclear;错误显示：? Error using = horzcatAll matrices on a row in the bracketed expression must have thesame number of rows.Error in = xnn at 41c2=ca3 cdd3soft cdd2soft cdd1soft;i以至于最后一个图形没法显示同时处理后的图形数据怎么调出来，以及第一个数据应该是【0，0】，但是处理后的不是。请问是我选择的分解重构函数部队还是其他原因请问如何修改需要说明的是：收敛表格的两组数据为：00.0020.4940.9451.0361.0371.48101.49111.86141.96191.78262.75282.78312.71342.54392.85472.71492.71545.43615.47685.47755.20请教程序问题，心电信号去噪的matlab程序想把这个程序。改为只用软硬阈值对比的心电信号去噪分析，该怎么办？ 求高手。程序见下：%信号小波分解%基于Haar小波%ecg=fopen(100.dat,r);N=1201;