硬汉嵌入式论坛

 找回密码
 立即注册
查看: 3914|回复: 3
收起左侧

[DSP] 【安富莱DSP教程】第30章 复数FFT的实现

[复制链接]

740

主题

1326

回帖

3546

积分

管理员

春暖花开

Rank: 9Rank: 9Rank: 9

积分
3546
QQ
发表于 2015-4-15 10:26:36 | 显示全部楼层 |阅读模式
特别说明:完整45期数字信号处理教程,原创高性能示波器代码全开源地址:链接
第30章 复数FFT的实现


    本章主要讲解复数FFT的浮点和定点Q31,Q15的实现。
    本章节使用的复数FFT函数来自ARM官方库的TransformFunctions部分
    30.1 复数FFT
    30.2 复数FFT-基2算法
    30.3 复数FFT-基4算法
    30.4 总结

30.1 复数FFT

30.1.1 描述

    当前复数FFT函数支持三种数据类型,分别是浮点,Q31和Q15。这些FFT函数有一个共同的特点,就是用于输入信号的缓冲,在转化结束后用来存储输出结果。这样做的好处是节省了RAM空间,不需要为输入和输出结果分别设置缓存。由于是复数FFT,所以输入和输出缓存要存储实部和虚部。存储顺序如下:{real[0], imag[0], real[1], imag[1],………………} ,在使用中切记不要搞错。

30.1.2 浮点

    浮点复数FFT使用了一个混合基数算法,通过多个基8与单个基2或基4算法实现。根据需要,该算法支持的长度[16,32,64,...,4096]和每个长度使用不同的旋转因子表。
    浮点复数FFT使用了标准的FFT定义,FFT正变换的输出结果会被放大fftLen倍数,计算FFT逆变换的时候会缩小到1/fftLen。这样就与教科书中的定义一致了。
    定义好的旋转因子和位反转表已经在头文件arm_const_structs.h中定义好了,调用浮点FFT函数arm_cfft_f32时,包含相应的头文件即可。比如:
       arm_cfft_f32(arm_cfft_sR_f32_len64, pSrc, 1, 1)
    上式就是计算一个64点的FFT逆变换包括位反转。数据结构arm_cfft_sR_f32_len64可以认为是常数,计算的过程中是不能修改的。同样是这种数据结构还能用于混合基的FFT正变换和逆变换。
    早期发布的浮点复数FFT函数版本包含基2和基4两种方法实现的,但是不推荐大家再使用了。现在全部用arm_cfft_f32代替了。

30.1.3 定点Q31和Q15

    定点库提供了基2和基4两种算法,基2算法支持的数据长度[16, 32, 64, ..., 4096],基4算法支持的数据长度[16, 32, 64, ..., 4096]。一般情况下,建议使用基4算法,基4算法比基2算法执行速度要快一些。
    为了防止计算结果溢出,定点FFT每个蝶形运算的结果都要做放缩处理。对于基2算法,每次蝶形运算的结果要做0.5倍的放缩。对于基4算法,每次蝶形运算的结果要做0.25倍的放缩。FFT逆变换也要做相同的处理。相对于标准教科书式的FFT定义,定点FFT的计算结果放缩了1/fftLen(数据)倍。定点FFT的逆变也要做放缩处理,但是跟教科书式的FFT定义是相符的。
    每个FFT变换都需要一个单独的结构体,但结构体中的旋转因子和位反转表可以被重新使用。
    每个数据类型都有一个相关的初始化函数,初始化函数主要完成如下操作:
l 初始化结构体成员。
l 初始化旋转因子和位反转表指针。
    使用初始化函数是可选的。尽管如此,如果使用了初始化函数,那么结构体不能放在const data段,如果要放在const data段,应当按照如下方法进行初始化:
