數(shù)字信號(hào)處理課程項(xiàng)目

1、數(shù)字信號(hào)處理B課程項(xiàng)目實(shí)施報(bào)告題 目： 數(shù)字音效處理器 . 目錄一、設(shè)計(jì)要求：3二、盲信號(hào)處理3（1）Fast ICA算法簡(jiǎn)介3（2）遇到問題5三、去噪5（1）遇到問題5（2）維納濾波器6四、性能與結(jié)果分析71維納濾波7（1）去噪時(shí)域圖7（2） 去噪頻域圖8（3） 頻域瀑布圖82盲信號(hào)分離9（1）分離男聲9（2） 分離女聲93用戶界面10五、源程序代碼11（1）FastICA盲信號(hào)分離代碼：11（2）維納濾波去噪：13（3） 瀑布圖代碼：15（4） 主界面代碼：16七、小組分工17一、設(shè)計(jì)要求：設(shè)計(jì)一個(gè)基于MATLAB的數(shù)字音頻分析與處理系統(tǒng)，能夠?qū)崿F(xiàn)對(duì)數(shù)字音頻的測(cè)試分析與處理。要求：（1）

2、輸入語音信號(hào)源為實(shí)際環(huán)境采集語音；（2） 至少實(shí)現(xiàn)1種音效測(cè)試分析功能（頻率響應(yīng)，瀑布頻譜圖，相位響應(yīng)曲線，抗阻曲線，諧波失真曲線，互調(diào)失真曲線，音壓曲線等）；并至少實(shí)現(xiàn)1中語音處理功能（聲紋識(shí)別，關(guān)鍵詞識(shí)別，語音去噪，聲源分離等）；并分別對(duì)其進(jìn)行實(shí)際運(yùn)用。（3） 用人機(jī)交互界面操控、揚(yáng)聲器/耳機(jī)輸出音效。二、盲信號(hào)處理（1）Fast ICA算法簡(jiǎn)介快速ICA算法，由于采用了定點(diǎn)迭代的思想，所以又稱固定點(diǎn)(Fixed-Point)算法或定點(diǎn)算法，是一種快速尋優(yōu)神經(jīng)網(wǎng)絡(luò)算法。與普通的神經(jīng)網(wǎng)絡(luò)算法不同的是這種算法采用了牛頓迭代的思想對(duì)采樣點(diǎn)采用批處理的方式，在每一步迭代中有大量的樣本數(shù)據(jù)參與運(yùn)算。

3、基于目標(biāo)函數(shù)的選擇不同，從而產(chǎn)生了FastICA算法的基于負(fù)熵最大、峭度、互信息最小、最大似然估計(jì)等的多種形式。FastICA算法的優(yōu)點(diǎn)如下：1.收斂速度快，F(xiàn)astICA收斂速度為2次以上，普通的ICA算法收斂?jī)H為1次。2.能利用任何的非線性函數(shù)直接找出任何非高斯分布的獨(dú)立分量，而對(duì)于其他的算法來說，概率密度函數(shù)的估計(jì)不得不首先進(jìn)行。3.獨(dú)立分量可被逐個(gè)估計(jì)出來，類似于做投影追蹤，這在僅需要估計(jì)幾個(gè)（不是全部）獨(dú)立分量的情況下，能減小計(jì)算量。4.Fast-ICA算法有許多神經(jīng)算法的優(yōu)點(diǎn)，它是并行、分布式的且計(jì)算簡(jiǎn)單，內(nèi)存要求很少，它的性能能夠通過選擇適當(dāng)?shù)姆蔷€性函數(shù)來最佳化。從Fast-IC

4、A算法流程圖可以看出，通常先對(duì)數(shù)據(jù)進(jìn)行初步的白化或球化處理，白化處理可去除各觀測(cè)信號(hào)之間的相關(guān)性，從而簡(jiǎn)化了后續(xù)獨(dú)立分量的提取過程。這是由于一般情況下，所獲得的數(shù)據(jù)都具有相關(guān)性。通過去均值和白化等預(yù)處理在對(duì)其進(jìn)行真正的算法分離。通常情況下，數(shù)據(jù)進(jìn)行白化處理與不對(duì)數(shù)據(jù)進(jìn)行白化處理相比，算法的收斂性比較好。概括來說，F(xiàn)ast-ICA算法主要包含兩部分:預(yù)處理部分和核心算法部分。預(yù)處理部分主要包括去均值化和白化處理；核心算法部分是基于擬牛頓算法，具體的實(shí)現(xiàn)見附錄程序。 Fast-ICA算法流程圖:觀測(cè)信號(hào)x去均值對(duì)去均值后的觀測(cè)信號(hào)x進(jìn)行白化基于負(fù)熵最大化理論確定目標(biāo)函數(shù)應(yīng)用牛頓迭代算法求

