The Museum of HP Calculators
The Museum of HP Calculators
This program is Copyright © 2005 by Jean-Marc Baillard and is used here by permission.
This program is supplied without representation or warranty of any kind. Jean-Marc Baillard and The Museum of HP Calculators therefore assume no responsibility and shall have no liability, consequential or otherwise, of any kind arising from the use of this program material or any part thereof.
-A loxodromic curve is a line on the surface of the Earth which cuts
all meridians at the same angle µ ( the azimuth )
-The 2 following programs compute D = the length of these curves and
the azimuth µ
1°) Spherical Model of the Earth
2°) Ellipsoidal model of the Earth
-The azimuths are measured positively clockwise from North
-And the longitudes are measured positively Eastwards from the meridian
of Greenwich
1°) Spherical Model of the Earth
Formulae: L2
- L1 = ( Tan µ ). Ln [ ( Tan ( 45° + b2/2
) ) / (Tan ( 45° + b1/2 ) ) ]
where Li = longitudes ; bi = latitudes
D = R.( b2 - b1 )/Cos µ
R = mean radius of the Earth = 6371 km
-If cos µ = 0 we have D = R.( L2
- L1 ).cos b where b = b1
= b2
Data Registers: /
Flags: /
Subroutines: /
-Synthetic registers M & N may be replaced by any standard registers.
01 LBL "LOX"
02 DEG
03 HR
04 STO M
05 STO N
06 X<> T
07 HMS-
08 HR
09 360
10 MOD
11 PI
12 R-D
13 X>Y?
14 CLX
15 ST+ X
16 -
17 X<>Y
18 HR
19 ST- M
20 2
21 ST/ T
22 /
23 45
24 ST+ T
25 +
26 TAN
27 R^
28 TAN
29 X<>Y
30 /
31 LN
32 R-D
33 R-P
34 X<>Y
35 HMS
36 RCL M
37 LASTX
38 COS
39 X#0?
40 GTO 00
41 +
42 X<> N
43 COS
44 R^
45 *
46 1
47 LBL 00
48 /
49 ABS
50 D-R
51 6371
52 *
53 CLA
54 END
( 78 bytes / SIZE 000 )
| STACK | INPUTS | OUTPUTS |
| T | L1 ( ° ' " ) | / |
| Z | b1 ( ° ' " ) | / |
| Y | L2 ( ° ' " ) | µ ( ° ' " ) |
| X | b2 ( ° ' " ) | D ( km ) |
Example1:
U.S. Naval Observatory at Washington (D.C.)
L1 = 77°03'56"0 W = -77°03'56"0 ;
b1 = 38°55'17"2 N = +38°55'17"2
The Observatoire de Paris:
L2 = 2°20'13"8 E = +2°20'13"8
; b2 = 48°50'11"2 N = +48°50'11"2
-77.03560 ENTER^
38.55172 ENTER^
2.20138 ENTER^
48.50112 XEQ "LOX" >>>>
D = 6436.5499 km X<>Y
µ = 80°08'14" ( in 4 seconds )
Example2:
0 ENTER^
30 ENTER^
120 ENTER^
30 R/S
>>>> D = 11555.7158 km X<>Y
µ = 90°
2°) Ellipsoidal Model of the Earth
-Higher accuracy is obtained if we take the Earth's flattening into account.
-Equatorial radius of the Earth = a = 6378.137 km
-Eccentricity of the ellipsoid = e = ( 0.006694385 )1/2
which corresponds to a flattening f = 1/298.257
( e2 = 2.f - f 2 )
Formulae:
L2 - L1 = ( Tan µ
).{ Ln [ ( Tan ( 45° + b2/2 ) ) / (Tan ( 45° +
b1/2 ) ) ] - (e/2).Ln [ ( 1 + e.sin b2 ).( 1 - e.sin
b1 )/( 1 - e.sin b2 )/( 1 + e.sin b1 )
] }
D = [ a.( 1 - e2 )/ cos µ ] §b1b2
( 1 - e2 sin2 t ) -3/2 dt
( § = Integral )
-One could use any integration formula but "LOX2" evaluates this integral
by a series expansion up to the terms in e6
( the first neglected term is proportional to e8
)
D
~ [ a.( 1 - e2 )/ cos µ ] [ ( 1 + 3e2/4 + 45e4/64
+ 175e6/256 ).( b2 - b1 ) - ( 3e2/8
+ 15e4/32 + 525e6/1024 ).( sin 2b2 - sin
2b1 )
+ ( 15e4/256 + 105e6/1024 ).( sin 4b2
- sin 4b1 ) - ( 35e6/3072 ).( sin 6b2
- sin 6b1 ) ]
-However, we cannot use this formula if cos µ = 0 , in this case ( i-e if b1 = b2 ) we have:
D = a.(
cos b ).( L2 - L1 ) / ( 1 - e2 sin2
b ) -1/2 the loxodromic line is the parallel
of latitude b = b1 = b2
Data Registers: R00 = e = 0.0066943851/2
; R01 = b1 ; R02 = b2 ;
R03 = L2 - L1 ; R04 = µ
( in degrees )
Flags: /
Subroutines: /
01 LBL "LOX2"
02 DEG
03 HR
04 STO 02
05 X<> T
06 HMS-
07 HR
08 360
09 MOD
10 PI
11 R-D
12 X>Y?
13 CLX
14 ST+ X
15 -
16 STO 03
17 X<>Y
18 HR
19 STO 01
20 2
21 ST/ T
22 /
23 45
24 ST+ T
25 +
26 TAN
27 R^
28 TAN
29 /
30 LN
31 STO 04
32 1
33 RCL 01
34 SIN
35 .006694385
36 SQRT
37 STO 00
38 *
39 ST+ Y
40 1
41 -
42 /
43 1
44 RCL 00
45 RCL 02
46 SIN
47 *
48 ST- Y
49 1
50 +
51 /
52 *
53 CHS
54 SQRT
55 LN
56 RCL 00
57 *
58 RCL 04
59 -
60 RCL 03
61 X<>Y
62 R-D
63 R-P
64 X<>Y
65 STO 04
66 RCL 01
67 6
68 *
69 SIN
70 RCL 02
71 6
72 *
73 SIN
74 -
75 3
76 *
77 RCL 02
78 4
79 *
80 SIN
81 RCL 01
82 4
83 *
84 SIN
85 -
86 2657
87 *
88 +
89 RCL 02
90 ST+ X
91 SIN
92 RCL 01
93 ST+ X
94 SIN
95 -
96 2531555
97 *
98 -
99 RCL 02
100 RCL 01
101 -
102 D-R
103 1005052504
104 *
105 +
106 E9
107 /
108 1
109 RCL 00
110 X^2
111 -
112 *
113 RCL 04
114 COS
115 X#0?
116 GTO 00
117 RCL 03
118 D-R
119 RCL 01
120 COS
121 *
122 1
123 RCL 01
124 SIN
125 RCL 00
126 *
127 X^2
128 -
129 SQRT
130 LBL 00
131 /
132 RCL 04
133 HMS
134 X<>Y
135 ABS
136 6378.137
137 *
138 END
( 194 bytes / SIZE 005 )
| STACK | INPUTS | OUTPUTS |
| T | L1 ( ° ' " ) | / |
| Z | b1 ( ° ' " ) | / |
| Y | L2 ( ° ' " ) | µ ( ° ' " ) |
| X | b2 ( ° ' " ) | D ( km ) |
Example1:
-77.03560 ENTER^
38.55172 ENTER^
2.20138 ENTER^
48.50112 XEQ "LOX2" >>>>
D = 6453.389979 km X<>Y
µ = 80°10'15"31 ( in 12 seconds )
the exact value is D = 6453.389986 km
-Compare this value with the length of the corresponding geodesic:
6181.621794 km
-Rhumb-sailing is easier but slower than orthodromy.
Example2:
0 ENTER^
30 ENTER^
120 ENTER^
30 R/S
>>>> D = 11578.35364 km X<>Y
µ = 90°
-Due to roundoff errors, there is a loss of significant digits if the
2 latitudes are very close but not equal.
-For instance: ( 0 ; 30° ) ( 120 ; 30°00'01" )
gives D = 11578.956 km instead of 11578.340
( µ = 89°59'59"45 is correctly computed, however
)
-Otherwise, the accuracy is of the order of a few centimeters.
-If you don't want to calculate D , simply replace lines 66 to 137
by a single HMS.
-In order to avoid a DATA ERROR & an OUT OF RANGE if the
latitudes are +/- 90°
key in 89°59'59"9999 instead of 90° , the difference
is negligible.
Go back to the HP-41 software library
Go back to the general software library
Go
back to the main exhibit hall