滑轨声源定向的理论模型与参数估计

通过Chirp信号的接受来确定声源的位置和方向，通过实验来研究它的精度和那些因素有关系？

一、实验的设置：

实验分为声音的发送和接收装置。下面依次介绍发送和接收实验装置的设置。

1. 接收双声道麦克风的布局

使用了两个驻极体麦克风来接收声源的声音。它们的固定方式如下图所示。
左右两个麦克风的距离为26.5厘米，距离地面的高度为15厘米。

接收声音的左右声道麦克风的布局

2. 发送声音的扬声器和滑轨

固定在滑轨上的声源采用8Ω\OmegaΩ的动圈式扬声器。它的尺寸和固定的高度见下图所示。它外部固定一个与其尺寸相当的原来的铝质灯罩，可以减少扬声器外部声音断路，提高发送声音的效率。

固定扬声器的滑轨可以做维运动，运动的范围大约0.7米。驱动器可以通过WiFi接收到运动距离的指令。

实验平台

定义滑轨的起始位置和结束位置：
起始位置，近端：带有步进电机驱动、齿轮传动的一端成为起始端；
结束位置，远端：相对应的另外一段定义为结束位置；

控制命令的功能：
rccw(nstep) : 控制滑轨从起始位置运动到结束位置，nstep:移动的步数，取值范围（0~14000）。由于有惯性，如果分布运行时，总的运动步数减小。比如：移动100步，每次移动115，总共11500,基本上可以从起始端移动到结束端；

rcccw(nstep): 控制滑轨从结束端移动到起始端。

3. 驱动声源喇叭的功放

驱动声源扬声器的音频功率放大器是一款双声道的音箱音频放大器。适配电源为12V，BTL输出。前端带有音量、左右平衡、高音、低音调解旋钮。

实验中只是用其中一路来驱动滑轨上的扬声器。

驱动声源扬声器的音频功率放大器

4. 数据发送和接收

下图是发送和接收声音数据的带有单片机的实验面包板。它的原理在如下相关博文中进行了介绍：

基于STM32F103F ADDA板制作
声音定位硬件电路

发送和接收声音的数据采集单片机实验面包板

二、实验数据采集

1. 实验的方案

设置麦克阵列在相距滑轨不同的距离下，控制扬声器从左移动到右，每个位置下发出Chirp声音。接收并存储该声音数据。

麦克与滑轨之间的距离

1. Chirp声音参数：频率250~2000Hz;时间长度0.2048秒；采样频率10kHz；
2. 滑轨移动参数：方向是从远端移动到近端。移动100次，控制指令中移动步长为115。
3. 麦克风与滑轨距离：距离分别设置在90厘米，50厘米，25厘米。总共采集三组数据。
4. 采集数据格式每个位置数据包括有100组数据，每组数据为(d1,d2, sd, rd1, rd2)
其中：d1,d2：是按照27.7摄氏度下声速计算出来的左右两个声道距离声源的距离；
sd:发送信号的原始数据； rd1,rd2：是接受到的左右两个声道的原始声音数据。数据的长度为2048。

d1,d2数据处理算法代码：

AMBIENT_TEMPERATURE = 27.7def soundspeed(temperature = 25):return 331.4+0.6*temperature
def procdatafftinterpolation(send, rece, interptime):sm = mean(send)rm = mean(rece)send = [d - sm for d in send]rece = [d - rm for d in rece]datalen = len(send)sendspace = zeros(datalen*2)sendspace[0:datalen] = sendrecespace = zeros(datalen*2)recespace[0:datalen] = recetempfft = fft.fft(sendspace) * conj(fft.fft(recespace))expandfft = zeros(datalen*2*interptime, dtype=complex64)expandfft[0:datalen] = tempfft[0:datalen]expandfft[-datalen:] = tempfft[datalen:]cor = list(fft.ifft(expandfft).real)maxpos = cor.index(max(cor))return (4096 * interptime - maxpos) / interptime * soundspeed(AMBIENT_TEMPERATURE) * 0.1e-3