*arm_cfft_radix2_instance_q31 S = {fftLen, ifftFlag, bitReverseFlag, pTwiddle, pBitRevTable, twidCoefModifier, bitRevFactor};   
*arm_cfft_radix2_instance_q15 S = {fftLen, ifftFlag, bitReverseFlag, pTwiddle, pBitRevTable, twidCoefModifier, bitRevFactor};   
*arm_cfft_radix4_instance_q31 S = {fftLen, ifftFlag, bitReverseFlag, pTwiddle, pBitRevTable, twidCoefModifier, bitRevFactor};   
*arm_cfft_radix4_instance_q15 S = {fftLen, ifftFlag, bitReverseFlag, pTwiddle, pBitRevTable, twidCoefModifier, bitRevFactor};   
*arm_cfft_instance_f32 S = {fftLen, pTwiddle, pBitRevTable, bitRevLength};
    Q15和Q31 FFT使用了一个比较大的旋转因数和位反转表。这个表定义了一个最大长度的变换,表的子集可以用于短的变化。

30.1.4 arm_cfft_f32

函数定义如下:
    void arm_cfft_f32(
         const arm_cfft_instance_f32 * S,
         float32_t * p1,
         uint8_t ifftFlag,
         uint8_t bitReverseFlag)
参数定义:
     [in]   *S    points to an instance of the floating-point CFFT structure.  
     [in, out] *p1   points to the complex data buffer of size <code>2*fftLen</code>. Processing
  occurs in-place.  
     [in] ifftFlag         flag that selects forward (ifftFlag=0) or inverse (ifftFlag=1) transform.  
     [in] bitReverseFlag flag that enables (bitReverseFlag=1) or disables (bitReverseFlag=0) bit
reversal of output.  
注意事项:
结构const arm_cfft_instance_f32的定义如下(在文件arm_math.h文件):
      typedef struct
      {
          uint16_t fftLen;                  
          const float32_t *pTwiddle;         
          const uint16_t *pBitRevTable;      
          uint16_t bitRevLength;            
      } arm_cfft_instance_f32;

    下面通过在开发板上运行这个函数并计算幅频相应,然后再与Matlab计算的结果做对比。
  1. #define TEST_LENGTH_SAMPLES 2048
  2. /* 输入和输出缓冲 */
  3. static float32_t testOutput[TEST_LENGTH_SAMPLES/2];
  4. static float32_t testInput_f32_10khz[TEST_LENGTH_SAMPLES];
  5. /* 变量 */
  6. uint32_t fftSize = 1024;
  7. uint32_t ifftFlag = 0;
  8. uint32_t doBitReverse = 1;
  9. /*
  10. *********************************************************************************************************
  11. *    函 数 名: arm_cfft_f32_app
  12. *    功能说明: 调用函数arm_cfft_f32_app计算幅频。
  13. *    形    参:无
  14. *    返 回 值: 无
  15. *********************************************************************************************************
  16. */
  17. void arm_cfft_f32_app(void)
  18. {
  19. uint16_t i;
  20. /* 按照实部,虚部,实部,虚部..... 的顺序存储数据 */
  21. for(i=0; i<1024; i++)
  22. {
  23. /* 虚部全部置零 */
  24. testInput_f32_10khz[i*2+1] = 0;
  25. /* 50Hz正弦波,采样率1KHz ,作为实部 */
  26. testInput_f32_10khz[i*2] = arm_sin_f32(2*3.1415926f*50*i/1000);
  27. }
  28. /* CFFT变换 */
  29. arm_cfft_f32(&arm_cfft_sR_f32_len1024, testInput_f32_10khz, ifftFlag, doBitReverse);
  30. /* 求解模值  */
  31. arm_cmplx_mag_f32(testInput_f32_10khz, testOutput, fftSize);
  32. /* 串口打印求解的模值 */
  33. for(i=0; i<1024; i++)
  34. {
  35. printf("%f\r\n", testOutput[i]);
  36. }
  37. }
复制代码

运行如上函数可以通过串口打印出计算的模值,下面我们就通过Matlab计算的模值跟arm_cfft_f32计算的模值做对比。
    对比前需要先将串口打印出的数据加载到Matlab中,并给这个数组起名sampledata,加载方法在前面的教程中已经讲解过,这里不做赘述了。Matlab中运行的代码如下:
Fs = 1000;                % 采样率
N  = 1024;               % 采样点数
n  = 0:N-1;               % 采样序列
t  = 0:1/Fs:1-1/Fs;        % 时间序列
f = n * Fs / N;             %真实的频率

