【STM32H7的DSP教程】第30章 STM32H7復數浮點FFT（支持單精度和雙精度）

本文轉載自查看原文 2021-06-11 16:59 1298 第2版DSP和CMSIS-NN神經網絡教程

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

第30章 STM32H7復數浮點FFT（支持單精度和雙精度）

本章主要講解復數浮點FTT，支持單精度和雙精度。

30.1 初學者重要提示

30.2 復數浮點FFT 說明

30.3 單精度函數arm_cfft_f32的使用（含幅頻和相頻）

30.4 雙精度函數arm_cfft_f64的使用（含幅頻和相頻）

30.5 實驗例程說明(MDK)

30.6 實驗例程說明(IAR)

30.7 總結

30.1 初學者重要提示

新版DSP庫浮點FFT推薦使用混合基函數arm_cfft_f32，而基2函數arm_cfft_radix2_f32和基4函數arm_cfft_radix4_f32將廢棄。ARM說明如下：

Earlier releases of the library provided separate radix-2 and radix-4 algorithms that operated on floating-point data.  These functions are still provided but are deprecated.  The older functions are slower and less general than the new functions.
DSP庫的早期發行版提供了單獨的radix-2和radix-4對浮點數據進行運算的算法。 這些功能仍然提供，但已棄用。 相比新版函數，老版的功能較慢且通用性較低

30.2 復數浮點FFT說明

30.2.1 功能描述

當前復數FFT函數支持三種數據類型，分別是浮點，定點Q31和Q15。這些FFT函數有一個共同的特點，就是用於輸入信號的緩沖，在轉化結束后用來存儲輸出結果。這樣做的好處是節省了RAM空間，不需要為輸入和輸出結果分別設置緩存。由於是復數FFT，所以輸入和輸出緩存要存儲實部和虛部。存儲順序如下：{real[0], imag[0], real[1], imag[1],………………} ，在使用中切記不要搞錯。

30.2.2 浮點FFT

浮點復數FFT使用了一個混合基數算法，通過多個基8與單個基2或基4算法實現。根據需要，該算法支持的長度[16，32，64，...，4096]和每個長度使用不同的旋轉因子表。

浮點復數FFT使用了標准的FFT定義，FFT正變換的輸出結果會被放大fftLen倍數，計算FFT逆變換的時候會縮小到1/fftLen。這樣就與教科書中的定義一致了。

定義好的旋轉因子和位反轉表已經在頭文件arm_const_structs.h中定義好了，調用浮點FFT函數arm_cfft_f32時，包含相應的頭文件即可。比如：

arm_cfft_f32(arm_cfft_sR_f32_len64, pSrc, 1, 1)

上式就是計算一個64點的FFT逆變換包括位反轉。數據結構arm_cfft_sR_f32_len64可以認為是常數，計算的過程中是不能修改的。同樣是這種數據結構還能用於混合基的FFT正變換和逆變換。

早期發布的浮點復數FFT函數版本包含基2和基4兩種方法實現的，但是不推薦大家再使用。現在全部用arm_cfft_f32代替了。

30.3 單精度函數arm_cfft_f32的使用(含幅頻和相頻)

30.3.1 函數說明

函數原型：

void arm_cfft_f32(
  const arm_cfft_instance_f32 * S,
        float32_t * p1,
        uint8_t ifftFlag,
        uint8_t bitReverseFlag)

函數描述：

這個函數用於單精度浮點復數FFT。

函數參數：

1、第1個參數是封裝好的浮點FFT例化，支持的參數如下：

arm_cfft_sR_f32_len16，16點FFT
arm_cfft_sR_f32_len32，32點FFT
arm_cfft_sR_f32_len64，64點FFT
arm_cfft_sR_f32_len128，128點FFT
arm_cfft_sR_f32_len256，256點FFT
arm_cfft_sR_f32_len512，512點FFT
arm_cfft_sR_f32_len1024，1024點FFT
arm_cfft_sR_f32_len2048，2048點FFT
arm_cfft_sR_f32_len4096，4096點FFT

