Go Premium for a chance to win a PS4. Enter to Win

x
Solved

# VB Solution to return the angle, in degrees, from a given lat, lon to current position.

Posted on 2006-11-10
Medium Priority
799 Views
Hi,

Basically the title says it all. What I need is a routine which will return the angle, in degrees, from a given point (given as Lat,Lon) to my current position (also given as Lat,Lon), using VB6.

Thanks

Dave
0
Question by:dmarch100
• 5
• 5
• 3

LVL 22

Assisted Solution

danaseaman earned 750 total points
ID: 17919661
Here are 2 Functions. 1 for Bearing and another for GeoDistBetween.

'To call Bearing:

Heading = Bearing(LatDest - LatOrig, LongDest - LongOrig)

' ----------------------------------------------------------------------
' Inverse longitudinal tangent function to return angles in degrees.
' This function uses dual arguments and returns an angle in the range
' from 0 to 360 degrees relative to True North.

Public Function Bearing(ArgY, ArgX) As Single
'  Check for special coordinate zero value cases
On Error GoTo Bearing_Err
100   If ArgY = 0 And ArgX >= 0 Then
102      Bearing = 90
104   ElseIf ArgY > 0 And ArgX = 0 Then
106      Bearing = 180
108   ElseIf ArgY = 0 And ArgX < 0 Then
110      Bearing = 270
112   ElseIf ArgY < 0 And ArgX = 0 Then
114      Bearing = 0
Else
116      Bearing = 90 + 45 * Atn(ArgY / ArgX) / Atn1
118      If ArgX < 0 Then
120         Bearing = Bearing + 180
End If
122      If Bearing < 0 Then
124         Bearing = Bearing + 360
End If
End If

'<EhFooter>
Exit Function

Bearing_Err:
Select Case ErrMsgBox("Cyberime2002.modAstro.Bearing")
Case vbAbort
Screen.MousePointer = vbDefault
Exit Function
Case vbRetry
Resume
Case vbIgnore
Resume Next
End Select
'</EhFooter>
End Function

'--------

Private Function GeoDistBetween(Long1 As Single, Lat1 As Single, Long2 As Single, Lat2 As Single) As Single
' V2.0
'
' Compute the gedesic distance between two Earth surface
' coordinates to an accuracy of about Â±50 meters.  To get
' this degree of accuracy, this function takes into account
' the spheroidal flattening factor of the Earth rather than
' assuming that the Earth is a perfect sphere.
'
' The arguments (Lng1, Lat1, Lng2, Lat2) are the longitudes
' and latitudes of the two locations expressed in degrees.
On Error GoTo GeoDistBetween_Err

' Positive longitude is west and negative is east.
'
' ===========================================================
' Function written by Jay Tanner - Based on the mathematical
' method of H. Andoyer.
'
' References:
'
' International Earth Rotation Service,
' Annual report for 1996 (Observatorie de Paris, 1997)
'
' Annuaire du Bureau des Longitudes pour 1950 (Paris) pg.145
'
' ==========================================================
' Flattening factor of the earth geoid
Const ff As Double = 1 / 298.257

' Auxiliary working variables
Dim c  As Double
Dim D  As Double
Dim F  As Double
Dim g  As Double
Dim H1 As Double
Dim H2 As Double
Dim L  As Double
Dim O  As Double
Dim r  As Double
Dim s  As Double
Dim W1 As Double
Dim W2 As Double
Dim SG As Double
Dim cG As Double
Dim SF As Double
Dim CF As Double
Dim sl As Double
Dim CL As Double

' Compute auxiliary angles
100   F = (Lat1 + Lat2) / 2
102   g = (Lat1 - Lat2) / 2
104   L = (Long1 - Long2) / 2
106   If g + L = 0 Then
108      MsgBox "Warning." & vbCrLf & "Source and Destin are equal."
Exit Function
End If
' Compute sines and cosines of auxiliary angles
110   SG = Sin(g * Atn(1) / 45)
112   cG = Cos(g * Atn(1) / 45)
114   SF = Sin(F * Atn(1) / 45)
116   CF = Cos(F * Atn(1) / 45)
118   sl = Sin(L * Atn(1) / 45)
120   CL = Cos(L * Atn(1) / 45)

122   W1 = SG * CL: W1 = W1 * W1
124   W2 = CF * sl: W2 = W2 * W2
126   s = W1 + W2

128   W1 = cG * CL: W1 = W1 * W1
130   W2 = SF * sl: W2 = W2 * W2
132   c = W1 + W2

134   O = Atn(Sqr(s / c))
136   r = Sqr(s * c) / O
138   D = 2 * O * 6378.14

