這次博文主要內(nèi)容是對離散傅里葉變換和希爾伯特變換進行簡單介紹。在MATLAB中仿真出一些例子，讓我們學(xué)習(xí)時對理論知識有一個形象的理性認識，這樣理解起來比較深刻。

    內(nèi)容概要：利用離散傅里葉變換將加噪的調(diào)制信號變換到頻域，用去除高頻的高斯白噪聲干擾的方法進行降噪。然后利用希爾伯特變換求得調(diào)制信號的解析信號，根據(jù)解析式得到調(diào)制信號的瞬時參數(shù)：瞬時幅度，瞬時頻率，瞬時相位。最后用自相關(guān)函數(shù)檢測調(diào)制信號的碼元速率。

一．理論分析（理論部分這里博主略去，很多公式不方便錄入，學(xué)習(xí)者自己找，很多）

1.離散傅里葉變換（Discrete Fourier Transform，縮寫為DFT），是傅里葉變換在時域和頻域上都呈離散的形式，將信號的時域采樣變換為其DTFT的頻域采樣。在形式上，變換兩端（時域和頻域上）的序列是有限長的，而實際上這兩組序列都應(yīng)當被認為是離散周期信號的主值序列。即使對有限長的離散信號作DFT，也應(yīng)當將其看作其周期延拓的變換。在實際應(yīng)用中通常采用快速傅里葉變換計算DFT。

2.希爾伯特變換

3.自相關(guān)函數(shù)

二．針對問題和方法

1.信號降噪：現(xiàn)實中的信號一般都受到噪聲的干擾，而導(dǎo)致在提取信號特征信息的造成誤差，所以在處理信號之前先進行降噪。一個含噪聲的一維信號模型可表示為如下形式：

S(k)=f(k)+e(k)

S(k)為含噪信號，f(k)是有用信號，e(k)是噪聲信號。這里假定噪聲是高斯白噪聲，頻譜一般分布為整個頻域。

20db信噪比4psk的頻譜

實際工程中f(k)通常為低頻信號或者頻譜范圍分布有限的信號。因此，通過離散傅里葉變換得到信號的頻譜以后，可以把有用的頻譜之外的噪聲頻譜去除，然后再經(jīng)過離散傅里葉逆變換就可以得到降噪的信號。

20db信噪比4psk去除一段噪聲的頻譜

2.提取調(diào)制信號瞬時參數(shù)：利用信號的解析函數(shù)的特性計算調(diào)制信號的特征。

用原信號和其希爾伯特變換表示的復(fù)函數(shù)稱作這個原實信號的解析函數(shù)

z(tn)=u(n)+jv(n)

采用頻域的方法計算采樣信號u(n)的希爾伯特變換v(n)。

3.利用自相關(guān)函數(shù)求調(diào)制信號碼元速率

         調(diào)制信號在相位突變的時候，由解析函數(shù)求得的瞬時頻率在碼元跳變處會頻率突變。利用這一點來檢測調(diào)制信號的碼元速率。先把瞬時頻率信號減去瞬時頻率的平均值，使得頻率中心盡量歸零，接著對處理后的瞬時頻率信號進行自相關(guān)，自相關(guān)后的信號每隔一段采樣點就會出現(xiàn)突變，這些突變的最小間隔時間就是碼元速率Rs。

三.MATLAB仿真

調(diào)制信號的瞬時參數(shù)在很多方面都有重要應(yīng)用，下面針對數(shù)字相位調(diào)制信號4psk（加高斯白噪聲），首先進行降噪，然后根據(jù)解析信號的求得瞬時參數(shù)。仿真如下：

         仿真信號s(n)采用4psk離散調(diào)制信號，載波選擇fs=1kHz，碼元速率fd=100Hz，采樣頻率選擇fs=10kHz，疊加信噪比20db的高斯白噪聲。

信噪比為20db的4psk信號（未降噪）瞬時參數(shù)提取結(jié)果：

信噪比為20db的4psk信號（降噪后）瞬時參數(shù)提取結(jié)果：

由這兩幅仿真圖看出，信號降噪前后有明顯的效果。由解析信號提取的瞬時幅度在1附近變化，瞬時頻率在1kHz附近變化，在相位突變處瞬時頻率發(fā)生劇烈變化。調(diào)制信號相位對應(yīng)的碼元xi=[0 1 2 3 3 2 1 0 2 3 0 1]

瞬時頻率信號歸零后的自相關(guān)的絕對值：

   由圖可得，基本每隔Nt=100個采樣點就會出現(xiàn)一個自相關(guān)函數(shù)值的劇烈增大，而信號的采樣頻率fs=10000Hz，所以碼元速率Rs=Nt/fs=0.01b/s.調(diào)制的碼元頻率fd=100Hz，理論上碼元速率是1/fd,即0.01b/s。碼元速率的檢測值和理論值接近。

主m程序

clear all

M=4;

fc=1000;

fd=100;

fs=10000;

xi=[0 1 2 3 3 2 1 0 2 3 0 1];

[u,t]=dmod(xi,fc,fd,fs,'psk',M);

sn0=awgn(u,20,'measured','db');%?ó??éù

sn=fftlvbo(sn0);

z=hilbert(sn);

fai1=angle(z);

p=FXXXWTQ3(z,fc,fs);

fi=f(z,fc,fs);

z=abs(z);

figure(1)

subplot(2,2,1);

plot(t,sn);

title('?ó??μ÷??D?o?sn=f(t)')

subplot(2,2,2);

plot(t,z);

title('μ÷??D?o??2ê±·ù?è|z|')

subplot(2,2,3);

plot(t,fi);

title('μ÷??D?o??2ê±?μ?êfi')

subplot(2,2,4);

dt=1/fs;

Np=length(p);

tp=[1:Np]*dt;

plot(tp,p);

title('μ÷??D?o??2ê±?à??p￡¨?à??è￥??μt￡?')

fi1=fi-mean(fi);

r = xcorr(fi1);

figure(3)

plot(abs(r))

降噪m程序

function u=fftlvbo(sn)

f=fft(sn);

L=length(f);

for i=270:930

    if(abs(f(i))<20)

        f(i)=0;

    end

end

u=ifft(f);

figure(2)

plot(abs(f))

解析信號構(gòu)建m文件

U=fft(u,N);%à?ó??￡??2?ì?±???￡?′ó?μóò±???μ?μ?v￡¨t￡?

k=1:N/2;%V(f)=-j*sgn(f)*U(f)

V(k)=-j*U(k);

k=(N+2)/2:N;

V(k)=j*U(k);

v=ifft(V,N);%?μà?ò???±???μ?μ?v￡¨t￡?

z=u+j*v;%?a??D?o?