16

Please refer to the example and corresponding image.

I would like to achieve the following: provide two locations (lat/lng), which are shown below as A and B. From this, a virtual line would be drawn and then the distance between this line and C would be calculated (in any measurement).

drawing

I have achieved this currently in Google Maps API v3 but would want to also be able to perform this behind the scenes in my language of choice. Any tips/ideas would be greatly appreciated!

whuber
  • 69,783
  • 15
  • 186
  • 281
Prisoner
  • 283
  • 1
  • 3
  • 10
  • @Kirk, No, AB is just a straight line – Prisoner Jun 23 '11 at 16:14
  • @Michael, thats an interesting point. I will have to have a look into it! – Prisoner Jun 23 '11 at 16:15
  • @Prisoner @Kirk Literally, a "straight line" will pass beneath the earth's surface. In general, its radial projection back onto the surface will indeed be a segment of a great circle (using a spherical earth model). – whuber Jun 23 '11 at 16:26
  • @whuber, For the distance that are being used (20,30 meters) this shouldn't be too much a problem. Correct? – Prisoner Jun 23 '11 at 16:27
  • 1
    @Prisoner That is an extremely useful extra piece of information! Yes, you are correct. You still have to compensate for the fact that using (lat,lon) differentially distorts east-west distances compared to north-south. As @Jose advises, project the coordinates. This can be as simple as pre-multiplying the longitudes by the cosine of the average latitude and then pretending you're on a Euclidean plane. – whuber Jun 23 '11 at 16:32
  • I cannot see the image, however, I can proide some simple code in python which takes A and B (two points representing a straight line), C, which is the offset point and returns d, a point perpendicular to the point. In it, it contains the shortest distance which is the measurable you want, I believe? – Hairy Jun 24 '11 at 08:42
  • Sorry to bring this thread back but I found some issues with it. If using the a point where a lat/long has the same value as one of the line segment points it will not work. Some thing if a line is due north due east/west or due south as it relies on the cross products which would be 0. I have taken this code and implemented in PHP Function IntersectionDAN ($startLat, $startLon, $endLat, $endLon, $pointlat, $pointlon) { //Find Instersection Point $startLat = deg2rad($startLat); $startLon = deg2rad($startLon); $endLat = deg2rad($endLat); $endLon = deg2rad($endLon); $pointlat = deg2rad($pointlat – dan Jun 13 '14 at 17:06

4 Answers4

11

Maybe I'm making it too complicated, but what you want is the distance from a point to a line. That is the distance from a point along AB that links AB with C with a line orthogonal to AB. This vector perpendicular to AB is given by

v=[x2-x1, -(y2-y1)] # Point A is [x1,y1] Point B is [x2,y2]

(I have used the square brackets to define a vector, or a two-element array). The distance between C [xp, yp] and point A is

u=[x1-xp, y1-xp]

The distance between the line and C is just the projection of u on to v. If we assume mod(v) = 1 (just normalise it), then

distance = u*v = abs( (x2-x1)*(y1-yp) - (x1-xp)*(y2-y1) )

The only complication is that you probably want to make sure your coordinates are not WGS84 lat/log pairs, but projected (or use geodetic coordinates). You can use OGR or Proj4 for this.

Glorfindel
  • 1,096
  • 2
  • 9
  • 14
Jose
  • 3,332
  • 20
  • 21
  • 3
    +several million pseudo-points for not using trigonometric functions, by the way. All too many people pull out ArcTan when they should be looking at this: http://en.wikipedia.org/wiki/Dot_product – Herb Jun 23 '11 at 16:11
  • @Jose, thanks for the reply! I'm using the lat/long from google maps API. The maths part is quite new to me so I'll give it a go and see what I can come up with. Any tips with the maths? E.g. [x2-x1, -(y2-y1)], what does that translate to? – Prisoner Jun 23 '11 at 16:21
  • I have added a short edit for this. Basically, it's an array notation, but if you store your coordinates in variables x1, x2, y1, y2 and xp, yp, you only need to write the right hand side of the last equation I've provided. This is pretty much valid C, Java, JS, Python etc code :) – Jose Jun 23 '11 at 16:26
  • 2
    @Jose You are computing the distance from C to the line AB. Based on the figure, I believe the OP wants the distance from C to the line segment AB. This requires extra work to check whether the projection of C onto line AB lies between A and B or not. In the latter case, use the shorter of the two lengths CA and CB. – whuber Jun 23 '11 at 16:29
  • @whuber You are right, but he only mentions line in his post :) Do I still get away with it? :) – Jose Jun 23 '11 at 16:38
  • @Jose Now that we have clarified what your solution actually does, +1. – whuber Jun 23 '11 at 16:42
  • What is the difference between line and line segment? – Prisoner Jun 23 '11 at 17:08
  • 2
    @Prisoner The main difference is that a line extends forever (it is only defined by a direction vector and a point, or by two points), whereas a segment between A and B is the bit of the the infinite line that goes between A and B (it has a finite length) – Jose Jun 23 '11 at 17:25
