Calculating coordinates of points from just distances, no bearings

The task I'm giving myself is about enabling people with little equipment to produce a map of a forest, or a garden with an important presence of trees.

I participate to openstreetmap and I plan to (enable people to) contribute tree locations for botanic gardens in areas where Bing aerial coverage is not available, or where it is useless because tree canopy is too thick to allow for identification of individual trees.

My hope is that there is a QGIS plug-in (and if it does not exist, to develop one myself) that accepts an amount of reference points distributed around the area of interest, a matrix of distances among points (referenced or not), and populates a vector(points) layer with the calculated 3D coordinates.

I have been looking at trilateration, and I found this 22 years old article.

I suppose the math still holds.

I do not want to work with angles, only distances.

not yet an answer, but I have been working at it, and I have a couple more details about the issue.

• I need my input coordinates in metres (or feet, or millimetres), not in degrees, because I need an isometric coordinate system.
• The math in that article still holds, and numpy helps a lot, specifically, the article suggests solving the overdefined problem by least squares, and numpy contains the function np.linalg.lstsq.
• I gave up hope to work in 3D and I am going to assume earth is flat, but really flat like the Netherlands, or at least that the digital device I will acquire to measure distances is able to split distances in their horizontal and vertical components.
• I have a working prototype, but I still need to discover how to create a point feature using the coordinates in the coordinate system of the layer.

the code for the coordinates calculation is not at all long, while the difficult part was about choosing which point I can reference given the points I have already calculated:

def find_point_coordinates(points, distances, point_id):
import numpy as np
from numpy.core.umath_tests import matrix_multiply

connected_to = [id for id in sorted(distances[point_id]) if points[id].get('coordinates')]
connected_matrix = np.array([points[id]['coordinates'] for id in connected_to])
A = connected_matrix[1:,] - connected_matrix[0,]
A = 1.0 * A  # make sure we work with floating point values

## squared distances vector, beacon_i to first beacon for which we have distances
D_i1_2 = matrix_multiply(A*A, [,])
## distances of targeted point from used reference points
dfb_sel = np.array([distances[point_id][ref_id] for ref_id in connected_to])
r2 = dfb_sel * dfb_sel

rhs = ((r2 - r2[1:]).reshape(D_i1_2.shape) + D_i1_2) / 2.0
r1, r2, r3, r4 = np.linalg.lstsq(A, rhs.reshape(rhs.shape[:1]))
return connected_matrix[0,] + r1
• you are essentially solving the same problem as a GPS receiver. You know the distances to multiple points, and you know the exact locations of those points. I suggest you read up on how a GPS receiver derives its location – tomfumb Jan 3 '17 at 17:42
• Surveying Calculation plug-in can solve arc section (using hundred year old formulas) in 2D. For a 3D solution zenith angles are necessary. – Zoltan Jan 3 '17 at 17:51
• I am reading about it, and writing code that I can test, and understanding details of the problem. As of now, I think that I will look for ways to work in 2D. The question would be if a "Digital Handheld Laser Distance Meter" will give me decent horizontal estimations given a measurement taken at an angle. – mariotomo Jan 3 '17 at 21:59
• in my case, I needed transf = QgsCoordinateTransform(QgsCoordinateReferenceSystem(3117),QgsCoordinateReferenceSystem(4326)) and then transf.transform(QgsPoint(point, point)) – mariotomo Jan 4 '17 at 20:16
• @Zoltan, I have installed the surveying calculation plugin, it looks a lot more powerful than what I need, I will reserve some time to study it. thank you for the hint. – mariotomo Jan 5 '17 at 0:33