2. 实验数据

每个位置的数据中，d1,d2数量比较小，下面分别给出了三个位置下，d1,d2的波形和数据。原始的sd,rd1,rd2由于维度（2048维）比较大，可以在下载资源里进行下载。

1. 在一米处的采集数据

在1米处采集到的D1,D2数据
在实际处理中，估计相关位置进行插值，插值倍数取20。

d1=[1.2946344, 1.2946344, 1.2894141, 1.2859339, 1.2824537, 1.2789735, 1.2754933, 1.2720131, 1.2685329, 1.2650527, 1.2615725, 1.2580923, 1.2546121, 1.2511319, 1.2476517, 1.2459116, 1.2424314, 1.2406913, 1.2372111, 1.235471, 1.2319908, 1.2302507, 1.2267705, 1.2250304, 1.2232903, 1.2198101, 1.21807, 1.2163299, 1.2145898, 1.2128497, 1.2111096, 1.2093695, 1.2076294, 1.2058893, 1.2041492, 1.200669, 1.1989289, 1.1971888, 1.1954487, 1.1937086, 1.1902284, 1.1884883, 1.1867482, 1.1850081, 1.183268, 1.183268, 1.1815279, 1.1797878, 1.1850081, 1.1763076, 1.1745675, 1.1745675, 1.1728274, 1.1728274, 1.1710873, 1.1710873, 1.1693472, 1.1693472, 1.1693472, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.165867, 1.165867, 1.165867, 1.165867, 1.165867, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1676071, 1.1693472, 1.1693472, 1.1710873, 1.1710873, 1.1728274, 1.1728274, 1.1745675, 1.1763076, 1.1780477, 1.1797878, 1.183268, 1.1867482, 1.1902284, 1.1919685]
d2=[1.1971888, 1.1989289, 1.1954487, 1.1937086, 1.1919685, 1.1902284, 1.1884883, 1.1867482, 1.1850081, 1.1815279, 1.1797878, 1.1797878, 1.1780477, 1.1763076, 1.1745675, 1.1728274, 1.1728274, 1.1728274, 1.1710873, 1.1693472, 1.1693472, 1.1693472, 1.1676071, 1.165867, 1.1676071, 1.165867, 1.165867, 1.1641269, 1.1641269, 1.1641269, 1.1641269, 1.1623868, 1.1623868, 1.1623868, 1.1623868, 1.1623868, 1.1623868, 1.1623868, 1.1623868, 1.1623868, 1.1641269, 1.1641269, 1.1641269, 1.165867, 1.165867, 1.1676071, 1.1676071, 1.1693472, 1.1763076, 1.1710873, 1.1710873, 1.1745675, 1.1763076, 1.1780477, 1.1797878, 1.1815279, 1.183268, 1.1850081, 1.1867482, 1.1884883, 1.1902284, 1.1902284, 1.1937086, 1.1937086, 1.1954487, 1.1971888, 1.1989289, 1.200669, 1.200669, 1.2024091, 1.2041492, 1.2058893, 1.2093695, 1.2111096, 1.2128497, 1.2145898, 1.21807, 1.21807, 1.2215502, 1.2232903, 1.2250304, 1.2267705, 1.2285106, 1.2302507, 1.2302507, 1.2319908, 1.235471, 1.2372111, 1.2406913, 1.2441715, 1.2476517, 1.2511319, 1.2546121, 1.2580923, 1.2598324, 1.2633126, 1.2633126, 1.2667928, 1.270273, 1.2737532]

2. 在50厘米处采集的数据

在50厘米处采集到的D1,D2数据

