全相位fft程序

㈠ FFT的使用方法

一.調用方法X=FFT(x)；X=FFT(x，N)；
x=IFFT(X);x=IFFT(X,N)用MATLAB進行譜分析時注意：（1）函數FFT返回值的數據結構具有對稱性。例：N=8;n=0:N-1;xn=[4 3 2 6 7 8 9 0];Xk=fft(xn)→Xk =
39.0000 -10.7782 + 6.2929i 0 - 5.0000i 4.7782 - 7.7071i 5.0000 4.7782 + 7.7071i 0 + 5.0000i -10.7782 - 6.2929iXk與xn的維數相同，共有8個元素。Xk的第一個數對應於直流分量，即頻率值為0。（2）做FFT分析時，幅值大小與FFT選擇的點數有關，但不影響分析結果。在IFFT時已經做了處理。要得到真實的振幅值的大小，只要將得到的變換後結果乘以2除以N即可。二.FFT應用舉例例1：x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t)。采樣頻率fs=100Hz，分別繪制N=128、1024點幅頻圖。clf;fs=100;N=128; %采樣頻率和數據點數n=0:N-1;t=n/fs; %時間序列x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t); %信號y=fft(x,N); %對信號進行快速Fourier變換mag=abs(y); %求得Fourier變換後的振幅f=n*fs/N; %頻率序列subplot(2,2,1),plot(f,mag); %繪出隨頻率變化的振幅xlabel('頻率/Hz');ylabel('振幅');title('N=128');grid on;subplot(2,2,2),plot(f(1:N/2),mag(1:N/2)); %繪出Nyquist頻率之前隨頻率變化的振幅xlabel('頻率/Hz');ylabel('振幅');title('N=128');grid on;%對信號采樣數據為1024點的處理fs=100;N=1024;n=0:N-1;t=n/fs;x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t); %信號y=fft(x,N); %對信號進行快速Fourier變換mag=abs(y); %求取Fourier變換的振幅f=n*fs/N;subplot(2,2,3),plot(f,mag); %繪出隨頻率變化的振幅xlabel('頻率/Hz');ylabel('振幅');title('N=1024');grid on;subplot(2,2,4)plot(f(1:N/2),mag(1:N/2)); %繪出Nyquist頻率之前隨頻率變化的振幅xlabel('頻率/Hz');ylabel('振幅');title('N=1024');grid on;運行結果：
fs=100Hz，Nyquist頻率為fs/2=50Hz。整個頻譜圖是以Nyquist頻率為對稱軸的。並且可以明顯識別出信號中含有兩種頻率成分：15Hz和40Hz。由此可以知道FFT變換數據的對稱性。因此用FFT對信號做譜分析，只需考察0~Nyquist頻率范圍內的福頻特性。若沒有給出采樣頻率和采樣間隔，則分析通常對歸一化頻率0~1進行。另外，振幅的大小與所用采樣點數有關，採用128點和1024點的相同頻率的振幅是有不同的表現值，但在同一幅圖中，40Hz與15Hz振動幅值之比均為4：1，與真實振幅0.5：2是一致的。為了與真實振幅對應，需要將變換後結果乘以2除以N。例2：x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t),fs=100Hz,繪制：（1）數據個數N=32，FFT所用的采樣點數NFFT=32；（2）N=32，NFFT=128；（3）N=136，NFFT=128；（4）N=136，NFFT=512。clf;fs=100; %采樣頻率Ndata=32; %數據長度N=32; %FFT的數據長度n=0:Ndata-1;t=n/fs; %數據對應的時間序列x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t); %時間域信號y=fft(x,N); %信號的Fourier變換mag=abs(y); %求取振幅f=(0:N-1)*fs/N; %真實頻率subplot(2,2,1),plot(f(1:N/2),mag(1:N/2)*2/N); %繪出Nyquist頻率之前的振幅xlabel('頻率/Hz');ylabel('振幅');title('Ndata=32 Nfft=32');grid on;
Ndata=32; %數據個數N=128; %FFT採用的數據長度n=0:Ndata-1;t=n/fs; %時間序列x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t);y=fft(x,N);mag=abs(y);f=(0:N-1)*fs/N; %真實頻率subplot(2,2,2),plot(f(1:N/2),mag(1:N/2)*2/N); %繪出Nyquist頻率之前的振幅xlabel('頻率/Hz');ylabel('振幅');title('Ndata=32 Nfft=128');grid on;
Ndata=136; %數據個數N=128; %FFT採用的數據個數n=0:Ndata-1;t=n/fs; %時間序列x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t);y=fft(x,N);mag=abs(y);f=(0:N-1)*fs/N; %真實頻率subplot(2,2,3),plot(f(1:N/2),mag(1:N/2)*2/N); %繪出Nyquist頻率之前的振幅xlabel('頻率/Hz');ylabel('振幅');title('Ndata=136 Nfft=128');grid on;
Ndata=136; %數據個數N=512; %FFT所用的數據個數n=0:Ndata-1;t=n/fs; %時間序列x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t);y=fft(x,N);mag=abs(y);f=(0:N-1)*fs/N; %真實頻率subplot(2,2,4),plot(f(1:N/2),mag(1:N/2)*2/N); %繪出Nyquist頻率之前的振幅xlabel('頻率/Hz');ylabel('振幅');title('Ndata=136 Nfft=512');grid on;
結論：（1）當數據個數和FFT採用的數據個數均為32時，頻率解析度較低，但沒有由於添零而導致的其他頻率成分。（2）由於在時間域內信號加零，致使振幅譜中出現很多其他成分，這是加零造成的。其振幅由於加了多個零而明顯減小。（3）FFT程序將數據截斷，這時解析度較高。（4）也是在數據的末尾補零，但由於含有信號的數據個數足夠多，FFT振幅譜也基本不受影響。對信號進行頻譜分析時，數據樣本應有足夠的長度，一般FFT程序中所用數據點數與原含有信號數據點數相同，這樣的頻譜圖具有較高的質量，可減小因補零或截斷而產生的影響。
例3：x=cos(2*pi*0.24*n)+cos(2*pi*0.26*n)
（1）數據點過少，幾乎無法看出有關信號頻譜的詳細信息；（2）中間的圖是將x(n)補90個零，幅度頻譜的數據相當密，稱為高密度頻譜圖。但從圖中很難看出信號的頻譜成分。（3）信號的有效數據很長，可以清楚地看出信號的頻率成分，一個是0.24Hz，一個是0.26Hz，稱為高解析度頻譜。可見，采樣數據過少，運用FFT變換不能分辨出其中的頻率成分。添加零後可增加頻譜中的數據個數，譜的密度增高了，但仍不能分辨其中的頻率成分，即譜的解析度沒有提高。只有數據點數足夠多時才能分辨其中的頻率成分。

