STM32
直播中

杨火亭

7年用户 989经验值
擅长:控制/MCU
私信 关注
[问答]

如何利用arm-DSP库进行FFT计算获得信号的频谱/幅值及相位?

如何利用arm-DSP库进行FFT计算获得信号的频谱/幅值及相位?

回帖(1)

程家文

2021-11-22 11:27:40
代码终于要来了,那么咱们就一步一步来。核心代码就三个文件——data.c data.h ffttest.c,其中data.c data.h为数据的准备部分,数据当然来自Matlab。ffttest.c里才是干货。(没看过上的自行翻阅)

一.data.c/data.h数据准备
data.h


#ifndef _DATA_H
#define _DATA_H
#include "arm_math.h"
#include  
#define SAMPLE_N    64
#define SAMPLE_FS   3200
#define FFT_DPI     (SAMPLE_FS / SAMPLE_N);


extern float32_t data2handle[SAMPLE_N * 2];
extern float32_t FFT_Output[SAMPLE_N];


extern float32_t xiebo[SAMPLE_N];
extern float32_t phaseA[SAMPLE_N];
extern float32_t phaseB[SAMPLE_N];
extern float32_t phaseC[SAMPLE_N];


#endif


data.c


#include "data.h"


float32_t data2handle[SAMPLE_N * 2] = {0};
float32_t FFT_Output[SAMPLE_N]      = {0};


float32_t xiebo[SAMPLE_N]  = {0,12.0830468476080,23.0625562988315,31.9767233645297,38.1266205974219,41.1592650408757,41.1008388582348,38.3357763183915,33.5355339059327,27.5483591934343,21.2671323190531,15.4955017031488,10.8326010981529,7.59361568157981,5.77782765517063,5.08838180836849,5,4.86346545835348,4.03002514886167,1.97578767574228,-1.59380577304000,-6.67628905966528,-12.9524361960277,-19.8182546598069,-26.4644660940673,-31.9918434767550,-35.5451365280388,-36.4452976726157,-34.2997862737710,-29.0738765919851,-21.1116530786703,-11.1028754443124,-1.40834381901946e-14,11.1028754443124,21.1116530786703,29.0738765919851,34.2997862737710,36.4452976726157,35.5451365280388,31.9918434767550,26.4644660940672,19.8182546598069,12.9524361960277,6.67628905966530,1.59380577303998,-1.97578767574228,-4.03002514886168,-4.86346545835348,-5,-5.08838180836849,-5.77782765517063,-7.59361568157980,-10.8326010981529,-15.4955017031488,-21.2671323190531,-27.5483591934343,-33.5355339059328,-38.3357763183915,-41.1008388582348,-41.1592650408757,-38.1266205974219,-31.9767233645297,-23.0625562988315,-12.0830468476080};


float32_t phaseA[SAMPLE_N] = {0,0.490085701647803,0.975451610080641,1.45142338627231,1.91341716182545,2.35698368412999,2.77785116509801,3.17196642081823,3.53553390593274,3.86505226681368,4.15734806151273,4.40960632174178,4.61939766255643,4.78470167866105,4.90392640201615,4.97592363336099,5,4.97592363336099,4.90392640201615,4.78470167866104,4.61939766255643,4.40960632174178,4.15734806151273,3.86505226681368,3.53553390593274,3.17196642081823,2.77785116509801,2.35698368412999,1.91341716182545,1.45142338627231,0.975451610080641,0.490085701647802,6.12323399573677e-16,-0.490085701647803,-0.975451610080642,-1.45142338627231,-1.91341716182545,-2.35698368412999,-2.77785116509801,-3.17196642081823,-3.53553390593274,-3.86505226681369,-4.15734806151273,-4.40960632174178,-4.61939766255643,-4.78470167866104,-4.90392640201615,-4.97592363336099,-5,-4.97592363336098,-4.90392640201615,-4.78470167866105,-4.61939766255643,-4.40960632174178,-4.15734806151273,-3.86505226681368,-3.53553390593274,-3.17196642081823,-2.77785116509801,-2.35698368412999,-1.91341716182545,-1.45142338627231,-0.975451610080639,-0.490085701647803};
float32_t phaseB[SAMPLE_N] = {-4.33012701892219,-4.55431912460588,-4.73465064747553,-4.86938489638667,-4.95722430686905,-4.99732293738183,-4.98929461619302,-4.93321666042440,-4.82962913144534,-4.67952963378663,-4.48436370766344,-4.24601090763289,-3.96676670145618,-3.64932036348918,-3.29672907550034,-2.91238848433901,-2.50000000000000,-2.06353514902197,-1.60719732651581,-1.13538131517186,-0.652630961100256,-0.163595414108879,0.327015646150717,0.814477366972946,1.29409522551260,1.76125023960617,2.21144345109501,2.64033925325184,3.04380714504361,3.41796151011436,3.75919903739489,4.06423342295808,4.33012701892219,4.55431912460588,4.73465064747553,4.86938489638667,4.95722430686905,4.99732293738183,4.98929461619302,4.93321666042440,4.82962913144534,4.67952963378663,4.48436370766344,4.24601090763289,3.96676670145617,3.64932036348918,3.29672907550034,2.91238848433901,2.50000000000000,2.06353514902197,1.60719732651581,1.13538131517187,0.652630961100256,0.163595414108880,-0.327015646150719,-0.814477366972946,-1.29409522551261,-1.76125023960617,-2.21144345109501,-2.64033925325184,-3.04380714504361,-3.41796151011436,-3.75919903739489,-4.06423342295808};
float32_t phaseC[SAMPLE_N] = {4.33012701892219,4.06423342295808,3.75919903739489,3.41796151011436,3.04380714504360,2.64033925325184,2.21144345109501,1.76125023960617,1.29409522551261,0.814477366972945,0.327015646150716,-0.163595414108881,-0.652630961100257,-1.13538131517187,-1.60719732651581,-2.06353514902197,-2.50000000000000,-2.91238848433901,-3.29672907550034,-3.64932036348918,-3.96676670145617,-4.24601090763289,-4.48436370766344,-4.67952963378663,-4.82962913144534,-4.93321666042440,-4.98929461619302,-4.99732293738183,-4.95722430686905,-4.86938489638667,-4.73465064747553,-4.55431912460588,-4.33012701892220,-4.06423342295808,-3.75919903739489,-3.41796151011436,-3.04380714504360,-2.64033925325184,-2.21144345109501,-1.76125023960617,-1.29409522551260,-0.814477366972945,-0.327015646150718,0.163595414108876,0.652630961100260,1.13538131517187,1.60719732651581,2.06353514902197,2.50000000000000,2.91238848433901,3.29672907550035,3.64932036348918,3.96676670145617,4.24601090763289,4.48436370766344,4.67952963378663,4.82962913144534,4.93321666042440,4.98929461619302,4.99732293738183,4.95722430686905,4.86938489638667,4.73465064747553,4.55431912460588};


