隨機信號估計演算法

A. 鏁板瓧淇″彿澶勭悊涓濡備綍鐢婚浂鐐瑰拰鏋佺偣鍒嗗竷鍥撅紵

鍏堢敾涓涓澶鍧愭爣緋錛岀劧鍚庢眰鍑浼犻掑嚱鏁G (x)鐨勯浂鐐瑰拰鏋佺偣錛屾爣鍦ㄥ潗鏍囩郴涓鍗沖彲錛岄浂鐐逛負鍒嗗瓙涓洪浂鐨勭偣錛屾瀬鐐逛負鍒嗘瘝涓洪浂鐨勭偣錛

鍦棰戝煙涓鎻忚堪淇″彿鐗規х殑涓縐嶅垎鏋愭柟娉曪紝涓嶄粎鍙鐢ㄤ簬紜瀹氭т俊鍙鳳紝涔熷彲鐢ㄤ簬闅忔満鎬т俊鍙楓傛墍璋撶『瀹氭т俊鍙峰彲鐢ㄦ棦瀹氱殑鏃墮棿鍑芥暟鏉ヨ〃紺猴紝瀹冨湪浠諱綍鏃跺埢鐨勫兼槸紜瀹氱殑錛涢殢鏈轟俊鍙峰垯涓嶅叿鏈夎繖鏍風殑鐗規э紝瀹冨湪鏌愪竴鏃跺埢鐨勫兼槸闅忔満鐨勩

鍥犳わ紝闅忔満淇″彿澶勭悊鍙鑳芥牴鎹闅忔満榪囩▼鐞嗚猴紝鍒╃敤緇熻℃柟娉曟潵榪涜屽垎鏋愬拰澶勭悊錛屽傜粡甯稿埄鐢ㄥ潎鍊箋佸潎鏂瑰箋佹柟宸銆佺浉鍏沖嚱鏁般鍔熺巼璋瀵嗗害鍑芥暟絳夌粺璁￠噺鏉ユ弿榪伴殢鏈鴻繃紼嬬殑鐗瑰緛鎴栭殢鏈轟俊鍙風殑鐗規с

鎵╁睍璧勬枡錛

鍦ㄥ疄闄呬腑瑙傛祴鍒扮殑鏁版嵁鏄鏈夐檺鐨勩傝繖灝遍渶瑕佸埄鐢ㄤ竴浜涗及璁＄殑鏂規硶錛屾牴鎹鏈夐檺鐨勫疄嫻嬫暟鎹浼拌″嚭鏁翠釜淇″彿鐨勫姛鐜囪氨銆傞拡瀵逛笉鍚岀殑瑕佹眰錛屽傚噺灝忚氨鍒嗘瀽鐨勫嚲宸錛屽噺灝忓瑰櫔澹扮殑鐏墊晱紼嬪害錛屾彁楂樿氨鍒嗚鯨鐜囩瓑銆傚凡鎻愬嚭璁稿氫笉鍚岀殑璋變及璁℃柟娉曘

鍦ㄧ嚎鎬т及璁℃柟娉曚腑錛屾湁鍛ㄦ湡鍥炬硶錛岀浉鍏蟲硶鍜屽崗鏂瑰樊娉曪紱鍦ㄩ潪綰挎т及璁℃柟娉曚腑,鏈鏈澶т技鐒舵硶,鏈澶х喌娉曪紝鑷鍥炲綊婊戝姩騫沖潎淇″彿妯″瀷娉曠瓑銆傝氨鍒嗘瀽鍜岃氨浼拌′粛鍦ㄧ爺絀跺拰鍙戝睍涓銆

鏁板瓧淇″彿澶勭悊涓庢ā鎷熶俊鍙峰勭悊閮芥槸淇″彿澶勭悊鐨勫瓙闆嗭紝鑼冪暣鍧囧睘浜庝俊鍙峰勭悊錛屾墍璋"淇″彿澶勭悊"錛屽氨鏄瑕佹妸璁板綍鍦ㄦ煇縐嶅獟浣撲笂鐨勪俊鍙瘋繘琛屽勭悊錛屼互渚挎娊鍙栧嚭鏈夌敤淇℃伅鐨勮繃紼嬨

瀹冩槸瀵逛俊鍙瘋繘琛屾彁鍙栥佸彉鎹銆佸垎鏋愩佺患鍚堢瓑澶勭悊榪囩▼鐨勭粺縐般備絾鏁板瓧淇″彿澶勭悊浠ュ強妯℃嫙淇″彿澶勭悊鎵澶勭悊鐨勫硅薄涓嶄竴鑷達紝鍥犳ゅ勭悊鐨勫叿浣撴祦紼嬩篃涓嶅敖鐩稿悓錛屼絾鐩鐨勯兘鏄涓轟簡鎻愬彇鍑烘湁鐢ㄧ殑淇℃伅銆