x = sin(2*pi*50*t) ;  %原始信号
y = fft(x, N);        %对原始信号做FFT变换

subplot(2,1,2);
plot(f, abs(y));       %绘制幅频相应曲线
title('Matlab计算结果');
xlabel('频率');
ylabel('幅度');

subplot(2,1,1);
plot(f,  sampledata);       %绘制幅频相应曲线
title('复数FFT计算结果');
xlabel('频率');
ylabel('幅度');

Matlab运行结果如下:

从上面的对比结果中可以看出,Matlab和函数arm_cfft_f32计算的结果基本是一直的。
30.1.png
努力打造安富莱高质量微信公众号:点击扫描图片关注
回复

使用道具 举报

740

主题

1326

回帖

3546

积分

管理员

春暖花开

Rank: 9Rank: 9Rank: 9

积分
3546
QQ
 楼主| 发表于 2015-4-15 10:37:40 | 显示全部楼层
30.2 复数FFT基2算法

30.2.1 arm_cfft_radix2_f32

    此函数已经不推荐使用,后面的版本会被删除,故不做介绍。

30.2.2 arm_cfft_radix2_q31

函数定义如下:
    void arm_cfft_radix2_q31(
    const arm_cfft_radix2_instance_q31 * S,
    q31_t * pSrc)
参数定义:
    [in]      *S    points to an instance of the fixed-point CFFT/CIFFT structure.  
    [in, out] *pSrc points to the complex data buffer of size <code>2*fftLen</code>. Processing occurs in-place.
注意事项:
结构const arm_cfft_radix2_instance_q31的定义如下(在文件arm_math.h文件):
      typedef struct
      {
          uint16_t fftLen;               
          uint8_t ifftFlag;               
          uint8_t bitReverseFlag;         
          q31_t *pTwiddle;                 
          uint16_t *pBitRevTable;         
          uint16_t twidCoefModifier;      
          uint16_t bitRevFactor;         
      } arm_cfft_radix2_instance_q31;

    下面通过在开发板上运行这个函数并计算幅频相应,然后再与Matlab计算的结果做对比。
  1. q31_t testInput_radix2_q31_50hz[TEST_LENGTH_SAMPLES];
  2. q31_t testOutputQ31[TEST_LENGTH_SAMPLES/2];
  3. /*
  4. *********************************************************************************************************
  5. *    函 数 名: arm_cfft_radix2_q31_app
  6. *    功能说明: 调用函数arm_cfft_radix2_q31计算幅频。
  7. *    形    参:无
  8. *    返 回 值: 无
  9. *********************************************************************************************************
  10. */
  11. void arm_cfft_radix2_q31_app(void)
  12. {
  13. uint16_t i,j;
  14. arm_cfft_radix2_instance_q31 S;
  15. fftSize = 1024;
  16.     ifftFlag = 0;
  17.     doBitReverse = 1;
  18. /* 初始化结构S */
  19. arm_cfft_radix2_init_q31(&S, fftSize, ifftFlag, doBitReverse);
  20. /* 按照实部,虚部,实部,虚部..... 的顺序存储数据 */
  21. for(i=0; i<1024; i++)
  22. {
  23. testInput_radix2_q31_50hz[i*2+1] = 0;
  24. /* 51.2Hz正弦波,采样率1024Hz。
  25.    arm_sin_q31输入参数的范围0-2^31, 这里每20次为一个完整的正弦波,
  26.    2^31 / 20 = 107374182.4
  27. */
  28. j = i % 20;
  29.       testInput_radix2_q31_50hz[i*2] = arm_sin_q31(107374182*j);
  30.       printf("%drn", testInput_radix2_q31_50hz[i*2]);
  31. }
  32. /* 输出结果分割线 */
  33. printf("**************************************************rn");
  34. printf("**************************************************rn");
  35. /* 计算CFFT */
  36. arm_cfft_radix2_q31(&S, testInput_radix2_q31_50hz);
  37. /* 计算模值 */
  38. arm_cmplx_mag_q31(testInput_radix2_q31_50hz, testOutputQ31, fftSize);
  39. /* 串口打印求解的模值 */
  40. for(i=0; i<1024; i++)
  41. {
  42. printf("%drn", testOutputQ31[i]);
  43. }
  44. }
