FFT

2021-01-14 創客小巴
最近在項目中需要用到FFT，之前對於FFT也只是有一個模糊的印象也並不清楚他的具體物理意義，之前幾次想學習都被擱置了，現在項目需要又從新學習，在此把我收穫的和大家分享一下：


1- FFT簡介

FFT是一種DFT的高效算法，稱為快速傅立葉變換（fast Fourier transform）。傅立葉變換是時域--頻域變換分析中最基本的方法之一。可以將一個信號變換到頻域。有些信號在時域上很難看出什麼特徵，不利於分析，但是如果變換到頻域之後，就很容易看出信號的特徵了。這就是很多信號分析採用FFT變換的原因。FFT也可以將一個信號的頻譜提取出來，常應用於頻譜分析上。

2 -採樣定理

採樣頻率必須是信號的最高頻率的兩倍及其以上，才能保證被採樣的信號不失真。

3- FFT的物理意義

現假設我們需要對一個信號進行採樣然後做FFT分析，設定其採樣頻率為Fs，信號的頻率F，採樣點數為N。那麼FFT之後結果其實就是一個為N點的複數。每一個點就對應著一個頻率點。該點的模值，就是該頻率值下的幅度特性。

假設經過FFT之後得到的某點n，使用複數表示為a+bi，則該點的參數如下：


模值為A(n)                                                                                                    

相位為P(n) = atan2(b,    a)

頻率表達式為：Fn = (n - 1) * Fs / N

幅度(n ¹ 1) = A(n) / (N / 2)

幅度(n =    1) = A(n) / N     (直流分量0Hz)

解析度 = Fs / N

 

例：


假設現在有一個信號由三個波形組成，分別是幅度為2的直流信號、幅度為3頻率為50Hz相位為-30度的正弦波、幅度為1.5頻率為75Hz相位為90度的正弦波。使用數學表達式表示為如下所示：

S=Adc+A1*cos(2*pi*F1*t+pi*P1/180)+A2*cos(2*pi*F2*t+pi*P2/180);


假設我們以256Hz採樣頻率對該信號進行採樣，採樣點數為256點，即N = 256。由上面總結的公式可是某點n對應的頻率為 Fn = (n - 1) * 256 / 256 = n-1，我們使用數學模型模擬的信號頻率分別為0Hz、50Hz、75Hz，由於採樣的點是從1開始的，因此對應FFT運算之後的數據應該是1點、51點、76點，只要分析以上三點的數據即可知道對應波形的幅度、相位等信息。

1點： 512+0i
2點： -2.6195E-14 - 1.4162E-13i
3點： -2.8586E-14 - 1.1898E-13i

…
50點：-6.2076E-13 - 2.1713E-12i
51點：332.55 - 192i
52點：-1.6707E-12 - 1.5241E-12i

…
75點：-2.2199E-13 -1.0076E-12i
76點：3.4315E-12 + 192i
77點：-3.0263E-14 +7.5609E-13i

 

由以上數據可以分析出1點、51點、76點對應的模值分別為512、384、192，因此按照以上公式可以得出n=1時頻率為0Hz點對應的為直流分量，該點的幅度為：A(1) = 512 / N = 2，51點對應的幅度為：A(51)= 384 / (N / 2)= 3，76點對應的幅度為：A(76)= 192 / (N / 2)= 1.5。

相位計算：對於直流信號來說無相位可言，51點對應的相位為：atan(-192, 332.55) = -0.5236，計算的結果為弧度，轉換成角度為：180*（-0.5236）/ p = -30.0001度，76點對應的相位為：atan(192,3.4315E-12) = 1.5708，換算成角度為：90.0002，由此可見，分析後的數據和我們設定的數學模型是相符合的。Matlab仿真結果如下：



Matlab仿真代碼如下： 


close    all;          %先關閉所有圖片

Adc=2;              %直流分量幅度

A1=3;               %頻率F1信號的幅度

A2=1.5;             %頻率F2信號的幅度

F1=50;              %信號1頻率(Hz)

F2=75;              %信號2頻率(Hz)

Fs=256;             %採樣頻率(Hz)

P1=-30;             %信號1相位(度)

P2=90;             %信號相位(度)

N=256;              %採樣點數