无非就是准备了几个数组:
1.float32_t data2handle[SAMPLE_N * 2],这个用来装处理后的采样数据,所谓处理就是之气那说的虚部填上0,所以对于SAMPLE_N 采样点,这个数组的长度为SAMPLE_N *2;
2.float32_t FFT_Output[SAMPLE_N],这个就是FFT输出的频谱数组,咋看之前都讲过了;
3.float32_t xiebo[SAMPLE_N],就是那个50Hz、100Hz、200Hz叠加后的信号;
4.float32_t phaseA、B、C[SAMPLE_N],就是标准三相电的信号;


二.ffttest.c数据处理
ffttest.c


#include "data.h"
#include "arm_const_structs.h"
#include "app_serial.h"
#include "os.h"


void Create_data2handle(float32_t *p)
{
    for(int i = 0;i < SAMPLE_N; i++)
    {
        data2handle[2 * i]     = p;
        data2handle[2 * i + 1] = 0;
    }
}


void CalXiebo(void)
{
   
    Create_data2handle(xiebo);
    arm_cfft_f32(&arm_cfft_sR_f32_len64, data2handle, 0, 1);
    arm_cmplx_mag_f32(data2handle, FFT_Output, SAMPLE_N);


    App_SerPrintf("/********************************************************************************/n");
    for(int i = 0,j = 0; i < SAMPLE_N;i++,j++)
    {
        App_SerPrintf("%.8ft",FFT_Output);
        if(j >= 5)
        {
            App_SerPrintf("n");
            j = -1;
        }
    }
    App_SerPrintf("n/********************************************************************************/n");   
}


void CalPhase(float32_t *phase,char phase_name)
{
    int Index_50Hz = 50 / FFT_DPI;
    float32_t Real,Imaginary,hudu,jiaodu;
    Create_data2handle(phase);
    arm_cfft_f32(&arm_cfft_sR_f32_len64, data2handle, 0, 1);
    arm_cmplx_mag_f32(data2handle, FFT_Output, SAMPLE_N);


    Real      = data2handle[Index_50Hz * 2];
    Imaginary = data2handle[Index_50Hz * 2 + 1];


    hudu = atan2(Imaginary , Real);
    jiaodu = hudu * 180.0 / 3.1415926;
    App_SerPrintf("%c相相位:%f°n",phase_name,jiaodu);


}


void FFTtest(void)
{
    CalXiebo();
    CalPhase(phaseA,'A');
    CalPhase(phaseB,'B');
    CalPhase(phaseC,'C');
}


三.万事俱备,只差测试






怎么样。。。在串口助手里打印出来还算能看吧。。
将输出的FFT_Output复制到Matlab里画图




可以看到,在分辨率50Hz的情况下,0代表的是直流分量,1为50Hz分量,2为100Hz分量,4为200Hz分量。之前50、100、200设置的幅值分别是5、30、15,那么用上一篇中所说的计算实际幅值的办法,160/32 = 5;960/32=30;480/32=15,完全OK。
那么再来看相位,A相-90°,B相150°,C相30°,向量图如下图,凑活看吧。跟设置的相位一毛一样。




所以给定一个信号,用FFT就能分解出你想看的分量的全部信息,即频率、幅值、相位,可以说非常给力了。
举报

更多回帖

×
20
完善资料,
赚取积分