aboutsummaryrefslogtreecommitdiffstats
path: root/libtellurian_sea_level_radians.c
blob: ec03985d7c55eb05c5882005dcd50c75e575777c (plain) (blame)
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
/* See LICENSE file for copyright and license details. */
#include "common.h"
#ifndef TEST


double
libtellurian_sea_level_radians(double latitude)
{
	double x = cos(latitude) * LIBTELLURIAN_EQUATORIAL_RADIUS;
	double y = sin(latitude) * LIBTELLURIAN_POLAR_RADIUS;
	return sqrt(fma(x, x, y * y));
}


#else


int
main(void)
{
	double x = pow(LIBTELLURIAN_EQUATORIAL_RADIUS, 4.0) + pow(LIBTELLURIAN_POLAR_RADIUS, 4.0);
	x = sqrt(x / (pow(LIBTELLURIAN_EQUATORIAL_RADIUS, 2.0) + pow(LIBTELLURIAN_POLAR_RADIUS, 2.0)));
	ASSERT(libtellurian_sea_level_radians(0) == LIBTELLURIAN_EQUATORIAL_RADIUS);
	ASSERT(libtellurian_sea_level_radians(D90) == LIBTELLURIAN_POLAR_RADIUS);
	ASSERT(libtellurian_sea_level_radians(D180) == LIBTELLURIAN_EQUATORIAL_RADIUS);
	ASSERT(libtellurian_sea_level_radians(-D90) == LIBTELLURIAN_POLAR_RADIUS);
	ASSERT(libtellurian_sea_level_radians(D45) < LIBTELLURIAN_EQUATORIAL_RADIUS);
	ASSERT(libtellurian_sea_level_radians(D45) > LIBTELLURIAN_POLAR_RADIUS);
	ASSERT(libtellurian_sea_level_radians(-D45) == libtellurian_sea_level_radians(D45));
	ASSERT(fabs(libtellurian_sea_level_radians(D45) / x - 1.0) < 1.0e-12);
	return 0;
}


#endif