快速傅里叶变换(FFT)是一种数字信号处理中常用的技术,用于将快速序列转换为频域表示。在嵌入式系统中,如基于STM32的微控制器,实现FFT可以帮助解决信号处理的需求,例如声音处理、图像处理等。本文将介绍基于STM32的离散傅里叶变换的原理、实现方法和应用。: N$ e' P) d0 Q; r* c
( Y! g: x9 L* k( _ r
0 f) z0 Y, n8 u* W& h+ ?6 @. Q: }% c# n" v
FFT是一种将时域序列转换为频域表示的技术,它将一个序列的N个采样点映射到频域中N个频率分量。其数学表达式如下:3 B; M5 Z M% n+ A4 g0 r
4 U$ b1 ~4 n6 X2 m; _/ j$ b& b
5 Q( n% p6 C% n* w
0 X! l1 p1 I9 W* S |3 b i其中,x(n) 是输入序列,X(k) 是输出的频域表示。4 X& Z8 Q0 t9 u. @5 \
6 T$ K( x+ y4 [* R4 S% U# ?
准备工作:9 D `4 K, w2 e( t
7 w# L0 w9 X5 I) Z4 U5 g
, f: W! d9 o9 f& O( E/ p) D
1 X; k+ B1 |8 ?; q- p& D3 h
/ J9 L7 O; _% ^4 |+ O! Q$ `
0 C" L. J) V1 z6 e& i5 v
; S/ J3 Q6 O2 F" o2 @. \# z$ _& ?
Keil中的DSP库(Digital Signal Processing Library,数字信号处理库)是针对ARM Cortex-M处理器系列的一组软件库,用于提供各种数字信号处理功能的支持。这些库提供了一系列优化过的算法,可以帮助开发人员在嵌入式系统中高效地实现音频处理、图像处理、通信系统等各种信号处理应用。2 [( s( K3 Z4 g3 y
! Q3 s0 o, y$ _5 u' y. u* l
因此我们需要在Keil中安装我们的DSP库。; f# S! [ V, a# G# i8 G
- #include "arm_math.h" // 包含DSP库
复制代码
* P: C. U) j4 y( G7 x, B首先包含我们的DSP库。4 R& H. ^% g# Z" Z4 o v% {6 r
* @1 D' f0 U2 F# g; f
( m: `3 j/ p1 [+ d- #define FFT_LENGTH 1004 x8 L2 _4 U' j( B- M; r
- // 输入序列1 X0 H( Y, Y2 r. X# _5 y
- float32_t inputSignal[FFT_LENGTH*2];: u7 J- h( L3 X T/ ]8 G
) T$ t8 \# j9 ~0 E/ E2 L0 r ]- E- // 输出序列,存储变换后的结果
# Z8 r$ D2 C- H5 {- |4 Q+ N - float32_t outputSignal[FFT_LENGTH];
复制代码
# w, d, R' [! \定义FFT的的输入和输出数组还有数组长度
; w, M' u( U8 @; ]- arm_status status;! ~* M$ E' s8 I. u4 b# F
- arm_cfft_radix4_instance_f32 fft_inst;. h& V+ W( @1 M1 Q. z: P( O
- status = arm_cfft_radix4_init_f32(&fft_inst, FFT_LENGTH,0,1);
复制代码- void arm_cfft_radix4_init_f32(, e9 D4 I0 b2 D S
- arm_cfft_radix4_instance_f32 * S,$ M" h5 a- w7 m$ A7 ], t) Z# {
- uint16_t fftLen,+ O- ]; B- k0 C, g8 E2 B8 `
- uint8_t ifftFlag,
5 d* g5 R N8 S" ^ - uint8_t bitReverseFlag& `4 q& q- ^/ p* J2 \
- );
复制代码
2 g, L6 \2 l* w( {6 Q定义一个状态变量用来显示FFT的初始化是否成功。: k0 M% L* f0 D& Q
定义一个FFT的配置变量。
8 X4 n) A/ w0 I+ n! R, i5 C初始化FFT。
& ]! F: g- i' N lS:指向 arm_cfft_radix4_instance_f32 结构体的指针,该结构体定义了 FFT 实例的状态信息。
5 c$ |9 n8 S- }- |0 m3 f) Y9 b E4 p* f4 m
fftLen:FFT 的长度。
$ ]7 u' g- h' @# X' x, l9 {3 @. E6 a. ?' N/ T" ]0 A
ifftFlag:指定是否进行逆变换。如果为 1,则表示初始化的是逆变换的 FFT;如果为 0,则表示初始化的是正变换的 FFT。
: |# a" h) u; F" L4 [. m; E( H3 F. I& W: Z$ t
bitReverseFlag:指定是否进行比特翻转。如果为 1,则表示进行比特翻转;如果为 0,则表示不进行比特翻转。
: w. D4 \2 v+ N! Y; g
& u' L' }: c: g* ~在FFT算法中,比特(bit)反转是一种关键的步骤,用于将输入数据重新排列为正确的顺序,以便在后续的计算中进行有效处理。
& g! P% m F0 D* }- S% D) \- D5 |/ @* |; C
当进行快速傅立叶变换时,算法要求输入数据的顺序是按照特定的方式排列的。特别是在使用基于分治法的算法(如Cooley-Tukey算法)时,输入数据的顺序必须满足按照一定规律的排列。9 V( [6 O# z- w7 z) [9 P0 I
3 E* r. A1 F, Y$ b2 {3 }0 q
在实际的FFT实现中,最常见的方式是通过比特反转来重新排列输入数据。比特反转就是将输入数据的比特位(二进制位)的顺序进行颠倒。这是因为在FFT算法中,数据会被分组,并按照一定规则进行反转,以便在每个阶段的运算中,数据可以正确地与其它组合进行配对。
1 E4 h' ~: ?8 T/ a" ^9 I' N
, g1 C4 b# m5 K8 b% ^* A7 f2 t9 ^ ~举个简单的例子,假设有一个长度为8的数据序列,按照0到7的顺序排列:
' C* x8 ^% R/ f# @, T% o0 F# N
& H7 U8 Q& {. y0 b5 Q2 }$ L0 1 2 3 4 5 6 76 n8 A' s0 W+ v1 v' u
" G- `1 \# F3 P# X$ ~4 l- B. S9 R
在进行FFT时,需要按照一定规则重新排列这些数据。比特反转操作将会对这个数据序列进行如下的重新排列:
! N" W/ l/ B: j( [' c% N( U! `
/ F! [3 L' p( }( Q" {0 o( q( V0 4 2 6 1 5 3 7 L( o( ` h9 |0 S' s
8 Q1 N, h/ f' }7 n' C1 t在FFT算法的每个阶段中,这种重新排列都会使得数据正确地与其它组合进行配对,从而实现快速傅立叶变换的计算。0 }) {7 I6 G6 j6 ]* X
# C, H' `- Q% P进行FFT并转换为模值
4 L0 X9 }2 A5 B7 q `3 K; s, |- Z- arm_cfft_radix4_f32(&fft_inst,inputSignal); //FFT计算+ N3 K. H9 U- i1 _) V9 n: \
- arm_cmplx_mag_f32(inputSignal,outputSignal,FFT_LENGTH); //取模得幅值
复制代码
! E4 e2 x) C3 [1 g对输入数组进行FFT变换,并将FFT的结果转化为模值。
. R: x) y8 J/ q9 @ N$ d: X
/ N; l) J. E6 N7 u8 j1 P! V测试# W0 j. {. b, x" W& Q, O
我们进行一个简单的测试
# K/ |6 ~5 H0 Z- #define FFT_SIZE 1024
$ r" F) ~' d; x - #define SAMPLE_RATE 1000
) R' |5 {1 V7 b* y* t - #define NUM_SAMPLES 1000+ O( k% ^" E p" u
- #define FREQ_OF_INTEREST 100$ T, `3 M# z3 p2 B$ @, V
- for (int i = 0; i < NUM_SAMPLES; i++) {
6 A4 J) S& | ~9 D - float32_t t = (float32_t)i / SAMPLE_RATE;
4 Q0 u3 u6 ]9 I0 I( z% [3 } - float32_t sin_value = sinf(2 * PI * FREQ_OF_INTEREST * t); // 计算正弦波值; [: g$ s: R( i
- inputSignal[i * 2] = sin_value; // 实部5 {/ E/ U1 W% _' @9 g
- inputSignal[i * 2 + 1] = 0; // 虚部
c _- J& {; R1 p - }
复制代码 & h+ A8 \9 a x) M. }1 B
一千个点的采样值,频率假设为100HZ作为输入信号。
; f6 Q$ `+ x) K9 z* H- for (int i = 0; i < FFT_SIZE; i++) {
% |* h4 i2 u2 \ - // 计算复数的模值% |8 W+ T }+ `0 m; ^
- float32_t real = inputSignal[2 * i];
+ {3 i \2 C" s! `+ o0 v - float32_t imag = inputSignal[2 * i + 1];# J+ g1 y5 x( X v) c5 `( `4 m2 V; W# d
- float32_t magnitude = sqrtf(real * real + imag * imag);
+ _; u3 m2 c3 I2 V* x - 2 [) }" x8 W* [6 U# K1 ]
- // 打印每个频率分量的模值
1 C3 K, H4 T% y+ a. C" Z, j - printf("Magnitude: %f\n", magnitude);2 o3 u1 ~- a3 T9 b# P
- }
复制代码 ! r' X1 V& U4 ~1 q0 |* L
进行傅里叶变换后打印模值。
2 O8 K. T; L( e4 k! \% U. @
9 i/ I3 l: i0 T& M' G5 m9 c1 }4 B/ I+ X
1 f( [4 X' T5 ^; v2 h( i
1 I1 J4 M( [- M可以看到傅里叶变换执行成功。
# x5 q2 H* M" Z! y' B- F6 c- for (int i = 0; i < NUM_SAMPLES; i++) {# k+ L8 h1 }3 ]4 ?7 s, ?6 V$ n# j
- float32_t t = (float32_t)i / SAMPLE_RATE;
6 q1 l. [" R O+ G - 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); // 计算正弦波值
9 O4 s! l( j; k# K - inputSignal[i * 2] = sin_value; // 实部5 t1 z5 I+ \7 l5 \+ d L! v
- inputSignal[i * 2 + 1] = 0; // 虚部
6 W2 j+ D* B5 f2 d* g+ c! L - }
复制代码
- u: |; x3 T8 Q; Y* o) m我们将信号制作成100HZ+200HZ+300HZ的信号。
# s* G6 D) X0 F$ O3 d& c
3 r. B6 n- A. z2 r) K7 ]+ {: C
* {: t+ T, Z8 ^( J
4 a0 {: q7 u O& x, q8 l7 o+ M) ^: R; r$ ^
4 B' F) g" Y+ ]6 P; V
转载自:电路小白 q& @% W7 S. o) l7 d3 T9 z
如有侵权请联系删除4 ^) P7 m. W8 ?4 ^: l1 ]
1 C# h8 W, c- p7 f' E0 O7 g2 S0 ]
|
这个FFT不错,学习参考一下