基于STM32的FFT频谱分析+波形识别

65 阅读 0 评论 43 点赞

我是靠谱客的博主酷炫大门，最近开发中收集的这篇文章主要介绍基于STM32的FFT频谱分析+波形识别，觉得挺不错的，现在分享给大家，希望可以做个参考。

概述

一.硬件部分

信号发生器，正点原子精英板，3.5’TFTLCD，两根杜邦线（接PC1和GND）

二.基本思路

1.使用ADC采集音频信号

2.使用官方提供的FFT函数（1024点）对采集到的信号进行处理

3.量化、频谱图显示

采样频率：Fs = 2400Hz（触摸版本可以根据实际情况调整）

样本数量：NPT = 1024

三.程序编写

1.ADC采样
这次用到的是ADC的DMA传输，这可以减少对程序的占用。这里着重关注这几行代码：

 ADC_InitStructure.ADC_Mode = ADC_Mode_Independent; //ADC1 工作在独立模式
 ADC_InitStructure.ADC_ScanConvMode =DISABLE; //模数转换工作在非扫描模式
 ADC_InitStructure.ADC_ContinuousConvMode =DISABLE; //模数转换工作在不连续转换模式
 ADC_InitStructure.ADC_ExternalTrigConv = ADC_ExternalTrigConv_T1_CC1; //Timer1触发转换开启(定时器T1的CC1通道，控制采样频率)

假设有ADC1有四个通道需要转换（本例中只有一个通道）：

	连续转换	单次转换
扫描模式	开始转换CH0、CH1、CH2、CH3转换完成后，继续转换CH0	开始转换CH0、CH1、CH2、CH3转换完成后停止，等待ADC的下一次启动
非扫描模式	开始转换CH0，转换完成后，继续转换CH0	转换完成后停止，等待ADC的下一次启动

我的代码就是按这个写的，不这样写应该也可以，不过我认为还是要多注意细节，不然根本不知道问题出在哪儿。

 DMA_InitStructure.DMA_PeripheralBaseAddr = (u32)&ADC1->DR; //DMA 外设 ADC 基地址
 DMA_InitStructure.DMA_MemoryBaseAddr = (u32)&ADC_Value; //DMA 内存基地址

DMA传输注意这两句，存放地址。

2.FFT频谱分析
先下载一个FFT官方库，添加这几个文件
在这里插入图片描述
主要就是这几个函数

for(i=0;i<NPT;i++)
  {
    lBufInArray[i]=ADC_Value[i]<<16;
  }
  cr4_fft_1024_stm32(lBufOutArray, lBufInArray, NPT);    
  GetPowerMag();
/******************************************************************
函数名称:GetPowerMag()
函数功能:计算各次谐波幅值  [short 的范围，是-32767 到 32767 。也就是 -(2^15 - 1)到(2^15 - 1)。]
参数说明:
备　　注:先将lBufOutArray分解成实部(X)和虚部(Y)，然后计算幅值(sqrt(X*X+Y*Y)
*******************************************************************/
void GetPowerMag(void)
{
    signed short lX,lY;                                                  //算频率的话Fn=i*Fs/NPT  //由于此处i是从0开始的，所以不需要再减1
    float X,Y,Mag;                                                      
    unsigned short i;
    for(i=0; i<NPT/2; i++)                                                  //经过FFT后，每个频率点处的真实幅值  A0=lBufOutArray[0]/NPT
    {                                                                       //                                 Ai=lBufOutArray[i]*2/NPT
        lX  = (lBufOutArray[i] << 16) >> 16;  //lX  = lBufOutArray[i];
        lY  = (lBufOutArray[i] >> 16);
                                    
        X = NPT * ((float)lX) / 32768;//除以32768再乘65536是为了符合浮点数计算规律,不管他
        Y = NPT * ((float)lY) / 32768;
        Mag = sqrt(X * X + Y * Y) / NPT;
        if(i == 0)
            lBufMagArray[i] = (unsigned long)(Mag * 32768);   //0Hz是直流分量，直流分量不需要乘以2
        else
            lBufMagArray[i] = (unsigned long)(Mag * 65536);
    }
}

需要说明的是：按照FFT官方库的说明，lBufOutArray和lBufInArray都必须是32位的数据类型，其中高16位存储实部，低16位存储虚部。因为信号发生器输出的不可能是虚数，所以对于lBufInArray来说，低16位存储的虚部总是为0。上面的代码，主要就是计算各次谐波的幅值，就是把虚部和实部取出来平方开根号。

void lcd_show_fft(unsigned int *p)
{
   unsigned int *pp = p+1;             //p+1相当于我直接把0HZ部分滤掉了
   unsigned int i = 0;
   for(i = 0;i<480;i++)
   {
      LCD_Fill(0,        i, *pp*0.11, (i+1), WHITE);     //有效部分白色       
      LCD_Fill(*pp*0.11, i, 270,       (i+1), BLACK);   //其他就黑色
      pp++;
   }
 }

显示的函数也没什么可讲的，注意不要用画直线的函数，刷屏很慢。

3.波形识别
先说说我的思路，这个思路肯定不是最好的，我们先来看看各个波形的特点：
在这里插入图片描述

在这里插入图片描述

