CMSIS 的反向 FFT 是错误的
Inverse FFT with CMSIS is wrong
我正在尝试对信号执行 FFT,并使用生成的数据通过 IFFT 检索原始样本。我在带 M3 的 STM32 上使用 CMSIS DSP 库。
我的问题是了解 FFT 发生的缩放,以及如何获得正确的 IFFT。目前,IFFT 产生与输入相似的波形,但点缩放比例在原始波形的 120x-140x 之间。这仅仅是 q15 精度误差的结果吗?我是否也将 IFFT 结果缩放了 7 位?我的代码如下
文档还提到“对于 RIFFT,源缓冲区的长度必须至少为 fftLenReal + 2。最后两个元素必须等于 RFFT 生成的元素:(pSrc[0] - pSrc[ 1]) >> 1 和 0"。这个是来做什么的?将这些操作分别应用于 FFT_SIZE2 - 2 和 FFT_SIZE2 - 1 根本没有改变 IFFT 的结果。
//128 point FFT
#define FFT_SIZE 128
arm_rfft_instance_q15 fft_instance;
arm_rfft_instance_q15 ifft_instance;
//time domain signal buffers
float32_t sinetbl_in[FFT_SIZE];
float32_t sinetbl_out[FFT_SIZE];
//a copy for comparison after RFFT since function modifies input buffer
volatile q15_t fft_in_buf_cpy[FFT_SIZE];
q15_t fft_in_buf[FFT_SIZE];
//output for FFT, RFFT provides real and complex data organized as re[0], im[0], re[1], im[1]
q15_t fft_out_buf[FFT_SIZE*2];
q15_t fft_out_buf_mag[FFT_SIZE*2];
//inverse fft buffer result
q15_t ifft_out_buf[FFT_SIZE];
//generate 1kHz sinewave with a sample frequency of 8kHz for 128 samples, amplitude is 1
for(int i = 0; i < FFT_SIZE; ++i){
sinetbl_in[i] = arm_sin_f32(2*3.14*1000 *i/8000);
sinetbl_out[i] = 0;
}
//convert buffer to q15 (not enough flash to use f32 fft functions)
arm_float_to_q15(sinetbl_in, fft_in_buf, FFT_SIZE);
memcpy(fft_in_buf_cpy, fft_in_buf, FFT_SIZE*2);
//perform RFFT
arm_rfft_init_q15(&fft_instance, FFT_SIZE, 0, 1);
arm_rfft_q15(&fft_instance, fft_in_buf, fft_out_buf);
//calculate magnitude, skip 1st real and img numbers as they are DC and both real
arm_cmplx_mag_q15(fft_out_buf + 2, fft_out_buf_mag + 1, FFT_SIZE/2-1);
//weird operations described by documentation, does not change results
//fft_out_buf[FFT_SIZE*2 - 2] = (fft_out_buf[0] - fft_out_buf[1]) >> 1;
//fft_out_buf[FFT_SIZE*2 - 1] = 0;
//perform inverse FFT
arm_rfft_init_q15(&ifft_instance, FFT_SIZE, 1, 1);
arm_rfft_q15(&ifft_instance, fft_out_buf, ifft_out_buf);
//closest approximation to get to original scaling
//arm_shift_q15(ifft_out_buf, 7, ifft_out_buf, FFT_SIZE);
//convert back to float for comparison with input
arm_q15_to_float(ifft_out_buf, sinetbl_out, FFT_SIZE);
我觉得我用精确的评论回答了我自己的问题,但我想确定一下。我做这个 FFT 的东西对吗?
提前致谢
正如 Cris 所指出的,一些图书馆跳过了规范化过程。 CMSIS DSP 是这些库之一,因为它的目的是要快。对于 CMSIS,根据 FFT 大小,您必须将数据左移一定量才能返回到原始范围。在我的情况下,FFT 大小为 128 并且幅度计算也是我最初推测的 7。
我正在尝试对信号执行 FFT,并使用生成的数据通过 IFFT 检索原始样本。我在带 M3 的 STM32 上使用 CMSIS DSP 库。
我的问题是了解 FFT 发生的缩放,以及如何获得正确的 IFFT。目前,IFFT 产生与输入相似的波形,但点缩放比例在原始波形的 120x-140x 之间。这仅仅是 q15 精度误差的结果吗?我是否也将 IFFT 结果缩放了 7 位?我的代码如下
文档还提到“对于 RIFFT,源缓冲区的长度必须至少为 fftLenReal + 2。最后两个元素必须等于 RFFT 生成的元素:(pSrc[0] - pSrc[ 1]) >> 1 和 0"。这个是来做什么的?将这些操作分别应用于 FFT_SIZE2 - 2 和 FFT_SIZE2 - 1 根本没有改变 IFFT 的结果。
//128 point FFT
#define FFT_SIZE 128
arm_rfft_instance_q15 fft_instance;
arm_rfft_instance_q15 ifft_instance;
//time domain signal buffers
float32_t sinetbl_in[FFT_SIZE];
float32_t sinetbl_out[FFT_SIZE];
//a copy for comparison after RFFT since function modifies input buffer
volatile q15_t fft_in_buf_cpy[FFT_SIZE];
q15_t fft_in_buf[FFT_SIZE];
//output for FFT, RFFT provides real and complex data organized as re[0], im[0], re[1], im[1]
q15_t fft_out_buf[FFT_SIZE*2];
q15_t fft_out_buf_mag[FFT_SIZE*2];
//inverse fft buffer result
q15_t ifft_out_buf[FFT_SIZE];
//generate 1kHz sinewave with a sample frequency of 8kHz for 128 samples, amplitude is 1
for(int i = 0; i < FFT_SIZE; ++i){
sinetbl_in[i] = arm_sin_f32(2*3.14*1000 *i/8000);
sinetbl_out[i] = 0;
}
//convert buffer to q15 (not enough flash to use f32 fft functions)
arm_float_to_q15(sinetbl_in, fft_in_buf, FFT_SIZE);
memcpy(fft_in_buf_cpy, fft_in_buf, FFT_SIZE*2);
//perform RFFT
arm_rfft_init_q15(&fft_instance, FFT_SIZE, 0, 1);
arm_rfft_q15(&fft_instance, fft_in_buf, fft_out_buf);
//calculate magnitude, skip 1st real and img numbers as they are DC and both real
arm_cmplx_mag_q15(fft_out_buf + 2, fft_out_buf_mag + 1, FFT_SIZE/2-1);
//weird operations described by documentation, does not change results
//fft_out_buf[FFT_SIZE*2 - 2] = (fft_out_buf[0] - fft_out_buf[1]) >> 1;
//fft_out_buf[FFT_SIZE*2 - 1] = 0;
//perform inverse FFT
arm_rfft_init_q15(&ifft_instance, FFT_SIZE, 1, 1);
arm_rfft_q15(&ifft_instance, fft_out_buf, ifft_out_buf);
//closest approximation to get to original scaling
//arm_shift_q15(ifft_out_buf, 7, ifft_out_buf, FFT_SIZE);
//convert back to float for comparison with input
arm_q15_to_float(ifft_out_buf, sinetbl_out, FFT_SIZE);
我觉得我用精确的评论回答了我自己的问题,但我想确定一下。我做这个 FFT 的东西对吗? 提前致谢
正如 Cris 所指出的,一些图书馆跳过了规范化过程。 CMSIS DSP 是这些库之一,因为它的目的是要快。对于 CMSIS,根据 FFT 大小,您必须将数据左移一定量才能返回到原始范围。在我的情况下,FFT 大小为 128 并且幅度计算也是我最初推测的 7。