1 音频信号频谱分析简介
1.1 设计背景
在当今数字化时代的背景下，由于DSP具有处理速度快，功耗低，性能好，存储容量大等特点，因此已被广泛应用于通信、计算机、消费类电子产品等多个领域中。
而对于音频信号的分析来说，单纯只对其时域进行分析，有时候往往会不尽如人意，而如果换一种角度，从频域上对其进行分析，有时可能会有惊人的效果。如在一段嘈杂的音频中对其频域进行分析，便有可能从中确定声音的来源，再通过对其进行相应的滤波处理便能够很好地提取到对我们有用的信息。
因此通过DSP对音频信号进行频谱分析将会具有很好的实际意义和应用前景。
1.2 设计方法
本次课程设计将采用AIC23音频处理芯片完成对音频信号的采集，同时通过TMS320F28335对采集进来的音频信号进行FFT处理并将结果显示于CCS中从而完成对音频信号的频谱分析。
1.3 设计要求
通过A/D采集音频信号(fmax<10kHz)，编写DSP的FFT处理程序（自定频谱分辨力），获得幅频特性之后，在点阵液晶中大致显示出幅频图。并在液晶中用文字显示频率幅值前三的频率值。
1、DSP与A/D接口电路的原理图绘制；
2、DSP控制A/D的程序编写与调试；
3、编写DSP的FFT处理程序；
4、控制点阵液晶，实现绘图功能，将幅频图显示出来
5、按要求编写课程设计报告书，正确、完整的阐述设计和实验结果。
6、在报告中绘制程序的流程图，并文字说明。
2 硬件设计
2.1 TMS320F28335最小系统
TMS320F28335最小系统是指能使DSP正常工作的最简电路，如图2.1所示，其中大体包括以下四个部分：
- 主芯片；
- 晶振电路模块：外接30MHz的晶振，内部可通过锁相环PLL将其倍频至150MHz作为系统的时钟频率；
- 复位电路模块：复位时低电平有效。上电后，由于电容两端电压不能突变的缘故，所以XRS端会在很短的一段时间内保持低电平，再被拉到高电平。复位后，DSP会先响应复位中断信号进入到中断向量表中，再通过中断向量表的指引执行BootLoader程序，完成对DSP的初始化工作；
- JTAG电路模块：用于下载程序和对程序进行仿真。
图2.1 TMS320F28335最小系统
2.2 音频编解码模块
2.2.1 TLV320AIC23芯片介绍
TLV320AIC23是一个高性能的多媒体数字语音编解码器，它的内部ADC和DAC转换模块带有完整的数字滤波器。详细指标如下：
①在 ADC 采用 96kHz 采样率时噪音 90-dB；
②在 DAC 采用 96kHz 采样率时噪音 100-dB；
③ 1.42V – 3.6V 核心数字电压：兼容 TI F28x DSP 内核电压；
④ 2.7V – 3.6V 缓冲器和模拟：兼容 TI F28x DSP 内核电压；
⑤支持 8kHz – 96kHz 的采样频率；
⑥软件控制通过 TI IIC 接口；
⑦音频数据输入输出通过 TI McBSP 接口。
2.2.2 TLV320AIC23硬件连接
如图2.2.2所示，LINE_IN2为音频输出口，可接耳机或音响播放音频；LINE_IN3为语音输入口，可接话筒。同时F28335通过IIC接口和McBSP接口对TLV320AIC23进行控制。
图2.2.2 TLV320AIC32硬件连接
3 软件设计
3.1 音频编解码模块初始化
如表3.1.1所示，根据芯片的Register Map，可通过IIC协议修改芯片内对应寄存器的值。
表3.1.1 TLV320AIC23芯片Register Map
|
ADDRESS |
REGISTER |
|
0000000 |
Left line input channel volume control |
|
0000001 |
Right line input channel volume control |
|
0000010 |
Left channel headphone volume control |
|
0000011 |
Right channel headphone volume control |
|
0000100 |
Analog audio path control |
|
0000101 |
Digital audio path control |
|
0000110 |
Power down control |
|
0000111 |
Digital audio interface format |
|
0001000 |
Sample rate control |
|
0001001 |
Digital interface activation |
|
0001111 |
Reset register |
如表3.1.2所示，再根据芯片手册所提供的对”Sample rate control”的设定，即可将采样速率设置为8kHz。至此，音频编解码模块的初始化工作便已基本结束。
表3.1.2 TLV320AIC23芯片Sample rate control
|
SAMPLING RATE |
FILTER TYPE |
SAMPLING-RATE CONTROL SETTINGS |
|||||
|
ADC |
DAC |
SR3 |
SR2 |
SR1 |
SR0 |
BOSR |
|
|
96 |
96 |
3 |
0 |
1 |
1 |
1 |
0 |
|
88.2 |
88.2 |
2 |
1 |
1 |
1 |
1 |
1 |
|
48 |
48 |
0 |
0 |
0 |
0 |
0 |
0 |
|
44.1 |
44.1 |
1 |
1 |
0 |
0 |
0 |
1 |
|
32 |
32 |
0 |
0 |
1 |
1 |
0 |
0 |
|
8.021 |
8.021 |
1 |
1 |
0 |
1 |
1 |
1 |
|
8 |
8 |
0 |
0 |
0 |
1 |
1 |
1 |
|
48 |
8 |
0 |
0 |
0 |
0 |
1 |
0 |
|
44.1 |
8.021 |
1 |
1 |
0 |
0 |
1 |
1 |
|
8 |
48 |
0 |
0 |
0 |
1 |
0 |
0 |
|
8.021 |
44.1 |
1 |
1 |
0 |
1 |
0 |
1 |
3.2 音频编解码模块中断程序
由于本次FFT采用的点数是256，点数较多，运算量较大，因此不可能每次进入中断后都对音频进行FFT的处理，因此为了解决这个问题，我通过将每次采集到的音频都用FIFO的方法放入一个数组中储存起来（也即是队列），再判断新增音频的数量（也即是进入中断的次数），如果新增音频的数量等于FFT点数则对整个数组进行FFT处理并将记录新增音频数的变量清0；如果不是，则退出中断，如图3.2所示。这样便能避免每新增一个音频点便进行一次FFT计算这样低效率的冗余操作了。
图3.2 音频编解码模块中断程序流程图
下面对FFT程序进行介绍。注：FFT程序是先在Visual Studio中利用图形化界面调试好，再移植到CCS中的。
3.2.1 基本复数运算的定义
由于FFT程序涉及到大量的复数运算，因此要想进行FFT，首先最基本的便是定义基本的复数运算，如程序3.2.1所示，这里我们分别定义了复数的加法运算、减法运算和乘法运算。
1 struct complex c_add(struct complex a, struct complex b)//复数加法 2 { 3 struct complex c; 4 c.imag = a.imag + b.imag; 5 c.real = a.real + b.real; 6 return c; 7 } 8 struct complex c_sub(struct complex a, struct complex b)//复数减法 9 { 10 struct complex c; 11 c.imag = a.imag - b.imag; 12 c.real = a.real - b.real; 13 return c; 14 } 15 struct complex c_mul(struct complex a, struct complex b)//复数乘法 16 { 17 struct complex c; 18 c.real = a.real*b.real - a.imag*b.imag; 19 c.imag = a.real*b.imag + a.imag*b.real; 20 return c; 21 }
程序3.2.1 基本复数运算的定义程序
3.2.2 旋转因子的定义
旋转因子的数学表达式如下，根据该表达式便可写出相应函数，如程序3.2.2所示。
$W_{N}^{kn}=e^{j\left( -\frac{2\pi}{N} \right) kn}$
1 struct complex get_WN_kn(double n, double N, double k)//获得旋转因子 2 { 3 int i = 0; 4 struct complex ex, W; 5 ex.real = 1.0*cos(2 * PI / N); 6 ex.imag = -1.0*sin(2 * PI / N); 7 W.real = 1; W.imag = 0; 8 for (i = 0; i < k*n; i++) 9 W = c_mul(ex, W); 10 return W; 11 }
程序3.2.2 旋转因子的定义程序
3.2.3 码位颠倒函数
由于在FFT运算流图中，第一步便是把输入序列以倒序的形式输入，因此如程序3.2.3所示，该函数的大致思路是将原码通过右移的方法一位一位地送到倒序的码位上，最终实现码位颠倒，如图3.2.3所示。
图3.2.3 码位颠倒示意图
1 void change(struct complex x[]) 2 { 3 int i = 0, j = 0, k = 0; 4 double t; 5 struct complex temp; 6 for (i = 0; i<fft_N; i++) 7 { 8 k = i; j = 0; 9 t = (log(fft_N) / log(2)); 10 while ((t--)>0) //利用按位与以及循环实现码位颠倒 11 { 12 j = j << 1; 13 j |= (k & 1); 14 k = k >> 1; 15 } 16 if (j>i) //将x(n)的码位互换 17 { 18 temp = x[i]; 19 x[i] = x[j]; 20 x[j] = temp; 21 } 22 } 23 }
程序3.2.3 码位颠倒函数程序
3.2.4 FFT之蝶形运算
首先调用码位颠倒函数，再根据蝶形运算的方法和流程进行FFT的最终运算，如图3.2.4和程序3.2.4所示。
图3.2.4 蝶形运算示意图
1 void fft(struct complex x[], int N) 2 { 3 change(x); 4 struct complex W, up, down, product; 5 int i, l, j, k; 6 for (i = 0; i<log(fft_N) / log(2); i++) /*一级蝶形运算*/ 7 { 8 l = 1 << i; 9 for (j = 0; j<fft_N; j += 2 * l) /*一组蝶形运算*/ 10 { 11 for (k = 0; k<l; k++) /*一个蝶形运算*/ 12 { 13 W = get_WN_kn(1, fft_N, k * fft_N / (2 * l)); 14 product = c_mul(x[j + k + l], W); 15 up = c_add(x[j + k], product); 16 down = c_sub(x[j + k], product); 17 x[j + k] = up; 18 x[j + k + l] = down; 19 } 20 } 21 } 22 }
程序3.2.4 FFT之蝶形运算程序
3.2.5 IFFT函数
该函数在本次课设中其实并没有起到实际的作用，但是在调试过程中却能用于检验FFT函数的正确性，结果如图3.2.5所示。
图3.2.5 FFT&IFFT函数在VS上的运行结果
1 void ifft(struct complex x[], int N) 2 { 3 int i; 4 for (i = 0; i<fft_N; i++) 5 x[i].imag = -x[i].imag; 6 fft(x, fft_N); 7 for (i = 0; i<fft_N; i++) 8 { 9 x[i].real = x[i].real*1.0 / fft_N; 10 x[i].imag = -x[i].imag*1.0 / fft_N; 11 } 12 }
程序3.2.5.2 IFFT函数程序
4 运行结果
通过XDS100v2仿真器运行程序，并使用CCS中的画图工具，对经过FFT处理后的数组数据进行绘制便可得到以下结果。
4.1 “1000Hz+2500Hz”正弦波运行结果
如图4.1所示，该段音频是通过Matlab编写函数来产生的，Matlab函数请见附录A。
图4.1 “1000Hz+2500Hz”正弦波运行结果
4.2 “1500Hz+2000Hz”正弦波运行结果
如图4.2所示，该段音频是通过Matlab编写函数来产生的，Matlab函数请见附录A。
图4.2 “1500Hz+2000Hz”正弦波运行结果
4.3 “1700Hz+1800Hz”正弦波运行结果
如图4.3所示，该段音频是通过Matlab编写函数来产生的，Matlab函数请见附录A。
图4.3 “1700Hz+1800Hz”正弦波运行结果
4.4 “光辉岁月”音乐运行结果
如图4.4所示，光辉岁月这首歌的频率大多集中在342Hz这个频率左右。
图4.4 “光辉岁月”音乐运行结果
4.5 “Opera 2″音乐运行结果
如图4.5所示，Opera 2这首歌的最高频率可以达到2090Hz。
图4.5 “Opera 2″音乐运行结果
4.6 “相爱很难”音乐运行结果
如图4.6所示，相爱很难这首歌包含了两个主要频率，男音部分的频率大概是228Hz左右，女音部分的频率大概是532Hz左右。
图4.6 “相爱很难”音乐运行结果
附录A Matlab源程序代码
f1=600;%正弦频率 f2=200;%正弦频率 fs=22050;%采样频率 T=5;%总时长 x=1:T*fs; y1=sin(f1*2*pi*x/fs); y2=sin(f2*2*pi*x/fs); y=y1; sound(y,fs)
附录B DSP源程序代码
1 #include "math.h" 2 #define PI 3.1415926 3 #define fft_N 256 4 struct complex 5 { 6 double real; 7 double imag; 8 }; 9 struct complex c_mul(struct complex a, struct complex b); 10 void change(struct complex x[]); 11 struct complex get_WN_kn(double r, double N, double k); 12 struct complex c_add(struct complex a, struct complex b); 13 struct complex c_sub(struct complex a, struct complex b); 14 void fft(struct complex x[], int N); 15 void ifft(struct complex x[], int N); 16 void fft1(struct complex x[], struct complex y[]); 17 double abs2(struct complex a); 18 double result[fft_N]; 19 struct complex x[fft_N]; 20 unsigned int i, j; 21 struct complex fft_complex[fft_N]; 22 struct complex fft_result[fft_N]; 23 double fft_double[fft_N]; 24 #include "DSP2833x_Device.h" // DSP2833x Headerfile Include File 25 #include "DSP2833x_Examples.h" // DSP2833x Examples Include File 26 #include "leds.h" 27 #include "time.h" 28 #include "uart.h" 29 #include "stdio.h" 30 void I2CA_Init(void); 31 Uint16 AIC23Write(int Address,int Data); 32 void Delay(int time); 33 void Mydelay(Uint32 k); 34 interrupt void ISRMcbspSend(); 35 #define LuYin GpioDataRegs.GPADAT.bit.GPIO12 36 #define LuYin_ST GpioDataRegs.GPADAT.bit.GPIO13 37 #define BoYin GpioDataRegs.GPADAT.bit.GPIO14 38 Uint16 temp,i,j,m=0; 39 Uint32 in_H,in_L; 40 Uint16 out_H,out_L; 41 Uint32 in[fft_N]; 42 Uint32 out[fft_N]; 43 void main() 44 { 45 InitSysCtrl(); 46 InitPieCtrl(); 47 IER = 0x0000; 48 IFR = 0x0000; 49 InitPieVectTable(); 50 51 LED_Init(); 52 TIM0_Init(150,200000);//200ms 53 UARTa_Init(4800); 54 55 InitMcbspaGpio(); //zq 56 InitI2CGpio(); 57 58 I2CA_Init(); 59 60 AIC23Write(0x00,0x00); 61 Delay(100); 62 AIC23Write(0x02,0x00); 63 Delay(100); 64 AIC23Write(0x04,0x7f); 65 Delay(100); 66 AIC23Write(0x06,0x7f); 67 Delay(100); 68 AIC23Write(0x08,0x14); 69 Delay(100); 70 AIC23Write(0x0A,0x00); 71 Delay(100); 72 AIC23Write(0x0C,0x00); 73 Delay(100); 74 AIC23Write(0x0E,0x43); 75 Delay(100); 76 AIC23Write(0x10,0x23); 77 Delay(100); 78 AIC23Write(0x12,0x01); 79 Delay(100); //AIC23Init 80 81 InitMcbspa(); // Initalize the Mcbsp-A 82 83 EALLOW; // This is needed to write to EALLOW protected registers 84 PieVectTable.MRINTA = &ISRMcbspSend; 85 EDIS; // This is needed to disable write to EALLOW protected registers 86 87 PieCtrlRegs.PIECTRL.bit.ENPIE = 1; // Enable the PIE block 88 PieCtrlRegs.PIEIER6.bit.INTx5=1; // Enable PIE Group 6, INT 5 89 IER |= M_INT6; // Enable CPU INT6 90 91 EINT; // Enable Global interrupt INTM 92 ERTM; // Enable Global realtime interrupt DBGM 93 94 for (i = 0; i <= fft_N - 1; i++)//产生待测序列，模拟在一个周期内采样的点 95 { 96 x[i].real = 100* sin(10 * i / (1.0*fft_N)*2.0*PI); 97 x[i].imag = 0; 98 } 99 fft(x, fft_N); 100 for (i = 0; i < fft_N; i++) 101 { 102 result[i] = abs2(x[i]); 103 } 104 while(1) 105 { 106 } 107 } 108 109 void I2CA_Init(void) 110 { 111 // Initialize I2C 112 I2caRegs.I2CSAR = 0x001A; // Slave address - EEPROM control code 113 114 #if (CPU_FRQ_150MHZ) // Default - For 150MHz SYSCLKOUT 115 I2caRegs.I2CPSC.all = 14; // Prescaler - need 7-12 Mhz on module clk (150/15 = 10MHz) 116 #endif 117 #if (CPU_FRQ_100MHZ) // For 100 MHz SYSCLKOUT 118 I2caRegs.I2CPSC.all = 9; // Prescaler - need 7-12 Mhz on module clk (100/10 = 10MHz) 119 #endif 120 121 I2caRegs.I2CCLKL = 100; // NOTE: must be non zero 122 I2caRegs.I2CCLKH = 100; // NOTE: must be non zero 123 I2caRegs.I2CIER.all = 0x24; // Enable SCD & ARDY interrupts 124 125 // I2caRegs.I2CMDR.all = 0x0020; // Take I2C out of reset 126 I2caRegs.I2CMDR.all = 0x0420; // Take I2C out of reset //zq 127 // Stop I2C when suspended 128 129 I2caRegs.I2CFFTX.all = 0x6000; // Enable FIFO mode and TXFIFO 130 I2caRegs.I2CFFRX.all = 0x2040; // Enable RXFIFO, clear RXFFINT, 131 132 return; 133 } 134 135 Uint16 AIC23Write(int Address,int Data) 136 { 137 138 139 if (I2caRegs.I2CMDR.bit.STP == 1) 140 { 141 return I2C_STP_NOT_READY_ERROR; 142 } 143 144 // Setup slave address 145 I2caRegs.I2CSAR = 0x1A; 146 147 // Check if bus busy 148 if (I2caRegs.I2CSTR.bit.BB == 1) 149 { 150 return I2C_BUS_BUSY_ERROR; 151 } 152 153 // Setup number of bytes to send 154 // MsgBuffer + Address 155 I2caRegs.I2CCNT = 2; 156 I2caRegs.I2CDXR = Address; 157 I2caRegs.I2CDXR = Data; 158 // Send start as master transmitter 159 I2caRegs.I2CMDR.all = 0x6E20; 160 return I2C_SUCCESS; 161 162 } 163 164 void Delay(int time) 165 { 166 int i,j,k=0; 167 for(i=0;i<time;i++) 168 for(j=0;j<1024;j++) 169 k++; 170 } 171 172 void Mydelay(Uint32 k) 173 { 174 while(k--); 175 } 176 177 interrupt void ISRMcbspSend(void) 178 { 179 int temp1,temp2; 180 181 temp1=McbspaRegs.DRR1.all+100; 182 temp2=McbspaRegs.DRR2.all; 183 184 in_L=McbspaRegs.DRR1.all+100; 185 in_H=McbspaRegs.DRR2.all; 186 for(j=0;j<fft_N-1;j++) 187 { 188 in[j]=in[j+1]; 189 } 190 in[fft_N-1]=(in_H<<16)+in_L; 191 for(j=0;j<fft_N;j++) 192 { 193 fft_complex[j].real=in[j]; 194 fft_complex[j].imag=0; 195 } 196 197 if(m >= fft_N) 198 { 199 fft1(fft_complex, fft_result); 200 for(j=0;j<fft_N;j++) 201 { 202 fft_double[j]=abs2(fft_result[j]); 203 } 204 m=0; 205 } 206 m++; 207 McbspaRegs.DXR1.all = temp1; //放音 208 McbspaRegs.DXR2.all = temp2; 209 210 PieCtrlRegs.PIEACK.all = 0x0020; 211 // PieCtrlRegs.PIEIFR6.bit.INTx5 = 0; 212 // ERTM; 213 } 214 215 void ifft(struct complex x[], int N) 216 { 217 int i; 218 for (i = 0; i<fft_N; i++) 219 { 220 x[i].imag = -x[i].imag; 221 } 222 fft(x, fft_N); 223 224 for (i = 0; i<fft_N; i++) 225 { 226 x[i].real = x[i].real*1.0 / fft_N; 227 x[i].imag = -x[i].imag*1.0 / fft_N; 228 } 229 } 230 void fft(struct complex x[], int N) 231 { 232 change(x); 233 struct complex W; 234 struct complex up, down, product; 235 int i, l, j, k; 236 for (i = 0; i<log(fft_N) / log(2); i++) /*一级蝶形运算 stage */ 237 { 238 l = 1 << i; 239 for (j = 0; j<fft_N; j += 2 * l) /*一组蝶形运算 group,每个group的蝶形因子乘数不同*/ 240 { 241 for (k = 0; k<l; k++) /*一个蝶形运算 每个group内的蝶形运算*/ 242 { 243 W = get_WN_kn(1, fft_N, k * fft_N / (2 * l)); 244 product = c_mul(x[j + k + l], W); 245 up = c_add(x[j + k], product); 246 down = c_sub(x[j + k], product); 247 x[j + k] = up; 248 x[j + k + l] = down; 249 } 250 } 251 } 252 } 253 void fft1(struct complex x[], struct complex y[]) 254 { 255 int n; 256 for (n = 0; n < fft_N; n++) 257 { 258 y[n] = x[n]; 259 } 260 change(y); 261 struct complex W; 262 struct complex up, down, product; 263 int i, l, j, k; 264 for (i = 0; i<log(fft_N) / log(2); i++) /*一级蝶形运算 stage */ 265 { 266 l = 1 << i; 267 for (j = 0; j<fft_N; j += 2 * l) /*一组蝶形运算 group,每个group的蝶形因子乘数不同*/ 268 { 269 for (k = 0; k<l; k++) /*一个蝶形运算 每个group内的蝶形运算*/ 270 { 271 W = get_WN_kn(1, fft_N, k * fft_N / (2 * l)); 272 product = c_mul(y[j + k + l], W); 273 up = c_add(y[j + k], product); 274 down = c_sub(y[j + k], product); 275 y[j + k] = up; 276 y[j + k + l] = down; 277 } 278 } 279 } 280 } 281 void change(struct complex x[]) 282 { 283 int i = 0, j = 0, k = 0; 284 double t; 285 struct complex temp; 286 for (i = 0; i<fft_N; i++) 287 { 288 k = i; j = 0; 289 t = (log(fft_N) / log(2)); 290 while ((t--)>0) //利用按位与以及循环实现码位颠倒 291 { 292 j = j << 1; 293 j |= (k & 1); 294 k = k >> 1; 295 } 296 if (j>i) //将x(n)的码位互换 297 { 298 temp = x[i]; 299 x[i] = x[j]; 300 x[j] = temp; 301 } 302 } 303 } 304 struct complex c_mul(struct complex a, struct complex b)//复数乘法 305 { 306 struct complex c; 307 c.real = a.real*b.real - a.imag*b.imag; 308 c.imag = a.real*b.imag + a.imag*b.real; 309 return c; 310 } 311 struct complex c_add(struct complex a, struct complex b)//复数加法 312 { 313 struct complex c; 314 c.imag = a.imag + b.imag; 315 c.real = a.real + b.real; 316 return c; 317 } 318 struct complex c_sub(struct complex a, struct complex b)//复数减法 319 { 320 struct complex c; 321 c.imag = a.imag - b.imag; 322 c.real = a.real - b.real; 323 return c; 324 } 325 struct complex get_WN_kn(double n, double N, double k)//获得旋转因子 326 { 327 int i = 0; 328 struct complex ex; 329 struct complex W; 330 ex.real = 1.0*cos(2 * PI / N); 331 ex.imag = -1.0*sin(2 * PI / N); 332 W.real = 1; 333 W.imag = 0; 334 for (i = 0; i < k*n; i++) 335 { 336 W = c_mul(ex, W); 337 } 338 return W; 339 } 340 double abs2(struct complex a) 341 { 342 double y; 343 y = a.real * a.real + a.imag * a.imag; 344 return y; 345 }
原文地址：http://www.cnblogs.com/Kelvin-Wu/p/16795037.html
