統(tǒng)計信號處理實驗報告一---沈昊骍-04008409

1、統(tǒng)計信號處理實驗報告（實驗一）姓名： 沈昊骍學(xué)號：04008409 時間： 2011.11.18統(tǒng)計信號處理實驗一實驗?zāi)康模?、掌握噪聲中信號檢測的方法； 2、熟悉Matlab的使用； 3、掌握用計算機進行數(shù)據(jù)分析的方法；實驗內(nèi)容：假設(shè)信號為波形如下圖所示在有信號到達時接收到的信號為，在沒有信號到達時接收到的信號為。其中是均值為零、方差為的高斯白噪聲。對接受到的信號分別在上進行取樣，得到觀測序列。利用似然比檢測方法，對信號是否到達進行檢測；假設(shè)有信號到達的概率P(H1)=0.6，沒有信號到達的概率P(H0)=0.4，。利用Bayes檢測方法，對信號是否到達進行檢測；通過計算機產(chǎn)生的仿真數(shù)據(jù)，對

2、兩種方法的檢測概率、誤警概率、漏警概率和Bayes風(fēng)險進行仿真計算；改變判決的門限，觀察檢測方法的、和Bayes風(fēng)險的變化；改變噪聲的方差，觀察檢測方法的、和Bayes風(fēng)險的變化；將信號取樣間隔減小一倍(相應(yīng)的取樣點數(shù)增加一倍)，觀察似然比檢測方法的、和Bayes風(fēng)險的變化；根據(jù)設(shè)計一個離散匹配濾波器，并觀察經(jīng)過該濾波器以后的輸出。實驗要求：設(shè)計仿真計算的Matlab程序，給出軟件清單完成實驗報告，對實驗過程進行描述，并給出實驗結(jié)果，對實驗數(shù)據(jù)進行分析，給出結(jié)論。實驗過程：1）首先產(chǎn)生信號s(t),n(t),x(t)，即s(i),n(i),x(i), 其中i=0，1， 200；2）根據(jù)定義似然

3、比函數(shù)：，門限；如果，則判定；否則，判定。這就是似然比檢測準(zhǔn)則。假設(shè)似然比為1，根據(jù)似然比檢測準(zhǔn)則：兩邊取對數(shù)后得：由此對信號是否到達進行檢測；3） Bayes判決準(zhǔn)則如下：準(zhǔn)則或風(fēng)險函數(shù)：其中的諸系數(shù)是根據(jù)實際需要設(shè)定的風(fēng)險系數(shù)。時判，否則判。假設(shè)有信號到達的概率P(H1)=0.6，沒有信號到達的概率P(H0)=0.4，。由此計算判決門限為（2*0.4）/（1*0.6）=4/3。由此對信號是否到達進行檢測；4)根據(jù)蒙特卡洛仿真方法分別對以上兩種方法下的檢測概率、誤警概率、漏警概率和Bayes風(fēng)險進行仿真計算：共做M=10000次統(tǒng)計：在x(t)=s(t)+n(t)的情況下，每次出現(xiàn)There

4、 is a signal 則檢測到信號的次數(shù)n1加1； There is not a signal則未檢測到信號的次數(shù)n0加1；在x(t)= n(t)的情況下，每次出現(xiàn)There is a signal 則檢測到信號的次數(shù)n2加1；其中：檢測概率=n1/M;漏警概率=n0/M;誤警概率=n2/M;Bayes風(fēng)險系數(shù)r= c00*(1-pf)+c10*pf+c01*pm+c11*pd.5) 用同（4）的方法，通過改變判決的門限，觀察檢測方法的、和Bayes風(fēng)險的變化； 6) 用同（4）的方法，通過改變噪聲的方差，觀察檢測方法的、和Bayes風(fēng)險的變化；7) 通過改變是s(t)的取樣間隔（由1變?yōu)?/p>

