15 double precision tstdat(20, 12)
16 common /tstcom/ tstdat
17 data (tstdat(1,j), j = 1,12) /
18 + 35.60777d0,-139.44815d0,111.098748429560326d0,
19 + -11.17491d0,-69.95921d0,129.289270889708762d0,
20 + 8935244.5604818305d0,80.50729714281974d0,6273170.2055303837d0,
21 + 0.16606318447386067d0,0.16479116945612937d0,
22 + 12841384694976.432d0 /
23 data (tstdat(2,j), j = 1,12) /
24 + 55.52454d0,106.05087d0,22.020059880982801d0,
25 + 77.03196d0,197.18234d0,109.112041110671519d0,
26 + 4105086.1713924406d0,36.892740690445894d0,
27 + 3828869.3344387607d0,
28 + 0.80076349608092607d0,0.80101006984201008d0,
29 + 61674961290615.615d0 /
30 data (tstdat(3,j), j = 1,12) /
31 + -21.97856d0,142.59065d0,-32.44456876433189d0,
32 + 41.84138d0,98.56635d0,-41.84359951440466d0,
33 + 8394328.894657671d0,75.62930491011522d0,6161154.5773110616d0,
34 + 0.24816339233950381d0,0.24930251203627892d0,
35 + -6637997720646.717d0 /
36 data (tstdat(4,j), j = 1,12) /
37 + -66.99028d0,112.2363d0,173.73491240878403d0,
38 + -12.70631d0,285.90344d0,2.512956620913668d0,
39 + 11150344.2312080241d0,100.278634181155759d0,
40 + 6289939.5670446687d0,
41 + -0.17199490274700385d0,-0.17722569526345708d0,
42 + -121287239862139.744d0 /
43 data (tstdat(5,j), j = 1,12) /
44 + -17.42761d0,173.34268d0,-159.033557661192928d0,
45 + -15.84784d0,5.93557d0,-20.787484651536988d0,
46 + 16076603.1631180673d0,144.640108810286253d0,
47 + 3732902.1583877189d0,
48 + -0.81273638700070476d0,-0.81299800519154474d0,
49 + 97825992354058.708d0 /
50 data (tstdat(6,j), j = 1,12) /
51 + 32.84994d0,48.28919d0,150.492927788121982d0,
52 + -56.28556d0,202.29132d0,48.113449399816759d0,
53 + 16727068.9438164461d0,150.565799985466607d0,
54 + 3147838.1910180939d0,
55 + -0.87334918086923126d0,-0.86505036767110637d0,
56 + -72445258525585.010d0 /
57 data (tstdat(7,j), j = 1,12) /
58 + 6.96833d0,52.74123d0,92.581585386317712d0,
59 + -7.39675d0,206.17291d0,90.721692165923907d0,
60 + 17102477.2496958388d0,154.147366239113561d0,
61 + 2772035.6169917581d0,
62 + -0.89991282520302447d0,-0.89986892177110739d0,
63 + -1311796973197.995d0 /
64 data (tstdat(8,j), j = 1,12) /
65 + -50.56724d0,-16.30485d0,-105.439679907590164d0,
66 + -33.56571d0,-94.97412d0,-47.348547835650331d0,
67 + 6455670.5118668696d0,58.083719495371259d0,
68 + 5409150.7979815838d0,
69 + 0.53053508035997263d0,0.52988722644436602d0,
70 + 41071447902810.047d0 /
71 data (tstdat(9,j), j = 1,12) /
72 + -58.93002d0,-8.90775d0,140.965397902500679d0,
73 + -8.91104d0,133.13503d0,19.255429433416599d0,
74 + 11756066.0219864627d0,105.755691241406877d0,
75 + 6151101.2270708536d0,
76 + -0.26548622269867183d0,-0.27068483874510741d0,
77 + -86143460552774.735d0 /
78 data (tstdat(10,j), j = 1,12) /
79 + -68.82867d0,-74.28391d0,93.774347763114881d0,
80 + -50.63005d0,-8.36685d0,34.65564085411343d0,
81 + 3956936.926063544d0,35.572254987389284d0,3708890.9544062657d0,
82 + 0.81443963736383502d0,0.81420859815358342d0,
83 + -41845309450093.787d0 /
84 data (tstdat(11,j), j = 1,12) /
85 + -10.62672d0,-32.0898d0,-86.426713286747751d0,
86 + 5.883d0,-134.31681d0,-80.473780971034875d0,
87 + 11470869.3864563009d0,103.387395634504061d0,
88 + 6184411.6622659713d0,
89 + -0.23138683500430237d0,-0.23155097622286792d0,
90 + 4198803992123.548d0 /
91 data (tstdat(12,j), j = 1,12) /
92 + -21.76221d0,166.90563d0,29.319421206936428d0,
93 + 48.72884d0,213.97627d0,43.508671946410168d0,
94 + 9098627.3986554915d0,81.963476716121964d0,
95 + 6299240.9166992283d0,
96 + 0.13965943368590333d0,0.14152969707656796d0,
97 + 10024709850277.476d0 /
98 data (tstdat(13,j), j = 1,12) /
99 + -19.79938d0,-174.47484d0,71.167275780171533d0,
100 + -11.99349d0,-154.35109d0,65.589099775199228d0,
101 + 2319004.8601169389d0,20.896611684802389d0,
102 + 2267960.8703918325d0,
103 + 0.93427001867125849d0,0.93424887135032789d0,
104 + -3935477535005.785d0 /
105 data (tstdat(14,j), j = 1,12) /
106 + -11.95887d0,-116.94513d0,92.712619830452549d0,
107 + 4.57352d0,7.16501d0,78.64960934409585d0,
108 + 13834722.5801401374d0,124.688684161089762d0,
109 + 5228093.177931598d0,
110 + -0.56879356755666463d0,-0.56918731952397221d0,
111 + -9919582785894.853d0 /
112 data (tstdat(15,j), j = 1,12) /
113 + -87.85331d0,85.66836d0,-65.120313040242748d0,
114 + 66.48646d0,16.09921d0,-4.888658719272296d0,
115 + 17286615.3147144645d0,155.58592449699137d0,
116 + 2635887.4729110181d0,
117 + -0.90697975771398578d0,-0.91095608883042767d0,
118 + 42667211366919.534d0 /
119 data (tstdat(16,j), j = 1,12) /
120 + 1.74708d0,128.32011d0,-101.584843631173858d0,
121 + -11.16617d0,11.87109d0,-86.325793296437476d0,
122 + 12942901.1241347408d0,116.650512484301857d0,
123 + 5682744.8413270572d0,
124 + -0.44857868222697644d0,-0.44824490340007729d0,
125 + 10763055294345.653d0 /
126 data (tstdat(17,j), j = 1,12) /
127 + -25.72959d0,-144.90758d0,-153.647468693117198d0,
128 + -57.70581d0,-269.17879d0,-48.343983158876487d0,
129 + 9413446.7452453107d0,84.664533838404295d0,
130 + 6356176.6898881281d0,
131 + 0.09492245755254703d0,0.09737058264766572d0,
132 + 74515122850712.444d0 /
133 data (tstdat(18,j), j = 1,12) /
134 + -41.22777d0,122.32875d0,14.285113402275739d0,
135 + -7.57291d0,130.37946d0,10.805303085187369d0,
136 + 3812686.035106021d0,34.34330804743883d0,3588703.8812128856d0,
137 + 0.82605222593217889d0,0.82572158200920196d0,
138 + -2456961531057.857d0 /
139 data (tstdat(19,j), j = 1,12) /
140 + 11.01307d0,138.25278d0,79.43682622782374d0,
141 + 6.62726d0,247.05981d0,103.708090215522657d0,
142 + 11911190.819018408d0,107.341669954114577d0,
143 + 6070904.722786735d0,
144 + -0.29767608923657404d0,-0.29785143390252321d0,
145 + 17121631423099.696d0 /
146 data (tstdat(20,j), j = 1,12) /
147 + -29.47124d0,95.14681d0,-163.779130441688382d0,
148 + -27.46601d0,-69.15955d0,-15.909335945554969d0,
149 + 13487015.8381145492d0,121.294026715742277d0,
150 + 5481428.9945736388d0,
151 + -0.51527225545373252d0,-0.51556587964721788d0,
152 + 104679964020340.318d0 /
155 integer function assert(x, y, d)
156 double precision x, y, d
158 if (abs(x - y) .le. d)
then
163 10
format(1x,
'assert fails: ',
164 + g14.7,
' != ', g14.7,
' +/- ', g10.3)
170 integer function chknan(x)
182 integer function tstinv()
183 double precision tstdat(20, 12)
184 common /tstcom/ tstdat
185 double precision lat1, lon1, azi1, lat2, lon2, azi2,
186 + s12, a12, m12, MM12, MM21, SS12
187 double precision azi1a, azi2a, s12a, a12a,
188 + m12a, MM12a, MM21a, SS12a
189 double precision a, f
190 integer r, assert, i, omask
191 include
'geodesic.inc'
195 f = 1/298.257223563d0
196 omask = 1 + 2 + 4 + 8
212 call invers(a, f, lat1, lon1, lat2, lon2,
213 + s12a, azi1a, azi2a, omask, a12a, m12a, mm12a, mm21a, ss12a)
214 r = r + assert(azi1, azi1a, 1d-13)
215 r = r + assert(azi2, azi2a, 1d-13)
216 r = r + assert(s12, s12a, 1d-8)
217 r = r + assert(a12, a12a, 1d-13)
218 r = r + assert(m12, m12a, 1d-8)
219 r = r + assert(mm12, mm12a, 1d-15)
220 r = r + assert(mm21, mm21a, 1d-15)
221 r = r + assert(ss12, ss12a, 0.1d0)
228 integer function tstdir()
229 double precision tstdat(20, 12)
230 common /tstcom/ tstdat
231 double precision lat1, lon1, azi1, lat2, lon2, azi2,
232 + s12, a12, m12, MM12, MM21, SS12
233 double precision lat2a, lon2a, azi2a, a12a,
234 + m12a, MM12a, MM21a, SS12a
235 double precision a, f
236 integer r, assert, i, omask, flags
237 include
'geodesic.inc'
241 f = 1/298.257223563d0
242 omask = 1 + 2 + 4 + 8
259 call direct(a, f, lat1, lon1, azi1, s12, flags,
260 + lat2a, lon2a, azi2a, omask, a12a, m12a, mm12a, mm21a, ss12a)
261 r = r + assert(lat2, lat2a, 1d-13)
262 r = r + assert(lon2, lon2a, 1d-13)
263 r = r + assert(azi2, azi2a, 1d-13)
264 r = r + assert(a12, a12a, 1d-13)
265 r = r + assert(m12, m12a, 1d-8)
266 r = r + assert(mm12, mm12a, 1d-15)
267 r = r + assert(mm21, mm21a, 1d-15)
268 r = r + assert(ss12, ss12a, 0.1d0)
275 integer function tstarc()
276 double precision tstdat(20, 12)
277 common /tstcom/ tstdat
278 double precision lat1, lon1, azi1, lat2, lon2, azi2,
279 + s12, a12, m12, MM12, MM21, SS12
280 double precision lat2a, lon2a, azi2a, s12a,
281 + m12a, MM12a, MM21a, SS12a
282 double precision a, f
283 integer r, assert, i, omask, flags
284 include
'geodesic.inc'
288 f = 1/298.257223563d0
289 omask = 1 + 2 + 4 + 8
306 call direct(a, f, lat1, lon1, azi1, a12, flags,
307 + lat2a, lon2a, azi2a, omask, s12a, m12a, mm12a, mm21a, ss12a)
308 r = r + assert(lat2, lat2a, 1d-13)
309 r = r + assert(lon2, lon2a, 1d-13)
310 r = r + assert(azi2, azi2a, 1d-13)
311 r = r + assert(s12, s12a, 1d-8)
312 r = r + assert(m12, m12a, 1d-8)
313 r = r + assert(mm12, mm12a, 1d-15)
314 r = r + assert(mm21, mm21a, 1d-15)
315 r = r + assert(ss12, ss12a, 0.1d0)
322 integer function tstg0()
323 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
324 double precision a, f
325 integer r, assert, omask
326 include
'geodesic.inc'
330 f = 1/298.257223563d0
333 call invers(a, f, 40.6d0, -73.8d0, 49.01666667d0, 2.55d0,
334 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
335 r = r + assert(azi1, 53.47022d0, 0.5d-5)
336 r = r + assert(azi2, 111.59367d0, 0.5d-5)
337 r = r + assert(s12, 5853226d0, 0.5d0)
343 integer function tstg1()
344 double precision lat2, lon2, azi2, a12, m12, MM12, MM21, SS12
345 double precision a, f
346 integer r, assert, omask, flags
347 include
'geodesic.inc'
351 f = 1/298.257223563d0
355 call direct(a, f, 40.63972222d0, -73.77888889d0, 53.5d0, 5850d3,
356 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
357 r = r + assert(lat2, 49.01467d0, 0.5d-5)
358 r = r + assert(lon2, 2.56106d0, 0.5d-5)
359 r = r + assert(azi2, 111.62947d0, 0.5d-5)
365 integer function tstg2()
367 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
368 double precision a, f
369 integer r, assert, omask
370 include
'geodesic.inc'
376 call invers(a, f, 0.07476d0, 0d0, -0.07476d0, 180d0,
377 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
378 r = r + assert(azi1, 90.00078d0, 0.5d-5)
379 r = r + assert(azi2, 90.00078d0, 0.5d-5)
380 r = r + assert(s12, 20106193d0, 0.5d0)
381 call invers(a, f, 0.1d0, 0d0, -0.1d0, 180d0,
382 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
383 r = r + assert(azi1, 90.00105d0, 0.5d-5)
384 r = r + assert(azi2, 90.00105d0, 0.5d-5)
385 r = r + assert(s12, 20106193d0, 0.5d0)
391 integer function tstg4()
393 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
394 double precision a, f
395 integer r, assert, omask
396 include
'geodesic.inc'
400 f = 1/298.257223563d0
404 + 36.493349428792d0, 0d0, 36.49334942879201d0, .0000008d0,
405 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
406 r = r + assert(s12, 0.072d0, 0.5d-3)
412 integer function tstg5()
413 double precision lat2, lon2, azi2, a12, m12, MM12, MM21, SS12
414 double precision a, f
415 integer r, assert, omask, flags
416 include
'geodesic.inc'
420 f = 1/298.257223563d0
424 call direct(a, f, 0.01777745589997d0, 30d0, 0d0, 10d6,
425 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
426 if (lon2 .lt. 0)
then
427 r = r + assert(lon2, -150d0, 0.5d-5)
428 r = r + assert(abs(azi2), 180d0, 0.5d-5)
430 r = r + assert(lon2, 30d0, 0.5d-5)
431 r = r + assert(azi2, 0d0, 0.5d-5)
438 integer function tstg6()
441 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
442 double precision a, f
443 integer r, assert, omask
444 include
'geodesic.inc'
448 f = 1/298.257223563d0
451 call invers(a, f, 88.202499451857d0, 0d0,
452 + -88.202499451857d0, 179.981022032992859592d0,
453 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
454 r = r + assert(s12, 20003898.214d0, 0.5d-3)
455 call invers(a, f, 89.262080389218d0, 0d0,
456 + -89.262080389218d0, 179.992207982775375662d0,
457 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
458 r = r + assert(s12, 20003925.854d0, 0.5d-3)
459 call invers(a, f, 89.333123580033d0, 0d0,
460 + -89.333123580032997687d0, 179.99295812360148422d0,
461 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
462 r = r + assert(s12, 20003926.881d0, 0.5d-3)
468 integer function tstg9()
470 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
471 double precision a, f
472 integer r, assert, omask
473 include
'geodesic.inc'
477 f = 1/298.257223563d0
480 call invers(a, f, 56.320923501171d0, 0d0,
481 + -56.320923501171d0, 179.664747671772880215d0,
482 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
483 r = r + assert(s12, 19993558.287d0, 0.5d-3)
489 integer function tstg10()
492 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
493 double precision a, f
494 integer r, assert, omask
495 include
'geodesic.inc'
499 f = 1/298.257223563d0
502 call invers(a, f, 52.784459512564d0, 0d0,
503 + -52.784459512563990912d0, 179.634407464943777557d0,
504 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
505 r = r + assert(s12, 19991596.095d0, 0.5d-3)
511 integer function tstg11()
514 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
515 double precision a, f
516 integer r, assert, omask
517 include
'geodesic.inc'
521 f = 1/298.257223563d0
524 call invers(a, f, 48.522876735459d0, 0d0,
525 + -48.52287673545898293d0, 179.599720456223079643d0,
526 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
527 r = r + assert(s12, 19989144.774d0, 0.5d-3)
533 integer function tstg12()
537 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
538 double precision a, f
539 integer r, assert, omask
540 include
'geodesic.inc'
546 call invers(a, f, 0d0, 0d0, -10d0, 160d0,
547 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
548 r = r + assert(azi1, 120.27d0, 1d-2)
549 r = r + assert(azi2, 105.15d0, 1d-2)
550 r = r + assert(s12, 266.7d0, 1d-1)
556 integer function tstg14()
558 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
559 double precision a, f, LatFix
560 integer r, chknan, omask
561 include
'geodesic.inc'
565 f = 1/298.257223563d0
568 call invers(a, f, 0d0, 0d0, 1d0, latfix(91d0),
569 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
578 integer function tstg15()
581 double precision lat2, lon2, azi2, a12, m12, MM12, MM21, SS12
582 double precision a, f
583 integer r, assert, omask, flags
584 include
'geodesic.inc'
591 call direct(a, f, 1d0, 2d0, 3d0, 4d0,
592 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
593 r = r + assert(ss12, 23700d0, 0.5d0)
599 integer function tstg17()
601 double precision lat2, lon2, azi2, a12, m12, MM12, MM21, SS12
602 double precision a, f
603 integer r, assert, omask, flags
604 include
'geodesic.inc'
608 f = 1/298.257223563d0
612 call direct(a, f, 40d0, -75d0, -10d0, 2d7,
613 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
614 r = r + assert(lat2, -39d0, 1d0)
615 r = r + assert(lon2, -254d0, 1d0)
616 r = r + assert(azi2, -170d0, 1d0)
618 call direct(a, f, 40d0, -75d0, -10d0, 2d7,
619 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
620 r = r + assert(lat2, -39d0, 1d0)
621 r = r + assert(lon2, 105d0, 1d0)
622 r = r + assert(azi2, -170d0, 1d0)
628 integer function tstg26()
630 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
631 double precision a, f
632 integer r, assert, omask
633 include
'geodesic.inc'
639 call invers(a, f, 1d0, 2d0, 3d0, 4d0,
640 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
641 r = r + assert(ss12, 49911046115.0d0, 0.5d0)
647 integer function tstg28()
649 double precision lat2, lon2, azi2, a12, m12, MM12, MM21, SS12
650 double precision a, f
651 integer r, assert, omask, flags
652 include
'geodesic.inc'
656 omask = 1 + 2 + 4 + 8
659 call direct(a, f, 1d0, 2d0, 10d0, 5d6,
660 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
661 r = r + assert(a12, 48.55570690d0, 0.5d-8)
667 integer function tstg33()
671 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
672 double precision a, f
673 integer r, assert, omask
674 include
'geodesic.inc'
678 f = 1/298.257223563d0
681 call invers(a, f, 0d0, 0d0, 0d0, 179d0,
682 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
683 r = r + assert(azi1, 90.00000d0, 0.5d-5)
684 r = r + assert(azi2, 90.00000d0, 0.5d-5)
685 r = r + assert(s12, 19926189d0, 0.5d0)
686 call invers(a, f, 0d0, 0d0, 0d0, 179.5d0,
687 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
688 r = r + assert(azi1, 55.96650d0, 0.5d-5)
689 r = r + assert(azi2, 124.03350d0, 0.5d-5)
690 r = r + assert(s12, 19980862d0, 0.5d0)
691 call invers(a, f, 0d0, 0d0, 0d0, 180d0,
692 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
693 r = r + assert(azi1, 0.00000d0, 0.5d-5)
694 r = r + assert(abs(azi2), 180.00000d0, 0.5d-5)
695 r = r + assert(s12, 20003931d0, 0.5d0)
696 call invers(a, f, 0d0, 0d0, 1d0, 180d0,
697 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
698 r = r + assert(azi1, 0.00000d0, 0.5d-5)
699 r = r + assert(abs(azi2), 180.00000d0, 0.5d-5)
700 r = r + assert(s12, 19893357d0, 0.5d0)
703 call invers(a, f, 0d0, 0d0, 0d0, 179d0,
704 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
705 r = r + assert(azi1, 90.00000d0, 0.5d-5)
706 r = r + assert(azi2, 90.00000d0, 0.5d-5)
707 r = r + assert(s12, 19994492d0, 0.5d0)
708 call invers(a, f, 0d0, 0d0, 0d0, 180d0,
709 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
710 r = r + assert(azi1, 0.00000d0, 0.5d-5)
711 r = r + assert(abs(azi2), 180.00000d0, 0.5d-5)
712 r = r + assert(s12, 20106193d0, 0.5d0)
713 call invers(a, f, 0d0, 0d0, 1d0, 180d0,
714 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
715 r = r + assert(azi1, 0.00000d0, 0.5d-5)
716 r = r + assert(abs(azi2), 180.00000d0, 0.5d-5)
717 r = r + assert(s12, 19994492d0, 0.5d0)
720 call invers(a, f, 0d0, 0d0, 0d0, 179d0,
721 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
722 r = r + assert(azi1, 90.00000d0, 0.5d-5)
723 r = r + assert(azi2, 90.00000d0, 0.5d-5)
724 r = r + assert(s12, 19994492d0, 0.5d0)
725 call invers(a, f, 0d0, 0d0, 0d0, 180d0,
726 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
727 r = r + assert(azi1, 90.00000d0, 0.5d-5)
728 r = r + assert(azi2, 90.00000d0, 0.5d-5)
729 r = r + assert(s12, 20106193d0, 0.5d0)
730 call invers(a, f, 0d0, 0d0, 0.5d0, 180d0,
731 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
732 r = r + assert(azi1, 33.02493d0, 0.5d-5)
733 r = r + assert(azi2, 146.97364d0, 0.5d-5)
734 r = r + assert(s12, 20082617d0, 0.5d0)
735 call invers(a, f, 0d0, 0d0, 1d0, 180d0,
736 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
737 r = r + assert(azi1, 0.00000d0, 0.5d-5)
738 r = r + assert(abs(azi2), 180.00000d0, 0.5d-5)
739 r = r + assert(s12, 20027270d0, 0.5d0)
745 integer function tstg55()
748 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
749 double precision a, f
750 integer r, chknan, omask
751 include
'geodesic.inc'
755 f = 1/298.257223563d0
758 call invers(a, f, 91d0, 0d0, 0d0, 90d0,
759 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
763 call invers(a, f, 91d0, 0d0, 90d0, 9d0,
764 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
773 integer function tstg59()
775 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
776 double precision a, f
777 integer r, assert, omask
778 include
'geodesic.inc'
782 f = 1/298.257223563d0
785 call invers(a, f, 5d0, 0.00000000000001d0, 10d0, 180d0,
786 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
787 r = r + assert(azi1, 0.000000000000035d0, 1.5d-14)
788 r = r + assert(azi2, 179.99999999999996d0, 1.5d-14)
789 r = r + assert(s12, 18345191.174332713d0, 5d-9)
795 integer function tstg61()
797 double precision lat2, lon2, azi2, a12, m12, MM12, MM21, SS12
798 double precision a, f
799 integer r, assert, omask, flags
800 include
'geodesic.inc'
804 f = 1/298.257223563d0
808 call direct(a, f, 45d0, 0d0, -0.000000000000000003d0, 1d7,
809 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
810 r = r + assert(lat2, 45.30632d0, 0.5d-5)
811 r = r + assert(lon2, -180d0, 0.5d-5)
812 r = r + assert(abs(azi2), 180d0, 0.5d-5)
818 integer function tstg73()
824 double precision lat2, lon2, azi2, a12, m12, MM12, MM21, SS12
825 double precision a, f
826 integer r, assert, omask, flags
827 include
'geodesic.inc'
831 f = 1/298.257223563d0
835 call direct(a, f, 90d0, 10d0, 180d0, -1d6,
836 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
837 r = r + assert(lat2, 81.04623d0, 0.5d-5)
838 r = r + assert(lon2, -170d0, 0.5d-5)
839 r = r + assert(azi2, 0d0, 0d0)
840 r = r + assert(sign(1d0, azi2), 1d0, 0d0)
846 integer function tstg74()
849 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
850 double precision a, f
851 integer r, assert, omask
852 include
'geodesic.inc'
856 f = 1/298.257223563d0
857 omask = 1 + 2 + 4 + 8
859 call invers(a, f, 54.1589d0, 15.3872d0, 54.1591d0, 15.3877d0,
860 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
861 r = r + assert(azi1, 55.723110355d0, 5d-9)
862 r = r + assert(azi2, 55.723515675d0, 5d-9)
863 r = r + assert(s12, 39.527686385d0, 5d-9)
864 r = r + assert(a12, 0.000355495d0, 5d-9)
865 r = r + assert(m12, 39.527686385d0, 5d-9)
866 r = r + assert(mm12, 0.999999995d0, 5d-9)
867 r = r + assert(mm21, 0.999999995d0, 5d-9)
868 r = r + assert(ss12, 286698586.30197d0, 5d-4)
874 integer function tstg76()
877 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
878 double precision a, f
879 integer r, assert, omask
880 include
'geodesic.inc'
884 f = 1/298.257223563d0
888 + -(41+19/60d0), 174+49/60d0, 40+58/60d0, -(5+30/60d0),
889 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
890 r = r + assert(azi1, 160.39137649664d0, 0.5d-11)
891 r = r + assert(azi2, 19.50042925176d0, 0.5d-11)
892 r = r + assert(s12, 19960543.857179d0, 0.5d-6)
898 integer function tstg78()
900 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
901 double precision a, f
902 integer r, assert, omask
903 include
'geodesic.inc'
907 f = 1/298.257223563d0
910 call invers(a, f, 27.2d0, 0d0, -27.1d0, 179.5d0,
911 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
912 r = r + assert(azi1, 45.82468716758d0, 0.5d-11)
913 r = r + assert(azi2, 134.22776532670d0, 0.5d-11)
914 r = r + assert(s12, 19974354.765767d0, 0.5d-6)
920 integer function tstg80()
923 double precision azi1, azi2, s12, a12, m12, MM12, MM21, SS12
924 double precision a, f
925 integer r, assert, omask
926 include
'geodesic.inc'
930 f = 1/298.257223563d0
934 call invers(a, f, 0d0, 0d0, 0d0, 90d0,
935 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
936 r = r + assert(mm12, -0.00528427534d0, 0.5d-10)
937 r = r + assert(mm21, -0.00528427534d0, 0.5d-10)
939 call invers(a, f, 0d0, 0d0, 1d-6, 1d-6,
940 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
941 r = r + assert(mm12, 1d0, 0.5d-10)
942 r = r + assert(mm21, 1d0, 0.5d-10)
945 call invers(a, f, 20.001d0, 0d0, 20.001d0, 0d0,
946 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
947 r = r + assert(a12, 0d0, 1d-13)
948 r = r + assert(s12, 0d0, 1d-8)
949 r = r + assert(azi1, 180d0, 1d-13)
950 r = r + assert(azi2, 180d0, 1d-13)
951 r = r + assert(m12, 0d0, 1d-8)
952 r = r + assert(mm12, 1d0, 1d-15)
953 r = r + assert(mm21, 1d0, 1d-15)
954 r = r + assert(ss12, 0d0, 1d-10)
956 call invers(a, f, 90d0, 0d0, 90d0, 180d0,
957 + s12, azi1, azi2, omask, a12, m12, mm12, mm21, ss12)
958 r = r + assert(a12, 0d0, 1d-13)
959 r = r + assert(s12, 0d0, 1d-8)
960 r = r + assert(azi1, 0d0, 1d-13)
961 r = r + assert(azi2, 180d0, 1d-13)
962 r = r + assert(m12, 0d0, 1d-8)
963 r = r + assert(mm12, 1d0, 1d-15)
964 r = r + assert(mm21, 1d0, 1d-15)
965 r = r + assert(ss12, 127516405431022d0, 0.5d0)
971 integer function tstg84()
974 double precision lat2, lon2, azi2, a12, m12, MM12, MM21, SS12
975 double precision a, f, nan, inf, LatFix
976 integer r, assert, chknan, omask, flags
977 include
'geodesic.inc'
981 f = 1/298.257223563d0
984 inf = 1d0/latfix(0d0)
987 call direct(a, f, 0d0, 0d0, 90d0, inf,
988 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
992 call direct(a, f, 0d0, 0d0, 90d0, nan,
993 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
997 call direct(a, f, 0d0, 0d0, inf, 1000d0,
998 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
1000 r = r + chknan(lon2)
1001 r = r + chknan(azi2)
1002 call direct(a, f, 0d0, 0d0, nan, 1000d0,
1003 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
1004 r = r + chknan(lat2)
1005 r = r + chknan(lon2)
1006 r = r + chknan(azi2)
1007 call direct(a, f, 0d0, inf, 90d0, 1000d0,
1008 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
1009 r = r + assert(lat2, 0d0, 0d0)
1010 r = r + chknan(lon2)
1011 r = r + assert(azi2, 90d0, 0d0)
1012 call direct(a, f, 0d0, nan, 90d0, 1000d0,
1013 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
1014 r = r + assert(lat2, 0d0, 0d0)
1015 r = r + chknan(lon2)
1016 r = r + assert(azi2, 90d0, 0d0)
1017 call direct(a, f, inf, 0d0, 90d0, 1000d0,
1018 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
1019 r = r + chknan(lat2)
1020 r = r + chknan(lon2)
1021 r = r + chknan(azi2)
1022 call direct(a, f, nan, 0d0, 90d0, 1000d0,
1023 + flags, lat2, lon2, azi2, omask, a12, m12, mm12, mm21, ss12)
1024 r = r + chknan(lat2)
1025 r = r + chknan(lon2)
1026 r = r + chknan(azi2)
1032 integer function tstp0()
1034 double precision lata(4), lona(4)
1035 data lata / 89d0, 89d0, 89d0, 89d0 /
1036 data lona / 0d0, 90d0, 180d0, 270d0 /
1037 double precision latb(4), lonb(4)
1038 data latb / -89d0, -89d0, -89d0, -89d0 /
1039 data lonb / 0d0, 90d0, 180d0, 270d0 /
1040 double precision latc(4), lonc(4)
1041 data latc / 0d0, -1d0, 0d0, 1d0 /
1042 data lonc / -1d0, 0d0, 1d0, 0d0 /
1043 double precision latd(3), lond(3)
1044 data latd / 90d0, 0d0, 0d0 /
1045 data lond / 0d0, 0d0, 90d0 /
1046 double precision a, f, AA, PP
1048 include
'geodesic.inc'
1052 f = 1/298.257223563d0
1055 call area(a, f, lata, lona, 4, aa, pp)
1056 r = r + assert(pp, 631819.8745d0, 1d-4)
1057 r = r + assert(aa, 24952305678.0d0, 1d0)
1059 call area(a, f, latb, lonb, 4, aa, pp)
1060 r = r + assert(pp, 631819.8745d0, 1d-4)
1061 r = r + assert(aa, -24952305678.0d0, 1d0)
1063 call area(a, f, latc, lonc, 4, aa, pp)
1064 r = r + assert(pp, 627598.2731d0, 1d-4)
1065 r = r + assert(aa, 24619419146.0d0, 1d0)
1067 call area(a, f, latd, lond, 3, aa, pp)
1068 r = r + assert(pp, 30022685d0, 1d0)
1069 r = r + assert(aa, 63758202715511.0d0, 1d0)
1075 integer function tstp5()
1077 double precision lat(3), lon(3)
1078 data lat / 89d0, 89d0, 89d0 /
1079 data lon / 0.1d0, 90.1d0, -179.9d0 /
1080 double precision a, f, AA, PP
1082 include
'geodesic.inc'
1086 f = 1/298.257223563d0
1089 call area(a, f, lat, lon, 3, aa, pp)
1090 r = r + assert(pp, 539297d0, 1d0)
1091 r = r + assert(aa, 12476152838.5d0, 1d0)
1097 integer function tstp6()
1099 double precision lata(3), lona(3)
1100 data lata / 9d0, 9d0, 9d0 /
1101 data lona / -0.00000000000001d0, 180d0, 0d0 /
1102 double precision latb(3), lonb(3)
1103 data latb / 9d0, 9d0, 9d0 /
1104 data lonb / 0.00000000000001d0, 0d0, 180d0 /
1105 double precision latc(3), lonc(3)
1106 data latc / 9d0, 9d0, 9d0 /
1107 data lonc / 0.00000000000001d0, 180d0, 0d0 /
1108 double precision latd(3), lond(3)
1109 data latd / 9d0, 9d0, 9d0 /
1110 data lond / -0.00000000000001d0, 0d0, 180d0 /
1111 double precision a, f, AA, PP
1113 include
'geodesic.inc'
1117 f = 1/298.257223563d0
1120 call area(a, f, lata, lona, 3, aa, pp)
1121 r = r + assert(pp, 36026861d0, 1d0)
1122 r = r + assert(aa, 0d0, 1d0)
1128 integer function tstp12()
1130 double precision lat(2), lon(2)
1131 data lat / 66.562222222d0, 66.562222222d0 /
1132 data lon / 0d0, 180d0 /
1133 double precision a, f, AA, PP
1135 include
'geodesic.inc'
1139 f = 1/298.257223563d0
1142 call area(a, f, lat, lon, 2, aa, pp)
1143 r = r + assert(pp, 10465729d0, 1d0)
1144 r = r + assert(aa, 0d0, 1d0)
1150 integer function tstp13()
1152 double precision lat(6), lon(6)
1153 data lat / 89d0, 89d0, 89d0, 89d0, 89d0, 89d0 /
1154 data lon / -360d0, -240d0, -120d0, 0d0, 120d0, 240d0 /
1155 double precision a, f, AA, PP
1157 include
'geodesic.inc'
1161 f = 1/298.257223563d0
1164 call area(a, f, lat, lon, 6, aa, pp)
1165 r = r + assert(pp, 1160741d0, 1d0)
1166 r = r + assert(aa, 32415230256.0d0, 1d0)
1172 integer function tstp15()
1176 double precision lat(3), lon(3)
1177 data lat / 2d0, 1d0, 3d0 /
1178 data lon / 1d0, 2d0, 3d0 /
1179 double precision a, f, AA, PP
1181 include
'geodesic.inc'
1185 f = 1/298.257223563d0
1188 call area(a, f, lat, lon, 3, aa, pp)
1189 r = r + assert(aa, 18454562325.45119d0, 1d0)
1192 call area(a, f, lon, lat, 3, aa, pp)
1193 r = r + assert(aa, -18454562325.45119d0, 1d0)
1199 integer function tstp19()
1202 double precision lat(1), lon(1)
1205 double precision a, f, AA, PP
1207 include
'geodesic.inc'
1211 f = 1/298.257223563d0
1214 call area(a, f, lat, lon, 1, aa, pp)
1215 r = r + assert(aa, 0d0, 0d0)
1216 r = r + assert(pp, 0d0, 0d0)
1222 integer function tstp21()
1225 double precision lat(12), lon(12), lonr(12)
1226 data lat / 12*45d0 /
1227 data lon / 60d0, 180d0, -60d0,
1228 + 60d0, 180d0, -60d0,
1229 + 60d0, 180d0, -60d0,
1230 + 60d0, 180d0, -60d0 /
1231 data lonr / -60d0, 180d0, 60d0,
1232 + -60d0, 180d0, 60d0,
1233 + -60d0, 180d0, 60d0,
1234 + -60d0, 180d0, 60d0 /
1235 double precision a, f, AA, PP, AA1
1237 include
'geodesic.inc'
1241 f = 1/298.257223563d0
1244 aa1 = 39433884866571.4277d0
1247 call area(a, f, lat, lon, 3*i, aa, pp)
1248 r = r + assert(aa, aa1*i, 0.5d0)
1249 call area(a, f, lat, lonr, 3*i, aa, pp)
1250 r = r + assert(aa, -aa1*i, 0.5d0)
1259 integer tstinv, tstdir, tstarc,
1260 + tstg0, tstg1, tstg2, tstg5, tstg6, tstg9, tstg10, tstg11,
1261 + tstg12, tstg14, tstg15, tstg17, tstg26, tstg28, tstg33,
1262 + tstg55, tstg59, tstg61, tstg73, tstg74, tstg76, tstg78,
1264 + tstp0, tstp5, tstp6, tstp12, tstp13, tstp15, tstp19, tstp21
1270 print *,
'tstinv fail:', i
1275 print *,
'tstdir fail:', i
1280 print *,
'tstarc fail:', i
1285 print *,
'tstg0 fail:', i
1290 print *,
'tstg1 fail:', i
1295 print *,
'tstg2 fail:', i
1300 print *,
'tstg5 fail:', i
1305 print *,
'tstg6 fail:', i
1310 print *,
'tstg9 fail:', i
1315 print *,
'tstg10 fail:', i
1320 print *,
'tstg11 fail:', i
1325 print *,
'tstg12 fail:', i
1330 print *,
'tstg14 fail:', i
1335 print *,
'tstg15 fail:', i
1340 print *,
'tstg17 fail:', i
1345 print *,
'tstg26 fail:', i
1350 print *,
'tstg28 fail:', i
1355 print *,
'tstg33 fail:', i
1360 print *,
'tstg55 fail:', i
1365 print *,
'tstg59 fail:', i
1370 print *,
'tstg61 fail:', i
1375 print *,
'tstg73 fail:', i
1380 print *,
'tstg74 fail:', i
1385 print *,
'tstg76 fail:', i
1390 print *,
'tstg78 fail:', i
1395 print *,
'tstg80 fail:', i
1400 print *,
'tstg84 fail:', i
1405 print *,
'tstp0 fail:', i
1410 print *,
'tstp5 fail:', i
1415 print *,
'tstp6 fail:', i
1420 print *,
'tstp12 fail:', i
1425 print *,
'tstp13 fail:', i
1430 print *,
'tstp15 fail:', i
1435 print *,
'tstp19 fail:', i
1440 print *,
'tstp21 fail:', i