[Webinar] Streamline your web hosting managementRegister Today

x
  • Status: Solved
  • Priority: Medium
  • Security: Public
  • Views: 3915
  • Last Modified:

PHP function to calculate bearing between latitudes and longitudes

I need to be able to calculate the bearing between two latitude and longitude pairs, in PHP. They are decimal degrees, and since we are in Australia, our points we are working with are negative latitude, positive longitude. eg. -29.54546 Latitude and 149.68756 longitude.

I have tried this:
$bearingradians = atan2(asin($long1-$long2)*cos($lat2),
cos($lat1)*sin($lat2) - sin($lat1)*cos($lat2)*cos($long1-$long2));
$bearingdegrees = rad2deg($bearing);

This seems to work, but on closer inspection, it gets east and west mixed up! ie: When I know that a point should be west of a certain point, it comes up with say, 88degrees, which is obviously east! I am putting the lats and longs around the right way, so thats not the problem. I am trying to make a function which will take the four numbers and spit out the bearing between them, like this:

function bearing($lat1, $long1, $lat2, $long2){

$bearingradians = atan2(asin($long1-$long2)*cos($lat2),
cos($lat1)*sin($lat2) - sin($lat1)*cos($lat2)*cos($long1-$long2));
$bearingdegrees = rad2deg($bearing);

return $bearingdegrees
}

Can anyone see what im doing wrong mathematically? I need this working asap...
0
Dr_Snapid
Asked:
Dr_Snapid
  • 10
  • 6
1 Solution
 
catknowsCommented:
My first impression is that the unit is wrong:


1. function bearing($lat1, $long1, $lat2, $long2){
2. $bearingradians = atan2(asin($long1-$long2)*cos($lat2),
3. cos($lat1)*sin($lat2) - sin($lat1)*cos($lat2)*cos($long1-$long2));
4. $bearingdegrees = rad2deg($bearing);

5. return $bearingdegrees
6. }

On line 3, the first term and second term are not of the same unit. You cannot do subtraction between the two term.

If you can tell me how you derive this formular, you can help you right away.

0
 
catknowsCommented:
Continuing above:

I don't know how you derive your formula, but if I use formula #6 on of the below link,

http://mathworld.wolfram.com/TrigonometricAdditionFormulas.html

This is what I get:

function b2($y1, $x1, $y2, $x2)
{

  $bearingradians = atan2( $y1*$x2 - $y2*$x1, $x1*$x2 + $y1*$y2 );

  $bearingdegrees = rad2deg($bearingradians);

  return $bearingdegrees;
}
0
 
catknowsCommented:
Continuing above:

Just a note: The function I wrote takes values that are the coordinates of the longitude and latitude pairs.  The bearing angle is the angle on the plane that is tangent to the surface of the earth, assuming the two point are close enough.  Computing the angle bearing relative to the center of the earth would involve a little more.

Cheers.

Long
0
Never miss a deadline with monday.com

The revolutionary project management tool is here!   Plan visually with a single glance and make sure your projects get done.

 
Dr_SnapidAuthor Commented:
Cat knows, heres the code I tried, it didnt work.


function getbearing($y1, $x1, $y2, $x2)
{

  $bearingradians = atan2( $y1*$x2 - $y2*$x1, $x1*$x2 + $y1*$y2 );

  $bearingdegrees = rad2deg($bearingradians);

  return $bearingdegrees;
}


#test the function
// Grafton NSW
$lat1 = -29.6844;
$long1 = 152.9325;
//Coffs Harbour NSW, about 170 degrees from grafton
$lat2 = -30.31319;
$long2 = 153.13753;

echo "The bearing from Grafton to Coffs Harbour is ".getbearing($long1,$lat1,$long2,$lat2);
-------------------------------------------------------------------------------------------------------------
This returns -0.212214774169, which is wrong.... Should be about 170 or so. Can you see anything wrong?
0
 
catknowsCommented:
Oh, I got it.

By the way, if you can tell me from what source you derive your formula that would confirm my finding, too.  Thanks.

My equation compute the angle formed by the vectors having the coordinates of x1,y1, and x2,y2 that is why it does not come out right.

But hang on. Will get back to you.
0
 
Dr_SnapidAuthor Commented:
It came from Dr. Math, and now I cant find the URL, will post as soon as i find it.
0
 
catknowsCommented:
Ok,
Try this...

function getbearing($y1, $x1, $y2, $x2)
{

  $bearingradians = atan2( tan($y2-$y1), tan($x2-$x1)  );

  $bearingdegrees = rad2deg($bearingradians);

  return $bearingdegrees;
}

This is the angle originated from the destination relative to the source.
0
 
Dr_SnapidAuthor Commented:
Jeez we're close now, worked perfectly for grafton to Coffs Harbour, but then was -108.95732909deg for Grafton to Armidale, which should be 233.2deg.

