快速傅里叶变换(FFT)是一种数字信号处理中常用的技术,用于将快速序列转换为频域表示。在嵌入式系统中,如基于STM32的微控制器,实现FFT可以帮助解决信号处理的需求,例如声音处理、图像处理等。本文将介绍基于STM32的离散傅里叶变换的原理、实现方法和应用。. e |2 ^& c* E( L! b
q) y# Y# p0 T! @; T" u. H/ c) q
+ L ~: U( e$ E5 o" ~
( o. |, s2 a$ M1 d- n
FFT是一种将时域序列转换为频域表示的技术,它将一个序列的N个采样点映射到频域中N个频率分量。其数学表达式如下:! o: V: w Z5 x5 r
$ J9 U. `7 c8 i. e( W6 B% X X# p
* h: N' a; e8 q. M: _% E
6 \6 \8 q+ w- ~0 }
其中,x(n) 是输入序列,X(k) 是输出的频域表示。
4 y4 d! D6 P. m( J1 e# A8 d. Z3 D0 r n; B
准备工作:
9 L5 ^2 F: f c0 [/ A' r0 ]7 B2 l
* J9 J0 u8 l4 W) c4 z8 b7 Q* I" C4 M9 r# d% }8 C! t% V
. y/ C1 g. G0 _
3 t6 B U4 @0 M; Y+ u! i
2 l8 q8 |6 x2 r
' P* U% H, V/ q( m, M/ b6 j! H/ H( b
Keil中的DSP库(Digital Signal Processing Library,数字信号处理库)是针对ARM Cortex-M处理器系列的一组软件库,用于提供各种数字信号处理功能的支持。这些库提供了一系列优化过的算法,可以帮助开发人员在嵌入式系统中高效地实现音频处理、图像处理、通信系统等各种信号处理应用。+ f! \ ~8 w6 @% ?0 g5 Q& H$ F
9 H1 y( @5 h# J
因此我们需要在Keil中安装我们的DSP库。' W7 D0 P$ i/ ] W. W1 r
- #include "arm_math.h" // 包含DSP库
复制代码 " E; B Z8 e( ^7 m4 \
首先包含我们的DSP库。# U3 \& l3 o7 t5 V3 C& p3 H
$ v. q! f/ t3 i5 i5 }. K P
: j" T: @, F) _5 W( R- ~/ d- #define FFT_LENGTH 100
2 y: ~4 |5 R& |# i; e/ A$ c9 F - // 输入序列( I! L! G8 n- M0 |
- float32_t inputSignal[FFT_LENGTH*2];9 z' k+ ^3 u) w( ? t+ K
6 ~9 H. B* a% v# z; p- // 输出序列,存储变换后的结果+ q9 S4 y+ ?: `. z" O% E
- float32_t outputSignal[FFT_LENGTH];
复制代码
5 e1 R/ Y1 B: u- R" r2 {定义FFT的的输入和输出数组还有数组长度/ i1 w2 w" U1 L4 w
- arm_status status;" u4 m5 n( R0 w( O4 e7 I( W) @4 D
- arm_cfft_radix4_instance_f32 fft_inst;8 _- B% R8 _# |( y5 w3 ^1 H, \: {' G& z3 b
- status = arm_cfft_radix4_init_f32(&fft_inst, FFT_LENGTH,0,1);
复制代码- void arm_cfft_radix4_init_f32(( K, S3 t2 V2 T1 e* G
- arm_cfft_radix4_instance_f32 * S,
$ C1 _$ v: B# p8 A - uint16_t fftLen,
! G& n G. T# ?( j - uint8_t ifftFlag,! F1 k6 _' N& [6 V- [4 I
- uint8_t bitReverseFlag' X1 ]( l u/ j0 f1 z& a6 @ S
- );
复制代码 ; u2 M/ s& V0 X) t3 p ~2 G
定义一个状态变量用来显示FFT的初始化是否成功。
/ P6 O, ^: f8 V! z6 u# g9 }7 M定义一个FFT的配置变量。4 s9 X9 N) W3 u. x0 Y* N
初始化FFT。8 N% }0 J0 L7 k! q4 c& R6 D! V
S:指向 arm_cfft_radix4_instance_f32 结构体的指针,该结构体定义了 FFT 实例的状态信息。4 f* c9 \* |% S& p2 u" A
( X( V s) w! o) Q& F* I! U% x, Q
fftLen:FFT 的长度。( }" u$ v* D/ o
+ n7 k6 w" ?# @7 b7 HifftFlag:指定是否进行逆变换。如果为 1,则表示初始化的是逆变换的 FFT;如果为 0,则表示初始化的是正变换的 FFT。
! K- x6 }: t6 o) S3 F J, R2 a _& K6 M, s& t( d$ F) q
bitReverseFlag:指定是否进行比特翻转。如果为 1,则表示进行比特翻转;如果为 0,则表示不进行比特翻转。" b2 |0 f% Y% s8 A: S
. C# `/ F. K9 y* T$ p- D% @5 B5 E
在FFT算法中,比特(bit)反转是一种关键的步骤,用于将输入数据重新排列为正确的顺序,以便在后续的计算中进行有效处理。
/ Z; ?( E8 Z3 b, C, v( n) y
* x7 g% o5 D# z- L2 M当进行快速傅立叶变换时,算法要求输入数据的顺序是按照特定的方式排列的。特别是在使用基于分治法的算法(如Cooley-Tukey算法)时,输入数据的顺序必须满足按照一定规律的排列。
$ n* F* B7 l0 J7 Z2 c3 t7 k% Y9 n/ M8 [% M9 V% {
在实际的FFT实现中,最常见的方式是通过比特反转来重新排列输入数据。比特反转就是将输入数据的比特位(二进制位)的顺序进行颠倒。这是因为在FFT算法中,数据会被分组,并按照一定规则进行反转,以便在每个阶段的运算中,数据可以正确地与其它组合进行配对。
$ V- `# E5 E2 K; D; ^/ D; ^, h) T! _8 ^( h- ]: a y, K7 d
举个简单的例子,假设有一个长度为8的数据序列,按照0到7的顺序排列:
& @$ ~8 l+ g1 J
* y8 ~. j, w! \0 1 2 3 4 5 6 73 X( ~: ?& B2 |3 u% S) M' d. F
9 q9 ^% i( l: G o7 `在进行FFT时,需要按照一定规则重新排列这些数据。比特反转操作将会对这个数据序列进行如下的重新排列:# S2 v' n. d& ?; d% r$ A8 Y& c
! v& @. ?; X/ K9 a; ], w0 4 2 6 1 5 3 74 y: B8 l) \0 O% `
s; @( z2 V' Z( K E& v. Q
在FFT算法的每个阶段中,这种重新排列都会使得数据正确地与其它组合进行配对,从而实现快速傅立叶变换的计算。 j/ ?" T, p2 e: \; e
( {$ B3 ?! x) _- \进行FFT并转换为模值
9 I8 o6 [. y, B! N' a- arm_cfft_radix4_f32(&fft_inst,inputSignal); //FFT计算/ {" W: J; [/ Y& \5 g
- arm_cmplx_mag_f32(inputSignal,outputSignal,FFT_LENGTH); //取模得幅值
复制代码 - n0 k# D K6 S
对输入数组进行FFT变换,并将FFT的结果转化为模值。
8 b/ b' s H# V- L9 R, D' b4 n7 `) m+ a3 r- X6 A' C
测试
0 J" [5 w* r6 ]' U' N: \我们进行一个简单的测试' ^3 J( Q- D! k. j
- #define FFT_SIZE 1024
" S+ p& H( {' p' c - #define SAMPLE_RATE 1000
, `2 e; ]% N8 k7 z) y2 z2 ~# A - #define NUM_SAMPLES 10006 ]7 a+ o5 K4 a9 |
- #define FREQ_OF_INTEREST 100
5 k" q& L/ _7 Y! k. D. |! I - for (int i = 0; i < NUM_SAMPLES; i++) {+ f! h. s, d1 C/ T" z
- float32_t t = (float32_t)i / SAMPLE_RATE;- B( a3 A4 y: Z
- float32_t sin_value = sinf(2 * PI * FREQ_OF_INTEREST * t); // 计算正弦波值$ Y; M: m" M7 R8 G
- inputSignal[i * 2] = sin_value; // 实部0 N( W& m! o9 c
- inputSignal[i * 2 + 1] = 0; // 虚部' ~# M" C( w( I9 D
- }
复制代码
4 I. n6 m2 ^8 H# `4 O, o$ i, z一千个点的采样值,频率假设为100HZ作为输入信号。0 G( s6 ~" a' t
- for (int i = 0; i < FFT_SIZE; i++) {
* V. B& t. \5 k% R2 t/ P - // 计算复数的模值9 K A% x2 x1 Q- V0 F4 }
- float32_t real = inputSignal[2 * i];* x ^* A6 D. N
- float32_t imag = inputSignal[2 * i + 1];( O& p( i# W" @2 r, F
- float32_t magnitude = sqrtf(real * real + imag * imag);3 u# d' S, G& l3 z
- 3 U- D. v1 [( Z& [# [. n" t1 x
- // 打印每个频率分量的模值
' P1 g3 Q8 w( Z$ p0 z& H, U - printf("Magnitude: %f\n", magnitude);8 Z% }* z' [: _3 h4 K1 E% s$ q
- }
复制代码
- N. O0 U+ T. i: n" m进行傅里叶变换后打印模值。
% f# Q/ h! _* ?7 {9 a
/ b$ v" a/ G" ]
( \7 S9 V" {( |4 N+ i( d/ q5 h: g
f7 s( q0 H& E: q, X4 X4 }
可以看到傅里叶变换执行成功。3 E8 p" b3 e K- a2 p9 A$ Q+ ]* |
- for (int i = 0; i < NUM_SAMPLES; i++) {7 K: `0 m- t- r5 q
- float32_t t = (float32_t)i / SAMPLE_RATE;
. R$ ~* m* @' n! G3 s) L( b - 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); // 计算正弦波值$ ]" L, D3 _ v* N
- inputSignal[i * 2] = sin_value; // 实部
/ ~7 y H- Z- r, X - inputSignal[i * 2 + 1] = 0; // 虚部
5 q* p6 u2 {( v - }
复制代码
) O5 s/ H+ S0 f5 Y/ e我们将信号制作成100HZ+200HZ+300HZ的信号。
4 @5 P0 J; }" u3 o) F$ D7 ^' `& L3 S7 T
* v1 L$ E' ?+ q- u& F3 u- x- Q
7 E3 i9 {- m3 e) d3 u9 ^
' h% V) U3 @! c$ M3 Z% p! c& q3 O+ h6 s/ A
转载自:电路小白
- o+ Q3 i. j9 I' W1 _& h1 k如有侵权请联系删除
5 H' g ~9 H: i# o
5 x: M0 h) k, W5 U$ h/ u |
这个FFT不错,学习参考一下