复制代码
运行如上函数可以通过串口打印出原始信号和计算的模值,下面我们就通过Matlab计算的模值跟arm_cfft_radix2_q31计算的模值做对比。
    对比前需要先将串口打印出的数据加载到Matlab中,原始信号起名signal,函数arm_cmplx_mag_q31计算的模值起名sampledata,加载方法在前面的教程中已经讲解过,这里不做赘述了。Matlab中运行的代码如下:
Fs = 1024;              % 采样率
N  = 1024;             % 采样点数
n  = 0:N-1;             % 采样序列
f = n * Fs / N;           %真实的频率

y = fft(signal, N);        %对原始信号做FFT变换

subplot(2,1,2);
plot(f, abs(y));           %绘制幅频相应曲线
title('Matlab计算结果');
xlabel('频率');
ylabel('幅度');

subplot(2,1,1);
plot(f,  sampledata);      %绘制幅频相应曲线
title('复数FFT计算结果');
xlabel('频率');
ylabel('幅度');

Matlab运行结果如下:
30.2.png

从上面的对比结果中可以看出,Matlab和函数arm_cfft_radix2_q31计算的频率点基本是一致的,而幅值大小不一样是因为调用函数arm_cmplx_mag_q31对数据结果做了移位处理。

30.2.3 arm_cfft_radix2_q15

函数定义:
    void arm_cfft_radix2_q15(
        const arm_cfft_radix2_instance_q15 * S,
        q15_t * pSrc)
参数定义:
    [in]      *S    points to an instance of the fixed-point CFFT/CIFFT structure.  
    [in, out] *pSrc points to the complex data buffer of size <code>2*fftLen</code>. Processing occurs in-place.  
注意事项:
结构const arm_cfft_radix2_instance_q15的定义如下(在文件arm_math.h文件):
      typedef struct
      {
          uint16_t fftLen;               
          uint8_t ifftFlag;               
          uint8_t bitReverseFlag;         
          q15_t *pTwiddle;                    
      uint16_t *pBitRevTable;         
      uint16_t twidCoefModifier;      
          uint16_t bitRevFactor;         
      } arm_cfft_radix2_instance_q15;

    下面通过在开发板上运行这个函数并计算幅频相应,然后再与Matlab计算的结果做对比。
  1. q15_t testInput_radix2_q15_50hz[TEST_LENGTH_SAMPLES];
  2. q15_t testOutputQ15[TEST_LENGTH_SAMPLES/2];
  3. /*
  4. *********************************************************************************************************
  5. *    函 数 名: arm_cfft_radix2_q15_app
  6. *    功能说明: 调用函数arm_cfft_radix2_q15计算幅频。
  7. *    形    参:无
  8. *    返 回 值: 无
  9. *********************************************************************************************************
  10. */
  11. void arm_cfft_radix2_q15_app(void)
  12. {
  13. uint16_t i,j;
  14. arm_cfft_radix2_instance_q15 S;
  15. fftSize = 1024;
  16.     ifftFlag = 0;
  17.     doBitReverse = 1;
  18. /* 初始化结构S */
  19. arm_cfft_radix2_init_q15(&S, fftSize, ifftFlag, doBitReverse);
  20. /* 按照实部,虚部,实部,虚部..... 的顺序存储数据 */
  21. for(i=0; i<1024; i++)
  22. {
  23. /* 虚部全部置0 */
  24. testInput_radix2_q15_50hz[i*2+1] = 0;
  25. /* 51.2Hz正弦波,采样率1024Hz。
  26.    arm_sin_q15输入参数的范围[0, 32768), 这里每20次为一个完整的正弦波,
  27.    32768 / 20 = 1638.4
  28. */
  29. j = i % 20;
  30.       testInput_radix2_q15_50hz[i*2] = arm_sin_q15(1638*j);
  31.       printf("%drn", testInput_radix2_q15_50hz[i*2]);
  32. }
  33. /* 输出结果分割线 */
  34. printf("**************************************************rn");
  35. printf("**************************************************rn");
  36. /* 计算CFFT */
  37. arm_cfft_radix2_q15(&S, testInput_radix2_q15_50hz);
  38. /* 计算模值 */
  39. arm_cmplx_mag_q15(testInput_radix2_q15_50hz, testOutputQ15, fftSize);
  40. /* 串口打印求解的模值 */
  41. for(i=0; i<1024; i++)
  42. {
  43. printf("%drn", testOutputQ15[i]);
  44. }
  45. }
