diff options
Diffstat (limited to '')
50 files changed, 1564 insertions, 92 deletions
@@ -12,3 +12,5 @@ *.gcov *.gcno *.gcda +*.to +*.test @@ -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 @@ -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. @@ -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 @@ -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 @@ -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 |