完整版教程下載地址：http://www.armbbs.cn/forum.php?mod=viewthread&tid=94547

第33章 STM32F407不限制點數FFT實現

本章主要講解不限制點數FFT的實現。

33.1 初學者重要提示

33.2 不限制點數FFT移植

33.3 不限制點數FFT應用說明

33.4 實驗例程說明(MDK)

33.5 實驗例程說明(IAR)

33.6 總結

33.1 初學者重要提示

1. 由於ARM DSP庫限制最大隻能4096點，有點不夠用，所以整理了個更大點數的。不限制點數，滿足2^n即可，n最小值4， 即16個點的FFT，而最大值不限。
2. 此FFT演算法沒有再使用ARM DSP庫，重新實現了一個。

33.2 不限制點數FFT移植

33.2.1 移植FFT相關檔案

移植下面兩個檔案fft.c和FFTInc.h到工程：

33.2.2 新增路徑

新增路徑，大家根據自己的工程來設定即可：

33.3 不限制點數FFT應用說明

33.3.1 支援的點數範圍

支援最小16點FFT，最大值不限，但需滿足2^n。

33.3.2 函式InitTableFFT

函式原型：

/*
*********************************************************************************************************
*    函 數 名: Int_FFT_TAB
*    功能說明: 正弦和餘弦表
*    形    參: 點數
*    返 回 值: 無
*********************************************************************************************************
 
*/
#if (MAX_FFT_N != 8192) && (MAX_FFT_N != 16384)
float32_t   costab[MAX_FFT_N/2];
float32_t   sintab[MAX_FFT_N/2];
void InitTableFFT(uint32_t n)
{
    uint32_t i;

/* 正常使用下面獲取cos和sin值 */
#if 1
    for (i = 0; i < n; i ++ )
    {
         sintab[ i ]=  sin( 2 * PI * i / MAX_FFT_N ); 
        costab[ i ]
 
=  cos( 2 * PI * i / MAX_FFT_N ); 
    }
    
/* 打印出來是方便整理cos值和sin值陣列，將其放到內部Flash，從而節省RAM */
#else
    printf("const float32_t sintab[%d] = {\r\n", n);
    for (i = 0; i < n; i ++ )
    {
         sintab[ i ]=  sin( 2 * PI * i / MAX_FFT_N ); 
        printf("%.11ff,\r\n", sintab[ i ]);
    }    
    printf("};\r\n");
    
    printf("const float32_t costab[%d] = {\r\n", n);
    for (i = 0; i < n; i ++ )
    {
         sintab[ i ]=  cos( 2 * PI * i / MAX_FFT_N ); 
        printf("%.11ff,\r\n", sintab[ i ]);
    }    
    printf("};\r\n");
#endif
}
#endif

函式描述：

這個函式用於FFT計算過程中需要用到的正弦和餘弦表。對於8192點和16384點已經專門製作了數值表，存到內部Flash，其它點數繼續使用的RAM空間，大家可以根據所使用晶片的RAM和Flash大小，選擇正弦和餘弦值存到RAM還是Flash。

函式引數：

• 第1個引數是FFT點數。

33.3.3 函式cfft

函式原型：

void cfft(struct compx *_ptr, uint32_t FFT_N )

