1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
|
/* See LICENSE file for copyright and license details. */
#include "common.h"
#define R LIBTELLURIAN_AUTHALIC_RADIUS
#ifndef TEST
const double libtellurian_coarse_distance_radius = R;
double
libtellurian_coarse_distance_radians(double latitude1, double longitude1,
double latitude2, double longitude2)
{
double h = fma(haversin(longitude2 - longitude1),
cos(latitude1) * cos(latitude2),
haversin(latitude2 - latitude1));
return 2.0 * R * asin(sqrt(h));
}
#else
static int
check(double a, double b, double c, double d, double expect)
{
double got = libtellurian_coarse_distance_radians(a, b, c, d);
double rgot = libtellurian_coarse_distance_radians(c, d, a, b);
double err = expect ? fabs(got / (R * expect) - 1.0) : got;
return err <= 1.0e-8 && got == rgot;
}
int
main(void)
{
ASSERT(check(0, 0, 0, 0, 0));
ASSERT(check(0, 0, D90, 0, M_PI / 2.0));
ASSERT(check(0, 0, 0, D90, M_PI / 2.0));
ASSERT(check(0, 0, D180, 0, M_PI));
ASSERT(check(0, 0, 0, D180, M_PI));
ASSERT(check(0, 0, D45, 0, M_PI / 4.0));
ASSERT(check(0, 0, 0, D45, M_PI / 4.0));
ASSERT(check(-D90, 0, D90, 0, M_PI));
ASSERT(check(0, -D90, 0, D90, M_PI));
ASSERT(check(0, 0, D60, D90, M_PI / 2.0));
ASSERT(check(0, 0, D30, D90, M_PI / 2.0));
ASSERT(check(0, 0, D90, D90, M_PI / 2.0));
ASSERT(check(0, 0, D45, D45, M_PI / 3.0));
ASSERT(check(-D45, -D45, D45, D45, M_PI / 1.5));
ASSERT(check(D45, -D45, D45, D45, M_PI / 3.0));
ASSERT(check(-D45, D45, D45, D45, M_PI / 2.0));
ASSERT(libtellurian_coarse_distance_radius == R);
return 0;
}
#endif
|