140   H1 = (3 * r - 1) / (2 * c)
142   H2 = (3 * r + 1) / (2 * s)

' Compute the angle between the points on a synthetic
' geodesic sphereoid connecting the two points.
144   W1 = SF * cG: W1 = W1 * W1 * H1 * ff + 1
146   W2 = CF * SG: W2 = W2 * W2 * H2 * ff

' Return the distance between the given locations in
' Kilometers
148   GeoDistBetween = D * (W1 - W2)

'<EhFooter>
Exit Function

GeoDistBetween_Err:
Select Case ErrMsgBox("Cyberime2002.frmWorld.GeoDistBetween")
Case vbAbort
Screen.MousePointer = vbDefault
Exit Function
Case vbRetry
Resume
Case vbIgnore
Resume Next
End Select
'</EhFooter>
End Function

Public Function ErrMsgBox(Msg As String) As Integer
'*****WARNING*****
'Do not put error handling in this sub.
'Doing so will clear the Err Object and Erl

100   Msg = "Error Number " & Err.Number & vbCrLf & _
Err.Description & vbCrLf & _
"in " & Msg & vbCrLf & _
"at line " & Erl
102   ErrMsgBox = MsgBox(Msg, vbAbortRetryIgnore + vbCritical, "Error details")

End Function
0

Author Comment

ID: 17920461
Hi danaseaman

Thanks for your answer but how does the 'Bearing' function relate to my current position in any way? I need to know the angle in degrees from a given point relative to my current location

Dave
0

LVL 22

Expert Comment

ID: 17922523
You call it with 4 parameters. Latitude and Longitude of both source and destination.
It will return the angle relative to the source Lat/Long.

Heading = Bearing(LatDest - LatOrig, LongDest - LongOrig)
0

Author Comment

ID: 17924829
Hi danaseaman,

Sorry for not noting the call procedure as clearly shown in your original post.

However, alsthough on the right track it doesn't actually work. If we concentrate on the Bearing function only at line 116 the Atn1 variable is not defined or used and subsequently kicks out a Division by Zero error. With that variable removed from the equation it works but obviously returns an incorrect radial.  My test data is Glasgow airport as the current position and Heathrow as the destination....

Heading = Bearing(51.28 - 55.52, 0.27 - 4.25)

Dave
0

LVL 22

Expert Comment

ID: 17925527
Oops, forgot the constant.

Public Const Atn1     As Double = PI / 4

0

LVL 22

Expert Comment

ID: 17926075
Public Const PI       As Double = 3.14159265358979
Public Const PI2      As Double = PI * 2
Public Const Atn1     As Double = PI / 4
0

LVL 28

Expert Comment

ID: 17927525
Which bearing do you need - rhumb line or great circle?
0

LVL 28

Accepted Solution

Ark earned 750 total points
ID: 17927719
Public Enum bTypes
RhumbLine
GreatCircle
End Type

Public Type POSITION
lat As Double
lon As Double
End Type

Const PI As Double = 3.14159265358979
Const e As Double = 0.081813333

Public Function Azimuth(ptFrom As POSITION, ptTo As POSITION, Optional bType As bTypes = RhumbLine) As Double
If bType = GreatCircle Then
Azimuth = AzimuthGC(ptFrom, ptTo)
Else
Azimuth = AzimuthRL(ptFrom, ptTo)
End If
End Function

Public Function Distance(ptFrom As POSITION, ptTo As POSITION, Optional bType As bTypes = RhumbLine) As Double
If bType = GreatCircle Then
Distance = DistanceGC(ptFrom, ptTo)
Else
Distance = DistanceRL(ptFrom, ptTo)
End If
End Function

' Great Circle staff
Private Function AzimuthGC(ptFrom As POSITION, ptTo As POSITION) As Double
AzimuthGC = NormalizeAngle(Deg(Atan2(Sind(ptTo.lon - ptFrom.lon), Tand(ptTo.lat) * Cosd(ptFrom.lat) - Sind(ptFrom.lat) * Cosd(ptTo.lon - ptFrom.lon))))
End Function

Private Function DistanceGC(ptFrom As POSITION, ptTo As POSITION) As Double
DistanceGC = Abs(Arccosd(Sind(ptFrom.lat) * Sind(ptTo.lat) + Cosd(ptFrom.lat) * Cosd(ptTo.lat) * Cosd(ptTo.lon - ptFrom.lon)) * 60)
End Function

' Rhumb Line staff
Private Function AzimuthRL(ptFrom As POSITION, ptTo As POSITION) As Double
Dim dL As Double, dR As Double
dL = ptTo.lon - ptFrom.lon
If dL > 180 Then
dL = dL - 360
ElseIf dL < -180 Then
dL = dL + 360
End If
dR = MCH(ptTo.lat) - MCH(ptFrom.lat)
AzimuthRL = NormalizeAngle(Deg(Atan2(dL * 60, dR)))
End Function