5、目標(biāo)函數(shù)的最大值，即分離矩陣w收斂？由獲得信號(hào)源是否在所有等方差的隨機(jī)變量中，高斯變量的熵最大，因而我們可以利用熵來度量非高斯性，常用熵的修正形式，即負(fù)熵。在分離過程中，可通過對(duì)分離結(jié)果的非高斯性度量來表示分離結(jié)果間的相互獨(dú)立性，當(dāng)非高斯性度量達(dá)到最大時(shí)，則表明已完成對(duì)各獨(dú)立分量的分離。負(fù)熵的定義： ，H為隨機(jī)變量的微分熵當(dāng)Y具有高斯分布時(shí)Ng（Y）=0，Y的非高斯性越強(qiáng)，其微分熵越小，所以Ng（Y）可以作為隨機(jī)變量Y非高斯性的測(cè)度。快速ICA學(xué)習(xí)規(guī)則是找一個(gè)方向以便具有最大的非高斯性。FastICA算法的迭代公式：（2）遇到問題對(duì)于聲源分離的理解一開始我們將其理解為經(jīng)典的“雞尾酒會(huì)效應(yīng)”，經(jīng)

6、過大量的數(shù)據(jù)查詢和整理，我們發(fā)現(xiàn)目前情況下，我們只能用機(jī)器實(shí)現(xiàn)分離簡(jiǎn)單疊加的盲信號(hào)，或者說是以特定形式混合在一起的信號(hào)，不能分離一些復(fù)雜的盲信號(hào)。在近幾年的有效方法中，F(xiàn)astICA算法可以有效的分離簡(jiǎn)單疊加的盲信號(hào)，并且可以用matlab編程的方式得以實(shí)現(xiàn)。于是，我們使用matlab中的FastICA算法進(jìn)行了簡(jiǎn)單疊加的盲信號(hào)分離。三、去噪（1）遇到問題在本次課程項(xiàng)目中，我們一開始選用了最基礎(chǔ)的低通濾波器，并花了一些時(shí)間做出了巴特沃斯低通濾波器進(jìn)行測(cè)試。但是除了語音變得更沉悶一些了對(duì)噪音并沒有什么實(shí)質(zhì)性的消除效果。因?yàn)槲覀兘M要消除的噪音是雪花噪音。我們嘗試了維納濾波器，維納濾波器可以將大部分

7、噪音都濾除，只留下需要的人聲信號(hào)以及幾乎聽不見的混在人聲特定地方的噪音。（2）維納濾波器維納（Wiener）是用來解決從噪聲中提取信號(hào)的一種過濾（或?yàn)V波）方法。這種線性濾波問題，可以看作是一種估計(jì)問題或一種線性估計(jì)問題。 一個(gè)線性系統(tǒng)，如果它的單位樣本響應(yīng)為，當(dāng)輸入一個(gè)隨機(jī)信號(hào)，且 （1）其中表示信號(hào)，)表示噪聲，則輸出為 （2）我們希望通過線性系統(tǒng)后得到的盡量接近于，因此稱為的估計(jì)值，用表示，即 （3） 則維納濾波器的輸入輸出關(guān)系可用下面圖1表示。圖1實(shí)際上，式（2）所示的卷積形式可以理解為從當(dāng)前和過去的觀察值，來估計(jì)信號(hào)的當(dāng)前值。因此，用進(jìn)行過濾問題實(shí)際上是一種統(tǒng)計(jì)估計(jì)問題。一般地，從當(dāng)前

8、的和過去的觀察值，估計(jì)當(dāng)前的信號(hào)值成為過濾或?yàn)V波；從過去的觀察值，估計(jì)當(dāng)前的或者將來的信號(hào)值稱為外推或預(yù)測(cè)；從過去的觀察值，估計(jì)過去的信號(hào)值稱為平滑或內(nèi)插。因此維納濾波器又常常被稱為最佳線性過濾與預(yù)測(cè)或線性最優(yōu)估計(jì)。這里所謂的最佳與最優(yōu)是以最小均方誤差為準(zhǔn)則的。如果我們分別以與表示信號(hào)的真實(shí)值與估計(jì)值，而用表示他們之間的誤差，即 （4）顯然可能是正值，也可能是負(fù)值，并且它是一個(gè)隨機(jī)變量。因此，用它的均方誤差來表達(dá)誤差是合理的，所謂均方誤差最小即它的平方的統(tǒng)計(jì)期望最?。?（5）四、性能與結(jié)果分析1維納濾波（1）去噪時(shí)域圖從無語音段提取出的噪聲波形維納濾波后波形 信噪比=14.8306dB帶噪語音

