【STM32F429的DSP教程】第30章 STM32F429複數浮點FFT（支援單精度和雙精度）
阿新 • • 發佈：2021-06-11
完整版教程下載地址：http://www.armbbs.cn/forum.php?mod=viewthread&tid=94547
第30章 STM32F429複數浮點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計算的做對比。
對比前需要先將串列埠打印出的資料載入到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_f64計算的結果基本是一直的。幅頻響應求出的幅值和相頻響應中的求出的初始相角都是沒問題的。
30.5 實驗例程說明（MDK）
配套例子：
V6-220_複數浮點FTT(支援單精度和雙精度)
實驗目的：
- 學習複數浮點FFT，支援單精度浮點和雙精度浮點
實驗內容：
- 啟動一個自動重灌軟體定時器，每100ms翻轉一次LED2。
- 按下按鍵K1，串列埠列印1024點複數單精度FFT的幅頻響應和相頻響應。
- 按下按鍵K2，串列埠列印1024點複數雙精度FFT的幅頻響應和相頻響應。
使用AC6注意事項
特別注意附件章節C的問題
上電後串列埠列印的資訊：
波特率 115200，資料位 8，奇偶校驗位無，停止位 1。
RTT方式列印資訊：
程式設計：
系統棧大小分配：
硬體外設初始化
硬體外設的初始化是在 bsp.c 檔案實現：
/* ********************************************************************************************************* * 函 數 名: bsp_Init * 功能說明: 初始化所有的硬體裝置。該函式配置CPU暫存器和外設的暫存器並初始化一些全域性變數。只需要呼叫一次 * 形 參：無 * 返 回 值: 無 ********************************************************************************************************* */ void bsp_Init(void) { /* STM32F429 HAL 庫初始化，此時系統用的還是F429自帶的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_InitExtIO(); /* 初始化擴充套件IO */ bsp_InitLed(); /* 初始化LED */ }
主功能：
主程式實現如下操作：
- 啟動一個自動重灌軟體定時器，每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）
配套例子：
V6-220_複數浮點FTT(支援單精度和雙精度)
實驗目的：
- 學習複數浮點FFT，支援單精度浮點和雙精度浮點
實驗內容：
- 啟動一個自動重灌軟體定時器，每100ms翻轉一次LED2。
- 按下按鍵K1，串列埠列印1024點複數單精度FFT的幅頻響應和相頻響應。
- 按下按鍵K2，串列埠列印1024點複數雙精度FFT的幅頻響應和相頻響應。
使用AC6注意事項
特別注意附件章節C的問題
上電後串列埠列印的資訊：
波特率 115200，資料位 8，奇偶校驗位無，停止位 1。
RTT方式列印資訊：
程式設計：
系統棧大小分配：
硬體外設初始化
硬體外設的初始化是在 bsp.c 檔案實現：
/* ********************************************************************************************************* * 函 數 名: bsp_Init * 功能說明: 初始化所有的硬體裝置。該函式配置CPU暫存器和外設的暫存器並初始化一些全域性變數。只需要呼叫一次 * 形 參：無 * 返 回 值: 無 ********************************************************************************************************* */ void bsp_Init(void) { /* STM32F429 HAL 庫初始化，此時系統用的還是F429自帶的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_InitExtIO(); /* 初始化擴充套件IO */ bsp_InitLed(); /* 初始化LED */ }
主功能：
主程式實現如下操作：
- 啟動一個自動重灌軟體定時器，每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.7 總結
本章節設計到FFT實現，有興趣的可以深入瞭解原始碼的實現。