2、第2個參數是復數地址，存儲順序是實部，虛部，實部，虛部，依次類推。

3、第3個參數用於設置正變換和逆變換，ifftFlag=0表示正變換，ifftFlag=1表示逆變換。

4、第4個參數用於設置輸出位反轉，bitReverseFlag=1表示使能，bitReverseFlag=0表示禁止。

30.3.2 使用舉例並和Matlab比較

下面通過在開發板上運行這個函數並計算幅頻相應，然后再與Matlab計算的結果做對比。

/*
*********************************************************************************************************
*    函 數 名: arm_cfft_f32_app
*    功能說明: 調用函數arm_cfft_f32計算幅頻和相頻
*    形    參：無
*    返 回 值: 無
*********************************************************************************************************
*/
static void arm_cfft_f32_app(void)
{
    uint16_t i;
    
    ifftFlag = 0; 
    doBitReverse = 1; 
    
    /* 按照實部，虛部，實部，虛部..... 的順序存儲數據 */
    for(i=0; i<TEST_LENGTH_SAMPLES; i++)
    {
        /* 波形是由直流分量，50Hz正弦波組成，波形采樣率1024，初始相位60° */
        testInput_f32[i*2] = 1 + cos(2*3.1415926f*50*i/1024 + 3.1415926f/3);
        testInput_f32[i*2+1] = 0;
    }
    
    /* CFFT變換 */ 
    arm_cfft_f32(&arm_cfft_sR_f32_len1024, testInput_f32, ifftFlag, doBitReverse);

    /* 求解模值  */ 
    arm_cmplx_mag_f32(testInput_f32, testOutput_f32, TEST_LENGTH_SAMPLES);
    

    printf("=========================================\r\n");    
    
    /* 求相頻 */
    PowerPhaseRadians_f32(testInput_f32, Phase_f32, TEST_LENGTH_SAMPLES, 0.5f);
    
    /* 串口打印求解的模值 */
    for(i=0; i<TEST_LENGTH_SAMPLES; i++)
    {
        printf("%f, %f\r\n", testOutput_f32[i], Phase_f32[i]);
    }    
}

運行函數arm_cfft_f32_app可以通過串口打印出計算的模值和相角，下面我們就通過Matlab計算的模值和相角跟arm_cfft_f32計算的做對比。

對比前需要先將串口打印出的數據加載到Matlab中，並給這個數組起名sampledata，加載方法在前面的教程的第13章13.6小結已經講解，這里不做贅述了。Matlab中運行的代碼如下：：

Fs = 1024;               % 采樣率
N  = 1024;               % 采樣點數
n  = 0:N-1;              % 采樣序列
t  = 0:1/Fs:1-1/Fs;      % 時間序列
f = n * Fs / N;          %真實的頻率

%波形是由直流分量，50Hz正弦波正弦波組成
x = 1 + cos(2*pi*50*t + pi/3)   ;  
y = fft(x, N);               %對原始信號做FFT變換
Mag = abs(y);

subplot(2,2,1);
plot(f, Mag); 
title('Matlab計算幅頻響應');
xlabel('頻率');
ylabel('賦值');

subplot(2,2,2);
realvalue = real(y);
imagvalue = imag(y);
plot(f, atan2(imagvalue, realvalue)*180/pi.*(Mag>=200)); 
title('Matlab計算相頻響應');
xlabel('頻率');
ylabel('相角');

subplot(2,2,3);
plot(f, sampledata1);  %繪制STM32計算的幅頻相應
title('STM32計算幅頻響應');
xlabel('頻率');
ylabel('賦值');

subplot(2,2,4);
plot(f, sampledata2);   %繪制STM32計算的相頻相應
title('STM32計算相頻響應');
xlabel('頻率');
ylabel('相角');

運行Matlab后的輸出結果如下：