t=[0:1/Fs:N/Fs];      %採樣時刻

 

%信號

S=Adc+A1*cos(2*pi*F1*t+pi*P1/180)+A2*cos(2*pi*F2*t+pi*P2/180);

%顯示原始信號

plot(S);

title('原始信號');

 

figure;

Y =    fft(S,N);      %做FFT變換

Ayy =    (abs(Y));    %取模

plot(Ayy(1:N));       %顯示原始的FFT模值結果

title('FFT 模值');

 

figure;

Ayy=Ayy/(N/2);              %換算成實際的幅度

Ayy(1)=Ayy(1)/2;

F=([1:N]-1)*Fs/N;              %換算成實際的頻率值

plot(F(1:N/2),Ayy(1:N/2));   %顯示換算後的FFT模值結果

title('幅度-頻率曲線圖');

figure;

Pyy=[1:N/2];

for    i=1:N/2

Pyy(i)=phase(Y(i));             %計算相位

Pyy(i)=Pyy(i)*180/pi;           %換算為角度

end;

plot(F(1:N/2),Pyy(1:N/2));      %顯示相位圖

title('相位-頻率曲線圖');


4- FFT的算法實現

本節將使用STM32官方提供的DSP庫運行FFT算法，使用DSP庫的好處在於FFT的運算時間很快，在72Mhz的主頻下運行256點的FFT僅僅只需要零點幾毫秒，能基本滿足我們的測試FFT的需求。

1.  void Init_Single(void)  

2.  {  

3.      unsigned short i;  

4.      float fx;  

5.        

6.      for(i = 0; i < N; i++)  

7.      {  

8.          fx = 1500 * sin(PI2 * i * 350.0 / Fs) +  

9.               2700 * sin(PI2 * i * 8400.0 / Fs) +  

10.             4000 * sin(PI2 * i * 18725.0 / Fs);  

11.        In_Array[i] = ((signed short)fx) << 16;  

12.    }  

13.}  

14.  

15.void Get_Single_Message()  

16.{  

17.    signed short lX,lY;  

18.    float X,Y,Mag;  

19.    unsigned short i;  

20.    unsigned int f = 0;  

21.      

22.    SEGGER_RTT_printf(0, "Num         f(Hz)         A         X         Y\n");  

23.      

24.    for(i = 0; i < N/2; i++)  

25.    {  

26.        lX  = (Out_Array[i] << 16) >> 16;  

27.        lY  = (Out_Array[i] >> 16);  

28.        X = N * ((float)lX) / 32768;  

29.        Y = N * ((float)lY) / 32768;  

30.        Mag = sqrt(X * X + Y * Y) / (N / 2);  //某點的幅值=該點的模值/（N/2）  

31.        if(i == 0)  

32.            MagArray[i] = (unsigned long)(Mag * 32768);  

33.        else  

34.        {  

35.       MagArray[i] = (unsigned long)(Mag * 32768);  

36.       f += 175;  

37.    }  

38.                  

39.    SEGGER_RTT_printf(0, " %d    %d    %d     %d     %d\n", i, f, MagArray[i], Mag, Y);  

40.    }  

41.          

42.    SEGGER_RTT_printf(0, "\n");  

43.}  


以上若有不嚴謹之處還請各位指出，謝謝！