/*
*********************************************************************************************************
*    函 數 名: cfft
*    功能說明: 對輸入的複數組進行快速傅立葉變換（FFT）
*    形    參: *_ptr 複數結構體組的首地址指標struct型 
*             FFT_N 表示點數
*    返 回 值: 無
*********************************************************************************************************
*/
void cfft(struct compx *_ptr, uint32_t FFT_N )
{
    float32_t TempReal1, TempImag1, TempReal2, TempImag2;
    uint32_t k,i,j,z;
    uint32_t Butterfly_NoPerColumn;                    /* 每級蝶形的蝶形組數 */
    uint32_t Butterfly_NoOfGroup;                    /* 蝶形組的第幾組 */
    uint32_t Butterfly_NoPerGroup;                    /* 蝶形組的第幾個蝶形 */
    uint32_t ButterflyIndex1,ButterflyIndex2,P,J;
    uint32_t L;
    uint32_t M;

    z=FFT_N/2;    /* 變址運算，即把自然順序變成倒位序，採用雷德演算法 */
    
    for(i=0,j=0;i<FFT_N-1;i++)        
    {
        /* 
          如果i<j,即進行變址 i=j說明是它本身，i>j說明前面已經變換過了，不許再變化，注意這裡一般是實數
 有虛數部分 設定成結合體 
        */
        if(i<j)                                        
        {                                            
            TempReal1  = _ptr[j].real;               
            _ptr[j].real= _ptr[i].real;
            _ptr[i].real= TempReal1;
        }
         
        k=z;         /*求j的下一個倒位序 */
        
        while(k<=j)  /* 如果k<=j,表示j的最高位為1 */  
        {           
            j=j-k; /* 把最高位變成0 */
            k=k/2; /* k/2，比較次高位，依次類推，逐個比較，直到某個位為0，通過下面那句j=j+k使其變為1 */
        }
        
        j=j+k;     /* 求下一個反序號，如果是0，則把0改為1 */
    }
    
/* 第L級蝶形(M)第Butterfly_NoOfGroup組(Butterfly_NoPerColumn)第J個蝶形(Butterfly_NoPerGroup)****** */
/* 蝶形的組數以2的倍數遞減Butterfly_NoPerColumn，每組中蝶形的個數以2的倍數遞增Butterfly_NoPerGroup */
/* 在計算蝶形時，每L列的蝶形組數,一共有M列，每組蝶形中蝶形的個數,蝶形的階數(0,1,2.....M-1) */
    Butterfly_NoPerColumn = FFT_N;                             
    Butterfly_NoPerGroup = 1;    
    M = log2(FFT_N);
    for ( L = 0;L < M; L++ )                                 
    {
        Butterfly_NoPerColumn >>= 1;        /* 蝶形組數 假如N=8，則(4,2,1) */
        
        //第L級蝶形 第Butterfly_NoOfGroup組    （0,1，....Butterfly_NoOfGroup-1）                    
        for ( Butterfly_NoOfGroup = 0;Butterfly_NoOfGroup < Butterfly_NoPerColumn;Butterfly_NoOfGroup++ )
        {  
            /* 第 Butterfly_NoOfGroup 組中的第J個 */
            for ( J = 0;J < Butterfly_NoPerGroup;J ++ )        
            {   /* 第 ButterflyIndex1 和第 ButterflyIndex2 個元素作蝶形運算,WNC */
/* (0,2,4,6)(0,1,4,5)(0,1,2,3) */
/* 兩個要做蝶形運算的數相距Butterfly_NoPerGroup，ge=1,2,4 */                        
/* 這裡相當於P=J*2^(M-L),做了一個換算下標都是N (0,0,0,0)(0,2,0,2)(0,1,2,3) */
                ButterflyIndex1 = ( ( Butterfly_NoOfGroup * Butterfly_NoPerGroup ) << 1 ) + J;
                ButterflyIndex2 = ButterflyIndex1 + Butterfly_NoPerGroup;
                P = J * Butterfly_NoPerColumn;                
                
                /* 計算和轉換因子乘積 */
                TempReal2 = _ptr[ButterflyIndex2].real * costab[ P ] +  _ptr[ButterflyIndex2].imag * 
sintab[ P ];
                TempImag2 = _ptr[ButterflyIndex2].imag * costab[ P ] -  _ptr[ButterflyIndex2].real * 
sintab[ P ] ;
                TempReal1 = _ptr[ButterflyIndex1].real;
                TempImag1 = _ptr[ButterflyIndex1].imag;
                
                /* 蝶形運算 */
                _ptr[ButterflyIndex1].real = TempReal1 + TempReal2;    
                _ptr[ButterflyIndex1].imag = TempImag1 + TempImag2;
                _ptr[ButterflyIndex2].real = TempReal1 - TempReal2;
                _ptr[ButterflyIndex2].imag = TempImag1 - TempImag2;
            }
        } 
        
        Butterfly_NoPerGroup<<=1; /* 一組中蝶形的個數(1,2,4) */
    }
}

函式描述：

這個函式用於複數FFT變換。

函式引數：

• 第1個引數是複數格式。
• 第2個引數是FFT點數，最小值16，最大值不限，滿足滿足2^n即可。

33.3.4 FFT幅頻響應舉例

下面通過函式cff將正弦波做8192點FFT。