5、0.5），將取樣間隔減小一倍(相應(yīng)的取樣點數(shù)增加一倍)，n(t)也變?yōu)?00個元素的矩陣，然后再來觀察似然比檢測方法的、和Bayes風(fēng)險的變化；8）設(shè)計匹配濾波器h(t)=c*s(T-t)，通過使待檢測信號x(t)經(jīng)過匹配濾波器，即和h(t)進行卷積，得到濾波以后的輸出X(t)。實驗結(jié)果及分析：1）利用似然比和Bayes兩種檢測方法進行信號檢測：ans =最大似然,檢測到有信號到達的次數(shù)為：out1 = 148ans =貝艾斯檢測,檢測到有信號到達的次數(shù)為：out2 = 132分析：由圖像可知，x(t)受高斯隨機噪聲影響比較嚴(yán)重，兩種方法都檢測到信號。2）根據(jù)蒙特卡洛仿真方法: pd1 pm1

6、 pf1 -似然比檢測方法ans = 0.7954 0.2046 0.2126 pd2 pm2 pf2 r-Bayes檢測方法ans = 0.7472 0.2528 0.1645 0.5818分析：似然比檢測方法比門限值為4/3的Bayes檢測方法檢測結(jié)果好一些。3）改變門限值gama=0.7,1,1.1,1.3,1.6 pd1-似然比檢測方法pd1 = 0.8502 0.7906 0.7693 0.7423 0.6974 pm1pm1 = 0.1498 0.2094 0.2307 0.2577 0.3026 pf1pf1 = 0.2693 0.2076 0.1949 0.1613 0.136

7、2 pd2-Bayes檢測方法pd2 = 0.7294 0.6548 0.6181 0.6011 0.5445 pm2pm2 =0.2706 0.3452 0.3819 0.3989 0.4555 pf2pf2 = 0.1589 0.1031 0.0983 0.0815 0.0622 rr = 0.5884 0.5514 0.5785 0.5619 0.5799分析：隨著門限值gama的增大，兩種檢測方法的檢測概率pd都在變小，漏警概率pm增大，而虛警概率pf卻隨門限增大在減小，風(fēng)險系數(shù)r在增大。4）改變方差sigma=4,16,25,36,64 pd1-似然比檢測方法pd1 = 0.9804

8、 0.8409 0.7975 0.7447 0.7015 pm1pm1 = 0.0196 0.1591 0.2025 0.2553 0.2985 pf1pf1 = 0.0239 0.1528 0.2138 0.2503 0.3164 pd2-Bayes檢測方法pd2 = 0.9738 0.8124 0.7374 0.6853 0.5975 pm2pm2 = 0.0262 0.1876 0.2626 0.3147 0.4025 pf2pf2 = 0.0159 0.1204 0.1661 0.1879 0.2210 rr = 0.0580 0.4284 0.5948 0.6905 0.8445 分

9、析：隨著方差sigma的增大，兩種檢測的檢測概率pd都在減小，漏警概率pm和虛警概率pf都在增大，風(fēng)險系數(shù)r快速增大，檢測情況變差，可靠性降低。5）將信號取樣間隔減小一倍(相應(yīng)的取樣點數(shù)增加一倍)： pd1 pm1 pf1-似然比檢測方法ans = 0.8718 0.1282 0.1286 pd2 pm2 pf2 r-Bayes檢測方法ans = 0.8526 0.1474 0.0999 0.3472分析：兩種檢測對比可以看出，信號取樣間隔越小，檢測得到的結(jié)果越好，檢測概率越高，虛警漏警概率越低，風(fēng)險系數(shù)越小。6）根據(jù)s(t)設(shè)計一個離散匹配濾波器h(n):分析：在t0=200處的取值進行觀察

10、；當(dāng)輸入信號中有有用信號時，系統(tǒng)的輸出值可以達到100以上；沒有信號時，系統(tǒng)輸出小于30。這說明通過匹配濾波器后，信號中的有用信號分量得到了加強，信噪比得到了提高。MATLAB程序清單：homework1%To generate signals that should be testedt1=0:49;t2=50:149;t3=150:199;s1=t1/50;s2=-t2/50+2;s3=t3/50-4;s=s1 s2 s3; %To generate signalsfor m=1:200n1=5.*randn(200,1); %To generate noisex=s+n1;subplot

