Plaatsbepaling met signaalbronnen

[javi]

Ik ben op zoek naar een formule voor de volgende situatie. Het gaat om een plaatsbepaling van een punt ten opzichte van drie audio signaalbronnen waarvan alleen de relatieve onderlinge intensiteit bepaald kan worden, en niet de richting. Dit is terug te brengen tot het volgende goniometrische probleem. Ik heb verschillende oplossingsrichtingen geprobeerd, maar ben te weinig thuis in goniometrie om er uit te komen.

Van een willekeurige driehoek ABC weet ik van ieder punt A, B en C de (x,y) coordinaten. Ergens op het vlak ligt punt P. Van de lijnstukken AP, BP en CP weet ik de onderlinge lengte verhoudingen. Dus AP:BP:CP is bekend, maar niet de daadwerkelijke lengte van de lijnstukken, en ook niet de hoeken die de lijnstukken maken ten opzichte van de zijden van de driehoek. Twee vragen: 1) zijn op basis van deze gegevens de (x,y) coordinaten van P altijd te berekenen? Mijn intuitie zegt van wel, maar ik kan me bijvoorbeeld voorstellen dat dat alleen eenduidig kan wanneer P binnen de driehoek ligt. 2) Zo ja, wat is dan de formule? Hoe druk ik P uit als functie van de coordinaten van de punten A, B en C en de verhouding tussen de lijnstukken AP, BP en CP?



Alvast bedankt voor de hulp!

[David]

Stel dat je de coordinaten van P, zeg (xp, yp) kent. Dan kan je met de stelling van pythagoras de coordinaten uitdrukken in de verhouding. Met de stelling van Pythagoras kan je ook de afstanden tussen de bekende punten bepalen.

[javi]

Dank voor de razendsnelle reactie! Uiteraard heb ik geprobeerd met pythagoras hier uit te komen. Daarvoor heb ik wel driehoeken nodig met een rechte hoek, en dat is juist iets waar ik niet goed uit kwam. Welke hulppunten zou ik moeten definieren om geschikte rechthoekige driehoeken te krijgen, wat zou een geschikte strategie zijn daarvoor?

[David]

Dit is wat je kan doen. Hier zijn twee punten, A en B, met bekende coördinaten. Ik heb een voorbeeld gekozen zodat het wellicht intuïtiever wordt. Je kan nu een hulppunt a kiezen met een x-coordinaat van een punt en het y - coordinaat van het andere punt. Dan krijg je een rechthoekige driehoek. (zo kan je twee hulppunten H kiezen). Vervolgens kan je afstand AH bepalen en afstand BH en zo AB bepalen m.b.v. de stelling van Pythagoras. Lukt het zo?

[javi]

Dank je! Op deze wijze met Pythagoras de lengtes van AB, BC en CA bepalen is voor mij geen probleem. Het probleem is echter: hoe kom ik vanuit deze gegeven waardes (bekend zijn dus: de xy coordinaten van A, B, C, daarmee inderdaad de lengtes van lijnstukken AB, BC en CA, en de verhoudingen tussen de lengtes van lijnstukken AP, BP en CP) tot de xy coordinaat van P?
Dus bijvoorbeeld: A(0,0), B(4,0), C(1,3) en AP:BP:CP=2:3:5 Hoe bereken ik dan de coordinaten van P?

[arie]

In feite is de sleutel tot de oplossing geleverd door David.
De rest is een kwestie van stug doorrekenen.

In het algemeen: neem
A = (ax, ay)
B = (bx, by)
C = (cx, cy)
P = (x, y)
en
|PA| : |PB| : |PC| = lva : lvb : lvc
waarbij lva, lvb en lvc de lengteverhoudingen van P tot resp A, B en C zijn.

Dan geldt:
lvb * |PA| = lva * |PB|
ofwel
lvb^2 * |PA|^2 = lva^2 * |PB|^2
Gebruik hier de stelling van Pythagoras:
lvb^2 * [ (x - ax)^2 + (y - ay)^2 ] = lva^2 * [ (x - bx)^2 + (y - by)^2 ]
werk dit netjes uit naar de vorm:
p1*x^2 + q1*y^2 + r1*x + s1*y + t1 = 0
(merk op dat altijd p1=q1)
Ik kom dan uit op:
p1=lvb^2 - lva^2;
q1=lvb^2 - lva^2;
r1=2*bx*lva^2 - 2*ax*lvb^2;
s1=2*by*lva^2 - 2*ay*lvb^2;
t1=(ax^2+ay^2)*lvb^2 - (bx^2+by^2)*lva^2;