正弦波：只有基波分量，基本无谐波分量，没什么好说的。
方波：除了基波，还有3，5，7次谐波分量，且3次谐波分量为基波分量的1/3.
三角波：除了基波，还有3，5，7次谐波分量，但3次谐波分量为基波分量的1/9.
锯齿波：除了基波，还有2，3，4次谐波分量.

/***********************************************
找最大值，次大值……对应的频率，分析波形
*************************************************/
void select_max(float *f,float *a)
{
   int i,j;
   float k,k1,m;
    float aMax =0.0,aSecondMax = 0.0,aThirdMax = 0.0,aFourthMax=0.0;
    float fMax =0.0,fSecondMax = 0.0,fThirdMax = 0.0,fFourthMax=0.0;
   int nMax=0,nSecondMax=0,nThirdMax=0,nFourthMax=0;
   for ( i = 1; i < NPT/2; i++)//i必须是1，是0的话，会把直流分量加进去！！！！
    {
        if (a[i]>aMax)
        {
            aMax = a[i]; 
       nMax=i;
       fMax=f[nMax];
        }
    }
  for ( i=1; i < NPT/2; i++)
    {
    if (nMax == i)
    {
      continue;//跳过原来最大值的下标，直接开始i+1的循环
    }
        if (a[i]>aSecondMax&&a[i]>a[i+1]&&a[i]>a[i-1])
        {
            aSecondMax = a[i]; 
       nSecondMax=i;
       fSecondMax=f[nSecondMax];
        }
    }
  for ( i=1; i < NPT/2; i++)
    {
    if (nMax == i||nSecondMax==i)
    {
      continue;//跳过原来最大值的下标，直接开始i+1的循环
    }
        if (a[i]>aThirdMax&&a[i]>a[i+1]&&a[i]>a[i-1])
        {
            aThirdMax = a[i]; 
       nThirdMax=i;
       fThirdMax=f[nThirdMax];
        }
    }
  for ( i=1; i < NPT/2; i++)
    {
    if (nMax == i||nSecondMax==i||nThirdMax==i)
    {
      continue;//跳过原来最大值的下标，直接开始i+1的循环
    }
        if (a[i]>aFourthMax&&a[i]>a[i+1]&&a[i]>a[i-1])
        {
            aFourthMax = a[i]; 
       nFourthMax=i;
       fFourthMax=f[nFourthMax];
        }
    }
  k=fabs(2*fMax-fSecondMax);
  k1=fabs(3*fMax-fSecondMax);
  m=fabs((float)(aMax-3.0*aSecondMax)); 
  if(k<=5)
     LCD_ShowString(275,230,12*4,12,12,"JvChi  ");
  else if(k1<=5&&m<0.4) 
     LCD_ShowString(275,230,12*4,12,12,"Fang   ");
  else if(k1<=5&&m>=0.4)
     LCD_ShowString(275,230,12*4,12,12,"SanJiao");
  else LCD_ShowString(275,230,12*4,12,12,"Sin    ");
}

这里我们要做的，就是把各次谐波的频率和对应的的幅值提取出来，我采用了两种方法，并把它们结合起来使用。首先是求幅值的最大值、次大值和次次大值，网上也有人写过，我直接拿来用感觉并不好用，就自己写了一个。求出来之后还要和两侧的值比较，看他是否为极大值，是的话才保留。最后，按前面的思路进行比较，比较时的参数是我自己选的（适用于50Hz~200Hz，其他范围的可能也可以）。

4.触屏调整采样频率
我最初定的采样频率是2400Hz，因为要求的最大可以分析1kHz的谐波，采样频率要大于2倍的信号频率。通过计算，频率分辨率只有2.34Hz，这样导致有些频率和幅值测得不准确。补零也不能提高频率分辨率，就考虑实时的改变采样频率。

if(tp_dev.sta&TP_PRES_DOWN)   //触摸屏被按下
 { 
  if(tp_dev.x[0]>270&&tp_dev.x[0]<320&&tp_dev.y[0]>360&&tp_dev.y[0]<400)
  { 
    LCD_Fill(270,360,320,400,BLACK);
    delay_ms(200);
    LCD_Fill(270,360,320,400,YELLOW);
    ji_shu=ji_shu-200;
    TIM1_Int_Init(ji_shu-1,fen_pin-1);//只有放这里才能保证只运行一次，不然可能会卡死
    POINT_COLOR=BLACK; //画笔颜色
       BACK_COLOR=YELLOW;  //背景色 
       LCD_ShowString(284,365,16*2,16,24,"up");
    //i++;
  }
  if(tp_dev.x[0]>270&&tp_dev.x[0]<320&&tp_dev.y[0]>440&&tp_dev.y[0]<480)
  { 
    LCD_Fill(270,440,320,480,BLACK);
    delay_ms(200);
    LCD_Fill(270,440,320,480,YELLOW);
    ji_shu=ji_shu+200;
    TIM1_Int_Init(ji_shu-1,fen_pin-1);
    POINT_COLOR=BLACK; //画笔颜色
       BACK_COLOR=YELLOW;  //背景色 
       LCD_ShowString(272,442,16*4,16,24,"down");
    //i--;
  }
    //TIM1_Int_Init(ji_shu,fen_pin);//取消注释，触屏后会卡死
  Fs=2400000/ji_shu;