鍙傝冭祫鏂欐潵婧愶細鐧懼害鐧劇--鏁板瓧淇″彿澶勭悊

B. 信號功率的求法

能量有限平均功率為零的信號是能量信號，如單位沖擊信號；
能量無限，平均功率有限的信號是功率信號。
題中信號有一個常量10，平均功率不為零且有限，能量無限，屬於"功率信號「
功率為：10^2+4^2+2^2=120

C. 運用matlab進行隨機信號的功率譜密度估計模擬【急求】

好眼熟啊。。。通院08級的隨機信號實驗題目？
這份是我在網路上找的：
（1）.周期圖法： 思想：周期圖法是為了得到功率譜估值，先取信號序列的離散傅里葉變換，然後取其幅頻特性的平方並除以序列長度 N。由於序列 x(n)的離散傅里葉變換 X(k) 具有周期性，因而這種功率譜也具有周期性，常稱為周期圖。周期圖是信號功率譜的一個有偏估值；而且，當信號序列的長度增大到無窮時，估值的方差不趨於零。因此，隨著所取的信號序列長度的不同，所得到的周期圖也不同，這種現象稱為隨機起伏。由於隨機起伏大，使用周期圖不能得到比較穩定的估值。

程序：
%首先，生成輸入信號的程序為：
clear ; fs=20000; n=0:1/ fs: 0. 1; N=lengt h( n); W=2000*pi;%因方波頻率 F=1000HZ 所以角頻率 W=2000pi X1n=square( W*n) ;%方波信號 X2n=randn( 1, N);%白雜訊信號 xn=X1n+0. 2*X2n; subplot (3, 1, 1) plot (n, xn) ; xlabel( ' n' ) ylabel('輸入信號') %其次，開始用周期圖法進行估計； clear all fs=20000; n=0:1/fs:0.1; N=length(n); W=2000*pi; x1n=square(W*n); x2n=randn(1,N); xn=x1n+0.2*x2n; subplot(2,1,1) plot(n,xn); Nfft=256;N=256;%傅里葉變換的采樣點數256 Pxx=abs(fft(xn,Nfft).^2)/N; f=(0:length(Pxx)-1)*fs/length(Pxx); subplot(2,1,2), plot(f,10*log10(Pxx)),%轉成DB 單位

（2）.自相關函數法： 思想： 隨機信號 x(n)的相關函數是在時間域內描述隨機過程的重要特徵。自相關函數是隨機信號在不同時刻的值之間的依賴性的量度，是一個很有用的統計平均量，其定義為自相關函數 (1) 式中E【·】表示數學期望，*表示共軛值,m 為時間滯後數。在隨機信號處理中，自相關函數可以用來檢測淹沒在隨機雜訊干擾中的信號，隨機信號的自功率譜等於它的自相關函數的傅里葉變換。因此，通過自相關估計可求得信號的功率譜。利用計算機計算自相關估值有兩種方法。一種是直接方法，先計算出隨機信號和它的滯後序列的乘積，再取其平均值即得相關函數的估計值。另一種是間接方法，先用快速變換演算法計算隨機序列的功率譜密度，再作反變換計算出相關函數。