㈡ 求用C語言實現FFT變換的程序（見下面）

你好，這是我的回答，希望可以幫到你。

1）結果討論

一，如果對信號進行同樣點數N的FFT變換，采樣頻率fs越高，則可以分析越高頻的信號；與此同時，采樣頻率越低，對於低頻信號的頻譜解析度則越好。

二，假設采樣點不在正弦信號的波峰、波谷、以及0電壓處，頻譜則會產生泄露（leakage）。

三，對於同樣的采樣率fs，提高FFT的點數N，則可提高頻譜的解析度。

四，如果采樣頻率fs小於2倍信號頻率2*fs（奈圭斯特定理），則頻譜分析結果會出錯。

五，對於（二）中泄露現象，可以通過在信號後面補零點解決。

2）程序及註解如下

%清除命令窗口及變數
clc;
clear all;

%輸入f、N、T、是否補零（補幾個零）
f=input('Input frequency of the signal: f\n');
N=input('Input number of pointsl: N\n');
T=input('Input sampling time: T\n');
flag=input('Add zero too sampling signal or not? yes=1 no=0\n');
if(flag)
ZeroNum=input('Input nmber of zeros\n');
else
ZeroNum=0;
end

%生成信號，signal是原信號。signal為采樣信號。
fs=1/T;
t=0:0.00001:T*(N+ZeroNum-1);
signal=sin(2*pi*f*t);
t2=0:T:T*(N+ZeroNum-1);
signal2=sin(2*pi*f*t2);
if (flag)
signal2=[signal2 zeros(1, ZeroNum)];
end

%畫出原信號及采樣信號。
figure;
subplot(2,1,1);
plot(t,signal);
xlabel('Time(s)');
ylabel('Amplitude(volt)');
title('Singnal');
hold on;
subplot(2,1,1);
stem(t2,signal2,'r');
axis([0 T*(N+ZeroNum) -1 1]);

%作FFT變換，計算其幅值，歸一化處理，並畫出頻譜。
Y = fft(signal2,N);
Pyy = Y.* conj(Y) ;
Pyy=(Pyy/sum(Pyy))*2;
f=0:fs/(N-1):fs/2;4
subplot(2,1,2);
bar(f,Pyy(1:N/2));
xlabel('Frequency(Hz)');
ylabel('Amplitude');
title('Frequency compnents of signal');
axis([0 fs/2 0 ceil(max(Pyy))])
grid on;

祝你好運！
我可以幫助你，你先設置我最佳答案後，我網路Hii教你。

㈢ 用matlab編寫實現fft的程序。

