4

I'm trying to develop a function in R that will locate a point (A) on a map given distances between that point A and two other points (B and C), for which exact coordinates are known. The problem I'm coming up against is that my approach is sensitive to the shape of the triangle, and thus the number of conditions in my algorithm is higher than I think it needs to be, probably because my middle-school maths aren't good enough, so I'm looking for somebody to be really condescending, tell me how silly I'm being, and help me to a general approach to this problem.

For example, I'm looking for a single algorithm that can place point A in all of these circumstances:

  1. (x,y): B(1,2); C(4,1), distances: A to B = 6, A to C = 2
  2. (x,y): B(1,1); C(4,2), distances: A to B = 2, A to C = 3
  3. (x,y): B(1,1); C(4,1), distances: A to B = 2, A to C = 6

For the record, I'm aware that there are two possible points in space that satisfy these assumptions. I'm happy to have both.

Lastly, I'll point out that while this sounds like homework, it's actually ecological research, albeit performed by somebody apparently underqualified to be doing it.

  • Have you tried defining two circles, one each around B and C, using the distances as radii, then finding the intersection points of the two circles? There are tutorial web pages out there with the equations. – KevinMayall Sep 6 '14 at 0:20
1

KevinMayall's suggestion of intersecting circles is the easiest. A good theoretical treatment can be found at Wolfram Mathworld, which may get into more detail than you want to get this working (but is nice background for writing a paper about your methodology).

This StackOverflow answer has a broken link but outlines a basic approach.

Intersection of 2 circles

First calculate the distance d between the center of the circles. d = ||P1 - P0||.

  • If d > r0 + r1 then there are no solutions, the circles are separate.
  • If d < | r0 - r1| then there are no solutions because one circle is contained within the other.

  • If d = 0 and r0 = r1 then the circles are coincident and there are an infinite number of solutions.

Considering the two triangles P0P2P3 and P1P2P3 we can write

a2 + h2 = r02 and b2 + h2 = r12

Using d = a + b we can solve for a,

a = (r02 - r12 + d2 ) / (2 d)

It can be readily shown that this reduces to r0 when the two circles touch at one point, ie: d = r0 + r1

Solve for h by substituting a into the first equation, h2 = r02 - a2

So

P2 = P0 + a ( P1 - P0 ) / d

And finally, P3 = (x3,y3) in terms of P0 = (x0,y0), P1 = (x1,y1) and P2 = (x2,y2), is

x3 = x2 +- h ( y1 - y0 ) / d

y3 = y2 -+ h ( x1 - x0 ) / d

One important note, however, is that the mathematics here is based on Euclidean distances. The effect of curvature of the earth is not included. This is OK for small distances (< 20km, see this GIS FAQ, scroll to Q5.1) but if you're working on large analysis it begins to be inaccurate.

2

The gIntersection function from the rgeos library might help you. See the commented code below.

library(sp)  
library(rgeos)

pointB <- SpatialPoints(cbind(1,1))
pointC <- SpatialPoints(cbind(4,2))
distanceA2B <- 2
distanceA2C <- 3

# create the circle polygons around the points with the distances
polyB <- gBuffer(pointB, width = distanceA2B)
polyC <- gBuffer(pointC, width = distanceA2C)

# extract the feature coordinates of the polygon - we need to intersect lines, not polygons
coordsB <- lapply(polyB@polygons, function(x) {x@Polygons[[1]]@coords})
coordsC <- lapply(polyC@polygons, function(x) {x@Polygons[[1]]@coords})

# create circles as lines
lineB <- SpatialLines(list(Lines(Line(coordsB[[1]]), "B")))
lineC <- SpatialLines(list(Lines(Line(coordsC[[1]]), "C")))

# find intersections
if (is.null(gIntersection(lineB, lineC))) coords <- NULL else coords <- coordinates(gIntersection(lineB, lineC))
coords

plot(lineC)
lines(lineB)
points(coords, col="red", pch=16)
points(pointB, pch=4); text(coordinates(pointB) + .4, "Point B")
points(pointC, pch=4); text(coordinates(pointC) + .4, "Point C")

enter image description here

(inspired by https://stackoverflow.com/a/21648987)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.