diff --git a/CHANGELOG.md b/CHANGELOG.md index e5912459..ef986776 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,7 +7,7 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.1.0/), [Semantic Versioning](https://semver.org/spec/v2.0.0.html). -## [1.2.0-rc5] - 2025-01-01 +## [Unreleased] Release candidate for the next feature release, expected around 1 February 2025. @@ -80,6 +80,9 @@ Release candidate for the next feature release, expected around 1 February 2025. - #98: Added `gcrs_to_tod()` / `tod_to_gcrs()` and `gcrs_to_mod()` / `mod_to_gcrs()` vector conversion functions for convenience. + + - #106: New example files under `examples/` demonstrating the recommended approach for using SuperNOVAS to calculate + positions for various types of object. - Added various `object` initializer macros in `novas.h` for the major planets, Sun, Moon, and barycenters, e.g. `NOVAS_EARTH_INIT` or `NOVAS_SSB_INIT`. These wrap the parametric `NOVAS_PLANET_INIT(num, name)` macro, and can be @@ -112,6 +115,9 @@ Release candidate for the next feature release, expected around 1 February 2025. - #97: Updated `NOVAS_PLANETS`, `NOVAS_PLANET_NAMES_INIT`, and `NOVAS_RMASS_INIT` macros to include the added planet constants. + - #106: The old (legacy) NOVAS C example has been removed. Instead a new set of examples are provided, which are + better suited for SuperNOVAS. + - `make check` now runs both static analysis by cppcheck (new `analysis` target) and regression tests (`test` target), in closer conformance to GNU Makefile standards. diff --git a/examples/example-calceph b/examples/example-calceph new file mode 100755 index 00000000..87c613c5 Binary files /dev/null and b/examples/example-calceph differ diff --git a/examples/example-calceph.c b/examples/example-calceph.c new file mode 100644 index 00000000..cf509d5e --- /dev/null +++ b/examples/example-calceph.c @@ -0,0 +1,179 @@ +/** + * @file + * + * @date Created on Jan 9, 2025 + * @author Attila Kovacs + * + * Example file for using the SuperNOVAS C/C++ library for determining positions for + * Solary-system sources, with the CALCEPH C library providing access to ephemeris + * files. + * + * You will need access to the CALCEPH library (unversioned `libcalceph.so` or else + * `libcalceph.a`) and C headers (`calceph.h`), and the SuperNOVAS + * `libsolsys-calceph.so` (or `libsoslsys-calceph.a`) module. + * + * Link with: + * + * ``` + * -lsupernovas -lsolsys-calceph -lcalceph + * ``` + */ + +#include +#include +#include + +#include ///< SuperNOVAS functions and definitions +#include ///< CALCEPH adapter functions to SuperNOVAS + + +// Below are some Earth orientation values. Here we define them as constants, but they may +// of course be variables. They should be set to the appropriate values for the time +// of observation based on the IERS Bulletins... + +#define LEAP_SECONDS 37 ///< [s] current leap seconds +#define DUT1 0.114 ///< [s] current UT1 - UTC time difference from IERS Bulletin A +#define POLAR_DX 230.0 ///< [mas] Earth polar offset x, e.g. from IERS Bulletin A. +#define POLAR_DY -62.0 ///< [mas] Earth polar offset y, e.g. from IERS Bulletin A. + +int main() { + // SuperNOVAS aariables used for the calculations -------------------------> + cat_entry star = CAT_ENTRY_INIT; // catalog information about a sidereal source + object source; // observed source + observer obs; // observer location + novas_timespec obs_time; // astrometric time of observation + novas_frame obs_frame; // observing frame defined for observing time and location + enum novas_accuracy accuracy; // NOVAS_FULL_ACCURACY or NOVAS_REDUCED_ACCURACY + sky_pos apparent; // calculated precise observed (apparent) position of source + + // Calculated quantities -------------------------------------------------> + double az, el; // calculated azimuth and elevation at observing site + + // Intermediate variables we'll use --------------------------------------> + t_calcephbin *de440; // CALCEPH ephemeris binary + struct timespec unix_time; // Standard precision UNIX time structure + + + // We'll print debugging messages and error traces... + novas_debug(NOVAS_DEBUG_ON); + + // ------------------------------------------------------------------------- + // We'll use the CALCEPH library to provide ephemeris data + + // First open one or more ephemeris files with CALCEPH to use + // E.g. the DE440 (short-term) ephemeris data from JPL. + de440 = calceph_open("path/to/de440s.bsp"); + if(!de440) { + fprintf(stderr, "ERROR! could not open ephemeris data\n"); + return 1; + } + + // Make de440 provide ephemeris data for the major planets. + novas_use_calceph_planets(de440); + + // We could specify to use a calceph ephemeris binary for generic + // Solar-systems sources also (including planets too if + // novas_use_calceph_planets() is not called separately + // + // E.g. Jovian satellites: + // t_calcephbin jovian = calceph_open("/path/to/jup365.bsp"); + // novas_use_calceph(jovian); + + // Since we have an ephemeris provider for major planets, we can unlock the + // ultimate accuracy of SuperNOVAS. + accuracy = NOVAS_FULL_ACCURACY; // sub-uas precision + + + // ------------------------------------------------------------------------- + // Define a Solar-system source + + // To define a major planet (or Sun, Moon, SSB, or EMB): + if(make_planet(NOVAS_MARS, &source) != 0) { + fprintf(stderr, "ERROR! defining planet.\n"); + return 1; + } + + // ... Or, to define a minor body, such as an asteroid or satellite + // with a name and ID number. + /* + if(make_ephem_object("Io", 501, &source) != 0) { + fprintf(stderr, "ERROR! defining ephemeris body.\n"); + return 1; + } + */ + + /* + // If the object uses CALCEPH IDs instead of NAIF, then + novas_calceph_use_ids(NOVAS_ID_CALCEPH); + */ + + + // ------------------------------------------------------------------------- + // Define observer somewhere on Earth (we can also define observers in Earth + // or Sun orbit, at the geocenter or at the Solary-system barycenter...) + + // Specify the location we are observing from + // 50.7374 deg N, 7.0982 deg E, 60m elevation + if(make_observer_on_surface(50.7374, 7.0982, 60.0, 0.0, 0.0, &obs) != 0) { + fprintf(stderr, "ERROR! defining Earth-based observer location.\n"); + return 1; + } + + + // ------------------------------------------------------------------------- + // Set the astrometric time of observation... + + // Get the current system time, with up to nanosecond resolution... + clock_gettime(CLOCK_REALTIME, &unix_time); + + // Set the time of observation to the precise UTC-based UNIX time + if(novas_set_unix_time(unix_time.tv_sec, unix_time.tv_nsec, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; + } + + // ... Or you could set a time explicily in any known timescale. + /* + // Let's set a TDB-based time for the start of the J2000 epoch exactly... + if(novas_set_time(NOVAS_TDB, NOVAS_JD_J2000, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; + } + */ + + // ------------------------------------------------------------------------- + // Initialize the observing frame with the given observing and Earth + // orientation patameters. + // + if(novas_make_frame(accuracy, &obs, &obs_time, POLAR_DX, POLAR_DY, &obs_frame) != 0) { + fprintf(stderr, "ERROR! failed to define observing frame.\n"); + return 1; + } + + + // ------------------------------------------------------------------------- + // Calculate the precise apparent position (here in CIRS coordinates) + if(novas_sky_pos(&source, &obs_frame, NOVAS_CIRS, &apparent) != 0) { + fprintf(stderr, "ERROR! failed to calculate apparent position.\n"); + return 1; + } + + // Let's print the apparent position + printf(" RA = %.9f h, Dec = %.9f deg, rad_vel = %.6f km/s\n", apparent.ra, apparent.dec, apparent.rv); + + + // ------------------------------------------------------------------------- + // Convert the apparent position in CIRS on sky to horizontal coordinates + // We'll use a standard (fixed) atmospheric model to estimate an optical refraction + // (You might use other refraction models, or NULL to ignore refraction corrections) + if(novas_app_to_hor(&obs_frame, NOVAS_CIRS, apparent.ra, apparent.dec, novas_standard_refraction, &az, &el) != 0) { + fprintf(stderr, "ERROR! failed to calculate azimuth / elevation.\n"); + return 1; + } + + // Let's print the calculated azimuth and elevation + printf(" Az = %.6f deg, El = %.6f deg\n", az, el); + + return 0; +} + diff --git a/examples/example-cspice b/examples/example-cspice new file mode 100755 index 00000000..c333358c Binary files /dev/null and b/examples/example-cspice differ diff --git a/examples/example-cspice.c b/examples/example-cspice.c new file mode 100644 index 00000000..cafc9a41 --- /dev/null +++ b/examples/example-cspice.c @@ -0,0 +1,175 @@ +/** + * @file + * + * @date Created on Jan 9, 2025 + * @author Attila Kovacs + * + * Example file for using the SuperNOVAS C/C++ library for determining positions for + * Solary-system sources, with the NAIF CSPICE toolkit providing access to ephemeris + * files. + * + * You will need access to the NAIF CSPICE library (unversioned `libcspice.so` or else + * `libcspice.a`) and C headers (`cspice/*.h`), and the SuperNOVAS `libsolsys-cspice.so` + * (or `libsoslsys-cspice.a`) module. + * + * To compile CSPICE as a shared (.so) library, you may want to check out the GitHub + * repository: + * + * - https://github.com/Smithsonian/cspice-sharedlib + * + * Link with: + * + * ``` + * -lsupernovas -lsolsys-cspice -lcspice + * ``` + */ + +#include +#include +#include + +#include ///< SuperNOVAS functions and definitions +#include ///< CSPICE adapter functions to SuperNOVAS + + +// Below are some Earth orientation values. Here we define them as constants, but they may +// of course be variables. They should be set to the appropriate values for the time +// of observation based on the IERS Bulletins... + +#define LEAP_SECONDS 37 ///< [s] current leap seconds +#define DUT1 0.114 ///< [s] current UT1 - UTC time difference from IERS Bulletin A +#define POLAR_DX 230.0 ///< [mas] Earth polar offset x, e.g. from IERS Bulletin A. +#define POLAR_DY -62.0 ///< [mas] Earth polar offset y, e.g. from IERS Bulletin A. + +int main() { + // SuperNOVAS aariables used for the calculations -------------------------> + cat_entry star = CAT_ENTRY_INIT; // catalog information about a sidereal source + object source; // observed source + observer obs; // observer location + novas_timespec obs_time; // astrometric time of observation + novas_frame obs_frame; // observing frame defined for observing time and location + enum novas_accuracy accuracy; // NOVAS_FULL_ACCURACY or NOVAS_REDUCED_ACCURACY + sky_pos apparent; // calculated precise observed (apparent) position of source + + // Calculated quantities -------------------------------------------------> + double az, el; // calculated azimuth and elevation at observing site + + // Intermediate variables we'll use --------------------------------------> + struct timespec unix_time; // Standard precision UNIX time structure + + + // We'll print debugging messages and error traces... + novas_debug(NOVAS_DEBUG_ON); + + + // ------------------------------------------------------------------------- + // We'll use the NAIF CSPICE Toolkit to provide ephemeris data + + // Open one or more ephemeris files to use...' + // E.g. the DE440 (short-term) ephemeris data from JPL. + if(cspice_add_kernel("path/to/de440s.bsp") != 0) { + fprintf(stderr, "ERROR! could not open ephemeris data\n"); + return 1; + } + + // ... You can open multiple NAIF kernels + // E.g. add Jovian satellites... + // cspice_add_kernel("path/to/jup365.bsp"); + + // Now we can use the loaded ephemeris files for Solar-system objects. + // (major planets and minor bodies alike). + novas_use_cspice(); + + // And, since we have an ephemeris provider for major planets, we can unlock + // the ultimate accuracy of SuperNOVAS. + accuracy = NOVAS_FULL_ACCURACY; // sub-uas precision + + + // ------------------------------------------------------------------------- + // Define a Solar-system source + + // To define a major planet (or Sun, Moon, SSB, or EMB): + if(make_planet(NOVAS_MARS, &source) != 0) { + fprintf(stderr, "ERROR! defining planet.\n"); + return 1; + } + + // ... Or, to define a minor body, such as an asteroid or satellite + // with a name and NAIF ID. + /* + if(make_ephem_object("Io", 501, &source) != 0) { + fprintf(stderr, "ERROR! defining ephemeris body.\n"); + return 1; + } + */ + + + // ------------------------------------------------------------------------- + // Define observer somewhere on Earth (we can also define observers in Earth + // or Sun orbit, at the geocenter or at the Solary-system barycenter...) + + // Specify the location we are observing from + // 50.7374 deg N, 7.0982 deg E, 60m elevation + if(make_observer_on_surface(50.7374, 7.0982, 60.0, 0.0, 0.0, &obs) != 0) { + fprintf(stderr, "ERROR! defining Earth-based observer location.\n"); + return 1; + } + + + // ------------------------------------------------------------------------- + // Set the astrometric time of observation... + + // Get the current system time, with up to nanosecond resolution... + clock_gettime(CLOCK_REALTIME, &unix_time); + + // Set the time of observation to the precise UTC-based UNIX time + if(novas_set_unix_time(unix_time.tv_sec, unix_time.tv_nsec, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; + } + + // ... Or you could set a time explicily in any known timescale. + /* + // Let's set a TDB-based time for the start of the J2000 epoch exactly... + if(novas_set_time(NOVAS_TDB, NOVAS_JD_J2000, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; + } + */ + + // ------------------------------------------------------------------------- + // Initialize the observing frame with the given observing and Earth + // orientation patameters. + // + if(novas_make_frame(accuracy, &obs, &obs_time, POLAR_DX, POLAR_DY, &obs_frame) != 0) { + fprintf(stderr, "ERROR! failed to define observing frame.\n"); + return 1; + } + + + // ------------------------------------------------------------------------- + // Calculate the precise apparent position (here in CIRS coordinates) + if(novas_sky_pos(&source, &obs_frame, NOVAS_CIRS, &apparent) != 0) { + fprintf(stderr, "ERROR! failed to calculate apparent position.\n"); + return 1; + } + + // Let's print the apparent position + printf(" RA = %.9f h, Dec = %.9f deg, rad_vel = %.6f km/s\n", apparent.ra, apparent.dec, apparent.rv); + + + // ------------------------------------------------------------------------- + // Convert the apparent position in CIRS on sky to horizontal coordinates + // We'll use a standard (fixed) atmospheric model to estimate an optical refraction + // (You might use other refraction models, or NULL to ignore refraction corrections) + if(novas_app_to_hor(&obs_frame, NOVAS_CIRS, apparent.ra, apparent.dec, novas_standard_refraction, &az, &el) != 0) { + fprintf(stderr, "ERROR! failed to calculate azimuth / elevation.\n"); + return 1; + } + + // Let's print the calculated azimuth and elevation + printf(" Az = %.6f deg, El = %.6f deg\n", az, el); + + return 0; +} + diff --git a/examples/example-high-z b/examples/example-high-z new file mode 100755 index 00000000..1ceb3570 Binary files /dev/null and b/examples/example-high-z differ diff --git a/examples/example-high-z.c b/examples/example-high-z.c new file mode 100644 index 00000000..6f867540 --- /dev/null +++ b/examples/example-high-z.c @@ -0,0 +1,168 @@ +/** + * @file + * + * @date Created on Jan 9, 2025 + * @author Attila Kovacs + * + * Example file for using the SuperNOVAS C/C++ library for determining positions for + * distant galaxies and quasars, or other high-redshift objects. + * + * It's the same recipe as `example-star.c`, except that we define the object of + * interest a little differently. + * + * Link with: + * + * ``` + * -lsupernovas + * ``` + * + */ + +#include +#include +#include + +#include ///< SuperNOVAS functions and definitions + +// Below are some Earth orientation values. Here we define them as constants, but they may +// of course be variables. They should be set to the appropriate values for the time +// of observation based on the IERS Bulletins... + +#define LEAP_SECONDS 37 ///< [s] current leap seconds +#define DUT1 0.114 ///< [s] current UT1 - UTC time difference from IERS Bulletin A +#define POLAR_DX 230.0 ///< [mas] Earth polar offset x, e.g. from IERS Bulletin A. +#define POLAR_DY -62.0 ///< [mas] Earth polar offset y, e.g. from IERS Bulletin A. + +int main() { + // SuperNOVAS aariables used for the calculations -------------------------> + object source; // a celestial object: sidereal, planet, ephemeris or orbital source + observer obs; // observer location + novas_timespec obs_time; // astrometric time of observation + novas_frame obs_frame; // observing frame defined for observing time and location + enum novas_accuracy accuracy; // NOVAS_FULL_ACCURACY or NOVAS_REDUCED_ACCURACY + sky_pos apparent; // calculated precise observed (apparent) position of source + + // Calculated quantities -------------------------------------------------> + double az, el; // calculated azimuth and elevation at observing site + + // Intermediate variables we'll use --------------------------------------> + struct timespec unix_time; // Standard precision UNIX time structure + + + // We'll print debugging messages and error traces... + novas_debug(NOVAS_DEBUG_ON); + + + // ------------------------------------------------------------------------- + // Define a high-z source. + + // 12h29m6.6997s +2d3m8.598s (ICRS) z=0.158339 + if(make_redshifted_object("3c273", 12.4851944, 2.0523883, 0.158339, &source) != 0) { + fprintf(stderr, "ERROR! defining cat_entry.\n"); + return 1; + } + + // If we did not use ICRS catalog coordinates, we would have to convert them to ICRS... + + /* E.g. change B1950 to the J2000 (FK5) system... + if(transform_cat(CHANGE_EPOCH, NOVAS_JD_B1950, &source.star, NOVAS_JD_J2000, "FK5", &source.star) != 0) { + fprintf(stderr, "ERROR! converting B1950 catalog coordinates to J2000.\n"); + return 1; + } */ + + /* Then convert J2000 coordinates to ICRS (also in place). Here the dates don't matter... + if(transform_cat(CHANGE_J2000_TO_ICRS, 0.0, &source.star, 0.0, "ICRS", &source.star) != 0) { + fprintf(stderr, "ERROR! converting J2000 catalog coordinates to ICRS.\n"); + return 1; + } */ + + + // ------------------------------------------------------------------------- + // Define observer somewhere on Earth (we can also define observers in Earth + // or Sun orbit, at the geocenter or at the Solary-system barycenter...) + + // Specify the location we are observing from + // 50.7374 deg N, 7.0982 deg E, 60m elevation + if(make_observer_on_surface(50.7374, 7.0982, 60.0, 0.0, 0.0, &obs) != 0) { + fprintf(stderr, "ERROR! defining Earth-based observer location.\n"); + return 1; + } + + + // ------------------------------------------------------------------------- + // Set the astrometric time of observation... + + // Get the current system time, with up to nanosecond resolution... + clock_gettime(CLOCK_REALTIME, &unix_time); + + // Set the time of observation to the precise UTC-based UNIX time + if(novas_set_unix_time(unix_time.tv_sec, unix_time.tv_nsec, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; + } + + // ... Or you could set a time explicily in any known timescale. + /* + // Let's set a TDB-based time for the start of the J2000 epoch exactly... + if(novas_set_time(NOVAS_TDB, NOVAS_JD_J2000, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; + } + */ + + + // ------------------------------------------------------------------------- + // You might want to set a provider for precise planet positions so we might + // calculate Earth, Sun and major planet positions accurately. If an planet + // provider is configured, we can unlock the ultimate (sub-uas) accuracy of + // SuperNOVAS. + // + // There are many ways to set a provider of planet positions. For example, + // you may use CALCEPH: + // + // t_calcephbin *planets = calceph_open("path/to/de440s.bsp"); + // novas_use_calceph(planets); + // + // accuracy = NOVAS_FULL_ACCURACY; // sub-uas precision + + // Without a planet provider, we are stuck with reduced (mas) precisions + // only... + accuracy = NOVAS_REDUCED_ACCURACY; // mas-level precision, typically + + + // ------------------------------------------------------------------------- + // Initialize the observing frame with the given observing and Earth + // orientation patameters. + // + if(novas_make_frame(accuracy, &obs, &obs_time, POLAR_DX, POLAR_DY, &obs_frame) != 0) { + fprintf(stderr, "ERROR! failed to define observing frame.\n"); + return 1; + } + + + // ------------------------------------------------------------------------- + // Calculate the precise apparent position (here in CIRS coordinates) + if(novas_sky_pos(&source, &obs_frame, NOVAS_CIRS, &apparent) != 0) { + fprintf(stderr, "ERROR! failed to calculate apparent position.\n"); + return 1; + } + + // Let's print the apparent position + printf(" RA = %.9f h, Dec = %.9f deg, z_obs = %.9f\n", apparent.ra, apparent.dec, novas_v2z(apparent.rv)); + + + // ------------------------------------------------------------------------- + // Convert the apparent position in CIRS on sky to horizontal coordinates + // We'll use a standard (fixed) atmospheric model to estimate an optical refraction + // (You might use other refraction models, or NULL to ignore refraction corrections) + if(novas_app_to_hor(&obs_frame, NOVAS_CIRS, apparent.ra, apparent.dec, novas_standard_refraction, &az, &el) != 0) { + fprintf(stderr, "ERROR! failed to calculate azimuth / elevation.\n"); + return 1; + } + + // Let's print the calculated azimuth and elevation + printf(" Az = %.6f deg, El = %.6f deg\n", az, el); + + return 0; +} + diff --git a/examples/example-orbital b/examples/example-orbital new file mode 100755 index 00000000..bd589fe2 Binary files /dev/null and b/examples/example-orbital differ diff --git a/examples/example-orbital.c b/examples/example-orbital.c new file mode 100644 index 00000000..a4559492 --- /dev/null +++ b/examples/example-orbital.c @@ -0,0 +1,184 @@ +/** + * @file + * + * @date Created on Jan 9, 2025 + * @author Attila Kovacs + * + * Example file for using the SuperNOVAS C/C++ library for determining positions for + * Solar-system objects define through a set of orbital parameters. + * + * For example, the IAU Minor Planet Center (MPC) publishes current orbital + * parameters for known asteroids, comets, and near-Earth objects. While orbitals are + * not super precise in general, they can provide sufficienly accurate positions on + * the arcsecond level (or below), and may be the best/only source of position data + * for newly discovered objects. + * + * See https://minorplanetcenter.net/data + * + * Link with + * + * ``` + * -lsupernovas + * ``` + */ + +#include +#include +#include + +#include ///< SuperNOVAS functions and definitions + + +// Below are some Earth orientation values. Here we define them as constants, but they may +// of course be variables. They should be set to the appropriate values for the time +// of observation based on the IERS Bulletins... + +#define LEAP_SECONDS 37 ///< [s] current leap seconds +#define DUT1 0.114 ///< [s] current UT1 - UTC time difference from IERS Bulletin A +#define POLAR_DX 230.0 ///< [mas] Earth polar offset x, e.g. from IERS Bulletin A. +#define POLAR_DY -62.0 ///< [mas] Earth polar offset y, e.g. from IERS Bulletin A. + +int main() { + // SuperNOVAS aariables used for the calculations -------------------------> + novas_orbital orbit = NOVAS_ORBIT_INIT; // Orbital parameters + object source; // a celestial object: sidereal, planet, ephemeris or orbital source + observer obs; // observer location + novas_timespec obs_time; // astrometric time of observation + novas_frame obs_frame; // observing frame defined for observing time and location + enum novas_accuracy accuracy; // NOVAS_FULL_ACCURACY or NOVAS_REDUCED_ACCURACY + sky_pos apparent; // calculated precise observed (apparent) position of source + + // Calculated quantities -------------------------------------------------> + double az, el; // calculated azimuth and elevation at observing site + + // Intermediate variables we'll use --------------------------------------> + struct timespec unix_time; // Standard precision UNIX time structure + + + // We'll print debugging messages and error traces... + novas_debug(NOVAS_DEBUG_ON); + + + // Orbitals assume Keplerian motion, and are never going to be accurate much below the + // tens of arcsec level even for the most current MPC orbits. Orbitals for planetary + // satellites are even less precise. So, with orbitals, there is no point on pressing + // for ultra-high (sub-uas level) accuracy... + accuracy = NOVAS_REDUCED_ACCURACY; // mas-level precision, typically + + + // ------------------------------------------------------------------------- + // Define a sidereal source + + // Orbital Parameters for the asteroid Ceres from the Minor Planet Center + // (MPC) at JD 2460600.5 + orbit.jd_tdb = 2460600.5; // [day] TDB date + orbit.a = 2.7666197; // [AU] + orbit.e = 0.079184; + orbit.i = 10.5879; // [deg] + orbit.omega = 73.28579; // [deg] + orbit.Omega = 80.25414; // [deg] + orbit.M0 = 145.84905; // [deg] + orbit.n = 0.21418047; // [deg/day] + + // Define Ceres as the observed object (we can use whatever ID numbering + // system here, since it's irrelevant to SuperNOVAS in this context). + make_orbital_object("Ceres", 2000001, &orbit, &source); + + + // ... Or, you could define orbitals for a satellite instead: + /* + // E.g. Callisto's orbital parameters from JPL Horizons + // https://ssd.jpl.nasa.gov/sats/elem/sep.html + // 1882700. 0.007 43.8 87.4 0.3 309.1 16.690440 277.921 577.264 268.7 64.8 + orbit.system.center = NOVAS_JUPITER; + novas_set_orbsys_pole(NOVAS_GCRS, 268.7 / 15.0, 64.8, &orbit->system); + + orbit.jd_tdb = NOVAS_JD_J2000; + orbit.a = 1882700.0 * 1e3 / AU; + orbit.e = 0.007; + orbit.omega = 43.8; + orbit.M0 = 87.4; + orbit.i = 0.3; + orbit.Omega = 309.1; + orbit.n = TWOPI / 16.690440; + orbit.apsis_period = 277.921 * 365.25; + orbit.node_period = 577.264 * 365.25; + + // Set Callisto as the observed object + make_orbital_object("Callisto", 501, &orbit, &source); + */ + + + // ------------------------------------------------------------------------- + // Define observer somewhere on Earth (we can also define observers in Earth + // or Sun orbit, at the geocenter or at the Solary-system barycenter...) + + // Specify the location we are observing from + // 50.7374 deg N, 7.0982 deg E, 60m elevation + if(make_observer_on_surface(50.7374, 7.0982, 60.0, 0.0, 0.0, &obs) != 0) { + fprintf(stderr, "ERROR! defining Earth-based observer location.\n"); + return 1; + } + + + // ------------------------------------------------------------------------- + // Set the astrometric time of observation... + + // Get the current system time, with up to nanosecond resolution... + clock_gettime(CLOCK_REALTIME, &unix_time); + + // Set the time of observation to the precise UTC-based UNIX time + // (We can set astromtric time using an other time measure also...) + if(novas_set_unix_time(unix_time.tv_sec, unix_time.tv_nsec, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; + } + + + // ------------------------------------------------------------------------- + // You might want to set a provider for precise planet positions so we might + // calculate Earth, Sun and major planet positions accurately. It is needed + // if you have orbitals defined around a major planet. + // + // There are many ways to set a provider of planet positions. For example, + // you may use CALCEPH: + // + // t_calcephbin *planets = calceph_open("path/to/de440s.bsp"); + // novas_use_calceph(planets); + + + // ------------------------------------------------------------------------- + // Initialize the observing frame with the given observing and Earth + // orientation patameters. + // + if(novas_make_frame(accuracy, &obs, &obs_time, POLAR_DX, POLAR_DY, &obs_frame) != 0) { + fprintf(stderr, "ERROR! failed to define observing frame.\n"); + return 1; + } + + // ------------------------------------------------------------------------- + // Calculate the precise apparent position (here in CIRS coordinates) + if(novas_sky_pos(&source, &obs_frame, NOVAS_CIRS, &apparent) != 0) { + fprintf(stderr, "ERROR! failed to calculate apparent position.\n"); + return 1; + } + + // Let's print the apparent position + printf(" RA = %.9f h, Dec = %.9f deg, rad_vel = %.6f km/s\n", apparent.ra, apparent.dec, apparent.rv); + + + // ------------------------------------------------------------------------- + // Convert the apparent position in CIRS on sky to horizontal coordinates + // We'll use a standard (fixed) atmospheric model to estimate an optical refraction + // (You might use other refraction models, or NULL to ignore refraction corrections) + if(novas_app_to_hor(&obs_frame, NOVAS_CIRS, apparent.ra, apparent.dec, novas_standard_refraction, &az, &el) != 0) { + fprintf(stderr, "ERROR! failed to calculate azimuth / elevation.\n"); + return 1; + } + + // Let's print the calculated azimuth and elevation + printf(" Az = %.6f deg, El = %.6f deg\n", az, el); + + return 0; +} + diff --git a/examples/example-star.c b/examples/example-star.c index 7c811827..7deab73f 100644 --- a/examples/example-star.c +++ b/examples/example-star.c @@ -1,249 +1,174 @@ -/* - Naval Observatory Vector Astrometry Software (NOVAS) - C Edition, Version 3.1 - - example.c: Examples of NOVAS calculations - - U. S. Naval Observatory - Astronomical Applications Dept. - Washington, DC - http://www.usno.navy.mil/USNO/astronomical-applications +/** + * @file + * + * @date Created on Jan 9, 2025 + * @author Attila Kovacs + * + * Example file for using the SuperNOVAS C/C++ library for determining positions for + * nearby (non-high-z) sidereal sources, such as a star. + * + * Link with + * + * ``` + * -lsupernovas + * ``` */ #include #include -#include +#include -#include "novas.h" -#include "eph_manager.h" // For solsys1.c +#include ///< SuperNOVAS functions and definitions -int main(void) { - /* - NOVAS 3.1 Example Calculations - - See Chapter 3 of User's Guide for explanation. - - Written for use with solsys version 1. - - To adapt for use with solsys version 2, see comments throughout file. - Assumes JPL ephemeris file "JPLEPH" located in same directory as - application. - */ - - const double latitude = 42.0; ///< [deg] geodetic latitude of observer - const double longitude = -70; ///< [deg] geodetic longitude of observer - const double height = 0.0; ///< [m] observer's altitude above sea level - const double temperature = 10.0; ///< [C] Ambient temperature at observing location (for refraction correction) - const double pressure = 1010.0; ///< [mbar] Atmospheric pressure at observing location (for refraction correction) - - const int leap_secs = 33; ///< [s] Leap seconds (UTC - TAI; see IERS Bulletins) - const double x_pole = -0.002; ///< [arcsec] Celestial pole offset in x (see IERS Bulletins) - const double y_pole = +0.529; ///< [arcsec] Celestial pole offset in y (see IERS Bulletins) - const double ut1_utc = -0.387845; ///< [s] UT1 - UTC time difference (see IERS Bulletins) - - // Various variables we will calculate / populate.... - double jd_beg, jd_end, jd_utc, jd_tt, jd_ut1, jd_tdb, delta_t, ra, dec, dis, rat, dect, dist, zd, az, rar, decr, gast, last, theta, jd[2], - pos[3], vel[3], pose[3], elon, elat, r, lon_rad, lat_rad, sin_lon, cos_lon, sin_lat, cos_lat, vter[3], vcel[3]; - - - observer obs_loc; - on_surface *geo_loc; - cat_entry star; - object moon, mars; - sky_pos t_place; - - int error = 0; - short de_num = 0; - - // Make structures of type 'on_surface' and 'observer-on-surface' containing - // the observer's position and weather (latitude, longitude, height, - //temperature, and atmospheric pressure) - make_observer_on_surface(latitude, longitude, height, temperature, pressure, &obs_loc); - geo_loc = &obs_loc.on_surf; - - // Make a structure of type 'cat_entry' containing the ICRS position - // and motion of star FK6 1307. - make_cat_entry("GMB 1830", "FK6", 1307, 11.88299133, 37.71867646, 4003.27, -5815.07, 109.21, -98.8, &star); - - if((error = make_object(NOVAS_PLANET, NOVAS_MOON, "Moon", NULL, &moon)) != 0) { - printf("Error %d from make_object (Moon)\n", error); - return (error); - } +// Below are some Earth orientation values. Here we define them as constants, but they may +// of course be variables. They should be set to the appropriate values for the time +// of observation based on the IERS Bulletins... - if((error = make_object(NOVAS_PLANET, NOVAS_MARS, "Mars", NULL, &mars)) != 0) { - printf("Error %d from make_object (Mars)\n", error); - return (error); - } - - // Open the JPL binary ephemeris file, here named "JPLEPH". - error = ephem_open("JPLEPH", &jd_beg, &jd_end, &de_num); - if(error != 0) { - if(error == 1) printf("JPL ephemeris file not found.\n"); - else printf("Error reading JPL ephemeris file header.\n"); - return (error); - } - else { - printf("JPL ephemeris DE%d open. Start JD = %10.2f End JD = %10.2f\n", de_num, jd_beg, jd_end); - printf("\n"); - } +#define LEAP_SECONDS 37 ///< [s] current leap seconds +#define DUT1 0.114 ///< [s] current UT1 - UTC time difference from IERS Bulletin A +#define POLAR_DX 230.0 ///< [mas] Earth polar offset x, e.g. from IERS Bulletin A. +#define POLAR_DY -62.0 ///< [mas] Earth polar offset y, e.g. from IERS Bulletin A. - // Set the planet calculator to the one using eph_manager - set_planet_provider(planet_eph_manager); - set_planet_provider_hp(planet_eph_manager_hp); +int main() { + // SuperNOVAS aariables used for the calculations -------------------------> + cat_entry star = CAT_ENTRY_INIT; // catalog information about a sidereal source + object source; // a celestial object: sidereal, planet, ephemeris or orbital source + observer obs; // observer location + novas_timespec obs_time; // astrometric time of observation + novas_frame obs_frame; // observing frame defined for observing time and location + enum novas_accuracy accuracy; // NOVAS_FULL_ACCURACY or NOVAS_REDUCED_ACCURACY + sky_pos apparent; // calculated precise observed (apparent) position of source - // Write banner. - printf("NOVAS Sample Calculations\n"); - printf("-------------------------\n"); - printf("\n"); + // Calculated quantities -------------------------------------------------> + double az, el; // calculated azimuth and elevation at observing site - // Write assumed longitude, latitude, height (ITRS = WGS-84). - printf("Geodetic location:\n"); - printf("%15.10f %15.10f %15.10f\n\n", geo_loc->longitude, geo_loc->latitude, geo_loc->height); + // Intermediate variables we'll use --------------------------------------> + struct timespec unix_time; // Standard precision UNIX time structure - // Establish time arguments. - jd_utc = julian_date(2008, 4, 24, 10.605); - jd_tt = jd_utc + ((double) leap_secs + 32.184) / 86400.0; - jd_ut1 = jd_utc + ut1_utc / 86400.0; - delta_t = get_ut1_to_tt(leap_secs, ut1_utc); - jd_tdb = jd_tt; /* Approximation good to 0.0017 seconds. */ + // We'll print debugging messages and error traces... + novas_debug(NOVAS_DEBUG_ON); - printf("TT and UT1 Julian Dates and Delta-T:\n"); - printf("%15.6f %15.6f %16.11f\n", jd_tt, jd_ut1, delta_t); - printf("\n"); - // Apparent and topocentric place of star FK6 1307 = GMB 1830. - error = app_star(jd_tt, &star, NOVAS_FULL_ACCURACY, &ra, &dec); - if(error != 0) { - printf("Error %d from app_star.\n", error); - return (error); - } + // ------------------------------------------------------------------------- + // Define a sidereal source - error = topo_star(jd_tt, delta_t, &star, &obs_loc.on_surf, NOVAS_FULL_ACCURACY, &rat, &dect); - if(error != 0) { - printf("Error %d from topo_star.\n", error); - return (error); + // Let's assume we have B1950 (FK4) coordinates... + // 16h26m20.1918s, -26d19m23.138s (B1950), proper motion -12.11, -23.30 mas/year, + // parallax 5.89 mas, radial velocity -3.4 km/s. + if(make_cat_entry("Antares", "FK4", 1, 16.43894213, -26.323094, -12.11, -23.30, 5.89, -3.4, &star) != 0) { + fprintf(stderr, "ERROR! defining cat_entry.\n"); + return 1; } - printf("FK6 1307 geocentric and topocentric positions:\n"); - printf("%15.10f %15.10f\n", ra, dec); - printf("%15.10f %15.10f\n", rat, dect); - printf("\n"); - - // Apparent and topocentric place of the Moon. - error = app_planet(jd_tt, &moon, NOVAS_FULL_ACCURACY, &ra, &dec, &dis); - if(error != 0) { - printf("Error %d from app_planet.", error); - return (error); + // First change the catalog coordinates (in place) to the J2000 (FK5) system... + if(transform_cat(CHANGE_EPOCH, NOVAS_JD_B1950, &star, NOVAS_JD_J2000, "FK5", &star) != 0) { + fprintf(stderr, "ERROR! converting B1950 catalog coordinates to J2000.\n"); + return 1; } - error = topo_planet(jd_tt, &moon, delta_t, &obs_loc.on_surf, NOVAS_FULL_ACCURACY, &rat, &dect, &dist); - if(error != 0) { - printf("Error %d from topo_planet.", error); - return (error); + // Then convert J2000 coordinates to ICRS (also in place). Here the dates don't matter... + if(transform_cat(CHANGE_J2000_TO_ICRS, 0.0, &star, 0.0, "ICRS", &star) != 0) { + fprintf(stderr, "ERROR! converting J2000 catalog coordinates to ICRS.\n"); + return 1; } - printf("Moon geocentric and topocentric positions:\n"); - printf("%15.10f %15.10f %15.12f\n", ra, dec, dis); - printf("%15.10f %15.10f %15.12f\n", rat, dect, dist); - // Topocentric (True of Date) place of the Moon using function 'place'-- should be - // same as above. - error = place(jd_tt, &moon, &obs_loc, delta_t, NOVAS_TOD, NOVAS_FULL_ACCURACY, &t_place); - if(error != 0) { - printf("Error %d from place.", error); - return (error); + // ------------------------------------------------------------------------- + // Wrap the sidereal souce into an object structure... + if(make_cat_object(&star, &source) != 0) { + fprintf(stderr, "ERROR! configuring observed object\n"); + return 1; } - printf("%15.10f %15.10f %15.12f\n", t_place.ra, t_place.dec, t_place.dis); - printf("\n"); - // Position of the Moon in local horizon coordinates. (Polar motion ..ignored here.) - equ2hor(jd_ut1, delta_t, NOVAS_FULL_ACCURACY, 0.0, 0.0, &obs_loc.on_surf, rat, dect, NOVAS_STANDARD_ATMOSPHERE, // - &zd, &az, &rar, &decr); + // ------------------------------------------------------------------------- + // Define observer somewhere on Earth (we can also define observers in Earth + // or Sun orbit, at the geocenter or at the Solary-system barycenter...) - printf("Moon zenith distance and azimuth:\n"); - printf("%15.10f %15.10f\n", zd, az); - printf("\n"); - - // Greenwich and local apparent sidereal time and Earth Rotation Angle. - error = sidereal_time(jd_ut1, 0.0, delta_t, NOVAS_TRUE_EQUINOX, EROT_ERA, NOVAS_FULL_ACCURACY, &gast); - if(error != 0) { - printf("Error %d from sidereal_time.", error); - return (error); + // Specify the location we are observing from + // 50.7374 deg N, 7.0982 deg E, 60m elevation + if(make_observer_on_surface(50.7374, 7.0982, 60.0, 0.0, 0.0, &obs) != 0) { + fprintf(stderr, "ERROR! defining Earth-based observer location.\n"); + return 1; } - last = remainder(gast + geo_loc->longitude / 15.0, 24.0); - if(last < 0.0) last += 24.0; - - theta = era(jd_ut1, 0.0); - printf("Greenwich and local sidereal time and Earth Rotation Angle:\n"); - printf("%16.11f %16.11f %15.10f\n", gast, last, theta); - printf("\n"); + // ------------------------------------------------------------------------- + // Set the astrometric time of observation... - // Heliocentric position of Mars in BCRS. + // Get the current system time, with up to nanosecond resolution... + clock_gettime(CLOCK_REALTIME, &unix_time); - // TDB ~ TT approximation could lead to error of ~50 m in position of Mars. - jd[0] = jd_tdb; - jd[1] = 0.0; - - error = ephemeris(jd, &mars, NOVAS_HELIOCENTER, NOVAS_FULL_ACCURACY, pos, vel); - if(error != 0) { - printf("Error %d from ephemeris (Mars).", error); - return (error); - } - - error = equ2ecl_vec(NOVAS_JD_J2000, NOVAS_MEAN_EQUATOR, NOVAS_FULL_ACCURACY, pos, pose); - if(error != 0) { - printf("Error %d from equ2ecl_vec.", error); - return (error); + // Set the time of observation to the precise UTC-based UNIX time + if(novas_set_unix_time(unix_time.tv_sec, unix_time.tv_nsec, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; } - if((error = vector2radec(pose, &elon, &elat)) != 0) { - printf("Error %d from vector2radec.", error); - return (error); + // ... Or you could set a time explicily in any known timescale. + /* + // Let's set a TDB-based time for the start of the J2000 epoch exactly... + if(novas_set_time(NOVAS_TDB, NOVAS_JD_J2000, LEAP_SECONDS, DUT1, &obs_time) != 0) { + fprintf(stderr, "ERROR! failed to set time of observation.\n"); + return 1; } - elon *= 15.0; - - r = sqrt(pose[0] * pose[0] + pose[1] * pose[1] + pose[2] * pose[2]); - - printf("Mars heliocentric ecliptic longitude and latitude and radius vector:\n"); - printf("%15.10f %15.10f %15.12f\n", elon, elat, r); - printf("\n"); - - // Terrestrial to celestial transformation. - lon_rad = geo_loc->longitude * DEG2RAD; - lat_rad = geo_loc->latitude * DEG2RAD; - sin_lon = sin(lon_rad); - cos_lon = cos(lon_rad); - sin_lat = sin(lat_rad); - cos_lat = cos(lat_rad); - - // Form vector toward local zenith (orthogonal to ellipsoid) in ITRS. - vter[0] = cos_lat * cos_lon; - vter[1] = cos_lat * sin_lon; - vter[2] = sin_lat; - - // Transform vector to GCRS. - error = ter2cel(jd_ut1, 0.0, delta_t, EROT_ERA, NOVAS_FULL_ACCURACY, NOVAS_REFERENCE_CLASS, x_pole, y_pole, vter, vcel); - if(error != 0) { - printf("Error %d from ter2cel.", error); - return (error); + */ + + + // ------------------------------------------------------------------------- + // You might want to set a provider for precise planet positions so we might + // calculate Earth, Sun and major planet positions accurately. If an planet + // provider is configured, we can unlock the ultimate (sub-uas) accuracy of + // SuperNOVAS. + // + // There are many ways to set a provider of planet positions. For example, + // you may use CALCEPH: + // + // t_calcephbin *planets = calceph_open("path/to/de440s.bsp"); + // novas_use_calceph(planets); + // + // accuracy = NOVAS_FULL_ACCURACY; // sub-uas precision + + // Without a planet provider, we are stuck with reduced (mas) precisions + // only... + accuracy = NOVAS_REDUCED_ACCURACY; // mas-level precision, typically + + + // ------------------------------------------------------------------------- + // Initialize the observing frame with the given observing and Earth + // orientation patameters. + // + if(novas_make_frame(accuracy, &obs, &obs_time, POLAR_DX, POLAR_DY, &obs_frame) != 0) { + fprintf(stderr, "ERROR! failed to define observing frame.\n"); + return 1; } - error = vector2radec(vcel, &ra, &dec); - if(error != 0) { - printf("Error %d from vector2radec.", error); - return (error); + + // ------------------------------------------------------------------------- + // Calculate the precise apparent position (here in CIRS coordinates) + if(novas_sky_pos(&source, &obs_frame, NOVAS_CIRS, &apparent) != 0) { + fprintf(stderr, "ERROR! failed to calculate apparent position.\n"); + return 1; } - printf("Direction of zenith vector (RA & Dec) in GCRS:\n"); - printf("%15.10f %15.10f\n", ra, dec); - printf("\n"); + // Let's print the apparent position + printf(" RA = %.9f h, Dec = %.9f deg, rad_vel = %.6f km/s\n", apparent.ra, apparent.dec, apparent.rv); - ephem_close(); /* remove this line for use with solsys version 2 */ - return (0); + // ------------------------------------------------------------------------- + // Convert the apparent position in CIRS on sky to horizontal coordinates + // We'll use a standard (fixed) atmospheric model to estimate an optical refraction + // (You might use other refraction models, or NULL to ignore refraction corrections) + if(novas_app_to_hor(&obs_frame, NOVAS_CIRS, apparent.ra, apparent.dec, novas_standard_refraction, &az, &el) != 0) { + fprintf(stderr, "ERROR! failed to calculate azimuth / elevation.\n"); + return 1; + } + + // Let's print the calculated azimuth and elevation + printf(" Az = %.6f deg, El = %.6f deg\n", az, el); + + return 0; } + diff --git a/examples/example-usno.txt b/examples/example-usno.txt deleted file mode 100644 index 682aded3..00000000 --- a/examples/example-usno.txt +++ /dev/null @@ -1,32 +0,0 @@ -JPL ephemeris DE405 open. Start JD = 2305424.50 End JD = 2525008.50 - -NOVAS Sample Calculations -------------------------- - -Geodetic location: - -70.0000000000 42.0000000000 0.0000000000 - -TT and UT1 Julian Dates and Delta-T: - 2454580.942629 2454580.941871 65.57184500000 - -FK6 1307 geocentric and topocentric positions: - 11.8915509892 37.6586357955 - 11.8915479153 37.6586695456 - -Moon geocentric and topocentric positions: - 17.1390774264 -27.5374448869 0.002710296515 - 17.1031967646 -28.2902502967 0.002703785126 - 17.1031967646 -28.2902502967 0.002703785126 - -Moon zenith distance and azimuth: - 81.6891016502 219.2708903405 - -Greenwich and local sidereal time and Earth Rotation Angle: - 0.79362134148 20.12695467481 11.7956158462 - -Mars heliocentric ecliptic longitude and latitude and radius vector: - 148.0032235906 1.8288284075 1.664218258879 - -Direction of zenith vector (RA & Dec) in GCRS: - 20.1221838608 41.9769823554 -