程序： n=0:1/fs:1; N=length(n); W=2000*pi; x1n=square(W*n); x2n=randn(1,N); xn=x1n+0.2*x2n; figure(2) subplot(3,1,1) plot(n,xn);%輸入信號 axis([0 0.01 -5 5]) m =-100:100 [r,lag]=xcorr(xn,100,'biased')%求XN 的自相關函數R，biased 為有偏估計,lag 為R 的序列號 subplot(3,1,2) hndl=stem(m,r);%繪制離散圖，分布點從-100—+100 set(hndl,'Marker','.') set(hndl,'MarkerSize',2); ylabel('自相關函數R(m)') %利用間接法計算功率譜 k=0:1000;%取1000 個點 w=(pi/500)*k; M=k/500; X=r*(exp(-li*pi/500).^(m'*k));%對R 求傅里葉變換,li 就是j magX=abs(X); subplot(3,1,3) plot(M,10*log10(magX)); xlabel('功率譜的改進直接法估計')

（3）.自協方差法： 思想：在實際中,有時用隨機信號在兩個不同時刻t 1 ,t 2 的取值X(t 1 )和X(t 2 )之間的二階混合中心矩來描述隨機信號 X(t)在任意兩個時刻取值起伏變化的相依程度, 即自協方差。自協方差函數與自相關函數描述的特性是一致的，所以其原理與相關函數法近似。

clear all; fs=20000; n=0:1/fs:0.1; P=2000*pi; y=square(P*n); xn=y+0.2*randn(size(n));%繪制信號波形 figure(3) subplot(2,1,1) plot(n,xn) xlabel('時間（s）') ylabel('幅度') title('y+randn(size(n))') ymax_xn=max(xn)+0.2; ymin_xn=min(xn)-0.2; axis([0 0.3 ymin_xn ymax_xn]) %使用協方差法估計序列功率譜 p=floor(length(xn)/3)+1; nfft=1024; [xpsd,f]=pcov(xn,p,nfft,fs,'half'); %繪制功率譜估計 pmax=max(xpsd); xpsd=xpsd/pmax; xpsd=10*log10(xpsd+0.000001); subplot(2,1,2) plot(f,xpsd) title('基於協方差的功率譜估計') ylabel('功率譜估計（db）') xlabel('頻率') grid on; ymin=min(xpsd)-2; ymax=max(xpsd)+2; axis([0 fs/2 ymin ymax])

（4）.最大熵法 思想：上網查的最大熵法原理是取一組時間序列，使其自相關函數與一組已知數據的自相關函數相同，同時使已知自相關函數以外的部分的隨機性最強，以所取時間序列的譜作為已知數據的譜估值。它等效於根據使隨機過程的熵為最大的原則,利用 N 個已知的自相關函數值來外推其他未知的自相關函數值所得到的功率譜。最大熵法功率譜估值是一種可獲得高解析度的非線性譜估值方法，特別適用於數據長度較短的情況。最大熵法譜估值對未知數據的假定 ，一個平穩的隨機序列，可以用周期圖法對其功率譜進行估值。這種估值方法隱含著假定未知數據是已知數據的周期性重復。現有的線性譜估計方法是假定未知數據的自相關函數值為零，這種人為假定帶來的誤差較大。最大熵法是利用已知的自相關函數值來外推未知的自相關函數值，去除了對未知數據的人為假定，從而使譜估計的結果更為合理。熵在資訊理論中是信息的度量，事件越不確定，其信息量越大，熵也越大。對於上述問題來說，對隨機過程的未知的自相關函數值，除了從已知的自相關函數值得到有關它的信息以外，沒有其他的先驗知識。因而，在外推時，不希望加以其他任何新的限制，亦即使之「最不確定」。換言之，就是使隨機過程的熵最大。

程序： fs=20000; n=0:1/fs:1; N=length(n); W=2000*pi; x1n=square(W*n); x2n=randn(1,N); xn=x1n+0.2*x2n; figure(5) subplot(3,1,1) plot(n,xn); asis([0 0.01 -5 5]) Nfft=256;%分段長度256 [Pxx,f]=pmem(xn,14,Nfft,fs);%調用最大熵函數pmem，濾波器階數14 subplot(2,1,2), plot(f,10*log10(Pxx)), title(' 最大熵法，濾波器14'),xlabel('頻率'),ylabel('功率譜db');

（5）.最大似然法： 思想：最大似然法原理是讓信號通過一個濾波器，選擇濾波器的參數使所關心的頻率的正弦波信號能夠不失真地通過，同時，使所有其他頻率的正弦波通過這個濾波器後輸出的均方值最小。在這個條件下，信號經過這個濾波器後輸出的均方值就作為其最大似然法功率譜估值。可以證明,如果信號x 是由一個確定性信號S 加上一個高斯白雜訊n 所組成，則上述濾波器的輸出是信號S 的最大似然估值。如果n 不是高斯雜訊，則上述濾波器的輸出是信號S 的最小方差的線性的無偏估值，而且它能得到比使用固定的窗口函數的周期圖法更高的解析度。

程序： fs=20000; n=0:1/fs:1; N=length(n); W=2000*pi; x1n=square(W*n); x2n=randn(1,N); xn=x1n+0.2*x2n; figure(4) subplot(3,1,1) plot(n,xn); axis([0 0.01 -5 5]) %估計自相關函數 m=-500:500; [r,lag]=xcorr(xn,500,'biased'); R=[r(501) r(502) r(503) r(504); r(500) r(501) r(502) r(503); r(499) r(500) r(501) r(502); r(498) r(499) r(500) r(501)]; [V,D]=eig(R); V3=[V(1,3),V(2,3),V(3,3),V(4,3)].'; V3=[V(1,4),V(2,4),V(3,4),V(4,4)].'; p=0:3; wm=[0:0.002*pi:2*pi]; B=[(exp(-li)).^(wm'*p)];%li 就是虛數單位j A=B; %最小方差功率譜估計 z=A*inv(R)*A'; Z=diag(z'); pmv=1./Z; subplot(2,1,2) plot(wm/pi,pmv); title('基於最大似然的功率譜估計') ylabel('功率譜幅度（db）') xlabel('角度頻率w/pi')
自己分下段就能運行了。