8

I am working with Latitude / Longitude coordinates in a google map.

I have two lines :

  • Line A : 48.31508162629726, -2.591741396838972 to 48.40216156645915, -2.2218462112093404
  • Line B : 48.383816077371215, -2.274292940053768 to 48.66103546935337, -1.7066197241571377

I then use the following formula to find the point where they cross.

var XAsum = A.LngStart - A.LngEnd;
var XBsum = B.LngStart - B.LngEnd;
var YAsum = A.LatStart - A.LatEnd;
var YBsum = B.LatStart - B.LatEnd;

var LineDenominator = XAsum * YBsum - YAsum * XBsum;
if(LineDenominator == 0.0)
    return false;

var a = A.LngStart * A.LatEnd - A.LatStart * A.LngEnd;
var b = B.LngStart * B.LatEnd - B.LatStart * B.LngEnd;

var x = (a * XBsum - b * XAsum) / LineDenominator;
var y = (a * YBsum - b * YAsum) / LineDenominator;

This tells me that the lines do indeed cross and returns the x and y values.

However, when I plot the returned point, it is offset (not much) from the real intersection.

Is there a better and just as fast algorithm I could use which will return me the correct intersection point ?

It needs to be fast as I am iterating over a large number of lines (~1000).

EDIT : Note this is giving me an error offset of around 7.5 meters

Simon
  • 2,208
  • 4
  • 32
  • 47
  • What language are you programming in? Are there any libraries you can use (or are using) to help you solve GIS/spatial-analysis problems? – Jeff B Jul 25 '12 at 16:39
  • Did you solve your problem? I have the same challenge to find the intersection of two lines (arcs) on google maps. The Cartesian intersection produces errors on the map. For each arc I have a pair of latlons and I need to find their intersection on google maps. – Ash Apr 20 '15 at 01:34
  • No, I had to live with the error. Sorry – Simon Apr 20 '15 at 05:55
  • Try with negative lat (south part or earth). It don't work at all. I'll have a look at this (I'm in Brasil) – Peter Dec 17 '22 at 12:48
  • It seems when LineDenominator = 0, line are stricly paralell. But even if they are not, the question is: do you want to get intersection of lines or of segment? Because here, you seem to het intersection of the line and not of the segment. So in some case you get an intersection coord which is far from the segments. – Peter Dec 17 '22 at 13:01

2 Answers2

4

I'm assuming the algorithm you're using is the one for finding line intersections on a Cartesian coordinate system (i.e. a plane). Unfortunately, the earth is not a plane (or even a sphere) so using that algorithm will introduce error. Google Maps uses an ellipsoid (specifically WGS84) to approximate the surface of the earth, so you'll need an algorithm for finding the intersections of arcs on an ellipsoid.

This page might contains some helpful information: http://mathhelpforum.com/calculus/90196-point-intersection-two-lines.html

Jeff B
  • 8,572
  • 17
  • 61
  • 140
  • thank you for taking the trouble to answer. You are right in your assumptions. I looked at the page you mention, but I am not sufficiently advanced (!) in maths to be able to translate what is indicated to my use. Would you have any further suggestions ? – Simon Jul 25 '12 at 18:28
  • Any particular language you care to use? – Jeff B Jul 26 '12 at 00:57
  • Didn't have time to do much research, but this might help? http://geotools.org/ Or maybe this? http://geospatialmethods.org/spheres/ – Jeff B Jul 26 '12 at 17:40
0

After searching and trying to use your code, I discover the main problem came from the fact you're not checking for intersection of segment. In many case the coord of the intersection are outside of the segment, meanig you're checking for intersection of line.

Here is a piece of PHP code. Assuming you have 2 lines in GPS Coord (x1y1-x2y2 for the first one and x3y3-x4y4 for the second one)

$x1 = deg2rad($x1);
$y1 = deg2rad($y1);

$x2 = deg2rad($x2);
$y2 = deg2rad($y2); 

$x3 = deg2rad($x3);
$y3 = deg2rad($y3);

$x4 = deg2rad($x4);
$y4 = deg2rad($y4);     

$denom = (($y4 - $y3) * ($x2 - $x1)) - (($x4 - $x3) * ($y2 - $y1));

// Id denom = 0 -> parallele line
if ($denom != 0)
{
    $ua = (($x4 - $x3) * ($y1 - $y3) - ($y4 - $y3) * ($x1 - $x3))/$denom;
    $ub = (($x2 - $x1) * ($y1 - $y3) - ($y2 - $y1) * ($x1 - $x3))/$denom;

    if ($ua >= 0 && $ua <= 1 && $ub >= 0 && $ub <= 1)
    {
        $intersection_x = ($x1 + $ua*($x2 - $x1));
        $intersection_y = ($y1 + $ua*($y2 - $y1));
        echo "<br><b>Intersection: ".round(rad2deg($intersection_x),6).",".round(rad2deg($intersection_y),6)."</b>";
    }
    else
    {
        echo "<br>No intersection";
    }
}
Peter
  • 1,247
  • 19
  • 33