Doe ditzelfde voor P, A en C:
lvc * |PA| = lva * |PC|
ofwel
lvc^2 * |PA|^2 = lva^2 * |PC|^2
en opnieuw de stelling van Pythagoras:
lvc^2 * [ (x - ax)^2 + (y - ay)^2 ] = lva^2 * [ (x - cx)^2 + (y - cy)^2 ]
werk dit uit naar de vorm:
p2*x^2 + q2*y^2 + r2*x + s2*y + t2 = 0
(merk op dat ook hier weer altijd p2=q2)
Ik kom dan uit op:
p2=lvc^2 - lva^2;
q2=lvc^2 - lva^2;
r2=2*cx*lva^2 - 2*ax*lvc^2;
s2=2*cy*lva^2 - 2*ay*lvc^2;
t2=(ax^2+ay^2)*lvc^2 - (cx^2+cy^2)*lva^2;

Vermenigvuldig dan je eerste vergelijking met p2 en de tweede met p1:
p2*p1*x^2 + p2*p1*y^2 + p2*r1*x + p2*s1*y + p2*t1 = 0
p2*p1*x^2 + p2*p1*y^2 + p1*r2*x + p1*s2*y + p2*t2 = 0
trek deze 2 van elkaar af en je krijgt een vergelijking in de vorm:
lx*x + ly*y + lc = 0
waarbij
lx = p2*r1 - p1*r2;
ly = p2*s1 - p1*s2;
lc = p2*t1 - p1*t2;

werk dit om naar
y = la*x + lb
dan is
la = -lx/ly;
lb = -lc/ly;
(noot: in het zeldzame geval dat ly=0 weet je x, en kan je vervolgens y uit 1 van bovenstaande kwadratische vergelijkingen bepalen).

Substitueer deze vergelijking (y = la*x + lb) in
p1*x^2 + q1*y^2 + r1*x + s1*y + t1 = 0
en je krijgt een tweedegraadsfunctie in x:
a*x^2 + b*x + c = 0
met
a = p1 + q1*la^2;
b = 2*la*lb*q1 + s1*la + r1;
c = q1*lb^2 + s1*lb + t1;

en los x hieruit op met de abc-formule.

Als je 2 oplossingen vindt, x1 en x2, zijn de bijbehorende waarden voor y
te vinden via y = la*x + lb, dus:
y1 = la*x1 + lb
y2 = la*x2 + lb

Bovenstaande formules kan je in veel programmeertalen vrijwel direct verwerken.


Voorbeeld:
A = (2, 1)
B = (8, 2)
C = (4, 5)

en |PA| : |PB| : |PC| = 15.2 : 13.3 : 10.0
dus
lva=15.2;
lvb=13.3;
lvc=10.0;

levert 2 mogelijkheden voor P:
P ~= (6.3143, 13.1798)
of
P ~= (5.0157, 2.8902)

[javi]

Zeer veel dank Arie, ik heb het probleem bijna rond! Ik heb het geheel vandaag in een javascript gezet, dat werkt nu in alle gevallen waarbij de verhouding tussen AP en BP niet 1:1 is. Is die verhouding wel 1:1, dan wordt ergens in het script door 0 gedeeld en slaat het script 'op hol'. Ik vermoed dat dit te maken heeft met de volgende passage in de oplossing:

(noot: in het zeldzame geval dat ly=0 weet je x, en kan je vervolgens y uit 1 van bovenstaande kwadratische vergelijkingen bepalen).

Ik heb dat geprobeerd op te maken uit

lx*x + ly*y + lc = 0

want als ly=0 dan volgt x=-lc/lx, en y is in dat geval inderdaad te bepalen met behulp van een van de eerdere kwadratische vergelijkingen.

Maar in het geval de verhouding AP:BP 1:1 is, dan is ook lx gelijk aan 0. Daarmee deel ik nog steeds door 0 om x te krijgen en loopt mijn javascript vast.

Voorlopige oplossing in mijn script is de verhouding AP:BP checken op 1:1, en dan een heel klein verschil maken (door iets als 0.0000001 bij AP op te tellen ). Dat werkt uiteraard, de P coordinaten worden daarmee in ieder geval benaderd, maar het is wiskundig natuurlijk nogal twijfelachtig!

