Terrestrial Geodesic Distance(2) for the HP-41
Overview
0°) Method
1°) With Andoyer's formula of order
2
2°) With Sodano's Formula - order
3
a) Sodano's Formula ( Ellipsoid
of Revolution )
b) Triaxial Ellipsoid
3°) With Andoyer's formula of order 3
a) Andoyer's Formula ( Ellipsoid
of Revolution )
b) Triaxial Ellipsoid
4°) With Vincenty's Formula
5°) Minimizing the Sum of 2 Distances
6°) A Slightly Different Method
a) With Andoyer's formula of order
2
b) With Andoyer's formula of order
3
c) Almost With Andoyer's formula
of order 3
d) Geodetic or Geocentric Longitudes
& Latitudes ?
-These programs evaluate the geodesic distance between 2 points M(L,B)
& N(L',B') on a triaxial ellipsoid model of the Earth.
-The latitudes L and the longitudes B are supposed geodetic.
a = 6378.172 km
-The semi-axes are: b = 6378.102 km
and the longitude of the equatorial major-axis = 14°92911 West
c = 6356.752314 km
-So, the mean radius of the equator and the flattening correspond to
WGS84.
0°) Method
-First, the distances of the 2 points M & N are computed on the great ellipses passing by M & N:
-For the triaxial model, it gives, say A
-For the WGS84 ellipsoid, it gives A'
-Then, the geodesic distance d is calculated on the WGS84 ellipsoid of revolution.
-Since the flattening of the Earth equator is very small, the geodesic distance D on the triaxial ellipsoid is very well approximated by the simple relation:
D = (A/A') d
1°) With Andoyer's Formula of order 2
"TGD" determines the elements of the great ellipse (E) passing by M & N as follows:
With OM / OM = ( x , y , z ) & ON / ON = ( x' , y' ,z' ) ( the unit vectors defined by OM & ON ) we have
OP = ( Sin W ) / ( A2 + B2 + C2 )1/2 for a point P on the ellipsoid and with µ = ( OM , OP )
A = [ x Sin ( W - µ ) + x' Sin
µ ] / a
B = [ y Sin ( W - µ ) + y' Sin
µ ] / b
W = the angle ( OM , ON )
C = [ z Sin ( W - µ ) + z' Sin
µ ] / c
N P
/ /
/ /
/ / µ
O//-----------------------------M
>>> The value of µ for the axes of the great ellipse (E) is found by equating to zero the derivative of OP with respect to µ
-After some calculations, it yields: ( Sin 2µ ) / ( Cos 2µ ) = P / Q with
P = -2 [ ( x Sin W ) ( x Cos W - x' )
/ a2 + ( y Sin W ) ( y Cos W - y' ) / b2 + ( z Sin
W ) ( z Cos W - z' ) / c2 ]
Q = [ x2 Sin2 W - (
x Cos W - x' )2 ] / a2 + [ y2 Sin2
W - ( y Cos W - y' )2 ] / b2 + [ z2 Sin2
W - ( z Cos W - z' )2 ] / c2
>>> However, it's simpler & faster to compute
A2 + B2 + C2 = T Sin2 ( W - µ ) + S Sin2 µ + 2 U Sin µ Sin ( W - µ )
where S = x'2 / a2 + y'2 / b2 + z'2 / c2 , T = x2 / a2 + y2 / b2 + z2 / c2 , U = xx' / a2 + yy' / b2 + zz' / c2
and P = 2 U - 2 T Cos W ,
Q = T Sin W + ( 2 U Cos W - S - T Cos2 W ) / Sin W
-Then, the geodesic distances - on the great ellipses corresponding
to WGS84 & the triaxial ellipsoid - between M & N ( say A' &
A )
are computed by Andoyer's formula on order 2.
-Finally, the geodesic distance d on the WGS84 ellipsoid is calculated
by 2nd order Andoyer's formula
and the geodesic distance D on the triaxial ellipsoid is approximated
by D = (A/A') d
Andoyer's Formula of order 2:
-With H = ( L' - L )/2 ; F = ( B + B' )/2 ; G = ( B' - B )/2
S = sin2 G cos2 H + cos2 F sin2 H ; µ = Arcsin S1/2 ; R = ( S.(1-S) )1/2
we define A = ( sin2 G cos2
F )/S + ( sin2 F cos2 G )/( 1-S )
and B
= ( sin2 G cos2 F )/S - ( sin2 F
cos2 G )/( 1-S )
-The geodesic distance on the biaxial ellipsoid is then: d = 2a.µ ( 1 + eps )
where eps = (f/2) ( -A - 3B R/µ
)
+ (f2/4) { { - ( µ/R + 6R/µ
).B + [ -15/4 + ( 2S - 1 ) ( µ/R + (15/4) R/µ ) ] A + 4 - (
2S - 1 ) µ/R } A
- [ (15/2) ( 2S - 1 ).B R/µ - ( µ/R + 6R/µ ) ]
B
}
and f = flattening = 1/298.25722
Data Registers: R00 to
R15: temp
Flags: /
Subroutines: /
-Line 38 is a three-byte GTO 14
01 LBL "TGD"
02 DEG 03 HR 04 STO 15 05 RDN 06 HR 07 14.929 08 STO 00 09 + 10 STO 14 11 RDN 12 HR 13 STO 13 14 X<>Y 15 HR 16 RCL 00 17 + 18 STO 12 19 6378.172 20 STO 01 21 .07 22 - 23 STO 02 24 21.349686 25 - 26 STO 03 27 XEQ 01 28 STO 00 29 RCL 01 30 RCL 02 31 + 32 2 33 / 34 STO 01 35 STO 02 36 XEQ 01 37 ST/ 00 38 GTO 14 39 LBL 01 40 RCL 12 41 RCL 13 42 XEQ 02 43 RCL 08 44 STO 05 45 RCL 09 46 STO 06 47 RCL 10 48 STO 07 49 RCL 14 50 RCL 15 51 XEQ 02 52 RCL 05 53 RCL 08 54 * 55 RCL 06 |
56 RCL 09
57 * 58 RCL 07 59 RCL 10 60 * 61 + 62 + 63 STO 04 64 RCL 01 65 ST/ 05 66 ST/ 08 67 RCL 02 68 ST/ 06 69 ST/ 09 70 RCL 03 71 ST/ 07 72 ST/ 10 73 RCL 05 74 X^2 75 RCL 06 76 X^2 77 RCL 07 78 X^2 79 + 80 + 81 STO 11 82 RCL 08 83 X^2 84 RCL 09 85 X^2 86 RCL 10 87 X^2 88 + 89 + 90 X<> 10 91 RCL 07 92 * 93 RCL 06 94 RCL 09 95 * 96 RCL 05 97 RCL 08 98 * 99 + 100 + 101 STO 09 102 ST+ 09 103 RCL 04 104 RCL 11 105 * 106 - 107 ST+ X 108 RCL 04 109 RCL 09 110 * |
111 RCL 10
112 - 113 RCL 11 114 RCL 04 115 X^2 116 * 117 - 118 RCL 11 119 RCL 04 120 ACOS 121 STO 04 122 SIN 123 STO 08 124 ST/ Z 125 * 126 + 127 X#0? 128 / 129 ATAN 130 2 131 / 132 90 133 STO 06 134 X<>Y 135 X<0? 136 + 137 STO 05 138 XEQ 03 139 X<> 06 140 RCL 05 141 + 142 XEQ 03 143 STO 07 144 RCL 05 145 CHS 146 ST+ 04 147 1 148 P-R 149 RCL 07 150 RCL 06 151 / 152 X^2 153 STO 09 154 * 155 R-P 156 X<> 04 157 1 158 P-R 159 RCL 09 160 * 161 R-P 162 RDN 163 ST- Z 164 + 165 COS |
166 STO 05
167 X^2 168 ST+ X 169 X<>Y 170 ABS 171 ENTER 172 D-R 173 STO 10 174 SIGN 175 P-R 176 X<>Y 177 RCL 10 178 / 179 STO 11 180 * 181 ST* Y 182 - 183 15 184 * 185 1 186 + 187 4 188 / 189 RCL 07 190 CHS 191 RCL 06 192 ST+ Y 193 ST+ X 194 / 195 STO 09 196 * 197 RCL 05 198 3 199 * 200 RCL 11 201 * 202 - 203 RCL 09 204 ST* Y 205 - 206 RCL 10 207 ST* Y 208 + 209 RCL 06 210 * 211 RTN 212 LBL 02 213 1 214 P-R 215 X<>Y 216 RCL 03 217 X^2 218 * 219 STO 10 220 X^2 |
221 X<> Z
222 X<>Y 223 P-R 224 RCL 01 225 X^2 226 * 227 STO 08 228 X^2 229 X<>Y 230 RCL 02 231 X^2 232 * 233 STO 09 234 X^2 235 + 236 + 237 SQRT 238 ST/ 08 239 ST/ 09 240 ST/ 10 241 RTN 242 LBL 03 243 ENTER 244 SIN 245 STO 07 246 X^2 247 RCL 10 248 * 249 X<>Y 250 RCL 04 251 - 252 SIN 253 ST* 07 254 X^2 255 RCL 11 256 * 257 + 258 RCL 07 259 RCL 09 260 * 261 - 262 SQRT 263 RCL 08 264 X<>Y 265 / 266 RTN 267 LBL 14 268 RCL 15 269 RCL 13 270 RCL 14 271 RCL 12 272 - 273 2 274 / 275 SIN |
276 X^2
277 RDN 278 + 279 2 280 / 281 ST- Y 282 1 283 P-R 284 X^2 285 STO 08 286 X<>Y 287 X^2 288 STO 07 289 RDN 290 X<>Y 291 1 292 P-R 293 X^2 294 ST* 07 295 RDN 296 X^2 297 ST* 08 298 STO T 299 - 300 * 301 + 302 ST/ 08 303 STO 09 304 SQRT 305 ASIN 306 D-R 307 RCL 09 308 1 309 RCL 09 310 - 311 ST/ 07 312 ST- 09 313 * 314 SQRT 315 / 316 STO 05 317 CLX 318 RCL 08 319 RCL 07 320 ST- 08 321 + 322 STO 07 323 RCL 09 324 RCL 05 325 ST* 09 326 / 327 7.5 328 * 329 STO 10 |
330 LASTX
331 - 332 2 333 / 334 RCL 09 335 + 336 * 337 4 338 + 339 RCL 09 340 - 341 6 342 RCL 05 343 ST/ Y 344 + 345 RCL 08 346 * 347 STO Z 348 - 349 RCL 07 350 * 351 + 352 RCL 08 353 X^2 354 RCL 10 355 * 356 - 357 RCL 03 358 CHS 359 RCL 01 360 ST+ Y 361 ST+ X 362 STO 09 363 / 364 STO 10 365 * 366 RCL 07 367 - 368 RCL 08 369 3 370 * 371 RCL 05 372 / 373 - 374 RCL 10 375 * 376 * 377 + 378 RCL 09 379 * 380 STO Y 381 RCL 00 382 * 383 END |
( 471 bytes / SIZE 016 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | d ( km ) |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial
ellipsoid
and d is the geodesic
distance on the ellipsoid of revolution ( WGS 84 )
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12138.65748 km
---Execution time = 47s---
X<>Y d = 12138.68425 km
-The exact results are D = 12138.65757 & d =
12138.68433 km ( rounded to 5 decimals )
Notes:
-Longitudes are geodetic and positive East.
-If the points are not nearly antipodal, the precision is of the order
of one meter.
-With (L,b,L',b') = (70°,0°,270°,10°) "TGD" gives a distance of 17554.6159 instead of 17554.6140 km
-Here are a few results about the error of the results:
If | L' - L | < 120° rms error ~ 6 cm
and max error ~ 26 cm
If | L' - L | < 140° rms error ~ 7 cm
and max error ~ 38 cm
If | L' - L | < 150° rms error ~ 8 cm
and max error ~ 70 cm
If | L' - L | < 160° rms error ~ 11 cm and
max error ~ 168 cm
-A better accuracy may be obtained by the programs listed in "Geodesics
on a triaxial ellipsoid for the HP-41" but the program above is faster
!
2°) With Sodano's Formula of order 3
a) Sodano's Formula ( Ellipsoid
of Revolution )
-"SOD" uses a 3rd order formula in the flattening to calculate the geodesic distance on an ellipsoid of revolution ( WGS84 ) - cf references [2] & [3]
-I've rewritten the formula as follows, after multiplying both sides by (1-f)
Dist = a { µ + f [ p Sin µ - (m/2) ( µ + Sin µ Cos µ )
+ f2 [ p ((m-1)/2) µ2 Csc µ + m ((1-m)/2) µ2 Cot µ + ((m2-8.p2)/16) Sin µ Cos µ + (m2/16) µ - (m2/8) Sin µ Cos3 µ + (mp/2) Sin µ Cos2 µ ]
+ f3
[ µ2 Csc µ (p/4) (-2+5m-3m2) + µ3
Csc µ Cot µ (p/2) (1-4m+3m2) - µ2
Cot µ (m/4) (-2+5m-3m2) + µ3 (m/6) (-1+2m-m2)
+ µ3 Cot2 µ (m/2) (-1+3m-2m2)
+ Sinµ Cos µ (-8 p2 + m2 - 8 p2m
- m3/2 )/16 + µ (8 p2 + m2
- 8 p2m - m3/2 )/16
+ µ3 Csc2 µ (p2/2) (1-m) +
Sin µ Cos3 µ (m2/16) (-2+m) + µ
Cos2 µ (m2/4) (1-m) + Sin µ Cos2
µ (pm/2) (1+m) + p2m Sin3 µ Cos µ
+ µ Cos µ (3pm/4) (-1+m) + (pm2/2) Sin5
µ + Sin µ (p/2) (p2-m2) + (m3/12)
Sin3 µ Cos3 µ - (2p3/3) Sin3
µ ] }
where µ = 2 Arc Sin [ sin2 ((B'-B)/2) cos2 ((L'-L)/2) + cos2 ((B+B')/2) sin2 ((L'-L)/2) ]1/2 L , L' = longitudes & B , B' = parametric latitudes
Tan B = (1-f) Tan l where l = geodetic latitude
p = Sin B Sin B' m = 1
- [ Cos B Cos B' Sin (L'-L) ]2 / Sin2 µ
, a = 6378.137 km , b = 6356.752314 km
( WGS84 )
Data Registers: R00 thru R13: temp
Flags: /
Subroutines: /
01 LBL "SOD"
02 DEG 03 HR 04 STO 06 05 X<> T 06 HMS- 07 HR 08 STO 07 09 RDN 10 HR 11 298.25722 12 1/X 13 STO 04 14 SIGN 15 LASTX 16 - 17 STO 03 18 P-R 19 LASTX 20 / 21 R-P 22 X<>Y 23 STO 05 24 1 25 P-R 26 STO 01 27 X<>Y 28 STO 02 29 RCL 06 30 RCL 03 31 P-R 32 LASTX 33 / 34 R-P 35 X<>Y 36 STO 06 37 1 38 P-R 39 STO 08 40 X<>Y 41 STO 09 42 RCL 02 43 * 44 STO 03 45 RCL 01 46 ST* 09 47 RCL 08 48 ST* 02 |
49 *
50 STO 10 51 RCL 07 52 2 53 / 54 SIN 55 X^2 56 * 57 RCL 06 58 RCL 05 59 - 60 2 61 / 62 SIN 63 X^2 64 + 65 SQRT 66 ASIN 67 ST+ X 68 D-R 69 STO 13 70 LASTX 71 1 72 P-R 73 STO 11 74 X<>Y 75 STO 12 76 RCL 07 77 SIN 78 RCL 10 79 * 80 X<>Y 81 / 82 X^2 83 STO 10 84 SIGN 85 LASTX 86 - 87 STO 05 88 3 89 * 90 4 91 - 92 RCL 05 93 * 94 STO 02 95 RCL 03 96 ST* Y |
97 +
98 RCL 11 99 RCL 12 100 / 101 STO 08 102 * 103 RCL 12 104 ST+ X 105 / 106 RCL 10 107 X^2 108 RCL 05 109 * 110 6 111 / 112 - 113 3 114 RCL 05 115 ST+ X 116 - 117 RCL 05 118 X^2 119 STO 00 120 * 121 RCL 05 122 - 123 RCL 08 124 X^2 125 * 126 2 127 / 128 + 129 RCL 10 130 RCL 03 131 X^2 132 * 133 RCL 12 134 X^2 135 STO 09 136 ST+ X 137 / 138 + 139 RCL 13 140 * 141 RCL 02 142 RCL 05 143 - 144 2 |
145 +
146 4 147 / 148 STO 02 149 RCL 05 150 * 151 RCL 08 152 * 153 + 154 RCL 02 155 RCL 03 156 * 157 RCL 12 158 / 159 - 160 RCL 13 161 X^2 162 STO 07 163 * 164 RCL 00 165 2 166 / 167 RCL 03 168 X^2 169 8 170 * 171 STO 06 172 + 173 RCL 05 174 * 175 RCL 00 176 - 177 STO 02 178 RCL 06 179 ST- 02 180 + 181 RCL 11 182 RCL 12 183 * 184 STO 01 185 * 186 RCL 02 187 RCL 13 188 * 189 + 190 16 191 STO 02 192 / |
193 -
194 RCL 10 195 RCL 13 196 * 197 4 198 * 199 RCL 10 200 RCL 01 201 ST* Y 202 + 203 - 204 RCL 00 205 * 206 RCL 11 207 X^2 208 * 209 RCL 02 210 / 211 + 212 RCL 05 213 RCL 01 214 ST* Y 215 + 216 RCL 10 217 RCL 13 218 * 219 1.5 220 STO 08 221 * 222 - 223 2 224 / 225 RCL 03 226 RCL 09 227 * 228 RCL 12 229 * 230 + 231 RCL 03 232 * 233 RCL 05 234 * 235 RCL 11 236 * 237 + 238 RCL 09 239 X^2 |
240 RCL 00
241 ST* Y 242 - 243 RCL 03 244 X^2 245 + 246 RCL 03 247 * 248 RCL 12 249 * 250 RCL 01 251 RCL 05 252 * 253 X^2 254 LASTX 255 * 256 6 257 / 258 + 259 2 260 / 261 + 262 RCL 03 263 RCL 12 264 * 265 X^2 266 LASTX 267 * 268 RCL 08 269 / 270 - 271 RCL 04 272 * 273 RCL 05 274 RCL 11 275 * 276 RCL 03 277 - 278 RCL 10 279 * 280 RCL 07 281 * 282 RCL 12 283 ST+ X 284 / 285 + 286 RCL 13 |
287 RCL 01
288 RCL 11 289 X^2 290 * 291 ST+ X 292 - 293 RCL 00 294 ST* Y 295 RCL 06 296 - 297 RCL 01 298 * 299 + 300 RCL 02 301 / 302 + 303 RCL 03 304 RCL 05 305 * 306 RCL 01 307 * 308 RCL 11 309 * 310 2 311 / 312 + 313 RCL 04 314 * 315 RCL 03 316 RCL 12 317 * 318 + 319 RCL 01 320 RCL 13 321 + 322 RCL 05 323 * 324 2 325 / 326 - 327 RCL 04 328 * 329 RCL 13 330 + 331 6378.137 332 * 333 END |
( 374 bytes / SIZE 014 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | / |
X | b' ( ° ' " ) | D ( km ) |
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "SOD"
>>>> D = 12138.68432 km
---Execution time = 19s---
Notes:
-If | L' - L | < 120° , the RMS error is about 1.4
mm and the maximum error about 9 mm
-If | L' - L | < 170° , -----------------------
12 mm ----------------------------- 22 cm
-If | L' - L | < 175° , -----------------------
5 cm ----------------------------- 2 m
-When Vincenty's iterative method requires more than 3 iterations, even
a 3rd order formula may produce an error of 2 meters:
-For example, with Long1 = 0° , Lat1 = 2° &
Long2 = 175° , Lat2 = -4°
"SOD" gives 19434.163468 km instead of 19434.161314 km returned by Vincenty's formula !
-These estimations were obtained with free42.
-With an HP41 - which works with 10 digits - the errors may be slightly
larger...
b) Triaxial Ellipsoid
-This variant employs Sodano's formula for the last calculation and
Andoyer's formula of order 2 to estimate (A/A')
Data Registers: R00 thru R17: temp
Flags: /
Subroutines: /
-Line 38 is a three-byte GTO 14
01 LBL "TGD"
02 DEG 03 HR 04 STO 15 05 RDN 06 HR 07 14.92911 08 STO 00 09 + 10 STO 14 11 RDN 12 HR 13 STO 13 14 X<>Y 15 HR 16 RCL 00 17 + 18 STO 12 19 6378.172 20 STO 01 21 .07 22 - 23 STO 02 24 21.349686 25 - 26 STO 03 27 XEQ 01 28 STO 17 29 RCL 01 30 RCL 02 31 + 32 2 33 / 34 STO 01 35 STO 02 36 XEQ 01 37 ST/ 17 38 GTO 14 39 LBL 01 40 RCL 12 41 RCL 13 42 XEQ 02 43 RCL 08 44 STO 05 45 RCL 09 46 STO 06 47 RCL 10 48 STO 07 49 RCL 14 50 RCL 15 51 XEQ 02 52 RCL 05 53 RCL 08 54 * 55 RCL 06 56 RCL 09 57 * 58 RCL 07 59 RCL 10 60 * 61 + 62 + 63 STO 04 64 RCL 01 65 ST/ 05 66 ST/ 08 67 RCL 02 68 ST/ 06 69 ST/ 09 70 RCL 03 71 ST/ 07 72 ST/ 10 73 RCL 05 74 X^2 75 RCL 06 76 X^2 77 RCL 07 78 X^2 79 + 80 + 81 STO 11 82 RCL 08 83 X^2 84 RCL 09 85 X^2 86 RCL 10 |
87 X^2
88 + 89 + 90 X<> 10 91 RCL 07 92 * 93 RCL 06 94 RCL 09 95 * 96 RCL 05 97 RCL 08 98 * 99 + 100 + 101 STO 09 102 ST+ 09 103 RCL 04 104 RCL 11 105 * 106 - 107 ST+ X 108 RCL 04 109 RCL 09 110 * 111 RCL 10 112 - 113 RCL 11 114 RCL 04 115 X^2 116 * 117 - 118 RCL 11 119 RCL 04 120 ACOS 121 STO 04 122 SIN 123 STO 08 124 ST/ Z 125 * 126 + 127 X#0? 128 / 129 ATAN 130 2 131 / 132 90 133 STO 06 134 X<>Y 135 X<0? 136 + 137 STO 05 138 XEQ 03 139 X<> 06 140 RCL 05 141 + 142 XEQ 03 143 STO 07 144 RCL 05 145 CHS 146 ST+ 04 147 1 148 P-R 149 RCL 07 150 RCL 06 151 / 152 X^2 153 STO 09 154 * 155 R-P 156 X<> 04 157 1 158 P-R 159 RCL 09 160 * 161 R-P 162 RDN 163 ST- Z 164 + 165 COS 166 STO 05 167 X^2 168 ST+ X 169 X<>Y 170 ABS 171 ENTER 172 D-R |
173 STO 10
174 SIGN 175 P-R 176 X<>Y 177 RCL 10 178 / 179 STO 11 180 * 181 ST* Y 182 - 183 15 184 * 185 1 186 + 187 4 188 / 189 RCL 07 190 CHS 191 RCL 06 192 ST+ Y 193 ST+ X 194 / 195 STO 09 196 * 197 RCL 05 198 3 199 * 200 RCL 11 201 * 202 - 203 RCL 09 204 ST* Y 205 - 206 RCL 10 207 ST* Y 208 + 209 RCL 06 210 * 211 RTN 212 LBL 02 213 1 214 P-R 215 X<>Y 216 RCL 03 217 X^2 218 * 219 STO 10 220 X^2 221 X<> Z 222 X<>Y 223 P-R 224 RCL 01 225 X^2 226 * 227 STO 08 228 X^2 229 X<>Y 230 RCL 02 231 X^2 232 * 233 STO 09 234 X^2 235 + 236 + 237 SQRT 238 ST/ 08 239 ST/ 09 240 ST/ 10 241 RTN 242 LBL 03 243 ENTER 244 SIN 245 STO 07 246 X^2 247 RCL 10 248 * 249 X<>Y 250 RCL 04 251 - 252 SIN 253 ST* 07 254 X^2 255 RCL 11 256 * 257 + 258 RCL 07 |
259 RCL 09
260 * 261 - 262 SQRT 263 RCL 08 264 X<>Y 265 / 266 RTN 267 LBL 14 268 RCL 14 269 RCL 12 270 - 271 STO 07 272 RCL 13 273 RCL 03 274 CHS 275 RCL 01 276 STO 16 277 ST+ Y 278 / 279 STO 04 280 SIGN 281 LASTX 282 - 283 STO 03 284 P-R 285 LASTX 286 / 287 R-P 288 X<>Y 289 STO 05 290 1 291 P-R 292 STO 01 293 X<>Y 294 STO 02 295 RCL 15 296 RCL 03 297 P-R 298 LASTX 299 / 300 R-P 301 X<>Y 302 STO 06 303 1 304 P-R 305 STO 08 306 X<>Y 307 STO 09 308 RCL 02 309 * 310 STO 03 311 RCL 01 312 ST* 09 313 RCL 08 314 ST* 02 315 * 316 STO 10 317 RCL 07 318 2 319 STO 14 320 / 321 SIN 322 X^2 323 * 324 RCL 06 325 RCL 05 326 - 327 RCL 14 328 / 329 SIN 330 X^2 331 + 332 SQRT 333 ASIN 334 ST+ X 335 D-R 336 STO 13 337 LASTX 338 1 339 P-R 340 STO 11 341 X<>Y 342 STO 12 343 RCL 07 344 SIN |
345 RCL 10
346 * 347 X<>Y 348 / 349 X^2 350 STO 10 351 SIGN 352 LASTX 353 - 354 STO 05 355 3 356 * 357 4 358 - 359 RCL 05 360 * 361 STO 02 362 RCL 03 363 ST* Y 364 + 365 RCL 11 366 RCL 12 367 / 368 STO 08 369 * 370 RCL 12 371 ST+ X 372 / 373 RCL 10 374 X^2 375 RCL 05 376 * 377 6 378 / 379 - 380 3 381 RCL 05 382 ST+ X 383 - 384 RCL 05 385 X^2 386 STO 00 387 * 388 RCL 05 389 - 390 RCL 08 391 X^2 392 * 393 RCL 14 394 / 395 + 396 RCL 10 397 RCL 03 398 X^2 399 STO 15 400 * 401 RCL 12 402 X^2 403 STO 09 404 ST+ X 405 / 406 + 407 RCL 13 408 * 409 RCL 02 410 RCL 05 411 - 412 RCL 14 413 + 414 4 415 / 416 STO 02 417 RCL 05 418 * 419 RCL 08 420 * 421 + 422 RCL 02 423 RCL 03 424 * 425 RCL 12 426 / 427 - 428 RCL 13 429 X^2 430 STO 07 |
431 *
432 RCL 00 433 RCL 14 434 / 435 RCL 15 436 8 437 * 438 STO 06 439 + 440 RCL 05 441 * 442 RCL 00 443 - 444 STO 02 445 RCL 06 446 ST- 02 447 + 448 RCL 11 449 RCL 12 450 * 451 STO 01 452 * 453 RCL 02 454 RCL 13 455 * 456 + 457 16 458 STO 02 459 / 460 - 461 RCL 10 462 RCL 13 463 * 464 4 465 * 466 RCL 10 467 RCL 01 468 ST* Y 469 + 470 - 471 RCL 00 472 * 473 RCL 11 474 X^2 475 * 476 RCL 02 477 / 478 + 479 RCL 05 480 RCL 01 481 ST* Y 482 + 483 RCL 10 484 RCL 13 485 * 486 1.5 487 STO 08 488 * 489 - 490 RCL 14 491 / 492 RCL 03 493 RCL 09 494 * 495 RCL 12 496 * 497 + 498 RCL 03 499 * 500 RCL 05 501 * 502 RCL 11 503 * 504 + 505 RCL 09 506 X^2 507 RCL 00 508 ST* Y 509 - 510 RCL 15 511 + 512 RCL 03 513 * 514 RCL 12 515 * 516 RCL 01 |
517 RCL 05
518 * 519 X^2 520 LASTX 521 * 522 6 523 / 524 + 525 RCL 14 526 / 527 + 528 RCL 03 529 RCL 12 530 * 531 X^2 532 LASTX 533 * 534 RCL 08 535 / 536 - 537 RCL 04 538 * 539 RCL 05 540 RCL 11 541 * 542 RCL 03 543 - 544 RCL 07 545 * 546 RCL 10 547 * 548 RCL 12 549 ST+ X 550 / 551 + 552 RCL 13 553 RCL 01 554 RCL 11 555 X^2 556 * 557 ST+ X 558 - 559 RCL 00 560 ST* Y 561 RCL 06 562 - 563 RCL 01 564 * 565 + 566 RCL 02 567 / 568 + 569 RCL 01 570 RCL 03 571 * 572 RCL 05 573 * 574 RCL 11 575 * 576 RCL 14 577 / 578 + 579 RCL 04 580 * 581 RCL 03 582 RCL 12 583 * 584 + 585 RCL 01 586 RCL 13 587 + 588 RCL 05 589 * 590 RCL 14 591 / 592 - 593 RCL 04 594 * 595 RCL 13 596 + 597 RCL 16 598 * 599 STO Y 600 RCL 17 601 * 602 END |
( 699 bytes / SIZE 018 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | d ( km ) |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial
ellipsoid
and d is the geodesic
distance on the ellipsoid of revolution ( WGS 84 )
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12138.65755 km
---Execution time = 59s---
X<>Y d = 12138.68432 km
3°) With Andoyer's Formula of order 3
a) Andoyer's Formula ( Ellipsoid
of Revolution )
-Sodano's formula may be re-written as an Andoyer's formula.
( In reference [1] - page 87 - Paul D. Thomas gives the key to
perform this transformation )
With L , L' = longitudes & B , B' = parametric latitudes a = 6378.137 km , b = 6356.752314 km ( WGS84 ) , f = (a-b) / a
Tan B = (1-f) Tan l where l = geodetic latitude
H = ( L' - L )/2 ; F = ( B + B' )/2 ; G = ( B' - B )/2
S = sin2 G cos2 H + cos2 F sin2 H , µ = 2 Arc Sin S1/2
we define A = ( sin2 G cos2
F )/S + ( sin2 F cos2 G )/( 1-S )
and B
= ( sin2 G cos2 F )/S - ( sin2 F
cos2 G )/( 1-S )
-After "some" calculations, it yields:
Dist = a { µ - (f/2) ( A.µ + B Sin µ )
+ (f2/4) [ A µ2 Cot µ + B µ2 / Sin µ + A2 ( µ + Sin µ Cos µ - 4 µ2 Cot µ ) / 4 - A.B µ2 / Sin µ - ( B2/2 ) Sin µ Cos µ ]
+ (f3/8) [ A ( 2 µ2 Cot µ - 2 µ3
Cot2 µ - 8 µ3 / 6 ) + A2 (
-5 µ2 Cot µ + 5 µ3 Cot2
µ + 8 µ3 / 3 + (1/2) Sin µ Cos µ + µ
/2 )
+ A3 ( 3 µ2 Cot µ - 3 µ3
Cot2 µ - 8 µ3 / 6 + 2 Sin µ Cos3
µ - (9/4) Sin µ Cos µ - µ / 4 + 2 Sin3
µ Cos3 µ + 2 Sin5 µ Cos µ
)
+ B ( 2 µ2 / Sin µ - 2 µ3 Cos µ
/ Sin2 µ ) + A.B ( -5 µ2 / Sin µ
+ 6 µ3 Cos µ / Sin2 µ + µ
Cos µ )
+ A2 B ( 3 µ2 / Sin µ - 4 µ3
Cos µ / Sin2 µ - 2 Sin3 µ
Cos2 µ - µ Cos µ - 2 Sin5 µ
+ 2 Sin µ - (3/2) Sin µ Cos2 µ )
+ B2 ( - Sin µ Cos µ + µ3 / Sin2
µ + µ ) + A B2 ( (1/2) Sin µ Cos µ
- µ - µ3 / Sin2 µ )
+ B3 ( (2/3) Sin3 µ - (1/2) Sin µ ) ]
}
Data Registers: R00 thru R15: temp
Flags: /
Subroutines: /
01 LBL "AND"
02 DEG 03 HR 04 STO 01 05 X<> T 06 HMS- 07 HR 08 2 09 / 10 SIN 11 X^2 12 STO 00 13 X<>Y 14 HR 15 1 16 P-R 17 LASTX 18 298.25722 19 STO 04 20 ST+ 04 21 1/X 22 - 23 STO 05 24 / 25 R-P 26 X<> 01 27 1 28 P-R 29 RCL 05 30 / 31 R-P 32 RDN 33 + 34 2 35 / 36 ST- Y 37 1 38 P-R 39 X^2 40 STO 02 41 X<>Y 42 X^2 43 STO 01 |
44 RDN
45 X<>Y 46 1 47 P-R 48 X^2 49 ST* 01 50 RDN 51 X^2 52 ST* 02 53 STO Z 54 - 55 RCL 00 56 * 57 + 58 ST/ 02 59 STO 05 60 SQRT 61 ASIN 62 ST+ X 63 STO 03 64 D-R 65 STO 00 66 1 67 RCL 05 68 - 69 ST/ 01 70 RCL 02 71 RCL 01 72 ST- 02 73 + 74 STO 01 75 RCL 03 76 1 77 P-R 78 STO 06 79 X<>Y 80 STO 05 81 / 82 STO 07 83 ENTER 84 X^2 85 RCL 00 86 * |
87 -
88 RCL 00 89 X^2 90 STO 09 91 * 92 STO 03 93 3 94 * 95 RCL 00 96 RCL 09 97 * 98 STO 10 99 .75 100 STO 15 101 / 102 STO 13 103 - 104 RCL 05 105 RCL 06 106 * 107 STO 11 108 LASTX 109 X^2 110 * 111 RCL 11 112 ENTER 113 X^2 114 STO 08 115 * 116 + 117 RCL 05 118 ST* 08 119 X^2 120 X^2 121 STO 12 122 RCL 11 123 * 124 + 125 ST+ X 126 + 127 RCL 11 128 3 129 * |
130 RCL 15
131 * 132 - 133 RCL 00 134 4 135 / 136 - 137 RCL 01 138 * 139 RCL 03 140 5 141 * 142 - 143 RCL 13 144 ST+ X 145 + 146 RCL 00 147 RCL 11 148 + 149 2 150 / 151 + 152 RCL 05 153 RCL 05 154 RCL 12 155 * 156 - 157 RCL 08 158 - 159 ST+ X 160 RCL 10 161 RCL 05 162 X^2 163 / 164 STO 12 165 RCL 06 166 * 167 STO 13 168 4 169 * 170 - 171 RCL 09 172 RCL 05 |
173 /
174 STO 14 175 3 176 * 177 + 178 RCL 00 179 RCL 06 180 * 181 STO 08 182 - 183 RCL 06 184 RCL 11 185 * 186 RCL 15 187 ST+ X 188 STO 06 189 * 190 - 191 RCL 02 192 * 193 + 194 RCL 01 195 * 196 RCL 03 197 ST+ X 198 + 199 RCL 10 200 RCL 15 201 / 202 - 203 RCL 13 204 6 205 * 206 RCL 14 207 5 208 * 209 - 210 RCL 08 211 + 212 RCL 02 213 * 214 + 215 RCL 01 |
216 *
217 RCL 11 218 2 219 / 220 RCL 00 221 - 222 RCL 12 223 - 224 RCL 01 225 * 226 RCL 05 227 X^2 228 RCL 06 229 / 230 .5 231 - 232 RCL 05 233 * 234 RCL 02 235 * 236 + 237 RCL 11 238 - 239 RCL 12 240 + 241 RCL 00 242 + 243 RCL 02 244 * 245 RCL 14 246 RCL 13 247 - 248 ST+ X 249 + 250 RCL 02 251 * 252 + 253 RCL 04 254 / 255 RCL 00 256 RCL 11 257 + 258 4 |
259 /
260 RCL 07 261 RCL 09 262 * 263 STO 13 264 - 265 RCL 01 266 * 267 RCL 02 268 ST* 05 269 ST* 11 270 RCL 14 271 * 272 - 273 RCL 13 274 + 275 RCL 01 276 * 277 + 278 RCL 14 279 RCL 11 280 2 281 / 282 - 283 RCL 02 284 * 285 + 286 RCL 04 287 / 288 RCL 00 289 RCL 01 290 * 291 - 292 RCL 05 293 - 294 RCL 04 295 / 296 RCL 00 297 + 298 6378.137 299 * 300 END |
( 346 bytes / SIZE 015 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | / |
X | b' ( ° ' " ) | D ( km ) |
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "AND"
>>>> D = 12138.68432 km
---Execution time = 16s---
Notes:
-"AND" is slightly shorter & faster than "SOD"
-The formulas are quite equivalent, so the differences in the results
are only due to roundoff-errors.
-It is certainly possible to express the formula without calculating the parametric latitudes, but I've not found how.
-In order to calculate the parametric latitudes, we can also use TAN
& ATAN instead of P-R & R-P
-But the HP41 gives TAN (-90°) = +9.999999999 E99 instead of -9.999999999
E99
-So, we must use a trick to take this case into account:
-We can replace lines 01 to 32 by
01 LBL "AND"
02 DEG 03 HR 04 X<> T 05 HMS- 06 HR 07 2 08 / 09 SIN 10 X^2 11 STO 00 12 RDN 13 HR 14 X<Y? 15 X<>Y 16 TAN 17 PI 18 SIGN 19 298.25722 20 STO 04 21 ST+ 04 22 1/X 23 - 24 STO 05 25 * 26 ATAN 27 X<>Y 28 SIGN 29 LASTX 30 ABS 31 TAN 32 RCL 05 33 * 34 * 35 ATAN |
-It costs 2 extra bytes but it saves almost 1 second with an HP41 !
-However, this method cannot be used with free42 because 90 TAN gives
OUT OF RANGE
b) Triaxial Ellipsoid
-The following "TGD" employs Andoyer's formula of order 3 for the last
calculation and Andoyer's formula of order 2 to estimate (A/A')
Data Registers: R00 thru R17: temp
Flags: /
Subroutines: /
-Line 38 is a three-byte GTO 14
01 LBL "TGD"
02 DEG 03 HR 04 STO 15 05 RDN 06 HR 07 14.92911 08 STO 00 09 + 10 STO 14 11 RDN 12 HR 13 STO 13 14 X<>Y 15 HR 16 RCL 00 17 + 18 STO 12 19 6378.172 20 STO 01 21 .07 22 - 23 STO 02 24 21.349686 25 - 26 STO 03 27 XEQ 01 28 STO 17 29 RCL 01 30 RCL 02 31 + 32 2 33 / 34 STO 01 35 STO 02 36 XEQ 01 37 ST/ 17 38 GTO 14 39 LBL 01 40 RCL 12 41 RCL 13 42 XEQ 02 43 RCL 08 44 STO 05 45 RCL 09 46 STO 06 47 RCL 10 48 STO 07 49 RCL 14 50 RCL 15 51 XEQ 02 52 RCL 05 53 RCL 08 54 * 55 RCL 06 56 RCL 09 57 * 58 RCL 07 59 RCL 10 60 * 61 + 62 + 63 STO 04 64 RCL 01 65 ST/ 05 66 ST/ 08 67 RCL 02 68 ST/ 06 69 ST/ 09 70 RCL 03 71 ST/ 07 72 ST/ 10 73 RCL 05 74 X^2 75 RCL 06 76 X^2 77 RCL 07 78 X^2 79 + 80 + |
81 STO 11
82 RCL 08 83 X^2 84 RCL 09 85 X^2 86 RCL 10 87 X^2 88 + 89 + 90 X<> 10 91 RCL 07 92 * 93 RCL 06 94 RCL 09 95 * 96 RCL 05 97 RCL 08 98 * 99 + 100 + 101 STO 09 102 ST+ 09 103 RCL 04 104 RCL 11 105 * 106 - 107 ST+ X 108 RCL 04 109 RCL 09 110 * 111 RCL 10 112 - 113 RCL 11 114 RCL 04 115 X^2 116 * 117 - 118 RCL 11 119 RCL 04 120 ACOS 121 STO 04 122 SIN 123 STO 08 124 ST/ Z 125 * 126 + 127 X#0? 128 / 129 ATAN 130 2 131 / 132 90 133 STO 06 134 X<>Y 135 X<0? 136 + 137 STO 05 138 XEQ 03 139 X<> 06 140 RCL 05 141 + 142 XEQ 03 143 STO 07 144 RCL 05 145 CHS 146 ST+ 04 147 1 148 P-R 149 RCL 07 150 RCL 06 151 / 152 STO 09 153 * 154 R-P 155 X<> 04 156 1 157 P-R 158 RCL 09 159 * 160 R-P |
161 RDN
162 ST- Z 163 + 164 COS 165 STO 05 166 X^2 167 ST+ X 168 CHS 169 X<>Y 170 ABS 171 ENTER 172 D-R 173 STO 10 174 SIGN 175 P-R 176 X<>Y 177 ST* 05 178 * 179 ST* Y 180 + 181 RCL 10 182 + 183 4 184 / 185 RCL 07 186 CHS 187 RCL 06 188 ST+ Y 189 ST+ X 190 / 191 STO 09 192 * 193 RCL 05 194 - 195 RCL 10 196 - 197 RCL 09 198 * 199 RCL 10 200 + 201 RCL 06 202 * 203 RTN 204 LBL 02 205 1 206 P-R 207 X<>Y 208 RCL 03 209 X^2 210 * 211 STO 10 212 X^2 213 X<> Z 214 X<>Y 215 P-R 216 RCL 01 217 X^2 218 * 219 STO 08 220 X^2 221 X<>Y 222 RCL 02 223 X^2 224 * 225 STO 09 226 X^2 227 + 228 + 229 SQRT 230 ST/ 08 231 ST/ 09 232 ST/ 10 233 RTN 234 LBL 03 235 ENTER 236 SIN 237 STO 07 238 X^2 239 RCL 10 240 * |
241 X<>Y
242 RCL 04 243 - 244 SIN 245 ST* 07 246 X^2 247 RCL 11 248 * 249 + 250 RCL 07 251 RCL 09 252 * 253 - 254 SQRT 255 RCL 08 256 X<>Y 257 / 258 RTN 259 LBL 14 260 RCL 14 261 RCL 12 262 - 263 2 264 / 265 SIN 266 X^2 267 STO 00 268 RCL 13 269 1 270 P-R 271 RCL 03 272 CHS 273 RCL 01 274 STO 16 275 ST+ Y 276 / 277 STO 04 278 SIGN 279 LASTX 280 - 281 STO 05 282 / 283 R-P 284 X<> 15 285 1 286 P-R 287 RCL 05 288 / 289 R-P 290 RDN 291 + 292 2 293 ST/ 04 294 / 295 ST- Y 296 1 297 P-R 298 X^2 299 STO 02 300 X<>Y 301 X^2 302 STO 01 303 RDN 304 X<>Y 305 1 306 P-R 307 X^2 308 ST* 01 309 RDN 310 X^2 311 ST* 02 312 STO Z 313 - 314 RCL 00 315 * 316 + 317 ST/ 02 318 STO 05 319 SQRT 320 ASIN |
321 ST+ X
322 STO 03 323 D-R 324 STO 00 325 1 326 RCL 05 327 - 328 ST/ 01 329 RCL 02 330 RCL 01 331 ST- 02 332 + 333 STO 01 334 RCL 03 335 1 336 P-R 337 STO 06 338 X<>Y 339 STO 05 340 / 341 STO 07 342 ENTER 343 X^2 344 RCL 00 345 * 346 - 347 RCL 00 348 X^2 349 STO 09 350 * 351 STO 03 352 3 353 * 354 RCL 00 355 RCL 09 356 * 357 STO 10 358 .75 359 STO 15 360 / 361 STO 13 362 - 363 RCL 05 364 RCL 06 365 * 366 STO 11 367 LASTX 368 X^2 369 * 370 RCL 11 371 ENTER 372 X^2 373 STO 08 374 * 375 + 376 RCL 05 377 ST* 08 378 X^2 379 X^2 380 STO 12 381 RCL 11 382 * 383 + 384 ST+ X 385 + 386 RCL 11 387 3 388 * 389 RCL 15 390 * 391 - 392 RCL 00 393 4 394 / 395 - 396 RCL 01 397 * 398 RCL 03 399 5 400 * |
401 -
402 RCL 13 403 ST+ X 404 + 405 RCL 00 406 RCL 11 407 + 408 2 409 / 410 + 411 RCL 05 412 RCL 05 413 RCL 12 414 * 415 - 416 RCL 08 417 - 418 ST+ X 419 RCL 10 420 RCL 05 421 X^2 422 / 423 STO 12 424 RCL 06 425 * 426 STO 13 427 4 428 * 429 - 430 RCL 09 431 RCL 05 432 / 433 STO 14 434 3 435 * 436 + 437 RCL 00 438 RCL 06 439 * 440 STO 08 441 - 442 RCL 06 443 RCL 11 444 * 445 RCL 15 446 ST+ X 447 STO 06 448 * 449 - 450 RCL 02 451 * 452 + 453 RCL 01 454 * 455 RCL 03 456 ST+ X 457 + 458 RCL 10 459 RCL 15 460 / 461 - 462 RCL 13 463 6 464 * 465 RCL 14 466 5 467 * 468 - 469 RCL 08 470 + 471 RCL 02 472 * 473 + 474 RCL 01 475 * 476 RCL 11 477 2 478 / 479 RCL 00 480 - 481 RCL 12 |
482 -
483 RCL 01 484 * 485 RCL 05 486 X^2 487 RCL 06 488 / 489 .5 490 - 491 RCL 05 492 * 493 RCL 02 494 * 495 + 496 RCL 11 497 - 498 RCL 12 499 + 500 RCL 00 501 + 502 RCL 02 503 * 504 RCL 14 505 RCL 13 506 - 507 ST+ X 508 + 509 RCL 02 510 * 511 + 512 RCL 04 513 * 514 RCL 00 515 RCL 11 516 + 517 4 518 / 519 RCL 07 520 RCL 09 521 * 522 STO 13 523 - 524 RCL 01 525 * 526 RCL 02 527 ST* 05 528 ST* 11 529 RCL 14 530 * 531 - 532 RCL 13 533 + 534 RCL 01 535 * 536 + 537 RCL 14 538 RCL 11 539 2 540 / 541 - 542 RCL 02 543 * 544 + 545 RCL 04 546 * 547 RCL 00 548 RCL 01 549 * 550 - 551 RCL 05 552 - 553 RCL 04 554 * 555 RCL 00 556 + 557 RCL 16 558 * 559 STO Y 560 RCL 17 561 * 562 END |
( 662 bytes / SIZE 018 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | d ( km ) |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial
ellipsoid
and d is the geodesic
distance on the ellipsoid of revolution ( WGS 84 )
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12138.65755 km
---Execution time = 56s---
X<>Y d = 12138.68432 km
Note:
-If you want to use TAN & ATAN instead of P-R
& R-P to calculate the parametric latitudes, replace lines 268
to 290 by
268 RCL 13
269 RCL 15 270 X<Y? 271 X<>Y 272 TAN 273 RCL 03 274 CHS 275 RCL 01 276 STO 16 277 ST+ Y 278 / 279 STO 04 280 SIGN 281 LASTX 282 - 283 STO 05 284 * 285 ATAN 286 X<>Y 287 SIGN 288 LASTX 289 ABS 290 TAN 291 RCL 05 292 * 293 * 294 ATAN |
-Just 3 extra bytes, but almost 1 second saved.
4°) With Vincenty's Formula
-Vincenty's formula is probably the best method for an ellipsoid of
revolution.
-The following version of "TGD" - renamed "TGDV"- employs this
method 3 times.
Data Registers: R00 thru R21: temp
Flags: /
Subroutines: /
-Line 47 is a three-byte GTO 14
01 LBL "TGDV"
02 DEG 03 HR 04 STO 21 05 RDN 06 HR 07 14.92911 08 STO 00 09 + 10 STO 20 11 RDN 12 HR 13 STO 19 14 X<>Y 15 HR 16 RCL 00 17 + 18 STO 18 19 6378.172 20 STO 01 21 .07 22 - 23 STO 02 24 6356.752314 25 STO 03 26 XEQ 01 27 STO 17 28 RCL 01 29 RCL 02 30 + 31 2 32 / 33 STO 01 34 STO 02 35 XEQ 01 36 ST/ 17 37 RCL 01 38 STO 06 39 RCL 03 40 STO 07 41 RCL 21 42 RCL 19 43 RCL 20 44 RCL 18 45 XEQ 13 46 ENTER 47 GTO 14 48 LBL 01 49 RCL 18 50 RCL 19 51 XEQ 02 52 RCL 08 53 STO 05 54 RCL 09 55 STO 06 56 RCL 10 57 STO 07 58 RCL 20 59 RCL 21 60 XEQ 02 61 RCL 05 62 RCL 08 63 * 64 RCL 06 |
65 RCL 09
66 * 67 RCL 07 68 RCL 10 69 * 70 + 71 + 72 STO 00 73 ACOS 74 STO 04 75 RCL 01 76 ST/ 05 77 ST/ 08 78 RCL 02 79 ST/ 06 80 ST/ 09 81 RCL 03 82 ST/ 07 83 ST/ 10 84 RCL 05 85 X^2 86 RCL 06 87 X^2 88 RCL 07 89 X^2 90 + 91 + 92 STO 11 93 RCL 08 94 X^2 95 RCL 09 96 X^2 97 RCL 10 98 X^2 99 + 100 + 101 X<> 10 102 RCL 07 103 * 104 RCL 06 105 RCL 09 106 * 107 RCL 05 108 RCL 08 109 * 110 + 111 + 112 STO 09 113 RCL 00 114 * 115 ST+ X 116 RCL 10 117 - 118 RCL 11 119 RCL 00 120 X^2 121 * 122 - 123 RCL 11 124 RCL 04 125 SIN 126 STO 08 127 ST/ Z 128 * |
129 +
130 RCL 09 131 RCL 11 132 RCL 00 133 * 134 - 135 ST+ X 136 X<>Y 137 X#0? 138 / 139 ATAN 140 2 141 / 142 90 143 STO 00 144 X<>Y 145 X<0? 146 + 147 STO 05 148 XEQ 03 149 STO 06 150 RCL 05 151 RCL 00 152 + 153 XEQ 03 154 STO 07 155 RCL 06 156 X>Y? 157 GTO 00 158 X<> 07 159 STO 06 160 RCL 00 161 ST+ 05 162 LBL 00 163 PI 164 R-D 165 RCL 05 166 X=Y? 167 CLX 168 STO 05 169 CHS 170 1 171 P-R 172 RCL 07 173 RCL 06 174 / 175 X^2 176 STO 09 177 * 178 R-P 179 CLX 180 RCL 04 181 RCL 05 182 - 183 1 184 P-R 185 RCL 09 186 * 187 R-P 188 CLX 189 ENTER 190 LBL 13 191 - 192 360 |
193 MOD
194 PI 195 R-D 196 STO 16 197 X>Y? 198 CLX 199 ST+ X 200 - 201 STO 11 202 STO 04 203 CLX 204 RCL 07 205 CHS 206 RCL 06 207 ST+ Y 208 / 209 STO 08 210 SIGN 211 LASTX 212 - 213 STO 15 214 P-R 215 LASTX 216 / 217 R-P 218 RDN 219 STO 09 220 CLX 221 RCL 15 222 P-R 223 LASTX 224 / 225 R-P 226 X<>Y 227 STO 10 228 LBL 12 229 RCL 16 230 RCL 04 231 X>Y? 232 ACOS 233 RCL 10 234 COS 235 STO 05 236 STO 12 237 RCL 04 238 SIN 239 * 240 STO 14 241 RCL 09 242 COS 243 ST* 05 244 ST* 14 245 RCL 10 246 SIN 247 STO 13 248 * 249 RCL 09 250 SIN 251 ST* 13 252 RCL 12 253 * 254 RCL 04 255 COS 256 ST* 05 |
257 *
258 - 259 R-P 260 RCL 05 261 RCL 13 262 + 263 R-P 264 X<> 14 265 X<>Y 266 STO 05 267 SIN 268 X#0? 269 / 270 ASIN 271 STO 12 272 COS 273 X^2 274 RCL 05 275 COS 276 RCL 13 277 ENTER 278 + 279 R^ 280 X=0? 281 SIGN 282 / 283 - 284 STO 13 285 CLX 286 3 287 * 288 4 289 X<>Y 290 - 291 RCL 08 292 * 293 4 294 + 295 * 296 RCL 08 297 * 298 16 299 / 300 RCL 05 301 SIN 302 RCL 13 303 RCL Z 304 * 305 * 306 R-D 307 RCL 05 308 + 309 RCL 12 310 SIN 311 * 312 RCL 08 313 * 314 1 315 R^ 316 - 317 * 318 RCL 11 319 + 320 ENTER |
321 X<> 04
322 - 323 ABS 324 E-7 325 X<Y? 326 GTO 12 327 2 328 RCL 08 329 ST- Y 330 * 331 RCL 12 332 COS 333 RCL 15 334 / 335 X^2 336 * 337 12 338 RCL Y 339 5 340 * 341 - 342 * 343 64 344 - 345 * 346 256 347 ST- Y 348 / 349 STO 14 350 CLX 351 37 352 * 353 64 354 - 355 * 356 512 357 / 358 4 359 1/X 360 + 361 * 362 RCL 13 363 X^2 364 ST+ X 365 1 366 - 367 RCL 05 368 COS 369 4 370 / 371 * 372 * 373 RCL 13 374 + 375 * 376 RCL 05 377 SIN 378 * 379 RCL 05 380 D-R 381 - 382 RCL 14 383 * 384 RCL 15 |
385 *
386 RCL 06 387 * 388 RTN 389 LBL 02 390 1 391 P-R 392 X<>Y 393 RCL 03 394 X^2 395 * 396 STO 10 397 X^2 398 STO 04 399 RDN 400 P-R 401 RCL 01 402 X^2 403 * 404 STO 08 405 X^2 406 ST+ 04 407 X<>Y 408 RCL 02 409 X^2 410 * 411 STO 09 412 X^2 413 RCL 04 414 + 415 SQRT 416 ST/ 08 417 ST/ 09 418 ST/ 10 419 RTN 420 LBL 03 421 ENTER 422 SIN 423 STO 07 424 X^2 425 RCL 10 426 * 427 X<>Y 428 RCL 04 429 - 430 SIN 431 ST* 07 432 X^2 433 RCL 11 434 * 435 + 436 RCL 07 437 RCL 09 438 ST+ X 439 * 440 - 441 SQRT 442 RCL 08 443 X<>Y 444 / 445 RTN 446 LBL 14 447 RCL 17 448 * 449 END |
( 564 bytes / SIZE 022 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | d ( km ) |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial
ellipsoid
and d is the geodesic
distance on the ellipsoid of revolution ( WGS 84 )
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGDV"
>>>> D = 12138.65757 km
---Execution time = 119s---
X<>Y d = 12138.68433 km
Notes:
-This version gives the best accuarcy, provided the points are not (
almost ) antipodal.
-For nearly antipodal points, the routine stops and the HP41 displays
DATA ERROR ( line 235 )
5°) Minimizing the Sum of 2 Distances
-No one of the programs above works well for nearly antipodal points.
-The following one uses the golden search to minimize the sum of 2
geodesic distances:
Let M(L,b) & N(L',b') , "D+D" minimizes the sum D = Dist(M,P) + Dist(P,N)
where P(L",b") with L" = the longitude of the point on a sphere in the middle of the geodesic MN on the sphere.
-The program starts with b" = -90° and +90° and
employs the golden search to find the minimum D-value.
Data Registers: R00 thru R35: temp
Flags: /
Subroutine: "TGD" ( cf paragraph 2-b
or 3-b above, using 3rd order formulas )
01 LBL "D+D"
02 STO 26 03 RDN 04 STO 25 05 RDN 06 STO 24 07 X<>Y 08 STO 23 09 R^ 10 X<>Y 11 HMS- 12 HR 13 RCL 26 14 HR 15 COS 16 P-R 17 RCL 24 18 HR 19 COS 20 + 21 R-P 22 CLX |
23 RCL 23
24 HR 25 + 26 HMS 27 STO 27 28 GTO 00 29 LBL 12 30 STO 35 31 RCL 27 32 X<>Y 33 RCL 25 34 RCL 26 35 XEQ "TGD" 36 STO 31 37 RCL 23 38 RCL 24 39 RCL 27 40 RCL 35 41 XEQ "TGD" 42 RCL 31 43 + 44 RTN |
45 LBL 00
46 90 47 ENTER 48 STO 29 49 CHS 50 STO 28 51 - 52 5 53 SQRT 54 1 55 + 56 2 57 / 58 STO 34 59 / 60 RCL 28 61 HR 62 + 63 HMS 64 STO 32 65 XEQ 12 66 STO 33 |
67 LBL 10
68 RCL 32 69 RCL 28 70 HMS- 71 HR 72 RCL 34 73 / 74 RCL 28 75 HR 76 + 77 HMS 78 STO 30 79 XEQ 12 80 STO 31 81 VIEW 31 82 RCL 33 83 X<Y? 84 GTO 02 85 RCL 30 86 X<> 32 87 STO 29 88 RCL 31 |
89 STO 33
90 GTO 04 91 LBL 02 92 RCL 29 93 STO 28 94 RCL 30 95 STO 29 96 LBL 04 97 RCL 13 98 PI 99 R-D 100 / 101 RCL 28 102 RCL 30 103 HMS- 104 HR 105 ABS 106 X>Y? 107 GTO 10 108 RCL 30 109 RCL 31 110 CLD 111 END |
( 181 bytes / SIZE 036 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | b" ( ° ' " ) |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial
ellipsoid
and b" is the latitude of
the "midpoint"
Example1: With (L,b) = (0°,0°) & (L',b') = (179.51°,0°)
0
ENTER^ ENTER^
179.51 ENTER^
0
XEQ "D+D" >>>> D = 20001.89955 km
X<>Y b" ~ -75°38'
Example2: With (L,b) = (-40°,-40°) & (L',b') = (120°,50°)
40 CHS ENTER^ ENTER^
120
ENTER^
50
XEQ "D+D" >>>> D = 18095.49447 km
X<>Y b" ~ 24°58'
Notes:
-The HP41 displays the successive approximations.
-The angle b" is not computed accurately: it's not necessary
to get a good precision for the distance.
-You can press XEQ 10 to make another loop, or - for instance
- replace lines 98-99 by a larger number.
-This program is very slow with a true HP41, so a good emulator is not
superfluous !
-With free42, you get the results 20001.89954 & 18095.49449
almost instantaneously...
-This program may also be used for very long geodesics when the points
are not almost antipodal,
but when Andoyer's or Sodano's formulas give an error of a few
meters.
-Since formulas of order 3 is very accurate when the distance is about
10000 km, the sum will remain very accurate too !
Variant:
-The following routine may be used to find the minimum distance from
several values of b"
01 LBL "D+D"
02 STO 26 03 RDN 04 STO 25 05 RDN 06 STO 24 07 X<>Y 08 STO 23 09 R^ 10 X<>Y 11 HMS- 12 HR 13 RCL 26 |
14 HR
15 COS 16 P-R 17 RCL 24 18 HR 19 COS 20 + 21 R-P 22 RCL 26 23 HR 24 SIN 25 RCL 24 26 HR |
27 SIN
28 + 29 X<>Y 30 R-P 31 RDN 32 HMS 33 STO 28 34 X<>Y 35 RCL 23 36 HR 37 + 38 HMS 39 STO 27 |
40 X<>Y
41 " B0=" 42 ARCL X 43 PROMPT 44 LBL 10 45 " B=?" 46 PROMPT 47 STO 29 48 RCL 27 49 X<>Y 50 RCL 25 51 RCL 26 52 XEQ "TGD" |
53 STO 30
54 RCL 23 55 RCL 24 56 RCL 27 57 RCL 29 58 XEQ "TGD" 59 RCL 30 60 + 61 " D=" 62 ARCL X 63 PROMPT 64 GTO 10 65 END |
( 116 bytes )
STACK | INPUTS |
T | L ( ° ' " ) |
Z | b ( ° ' " ) |
Y | L' ( ° ' " ) |
X | b' ( ° ' " ) |
Example: With (L,b) = (-40°,-40°)
& (L',b') = (120°,50°)
40 CHS ENTER^ ENTER^
120
ENTER^
50
XEQ "D+D" >>>> B0=24.1727 i-e
24°17'27"
R/S >>>> B=? if you choose for instance: b" = 12° 17° 22° 27° 32°
12 R/S >>>> D=18114.88342
R/S >>>> B=?
17 R/S >>>> D=18102.74082
R/S >>>> B=?
22 R/S >>>> D=18096.49637
R/S >>>> B=?
27 R/S >>>> D=18095.96331
R/S >>>> B=?
32 R/S >>>> D=18101.15912
-Store these 5 results in R01 thru R05 and if you use the programs listed in "Extrema for the HP41" paragraph 1-d)
5 ENTER^
22 XEQ "DEX5" >>>> b" =
24.96999 = 24°58'12"
X<>Y D = 18095.49452 km
-With a true HP41, it will be less slow than the first "D+D"
-However, it will not be always easy to estimate the proper - and small
- interval for b" for nearly antipodal points...
6°) A Slightly Different Method
a) With Andoyer's Formula
of Order 2
-In the 3 following programs, intead of calculating D ~
(A/A') d , we compute D ~ A + ( d - A' )
-In other words, the part of the formula giving the great elliptic
arc distance on the ellipsoid of revolution in the WGS84 geodesic distance
is replaced by the great elliptic arc distance on the triaxial
ellipsoid.
-The precision is almost the same, but the routines are faster !
Formula:
d - A' = (f2/4) µ [ (A2-A)
( 1 - µ Cot µ ) - B.(A-1) ( µ / Sin µ - (Sin µ)
/ µ ) ] where µ , A , B
are computed with parametric latitudes
Data Registers: R00 thru R15: temp
Flags: /
Subroutines: /
01 LBL "TGD"
02 DEG 03 HR 04 STO 15 05 RDN 06 HR 07 14.92911 08 STO 00 09 + 10 STO 14 11 RDN 12 HR 13 STO 13 14 X<>Y 15 HR 16 RCL 00 17 + 18 STO 12 19 6378.172 20 STO 01 21 .07 22 - 23 STO 02 24 21.349686 25 - 26 STO 03 27 RCL 12 28 RCL 13 29 XEQ 02 30 RCL 08 31 STO 05 32 RCL 09 33 STO 06 34 RCL 10 35 STO 07 36 RCL 14 37 RCL 15 38 XEQ 02 39 RCL 05 40 RCL 08 41 * 42 RCL 06 43 RCL 09 44 * 45 RCL 07 46 RCL 10 47 * 48 + 49 + 50 STO 04 51 RCL 01 |
52 ST/ 05
53 ST/ 08 54 RCL 02 55 ST/ 06 56 ST/ 09 57 RCL 03 58 ST/ 07 59 ST/ 10 60 RCL 05 61 X^2 62 RCL 06 63 X^2 64 RCL 07 65 X^2 66 + 67 + 68 STO 11 69 RCL 08 70 X^2 71 RCL 09 72 X^2 73 RCL 10 74 X^2 75 + 76 + 77 X<> 10 78 RCL 07 79 * 80 RCL 06 81 RCL 09 82 * 83 RCL 05 84 RCL 08 85 * 86 + 87 + 88 STO 09 89 ST+ 09 90 RCL 04 91 RCL 11 92 * 93 - 94 ST+ X 95 RCL 04 96 RCL 09 97 * 98 RCL 10 99 - 100 RCL 11 101 RCL 04 102 X^2 |
103 *
104 - 105 RCL 11 106 RCL 04 107 ACOS 108 STO 04 109 SIN 110 STO 08 111 ST/ Z 112 * 113 + 114 X#0? 115 / 116 ATAN 117 2 118 / 119 90 120 STO 06 121 X<>Y 122 X<0? 123 + 124 STO 05 125 XEQ 03 126 X<> 06 127 RCL 05 128 + 129 XEQ 03 130 STO 07 131 RCL 05 132 CHS 133 ST+ 04 134 1 135 P-R 136 RCL 07 137 RCL 06 138 / 139 STO 09 140 * 141 R-P 142 X<> 04 143 1 144 P-R 145 RCL 09 146 * 147 R-P 148 RDN 149 ST- Z 150 + 151 COS 152 STO 05 153 X^2 |
154 ST+ X
155 CHS 156 X<>Y 157 ABS 158 ENTER 159 D-R 160 STO 10 161 SIGN 162 P-R 163 X<>Y 164 ST* 05 165 * 166 ST* Y 167 + 168 RCL 10 169 + 170 4 171 / 172 RCL 07 173 CHS 174 RCL 06 175 ST+ Y 176 ST+ X 177 / 178 STO 09 179 * 180 RCL 05 181 - 182 RCL 10 183 - 184 RCL 09 185 * 186 RCL 10 187 + 188 ST* 06 189 GTO 14 190 LBL 02 191 1 192 P-R 193 X<>Y 194 RCL 03 195 X^2 196 * 197 STO 10 198 X^2 199 X<> Z 200 X<>Y 201 P-R 202 RCL 01 203 X^2 204 * |
205 STO 08
206 X^2 207 X<>Y 208 RCL 02 209 X^2 210 * 211 STO 09 212 X^2 213 + 214 + 215 SQRT 216 ST/ 08 217 ST/ 09 218 ST/ 10 219 RTN 220 LBL 03 221 ENTER 222 SIN 223 STO 07 224 X^2 225 RCL 10 226 * 227 X<>Y 228 RCL 04 229 - 230 SIN 231 ST* 07 232 X^2 233 RCL 11 234 * 235 + 236 RCL 07 237 RCL 09 238 * 239 - 240 SQRT 241 RCL 08 242 X<>Y 243 / 244 RTN 245 LBL 14 246 RCL 01 247 RCL 02 248 + 249 2 250 / 251 STO 11 252 RCL 14 253 RCL 12 254 - 255 2 |
256 /
257 SIN 258 X^2 259 STO 10 260 RCL 13 261 1 262 P-R 263 RCL 03 264 CHS 265 RCL 11 266 ST+ Y 267 / 268 STO 04 269 SIGN 270 LASTX 271 - 272 STO 05 273 / 274 R-P 275 X<> 15 276 1 277 P-R 278 RCL 05 279 / 280 R-P 281 RDN 282 + 283 2 284 ST/ 04 285 / 286 ST- Y 287 1 288 P-R 289 X^2 290 STO 09 291 X<>Y 292 X^2 293 STO 08 294 RDN 295 X<>Y 296 1 297 P-R 298 X^2 299 ST* 08 300 RDN 301 X^2 302 ST* 09 303 STO Z 304 - 305 RCL 10 306 * |
307 +
308 ST/ 09 309 STO 05 310 SQRT 311 ASIN 312 ST+ X 313 STO 10 314 D-R 315 STO 07 316 1 317 RCL 05 318 - 319 ST/ 08 320 RCL 09 321 RCL 08 322 ST- 09 323 + 324 STO 08 325 RCL 10 326 1 327 P-R 328 X<>Y 329 STO 05 330 / 331 RCL 07 332 * 333 RCL 08 334 ST* Y 335 - 336 RCL 05 337 RCL 07 338 / 339 ENTER 340 1/X 341 - 342 RCL 09 343 * 344 - 345 1 346 RCL 08 347 - 348 * 349 RCL 04 350 X^2 351 * 352 RCL 07 353 * 354 RCL 11 355 * 356 RCL 06 357 + 358 END |
( 437 bytes / SIZE 016 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | / |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial ellipsoid
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12138.65756 km
---Execution time = 28s---
Notes:
-The root-mean-square error is about 3 cm if 80° <
| L' - L | < 120° and the maximum error is about 16 cm
-The root-mean-square error is about 6 cm if 110° < | L' - L
| < 150° and the maximum error is about 78 cm
-When the difference in the longitudes is 165°, the error can reach
about 4 meters.
-Registers L & R06 contain the great elliptic arc distance.
b) With Andoyer's Formula
of Order 3
-In this variant, the WGS84 geodesic distance is computed with Andoyer's
formula of order 3 ( the other distances with formulas of order 2 )
Data Registers: R00 thru R17: temp
Flags: /
Subroutines: /
01 LBL "TGD"
02 DEG 03 HR 04 STO 15 05 RDN 06 HR 07 14.92911 08 STO 00 09 + 10 STO 14 11 RDN 12 HR 13 STO 13 14 X<>Y 15 HR 16 RCL 00 17 + 18 STO 12 19 6378.172 20 STO 01 21 .07 22 - 23 STO 02 24 21.349686 25 - 26 STO 03 27 RCL 12 28 RCL 13 29 XEQ 02 30 RCL 08 31 STO 05 32 RCL 09 33 STO 06 34 RCL 10 35 STO 07 36 RCL 14 37 RCL 15 38 XEQ 02 39 RCL 05 40 RCL 08 41 * 42 RCL 06 43 RCL 09 44 * 45 RCL 07 46 RCL 10 47 * 48 + 49 + 50 STO 04 51 RCL 01 52 ST/ 05 53 ST/ 08 54 RCL 02 55 ST/ 06 56 ST/ 09 57 RCL 03 58 ST/ 07 59 ST/ 10 60 RCL 05 61 X^2 62 RCL 06 63 X^2 64 RCL 07 65 X^2 66 + 67 + 68 STO 11 69 RCL 08 70 X^2 71 RCL 09 72 X^2 73 RCL 10 74 X^2 75 + 76 + 77 X<> 10 |
78 RCL 07
79 * 80 RCL 06 81 RCL 09 82 * 83 RCL 05 84 RCL 08 85 * 86 + 87 + 88 STO 09 89 ST+ 09 90 RCL 04 91 RCL 11 92 * 93 - 94 ST+ X 95 RCL 04 96 RCL 09 97 * 98 RCL 10 99 - 100 RCL 11 101 RCL 04 102 X^2 103 * 104 - 105 RCL 11 106 RCL 04 107 ACOS 108 STO 04 109 SIN 110 STO 08 111 ST/ Z 112 * 113 + 114 X#0? 115 / 116 ATAN 117 2 118 / 119 90 120 STO 06 121 X<>Y 122 X<0? 123 + 124 STO 05 125 XEQ 03 126 X<> 06 127 RCL 05 128 + 129 XEQ 03 130 STO 07 131 RCL 05 132 CHS 133 ST+ 04 134 1 135 P-R 136 RCL 07 137 RCL 06 138 / 139 STO 09 140 * 141 R-P 142 X<> 04 143 1 144 P-R 145 RCL 09 146 * 147 R-P 148 RDN 149 ST- Z 150 + 151 COS 152 STO 05 153 X^2 154 ST+ X |
155 CHS
156 X<>Y 157 ABS 158 ENTER 159 D-R 160 STO 10 161 SIGN 162 P-R 163 X<>Y 164 ST* 05 165 * 166 ST* Y 167 + 168 RCL 10 169 + 170 4 171 / 172 RCL 07 173 CHS 174 RCL 06 175 ST+ Y 176 ST+ X 177 / 178 STO 09 179 * 180 RCL 05 181 - 182 RCL 10 183 - 184 RCL 09 185 * 186 RCL 10 187 + 188 RCL 06 189 * 190 STO 17 191 GTO 14 192 LBL 02 193 1 194 P-R 195 X<>Y 196 RCL 03 197 X^2 198 * 199 STO 10 200 X^2 201 X<> Z 202 X<>Y 203 P-R 204 RCL 01 205 X^2 206 * 207 STO 08 208 X^2 209 X<>Y 210 RCL 02 211 X^2 212 * 213 STO 09 214 X^2 215 + 216 + 217 SQRT 218 ST/ 08 219 ST/ 09 220 ST/ 10 221 RTN 222 LBL 03 223 ENTER 224 SIN 225 STO 07 226 X^2 227 RCL 10 228 * 229 X<>Y 230 RCL 04 231 - |
232 SIN
233 ST* 07 234 X^2 235 RCL 11 236 * 237 + 238 RCL 07 239 RCL 09 240 * 241 - 242 SQRT 243 RCL 08 244 X<>Y 245 / 246 RTN 247 LBL 14 248 RCL 01 249 RCL 02 250 + 251 2 252 / 253 STO 16 254 RCL 14 255 RCL 12 256 - 257 2 258 / 259 SIN 260 X^2 261 STO 10 262 RCL 13 263 1 264 P-R 265 RCL 03 266 CHS 267 RCL 16 268 ST+ Y 269 / 270 STO 04 271 SIGN 272 LASTX 273 - 274 STO 05 275 / 276 R-P 277 X<> 15 278 1 279 P-R 280 RCL 05 281 / 282 R-P 283 RDN 284 + 285 2 286 ST/ 04 287 / 288 ST- Y 289 1 290 P-R 291 X^2 292 STO 02 293 X<>Y 294 X^2 295 STO 01 296 RDN 297 X<>Y 298 1 299 P-R 300 X^2 301 ST* 01 302 RDN 303 X^2 304 ST* 02 305 STO Z 306 - 307 RCL 10 308 * |
309 +
310 ST/ 02 311 STO 05 312 SQRT 313 ASIN 314 ST+ X 315 STO 03 316 D-R 317 STO 00 318 RCL 02 319 RCL 01 320 1 321 RCL 05 322 - 323 / 324 ST- 02 325 + 326 STO 01 327 RCL 03 328 1 329 P-R 330 STO 06 331 X<>Y 332 STO 05 333 / 334 STO 07 335 ENTER 336 X^2 337 RCL 00 338 * 339 - 340 RCL 00 341 X^2 342 STO 09 343 * 344 STO 03 345 3 346 * 347 RCL 00 348 RCL 09 349 * 350 STO 10 351 .75 352 STO 15 353 / 354 STO 13 355 - 356 RCL 05 357 RCL 06 358 * 359 STO 11 360 LASTX 361 X^2 362 * 363 RCL 11 364 ENTER 365 X^2 366 STO 08 367 * 368 + 369 RCL 05 370 ST* 08 371 X^2 372 X^2 373 STO 12 374 RCL 11 375 * 376 + 377 ST+ X 378 + 379 RCL 11 380 3 381 * 382 RCL 15 383 * 384 - 385 RCL 00 |
386 4
387 / 388 - 389 RCL 01 390 * 391 RCL 03 392 5 393 * 394 - 395 RCL 13 396 ST+ X 397 + 398 RCL 00 399 RCL 11 400 + 401 2 402 / 403 + 404 RCL 05 405 RCL 05 406 RCL 12 407 * 408 - 409 RCL 08 410 - 411 ST+ X 412 RCL 10 413 RCL 05 414 X^2 415 / 416 STO 12 417 RCL 06 418 * 419 STO 13 420 4 421 * 422 - 423 RCL 09 424 RCL 05 425 / 426 STO 14 427 3 428 * 429 + 430 RCL 00 431 RCL 06 432 * 433 STO 08 434 - 435 RCL 06 436 RCL 11 437 * 438 RCL 15 439 ST+ X 440 STO 06 441 * 442 - 443 RCL 02 444 * 445 + 446 RCL 01 447 * 448 RCL 03 449 ST+ X 450 + 451 RCL 10 452 RCL 15 453 / 454 - 455 RCL 13 456 6 457 * 458 RCL 14 459 5 460 * 461 - |
462 RCL 08
463 + 464 RCL 02 465 * 466 + 467 RCL 01 468 * 469 RCL 11 470 2 471 / 472 RCL 00 473 - 474 RCL 12 475 - 476 RCL 01 477 * 478 RCL 05 479 X^2 480 RCL 06 481 / 482 .5 483 - 484 RCL 05 485 * 486 RCL 02 487 * 488 + 489 RCL 11 490 - 491 RCL 12 492 + 493 RCL 00 494 + 495 RCL 02 496 * 497 RCL 14 498 RCL 13 499 - 500 ST+ X 501 + 502 RCL 02 503 * 504 + 505 RCL 04 506 * 507 RCL 05 508 RCL 00 509 / 510 ENTER 511 1/X 512 - 513 RCL 02 514 * 515 RCL 00 516 RCL 07 517 * 518 1 519 - 520 RCL 01 521 * 522 - 523 RCL 01 524 1 525 - 526 * 527 RCL 00 528 * 529 + 530 RCL 04 531 X^2 532 * 533 RCL 16 534 * 535 RCL 17 536 + 537 END |
( 628 bytes / SIZE 018 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | / |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial ellipsoid
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12138.65757 km
---Execution time = 35s---
Notes:
-The root-mean-square error is about 2 cm if 80° <
| L' - L | < 120° and the maximum error is about 13 cm
-The root-mean-square error is about 4 cm if 110° < | L' - L
| < 150° and the maximum error is about 34 cm
-When the difference in the longitudes is 165°, the error can reach
about 75 cm.
-Registers L & R17 contain the great elliptic arc distance.
c) Almost With Andoyer's Formula
of Order 3
-We can save about 100 bytes if we only use the larger terms in the
3rd order formula:
-The precision for very long geodesics are of the same order as the
program above.
Data Registers: R00 thru R15: temp
Flags: /
Subroutines: /
01 LBL "TGD"
02 DEG 03 HR 04 STO 15 05 RDN 06 HR 07 14.92911 08 STO 00 09 + 10 STO 14 11 RDN 12 HR 13 STO 13 14 X<>Y 15 HR 16 RCL 00 17 + 18 STO 12 19 6378.172 20 STO 01 21 .07 22 - 23 STO 02 24 21.349686 25 - 26 STO 03 27 RCL 12 28 RCL 13 29 XEQ 02 30 RCL 08 31 STO 05 32 RCL 09 33 STO 06 34 RCL 10 35 STO 07 36 RCL 14 37 RCL 15 38 XEQ 02 39 RCL 05 40 RCL 08 41 * 42 RCL 06 43 RCL 09 44 * 45 RCL 07 46 RCL 10 47 * 48 + 49 + 50 STO 04 51 RCL 01 52 ST/ 05 53 ST/ 08 54 RCL 02 55 ST/ 06 56 ST/ 09 57 RCL 03 58 ST/ 07 59 ST/ 10 60 RCL 05 61 X^2 62 RCL 06 63 X^2 64 RCL 07 |
65 X^2
66 + 67 + 68 STO 11 69 RCL 08 70 X^2 71 RCL 09 72 X^2 73 RCL 10 74 X^2 75 + 76 + 77 X<> 10 78 RCL 07 79 * 80 RCL 06 81 RCL 09 82 * 83 RCL 05 84 RCL 08 85 * 86 + 87 + 88 STO 09 89 ST+ 09 90 RCL 04 91 RCL 11 92 * 93 - 94 ST+ X 95 RCL 04 96 RCL 09 97 * 98 RCL 10 99 - 100 RCL 11 101 RCL 04 102 X^2 103 * 104 - 105 RCL 11 106 RCL 04 107 ACOS 108 STO 04 109 SIN 110 STO 08 111 ST/ Z 112 * 113 + 114 X#0? 115 / 116 ATAN 117 2 118 / 119 90 120 STO 06 121 X<>Y 122 X<0? 123 + 124 STO 05 125 XEQ 03 126 X<> 06 127 RCL 05 128 + |
129 XEQ 03
130 STO 07 131 RCL 05 132 CHS 133 ST+ 04 134 1 135 P-R 136 RCL 07 137 RCL 06 138 / 139 STO 09 140 * 141 R-P 142 X<> 04 143 1 144 P-R 145 RCL 09 146 * 147 R-P 148 RDN 149 ST- Z 150 + 151 COS 152 STO 05 153 X^2 154 ST+ X 155 CHS 156 X<>Y 157 ABS 158 ENTER 159 D-R 160 STO 10 161 SIGN 162 P-R 163 X<>Y 164 ST* 05 165 * 166 ST* Y 167 + 168 RCL 10 169 + 170 4 171 / 172 RCL 07 173 CHS 174 RCL 06 175 ST+ Y 176 ST+ X 177 / 178 STO 09 179 * 180 RCL 05 181 - 182 RCL 10 183 - 184 RCL 09 185 * 186 RCL 10 187 + 188 ST* 06 189 GTO 14 190 LBL 02 191 1 192 P-R |
193 X<>Y
194 RCL 03 195 X^2 196 * 197 STO 10 198 X^2 199 X<> Z 200 X<>Y 201 P-R 202 RCL 01 203 X^2 204 * 205 STO 08 206 X^2 207 X<>Y 208 RCL 02 209 X^2 210 * 211 STO 09 212 X^2 213 + 214 + 215 SQRT 216 ST/ 08 217 ST/ 09 218 ST/ 10 219 RTN 220 LBL 03 221 ENTER 222 SIN 223 STO 07 224 X^2 225 RCL 10 226 * 227 X<>Y 228 RCL 04 229 - 230 SIN 231 ST* 07 232 X^2 233 RCL 11 234 * 235 + 236 RCL 07 237 RCL 09 238 * 239 - 240 SQRT 241 RCL 08 242 X<>Y 243 / 244 RTN 245 LBL 14 246 RCL 01 247 RCL 02 248 + 249 2 250 / 251 STO 11 252 RCL 14 253 RCL 12 254 - 255 2 256 / |
257 SIN
258 X^2 259 STO 10 260 RCL 13 261 1 262 P-R 263 RCL 03 264 CHS 265 RCL 11 266 ST+ Y 267 / 268 STO 04 269 SIGN 270 LASTX 271 - 272 STO 05 273 / 274 R-P 275 X<> 15 276 1 277 P-R 278 RCL 05 279 / 280 R-P 281 RDN 282 + 283 2 284 ST/ 04 285 / 286 ST- Y 287 1 288 P-R 289 X^2 290 STO 09 291 X<>Y 292 X^2 293 STO 08 294 RDN 295 X<>Y 296 1 297 P-R 298 X^2 299 ST* 08 300 RDN 301 X^2 302 ST* 09 303 STO Z 304 - 305 RCL 10 306 * 307 + 308 ST/ 09 309 STO 05 310 SQRT 311 ASIN 312 ST+ X 313 STO 10 314 D-R 315 STO 07 316 RCL 09 317 RCL 08 318 1 319 RCL 05 320 - |
321 /
322 ST- 09 323 + 324 STO 08 325 RCL 10 326 1 327 P-R 328 STO 12 329 X<>Y 330 STO 05 331 / 332 STO 13 333 ENTER 334 X^2 335 RCL 07 336 * 337 - 338 RCL 07 339 X^2 340 STO 14 341 * 342 STO 10 343 3 344 * 345 RCL 07 346 RCL 14 347 * 348 STO 00 349 .75 350 / 351 STO 15 352 - 353 RCL 08 354 * 355 RCL 10 356 5 357 * 358 - 359 RCL 15 360 ST+ X 361 + 362 RCL 14 363 RCL 05 364 / 365 STO 14 366 3 367 * 368 RCL 00 369 RCL 05 370 X^2 371 / 372 STO 00 373 RCL 12 374 * 375 STO 12 376 4 377 * 378 - 379 RCL 09 380 * 381 + 382 RCL 08 383 * |
384 RCL 10
385 ST+ X 386 + 387 RCL 15 388 - 389 RCL 12 390 6 391 * 392 RCL 14 393 5 394 * 395 - 396 RCL 09 397 * 398 + 399 RCL 08 400 * 401 RCL 00 402 RCL 00 403 LASTX 404 * 405 - 406 RCL 09 407 * 408 RCL 14 409 RCL 12 410 - 411 ST+ X 412 + 413 RCL 09 414 * 415 + 416 RCL 04 417 * 418 RCL 05 419 RCL 07 420 ST* 13 421 / 422 ENTER 423 1/X 424 - 425 RCL 09 426 * 427 RCL 13 428 RCL 08 429 ST* Y 430 - 431 - 432 RCL 08 433 1 434 - 435 * 436 RCL 07 437 * 438 + 439 RCL 04 440 X^2 441 * 442 RCL 11 443 * 444 RCL 06 445 + 446 END |
( 530 bytes / SIZE 016 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | / |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial ellipsoid
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12138.65756 km
---Execution time = 31s---
Notes:
-The root-mean-square error is about 3 cm if 80° <
| L' - L | < 120° and the maximum error is about 15 cm
-The root-mean-square error is about 6 cm if 110° < | L' - L
| < 150° and the maximum error is about 34 cm
-When the difference in the longitudes is 165°, the error can reach
about 78 cm.
-I've found these results with free42binary and with latitudes between
-87° & +87°
-Registers L & R06 contain the great elliptic arc distance.
d)
Geodetic or Geocentric Longitudes & Latitudes ?
-This program employs the same formulas as above, but you can choose
the semi-axes of the Earth, and geodetic or geocentric longitudes and latitudes.
Data Registers: (•) R00 = L0 ( Registers R00 thru R04 are to be initialized IF Flag F00 is set )
(•) R01 = a
R04 to R16: temp
(•) R02 = b
(•) R03 = c
Flags: F00-F03-F04
If F00 is clear, "TGD" uses a = 6378.172 km , b = 6378.102
km , c = 6356.752314 km , L0 = -14°92911
If F00 is set , initialize R00 to R03
If F04 is clear & F03 is clear , longitudes and latitudes are
geodetic
If F04 is clear & F03 is set , longitudes are geocentric but
latitudes are geocentric
If F04 is set, longitudes and latitudes are geocentric.
Subroutines: /
-Line 190 is a three-byte GTO 14
01 LBL "TGD"
02 DEG 03 HR 04 STO 15 05 RDN 06 HR 07 STO 14 08 RDN 09 HR 10 STO 13 11 X<>Y 12 HR 13 STO 12 14 FS? 00 15 GTO 00 16 14.92911 17 CHS 18 STO 00 19 6378.172 20 STO 01 21 .07 22 - 23 STO 02 24 21.349686 25 - 26 STO 03 27 LBL 00 28 RCL 13 29 RCL 12 30 XEQ 02 31 RCL 08 32 STO 05 33 RCL 09 34 STO 06 35 RCL 10 36 STO 07 37 RCL 15 38 RCL 14 39 XEQ 02 40 RCL 05 41 RCL 08 42 * 43 RCL 06 44 RCL 09 45 * 46 RCL 07 47 RCL 10 48 * 49 + 50 + 51 STO 04 52 RCL 01 53 ST/ 05 54 ST/ 08 55 RCL 02 56 ST/ 06 57 ST/ 09 58 RCL 03 59 ST/ 07 60 ST/ 10 61 RCL 05 62 X^2 63 RCL 06 64 X^2 65 RCL 07 66 X^2 67 + 68 + 69 STO 11 70 RCL 08 71 X^2 |
72 RCL 09
73 X^2 74 RCL 10 75 X^2 76 + 77 + 78 X<> 10 79 RCL 07 80 * 81 RCL 06 82 RCL 09 83 * 84 RCL 05 85 RCL 08 86 * 87 + 88 + 89 STO 09 90 ST+ 09 91 RCL 04 92 RCL 11 93 * 94 - 95 ST+ X 96 RCL 04 97 RCL 09 98 * 99 RCL 10 100 - 101 RCL 11 102 RCL 04 103 X^2 104 * 105 - 106 RCL 11 107 RCL 04 108 ACOS 109 STO 04 110 SIN 111 STO 08 112 ST/ Z 113 * 114 + 115 X#0? 116 / 117 ATAN 118 2 119 / 120 90 121 STO 06 122 X<>Y 123 X<0? 124 + 125 STO 05 126 XEQ 03 127 X<> 06 128 RCL 05 129 + 130 XEQ 03 131 STO 07 132 RCL 05 133 CHS 134 ST+ 04 135 1 136 P-R 137 RCL 07 138 RCL 06 139 / 140 STO 09 141 * 142 R-P |
143 X<> 04
144 1 145 P-R 146 RCL 09 147 * 148 R-P 149 RDN 150 ST- Z 151 + 152 COS 153 STO 05 154 X^2 155 ST+ X 156 CHS 157 X<>Y 158 ABS 159 ENTER 160 D-R 161 STO 10 162 SIGN 163 P-R 164 X<>Y 165 ST* 05 166 * 167 ST* Y 168 + 169 RCL 10 170 + 171 4 172 / 173 RCL 07 174 CHS 175 RCL 06 176 ST+ Y 177 ST+ X 178 / 179 STO 09 180 * 181 RCL 05 182 - 183 RCL 10 184 - 185 RCL 09 186 * 187 RCL 10 188 + 189 ST* 06 190 GTO 14 191 LBL 02 192 RCL 00 193 - 194 FS? 04 195 GTO 07 196 FS? 03 197 GTO 06 198 X<>Y 199 1 200 P-R 201 X<>Y 202 RCL 03 203 X^2 204 * 205 STO 10 206 X^2 207 X<> Z 208 X<>Y 209 P-R 210 RCL 01 211 X^2 212 * 213 STO 08 |
214 X^2
215 X<>Y 216 RCL 02 217 X^2 218 * 219 STO 09 220 X^2 221 + 222 + 223 SQRT 224 ST/ 08 225 ST/ 09 226 ST/ 10 227 RTN 228 LBL 06 229 X<>Y 230 STO 11 231 X<>Y 232 STO 04 233 1 234 P-R 235 X^2 236 X<>Y 237 RCL 01 238 RCL 02 239 / 240 X^2 241 * 242 X^2 243 + 244 SQRT 245 RCL 03 246 RCL 01 247 / 248 X^2 249 * 250 RCL 11 251 1 252 P-R 253 RCL Z 254 / 255 R-P 256 X<> 04 257 LBL 07 258 1 259 X<>Y 260 RDN 261 P-R 262 R^ 263 X<>Y 264 P-R 265 STO 08 266 RDN 267 STO 09 268 X<>Y 269 STO 10 270 RTN 271 LBL 03 272 ENTER 273 SIN 274 STO 07 275 X^2 276 RCL 10 277 * 278 X<>Y 279 RCL 04 280 - 281 SIN 282 ST* 07 283 X^2 284 RCL 11 |
285 *
286 + 287 RCL 07 288 RCL 09 289 * 290 - 291 SQRT 292 RCL 08 293 X<>Y 294 / 295 RTN 296 LBL 14 297 RCL 14 298 RCL 12 299 - 300 2 301 / 302 SIN 303 X^2 304 STO 10 305 RCL 01 306 RCL 02 307 + 308 2 309 / 310 STO 16 311 RCL 13 312 1 313 P-R 314 RCL 03 315 CHS 316 R^ 317 ST+ Y 318 / 319 STO 04 320 SIGN 321 LASTX 322 - 323 STO 05 324 / 325 R-P 326 X<> 15 327 1 328 P-R 329 RCL 05 330 / 331 R-P 332 RDN 333 + 334 2 335 ST/ 04 336 / 337 ST- Y 338 1 339 P-R 340 X^2 341 STO 09 342 X<>Y 343 X^2 344 STO 08 345 RDN 346 X<>Y 347 1 348 P-R 349 X^2 350 ST* 08 351 RDN 352 X^2 353 ST* 09 354 STO Z 355 - |
356 RCL 10
357 * 358 + 359 ST/ 09 360 STO 05 361 SQRT 362 ASIN 363 ST+ X 364 STO 10 365 D-R 366 STO 07 367 RCL 09 368 RCL 08 369 1 370 RCL 05 371 - 372 / 373 ST- 09 374 + 375 STO 08 376 RCL 10 377 1 378 P-R 379 STO 12 380 X<>Y 381 STO 05 382 / 383 STO 13 384 ENTER 385 X^2 386 RCL 07 387 * 388 - 389 RCL 07 390 X^2 391 STO 14 392 * 393 STO 10 394 3 395 * 396 RCL 07 397 RCL 14 398 * 399 STO 11 400 .75 401 / 402 STO 15 403 - 404 RCL 08 405 * 406 RCL 10 407 5 408 * 409 - 410 RCL 15 411 ST+ X 412 + 413 RCL 14 414 RCL 05 415 / 416 STO 14 417 3 418 * 419 RCL 11 420 RCL 05 421 X^2 422 / 423 STO 11 424 RCL 12 425 * 426 STO 12 |
427 4
428 * 429 - 430 RCL 09 431 * 432 + 433 RCL 08 434 * 435 RCL 10 436 ST+ X 437 + 438 RCL 15 439 - 440 RCL 12 441 6 442 * 443 RCL 14 444 5 445 * 446 - 447 RCL 09 448 * 449 + 450 RCL 08 451 * 452 RCL 11 453 RCL 11 454 LASTX 455 * 456 - 457 RCL 09 458 * 459 RCL 14 460 RCL 12 461 - 462 ST+ X 463 + 464 RCL 09 465 * 466 + 467 RCL 04 468 * 469 RCL 05 470 RCL 07 471 ST* 13 472 / 473 ENTER 474 1/X 475 - 476 RCL 09 477 * 478 RCL 13 479 RCL 08 480 ST* Y 481 - 482 - 483 RCL 08 484 1 485 - 486 * 487 RCL 07 488 * 489 + 490 RCL 04 491 X^2 492 * 493 RCL 16 494 * 495 RCL 06 496 + 497 END |
( 592 bytes
/ SIZE 017 )
STACK | INPUTS | OUTPUTS |
T | L ( ° ' " ) | / |
Z | b ( ° ' " ) | / |
Y | L' ( ° ' " ) | / |
X | b' ( ° ' " ) | D ( km ) |
Where D is the geodesic distance on the triaxial ellipsoid
Example: Calculate the geodesic distance
between the Sydney Observatory ( Long = 151°12'17"8 E , Lat = 33°51'41"1
S )
and the Palomar Observatory ( Long = 116°51'50"4 W , Lat = 33°21'22"4
N )
• If CF 00 CF 03 CF 04 geodetic longitudes & latitudes
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12138.65756 km
LASTX d = 12138.65876 km = R06 = Great elliptic arc distance
• If CF 00 SF 03 CF 04 geocentric longitudes & geodetic latitudes
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12138.70239 km
LASTX d = 12138.70359 km = R06 = Great elliptic arc distance
• If CF 00 SF 04 geocentric longitudes & latitudes
151.12178 ENTER^
-33.51411 ENTER^
-116.51504 ENTER^
33.21224 XEQ "TGD"
>>>> D = 12157.22639 km
LASTX d = 12157.22759 km = R06 = Great elliptic arc distance
Notes:
-The latitudes given on the web are always geodetic, but I don't
really know if the longitudes are geocentric or geodetic !
-It's the same for an ellipsoid of revolution, but not for a triaxial
ellipsoid...
References:
[1] Paul D. Thomas - "Spheroidal Geodesics, Reference Systems
& Local Geometry" - US Naval Oceanographic Office.
[2] Richard H. Rapp - "Geometric Geodesy, Part II - The Ohio
State University
[3] Emanuel M. Sodano & Thelma A. Robinson - "Direct and
Inverse Solutions of Geodesics" - Army Map Service, Technical Report n°7
-Many thanks to Professor Richard H. Rapp who sent me reference [3]
to find the typo in Sodano's formula !