分享一下18,19届声音信标软件部分处理的思路,该部分只是本人对声音处理部分的思考,如发现错误,欢迎大家在评论区或私信指正。
一、比赛规则![](https://i-blog.csdnimg.cn/direct/2c05ce5cfbe44b7bb6fcd45e35df4390.png)
越野组如果在室内比赛, 最终退化成声音信标组; 不难看出声音处理在信标组的两届智能车比赛中发挥着至关重要的作用。
硬件:
(1)主控:Infineon微控制器cyt4bb7
(2)IMU: icm963
(3)gps:TAU1201
(4)车模:R车
(5)无线通信模块:4C39
...
任务:
使用硅麦对发出声音的信标进行定位。
二、TDOA声音定位算法
在小车运动中,如何求出声源角度是一个至关重要的问题, 由于 Chirp 信号的声音长度在 0.2048 秒,其频率在 250-2000HZ。于是我们开启一个 100us 的中断,以10KHZ 来采集声音信号。并使用一个长度为 2048 的数组来存储声音信号,这样每采集 2048 次,就使用一次完整的声音信号用于声源角度的解算, 对于角度的解算,我们首先要做的是正确求出硅麦所收集的声音信号之间的时延。而对于时延的获取,最为广泛使用的就是进行互相关运算。互相关运算可看作为对两组信号在时域进行卷积运算,而对如此长的数组进行卷积运算是非常耗费时间的。而信号在频域上相乘等于在时域上进行卷积。为加计算速度,我们对信号进行 FFT 运算使其变化到频域,再在频域上进行相乘处理,再进行 IFFT 运算使其变换到时域,测量其峰值下标即可获得时延。
互相关算法通俗解释以及实现
互相关运算实际上就是计算在每一个时刻两个波形的相似程度,最相似的时候计算出来的值也最大,因此我们只需要找到最大值所对应的时刻,也就得到了时延差。例如对两个麦克风数据进行互相关,得到的时延差则表示声音到达两个麦克风的时间差。
互相关计算的具体实现如下:
麦克风0的信号用A表示,麦克风1的信号用B表示。A、B信号长度为2048个16位数据。由于iar不支持复数的表达,则用奇数代替虚部,偶数代替实部。再将实部剪掉对应的直流偏量(不减影响也不大)。
a= inputSignal_1[2*i] = ((float)g_adc3Data[AdcBuffIndex][i]*3.3/4096)-1.6;
inputSignal_1[2*i + 1]= 0.00;
b=inputSignal_2[2*i] = ((float)g_adc3Data[AdcBuffIndex][i]*3.3/4096)-1.6;//将实数信号变成交错复数的形式
inputSignal_2[2*i + 1]= 0.00;
再限幅。
if(inputSignal_1[2*i] >= 3.30)
{
inputSignal_1[2*i] = 3.30;
}
if(inputSignal_2[2*i] >= 3.30)
{
inputSignal_2[2*i] = 3.30;
}
对A与B信号进行补零得到A2与B2信号(补零实际上也就是后面增加2048个数据,数据全部是0而已),使得补零后的A2与B2信号长度等于A与B长度之和。
for(i = ADC_DATA_LEN; i < ADC_DATA_LEN*2; i++) //补零
{
inputSignal_1[2*i] = 0;
inputSignal_1[2*i + 1]= 0;
inputSignal_2[2*i] = 0;
inputSignal_2[2*i + 1]= 0;
}
对A2与B2都进行FFT运算。
arm_cfft_instance_f32 arm_cfft_instance_f32_len_1024; // 定义FFT对象
arm_cfft_init_f32(&arm_cfft_instance_f32_len_1024, ADC_DATA_LEN*2); //是否乘2?// 初始化FFT对象 赋予计算长度//如果你需要更高的频率分辨率,你可以增加 FFT 的长度
arm_cfft_f32(&arm_cfft_instance_f32_len_1024 , inputSignal_1 , 0 , 1); // 32位浮点FFT运算
arm_cfft_f32(&arm_cfft_instance_f32_len_1024 , inputSignal_2 , 0 , 1); // 32位浮点FFT运算
对A2信号取共轭得到A2_C。
arm_cmplx_conj_f32(inputSignal_2, Conjugate_date,FFT_SIZE); //共轭
再求A2_C与B2的乘积。
arm_cmplx_mult_cmplx_f32(inputSignal_1,Conjugate_date,Correlation_data,FFT_SIZE);//相乘
再对求出的乘积求逆FFT运算。
arm_cfft_f32(&arm_cfft_instance_f32_len_1024 , Correlation_data , 1 , 1); // 32位浮点FFT运算
继续对逆FFT后的结果进行位置调换,把数据前面一半与后面一半的位置进行交换(位置0与位置2048进行交换以此类推),如上图所示。
再取模。
arm_cmplx_mag_f32(Correlation_data,Correlation_mag_data,FFT_SIZE);//取模
最后在交换顺序后得到的结果中找到最大值(幅值)所对应的下标,用2048减去这个下标变最终得到了时延差,并且具有正负,正负可以很便于我们判断声源距离哪个麦克风更近一些。正数靠近A,负数靠近B。
arm_max_f32(Correlation_mag_data,2*FFT_SIZE, &MaxValue, &MaxValue_subscript);//找最大值
声音信号的滤波的处理 在声源角度的解算中,噪声的影响对于角度的解算的准确率有着极大的影响,列如车模运动时电机产生的噪音,车模本身运动时所产生的噪音等等。
在硬件上,为了减少噪声带来的影响,我们选择了将硅麦支架增高,从而尽量避免电机的噪声影响,并且硅麦支架的提高也可以更好的采集声音信号, 并在硅麦上套个防噪棉,避免风声影响。并通过改变结构来增加车模运动过程 中的稳定性,使其产生的噪声更小。 这些做法相比于在软件上进行滤波处理,本质上减少噪音的产生往往更加有效。
软件处理由于我们信号进行了 FFT 处理使其变换到频域,而滤波处理可直接在频域 上进行操作,因 Chirp 信号的声音频率在 250HZ-2000HZ。我们可以在频域上直接去除无关的频率片段,只保留真正的 Chirp 信号片段,这样在计算时可以滤掉 大部分的无关声音。使角度解算的准确率得到很大的提升。
for( i = 0; i < 250; i++) //理想带通滤波,保留200-2000Hz频率
{//只需对250-2000hz范围内信号做互相关 (10000hz/2048=4.88 //250/4.88=51 51/2=25)//相当于带通滤波
inputSignal_1[2*i] = 0.00;
inputSignal_1[2*i + 1] = 0.00;
inputSignal_2[2*i] = 0.00;
inputSignal_2[2*i + 1] = 0.00;
}
for(i = 2001; i < ADC_DATA_LEN*2; i++)//i为频率与长度对应的索引,索引=频率/ (采样率 / FFT长度)
{//只需对250-2000hz范围内信号做互相关 (10000hz/2048=4.88 //2000/4.88=409 409/2=204)//相当于带通滤波
inputSignal_1[2*i] = 0.00;
inputSignal_1[2*i + 1] = 0.00;
inputSignal_2[2*i] = 0.00;
inputSignal_2[2*i + 1] = 0.00;
}
另一个滤波方法就是广义互相关的使用,其公式为: 也就是在普通互相关时乘上一个加权函数,其中为加权函数。而加权函 数有 PHAT,ROTH,SCOT 等。我们使用了PHAT 函数进行处理。 而在实际测试中,我们将带 PHAT 和不带 PHAT 运算后的数组传入上位机看形,当信噪比较低时,其波形图并没有太大的区别。当信噪比较高的时候, 加权函数的作用才有一点作用。
/*path滤波*/
Correlation_data[2*i] = Correlation_data[2*i] * (float32_t)(1.00/(sqrt(Correlation_data[2*i] *Correlation_data[2*i] +Correlation_data[2*i+1] * Correlation_data[2*i+1])));
Correlation_data[2*i+1] = Correlation_data[2*i+1] * (float32_t)(1.00/(sqrt(Correlation_data[2*i] * Correlation_data[2*i] + Correlation_data[2*i+1] * Correlation_data[2*i+1])));
该算法解算一次完整的声音角度仅需25ms,快速结算可达到10ms以内,大家可以通过精简代码,优化算法来减少运算时间,时间越短,结算的角度越准,最后祝大家的小车都能吹到国赛的风。