相關焦點

  • 每日函數——fft
    fft快速傅立葉變換語法Y  = fft(X)Y  = fft(X
  • 通俗易懂的講解FFT的讓你快速了解FFT
    =1; //經過觀察，發現每級蝶形運算需要N/2次運算，共運算N/2*log2N  次   unsigned char fft_counter=0; //在此要進行補2   N必須是2^n   在此略 //蝶形級數  (L級) L=log2(N);         //計算每級蝶形計算的次數(這裡只是一個初始值)
  • matlab下實現FFT信號分析
    利用matlab做頻譜分析前我們需要了解分析過程中的一些基礎知識，matlab中的 fft 函數用法、fftshift 函數的用法函數 1  fft ：作用：快速傅立葉變換。語法：Y = fft(X)Y = fft(X,n)Y = fft(X,n,dim)語法：Y = fft(X) 用快速傅立葉變換 (FFT) 算法計算 X 的離散傅立葉變換 (DFT)。
  • 用FPGA實現FFT算法(圖)
    快速傅立葉變換(fast fourier transformation，簡稱fft)使dft運算效率提高1～2個數量級。其原因是當n較大時，對dft進行了基4和基2分解運算。fft算法除了必需的數據存儲器ram和旋轉因子rom外，仍需較複雜的運算和控制電路單元，即使現在,實現長點數的fft仍然是很困難。
  • 史上最簡單的FFT(快速傅立葉變換)
    asin(1.0)*2; typedef complex<double> cp; cp tmp[N<<1],a[N<<1],b[N<<1]; int n,m;cp c(int n,int k){ return cp(cos(2*k*PI/n),sin(2*k*PI/n)); } void fft
  • 數字與信號處理實驗3 用FFT進行譜分析
    處，因而用到fftshift語句。移位後頻率範圍為範圍的值，若要處於對稱的頻率範圍，則需採用fftshift(X)。       2.有限長序列x(n)的頻譜^n1; X1=fft(x1);                     %長度為N1的序列及其FFT                            N2=2*N1; n2=0 : N2-1;                     %數據長度加倍為N2                            x2=0.5.
  • 使用STM32 的DSP庫進行FFT變換
    * 使用三角函數生成採樣點，供FFT計算* 進行FFT測試時，按下面順序調用函數即可：* dsp_asm_init();* dsp_asm_test();*/#include "stm32f10x.h"#include "dsp_asm.h"#include "stm32_dsp.h"#include "table_fft.h"
  • 快速傅立葉變換FFT在MATLAB中的實現
    (y); %採集信號的長度t=(0:1:N-1)*T; %定義整個採集時間點t=t';  %轉置成列向量figureplot(t,y)xlabel('時間')ylabel('信號值')title('時域信號')2. fft
  • MATLAB實驗——FFT變換
    實驗代碼：clc;clear;alex=imread('E:\SEO\公眾號\0最菜程序猿\圖片\1.jpg');alex1=rgb2gray(alex);alex2=im2double(alex1);alex=fft2(alex2);subplot(2,1,1);imshow(alex1);alex1=fftshift(alex);alex2
  • 用C語言實現FFT算法
    /*****************fft programe*********************/#include typedef.h #include math.h本文引用地址：http://www.eepw.com.cn/article/150637.htmstruct compx
  • 基於MSP430系列微控制器的FFT算法實現
    由於MSP430系列微控制器的開發軟體不支持複數運算，這裡複數運算需要分解成實部和虛部分別來完成，下面的函數「fft_2sin」用來實現蝶形運算。在前面給出的函數「fft_2 sin」中需要通過三角運算分別完成相位因子實部和虛部的計算。三角函數計算需要花費大量的時間，但是在分析的數據點數量給定以後可以首先完成相位因子的計算，將計算值存儲在一個數據表中，通過查表的方法代替三角函數計算。修改後的基2的FFT算法函數如下。
  • PyTorch中的傅立葉卷積:通過FFT計算大核卷積的數學原理和代碼
    從概念上講，此功能的內部工作原理是：def fft_conv( signal: Tensor, kernel: Tensor, bias: Tensor = None, padding: int = 0, ) -> Tensor: # 1.
  • FFT 的物理意義
    度)P2=90; %信號相位(度)N=256; %採樣點數t=[0:1/Fs:N/Fs]; %採樣時刻%信號S=Adc+A1*cos(2*pi*F1*t+pi*P1/180)+A2*cos(2*pi*F2*t+pi*P2/180);%顯示原始信號plot(S);title(原始信號);figure;Y = fft
  • 基於DSP的FFT算法實現
    2、基於DSP的FFT算法實現  用C語言實現FFT算法/*****************fft programe*********************/#include "typedef.h" #include "math.h" struct compx EE(struct compx b1,struct compx b2){
  • FFT相關原理及使用注意事項
    在網絡上簡單搜索一下例程，在Matlab軟體中簡單敲入fft()，即可做信號分析。而ZLG立功科技-致遠電子的高性能數據挖掘性示波器，FFT分析的樣本數可達4Mpts，這使得示波器可以在最高採樣率下，採樣更長時間的波形。這樣在FFT後，數據的奈奎斯特區就相當寬，而頻率解析度又相當窄，非常適用信號分析與噪聲定位。