從上面的對比結果中可以看出，Matlab和函數arm_cfft_f32計算的結果基本是一直的。幅頻響應求出的幅值和相頻響應中的求出的初始相角都是沒問題的。

30.4 雙精度函數arm_cfft_f64的使用(含幅頻和相頻)

30.4.1 函數說明

函數原型：

void arm_cfft_f64(
  const arm_cfft_instance_f64 * S,
        float64_t * p1,
        uint8_t ifftFlag,
        uint8_t bitReverseFlag)

函數描述：

這個函數用於雙精度浮點復數FFT。

函數參數：

1、第1個參數是封裝好的浮點FFT例化，支持的參數如下：

arm_cfft_sR_f64_len16，16點FFT
arm_cfft_sR_f64_len32，32點FFT
arm_cfft_sR_f64_len64，64點FFT
arm_cfft_sR_f64_len128，128點FFT
arm_cfft_sR_f64_len256，256點FFT
arm_cfft_sR_f64_len512，512點FFT
arm_cfft_sR_f64_len1024，1024點FFT
arm_cfft_sR_f64_len2048，2048點FFT
arm_cfft_sR_f64_len4096，4096點FFT

2、第2個參數是復數地址，存儲順序是實部，虛部，實部，虛部，依次類推。

3、第3個參數用於設置正變換和逆變換，ifftFlag=0表示正變換，ifftFlag=1表示逆變換。

4、第4個參數用於設置輸出位反轉，bitReverseFlag=1表示使能，bitReverseFlag=0表示禁止。

30.4.2 使用舉例並和Matlab比較

下面通過在開發板上運行這個函數並計算幅頻相應，然后再與Matlab計算的結果做對比。

/*
*********************************************************************************************************
*    函 數 名: arm_cfft_f64_app
*    功能說明: 調用函數arm_cfft_f64計算幅頻和相頻
*    形    參：無
*    返 回 值: 無
*********************************************************************************************************
*/
static void arm_cfft_f64_app(void)
{
    uint16_t i;
    float64_t lX,lY;
    
    ifftFlag = 0; 
    doBitReverse = 1; 
    
    /* 按照實部，虛部，實部，虛部..... 的順序存儲數據 */
    for(i=0; i<TEST_LENGTH_SAMPLES; i++)
    {
        /* 波形是由直流分量，50Hz正弦波組成，波形采樣率1024，初始相位60° */
        testInput_f64[i*2] = 1 + cos(2*3.1415926*50*i/1024 + 3.1415926/3);
        testInput_f64[i*2+1] = 0;
    }
    
    /* CFFT變換 */ 
    arm_cfft_f64(&arm_cfft_sR_f64_len1024, testInput_f64, ifftFlag, doBitReverse);

    /* 求解模值  */ 
    for (i =0; i < TEST_LENGTH_SAMPLES; i++)
    {
         lX = testInput_f64[2*i];            /* 實部*/
        lY = testInput_f64[2*i+1];          /* 虛部 */  
        testOutput_f64[i] = sqrt(lX*lX+ lY*lY);   /* 求模 */
    }
    
    printf("=========================================\r\n");    
    
    /* 求相頻 */
    PowerPhaseRadians_f64(testInput_f64, Phase_f64, TEST_LENGTH_SAMPLES, 0.5);
    
    
    /* 串口打印求解的模值 */
    for(i=0; i<TEST_LENGTH_SAMPLES; i++)
    {
        printf("%.11f, %.11f\r\n", testOutput_f64[i], Phase_f64[i]);
    }    
    
}

運行函數arm_cfft_f64_app可以通過串口打印出計算的模值和相角，下面我們就通過Matlab計算的模值和相角跟arm_cfft_f64計算的做對比。

Fs = 1024;               % 采樣率
N  = 1024;               % 采樣點數
n  = 0:N-1;              % 采樣序列
t  = 0:1/Fs:1-1/Fs;      % 時間序列
f = n * Fs / N;          %真實的頻率

