快速傅里叶变换(FFT)是一种数字信号处理中常用的技术,用于将快速序列转换为频域表示。在嵌入式系统中,如基于STM32的微控制器,实现FFT可以帮助解决信号处理的需求,例如声音处理、图像处理等。本文将介绍基于STM32的离散傅里叶变换的原理、实现方法和应用。
! t1 W+ s. x: V2 b
~! |7 Z9 k& `& m0 I @
) B6 T) [4 ^8 y- x. x: D7 J% W6 u* P& g# Z3 @! p+ l5 K9 O
FFT是一种将时域序列转换为频域表示的技术,它将一个序列的N个采样点映射到频域中N个频率分量。其数学表达式如下:" p, W, o4 ?; c$ K0 h* l
( g1 J: [/ v$ T! Q
1 _$ h9 H* U' r) Z
7 x2 b, a2 w! \) U n( Y- v/ A其中,x(n) 是输入序列,X(k) 是输出的频域表示。
2 F3 m2 p. X5 B+ Z$ z7 f
8 X* B1 s) J9 ]准备工作:
4 U9 e: w& F# a) F5 b, y( L/ K# z1 T3 M% B$ u
* e% y D2 F+ \ r$ {
0 s/ O S8 T. q" I) o
+ F4 E: S) e4 o. S. m0 n
# [( ^$ k0 w" c0 l
; F+ C* ~- R: x5 |# v
Keil中的DSP库(Digital Signal Processing Library,数字信号处理库)是针对ARM Cortex-M处理器系列的一组软件库,用于提供各种数字信号处理功能的支持。这些库提供了一系列优化过的算法,可以帮助开发人员在嵌入式系统中高效地实现音频处理、图像处理、通信系统等各种信号处理应用。& M# k% ]# W( N
- D; Y5 G- G+ w因此我们需要在Keil中安装我们的DSP库。) u4 Q7 y- k5 _) [" `! D% K, a+ H
- #include "arm_math.h" // 包含DSP库
复制代码 * |, J6 t0 F# W2 C% u* S0 M4 w; w
首先包含我们的DSP库。
. l+ v1 C' a2 `5 H# Y3 k6 O% l5 O8 g9 p* u
) ^, r8 r1 y: ~
- #define FFT_LENGTH 100' b# ]7 ?& b; c" y* [; p
- // 输入序列
3 D; D9 t6 R$ Q [' \ i - float32_t inputSignal[FFT_LENGTH*2];' {8 W4 m% D6 I2 {2 a# g- t$ j
% W. w: K: [/ i3 I3 l3 V- // 输出序列,存储变换后的结果
& J1 |, Z( |3 p1 u: v2 U) O - float32_t outputSignal[FFT_LENGTH];
复制代码
& a H) [0 Y9 x& R, Y( \, ^1 `定义FFT的的输入和输出数组还有数组长度" |" M4 Y$ e4 q- K5 f2 s
- arm_status status;
9 b( i7 t# N- v1 B0 D. [- U - arm_cfft_radix4_instance_f32 fft_inst;
, n) d2 k/ s) `) Z" D3 R y - status = arm_cfft_radix4_init_f32(&fft_inst, FFT_LENGTH,0,1);
复制代码- void arm_cfft_radix4_init_f32() g: _" v) D$ e' |, v/ R( ~% J
- arm_cfft_radix4_instance_f32 * S,7 n6 r3 T% M8 }$ D) M2 ~% k) b" }; B5 B
- uint16_t fftLen,$ a, m1 ?- s3 D* O& P
- uint8_t ifftFlag,8 X, B! y( ?. R
- uint8_t bitReverseFlag
- l5 t* N B- V7 z - );
复制代码
9 p5 e/ ~. g' y W! U& v/ N定义一个状态变量用来显示FFT的初始化是否成功。; A# K& \/ V1 f+ M! q7 ?% z9 n
定义一个FFT的配置变量。
) ~6 {' y- k; e X初始化FFT。; q7 t+ N* a* h2 R
S:指向 arm_cfft_radix4_instance_f32 结构体的指针,该结构体定义了 FFT 实例的状态信息。% P5 @: @5 m* S* Y8 C
1 U& P% E. o, `8 e; C
fftLen:FFT 的长度。$ T9 ]0 |$ P9 A( {$ }2 {8 y9 D. ^. V
. H2 V2 D, F$ B( X# Y$ Z+ D$ D2 l
ifftFlag:指定是否进行逆变换。如果为 1,则表示初始化的是逆变换的 FFT;如果为 0,则表示初始化的是正变换的 FFT。
1 l5 e3 X4 D8 ~. m, \6 L
, Q5 U A2 U3 {. IbitReverseFlag:指定是否进行比特翻转。如果为 1,则表示进行比特翻转;如果为 0,则表示不进行比特翻转。
' [, s) S' K4 v/ Q+ A9 @! F" {8 n& y6 W s! B9 n) T4 B t
在FFT算法中,比特(bit)反转是一种关键的步骤,用于将输入数据重新排列为正确的顺序,以便在后续的计算中进行有效处理。' S. @7 g" O. N, {& M) G6 S
9 o9 Q) [* ]+ K7 R& r当进行快速傅立叶变换时,算法要求输入数据的顺序是按照特定的方式排列的。特别是在使用基于分治法的算法(如Cooley-Tukey算法)时,输入数据的顺序必须满足按照一定规律的排列。
4 E0 X8 ~$ Y7 ^6 ^4 o' |) U9 H6 B
) o0 f; n$ j' r8 W+ Y& \0 ^在实际的FFT实现中,最常见的方式是通过比特反转来重新排列输入数据。比特反转就是将输入数据的比特位(二进制位)的顺序进行颠倒。这是因为在FFT算法中,数据会被分组,并按照一定规则进行反转,以便在每个阶段的运算中,数据可以正确地与其它组合进行配对。! u+ a5 z+ n x1 l7 r3 @+ { N; A
* N. Y! i4 _ ^8 i
举个简单的例子,假设有一个长度为8的数据序列,按照0到7的顺序排列:: {- v2 E8 i) X
" x3 ?6 |, c! v, b$ k0 1 2 3 4 5 6 7
+ h# u) _# M5 m3 O& N
: R6 N9 |* c6 Z3 |) d在进行FFT时,需要按照一定规则重新排列这些数据。比特反转操作将会对这个数据序列进行如下的重新排列:' G+ N. M9 p" H2 Q, U/ K5 \
" Q- k2 q. }( ?/ I2 }8 W7 U- `0 4 2 6 1 5 3 7
& i0 J6 P* K% }1 s
; O3 Q- | T% G4 }$ A( r4 V f. `5 f在FFT算法的每个阶段中,这种重新排列都会使得数据正确地与其它组合进行配对,从而实现快速傅立叶变换的计算。
" W0 J! D4 D# F% X" D, U% [: k
6 g8 U, p& p1 K9 y9 n) ~进行FFT并转换为模值
* |, U7 D& T- I5 A. ]- |& o- arm_cfft_radix4_f32(&fft_inst,inputSignal); //FFT计算
6 r; ^ o% _1 u# A- {* y1 ] - arm_cmplx_mag_f32(inputSignal,outputSignal,FFT_LENGTH); //取模得幅值
复制代码
" K- ~6 ]0 I$ D, _, K# ]对输入数组进行FFT变换,并将FFT的结果转化为模值。/ t! G: S( E/ x: r1 ]
; X$ J7 u1 Y) [3 Y v测试
- \5 ?; z7 D8 f% h7 c我们进行一个简单的测试9 e* `' v# M0 Q
- #define FFT_SIZE 1024, n& b* T% Q/ f2 t
- #define SAMPLE_RATE 1000' m3 c0 P: ]+ l# _- ^
- #define NUM_SAMPLES 1000
! u5 W9 p7 ^/ H J/ u - #define FREQ_OF_INTEREST 100- U8 H# Q8 [: z9 `. C
- for (int i = 0; i < NUM_SAMPLES; i++) {4 U4 N7 x* {& O- ^: Y
- float32_t t = (float32_t)i / SAMPLE_RATE;
- |" L2 _. o" d# M" V2 X8 O - float32_t sin_value = sinf(2 * PI * FREQ_OF_INTEREST * t); // 计算正弦波值
; y7 T* T/ R* r m Z& b L5 Q% T - inputSignal[i * 2] = sin_value; // 实部2 e+ Y8 G, z% r. D
- inputSignal[i * 2 + 1] = 0; // 虚部
8 K, J+ l/ f+ K; T3 b - }
复制代码 , L/ M* Z! \% r7 Y1 l; ?
一千个点的采样值,频率假设为100HZ作为输入信号。
6 {5 I+ F* V; F" J3 V1 C. ?, K8 N( R- for (int i = 0; i < FFT_SIZE; i++) {
- l. B2 E1 Z3 J9 Y- ^ - // 计算复数的模值
1 ^, E$ }! y M- k, C - float32_t real = inputSignal[2 * i];
_% ?7 V) s4 ]/ y - float32_t imag = inputSignal[2 * i + 1];# }: A$ T' w; k' p I \1 d
- float32_t magnitude = sqrtf(real * real + imag * imag);' I6 n; \+ y) g9 v
- 6 _1 [3 {, y, d' C1 h
- // 打印每个频率分量的模值 F" d) S7 P! t8 L1 e0 b2 N
- printf("Magnitude: %f\n", magnitude);. B5 Z5 `% {# \( d. N, R# P& o7 R
- }
复制代码
2 t+ o, D+ o5 J3 T进行傅里叶变换后打印模值。* C, h0 N; l& c4 L9 |) ~% ^
& L% F" N; ]5 U. A% ~, S
; o! n: z$ r) H
5 J6 L7 V- e0 j6 @& R
可以看到傅里叶变换执行成功。7 z4 y J. T1 _; a. `5 ^
- for (int i = 0; i < NUM_SAMPLES; i++) {
/ M' Q) W2 e( o - float32_t t = (float32_t)i / SAMPLE_RATE;5 \ A# t* t) B6 T
- 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); // 计算正弦波值
4 Q; {6 \, v8 I, U: u - inputSignal[i * 2] = sin_value; // 实部
# a; r+ P2 M6 f& Z7 ] - inputSignal[i * 2 + 1] = 0; // 虚部7 i4 J; g. h5 [& D' M
- }
复制代码
% o: S& T5 \8 c/ C M2 f% d: o" @我们将信号制作成100HZ+200HZ+300HZ的信号。
1 g# A* x- x- v4 X3 _$ Q4 ~5 b
) |8 l# A! H" P' q. q! [0 s2 C: N
( F( O0 [* u0 t9 u/ I8 N3 N8 k0 _" H, J- p4 u( l" O) Q
& j% N2 l. p" K+ r3 ]1 A
) Y3 C C0 O4 E( d: Z% I转载自:电路小白9 u( O: w9 g4 {# \" [4 v) {
如有侵权请联系删除, I$ `% ~& ~! R1 N6 z* p
8 W+ ^7 _1 z9 E9 j
|
这个FFT不错,学习参考一下