d1=[0.9883767999999998, 0.9883767999999998, 0.9831564999999999, 0.9779361999999999, 0.9709758000000001, 0.9657555000000001, 0.9622753, 0.957055, 0.9518347, 0.9466143999999999, 0.9413940999999999, 0.9361738, 0.9309535, 0.9257332, 0.9222529999999999, 0.9187728, 0.9135525, 0.9083322, 0.9048520000000001, 0.8996317000000001, 0.8961515, 0.8926713000000001, 0.8891911, 0.8857108999999999, 0.8804905999999999, 0.8770104, 0.8735302, 0.87005, 0.8665698, 0.8648297, 0.8613495, 0.8596094, 0.8561292, 0.8543891, 0.8509089000000001, 0.8474287, 0.8439484999999999, 0.8422084, 0.8387281999999999, 0.835248, 0.8335079, 0.8317678000000001, 0.8300277, 0.8282876, 0.8265475000000001, 0.8230673, 0.8213272000000001, 0.8195871, 0.817847, 0.8143668, 0.8126267, 0.8108866, 0.8091465, 0.8074064, 0.8056663, 0.8039262, 0.8021861, 0.800446, 0.800446, 0.7987059, 0.7987059, 0.7969658, 0.7969658, 0.7969658, 0.7969658, 0.7969658, 0.7969658, 0.7969658, 0.7969658, 0.7969658, 0.7987059, 0.7987059, 0.7987059, 0.7987059, 0.7987059, 0.800446, 0.7987059, 0.800446, 0.8021861, 0.8021861, 0.8039262, 0.8039262, 0.8056663, 0.8074064, 0.8108866, 0.8126267, 0.8143668, 0.8161069, 0.817847, 0.8195871, 0.8213272000000001, 0.8230673, 0.8248073999999999, 0.8265475000000001, 0.8282876, 0.8300277, 0.8335079, 0.835248, 0.8404683, 0.8422084]
d2=[0.835248, 0.8335079, 0.8317678000000001, 0.8282876, 0.8265475000000001, 0.8248073999999999, 0.8230673, 0.8195871, 0.8161069, 0.8126267, 0.8091465, 0.8074064, 0.8039262, 0.800446, 0.7987059, 0.7969658, 0.7952257, 0.7934856, 0.7917455, 0.7900054, 0.7882653, 0.7865252, 0.7847851, 0.783045, 0.7813049, 0.7813049, 0.7795648, 0.7795648, 0.7795648, 0.7778247, 0.7778247, 0.7778247, 0.7778247, 0.7760846, 0.7778247, 0.7778247, 0.7778247, 0.7795648, 0.7795648, 0.7813049, 0.7813049, 0.7813049, 0.7813049, 0.7813049, 0.783045, 0.7847851, 0.7847851, 0.7865252, 0.7882653, 0.7900054, 0.7934856, 0.7952257, 0.7987059, 0.8021861, 0.8056663, 0.8108866, 0.8143668, 0.8161069, 0.8195871, 0.8213272000000001, 0.8248073999999999, 0.8265475000000001, 0.8282876, 0.8317678000000001, 0.835248, 0.8387281999999999, 0.8422084, 0.8456886000000001, 0.8491687999999998, 0.8509089000000001, 0.852649, 0.8561292, 0.8596094, 0.8613495, 0.8665698, 0.87005, 0.8735302, 0.8770104, 0.8804905999999999, 0.8839708, 0.8874510000000001, 0.8909312, 0.8944114, 0.8996317000000001, 0.9048520000000001, 0.9083322, 0.9118124000000001, 0.9152926, 0.9187728, 0.9222529999999999, 0.9274732999999998, 0.9344337, 0.939654, 0.9466143999999999, 0.9518347, 0.957055, 0.9605352000000001, 0.9657555000000001, 0.9727159, 0.9761961]

3. 在25厘米处采集的数据

在25厘米处接收到的D1,D2的数据