9、 信噪比=4.8483dB（2） 去噪頻域圖帶噪語音頻譜噪聲頻譜處理后的語音頻譜（3） 頻域瀑布圖2盲信號(hào)分離（1）分離男聲時(shí)域圖頻域圖（2） 分離女聲時(shí)域圖頻域圖3用戶界面五、源程序代碼（1）FastICA盲信號(hào)分離代碼：x,fs=audioread('D:shuzixinhaogongchengXMhunhe1.wav');y,fs=audioread('D:shuzixinhaogongchengXMhunhe2.wav');%sound(x,fs);  % /播放原語音信號(hào)N1=length(x);  %

10、 /求出語音信號(hào)的長(zhǎng)度N2=length(y);t1=(0:N1-1)/fs;    % /語音信號(hào)的長(zhǎng)度除以采樣頻率t2=(0:N2-1)/fs;w1=fs*linspace(0,1,N1); w2=fs*linspace(0,1,N2); x1=abs(fft(x,N1);  y1=abs(fft(y,N2);  x=x(:,1);y=y(:,1);x=x'y=y'MixedS = x;y;MixedS_bak=MixedS; 

11、;  % 將混合后的數(shù)據(jù)備份，以便在恢復(fù)時(shí)直接調(diào)用%  中心化  %MixedS_mean=zeros(2,1);for i=1:2    MixedS_mean(i)=mean(MixedS(i,:);endfor i=1:2    for j=1:size(MixedS,2)        MixedS(i,j)=MixedS(i,j

12、)-MixedS_mean(i);    endend%  白化  %MixedS_cov=cov(MixedS');           % cov為求協(xié)方差的函數(shù)E,D=eig(MixedS_cov);            % 對(duì)信號(hào)矩陣的協(xié)方差函數(shù)進(jìn)行特

13、征值分解Q=inv(sqrt(D)*(E)'           % Q為白化矩陣MixedS_white=Q*MixedS;     % MixedS_white為白化后的信號(hào)矩陣IsI=cov(MixedS_white');   % IsI應(yīng)為單位陣  %FASTICA算法  %X=MixedS_white; 

14、;        % 以下算法將對(duì)X進(jìn)行操作VariableNum,SampleNum=size(X);numofIC=VariableNum;   % 在此應(yīng)用中，獨(dú)立元個(gè)數(shù)等于變量個(gè)數(shù)B=zeros(numofIC,VariableNum);  % 初始化列向量w的寄存矩陣,B=b1  b2  .   bdfor r=1:numofIC   &

15、#160;i=1;maxIterationsNum=100;                % 設(shè)置最大迭代次數(shù)（即對(duì)于每個(gè)獨(dú)立分量而言迭代均不超過此次數(shù)）    b=rand(numofIC,1)-.5;                &

16、#160;  % 隨機(jī)設(shè)置b初值    b=b/norm(b);                           % 對(duì)b標(biāo)準(zhǔn)化 norm(b):向量元素平方和開根號(hào)    while i<

17、=maxIterationsNum+1        bOld=b;                                  t=X'*b;  

18、;      g=t.*exp(-t.2/2);        dg=(1-t.2).*exp(-t.2/2);        b=X*g/SampleNum-mean(dg)*b;        b=b-B*B'*b;      

19、           % 對(duì)b正交化        b=b/norm(b);         if abs(abs(b'*bOld)-1)<1e-9    % 如果收斂，則     &#

20、160;       B(:,r)=b;          % 保存所得向量b             break;         end     

21、   i=i+1;            endendif i = maxIterationsNum+1          % 循環(huán)結(jié)束處理      fprintf('n第%d分量在%d次迭代內(nèi)并不收斂。', r,maxI

22、terationsNum);      return;end%  ICA計(jì)算的數(shù)據(jù)復(fù)原并構(gòu)圖  %ICAedS=B'*Q*MixedS_bak;          % 計(jì)算ICA后的矩陣% 將混合矩陣重新排列并輸出sound(ICAedS(1,:),fs);Ws=sum(ICAedS(1,:).2);%計(jì)算信噪比Wn=sum(ICAedS(2,:).2);A=B'

23、;*Q;a11=A(1,1);a12=A(1,2);a21=A(2,1);a22=A(2,2);SNRy1=10*(a112*Ws)/(a122*Wn)（2）維納濾波去噪：if (nargin<3 | isstruct(IS)   % 如果輸入?yún)?shù)小于3個(gè)或IS是結(jié)構(gòu)數(shù)據(jù)   IS=.25; endW=fix(.025*fs);             

24、;  % 幀長(zhǎng)為25msSP=.4;                                      % 幀移比例取40%(10ms) wnd=hamming(

25、W);                % 設(shè)置窗函數(shù)% 如果輸入?yún)?shù)大于或等于3個(gè)并IS是結(jié)構(gòu)數(shù)據(jù)(為了兼容其他程序)if (nargin>=3 && isstruct(IS)    SP=IS.shiftsize/W;    nfft=IS.nfft;   

26、; wnd=IS.window;    if isfield(IS,'IS')        IS=IS.IS     else        IS=.25;    endendpre_emph=0;signal=filter(1 -pre_emph,1,sign

27、al);      % 預(yù)加重NIS=fix(IS*fs-W)/(SP*W) +1);           % 計(jì)算無話段幀數(shù)y=segment(signal,W,SP,wnd);            % 分幀 Y=fft(y);     

28、;                              % FFTYPhase=angle(Y(1:fix(end/2)+1,:);         % 帶噪語音的相位角 Y=abs(

29、Y(1:fix(end/2)+1,:);            % 取正頻率譜值numberOfFrames=size(Y,2);             % 計(jì)算總幀數(shù)FreqResol=size(Y,1);          

30、     % 計(jì)算頻譜中的譜線數(shù)N=mean(Y(:,1:NIS)')'                % 計(jì)算無話段噪聲平均譜值 LambdaD=mean(Y(:,1:NIS)').2)'          % 初始噪聲功率譜方差alpha=.99;

31、60;                    % 設(shè)置平滑系數(shù)fn=numberOfFrames;miniL=5;                        

32、;            % 設(shè)置miniLvoiceseg,vosl,SF,Ef=pitch_vad1(y,fn,T1,miniL); %端點(diǎn)檢測(cè) NoiseCounter=0;                 % 初始化NoiseCounterNoiseLength

33、=9;                  % 設(shè)置噪聲平滑區(qū)間長(zhǎng)度G=ones(size(N);                % 初始化譜估計(jì)器Gamma=G;X=zeros(size(Y);    

34、;            % 初始化Xh=waitbar(0,'Wait.');            % 設(shè)置運(yùn)行進(jìn)度條圖 for i=1:numberOfFrames    SpeechFlag=SF(i);    if

35、 i<=NIS                               % 若i<=NIS在前導(dǎo)無聲(噪聲)段        SpeechFlag=0;  

36、;      NoiseCounter=100;    %else                   % i>NIS判斷是否為有話幀        %NoiseFlag, SpeechFlag,

37、0;NoiseCounter, Dist=vad(Y(:,i),N,NoiseCounter);     end    if SpeechFlag=0              % 在無話段中平滑更新噪聲譜值        N=(NoiseLength*N+Y

38、(:,i)/(NoiseLength+1);         LambdaD=(NoiseLength*LambdaD+(Y(:,i).2)./(1+NoiseLength);%更新和平滑噪聲方差    end     gammaNew=(Y(:,i).2)./LambdaD;          % 計(jì)算后驗(yàn)

39、信噪比    xi=alpha*(G.2).*Gamma+(1-alpha).*max(gammaNew-1,0); % 計(jì)算先驗(yàn)信噪比    Gamma=gammaNew;    G=(xi./(xi+1);         % 計(jì)算維納濾波器的譜估計(jì)器    X(:,i)=G.*Y(:,i);  

40、0;       % 維納濾波后的幅值%顯示運(yùn)行進(jìn)度條圖    waitbar(i/numberOfFrames,h,num2str(fix(100*i/numberOfFrames);endclose(h);                       &#

41、160;           % 關(guān)閉運(yùn)行進(jìn)度條圖output=OverlapAdd2(X,YPhase,W,SP*W);        % 語音合成output=filter(1,1 -pre_emph,output);      % 消除預(yù)加重影響output=output/max(abs(output)

42、; （3） 瀑布圖代碼：b,Fs=audioread('D:shuzixinhaogongchengXMxuehua.wav');b=b(:,1);t=size(b)/44100;f=0:1:22049;R=1;G=0;B=0.5;figurefor i=1:0.2:t/2; z=b(i*44100-44099):(i*44100); z=abs(fft(z,44100)/22050; z=z' f,i=meshgrid(f,i); plot3(f,i,z(1:22050),'color',R G B); hold on; G=G+0.035;endR=0.5;G=1;B=0;for i=t/2:0