复制代码
运行如上函数可以通过串口打印出原始信号和计算的模值,下面我们就通过Matlab计算的模值跟arm_cfft_radix2_q15计算的模值做对比。
    对比前需要先将串口打印出的数据加载到Matlab中,原始信号起名signal,函数arm_cmplx_mag_q15计算的模值起名sampledata,加载方法在前面的教程中已经讲解过,这里不做赘述了。Matlab中运行的代码如下:
Fs = 1024;                     % 采样率
N  = 1024;                     % 采样点数
n  = 0:N-1;                     % 采样序列
f = n * Fs / N;                %真实的频率

y = fft(signal, N);           %对原始信号做FFT变换

subplot(2,1,2);
plot(f, abs(y));           %绘制幅频相应曲线
title('Matlab计算结果');
xlabel('频率');
ylabel('幅度');

subplot(2,1,1);
plot(f,  sampledata);      %绘制幅频相应曲线
title('复数FFT计算结果');
xlabel('频率');
ylabel('幅度');

Matlab运行结果如下:
30.3.png

从上面的对比结果中可以看出,Matlab和函数arm_cfft_radix2_q15计算的频率点基本是一致的,而幅值大小不一样是因为调用函数arm_cmplx_mag_q15对数据结果做了移位处理。
努力打造安富莱高质量微信公众号:点击扫描图片关注
回复

使用道具 举报

740

主题

1326

回帖

3546

积分

管理员

春暖花开

Rank: 9Rank: 9Rank: 9

积分
3546
QQ
 楼主| 发表于 2015-4-15 10:48:12 | 显示全部楼层
30.3 复数FFT基4算法
    如果希望直接调用FFT程序计算IFFT,可以用下面的方法:

30.3.1 arm_cfft_radix4_f32

    此函数已经不推荐使用,后面的版本会被删除,故不做介绍。

30.3.2 arm_cfft_radix4_q31

函数定义如下:
    void arm_cfft_radix4_q31(
    const arm_cfft_radix4_instance_q31 * S,
    q31_t * pSrc)
参数定义:
    [in]      *S    points to an instance of the fixed-point CFFT/CIFFT structure.  
    [in, out] *pSrc points to the complex data buffer of size <code>2*fftLen</code>. Processing occurs in-place.
注意事项:
1. 结构const arm_cfft_radix4_instance_q31的定义如下(在文件arm_math.h文件):
      typedef struct
      {
          uint16_t fftLen;               
          uint8_t ifftFlag;               
          uint8_t bitReverseFlag;         
          q31_t *pTwiddle;                 
          uint16_t *pBitRevTable;         
          uint16_t twidCoefModifier;      
          uint16_t bitRevFactor;         
      } arm_cfft_radix4_instance_q31;
2. 为了防止数据饱和,每次蝶形运行的结果都要除以2,故不同的长度的FFT运算的最终结果输出格式不同。具体信息如下:
    Q31 CFFT