7
def get_perp( X1, Y1, X2, Y2, X3, Y3):
    """************************************************************************************************ 
    Purpose - X1,Y1,X2,Y2 = Two points representing the ends of the line segment
              X3,Y3 = The offset point 
    'Returns - X4,Y4 = Returns the Point on the line perpendicular to the offset or None if no such
                        point exists
    '************************************************************************************************ """
    XX = X2 - X1 
    YY = Y2 - Y1 
    ShortestLength = ((XX * (X3 - X1)) + (YY * (Y3 - Y1))) / ((XX * XX) + (YY * YY)) 
    X4 = X1 + XX * ShortestLength 
    Y4 = Y1 + YY * ShortestLength
    if X4 < X2 and X4 > X1 and Y4 < Y2 and Y4 > Y1:
        return X4,Y4
    return None

The shortest length is the distance you require, unless I am mistaken?

Hairy
  • 4,843
  • 1
  • 26
  • 45
  • Yes, I am looking for the shortest distance from C to the line segment. Is this what this math computes? – Prisoner Jun 24 '11 at 09:30
  • this gives you the point on the line and the shortest distance. – Hairy Jun 24 '11 at 09:35
  • Excuse the ignorance, but what format is X1-Y3? As mentioned I'm using lat and longs but I assume this is not what I need! – Prisoner Jun 24 '11 at 10:18
  • x and y are your points: x1,y1 and x2,y2 is the line, x3,y3 is the offset point. Lat longs are simply x,y pairs? – Hairy Jun 24 '11 at 10:36
  • I realise they are x and y, I just wasnt sure if they needed to be converted to anything. I'll give it a go and regardless, thanks a lot for the answer! – Prisoner Jun 24 '11 at 10:46
  • No, they should work fine; effectively, you're just operating on a grid. – Hairy Jun 24 '11 at 10:55
  • 1
    It did indeed work fine, I passed the three points (A,B,C) in the following: http://i.imgur.com/bK9oB.jpg and it came back with the lat/lng of X. Great work! – Prisoner Jun 24 '11 at 11:59
  • Whuber may disagree with me, but things like this do belong out of the GIS and in Math where it works better. I have quite a few little utility scripts like this which work well and I am only glad this wored for you. – Hairy Jun 24 '11 at 12:01
  • 1
    @Hairy, One last thing, how would I go about modifying this to go to the nearest point (not just line) so if I had put it past the point of joing a line, how could I get it to check the distance to the point? – Prisoner Jun 24 '11 at 14:47
  • 1
    @Hairy Good start, but it seems that too often this code returns None when a legitimate solution exists. The problem is that the last conditional assumes X1 < X2 and Y1 < Y2, which cannot always be assured. A better test of betweenness is needed. – whuber Jun 26 '11 at 18:11
  • I disagree whuber. If you have produced a point,. which is on a line segment between x1,y1 and x2,y2, then this has to be right,. If if X4 > X2, or < X1 then it isn't on the segment between X1 and X2. I assumed the object was to create a point on this segment which means that last conditional is essential – Hairy Jun 27 '11 at 09:37
  • @whuber, can I ask how I have been marked down 15 points on this post without it being noted? – Hairy Jun 27 '11 at 09:38
  • This was accepted as the answer and I appear to have lost those points; cna they be reinsted please? – Hairy Jun 27 '11 at 09:39
  • I had actually made those changed with the condition and it works as expected, the reason I had unmarked it as an answer was due to the reason stated above: if the point was within the correct bound it returned nothing (on some occasions) and if I removed the check it returned correctly for the line not the line segment. Anyway, the check above has been changed to look like: if(($pointCrossSectionLat < $endLat && $pointCrossSectionLat > $startLat) || ($pointCrossSectionLng < $endLng && $pointCrossSectionLng > $startLng)) – Prisoner Jun 27 '11 at 11:28
  • AH, ok. But surely, if the x4 is within x1 and x2 it wouldn't have failed? I am slightly confused as to how this logic could fail? – Hairy Jun 27 '11 at 11:34
  • Sorry whuber/Prisoner, I ahve understood what you are saying ehre. Yes, the presumption us that it is operating where X1 is greater than X2 and that the normalisation of the data happens pre processing i.e. X1 could be -/+ X2 etc, so as to simplify the distance calculations.; in short, yes, it could be written better and as a standalone function, but I put it here for Prisoner to work with as a start in the first instance. I pre process my data in order the conditional is never false, but is there as a matter of completeness. – Hairy Jun 27 '11 at 11:42
  • 1
    @Hairy It sounds like this exchange between you and @prisoner has been productive. I would like to emphasize that I have had nothing to do with (or even any control over) any changes in voting or points that may have occurred and that my comment was intended only to help you improve your reply. – whuber Jun 27 '11 at 13:18
  • I appreciate that @Whuber and apologise if anything unpleasant was inferred. I didn't understand the statement, as I pre process. Thanks for your comments. – Hairy Jun 28 '11 at 06:57
  • Oh well, @Hairy you have your points and I have a working function. Everyone wins! – Prisoner Jun 28 '11 at 08:00