11、(211); %To generate figureplot(s);xlabel(t);ylabel(s);title(信號 s(t);subplot(212);plot(x);xlabel(t);ylabel(x);title(有噪聲x(t);%(1)Using the likelihood ratio(1) to test the signalss1=x.*s;s2=s.*s;s1=sum(s1);s2=sum(s2);if s10.5*s2 a(m)=1; %檢測到信號else a(m)=0; %未檢測到信號endend%(2)Using Bayes method to test the

12、 signalsfor m=1:200c10=2;c01=1;p1=0.6;p0=0.4;r=(c10*p0)/(c01*p1);r=log(r);n2=5.*randn(200,1); %To generate noisex=s+n2;s1=x.*s;s2=s.*s;s1=sum(s1);s2=sum(s2);if s10.5*s2+25*r b(m)=1; %有信號產(chǎn)生Method(2):檢測到信號else b(m)=0; %無信號產(chǎn)生Method(2):未檢測到信號endend最大似然,檢測到有信號到達的次數(shù)為：out1=sum(a)貝艾斯檢測,檢測到有信號到達的次數(shù)為：out2=sum

13、(b) homework1-1%(1)Using the likelihood ratio(1) to calculate pd, pf, pm and rt1=0:49;t2=50:149;t3=150:199;s1=t1/50;s2=-t2/50+2;s3=t3/50-4;s=s1 s2 s3;n1=0; %在x(t)=s(t)+n(t)的情況下，檢測到信號的次數(shù)n0=0; %在x(t)=s(t)+n(t)的情況下，未檢測到信號的次數(shù)n2=0; %在x(t)= n(t)的情況下，檢測到信號的次數(shù)M=10000;for i=1:Mn=5.*randn(200,1); %To generate

14、 noisex2=s+n; %To generate signals that to be testedx1=n;s1=x2.*s;s2=s.*s;s3=x1.*s;s1=sum(s1);s2=sum(s2);s3=sum(s3);if s10.5*s2 n1=n1+1;else n0=n0+1;endif s30.5*s2 n2=n2+1;endendpd1=n1/M; %Using Monte Carlo method to calculate possiblilities 檢測概率pm1=n0/M; %漏警概率pf1=n2/M; %誤警概率%(2)Using Bayes method t

15、o calculate pd, pf, pm and rc10=2;c01=1;c00=0;c11=0;p1=0.6;p0=0.4;gama=(c10*p0)/(c01*p1);gama=log(gama);n1=0;n0=0;n2=0;M=10000;for i=1:Mn=5.*randn(200,1); %To generate noisex2=s+n; %To generate signals that to be testedx1=n;s1=x2.*s;s2=s.*s;s3=x1.*s;s1=sum(s1);s2=sum(s2);s3=sum(s3);if s10.5*s2+25*ga

16、ma n1=n1+1;else n0=n0+1;endif s30.5*s2+25*gama n2=n2+1;endendpd2=n1/M; %Using Monte Carlo method to calculate possibilitiespm2=n0/M;pf2=n2/M;r=c00*(1-pf2)+c10*pf2+c01*pm2+c11*pd2; %Calculate risk indexhomework1-2%(1)(Changing gama)-Using the likelihood ratio(1) to calculate pd, pf, pm and rt1=0:49;t

17、2=50:149;t3=150:199;s1=t1/50;s2=-t2/50+2;s3=t3/50-4;s=s1 s2 s3;gama=0.7,1,1.1,1.3,1.6;for j=1:5n1=0;n0=0;n2=0;M=10000;for i=1:Mn=5.*randn(200,1); %To generate noisex2=s+n; %To generate signals that to be testedx1=n;s1=x2.*s;s2=s.*s;s3=x1.*s;s1=sum(s1);s2=sum(s2);s3=sum(s3);if s10.5*s2+25*log(gama(j)

18、 n1=n1+1;else n0=n0+1;endif s30.5*s2+25*log(gama(j) n2=n2+1;endendpd1(j)=n1/M; %Using Monte Carlo method to calculate possibilitiespm1(j)=n0/M;pf1(j)=n2/M;end%(2)(Changing gama)-Using Bayes method to calculate pd, pf, pm and rc10=2;c01=1;c00=0;c11=0;p1=0.6;p0=0.4;gama=0.7,1,1.1,1.3,1.6;gama=(c10/c01

19、).*gama;for j=1:5n1=0;n0=0;n2=0;M=10000;for i=1:Mn=5.*randn(200,1); %To generate noisex2=s+n;x1=n;s1=x2.*s;s2=s.*s;s3=x1.*s;s1=sum(s1);s2=sum(s2);s3=sum(s3);if s10.5*s2+25*log(gama(j) n1=n1+1;else n0=n0+1;endif s30.5*s2+25*log(gama(j) n2=n2+1;endendpd2(j)=n1/M; %Using Monte Carlo method to calculate

20、 possibilitiespm2(j)=n0/M;pf2(j)=n2/M;r(j)=c00*(1-pf2(j)+c10*pf2(j)+c01*pm2(j)+c11*pd2(j);%Calculate risk indexendhomework1-3%(1)(Changing sigma)-Using the likelihood ratio(1) to calculate pd, pf, pm and rt1=0:49;t2=50:149;t3=150:199;s1=t1/50;s2=-t2/50+2;s3=t3/50-4;s=s1 s2 s3;sigma=4,16,25,36,64;for

21、 j=1:5n1=0;n0=0;n2=0;M=10000;for i=1:Mn=sqrt(sigma(j).*randn(200,1); %To generate noisex2=s+n; %To generate signals that to be testedx1=n;s1=x2.*s;s2=s.*s;s3=x1.*s;s1=sum(s1);s2=sum(s2);s3=sum(s3);if s10.5*s2 n1=n1+1;else n0=n0+1;endif s30.5*s2 n2=n2+1;endendpd1(j)=n1/M; %Using Monte Carlo method to

22、 calculate possibilitiespm1(j)=n0/M;pf1(j)=n2/M;end%(2)(Changing sigma)-Using Bayes method to calculate pd, pf, pm and rc10=2;c01=1;c00=0;c11=0;p1=0.6;p0=0.4;gama=(c10*p0)/(c01*p1);gama=log(gama);sigma=4,16,25,36,64;for j=1:5n1=0;n0=0;n2=0;M=10000;for i=1:Mn=sqrt(sigma(j).*randn(200,1); %To generate

23、 noisex2=s+n; %To generate signals that to be testedx1=n;s1=x2.*s;s2=s.*s;s3=x1.*s;s1=sum(s1);s2=sum(s2);s3=sum(s3);if s10.5*s2+sigma(j)*gama n1=n1+1;else n0=n0+1;endif s30.5*s2+sigma(j)*gama n2=n2+1;endendpd2(j)=n1/M; %Using Monte Carlo method to calculate possibilitiespm2(j)=n0/M;pf2(j)=n2/M;r(j)=

24、c00*(1-pf2(j)+c10*pf2(j)+c01*pm2(j)+c11*pd2(j); %Calculate risk indexendhomework1-4%(1)(Increasing sampling points)-Using the likelihood ratio(1) to calculate pd, pf, pm and rt1=0:0.5:49.5;t2=50:0.5:149.5;t3=150:0.5:199.5;s1=t1/50;s2=-t2/50+2;s3=t3/50-4;s=s1 s2 s3;n1=0;n0=0;n2=0;M=10000;for i=1:Mn=5

25、.*randn(400,1); %To generate noisex2=s+n; %To generate signals that to be testedx1=n;s1=x2.*s;s2=s.*s;s3=x1.*s;s1=sum(s1);s2=sum(s2);s3=sum(s3);if s10.5*s2 n1=n1+1;else n0=n0+1;endif s30.5*s2 n2=n2+1;endendpd1=n1/M; %Using Monte Carlo method to calculate possibilitiespm1=n0/M;pf1=n2/M;%(2)(Increasin

26、g sampling points)-Using Bayes method to calculate pd, pf, pm and rc10=2;c01=1;c00=0;c11=0;p1=0.6;p0=0.4;gama=(c10*p0)/(c01*p1);gama=log(gama);n1=0;n0=0;n2=0;M=10000;for i=1:Mn=5.*randn(400,1); %To generate noisex2=s+n; %To generate signals that to be testedx1=n;s1=x2.*s;s2=s.*s;s3=x1.*s;s1=sum(s1);s2=sum(s2);s3=sum(s3);if s10.5*s2+25*gama n1=n1+1