Private Function DistanceRL(ptFrom As POSITION, ptTo As POSITION) As Double
Dim dF As Double, dL As Double, dW As Double, fAver As Double
dF = ptTo.lat - ptFrom.lat
dL = ptTo.lon - ptFrom.lon
fAver = (ptTo.lat + ptFrom.lat) / 2
dW = dL * Cosd(fAver)
DistanceRL = Sqr(dF * dF + dW * dW)
End Function

Private Function Rad(ByVal x As Double) As Double
Rad = x * PI / 180
End Function

Private Function Tand(ByVal x As Double) As Double
End Function

Private Function Deg(ByVal x As Double) As Double
Deg = x * 180 / PI
End Function

'Sin with argument in degrees
Private Function Sind(ByVal x As Double) As Double
End Function

'Cos with argument in degrees
Private Function Cosd(ByVal x As Double) As Double
End Function

Private Function Atan2(ByVal Y As Double, ByVal x As Double) As Double
If x = 0 Then
Atan2 = IIf(Y < 0, 3 * PI / 2, PI / 2)
Else
Atan2 = Atn(Y / x)
If x < 0 Then
If Y < 0 Then Atan2 = Atan2 + PI Else Atan2 = Atan2 - PI
End If
End If
End Function

'ArcCos in degrees
Private Function Arccosd(ByVal x As Double) As Double
On Error Resume Next
If x = 1 Then
Arccosd = 0
ElseIf x = -1 Then
Arccosd = 180
Else
Arccosd = Deg(Atan2(-x, Sqr(1 - x * x))) + 90
End If
On Error GoTo 0
End Function

'ArcSin in degrees
Private Function Arcsind(ByVal x As Double) As Double
If x = 1 Then
Arcsind = 90
ElseIf x = -1 Then
Arcsind = 270
Else
Arcsind = Deg(Atan2(x, Sqr(-x * x + 1)))
End If
End Function

Private Function NormalizeAngle(x As Double) As Double
NormalizeAngle = x
If x < 0 Then NormalizeAngle = x + 360
If x > 360 Then NormalizeAngle = x - 360
End Function

'Meridional parts
Private Function MCH(ByVal x As Double) As Double
Dim a As Double
a = 45 * 60 / Atn(1)
On Error Resume Next
MCH = a * Log(Tan(Rad(45) + Rad(x) / 2) * ((1 - e * Sind(x)) / (1 + e * Sind(x))) ^ (e / 2))
End Function

'=====Using=========
Private Sub Command1_Click()
Dim ptCurrent As POSITION, ptDestination As POSITION
'Glasgow
ptCurrent.lat = 55 + 50 / 60
ptCurrent.lon = -1 * (4 + 15 / 60)
'Heathrow
ptDestination.lat = 51 + 28 / 60 + 11 / 3600
ptDestination.lon = -1 * (0 + 27 / 60 + 5 / 3600)
'Calculation
Dim A_GC As Double, A_RL As Double
Dim D_GC As Double, D_RL As Double
A_GC = Azimuth(ptCurrent, ptDestination, GreatCircle)
A_RL = Azimuth(ptCurrent, ptDestination, RhumbLine)
D_GC = Distance(ptCurrent, ptDestination, GreatCircle)
D_RL = Distance(ptCurrent, ptDestination, RhumbLine)
Dim sMsg As String
sMsg = "From Glasgow to Heathrow:" & vbCrLf & "   Rhumb Line:" & vbCrLf
sMsg = sMsg & "Azimuth: " & Format(A_RL, "000.00 deg.") & vbCrLf & "Distance: " & Format(D_RL, "0000.00 nm.") & vbCrLf
sMsg = sMsg & "    Great Gircle: " & vbCrLf
sMsg = sMsg & "Azimuth: " & Format(A_GC, "000.00 deg.") & vbCrLf & "Distance: " & Format(D_GC, "0000.00 nm.")
MsgBox sMsg, vbInformation
End Sub
0

LVL 28

Expert Comment

ID: 17927757
Oops, sorry,
DistanceRL = Sqr(dF * dF + dW * dW) 'is in degrees
'Should be
DistanceRL = 60 * Sqr(dF * dF + dW * dW) 'in nautical miles

0

LVL 28

Expert Comment

ID: 17927838
PS
You can check it here: http://www.movable-type.co.uk/scripts/LatLong.html
Note that small difference in Rhumb Line calculations between above page and my code depends on spheroidical Earth calculation which I used(unlike above page used sphere Earth)
0

