请选择 进入手机版 | 继续访问电脑版
查看: 183|回复: 9

[原创] 有趣|如何用低成本MCU实现音乐频谱显示

[复制链接]
  • TA的每日心情
    开心
    2022-11-30 13:07
  • 签到天数: 209 天

    [LV.7]常住居民III

    2733

    主题

    5382

    帖子

    0

    管理员

    Rank: 9Rank: 9Rank: 9

    积分
    25599
    最后登录
    2022-12-7
    发表于 2022-9-22 10:15:30 | 显示全部楼层 |阅读模式
    有趣|如何用低成本MCU实现音乐频谱显示

    ♦效果展示
    下面是实际频谱显示的测试效果。


    最近在做一个有趣的小项目,其中有一小部分的内容的是使用FFT做音乐频谱显示。于是就有了下面这个音乐频谱显示的低成本方案,话不多说看看低成本MCU如何实现FFT音乐频谱显示吧。

    ♦音频采集硬件电路
    音频采集的硬件电路比较简单,主要的器件就是麦克风和LM358运放。
    11.png
    图中电路R5可调电阻的作用是来调节运放的增益。R4的作用的是给运放一个VDD*R4/(R3+R4) 的直流偏置,这里加直流偏置是由于ADC只能采集正电压值,为了不丢失负电压的音频信号,给信号整体加了一个直流偏置。


    但是这个图还有一个小问题,运放的输出端加了一个电容C2,C2会把直流偏置给隔掉。在设计时,这个电容可以去掉。


    下图是按照上图搭建的音频采集电路的输出信号,图中波动信号是施加的外部音频,是我们需要做音乐频谱显示需要的信号。该信号有一个2.3v的直流偏置,在后续处理时需要减去这个偏置。
    12.png
    ♦CTimer+ADC+DMA 音频信号采集
    为了呼应标题,我们选择的MCU是LPC845,这是NXP的一款低成本的MCU。考虑到我们平常听的音乐频率大都低于5kHz,在软件设计时设置ADC采样频率为10kHz。不要问为什么,问就是采样定理。

    LPC845的ADC有8个触发源,我们使用CTiimer match3来触发ADC,将寄存器SEQA_CTRL的bit 14:12设置为0x5。CTimer match 3的输出频率为10kHz。

    为了确保我们采集数据的实时性,DMA建议配置成双buffer模式,以防止采样的数据被覆盖掉。
    ♦FFT音频信号处理
    在DMA搬运ADC采样值时,使用了双buffer来搬,ADC采样值需要减去一个2.3V的直流偏置。Samples[]数组用于FFT计算。
    1.    //Calculate the FFT input buffer
    2.    if(g_DmaTransferDoneFlag_A == true)
    3.    {
    4.            for (i=0; i<128; i++)
    5.            {
    6.                   Samples[i] =(int16_t)(((g_AdcConvResult_A[i] & 0xfff0) >> 4) - 2979);//substract the 2.3v offset in the Amplifier output
    7.             }
    8.           g_DmaTransferDoneFlag_A = false;

    9.     }
    10.    else if(g_DmaTransferDoneFlag_B == true)
    11.    {
    12.            for (i=0; i<128; i++)
    13.            {
    14.                   Samples[i] =(int16_t)(((g_AdcConvResult_B[i] & 0xfff0) >> 4) - 2979);//substract the 2.3v offset in the Amplifier output
    15.            }
    16.            g_DmaTransferDoneFlag_B = false;
    17.    }
    复制代码
    根据FFT算法的原理,在进行FFT计算之前,还需要将ADC的采样值Samples[]乘上一个窗函数,这里我们使用的汉宁窗函数,由于篇幅限制,具体原理可以去查看FFT算法相关的资料。
    1.   //If 'Window' isn't rectangular, apply window
    2.     if(Window == Triangular){
    3.         //Apply a triangular window to the data.
    4.         for(Cnt = 0; Cnt<Len; Cnt++){
    5.             if(Cnt<(Len/2)) Samples[Cnt] = ((int32_t)Samples[Cnt]*Cnt)>>L2Len;
    6.             else Samples[Cnt] = ((int32_t)Samples[Cnt]*((Len/2)-Cnt))>>L2Len;
    7.         }
    8.     }
    9.     else if(Window == Hann){
    10.         //Use the cosine window wavetable to apply a Hann windowing function to the samples
    11.         for(Cnt = 0; Cnt<Len; Cnt++){
    12.             Index = (Cnt*CWLen)>>L2Len;
    13.             Samples[Cnt] = ((int32_t)Samples[Cnt]*(int32_t)CosWindow[Index])>>(CWBD);
    14.         }
    15.     }
    复制代码
    前面说了这么多,FFT算法才是实现音乐频谱显示的关键部分(其实上边每一步都缺一不可)。

    我在网上找了好多FFT算法的资料,大家在做频谱显示时,用到最多的就是CMSIS DSP的算法库。于是乎,采用CMSIS DSP的库貌似是首选。

    但是不用不知道,一用才发现,由于CMSIS DSP的库使用的是查表的方式,我的64K Flash的LPC845轻轻松松就被撑爆了。没办法,只能改用其他方案。经过不懈的查阅资料,在GitHub找到一份FFT算法的代码,这个代码写的非常简洁,而且用起来很好用,感谢发布者pyrohaz,下面是FFT代码的一部分。
    1. /*
    2.   FIX_MPY() - fixed-point multiplication & scaling.
    3.   Substitute inline assembly for hardware-specific
    4.   optimization suited to a particluar DSP processor.
    5.   Scaling ensures that result remains 16-bit.
    6. */
    7. inline short FIX_MPY(short a, short b)
    8. {
    9.     /* shift right one less bit (i.e. 15-1) */
    10.     int c = ((int)a * (int)b) >> 14;
    11.     /* last bit shifted out = rounding-bit */
    12.     b = c & 0x01;
    13.     /* last shift + rounding bit */
    14.     a = (c >> 1) + b;
    15.     return a;
    16. }
    复制代码
    fix_fft(short fr[], short fi[], short m, short inverse)函数,FFT计算函数
    1. int fix_fft(short fr[], short fi[], short m, short inverse)
    2. {
    3.     int mr, nn, i, j, l, k, istep, n, scale, shift;
    4.     short qr, qi, tr, ti, wr, wi;

    5.     n = 1 << m;

    6.     /* max FFT size = N_WAVE */
    7.     if (n > N_WAVE)
    8.         return -1;

    9.     mr = 0;
    10.     nn = n - 1;
    11.     scale = 0;

    12.     /* decimation in time - re-order data */
    13.     for (m=1; m<=nn; ++m) {
    14.         l = n;
    15.         do {
    16.             l >>= 1;
    17.         } while (mr+l > nn);
    18.         mr = (mr & (l-1)) + l;

    19.         if (mr <= m)
    20.             continue;
    21.         tr = fr[m];
    22.         fr[m] = fr[mr];
    23.         fr[mr] = tr;
    24.         ti = fi[m];
    25.         fi[m] = fi[mr];
    26.         fi[mr] = ti;
    27.     }      
    复制代码
    接 fix_fft(short fr[], short fi[], short m, short inverse)函数
    1.   l = 1;
    2.     k = LOG2_N_WAVE-1;
    3.     while (l < n) {
    4.         if (inverse) {
    5.             /* variable scaling, depending upon data */
    6.             shift = 0;
    7.             for (i=0; i<n; ++i) {
    8.                 j = fr[i];
    9.                 if (j < 0)
    10.                     j = -j;
    11.                 m = fi[i];
    12.                 if (m < 0)
    13.                     m = -m;
    14.                 if (j > 16383 || m > 16383) {
    15.                     shift = 1;
    16.                     break;
    17.                 }
    18.             }
    19.             if (shift)
    20.                 ++scale;
    21.         } else {
    22.             /*
    23.                 fixed scaling, for proper normalization --
    24.                 there will be log2(n) passes, so this results
    25.                 in an overall factor of 1/n, distributed to
    26.                 maximize arithmetic accuracy.
    27.             */
    28.             shift = 1;
    29.         }
    复制代码
    接fix_fftr(short f[], int m, int inverse)函数
    1.   /*
    2.             it may not be obvious, but the shift will be
    3.             performed on each data point exactly once,
    4.             during this pass.
    5.         */
    6.         istep = l << 1;
    7.         for (m=0; m<l; ++m) {
    8.             j = m << k;
    9.             /* 0 <= j < N_WAVE/2 */
    10.             wr =  Sinewave[j+N_WAVE/4];
    11.             wi = -Sinewave[j];
    12.             if (inverse)
    13.                 wi = -wi;
    14.             if (shift) {
    15.                 wr >>= 1;
    16.                 wi >>= 1;
    17.             }
    18.             for (i=m; i<n; i+=istep) {
    19.                 j = i + l;
    20.                 tr = FIX_MPY(wr,fr[j]) - FIX_MPY(wi,fi[j]);
    21.                 ti = FIX_MPY(wr,fi[j]) + FIX_MPY(wi,fr[j]);
    22.                 qr = fr[i];
    23.                 qi = fi[i];
    24.                 if (shift) {
    25.                     qr >>= 1;
    26.                     qi >>= 1;
    27.                 }
    28.                 fr[j] = qr - tr;
    29.                 fi[j] = qi - ti;
    30.                 fr[i] = qr + tr;
    31.                 fi[i] = qi + ti;
    32.             }
    33.         }
    34.         --k;
    35.         l = istep;
    36.     }
    37.     return scale;
    38. }
    复制代码
    1. /*
    2.   fix_fftr() - forward/inverse FFT on array of real numbers.
    3.   Real FFT/iFFT using half-size complex FFT by distributing
    4.   even/odd samples into real/imaginary arrays respectively.
    5.   In order to save data space (i.e. to avoid two arrays, one
    6.   for real, one for imaginary samples), we proceed in the
    7.   following two steps: a) samples are rearranged in the real
    8.   array so that all even samples are in places 0-(N/2-1) and
    9.   all imaginary samples in places (N/2)-(N-1), and b) fix_fft
    10.   is called with fr and fi pointing to index 0 and index N/2
    11.   respectively in the original array. The above guarantees
    12.   that fix_fft "sees" consecutive real samples as alternating
    13.   real and imaginary samples in the complex array.
    14. */
    15. int fix_fftr(short f[], int m, int inverse)
    16. {
    17.     int i, N = 1<<(m-1), scale = 0;
    18.     short tt, *fr=f, *fi=&f[N];

    19.     if (inverse)
    20.         scale = fix_fft(fi, fr, m-1, inverse);
    21.     for (i=1; i<N; i+=2) {
    22.         tt = f[N+i-1];
    23.         f[N+i-1] = f[i];
    24.         f[i] = tt;
    25.     }
    26.     if (! inverse)
    27.         scale = fix_fft(fi, fr, m-1, inverse);
    28.     return scale;
    29. }
    复制代码
    int fix_fft(short fr[], short fi[], short m, short inverse) 是FFT算法的计算函数,fr[]是ADC采集到信号值的实部,fi[]是ADC采集到信号值的虚部。经过fix_fft函数处理之后,fr[]是FFT计算所得实部,fi[]是计算所得的虚部。

    我们最终要显示的音乐频谱其实是FFT频域中音频的幅值,幅值的计算是实部的平方+虚部的平方开根号。下面是具体的幅值计算部分的代码,每一个幅值点对应OLED的一列像素点。
    1. //Calculate the magnitude
    2.     for(Cnt = 0; Cnt<Len; Cnt++){
    3.         //Square the real components
    4.         R = (int32_t)Samples[Cnt]*(int32_t)Samples[Cnt];

    5.         //Square the imaginary components
    6.         I = (int32_t)Imag[Cnt]*(int32_t)Imag[Cnt];

    7.         //Calculate the magnitude and sum into BufSum. As there
    8.         //are more samples than columns, sum the magnitude of
    9.         //the bins. E.g. for 128 samples and 64 pixels on the LCD,
    10.         //each column should be the sum of 2 bins (128/64). If there
    11.         //are 512 samples and 64 pixels on the LCD, each column should
    12.         //be the sum of 8 bins (512/64)
    13.         BufSum += isqrt(R + I);

    14.         //Calculate the index for the column
    15.         Index = Cnt;

    16.         //Reset BufSum once the index has changed
    17.         if(Index != IndO){
    18.             //Write buffer sum to the column. Scale BufSum dependent on
    19.             //input amplitude range. If 'ColumnFilter' is more than 0,
    20.             //low pass filter the value to smooth out erratic changes.
    21.             if(Decibels){
    22.                 Col[Index] += (int32_t)(20.0f*log10f(BufSum)-Col[Index])>>ColumnFilter; //calculate the DB
    23.             }
    24.             else{
    25.                 Col[Index] += (BufSum-Col[Index])>>ColumnFilter;  //calculate the amplitude
    26.             }

    27.             //Limit maximum column value
    28.             if(Col[Index] >= YPix-9) Col[Index] = YPix-10;

    29.             IndO = Index;
    30.             BufSum = 0;
    31.         }
    32.     }
    复制代码







    签到www.nxpic.org.cn
    回复

    使用道具 举报

  • TA的每日心情
    开心
    2022-5-26 10:41
  • 签到天数: 40 天

    [LV.5]常住居民I

    39

    主题

    385

    帖子

    0

    金牌会员

    Rank: 6Rank: 6

    积分
    1131

    热心会员

    最后登录
    2022-12-7
    发表于 2022-9-22 10:30:26 | 显示全部楼层
    这个看着不错
    签到
    回复 支持 反对

    使用道具 举报

  • TA的每日心情
    奋斗
    5 小时前
  • 签到天数: 1472 天

    [LV.10]以坛为家III

    37

    主题

    9075

    帖子

    3

    版主

    Rank: 7Rank: 7Rank: 7

    积分
    12549
    最后登录
    2022-12-8
    发表于 2022-9-22 10:33:13 | 显示全部楼层
    有机会试试
    该会员没有填写今日想说内容.
    回复 支持 反对

    使用道具 举报

  • TA的每日心情
    奋斗
    昨天 08:17
  • 签到天数: 1057 天

    [LV.10]以坛为家III

    200

    主题

    1万

    帖子

    64

    超级版主

    Rank: 8Rank: 8

    积分
    59585
    最后登录
    2022-12-7
    发表于 2022-9-22 10:57:55 | 显示全部楼层
    你们能看清到视频嘛?
    该会员没有填写今日想说内容.
    回复 支持 反对

    使用道具 举报

  • TA的每日心情
    开心
    2022-11-30 13:07
  • 签到天数: 209 天

    [LV.7]常住居民III

    2733

    主题

    5382

    帖子

    0

    管理员

    Rank: 9Rank: 9Rank: 9

    积分
    25599
    最后登录
    2022-12-7
     楼主| 发表于 2022-9-22 12:39:46 | 显示全部楼层
    stm1024 发表于 2022-9-22 10:57
    你们能看清到视频嘛?

    12.png
    可以吧
    签到www.nxpic.org.cn
    回复 支持 反对

    使用道具 举报

  • TA的每日心情
    奋斗
    昨天 09:05
  • 签到天数: 1357 天

    [LV.10]以坛为家III

    64

    主题

    5848

    帖子

    1

    金牌会员

    Rank: 6Rank: 6

    积分
    10298
    最后登录
    2022-12-7
    发表于 2022-9-22 14:03:51 | 显示全部楼层
    很不错。。。
    回复

    使用道具 举报

  • TA的每日心情
    开心
    2019-3-5 08:47
  • 签到天数: 1 天

    [LV.1]初来乍到

    69

    主题

    1942

    帖子

    2

    金牌会员

    Rank: 6Rank: 6

    积分
    5962
    最后登录
    2022-12-6
    发表于 2022-9-22 14:33:20 | 显示全部楼层

    换个计算机试试,看不到视频,漆黑一片。
    加油哦
    回复 支持 反对

    使用道具 举报

  • TA的每日心情
    开心
    昨天 20:49
  • 签到天数: 786 天

    [LV.10]以坛为家III

    20

    主题

    9038

    帖子

    1

    金牌会员

    Rank: 6Rank: 6

    积分
    9520
    最后登录
    2022-12-7
    发表于 2022-9-22 15:32:03 | 显示全部楼层
    视频看不到
    跟着日天混,三天饱九顿!
    回复 支持 反对

    使用道具 举报

  • TA的每日心情
    奋斗
    昨天 11:58
  • 签到天数: 1577 天

    [LV.Master]伴坛终老

    17

    主题

    3884

    帖子

    5

    金牌会员

    Rank: 6Rank: 6

    积分
    7850
    最后登录
    2022-12-7
    发表于 2022-9-22 20:43:56 | 显示全部楼层
    收藏了。。
    该会员没有填写今日想说内容.
    回复

    使用道具 举报

  • TA的每日心情
    慵懒
    昨天 00:19
  • 签到天数: 856 天

    [LV.10]以坛为家III

    16

    主题

    3181

    帖子

    0

    金牌会员

    Rank: 6Rank: 6

    积分
    5499

    活跃会员

    最后登录
    2022-12-8
    发表于 2022-9-23 12:31:38 | 显示全部楼层
    做一个成品当转盘一等奖
    该会员没有填写今日想说内容.
    回复 支持 反对

    使用道具 举报

    您需要登录后才可以回帖 注册/登录

    本版积分规则

    关闭

    站长推荐上一条 /3 下一条

    Archiver|手机版|小黑屋|恩智浦技术社区

    GMT+8, 2022-12-8 05:01 , Processed in 0.107456 second(s), 33 queries , MemCache On.

    Powered by Discuz! X3.4

    Copyright © 2001-2021, Tencent Cloud.

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