De vraag die nu nog open staat is: hoe reken ik in dit specifieke geval (wanneer ly=0) x uit? Zie ik iets over het hoofd?

[arie]

Normaliter zoek je de 0, 1 of 2 snijpunten van de krommen:
kromme K1 = p1*x^2 + q1*y^2 + r1*x + s1*y + t1 = 0
en
kromme K2 = p2*x^2 + q2*y^2 + r2*x + s2*y + t2 = 0

waarbij
p1 = q1 = lvb^2 - lva^2
en
p2 = q2 = lvc^2 - lva^2

Maar:

[1] als lvb = lva dan is p1 = q1 = 0 en is K1 een rechte lijn die je met K2 moet snijden (om precies te zijn: nu is K1 = de middelloodlijn van lijnstuk AB).
K1 is in dit geval r1*x + s1*y + t1 = 0, die je kan herschrijven als ofwel y = a*x + b ofwel x = c (indien s1 = 0), substitueer dat in K2 en los dan x of y op.
(noot: in dit geval vermenigvuldig je in de algemene oplossing alle parameters van K2 met nul, waardoor verderop in de abc-formule a = 0 wordt; je tweede "deling door nul"-error komt door delen door a=0 in de abc-formule).

[2] als lvc = lva dan is p2 = q2 = 0 en is K2 een rechte lijn die je met K1 moet snijden, analoog aan situatie [1]

[3] als lva = lvb = lvc dan is P het driehoekscentrum van driehoek ABC = het snijpunt van de middelloodlijnen van de zijden van ABC = het middelpunt van de omschreven cirkel.
Dit zou je ook via de methode van situatie [1] moeten kunnen oplossen

Lukt het je om deze situaties af te vangen en op te lossen?


PS: moet je ook nog controleren of de 3 punten A, B en C echt alle drie verschillend zijn en bovendien niet alle 3 op een rechte lijn liggen (en dan dus geen driehoek vormen)?

[javi]

Ik heb nu alles werkend in het javascript, echter wel met een kleine 'cheat' om sommige instanties van lva=lvb af te vangen.

Bovenstaande situatie [3] (lva=lvb=lvc) los ik nu eenvoudig op door de twee lijnen als rechte lijnen te benaderen, er is tenslotte maar 1 oplossing mogelijk in dat geval. Daarbij check ik nog of AB parallel loopt aan de x-as, want dan loopt de middelloodlijn parallel aan de y-as en moet ik uitgaan van x=c, of niet, in dat geval ga ik uit van y=ax+b

De situatie lva=lvb, maar beide niet gelijk aan lvc, heb ik eerst netjes geprobeerd op te lossen, inclusief de verschillende oplossingsrichtingen bij lijnstuk AB parallel aan de x-as of niet. Dat lukte niet, ik vermoed (achteraf) doordat ik soms per ongeluk testcoordinaten en testverhoudingen heb gebruikt die geen oplossingsmogelijkheid hebben en daarmee voortdurend letterlijk 0 op het rekest kreeg, dat gaf verwarring. De uiteindelijke 'cheat' daarvoor werkt als volgt: ik heb in het begin van het javascript een check op deze situatie gemaakt (lva=lvb && lva!=lvc). Als die zich voordoet, schuif ik simpelweg de coordinaten door: punt A krijgt de coordinaten van B, B van C, en C van A. Daarmee is automatisch lva niet meer gelijk aan lvb, dat werkte al, en de uitkomst voor P is dezelfde. Wie niet slim is moet af en toe een list bedenken

Ik heb mijn javascript nu aan de praat, zeer veel dank voor alle hulp!

[arie]

Mooi dat het werkt.
Het snijden van een lijn met een kromme had je in feite ook al gedaan in de algemene oplossing. Dat zou in de speciale gevallen dus ook moeten werken. Ik denk ook dat jouw testgevallen geen oplossingen (= snijpunten) hadden waardoor je in de problemen kwam.
Verhoudingen die allemaal dicht bij 1 liggen zouden wel moeten werken, bv 1.2 : 1.2 : 1.0.
P ligt dan dicht bij het driehoekscentrum (in ieder geval 1 van je P's).
Bedenk wel dat cheats de duidelijkheid en onderhoudsvriendelijkheid van je software doorgaans verminderen. Maar de afweging of dat wel of niet OK is ligt uiteindelijk bij jou.