%波形是由直流分量，50Hz正弦波正弦波組成
x = 1 + cos(2*pi*50*t + pi/3)   ;  
y = fft(x, N);               %對原始信號做FFT變換
Mag = abs(y);

subplot(2,2,1);
plot(f, Mag); 
title('Matlab計算幅頻響應');
xlabel('頻率');
ylabel('賦值');

subplot(2,2,2);
realvalue = real(y);
imagvalue = imag(y);
plot(f, atan2(imagvalue, realvalue)*180/pi.*(Mag>=200)); 
title('Matlab計算相頻響應');
xlabel('頻率');
ylabel('相角');

subplot(2,2,3);
plot(f, sampledata1);  %繪制STM32計算的幅頻相應
title('STM32計算幅頻響應');
xlabel('頻率');
ylabel('賦值');

subplot(2,2,4);
plot(f, sampledata2);   %繪制STM32計算的相頻相應
title('STM32計算相頻響應');
xlabel('頻率');
ylabel('相角');

運行Matlab后的輸出結果如下：

從上面的對比結果中可以看出，Matlab和函數arm_cfft_f64計算的結果基本是一直的。幅頻響應求出的幅值和相頻響應中的求出的初始相角都是沒問題的。

30.5 實驗例程說明（MDK）

配套例子：

V7-220_復數浮點FTT(支持單精度和雙精度)

實驗目的：

學習復數浮點FFT，支持單精度浮點和雙精度浮點

實驗內容：

啟動一個自動重裝軟件定時器，每100ms翻轉一次LED2。
按下按鍵K1，串口打印1024點復數單精度FFT的幅頻響應和相頻響應。
按下按鍵K2，串口打印1024點復數雙精度FFT的幅頻響應和相頻響應。

使用AC6注意事項

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

上電后串口打印的信息：

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

RTT方式打印信息：

程序設計：

系統棧大小分配：

RAM空間用的DTCM：

硬件外設初始化

硬件外設的初始化是在 bsp.c 文件實現：

/*
*********************************************************************************************************
*    函 數 名: bsp_Init
*    功能說明: 初始化所有的硬件設備。該函數配置CPU寄存器和外設的寄存器並初始化一些全局變量。只需要調用一次
*    形    參：無
*    返 回 值: 無
*********************************************************************************************************
*/
void bsp_Init(void)
{
    /* 配置MPU */
    MPU_Config();
    
    /* 使能L1 Cache */
    CPU_CACHE_Enable();

    /* 
       STM32H7xx HAL 庫初始化，此時系統用的還是H7自帶的64MHz，HSI時鍾:
       - 調用函數HAL_InitTick，初始化滴答時鍾中斷1ms。
       - 設置NVIC優先級分組為4。
     */
    HAL_Init();

    /* 
       配置系統時鍾到400MHz
       - 切換使用HSE。
       - 此函數會更新全局變量SystemCoreClock，並重新配置HAL_InitTick。
    */
    SystemClock_Config();

    /* 
       Event Recorder：
       - 可用於代碼執行時間測量，MDK5.25及其以上版本才支持，IAR不支持。
       - 默認不開啟，如果要使能此選項，務必看V7開發板用戶手冊第8章
    */    
#if Enable_EventRecorder == 1  
    /* 初始化EventRecorder並開啟 */
    EventRecorderInitialize(EventRecordAll, 1U);
    EventRecorderStart();
#endif
    
    bsp_InitKey();        /* 按鍵初始化，要放在滴答定時器之前，因為按鈕檢測是通過滴答定時器掃描 */
    bsp_InitTimer();      /* 初始化滴答定時器 */
    bsp_InitUart();    /* 初始化串口 */
    bsp_InitExtIO();    /* 初始化FMC總線74HC574擴展IO. 必須在 bsp_InitLed()前執行 */    
    bsp_InitLed();        /* 初始化LED */    
}

MPU配置和Cache配置：