/*
*********************************************************************************************************
*    函 數 名: cfft_f32_mag
*    功能說明: 計算幅頻
*    形    參：無
*    返 回 值: 無
*********************************************************************************************************
*/
static void cfft_f32_mag(void)
{
    uint32_t i;
    
    /* 計算一批sin，cos係數 */
#if (MAX_FFT_N != 8192) && (MAX_FFT_N != 16384)
    InitTableFFT(MAX_FFT_N);
#endif
    
    for(i=0; i<MAX_FFT_N; i++)
    {
        /* 波形是由直流分量，500Hz正弦波組成，波形取樣率MAX_FFT_N，初始相位60° */
        s[i].real = 1 + cos(2*3.1415926f*500*i/MAX_FFT_N + 3.1415926f/3);        
        s[i].imag = 0;
    }
    
    /* MAX_FFT_N點快速FFT */ 
    cfft(s, MAX_FFT_N);
    
    /* 計算幅頻 */ 
    for(i=0; i<MAX_FFT_N; i++)
    {
        arm_sqrt_f32((float32_t)(s[i].real *s[i].real+ s[i].imag*s[i].imag ), &s[i].real); 
    }
    
    /* 串列埠列印求解的幅頻 */
    for(i=0; i<MAX_FFT_N; i++)
    {
        printf("%f\r\n", s[i].real);
    }
}

執行函式cfft_f32_mag可以通過串列埠列印FFT結果：

從上面的結果中可以出直流分量和正弦波幅值都是沒有問題的。

33.3.5 FFT相頻響應舉例

下面通過函式cff將正弦波做8192點FFT變換：

/*
*********************************************************************************************************
*    函 數 名: PowerPhaseRadians_f32
*    功能說明: 求相位
*    形    參：_usFFTPoints  複數個數，每個複數是兩個float32_t數值
*             _uiCmpValue  比較值，需要求出相位的數值
*    返 回 值: 無
*********************************************************************************************************
*/
void PowerPhaseRadians_f32(uint16_t _usFFTPoints, float32_t _uiCmpValue)        
{
    float32_t lX, lY;
    uint32_t i;
    float32_t phase;
    float32_t mag;
    
    
    for (i=0; i <_usFFTPoints; i++)
    {
        lX= s[i].real;       /* 實部 */
        lY= s[i].imag;    /* 虛部 */ 
        
         phase = atan2f(lY, lX);                            /* atan2求解的結果範圍是(-pi, pi], 弧度制 */
        arm_sqrt_f32((float32_t)(lX*lX+ lY*lY), &mag);   /* 求模 */
        
        if(_uiCmpValue > mag)
        {
            s[i].real = 0;            
        }
        else
        {
            s[i].real= phase* 180.0f/3.1415926f;   /* 將求解的結果由弧度轉換為角度 */
        }
    }
}

/*
*********************************************************************************************************
*    函 數 名: cfft_f32_phase
*    功能說明: 計算相頻
*    形    參：無
*    返 回 值: 無
*********************************************************************************************************
*/
static void cfft_f32_phase(void)
{
    uint32_t i;
    

    /* 計算一批sin，cos係數 */
#if (MAX_FFT_N != 8192) && (MAX_FFT_N != 16384)
    InitTableFFT(MAX_FFT_N);
#endif
    
    for(i=0; i<MAX_FFT_N; i++)
    {
        /* 波形是由直流分量，500Hz正弦波組成，波形取樣率MAX_FFT_N，初始相位60° */
        s[i].real = 1 + cos(2*3.1415926f*500*i/MAX_FFT_N + 3.1415926f/3);        
        s[i].imag = 0;
    }
    
    /* MAX_FFT_N點快速FFT */ 
    cfft(s, MAX_FFT_N);
    
    /* 求相頻 */
    PowerPhaseRadians_f32(MAX_FFT_N, 0.5f);
    
    /* 串列埠列印求解相頻 */
    for(i=0; i<MAX_FFT_N; i++)
    {
        printf("%f\r\n", s[i].real);
    }
    
}

執行函式cfft_f32_phase可以通過串列埠列印FFT結果：

從上面的結果中可以出計算的初始相位是沒有問題的。

33.4 實驗例程說明（MDK）

配套例子：

V5-223_不限制點數FFT實現

