Bearing Calculation

Date: 09/01/97 at 00:39:54
From: Hrvoje Vojnic
Subject: Bearing calculation

Hello!

I have a question regarding the calculation of bearings.

Given two cities of geographic coordinates (xA,yA) and (xB,yB) is 
there a formula which calculates the bearing from city A to city B?  
If so, would you please write it, along with any specifics which it 
may include (angle in radians or degrees, distance in miles or 
kilometers, etc.)

Thank you in advance!

Best regards,

Hrvoje

Date: 09/01/97 at 14:56:04
From: Doctor Jerry
Subject: Re: Bearing calculation

Dear Hrvoje,

Do you mean the angle made by the great circle route with the meridian 
through city A?  Please let me know and I'll try to answer your 
question.

-Doctor Jerry,  The Math Forum
 Check out our web site!  http://mathforum.org/dr.math/   

Date: 09/04/97 at 01:49:03
From: Hrvoje Vojnic
Subject: Re: Bearing calculation

Hello again, Dr. Jerry!

Regarding my question, (bearing from city A to city B, given their 
geographic coordinates) I apologize for not being precise. It is 
exactly like you understood my question, I'm interested how to 
calculate the angle between the great circle route from city A to 
city B and the meridian of city A. That is, the bearing at city B 
from city A. I know how to calculate the great circle distance, but 
bearing is a puzzle.

Thank you!

Hrvoje

Date: 09/04/97 at 12:17:30
From: Doctor Jerry
Subject: Re: Bearing calculation

Hi Hrvoje,

It's easiest for me to explain in terms of an example.  

I'll use spherical coordinates, which are closely related to the 
geographic coordinates of latitude and longitude, and I'll work on a 
sphere with radius 1. I think the solution procedure will be clear for 
the general case. I don't know that the solution I'll explain is the 
simplest. If one were to read a book written for navigators, I presume 
that relatively simple ways to calculate the bearing of a great circle 
route would be explained.  However, from my experience in looking at 
such books, to understand the procedures I would have to read more of 
the book than I want.  

I'll consider points S (starting point) and E (ending point); 
specifically, S has spherical coordinates (1,pi/4,pi/6) and 
E has spherical coordinates (1,pi/6,pi/3).  These correspond to 
(on the unit sphere) points with latitudes 45 deg N and 60 deg N 
and longitudes differing by 30 deg.

I'll use the spherical to rectangular transformation

r(phi,theta)={r*sin(phi)*cos(theta),r*sin(phi)*sin(theta),r*cos(phi)}.

Start by calculating the rectangular coordinates of S and E; call 
these P1=r(pi/4,pi/6) and P2=r(pi/6,pi/3). Think about the unit 
vector u joining the origin to the point (x1,y1,z1). The way I wrote 
the r function, using curly brackets { }, already puts us in the 
vector format. 

Form the unit vector v lying in the plane through the origin, P1, 
and P2 (the intersection of this plane with the sphere is the great 
circle we're working with), as follows. I'll give part of the 
Mathematica code I used to check all of this. The code is a litte 
peculiar, but I think you can follow it with no difficulty.

r[phi_,theta_]:={Sin[phi]*Cos[theta],Sin[phi]*Sin[theta],Cos[phi]} 

P1=r[Pi/4,Pi/6]; P2=r[Pi/6,Pi/3];

u=P1;

The vector w (below) is the vector along u onto which the vector P2 
projects.  The expression P1.P1 is dot product. Lgth is the length 
function.

w=P2-((P1.P2)/(P1.P1))P1

v=w/Lgth[w];

The great circle route through S and E is the parametric curve

r1[t_]:=Cos[t]*u+Sin[t]*v

for 0 <= t<= 2*pi.

Note that r1[0] is the point S.  To calculate the bearing I'll 
calculate the tangent vector r1'[0] and dot it with the north-pointing 
unit vector m through S.  Differentiating r1,

r1'(t)=-sin(t) + cos(t) v.

Hence, r1'(0)=v. Also,

m={-sqrt(3)/(2*sqrt(2)), -1/(2*sqrt(2)), 1/sqrt(2)}

This vector can be calculated by forming the negative of the tangent 
vector to the curve with parametric equation r(phi,pi/6), and 
evaluating at phi=pi/4.  The bearing is then 

b1=ArcCos[(m.v)/(Lgth[m]*Lgth[v])]

In particular, 

b1=ArcCos[3*sqrt(3)/(8*sqrt(5))+sqrt(15)/8]=0.684719 or 39.23 degrees.

Although this procedure probably will require some special cases for 
points in different hemispheres and other quirks of spherical or 
geographic coordinates, it is fairly straightforward.  

I checked all of this by doing it a different way.  I assumed the 
equation of the plane through the center of the earth, S, and E had 
the form A*x+B*y+z=0 and, with the rectangular coordinates of S and E, 
calculated A and B. This gave the plane -2y+z=0.  I solved this 
simultaneously with x^2+y^2+z^2=1, thinking of z as fixed.

This gave me parametric equation (parameter z) for the great circle.  
I found

r2[z_]:={Sqrt[4-5z^2]/2,z/2,z}

I calculated r2'[0] and dotted it with the meridian vector as before 
and got the same bearing.

-Doctor Jerry,  The Math Forum
 Check out our web site!  http://mathforum.org/dr.math/