diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index d5e3222f..af131680 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -33,8 +33,9 @@ jobs: - name: Generate coverage data run: make coverage - - name: Upload coverage to Codecov + - name: Upload coverage to Codecov.io uses: codecov/codecov-action@v4 + continue-on-error: true with: fail_ci_if_error: false @@ -42,4 +43,8 @@ jobs: name: codecov token: ${{ secrets.CODECOV_TOKEN }} verbose: true - + + - name: Upload coverage to Coveralls.io + continue-on-error: true + uses: coverallsapp/github-action@v2 + diff --git a/CHANGELOG.md b/CHANGELOG.md index 495dd3a8..380ed7d6 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,15 +2,15 @@ ## [Unreleased] -Changes expected for the next release... +Changes planned/expected for the next release... ### Added - - Add debug mode, which can be activated via `novas_debug(1)`. When debug mode is enabled, all errors are printed to - the console with a trace (the trace does not contain line numbers but provides information on the functions and + - Add debug mode, which can be activated via `novas_debug(1)`. When debug mode is enabled, all errors will be printed + to the console with a trace (the trace does not contain line numbers but provides information on the functions and locations within them where the error occurred. - + ## [1.0.0] - 2024-03-01 @@ -18,38 +18,44 @@ This is the initial release of the SuperNOVAS library. ### Fixed - - Fixed the [sidereal_time bug](https://aa.usno.navy.mil/software/novas_faq), whereby the `sidereal_time()` + - Fixes the [sidereal_time bug](https://aa.usno.navy.mil/software/novas_faq), whereby the `sidereal_time()` function had an incorrect unit cast. This is a known issue of NOVAS C 3.1. - - Fixed the [ephem_close bug](https://aa.usno.navy.mil/software/novas_faq), whereby `ephem_close()` in + - Fixes the [ephem_close bug](https://aa.usno.navy.mil/software/novas_faq), whereby `ephem_close()` in `ephem_manager.c` did not reset the `EPHFILE` pointer to NULL. This is a known issue of NOVAS C 3.1. - - Fixed antedating velocities and distances for light travel time in `ephemeris()`. When getting positions and + - Fixes antedating velocities and distances for light travel time in `ephemeris()`. When getting positions and velocities for Solar-system sources, it is important to use the values from the time light originated from the observed body rather than at the time that light arrives to the observer. This correction was done properly for positions, but not for velocities or distances, resulting in incorrect observed radial velocities or apparent distances being reported for spectroscopic observations or for angular-physical size conversions. - - Fixed bug in `ira_equinox()` which may return the result for the wrong type of equinox (mean vs. true) if the the + - Fixes bug in `ira_equinox()` which may return the result for the wrong type of equinox (mean vs. true) if the the `equinox` argument was changing from 1 to 0, and back to 1 again with the date being held the same. This affected routines downstream also, such as `sidereal_time()`. - - Fixed accuracy switching bug in `cio_basis()`, `cio_location()`, `ecl2equ`, `equ2ecl_vec()`, `ecl2equ_vec()`, + - Fixes accuracy switching bug in `cio_basis()`, `cio_location()`, `ecl2equ`, `equ2ecl_vec()`, `ecl2equ_vec()`, `geo_posvel()`, `place()`, and `sidereal_time()`. All these functions returned a cached value for the other accuracy if the other input parameters are the same as a prior call, except the accuracy. - - Fixed multiple bugs in using cached values in `cio_basis()` with alternating CIO location reference systems. + - Fixes multiple bugs related to using cached values in `cio_basis()` with alternating CIO location reference + systems. - - Fixed bug in `equ2ecl_vec()` and `ecl2equ_vec()` whereby a query with `coord_sys = 2` (GCRS) has overwritten the + - Fixes bug in `equ2ecl_vec()` and `ecl2equ_vec()` whereby a query with `coord_sys = 2` (GCRS) has overwritten the cached mean obliquity value for `coord_sys = 0` (mean equinox of date). As a result, a subsequent call with `coord_sys = 0` and the same date as before would return the results GCRS coordinates instead of the requested mean equinox of date coordinates. - The use of `fmod()` in NOVAS C 3.1 led to the wrong results when the numerator was negative and the result was not turned into a proper remainder. This affected the calculation of the mean anomaly in `solsys3.c` (line 261) - and the fundamental arguments calculted in `fund_args()` and `ee_ct()` for dates prior to J2000. Less + and the fundamental arguments calculated in `fund_args()` and `ee_ct()` for dates prior to J2000. Less critically, it also was the reason `cal_date()` did not work for negative JD values. + - Fixes `aberration()` returning NAN vectors if the `ve` argument is 0. It now returns the unmodified input + vector appropriately instead. + + - Fixes potential string overflows and eliminates associated compiler warnings. + ### Added @@ -78,7 +84,7 @@ This is the initial release of the SuperNOVAS library. approximation via `set_nutation_lp_provider()`. For example, the user may want to use the `iau2000b()` model instead or some custom algorithm instead. - - New intutitive XYZ coordinate coversion functions: + - New intutitive XYZ coordinate conversion functions: * for GCRS - CIRS - ITRS (IAU 2000 standard): `gcrs_to_cirs()`, `cirs_to_itrs()`, and `itrs_to_cirs()`, `cirs_to_gcrs()`. * for GCRS - J2000 - TOD - ITRS (old methodology): `gcrs_to_j2000()`, `j2000_to_tod()`, `tod_to_itrs()`, and @@ -89,6 +95,8 @@ This is the initial release of the SuperNOVAS library. version of the existing `equ2hor()` for converting from TOD to to local horizontal (old methodology), the `cirs_to_itrs()` followed by `itrs_to_hor()` does the same from CIRS (new IAU standard methodology), and had no equivalent in NOVAS C 3.1. + + - New `ecl2equ()` for converting ecliptic coordinates to equatorial, complementing existing `equ2ecl()`. - New `gal2equ()` for converting galactic coordinates to ICRS equatorial, complementing existing `equ2gal()`. @@ -98,7 +106,7 @@ This is the initial release of the SuperNOVAS library. - New convenience functions to wrap `place()` for simpler specific use: `place_star()`, `place_icrs()`, `place_gcrs()`, `place_cirs()`, and `place_tod()`. - - New `radec_star()` and `radec_planet()` methods as the common point for all existing methods sych as `astro_star()` + - New `radec_star()` and `radec_planet()` methods as the common point for all existing methods such as `astro_star()` `local_star()`, `topo_planet()` etc. - New time conversion utilities `tt2tdb()` and `get_ut1_to_tt()` make it simpler to convert between UT1, TT, and TDB @@ -107,8 +115,9 @@ This is the initial release of the SuperNOVAS library. - Co-existing `solarsystem()` variants. It is possible to use the different `solarsystem()` implementations provided by `solsys1.c`, `solsys2.c`, `solsys3.c` and/or `solsys-ephem.c` side-by-side, as they define their functionalities with distinct, non-conflicting names, e.g. `earth_sun_calc()` vs `planet_jplint()` vs - `planet_ephem_manager()` vs `planet_ephem_provider()`. See the section on [Building and installation](#installation) - further above on including a selection of these in your library build.) + `planet_ephem_manager()` vs `planet_ephem_provider()`. See the section on + [Building and installation](#installation) further above on including a selection of these in your library + build.) - New `novas_case_sensitive(int)` method to enable (or disable) case-sensitive processing of object names. (By default NOVAS object names were converted to upper-case, making them effectively case-insensitive.) @@ -131,32 +140,50 @@ This is the initial release of the SuperNOVAS library. - All erroneous returns now set `errno` so that users can track the source of the error in the standard C way and use functions such as `perror()` ans `strerror()` to print human-readable error messages. - - Output values supplied via pointers are set to clearly invalid values in case of erroneous returns, such as - `NAN` so that even if the caller forgets to check the error code, it becomes obvious that the values returned are - should not be used as if they were valid + - Many output values supplied via pointers are set to clearly invalid values in case of erroneous returns, such as + `NAN` so that even if the caller forgets to check the error code, it becomes obvious that the values returned + should not be used as if they were valid. (No more sneaky silent errors.) - Many SuperNOVAS functions allow `NULL` arguments, both for optional input values as well as outputs that are not required. See the [API Documentation](https://smithsonian.github.io/SuperNOVAS.home/apidoc/html/) for specifics). - This eliminates the need to declare dummy variables in your application code. + This eliminates the need to declare dummy variables in your application code for quantities you do not require. - All SuperNOVAS functions that take an input vector to produce an output vector allow the output vector argument be the same as the input vector argument. For example, `frame_time(pos, J2000_TO_ICRS, pos)` using the same `pos` vector both as the input and the output. In this case the `pos` vector is modified in place by the call. - This can greatly simplify usage, and eliminate extraneous declarations, when intermediates are not required. + This can greatly simplify usage, and can eliminate extraneous declarations, when intermediates are not required. - - SuperNOVAS prototypes declare function pointer arguments as `const` whenever the function does not modify the - data content being pointed at. This supports better programming practices that generally aim to avoid unintended - data modifications. - - - Source names and catalog names can both be up to 64 bytes (including termination), up from 51 and 4 respectively - NOVAS C, while keeping `struct` layouts the same thanks to alignment. + - SuperNOVAS declares function pointer arguments as `const` whenever the function does not modify the data content + being referenced. This supports better programming practices that generally aim to avoid unintended data + modifications. + + - Catalog names can be up to 6 bytes (including termination), up from 4 in NOVAS C, while keeping `struct` layouts + the same as NOVAS C thanks to alignment, thus allowing cross-compatible binary exchage of `cat_entry` records + with NOVAS C 3.1. - Object ID numbers are `long` instead of `short` to accommodate NAIF IDs, which require minimum 32-bit integers. + + - `cel2ter()` and `tel2cel()` can now process 'option'/'class' = 1 (`NOVAS_REFERENCE_CLASS`) regardless of the + methodology (`EROT_ERA` or `EROT_GST`) used to input or output coordinates in GCRS. + + - `make_object()` ignored the specified number argument for sidereal sources (set to 0), but we set it to the + specified value assuming the caller provided it for a reason. (It does not change the separate `starnumber` value + that is included in the `star` argument however) + + - `sun_eph()` in `solsysl3.c` evaluates the series in reverse order compared to NOVAS C 3.1, accumulating the least + significant terms first, and thus resulting in higher precision result in the end. + + - Changed `vector2radec()` to return NAN values if the input is a null-vector (i.e. all components are zero). + + - Changed the standard atmospheric model for (optical) refraction calculation to include a simple model for the + annual average temperature at the site (based on latitude and elevation). This results is a slightly more educated + guess of the actual refraction than the global fixed temperature of 10 °C assumed by NOVAC C 3.1 regardless of + observing location. ### Deprecated - - `novascon.h` / `novascon.c`: These definitions of constants was troublesome fow two reasons: (1) They were + - `novascon.h` / `novascon.c`: These definitions of constants was troublesome for two reasons: (1) They were primarily meant for use internally within the library itself. As the library clearly defines in what units input and output quantities are expressed, the user code can apply its own appropriate conversions that need not match the internal system used by the library. Hence exposing these constants to users was half baked. (2) The naming of @@ -167,21 +194,21 @@ This is the initial release of the SuperNOVAS library. - `equ2hor()`: It's name does not make it clear that this function is suitable only for converting TOD (old methodology) to horizontal but not CIRS to horizontal (IAU 2000 standard). You should use the equivalent but more - specific `tod_to_itrs()`, or else the newly added `cirs_to_itrs()`, followed by `itrs_to_hor()` instead. + specific `tod_to_itrs()` or the newly added `cirs_to_itrs()`, followed by `itrs_to_hor()` instead. - `cel2ter()` / `ter2cel()`: These function can be somewhat confusing to use. You are likely better off with - `tod_to_itrs()` and `cirs_to_itrs()` instead. + `tod_to_itrs()` and `cirs_to_itrs()` instead, and possibly followed by further conversions if desired. - `app_star()`, `app_planet()`, `topo_star()` and `topo_planet()`: These use the old (pre IAU 2000) methodology, which isn't clear from their naming. Use `place()` or `place_star()` with `NOVAS_TOD` or `NOVAS_CIRS` as the system instead, as appropriate. - - `readeph()`: prone to memtoy leaks, and not flexible with its origin (e.g. barycenter vs heliocenter). Instead, use - a similar `novas_ephem_provider` implementation and `set_ephem_provider()` instead for a more flexible and less - troublesome equivalent, which also does not need to be baked into the library but can be configured at runtime. + - `readeph()`: prone to memory leaks, and not flexible with its origin (necessarily at the barycenter). Instead, use + a similar `novas_ephem_provider` implementation and `set_ephem_provider()` for a more flexible and less + troublesome equivalent, which also does not need to be baked into the library and can be configured at runtime. - - `tdb2tt()`. Use `tt2tdb()` instead. It's both more intuitive (returning the time difference as a double) and faster - to calculate, not to mention that it implements the more standard approach. + - `tdb2tt()`. Use `tt2tdb()` instead. It's both more intuitive to use (returning the time difference as a double) and + faster to calculate, not to mention that it implements the more standard approach. diff --git a/CONTRIBUTING.md b/CONTRIBUTING.md index 538d3267..a41aabf3 100644 --- a/CONTRIBUTING.md +++ b/CONTRIBUTING.md @@ -31,14 +31,14 @@ upstream every once in a while, and merging them into your development branch. D tests to `test/src/test.c` as appropriate, and if you add new source files, you may need to adjust `test/Makefile` so that `gcov` coverage is generated for your new source files. You should aim for 100% diff coverage. When pushing changes to your fork, you can get a coverage report by checking - the Github Actions result of your commit (click the Codecov link), and you can analyze what line(s) of code need to - have tests added. Try to create tests that are simple but meaningful (i.e. check for valid results, rather than just - confirm existing behavior), and try to cover as many realistic scenarios as appropriate. Write lots of tests if you - need to. + the Github Actions result of your commit (click the Coveralls link), and you can analyze what line(s) of code need + to have tests added. Try to create tests that are simple but meaningful (i.e. check for valid results, rather than + just confirm existing behavior), and try to cover as many realistic scenarios as appropriate. Write lots of tests + if you need to. 4. __Pull Request__. Once you feel your work can be integrated, create a pull request from your fork/branch. You can do that easily from the github page of your fork/branch directly. In the pull request, provide a concise description -of what you added or changed. Your pull request will be reviwed. You may get some feedback at this point, and maybe +of what you added or changed. Your pull request will be reviewed. You may get some feedback at this point, and maybe there will be discussions about possible improvements or regressions etc. It's a good thing too, and your changes will likely end up with added polish as a result. You can be all the more proud of it in the end! diff --git a/README.md b/README.md index 3fb9d234..fb07299e 100644 --- a/README.md +++ b/README.md @@ -2,6 +2,9 @@ ![Test](https://github.com/Smithsonian/SuperNOVAS/actions/workflows/test.yml/badge.svg) ![Static Analysis](https://github.com/Smithsonian/SuperNOVAS/actions/workflows/check.yml/badge.svg) ![API documentation](https://github.com/Smithsonian/SuperNOVAS/actions/workflows/dox.yml/badge.svg) + + ![Coverage Status](https://coveralls.io/repos/github/Smithsonian/SuperNOVAS/badge.svg?branch=main) + @@ -11,17 +14,18 @@
# SuperNOVAS -The Naval Observatory NOVAS C astrometry library, made better. +The NOVAS C astrometry library, made better. [SuperNOVAS](https://github.com/Smithsonian/SuperNOVAS/) is a positional astronomy library for the for the C programming language, providing high-precision astrometry such as one might need for running an observatory or a -precise planetarium program. Its source code is compatible with the C90 standard, and hence should be suitable for -many older platforms also. It is light-weight and easy to use, with full support for the IAU 2000/2006 standards -for sub-microarcsecond position calculations. +precise planetarium program. It is a fork of the Naval Observatory Vector Astrometry Software +([NOVAS](https://aa.usno.navy.mil/software/novas_info)), with the aim of making it more user-friendly and +easier to use overall. + +SuperNOVAS is entirely free to use without any licensing restrictions. Its source code is compatible with the C90 +standard, and hence should be suitable for many older platforms also. It is light-weight and easy to use, with full +support for the IAU 2000/2006 standards for sub-microarcsecond position calculations. -SuperNOVAS is a fork of the Naval Observatory Vector Astrometry Software -([NOVAS](https://aa.usno.navy.mil/software/novas_info)), with the overall aim of making it more user-friendly and -easier to use. It is entirely free to use without any licensing restrictions. # Table of Contents @@ -32,9 +36,9 @@ easier to use. It is entirely free to use without any licensing restrictions. - [Fixed NOVAS C 3.1 issues](#fixed-issues) - [Building and installation](#installation) - [Building your application with SuperNOVAS](#building-your-application) - - [Basic usage examples](#examples) - - [SuperNOVAS specific features](#supernovas-features) + - [Example usage](#examples) - [Notes on precision](#precision) + - [SuperNOVAS specific features](#supernovas-features) - [External Solar-system ephemeris data or services](#solarsystem) - [Release schedule](#release-schedule) @@ -55,7 +59,7 @@ The primary goals of SuperNOVAS is to improve on the stock NOVAS C library by: - Improving the ease of use by using `enum`s instead of integer constants, which also allows for some checking of use during compilations (such as using the incorrect `enum` type). - Improving [API documentation](https://smithsonian.github.io/SuperNOVAS.home/apidoc/html/) with - [Doxygen](https://www.doxygen.nl/) to provide browseable cross-referenced API docs. + [Doxygen](https://www.doxygen.nl/) to provide browsable cross-referenced API docs. - Streamlining calculations where possible to improve performance. - Adding `const` modifier to prototype arguments where appropriate - Checking arguments and setting `errno` as appropriate (and returning -1 unless another appropriate error code was @@ -73,8 +77,8 @@ may have written for NOVAS C. SuperNOVAS is currently based on NOVAS C version 3.1. We plan to rebase SuperNOVAS to the latest upstream release of the NOVAS C library, if new releases become available. -SuperNOVAS is maintained by Attila Kovacs at the Center for Astrophysics \| Harvard and Smithsonian, and it is -available through the [Smithsonian/SuperNOVAS](https://github.com/Smithsonian/SuperNOVAS) repo on GitHub. +SuperNOVAS is maintained by Attila Kovacs at the Center for Astrophysics \| Harvard & Smithsonian, and it is +available through the [Smithsonian/SuperNOVAS](https://github.com/Smithsonian/SuperNOVAS) repository on GitHub. Outside contributions are very welcome. See [how you can contribute](https://github.com/Smithsonian/SuperNOVAS/CONTRIBUTING.md) to make SuperNOVAS even better. @@ -89,7 +93,8 @@ Here are some links to SuperNOVAS related content online: - [API Documentation](https://smithsonian.github.io/SuperNOVAS.home/apidoc/html/) - [Project site](https://github.com/Smithsonian/SuperNOVAS/) on GitHUB. - - [SuperNOVAS home page](https://smithsonian.github.io/SuperNOVAS.home) page on github.io. + - [SuperNOVAS home page](https://smithsonian.github.io/SuperNOVAS.home) page on github.io. + - [Regression Test Coverage](https://coveralls.io/github/Smithsonian/SuperNOVAS?branch=main) at Coveralls.io. - [How to Contribute](https://github.com/Smithsonian/SuperNOVAS/CONTRIBUTING.md) guide - [NOVAS](https://aa.usno.navy.mil/software/novas_info) home page at the US Naval Observatory. - [SPICE toolkit](https://naif.jpl.nasa.gov/naif/toolkit.html) for integrating Solar-system ephemeris @@ -103,15 +108,15 @@ https://github.com/attipaci/attipaci.gitgub.io ## Compatibility with NOVAS C 3.1 SuperNOVAS strives to maintain API compatibility with the upstream NOVAS C 3.1 library, but not binary -compatilibility. In practical terms it means that you cannot simply drop-in replace your compiled objects (e.g. +compatibility. In practical terms it means that you cannot simply drop-in replace your compiled objects (e.g. `novas.o`), or the static (e.g. `novas.a`) or shared (e.g. `novas.so`) libraries, from NOVAS C 3.1 with that from SuperNOVAS. Instead, you will need to (re)compile and or (re)link your application with the SuperNOVAS versions of these. This is because some function signatures have changed, e.g. to use an `enum` argument instead of the nondescript -`short int` argument of NOVAS C 3.1, or because we added a return value to a functions that was declared `void` -in NOVAS C 3.1. We also changed the `object` structure to contain a `long` ID number instead of `short` to -accommodate JPL NAIF values, which require a 32-bit width. +`short int` argument of NOVAS C 3.1, or because we added a return value to a function that was declared `void` +in NOVAS C 3.1. We also changed the `object` structure to contain a `long` ID number instead of `short` to +accommodate JPL NAIF codes, and for which 16-bit storage is insufficient. ----------------------------------------------------------------------------- @@ -141,7 +146,8 @@ provided by SuperNOVAS over the upstream NOVAS C 3.1 code: `geo_posvel()`, `place()`, and `sidereal_time()`. All these functions returned a cached value for the other accuracy if the other input parameters are the same as a prior call, except the accuracy. - - Fix multiple bugs in using cached values in `cio_basis()` with alternating CIO location reference systems. + - Fixes multiple bugs related to using cached values in `cio_basis()` with alternating CIO location reference + systems. - Fix bug in `equ2ecl_vec()` and `ecl2equ_vec()` whereby a query with `coord_sys = 2` (GCRS) has overwritten the cached mean obliquity value for `coord_sys = 0` (mean equinox of date). As a result, a subsequent call with @@ -150,9 +156,13 @@ provided by SuperNOVAS over the upstream NOVAS C 3.1 code: - The use of `fmod()` in NOVAS C 3.1 led to the wrong results when the numerator was negative and the result was not turned into a proper remainder. This affected the calculation of the mean anomaly in `solsys3.c` (line 261) - and the fundamental arguments calculted in `fund_args()` and `ee_ct()` for dates prior to J2000. Less + and the fundamental arguments calculated in `fund_args()` and `ee_ct()` for dates prior to J2000. Less critically, it also was the reason `cal_date()` did not work for negative JD values. + - Fixes `aberration()` returning NAN vectors if the `ve` argument is 0. It now returns the unmodified input + vector appropriately instead. + + - Fixes potential string overflows and eliminates associated compiler warnings. ----------------------------------------------------------------------------- @@ -160,7 +170,7 @@ provided by SuperNOVAS over the upstream NOVAS C 3.1 code: ## Building and installation -The SuperNOVAS distibution contains a `Makefile` for GNU make, which is suitable for compiling the library (as well as +The SuperNOVAS distribution contains a `Makefile` for GNU make, which is suitable for compiling the library (as well as local documentation, and tests, etc.) on POSIX systems such as Linux, BSD, MacOS X, or Cygwin or WSL on Windows. Before compiling the library take a look a `config.mk` and edit it as necessary for your needs, such as: @@ -241,14 +251,35 @@ being unset in `config.mk`). -## Basic usage examples +## Example usage + - [Note on alternative methodologies](#methodologies) - [Calculating positions for a sidereal source](#sidereal-example) - [Calculating positions for a Solar-system source](#solsys-example) - [Reduced accuracy shortcuts](#accuracy-notes) - [Performance considerations](#performance-note) - - [Note on alternative methodologies](#methodologies) + + + +### Note on alternative methodologies + +The IAU 2000 and 2006 resolutions have completely overhauled the system of astronomical coordinate transformations +to enable higher precision astrometry. (Super)NOVAS supports coordinate calculations both in the old (pre IAU 2000) +ways, and in the new IAU standard method. Here is an overview of how the old and new methods define some of the +terms differently: + + | Concept | Old standard | New IAU standard | + | -------------------------- | ----------------------------- | ------------------------------------------------- | + | Catalog coordinate system | J2000 or B1950 | International Celestial Reference System (ICRS) | + | Dynamical system | True of Date (TOD) | Celestial Intermediate Reference System (CIRS) | + | Dynamical R.A. origin | true equinox of date | Celestial Intermediate Origin (CIO) | + | Precession, nutation, bias | separate, no tidal terms | IAU 2006 precession/nutation model | + | Celestial Pole offsets | dψ, dε | _dx_, _dy_ | + | Earth rotation measure | Greenwich Sidereal Time (GST) | Earth Rotation Angle (ERA) | + | Fixed Earth System | WGS84 | International Terrestrial Reference System (ITRS) | +See the various enums and constants defined in `novas.h`, as well as the descriptions on the various NOVAS routines +on how they are appropriate for the old and new methodologies respectively. ### Calculating positions for a sidereal source @@ -292,8 +323,10 @@ Next, we define the location where we observe from. Here we can (but don't have make_observer_on_surface(50.7374, 7.0982, 60.0, 0.0, 0.0, &obs); ``` -We also need to set the time of observation. Out clocks usually measure UTC, but for NOVAS we usually need time -measured based on Terrestrial Time (TT) or Barycentric Time (TDB) or UT1. So +We also need to set the time of observation. Our clocks usually measure UTC, but for NOVAS we usually need time +measured based on Terrestrial Time (TT) or Barycentric Time (TDB) or UT1. Typically you will have to provide +NOVAS with the TT - UT1 time difference, which can be calculated from the current leap seconds and the UT1 - UTC +time difference (a.k.a. DUT1): ```c // The current value for the leap seconds (UTC - TAI) @@ -332,7 +365,7 @@ coordinate system and for the specified observing location). We also get radial distance (e.g. for apparent-to-physical size conversion): ```c - sky_pos pos; // We'll return the observable positions in this structure + sky_pos pos; // We'll return the observable positions (in CIRS) in this structure // Calculate the apparent (CIRS) topocentric positions for the above configuration int status = place_star(jd_tt, &source, &obs, ut1_to_tt, NOVAS_CIRS, NOVAS_FULL_ACCURACY, &pos); @@ -344,15 +377,15 @@ distance (e.g. for apparent-to-physical size conversion): } ``` -Finally, we may want to calculate the astrometric zenith distance (= 90° - azimuth) and elevation angles of the source -at the specified observing location (without refraction correction): +Finally, we may want to calculate the astrometric azimuth and zenith distance (= 90° - azimuth) angles of the +source at the specified observing location (without refraction correction): ```c double itrs[3]; // ITRS position vector of source to populate double az, zd; // [deg] local azimuth and zenith distance angles to populate // Convert CIRS to Earth-fixed ITRS using the pole offsets. - cirs_to_itrs(jd_tt, ut1_to_tt, NOVAS_FULL_ACCURACY, dx, dy, itrs); + cirs_to_itrs(jd_tt, 0.0, ut1_to_tt, NOVAS_FULL_ACCURACY, dx, dy, pos.r_hat, itrs); // Finally convert ITRS to local horizontal coordinates at the observing site itrs_to_hor(itrs, &obs.on_surface, &az, &zd); @@ -378,16 +411,16 @@ want, e.g.: Solar-system sources work similarly to the above with a few important differences. First, You will have to provide one or more functions to obtain the barycentric ICRS positions for your Solar-system -source of interest for the specific Barycentric Dynamical Time (TDB) of observation. See Section on integrating +source(s) of interest for the specific Barycentric Dynamical Time (TDB) of observation. See section on integrating [External Solar-system ephemeris data or services](#solarsystem) with SuperNOVAS. You can specify the functions that will handle the respective ephemeris data at runtime before making the NOVAS calls that need them, e.g.: ```c // Set the function to use for regular precision planet position calculations - set_planet_provider(my_planet_calcular_function); + set_planet_provider(my_planet_function); // Set the function for high-precision planet position calculations - set_planet_provider(my_very_precise_planet_calculator_function); + set_planet_provider(my_very_precise_planet_function); // Set the function to use for calculating all sorts of solar-system bodies set_ephem_provider(my_ephemeris_provider_function); @@ -410,11 +443,11 @@ more generic ephemeris handling via a user-provided `novas_ephem_provider`. E.g. object mars, ceres; // Hold data on solar-system bodies. // Mars will be handled by hte planet calculator function - make_object(NOVAS_PLANET_MARS, NOVAS_MARS, "Mars", NULL, &mars); + make_planet(NOVAS_MARS, &mars); // Ceres will be handled by the generic ephemeris reader, which say uses the // NAIF ID of 2000001 - make_object(NOVAS_EPHEM_OBJECT, 2000001, "Ceres", NULL, &ceres); + make_ephem_object("Ceres", 2000001, &ceres); ``` Other than that, it's the same spiel as before, except using the appropriate `place()` for generic celestial @@ -444,36 +477,57 @@ When one does not need positions at the microarcsecond level, some shortcuts can #### Performance considerations Some of the calculations involved can be expensive from a computational perspective. For the most typical use case -however, NOVAS (and SuperNOVAS) has a trick up its sleve: it caches the last result of intensive calculations so they +however, NOVAS (and SuperNOVAS) has a trick up its sleeve: it caches the last result of intensive calculations so they may be re-used if the call is made with the same environmental parameters again (such as JD time and accuracy). Therefore, when calculating positions for a large number of sources at different times: - It is best to iterate over the sources while keeping the time fixed in the inner loop. - - You probably want to stick to one accuracy morde (`NOVAS_FULL_ACCURACY` or `NOVAS_REDUCED_ACCURACY`) to prevent + - You probably want to stick to one accuracy mode (`NOVAS_FULL_ACCURACY` or `NOVAS_REDUCED_ACCURACY`) to prevent re-calculating the same quantities repeatedly to alternating precision. - If super-high accuracy is not required `NOVAS_REDUCED_ACCURACY` mode offers much faster calculations, in general. - - -### Note on alternative methodologies +----------------------------------------------------------------------------- -The IAU 2000 and 2006 resolutions have completely overhauled the system of astronomical coordinate transformations -to enable higher precision astrometry. (Super)NOVAS supports coordinate calculations both in the old (pre IAU 2000) -ways, and in the new IAU standard method. + +## Notes on precision - | Concept | Old standard | New IAU standard | - | -------------------------- | ----------------------------- | ------------------------------------------------- | - | Catalog coordinate system | J2000 or B1950 | International Celestial Reference System (ICRS) | - | Dynamical system | True of Date (TOD) | Celestial Intermediate Reference System (CIRS) | - | Dynamical R.A. origin | true equinox of date | Celestial Intermediate Origin (CIO) | - | Precession, nutation, bias | separate, no tidal terms | IAU 2006 precession/nutation model | - | Celestial Pole offsets | dψ, dε | _dx_, _dy_ | - | Earth rotation measure | Greenwich Sidereal Time (GST) | Earth Rotation Angle (ERA) | - | Fixed Earth System | WGS84 | International Terrestrial Reference System (ITRS) | +The SuperNOVAS library is in principle capable of calculating positions to sub-microarcsecond, and velocities to mm/s +precision for all types of celestial sources. However, there are certain prerequisites and practical considerations +before that level of accuracy is reached. + + + 1. __IAU 2000/2006 conventions__: High precision calculations will generally require that you use SuperNOVAS with the + new IAU standard quantities and methods. The old ways were simply not suited for precision much below the + milliarcsecond level. + + 2. __Earth's polar motion__: Calculating precise positions for any Earth-based observations requires precise + knowledge of Earth orientation at the time of observation. The pole is subject to predictable precession and + nutation, but also small irregular variations in the orientation of the rotational axis and the rotation period + (a.k.a polar wobble). The [IERS Bulletins](https://www.iers.org/IERS/EN/Publications/Bulletins/bulletins.html) + provide up-to-date measurements, historical data, and near-term projections for the polar offsets and the UT1-UTC + (DUT1) time difference and leap-seconds (UTC-TAI). In SuperNOVAS you can use `cel_pole()` and `get_ut1_to_tt()` + functions to apply / use the published values from these to improve the astrometic precision of Earth-orientation + based coordinate calculations. Without setting and using the actual polar offset values for the time of + observation, positions for Earth-based observations will be accurate at the arcsecond level only. -See the various enums and constands defined in `novas.h`, as well as the descriptions on the various NOVAS routines -on how they are appropriate for the old and new methodologies respectively. + 3. __Solar-system sources__: Precise calculations for Solar-system sources requires precise ephemeris data for both + the target object as well as for Earth, and the Sun vs the Solar-system barycenter. For the highest precision + calculations you also need positions for all major planets to calculate gravitational deflection precisely. By + default SuperNOVAS can only provide approximate positions for the Earth and Sun (see `earth_sun_calc()` in + `solsys3.c`), but certainly not at the sub-microarcsecond level, and not for other solar-system sources. You will + need to provide a way to interface SuperNOVAS with a suitable ephemeris source (such as the CSPICE toolkit from + JPL) if you want to use it to obtain precise positions for Solar-system bodies. See the + [section below](#solarsystem) for more information how you can do that. + + 4. __Refraction__: Ground based observations are also subject to atmospheric refraction. SuperNOVAS offers the + option to include _optical_ refraction corrections in `equ2hor()` either for a standard atmosphere or more + precisely using the weather parameters defined in the `on_surface` data structure that specifies the observer + locations. Note, that refraction at radio wavelengths is notably different from the included optical model. In + either case you may want to skip the refraction corrections offered in this library, and instead implement your + own as appropriate (or not at all). + + ----------------------------------------------------------------------------- @@ -490,6 +544,13 @@ on how they are appropriate for the old and new methodologies respectively. duplicate arguments) and will return -1 (with `errno` set, usually to `EINVAL`) if the arguments supplied are invalid (unless the NOVAS C API already defined a different return value for specific cases. If so, the NOVAS C error code is returned for compatibility). + + - All erroneous returns now set `errno` so that users can track the source of the error in the standard C way and + use functions such as `perror()` ans `strerror()` to print human-readable error messages. + + - Many output values supplied via pointers are set to clearly invalid values in case of erroneous returns, such as + `NAN` so that even if the caller forgets to check the error code, it becomes obvious that the values returned + should not be used as if they were valid. (No more sneaky silent errors.) - Many SuperNOVAS functions allow `NULL` arguments, both for optional input values as well as outputs that are not required. See the [API Documentation](https://smithsonian.github.io/SuperNOVAS.home/apidoc/html/) for specifics). @@ -507,7 +568,7 @@ on how they are appropriate for the old and new methodologies respectively. - Source names and catalog names can both be up to 64 bytes (including termination), up from 51 and 4 respectively NOVAS C, while keeping `struct` layouts the same thanks to alignment. - - Runtime configurability: + - Runtime configuration: * The planet position calculator function used by `ephemeris` can be set at runtime via `set_planet_provider()`, and `set_planet_provider_hp` (for high precision calculations). Similarly, if `planet_ephem_provider()` or @@ -523,7 +584,7 @@ on how they are appropriate for the old and new methodologies respectively. approximation via `set_nutation_lp_provider()`. For example, the user may want to use the `iau2000b()` model or some custom algorithm instead. - - New intutitive XYZ coordinate coversion functions: + - New intuitive XYZ coordinate conversion functions: * for GCRS - CIRS - ITRS (IAU 2000 standard): `gcrs_to_cirs()`, `cirs_to_itrs()`, and `itrs_to_cirs()`, `cirs_to_gcrs()`. * for GCRS - J2000 - TOD - ITRS (old methodology): `gcrs_to_j2000()`, `j2000_to_tod()`, `tod_to_itrs()`, and @@ -535,6 +596,8 @@ on how they are appropriate for the old and new methodologies respectively. `cirs_to_itrs()` followed by `itrs_to_hor()` does the same from CIRS (new IAU standard methodology), and had no equivalent in NOVAS C 3.1. + - New `ecl2equ()` for converting ecliptic coordinates to equatorial, complementing existing `equ2ecl()`. + - New `gal2equ()` for converting galactic coordinates to ICRS equatorial, complementing existing `equ2gal()`. - New `refract_astro()` function that complements the existing `refract()` but takes an unrefracted (astrometric) @@ -542,6 +605,9 @@ on how they are appropriate for the old and new methodologies respectively. - New convenience functions to wrap `place()` for simpler specific use: `place_star()`, `place_icrs()`, `place_gcrs()`, `place_cirs()`, and `place_tod()`. + + - New `radec_star()` and `radec_planet()` methods as the common point for all existing methods such as `astro_star()` + `local_star()`, `topo_planet()` etc. - New time conversion utilities `tt2tdb()` and `get_ut1_to_tt()` make it simpler to convert between UT1, TT, and TDB time scales, and to supply `ut1_to_tt` arguments to `place()` or topocentric calculations. @@ -549,56 +615,27 @@ on how they are appropriate for the old and new methodologies respectively. - Co-existing `solarsystem()` variants. It is possible to use the different `solarsystem()` implementations provided by `solsys1.c`, `solsys2.c`, `solsys3.c` and/or `solsys-ephem.c` side-by-side, as they define their functionalities with distinct, non-conflicting names, e.g. `earth_sun_calc()` vs `planet_jplint()` vs - `planet_ephem_manager()` vs `planet_ephem_provider()`. See the section on [Building and installation](#installation) - further above on including a selection of these in your library build.) + `planet_ephem_manager()` vs `planet_ephem_provider()`. See the section on + [Building and installation](#installation) further above on including a selection of these in your library + build.) - New `novas_case_sensitive(int)` method to enable (or disable) case-sensitive processing of object names. (By default NOVAS object names were converted to upper-case, making them effectively case-insensitive.) - New `make_planet()` and `make_ephem_object()` to make it simpler to configure Solar-system objects. - - - ------------------------------------------------------------------------------ - - -## Notes on precision - -The SuperNOVAS library is in principle capable of calculating positions to sub-microarcsecond, and velocities to mm/s -precision for all types of celestial sources. However, there are certain pre-requisites and practical considerations -before that level of accuracy is reached. - - - 1. __IAU 2000/2006 conventions__: High precision calculations will generally require that you use SuperNOVAS with the - new IAU standard quantities and methods. The old ways were simply not suited for precisions much below the - milliarcsecond level. - - 2. __Earth's polar motion__: Calculating precise positions for any Earth-based observations requires precise - knowledge of Earth orientation at the time of observation. The pole is subject to predictable precession and - nutation, but also small irregular variations in the orientation of the rotational axis and the rotation period - (a.k.a polar wobble). The [IERS Bulletins](https://www.iers.org/IERS/EN/Publications/Bulletins/bulletins.html) - provide up-to-date measurements, historical data, and near-term projections for the polar offsets and the UT1-UTC - (DUT1) time difference and leap-seconds (UTC-TAI). In SuperNOVAS you can use `cel_pole()` and `get_ut1_to_tt()` - functions to apply / use the published values from these to improve the astrometic precision of Earth-orientation - based coordinate calculations. Without setting and using the actual polar offset values for the time of - observation, positions for Earth-based observations will be accurate at the arcsecond level only. - 3. __Solar-system sources__: Precise calculations for Solar-system sources requires precise ephemeris data for both - the target object as well as for Earth, and the Sun vs the Solar-system barycenter. For the highest precision - calculations you also need positions for all major planets to calculate gravitational deflection precisely. By - default SuperNOVAS can only provide approximate positions for the Earth and Sun (see `earth_sun_calc()` in - `solsys3.c`), but certainly not at the sub-microarcsecond level, and not for other solar-system sources. You will - need to provide a way to interface SuperNOVAS with a suitable ephemeris source (such as the CSPICE toolkit from - JPL) if you want to use it to obtain precise positions for Solar-system bodies. See the - [section below](#solarsystem) for more information how you can do that. - - 4. __Refraction__: Ground based observations are also subject to atmospheric refraction. SuperNOVAS offers the - option to include _optical_ refraction corrections in `equ2hor()` either for a standard atmosphere or more - precisely using the weather parameters defined in the `on_surface` data structure that specifies the observer - locations. Note, that refraction at radio wavelengths is notably different from the included optical model. In - either case you may want to skip the refraction corrections offered in this library, and instead implement your - own as appropriate (or not at all). + - `cel2ter()` and `tel2cel()` can now process 'option'/'class' = 1 (`NOVAS_REFERENCE_CLASS`) regardless of the + methodology (`EROT_ERA` or `EROT_GST`) used to input or output coordinates in GCRS. + - Changed `make_object()` retains the specified number argument (which can be different from the `starnumber` value + in the supplied `cat_entry` structure). + + - Changed the standard atmospheric model for (optical) refraction calculation to include a simple model for the + annual average temperature at the site (based on latitude and elevation). This results is a slightly more educated + guess of the actual refraction than the global fixed temperature of 10 °C assumed by NOVAC C 3.1 regardless of + observing location. + + ----------------------------------------------------------------------------- @@ -624,7 +661,7 @@ Possibly the most universal way to integrate ephemeris data with SuperNOVAS is t `novas_ephem_provider`, e.g.: ```c - int my_ephem_reader(long id, const char *name, double jd_tdb_high, double jd_tdb_low, + int my_ephem_reader(const char *name, long id, double jd_tdb_high, double jd_tdb_low, enum novas_origin *origin, double *pos, double *vel) { // Your custom ephemeris reader implementation here ... @@ -718,7 +755,7 @@ your planetary ephemeris provider at runtime via: Integrating JPL ephemeris data this way can be arduous. You will need to compile and link FORTRAN with C (not the end of the world), but you may also have to modify `jplint.f` to work with the version of `pleph.f` that you will be using. -Unless you already have code that relies on this method, you are probably better off chosing one of the other ways +Unless you already have code that relies on this method, you are probably better off choosing one of the other ways for integrating planetary ephemeris data with SuperNOVAS. @@ -767,7 +804,7 @@ API -- in line with the desire to keep bug-fix releases fully backwards compatib In the month(s) preceding releases one or more _release candidates_ (e.g. `1.0.1-rc3`) will be available on github briefly, under [Releases](https://github.com/Smithsonian/SuperNOVAS/releases), so that changes can be tested by adopters before the releases are finalized. Please use due diligence to test such release candidates with your code -when they become available to avoid unexpected suprises when the finalized release is published. Release candidates +when they become available to avoid unexpected surprises when the finalized release is published. Release candidates are typically available for one week only before they are superseded either by another, or by the finalized release. diff --git a/include/novas.h b/include/novas.h index 80d5edaa..8d3a1e92 100644 --- a/include/novas.h +++ b/include/novas.h @@ -53,11 +53,12 @@ #define SUPERNOVAS_MAJOR_VERSION 1 ///< API major version #define SUPERNOVAS_MINOR_VERSION 0 ///< API minor version #define SUPERNOVAS_SUBVERSION 0 ///< Integer sub version of the release -#define SUPERNOVAS_RELEASE_STRING "" ///< Additional release information in version, e.g. "-1", or "-rc1". +#define SUPERNOVAS_RELEASE_STRING "-rc2" ///< Additional release information in version, e.g. "-1", or "-rc1". /// The version string for this library -#define SUPERNOVAS_VERSION_STRING #SUPERNOVAS_MAJOR_VERSION "." #SUPERNOVAS_MINOR_VERSION "." \ - #SUPERNOVAS_SUBVERSION SUPERNOVAS_RELEASE_STRING +#define SUPERNOVAS_VERSION_STRING #SUPERNOVAS_MAJOR_VERSION "." #SUPERNOVAS_MINOR_VERSION \ + (#SUPERNOVAS_SUBVERSION ? "." #SUPERNOVAS_SUBVERSION : "") \ + SUPERNOVAS_RELEASE_STRING #define NOVAS_MAJOR_VERSION 3 ///< Major version of NOVAS on which this library is based #define NOVAS_MINOR_VERSION 1 ///< Minor version of NOVAS on which this library is based @@ -275,7 +276,7 @@ enum novas_reference_system { enum novas_equator_type { NOVAS_MEAN_EQUATOR = 0, ///< Mean equator without nutation (pre IAU 2006 system). NOVAS_TRUE_EQUATOR, ///< True equator (pre IAU 2006 system). - NOVAS_GCRS_EQUATOR ///< Geocentric Celestial Reference system (GCRS). + NOVAS_GCRS_EQUATOR ///< Geocentric Celestial Reference System (GCRS). }; /** @@ -334,16 +335,16 @@ enum novas_earth_rotation_measure { /// 2006 standard) EROT_ERA = 0, - /// Use GST as the rotation measure, relative to the true equinox (before IAU 20006 standard) + /// Use GST as the rotation measure, relative to the true equinox (pre IAU 20006 standard) EROT_GST }; /** * Constants for ter2cel() and cel2ter() */ -enum novas_celestial_type { - CELESTIAL_GCRS = 0, ///< Celestial coordinates are in GCRS - CELESTIAL_APPARENT ///< Celestial coordinates are apparent values (CIRS or TOD) +enum novas_equatorial_class { + NOVAS_REFERENCE_CLASS = 0, ///< Celestial coordinates are in GCRS + NOVAS_DYNAMICAL_CLASS ///< Celestial coordinates are apparent values (CIRS or TOD) }; /** @@ -496,15 +497,16 @@ typedef struct { double Omega; ///< [rad] mean longitude of the Moon's ascending node. } novas_delaunay_args; -// These sit next to 64-bit values in structures, which means the structure is aligned to 64-bytes. So we -// might as well define names to contain up to 64 bytes, including termination. -#define SIZE_OF_OBJ_NAME 64 ///< Maximum bytes in object names including string termination. -#define SIZE_OF_CAT_NAME 64 ///< Maximum bytes in catalog IDs including string termination. + +#define SIZE_OF_OBJ_NAME 50 ///< Maximum bytes in object names including string termination. +#define SIZE_OF_CAT_NAME 6 ///< Maximum bytes in catalog IDs including string termination. /** - * Basic astrometric data for any celestial object - * located outside the solar system; the catalog - * data for a star + * Basic astrometric data for any sidereal object located outside the solar system. + * + * Note, that despite the slightly expanded catalog name, this has the same memory footprint + * as the original NOVAS C version, allowing for cross-compatible binary exchange (I/O) of + * these structures. * */ typedef struct { @@ -520,7 +522,11 @@ typedef struct { } cat_entry; /** - * Celestial object of interest + * Celestial object of interest. + * + * Note, the memory footprint is different from NOVAS C due to the use of the enum vs short 'type' + * and the long vs. short 'number' values -- hence it is not cross-compatible for binary data + * exchange with NOVAS C 3.1. */ typedef struct { enum novas_object_type type; ///< NOVAS object type @@ -648,12 +654,12 @@ short sidereal_time(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enu double era(double jd_ut1_high, double jd_ut1_low); -short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_earth_rotation_measure method, - enum novas_accuracy accuracy, enum novas_celestial_type class, double xp, double yp, const double *vec1, +short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_earth_rotation_measure erot, + enum novas_accuracy accuracy, enum novas_equatorial_class class, double xp, double yp, const double *vec1, double *vec2); -short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_earth_rotation_measure method, - enum novas_accuracy accuracy, enum novas_celestial_type class, double xp, double yp, const double *vec1, +short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_earth_rotation_measure erot, + enum novas_accuracy accuracy, enum novas_equatorial_class class, double xp, double yp, const double *vec1, double *vec2); int spin(double angle, const double *pos1, double *pos2); @@ -729,7 +735,7 @@ short cio_array(double jd_tdb, long n_pts, ra_of_cio *cio); double ira_equinox(double jd_tdb, enum novas_equinox_type equinox, enum novas_accuracy accuracy); -short ephemeris(const double jd_tdb[2], const object *body, enum novas_origin origin, enum novas_accuracy accuracy, +short ephemeris(const double *jd_tdb, const object *body, enum novas_origin origin, enum novas_accuracy accuracy, double *pos, double *vel); int transform_hip(const cat_entry *hipparcos, cat_entry *hip_2000); @@ -773,7 +779,7 @@ void novas_case_sensitive(int value); int make_planet(enum novas_planet num, object *planet); -int make_ephem_body(const char *name, long num, object *body); +int make_ephem_object(const char *name, long num, object *body); int set_cio_locator_file(const char *filename); @@ -807,15 +813,18 @@ double get_ut1_to_tt(int leap_seconds, double dut1); double get_utc_to_tt(int leap_seconds); +int ecl2equ(double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, double elon, double elat, + double *ra, double *dec); + int gal2equ(double glon, double glat, double *ra, double *dec); // GCRS - CIRS - ITRS conversions int gcrs_to_cirs(double jd_tt, enum novas_accuracy accuracy, const double *in, double *out); -int cirs_to_itrs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, +int cirs_to_itrs(double jd_tt_high, double jd_tt_low, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, double *vec2); -int itrs_to_cirs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, +int itrs_to_cirs(double jd_tt_high, double jd_tt_low, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, double *vec2); int cirs_to_gcrs(double jd_tt, enum novas_accuracy accuracy, const double *in, double *out); @@ -825,10 +834,10 @@ int gcrs_to_j2000(const double *in, double *out); int j2000_to_tod(double jd_tt, enum novas_accuracy accuracy, const double *in, double *out); -int tod_to_itrs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, +int tod_to_itrs(double jd_tt_high, double jd_tt_low, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, double *vec2); -int itrs_to_tod(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, +int itrs_to_tod(double jd_tt_high, double jd_tt_low, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, double *vec2); int tod_to_j2000(double jd_tt, enum novas_accuracy accuracy, const double *in, double *out); diff --git a/include/solarsystem.h b/include/solarsystem.h index 51158fa7..1e0e6f3f 100644 --- a/include/solarsystem.h +++ b/include/solarsystem.h @@ -133,7 +133,7 @@ typedef short (*novas_planet_provider_hp)(const double jd_tdb[2], enum novas_pla * @since 1.0 * @author Attila Kovacs */ -typedef int (*novas_ephem_provider)(long id, const char *name, double jd_tdb_high, double jd_tdb_low, enum novas_origin *origin, double *pos, double *vel); +typedef int (*novas_ephem_provider)(const char *name, long id, double jd_tdb_high, double jd_tdb_low, enum novas_origin *origin, double *pos, double *vel); @@ -268,6 +268,8 @@ short earth_sun_calc(double jd_tdb, enum novas_planet body, enum novas_origin or short earth_sun_calc_hp(const double jd_tdb[2], enum novas_planet body, enum novas_origin origin, double *position, double *velocity); +void enable_earth_sun_hp(int value); + short planet_ephem_provider(double jd_tdb, enum novas_planet body, enum novas_origin origin, double *position, double *velocity); short planet_ephem_provider_hp(const double jd_tdb[2], enum novas_planet body, enum novas_origin origin, double *position, double *velocity); diff --git a/src/novas.c b/src/novas.c index aa7de1cb..5e3ff5fc 100644 --- a/src/novas.c +++ b/src/novas.c @@ -246,7 +246,8 @@ int j2000_to_tod(double jd_tdb, enum novas_accuracy accuracy, const double *in, * coordinate frame. * @param[out] out Output position or velocity vector in rectangular equatorial coordinates at * J2000. It can be the same vector as the input. - * @return 0 if successful, or -1 if either of the vector arguments is NULL. + * @return 0 if successful, or -1 if either of the vector arguments is NULL or the + * 'accuracy' is invalid. * * @sa j2000_to_tod() * @sa tod_to_gcrs() @@ -260,6 +261,11 @@ int tod_to_j2000(double jd_tdb, enum novas_accuracy accuracy, const double *in, return -1; } + if(accuracy != NOVAS_FULL_ACCURACY && accuracy != NOVAS_REDUCED_ACCURACY) { + errno = EINVAL; + return -1; + } + nutation(jd_tdb, NUTATE_TRUE_TO_MEAN, accuracy, in, out); precession(jd_tdb, out, JD_J2000, out); @@ -361,8 +367,9 @@ static int tod_to_gcrs(double jd_tdb, enum novas_accuracy accuracy, const double * @param in GCRS Input (x, y, z) position or velocity vector * @param[out] out Output position or velocity 3-vector in the True equinox of Date coordinate * frame. It can be the same vector as the input. - * @return 0 if successful, or -1 if either of the vector arguments is NULL, or - * an error from cio_location(), or else 10 + the error from cio_basis(). + * @return 0 if successful, or -1 if either of the vector arguments is NULL or the + * accuracy is invalid, or an error from cio_location(), or + * else 10 + the error from cio_basis(). * * @sa gcrs_to_tod() * @sa cirs_to_gcrs() @@ -410,8 +417,9 @@ int gcrs_to_cirs(double jd_tdb, enum novas_accuracy accuracy, const double *in, * @param in CIRS Input (x, y, z) position or velocity vector * @param[out] out Output position or velocity 3-vector in the GCRS coordinate frame. * It can be the same vector as the input. - * @return 0 if successful, or -1 if either of the vector arguments is NULL, or - * an error from cio_location(), or else 10 + the error from cio_basis(). + * @return 0 if successful, or -1 if either of the vector arguments is NULL + * or the accuracy is invalid, or an error from cio_location(), or else + * 10 + the error from cio_basis(). * * @sa tod_to_gcrs() * @sa gcrs_to_cirs() @@ -659,11 +667,13 @@ int place_tod(double jd_tt, const object *source, enum novas_accuracy accuracy, * @param sys Coordinate reference system in which to produce output values * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param[out] ra [h] Topocentric right ascension in hours, referred to true equator and - * equinox of date 'jd_tt'. (It may be NULL if not required) + * equinox of date 'jd_tt' or NAN when returning with an error code. + * (It may be NULL if not required) * @param[out] dec [deg] Topocentric declination in degrees, referred to true equator and - * equinox of date 'jd_tt'. (It may be NULL if not required) - * @param[out] rv [AU/day] radial velocity relative ot observer. (It may be NULL if not - * required) + * equinox of date 'jd_tt' or NAN when returning with an error code. + * (It may be NULL if not required) + * @param[out] rv [AU/day] radial velocity relative ot observer, or NAN when returning with + * an error code. (It may be NULL if not required) * @return 0 if successful, -1 if a required pointer argument is NULL, or else * 20 + the error code from place_star(). * @@ -675,14 +685,22 @@ int place_tod(double jd_tt, const object *source, enum novas_accuracy accuracy, int radec_star(double jd_tt, const cat_entry *star, const observer *obs, double ut1_to_tt, enum novas_reference_system sys, enum novas_accuracy accuracy, double *ra, double *dec, double *rv) { sky_pos output = { }; - int error = place_star(jd_tt, star, obs, ut1_to_tt, sys, accuracy, &output); + int error; - if(ra) *ra = error ? NAN : output.ra; - if(dec) *dec = error ? NAN : output.dec; - if(rv) *rv = error ? NAN : output.rv; + // Default return values in case of error. + if(ra) *ra = NAN; + if(dec) *dec = NAN; + if(rv) *rv = NAN; + error = place_star(jd_tt, star, obs, ut1_to_tt, sys, accuracy, &output); prop_error(error, 20); + if(ra) *ra = output.ra; + if(dec) *dec = output.dec; + if(rv) *rv = output.rv; + + + return 0; } @@ -706,14 +724,16 @@ int radec_star(double jd_tt, const cat_entry *star, const observer *obs, double * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param[out] ra [h] Topocentric apparent right ascension in hours, referred to the - * true equator and equinox of date. (It may be NULL if not required) + * true equator and equinox of date, or NAN when returning with an error + * code. (It may be NULL if not required) * @param[out] dec [deg] Topocentric apparent declination in degrees referred to the - * true equator and equinox of date. (It may be NULL if not required) - * @param[out] dis [AU] True distance from Earth to the body at 'jd_tt' in AU (may be - * NULL if not needed). - * @param[out] rv [AU/day] radial velocity relative ot observer. (It may be NULL if - * not required) - * @return 0 if successful, or -1 if the object argument is NULL, or else 1 if + * true equator and equinox of date, or NAN when returning with an error + * code. (It may be NULL if not required) + * @param[out] dis [AU] True distance from Earth to the body at 'jd_tt' in AU, or NAN when + * returning with an error code. (It may be NULL if not needed). + * @param[out] rv [AU/day] radial velocity relative ot observer, or NAN when returning with + * an error code. (It may be NULL if not required) + * @return 0 if successful, or -1 if the object argument is NULL or if * the value of 'where' in structure 'location' is invalid, or 10 + the * error code from place(). * @@ -727,20 +747,25 @@ int radec_planet(double jd_tt, const object *ss_body, const observer *obs, doubl sky_pos output = { }; int error; + // Default return values in case of error. + if(ra) *ra = NAN; + if(dec) *dec = NAN; + if(dis) *dis = NAN; + if(rv) *rv = NAN; + if(ss_body->type != NOVAS_PLANET && ss_body->type != NOVAS_EPHEM_OBJECT) { errno = EINVAL; - return 1; + return -1; } error = place(jd_tt, ss_body, obs, ut1_to_tt, sys, accuracy, &output); - - if(ra) *ra = error ? NAN : output.ra; - if(dec) *dec = error ? NAN : output.dec; - if(dis) *dis = error ? NAN : output.dis; - if(rv) *rv = error ? NAN : output.rv; - prop_error(error, 10); + if(ra) *ra = output.ra; + if(dec) *dec = output.dec; + if(dis) *dis = output.dis; + if(rv) *rv = output.rv; + return 0; } @@ -1160,8 +1185,8 @@ short local_planet(double jd_tt, const object *ss_body, double ut1_to_tt, const * @param dec [deg] Apparent (TOD) declination in degrees, referred to true equator * and equinox of date. * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) - * @param[out] ira [h] ICRS right ascension in hours. - * @param[out] idec [deg] ICRS declination in degrees. + * @param[out] ira [h] ICRS right ascension in hours, or NAN when returning with an error code. + * @param[out] idec [deg] ICRS declination in degrees, or NAN when returning with an error code. * @return 0 if successful; -1 if the supplied output pointers are NULL, * 1 if the iterative process did not converge after 30 iterations, or an * error from vector2radec(), or else > 10 + an error from app_star(). @@ -1289,8 +1314,13 @@ short place(double jd_tt, const object *target, const observer *location, double int i, error = 0; + if(!target) { + errno = EINVAL; + return -1; + } + // Check for invalid value of 'coord_sys' or 'accuracy'. - if((coord_sys < 0) || (coord_sys >= NOVAS_REFERENCE_SYSTEMS)) { + if(coord_sys < 0 || coord_sys >= NOVAS_REFERENCE_SYSTEMS) { errno = EINVAL; return 1; } @@ -1546,7 +1576,7 @@ int equ2gal(double ra, double dec, double *glon, double *glat) { * @author Attila Kovacs */ int gal2equ(double glon, double glat, double *ra, double *dec) { - double pos1[3], pos2[3], xyproj, g, coslat; + double pos1[3], pos2[3], xyproj, lon, coslat; // Rotation matrix A_g from Hipparcos documentation eq. 1.5.11. // AK: Transposed compared to NOVAS C 3.1 for dot product handling. @@ -1576,10 +1606,10 @@ int gal2equ(double glon, double glat, double *ra, double *dec) { // Decompose galactic vector into longitude and latitude. xyproj = sqrt(pos2[0] * pos2[0] + pos2[1] * pos2[1]); - if(xyproj > 0.0) g = atan2(pos2[1], pos2[0]); - else g = 0.0; + if(xyproj > 0.0) lon = atan2(pos2[1], pos2[0]); + else lon = 0.0; - *ra = g / HOURANGLE; + *ra = lon / HOURANGLE; if(*ra < 0.0) *ra += 24.0; *dec = atan2(pos2[2], xyproj) / DEGREE; @@ -1590,14 +1620,14 @@ int gal2equ(double glon, double glat, double *ra, double *dec) { /** * Convert right ascension and declination to ecliptic longitude and latitude. To convert * GCRS RA and dec to ecliptic coordinates (mean ecliptic and equinox of J2000.0), set - * 'coord_sys' to NOVAS_GCRS_EQUATOR (2); in this case the value of 'jd_tt' can be set to + * 'coord_sys' to NOVAS_GCRS_EQUATOR(2); in this case the value of 'jd_tt' can be set to * anything, since J2000.0 is assumed. Otherwise, all input coordinates are dynamical at * 'jd_tt'. * * @param jd_tt [day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' - * is NOVAS_GCRS_EQUATOR [2]) + * is NOVAS_GCRS_EQUATOR[2]) * @param coord_sys The astrometric reference system of the coordinates. If 'coord_sys' is - * NOVAS_GCRS_EQUATOR (2), the input GCRS coordinates are converted to + * NOVAS_GCRS_EQUATOR(2), the input GCRS coordinates are converted to * J2000 ecliptic coordinates. * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param ra [h] Right ascension in hours, referred to specified equator and equinox @@ -1610,22 +1640,29 @@ int gal2equ(double glon, double glat, double *ra, double *dec) { * equinox of date. * @return 0 if successful, or else 1 if the value of 'coord_sys' is invalid. * - * @sa ecl2equ_vec() + * @sa equ2ecl_vec() + * @sa ecl2equ() + * */ short equ2ecl(double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, double ra, double dec, double *elon, double *elat) { - double r, d, cosd, pos[3], xyproj; + double cosd, pos[3], xyproj; int error; + if(!elon || !elat) { + errno = EINVAL; + return -1; + } + // Form position vector in equatorial system from input coordinates. - r = ra * HOURANGLE; - d = dec * DEGREE; - cosd = cos(d); + ra *= HOURANGLE; + dec *= DEGREE; + cosd = cos(dec); - pos[0] = cosd * cos(r); - pos[1] = cosd * sin(r); - pos[2] = sin(d); + pos[0] = cosd * cos(ra); + pos[1] = cosd * sin(ra); + pos[2] = sin(dec); // Convert the vector from equatorial to ecliptic system. error = equ2ecl_vec(jd_tt, coord_sys, accuracy, pos, pos); @@ -1642,27 +1679,91 @@ short equ2ecl(double jd_tt, enum novas_equator_type coord_sys, enum novas_accura return 0; } +/** + * Convert ecliptic longitude and latitude to right ascension and declination. To convert + * GCRS ecliptic coordinates (mean ecliptic and equinox of J2000.0), set 'coord_sys' to + * NOVAS_GCRS_EQUATOR(2); in this case the value of 'jd_tt' can be set to anything, since + * J2000.0 is assumed. Otherwise, all input coordinates are dynamical at'jd_tt'. + * + * @param jd_tt [day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' + * is NOVAS_GCRS_EQUATOR[2]) + * @param coord_sys The astrometric reference system of the coordinates. If 'coord_sys' is + * NOVAS_GCRS_EQUATOR(2), the input GCRS coordinates are converted to + * J2000 ecliptic coordinates. + * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) + * @param elon [deg] Ecliptic longitude in degrees, referred to specified ecliptic and + * equinox of date. + * @param elat [deg] Ecliptic latitude in degrees, referred to specified ecliptic and + * equinox of date. + * @param[out] ra [h] Right ascension in hours, referred to specified equator and equinox + * of date. + * @param[out] dec [deg] Declination in degrees, referred to specified equator and equinox + * of date. + + * @return 0 if successful, or else 1 if the value of 'coord_sys' is invalid. + * + * @sa ecl2equ_vec() + * @sa equ2ecl() + * + * @since 1.0 + * @author Attila Kovacs + */ +int ecl2equ(double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, double elon, double elat, + double *ra, double *dec) { + + double coslat, pos[3], xyproj; + int error; + + if(!ra || !dec) { + errno = EINVAL; + return -1; + } + + // Form position vector in equatorial system from input coordinates. + elon *= DEGREE; + elat *= DEGREE; + coslat = cos(elat); + + pos[0] = coslat * cos(elon); + pos[1] = coslat * sin(elon); + pos[2] = sin(elat); + + // Convert the vector from equatorial to ecliptic system. + error = ecl2equ_vec(jd_tt, coord_sys, accuracy, pos, pos); + if(error) return error; + + // Decompose ecliptic vector into ecliptic longitude and latitude. + xyproj = sqrt(pos[0] * pos[0] + pos[1] * pos[1]); + + *ra = (xyproj > 0.0) ? atan2(pos[1], pos[0]) / HOURANGLE : 0.0; + if(*ra < 0.0) *ra += 24.0; + + *dec = atan2(pos[2], xyproj) / DEGREE; + + return 0; +} + /** * Converts an equatorial position vector to an ecliptic position vector. To convert * ICRS RA and dec to ecliptic coordinates (mean ecliptic and equinox of J2000.0), set - * 'coord_sys' to NOVAS_GCRS_EQUATOR (2); in this case the value of 'jd_tt' can be set + * 'coord_sys' to NOVAS_GCRS_EQUATOR(2); in this case the value of 'jd_tt' can be set * to anything, since J2000.0 is assumed. Otherwise, all input coordinates are dynamical * at 'jd_tt'. * * @param jd_tt [day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' - * is NOVAS_GCRS_EQUATOR [2]) + * is NOVAS_GCRS_EQUATOR[2]) * @param coord_sys The astrometric reference system type of the coordinates. * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param pos1 Position vector, referred to specified equator and equinox of date. * @param[out] pos2 Position vector, referred to specified ecliptic and equinox of date. * It can be the same vector as the input. If 'coord_sys' is - * NOVAS_GCRS_EQUATOR (2), the input GCRS coordinates are converted to + * NOVAS_GCRS_EQUATOR(2), the input GCRS coordinates are converted to * J2000 ecliptic coordinates. * @return 0 if successful, -1 if either vector argument is NULL or the accuracy * is invalid, or else 1 if the value of 'coord_sys' is invalid. * - * @sa ecl2equ_vec() * @sa equ2ecl() + * @sa ecl2equ_vec() */ short equ2ecl_vec(double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, const double *pos1, double *pos2) { @@ -1732,12 +1833,12 @@ short equ2ecl_vec(double jd_tt, enum novas_equator_type coord_sys, enum novas_ac /** * Converts an ecliptic position vector to an equatorial position vector. To convert * ecliptic coordinates (mean ecliptic and equinox of J2000.0) to GCRS RA and dec to, set - * 'coord_sys' to NOVAS_GCRS_EQUATOR (2); in this case the value of 'jd_tt' can be set to + * 'coord_sys' to NOVAS_GCRS_EQUATOR(2); in this case the value of 'jd_tt' can be set to * anything, since J2000.0 is assumed. Otherwise, all input coordinates are dynamical at * 'jd_tt'. * * @param jd_tt [day] Terrestrial Time (TT) based Julian date. (Unused if 'coord_sys' - * is NOVAS_GCRS_EQUATOR [2]) + * is NOVAS_GCRS_EQUATOR[2]) * @param coord_sys The astrometric reference system type of the coordinates * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param pos1 Position vector, referred to specified ecliptic and equinox of date. @@ -1746,6 +1847,7 @@ short equ2ecl_vec(double jd_tt, enum novas_equator_type coord_sys, enum novas_ac * @return 0 if successful, -1 if either vector argument is NULL or the accuracy * is invalid, or else 1 if the value of 'coord_sys' is invalid. * + * @sa ecl2equ() * @sa equ2ecl_vec() */ short ecl2equ_vec(double jd_tt, enum novas_equator_type coord_sys, enum novas_accuracy accuracy, const double *pos1, @@ -1833,6 +1935,7 @@ short ecl2equ_vec(double jd_tt, enum novas_equator_type coord_sys, enum novas_ac * @sa hor_to_itrs() * @sa cirs_to_itrs() * @sa tod_to_itrs() + * @sa refract_astro() * * @since 1.0 * @author Attila Kovacs @@ -1843,6 +1946,8 @@ int itrs_to_hor(const on_surface *location, const double *in, double *az, double double pn, pw, pz, proj; if(!location || !in) { + if(az) *az = NAN; + if(za) *za = NAN; errno = EINVAL; return -1; } @@ -1909,6 +2014,7 @@ int itrs_to_hor(const on_surface *location, const double *in, double *az, double * @sa itrs_to_hor() * @sa itrs_to_cirs() * @sa itrs_to_tod() + * @sa refract() * * @since 1.0 * @author Attila Kovacs @@ -2072,9 +2178,9 @@ int equ2hor(double jd_ut1, double ut1_to_tt, enum novas_accuracy accuracy, doubl // Rotate Earth-fixed orthonormal basis vectors to celestial system // (wrt equator and equinox of date). - ter2cel(jd_ut1, 0.0, ut1_to_tt, EROT_GST, accuracy, CELESTIAL_APPARENT, xp, yp, uze, uz); - ter2cel(jd_ut1, 0.0, ut1_to_tt, EROT_GST, accuracy, CELESTIAL_APPARENT, xp, yp, une, un); - ter2cel(jd_ut1, 0.0, ut1_to_tt, EROT_GST, accuracy, CELESTIAL_APPARENT, xp, yp, uwe, uw); + ter2cel(jd_ut1, 0.0, ut1_to_tt, EROT_GST, accuracy, NOVAS_DYNAMICAL_CLASS, xp, yp, uze, uz); + ter2cel(jd_ut1, 0.0, ut1_to_tt, EROT_GST, accuracy, NOVAS_DYNAMICAL_CLASS, xp, yp, une, un); + ter2cel(jd_ut1, 0.0, ut1_to_tt, EROT_GST, accuracy, NOVAS_DYNAMICAL_CLASS, xp, yp, uwe, uw); // Compute coordinates of object w.r.t orthonormal basis. @@ -2229,12 +2335,12 @@ short gcrs2equ(double jd_tt, enum novas_dynamical_type sys, enum novas_accuracy * @param ut1_to_tt [s] TT - UT1 Time difference in seconds * @param gst_type NOVAS_MEAN_EQUINOX (0) or NOVAS_TRUE_EQUINOX (1), depending on whether * wanting mean or apparent GST, respectively. - * @param method EROT_ERA (0) or EROT_GST (1), depending on whether to use GST relative + * @param erot EROT_ERA (0) or EROT_GST (1), depending on whether to use GST relative * to equinox of date (pre IAU 2006) or ERA relative to the CIO (IAU 2006 * standard). * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param[out] gst [h] Greenwich (mean or apparent) sidereal time, in hours [0:24]. (In case - * the returned error code is >1 the gst value will be set to NAN. + * the returned error code is >1 the gst value will be set to NAN.) * @return 0 if successful, or -1 if the 'gst' argument is NULL, 1 if 'accuracy' is * invalid 2 if 'method' is invalid, or else 10--30 with 10 + the error from * cio_location(). @@ -2246,7 +2352,7 @@ short gcrs2equ(double jd_tt, enum novas_dynamical_type sys, enum novas_accuracy * @sa get_ut1_to_tt() */ short sidereal_time(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_equinox_type gst_type, - enum novas_earth_rotation_measure method, enum novas_accuracy accuracy, double *gst) { + enum novas_earth_rotation_measure erot, enum novas_accuracy accuracy, double *gst) { double jd_ut, jd_tt, jd_tdb, t, theta, st, eqeq; @@ -2255,6 +2361,7 @@ short sidereal_time(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enu return -1; } + // Default return value *gst = NAN; if(accuracy != NOVAS_FULL_ACCURACY && accuracy != NOVAS_REDUCED_ACCURACY) { @@ -2275,8 +2382,8 @@ short sidereal_time(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enu // Compute the equation of the equinoxes if needed, depending upon the // input values of 'gst_type' and 'method'. If not needed, set to zero. - if(((gst_type == NOVAS_MEAN_EQUINOX) && (method == EROT_GST)) // GMST; CIO-TIO - || ((gst_type == NOVAS_TRUE_EQUINOX) && (method == EROT_ERA))) { // GAST; equinox + if(((gst_type == NOVAS_MEAN_EQUINOX) && (erot == EROT_ERA)) // GMST; CIO-TIO + || ((gst_type == NOVAS_TRUE_EQUINOX) && (erot == EROT_GST))) { // GAST; equinox static enum novas_accuracy acc_last = -1; static double jd_last = -1e100; static double ee; @@ -2297,8 +2404,8 @@ short sidereal_time(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enu // Compute Greenwich sidereal time depending upon input values of // method' and 'gst_type'. - switch(method) { - case (EROT_GST): { + switch(erot) { + case EROT_ERA: { // Use 'CIO-TIO-theta' method. See Circular 179, Section 6.5.4. const double ux[3] = { 1.0, 0.0, 0.0 }; double ra_cio, ha_eq, x[3], y[3], z[3], w1[3], w2[3], eq[3]; @@ -2331,7 +2438,7 @@ short sidereal_time(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enu return 0; } - case (EROT_ERA): + case EROT_GST: // Use equinox method. See Circular 179, Section 2.6.2. @@ -2402,10 +2509,10 @@ double era(double jd_ut1_high, double jd_ut1_low) { * * If both 'xp' and 'yp' are set to 0 no polar motion is included in the transformation. * - * @deprecated This function can be confusing to use. For the IAU standard methodology - * use itrs_to_cirs() instead, while for the old methodology use itrs_to_tod() - * instead. You can then follow these with other conversions to GCRS or - * J2000 as appropriate. + * @deprecated This function can be confusing to use due to the output coordinate system + * being specified by a combination of two options. Use itrs_to_cirs() or + * itrs_to_tod() instead. You can then follow these with other conversions to + * GCRS (or whatever else) as appropriate. * * REFERENCES: *
    @@ -2417,12 +2524,16 @@ double era(double jd_ut1_high, double jd_ut1_low) { * @param jd_ut1_high [day] High-order part of UT1 Julian date. * @param jd_ut1_low [day] Low-order part of UT1 Julian date. * @param ut1_to_tt [s] TT - UT1 Time difference in seconds - * @param method EROT_ERA (0) or EROT_GST (1), depending on whether to use GST relative + * @param erot EROT_ERA (0) or EROT_GST (1), depending on whether to use GST relative * to equinox of date (pre IAU 2006) or ERA relative to the CIO (IAU 2006 - * standard). + * standard) as the Earth rotation measure. The main effect of this option + * is that it selects the output coordinate system as CIRS or TOD if + * the output coordinate class is NOVAS_DYNAMICAL_CLASS. * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) - * @param class Output coordinate class CELESTIAL_GCRS (0, or any value other than 1) - * or CELESTIAL_APPARENT (1). + * @param class Output coordinate class NOVAS_REFERENCE_CLASS (0, or any value other than 1) + * or NOVAS_DYNAMICAL_CLASS (1). Use the former if the output coordinates are + * to be in the GCRS, and the latter if they are to be in CIRS or TOD (the 'erot' + * parameter selects which dynamical system to use for the output.) * @param xp [arcsec] Conventionally-defined X coordinate of celestial intermediate * pole with respect to ITRS pole, in arcseconds. * @param yp [arcsec] Conventionally-defined Y coordinate of celestial intermediate @@ -2431,7 +2542,8 @@ double era(double jd_ut1_high, double jd_ut1_low) { * referred to ITRS axes (terrestrial system) in the normal case * where 'option' is NOVAS_GCRS (0). * @param[out] vec2 Position vector, equatorial rectangular coordinates in the specified - * output system + * output system (GCRS if 'class' is NOVAS_REFERENCE_CLASS; + * or else either CIRS if 'erot' is EROT_ERA, or TOD if 'erot' is EROT_GST). * @return 0 if successful, -1 if either of the vector arguments is NULL, 1 if * 'accuracy' is invalid, 2 if 'method' is invalid 10--20, or else 10 + the * error from cio_location(), or 20 + error from cio_basis(). @@ -2443,8 +2555,8 @@ double era(double jd_ut1_high, double jd_ut1_low) { * @sa frame_tie() * @sa cel2ter() */ -short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_earth_rotation_measure method, - enum novas_accuracy accuracy, enum novas_celestial_type class, double xp, double yp, const double *vec1, +short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_earth_rotation_measure erot, + enum novas_accuracy accuracy, enum novas_equatorial_class class, double xp, double yp, const double *vec1, double *vec2) { double jd_ut1, jd_tt, jd_tdb, gast; @@ -2473,7 +2585,7 @@ short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova } else wobble(jd_tt, WOBBLE_ITRS_TO_PEF, xp, yp, vec1, vec2); - switch(method) { + switch(erot) { case (EROT_ERA): { // 'CIO-TIO-THETA' method. See second reference, eq. (3) and (4). @@ -2482,7 +2594,7 @@ short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova const double theta = era(jd_ut1_high, jd_ut1_low); spin(-theta, vec2, vec2); - if(class != CELESTIAL_APPARENT) { + if(class != NOVAS_DYNAMICAL_CLASS) { const int error = cirs_to_gcrs(jd_tdb, accuracy, vec2, vec2); prop_error(error, 10); } @@ -2493,7 +2605,7 @@ short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova sidereal_time(jd_ut1_high, jd_ut1_low, ut1_to_tt, NOVAS_TRUE_EQUINOX, EROT_GST, accuracy, &gast); spin(-15.0 * gast, vec2, vec2); - if(class != CELESTIAL_APPARENT) { + if(class != NOVAS_DYNAMICAL_CLASS) { tod_to_gcrs(jd_tdb, accuracy, vec2, vec2); } break; @@ -2513,6 +2625,11 @@ short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova * * If both 'xp' and 'yp' are set to 0 no polar motion is included in the transformation. * + * If extreme (sub-microarcsecond) accuracy is not required, you can use UT1-based Julian date + * instead of the TT-based Julian date and set the 'ut1_to_tt' argument to 0.0. and you can + * use UTC-based Julian date the same way.for arcsec-level precision also. + * + * * REFERENCES: *
      *
    1. Kaplan, G. H. et. al. (1989). Astron. Journ. 97, 1197-1210.
      2. @@ -2520,7 +2637,8 @@ short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova * XXV Joint Discussion 16. *
    * - * @param jd_tt [day] Terrestrial Time (TT) based Julian date. + * @param jd_tt_high [day] High-order part of Terrestrial Time (TT) based Julian date. + * @param jd_tt_low [day] Low-order part of Terrestrial Time (TT) based Julian date. * @param ut1_to_tt [s] TT - UT1 Time difference in seconds * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param xp [arcsec] Conventionally-defined X coordinate of celestial intermediate @@ -2542,17 +2660,22 @@ short ter2cel(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova * @since 1.0 * @author Attila Kovacs */ -int itrs_to_cirs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, double xp, +int itrs_to_cirs(double jd_tt_high, double jd_tt_low, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, double *vec2) { - return ter2cel(jd_tt, -ut1_to_tt, ut1_to_tt, EROT_ERA, accuracy, CELESTIAL_APPARENT, xp, yp, vec1, vec2); + return ter2cel(jd_tt_high, jd_tt_low - ut1_to_tt, ut1_to_tt, EROT_ERA, accuracy, NOVAS_DYNAMICAL_CLASS, xp, yp, vec1, vec2); } + /** * Rotates a position vector from the Earth-fixed ITRS frame to the dynamical True of Date * (TOD) frame of date (pre IAU 2000 method). * * If both 'xp' and 'yp' are set to 0 no polar motion is included in the transformation. * + * If extreme (sub-microarcsecond) accuracy is not required, you can use UT1-based Julian date + * instead of the TT-based Julian date and set the 'ut1_to_tt' argument to 0.0. and you can + * use UTC-based Julian date the same way.for arcsec-level precision also. + * * REFERENCES: *
      *
    1. Kaplan, G. H. et. al. (1989). Astron. Journ. 97, 1197-1210.
      2. @@ -2560,7 +2683,8 @@ int itrs_to_cirs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, d * XXV Joint Discussion 16. *
    * - * @param jd_tt [day] Terrestrial Time (TT) based Julian date. + * @param jd_tt_high [day] High-order part of Terrestrial Time (TT) based Julian date. + * @param jd_tt_low [day] Low-order part of Terrestrial Time (TT) based Julian date. * @param ut1_to_tt [s] TT - UT1 Time difference in seconds * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param xp [arcsec] Conventionally-defined X coordinate of celestial intermediate @@ -2582,11 +2706,12 @@ int itrs_to_cirs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, d * @since 1.0 * @author Attila Kovacs */ -int itrs_to_tod(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, double xp, +int itrs_to_tod(double jd_tt_high, double jd_tt_low, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, double *vec2) { - return ter2cel(jd_tt, -ut1_to_tt, ut1_to_tt, EROT_GST, accuracy, CELESTIAL_APPARENT, xp, yp, vec1, vec2); + return ter2cel(jd_tt_high, jd_tt_low - ut1_to_tt, ut1_to_tt, EROT_GST, accuracy, NOVAS_DYNAMICAL_CLASS, xp, yp, vec1, vec2); } + /** * Rotates a vector from the celestial to the terrestrial system. Specifically, it transforms * a vector in the GCRS, or the dynamcal CIRS or TOD frames to the ITRS (a rotating earth-fixed @@ -2598,10 +2723,10 @@ int itrs_to_tod(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, do * * If both 'xp' and 'yp' are set to 0 no polar motion is included in the transformation. * - * @deprecated This function can be confusing to use. For the IAU standard methodology - * use cirs_to_itrs() instead, while for the old methodology use tod_to_itrs() - * instead. You can then precede these with other conversions from GCRS or - * J2000 as appropriate. + * @deprecated This function can be confusing to use due to the input coordinate system + * being specified by a combination of two options. Use itrs_to_cirs() or + * itrs_to_tod() instead. You can then follow these with other conversions to + * GCRS (or whatever else) as appropriate. * * REFERENCES: *
      @@ -2615,20 +2740,25 @@ int itrs_to_tod(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, do * @param jd_ut1_high [day] High-order part of UT1 Julian date. * @param jd_ut1_low [day] Low-order part of UT1 Julian date. * @param ut1_to_tt [s] TT - UT1 Time difference in seconds - * @param method EROT_ERA (0) or EROT_GST (1), depending on whether to use GST relative to - * equinox of date (pre IAU 2006) or ERA relative to the CIO (IAU 2006 standard). + * @param erot EROT_ERA (0) or EROT_GST (1), depending on whether to use GST relative to + * equinox of date (pre IAU 2006) or ERA relative to the CIO (IAU 2006 standard) + * as the Earth rotation measure. The main effect of this option + * is that it specifies the input coordinate system as CIRS or TOD when + * the input coordinate class is NOVAS_DYNAMICAL_CLASS. * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) - * @param class Input coordinate class, CELESTIAL_GCRS (0, or any value other than 1) or - * CELESTIAL_APPARENT (1). + * @param class Input coordinate class, NOVAS_REFERENCE_CLASS (0) or NOVAS_DYNAMICAL_CLASS (1). + * Use the former if the input coordinates are in the GCRS, and the latter if they + * are CIRS or TOD (the 'erot' parameter selects which dynamical system the input is + * specified in.) * @param xp [arcsec] Conventionally-defined X coordinate of celestial intermediate * pole with respect to ITRS pole, in arcseconds. * @param yp [arcsec] Conventionally-defined Y coordinate of celestial intermediate * pole with respect to ITRS pole, in arcseconds. - * @param vec1 Position vector, geocentric equatorial rectangular coordinates in the - * specified input coordinate system. - * @param[out] vec2 Position vector, geocentric equatorial rectangular coordinates, - * referred to ITRS axes (terrestrial system). It can be the same vector as the - * input. + * @param vec1 Input position vector, geocentric equatorial rectangular coordinates in the + * specified input coordinate system (GCRS if 'class' is NOVAS_REFERENCE_CLASS; + * or else either CIRS if 'erot' is EROT_ERA, or TOD if 'erot' is EROT_GST). + * @param[out] vec2 ITRS position vector, geocentric equatorial rectangular coordinates + * (terrestrial system). It can be the same vector as the input. * @return 0 if successful, -1 if either of the vector arguments is NULL, 1 if 'accuracy' * is invalid, 2 if 'method' is invalid, or else 10 + the error from * cio_location(), or 20 + error from cio_basis(). @@ -2640,8 +2770,8 @@ int itrs_to_tod(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, do * @sa tod_to_itrs() * @sa ter2cel() */ -short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_earth_rotation_measure method, - enum novas_accuracy accuracy, enum novas_celestial_type class, double xp, double yp, const double *vec1, +short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum novas_earth_rotation_measure erot, + enum novas_accuracy accuracy, enum novas_equatorial_class class, double xp, double yp, const double *vec1, double *vec2) { double jd_ut1, jd_tt, jd_tdb, gast, theta; @@ -2663,10 +2793,10 @@ short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova // Compute the TDB Julian date corresponding to the input UT1 Julian date jd_tdb = jd_tt + tt2tdb(jd_tt) / DAY; - switch(method) { + switch(erot) { case (EROT_ERA): // IAU 2006 standard method - if(class != CELESTIAL_APPARENT) { + if(class != NOVAS_DYNAMICAL_CLASS) { // See second reference, eq. (3) and (4). int error = gcrs_to_cirs(jd_tt, accuracy, vec1, vec2); prop_error(error, 10); @@ -2681,7 +2811,7 @@ short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova case (EROT_GST): // Pre IAU 2006 method - if(class == CELESTIAL_APPARENT) { + if(class == NOVAS_DYNAMICAL_CLASS) { if(vec2 != vec1) memcpy(vec2, vec1, XYZ_VECTOR_SIZE); } else { @@ -2712,6 +2842,11 @@ short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova * * If both 'xp' and 'yp' are set to 0 no polar motion is included in the transformation. * + * If extreme (sub-microarcsecond) accuracy is not required, you can use UT1-based Julian date + * instead of the TT-based Julian date and set the 'ut1_to_tt' argument to 0.0. and you can + * use UTC-based Julian date the same way.for arcsec-level precision also. + * + * * REFERENCES: *
        *
      1. Kaplan, G. H. et. al. (1989). Astron. Journ. 97, 1197-1210.
        2. @@ -2719,7 +2854,8 @@ short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova * Joint Discussion 16. *
      * - * @param jd_tt [day] Terrestrial Time (TT) based Julian date. + * @param jd_tt_high [day] High-order part of Terrestrial Time (TT) based Julian date. + * @param jd_tt_low [day] Low-order part of Terrestrial Time (TT) based Julian date. * @param ut1_to_tt [s] TT - UT1 Time difference in seconds * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param xp [arcsec] Conventionally-defined X coordinate of celestial intermediate @@ -2742,17 +2878,22 @@ short cel2ter(double jd_ut1_high, double jd_ut1_low, double ut1_to_tt, enum nova * @since 1.0 * @author Attila Kovacs */ -int cirs_to_itrs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, double xp, +int cirs_to_itrs(double jd_tt_high, double jd_tt_low, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, double *vec2) { - return cel2ter(jd_tt, -ut1_to_tt, ut1_to_tt, EROT_ERA, accuracy, CELESTIAL_APPARENT, xp, yp, vec1, vec2); + return cel2ter(jd_tt_high, jd_tt_low - ut1_to_tt, ut1_to_tt, EROT_ERA, accuracy, NOVAS_DYNAMICAL_CLASS, xp, yp, vec1, vec2); } + /** * Rotates a position vector from the dynamical True of Date (TOD) frame of date the Earth-fixed * ITRS frame (pre IAU 2000 method). * * If both 'xp' and 'yp' are set to 0 no polar motion is included in the transformation. * + * If extreme (sub-microarcsecond) accuracy is not required, you can use UT1-based Julian date + * instead of the TT-based Julian date and set the 'ut1_to_tt' argument to 0.0. and you can + * use UTC-based Julian date the same way.for arcsec-level precision also. + * * REFERENCES: *
        *
      1. Kaplan, G. H. et. al. (1989). Astron. Journ. 97, 1197-1210.
        2. @@ -2760,8 +2901,9 @@ int cirs_to_itrs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, d * Joint Discussion 16. *
      * - * @param jd_tt [day] Terrestrial Time (TT) based Julian date. - * @param ut1_to_tt [s] TT - UT1 Time difference in seconds + * @param jd_tt_high [day] High-order part of Terrestrial Time (TT) based Julian date. + * @param jd_tt_low [day] Low-order part of Terrestrial Time (TT) based Julian date. + * @param ut1_to_tt [s] TT - UT1 Time difference in seconds. * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param xp [arcsec] Conventionally-defined X coordinate of celestial intermediate * pole with respect to ITRS pole, in arcseconds. @@ -2783,9 +2925,9 @@ int cirs_to_itrs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, d * @since 1.0 * @author Attila Kovacs */ -int tod_to_itrs(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, double xp, +int tod_to_itrs(double jd_tt_high, double jd_tt_low, double ut1_to_tt, enum novas_accuracy accuracy, double xp, double yp, const double *vec1, double *vec2) { - return cel2ter(jd_tt, -ut1_to_tt, ut1_to_tt, EROT_ERA, accuracy, CELESTIAL_APPARENT, xp, yp, vec1, vec2); + return cel2ter(jd_tt_high, jd_tt_low - ut1_to_tt, ut1_to_tt, EROT_ERA, accuracy, NOVAS_DYNAMICAL_CLASS, xp, yp, vec1, vec2); } /** @@ -2937,8 +3079,8 @@ int wobble(double jd_tt, enum novas_wobble_direction direction, double xp, doubl * and equinox of date, components in AU/day. (It may be NULL if * no velocity data is required). * - * @return 0 if successful, or -1 if the pos and vel output arguments are identical - * pointers. + * @return 0 if successful, or -1 if location is NULL or if the pos and vel output + * arguments are identical pointers. * * @sa make_on_surface() * @sa geo_posvel() @@ -2949,7 +3091,7 @@ int terra(const on_surface *location, double lst, double *pos, double *vel) { double ht_km; int j; - if(pos == vel) { + if(!location || pos == vel) { errno = EINVAL; return -1; } @@ -3399,8 +3541,8 @@ int frame_tie(const double *pos1, enum novas_frametie_direction direction, doubl *
    * * @param jd_tdb1 [day] Barycentric Dynamical Time (TDB) based Julian date of the first epoch. - * @param pos Position vector at first epoch. - * @param vel Velocity vector at first epoch. + * @param pos [AU] Position vector at first epoch. + * @param vel [AU/day] Velocity vector at first epoch. * @param jd_tdb2 [day] Barycentric Dynamical Time (TDB) based Julian date of the second epoch. * @param[out] pos2 Position vector at second epoch. It can be the same vector as the input. * @return 0 if successful, or -1 if any of the vector areguments is NULL. @@ -3451,14 +3593,13 @@ int bary2obs(const double *pos, const double *pos_obs, double *pos2, double *lig int j; if(!pos || !pos_obs || !pos2) { + if(lighttime) *lighttime = NAN; errno = EINVAL; return -1; } // Translate vector to geocentric coordinates. - for(j = 3; --j >= 0;) { - pos2[j] = pos[j] - pos_obs[j]; - } + for(j = 3; --j >= 0;) pos2[j] = pos[j] - pos_obs[j]; // Calculate length of vector in terms of light time. if(lighttime) *lighttime = vlen(pos2) / C_AUDAY; @@ -3587,7 +3728,7 @@ short geo_posvel(double jd_tt, double ut1_to_tt, enum novas_accuracy accuracy, c * geocenter), referred to ICRS axes, components in AU. * @param[out] vsb [AU/day] Velocity 3-vector of body, with respect to the Solar-system * barycenter, referred to ICRS axes, components in AU/day. - * @param[out] tlight [day] Calculated light time + * @param[out] tlight [day] Calculated light time, or NAN when returning with an error code. * * @return 0 if successful, -1 if any of the pointer arguments is NULL or if the * output vectors are the same or if they are the same as pos_obs, 1 if @@ -3611,6 +3752,7 @@ int light_time2(double jd_tdb, const object *body, const double *pos_obs, double return -1; } + // Default return value. *tlight = NAN; if(!body || !pos_obs || prel == pos_obs || vsb == pos_obs || prel == vsb) { @@ -3685,7 +3827,7 @@ short light_time(double jd_tdb, const object *body, const double *pos_obs, doubl } /** - * Teturns the difference in light-time, for a star, between the barycenter of the solar system + * Returns the difference in light-time, for a star, between the barycenter of the solar system * and the observer (or the geocenter). * * Alternatively, this function returns the light-time from the observer (or the geocenter) to @@ -3991,7 +4133,7 @@ int aberration(const double *pos, const double *vobs, double lighttime, double * vemag = vlen(vobs); if(!vemag) { - memcpy(pos2, pos, XYZ_VECTOR_SIZE); + if(pos2 != pos) memcpy(pos2, pos, XYZ_VECTOR_SIZE); return 0; } @@ -4082,6 +4224,7 @@ int rad_vel(const object *target, const double *pos, const double *vel, const do int i; if(!target || !pos || !vel || !vel_obs || !rv) { + if(rv) *rv = NAN; errno = EINVAL; return -1; } @@ -4435,6 +4578,8 @@ int set_nutation_lp_provider(novas_nutation_provider func) { */ int nutation_angles(double t, enum novas_accuracy accuracy, double *dpsi, double *deps) { if(!dpsi || !deps) { + if(dpsi) *dpsi = NAN; + if(deps) *deps = NAN; errno = EINVAL; return -1; } @@ -4510,7 +4655,8 @@ int fund_args(double t, novas_delaunay_args *a) { * * @param t [cy] Julian centuries since J2000 * @param planet Novas planet id, e.g. NOVAS_MARS. - * @return [rad] The approximate longitude of the planet in radians [-π:π]. + * @return [rad] The approximate longitude of the planet in radians [-π:π], + * or NAN if the `planet` id is out of range. * * @sa accum_prec() * @sa nutation_angles() @@ -4606,19 +4752,23 @@ double mean_obliq(double jd_tdb) { *
  1. Kaplan, G. H. et. al. (1989). Astron. Journ. 97, 1197-1210.
    2. *
* - * @param pos Position 3-vector, equatorial rectangular coordinates. - * @param ra [h] Right ascension in hours [0:24]. - * @param dec [deg] Declination in degrees [-90:90]. - * @return 0 if successful, -1 of any of the arguments are NULL, or - * 1 if all input components are 0 so 'ra' and 'dec' are indeterminate, - * or else 2 if both x and y are zero, but z is nonzero, and so 'ra' is - * indeterminate. + * @param pos Position 3-vector, equatorial rectangular coordinates. + * @param[out] ra [h] Right ascension in hours [0:24] or NAN if the position vector is NULL or a null-vector. + * @param[out] dec [deg] Declination in degrees [-90:90] or NAN if the position vector is NULL or a null-vector. + * @return 0 if successful, -1 of any of the arguments are NULL, or + * 1 if all input components are 0 so 'ra' and 'dec' are indeterminate, + * or else 2 if both x and y are zero, but z is nonzero, and so 'ra' is + * indeterminate. * * @sa radec2vector() */ short vector2radec(const double *pos, double *ra, double *dec) { double xyproj; + // Default return values. + if(ra) *ra = NAN; + if(dec) *dec = NAN; + if(!pos || !ra || !dec) { errno = EINVAL; return -1; @@ -4626,14 +4776,12 @@ short vector2radec(const double *pos, double *ra, double *dec) { xyproj = sqrt(pos[0] * pos[0] + pos[1] * pos[1]); if(xyproj == 0.0) { - *ra = 0.0; - if(pos[2] == 0) { - *dec = 0.0; errno = EINVAL; return 1; } + *ra = 0.0; *dec = (pos[2] < 0.0) ? -90.0 : 90.0; errno = EDOM; return 2; @@ -4808,7 +4956,8 @@ double get_ut1_to_tt(int leap_seconds, double dut1) { * * * @param jd_tdb [day] Barycentric Dynamic Time (TDB) based Julian date - * @param[out] jd_tt [day] Terrestrial Time (TT) based Julian date + * @param[out] jd_tt [day] Terrestrial Time (TT) based Julian date. (It may be NULL + * if not required) * @param[out] secdiff [s] Difference 'tdb_jd'-'tt_jd', in seconds. (It may be NULL if * not required) * @return 0 if successful, or -1 if the tt_jd pointer argument is NULL. @@ -4817,21 +4966,14 @@ double get_ut1_to_tt(int leap_seconds, double dut1) { */ int tdb2tt(double jd_tdb, double *jd_tt, double *secdiff) { const double t = (jd_tdb - JD_J2000) / JULIAN_CENTURY_DAYS; - double d; - - if(!jd_tt) { - errno = EINVAL; - return -1; - } // Expression given in USNO Circular 179, eq. 2.6. - d = 0.001657 * sin(628.3076 * t + 6.2401) + 0.000022 * sin(575.3385 * t + 4.2970) + const double d = 0.001657 * sin(628.3076 * t + 6.2401) + 0.000022 * sin(575.3385 * t + 4.2970) + 0.000014 * sin(1256.6152 * t + 6.1969) + 0.000005 * sin(606.9777 * t + 4.0212) + 0.000005 * sin(52.9691 * t + 0.4444) + 0.000002 * sin(21.3299 * t + 5.5431) + 0.000010 * t * sin(628.3076 * t + 4.2490); - *jd_tt = jd_tdb - d / DAY; - + if(jd_tt) *jd_tt = jd_tdb - d / DAY; if(secdiff) *secdiff = d; return 0; @@ -4874,7 +5016,7 @@ double tt2tdb(double jd_tt) { * @param jd_tt [day] Terrestrial Time (TT) based Julian date * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) * @param[out] ra_cio [h] Right ascension of the CIO, with respect to the true equinox of - * date, in hours (+ or -). + * date, in hours (+ or -), or NAN when returning with an error code. * @return 0 if successful, -1 if the output pointer argument is NULL, * 1 if 'accuracy' is invalid, 10--20: 10 + error code from cio_location(), * or else 20 + error from cio_basis() @@ -4891,6 +5033,7 @@ short cio_ra(double jd_tt, enum novas_accuracy accuracy, double *ra_cio) { return -1; } + // Default return value. *ra_cio = NAN; // Check for valid value of 'accuracy'. @@ -4978,11 +5121,13 @@ int set_cio_locator_file(const char *filename) { * * @param jd_tdb [day] Barycentric Dynamic Time (TDB) based Julian date * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1) - * @param[out] ra_cio [h] Right ascension of the CIO, in hours. + * @param[out] ra_cio [h] Right ascension of the CIO, in hours, or NAN if returning + * with an error. * @param[out] loc_type Pointer in which to return the reference system in which right * ascension is given, which is either CIO_VS_GCRS (1) if the * location was obtained via interpolation of the available data * file, or else CIO_VS_EQUINOX (2) if it was calculated locally. + * It is set to -1 if retrurning with an error. * * @return 0 if successful, -1 if one of the pointer arguments is NULL or the * accuracy is invalid. @@ -4997,6 +5142,10 @@ short cio_location(double jd_tdb, enum novas_accuracy accuracy, double *ra_cio, static double t_last = 0.0, ra_last = 0.0; static ra_of_cio cio[CIO_INTERP_POINTS]; + // Default return values... + if(ra_cio) *ra_cio = NAN; + if(loc_type) *loc_type = -1; + if(!ra_cio || !loc_type) { errno = EINVAL; return -1; @@ -5032,16 +5181,19 @@ short cio_location(double jd_tdb, enum novas_accuracy accuracy, double *ra_cio, *ra_cio *= ARCSEC / HOURANGLE; *loc_type = CIO_VS_GCRS; - t_last = jd_tdb; - ra_last = *ra_cio; - ref_sys_last = *loc_type; - return 0; } + else { + // Calculate the equation of origins. + *ra_cio = -1.0 * ira_equinox(jd_tdb, NOVAS_TRUE_EQUINOX, accuracy); + *loc_type = CIO_VS_EQUINOX; + } + + t_last = jd_tdb; + acc_last = accuracy; + ra_last = *ra_cio; + ref_sys_last = *loc_type; - // Calculate the equation of origins. - *ra_cio = -1.0 * ira_equinox(jd_tdb, NOVAS_TRUE_EQUINOX, accuracy); - *loc_type = CIO_VS_EQUINOX; return 0; } @@ -5267,7 +5419,7 @@ short cio_array(double jd_tdb, long n_pts, ra_of_cio *cio) { cache_count = N; } - if(n_pts > cache_count) { + if((index_rec - index_cache) + n_pts > cache_count) { errno = EOF; return 6; // Data requested beyond file... } @@ -5350,15 +5502,11 @@ double ira_equinox(double jd_tdb, enum novas_equinox_type equinox, enum novas_ac * @since 1.0 * @author Attila Kovacs * - * @sa get_ephem_reasder() + * @sa get_ephem_provider() * @sa ephemeris() * */ int set_ephem_provider(novas_ephem_provider func) { - if(!func) { - errno = EINVAL; - return -1; - } readeph2_call = func; return 0; } @@ -5391,30 +5539,31 @@ novas_ephem_provider get_ephem_provider() { * @param[out] pos [AU] Pointer to structure containing the designation of the body of interest * @param[out] vel [AU/day] Velocity vector of the body at 'jd_tdb'; equatorial rectangular * coordinates in AU/day referred to the ICRS. - * @return 0 if successful, -1 if the input object argument is NULL, or else - * 1 if 'origin' is invalid, 2 if cel_obj->type is invalid, + * @return 0 if successful, -1 if the 'jd_tdb' or input object argument is NULL, or + * else 1 if 'origin' is invalid, 2 if cel_obj->type is invalid, * 10 + the error code from solarsystem(), or 20 + the error code from * readeph(). * - * @sa ephem_open() * @sa set_planet_provider() * @sa set_planet_provider_hp() * @sa set_ephem_provider() - * @sa make_object() + * @sa ephem_open() + * @sa make_planet() + * @sa make_ephem_object() */ -short ephemeris(const double jd_tdb[2], const object *body, enum novas_origin origin, enum novas_accuracy accuracy, +short ephemeris(const double *jd_tdb, const object *body, enum novas_origin origin, enum novas_accuracy accuracy, double *pos, double *vel) { double posvel[6] = { }; int error = 0; - if(!body) { + if(!jd_tdb || !body) { errno = EINVAL; return -1; } // Check the value of 'origin'. - if((origin < 0) || (origin >= NOVAS_ORIGIN_TYPES)) { + if(origin < 0 || origin >= NOVAS_ORIGIN_TYPES) { errno = EINVAL; return 1; } @@ -5452,7 +5601,7 @@ short ephemeris(const double jd_tdb[2], const object *body, enum novas_origin or error = -1; if(readeph2_call) { // If there is a newstyle epehemeris access routine set, we will prefer it. - error = readeph2_call(body->number, body->name, jd_tdb[0], jd_tdb[1], &eph_origin, posvel, &posvel[3]); + error = readeph2_call(body->name, body->number, jd_tdb[0], jd_tdb[1], &eph_origin, posvel, &posvel[3]); } # ifdef DEFAULT_READEPH else { @@ -5510,13 +5659,13 @@ short ephemeris(const double jd_tdb[2], const object *body, enum novas_origin or * * @return 0 if successful, or -1 if either of the input pointer arguments is NULL. * - * @sa transform_cat() * @sa make_cat_entry() + * @sa NOVAS_JD_HIP */ int transform_hip(const cat_entry *hipparcos, cat_entry *hip_2000) { cat_entry scratch; - if(!hipparcos || !hip_2000) { + if(!hipparcos) { errno = EINVAL; return -1; } @@ -5530,9 +5679,7 @@ int transform_hip(const cat_entry *hipparcos, cat_entry *hip_2000) { scratch.ra /= 15.0; // Change the epoch of the Hipparcos data from J1991.25 to J2000.0. - transform_cat(1, NOVAS_JD_HIP, &scratch, JD_J2000, "HP2", hip_2000); - - return 0; + return transform_cat(1, NOVAS_JD_HIP, &scratch, JD_J2000, "HP2", hip_2000); } /** @@ -5576,15 +5723,14 @@ int transform_hip(const cat_entry *hipparcos, cat_entry *hip_2000) { * case the catalog name is inherited from the input. * @param[out] out The transformed catalog entry, with units as given in the struct * definition - * @return 0 if successful, -1 if any of the pointer arguments is NULL or - * if the 'option' is invalid, or else 2 if 'out_id' is too long. + * @return 0 if successful, -1 if either vector argument is NULL or if the + * 'option' is invalid, or else 2 if 'out_id' is too long. * * @sa transform_hip() * @sa make_cat_entry() - * @sa precession() - * @sa frame_tie() * @sa NOVAS_JD_J2000 * @sa NOVAS_JD_B1950 + * @sa NOVAS_JD_HIP */ short transform_cat(enum novas_transform_type option, double date_in, const cat_entry *in, double date_out, const char *out_id, cat_entry *out) { @@ -5746,17 +5892,18 @@ short transform_cat(enum novas_transform_type option, double date_in, const cat_ * Determines the angle of an object above or below the Earth's limb (horizon). The geometric * limb is computed, assuming the Earth to be an airless sphere (no refraction or oblateness * is included). The observer can be on or above the Earth. For an observer on the surface - * of the Earth, this function returns the approximate unrefracted altitude. + * of the Earth, this function returns the approximate unrefracted elevation. * * @param pos_obj [AU] Position 3-vector of observed object, with respect to origin at * geocenter, components in AU. * @param pos_obs [AU] Position 3-vector of observer, with respect to origin at * geocenter, components in AU. - * @param[out] limb_ang [deg] Angle of observed object above (+) or below (-) limb in degrees. - * It may be NULL if not required. + * @param[out] limb_ang [deg] Angle of observed object above (+) or below (-) limb in degrees, + * or NAN if reurning with an error. It may be NULL if not required. * @param[out] nadir_ang [deg] Nadir angle of observed object as a fraction of apparent radius * of limb: %lt;1.0 if below the limb; 1.0 on the limb; or >1.0 if - * above the limb. It may be NULL if not required. + * above the limb. Returns NAN in case of an error return. + * It may be NULL if not required. * * @return 0 if successful, or -1 if either of the input vectors is NULL. * @@ -5766,6 +5913,8 @@ int limb_angle(const double *pos_obj, const double *pos_obs, double *limb_ang, d double disobj, disobs, aprad, zdlim, coszd, zdobj; if(!pos_obj || !pos_obs) { + if(limb_ang) *limb_ang = NAN; + if(nadir_ang) *nadir_ang = NAN; errno = EINVAL; return -1; } @@ -5820,6 +5969,7 @@ int limb_angle(const double *pos_obj, const double *pos_obs, double *limb_ang, d * or 0.0 if the location is NULL or the option is invalid. * * @sa refract() + * @sa itrs_to_hor() * * @since 1.0 * @author Attila Kovacs @@ -5828,15 +5978,6 @@ double refract_astro(const on_surface *location, enum novas_refraction_model opt double refr = 0.0; int i; - if(!location) { - errno = EINVAL; - return 0.0; - } - if(option != NOVAS_STANDARD_ATMOSPHERE && option != NOVAS_WEATHER_AT_LOCATION) { - errno = EINVAL; - return 0.0; - } - for(i = 0; i < 30; i++) { double zd_obs = zd_astro - refr; refr = refract(location, option, zd_obs); @@ -5877,9 +6018,10 @@ double refract_astro(const on_surface *location, enum novas_refraction_model opt * @param zd_obs [deg] Observed (already refracted!) zenith distance through the * atmosphere. * @return [deg] the calculated optical refraction or 0.0 if the location is - * NULL or the option is invalid. + * NULL or the option is invalid or the 'zd_obs' is invalid (<90°). * * @sa refract_astro() + * @sa hor_to_itrs() */ double refract(const on_surface *location, enum novas_refraction_model option, double zd_obs) { // 's' is the approximate scale height of atmosphere in meters. @@ -5897,7 +6039,10 @@ double refract(const on_surface *location, enum novas_refraction_model option, d } // Compute refraction only for zenith distances between -0.1 and 90.1 degrees. - if((zd_obs < -0.1) || (zd_obs > 90.1)) return 0.0; + if((zd_obs < -0.1) || (zd_obs > 90.1)) { + if(zd_obs > 0) errno = EINVAL; + return 0.0; + } // If observed weather data are available, use them. Otherwise, use // crude estimates of average conditions. @@ -6107,16 +6252,14 @@ void novas_case_sensitive(int value) { * @param[out] cel_obj Pointer to the celestial object data structure to be populated. Used * only if 'type' is NOVAS_PLANET or NOVAS_EPHEM_OBJECT, otherwise * ignored and may be NULL. - * @return 0 if successful, or else 1 if 'type' is invalid, 2 if 'number' is out - * of range, 3 if cel_obj is NULL, 4 if star_data is NULL, or 5 if 'name' - * is too long. (The return values 3 and 4 are somewhat different than in - * the original NOVAS. However they are produced in very similar contexts, - * and for similar reasons). + * @return 0 if successful, or -1 if 'cel_obj' is NULL or when type is + * NOVAS_CATALOG_OBJECT and 'star' is NULL, or else 1 if 'type' is + * invalid, 2 if 'number' is out of legal range or 5 if 'name' is too long. * * @sa novas_case_sensitive() * @sa make_cat_entry() * @sa make_planet() - * @sa make_ephem_body() + * @sa make_ephem_object() * @sa place() * */ @@ -6124,13 +6267,13 @@ short make_object(enum novas_object_type type, long number, const char *name, co if(!cel_obj) { errno = EINVAL; - return 3; + return -1; } memset(cel_obj, 0, sizeof(*cel_obj)); // Set the object type. - if((type < 0) || (type >= NOVAS_OBJECT_TYPES)) { + if(type < 0 || type >= NOVAS_OBJECT_TYPES) { errno = EINVAL; return 1; } @@ -6138,18 +6281,11 @@ short make_object(enum novas_object_type type, long number, const char *name, co // Set the object number. if(type == NOVAS_PLANET) { - if((number < 0) || (number >= NOVAS_PLANETS)) { + if(number < 0 || number >= NOVAS_PLANETS) { errno = EINVAL; return 2; } } - else if(type == NOVAS_EPHEM_OBJECT) { - if(number <= 0) { - errno = EINVAL; - return 2; - } - } - else number = 0; cel_obj->number = number; @@ -6170,7 +6306,7 @@ short make_object(enum novas_object_type type, long number, const char *name, co if(type == NOVAS_CATALOG_OBJECT) { if(!star) { errno = EINVAL; - return 4; + return -1; } cel_obj->star = *star; } @@ -6186,7 +6322,7 @@ short make_object(enum novas_object_type type, long number, const char *name, co * @param[out] planet Pointer to structure to populate. * @return 0 if successful, or else -1 if the 'planet' pointer is NULL. * - * @sa make_ephem_body() + * @sa make_ephem_object() * @sa make_cat_entry() * @sa place() * @@ -6195,13 +6331,11 @@ short make_object(enum novas_object_type type, long number, const char *name, co */ int make_planet(enum novas_planet num, object *planet) { char *names[] = NOVAS_PLANET_NAMES_INIT; - - if(!planet || num < 0 || num >= NOVAS_PLANETS) { + if(num < 0 || num >= NOVAS_PLANETS) { errno = EINVAL; return -1; } - - return make_object(NOVAS_PLANET, num, names[num], NULL, planet); + return make_object(NOVAS_PLANET, num, names[num], NULL, planet) ? -1 : 0; } /** @@ -6213,7 +6347,7 @@ int make_planet(enum novas_planet num, object *planet) { * termination. * @param num Solar-system body ID number (e.g. NAIF) * @param[out] body Pointer to structure to populate. - * @return 0 if successful, or else -1 if the 'planet' pointer is NULL, or 5 if the name + * @return 0 if successful, or else -1 if the 'planet' pointer is NULL or the name * is too long. * * @@ -6224,13 +6358,8 @@ int make_planet(enum novas_planet num, object *planet) { * @since 1.0 * @author Attila Kovacs */ -int make_ephem_body(const char *name, long num, object *body) { - if(!name || !body) { - errno = EINVAL; - return -1; - } - - return make_object(NOVAS_EPHEM_OBJECT, num, name, NULL, body); +int make_ephem_object(const char *name, long num, object *body) { + return make_object(NOVAS_EPHEM_OBJECT, num, name, NULL, body) ? -1 : 0; } diff --git a/src/nutation.c b/src/nutation.c index a85b1460..513153ee 100644 --- a/src/nutation.c +++ b/src/nutation.c @@ -15,6 +15,7 @@ */ #include +#include #include "novas.h" @@ -4215,7 +4216,6 @@ int nu2000k(double jd_tt_high, double jd_tt_low, double *dpsi, double *deps) { int i; - // Compute fundamental arguments from Simon et al. (1994), in radians. fund_args(t, &a); diff --git a/src/solsys-ephem.c b/src/solsys-ephem.c index 0847f1aa..0e134d8a 100644 --- a/src/solsys-ephem.c +++ b/src/solsys-ephem.c @@ -41,7 +41,7 @@ short planet_ephem_provider_hp(const double jd_tdb[2], enum novas_planet body, enum novas_origin origin, double *position, double *velocity) { static const char *names[] = NOVAS_PLANET_NAMES_INIT; - novas_ephem_provider f; + novas_ephem_provider ephem_call; enum novas_origin o = NOVAS_BARYCENTER; int error; @@ -55,21 +55,21 @@ short planet_ephem_provider_hp(const double jd_tdb[2], enum novas_planet body, e return -1; } - f = get_ephem_provider(); - if(!f) { + ephem_call = get_ephem_provider(); + if(!ephem_call) { errno = ENOSYS; return -1; } - error = f(body, names[body], jd_tdb[0], jd_tdb[1], &o, position, velocity); + error = ephem_call(names[body], body, jd_tdb[0], jd_tdb[1], &o, position, velocity); if(error) return 2; if(o != origin) { double pos0[3], vel0[3]; int i; - int ref = (origin == NOVAS_BARYCENTER) ? NOVAS_SSB : NOVAS_SUN; + int ref = (o == NOVAS_BARYCENTER) ? NOVAS_SUN : NOVAS_SSB; - error = f(ref, names[origin], jd_tdb[0], jd_tdb[1], &o, pos0, vel0); + error = ephem_call(names[ref], ref, jd_tdb[0], jd_tdb[1], &o, pos0, vel0); if(error) return 2; for(i = 0; i < 3; i++) { diff --git a/src/solsys2.c b/src/solsys2.c index 9e636bed..1d3dbc6b 100644 --- a/src/solsys2.c +++ b/src/solsys2.c @@ -88,10 +88,6 @@ short planet_jplint(double jd_tdb, enum novas_planet body, enum novas_origin ori errno = EINVAL; return 1; } - else if((origin < 0) || (origin >= NOVAS_ORIGIN_TYPES)) { - errno = EINVAL; - return 1; - } // Select 'targ' according to the value of 'body'. if(body == NOVAS_SUN) targ = 11L; @@ -119,8 +115,8 @@ short planet_jplint(double jd_tdb, enum novas_planet body, enum novas_origin ori // Decompose 'posvel' into 'position' and 'velocity'. for(i = 3; --i >= 0; ) { - position[i] = posvel[i]; - velocity[i] = posvel[i + 3]; + if(position) position[i] = posvel[i]; + if(velocity) velocity[i] = posvel[i + 3]; } return 0; diff --git a/src/solsys3.c b/src/solsys3.c index ebe95cb8..2ef7c4b2 100644 --- a/src/solsys3.c +++ b/src/solsys3.c @@ -26,8 +26,8 @@ /// \endcond // Additional local function prototypes -void sun_eph(double jd, double *ra, double *dec, double *dis); -void earth_sun_enable_hp(int value); +int sun_eph(double jd, double *ra, double *dec, double *dis); + /** * Whether the high-precision call is allowed to return a low-precision result. If set to 0 @@ -54,7 +54,7 @@ static int allow_lp_for_hp = 0; * * @sa earth_sun_calc_hp() */ -void earth_sun_enable_hp(int value) { +void enable_earth_sun_hp(int value) { allow_lp_for_hp = (value != 0); } @@ -102,14 +102,14 @@ short earth_sun_calc(double jd_tdb, enum novas_planet body, enum novas_origin or // Explanatory Supplement (1992), p. 316) with angles in radians. These // data are used for barycenter computations only. - static const float pm[4] = { 1047.349, 3497.898, 22903.0, 19412.2 }; - static const float pa[4] = { 5.203363, 9.537070, 19.191264, 30.068963 }; - static const float pe[4] = { 0.048393, 0.054151, 0.047168, 0.008586 }; - static const float pj[4] = { 0.022782, 0.043362, 0.013437, 0.030878 }; - static const float po[4] = { 1.755036, 1.984702, 1.295556, 2.298977 }; - static const float pw[4] = { 0.257503, 1.613242, 2.983889, 0.784898 }; - static const float pl[4] = { 0.600470, 0.871693, 5.466933, 5.321160 }; - static const float pn[4] = { 1.450138e-3, 5.841727e-4, 2.047497e-4, 1.043891e-4 }; + static const double pm[4] = { 1047.349, 3497.898, 22903.0, 19412.2 }; + static const double pa[4] = { 5.203363, 9.537070, 19.191264, 30.068963 }; + static const double pe[4] = { 0.048393, 0.054151, 0.047168, 0.008586 }; + static const double pj[4] = { 0.022782, 0.043362, 0.013437, 0.030878 }; + static const double po[4] = { 1.755036, 1.984702, 1.295556, 2.298977 }; + static const double pw[4] = { 0.257503, 1.613242, 2.983889, 0.784898 }; + static const double pl[4] = { 0.600470, 0.871693, 5.466933, 5.321160 }; + static const double pn[4] = { 1.450138e-3, 5.841727e-4, 2.047497e-4, 1.043891e-4 }; // 'obl' is the obliquity of ecliptic at epoch J2000.0 in degrees. @@ -162,7 +162,7 @@ short earth_sun_calc(double jd_tdb, enum novas_planet body, enum novas_origin or } // Check if input Julian date is within range (within 3 centuries of J2000). - if((jd_tdb < 2340000.5) || (jd_tdb > 2560000.5)) { + if(jd_tdb < 2340000.5 || jd_tdb > 2560000.5) { errno = EDOM; return 1; } @@ -252,7 +252,7 @@ short earth_sun_calc(double jd_tdb, enum novas_planet body, enum novas_origin or /** * It may provide the position and velocity of the Earth and Sun, the same as - * solarsystem_earth_sun(), if earth_sun_enable_hp() is set to true (non-zero). Otherwise, + * solarsystem_earth_sun(), if enable_earth_sun_hp() is set to true (non-zero). Otherwise, * it will return with an error code of 3, indicating that high-precision calculations are * not provided by this implementation. * @@ -284,7 +284,7 @@ short earth_sun_calc(double jd_tdb, enum novas_planet body, enum novas_origin or * out of range, 2 if 'body' is invalid, or 3 if the high-precision * orbital data cannot be produced (default return value). * - * @sa earth_sun_enable_hp() + * @sa enable_earth_sun_hp() * @sa earth_sun_calc() * @sa set_planet_provider() * @sa solarsystem_hp() @@ -323,14 +323,15 @@ short earth_sun_calc_hp(const double jd_tdb[2], enum novas_planet body, enum nov * @param[out] ra [h] Right ascension referred to mean equator and equinox of date (hours). * @param[out] dec [deg] Declination referred to mean equator and equinox of date (degrees). * @param[out] dis [AU] Geocentric distance (AU). + * @return 0 if successful, or else -1 if any of the pointer arguments are NULL. * * @sa earth_sun_calc() */ -void sun_eph(double jd, double *ra, double *dec, double *dis) { +int sun_eph(double jd, double *ra, double *dec, double *dis) { struct sun_con { int l; int r; - float alpha; + double alpha; double nu; }; @@ -393,6 +394,11 @@ void sun_eph(double jd, double *ra, double *dec, double *dis) { { 10, 0, 1.50, 21463.25 }, // { 10, -9, 2.55, 157208.40 } }; + if(!ra || !dec || !dis) { + errno = EINVAL; + return -1; + } + // Define the time units 'u', measured in units of 10000 Julian years // from J2000.0, and 't', measured in Julian centuries from J2000.0. u = (jd - T0) / 3652500.0; @@ -427,7 +433,7 @@ void sun_eph(double jd, double *ra, double *dec, double *dis) { *dec = asin(sin(emean) * sin_lon) * RAD2DEG; - return; + return 0; } #if DEFAULT_SOLSYS == 3 diff --git a/test/Makefile b/test/Makefile index dd5f9a7a..bceff5ea 100644 --- a/test/Makefile +++ b/test/Makefile @@ -5,15 +5,18 @@ LDFLAGS += -fprofile-arcs -ftest-coverage -lm CFLAGS += -DDEFAULT_SOLSYS=3 -DDEFAULT_CIO_LOCATOR_FILE=\"cio_ra.bin\" -OBJECTS := $(subst obj/,,$(OBJECTS)) +#OBJECTS := $(subst obj/,,$(OBJECTS)) +OBJECTS := novas.o nutation.o solsys3.o .PHONY: run -run: clean-cov clean-data test test_cio_file test-super data ../cio_ra.bin - ./test +run: clean-cov clean-data test-compat test-cio_file test-super test-errors data ../cio_ra.bin + ./test-compat ln -s ../cio_ra.bin - ./test_cio_file + ./test-cio_file diff data reference ./test-super + rm -f cio_ra.bin + ./test-errors cov: cov.info genhtml $< -o $@ @@ -41,24 +44,18 @@ data: ../cio_ra.bin: make -C .. cio_ra.bin -test: test.o $(OBJECTS) | Makefile data +test-%: test-%.o $(OBJECTS) | Makefile data $(CC) -o $@ $^ $(LDFLAGS) -test_cio_file: test_cio_file.o $(OBJECTS) | Makefile data - $(CC) -o $@ $^ $(LDFLAGS) - -test-super: test-super.o $(OBJECTS) | Makefile data - $(CC) -o $@ $^ $(LDFLAGS) - -test%.o: src/test%.c | Makefile - $(CC) -c -o $@ -I../include $< +test-%.o: src/test-%.c | Makefile + $(CC) -c -o $@ -g -I../include $< %.o: ../src/%.c | Makefile $(CC) -c -o $@ $(CFLAGS) $< .PHONY: clean-test clean-test: - rm -f test test_cio_file test-super + rm -f test-* .PHONY: clean-cov clean-cov: diff --git a/test/novasc3.1/aberration.out b/test/novasc3.1/aberration.out new file mode 100644 index 00000000..5df27789 --- /dev/null +++ b/test/novasc3.1/aberration.out @@ -0,0 +1,120 @@ + +-10000.0 22+20 S2 O0 A0: 0.813836248 -0.469762512 0.342043481 0 0 0 + +-10000.0 22+20 S2 O1 A0: 0.813874074 -0.469679861 0.342066983 0 0 0 + +-10000.0 22+20 S2 O2 A0: 0.832499420 -0.447839966 0.326165725 0 0 0 + +-10000.0 16-20 S2 O0 A0: -0.469833971 -0.813804082 -0.342021863 0 0 0 + +-10000.0 16-20 S2 O1 A0: -0.469822163 -0.813810361 -0.342023146 0 0 0 + +-10000.0 16-20 S2 O2 A0: -0.423555119 -0.835114441 -0.350977111 0 0 0 + + 0.0 22+20 S2 O0 A0: 0.813759734 -0.469898540 0.342038678 0 0 0 + + 0.0 22+20 S2 O1 A0: 0.813722227 -0.469950306 0.342056789 0 0 0 + + 0.0 22+20 S2 O2 A0: 0.832358383 -0.448107854 0.326157743 0 0 0 + + 0.0 16-20 S2 O0 A0: -0.469916143 -0.813763333 -0.342005931 0 0 0 + + 0.0 16-20 S2 O1 A0: -0.469984792 -0.813729576 -0.341991918 0 0 0 + + 0.0 16-20 S2 O2 A0: -0.423734139 -0.835036390 -0.350946729 0 0 0 + + 10000.0 22+20 S2 O0 A0: 0.813812821 -0.469856085 0.341970689 0 0 0 + + 10000.0 22+20 S2 O1 A0: 0.813828043 -0.469865384 0.341921684 0 0 0 + + 10000.0 22+20 S2 O2 A0: 0.832454726 -0.448025345 0.326025183 0 0 0 + + 10000.0 16-20 S2 O0 A0: -0.469759945 -0.813839824 -0.342038499 0 0 0 + + 10000.0 16-20 S2 O1 A0: -0.469674751 -0.813881219 -0.342056997 0 0 0 + + 10000.0 16-20 S2 O2 A0: -0.423406001 -0.835177534 -0.351006902 0 0 0 + + 10000.0 22+20 S2 O0 A0: 0.813812821 -0.469856085 0.341970689 0 0 0 + + 10000.0 22+20 S2 O1 A0: 0.813828043 -0.469865384 0.341921684 0 0 0 + + 10000.0 22+20 S2 O2 A0: 0.832454726 -0.448025345 0.326025183 0 0 0 + + 10000.0 16-20 S2 O0 A0: -0.469759945 -0.813839824 -0.342038499 0 0 0 + + 10000.0 16-20 S2 O1 A0: -0.469674751 -0.813881219 -0.342056997 0 0 0 + + 10000.0 16-20 S2 O2 A0: -0.423406001 -0.835177534 -0.351006902 0 0 0 + + 10000.0 22+20 S2 O0 A0: 0.813812827 -0.469856070 0.341970694 0 0 0 + + 10000.0 22+20 S2 O1 A0: 0.813828009 -0.469865430 0.341921701 0 0 0 + + 10000.0 22+20 S2 O2 A0: 0.832454738 -0.448025317 0.326025191 0 0 0 + + 10000.0 16-20 S2 O0 A0: -0.469759944 -0.813839824 -0.342038499 0 0 0 + + 10000.0 16-20 S2 O1 A0: -0.469674776 -0.813881217 -0.342056969 0 0 0 + + 10000.0 16-20 S2 O2 A0: -0.423405998 -0.835177535 -0.351006902 0 0 0 + +-10000.0 22+20 S2 O0 A1: 0.813836 -0.469763 0.342043 0 0 0 + +-10000.0 22+20 S2 O1 A1: 0.813874 -0.469680 0.342067 0 0 0 + +-10000.0 22+20 S2 O2 A1: 0.832499 -0.447840 0.326166 0 0 0 + +-10000.0 16-20 S2 O0 A1: -0.469834 -0.813804 -0.342022 0 0 0 + +-10000.0 16-20 S2 O1 A1: -0.469822 -0.813810 -0.342023 0 0 0 + +-10000.0 16-20 S2 O2 A1: -0.423555 -0.835114 -0.350977 0 0 0 + + 0.0 22+20 S2 O0 A1: 0.813760 -0.469899 0.342039 0 0 0 + + 0.0 22+20 S2 O1 A1: 0.813722 -0.469950 0.342057 0 0 0 + + 0.0 22+20 S2 O2 A1: 0.832358 -0.448108 0.326158 0 0 0 + + 0.0 16-20 S2 O0 A1: -0.469916 -0.813763 -0.342006 0 0 0 + + 0.0 16-20 S2 O1 A1: -0.469985 -0.813730 -0.341992 0 0 0 + + 0.0 16-20 S2 O2 A1: -0.423734 -0.835036 -0.350947 0 0 0 + + 10000.0 22+20 S2 O0 A1: 0.813813 -0.469856 0.341971 0 0 0 + + 10000.0 22+20 S2 O1 A1: 0.813828 -0.469865 0.341922 0 0 0 + + 10000.0 22+20 S2 O2 A1: 0.832455 -0.448025 0.326025 0 0 0 + + 10000.0 16-20 S2 O0 A1: -0.469760 -0.813840 -0.342038 0 0 0 + + 10000.0 16-20 S2 O1 A1: -0.469675 -0.813881 -0.342057 0 0 0 + + 10000.0 16-20 S2 O2 A1: -0.423406 -0.835178 -0.351007 0 0 0 + + 10000.0 22+20 S2 O0 A1: 0.813813 -0.469856 0.341971 0 0 0 + + 10000.0 22+20 S2 O1 A1: 0.813828 -0.469865 0.341922 0 0 0 + + 10000.0 22+20 S2 O2 A1: 0.832455 -0.448025 0.326025 0 0 0 + + 10000.0 16-20 S2 O0 A1: -0.469760 -0.813840 -0.342038 0 0 0 + + 10000.0 16-20 S2 O1 A1: -0.469675 -0.813881 -0.342057 0 0 0 + + 10000.0 16-20 S2 O2 A1: -0.423406 -0.835178 -0.351007 0 0 0 + + 10000.0 22+20 S2 O0 A1: 0.813813 -0.469856 0.341971 0 0 0 + + 10000.0 22+20 S2 O1 A1: 0.813828 -0.469865 0.341922 0 0 0 + + 10000.0 22+20 S2 O2 A1: 0.832455 -0.448025 0.326025 0 0 0 + + 10000.0 16-20 S2 O0 A1: -0.469760 -0.813840 -0.342038 0 0 0 + + 10000.0 16-20 S2 O1 A1: -0.469675 -0.813881 -0.342057 0 0 0 + + 10000.0 16-20 S2 O2 A1: -0.423406 -0.835178 -0.351007 0 0 0 diff --git a/test/novasc3.1/ephemeris.out b/test/novasc3.1/ephemeris.out index a30130a8..4a0cde13 100644 --- a/test/novasc3.1/ephemeris.out +++ b/test/novasc3.1/ephemeris.out @@ -18,43 +18,43 @@ 10000.010 A0: ERROR 13 10000.010 A0: ERROR 13 10000.010 A0: ERROR 13 - -10000.000 A1: SUN 3 -0.000009 0.003773 0.001651 -0.009747 -0.002548 -0.000896 + -10000.000 A1: SUN 0 -0.000009 0.003773 0.001651 -0.009747 -0.002548 -0.000896 - -10000.000 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + -10000.000 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - -10000.000 A1: EARTH 3 0.810058 -0.553762 -0.240106 17.374036 21.766877 9.439017 + -10000.000 A1: EARTH 0 0.810058 -0.553762 -0.240106 17.374036 21.766877 9.439017 - -10000.000 A1: EARTH 3 0.810068 -0.557535 -0.241757 17.383783 21.769425 9.439914 + -10000.000 A1: EARTH 1 0.810068 -0.557535 -0.241757 17.383783 21.769425 9.439914 - 0.000 A1: SUN 3 -0.007138 -0.002644 -0.000921 0.009234 -0.011779 -0.005287 + 0.000 A1: SUN 0 -0.007138 -0.002644 -0.000921 0.009234 -0.011779 -0.005287 - 0.000 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + 0.000 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - 0.000 A1: EARTH 3 -0.184271 0.884782 0.383824 -29.784950 -5.028997 -2.180527 + 0.000 A1: EARTH 0 -0.184271 0.884782 0.383824 -29.784950 -5.028997 -2.180527 - 0.000 A1: EARTH 3 -0.177134 0.887425 0.384746 -29.794184 -5.017218 -2.175240 + 0.000 A1: EARTH 1 -0.177134 0.887425 0.384746 -29.794184 -5.017218 -2.175240 - 10000.000 A1: SUN 3 -0.000010 -0.004018 -0.001669 0.007845 0.006630 0.002675 + 10000.000 A1: SUN 0 -0.000010 -0.004018 -0.001669 0.007845 0.006630 0.002675 - 10000.000 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + 10000.000 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - 10000.000 A1: EARTH 3 -0.537235 -0.790493 -0.342586 24.755316 -14.602263 -6.330115 + 10000.000 A1: EARTH 0 -0.537235 -0.790493 -0.342586 24.755316 -14.602263 -6.330115 - 10000.000 A1: EARTH 3 -0.537226 -0.786475 -0.340917 24.747471 -14.608892 -6.332790 + 10000.000 A1: EARTH 1 -0.537226 -0.786475 -0.340917 24.747471 -14.608892 -6.332790 - 10000.000 A1: SUN 3 -0.000010 -0.004018 -0.001669 0.007845 0.006630 0.002675 + 10000.000 A1: SUN 0 -0.000010 -0.004018 -0.001669 0.007845 0.006630 0.002675 - 10000.000 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + 10000.000 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - 10000.000 A1: EARTH 3 -0.537235 -0.790493 -0.342586 24.755316 -14.602263 -6.330115 + 10000.000 A1: EARTH 0 -0.537235 -0.790493 -0.342586 24.755316 -14.602263 -6.330115 - 10000.000 A1: EARTH 3 -0.537226 -0.786475 -0.340917 24.747471 -14.608892 -6.332790 + 10000.000 A1: EARTH 1 -0.537226 -0.786475 -0.340917 24.747471 -14.608892 -6.332790 - 10000.010 A1: SUN 3 -0.000010 -0.004018 -0.001669 0.007844 0.006630 0.002675 + 10000.010 A1: SUN 0 -0.000010 -0.004018 -0.001669 0.007844 0.006630 0.002675 - 10000.010 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + 10000.010 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - 10000.010 A1: EARTH 3 -0.537093 -0.790577 -0.342623 24.757953 -14.598388 -6.328435 + 10000.010 A1: EARTH 0 -0.537093 -0.790577 -0.342623 24.757953 -14.598388 -6.328435 - 10000.010 A1: EARTH 3 -0.537083 -0.786560 -0.340954 24.750109 -14.605018 -6.331111 + 10000.010 A1: EARTH 1 -0.537083 -0.786560 -0.340954 24.750109 -14.605018 -6.331111 diff --git a/test/reference/aberration.out b/test/reference/aberration.out new file mode 100644 index 00000000..e0c3eded --- /dev/null +++ b/test/reference/aberration.out @@ -0,0 +1,120 @@ + +-10000.0 22+20 S2 O0 A0: 0.813836248 -0.469762512 0.342043481 1 1 1 + +-10000.0 22+20 S2 O1 A0: 0.813874074 -0.469679861 0.342066983 1 1 1 + +-10000.0 22+20 S2 O2 A0: 0.832499420 -0.447839966 0.326165725 1 1 1 + +-10000.0 16-20 S2 O0 A0: -0.469833971 -0.813804082 -0.342021863 1 1 1 + +-10000.0 16-20 S2 O1 A0: -0.469822163 -0.813810361 -0.342023146 1 1 1 + +-10000.0 16-20 S2 O2 A0: -0.423555119 -0.835114441 -0.350977111 1 1 1 + + 0.0 22+20 S2 O0 A0: 0.813759734 -0.469898540 0.342038678 1 1 1 + + 0.0 22+20 S2 O1 A0: 0.813722227 -0.469950306 0.342056789 1 1 1 + + 0.0 22+20 S2 O2 A0: 0.832358383 -0.448107854 0.326157743 1 1 1 + + 0.0 16-20 S2 O0 A0: -0.469916143 -0.813763333 -0.342005931 1 1 1 + + 0.0 16-20 S2 O1 A0: -0.469984792 -0.813729576 -0.341991918 1 1 1 + + 0.0 16-20 S2 O2 A0: -0.423734139 -0.835036390 -0.350946729 1 1 1 + + 10000.0 22+20 S2 O0 A0: 0.813812821 -0.469856085 0.341970689 1 1 1 + + 10000.0 22+20 S2 O1 A0: 0.813828043 -0.469865384 0.341921684 1 1 1 + + 10000.0 22+20 S2 O2 A0: 0.832454726 -0.448025345 0.326025183 1 1 1 + + 10000.0 16-20 S2 O0 A0: -0.469759945 -0.813839824 -0.342038499 1 1 1 + + 10000.0 16-20 S2 O1 A0: -0.469674751 -0.813881219 -0.342056997 1 1 1 + + 10000.0 16-20 S2 O2 A0: -0.423406001 -0.835177534 -0.351006902 1 1 1 + + 10000.0 22+20 S2 O0 A0: 0.813812821 -0.469856085 0.341970689 1 1 1 + + 10000.0 22+20 S2 O1 A0: 0.813828043 -0.469865384 0.341921684 1 1 1 + + 10000.0 22+20 S2 O2 A0: 0.832454726 -0.448025345 0.326025183 1 1 1 + + 10000.0 16-20 S2 O0 A0: -0.469759945 -0.813839824 -0.342038499 1 1 1 + + 10000.0 16-20 S2 O1 A0: -0.469674751 -0.813881219 -0.342056997 1 1 1 + + 10000.0 16-20 S2 O2 A0: -0.423406001 -0.835177534 -0.351006902 1 1 1 + + 10000.0 22+20 S2 O0 A0: 0.813812827 -0.469856070 0.341970694 1 1 1 + + 10000.0 22+20 S2 O1 A0: 0.813828009 -0.469865430 0.341921701 1 1 1 + + 10000.0 22+20 S2 O2 A0: 0.832454738 -0.448025317 0.326025191 1 1 1 + + 10000.0 16-20 S2 O0 A0: -0.469759944 -0.813839824 -0.342038499 1 1 1 + + 10000.0 16-20 S2 O1 A0: -0.469674776 -0.813881217 -0.342056969 1 1 1 + + 10000.0 16-20 S2 O2 A0: -0.423405998 -0.835177535 -0.351006902 1 1 1 + +-10000.0 22+20 S2 O0 A1: 0.813836 -0.469763 0.342043 1 1 1 + +-10000.0 22+20 S2 O1 A1: 0.813874 -0.469680 0.342067 1 1 1 + +-10000.0 22+20 S2 O2 A1: 0.832499 -0.447840 0.326166 1 1 1 + +-10000.0 16-20 S2 O0 A1: -0.469834 -0.813804 -0.342022 1 1 1 + +-10000.0 16-20 S2 O1 A1: -0.469822 -0.813810 -0.342023 1 1 1 + +-10000.0 16-20 S2 O2 A1: -0.423555 -0.835114 -0.350977 1 1 1 + + 0.0 22+20 S2 O0 A1: 0.813760 -0.469899 0.342039 1 1 1 + + 0.0 22+20 S2 O1 A1: 0.813722 -0.469950 0.342057 1 1 1 + + 0.0 22+20 S2 O2 A1: 0.832358 -0.448108 0.326158 1 1 1 + + 0.0 16-20 S2 O0 A1: -0.469916 -0.813763 -0.342006 1 1 1 + + 0.0 16-20 S2 O1 A1: -0.469985 -0.813730 -0.341992 1 1 1 + + 0.0 16-20 S2 O2 A1: -0.423734 -0.835036 -0.350947 1 1 1 + + 10000.0 22+20 S2 O0 A1: 0.813813 -0.469856 0.341971 1 1 1 + + 10000.0 22+20 S2 O1 A1: 0.813828 -0.469865 0.341922 1 1 1 + + 10000.0 22+20 S2 O2 A1: 0.832455 -0.448025 0.326025 1 1 1 + + 10000.0 16-20 S2 O0 A1: -0.469760 -0.813840 -0.342038 1 1 1 + + 10000.0 16-20 S2 O1 A1: -0.469675 -0.813881 -0.342057 1 1 1 + + 10000.0 16-20 S2 O2 A1: -0.423406 -0.835178 -0.351007 1 1 1 + + 10000.0 22+20 S2 O0 A1: 0.813813 -0.469856 0.341971 1 1 1 + + 10000.0 22+20 S2 O1 A1: 0.813828 -0.469865 0.341922 1 1 1 + + 10000.0 22+20 S2 O2 A1: 0.832455 -0.448025 0.326025 1 1 1 + + 10000.0 16-20 S2 O0 A1: -0.469760 -0.813840 -0.342038 1 1 1 + + 10000.0 16-20 S2 O1 A1: -0.469675 -0.813881 -0.342057 1 1 1 + + 10000.0 16-20 S2 O2 A1: -0.423406 -0.835178 -0.351007 1 1 1 + + 10000.0 22+20 S2 O0 A1: 0.813813 -0.469856 0.341971 1 1 1 + + 10000.0 22+20 S2 O1 A1: 0.813828 -0.469865 0.341922 1 1 1 + + 10000.0 22+20 S2 O2 A1: 0.832455 -0.448025 0.326025 1 1 1 + + 10000.0 16-20 S2 O0 A1: -0.469760 -0.813840 -0.342038 1 1 1 + + 10000.0 16-20 S2 O1 A1: -0.469675 -0.813881 -0.342057 1 1 1 + + 10000.0 16-20 S2 O2 A1: -0.423406 -0.835178 -0.351007 1 1 1 diff --git a/test/reference/ephemeris.out b/test/reference/ephemeris.out index 16ccbefe..4a0cde13 100644 --- a/test/reference/ephemeris.out +++ b/test/reference/ephemeris.out @@ -18,43 +18,43 @@ 10000.010 A0: ERROR 13 10000.010 A0: ERROR 13 10000.010 A0: ERROR 13 - -10000.000 A1: SUN 3 -0.000009 0.003773 0.001651 -0.009747 -0.002548 -0.000896 + -10000.000 A1: SUN 0 -0.000009 0.003773 0.001651 -0.009747 -0.002548 -0.000896 - -10000.000 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + -10000.000 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - -10000.000 A1: EARTH 3 0.810058 -0.553762 -0.240106 17.374036 21.766877 9.439017 + -10000.000 A1: EARTH 0 0.810058 -0.553762 -0.240106 17.374036 21.766877 9.439017 - -10000.000 A1: EARTH 3 0.810068 -0.557535 -0.241757 17.383783 21.769425 9.439914 + -10000.000 A1: EARTH 1 0.810068 -0.557535 -0.241757 17.383783 21.769425 9.439914 - 0.000 A1: SUN 3 -0.007138 -0.002644 -0.000921 0.009234 -0.011779 -0.005287 + 0.000 A1: SUN 0 -0.007138 -0.002644 -0.000921 0.009234 -0.011779 -0.005287 - 0.000 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + 0.000 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - 0.000 A1: EARTH 3 -0.184271 0.884782 0.383824 -29.784950 -5.028997 -2.180527 + 0.000 A1: EARTH 0 -0.184271 0.884782 0.383824 -29.784950 -5.028997 -2.180527 - 0.000 A1: EARTH 3 -0.177134 0.887425 0.384746 -29.794184 -5.017218 -2.175240 + 0.000 A1: EARTH 1 -0.177134 0.887425 0.384746 -29.794184 -5.017218 -2.175240 - 10000.000 A1: SUN 3 -0.000010 -0.004018 -0.001669 0.007845 0.006630 0.002675 + 10000.000 A1: SUN 0 -0.000010 -0.004018 -0.001669 0.007845 0.006630 0.002675 - 10000.000 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + 10000.000 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - 10000.000 A1: EARTH 3 -0.537235 -0.790493 -0.342586 24.755316 -14.602263 -6.330115 + 10000.000 A1: EARTH 0 -0.537235 -0.790493 -0.342586 24.755316 -14.602263 -6.330115 - 10000.000 A1: EARTH 3 -0.537226 -0.786475 -0.340917 24.747471 -14.608892 -6.332790 + 10000.000 A1: EARTH 1 -0.537226 -0.786475 -0.340917 24.747471 -14.608892 -6.332790 - 10000.000 A1: SUN 3 -0.000010 -0.004018 -0.001669 0.007845 0.006630 0.002675 + 10000.000 A1: SUN 0 -0.000010 -0.004018 -0.001669 0.007845 0.006630 0.002675 - 10000.000 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + 10000.000 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - 10000.000 A1: EARTH 3 -0.537235 -0.790493 -0.342586 24.755316 -14.602263 -6.330115 + 10000.000 A1: EARTH 0 -0.537235 -0.790493 -0.342586 24.755316 -14.602263 -6.330115 - 10000.000 A1: EARTH 3 -0.537226 -0.786475 -0.340917 24.747471 -14.608892 -6.332790 + 10000.000 A1: EARTH 1 -0.537226 -0.786475 -0.340917 24.747471 -14.608892 -6.332790 - 10000.010 A1: SUN 3 -0.000010 -0.004018 -0.001669 0.007844 0.006630 0.002675 + 10000.010 A1: SUN 0 -0.000010 -0.004018 -0.001669 0.007844 0.006630 0.002675 - 10000.010 A1: SUN 3 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 + 10000.010 A1: SUN 1 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 - 10000.010 A1: EARTH 3 -0.537093 -0.790577 -0.342623 24.757953 -14.598388 -6.328435 + 10000.010 A1: EARTH 0 -0.537093 -0.790577 -0.342623 24.757953 -14.598388 -6.328435 - 10000.010 A1: EARTH 3 -0.537083 -0.786560 -0.340954 24.750108 -14.605018 -6.331111 + 10000.010 A1: EARTH 1 -0.537083 -0.786560 -0.340954 24.750109 -14.605018 -6.331111 diff --git a/test/reference/make_object.out b/test/reference/make_object.out index 007df83b..bc639fef 100644 --- a/test/reference/make_object.out +++ b/test/reference/make_object.out @@ -1,3 +1,3 @@ EARTH 0 3 IO 1 501 -TEST 2 0 +TEST 2 1234567890 diff --git a/test/reference/nutation.out b/test/reference/nutation.out index ef4886a9..04162fda 100644 --- a/test/reference/nutation.out +++ b/test/reference/nutation.out @@ -1,8 +1,8 @@ -10000.000 A0: -nan -nan -nan - 0.000 A0: 0.801841322 -0.548228726 -0.237688365 - 10000.000 A0: -0.187719744 0.901076198 0.390925803 - 10000.000 A0: -0.529087777 -0.778583164 -0.337452783 - 10000.010 A0: -0.529087776 -0.778583166 -0.337452781 + 0.000 A0: 0.801841323 -0.548228725 -0.237688364 + 10000.000 A0: -0.187719741 0.901076198 0.390925803 + 10000.000 A0: -0.529087776 -0.778583165 -0.337452783 + 10000.010 A0: -0.529087775 -0.778583167 -0.337452781 -10000.000 A1: -0.528922 -0.778682 -0.337485 0.000 A1: 0.801841 -0.548229 -0.237688 10000.000 A1: -0.187720 0.901076 0.390926 diff --git a/test/src/test_cio_file.c b/test/src/test-cio_file.c similarity index 100% rename from test/src/test_cio_file.c rename to test/src/test-cio_file.c diff --git a/test/src/test.c b/test/src/test-compat.c similarity index 97% rename from test/src/test.c rename to test/src/test-compat.c index 9d150055..6afecbcc 100644 --- a/test/src/test.c +++ b/test/src/test-compat.c @@ -33,7 +33,10 @@ static int idx = -1; static char *header; -// cio_array + +static double vlen(double *pos) { + return sqrt(pos[0] * pos[0] + pos[1] * pos[1] + pos[2] * pos[2]); +} static void newline() { fprintf(fp, "\n%8.1f %-10s S%d O%d A%d: ", (tdb - J2000), source.name, source.type, obs.where, accuracy); @@ -270,8 +273,8 @@ static void test_ephemeris() { openfile("ephemeris"); if(is_ok(ephemeris(tdb2, &body[i], j, accuracy, pos1, vel1))) { - int j; - for(j = 0; j < 3; j++) vel1[j] *= 1e-3 * (1.4959787069098932e+11 / 86400.0); + int k; + for(k = 0; k < 3; k++) vel1[k] *= 1e-3 * (1.4959787069098932e+11 / 86400.0); fprintf(fp, "%-10s %d ", body[i].name, j); printvector(pos1); @@ -635,6 +638,23 @@ static void test_grav_def() { } } +static void test_aberration() { + double vo[3] = {}, v0[3] = {}, pos1[3]; + int i; + + // Calculate for sidereal sources only. + if (source.type != 2) return; + + for(i = 0; i < 3; i++) vo[i] = evel[i] + vobs[i]; + + openfile("aberration"); + aberration(pos0, vo, 0.0, pos1); + printunitvector(pos1); + + aberration(pos0, v0, 0.0, pos1); + for(i = 0; i < 3; i++) fprintf(fp, "%d ", fabs(pos0[i] - pos1[i]) < 1e-9 * vlen(pos0)); +} + static void test_place() { int i; @@ -883,6 +903,7 @@ static int test_source() { test_light_time(); test_grav_def(); test_place(); + test_aberration(); if(obs.where == 0) { test_astro_place(); diff --git a/test/src/test-errors.c b/test/src/test-errors.c new file mode 100644 index 00000000..ba74a8e6 --- /dev/null +++ b/test/src/test-errors.c @@ -0,0 +1,871 @@ +/** + * @file + * + * @date Created on Feb 19, 2024 + * @author Attila Kovacs + */ + + +#include +#include +#include +#include +#include + +#include "novas.h" + +static int check(const char *func, int exp, int error) { + if(error != exp) { + fprintf(stderr, "ERROR! %s: expected %d, got %d\n", func, exp, error); + return 1; + } + return 0; +} + +static int test_make_on_surface() { + on_surface loc; + if(check("make_on_surface", -1, make_on_surface(0.0, 0.0, 0.0, 0.0, 0.0, NULL))) return 1; + return 0; +} + +static int test_make_in_space() { + double p[3], v[3]; + in_space sp; + int n = 0; + + if(check("make_in_space", -1, make_in_space(p, v, NULL))) n++; + if(check("make_in_space:p", 0, make_in_space(NULL, v, &sp))) n++; + if(check("make_in_space:v", 0, make_in_space(p, NULL, &sp))) n++; + + return n; +} + +static int test_make_observer() { + in_space sp; + on_surface on; + observer obs; + int n = 0; + + if(check("make_observer:where", 1, make_observer(-1, &on, &sp, &obs))) n++; + if(check("make_observer", -1, make_observer(NOVAS_OBSERVER_AT_GEOCENTER, &on, &sp, NULL))) n++; + if(check("make_observer:on", -1, make_observer(NOVAS_OBSERVER_ON_EARTH, NULL, &sp, &obs))) n++; + if(check("make_observer:sp", -1, make_observer(NOVAS_OBSERVER_IN_EARTH_ORBIT, &on, NULL, &obs))) n++; + + return n; +} + +static int test_make_ephem_object() { + object o; + int n = 0; + + char longname[SIZE_OF_OBJ_NAME + 1]; + memset(longname, 'A', SIZE_OF_OBJ_NAME); + + if(check("make_ephem_object", -1, make_ephem_object("dummy", 1, NULL))) n++; + if(check("make_ephem_object:name", -1, make_ephem_object(longname, 1, &o))) n++; + + return n; +} + +static int test_make_planet() { + object o; + int n = 0; + + if(check("make_ephem_object:lo", -1, make_planet(-1, &o))) n++; + if(check("make_ephem_object:hi", -1, make_planet(NOVAS_PLANETS, &o))) n++; + + return n; +} + +static int test_make_cat_entry() { + cat_entry c; + int n = 0; + + char longname[SIZE_OF_OBJ_NAME + 1]; + char longcat[SIZE_OF_CAT_NAME + 1]; + + memset(longname, 'A', SIZE_OF_OBJ_NAME); + memset(longcat, 'A', SIZE_OF_CAT_NAME); + + if(check("make_cat_entry", -1, make_cat_entry("dummy", "cat", 1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, NULL))) n++; + if(check("make_cat_entry:name", 1, make_cat_entry(longname, "cat", 1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &c))) n++; + if(check("make_cat_entry:catname", 2, make_cat_entry("dummy", longcat, 1, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &c))) n++; + + return n; +} + +static int test_make_object() { + cat_entry s; + object o; + int n = 0; + + char longname[SIZE_OF_OBJ_NAME + 1]; + memset(longname, 'A', SIZE_OF_OBJ_NAME); + + if(check("make_object", -1, make_object(NOVAS_PLANET, 1, "dummy", &s, NULL))) n++; + if(check("make_object:star", -1, make_object(NOVAS_CATALOG_OBJECT, 1, "dummy", NULL, &o))) n++; + if(check("make_object:type", 1, make_object(-1, 1, "dummy", &s, &o))) n++; + if(check("make_object:pl:lo", 2, make_object(NOVAS_PLANET, -1, "dummy", NULL, &o))) n++; + if(check("make_object:pl:lo", 2, make_object(NOVAS_PLANET, NOVAS_PLANETS, "dummy", NULL, &o))) n++; + if(check("make_object:name", 5, make_object(NOVAS_PLANET, 1, longname, NULL, &o))) n++; + + return n; +} + +static int test_refract() { + on_surface o; + int n = 0; + + errno = 0; + double r = refract(NULL, NOVAS_STANDARD_ATMOSPHERE, 30.0); + if(check("refract_loc", 1, r == 0.0 && errno == EINVAL)) n++; + + errno = 0; + r = refract(&o, -1, 30.0); + if(check("refract_loc", 1, r == 0.0 && errno == EINVAL)) n++; + + errno = 0; + r = refract(&o, NOVAS_STANDARD_ATMOSPHERE, 90.11); + if(check("refract_loc", 1, r == 0.0 && errno == EINVAL)) n++; + + return n; +} + +static int test_limb_angle() { + double pos[3], a, b; + int n = 0; + + if(check("limb_angle:pos_obj", -1, limb_angle(NULL, pos, &a, &b))) n++; + if(check("limb_angle:pos_obs", -1, limb_angle(pos, NULL, &a, &b))) n++; + + return n; +} + + +static int test_transform_cat() { + cat_entry c, c1; + int n = 0; + + char longname[SIZE_OF_OBJ_NAME + 1]; + memset(longname, 'A', SIZE_OF_OBJ_NAME); + + if(check("transform_cat:in", -1, transform_cat(PRECESSION, NOVAS_JD_B1950, NULL, NOVAS_JD_J2000, "FK5", &c))) n++; + if(check("transform_cat:out", -1, transform_cat(PRECESSION, NOVAS_JD_B1950, &c, NOVAS_JD_J2000, "FK5", NULL))) n++; + if(check("transform_cat:option", -1, transform_cat(-1, NOVAS_JD_B1950, &c, NOVAS_JD_J2000, "FK5", &c1))) n++; + if(check("transform_cat:name", 2, transform_cat(PRECESSION, NOVAS_JD_B1950, &c, NOVAS_JD_J2000, longname, &c))) n++; + + return n; +} + +static int test_transform_hip() { + cat_entry c; + int n = 0; + + if(check("transform_hip:in", -1, transform_hip(NULL, &c))) n++; + if(check("transform_hip:in", -1, transform_hip(&c, NULL))) n++; + + return n; +} + +static int test_ephemeris() { + double p[3], v[3]; + double tdb[2] = { NOVAS_JD_J2000 }; + object mars; + int n = 0; + + if(check("ephemeris:body", -1, ephemeris(tdb, NULL, NOVAS_BARYCENTER, NOVAS_FULL_ACCURACY, p, v))) n++; + if(check("ephemeris:jd", -1, ephemeris(NULL, &mars, NOVAS_BARYCENTER, NOVAS_FULL_ACCURACY, p, v))) n++; + if(check("ephemeris:origin", 1, ephemeris(tdb, &mars, -1, NOVAS_FULL_ACCURACY, p, v))) n++; + + return n; +} + +static int test_j2000_to_tod() { + double p[3]; + int n = 0; + + if(check("j2000_to_tod:in", -1, j2000_to_tod(0.0, NOVAS_FULL_ACCURACY, NULL, p))) n++; + if(check("j2000_to_tod:out", -1, j2000_to_tod(0.0, NOVAS_FULL_ACCURACY, p, NULL))) n++; + if(check("j2000_to_tod:accuracy", -1, j2000_to_tod(0.0, -1, p, p))) n++; + + return n; +} + +static int test_tod_to_j2000() { + double p[3]; + int n = 0; + + if(check("tod_to_j2000:in", -1, tod_to_j2000(0.0, NOVAS_FULL_ACCURACY, NULL, p))) n++; + if(check("tod_to_j2000:out", -1, tod_to_j2000(0.0, NOVAS_FULL_ACCURACY, p, NULL))) n++; + if(check("tod_to_j2000:accuracy", -1, tod_to_j2000(0.0, -1, p, p))) n++; + + return n; +} + +static int test_gcrs_to_cirs() { + double p[3]; + int n = 0; + + if(check("gcrs_to_cirs:in", -1, gcrs_to_cirs(0.0, NOVAS_FULL_ACCURACY, NULL, p))) n++; + if(check("gcrs_to_cirs:out", -1, gcrs_to_cirs(0.0, NOVAS_FULL_ACCURACY, p, NULL))) n++; + if(check("gcrs_to_cirs:accuracy", -1, gcrs_to_cirs(0.0, -1, p, p))) n++; + + return n; +} + +static int test_cirs_to_gcrs() { + double p[3]; + int n = 0; + + if(check("cirs_to_gcrs:in", -1, cirs_to_gcrs(0.0, NOVAS_FULL_ACCURACY, NULL, p))) n++; + if(check("cirs_to_gcrs:out", -1, cirs_to_gcrs(0.0, NOVAS_FULL_ACCURACY, p, NULL))) n++; + if(check("cirs_to_gcrs:accuracy", -1, cirs_to_gcrs(0.0, -1, p, p))) n++; + + return n; +} + +static int test_set_planet_provider() { + if(check("set_planet_provider", -1, set_planet_provider(NULL))) return 1; + return 0; +} + +static int test_set_planet_provider_hp() { + if(check("set_planet_provider_hp", -1, set_planet_provider_hp(NULL))) return 1; + return 0; +} + +static int test_place_star() { + cat_entry c; + observer loc; + sky_pos pos; + int n = 0; + + if(check("place_star:in", -1, place_star(0.0, NULL, &loc, 0.0, NOVAS_GCRS, NOVAS_FULL_ACCURACY, &pos))) n++; + if(check("place_star:out", -1, place_star(0.0, &c, &loc, 0.0, NOVAS_GCRS, NOVAS_FULL_ACCURACY, NULL))) n++; + + return n; +} + +static int test_place() { + object o; + observer loc; + sky_pos pos; + int n = 0; + + if(check("place:object", -1, place(0.0, NULL, &loc, 0.0, NOVAS_GCRS, NOVAS_FULL_ACCURACY, &pos))) n++; + if(check("place:sys:lo", 1, place(0.0, &o, &loc, 0.0, -1, NOVAS_FULL_ACCURACY, &pos))) n++; + if(check("place:sys:hi", 1, place(0.0, &o, &loc, 0.0, NOVAS_REFERENCE_SYSTEMS, NOVAS_FULL_ACCURACY, &pos))) n++; + if(check("place:accuracy", 2, place(0.0, &o, &loc, 0.0, NOVAS_GCRS, -1, &pos))) n++; + + return n; +} + +static int test_radec_planet() { + object o; + observer loc; + double ra, dec, dis, rv; + + o.type = NOVAS_CATALOG_OBJECT; + + if(check("radec_planet", -1, radec_planet(0.0, &o, &loc, 0.0, NOVAS_GCRS, NOVAS_FULL_ACCURACY, &ra, &dec, &dis, &rv))) return 1; + return 0; +} + +static int test_mean_star() { + double x; + int n = 0; + + if(check("mean_star:ira", -1, mean_star(0.0, 0.0, 0.0, NOVAS_FULL_ACCURACY, NULL, &x))) n++; + if(check("mean_star:idec", -1, mean_star(0.0, 0.0, 0.0, NOVAS_FULL_ACCURACY, &x, NULL))) n++; + + return n; +} + +static int test_equ2gal() { + double x; + int n = 0; + + if(check("equ2gal:lon", -1, equ2gal(0.0, 0.0, NULL, &x))) n++; + if(check("equ2gal:lat", -1, equ2gal(0.0, 0.0, &x, NULL))) n++; + + return n; +} + +static int test_gal2equ() { + double x; + int n = 0; + + if(check("gal2equ:ra", -1, gal2equ(0.0, 0.0, NULL, &x))) n++; + if(check("gal2equ:dec", -1, gal2equ(0.0, 0.0, &x, NULL))) n++; + + return n; +} + +static int test_equ2ecl() { + double x; + int n = 0; + + if(check("equ2ecl:lon", -1, equ2ecl(0.0, NOVAS_GCRS_EQUATOR, NOVAS_FULL_ACCURACY, 0.0, 0.0, NULL, &x))) n++; + if(check("equ2ecl:lat", -1, equ2ecl(0.0, NOVAS_GCRS_EQUATOR, NOVAS_FULL_ACCURACY, 0.0, 0.0, &x, NULL))) n++; + + return n; +} + + +static int test_ecl2equ() { + double x; + int n = 0; + + if(check("ecl2equ:lon", -1, ecl2equ(0.0, NOVAS_GCRS_EQUATOR, NOVAS_FULL_ACCURACY, 0.0, 0.0, NULL, &x))) n++; + if(check("ecl2equ:lat", -1, ecl2equ(0.0, NOVAS_GCRS_EQUATOR, NOVAS_FULL_ACCURACY, 0.0, 0.0, &x, NULL))) n++; + + return n; +} + +static int test_equ2ecl_vec() { + double p[3]; + int n = 0; + + if(check("equ2ecl_vec:in", -1, equ2ecl_vec(0.0, NOVAS_MEAN_EQUATOR, NOVAS_FULL_ACCURACY, NULL, p))) n++; + if(check("equ2ecl_vec:out", -1, equ2ecl_vec(0.0, NOVAS_MEAN_EQUATOR, NOVAS_FULL_ACCURACY, p, NULL))) n++; + if(check("equ2ecl_vec:accuracy", -1, equ2ecl_vec(0.0, NOVAS_MEAN_EQUATOR, -1, p, p))) n++; + if(check("equ2ecl_vec:equator", 1, equ2ecl_vec(0.0, -1, NOVAS_FULL_ACCURACY, p, p))) n++; + + return n; +} + +static int test_ecl2equ_vec() { + double p[3]; + int n = 0; + + if(check("ecl2equ_vec:in", -1, ecl2equ_vec(0.0, NOVAS_MEAN_EQUATOR, NOVAS_FULL_ACCURACY, NULL, p))) n++; + if(check("ecl2equ_vec:out", -1, ecl2equ_vec(0.0, NOVAS_MEAN_EQUATOR, NOVAS_FULL_ACCURACY, p, NULL))) n++; + if(check("ecl2equ_vec:accuracy", -1, ecl2equ_vec(0.0, NOVAS_MEAN_EQUATOR, -1, p, p))) n++; + if(check("ecl2equ_vec:equator", 1, ecl2equ_vec(0.0, -1, NOVAS_FULL_ACCURACY, p, p))) n++; + + return n; +} + + +static int test_itrs_to_hor() { + on_surface loc; + double p[3], az, za; + int n = 0; + + if(check("itrs_to_hor:loc", -1, itrs_to_hor(NULL, p, &az, &za))) n++; + if(check("itrs_to_hor:in", -1, itrs_to_hor(&loc, NULL, &az, &za))) n++; + + return n; +} + +static int test_hor_to_itrs() { + on_surface loc; + double p[3]; + int n = 0; + + if(check("hor_to_itrs:loc", -1, hor_to_itrs(NULL, 0.0, 0.0, p))) n++; + if(check("hor_to_itrs:in", -1, hor_to_itrs(&loc, 0.0, 0.0, NULL))) n++; + + return n; +} + +static int test_equ2hor() { + on_surface loc; + double p[3], az, za, rar, decr; + int n = 0; + + if(check("equ2hor:loc", -1, equ2hor(0.0, 0.0, NOVAS_FULL_ACCURACY, 0.0, 0.0, NULL, 0.0, 0.0, NOVAS_STANDARD_ATMOSPHERE, &az, &za, &rar, &decr))) n++; + if(check("equ2hor:az", -1, equ2hor(0.0, 0.0, NOVAS_FULL_ACCURACY, 0.0, 0.0, &loc, 0.0, 0.0, NOVAS_STANDARD_ATMOSPHERE, NULL, &za, &rar, &decr))) n++; + if(check("equ2hor:zd", -1, equ2hor(0.0, 0.0, NOVAS_FULL_ACCURACY, 0.0, 0.0, &loc, 0.0, 0.0, NOVAS_STANDARD_ATMOSPHERE, &az, NULL, &rar, &decr))) n++; + + return n; +} + +static int test_gcrs2equ() { + double ra, dec; + int n = 0; + + if(check("gcrs2equ:ra", -1, gcrs2equ(0.0, NOVAS_DYNAMICAL_MOD, NOVAS_FULL_ACCURACY, 0.0, 0.0, NULL, &dec))) n++; + if(check("gcrs2equ:dec", -1, gcrs2equ(0.0, NOVAS_DYNAMICAL_MOD, NOVAS_FULL_ACCURACY, 0.0, 0.0, &ra, NULL))) n++; + if(check("gcrs2equ:sys", -1, gcrs2equ(0.0, -1, NOVAS_FULL_ACCURACY, 0.0, 0.0, &ra, &dec))) n++; + + return n; +} + +static int test_sidereal_time() { + double x; + int n = 0; + + if(check("sidereal_time:out", -1, sidereal_time(0.0, 0.0, 0.0, NOVAS_MEAN_EQUINOX, EROT_GST, NOVAS_FULL_ACCURACY, NULL))) n++; + if(check("sidereal_time:accuracy", 1, sidereal_time(0.0, 0.0, 0.0, NOVAS_MEAN_EQUINOX, EROT_GST, -1, &x))) n++; + if(check("sidereal_time:erot", 2, sidereal_time(0.0, 0.0, 0.0, NOVAS_MEAN_EQUINOX, -1, NOVAS_FULL_ACCURACY, &x))) n++; + + return n; +} + +static int test_ter2cel() { + double p[3]; + int n = 0; + + if(check("ter2cel:in", -1, ter2cel(0.0, 0.0, 0.0, EROT_GST, NOVAS_FULL_ACCURACY, NOVAS_DYNAMICAL_CLASS, 0.0, 0.0, NULL, p))) n++; + if(check("ter2cel:out", -1, ter2cel(0.0, 0.0, 0.0, EROT_GST, NOVAS_FULL_ACCURACY, NOVAS_DYNAMICAL_CLASS, 0.0, 0.0, p, NULL))) n++; + if(check("ter2cel:accuracy", 1, ter2cel(0.0, 0.0, 0.0, EROT_GST, -1, NOVAS_DYNAMICAL_CLASS, 0.0, 0.0, p, p))) n++; + if(check("ter2cel:erot", 2, ter2cel(0.0, 0.0, 0.0, -1, NOVAS_FULL_ACCURACY, NOVAS_DYNAMICAL_CLASS, 0.0, 0.0, p, p))) n++; + + return n; +} + +static int test_cel2ter() { + double p[3]; + int n = 0; + + if(check("cel2ter:in", -1, cel2ter(0.0, 0.0, 0.0, EROT_GST, NOVAS_FULL_ACCURACY, NOVAS_DYNAMICAL_CLASS, 0.0, 0.0, NULL, p))) n++; + if(check("cel2ter:out", -1, cel2ter(0.0, 0.0, 0.0, EROT_GST, NOVAS_FULL_ACCURACY, NOVAS_DYNAMICAL_CLASS, 0.0, 0.0, p, NULL))) n++; + if(check("cel2ter:accuracy", 1, cel2ter(0.0, 0.0, 0.0, EROT_GST, -1, NOVAS_DYNAMICAL_CLASS, 0.0, 0.0, p, p))) n++; + if(check("cel2ter:erot", 2, cel2ter(0.0, 0.0, 0.0, -1, NOVAS_FULL_ACCURACY, NOVAS_DYNAMICAL_CLASS, 0.0, 0.0, p, p))) n++; + + return n; +} + + +static int test_spin() { + double p[3]; + int n = 0; + + if(check("spin:in", -1, spin(0.0, NULL, p))) n++; + if(check("spin:out", -1, spin(0.0, p, NULL))) n++; + + return n; +} + +static int test_wobble() { + double p[3]; + int n = 0; + + if(check("wobble:in", -1, wobble(0.0, WOBBLE_ITRS_TO_PEF, 0.0, 0.0, NULL, p))) n++; + if(check("wobble:out", -1, wobble(0.0, WOBBLE_ITRS_TO_PEF, 0.0, 0.0, p, NULL))) n++; + + return n; +} + +static int test_terra() { + on_surface loc; + double p[3], v[3]; + int n = 0; + + if(check("terra:loc", -1, terra(NULL, 0.0, p, v))) n++; + if(check("terra:same", -1, terra(&loc, 0.0, p, p))) n++; + + return n; +} + +static int test_e_tilt() { + on_surface loc; + double p[3], v[3]; + int n = 0; + + if(check("e_tilt:accuracy", -1, e_tilt(0.0, -1, NULL, NULL, NULL, NULL, NULL))) n++; + + return n; +} + + + +static int test_cel_pole() { + on_surface loc; + double p[3], v[3]; + int n = 0; + + if(check("cel_pole:type", 1, cel_pole(0.0, -1, 0.0, 0.0))) n++; + + return n; +} + +static int test_frame_tie() { + double p[3]; + int n = 0; + + if(check("frame_tie:in", -1, frame_tie(NULL, 0, p))) n++; + if(check("frame_tie:out", -1, frame_tie(p, 0, NULL))) n++; + + return n; +} + +static int test_proper_motion() { + double p[3], v[3]; + int n = 0; + + if(check("frame_tie:p", -1, proper_motion(0.0, NULL, v, 1.0, p))) n++; + if(check("frame_tie:v", -1, proper_motion(0.0, p, NULL, 1.0, p))) n++; + if(check("frame_tie:out", -1, proper_motion(0.0, p, v, 1.0, NULL))) n++; + + return n; +} + +static int test_bary2obs() { + double p[3], po[3], out[3], lt; + int n = 0; + + if(check("bary2obs:pos", -1, bary2obs(NULL, po, out, <))) n++; + if(check("bary2obs:obs", -1, bary2obs(p, NULL, out, <))) n++; + if(check("bary2obs:out", -1, bary2obs(p, po, NULL, <))) n++; + + return n; +} + +static int test_geo_posvel() { + observer o; + double p[3], v[3]; + int n = 0; + + o.where = NOVAS_OBSERVER_ON_EARTH; + if(check("geo_posvel:loc", -1, geo_posvel(0.0, 0.0, NOVAS_FULL_ACCURACY, NULL, p, v))) n++; + if(check("geo_posvel:same", -1, geo_posvel(0.0, 0.0, NOVAS_FULL_ACCURACY, &o, p, p))) n++; + if(check("geo_posvel:accuracy", 1, geo_posvel(0.0, 0.0, -1, &o, p, v))) n++; + + o.where = -1; + if(check("geo_posvel:where", 2, geo_posvel(0.0, 0.0, NOVAS_FULL_ACCURACY, &o, p, v))) n++; + + return n; +} + +static int test_light_time2() { + object o; + double pos[3], p[3], v[3], t; + int n = 0; + + if(check("light_time2:tout", -1, light_time2(0.0, &o, pos, 0.0, NOVAS_FULL_ACCURACY, p, v, NULL))) n++; + if(check("light_time2:object", -1, light_time2(0.0, NULL, pos, 0.0, NOVAS_FULL_ACCURACY, p, v, &t))) n++; + if(check("light_time2:pos", -1, light_time2(0.0, &o, NULL, 0.0, NOVAS_FULL_ACCURACY, p, v, &t))) n++; + if(check("light_time2:same1", -1, light_time2(0.0, &o, pos, 0.0, NOVAS_FULL_ACCURACY, pos, v, &t))) n++; + if(check("light_time2:same2", -1, light_time2(0.0, &o, pos, 0.0, NOVAS_FULL_ACCURACY, p, pos, &t))) n++; + if(check("light_time2:same3", -1, light_time2(0.0, &o, pos, 0.0, NOVAS_FULL_ACCURACY, p, p, &t))) n++; + + return n; +} + +static int test_d_light() { + double p[3]; + int n = 0; + + if(check("d_light:1", 1, isnan(d_light(NULL, p)))) n++; + if(check("d_light:2", 1, isnan(d_light(p, NULL)))) n++; + + return n; +} + +static int test_cio_array() { + ra_of_cio x[5]; + int n = 0; + + if(check("cio_array:out", -1, cio_array(0.0, 5, NULL))) n++; + if(check("cio_array:n_pts:lo", 3, cio_array(0.0, 1, x))) n++; + if(check("cio_array:n_pts:hi", 3, cio_array(0.0, NOVAS_CIO_CACHE_SIZE + 1, x))) n++; + + set_cio_locator_file("blah"); + if(check("cio_array:file", 1, cio_array(0.0, 5, x))) n++; + + set_cio_locator_file("../cio_ra.bin"); + if(check("cio_array:beg", 2, cio_array(0.0, 5, x))) n++; + if(check("cio_array:end", 2, cio_array(1e20, 5, x))) n++; + + if(check("cio_array:corner:lo", 6, cio_array(2341952.6, 5, x))) n++; + if(check("cio_array:corner:hi", 6, cio_array(2561137.4, 5, x))) n++; + + return n; +} + +static int test_cio_basis() { + double x[3], y[3], z[3]; + int n = 0; + + if(check("cio_basis:x", -1, cio_basis(0.0, 0.0, CIO_VS_GCRS, NOVAS_FULL_ACCURACY, NULL, y, z))) n++; + if(check("cio_basis:y", -1, cio_basis(0.0, 0.0, CIO_VS_GCRS, NOVAS_FULL_ACCURACY, x, NULL, z))) n++; + if(check("cio_basis:z", -1, cio_basis(0.0, 0.0, CIO_VS_GCRS, NOVAS_FULL_ACCURACY, x, y, NULL))) n++; + if(check("cio_basis:accuracy", -1, cio_basis(0.0, 0.0, CIO_VS_GCRS, -1, x, y, z))) n++; + if(check("cio_basis:ref", 1, cio_basis(0.0, 0.0, -1, NOVAS_FULL_ACCURACY, x, y, z))) n++; + + return n; +} + +static int test_cio_location() { + double x; + short type; + int n = 0; + + if(check("cio_location:ra", -1, cio_location(0.0, NOVAS_FULL_ACCURACY, NULL, &type))) n++; + if(check("cio_location:type", -1, cio_location(0.0, NOVAS_FULL_ACCURACY, &x, NULL))) n++; + + return n; +} + +static int test_cio_ra() { + double x; + int n = 0; + + if(check("cio_location:ra", -1, cio_ra(0.0, NOVAS_FULL_ACCURACY, NULL))) n++; + if(check("cio_location:accuracy", 1, cio_ra(0.0, -1, &x))) n++; + + return n; +} + +static int test_starvectors() { + cat_entry star; + double p[3], v[3]; + int n = 0; + + if(check("starvectors:star", -1, starvectors(NULL, p, v))) n++; + if(check("starvectors:same", -1, starvectors(&star, p, p))) n++; + + return n; +} + +static int test_radec2vector() { + int n = 0; + + if(check("radec2vector", -1, radec2vector(0.0, 0.0, 1.0, NULL))) n++; + + return n; +} + +static int test_vector2radec() { + double p[3] = { }, ra, dec; + int n = 0; + + if(check("vector2radec:vec", -1, vector2radec(NULL, &ra, &dec))) n++; + + if(check("vector2radec:vec", 1, vector2radec(p, &ra, &dec))) n++; + + p[2] = 1.0; + if(check("vector2radec:vec", 2, vector2radec(p, &ra, &dec))) n++; + + return n; +} + +static int test_planet_lon() { + int n = 0; + + if(check("planet_lon", 1, isnan(planet_lon(0.0, -1)))) n++; + + return n; +} + +static int test_fund_args() { + int n = 0; + + if(check("find_args", -1, fund_args(0.0, NULL))) n++; + + return n; +} + +static int test_nutation_angles() { + double x; + int n = 0; + + if(check("nutation_angles:dpsi", -1, nutation_angles(0.0, NOVAS_FULL_ACCURACY, NULL, &x))) n++; + if(check("nutation_angles:deps", -1, nutation_angles(0.0, NOVAS_FULL_ACCURACY, &x, NULL))) n++; + + return n; +} + +static int test_set_nutation_lp_provider() { + int n = 0; + + if(check("set_nutation_lp_provider", -1, set_nutation_lp_provider(NULL))) n++; + + return n; +} + +static int test_nutation() { + double p[3]; + int n = 0; + + if(check("nutation:in", -1, nutation(0.0, NUTATE_MEAN_TO_TRUE, NOVAS_FULL_ACCURACY, NULL, p))) n++; + if(check("nutation:out", -1, nutation(0.0, NUTATE_MEAN_TO_TRUE, NOVAS_FULL_ACCURACY, p, NULL))) n++; + + return n; +} + +static int test_precession() { + double p[3]; + int n = 0; + + if(check("precesion:in", -1, precession(0.0, NULL, 1.0, p))) n++; + if(check("precesion:out", -1, precession(0.0, p, 1.0, NULL))) n++; + + return n; +} + +static int test_rad_vel() { + object o; + double p[3], v[3], vo[3], rv; + int n = 0; + + o.type = NOVAS_PLANET; + + if(check("rad_vel:object", -1, rad_vel(NULL, p, v, vo, 1.0, 1.0, 1.0, &rv))) n++; + if(check("rad_vel:pos", -1, rad_vel(&o, NULL, v, vo, 1.0, 1.0, 1.0, &rv))) n++; + if(check("rad_vel:vel", -1, rad_vel(&o, p, NULL, vo, 1.0, 1.0, 1.0, &rv))) n++; + if(check("rad_vel:vobs", -1, rad_vel(&o, p, v, NULL, 1.0, 1.0, 1.0, &rv))) n++; + if(check("rad_vel:out", -1, rad_vel(&o, p, v, vo, 1.0, 1.0, 1.0, NULL))) n++; + + o.type = -1; + if(check("rad_vel", -1, rad_vel(&o, p, v, vo, 1.0, 1.0, 1.0, &rv))) n++; + + return n; +} + +static int test_aberration() { + double p[3], v[3] = {}; + int n = 0; + + if(check("aberration:pos", -1, aberration(NULL, v, 0.0, p))) n++; + if(check("aberration:vel", -1, aberration(p, NULL, 0.0, p))) n++; + if(check("aberration:out", -1, aberration(p, v, 0.0, NULL))) n++; + + return n; +} + +static int test_grav_vec() { + double p[3], po[3], pb[3]; + int n = 0; + + if(check("grav_vec:pos", -1, grav_vec(NULL, po, pb, 1.0, p))) n++; + if(check("grav_vec:po", -1, grav_vec(p, NULL, pb, 1.0, p))) n++; + if(check("grav_vec:pb", -1, grav_vec(p, po, NULL, 1.0, p))) n++; + if(check("grav_vec:out", -1, grav_vec(p, po, pb, 1.0, NULL))) n++; + if(check("grav_vec:same", -1, grav_vec(p, po, pb, 1.0, po))) n++; + + return n; +} + +static int test_grav_def() { + double p[3], po[3], pb[3]; + int n = 0; + + if(check("grav_def:pos", -1, grav_def(0.0, NOVAS_OBSERVER_AT_GEOCENTER, NOVAS_FULL_ACCURACY, NULL, po, p))) n++; + if(check("grav_def:po", -1, grav_def(0.0, NOVAS_OBSERVER_AT_GEOCENTER, NOVAS_FULL_ACCURACY, p, NULL, p))) n++; + if(check("grav_def:out", -1, grav_def(0.0, NOVAS_OBSERVER_AT_GEOCENTER, NOVAS_FULL_ACCURACY, p, po, NULL))) n++; + if(check("grav_def:same", -1, grav_def(0.0, NOVAS_OBSERVER_AT_GEOCENTER, NOVAS_FULL_ACCURACY, p, po, po))) n++; + + return n; +} + +static int test_earth_sun_calc() { + double p[3], v[3]; + int n = 0; + + if(check("earth_sun_calc:pos", -1, earth_sun_calc(NOVAS_JD_J2000, NOVAS_SUN, NOVAS_BARYCENTER, NULL, v))) n++; + if(check("earth_sun_calc:vel", -1, earth_sun_calc(NOVAS_JD_J2000, NOVAS_SUN, NOVAS_BARYCENTER, p, NULL))) n++; + if(check("earth_sun_calc:tdb:lo", 1, earth_sun_calc(2340000.0, NOVAS_SUN, NOVAS_BARYCENTER, p, v))) n++; + if(check("earth_sun_calc:tdb:hi", 1, earth_sun_calc(2560001.0, NOVAS_SUN, NOVAS_BARYCENTER, p, v))) n++; + if(check("earth_sun_calc:number", 2, earth_sun_calc(NOVAS_JD_J2000, NOVAS_JUPITER, NOVAS_BARYCENTER, p, v))) n++; + + return n; +} + +static int test_sun_eph() { + extern int sun_eph(double jd, double *ra, double *dec, double *dis); + + double ra, dec, dis; + int n = 0; + + if(check("sun_eph", -1, sun_eph(NOVAS_JD_J2000, NULL, &dec, &dis))) n++; + if(check("sun_eph", -1, sun_eph(NOVAS_JD_J2000, &ra, NULL, &dis))) n++; + if(check("sun_eph", -1, sun_eph(NOVAS_JD_J2000, &ra, &dec, NULL))) n++; + + return n; +} + +int main() { + int n = 0; + + if(test_make_on_surface()) n++; + if(test_make_in_space()) n++; + if(test_make_observer()) n++; + + if(test_make_object()) n++; + if(test_make_ephem_object()) n++; + if(test_make_planet()) n++; + if(test_make_cat_entry()) n++; + if(test_transform_cat()) n++; + if(test_transform_hip()) n++; + + if(test_refract()) n++; + if(test_limb_angle()) n++; + + if(test_ephemeris()) n++; + + if(test_j2000_to_tod()) n++; + if(test_tod_to_j2000()) n++; + if(test_gcrs_to_cirs()) n++; + if(test_cirs_to_gcrs()) n++; + + if(test_set_planet_provider()) n++; + if(test_set_planet_provider_hp()) n++; + + if(test_place()) n++; + if(test_place_star()) n++; + if(test_radec_planet()) n++; + if(test_mean_star()) n++; + + if(test_equ2gal()) n++; + if(test_gal2equ()) n++; + + if(test_equ2ecl_vec()) n++; + if(test_ecl2equ_vec()) n++; + if(test_equ2ecl()) n++; + if(test_ecl2equ()) n++; + + if(test_itrs_to_hor()) n++; + if(test_hor_to_itrs()) n++; + + if(test_equ2hor()) n++; + if(test_gcrs2equ()) n++; + + if(test_sidereal_time()) n++; + if(test_ter2cel()) n++; + if(test_cel2ter()) n++; + + if(test_spin()) n++; + if(test_wobble()) n++; + if(test_terra()) n++; + if(test_e_tilt()) n++; + if(test_cel_pole()) n++; + if(test_frame_tie()) n++; + + if(test_proper_motion()) n++; + if(test_bary2obs()) n++; + if(test_geo_posvel()) n++; + + if(test_light_time2()) n++; + if(test_d_light()) n++; + + if(test_cio_array()) n++; + if(test_cio_basis()) n++; + if(test_cio_location()) n++; + if(test_cio_ra()) n++; + + if(test_starvectors()) n++; + if(test_radec2vector()) n++; + if(test_vector2radec()) n++; + + if(test_planet_lon()) n++; + if(test_fund_args()) n++; + if(test_nutation_angles()) n++; + if(test_set_nutation_lp_provider()) n++; + if(test_nutation()) n++; + if(test_precession()) n++; + if(test_rad_vel()) n++; + if(test_aberration()) n++; + if(test_grav_vec()) n++; + if(test_grav_def()) n++; + + if(test_earth_sun_calc()) n++; + if(test_sun_eph()) n++; + + if(n) fprintf(stderr, " -- FAILED %d tests\n", n); + else fprintf(stderr, " -- OK\n"); + + return n; +} diff --git a/test/src/test-super.c b/test/src/test-super.c index 3a6e7f5c..9faa2bfc 100644 --- a/test/src/test-super.c +++ b/test/src/test-super.c @@ -28,6 +28,30 @@ static double yp = -2.0; // Initialized quantities. static double pos0[3]; +static enum novas_origin ephem_origin; + +static short dummy_planet_hp(const double *jd_tdb, enum novas_planet body, enum novas_origin origin, double *position, double *velocity) { + memset(position, 0, 3 * sizeof(double)); + memset(velocity, 0, 3 * sizeof(double)); + position[0] = body % 10; + velocity[1] = 0.01 * (body % 10); + return 0; +} + +static short dummy_planet(double jd_tdb, enum novas_planet body, enum novas_origin origin, double *position, double *velocity) { + double tdb2[2] = { tdb }; + return dummy_planet_hp(tdb2, body, origin, position, velocity); +} + +static int dummy_ephem(const char *name, long id, double jd_tdb_high, double jd_tdb_low, enum novas_origin *origin, double *pos, double *vel) { + *origin = ephem_origin; + memset(pos, 0, 3 * sizeof(double)); + memset(vel, 0, 3 * sizeof(double)); + pos[0] = id % 100; + vel[1] = 0.01 * (id % 100); + return 0; +} + static int check_equal_pos(const double *posa, const double *posb, double tol) { int i; @@ -68,8 +92,8 @@ static int test_j2000_tod_j2000() { static int test_tod_itrs_tod() { double pos1[3]; - if(!is_ok("tod_to_itrs", tod_to_itrs(tdb, ut12tt, 0, xp, yp, pos0, pos1))) return 1; - if(!is_ok("itrs_to_tod", itrs_to_tod(tdb, ut12tt, 0, xp, yp, pos1, pos1))) return 1; + if(!is_ok("tod_to_itrs", tod_to_itrs(tdb, 0.0, ut12tt, 0, xp, yp, pos0, pos1))) return 1; + if(!is_ok("itrs_to_tod", itrs_to_tod(tdb, 0.0, ut12tt, 0, xp, yp, pos1, pos1))) return 1; if(!is_ok("tod_itrs_tod", check_equal_pos(pos0, pos1, 1e-9 * vlen(pos0)))) return 1; return 0; } @@ -87,8 +111,8 @@ static int test_gcrs_cirs_gcrs() { static int test_cirs_itrs_cirs() { double pos1[3]; - if(!is_ok("cirs_to_itrs", cirs_to_itrs(tdb, ut12tt, 0, xp, yp, pos0, pos1))) return 1; - if(!is_ok("itrs_to_cirs", itrs_to_cirs(tdb, ut12tt, 0, xp, yp, pos1, pos1))) return 1; + if(!is_ok("cirs_to_itrs", cirs_to_itrs(tdb, 0.0, ut12tt, 0, xp, yp, pos0, pos1))) return 1; + if(!is_ok("itrs_to_cirs", itrs_to_cirs(tdb, 0.0, ut12tt, 0, xp, yp, pos1, pos1))) return 1; if(!is_ok("cirs_itrs_cirs", check_equal_pos(pos0, pos1, 1e-9 * vlen(pos0)))) return 1; return 0; } @@ -114,8 +138,8 @@ static int test_tod_vs_cirs() { if(!is_ok("gcrs_to_j2000", gcrs_to_j2000(pos0, pos1))) return 1; if(!is_ok("j2000_to_tod", j2000_to_tod(tdb, 0, pos1, pos1))) return 1; - if(!is_ok("tod_to_itrs", tod_to_itrs(tdb, ut12tt, 0, xp, yp, pos1, pos1))) return 1; - if(!is_ok("itrs_to_cirs", itrs_to_cirs(tdb, ut12tt, 0, xp, yp, pos1, pos1))) return 1; + if(!is_ok("tod_to_itrs", tod_to_itrs(tdb, 0.0, ut12tt, 0, xp, yp, pos1, pos1))) return 1; + if(!is_ok("itrs_to_cirs", itrs_to_cirs(tdb, 0.0, ut12tt, 0, xp, yp, pos1, pos1))) return 1; if(!is_ok("cirs_to_gcrs", cirs_to_gcrs(tdb, 0, pos1, pos1))) return 1; if(!is_ok("tod_vs_cirs", check_equal_pos(pos0, pos1, 1e-9 * vlen(pos0)))) return 1; @@ -123,6 +147,20 @@ static int test_tod_vs_cirs() { return 0; } +static int test_equ_ecl() { + double ra, dec, elon, elat, pos1[3]; + + vector2radec(pos0, &ra, &dec); + if(!is_ok("equ2ecl", equ2ecl(tdb, NOVAS_GCRS_EQUATOR, NOVAS_FULL_ACCURACY, ra, dec, &elon, &elat))) return 1; + if(!is_ok("ecl2equ", ecl2equ(tdb, NOVAS_GCRS_EQUATOR, NOVAS_FULL_ACCURACY, elon, elat, &ra, &dec))) return 1; + radec2vector(ra, dec, vlen(pos0), pos1); + + if(!is_ok("equ_ecl_equ", check_equal_pos(pos0, pos1, 1e-9 * vlen(pos0)))) return 1; + + return 0; +} + + static int test_equ_gal() { double ra, dec, glon, glat, pos1[3]; @@ -131,7 +169,7 @@ static int test_equ_gal() { if(!is_ok("gal2equ", gal2equ(glon, glat, &ra, &dec))) return 1; radec2vector(ra, dec, vlen(pos0), pos1); - if(!is_ok("tod_vs_cirs", check_equal_pos(pos0, pos1, 1e-9 * vlen(pos0)))) return 1; + if(!is_ok("equ_gal_equ", check_equal_pos(pos0, pos1, 1e-9 * vlen(pos0)))) return 1; return 0; } @@ -243,6 +281,7 @@ static int test_source() { if(test_itrs_hor_itrs()) n++; + if(test_equ_ecl()) n++; if(test_equ_gal()) n++; if(test_place_star()) n++; @@ -269,17 +308,6 @@ static int test_make_planet() { return 0; } -static int test_make_ephem_body() { - object body; - - make_ephem_body("Ceres", 1000001, &body); - - if(!is_ok("make_ephem_body:type", body.type != NOVAS_EPHEM_OBJECT)) return 1; - if(!is_ok("make_ephem_body:number", body.number != 1000001)) return 1; - if(!is_ok("make_ephem_body:name", strcasecmp(body.name, "Ceres"))) return 1; - - return 0; -} static int test_precession() { double pos1[3], pos2[3]; @@ -295,6 +323,7 @@ static int test_precession() { } + static int test_radec_planet() { int i; object sun; @@ -322,12 +351,9 @@ static int test_observers() { double ps[3] = { 100.0, 30.0, 10.0 }, vs[3] = { 10.0 }; int n = 0; - if(test_make_planet()) n++; - if(test_make_ephem_body()) n++; if(test_precession()) n++; if(test_radec_planet()) n++; - make_observer_at_geocenter(&obs); n += test_source(); @@ -383,6 +409,30 @@ static int test_get_utc_to_tt() { return 0; } +static int test_nutation_lp_provider() { + double t = (tdb - NOVAS_JD_J2000) / 36525.0; + double de, dp, de0, dp0; + + int status = 1; + + + if(!is_ok("nutation_lp_provider:set_nutation_lp_provider", set_nutation_lp_provider(iau2000b))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("nutation_lp_provider:nutation_angles", nutation_angles(t, NOVAS_REDUCED_ACCURACY, &de, &dp))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("nutation_lp_provider:iau2000b", iau2000b(tdb, 0.0, &de0, &dp0))) goto cleanup; // @suppress("Goto statement used") + + de0 /= ASEC2RAD; + dp0 /= ASEC2RAD; + + if(!is_ok("nutation_lp_provider:check_de", fabs(de - de0) > 1e-7)) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("nutation_lp_provider:check_dp", fabs(dp - dp0) > 1e-7)) goto cleanup; // @suppress("Goto statement used") + + status = 0; + + cleanup: + + set_nutation_lp_provider(nu2000k); + return status; +} static int test_dates() { double offsets[] = {-10000.0, 0.0, 10000.0, 10000.0, 10000.01 }; @@ -390,6 +440,7 @@ static int test_dates() { if(test_get_ut1_to_tt()) n++; if(test_get_utc_to_tt()) n++; + if(test_nutation_lp_provider()) n++; for(i = 0; i < 5; i++) { tdb = J2000 + offsets[i]; @@ -442,13 +493,120 @@ static int test_case() { return 0; } +static int test_make_ephem_object() { + object body; + + make_ephem_object("Ceres", 1000001, &body); + + if(!is_ok("make_ephem_object:type", body.type != NOVAS_EPHEM_OBJECT)) return 1; + if(!is_ok("make_ephem_object:number", body.number != 1000001)) return 1; + if(!is_ok("make_ephem_object:name", strcasecmp(body.name, "Ceres"))) return 1; + + return 0; +} + +static int test_planet_provider() { + int status = 1; + object mars; + double p[3], v[3], p0[3], v0[3]; + double tdb2[2] = { tdb }; + + make_planet(NOVAS_MARS, &mars); + + if(!is_ok("planet_provider:set_planet_provider", set_planet_provider(dummy_planet))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:set_planet_provider_hp", set_planet_provider_hp(dummy_planet_hp))) goto cleanup; // @suppress("Goto statement used") + + if(!is_ok("planet_provider:ephemeris", ephemeris(tdb2, &mars, NOVAS_BARYCENTER, NOVAS_REDUCED_ACCURACY, p, v))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:control", dummy_planet(tdb, NOVAS_MARS, NOVAS_BARYCENTER, p0, v0))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:check_pos", check_equal_pos(p, p0, 1e-9 * vlen(p0)))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:check_vel", check_equal_pos(v, v0, 1e-9 * vlen(v0)))) goto cleanup; // @suppress("Goto statement used") + + if(!is_ok("planet_provider:ephemeris_hp", ephemeris(tdb2, &mars, NOVAS_BARYCENTER, NOVAS_FULL_ACCURACY, p, v))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:control_hp", dummy_planet_hp(tdb2, NOVAS_MARS, NOVAS_BARYCENTER, p0, v0))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:check_pos_hp", check_equal_pos(p, p0, 1e-9 * vlen(p0)))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:check_vel_hp", check_equal_pos(v, v0, 1e-9 * vlen(v0)))) goto cleanup; // @suppress("Goto statement used") + + status = 0; + + cleanup: + + set_planet_provider(earth_sun_calc); + set_planet_provider_hp(earth_sun_calc_hp); + return status; +} + + +static int test_ephem_provider() { + novas_ephem_provider prior = get_ephem_provider(); + + object body; + double p[3], v[3], p0[3], v0[3]; + double tdb2[2] = { tdb }; + int status = 1; + enum novas_origin o; + + make_ephem_object("Dummy", 1000001, &body); + + if(!is_ok("ephem_provider:set_ephem_provider", set_ephem_provider(dummy_ephem))) goto cleanup; // @suppress("Goto statement used") + + for(ephem_origin = 0; ephem_origin < 2; ephem_origin++) { + if(!is_ok("planet_provider:ephemeris", ephemeris(tdb2, &body, NOVAS_BARYCENTER, NOVAS_FULL_ACCURACY, p, v))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:control", dummy_ephem(body.name, body.number, tdb, 0.0, &o, p0, v0))) goto cleanup; // @suppress("Goto statement used") + if(o == NOVAS_BARYCENTER) { + if(!is_ok("planet_provider:check_pos", check_equal_pos(p, p0, 1e-9 * vlen(p0)))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:check_vel", check_equal_pos(v, v0, 1e-9 * vlen(v0)))) goto cleanup; // @suppress("Goto statement used") + } + + if(!is_ok("planet_provider:ephemeris", ephemeris(tdb2, &body, NOVAS_HELIOCENTER, NOVAS_FULL_ACCURACY, p, v))) goto cleanup; // @suppress("Goto statement used") + if(o == NOVAS_BARYCENTER) { + fprintf(stderr, " Expecing diffent A/B, twice:\n"); + if(!is_ok("planet_provider:check_pos", !check_equal_pos(p, p0, 1e-9 * vlen(p0)))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("planet_provider:check_vel", !check_equal_pos(v, v0, 1e-9 * vlen(v0)))) goto cleanup; // @suppress("Goto statement used") + fprintf(stderr, " OK.\n"); + } + } + + status = 0; + + cleanup: + + set_ephem_provider(prior); + return status; +} + + +static int test_enable_earth_sun_calc_hp() { + double tdb2[2] = { tdb }, p[3], v[3], p0[3], v0[3]; + int status = 1; + + enable_earth_sun_hp(1); + + if(!is_ok("enable_earth_sun_hp", earth_sun_calc(tdb, NOVAS_SUN, NOVAS_BARYCENTER, p0, v0))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("enable_earth_sun_hp", earth_sun_calc_hp(tdb2, NOVAS_SUN, NOVAS_BARYCENTER, p, v))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("enable_earth_sun_hp:check_pos", check_equal_pos(p, p0, 1e-9 * vlen(p0)))) goto cleanup; // @suppress("Goto statement used") + if(!is_ok("enable_earth_sun_hp:check_vel", check_equal_pos(v, v0, 1e-9 * vlen(v0)))) goto cleanup; // @suppress("Goto statement used") + + status = 0; + + cleanup: + + enable_earth_sun_hp(0); + return status; +} + int main() { int n = 0; make_object(NOVAS_CATALOG_OBJECT, 0, "None", NULL, &source); + if(test_make_planet()) n++; + if(test_make_ephem_object()) n++; if(test_refract_astro()) n++; if(test_case()) n++; + if(test_planet_provider()) n++; + if(test_ephem_provider()) n++; + if(test_enable_earth_sun_calc_hp()) n++; + n += test_dates(); return n;