2

I am currently attempting to solve a system of quadratic (and linear) systems that I have run into while trying to triangulate sound.

My hypothetical setup includes 3 sensors on a perfectly equilateral triangle, with one sensor located at $(0,0)$ and the other two located below it. (The specifics don't matter, as I am simply referring to the sensor locations using $a_1,a_2,a_3$ for the x-coordinates of the sensors, and $b_1,b_2,b_3$ for the y-coordinates of the sensors, with $r_1,r_2,r_3$ being the radii of the circles from each respective sensor to the sound point)

I am trying to specify equations for the x position of the sound, the y position of the sound, and finally the radius of the incident sensor to the sound (the sensor that picks up the sound wave first).

My equations are as follows:

$$(x - a_1)^2 + (y - b_1)^2 = r_1^2$$ $$(x - a_2)^2 + (y - b_2)^2 = r_2^2$$ $$(x - a_3)^2 + (y - b_3)^2 = r_3^2$$ $$r_3 = r_1 + (t_3 * \text{speed of sound})$$ $$r_2 = r_1 + (t_2 *\text{speed of sound})$$

In this example, I am assuming that the sound reaches sensor 1 first. I understand that a true solution requires 3 discrete solutions, one for each sensor being the "incident sensor". (assuming that there cannot be a scenario where sound perfectly reaches multiple sensors at the same time)

My known variables: $a_1,a_2,a_3,b_1,b_2,b_3,\text{speed of sound}, t_1,t_2,t_3$

My Unknown variables: $x,y,r_1,r_2,r_3$.

Now I understand that I can just substitute in the three linear equations, but that leaves me with three quadratic equations that I am unsure of how to solve and obtain a meaningful answer from.

I tried searching for revelant topics, and the closest I could come was this: https://math.stackexchange.com/a/187858/656339

Which has the same setup as I, but doesn't detail how to solve it.

Any help would be appreciated.

Rodrigo de Azevedo
  • 1
ReThink
  • 41
  • 4
  • What you’re trying to do is called “TDOA multilateration.” Technically speaking, triangulation involves angles, not arrival time differences. – amd Mar 21 '19 at 07:53
  • You can convert the system of circle equations into a system of linear equations: the intersections of a pair of circles lie on their radical axis, an equation for which you can obtain by subtracting one circle equation from another. However, since you don’t know $r_1$, what you’re really going to end up having to do is compute intersections of hyperbolas. – amd Mar 21 '19 at 07:55
  • The problem would be easier with a fourth sensor. – Claude Leibovici Mar 21 '19 at 09:16
  • @amd could you expound a little on how to go about doing that? Converting to a system of linear equations would be preferred – ReThink Mar 21 '19 at 13:16
  • Have a look at https://math.stackexchange.com/questions/2462039/triangulation-math/2462657#2462657 – Claude Leibovici Mar 21 '19 at 14:46
  • @ClaudeLeibovici I noticed you were the one to give an answer there as well, I was wondering if you could explain the method used in your answer a little more (I have not personally used linear regression before) – ReThink Mar 21 '19 at 18:13
  • @ReThink. You can reduce the problem to simple linear equations provided a fourth sensor to get $x,y,t_0$. In such a case, no regression at all for four sensors. It could be looking paradoxal but, with time data, the problem of three sensors is quite ugly. If you want me to elaborate, just tell. – Claude Leibovici Mar 22 '19 at 04:07
  • @ClaudeLeibovici Yes please! Four sensors works just as well for my purposes, I believe. – ReThink Mar 24 '19 at 19:15
  • Do you see how simple the problem becomes with a fourth sensor ? Just linear equations ! – Claude Leibovici Mar 25 '19 at 06:56
  • @ClaudeLeibovici, could you please point out to me how to do the formulation for an arbitrary number of microphones >4 ? Would adding more microphones improve the position estimation? – Loza Feb 01 '23 at 12:47

2 Answers2

3

After comments, let us work with $4$ sonsors. So we have $$(x - a_1)^2 + (y - b_1)^2 = c^2(t_1-\tau)^2\tag 1$$ $$(x - a_2)^2 + (y - b_2)^2 = c^2(t_2-\tau)^2\tag 2$$ $$(x - a_3)^2 + (y - b_3)^2 = c^2(t_3-\tau)^2\tag 3$$ $$(x - a_4)^2 + (y - b_4)^2 = c^2(t_4-\tau)^2\tag 4$$ where $c$ is the speed of sound and $\tau$ the time at which was produced the sound.

Now, subtract $(1)$ from $(2)$, $(3)$ and $(4)$ to get $$2(a_1-a_2)x+2(b_1-b_2)y+2c^2(t_2-t_1)\tau=(a_1^2+b_1^2-c^2t_1^2)-(a_2^2+b_2^2-c^2t_2^2)\tag 5$$ $$2(a_1-a_3)x+2(b_1-b_3)y+2c^2(t_3-t_1)\tau=(a_1^2+b_1^2-c^2t_1^2)-(a_3^2+b_3^2-c^2t_3^2)\tag 6$$ $$2(a_1-a_4)x+2(b_1-b_4)y+2c^2(t_4-t_1)\tau=(a_1^2+b_1^2-c^2t_1^2)-(a_4^2+b_4^2-c^2t_4^2)\tag 7$$

Define, for more simplicity, $$\alpha_i=2(a_1-a_i)\qquad \beta_i=2(b_1-b_i)\qquad \gamma_i=2c^2(t_i-t_1)$$ $$k_i=(a_1^2+b_1^2-c^2t_1^2)-(a_i^2+b_i^2-c^2t_i^2)$$ $(i=2,3,4)$ to make $$\alpha_2x+\beta_2y+\gamma_2 \tau=k_2\tag 8$$ $$\alpha_3x+\beta_3y+\gamma_3 \tau=k_3\tag 9$$ $$\alpha_4x+\beta_4y+\gamma_4 \tau=k_4\tag {10}$$ So, three linear equations in $(x,y,\tau)$ (easy to solve - use matrix calculations or successive elimination as I did below) and the explicit solutions are

$$\tau=\frac{k_4 (\alpha_3 \beta_2-\alpha_2 \beta_3)+k_3 (\alpha_2 \beta_4-\alpha_4 \beta_2)+k_2 (\alpha_4 \beta_3-\alpha_3 \beta_4) } {\alpha_4 (\beta_3 \gamma_2-\beta_2 \gamma_3)+\alpha_3 (\beta_2 \gamma_ 4-\beta_4 \gamma_2)+\alpha_2 (\beta_4 \gamma_3-\beta_3 \gamma_4) }\tag {11}$$ $$y=\frac{\alpha_3 (\gamma_2 \color{red}{\tau} -k_2)+\alpha_2 (k_3-\gamma_3 \color{red}{\tau} )}{\alpha_ 2 \beta_3-\alpha_3 \beta_2}\tag {12}$$ $$x=\frac{k_2-\beta_2 \color{red}{y}-\gamma_2 \color{red}{\tau} }{\alpha_2}\tag {13}$$

Claude Leibovici
  • 260,315
  • To clarify, is tau the time elapsed from the sound occurring relative to the first sensor that detects the sound wave? – ReThink Mar 25 '19 at 13:29
  • @ReThink. Suppose that the first sensor receives the information at 1:34.05, the second at 1:34.06, the third at 1:34.08 and the fourth at 1:34.08 and you find $\tau=1:31.01$; this is the time at which the sound was produced. – Claude Leibovici Mar 25 '19 at 17:16
  • Thanks a ton! One more question: When you define alpha, beta, gamma, and k, you use subscript '1', does this mean that these only hold if sensor #1 is the first to receive the sound? i.e., would I have to rearrange and replace with a subscript of 2,3,4 if the respective sensor is the first to receive the sound wave?

    ex: αi=2(a2−ai) if sensor #2 is the incident sensor? Or does your solution hold regardless of which sensor is the incident sensor?

     – ReThink Mar 25 '19 at 18:12
  • @ReThink. There is no incident sensor. You have sensors located at $x_i,y_i$ which at time $t_i$ receive an information. If you want to use this concept (which is dangerous at least to me), start defining $t_{min}$ and redefine the $t_i$ as $t_i-t_{min}$. – Claude Leibovici Mar 26 '19 at 05:36
  • I ran everything in Matlab just to visualize what's going on, and I understand now. Tau can be either an absolute or a relative time (it does not matter to the sensors) and t1 - t4 are offsets from this reference point (tau). The incident sensor does not matter, as that information is a part of the arrival times of the sound waves to the respective sensors, which is a fundamental part of t1-t4. Thank you very much for taking the time to explain this concept! It's fascinating that the problem can be simplified by just adding one more sensor. – ReThink Mar 27 '19 at 02:12
  • @ReThink. Good to know ! Now, where the problem starts to be very interesting is when $n>4$ because the $t_i$'s are all in error. In such a case, you need to use what I wrote in the linked page. – Claude Leibovici Mar 27 '19 at 05:12
  • Is there a way to account for practical error in these equations? I’m currently testing this physical setup, however slight error in measuring the t1-t4 seems to throw off x and y a great deal – ReThink Apr 17 '19 at 21:56
  • essentially anything further than ~1m or so is off by a fair amount – ReThink Apr 17 '19 at 22:39
  • forgot to include the tag: @Claude Leibovici – ReThink Apr 18 '19 at 03:09
0

Solving system of quadratic and indeed, general polynomial equations is possible with techniques like Buchberger's algorithm. See the first two chapters of the book by Cox et.al. There are also many numerical tools to assist you with this. Sympy in python is one alternative (the one I'd recommend) and there is also one I wrote in C# following the textbook cited.

Rohit Pandey
  • 6,803