d1=[0.8543891, 0.852649, 0.8474287, 0.8422084, 0.8369881, 0.8335079, 0.8265475000000001, 0.8195871, 0.8126267, 0.8074064, 0.8021861, 0.7952257, 0.7865252, 0.7778247, 0.7708643, 0.7673841, 0.7621638, 0.7552034, 0.7447628, 0.7360623, 0.730842, 0.7291019, 0.7256217, 0.7186613, 0.7082207, 0.6995202, 0.6925598, 0.6890796, 0.6855994, 0.6803791, 0.6734187, 0.6629781, 0.6542776, 0.6490573, 0.643837, 0.6368766, 0.6299162, 0.6229558, 0.6142553, 0.6055548, 0.5985944, 0.591634, 0.5829335, 0.5759731, 0.5707528, 0.5655325, 0.5620523, 0.556832, 0.5533518000000001, 0.5481315, 0.5411711, 0.5359508000000001, 0.5324706, 0.5289904, 0.5255102, 0.5237701, 0.5202899, 0.5185498, 0.5150696, 0.5115894, 0.5081092, 0.5063691, 0.5028889, 0.5011488000000001, 0.4994087, 0.4959285, 0.4959285, 0.4941883999999999, 0.49244830000000006, 0.49244830000000006, 0.4907082, 0.4907082, 0.4907082, 0.4907082, 0.4907082, 0.4907082, 0.4907082, 0.4907082, 0.49244830000000006, 0.49244830000000006, 0.4941883999999999, 0.4959285, 0.4994087, 0.5011488000000001, 0.504629, 0.5133295, 0.5115894, 0.5150696, 0.5168097, 0.5185498, 0.5202899, 0.52203, 0.5255102, 0.5289904, 0.5342107, 0.539431, 0.5429111999999999, 0.5481315, 0.5516116999999999, 0.5533518000000001]
d2=[0.6194756, 0.6194756, 0.6142553, 0.609035, 0.6038147, 0.5968543, 0.5898939, 0.5811934, 0.574233, 0.5690127, 0.5655325, 0.5603122, 0.5585721, 0.5550919, 0.5498716, 0.5446513000000001, 0.539431, 0.5342107, 0.5307305, 0.5289904, 0.5255102, 0.5237701, 0.5202899, 0.5185498, 0.5168097, 0.5150696, 0.5115894, 0.5098493, 0.5063691, 0.504629, 0.5028889, 0.5011488000000001, 0.4994087, 0.4976686, 0.4959285, 0.4959285, 0.4941883999999999, 0.4941883999999999, 0.4941883999999999, 0.4941883999999999, 0.4941883999999999, 0.4959285, 0.4959285, 0.4976686, 0.4994087, 0.5011488000000001, 0.5028889, 0.5028889, 0.5063691, 0.5081092, 0.5098493, 0.5115894, 0.5150696, 0.5185498, 0.52203, 0.5255102, 0.5307305, 0.5342107, 0.539431, 0.5429111999999999, 0.5463914, 0.5498716, 0.5550919, 0.5603122, 0.5655325, 0.5690127, 0.574233, 0.5794533, 0.5846736, 0.5898939, 0.5951142, 0.6020746, 0.609035, 0.6125152, 0.6212157, 0.6281761, 0.6333964, 0.6386167, 0.6455771, 0.6507974, 0.6594979, 0.6681984, 0.6734187, 0.678639, 0.6838593, 0.69604, 0.6995202, 0.7047405, 0.7117009, 0.7169212, 0.7221415, 0.7291019, 0.7360623, 0.7465029, 0.7534633, 0.7586836, 0.7639039, 0.7708643, 0.7795648, 0.7865252]

三、数据分析

1. 实验数据理论模型

根据勾股定理，声源在滑轨上的移动距离和接收声音的麦克风之间的距离为：L1=(x1−x0)2+y12L_1 = \sqrt {\left( {x_1 - x_0 } \right)^2 + y_1^2 }L1=(x1−x0)2+y12

滑块移动距离与声源和接收传感器之间的距离关系

下图就是在x0=50, y1=25, x1在[0,100]范围内，L1与x1之间的关系：

横向运动距离与直线距离之间的关系曲线

2. 实验数据参数估计方法

1. 分析模型已知和未知参数

从前面实验所得到的不同位置的d1,d2数据曲线来看，它们都大体呈现与上面理论模型相近的性质。只是在实验中粗略的知道y1,x0,x1的数据，实际上测量并不精确。再加上实验中声速也是一个变化量，所以直接使用原始数据来验证它是否符合理论模型就会有困难。

