今天测试的主要对象是NUCLEO-G431RB 的浮点计算能力，通过进行FFT侧面验证Nucleo-G431RB的计算能力。

$ \" Z$ v. E7 e$ D- x7 i

一、背景知识

* r" J5 ?4 d" n

FFT，是快速傅里叶变换（Fast Fourier Transform）的简称，即利用计算机计算离散傅里叶变换（DFT) 的高效、快速计算方法的统称；快速傅里叶变换是1965年由J.W.库利和T.W.图基提出的。采用这种算法能使计算机计算离散傅里叶变换所需要的乘法次数大为减少，特别是被变换的抽样点数N越多，FFT算法计算量的节省就越显著。
! L9 t; m$ ~" D* B1 e6 j  @# C

上图是有限长序列可以通过离散傅里叶变换(DFT）将其频域也离散化成有限长序列。但其计算量太大，很难实时地处理问题，因此引出了快速傅里叶变换(FFT)。1965年，Cooley和Tukey提出了计算离散傅里叶变换（DFT）的快速算法，将DFT的运算量减少了几个数量级。数字信号处理这门新兴学科也随FFT的出现和发展而迅速发展。根据对序列分解与选取方法的不同而产生了FFT的多种算法，基本算法是基2DIT和基2DIF。FFT在离散傅里叶反变换、线性卷积和线性相关等方面也有重要应用。

二、工程建立与代码移植

ST为此提供了一套静态库，分别为Keil、IAR、GCC这三个开发环境提供相关的数学运算库。
5 n6 N3 K' A0 A' u5 Z0 G" B

引用头文件位置在：G4固件库STM32Cube_FW_G4_Vx.x.x)\Drviers\CMSIS\DSP\Include 目录中；
静态库位置在：G4固件库(STM32Cube_FW_G4_Vx.x.x)\Drviers\CMSIS\DSP\Lib）目录中，根据开发环境选择不同的库文件，包括大小端系统环境选择，是否浮点的库。

其中arm_cortexM4b_math.lib代表设备为大端模式；其中arm_cortexM4bf_math.lib代表设备为浮点大端模式；其中arm_cortexM4l_math.lib代表设备为小端模式；其中arm_cortexM4lf_math.lib代表设备为浮点小端模式；

, Q- n; S& R% G" Z$ P  n  ]" n

其中libarm_cortexM4l_math.a 代表设备为小端模式；其中libarm_cortexM4lf_math.a 代表设备为浮点小端模式；

其中iar_cortexM4b_math.a代表设备为大端模式；其中iar_cortexM4bf_math.a代表设备为浮点大端模式；其中iar_cortexM4l_math.a代表设备为小端模式；其中iar_cortexM4lf_math.a代表设备为浮点小端模式； 此次体验期间适合STM32G431RB的免费的开发环境为STM32CubeIDE，该工具是在Eclipse 环境下使用的编译器为GCC。因此，适合搭建FFT测试内容是带浮点的库文件为：libarm_cortexM4lf_math.a

在STM32CubeIDE新建一个工程FFT_G431RB，打开FFT_G431RB.ioc后在CubeMX中设置主频170MHz，添加lpuart1设备使能，设置频率为115200，lpuart1用于连接STLINK的VCP作为计算结果输出。项目在CubeMX中构建好后生成相关代码。在项目中创建DSP\Include和Lib\GCC目录，将开发包文件中的.h头文件和.a静态库文件移植到项目中，并设置相关编译环境如下图：

添加编译静态库文件位置

添加文件路径

, S1 s/ \2 D# R0 I5 p

添加编译宏定义

在Main.c文件头部添加以下定义代码：

/* USER CODE BEGIN Includes */

#include <stdio.h>

#include <stdlib.h>

#include <inttypes.h>

#include "arm_math.h"

#include "arm_const_structs.h"

/* USER CODE END Includes */

……

/* USER CODE BEGIN PTD */

typedef unsigned int ee_u32;

typedef clock_t CORE_TICKS;

typedef ee_u32 secs_ret;

/* USER CODE END PTD */

……

# G3 d) o, p) q5 J# ]( g
) S4 R* S2 g/ _5 w1 I3 O

5 i+ x4 O+ x) S+ I8 V

/* USER CODE BEGIN PD */

#define EE_TICKS_PER_SEC 1000

//#define NPT 2048//4096

#define SysTick_Counter_Disable ((uint32_t)0xFFFFFFFE)

#define SysTick_Counter_Enable ((uint32_t)0x00000001)

#define SysTick_Counter_Clear ((uint32_t)0x00000000)

#define FFT_LENGTH        1024         //FFT长度，默认是1024点FFT

