'This is to calculate distance and bearing
'written by John Drew VK5DJ
'based on P149 of Amateur Radio Software by Morris published by RSGB
'It takes into account errors in latitude caused by elipticity of the earth and
'distance error caused by the ellipsoid nature of the earth.
'Accuracy likely to be within 0.2%, height not taken into account.
'To estimate distance/bearing for most radio applications this is about as accurate as you
'will get as terrain and propogation effects are a significant "noise" on accuracy
;-------------------------------------------------------------------------------
;**** Added by Fuse Configurator ****
; Use the Fuse Configurator plug-in to change these settings
Device = 18F452
Config_Start
OSC = HS ;HS oscillator
OSCS = OFF ;Oscillator system clock switch option is disabled (main oscillator is source)
PWRT = OFF ;PWRT disabled
BOR = On ;Brown-out Reset enabled
BORV = 20 ;VBOR set to 2.0V
WDT = OFF ;WDT disabled (control is placed on the SWDTEN bit)
WDTPS = 128 ;1:128
CCP2MUX = OFF ;CCP2 input/output is multiplexed with RB3
STVR = OFF ;Stack Full/Underflow will not cause RESET
LVP = OFF ;Low Voltage ICSP disabled
Debug = OFF ;Background Debugger disabled. RB6 and RB7 configured as general purpose I/O pins.
CP0 = OFF ;Block 0 (000200-001FFFh) not code protected
CP1 = OFF ;Block 1 (002000-003FFFh) not code protected
CP2 = OFF ;Block 2 (004000-005FFFh) not code protected
CP3 = OFF ;Block 3 (006000-007FFFh) not code protected
CPB = OFF ;Boot Block (000000-0001FFh) not code protected
CPD = OFF ;Data EEPROM not code protected
WRT0 = OFF ;Block 0 (000200-001FFFh) not write protected
WRT1 = OFF ;Block 1 (002000-003FFFh) not write protected
WRT2 = OFF ;Block 2 (004000-005FFFh) not write protected
WRT3 = OFF ;Block 3 (006000-007FFFh) not write protected
WRTC = OFF ;Configuration registers (300000-3000FFh) not write protected
WRTB = OFF ;Boot Block (000000-0001FFh) not write protected
WRTD = OFF ;Data EEPROM not write protected
EBTR0 = OFF ;Block 0 (000200-001FFFh) not protected from Table Reads executed in other blocks
EBTR1 = OFF ;Block 1 (002000-003FFFh) not protected from Table Reads executed in other blocks
EBTR2 = OFF ;Block 2 (004000-005FFFh) not protected from Table Reads executed in other blocks
EBTR3 = OFF ;Block 3 (006000-007FFFh) not protected from Table Reads executed in other blocks
EBTRB = OFF ;Boot Block (000000-0001FFh) not protected from Table Reads executed in other blocks
Config_End
;**** End of Fuse Configurator Settings ****
;-------------------------------------------------------------------------------
All_Digital =1
Declare Xtal 10
Declare LCD_Type 0 'text type
Declare LCD_DTPin PORTD.4 'assigns data lines to B4..7
Declare LCD_ENPin PORTE.1 'enable pin
Declare LCD_RSPin PORTE.0 'RS line pin
Declare LCD_Interface 4 '4 or 8 line interface
Declare LCD_Lines 2 'lines in the display
Declare Float_Rounding = On
Declare Float_Display_Type = LARGE
'Define some constants
Symbol PI = 3.14159265 'PI
Symbol PIby2 = PI / 2 'pi dividend by 2
Symbol TwoPI = PI * 2 'pi mulitplied by 2
Symbol DegToRad = PI/180 'multiply by this to convert degrees to radians
Symbol RadToDeg = 180/PI 'multiply by this to convert radians to degrees
Symbol MinRadiusEarth = 6356 'polar radius or semi minor axis
Symbol MaxRadiusEarth = 6378 'equatorial radius or semi major axis
Symbol MinbyMax = 0.9965506 'ratio of min radius / max radius
Symbol MeanRadius = 6367 '(Min + Max)/2
Symbol Eccentricity = 0.9931131 'measure of the ellipsoid minradius squared/maxradius squared
Dim HomeLong As Float 'set to home longitude in degrees -ve for West, +ve for East (input)
Dim HomeLat As Float 'set to home latitude in degrees -ve for South, +ve for North (input)
Dim RemoteLong As Float 'remote site long in degrees -ve for West, +ve for East (input)
Dim RemoteLat As Float 'remote site lat in degrees -ve for South, +ve for South (input
Dim Bearing As Float 'calculated bearing of the remote location from the home location (output)
Dim EccFactor As Float
Dim Azimuth As Float
Dim Temp1 As Float 'temp1,2,3 used as a means of overcoming any problems with trig functions
Dim Temp2 As Float 'in complex lines
Dim Temp3 As Float
Dim CO As Float
Dim SI As Float
Dim N1 As Float
Dim N2 As Float
Dim CA As Float
Dim AL As Float
Dim EP As Float
Dim Dist As Float
Dim T1 As Float
Dim T2 As Float
Dim TH As Float 'used as a parameter for Atan2 routine
Dim A2 As Float
Dim A1 As Float
Dim XI As Float
Dim Arctan2Flag As Bit
Dim Distance As Word 'distance between the two locations on a great circle (output)
GoTo Initialise 'jump over the sub routines
ArcTan2:
If Arctan2Flag=0 Then 'assign SI and CO while calculating central angle
SI = EccFactor * Sin(TH) 'TH is used to pass values to and from this routine
CO = Cos(TH)
EndIf
T1 = Abs(SI): T2 = Abs(CO) 'arctan2 routine
If T1<>0 Or T2<>0 Then
If T1>T2 Then
Temp1 = ATan (T2/T1)
TH = PIby2 - Temp1
EndIf
If T1<=T2 Then TH=ATan(T1/T2)
If CO<0 Then TH = PI - TH
If SI<0 Then TH = -TH
Else
TH=PIby2
EndIf
Return
LongLine1: 'just to separate out a long calculation from main code
Temp1 = DegToRad * (HomeLong - RemoteLong)
Temp1 = Cos(Temp1)
Temp2 = Cos(N1)
Temp3 = Cos(N2)
Temp1 = Temp1 * Temp2 * Temp3
Temp2 = Sin(N1)
Temp3 = Sin(N2)
Temp2 = Temp2 * Temp3
CO = Temp1 + Temp2
Temp1 = CO * CO
SI = Sqr (1- Temp1)
Return
LongLine2: 'just to separate out a long calculation from main code
Temp1 = DegToRad * (RemoteLong - HomeLong)
Temp1 = Sin (Temp1)
Temp2 = Cos (N2)
Temp3 = Cos (N1)
SI = Temp1 * Temp2 * Temp3
Temp1 = Sin (N1)
Temp2 = Cos (CA)
Temp1 = Temp1 * Temp2
Temp2 = Sin (N2)
CO = Temp2 - Temp1
Return
LongLine3: 'just to separate out a long calculation from main code
Temp1 = 1 - (SI * SI)
Temp2 = Abs Temp1
Temp3 = Sqr Temp2
CO = Temp3 * MinbyMax
Return
LongLine4: 'just to separate out a long calculation from main code
Temp1 = 1 - Eccentricity
Temp2 = Sin(AL)
EP = Temp2 * Temp2 * Temp1
Temp1 = Sqr(1-EP)
EccFactor = 1 / Temp1
Return
GetDistance:
Temp1 = A2 - A1
Temp2 = 1- EP/4
Temp1 = Temp1 * Temp2
Temp2 = Sin (A2 + A2)
Temp3 = Sin (A1 + A1)
Temp2 = Temp2 - Temp3
Temp2 = Temp2 * EP / 8
Dist = Temp1 - Temp2
Return
Initialise:
HomeLat = 51.59167 'home lat in decimal degrees
HomeLong = -1.246667 'home long in decimal degrees
RemoteLat = -50.00 'DX lat in decimal degrees
RemoteLong = 89.38333 'DX long in decimal degrees
Main:
TH = DegToRad * HomeLat
EccFactor = Eccentricity
Arctan2Flag=0 'get arctan2
GoSub ArcTan2
N1=TH 'adjusted latitude angle home
TH=DegToRad * RemoteLat
GoSub ArcTan2 'get atan2
N2=TH 'adjusted angle remote
GoSub LongLine1
Arctan2Flag=1
GoSub ArcTan2
CA=TH
GoSub LongLine2
GoSub ArcTan2
Azimuth=TH 'bearing of dist object from home in radians
If Azimuth<0 Then Azimuth=Azimuth+ TwoPI
If N1>=N2 Then
T1=N1 : N1=N2: N2=T1
EndIf
Temp1 = Sin (N1)
Temp2 = Sin (CA)
SI=Temp1 * Temp2
Temp1 = Cos(CA)
Temp2 = Sin(N1)
Temp1 = Temp1 * Temp2
Temp2 = Sin(N2)
CO=Temp2 - Temp1
GoSub ArcTan2
XI = TH : T1 = Sin(XI): T2 = Sin(XI+CA)
If Abs(T1)>Abs(T2) Then
SI = Sin(N1)/T1
Else
SI = Sin(N2)/T2
EndIf
GoSub LongLine3
GoSub ArcTan2
AL=TH
GoSub LongLine4
TH = XI
Arctan2Flag=0
GoSub ArcTan2
A1=TH
TH=XI+CA
GoSub ArcTan2
A2=TH
GoSub GetDistance
Distance = MaxRadiusEarth * Dist 'distance - remove numbers after decimal point
Bearing = Azimuth * RadToDeg 'bearing - convert to degrees
Print At 1,1, Dec5 Distance," km" 'using values in initialise the distance = 14122 km
Print At 2,1, Dec1 Bearing," deg" 'figures as above gives azimuth of 126.1 degrees
End;