數據Cache和指令Cache都開啟。配置了AXI SRAM區（本例子未用到AXI SRAM），FMC的擴展IO區。

/*
*********************************************************************************************************
*    函 數 名: MPU_Config
*    功能說明: 配置MPU
*    形    參: 無
*    返 回 值: 無
*********************************************************************************************************
*/
static void MPU_Config( void )
{
    MPU_Region_InitTypeDef MPU_InitStruct;

    /* 禁止 MPU */
    HAL_MPU_Disable();

    /* 配置AXI SRAM的MPU屬性為Write back, Read allocate，Write allocate */
    MPU_InitStruct.Enable           = MPU_REGION_ENABLE;
    MPU_InitStruct.BaseAddress      = 0x24000000;
    MPU_InitStruct.Size             = MPU_REGION_SIZE_512KB;
    MPU_InitStruct.AccessPermission = MPU_REGION_FULL_ACCESS;
    MPU_InitStruct.IsBufferable     = MPU_ACCESS_BUFFERABLE;
    MPU_InitStruct.IsCacheable      = MPU_ACCESS_CACHEABLE;
    MPU_InitStruct.IsShareable      = MPU_ACCESS_NOT_SHAREABLE;
    MPU_InitStruct.Number           = MPU_REGION_NUMBER0;
    MPU_InitStruct.TypeExtField     = MPU_TEX_LEVEL1;
    MPU_InitStruct.SubRegionDisable = 0x00;
    MPU_InitStruct.DisableExec      = MPU_INSTRUCTION_ACCESS_ENABLE;

    HAL_MPU_ConfigRegion(&MPU_InitStruct);
    
    
    /* 配置FMC擴展IO的MPU屬性為Device或者Strongly Ordered */
    MPU_InitStruct.Enable           = MPU_REGION_ENABLE;
    MPU_InitStruct.BaseAddress      = 0x60000000;
    MPU_InitStruct.Size             = ARM_MPU_REGION_SIZE_64KB;    
    MPU_InitStruct.AccessPermission = MPU_REGION_FULL_ACCESS;
    MPU_InitStruct.IsBufferable     = MPU_ACCESS_BUFFERABLE;
    MPU_InitStruct.IsCacheable      = MPU_ACCESS_NOT_CACHEABLE;    
    MPU_InitStruct.IsShareable      = MPU_ACCESS_NOT_SHAREABLE;
    MPU_InitStruct.Number           = MPU_REGION_NUMBER1;
    MPU_InitStruct.TypeExtField     = MPU_TEX_LEVEL0;
    MPU_InitStruct.SubRegionDisable = 0x00;
    MPU_InitStruct.DisableExec      = MPU_INSTRUCTION_ACCESS_ENABLE;
    
    HAL_MPU_ConfigRegion(&MPU_InitStruct);

    /*使能 MPU */
    HAL_MPU_Enable(MPU_PRIVILEGED_DEFAULT);
}

/*
*********************************************************************************************************
*    函 數 名: CPU_CACHE_Enable
*    功能說明: 使能L1 Cache
*    形    參: 無
*    返 回 值: 無
*********************************************************************************************************
*/
static void CPU_CACHE_Enable(void)
{
    /* 使能 I-Cache */
    SCB_EnableICache();

    /* 使能 D-Cache */
    SCB_EnableDCache();
}

主功能：

主程序實現如下操作：

啟動一個自動重裝軟件定時器，每100ms翻轉一次LED2。
按下按鍵K1，串口打印1024點復數單精度FFT的幅頻響應和相頻響應。
按下按鍵K2，串口打印1024點復數雙精度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鍵按下 */
                    arm_cfft_f32_app();
                    break;
                
                case KEY_DOWN_K2:            /* K2鍵按下 */
                    arm_cfft_f64_app();
                    break;
                
                    
                default:
                    /* 其它的鍵值不處理 */
                    break;
            }
        }

    }
}

30.6 實驗例程說明（IAR）