30.4.png
    Q31 CIFFT
   30.5.png
    下面通过在开发板上运行这个函数并计算幅频相应,然后再与Matlab计算的结果做对比。
  1. q31_t testInput_radix4_q31_50hz[TEST_LENGTH_SAMPLES];
  2. /*
  3. *********************************************************************************************************
  4. *    函 数 名: arm_cfft_radix4_q31_app
  5. *    功能说明: 调用函数arm_cfft_radix4_q31_app计算幅频。
  6. *    形    参:无
  7. *    返 回 值: 无
  8. *********************************************************************************************************
  9. */
  10. void arm_cfft_radix4_q31_app(void)
  11. {
  12. uint16_t i,j;
  13. arm_cfft_radix4_instance_q31 S;
  14. fftSize = 1024;
  15.     ifftFlag = 0;
  16.     doBitReverse = 1;
  17. /* 初始化结构S */
  18. arm_cfft_radix4_init_q31(&S, fftSize, ifftFlag, doBitReverse);
  19. /* 按照实部,虚部,实部,虚部..... 的顺序存储数据 */
  20. for(i=0; i<1024; i++)
  21. {
  22. testInput_radix4_q31_50hz[i*2+1] = 0;
  23. /* 51.2Hz正弦波,采样率1024Hz。
  24.    arm_sin_q31输入参数的范围0-2^31, 这里每20次为一个完整的正弦波,
  25.    2^31 / 20 = 107374182.4
  26. */
  27. j = i % 20;
  28.       testInput_radix4_q31_50hz[i*2] = arm_sin_q31(107374182*j);
  29.       printf("%drn", testInput_radix4_q31_50hz[i*2]);
  30. }
  31. /* 输出结果分割线 */
  32. printf("**************************************************rn");
  33. printf("**************************************************rn");
  34. /* 计算CFFT */
  35. arm_cfft_radix4_q31(&S, testInput_radix4_q31_50hz);
  36. /* 计算模值 */
  37. arm_cmplx_mag_q31(testInput_radix4_q31_50hz, testOutputQ31, fftSize);
  38. /* 串口打印求解的模值 */
  39. for(i=0; i<1024; i++)
  40. {
  41. printf("%drn", testOutputQ31[i]);
  42. }
  43. }
复制代码
运行如上函数可以通过串口打印出原始信号和计算的模值,下面我们就通过Matlab计算的模值跟arm_cfft_radix4_q31计算的模值做对比。
    对比前需要先将串口打印出的数据加载到Matlab中,原始信号起名signal,函数arm_cmplx_mag_q31计算的模值起名sampledata,加载
方法在前面的教程中已经讲解过,这里不做赘述了。Matlab中运行的代码如下:
Fs = 1024;              % 采样率
N  = 1024;             % 采样点数
n  = 0:N-1;             % 采样序列
f = n * Fs / N;           %真实的频率

y = fft(signal, N);        %对原始信号做FFT变换

subplot(2,1,2);
plot(f, abs(y));           %绘制幅频相应曲线
title('Matlab计算结果');
xlabel('频率');
ylabel('幅度');

subplot(2,1,1);
plot(f,  sampledata);      %绘制幅频相应曲线
title('复数FFT计算结果');
xlabel('频率');
ylabel('幅度');

Matlab运行结果如下:
30.6.png

从上面的对比结果中可以看出,Matlb和函数arm_cfft_radix4_q31计算的频率点基本是一致的,而幅值大小不一样是因为调用
函数arm_cmplx_mag_q31和arm_cfft_radix4_q31对数据结果做了移位处理。

30.3.3 arm_cfft_radix4_q15


函数定义:
    void arm_cfft_radix4_q15(
        const arm_cfft_radix4_instance_q15 * S,
        q15_t * pSrc)
参数定义:
    [in]      *S    points to an instance of the fixed-point CFFT/CIFFT structure.  
    [in, out] *pSrc points to the complex data buffer of size <code>2*fftLen</code>. Processing occurs in-place.  
注意事项:
1. 结构const arm_cfft_radix4_instance_q15的定义如下(在文件arm_math.h文件):
      typedef struct
      {
          uint16_t fftLen;               
          uint8_t ifftFlag;               
          uint8_t bitReverseFlag;         
          q15_t *pTwiddle;                    
     uint16_t *pBitRevTable;         
     uint16_t twidCoefModifier;      
          uint16_t bitRevFactor;         
      } arm_cfft_radix4_instance_q15;