Author Comment

ID: 17928673
Danaseaman and Ark,

Thanks for your efforts so far, however, there are still some problems somewhere and I'm confused.

Danaseaman, your code appears to work given some co-ordinates but is quite obvious it's failing on others given the result.

Ark, your code appears to be giving more accurate results until checking against the online calculator for which you supplied the link.

Example, using the co-ordinates mention above (Glasgow and Heathrow)

Danaseaman's code returns 136 degrees (after reversing the Dest/Cur co-ords)

Ark's code returns 152 degrees

The OnLine calculator returns 209 degrees.

I'll be more than happy to accept a couple of degree difference but the figures shown are a little to wide of the mark.

I will also be happy to split the points between you both should you both come back with a suitable answer

Thank you for your efforts so far

Dave
0

LVL 22

Expert Comment

ID: 17931743
One problem in testing is that longitudes you supplied for Glasgow and Heathrow should be negative.

'Glasgow
'55.870503  -4.445725

'Heathrow
'Lat 51.4833 Long -0.4500

Using online calculator with above values I get distance of 540.1 Km and bearing of 149Â°09&#8242;05&#8243;

Arks code with same values gives 151.62 deg. Rhumb Line or 150.03 deg. Great Circle
Distance looks close too, 299 nm(556 km)

I would go with his Arks solution.

0

LVL 28

Expert Comment

ID: 17934705
Thanks for points, glad I could help
Just checked for
'Glasgow
'55.870503  -4.445725

'Heathrow
'Lat 51.4833 Long -0.4500

Great Circle
=================
On-line calculator
A = 150Â°01&#8242;38&#8243;
D = 554.1 km

My program:
A = 150,027284473195 = 150Â°01&#8242;38.22&#8243;
D = 298,978134462858 nautical miles (ie. minutes of meridian) Equation 1 nm = 1852 m is average for meridional radius from equator to pole. For 6371 km radius (as given in on-line calculator page formulas) equation is 1 nm = 1853.24877740931 (PI * 6371000 / 180 / 60). In this case D = 554,080862165407 kM. If you need output in kilometers, you can use following function:

'Calculate distance via Great Circle in kilometers
Public Function DistanceGC_KM(ptFrom As POSITION, ptTo As POSITION) As Double
DistanceGC_KM = Rad(Abs(Arccosd(Sind(ptFrom.lat) * Sind(ptTo.lat) + Cosd(ptFrom.lat) * Cosd(ptTo.lat) * Cosd(ptTo.lon - ptFrom.lon)))) * 6371
End Function

'======================
Rhumb Line

On-line calculator:
A = 151Â°40&#8242;53&#8243;
D = 554.2 km

My program:
A = 151,624973888495 = 151Â°37&#8242;39.2&#8243;
D = 299,094954464381 = 554,297358690407 kM
Note, that for Rhumb Line Azimuth calculation I used Earth excentricity (more precessious), but for Distance I used approximate formula. More correct is to use already calculated azimuth:
D = dF/Cos(Az)

Private Function DistanceRL(ptFrom As POSITION, ptTo As POSITION) As Double
Dim dF As Double, dL As Double, Az As Double
dF = ptTo.lat - ptFrom.lat
dL = ptTo.lon - ptFrom.lon
If dL > 180 Then
dL = dL - 360
ElseIf dL < -180 Then
dL = dL + 360
End If
dR = MCH(ptTo.lat) - MCH(ptFrom.lat)
Az = Atan2(dL * 60, dR)
DistanceRL = dF / Cos(Az) * 60
End Function
0

## Featured Post

Question has a verified solution.

If you are experiencing a similar issue, please ask a related question

Introduction While answering a recent question (http://www.experts-exchange.com/Q_27402310.html) in the VB classic zone, I wrote some VB code in the (Office) VBA environment, rather than fire up my older PC.  I didn't post completely correct code oâ€¦
Have you ever wanted to restrict the users input in a textbox to numbers, and while doing that make sure that they can't 'cheat' by pasting in non-numeric text? Of course you can do that with code you write yourself but it's tedious and error-prone â€¦
Get people started with the utilization of class modules. Class modules can be a powerful tool in Microsoft Access. They allow you to create self-contained objects that encapsulate functionality. They can easily hide the complexity of a process fromâ€¦
This lesson covers basic error handling code in Microsoft Excel using VBA. This is the first lesson in a 3-part series that uses code to loop through an Excel spreadsheet in VBA and then fix errors, taking advantage of error handling code. This lâ€¦
###### Suggested Courses
Course of the Month6 days, 15 hours left to enroll