aboutsummaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
authorMattias Andrée <m@maandree.se>2024-10-20 17:28:46 +0200
committerMattias Andrée <m@maandree.se>2024-10-20 17:28:46 +0200
commit09e6857ae73249bc7433f2971dcf291c70e4c766 (patch)
treef07ef152c6372083ab87e2e4289c7ab8b1f2c1ad
parentThird commit (diff)
downloadlibtellurian-09e6857ae73249bc7433f2971dcf291c70e4c766.tar.gz
libtellurian-09e6857ae73249bc7433f2971dcf291c70e4c766.tar.bz2
libtellurian-09e6857ae73249bc7433f2971dcf291c70e4c766.tar.xz
Fourth commit
Signed-off-by: Mattias Andrée <m@maandree.se>
-rw-r--r--.gitignore2
-rw-r--r--Makefile77
-rw-r--r--README2
-rw-r--r--TODO5
-rw-r--r--common.h29
-rw-r--r--config.mk2
-rw-r--r--libtellurian.717
-rw-r--r--libtellurian.h92
-rw-r--r--libtellurian.h.0201
-rw-r--r--libtellurian_azimuth.c6
-rw-r--r--libtellurian_azimuth_radians.c42
-rw-r--r--libtellurian_azimuthal_radius.361
-rw-r--r--libtellurian_azimuthal_radius.c43
-rw-r--r--libtellurian_azimuthal_radius_radians.c61
-rw-r--r--libtellurian_coarse_distance.32
-rw-r--r--libtellurian_coarse_distance.c19
-rw-r--r--libtellurian_coarse_distance_radians.c44
-rw-r--r--libtellurian_distance.32
-rw-r--r--libtellurian_distance.c6
-rw-r--r--libtellurian_distance_radians.c6
-rw-r--r--libtellurian_effective_gravity.359
-rw-r--r--libtellurian_effective_gravity.c41
l---------libtellurian_effective_gravity_radians.31
-rw-r--r--libtellurian_effective_gravity_radians.c51
-rw-r--r--libtellurian_elevated_gravity.396
-rw-r--r--libtellurian_elevated_gravity.c31
-rw-r--r--libtellurian_elevated_gravity_radians.c36
-rw-r--r--libtellurian_end_point.32
-rw-r--r--libtellurian_end_point.c6
-rw-r--r--libtellurian_end_point_radians.c23
-rw-r--r--libtellurian_gaussian_radius.352
-rw-r--r--libtellurian_gaussian_radius.c31
-rw-r--r--libtellurian_gaussian_radius_radians.c45
-rw-r--r--libtellurian_meridan_radius.c10
l---------libtellurian_meridan_radius_radians.31
-rw-r--r--libtellurian_meridan_radius_radians.c14
-rw-r--r--libtellurian_meridian_radius.350
-rw-r--r--libtellurian_meridian_radius.c41
l---------libtellurian_meridian_radius_radians.31
-rw-r--r--libtellurian_meridian_radius_radians.c45
-rw-r--r--libtellurian_normal_gravity.353
-rw-r--r--libtellurian_normal_gravity.c31
-rw-r--r--libtellurian_normal_gravity_radians.c37
-rw-r--r--libtellurian_sea_level.34
-rw-r--r--libtellurian_sea_level.c31
-rw-r--r--libtellurian_sea_level_radians.c24
-rw-r--r--libtellurian_transverse_radius.350
-rw-r--r--libtellurian_transverse_radius.c31
-rw-r--r--libtellurian_transverse_radius_radians.c34
-rw-r--r--libtellurian_vincenty_inverse__.c6
50 files changed, 1564 insertions, 92 deletions
diff --git a/.gitignore b/.gitignore
index a071ed4..42c07ce 100644
--- a/.gitignore
+++ b/.gitignore
@@ -12,3 +12,5 @@
*.gcov
*.gcno
*.gcda
+*.to
+*.test
diff --git a/Makefile b/Makefile
index d79bd06..99e7f4d 100644
--- a/Makefile
+++ b/Makefile
@@ -29,10 +29,12 @@ OBJ_PUBLIC =\
libtellurian_end_point_radians.o\
libtellurian_normal_gravity.o\
libtellurian_normal_gravity_radians.o\
+ libtellurian_effective_gravity.o\
+ libtellurian_effective_gravity_radians.o\
libtellurian_elevated_gravity.o\
libtellurian_elevated_gravity_radians.o\
- libtellurian_meridan_radius.o\
- libtellurian_meridan_radius_radians.o\
+ libtellurian_meridian_radius.o\
+ libtellurian_meridian_radius_radians.o\
libtellurian_transverse_radius.o\
libtellurian_transverse_radius_radians.o\
libtellurian_azimuthal_radius.o\
@@ -48,12 +50,46 @@ HDR =\
libtellurian.h\
common.h
+MAN0 = libtellurian.h.0
+MAN3 = $(MAN3_FUNC) $(MAN3_CONST)
+MAN7 = libtellurian.7
+
+MAN3_FUNC = $(OBJ_PUBLIC:.o=.3)
+
+MAN3_CONST =\
+ LIBTELLURIAN_EQUATORIAL_RADIUS.3\
+ LIBTELLURIAN_POLAR_RADIUS.3\
+ LIBTELLURIAN_MEAN_RADIUS.3\
+ LIBTELLURIAN_VOLUMETRIC_RADIUS.3\
+ LIBTELLURIAN_AUTHALIC_RADIUS.3\
+ LIBTELLURIAN_RECTIFYING_RADIUS.3\
+ LIBTELLURIAN_NOMINAL_EQUATORIAL_RADIUS.3\
+ LIBTELLURIAN_NOMINAL_POLAR_RADIUS.3\
+ LIBTELLURIAN_NOMINAL_RADIUS.3\
+ LIBTELLURIAN_EQUATORIAL_CIRCUMFERENCE.3\
+ LIBTELLURIAN_POLAR_CIRCUMFERENCE.3\
+ LIBTELLURIAN_MEAN_CIRCUMFERENCE.3\
+ LIBTELLURIAN_VOLUMETRIC_CIRCUMFERENCE.3\
+ LIBTELLURIAN_AUTHALIC_CIRCUMFERENCE.3\
+ LIBTELLURIAN_MERIDIONAL_CIRCUMFERENCE.3\
+ LIBTELLURIAN_EQUATORIAL_GRAVITY.3\
+ LIBTELLURIAN_POLAR_GRAVITY.3\
+ LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY.3\
+ LIBTELLURIAN_NORMAL_POLAR_GRAVITY.3\
+ LIBTELLURIAN_MASS_OF_EARTH.3\
+ LIBTELLURIAN_ANGULAR_VELOCITY.3\
+ LIBTELLURIAN_GEOCENTRIC_GRAVITATIONAL_CONSTANT.3
+
LOBJ = $(OBJ:.o=.lo)
+TOBJ = $(OBJ:.o=.to)
+TEST = $(TOBJ:.to=.test)
-all: libtellurian.a libtellurian.$(LIBEXT)
+all: libtellurian.a libtellurian.$(LIBEXT) $(TEST)
$(OBJ): $(HDR)
$(LOBJ): $(HDR)
+$(TOBJ): $(HDR)
+$(TEST): libtellurian.a
.c.o:
$(CC) -c -o $@ $< $(CFLAGS) $(CPPFLAGS)
@@ -61,6 +97,12 @@ $(LOBJ): $(HDR)
.c.lo:
$(CC) -fPIC -c -o $@ $< $(CFLAGS) $(CPPFLAGS)
+.c.to:
+ $(CC) -DTEST -c -o $@ $< $(CFLAGS) $(CPPFLAGS)
+
+.to.test:
+ $(CC) -o $@ $< libtellurian.a $(LDFLAGS)
+
libtellurian.a: $(OBJ)
@rm -f -- $@
$(AR) rc $@ $(OBJ)
@@ -69,15 +111,36 @@ libtellurian.a: $(OBJ)
libtellurian.$(LIBEXT): $(LOBJ)
$(CC) $(LIBFLAGS) -o $@ $(LOBJ) $(LDFLAGS)
+check: $(TEST)
+ @for t in $(TEST); do\
+ printf '%s ' $(CHECK_PREFIX) ./$$t;\
+ printf '\n\033[31m';\
+ if ! $(CHECK_PREFIX) "./$$t"; then\
+ printf '\033[1m%s\n' 'TEST FAILED';\
+ fail=y;\
+ fi;\
+ printf '\033[m';\
+ done ; test -z "$$fail"
+
install: libtellurian.a libtellurian.$(LIBEXT)
mkdir -p -- "$(DESTDIR)$(PREFIX)/lib"
mkdir -p -- "$(DESTDIR)$(PREFIX)/include"
+ mkdir -p -- "$(DESTDIR)$(MANPREFIX)/man0"
+ mkdir -p -- "$(DESTDIR)$(MANPREFIX)/man3"
+ mkdir -p -- "$(DESTDIR)$(MANPREFIX)/man7"
cp -- libtellurian.a "$(DESTDIR)$(PREFIX)/lib/"
cp -- libtellurian.$(LIBEXT) "$(DESTDIR)$(PREFIX)/lib/libtellurian.$(LIBMINOREXT)"
$(FIX_INSTALL_NAME) "$(DESTDIR)$(PREFIX)/lib/libtellurian.$(LIBMINOREXT)"
ln -sf -- libtellurian.$(LIBMINOREXT) "$(DESTDIR)$(PREFIX)/lib/libtellurian.$(LIBMAJOREXT)"
ln -sf -- libtellurian.$(LIBMAJOREXT) "$(DESTDIR)$(PREFIX)/lib/libtellurian.$(LIBEXT)"
cp -- libtellurian.h "$(DESTDIR)$(PREFIX)/include/"
+ cp -P -- $(MAN0) "$(DESTDIR)$(MANPREFIX)/man0/"
+ cp -P -- $(MAN3_FUNC) "$(DESTDIR)$(MANPREFIX)/man3/"
+ set -e && for f in $(MAN3_CONST); do \
+ test ! -d "$(DESTDIR)$(MANPREFIX)/man3/$$f";\
+ ln -sf -- ../man0/libtellurian.h.0 "$(DESTDIR)$(MANPREFIX)/man3/$$f";\
+ done
+ cp -P -- $(MAN7) "$(DESTDIR)$(MANPREFIX)/man7/"
uninstall:
-rm -f -- "$(DESTDIR)$(PREFIX)/lib/libtellurian.a"
@@ -85,12 +148,16 @@ uninstall:
-rm -f -- "$(DESTDIR)$(PREFIX)/lib/libtellurian.$(LIBMINOREXT)"
-rm -f -- "$(DESTDIR)$(PREFIX)/lib/libtellurian.$(LIBEXT)"
-rm -f -- "$(DESTDIR)$(PREFIX)/include/libtellurian.h"
+ -cd -- "$(DESTDIR)$(MANPREFIX)/man0/" && rm -f -- $(MAN0)
+ -cd -- "$(DESTDIR)$(MANPREFIX)/man3/" && rm -f -- $(MAN3)
+ -cd -- "$(DESTDIR)$(MANPREFIX)/man7/" && rm -f -- $(MAN7)
clean:
-rm -f -- *.o *.a *.lo *.su *.so *.so.* *.dll *.dylib
-rm -f -- *.gch *.gcov *.gcno *.gcda *.$(LIBEXT)
+ -rm -f -- *.to *.test
.SUFFIXES:
-.SUFFIXES: .lo .o .c
+.SUFFIXES: .lo .o .c .to .test
-.PHONY: all install uninstall clean
+.PHONY: all check install uninstall clean
diff --git a/README b/README
index 7dc7dae..6916fa9 100644
--- a/README
+++ b/README
@@ -11,7 +11,7 @@ DESCRIPTION
and constants.
libtellurian always uses meters as the length unit,
- meters per square-second as the acceleration unit,
+ meters per square second as the acceleration unit,
and degrees for the angle unit, except when the
function name uses the suffix "_radians", in which
case the function uses radians instead of degrees.
diff --git a/TODO b/TODO
index 7974f98..0b30e62 100644
--- a/TODO
+++ b/TODO
@@ -1,4 +1,5 @@
-Add header man page
-Add const man pages
Finish libtellurian_end_point_radians
Finish libtellurian_vincenty_inverse__
+
+Add a method for creating a locally optimised spherical model of the earth,
+and add method for calculating distances on this model
diff --git a/common.h b/common.h
index 36ca06f..332e19a 100644
--- a/common.h
+++ b/common.h
@@ -5,7 +5,7 @@
#include <math.h>
-#if __GNUC__
+#if defined(__GNUC__) && !defined(__clang__)
# pragma GCC diagnostic ignored "-Wunsuffixed-float-constants"
#endif
@@ -17,3 +17,30 @@
void libtellurian_vincenty_inverse__(double latitude1, double longitude1, double latitude2, double longitude2,
double *distance_out, double *azimuth1_out, double *azimuth2_out);
+
+#ifdef TEST
+# include <stdio.h>
+# include <stdlib.h>
+# if __GNUC__
+# pragma GCC diagnostic ignored "-Wfloat-equal"
+# endif
+# define D180 M_PI
+# define D90 (M_PI / 2.0)
+# define D60 (M_PI / 3.0)
+# define D45 (M_PI / 4.0)
+# define D30 (M_PI / 6.0)
+# define ASSERT(ASSERTION)\
+ do {\
+ if (!(ASSERTION)) {\
+ fprintf(stderr, "Failed assertion at line %i: %s\n", __LINE__, #ASSERTION);\
+ exit(1);\
+ }\
+ } while (0)
+# if defined(__GNUC__)
+# define TO\
+DO_TEST __attribute__((__const__)) int main(void) { return 0; }
+# else
+# define TO\
+DO_TEST int main(void) { return 0; }
+# endif
+#endif
diff --git a/config.mk b/config.mk
index 037fe35..adb22b3 100644
--- a/config.mk
+++ b/config.mk
@@ -4,5 +4,5 @@ MANPREFIX = $(PREFIX)/share/man
CC = c99
CPPFLAGS = -D_DEFAULT_SOURCE -D_BSD_SOURCE -D_XOPEN_SOURCE=700 -D_GNU_SOURCE
-CFLAGS =
+CFLAGS = -O3
LDFLAGS = -lm
diff --git a/libtellurian.7 b/libtellurian.7
index 812fca8..0933aa9 100644
--- a/libtellurian.7
+++ b/libtellurian.7
@@ -3,7 +3,7 @@
libtellurian \- Geodesy library
.SH SYNPOSIS
-.ni
+.nf
#include <libtellurian.h>
.fi
.PP
@@ -18,7 +18,7 @@ and constants.
.PP
.B libtellurian
always uses meters as the length unit,
-meters per square-second as the acceleration unit,
+meters per square second as the acceleration unit,
and degrees for the angle unit, except when the
function name uses the suffix \(dq_radians\(dq, in which
case the function uses radians instead of degrees.
@@ -66,13 +66,18 @@ travelling for some distance in some direction.
Calculate the normal gravity for some point on
the Earth's ellipsoid.
.TP
+.BR libtellurian_effective_gravity (3)
+Calculates the effective gravity, at some point
+on Earth's ellipsoid, from the normal gravity
+at the same point.
+.TP
.BR libtellurian_elevated_gravity (3)
Adjust a gravity for some ellipsoidal height.
.PP
.I Radius of curvature calculation
.TP
-.BR libtellurian_meridan_radius (3)
-Calculate the Earth's meridan radius of
+.BR libtellurian_meridian_radius (3)
+Calculate the Earth's meridian radius of
curvature for some latitude.
.TP
.BR libtellurian_transverse_radius (3)
@@ -86,3 +91,7 @@ curvature for some latitude and azimuth.
.BR libtellurian_gaussian_radius (3)
Calculate the Earth's Gaussian (non-directional)
radius of curvature for some latitude.
+.PP
+See
+.BR libtellurian.h (0)
+for a listing of available constants.
diff --git a/libtellurian.h b/libtellurian.h
index 2113f29..609f61d 100644
--- a/libtellurian.h
+++ b/libtellurian.h
@@ -40,7 +40,7 @@
/**
* The Earth's rectifying radius, in meters
*/
-#define LIBTELLURIAN_RECTIFYING_RADIUS 6367449.14582 /* 2π⁻¹∫{0 → ½π} √(a² cos² φ + b² sin² φ) dφ ≅ ∛(½(√a³+√b³))² */
+#define LIBTELLURIAN_RECTIFYING_RADIUS 6367449.145823328 /* 2π⁻¹∫{0→½π} √(a²cos²φ + b²sin²φ) dφ ≅ 1.5(a+b) - ½√[(3a+b)(a+3b)] */
/**
* The Earth's nominal equatorial radius, in meters
@@ -63,9 +63,11 @@
#define LIBTELLURIAN_EQUATORIAL_CIRCUMFERENCE 40075016.68557849 /* 2bπ */
/**
- * The Earth's meridional circumference, in meters
+ * The circumference, in meters, of a sphere inscribed in the Earth's
+ * spheroid and intersecting with it's pole (the circumference of
+ * a circle with Earth's polar radius)
*/
-#define LIBTELLURIAN_MERIDIONAL_CIRCUMFERENCE 39940652.74224401 /* 2bπ */
+#define LIBTELLURIAN_POLAR_CIRCUMFERENCE 39940652.74224401 /* 2bπ */
/**
* The Earth's mean circumference (arithmetic mean), in meters
@@ -85,33 +87,49 @@
/**
* The Earth's rectifying circumference, in meters
*/
-#define LIBTELLURIAN_RECTIFYING_CIRCUMFERENCE 40007862.91722943 /* 4 ∫{0 → ½π} √(a² cos² φ + b² sin² φ) dφ */
+#define LIBTELLURIAN_MERIDIONAL_CIRCUMFERENCE 40007862.91722943 /* 4 ∫{0 → ½π} √(a² cos² φ + b² sin² φ) dφ */
/**
- * The Earth's gravity at the equator, in meters per square-second
+ * The Earth's gravity at the equator, in meters per square second
*/
-#define LIBTELLURIAN_EQUATORIAL_GRAVITY 9.7803253359
+#define LIBTELLURIAN_EQUATORIAL_GRAVITY 9.7803253359 /* 𝔾ₑ */
/**
- * The Earth's gravity at the poles, in meters per square-second
+ * The Earth's gravity at the poles, in meters per square second
*/
-#define LIBTELLURIAN_POLAR_GRAVITY 9.8321849378
+#define LIBTELLURIAN_POLAR_GRAVITY 9.8321849378 /* 𝔾ₚ */
/**
- * The Earth's normal gravity at the equator, in meters per square-second
+ * The Earth's normal gravity at the equator, in meters per square second
*/
-#define LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY 9.7803267715
+#define LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY 9.78379185547421 /* ½𝔾ₑ + √((½𝔾ₑ)² + aω²) */
/**
- * The Earth's normal gravity at the poles, in meters per square-second
+ * The Earth's normal gravity at the poles, in meters per square second
*/
-#define LIBTELLURIAN_NORMAL_POLAR_GRAVITY 9.8321863685
+#define LIBTELLURIAN_NORMAL_POLAR_GRAVITY LIBTELLURIAN_POLAR_GRAVITY
/**
* The Earth's mass, in kilograms
*/
#define LIBTELLURIAN_MASS_OF_EARTH 5.972168e24
+/**
+ * The Earth's nominal mean angular velocity, in radians per second
+ */
+#define LIBTELLURIAN_ANGULAR_VELOCITY 7.292115e-5 /* ω */
+
+/**
+ * The geocentric gravitational constant, in cubic meters per square second
+ *
+ * This is the (universal) gravitational constant (6.67430e-11) multiplied
+ * by the mass of the Earth, however this value more reliable that the
+ * gravitation constant and the mass of th Earth, and should thus be used
+ * instead of multiplying the universial gravitational constant with the
+ * the msas of the Earth
+ */
+#define LIBTELLURIAN_GEOCENTRIC_GRAVITATIONAL_CONSTANT 3.986004418e14
+
/**
* Calculate the distance of the nominal sea level (geocentric
@@ -311,7 +329,7 @@ void libtellurian_end_point_radians(double latitude1, double longitude1, double
* as an oblate spheroid
*
* @param latitude GPS latitude coordinate, in degrees
- * @return The normal gravity, in meters per square-second
+ * @return The normal gravity, in meters per square second
*/
LIBTELLURIAN_CONST__
double libtellurian_normal_gravity(double latitude);
@@ -322,22 +340,48 @@ double libtellurian_normal_gravity(double latitude);
* as an oblate spheroid
*
* @param latitude GPS latitude coordinate, in radians
- * @return The normal gravity, in meters per square-second
+ * @return The normal gravity, in meters per square second
*/
LIBTELLURIAN_CONST__
double libtellurian_normal_gravity_radians(double latitude);
/**
+ * Calculate the effective gravity based on an already
+ * calculated normal gravity, that is, the normal gravity
+ * subtracted by a centrifugal correction based on the
+ * normal gravity and the latitude it is calculated from
+ *
+ * @param gravity The normal gravity, in meters per square second
+ * @param latitude GPS latitude coordinate, in degrees
+ * @return The effective gravity, in meters per square second
+ */
+LIBTELLURIAN_CONST__
+double libtellurian_effective_gravity(double gravity, double latitude);
+
+/**
+ * Calculate the effective gravity based on an already
+ * calculated normal gravity, that is, the normal gravity
+ * subtracted by a centrifugal correction based on the
+ * normal gravity and the latitude it is calculated from
+ *
+ * @param gravity The normal gravity, in meters per square second
+ * @param latitude GPS latitude coordinate, in radians
+ * @return The effective gravity, in meters per square second
+ */
+LIBTELLURIAN_CONST__
+double libtellurian_effective_gravity_radians(double gravity, double latitude);
+
+/**
* Calculate the gravity adjusted for the elevation
* above the altitude where the gravity is measure
*
* Altitudes above circa 100000 meters is out of range
* for this function (that would be in outer space)
*
- * @param gravity The gravity at sea level, in meters per square-second
+ * @param gravity The gravity at sea level, in meters per square second
* @param latitude GPS latitude coordinate, in degrees
* @param altitude Elevation above the gravity's measurement point, in meters
- * @return The height-adjusted gravity, in meters per square-second
+ * @return The height-adjusted gravity, in meters per square second
*/
LIBTELLURIAN_CONST__
double libtellurian_elevated_gravity(double gravity, double latitude, double altitude);
@@ -349,31 +393,31 @@ double libtellurian_elevated_gravity(double gravity, double latitude, double alt
* Altitudes above circa 100000 meters is out of range
* for this function (that would be in outer space)
*
- * @param gravity The gravity at sea level, in meters per square-second
+ * @param gravity The gravity at sea level, in meters per square second
* @param latitude GPS latitude coordinate, in radians
* @param altitude Elevation above the gravity's measurement point, in meters
- * @return The height-adjusted gravity, in meters per square-second
+ * @return The height-adjusted gravity, in meters per square second
*/
LIBTELLURIAN_CONST__
double libtellurian_elevated_gravity_radians(double gravity, double latitude, double altitude);
/**
- * Calculate the Earth's meridan radius of curvature at some latitude
+ * Calculate the Earth's meridian radius of curvature at some latitude
*
* @param latitude GPS latitude coordinate, in degrees
- * @return The meridan radius of curvature, in meters
+ * @return The meridian radius of curvature, in meters
*/
LIBTELLURIAN_CONST__
-double libtellurian_meridan_radius(double latitude);
+double libtellurian_meridian_radius(double latitude);
/**
- * Calculate the Earth's meridan radius of curvature at some latitude
+ * Calculate the Earth's meridian radius of curvature at some latitude
*
* @param latitude GPS latitude coordinate, in radians
- * @return The meridan radius of curvature, in meters
+ * @return The meridian radius of curvature, in meters
*/
LIBTELLURIAN_CONST__
-double libtellurian_meridan_radius_radians(double latitude);
+double libtellurian_meridian_radius_radians(double latitude);
/**
* Calculate the Earth's transverse radius of curvature at some latitude
diff --git a/libtellurian.h.0 b/libtellurian.h.0
new file mode 100644
index 0000000..472d35c
--- /dev/null
+++ b/libtellurian.h.0
@@ -0,0 +1,201 @@
+.TH LIBTELLURIAN.H 0 libtellurian
+.SH NAME
+libtellurian.h \- Geodesy library header
+
+.SH SYNPOSIS
+.nf
+#include <libtellurian.h>
+
+#define LIBTELLURIAN_EQUATORIAL_RADIUS 6378137.0
+#define LIBTELLURIAN_POLAR_RADIUS 6356752.314245
+
+#define LIBTELLURIAN_MEAN_RADIUS /* derived value omitted */
+#define LIBTELLURIAN_VOLUMETRIC_RADIUS /* derived value omitted */
+#define LIBTELLURIAN_AUTHALIC_RADIUS /* derived value omitted */
+#define LIBTELLURIAN_RECTIFYING_RADIUS /* derived value omitted */
+
+#define LIBTELLURIAN_NOMINAL_EQUATORIAL_RADIUS 6378100.
+#define LIBTELLURIAN_NOMINAL_POLAR_RADIUS 6356800.
+#define LIBTELLURIAN_NOMINAL_RADIUS LIBTELLURIAN_NOMINAL_EQUATORIAL_RADIUS
+
+#define LIBTELLURIAN_EQUATORIAL_CIRCUMFERENCE /* derived value omitted */
+#define LIBTELLURIAN_POLAR_CIRCUMFERENCE /* derived value omitted */
+
+#define LIBTELLURIAN_MEAN_CIRCUMFERENCE /* derived value omitted */
+#define LIBTELLURIAN_VOLUMETRIC_CIRCUMFERENCE /* derived value omitted */
+#define LIBTELLURIAN_AUTHALIC_CIRCUMFERENCE /* derived value omitted */
+#define LIBTELLURIAN_MERIDIONAL_CIRCUMFERENCE /* derived value omitted */
+
+#define LIBTELLURIAN_EQUATORIAL_GRAVITY 9.7803253359
+#define LIBTELLURIAN_POLAR_GRAVITY 9.8321849378
+#define LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY /* derived value omitted */
+#define LIBTELLURIAN_NORMAL_POLAR_GRAVITY LIBTELLURIAN_POLAR_GRAVITY
+
+#define LIBTELLURIAN_MASS_OF_EARTH 5.972168e24
+#define LIBTELLURIAN_ANGULAR_VELOCITY 7.292115e-5
+#define LIBTELLURIAN_GEOCENTRIC_GRAVITATIONAL_CONSTANT 3.986004418e14
+.fi
+.PP
+Link with
+.I -ltellurian
+.IR -lm .
+
+.SH DESCRIPTION
+.B libtellurian.h
+provides a collection of constants useful for geodesy:
+.TP
+.B LIBTELLURIAN_EQUATORIAL_RADIUS
+The Earth's semimajor axes, in meters, which is the radius
+of the of the circle formed by the Earth's equator, when
+modelled as a spheroid.
+
+This value is a measurement and is subject to refinement.
+.TP
+.B LIBTELLURIAN_POLAR_RADIUS
+The Earth's semiminor axis, in meters, which is the distance
+the center of the Earth, when modelled as a spheroid, and
+either pole.
+
+This value is a measurement and is subject to refinement.
+.TP
+.B LIBTELLURIAN_MEAN_RADIUS
+The arithmetic mean of
+.I LIBTELLURIAN_EQUATORIAL_RADIUS
+and
+.IR LIBTELLURIAN_POLAR_RADIUS .
+.TP
+.B LIBTELLURIAN_VOLUMETRIC_RADIUS
+The geometric mean of
+.I LIBTELLURIAN_EQUATORIAL_RADIUS
+and
+.IR LIBTELLURIAN_POLAR_RADIUS ,
+which is the radius of a sphere with the same volume as
+the Earth.
+.TP
+.B LIBTELLURIAN_AUTHALIC_RADIUS
+The radius, in meters, of a sphere with the same area
+as the Earth.
+
+This is a value derived from
+.I LIBTELLURIAN_EQUATORIAL_RADIUS
+and
+.IR LIBTELLURIAN_POLAR_RADIUS .
+.TP
+.B LIBTELLURIAN_RECTIFYING_RADIUS
+The radius, in meters, of a sphere with the same
+circumference as the a great ellipse, on the Earth's
+ellipsoid, intersecting with both poles (the
+meridional circumference).
+
+This is a value derived from
+.I LIBTELLURIAN_EQUATORIAL_RADIUS
+and
+.IR LIBTELLURIAN_POLAR_RADIUS .
+.TP
+.B LIBTELLURIAN_NOMINAL_EQUATORIAL_RADIUS
+The Earth's nominal equatorial radius, in meters,
+which is a standardised approximate value of
+.IR LIBTELLURIAN_EQUATORIAL_RADIUS .
+.TP
+.B LIBTELLURIAN_NOMINAL_POLAR_RADIUS
+The Earth's nominal polar radius, in meters,
+which is a standardised approximate value of
+.IR LIBTELLURIAN_POLAR_RADIUS .
+.TP
+.B LIBTELLURIAN_NOMINAL_RADIUS
+The Earth's nominal radius, in meters, which is
+a standardised approximate value of
+.IR LIBTELLURIAN_EQUATORIAL_RADIUS ,
+and is a synonym for
+.IR LIBTELLURIAN_NOMINAL_EQUATORIAL_RADIUS .
+.TP
+.B LIBTELLURIAN_EQUATORIAL_CIRCUMFERENCE
+The circumference, in meters, of a circle with the radius
+.IR LIBTELLURIAN_EQUATORIAL_RADIUS .
+.TP
+.B LIBTELLURIAN_POLAR_CIRCUMFERENCE
+The circumference, in meters, of a circle with the radius
+.IR LIBTELLURIAN_POLAR_RADIUS .
+.TP
+.B LIBTELLURIAN_MEAN_CIRCUMFERENCE
+The circumference, in meters, of a circle with the radius
+.IR LIBTELLURIAN_MEAN_RADIUS .
+.TP
+.B LIBTELLURIAN_VOLUMETRIC_CIRCUMFERENCE
+The circumference, in meters, of a circle with the radius
+.IR LIBTELLURIAN_VOLUMETRIC_RADIUS .
+.TP
+.B LIBTELLURIAN_AUTHALIC_CIRCUMFERENCE
+The circumference, in meters, of a circle with the radius
+.IR LIBTELLURIAN_AUTHALIC_RADIUS .
+.TP
+.B LIBTELLURIAN_MERIDIONAL_CIRCUMFERENCE
+The circumference, in meters, of a circle with the radius
+.IR LIBTELLURIAN_RECTIFYING_RADIUS .
+.TP
+.B LIBTELLURIAN_EQUATORIAL_GRAVITY
+The Earth's gravity, in meters per square second, at the equator.
+
+This value is a measurement and is subject to refinement.
+.TP
+.B LIBTELLURIAN_POLAR_GRAVITY
+The Earth's gravity, in meters per square second, at the poles.
+
+This value is a measurement and is subject to refinement.
+.TP
+.B LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY
+The Earth's normal gravity, in meters per square second, at
+the equator. This is a gravity theoretical gravitational
+acceleration at the Earth's ellipsoid unaffected by rotation
+and assuming the Earth's is spheroidal.
+
+This value is derived from
+.IR LIBTELLURIAN_EQUATORIAL_GRAVITY ,
+.IR LIBTELLURIAN_EQUATORIAL_RADIUS ,
+and
+.IR LIBTELLURIAN_ANGULAR_VELOCITY .
+.TP
+.B LIBTELLURIAN_NORMAL_POLAR_GRAVITY
+The Earth's normal gravity, in meters per square second, at
+the poles. This is a gravity theoretical gravitational
+acceleration at the Earth's ellipsoid unaffected by rotation
+and assuming the Earth's is spheroidal.
+
+This is an alias for a measurement and is subject to refinement.
+.TP
+.B LIBTELLURIAN_MASS_OF_EARTH
+The mass of the Earth, in kilograms.
+
+This value is a measurement and is subject to refinement.
+.TP
+.B LIBTELLURIAN_ANGULAR_VELOCITY
+The Earth's nominal mean angular velocity, in radians per second.
+
+This value is a measurement and is subject to refinement.
+.TP
+.B LIBTELLURIAN_GEOCENTRIC_GRAVITATIONAL_CONSTANT
+The geocentric gravitational constant, in cubic meters per square
+second. This is the product of the (universal) gravitational
+constant and the Earth's mass, however this value is more reliable
+than both the mass of the Earth and the (even less reliable)
+universal gravitational constant.
+
+This value is a measurement and is subject to refinement.
+.PP
+.B libtellurian.h
+also defines a number of geodesy functions. See
+.BR libtellurian (7)
+for more information.
+
+.SH SEE ALSO
+.BR libtellurian (7),
+.BR libtellurian_sea_level (3),
+.BR libtellurian_coarse_distance (3),
+.BR libtellurian_distance (3),
+.BR libtellurian_end_point (3),
+.BR libtellurian_normal_gravity (3),
+.BR libtellurian_elevated_gravity (3),
+.BR libtellurian_meridian_radius (3),
+.BR libtellurian_transverse_radius (3),
+.BR libtellurian_azimuthal_radius (3),
+.BR libtellurian_gaussian_radius (3)
diff --git a/libtellurian_azimuth.c b/libtellurian_azimuth.c
index cc08037..4bfc743 100644
--- a/libtellurian_azimuth.c
+++ b/libtellurian_azimuth.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
void
@@ -17,3 +18,8 @@ libtellurian_azimuth(double latitude1, double longitude1,
if (azimuth2_out)
*azimuth2_out = degrees(*azimuth2_out);
}
+
+
+#else
+TODO_TEST
+#endif
diff --git a/libtellurian_azimuth_radians.c b/libtellurian_azimuth_radians.c
index fa729d6..18c954e 100644
--- a/libtellurian_azimuth_radians.c
+++ b/libtellurian_azimuth_radians.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
void
@@ -11,3 +12,44 @@ libtellurian_azimuth_radians(double latitude1, double longitude1,
latitude2, longitude2,
NULL, azimuth1_out, azimuth2_out);
}
+
+
+#else
+
+
+#if 1
+TODO_TEST
+#else
+
+static int
+check(double a, double b, double c, double d)
+{
+ double az1_ref = 1, az2_ref = 2, az1 = 3, az2 = 4, az1b = 5, az2b = 6;
+ double s;
+
+ s = libtellurian_distance_radians(a, b, c, d, &az1_ref, &az2_ref);
+ (void) s;
+
+ libtellurian_azimuth_radians(a, b, c, d, &az1, &az2);
+ libtellurian_azimuth_radians(a, b, c, d, &az1b, NULL);
+ libtellurian_azimuth_radians(a, b, c, d, NULL, &az2b);
+ libtellurian_azimuth_radians(a, b, c, d, NULL, NULL);
+ return az1 == az1_ref && az2 == az2_ref && az1 == az1b && az2 == az2b;
+}
+
+int
+main(void)
+{
+ ASSERT(check(1, 2, 3, 4));
+ ASSERT(check(2, 3, 4, 5));
+ ASSERT(check(3, 4, 5, 6));
+ ASSERT(check(7, 8, 7, 8));
+ ASSERT(check(0, 0, 0, 0));
+ ASSERT(check(D90, D90, D90, D90));
+ return 0;
+}
+
+#endif
+
+
+#endif
diff --git a/libtellurian_azimuthal_radius.3 b/libtellurian_azimuthal_radius.3
new file mode 100644
index 0000000..a1be8aa
--- /dev/null
+++ b/libtellurian_azimuthal_radius.3
@@ -0,0 +1,61 @@
+.TH LIBTELLURIAN_AZIMUTHAL_RADIUS 3 libtellurian
+.SH NAME
+libtellurian_azimuthal_radius \- Calculate a azimuthal radius of curvature
+
+.SH SYNPOSIS
+.nf
+#include <libtellurian.h>
+
+double libtellurian_azimuthal_radius(double \fIlatitude\fP, double \fIazimuth\fP);
+
+double libtellurian_azimuthal_radius_radians(double \fIlatitude\fP, double \fIazimuth\fP);
+.fi
+.PP
+Link with
+.I -ltellurian
+.IR -lm .
+
+.SH DESCRIPTION
+The
+.BR libtellurian_azimuthal_radius ()
+function calculates the Earth's azimuthal radius of
+curvature at a given
+.I latitude
+for a given
+.IR azimuth ,
+that is, the radius of curvature along the normal
+section in the direction of the
+.I azimuth
+at the specified
+.I latitude.
+.PP
+The
+.I latitude
+shall be specified according to GPS and in degrees,
+likewise, the
+.I azimuth
+shall be specified in degrees.
+.PP
+The
+.BR libtellurian_azimuthal_radius_radians ()
+function is identical to the
+.BR libtellurian_azimuthal_radius ()
+function except that
+.I latitude
+and
+.I azimuth
+shall be specified in radians.
+
+.SH RETURN VALUE
+The
+.BR libtellurian_azimuthal_radius ()
+and
+.BR libtellurian_azimuthal_radius_radians ()
+functions return the radius of curvature
+measured in meters.
+
+.SH ERRORS
+None.
+
+.SH SEE ALSO
+.BR libtellurian (7)
diff --git a/libtellurian_azimuthal_radius.c b/libtellurian_azimuthal_radius.c
index 3ad5841..335a351 100644
--- a/libtellurian_azimuthal_radius.c
+++ b/libtellurian_azimuthal_radius.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -9,3 +10,45 @@ libtellurian_azimuthal_radius(double latitude, double azimuth)
azimuth = radians(azimuth);
return libtellurian_azimuthal_radius_radians(latitude, azimuth);
}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_azimuthal_radius(0, 30) == libtellurian_azimuthal_radius_radians(0, D30));
+ ASSERT(approx(libtellurian_azimuthal_radius(180, 30), libtellurian_azimuthal_radius_radians(0, D30)));
+ ASSERT(approx(libtellurian_azimuthal_radius(360, 30), libtellurian_azimuthal_radius_radians(0, D30)));
+ ASSERT(approx(libtellurian_azimuthal_radius(90, 30), libtellurian_azimuthal_radius(-90, 30)));
+ ASSERT(approx(libtellurian_azimuthal_radius(90, 30), libtellurian_azimuthal_radius(270, 30)));
+ ASSERT(approx(libtellurian_azimuthal_radius(45, 30), libtellurian_azimuthal_radius(-45, 30)));
+ ASSERT(approx(libtellurian_azimuthal_radius(135, 30), libtellurian_azimuthal_radius(45, 30)));
+ ASSERT(approx(libtellurian_azimuthal_radius(-135, 30), libtellurian_azimuthal_radius(135, 30)));
+ ASSERT(libtellurian_azimuthal_radius(90, 30) == libtellurian_azimuthal_radius_radians(D90, D30));
+ ASSERT(libtellurian_azimuthal_radius(45, 30) == libtellurian_azimuthal_radius_radians(D45, D30));
+ ASSERT(libtellurian_azimuthal_radius(30, 30) == libtellurian_azimuthal_radius_radians(D30, D30));
+
+ ASSERT(libtellurian_azimuthal_radius(0, 60) == libtellurian_azimuthal_radius_radians(0, D60));
+ ASSERT(approx(libtellurian_azimuthal_radius(180, 60), libtellurian_azimuthal_radius_radians(0, D60)));
+ ASSERT(approx(libtellurian_azimuthal_radius(360, 60), libtellurian_azimuthal_radius_radians(0, D60)));
+ ASSERT(approx(libtellurian_azimuthal_radius(90, 60), libtellurian_azimuthal_radius(-90, 60)));
+ ASSERT(approx(libtellurian_azimuthal_radius(90, 60), libtellurian_azimuthal_radius(270, 60)));
+ ASSERT(approx(libtellurian_azimuthal_radius(45, 60), libtellurian_azimuthal_radius(-45, 60)));
+ ASSERT(approx(libtellurian_azimuthal_radius(135, 60), libtellurian_azimuthal_radius(45, 60)));
+ ASSERT(approx(libtellurian_azimuthal_radius(-135, 60), libtellurian_azimuthal_radius(135, 60)));
+ ASSERT(libtellurian_azimuthal_radius(90, 60) == libtellurian_azimuthal_radius_radians(D90, D60));
+ ASSERT(libtellurian_azimuthal_radius(45, 60) == libtellurian_azimuthal_radius_radians(D45, D60));
+ ASSERT(libtellurian_azimuthal_radius(30, 60) == libtellurian_azimuthal_radius_radians(D30, D60));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_azimuthal_radius_radians.c b/libtellurian_azimuthal_radius_radians.c
index 00671af..d58a903 100644
--- a/libtellurian_azimuthal_radius_radians.c
+++ b/libtellurian_azimuthal_radius_radians.c
@@ -1,11 +1,12 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
libtellurian_azimuthal_radius_radians(double latitude, double azimuth)
{
- double m = libtellurian_meridan_radius_radians(latitude);
+ double m = libtellurian_meridian_radius_radians(latitude);
double n = libtellurian_transverse_radius_radians(latitude);
double c2 = cos(azimuth) * cos(azimuth);
double s2 = sin(azimuth) * sin(azimuth);
@@ -13,3 +14,61 @@ libtellurian_azimuthal_radius_radians(double latitude, double azimuth)
double y = s2 / n;
return 1.0 / (x + y);
}
+
+
+#else
+
+
+static int
+approx(double a, double e)
+{
+ double err = fabs((a / e) - 1.0);
+ return err <= 1.0e-8;
+}
+
+
+int
+main(void)
+{
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, 0), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, D30), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, D45), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, D60), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, D90), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, D180), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, -D30), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, -D45), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, -D60), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, -D90), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, -D180), libtellurian_meridian_radius_radians(D90)));
+
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(0, 0), libtellurian_meridian_radius_radians(0)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D30, 0), libtellurian_meridian_radius_radians(D30)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D45, 0), libtellurian_meridian_radius_radians(D45)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D60, 0), libtellurian_meridian_radius_radians(D60)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D180, 0), libtellurian_meridian_radius_radians(D180)));
+
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(0, D180), libtellurian_meridian_radius_radians(0)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D30, D180), libtellurian_meridian_radius_radians(D30)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D45, D180), libtellurian_meridian_radius_radians(D45)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D60, D180), libtellurian_meridian_radius_radians(D60)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D180, D180), libtellurian_meridian_radius_radians(D180)));
+
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(0, D90), libtellurian_transverse_radius_radians(0)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D30, D90), libtellurian_transverse_radius_radians(D30)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D45, D90), libtellurian_transverse_radius_radians(D45)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D60, D90), libtellurian_transverse_radius_radians(D60)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, D90), libtellurian_transverse_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D180, D90), libtellurian_transverse_radius_radians(D180)));
+
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(0, -D90), libtellurian_transverse_radius_radians(0)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D30, -D90), libtellurian_transverse_radius_radians(D30)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D45, -D90), libtellurian_transverse_radius_radians(D45)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D60, -D90), libtellurian_transverse_radius_radians(D60)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D90, -D90), libtellurian_transverse_radius_radians(D90)));
+ ASSERT(approx(libtellurian_azimuthal_radius_radians(D180, -D90), libtellurian_transverse_radius_radians(D180)));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_coarse_distance.3 b/libtellurian_coarse_distance.3
index f1c308c..9a3973b 100644
--- a/libtellurian_coarse_distance.3
+++ b/libtellurian_coarse_distance.3
@@ -3,7 +3,7 @@
libtellurian_coarse_distance \- Calculate distance between two locations
.SH SYNPOSIS
-.ni
+.nf
#include <libtellurian.h>
double libtellurian_coarse_distance(double \fIlatitude1\fP, double \fIlongitude1\fP,
diff --git a/libtellurian_coarse_distance.c b/libtellurian_coarse_distance.c
index 8f056d1..3776855 100644
--- a/libtellurian_coarse_distance.c
+++ b/libtellurian_coarse_distance.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -12,3 +13,21 @@ libtellurian_coarse_distance(double latitude1, double longitude1,
longitude2 = radians(longitude2);
return libtellurian_coarse_distance_radians(latitude1, longitude1, latitude2, longitude2);
}
+
+
+#else
+
+
+int
+main(void)
+{
+ ASSERT(libtellurian_coarse_distance(0, 0, 0, 0) == libtellurian_coarse_distance_radians(0, 0, 0, 0));
+ ASSERT(libtellurian_coarse_distance(90, 0, 0, 0) == libtellurian_coarse_distance_radians(D90, 0, 0, 0));
+ ASSERT(libtellurian_coarse_distance(90, 180, 0, 0) == libtellurian_coarse_distance_radians(D90, D180, 0, 0));
+ ASSERT(libtellurian_coarse_distance(90, 180, -45, 0) == libtellurian_coarse_distance_radians(D90, D180, -D45, 0));
+ ASSERT(libtellurian_coarse_distance(90, 180, -45, -30) == libtellurian_coarse_distance_radians(D90, D180, -D45, -D30));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_coarse_distance_radians.c b/libtellurian_coarse_distance_radians.c
index eb9cac0..f22ad28 100644
--- a/libtellurian_coarse_distance_radians.c
+++ b/libtellurian_coarse_distance_radians.c
@@ -1,6 +1,10 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#define R LIBTELLURIAN_AUTHALIC_RADIUS
+
+#ifndef TEST
+
double
libtellurian_coarse_distance_radians(double latitude1, double longitude1,
@@ -9,5 +13,43 @@ libtellurian_coarse_distance_radians(double latitude1, double longitude1,
double h = fma(haversin(longitude2 - longitude1),
cos(latitude1) * cos(latitude2),
haversin(latitude2 - latitude1));
- return 2.0 * LIBTELLURIAN_AUTHALIC_RADIUS * asin(sqrt(h));
+ 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));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_distance.3 b/libtellurian_distance.3
index 47e0ee1..49181b3 100644
--- a/libtellurian_distance.3
+++ b/libtellurian_distance.3
@@ -3,7 +3,7 @@
libtellurian_distance \- Calculate distance between two locations
.SH SYNPOSIS
-.ni
+.nf
#include <libtellurian.h>
double libtellurian_distance(double \fIlatitude1\fP, double \fIlongitude1\fP,
diff --git a/libtellurian_distance.c b/libtellurian_distance.c
index 0a290d1..250dddc 100644
--- a/libtellurian_distance.c
+++ b/libtellurian_distance.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -19,3 +20,8 @@ libtellurian_distance(double latitude1, double longitude1,
*azimuth2_out = degrees(*azimuth2_out);
return r;
}
+
+
+#else
+TODO_TEST
+#endif
diff --git a/libtellurian_distance_radians.c b/libtellurian_distance_radians.c
index b2b2134..d7c5547 100644
--- a/libtellurian_distance_radians.c
+++ b/libtellurian_distance_radians.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -13,3 +14,8 @@ libtellurian_distance_radians(double latitude1, double longitude1,
&s, azimuth1_out, azimuth2_out);
return s;
}
+
+
+#else
+TODO_TEST
+#endif
diff --git a/libtellurian_effective_gravity.3 b/libtellurian_effective_gravity.3
new file mode 100644
index 0000000..d462964
--- /dev/null
+++ b/libtellurian_effective_gravity.3
@@ -0,0 +1,59 @@
+.TH LIBTELLURIAN_EFFECTIVE_GRAVITY 3 libtellurian
+.SH NAME
+libtellurian_effective_gravity \- Calculate effective gravity
+
+.SH SYNPOSIS
+.nf
+#include <libtellurian.h>
+
+double libtellurian_effective_gravity(double \fIgravity\fP, double \fIlatitude\fP);
+
+double libtellurian_effective_gravity_radians(double \fIgravity\fP, double \fIlatitude\fP);
+.fi
+.PP
+Link with
+.I -ltellurian
+.IR -lm .
+
+.SH DESCRIPTION
+The
+.BR libtellurian_effective_gravity ()
+function models the Earth as a spheroid and calculates
+the effective gravity for any point on a given
+.I latitude
+at the surface and the normal
+.I gravity
+at this point. This is an adjustment based on the
+centrifugal force at the given latitude due to the
+rotation of the Earth.
+.PP
+The
+.I latitude
+shall be specified according to GPS and in degrees,
+and the normal
+.I gravity
+shall be specified in meters per square second.
+.PP
+The
+.BR libtellurian_effective_gravity_radians ()
+function is identical to the
+.BR libtellurian_effective_gravity ()
+function except that
+.I latitude
+shall be specified in radians.
+
+.SH RETURN VALUE
+The
+.BR libtellurian_effective_gravity ()
+and
+.BR libtellurian_effective_gravity_radians ()
+functions return the effective gravity, at the specified
+latitude, measured in meters per square second.
+
+.SH ERRORS
+None.
+
+.SH SEE ALSO
+.BR libtellurian (7),
+.BR libtellurian_normal_gravity (3),
+.BR libtellurian_elevated_gravity (3)
diff --git a/libtellurian_effective_gravity.c b/libtellurian_effective_gravity.c
new file mode 100644
index 0000000..fe8270e
--- /dev/null
+++ b/libtellurian_effective_gravity.c
@@ -0,0 +1,41 @@
+/* See LICENSE file for copyright and license details. */
+#include "common.h"
+#ifndef TEST
+
+
+double
+libtellurian_effective_gravity(double gravity, double latitude)
+{
+ latitude = radians(latitude);
+ return libtellurian_effective_gravity_radians(gravity, latitude);
+}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_effective_gravity(9.1234, 0) == libtellurian_effective_gravity_radians(9.1234, 0));
+ ASSERT(approx(libtellurian_effective_gravity(9.1234, 180), libtellurian_effective_gravity_radians(9.1234, 0)));
+ ASSERT(approx(libtellurian_effective_gravity(9.1234, 360), libtellurian_effective_gravity_radians(9.1234, 0)));
+ ASSERT(approx(libtellurian_effective_gravity(9.1234, 90), libtellurian_effective_gravity(9.1234, -90)));
+ ASSERT(approx(libtellurian_effective_gravity(9.1234, 90), libtellurian_effective_gravity(9.1234, 270)));
+ ASSERT(approx(libtellurian_effective_gravity(9.1234, 45), libtellurian_effective_gravity(9.1234, -45)));
+ ASSERT(approx(libtellurian_effective_gravity(9.1234, 135), libtellurian_effective_gravity(9.1234, 45)));
+ ASSERT(approx(libtellurian_effective_gravity(9.1234, -135), libtellurian_effective_gravity(9.1234, 135)));
+ ASSERT(libtellurian_effective_gravity(9.1234, 90) == libtellurian_effective_gravity_radians(9.1234, D90));
+ ASSERT(libtellurian_effective_gravity(9.1234, 45) == libtellurian_effective_gravity_radians(9.1234, D45));
+ ASSERT(libtellurian_effective_gravity(9.1234, 30) == libtellurian_effective_gravity_radians(9.1234, D30));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_effective_gravity_radians.3 b/libtellurian_effective_gravity_radians.3
new file mode 120000
index 0000000..4e962d7
--- /dev/null
+++ b/libtellurian_effective_gravity_radians.3
@@ -0,0 +1 @@
+libtellurian_effective_gravity.3 \ No newline at end of file
diff --git a/libtellurian_effective_gravity_radians.c b/libtellurian_effective_gravity_radians.c
new file mode 100644
index 0000000..54e2b8e
--- /dev/null
+++ b/libtellurian_effective_gravity_radians.c
@@ -0,0 +1,51 @@
+/* See LICENSE file for copyright and license details. */
+#include "common.h"
+
+#define K (LIBTELLURIAN_EQUATORIAL_RADIUS * LIBTELLURIAN_ANGULAR_VELOCITY * LIBTELLURIAN_ANGULAR_VELOCITY)
+
+#ifndef TEST
+
+
+double
+libtellurian_effective_gravity_radians(double gravity, double latitude)
+{
+ return fma(-K / gravity, cos(latitude) * cos(latitude), gravity);
+}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_effective_gravity_radians(1, D90) == 1);
+ ASSERT(libtellurian_effective_gravity_radians(2, D90) == 2);
+ ASSERT(libtellurian_effective_gravity_radians(4, D90) == 4);
+ ASSERT(libtellurian_effective_gravity_radians(1, -D90) == 1);
+ ASSERT(libtellurian_effective_gravity_radians(2, -D90) == 2);
+ ASSERT(libtellurian_effective_gravity_radians(4, -D90) == 4);
+ ASSERT(libtellurian_effective_gravity_radians(1, 0) == fma(1.00, -K, 1.0));
+ ASSERT(libtellurian_effective_gravity_radians(2, 0) == fma(0.50, -K, 2.0));
+ ASSERT(libtellurian_effective_gravity_radians(4, 0) == fma(0.25, -K, 4.0));
+ ASSERT(libtellurian_effective_gravity_radians(1, D45) == fma(0.50, -K, 1.0));
+ ASSERT(libtellurian_effective_gravity_radians(2, D45) == fma(0.25, -K, 2.0));
+ ASSERT(libtellurian_effective_gravity_radians(1, D60) == fma(0.250, -K, 1.0));
+ ASSERT(libtellurian_effective_gravity_radians(2, D60) == fma(0.125, -K, 2.0));
+ ASSERT(libtellurian_effective_gravity_radians(4, D30) < 4);
+ ASSERT(libtellurian_effective_gravity_radians(4, D30) > libtellurian_effective_gravity_radians(4, 0));
+ ASSERT(approx(libtellurian_effective_gravity_radians(LIBTELLURIAN_NORMAL_POLAR_GRAVITY, D90),
+ LIBTELLURIAN_POLAR_GRAVITY));
+ ASSERT(approx(libtellurian_effective_gravity_radians(LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY, 0),
+ LIBTELLURIAN_EQUATORIAL_GRAVITY));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_elevated_gravity.3 b/libtellurian_elevated_gravity.3
new file mode 100644
index 0000000..fb728b9
--- /dev/null
+++ b/libtellurian_elevated_gravity.3
@@ -0,0 +1,96 @@
+.TH LIBTELLURIAN_ELEVATED_GRAVITY 3 libtellurian
+.SH NAME
+libtellurian_elevated_gravity \- Calculate gravity at some altitude
+
+.SH SYNPOSIS
+.nf
+#include <libtellurian.h>
+
+double libtellurian_elevated_gravity(double \fIgravity\fP, double \fIlatitude\fP, double \fIaltitude\fP);
+
+double libtellurian_elevated_gravity_radians(double \fIgravity\fP, double \fIlatitude\fP, double \fIaltitude\fP);
+.fi
+.PP
+Link with
+.I -ltellurian
+.IR -lm .
+
+.SH DESCRIPTION
+The
+.BR libtellurian_elevated_gravity ()
+function models the Earth as a spheroid and calculates
+the gravity at a given
+.I altitude
+and
+.IR latitude ,
+given the
+.I gravity
+at zero-altitude on the same
+.IR latitude .
+.PP
+The
+.I latitude
+shall be specified according to GPS and in degrees
+and the
+.I gravity
+shall be specified in meters per square second.
+However the
+.I altitude
+shall be the number of metres above a reference datum
+for which the input
+.I gravity
+applies. For the normal gravity, this would be the
+geocentric radius as calculated by the
+.BR libtellurian_sea_level (3)
+function, and so the
+.I altitude
+would be the altitude reported in raw GPS measurements,
+howver, for real world gravity, this would typically
+be the orthometric altitude.
+.PP
+The
+.BR libtellurian_elevated_gravity_radians ()
+function is identical to the
+.BR libtellurian_elevated_gravity ()
+function except that
+.I latitude
+shall be specified in radians.
+
+.SH RETURN VALUE
+The
+.BR libtellurian_elevated_gravity ()
+and
+.BR libtellurian_elevated_gravity_radians ()
+functions return the normal gravity, at the specified
+latitude, measured in meters per square second.
+
+.SH ERRORS
+None.
+
+.SH NOTES
+Some GPS systems, espcially those for hiking, surveying,
+and aviation, report orthometric height (height above
+sea level), rather than ellipsoidal height to make it
+easier for the user to interpret.
+.PP
+The
+.BR libtellurian_normal_gravity (3)
+function does not account for the Earth's rotation,
+however, the
+.BR libtellurian_elevated_gravity ()
+function does account for the Earth's rotation.
+Given a normal gravity, it should for first be corrected
+for the rotation of the Earth, using the
+.BR libtellurian_effective_gravity (3)
+function, before it is input to the
+.BR libtellurian_elevated_gravity ()
+function.
+.PP
+The adjustment this function creates is useful for
+common heights in aviation, does not work for outer
+space.
+
+.SH SEE ALSO
+.BR libtellurian (7),
+.BR libtellurian_effective_gravity (3),
+.BR libtellurian_normal_gravity (3)
diff --git a/libtellurian_elevated_gravity.c b/libtellurian_elevated_gravity.c
index 575c0e8..fcc2b4d 100644
--- a/libtellurian_elevated_gravity.c
+++ b/libtellurian_elevated_gravity.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -8,3 +9,33 @@ libtellurian_elevated_gravity(double gravity, double latitude, double altitude)
latitude = radians(latitude);
return libtellurian_elevated_gravity_radians(gravity, latitude, altitude);
}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_elevated_gravity(9.1234, 0, 123) == libtellurian_elevated_gravity_radians(9.1234, 0, 123));
+ ASSERT(approx(libtellurian_elevated_gravity(9.1234, 180, 123), libtellurian_elevated_gravity_radians(9.1234, 0, 123)));
+ ASSERT(approx(libtellurian_elevated_gravity(9.1234, 360, 123), libtellurian_elevated_gravity_radians(9.1234, 0, 123)));
+ ASSERT(approx(libtellurian_elevated_gravity(9.1234, 90, 123), libtellurian_elevated_gravity(9.1234, -90, 123)));
+ ASSERT(approx(libtellurian_elevated_gravity(9.1234, 90, 123), libtellurian_elevated_gravity(9.1234, 270, 123)));
+ ASSERT(approx(libtellurian_elevated_gravity(9.1234, 45, 123), libtellurian_elevated_gravity(9.1234, -45, 123)));
+ ASSERT(approx(libtellurian_elevated_gravity(9.1234, 135, 123), libtellurian_elevated_gravity(9.1234, 45, 123)));
+ ASSERT(approx(libtellurian_elevated_gravity(9.1234, -135, 123), libtellurian_elevated_gravity(9.1234, 135, 123)));
+ ASSERT(libtellurian_elevated_gravity(9.1234, 90, 123) == libtellurian_elevated_gravity_radians(9.1234, D90, 123));
+ ASSERT(libtellurian_elevated_gravity(9.1234, 45, 123) == libtellurian_elevated_gravity_radians(9.1234, D45, 123));
+ ASSERT(libtellurian_elevated_gravity(9.1234, 30, 123) == libtellurian_elevated_gravity_radians(9.1234, D30, 123));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_elevated_gravity_radians.c b/libtellurian_elevated_gravity_radians.c
index ffd1a5b..1c97f1a 100644
--- a/libtellurian_elevated_gravity_radians.c
+++ b/libtellurian_elevated_gravity_radians.c
@@ -1,14 +1,40 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
libtellurian_elevated_gravity_radians(double gravity, double latitude, double altitude)
{
- double k1 = -3.1570429877205807e-07;
- double k2 = 2.1026896504579084e-09;
- double k3 = -7.374516772941995e-14;
+ double a = LIBTELLURIAN_EQUATORIAL_RADIUS;
+ double b = LIBTELLURIAN_POLAR_RADIUS;
+ double omega = LIBTELLURIAN_ANGULAR_VELOCITY;
+ double GM = LIBTELLURIAN_GEOCENTRIC_GRAVITATIONAL_CONSTANT;
+ double f = 1.0 - b / a;
+ double m = (omega * a) * (omega * a) * (b / GM);
+ double neg_k1 = fma(-2.0, f + m, -2.0) / a;
+ double k2 = 4.0 * f / a;
+ double neg_k3 = -3.0 / (a * a);
double sin2_phi = sin(latitude) * sin(latitude);
- double f = fma(k3, altitude, fma(k2, sin2_phi, k1));
- return fma(f * altitude, gravity, gravity);
+ double s = fma(neg_k3, altitude, fma(k2, sin2_phi, neg_k1));
+ return fma(s * altitude, gravity, gravity);
}
+
+
+#else
+
+
+int
+main(void)
+{
+ ASSERT(libtellurian_elevated_gravity_radians(9.2, 0, 0) == 9.2);
+ ASSERT(libtellurian_elevated_gravity_radians(4.2, 0, 0) == 4.2);
+ ASSERT(libtellurian_elevated_gravity_radians(8.9, 2.0, 0) == 8.9);
+ ASSERT(libtellurian_elevated_gravity_radians(8.9, 2.0, 1.0) < 8.9);
+ ASSERT(libtellurian_elevated_gravity_radians(8.9, 2.0, 2.0) < libtellurian_elevated_gravity_radians(8.9, 2.0, 1.0));
+ ASSERT(libtellurian_elevated_gravity_radians(8.9, 3.0, 1.0) < libtellurian_elevated_gravity_radians(8.9, 2.0, 1.0));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_end_point.3 b/libtellurian_end_point.3
index 3a93fe9..15fe232 100644
--- a/libtellurian_end_point.3
+++ b/libtellurian_end_point.3
@@ -3,7 +3,7 @@
libtellurian_end_point \- Calculate travel end-point
.SH SYNPOSIS
-.ni
+.nf
#include <libtellurian.h>
void libtellurian_end_point(double \fIlatitude1\fP, double \fIlongitude1\fP,
diff --git a/libtellurian_end_point.c b/libtellurian_end_point.c
index fb2c8be..57af198 100644
--- a/libtellurian_end_point.c
+++ b/libtellurian_end_point.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
void
@@ -18,3 +19,8 @@ libtellurian_end_point(double latitude1, double longitude1, double azimuth1, dou
if (azimuth2_out)
*azimuth2_out = degrees(*azimuth2_out);
}
+
+
+#else
+TODO_TEST
+#endif
diff --git a/libtellurian_end_point_radians.c b/libtellurian_end_point_radians.c
index 79e759f..c512b59 100644
--- a/libtellurian_end_point_radians.c
+++ b/libtellurian_end_point_radians.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
/* TODO complete the implementation */
@@ -55,17 +56,18 @@ libtellurian_end_point_radians(double latitude1, double longitude1, double azimu
/* } repeat until sigma converges */
- if (latitude2_out || azimuth2_out)
+ if (latitude2_out || azimuth2_out) {
x0 = fma(cos_sigma * cos_alpha1, cos_u1, -sin_u1 * sin_alpha);
- if (latitude2_out) {
- y = fma(sin_sigma * cos_alpha1, cos_u1, sin_u1 * cos_sigma);
- x = sqrt(fma(sin_alpha, sin_alpha, x0 * x0)) * c;
- *latitude2_out = atan2(y, x);
- }
+ if (latitude2_out) {
+ y = fma(sin_sigma * cos_alpha1, cos_u1, sin_u1 * cos_sigma);
+ x = sqrt(fma(sin_alpha, sin_alpha, x0 * x0)) * c;
+ *latitude2_out = atan2(y, x);
+ }
- if (azimuth2_out)
- *azimuth2_out = atan2(sin_alpha, x0);
+ if (azimuth2_out)
+ *azimuth2_out = atan2(sin_alpha, x0);
+ }
if (longitude2_out) {
C = fma(fma(-0.75, cos2_alpha, 1.0), f, 1.0) * cos2_alpha * f / 4.0;
@@ -77,3 +79,8 @@ libtellurian_end_point_radians(double latitude1, double longitude1, double azimu
*longitude2_out = longitude1 + L;
}
}
+
+
+#else
+TODO_TEST
+#endif
diff --git a/libtellurian_gaussian_radius.3 b/libtellurian_gaussian_radius.3
new file mode 100644
index 0000000..aaa84dc
--- /dev/null
+++ b/libtellurian_gaussian_radius.3
@@ -0,0 +1,52 @@
+.TH LIBTELLURIAN_GAUSSIAN_RADIUS 3 libtellurian
+.SH NAME
+libtellurian_gaussian_radius \- Calculate a Gaussian radius of curvature
+
+.SH SYNPOSIS
+.nf
+#include <libtellurian.h>
+
+double libtellurian_gaussian_radius(double \fIlatitude\fP);
+
+double libtellurian_gaussian_radius_radians(double \fIlatitude\fP);
+.fi
+.PP
+Link with
+.I -ltellurian
+.IR -lm .
+
+.SH DESCRIPTION
+The
+.BR libtellurian_gaussian_radius ()
+function calculates the Earth's Gaussian radius of
+curvature at a given
+.IR latitude ,
+that is, the geometric mean of the meridian and
+the transverse radii of curvature at the given
+.IR latitude .
+.PP
+The
+.I latitude
+shall be specified according to GPS and in degrees.
+.PP
+The
+.BR libtellurian_gaussian_radius_radians ()
+function is identical to the
+.BR libtellurian_gaussian_radius ()
+function except that
+.I latitude
+shall be specified in radians.
+
+.SH RETURN VALUE
+The
+.BR libtellurian_gaussian_radius ()
+and
+.BR libtellurian_gaussian_radius_radians ()
+functions return the radius of curvature
+measured in meters.
+
+.SH ERRORS
+None.
+
+.SH SEE ALSO
+.BR libtellurian (7)
diff --git a/libtellurian_gaussian_radius.c b/libtellurian_gaussian_radius.c
index 24c9b54..3c1c96b 100644
--- a/libtellurian_gaussian_radius.c
+++ b/libtellurian_gaussian_radius.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -8,3 +9,33 @@ libtellurian_gaussian_radius(double latitude)
latitude = radians(latitude);
return libtellurian_gaussian_radius_radians(latitude);
}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_gaussian_radius(0) == libtellurian_gaussian_radius_radians(0));
+ ASSERT(approx(libtellurian_gaussian_radius(180), libtellurian_gaussian_radius_radians(0)));
+ ASSERT(approx(libtellurian_gaussian_radius(360), libtellurian_gaussian_radius_radians(0)));
+ ASSERT(approx(libtellurian_gaussian_radius(90), libtellurian_gaussian_radius(-90)));
+ ASSERT(approx(libtellurian_gaussian_radius(90), libtellurian_gaussian_radius(270)));
+ ASSERT(approx(libtellurian_gaussian_radius(45), libtellurian_gaussian_radius(-45)));
+ ASSERT(approx(libtellurian_gaussian_radius(135), libtellurian_gaussian_radius(45)));
+ ASSERT(approx(libtellurian_gaussian_radius(-135), libtellurian_gaussian_radius(135)));
+ ASSERT(libtellurian_gaussian_radius(90) == libtellurian_gaussian_radius_radians(D90));
+ ASSERT(libtellurian_gaussian_radius(45) == libtellurian_gaussian_radius_radians(D45));
+ ASSERT(libtellurian_gaussian_radius(30) == libtellurian_gaussian_radius_radians(D30));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_gaussian_radius_radians.c b/libtellurian_gaussian_radius_radians.c
index d61319a..02d0586 100644
--- a/libtellurian_gaussian_radius_radians.c
+++ b/libtellurian_gaussian_radius_radians.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -7,10 +8,48 @@ libtellurian_gaussian_radius_radians(double latitude)
{
double a = LIBTELLURIAN_EQUATORIAL_RADIUS;
double b = LIBTELLURIAN_POLAR_RADIUS;
- double f = 1.0 - b / a;
- double neg_e2 = (f - 2.0) * f;
- double neg_e2_plus_1 = fma(f - 2.0, f, 1.0);
+ double c = b / a;
+ double neg_e2 = fma(c, c, -1.0);
+ double neg_e2_plus_1 = c * c;
double s = sin(latitude);
double denom = fma(neg_e2, s * s, 1.0);
return a * sqrt(neg_e2_plus_1) / denom;
}
+
+
+#else
+
+
+static int
+check(double lat)
+{
+ double g = libtellurian_gaussian_radius_radians(lat);
+ double m = libtellurian_meridian_radius_radians(lat);
+ double n = libtellurian_transverse_radius_radians(lat);
+ double e = sqrt(m * n);
+ double err = fabs((g - e) / e);
+ return err <= 1.0e-8;
+}
+
+int
+main(void)
+{
+ ASSERT(check(0));
+ ASSERT(check(D30));
+ ASSERT(check(D45));
+ ASSERT(check(D60));
+ ASSERT(check(D90));
+ ASSERT(check(D180));
+ ASSERT(check(-D30));
+ ASSERT(check(-D45));
+ ASSERT(check(-D60));
+ ASSERT(check(-D90));
+ ASSERT(check(-D180));
+ ASSERT(check(1.0));
+ ASSERT(check(2.0));
+ ASSERT(check(3.0));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_meridan_radius.c b/libtellurian_meridan_radius.c
deleted file mode 100644
index 01c597f..0000000
--- a/libtellurian_meridan_radius.c
+++ /dev/null
@@ -1,10 +0,0 @@
-/* See LICENSE file for copyright and license details. */
-#include "common.h"
-
-
-double
-libtellurian_meridan_radius(double latitude)
-{
- latitude = radians(latitude);
- return libtellurian_meridan_radius_radians(latitude);
-}
diff --git a/libtellurian_meridan_radius_radians.3 b/libtellurian_meridan_radius_radians.3
deleted file mode 120000
index f8247e4..0000000
--- a/libtellurian_meridan_radius_radians.3
+++ /dev/null
@@ -1 +0,0 @@
-libtellurian_meridan_radius.3 \ No newline at end of file
diff --git a/libtellurian_meridan_radius_radians.c b/libtellurian_meridan_radius_radians.c
deleted file mode 100644
index 0a832e5..0000000
--- a/libtellurian_meridan_radius_radians.c
+++ /dev/null
@@ -1,14 +0,0 @@
-/* See LICENSE file for copyright and license details. */
-#include "common.h"
-
-
-double
-libtellurian_meridan_radius_radians(double latitude)
-{
- double a = LIBTELLURIAN_EQUATORIAL_RADIUS;
- double b = LIBTELLURIAN_POLAR_RADIUS;
- double f = 1.0 - b / a;
- double neg_e2 = (f - 2.0) * f;
- double s = sin(latitude);
- return fma(a, neg_e2, a) / pow(fma(neg_e2, s * s, 1.0), 1.5);
-}
diff --git a/libtellurian_meridian_radius.3 b/libtellurian_meridian_radius.3
new file mode 100644
index 0000000..a13b12b
--- /dev/null
+++ b/libtellurian_meridian_radius.3
@@ -0,0 +1,50 @@
+.TH LIBTELLURIAN_MERIDIAN_RADIUS 3 libtellurian
+.SH NAME
+libtellurian_meridian_radius \- Calculate a meridian radius of curvature
+
+.SH SYNPOSIS
+.nf
+#include <libtellurian.h>
+
+double libtellurian_meridian_radius(double \fIlatitude\fP);
+
+double libtellurian_meridian_radius_radians(double \fIlatitude\fP);
+.fi
+.PP
+Link with
+.I -ltellurian
+.IR -lm .
+
+.SH DESCRIPTION
+The
+.BR libtellurian_meridian_radius ()
+function calculates the Earth's meridian radius of
+curvature (meridional radius of curvature; the
+curvature in the north–south direction), at a given
+.IR latitude .
+.PP
+The
+.I latitude
+shall be specified according to GPS and in degrees.
+.PP
+The
+.BR libtellurian_meridian_radius_radians ()
+function is identical to the
+.BR libtellurian_meridian_radius ()
+function except that
+.I latitude
+shall be specified in radians.
+
+.SH RETURN VALUE
+The
+.BR libtellurian_meridian_radius ()
+and
+.BR libtellurian_meridian_radius_radians ()
+functions return the radius of curvature
+measured in meters.
+
+.SH ERRORS
+None.
+
+.SH SEE ALSO
+.BR libtellurian (7)
diff --git a/libtellurian_meridian_radius.c b/libtellurian_meridian_radius.c
new file mode 100644
index 0000000..f75862a
--- /dev/null
+++ b/libtellurian_meridian_radius.c
@@ -0,0 +1,41 @@
+/* See LICENSE file for copyright and license details. */
+#include "common.h"
+#ifndef TEST
+
+
+double
+libtellurian_meridian_radius(double latitude)
+{
+ latitude = radians(latitude);
+ return libtellurian_meridian_radius_radians(latitude);
+}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_meridian_radius(0) == libtellurian_meridian_radius_radians(0));
+ ASSERT(approx(libtellurian_meridian_radius(180), libtellurian_meridian_radius_radians(0)));
+ ASSERT(approx(libtellurian_meridian_radius(360), libtellurian_meridian_radius_radians(0)));
+ ASSERT(approx(libtellurian_meridian_radius(90), libtellurian_meridian_radius(-90)));
+ ASSERT(approx(libtellurian_meridian_radius(90), libtellurian_meridian_radius(270)));
+ ASSERT(approx(libtellurian_meridian_radius(45), libtellurian_meridian_radius(-45)));
+ ASSERT(approx(libtellurian_meridian_radius(135), libtellurian_meridian_radius(45)));
+ ASSERT(approx(libtellurian_meridian_radius(-135), libtellurian_meridian_radius(135)));
+ ASSERT(libtellurian_meridian_radius(90) == libtellurian_meridian_radius_radians(D90));
+ ASSERT(libtellurian_meridian_radius(45) == libtellurian_meridian_radius_radians(D45));
+ ASSERT(libtellurian_meridian_radius(30) == libtellurian_meridian_radius_radians(D30));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_meridian_radius_radians.3 b/libtellurian_meridian_radius_radians.3
new file mode 120000
index 0000000..1473698
--- /dev/null
+++ b/libtellurian_meridian_radius_radians.3
@@ -0,0 +1 @@
+libtellurian_meridian_radius.3 \ No newline at end of file
diff --git a/libtellurian_meridian_radius_radians.c b/libtellurian_meridian_radius_radians.c
new file mode 100644
index 0000000..f550289
--- /dev/null
+++ b/libtellurian_meridian_radius_radians.c
@@ -0,0 +1,45 @@
+/* See LICENSE file for copyright and license details. */
+#include "common.h"
+#ifndef TEST
+
+
+double
+libtellurian_meridian_radius_radians(double latitude)
+{
+ double a = LIBTELLURIAN_EQUATORIAL_RADIUS;
+ double b = LIBTELLURIAN_POLAR_RADIUS;
+ double neg_e2 = fma(b / a, b / a, -1.0);
+ double s = sin(latitude);
+ return fma(a, neg_e2, a) / pow(fma(neg_e2, s * s, 1.0), 1.5);
+}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ double x = 0.9933056200098025 * LIBTELLURIAN_EQUATORIAL_RADIUS;
+ ASSERT(fabs(libtellurian_meridian_radius_radians(0) - x) / x < 1.0e-8);
+ ASSERT(approx(libtellurian_meridian_radius_radians(D30), libtellurian_meridian_radius_radians(-D30)));
+ ASSERT(approx(libtellurian_meridian_radius_radians(D45), libtellurian_meridian_radius_radians(-D45)));
+ ASSERT(approx(libtellurian_meridian_radius_radians(D60), libtellurian_meridian_radius_radians(-D60)));
+ ASSERT(approx(libtellurian_meridian_radius_radians(D90), libtellurian_meridian_radius_radians(-D90)));
+ ASSERT(approx(libtellurian_meridian_radius_radians(D180), libtellurian_meridian_radius_radians(0)));
+ ASSERT(approx(libtellurian_meridian_radius_radians(-D180), libtellurian_meridian_radius_radians(0)));
+ ASSERT(approx(libtellurian_meridian_radius_radians(D90 + D180), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(approx(libtellurian_meridian_radius_radians(-D90), libtellurian_meridian_radius_radians(D90)));
+ ASSERT(fabs(libtellurian_meridian_radius_radians(D90) - 6399593.625758674) / 6399593.625758674 < 1.0e-8);
+ ASSERT(fabs(libtellurian_meridian_radius_radians(D30) - 6351377.103715289) / 6351377.103715289 < 1.0e-8);
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_normal_gravity.3 b/libtellurian_normal_gravity.3
new file mode 100644
index 0000000..cc3a6fa
--- /dev/null
+++ b/libtellurian_normal_gravity.3
@@ -0,0 +1,53 @@
+.TH LIBTELLURIAN_NORMAL_GRAVITY 3 libtellurian
+.SH NAME
+libtellurian_normal_gravity \- Calculate normal gravity
+
+.SH SYNPOSIS
+.nf
+#include <libtellurian.h>
+
+double libtellurian_normal_gravity(double \fIlatitude\fP);
+
+double libtellurian_normal_gravity_radians(double \fIlatitude\fP);
+.fi
+.PP
+Link with
+.I -ltellurian
+.IR -lm .
+
+.SH DESCRIPTION
+The
+.BR libtellurian_normal_gravity ()
+function models the Earth as a spheroid and calculates
+the normal (theoretical) gravity for any point on a given
+.I latitude
+at the surface. This acceleration does not account for
+the angular velocity of the Earth's rotation.
+.PP
+The
+.I latitude
+shall be specified according to GPS and in degrees.
+.PP
+The
+.BR libtellurian_normal_gravity_radians ()
+function is identical to the
+.BR libtellurian_normal_gravity ()
+function except that
+.I latitude
+shall be specified in radians.
+
+.SH RETURN VALUE
+The
+.BR libtellurian_normal_gravity ()
+and
+.BR libtellurian_normal_gravity_radians ()
+functions return the normal gravity, at the specified
+latitude, measured in meters per square second.
+
+.SH ERRORS
+None.
+
+.SH SEE ALSO
+.BR libtellurian (7),
+.BR libtellurian_effective_gravity (3),
+.BR libtellurian_elevated_gravity (3)
diff --git a/libtellurian_normal_gravity.c b/libtellurian_normal_gravity.c
index 27722b1..a0d3e35 100644
--- a/libtellurian_normal_gravity.c
+++ b/libtellurian_normal_gravity.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -8,3 +9,33 @@ libtellurian_normal_gravity(double latitude)
latitude = radians(latitude);
return libtellurian_normal_gravity_radians(latitude);
}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_normal_gravity(0) == libtellurian_normal_gravity_radians(0));
+ ASSERT(libtellurian_normal_gravity(180) == libtellurian_normal_gravity_radians(0));
+ ASSERT(libtellurian_normal_gravity(360) == libtellurian_normal_gravity_radians(0));
+ ASSERT(approx(libtellurian_normal_gravity(90), libtellurian_normal_gravity(-90)));
+ ASSERT(approx(libtellurian_normal_gravity(90), libtellurian_normal_gravity(270)));
+ ASSERT(approx(libtellurian_normal_gravity(45), libtellurian_normal_gravity(-45)));
+ ASSERT(approx(libtellurian_normal_gravity(135), libtellurian_normal_gravity(45)));
+ ASSERT(approx(libtellurian_normal_gravity(-135), libtellurian_normal_gravity(135)));
+ ASSERT(libtellurian_normal_gravity(90) == libtellurian_normal_gravity_radians(D90));
+ ASSERT(libtellurian_normal_gravity(45) == libtellurian_normal_gravity_radians(D45));
+ ASSERT(libtellurian_normal_gravity(30) == libtellurian_normal_gravity_radians(D30));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_normal_gravity_radians.c b/libtellurian_normal_gravity_radians.c
index 3fdb037..d1d55c8 100644
--- a/libtellurian_normal_gravity_radians.c
+++ b/libtellurian_normal_gravity_radians.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -8,11 +9,41 @@ libtellurian_normal_gravity_radians(double latitude)
double a = LIBTELLURIAN_EQUATORIAL_RADIUS;
double b = LIBTELLURIAN_POLAR_RADIUS;
double neg_e2 = fma(b / a, b / a, -1.0);
- double ag = a * LIBTELLURIAN_EQUATORIAL_GRAVITY;
- double bg = b * LIBTELLURIAN_POLAR_GRAVITY;
+ double ag = a * LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY;
+ double bg = b * LIBTELLURIAN_NORMAL_POLAR_GRAVITY;
double k = bg / ag - 1.0;
double sin2_phi = sin(latitude) * sin(latitude);
double num = fma(k, sin2_phi, 1.0);
double denom2 = fma(neg_e2, sin2_phi, 1.0);
- return LIBTELLURIAN_EQUATORIAL_GRAVITY * num / sqrt(denom2);
+ return LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY * num / sqrt(denom2);
}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a / b - 1.0) <= 1e-12;
+}
+
+int
+main(void)
+{
+ ASSERT(approx(libtellurian_normal_gravity_radians(D30), libtellurian_normal_gravity_radians(-D30)));
+ ASSERT(approx(libtellurian_normal_gravity_radians(D45), libtellurian_normal_gravity_radians(-D45)));
+ ASSERT(approx(libtellurian_normal_gravity_radians(D60), libtellurian_normal_gravity_radians(-D60)));
+ ASSERT(approx(libtellurian_normal_gravity_radians(D90), libtellurian_normal_gravity_radians(-D90)));
+ ASSERT(approx(libtellurian_normal_gravity_radians(D180), libtellurian_normal_gravity_radians(0)));
+ ASSERT(approx(libtellurian_normal_gravity_radians(-D180), libtellurian_normal_gravity_radians(0)));
+ ASSERT(approx(libtellurian_normal_gravity_radians(2.0), libtellurian_normal_gravity_radians(-2.0)));
+ ASSERT(approx(libtellurian_normal_gravity_radians(0), LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY));
+ ASSERT(approx(libtellurian_normal_gravity_radians(D90), LIBTELLURIAN_NORMAL_POLAR_GRAVITY));
+ ASSERT(libtellurian_normal_gravity_radians(D45) < LIBTELLURIAN_NORMAL_POLAR_GRAVITY);
+ ASSERT(libtellurian_normal_gravity_radians(D45) > LIBTELLURIAN_NORMAL_EQUATORIAL_GRAVITY);
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_sea_level.3 b/libtellurian_sea_level.3
index 80dc6ea..9401641 100644
--- a/libtellurian_sea_level.3
+++ b/libtellurian_sea_level.3
@@ -3,7 +3,7 @@
libtellurian_sea_level \- Calculate a geocentric radius
.SH SYNPOSIS
-.ni
+.nf
#include <libtellurian.h>
double libtellurian_sea_level(double \fIlatitude\fP);
@@ -19,7 +19,7 @@ Link with
The
.BR libtellurian_sea_level ()
function calculates the distance from the centre of
-the Earth the ellipsoid, for any point on a the given
+the Earth the ellipsoid, for any point on the given
.IR latitude ,
which givens the nominal sea level at this position,
which is close to the actual sea level, and is the
diff --git a/libtellurian_sea_level.c b/libtellurian_sea_level.c
index 76f817f..03893e2 100644
--- a/libtellurian_sea_level.c
+++ b/libtellurian_sea_level.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -8,3 +9,33 @@ libtellurian_sea_level(double latitude)
latitude = radians(latitude);
return libtellurian_sea_level_radians(latitude);
}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_sea_level(0) == libtellurian_sea_level_radians(0));
+ ASSERT(approx(libtellurian_sea_level(180), libtellurian_sea_level_radians(0)));
+ ASSERT(approx(libtellurian_sea_level(360), libtellurian_sea_level_radians(0)));
+ ASSERT(approx(libtellurian_sea_level(90), libtellurian_sea_level(-90)));
+ ASSERT(approx(libtellurian_sea_level(90), libtellurian_sea_level(270)));
+ ASSERT(approx(libtellurian_sea_level(45), libtellurian_sea_level(-45)));
+ ASSERT(approx(libtellurian_sea_level(135), libtellurian_sea_level(45)));
+ ASSERT(approx(libtellurian_sea_level(-135), libtellurian_sea_level(135)));
+ ASSERT(libtellurian_sea_level(90) == libtellurian_sea_level_radians(D90));
+ ASSERT(libtellurian_sea_level(45) == libtellurian_sea_level_radians(D45));
+ ASSERT(libtellurian_sea_level(30) == libtellurian_sea_level_radians(D30));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_sea_level_radians.c b/libtellurian_sea_level_radians.c
index a23b222..e536025 100644
--- a/libtellurian_sea_level_radians.c
+++ b/libtellurian_sea_level_radians.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -15,3 +16,26 @@ libtellurian_sea_level_radians(double latitude)
double denom = fma(x, c, y * s);
return sqrt(num / denom);
}
+
+
+#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
diff --git a/libtellurian_transverse_radius.3 b/libtellurian_transverse_radius.3
new file mode 100644
index 0000000..aa24be4
--- /dev/null
+++ b/libtellurian_transverse_radius.3
@@ -0,0 +1,50 @@
+.TH LIBTELLURIAN_TRANSVERSE_RADIUS 3 libtellurian
+.SH NAME
+libtellurian_transverse_radius \- Calculate a transverse radius of curvature
+
+.SH SYNPOSIS
+.nf
+#include <libtellurian.h>
+
+double libtellurian_transverse_radius(double \fIlatitude\fP);
+
+double libtellurian_transverse_radius_radians(double \fIlatitude\fP);
+.fi
+.PP
+Link with
+.I -ltellurian
+.IR -lm .
+
+.SH DESCRIPTION
+The
+.BR libtellurian_transverse_radius ()
+function calculates the Earth's transverse radius of
+curvature (prime-vertial radius of curvature; the
+curvature in the east–west direction), at a given
+.IR latitude .
+.PP
+The
+.I latitude
+shall be specified according to GPS and in degrees.
+.PP
+The
+.BR libtellurian_transverse_radius_radians ()
+function is identical to the
+.BR libtellurian_transverse_radius ()
+function except that
+.I latitude
+shall be specified in radians.
+
+.SH RETURN VALUE
+The
+.BR libtellurian_transverse_radius ()
+and
+.BR libtellurian_transverse_radius_radians ()
+functions return the radius of curvature
+measured in meters.
+
+.SH ERRORS
+None.
+
+.SH SEE ALSO
+.BR libtellurian (7)
diff --git a/libtellurian_transverse_radius.c b/libtellurian_transverse_radius.c
index d94ca86..15bc2f5 100644
--- a/libtellurian_transverse_radius.c
+++ b/libtellurian_transverse_radius.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -8,3 +9,33 @@ libtellurian_transverse_radius(double latitude)
latitude = radians(latitude);
return libtellurian_transverse_radius_radians(latitude);
}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_transverse_radius(0) == libtellurian_transverse_radius_radians(0));
+ ASSERT(approx(libtellurian_transverse_radius(180), libtellurian_transverse_radius_radians(0)));
+ ASSERT(approx(libtellurian_transverse_radius(360), libtellurian_transverse_radius_radians(0)));
+ ASSERT(approx(libtellurian_transverse_radius(90), libtellurian_transverse_radius(-90)));
+ ASSERT(approx(libtellurian_transverse_radius(90), libtellurian_transverse_radius(270)));
+ ASSERT(approx(libtellurian_transverse_radius(45), libtellurian_transverse_radius(-45)));
+ ASSERT(approx(libtellurian_transverse_radius(135), libtellurian_transverse_radius(45)));
+ ASSERT(approx(libtellurian_transverse_radius(-135), libtellurian_transverse_radius(135)));
+ ASSERT(libtellurian_transverse_radius(90) == libtellurian_transverse_radius_radians(D90));
+ ASSERT(libtellurian_transverse_radius(45) == libtellurian_transverse_radius_radians(D45));
+ ASSERT(libtellurian_transverse_radius(30) == libtellurian_transverse_radius_radians(D30));
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_transverse_radius_radians.c b/libtellurian_transverse_radius_radians.c
index dae8f04..c4f80aa 100644
--- a/libtellurian_transverse_radius_radians.c
+++ b/libtellurian_transverse_radius_radians.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
double
@@ -7,8 +8,35 @@ libtellurian_transverse_radius_radians(double latitude)
{
double a = LIBTELLURIAN_EQUATORIAL_RADIUS;
double b = LIBTELLURIAN_POLAR_RADIUS;
- double f = 1.0 - b / a;
- double neg_e2 = (f - 2.0) * f;
+ double neg_e2 = fma(b / a, b / a, -1.0);
double s = sin(latitude);
- return pow(fma(neg_e2, s * s, 1.0), -0.5);
+ return a * pow(fma(neg_e2, s * s, 1.0), -0.5);
}
+
+
+#else
+
+
+static int
+approx(double a, double b)
+{
+ return fabs(a - b) <= 1e-8 * (0.5 * (a + b));
+}
+
+int
+main(void)
+{
+ ASSERT(libtellurian_transverse_radius_radians(0) == LIBTELLURIAN_EQUATORIAL_RADIUS);
+ ASSERT(approx(libtellurian_transverse_radius_radians(D180), LIBTELLURIAN_EQUATORIAL_RADIUS));
+ ASSERT(approx(libtellurian_transverse_radius_radians(-D180), LIBTELLURIAN_EQUATORIAL_RADIUS));
+ ASSERT(approx(libtellurian_transverse_radius_radians(D30), libtellurian_transverse_radius_radians(-D30)));
+ ASSERT(approx(libtellurian_transverse_radius_radians(D45), libtellurian_transverse_radius_radians(-D45)));
+ ASSERT(approx(libtellurian_transverse_radius_radians(D60), libtellurian_transverse_radius_radians(-D60)));
+ ASSERT(approx(libtellurian_transverse_radius_radians(D90), libtellurian_transverse_radius_radians(-D90)));
+ ASSERT(fabs(libtellurian_transverse_radius_radians(D90) - 6399593.625758674) / 6399593.625758674 < 1.0e-8);
+ ASSERT(fabs(libtellurian_transverse_radius_radians(D30) - 6383480.917690154) / 6383480.917690154 < 1.0e-8);
+ return 0;
+}
+
+
+#endif
diff --git a/libtellurian_vincenty_inverse__.c b/libtellurian_vincenty_inverse__.c
index 746011f..c50752e 100644
--- a/libtellurian_vincenty_inverse__.c
+++ b/libtellurian_vincenty_inverse__.c
@@ -1,5 +1,6 @@
/* See LICENSE file for copyright and license details. */
#include "common.h"
+#ifndef TEST
/* TODO complete the implementation */
@@ -102,3 +103,8 @@ libtellurian_vincenty_inverse__(double latitude1, double longitude1, double lati
* https://en.wikipedia.org/wiki/Vincenty's_formulae#Nearly_antipodal_points
*/
}
+
+
+#else
+TODO_TEST
+#endif