So you can check:
Grafton -29.6844,152.9325
Coffs Harbour -30.31319,153.13753 (should be 164deg from Grafton, we got 164)
Armidale -30.51194,151.6675 (Should be 233deg from Grafton, we got -109)
Ballina -28.87398,153.56237 (Should be 34deg from Grafton, we got 34)
Tenterfield -29.05639, 152.01806 (should be 309deg, we got -60)

Looks like westerly directions are out.... if I add 360 degrees to them to get a positive value it is still out by a fair way...
You are good at this, what have we missed?
0
 
Dr_SnapidAuthor Commented:
If it helps, I chose those towns for being on opposite sides of grafton, Ballina is NE, Coffs is SE, Armidale is SW and Tenterfield is NW.
0
 
catknowsCommented:
Adding 360 is the correct way to get the number rightly orientated.  But since the numbers are still out, I can't think of any otherway.  I compute this purely on mathematically point of view, that is, assuming the earth is perfectly round, and angle are measured as though you are standing to the origin, which is Grafton, of the x-y coordinate.

Where do you get the measurements from?

The reason I am asking is that the formula I used is based on the tangent plane that touch the earth the origin where the angle is located.  If your measure is on the flat map, the angle might not come out accurate.  If you get your data out of a table, can you tell me where the source is?  If I can understand what the data you have really meant what they meant.
0
 
catknowsCommented:
I just add a line to the function to take care of the 360 degree wrap around:

>>
function getbearing($y1, $x1, $y2, $x2)
{

  $bearingradians = atan2( tan($y2-$y1), tan($x2-$x1)  );

  $bearingdegrees = rad2deg($bearingradians);
  if( $bearingdegrees < 0 )
    $bearingdegrees += 360;

  return $bearingdegrees;
}
>>

The result:

Coffs Harbour:
x2-x1=-0.62879
y2-y1=0.20502999999999
The bearing from Grafton to Coffs Harbour: 164.04275769504

Armidale:
x2-x1=-0.82754
y2-y1=-1.265
The bearing from Grafton to Armidale: 251.04267091014 [233]

Ballina:
x2-x1=0.81042
y2-y1=0.62986999999998
The bearing from Grafton to Ballina: 34.734391266122

Tenderfield:
x2-x1=0.62801
y2-y1=-0.91444000000001
The bearing from Grafton to Tenderfield: 299.21742521454 [309]

You can see that the larger the change in the latitude in the negative direction (more south), the worse it appears in the difference between your data and the formula results.

I think this is a matter of practical values, and theoretical values.

If you want to dig more into this, just let me know.

Cheers.
Got to go. Too hungry now.

0
 
Dr_SnapidAuthor Commented:
Hmmm I just had lunch. Check out the Chicken Fillet Burgers ant Country Fresh... Mmmmm!

So do you think there is a way to correct the function? It just seems to me there must be another thing weve missed... we are so close...

Someone has done it... look at this:

http://www.zeuswireless.com/support/lat-long.php

It works well. Just wish I could see his code...
0
 
catknowsCommented:
Ok, try this. The latest and greatest.

Cheers. By the way, no chicken for me. Thanks for the suggestion.


<?
function bearing( $lat1_d, $lon1_d, $lat2_d, $lon2_d )
{

   $lat1 = deg2rad($lat1_d);
   $lon1 = deg2rad($lon1_d);
   $lat2 = deg2rad($lat2_d);
   $lon2 = deg2rad($lon2_d);

   $L    = $lon2 - $lon1;

   $cosD = sin($lat1)*sin($lat2) + cos($lat1)*cos($lat2)*cos($L);
   $D    = acos($cosD);
   $cosC = (sin($lat2) - $cosD*sin($lat1))/(sin($D)*cos($lat1));
   
    $C = 180.0*acos($cosC)/pi();

    if( sin($L) < 0.0 )
        $C = 360.0 - $C;

    return $C;

}
?>
0
 
catknowsCommented:
Continue...

Your formula need to be fixed just a little bit then it will come out just fine, too.

function bearing( $lat1_d, $lon1_d, $lat2_d, $lon2_d )
{
   $lat1 = deg2rad($lat1_d);
   $long1 = deg2rad($lon1_d);
   $lat2 = deg2rad($lat2_d);
   $long2 = deg2rad($lon2_d);

$bearingradians = atan2(asin($long2-$long1)*cos($lat2),cos($lat1)*sin($lat2) - sin($lat1)*cos($lat2)*cos($long2-$long1));
$bearingdegrees = rad2deg($bearingradians);
$bearingdegrees = $bearingdegrees < 0? 360 + $bearingdegrees : $bearingdegrees;


return $bearingdegrees;
};

Cheers!
0
 
Dr_SnapidAuthor Commented:
Legend! I knew you could do it!
0
 
catknowsCommented:
Thank you! I am very glad that that helps.
0

Featured Post

Never miss a deadline with monday.com

The revolutionary project management tool is here!   Plan visually with a single glance and make sure your projects get done.

  • 10
  • 6
Tackle projects and never again get stuck behind a technical roadblock.
Join Now