下面先假设数据符合上面的理论模型，然后在对未知的x0, y1， x1进行估计，然后再验证模型的精确性。
在这里先假设以下两点是准确的：
（1）滑轨是直线的；滑轨移动距离ΔD\Delta DΔD是可以实现测量得到的，它等于：ΔD=0.79m\Delta D = 0.79mΔD=0.79m。
（2）每个数据点都是在起始点xstartx_{start}xstart和终止点xstopx_{stop}xstop之间的均匀采样。假设起始点xstart=0x_{start} = 0xstart=0，终止点xstop=ΔDx_{stop} = \Delta Dxstop=ΔD。
（3）采集数据的个数N=100，数据为：{Ln},n=0,1,...,N−1\left\{ {L_n } \right\},n = 0,1,...,N - 1{Ln},n=0,1,...,N−1
这个数据时通过声音传播速度计算得到的，它与真实距离之间相差一个比例aaa，即：dn=a⋅Lnd_n = a \cdot L_ndn=a⋅Ln

然后通过采集所得到的数据{Ln}\left\{ {L_n } \right\}{Ln}拟合下面的公式：

x1n=xstop−xstartN−1⋅n+xstart=ΔDN−1⋅nx_{1n} = {{x_{stop} - x_{start} } \over {N - 1}} \cdot n + x_{start} = {{\Delta D} \over {N - 1}} \cdot nx1n=N−1xstop−xstart⋅n+xstart=N−1ΔD⋅n

Lna=(x1n−x0)2+y12{{L_n } \over a} = \sqrt {\left( {x_{1n} - x_0 } \right)^2 + y_1^2 }aLn=(x1n−x0)2+y12

未知参数包括有：{a,x0,y1}\left\{ {a,x_0 ,y_1 } \right\}{a,x0,y1}。

2. 模型优化算法

模型优化算法采用python中的scipy.optimize优化工具完成模型参数估计。下面给出了相应的参考代码。

#------------------------------------------------------------
N = 100
xstart = 0
xstop = 0.7                    # as delta D
def distfunc(x, a, x0, y1):x1n = (xstop-xstart) * x / (N-1) + xstartreturn sqrt((x1n-x0)**2 + y1**2) / a
#------------------------------------------------------------
d1, d2 = dataall2lr(dataall2) # Read the L from the data file
x = linspace(0, N, N, endpoint=False)
fit_params, pcov = scipy.optimize.curve_fit(distfunc, x, d1)

3. 参数估计结果分析

1. 在1米处的数据进行参数估计
使用左右两个声道的距离完成估计的参数分别是：

左声道：a=0.8857, x0=0.5006, y1=1.0319
右声道：a=0.8958,x0=0.2341, y1=1.0445

两个声道的数据是同时测量的，测量的条件相似，它们的距离比例因子a的数值基本相同，大约在0.9左右。

两个声道的x0相差大约为0.266米，这与直接使用钢尺测量两个麦克风的中心位置距离0.265米在误差范围1mm内是相同的；

由于两个麦克风摆放位置是平行与滑轨，所以它们距离滑轨的距离y0应该是相同的，它们的结果反映了这一点。由于实际测两个麦克风与滑轨的距离为0.9日左右，所以将y1乘以距离因子a，则等于：D=a×y1=0.8857×1.0319=0.913D = a \times y_1 = 0.8857 \times 1.0319 = 0.913D=a×y1=0.8857×1.0319=0.913。这就与实际距离相吻合了。

两个声道数据和参数估计后拟合数据

2. 在0.5米处的数据处理
使用左右两个声道的距离完成估计的参数分别是：

左声道：a=0.8145, x0=0.4787, y1=0.6492
右声道：a=0.7907,x0=0.2326, y1=0.6177

左右两个声道原始数据和参数估计后的拟合数据

3. 在0.25米处的数据处理
使用左右两个声道的距离完成估计的参数分别是：

左声道：a=0.7180, x0=0.5241, y1=0.354
右声道：a=0.7039,x0=0.2644, y1=0.349