4

Being a bit averse to all this math as well, I would come at it from a different angle. I would make it an 'actual' line, rather than a virtual line, and then use existing tools.

If A and B share an attribute, you could connect them by drawing a line (Kosmo GIS has a tool that will create lines from points, and I believe there is also a QGIS plugin for this). Once you have the lines, a 'near' function on the 'C' point layer will give you the distance to the line. Let the software handle the math for you!

whuber
  • 69,783
  • 15
  • 186
  • 281
Darren Cope
  • 6,616
  • 2
  • 32
  • 47
  • Thanks for the comment but Hairy came up trumps on this one! – Prisoner Jun 24 '11 at 12:01
  • 1
    (+1) You make an excellent point. Computational geometry algorithms are notoriously tricky to get completely right in practice (as we can see from all the code offered so far, which is helpful and illustrative but doesn't yet fully work). Using a high-level GIS procedure often is a good way to assure you are getting the answer you expect and that it's right (provided you trust your GIS ;-). – whuber Jun 27 '11 at 13:22
1

If you were using java on android, its only one line with the library function

import static com.google.maps.android.PolyUtil.distanceToLine;

distanceToLine:

public static double distanceToLine(LatLng p, LatLng start,LatLng end)

Computes the distance on the sphere between the point p and the line segment start to end.

Parameters: p - the point to be measured

start - the beginning of the line segment

end - the end of the line segment

Returns: the distance in meters (assuming spherical earth)

Just add library to your

dependencies {
    compile 'com.google.maps.android:android-maps-utils:0.5+'
}
Prisoner
  • 283
  • 1
  • 3
  • 10
Indy
  • 11
  • 1