■ 音の出所はどこ？　[音源位置の同定] （JavaScript版）
　　　　　～　実際に連立方程式を解いて音源位置を求める　～

　「音の出所はどこ？」において、２次元平面上に設置された３個のマイクに音が到達する時間差から音源位置を推定できることが示されています。

　マイクの位置をM1、M2、M3とすると、音源位置はM1、M2を焦点とする双曲線1-2とM1、M3を焦点とする双曲線1-3の交点となります。　交点は最大８個ありますが、各マイクに音が到達する時刻の大小関係を考慮することでそれぞれの双曲線でどちらの曲線上の交点かが決まり、従ってその２曲線の交点として音源位置が通常１個求められます。
　しかし音源位置によってはその点を通る２つの双曲線（双曲線1-2と双曲線1-3）が２点で交わることもあります。

　次のアプリはマウスクリックされた音源位置を通るこれら２つの双曲線を描画します。　音源位置がマイク後方あるいはマイク近傍にあるときは、２つの双曲線が２点で交わることが確認できます。
　・図中の1, 2, 3（赤）がマイクです。
　・緑色の線はマイク 1、2 からの距離の差が一定の点の集まり（双曲線）です。

　また、ここでは音源(x0, y0)を未知として、ここから各マイクMi(xi, yi)に音が到達する時刻tiに関する後述の連立方程式を解いて音源位置を算出して表示します　－＞ ｸﾞﾗﾌｨｯｸ表示ｴﾘｱの左下隅に(x0, y0) = ・・・の形式で。
 　クリック位置X、Yが音源位置として求められることが分かります。
マイク位置

　・位置の単位はｍ、時刻、時間差は秒、図中のグリッド幅は１ｍです。
　・「双曲線 2-3 表示 On/Off」ボタンを押すと、音源位置を通り、M2、M3を焦点とする双曲線2-3の表示をOn/Offします。
　　双曲線1-2と双曲線1-3が２点で交わるとき、双曲線2-3もその２つの交点を通ります（理由は下記）。
　・「交点が2個のｴﾘｱ表示On/Off」ボタンを押すと、双曲線1-2と双曲線1-3が２点で交わる領域の表示をOn/Offします。

　教室内などで音源の位置をマイクを使用して推定する場合は、マイクを周囲の壁面などに設置することにより、マイク後方からの音の発生は除外することができます（候補点が２つある場合は室内の点を採用）。

　● 音源位置(x0, y0)に関する連立方程式とその解法
　各マイクの位置を(x1, y1), (x2, y2), (x3, y3)、音源位置を(x0, y0)、
　音の発生時刻をt0（未知）、各マイクに音が到達する時刻を t1, t2, t3 とする。

・マイク１，２に音が到達する時間差(t2-t1)から次式が成り立つ。
　　sqrt((x2-x0)^2+(y2-y0)^2) - sqrt((x1-x0)^2+(y1-y0)^2) = c(t2-t1)　…(式１)
　　ここで、sqrt : 平方根(ルート)
　　　　　　^2   : ２乗
　　　　　　c    : 音速（既知、気温により決まる）

　同様に、マイク２，３に音が到達する時間差(t3-t2)から次式が成り立つ。
　　sqrt((x3-x0)^2+(y3-y0)^2) - sqrt((x2-x0)^2+(y2-y0)^2) = c(t3-t2)　…(式２)

　しかし、ここでは音源位置を求めるために以下の手順を採用する。

・音源から各マイクまでの距離と音の到達時間より次式が成り立つ。
　　sqrt((x1-x0)^2 + (y1-y0)^2) = c(t1-t0)　…(式3.1)
　　sqrt((x2-x0)^2 + (y2-y0)^2) = c(t2-t0)　…(式3.2)
　　sqrt((x3-x0)^2 + (y3-y0)^2) = c(t3-t0)　…(式3.3)

・式3.1～式3.3を2乗すると、
　　(x1-x0)^2 + (y1-y0)^2 = [c(t1-t0)]^2　…(式4.1)
　　(x2-x0)^2 + (y2-y0)^2 = [c(t2-t0)]^2　…(式4.2)
　　(x3-x0)^2 + (y3-y0)^2 = [c(t3-t0)]^2　…(式4.3)

・式4.1から式4.2を引くと、
　　(x1-x0)^2 + (y1-y0)^2 - (x2-x0)^2 - (y2-y0)^2 = [c(t1-t0)]^2 - [c(t2-t0)]^2
　展開整理して、
　　(x1^2+y1^2-x2^2-y2^2)-2x0(x1-x2)-2y0(y1-y2) = c^2(t1^2-t2^2)-2c^2t0(t1-t2)　…(式5.1)

・同様に、式4.1から式4.3を引き、展開整理すると、
　　(x1^2+y1^2-x3^2-y3^2)-2x0(x1-x3)-2y0(y1-y3) = c^2(t1^2-t3^2)-2c^2t0(t1-t3)　…(式5.2)

・式5.1、式5.2を変形して、
　　-2(y1-y2)y0+2c^2(t1-t2)t0 = c^2(t1^2-t2^2)-(x1^2+y1^2-x2^2-y2^2)+2(x1-x2)x0　…(式5.1a)
　　-2(y1-y3)y0+2c^2(t1-t3)t0 = c^2(t1^2-t3^2)-(x1^2+y1^2-x3^2-y3^2)+2(x1-x3)x0　…(式5.2a)

・ここで、
　　a1=-2(y1-y2),  b1=2c^2(t1-t2),  c1=c^2(t1^2-t2^2)-(x1^2+y1^2-x2^2-y2^2),  d1=2(x1-x2)
　　a2=-2(y1-y3),  b2=2c^2(t1-t3),  c2=c^2(t1^2-t3^2)-(x1^2+y1^2-x3^2-y3^2),  d2=2(x1-x3)
　とおくと、式5.1a、式5.2aは
　　a1y0 + b1t0 = c1 + d1x0　…(式5.1b)
　　a2y0 + b2t0 = c2 + d2x0　…(式5.2b)
　となる。

・これをy0, t0に関する連立方程式と考えて解くと、
　　y0 = Ax0 + B　…(式6.1)
　　t0 = Cx0 + D　…(式6.2)
　　ここで、
　　　A = (b2d1-b1d2)/(a1b2-a2b1),   B = (b2c1-b1c2)/(a1b2-a2b1)
　　　C = (a1d2-a2d1)/(a1b2-a2b1),   D = (a1c2-a2c1)/(a1b2-a2b1)

・この y0、t0を式4.1に代入すると、x0に関する次の2次方程式を得る。
　　[x1-x0]^2 + [y1-(Ax0+B)]^2 = [c[t1-(Cx0+D)]]^2

　x0について整理すると、
　　Ex0^2 + Fx0 + G = 0
　　ここで、
　　　E = 1 + A^2 - c^2C^2　　　　　　　　（小文字cと大文字Cに注意）
　　　F = -2x1 - 2(y1-B)A + 2c^2(t1-D)C
　　　G = x1^2 + (y1-B)^2 - c^2(t1-D)^2

・この2次方程式を解くと音の発生位置x0が得られ、更に式6.1よりy0が求められる。
　　x0 = [-F ± sqrt(F^2 - 4EG)]/(2E)　　（但し、F^2 - 4EG ≧ 0 のとき）

（注）x0に関する2次方程式からは通常解が２個得られるが、このx0を式6.2に代入して
　　　得られるt0について、t0 < t1, t2, t3 を満たさないものは除外する。
音の出所はどこ？
音の出所はどこ？　[音源方向の推定と両耳間時間差]
ホーム