在這裡首先確定dft的物件為乙個有限長的離散非週期序列,這主要因為計算機處理的都是有限長的離散序列。如果你要處理的序列本身不是離散非週期的序列,可以通過擷取或者離散化等方法獲得所需的有限長的離散非週期序列。
對於有限長的離散非週期序列存在兩種計算n點dft的方法,一種方法是先將其通過dtft獲得序列的連續週期的頻譜,取其中乙個週期進行n點抽樣即可得到dft;另一種方法直接由 x(
n)計算n點的dft,但是這種方法要求 dft 變換的點數 n 大於等於原始時序序列的點數。下面將給出乙個具體的例子進行說明。
題目x
(n)=
r16(n
) ,求
n 分別取8,16,32,64時的離散傅利葉變換dft x(
k),並利用反變換求各個頻譜對應的時序序列,比較這些頻譜和序列。
第一種方法:直接計算
這種方法的原理在於首先將時序序列進行以 n 為週期的週期延拓,然後取出其中的乙個週期計算在 k=
0:n−
1 上的dtft(注意原始的是計算k∈
(−∞,
+∞) )。要如果計算dft的點數大於等於時序序列本身的點數,那麼不通過週期延拓直接計算就好,因為週期延拓之後後面的都是補上的0,沒有計算價值,因為計算完也是0。這時的計算公式如下 x(
k)=σ
n−1n
=0x(
n)e−
j2πn
nk=σ
n−1n
=0x(
n)wk
nn,k
∈[0,
n−1]
但是!這裡可以直接使用 x(
n)進行計算是因為 x(
n)以 n 點為週期進行週期延拓後面補的都是0,算不算這些點都沒有什麼關係,但是像這個問題,本身 x(
n)是乙個16點的矩形序列,如果做 8 點的dft時首先要以 8 為週期進行週期延拓,延拓後序列在時域上有混疊,所以不能直接將 x(
n)代入上式進行計算,應該帶入週期延拓之後時序序列的進行計算。但是這樣做明顯太麻煩了,這時就需要第二種方法了。
第二種方法:由 dtft 匯出 dft
這種方法是首先計算序列的 dtft,再取乙個週期進行 n 點抽樣。這種方法是符合dft原始定義的方法的一種匯出形式,無論是計算多少點數的dft變換均可以,而且在變換的過程中也不需要進行補零操作。
本**執行環境為 matlab 2014a
具體的計算思想可以參見從 dtft 的角度用矩陣相乘代替 for 迴圈進行計算。
clc;clear all;close all;
%% 繪製x(n)
x = ones(1,16);n=0:15;
subplot(5,2,[1,2]),stem(x);
%% x(n)的dtft
w=0:2*pi/3200:2*pi;
xdtft0=x*exp(-1
i*n'*w);
subplot(5,2,3),plot(w/pi,abs(xdtft0));
%% x(exp(jw))的idtft
syms ww
xdtft=x*exp(-1
i*n'*ww);
xidtft = int(xdtft*exp(1
i*ww'*n),ww,-pi,pi)/2/pi
subplot(5,2,4),stem(xidtft);
%% x(n) 的16點,32點,64點的 dft & idft
nn = [16,32,64];
for n = nn
nn = find(nn==n)
nnn=0:n-1;
k=0:n-1;
xn = [x,zeros(1,n- 16)];
xk=xn*exp(-1
j*2*pi/n).^(nnn'*k);
www = 0:2*pi/n:2*pi;
subplot(5,2,3+2*nn),stem(www(1,1:n)/pi,abs(xk));
xx=(xk*exp(j*2*pi/n).^(nnn'*k))/n;
subplot(5,2,4+2*nn),stem(abs(xx));
end
實驗結果如下圖所示
從圖中我們可以清楚的看到,dft的點數就是對原序列dtft乙個週期內頻譜的取樣點數。從頻率取樣定理可以知道,當頻域取樣點數大於等於時域實際序列的點數時,可完整恢復序列。但是這裡需要注意,這裡沒有計算序列的8點dft,該部分將在第二種方法中給出。
clc;clear all;close all;
%% 繪製x(n)
x = ones(1,16);n=0:15;
subplot(6,2,[1,2]),stem(x);
%% x(n)的dtft
w=0:2*pi/3200:2*pi;
xdtft0=x*exp(-1
i*n'*w);
subplot(6,2,3),plot(w/pi,abs(xdtft0));
%% x(exp(jw))的idtft
syms ww
xdtft=x*exp(-1
i*n'*ww);
xidtft = int(xdtft*exp(1
i*ww'*n),ww,-pi,pi)/2/pi
subplot(6,2,4),stem(xidtft);
%% x(n) 的16點,32點,64點的dft & idft
nn = [8,16,32,64]
for n = nn
hang = find(nn==n)
count = 1;
nn=0:n-1;
k=0:n-1;
for h = 1:round(3200/n):3200;
dft(hang,count) = xdtft0(h);
count = count+1;
endwww = 0:2*pi/n:2*pi;
subplot(6,2,3+2*hang),stem(www(1,1:n)/pi,abs(dft(hang,:)));
xx=(dft(hang,:)*exp(j*2*pi/n).^(nn'*k))/n;
subplot(6,2,4+2*hang),stem(abs(xx));
end
實驗結果如下圖所示
觀察上圖,圖中的第一行與第二行與圖1中的第一行、第二行影象相同;第三行為通過dtft的方式匯出的8點dft的幅頻特性和由8點dft反變換得到的時序序列,我們可以看到時序序列的幅值為2,而不是1,且乙個週期內只有8個點。產生這樣現象的原因就是之前說過的,頻域取樣的個數(即dft變換的點數)小於時序序列的點數造成的。
從上面可以看到,兩種方法中第一種方法(直接計算)適用於dft點數多餘時序序列點數的情況,而第二種計算方法(由 dtft 匯出 dft)可以不限制dft變換的點數。從理論上講,第二種方法似乎更好,但是考慮在實際使用中,dft的點數如果小於時序序列的點數,這樣的變換本身也是沒有意義的,所以一般也不會計算這樣的dft,所以了兩種方法均可以使用。但是為了保險,推薦使用第二種~
離散傅利葉變換的兩種實現方案
方案一 dftm sig ones 256,1 t 8 timeslice 8 timeslice f t nfft 256 sig cos 2 pi 1 t dft zeros nfft,nfft n 0 nfft 1 n的行向量,為1 n矩陣 k 0 nfft 1 k的行向量,為1 n矩陣 wn...
驗證離散傅利葉變換(DFT)的時域迴圈移位定理
離散傅利葉變換時域迴圈移位定理如下 接下來通過matlab程式設計進行驗證,為了便於給出任意移位點數與任意序列的驗證,此處將驗證 封裝為乙個dft time circshift函式,如下 function xk,yk,yk1 dft time circshift xn,m dft time circ...
離散化的兩種操作
離散化是程式設計中乙個常用的技巧,它可以有效的降低時間複雜度。其基本思想就是在眾多可能的情況中,只考慮需要用的值。離散化可以改進乙個低效的演算法,甚至實現根本不可能實現的演算法。要掌握這個思想,必須從大量的題目中理解此方法的特點。有些資料本身很大,自身無法作為陣列的下標儲存對應的屬性。如果這時只是需...