快速傅里叶变换(FFT)是一种数字信号处理中常用的技术,用于将快速序列转换为频域表示。在嵌入式系统中,如基于STM32的微控制器,实现FFT可以帮助解决信号处理的需求,例如声音处理、图像处理等。本文将介绍基于STM32的离散傅里叶变换的原理、实现方法和应用。
d1 U+ r, q' n' k. R G( C! { D. u' o3 G, P
: m" p7 V' i. H/ p+ c
3 l4 k0 B+ E8 U
FFT是一种将时域序列转换为频域表示的技术,它将一个序列的N个采样点映射到频域中N个频率分量。其数学表达式如下:
7 D ?9 c, A: Y8 y w: P. i7 D) H) p" i
4 x( f+ t0 V, W8 b
& k! X7 r% Q) \. u) ^其中,x(n) 是输入序列,X(k) 是输出的频域表示。
* y, x8 s! B2 p ~5 C' l$ I
3 e3 \- a5 x8 a7 Q3 d4 x& `准备工作:
6 B8 R u' b# k, B$ b" D! [, B. C& x6 P; D1 n
: i0 h, s; q) a" P$ ~5 R, W- T" n% N
) V0 f+ _% A) p/ e
9 w6 o& z. q/ Z' T! o9 Y0 @8 W; F! k, r
' c- ~0 B: X# `; u- Z5 z
# u1 d! ?* n" J4 ~/ DKeil中的DSP库(Digital Signal Processing Library,数字信号处理库)是针对ARM Cortex-M处理器系列的一组软件库,用于提供各种数字信号处理功能的支持。这些库提供了一系列优化过的算法,可以帮助开发人员在嵌入式系统中高效地实现音频处理、图像处理、通信系统等各种信号处理应用。$ [6 A; z% S) k- J! s
- I8 Y2 O* |- |9 B
因此我们需要在Keil中安装我们的DSP库。
5 m' Y0 ~9 n" h# D2 _+ q- #include "arm_math.h" // 包含DSP库
复制代码 / x# v. R- |! `9 J+ \/ E: Q: f! Z
首先包含我们的DSP库。( }6 U% |& d9 x& _3 i
, Z% o) a4 b& o2 X' w& h- B5 @9 Y
D7 ^% n s0 D- p/ {0 R- #define FFT_LENGTH 1007 e( y5 r6 ]; \& m. I2 s; b
- // 输入序列# Q+ n- M3 c% {/ X% k, `9 @. e1 L
- float32_t inputSignal[FFT_LENGTH*2];
& Y$ M6 t3 |1 s4 z+ _
$ O* A, h& c' o' h$ a) @- // 输出序列,存储变换后的结果/ X$ ]: S* z3 h/ F Y' b
- float32_t outputSignal[FFT_LENGTH];
复制代码
" b. N) m& ?' g3 H定义FFT的的输入和输出数组还有数组长度3 K5 V9 G% f I. Q
- arm_status status;
( g* _, s9 {) [6 d - arm_cfft_radix4_instance_f32 fft_inst;
& q! `2 X: @# g2 m - status = arm_cfft_radix4_init_f32(&fft_inst, FFT_LENGTH,0,1);
复制代码- void arm_cfft_radix4_init_f32(
+ ?) Z2 U6 {2 X7 }' v - arm_cfft_radix4_instance_f32 * S,$ g1 p0 H8 k E$ n
- uint16_t fftLen, i$ r! f+ _4 N/ I0 Q) k7 r0 f) }
- uint8_t ifftFlag,
* j3 ^% I- z$ n" k+ J: S* _ - uint8_t bitReverseFlag# r+ \: o: n/ X( ]# D) T
- );
复制代码
- g8 }5 ?4 K1 E定义一个状态变量用来显示FFT的初始化是否成功。
4 e; v' k, W! e1 h, M定义一个FFT的配置变量。$ N; i0 u$ S+ s
初始化FFT。5 j( Q0 T, z w3 `! b5 s
S:指向 arm_cfft_radix4_instance_f32 结构体的指针,该结构体定义了 FFT 实例的状态信息。
# @/ U& F- U" ?, F5 V6 i% c0 y) o/ l- @ y- f$ U& v/ r; d
fftLen:FFT 的长度。
. p" _1 ^( K( h2 F$ H- G( a; T, u$ d5 i5 s( I6 z) `) g
ifftFlag:指定是否进行逆变换。如果为 1,则表示初始化的是逆变换的 FFT;如果为 0,则表示初始化的是正变换的 FFT。3 Z, k! `6 u. d; G3 S
l1 Q* ^) X2 C X
bitReverseFlag:指定是否进行比特翻转。如果为 1,则表示进行比特翻转;如果为 0,则表示不进行比特翻转。
- ^, s; A. `3 L3 S
+ N9 R1 D; W& p8 i在FFT算法中,比特(bit)反转是一种关键的步骤,用于将输入数据重新排列为正确的顺序,以便在后续的计算中进行有效处理。
, H0 v! H }: q' w( R L5 W6 l0 f, }4 k7 u7 K) [. [( N
当进行快速傅立叶变换时,算法要求输入数据的顺序是按照特定的方式排列的。特别是在使用基于分治法的算法(如Cooley-Tukey算法)时,输入数据的顺序必须满足按照一定规律的排列。& F3 t! o Y. X( `
) c3 e& K5 `9 m: }$ E) M7 A- A
在实际的FFT实现中,最常见的方式是通过比特反转来重新排列输入数据。比特反转就是将输入数据的比特位(二进制位)的顺序进行颠倒。这是因为在FFT算法中,数据会被分组,并按照一定规则进行反转,以便在每个阶段的运算中,数据可以正确地与其它组合进行配对。* A2 X% t% G* ^5 _+ \
% j" W8 L+ N/ x) O6 g6 @举个简单的例子,假设有一个长度为8的数据序列,按照0到7的顺序排列:
" V& Q0 U9 d6 I" b& u7 t& M/ e* w5 Q- ^: L3 D+ b1 B2 |
0 1 2 3 4 5 6 7' V# d |' i* }7 G0 `
/ k0 k6 C" y2 f- t# L. L' A) R, C
在进行FFT时,需要按照一定规则重新排列这些数据。比特反转操作将会对这个数据序列进行如下的重新排列:
" c" b0 k6 s8 U* X+ V; {6 c8 F6 N; W/ g) U
0 4 2 6 1 5 3 7/ d! f L6 H6 \7 v. y) @ b9 h
3 R, v/ ?$ L7 L% M在FFT算法的每个阶段中,这种重新排列都会使得数据正确地与其它组合进行配对,从而实现快速傅立叶变换的计算。8 _5 y, I3 ~5 M4 q8 w
' P' l; Q0 c( O1 Z9 a6 N! F1 M' e进行FFT并转换为模值7 h2 Q3 a' _4 l& W' H1 z
- arm_cfft_radix4_f32(&fft_inst,inputSignal); //FFT计算
, m3 X: B `4 W6 q - arm_cmplx_mag_f32(inputSignal,outputSignal,FFT_LENGTH); //取模得幅值
复制代码
) h1 |, D4 v6 D对输入数组进行FFT变换,并将FFT的结果转化为模值。2 \5 `4 t r/ G) I
, I3 p$ C+ p m% ~测试" f- {% @$ n, A
我们进行一个简单的测试9 g( B& m: ~ k; `# s7 M
- #define FFT_SIZE 1024# s# w* V3 @( t( {, b* r5 m3 w
- #define SAMPLE_RATE 1000' i* f, a; [ l5 G4 s5 C! [
- #define NUM_SAMPLES 1000
7 Q u N0 d: C5 }7 h0 z# B - #define FREQ_OF_INTEREST 100: R3 x8 J; b* t1 U9 Z/ \
- for (int i = 0; i < NUM_SAMPLES; i++) {6 ^7 |0 C0 f- D* U1 ~) e
- float32_t t = (float32_t)i / SAMPLE_RATE;
' g2 ?- e" N) _% Q' h - float32_t sin_value = sinf(2 * PI * FREQ_OF_INTEREST * t); // 计算正弦波值3 B( M- l* q' |$ ]; O' i9 H- n
- inputSignal[i * 2] = sin_value; // 实部
4 S9 M; h- O+ i - inputSignal[i * 2 + 1] = 0; // 虚部! V4 ?! T% b2 }. c: y6 N
- }
复制代码
/ Z. \2 p% |# F2 e$ r/ W% P一千个点的采样值,频率假设为100HZ作为输入信号。4 t( d! [4 J, X
- for (int i = 0; i < FFT_SIZE; i++) {, k7 h, H! f/ `6 S, h
- // 计算复数的模值
5 ]5 R/ N( z5 j. b* o6 R4 O - float32_t real = inputSignal[2 * i];
E K6 L4 a( Z! O, [: g$ v - float32_t imag = inputSignal[2 * i + 1];
: U' [3 n( e: ^0 U% f/ u5 M# z l5 F6 |( I - float32_t magnitude = sqrtf(real * real + imag * imag);
+ g' X# N. C- E2 x& C - 3 N4 U% v6 x* [
- // 打印每个频率分量的模值5 h) m3 ~( @7 D( v7 A
- printf("Magnitude: %f\n", magnitude);
( N8 x q5 R/ g& y7 T - }
复制代码
3 Z8 _5 v( V4 r进行傅里叶变换后打印模值。6 }1 e5 f y7 z. S$ _
) H" r1 Q# ^! u. D/ q0 O
5 m! o1 G- x# ?# S3 O
: Y9 |6 j8 `3 T5 E! H3 V6 O可以看到傅里叶变换执行成功。
$ u) v' {0 q6 `; \' n' H- for (int i = 0; i < NUM_SAMPLES; i++) {
( t. {2 g8 g9 ~+ W G( n - float32_t t = (float32_t)i / SAMPLE_RATE;1 y3 t6 B, ^" V+ v ~& X0 n
- float32_t sin_value = sinf(2 * PI * FREQ_OF_INTEREST * t)+sinf(2 * PI * FREQ_OF_INTEREST * t*2)+sinf(3*2 * PI * FREQ_OF_INTEREST * t); // 计算正弦波值+ j: B+ ?, Q9 w+ L' ]) l
- inputSignal[i * 2] = sin_value; // 实部
' j) P, G3 Q. D, ~2 n - inputSignal[i * 2 + 1] = 0; // 虚部
! o# f) H* b6 Y1 F - }
复制代码 , K/ M9 N6 {4 j5 J9 H5 |
我们将信号制作成100HZ+200HZ+300HZ的信号。
/ f! W- j7 }# r- @, o
9 c9 A0 t. c7 e3 y
1 x m6 j# A* n3 u, v
6 e- l: x! _+ m" J" W- b9 R% U5 r5 H* y( B/ f: z
9 n' E, a8 p: M转载自:电路小白7 ~1 @7 y3 A8 U% _: y: [
如有侵权请联系删除
6 t" B1 x5 h4 a+ Z
0 X0 K* `. J& @0 E |
这个FFT不错,学习参考一下