function y=myditfft(x)
%本程序對輸入序列實現DIT-FFT基2演算法，點數取大於等於長度的2的冪次
%------------------------------------
% myditfft.c
%------------------------------------
m=nextpow2(x); %求的x長度對應的2的最低冪次m
N=2^m;
if length(x)<N
x=[x,zeros(1,N-length(x))]; %若的長度不是2的冪，補0到2的整數冪
end
nxd=bin2dec(fliplr(dec2bin([1:N]-1,m)))+1; %求1：2^m數列的倒序
y=x(nxd); %將倒序排列作為的初始值
for mm=1:m %將DFT做m次基2分解，從左到右，對每次分解作DFT運算
Nmr=2^mm;
u=1; %旋轉因子u初始化
WN=exp(-i*2*pi/Nmr); %本次分解的基本DFT因子WN＝exp(-i*2*pi/Nmr)
for j=1:Nmr/2 %本次跨越間隔內的各次碟形運算
for k=j:Nmr:N %本次碟形運算的跨越間隔為Nmr=2^mm
kp=k+Nmr/2; %確定碟形運算的對應單元下標
t=y(kp)*u; %碟形運算的乘積項
y(kp)=y(k)-t; %碟形運算的加法項
y(k)=y(k)+t;
end
u=u*WN; %修改旋轉因子，多乘一個基本DFT因子WN
end
end

㈣ 如何使用fft函數進行編程序和進行快速傅里葉逆變換

在圖象處理的廣泛應用領域中，傅立葉變換起著非常重要的作用，具體表現在包括圖象分析、圖象增強及圖象壓縮等方面。
fftshift的作用正是讓正半軸部分和負半軸部分的圖像分別關於各自的中心對稱。因為直接用fft得出的數據與頻率不是對應的，fftshift可以糾正過來。
假設f（x,y）是一個離散空間中的二維函數，則該函數的二維傅立葉變換的定義如下：

p=0,1…M-1 q=0,1…N-1 (1)

或 p=0,1…M-1 q=0,1…N-1 （2）

離散傅立葉反變換的定義如下：

m=0,1…M-1 n=0,1…N-1（3）
F（p,q）稱為f（m,n）的離散傅立葉變換系數。這個式子表明，函數f（m,n）可以用無數個不同頻率的復指數信號和表示，而在頻率（w1，w2）處的復指數信號的幅度和相位是F（w1，w2）。
2、MATLAB提供的快速傅立葉變換函數
（1）fft2
fft2函數用於計算二維快速傅立葉變換，其語法格式為：
B = fft2(I)
B = fft2(I)返回圖象I的二維fft變換矩陣，輸入圖象I和輸出圖象B大小相同。
例如，計算圖象的二維傅立葉變換，並顯示其幅值的結果，如圖所示，其命令格式如下
load imdemos saturn2
imshow(saturn2)
B = fftshift(fft2(saturn2));
imshow(log(abs(B)),[],'notruesize')
（2）fftshift
MATLAB提供的fftshift函數用於將變換後的圖象頻譜中心從矩陣的原點移到矩陣的中心，其語法格式為：
B = fftshift(I)
對於矩陣I，B = fftshift(I)將I的一、三象限和二、四象限進行互換。
（2）ifft2
ifft2函數用於計算圖象的二維傅立葉反變換，其語法格式為：
B = ifft2(I)
B = ifft2(A)返回圖象I的二維傅立葉反變換矩陣，輸入圖象I和輸出圖象B大小相同。其語法格式含義與fft2函數的語法格式相同，可以參考fft2函數的說明。
如果信號是二維的，用上面的函數即可！直接調用。
如果信號是一維的，給你下面的例子，你應該能明白！
clear
fs=100;N=128; %采樣頻率和數據點數
n=0:N-1;t=n/fs; %時間序列
x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t); %信號
y=fft(x,N); %對信號進行快速Fourier變換，逆變換函數為ifft
mag=abs(y); %求得Fourier變換後的振幅
f=n*fs/N; %頻率序列
subplot(2,2,1),plot(f,mag); %繪出隨頻率變化的振幅
xlabel('頻率/Hz');
ylabel('振幅');title('N=128');grid on;
subplot(2,2,2),plot(f(1:N/2),mag(1:N/2)); %繪出Nyquist頻率之前隨頻率變化的振幅
xlabel('頻率/Hz');
ylabel('振幅');title('N=128');grid on;
%對信號采樣數據為1024點的處理
fs=100;N=1024;n=0:N-1;t=n/fs;
x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t); %信號
y=fft(x,N); %對信號進行快速Fourier變換
mag=abs(y); %求取Fourier變換的振幅
f=n*fs/N;
subplot(2,2,3),plot(f,mag); %繪出隨頻率變化的振幅
xlabel('頻率/Hz');
ylabel('振幅');title('N=1024');grid on;
subplot(2,2,4)
plot(f(1:N/2),mag(1:N/2)); %繪出Nyquist頻率之前隨頻率變化的振幅
xlabel('頻率/Hz');
ylabel('振幅');title('N=1024');grid on;