實驗目的：

1. 學習不限制點數FFT。

實驗內容：

1. 啟動一個自動重灌軟體定時器，每100ms翻轉一次LED2。
2. 按下按鍵K1，串列埠列印8192點FFT的幅頻響應。
3. 按下按鍵K2，串列埠列印8192點FFT的相頻響應。

使用AC6注意事項

特別注意附件章節C的問題

上電後串列埠列印的資訊：

波特率 115200，資料位 8，奇偶校驗位無，停止位 1。

RTT方式列印資訊：

程式設計：

系統棧大小分配：

硬體外設初始化

硬體外設的初始化是在 bsp.c 檔案實現：

/*
*********************************************************************************************************
*    函 數 名: bsp_Init
*    功能說明: 初始化所有的硬體裝置。該函式配置CPU暫存器和外設的暫存器並初始化一些全域性變數。只需要呼叫一次
*    形    參：無
*    返 回 值: 無
*********************************************************************************************************
*/
void bsp_Init(void)
{
    /* 
       STM32F407 HAL 庫初始化，此時系統用的還是F407自帶的16MHz，HSI時鐘:
       - 呼叫函式HAL_InitTick，初始化滴答時鐘中斷1ms。
       - 設定NVIC優先順序分組為4。
     */
    HAL_Init();

    /* 
       配置系統時鐘到168MHz
       - 切換使用HSE。
       - 此函式會更新全域性變數SystemCoreClock，並重新配置HAL_InitTick。
    */
    SystemClock_Config();

    /* 
       Event Recorder：
       - 可用於程式碼執行時間測量，MDK5.25及其以上版本才支援，IAR不支援。
       - 預設不開啟，如果要使能此選項，務必看V5開發板使用者手冊第8章
    */    
#if Enable_EventRecorder == 1  
    /* 初始化EventRecorder並開啟 */
    EventRecorderInitialize(EventRecordAll, 1U);
    EventRecorderStart();
#endif
    
    bsp_InitKey();        /* 按鍵初始化，要放在滴答定時器之前，因為按鈕檢測是通過滴答定時器掃描 */
    bsp_InitTimer();      /* 初始化滴答定時器 */
    bsp_InitUart();    /* 初始化串列埠 */
    bsp_InitLed();        /* 初始化LED */        
}

主功能：

主程式實現如下操作：

• 啟動一個自動重灌軟體定時器，每100ms翻轉一次LED2。
• 按下按鍵K1，串列埠列印8192點FFT的幅頻響應。
• 按下按鍵K2，串列埠列印8192點FFT的相頻響應。

/*
*********************************************************************************************************
*    函 數 名: main
*    功能說明: c程式入口
*    形    參: 無
*    返 回 值: 錯誤程式碼(無需處理)
*********************************************************************************************************
*/
int main(void)
{
    uint8_t ucKeyCode;        /* 按鍵程式碼 */
    

    bsp_Init();        /* 硬體初始化 */
    PrintfLogo();    /* 列印例程資訊到串列埠1 */

    PrintfHelp();    /* 列印操作提示資訊 */
    

    bsp_StartAutoTimer(0, 100);    /* 啟動1個100ms的自動重灌的定時器 */

    /* 進入主程式迴圈體 */
    while (1)
    {
        bsp_Idle();        /* 這個函式在bsp.c檔案。使用者可以修改這個函式實現CPU休眠和喂狗 */
        

        if (bsp_CheckTimer(0))    /* 判斷定時器超時時間 */
        {
            /* 每隔100ms 進來一次 */
            bsp_LedToggle(4);    /* 翻轉LED2的狀態 */   
        }
        
        ucKeyCode = bsp_GetKey();    /* 讀取鍵值, 無鍵按下時返回 KEY_NONE = 0 */
        if (ucKeyCode != KEY_NONE)
        {
            switch (ucKeyCode)
            {
                case KEY_DOWN_K1:            /* K1鍵按下 */
                    cfft_f32_mag();
                    break;
                
                case KEY_DOWN_K2:            /* K2鍵按下 */
                    cfft_f32_phase();
                    break;
                
                    
                default:
                    /* 其它的鍵值不處理 */
                    break;
            }
        }

    }
}

33.5 實驗例程說明（IAR）

配套例子：

V5-223_不限制點數FFT實現