volatile __IO uint32_t Tick;

uint32_t ifftFlag = 0;

uint32_t doBitReverse = 1;

$ @2 q( r; c0 ?" ~& U8 n

float32_t fft_x[FFT_LENGTH * 2];

float32_t fft_input[FFT_LENGTH * 2];

float32_t fft_output[FFT_LENGTH]; //NPT/2

CORE_TICKS total_time;

/* USER CODE END PD */

……

& d+ a9 O. {1 r/ [- z' ]

/* USER CODE BEGIN PFP */

void fft_test();

void start_time(void);

void stop_time(void);

CORE_TICKS get_time(void);

secs_ret time_in_secs(CORE_TICKS ticks);

secs_ret time_in_subsecs(CORE_TICKS ticks);

/* USER CODE END PFP */

……

在main()函数中增加以下代码

……

  /* USER CODE BEGIN 2 */

  printf("NucleoG431RB Test FFT Speed!\n");

  printf("FFT Point Number:%d\n",FFT_LENGTH);

  printf("CLKFreq:%d\n", HAL_RCC_GetHCLKFreq());

  fft_test();

  /* USER CODE END 2 */

……

2 M  p( v7 D1 H9 Y# L$ F& d

添加相关函数定义

' y+ L/ C+ \' L( _) U, v
( X3 w1 l. o2 G7 n. A5 J# u; h

……

/* USER CODE BEGIN 4 */

void start_time(void) {

  //GETMYTIME(&start_time_val );

  Tick = 0;

  SysTick_Config(HAL_RCC_GetHCLKFreq() / 1000);

}

; p. ]& M1 }' L7 o1 M- ~9 B5 s

void stop_time(void) {

  //GETMYTIME(&stop_time_val );

  /* Stop the Timer and get the encoding time */

  SysTick->CTRL &= SysTick_Counter_Disable;

  /* Clear the SysTick Counter */

  SysTick->VAL = SysTick_Counter_Clear;

}

/ M4 P) B* u2 e: m' i5 t" J" P

CORE_TICKS get_time(void) {

  //CORE_TICKS elapsed = (CORE_TICKS)(MYTIMEDIFF(stop_time_val, start_time_val));

  CORE_TICKS elapsed = (CORE_TICKS) Tick;

  return elapsed;

}

secs_ret time_in_secs(CORE_TICKS ticks) {

  secs_ret retval = ((secs_ret) ticks) / (secs_ret) EE_TICKS_PER_SEC;

  return retval;

}

secs_ret time_in_subsecs(CORE_TICKS ticks) {

  secs_ret retval = ((secs_ret) ticks);

  return retval;

}

……

FFT主函数

……

void fft_test() {

  int i, j = 0;

  int Adc = 3; //直流分量幅度

  float32_t A1 = 7.0; //频率F1 信号的幅度

  float32_t A2 = 1.5; //频率F2信号的幅度

  float32_t A3 = 5.1; //频率F3信号的幅度

  float32_t F1 = 180.0; //信号1频率(Hz)

  float32_t F2 = 390.0; //信号2频率(Hz)

  float32_t F3 = 600.0; //信号3频率(Hz)

  float32_t P1 = 0.0;     //信号1相位(度)

  float32_t P2 = 0.0;     //信号2相位(度)

  float32_t P3 = 0.0;     //信号3相位(度)

  //float32_t Fs = (float32_t)FFT_LENGTH;  //采样频率(Hz)

  //float32_t N  = (float32_t)FFT_LENGTH*2;  //采样点数

  //生成信号序列

  for (i = 0; i < FFT_LENGTH; i++) {

    //采样时刻 i/Fs

    fft_x[2 * i] = 10

        + A1 * arm_sin_f32(2 * PI * i * F1 / FFT_LENGTH + PI * P1 / 180)

        + A2 * arm_sin_f32(2 * PI * i * F2 / FFT_LENGTH + PI * P2 / 180)

        + A3 * arm_sin_f32(2 * PI * i * F3 / FFT_LENGTH + PI * P3 / 180);

    //虚部全部为0

    fft_x[2 * i + 1] = 0;

  }

  start_time();

  for (j = 0; j < 1000; j++) {

    //arm_cfft_sR_f32_len1024，该变量即为"arm_const_structs.h"提供的配置变量，包含头文件后，直接调用即可。

    //memset(fft_output, 0x00, sizeof(float32_t)*FFT_LENGTH);

    memcpy(fft_input, fft_x, sizeof(float32_t)*FFT_LENGTH * 2);

    arm_cfft_f32(&arm_cfft_sR_f32_len1024, fft_input, ifftFlag,  doBitReverse);

    //把运算结果复数求模得幅值

    arm_cmplx_mag_f32(fft_input, fft_output, FFT_LENGTH);

  }

  stop_time();

  total_time = get_time();

  printf("Total time (sub_secs): %d\n", time_in_subsecs(total_time));

  for (i = 0; i < FFT_LENGTH; i++) {

    printf("%f\r\n", i, fft_output);

  }

}……

系统printf输出串口定义

8 v& G+ [+ F2 u; d% b

#ifdef __GNUC__

/* With GCC/RAISONANCE, small printf (option LD Linker->Libraries->Small printf

* set to 'Yes') calls __io_putchar() */

#define PUTCHAR_PROTOTYPE int __io_putchar(int ch)

#else

#define PUTCHAR_PROTOTYPE int fputc(int ch, FILE *f)

#endif /* __GNUC__ */

! n, t# I. w9 H! t

/**

* @brief Retargets the C library printf function to the USART.

* @param None

* @retval None

*/

PUTCHAR_PROTOTYPE {

  /* Place your implementation of fputc here */

  /* e.g. write a character to the EVAL_COM1 and Loop until the end of transmission */

  HAL_UART_Transmit(&hlpuart1, (uint8_t*) &ch, 1, 0xFFFF);

  return ch;

}

2 [: l5 m& N% V' b5 K

/* USER CODE END 4 */

▲ 左右移动，查看完整代码

编译后下载即可实现FFT相关计算。

三、FFT计算能力对比及数值校验

9 k; h- Y/ h, ^3 X

因计算单次FFT速度过快，Tick时钟无法记录，改用计算1000次的1024个点的FFT计算汇总时间，以便对比计算速度。
同样的代码移植到Nucleo-F446ZE的项目中，这里有一个小插曲，本想用Nucleo-F301R8的板子对比速度，但因F301R8的Flash、RAM比较小，无法支持带DSP静态库的应用代码，改用速度相近的F4板子Nucleo-F446进行功能验证和对比，以下是对比结果：

这两个个有FFT代码的Nucleo板输出的数值几乎一样，差别只在小数点后第4位、5位几乎可以忽略。其余的内容可以看出STM32G431RB在170MHz主频下浮点处理能力略低于180MHz主频下的STM32F446。不过G431RB具备数学运算加速器，在个别计算领域上可以向其靠拢。不过这也比F3系列的芯片有了很大的提高。
6 g9 k2 ]% W8 C: ]/ I1 M5 l数值图形展示及验证：

- W7 E6 E& n  P. I. @8 k

将输出的数据通过Excel绘图的方式展现全频域信息，符合预期。
- @' _& w' c3 c相关频点和幅值计算的结果如下图：

, P) y0 h! k" A$ z) K9 J

0频率，约10。
180频率，约7。
' W0 p9 ^& n. g# e$ [2 p390频率，约1.5。

600频率，约5.1。

与代码中生成信号序列的公式相关参数相近。
* |7 {  y, V( k6 fFFT计算及验证正确，单次FFT计算速度接近1毫秒，可用于DSP某些领域的数据处理。另外，相关数值计算还可通过上位机的Matelab、Python等工具进行计算并展示相关图形信息，用于验证FFT的处理结果对比。
6 s2 o& S/ m' \. M3 A1 {

很有用的分享

如果代码格式稍稍调整会更好些。

cyk128 发表于 2020-3-4 17:26
很有用的分享

& C3 Q7 {$ h5 r  m# @已修改

ldptest 发表于 2020-3-5 10:22
, r) D) o- H. ]& ?! v如果代码格式稍稍调整会更好些。

, J" u4 |0 v6 m+ f7 S6 ]好的，谢谢，已经调整

gaosmile 发表于 2020-3-5 15:13
好的，谢谢，已经调整

谢谢！

good things

平均 1.2ms ，计算 1024 点FFT？

楼主，向你请教一个问题，Keil和IAR编译浮点运算的时候，这2个编译器的编译时间差异大吗？我碰到一个奇怪的问题，我有一个文件需要使用大量的浮点运算，已经开启了FPU功能，但是在IAR编译的时候，编译很慢，如果不编译那个浮点运算的文件，只要几十秒就编译完了，但是加上需要大量计算浮点运算的那个文件，编译要20多分钟，编译时间相差10倍，但是在KEIL编译的时候，包括那个需要大量浮点计算的文件，整个工程编译也就几十秒，编译的很快，这个问题不知道楼主有没有碰到过？你在编译傅里叶计算的工程的时候，有没有感觉到这2个编译器编译时间存在差异的？有什么好的建议不？在线等，感谢感谢！！！！
% i+ z& K0 ], ]0 m+ {