快速傅里叶变换(FFT)是一种数字信号处理中常用的技术,用于将快速序列转换为频域表示。在嵌入式系统中,如基于STM32的微控制器,实现FFT可以帮助解决信号处理的需求,例如声音处理、图像处理等。本文将介绍基于STM32的离散傅里叶变换的原理、实现方法和应用。
1 M9 [% M. G; a6 E4 c+ a. d
/ _: U5 w: E- Y' \8 K
- o$ e, G- P% ^$ Z0 v: N' U
* A6 l3 I, m& A3 l6 LFFT是一种将时域序列转换为频域表示的技术,它将一个序列的N个采样点映射到频域中N个频率分量。其数学表达式如下:
( ^ O4 c! l* N0 h
: E C2 A: ~- @$ ^$ {4 e
4 c0 d9 `0 q/ o1 W# h" q
( b) t, J& b; w' W其中,x(n) 是输入序列,X(k) 是输出的频域表示。- a' A0 e& e' o
, v# c; m6 |0 | B" B l( \- C
准备工作:
9 s/ N9 q4 b7 ^4 ?' c' s7 v8 i9 T5 s5 F% ?0 f
7 r* ?. x* K3 K" ~) k3 m
# ]/ K" a1 K8 C/ z: C1 y
6 q- m* p+ A0 r+ A' S u# s
' Y/ c4 t' E; }# l
( O8 { \5 z4 V, uKeil中的DSP库(Digital Signal Processing Library,数字信号处理库)是针对ARM Cortex-M处理器系列的一组软件库,用于提供各种数字信号处理功能的支持。这些库提供了一系列优化过的算法,可以帮助开发人员在嵌入式系统中高效地实现音频处理、图像处理、通信系统等各种信号处理应用。" f; x; F U% p) U0 M
P+ h5 P) X3 Z5 |* p4 {因此我们需要在Keil中安装我们的DSP库。
, J; D& w- E" k! E2 Y6 R2 } K- #include "arm_math.h" // 包含DSP库
复制代码
% Y5 U# R" ?! Y4 E- Y3 {首先包含我们的DSP库。. m# G" {2 r+ g U6 O1 n
- | J6 G4 {8 X. N0 M
* p4 G9 I2 V% o, K
- #define FFT_LENGTH 100% z( K4 F* T2 n5 t4 P
- // 输入序列: V" W6 W& l* u9 {. }1 {0 P8 v" h
- float32_t inputSignal[FFT_LENGTH*2];3 ~+ c8 f( @5 k
H4 {& R+ o" D% d, c- // 输出序列,存储变换后的结果
5 w' \# j: C% w/ i' g9 I - float32_t outputSignal[FFT_LENGTH];
复制代码
. q5 @" Y; t4 ^4 @4 ]定义FFT的的输入和输出数组还有数组长度: [* b0 N+ @0 q# I s5 G
- arm_status status;( a6 W9 m F" {# W, D
- arm_cfft_radix4_instance_f32 fft_inst;3 o2 u7 S* p r
- status = arm_cfft_radix4_init_f32(&fft_inst, FFT_LENGTH,0,1);
复制代码- void arm_cfft_radix4_init_f32(( Z1 G( A. K6 L, w, E0 n/ M- p
- arm_cfft_radix4_instance_f32 * S,2 F3 l( _' ~8 i' a j8 [& k4 @
- uint16_t fftLen,; s6 L7 d, h2 q6 Y, v( r4 a: {4 I
- uint8_t ifftFlag,
$ A |3 j; l, u& ^) [9 m - uint8_t bitReverseFlag
6 ?' C; S2 k& E: B; ^4 H - );
复制代码
6 q! n5 [, u) A0 H定义一个状态变量用来显示FFT的初始化是否成功。
" b# a0 b" R6 D; O定义一个FFT的配置变量。
( m2 ]+ V: i+ x# q初始化FFT。
$ }) x5 S* C' D* E# z2 _4 C, Q# fS:指向 arm_cfft_radix4_instance_f32 结构体的指针,该结构体定义了 FFT 实例的状态信息。: P% |, [& S) g% n& a
: E# V& A' X/ G5 Q& t) H2 OfftLen:FFT 的长度。7 w# c* a. s' \. P. H v: @
& o: D4 v6 _: u: Z9 ^5 y0 W% GifftFlag:指定是否进行逆变换。如果为 1,则表示初始化的是逆变换的 FFT;如果为 0,则表示初始化的是正变换的 FFT。
6 @6 a+ b W9 F i
0 A! h* T+ N' c3 DbitReverseFlag:指定是否进行比特翻转。如果为 1,则表示进行比特翻转;如果为 0,则表示不进行比特翻转。1 w/ q% |* E2 g% o6 y: a1 m
: Y# M' C! L0 H0 m% X1 f2 |0 Z1 W
在FFT算法中,比特(bit)反转是一种关键的步骤,用于将输入数据重新排列为正确的顺序,以便在后续的计算中进行有效处理。
* R' ?4 \' e" q. {
: }3 O& Z2 C; {# W$ ^( N当进行快速傅立叶变换时,算法要求输入数据的顺序是按照特定的方式排列的。特别是在使用基于分治法的算法(如Cooley-Tukey算法)时,输入数据的顺序必须满足按照一定规律的排列。
% @1 P: {, J, e/ m# u+ v! M* U. H) q- i8 H/ R# b" E% L( h
在实际的FFT实现中,最常见的方式是通过比特反转来重新排列输入数据。比特反转就是将输入数据的比特位(二进制位)的顺序进行颠倒。这是因为在FFT算法中,数据会被分组,并按照一定规则进行反转,以便在每个阶段的运算中,数据可以正确地与其它组合进行配对。* W A2 k# V& f, ~$ ~+ C9 X
# K& [) L6 Q; A# r9 Z1 f
举个简单的例子,假设有一个长度为8的数据序列,按照0到7的顺序排列:# f1 _. K- H% N4 {
* V- h( X3 w) Z# X M6 I1 E- K
0 1 2 3 4 5 6 7
/ x6 B/ T9 D& O" l0 \& K5 n
- E6 N( r' l- J' x在进行FFT时,需要按照一定规则重新排列这些数据。比特反转操作将会对这个数据序列进行如下的重新排列:
3 P( h6 b" j* E& C1 F
4 l" f- `) S* \3 }, p$ E# R( E0 4 2 6 1 5 3 7
0 }: s1 H% Y- {4 G$ f7 H) R
0 X q8 ~7 n, a. w在FFT算法的每个阶段中,这种重新排列都会使得数据正确地与其它组合进行配对,从而实现快速傅立叶变换的计算。
: d: I) W7 z2 _4 \) U( y) \ y$ @ n, b! Z% ?7 B
进行FFT并转换为模值
5 T4 S1 ^% [ j4 c, {% ~- arm_cfft_radix4_f32(&fft_inst,inputSignal); //FFT计算) O! L& u) h5 c4 ~: k. U
- arm_cmplx_mag_f32(inputSignal,outputSignal,FFT_LENGTH); //取模得幅值
复制代码 : O+ Y( a6 C8 b1 y* I
对输入数组进行FFT变换,并将FFT的结果转化为模值。
* G/ T0 \7 {0 z! N+ ^* m" }& q" K! n! L7 H1 {" w6 F
测试
, L7 [6 c2 V* b我们进行一个简单的测试
( `8 f1 T- f, ~) s2 M! d* l- #define FFT_SIZE 1024
. X8 ?( R. ~. H+ P' ^1 w7 G) @ - #define SAMPLE_RATE 1000, T6 j3 L! Z, h& Z
- #define NUM_SAMPLES 1000
+ R& `* C+ p( i; M - #define FREQ_OF_INTEREST 100
2 Y: u$ k/ [6 }5 j+ ?6 w$ n - for (int i = 0; i < NUM_SAMPLES; i++) {$ A7 Z: V- |" f) r2 F8 w# @
- float32_t t = (float32_t)i / SAMPLE_RATE;$ m0 q( d# m3 D3 D0 u) Z# E$ G
- float32_t sin_value = sinf(2 * PI * FREQ_OF_INTEREST * t); // 计算正弦波值
" q+ Y6 B1 x- N0 P - inputSignal[i * 2] = sin_value; // 实部% [- \8 `3 W. z7 c' \
- inputSignal[i * 2 + 1] = 0; // 虚部
! y( F! W% K, j: r5 n4 i1 Y - }
复制代码 3 k0 |: N& b( X+ B9 | {* e3 \
一千个点的采样值,频率假设为100HZ作为输入信号。
: U1 ~# m$ R+ A! k& S1 d- for (int i = 0; i < FFT_SIZE; i++) {& N! p% H* L$ m5 Y9 Z
- // 计算复数的模值
: k* ~& O" a8 m9 E2 N - float32_t real = inputSignal[2 * i];, ]. k" `- t. N4 l
- float32_t imag = inputSignal[2 * i + 1];0 {# e; H7 A& c
- float32_t magnitude = sqrtf(real * real + imag * imag);
5 N$ g" W! [& ` ]. v5 u -
& l7 G) Q7 |2 n: a) l - // 打印每个频率分量的模值
7 w; l9 K5 x8 A! I6 M - printf("Magnitude: %f\n", magnitude);% ~; h3 n0 t+ W5 j P
- }
复制代码
2 U1 L: J6 z; r1 R3 V进行傅里叶变换后打印模值。
$ T! n0 N- t+ R
4 H# H, |5 {. h4 p# _- E, b* m7 D
1 M4 M8 m, y6 \8 c# I& U) j
2 q5 N! F9 c; L$ b7 F( R l
可以看到傅里叶变换执行成功。
! `; |6 K: f: ?* o% }1 b9 q- for (int i = 0; i < NUM_SAMPLES; i++) {9 w0 k% E' Y5 ~% L6 }! c
- float32_t t = (float32_t)i / SAMPLE_RATE;
8 Q7 `& w8 f( k$ Z4 j1 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); // 计算正弦波值
: Z$ y, O* e: H k+ g5 V- T - inputSignal[i * 2] = sin_value; // 实部
* w" z9 D* l, a+ K2 [. R - inputSignal[i * 2 + 1] = 0; // 虚部( K) G; j6 k- ?; ]) G# \! s
- }
复制代码 % r$ n. {9 }. S0 J) ~$ m( ]
我们将信号制作成100HZ+200HZ+300HZ的信号。
! A/ k8 V# g) w% c+ G/ A# G
( E4 u9 P$ g) |- `( M
$ y4 c% S. z2 I3 `9 M+ z! g0 `! J6 S# l, R* {. b i1 j5 J; u
& n7 K' _6 w) B% B2 s
2 Z5 Q4 P& y6 f9 s- Q
转载自:电路小白
( v' {; Q3 J# N8 a" f& N7 W+ j如有侵权请联系删除
9 h- a2 x2 O- y9 ]! f* W- {# ^- q8 N
! f% g: O* f/ m! w' O |
这个FFT不错,学习参考一下