左右两个声道原始数据和参数估计后的拟合数据

4. 声音测距因子

（1）实验结果中声音测距因子存在的问题

前面使用了基于理论模型参数估计的方法，从直线滑轨采集的数据中获得了实际滑轨和拾音麦克之间的相对位置（x0,y1），但是有一个参数：声音测距因子 a 在每组数据都各不相同。a 反映了通过声音延迟所获得的距离与实际距离之间的差别比例。

如果这个因子仅仅是由于利用空气声音传播速度、采样时间延迟等因素所造成的，那么因子a的数值有两点问题：
（1）因子a的数值不会这么大。空气声音传播速度会因为温度不同而变化，但是在室温条件下，这个变化范围不会超过5%。按照实际测量的温度相差十度来估算的误差范围；考虑到采样系统中AD，DA延迟，它们都在十几个微秒的范围，所造成的距离误差在几个毫米。
（2）在三个不同的位置所测量得到的距离因子 a 它们应该相同，但实际上，三个位置上的 a各不相同。

分析造成声音测距因子的来源应该是地面反射声波的影响。

（2）考虑地面反射波的影响

假设麦克风距离地面高度为H1H_1H1，扬声器距离地面的高度为H2H_2H2，它们之间的水平距离为L1L_1L1。

麦克与扬声器之间的直线距离：L2=L12+(H1−H2)2L_2 = \sqrt {L_1^2 + \left( {H_1 - H_2 } \right)^2 }L2=L12+(H1−H2)2

麦克与扬声器之间的位置关系

地面反射波在地面的反射点距离麦克的距离为m1m_1m1，距离扬声器的距离为m2m_2m2。那么：m1+m2=L1m_1 + m_2 = L_1m1+m2=L1
再根据反射波的入射角与出射角相同，则有：m1H1=m2H2{{m_1 } \over {H_1 }} = {{m_2 } \over {H_2 }}H1m1=H2m2

那么可以得到：m1=H1H1+H2⋅L1,m2=H2H1+H2⋅L1m_1 = {{H_1 } \over {H_1 + H_2 }} \cdot L_1 ,\,\,\,\,m_2 = {{H_2 } \over {H_1 + H_2 }} \cdot L_1m1=H1+H2H1⋅L1,m2=H1+H2H2⋅L1
反射距离L3L_3L3为：

L3=H12+(H1H1+H2L1)2+H22+(H2H1+H2L1)2L_3 = \sqrt {H_1^2 + \left( {{{H_1 } \over {H_1 + H_2 }}L_1 } \right)^2 } + \sqrt {H_2^2 + \left( {{{H_2 } \over {H_1 + H_2 }}L_1 } \right)^2 }L3=H12+(H1+H2H1L1)2+H22+(H1+H2H2L1)2

猜测实际声波传送的距离应该是L2L_2L2和L3L_3L3的某种加权平均。

根据前面实验中的数据：H1H_1H1=0.15m，H2H_2H2=0.12m。下面分别列些出L1L_1L1在0.9m，0.5m，0.25m的情况下，L2L_2L2，L3L_3L3的取值：

序号	L1L_1L1	L2L_2L2	L3L_3L3	L2,L3平均值
1	0.9	0.9005	0.9396	0.92
2	0.5	0.5009	0.5682	0.54
3	0.25	0.252	0.3680	0.31

从上面表格中可以看到在相距0.9,0.5,0.25米的距离上，直接传播距离L2，反射波传播距离L3都比L1大。在这里取L1与L3的比值来计算一下距离因子 a，可以得到三个距离下的距离因子a:

序号	L1L_1L1	L3L_3L3	a=L1L_1L1/L3L_3L3	数据拟合中的a
1	0.9	0.9396	0.9578	0.9
2	0.5	0.5682	0.88	0.8
3	0.25	0.368	0.679	0.7

在L1L_1L1等于0.25米是，对应的计算因子与实际数据拟合中的a相接近，但其它位置的a都相差较远。

具体的差别现在还需要进一步分析和探索。