2. 为了防止数据饱和,每次蝶形运行的结果都要除以2,故不同的长度的FFT运算的最终结果输出格式不同。具体信息如下:
    Q15 CFFT
       30.7.png

    Q15 CIFFT
       30.8.png

    下面通过在开发板上运行这个函数并计算幅频相应,然后再与Matlab计算的结果做对比。
  1. q15_t testInput_radix4_q15_50hz[TEST_LENGTH_SAMPLES];
  2. /*
  3. *********************************************************************************************************
  4. *    函 数 名: arm_cfft_radix4_q15_app
  5. *    功能说明: 调用函数arm_cfft_radix4_q15计算幅频。
  6. *    形    参:无
  7. *    返 回 值: 无
  8. *********************************************************************************************************
  9. */
  10. void arm_cfft_radix4_q15_app(void)
  11. {
  12. uint16_t i,j;
  13. arm_cfft_radix4_instance_q15 S;
  14. fftSize = 1024;
  15.     ifftFlag = 0;
  16.     doBitReverse = 1;
  17. /* 初始化结构S */
  18. arm_cfft_radix4_init_q15(&S, fftSize, ifftFlag, doBitReverse);
  19. /* 按照实部,虚部,实部,虚部..... 的顺序存储数据 */
  20. for(i=0; i<1024; i++)
  21. {
  22. /* 虚部全部置0 */
  23. testInput_radix4_q15_50hz[i*2+1] = 0;
  24. /* 51.2Hz正弦波,采样率1024Hz。
  25.    arm_sin_q15输入参数的范围[0, 32768), 这里每20次为一个完整的正弦波,
  26.    32768 / 20 = 1638.4
  27. */
  28. j = i % 20;
  29.       testInput_radix4_q15_50hz[i*2] = arm_sin_q15(1638*j);
  30.       printf("%drn", testInput_radix4_q15_50hz[i*2]);
  31. }
  32. /* 输出结果分割线 */
  33. printf("**************************************************rn");
  34. printf("**************************************************rn");
  35. /* 计算CFFT */
  36. arm_cfft_radix4_q15(&S, testInput_radix4_q15_50hz);
  37. /* 计算模值 */
  38. arm_cmplx_mag_q15(testInput_radix4_q15_50hz, testOutputQ15, fftSize);
  39. /* 串口打印求解的模值 */
  40. for(i=0; i<1024; i++)
  41. {
  42. printf("%drn", testOutputQ15[i]);
  43. }
  44. }
复制代码

运行如上函数可以通过串口打印出原始信号和计算的模值,下面我们就通过Matlab计算的模值跟arm_cfft_radix4_q15计算的模值做对比。
    对比前需要先将串口打印出的数据加载到Matlab中,原始信号起名signal,函数arm_cmplx_mag_q15计算的模值起名sampledata,加载
方法在前面的教程中已经讲解过,这里不做赘述了。Matlab中运行的代码如下:
Fs = 1024;                     % 采样率
N  = 1024;                     % 采样点数
n  = 0:N-1;                     % 采样序列
f = n * Fs / N;                %真实的频率

y = fft(signal, N);           %对原始信号做FFT变换

subplot(2,1,2);
plot(f, abs(y));           %绘制幅频相应曲线
title('Matlab计算结果');
xlabel('频率');
ylabel('幅度');

subplot(2,1,1);
plot(f,  sampledata);      %绘制幅频相应曲线
title('复数FFT计算结果');
xlabel('频率');
ylabel('幅度');

Matlab运行结果如下:
30.9.png

从上面的对比结果中可以看出,Matlab和函数arm_cfft_radix4_q15计算的频率点基本是一致的,而幅值大小不一样是因为调用函数arm_cmplx_mag_q15
arm_cfft_radix4_q15对数据结果做了移位处理。
努力打造安富莱高质量微信公众号:点击扫描图片关注
回复

使用道具 举报

740

主题

1326

回帖

3546

积分

管理员

春暖花开

Rank: 9Rank: 9Rank: 9

积分
3546
QQ
 楼主| 发表于 2015-4-15 10:54:25 | 显示全部楼层
30.4  总结


    本章节内容比较多,有兴趣的可以深入了解源码的实现。
努力打造安富莱高质量微信公众号:点击扫描图片关注
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

QQ|小黑屋|Archiver|手机版|硬汉嵌入式论坛

GMT+8, 2024-5-4 21:33 , Processed in 0.205524 second(s), 28 queries .

Powered by Discuz! X3.4 Licensed

Copyright © 2001-2023, Tencent Cloud.

快速回复 返回顶部 返回列表