實驗目的：

1. 學習不限制點數FFT。

實驗內容：

1. 啟動一個自動重灌軟體定時器，每100ms翻轉一次LED2。
2. 按下按鍵K1，串列埠列印8192點FFT的幅頻響應。
3. 按下按鍵K2，串列埠列印8192點FFT的相頻響應。

上電後串列埠列印的資訊：

波特率 115200，資料位 8，奇偶校驗位無，停止位 1。

RTT方式列印資訊：

程式設計：

系統棧大小分配：

硬體外設初始化

硬體外設的初始化是在 bsp.c 檔案實現：

/*
*********************************************************************************************************
*    函 數 名: bsp_Init
*    功能說明: 初始化所有的硬體裝置。該函式配置CPU暫存器和外設的暫存器並初始化一些全域性變數。只需要呼叫一次
*    形    參：無
*    返 回 值: 無
*********************************************************************************************************
*/
void bsp_Init(void)
{
    /* 
       STM32F407 HAL 庫初始化，此時系統用的還是F407自帶的16MHz，HSI時鐘:
       - 呼叫函式HAL_InitTick，初始化滴答時鐘中斷1ms。
       - 設定NVIC優先順序分組為4。
     */
    HAL_Init();

    /* 
       配置系統時鐘到168MHz
       - 切換使用HSE。
       - 此函式會更新全域性變數SystemCoreClock，並重新配置HAL_InitTick。
    */
    SystemClock_Config();

    /* 
       Event Recorder：
       - 可用於程式碼執行時間測量，MDK5.25及其以上版本才支援，IAR不支援。
       - 預設不開啟，如果要使能此選項，務必看V5開發板使用者手冊第8章
    */    
#if Enable_EventRecorder == 1  
    /* 初始化EventRecorder並開啟 */
    EventRecorderInitialize(EventRecordAll, 1U);
    EventRecorderStart();
#endif
    
    bsp_InitKey();        /* 按鍵初始化，要放在滴答定時器之前，因為按鈕檢測是通過滴答定時器掃描 */
    bsp_InitTimer();      /* 初始化滴答定時器 */
    bsp_InitUart();    /* 初始化串列埠 */
    bsp_InitLed();        /* 初始化LED */        
}

主功能：

主程式實現如下操作：

• 啟動一個自動重灌軟體定時器，每100ms翻轉一次LED2。
• 按下按鍵K1，串列埠列印8192點FFT的幅頻響應。
• 按下按鍵K2，串列埠列印8192點FFT的相頻響應。

/*
*********************************************************************************************************
*    函 數 名: main
*    功能說明: c程式入口
*    形    參: 無
*    返 回 值: 錯誤程式碼(無需處理)
*********************************************************************************************************
*/
int main(void)
{
    uint8_t ucKeyCode;        /* 按鍵程式碼 */
    

    bsp_Init();        /* 硬體初始化 */
    PrintfLogo();    /* 列印例程資訊到串列埠1 */

    PrintfHelp();    /* 列印操作提示資訊 */
    

    bsp_StartAutoTimer(0, 100);    /* 啟動1個100ms的自動重灌的定時器 */

    /* 進入主程式迴圈體 */
    while (1)
    {
        bsp_Idle();        /* 這個函式在bsp.c檔案。使用者可以修改這個函式實現CPU休眠和喂狗 */
        

        if (bsp_CheckTimer(0))    /* 判斷定時器超時時間 */
        {
            /* 每隔100ms 進來一次 */
            bsp_LedToggle(4);    /* 翻轉LED2的狀態 */   
        }
        
        ucKeyCode = bsp_GetKey();    /* 讀取鍵值, 無鍵按下時返回 KEY_NONE = 0 */
        if (ucKeyCode != KEY_NONE)
        {
            switch (ucKeyCode)
            {
                case KEY_DOWN_K1:            /* K1鍵按下 */
                    cfft_f32_mag();
                    break;
                
                case KEY_DOWN_K2:            /* K2鍵按下 */
                    cfft_f32_phase();
                    break;
                
                    
                default:
                    /* 其它的鍵值不處理 */
                    break;
            }
        }

    }
}

33.6 總結

本章節主要設計一個不限制點數的FFT功能，實際專案用到的地方也比較多，望初學著掌握。