@@ -81,19 +81,44 @@ static double PSI_COR = 0.0;
81
81
*/
82
82
static double EPS_COR = 0.0 ;
83
83
84
- static FILE * cio_file ; ///< Opened CIO locator data file, or NULL.
84
+ ///< Opened CIO locator data file, or NULL.
85
+ static FILE * cio_file ;
85
86
87
+ /// function to use for calculating positions for major planets
86
88
static novas_planet_calculator planetcalc = NULL ;
89
+
90
+ /// function to use for calculating positions for major planets with high precison
87
91
static novas_planet_calculator_hp planetcalc_hp = NULL ;
88
92
93
+ /// function to use for reading ephemeris data for all types of solar system sources
89
94
static novas_ephem_reader_func readeph2_call = NULL ;
90
95
96
+ /// Function to use for reduced-precision calculations. (The full IAU 2000A model is used always for high-precision calculations)
91
97
static novas_nutate_func nutate_lp = nu2000k ;
92
98
99
+ /**
100
+ * Calculates the length of a 3-vector
101
+ *
102
+ * @param v Pointer to a 3-component (x, y, z) vector. The argument cannot be NULL
103
+ * @return the length of the vector
104
+ *
105
+ * @sa vdot
106
+ * @sa vdist
107
+ */
93
108
static double vlen (const double * v ) {
94
109
return sqrt (v [0 ] * v [0 ] + v [1 ] * v [1 ] + v [2 ] * v [2 ]);
95
110
}
96
111
112
+ /**
113
+ * Calculates the distance between two 3-vectors.
114
+ *
115
+ * @param v1 Pointer to a 3-component (x, y, z) vector. The argument cannot be NULL
116
+ * @param v2 Pointer to another 3-component (x, y, z) vector. The argument cannot be NULL
117
+ * @return The distance between the two vectors
118
+ *
119
+ * @sa vlen()
120
+ * @sa vdot()
121
+ */
97
122
static double vdist (const double * v1 , const double * v2 ) {
98
123
double d2 = 0.0 ;
99
124
int i ;
@@ -104,6 +129,16 @@ static double vdist(const double *v1, const double *v2) {
104
129
return sqrt (d2 );
105
130
}
106
131
132
+ /**
133
+ * Calculates the dot product between two 3-vectors.
134
+ *
135
+ * @param v1 Pointer to a 3-component (x, y, z) vector. The argument cannot be NULL
136
+ * @param v2 Pointer to another 3-component (x, y, z) vector. The argument cannot be NULL
137
+ * @return The dot product between the two vectors.
138
+ *
139
+ * @sa vlen()
140
+ * @sa vdist()
141
+ */
107
142
static double vdot (const double * v1 , const double * v2 ) {
108
143
return (v1 [0 ] * v2 [0 ]) + (v1 [1 ] * v2 [1 ]) + (v1 [2 ] * v2 [2 ]);
109
144
}
@@ -125,10 +160,36 @@ static void tiny_rotate(const double *in, double ax, double ay, double az, doubl
125
160
out [2 ] = in [2 ] - 0.5 * (A [0 ] + A [1 ]) * in [2 ] - ay * in [0 ] + ax * in [1 ];
126
161
}
127
162
163
+ /**
164
+ * Checks if two Julian dates are equals under the precision that can be handled by this library.
165
+ * In practive two dates are considered equal if they agree within 10<sup>-8</sup> days (or about
166
+ * 1 ms) of each other.
167
+ *
168
+ *
169
+ * @param jd1 [day] a Julian date (in any time measure)
170
+ * @param jd2 [day] a Julian date in the same time measure as the first argument
171
+ * @return TRUE (1) if the two dates are effectively the same at the precision of comparison,
172
+ * or else FALSE (0) if they differ by more than the allowed tolerance.
173
+ */
128
174
static int time_equals (double jd1 , double jd2 ) {
129
175
return fabs (jd1 - jd2 ) <= 1.0e-8 ;
130
176
}
131
177
178
+ /**
179
+ * Transforms a rectangular equatorial (x, y, z) vector from J2000 coordinates to the
180
+ * True equinox Of Date (TOD) reference frame at the given epoch
181
+ *
182
+ * @param jd_tdb [day] Barycentric Dynamical Time (TDB) based Julian date that defines the output
183
+ * epoch. Typically it does not require much precision, and Julian dates
184
+ * in other time measures will be unlikely to affect the result
185
+ * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
186
+ * @param in Input (x, y, z) position or velocity vector in rectangular equatorial coordinates at J2000
187
+ * @param[out] out Output position or velocity 3-vector in the True equinox of Date coordinate frame.
188
+ * @return 0 if successful, or -1 if either of the vector arguments is NULL.
189
+ *
190
+ * @sa tod_to_j2000()
191
+ * @sa icrs_to_tod()
192
+ */
132
193
int j2000_to_tod (double jd_tdb , enum novas_accuracy accuracy , const double * in , double * out ) {
133
194
double v [3 ];
134
195
@@ -147,6 +208,21 @@ int j2000_to_tod(double jd_tdb, enum novas_accuracy accuracy, const double *in,
147
208
return 0 ;
148
209
}
149
210
211
+ /**
212
+ * Transforms a rectangular equatorial (x, y, z) vector from True equinox Of Date (TOD) reference frame at the
213
+ * given epoch to the J2000 coordinates.
214
+ *
215
+ * @param jd_tdb [day] Barycentric Dynamical Time (TDB) based Julian date that defines the input
216
+ * epoch. Typically it does not require much precision, and Julian dates
217
+ * in other time measures will be unlikely to affect the result
218
+ * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
219
+ * @param in Input (x, y, z) position or velocity 3-vector in the True equinox of Date coordinate frame.
220
+ * @param[out] out Output position or velocity vector in rectangular equatorial coordinates at J2000
221
+ * @return 0 if successful, or -1 if either of the vector arguments is NULL.
222
+ *
223
+ * @sa j2000_to_tod()
224
+ * @sa tod_to_icrs()
225
+ */
150
226
int tod_to_j2000 (double jd_tdb , enum novas_accuracy accuracy , const double * in , double * out ) {
151
227
double v [3 ];
152
228
@@ -166,13 +242,43 @@ int tod_to_j2000(double jd_tdb, enum novas_accuracy accuracy, const double *in,
166
242
return 0 ;
167
243
}
168
244
245
+ /**
246
+ * Transforms a rectangular equatorial (x, y, z) vector from the International Celestial
247
+ * Reference System (ICRS) to the True equinox Of Date (TOD) reference frame at the given epoch
248
+ *
249
+ * @param jd_tdb [day] Barycentric Dynamical Time (TDB) based Julian date that defines the output
250
+ * epoch. Typically it does not require much precision, and Julian dates
251
+ * in other time measures will be unlikely to affect the result
252
+ * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
253
+ * @param in ICRS Input (x, y, z) position or velocity vector
254
+ * @param[out] out Output position or velocity 3-vector in the True equinox of Date coordinate frame.
255
+ * @return 0 if successful, or -1 if either of the vector arguments is NULL.
256
+ *
257
+ * @sa tod_to_icrs()
258
+ * @sa j2000_to_tod()
259
+ */
169
260
int icrs_to_tod (double jd_tdb , enum novas_accuracy accuracy , const double * in , double * out ) {
170
261
double j2000 [3 ];
171
262
int error = frame_tie (in , TIE_ICRS_TO_J2000 , j2000 );
172
263
if (error ) return error ;
173
264
return j2000_to_tod (jd_tdb , accuracy , j2000 , out );
174
265
}
175
266
267
+ /**
268
+ * Transforms a rectangular equatorial (x, y, z) vector from True equinox Of Date (TOD) reference frame at the
269
+ * given epoch to the International Celestial Reference System(ICRS)
270
+ *
271
+ * @param jd_tdb [day] Barycentric Dynamical Time (TDB) based Julian date that defines the input
272
+ * epoch. Typically it does not require much precision, and Julian dates
273
+ * in other time measures will be unlikely to affect the result
274
+ * @param accuracy NOVAS_FULL_ACCURACY (0) or NOVAS_REDUCED_ACCURACY (1)
275
+ * @param in Input (x, y, z) position or velocity 3-vector in the True equinox of Date coordinate frame.
276
+ * @param[out] out Output ICRS position or velocity vector
277
+ * @return 0 if successful, or -1 if either of the vector arguments is NULL.
278
+ *
279
+ * @sa j2000_to_tod()
280
+ * @sa tod_to_icrs()
281
+ */
176
282
int tod_to_icrs (double jd_tdb , enum novas_accuracy accuracy , const double * in , double * out ) {
177
283
double j2000 [3 ];
178
284
int error = tod_to_j2000 (jd_tdb , accuracy , in , j2000 );
@@ -1693,7 +1799,7 @@ short sidereal_time(double jd_high, double jd_low, double ut1_to_tt, enum novas_
1693
1799
// AK: Fix for documented bug in NOVAS 3.1 --> 3.1.1
1694
1800
ha_eq -= (eqeq / 3600.0 );
1695
1801
1696
- ha_eq = fmod (ha_eq , DEG360 ) / 15.0 ;
1802
+ ha_eq = remainder (ha_eq , DEG360 ) / 15.0 ;
1697
1803
if (ha_eq < 0.0 ) ha_eq += DAY_HOURS ;
1698
1804
* gst = ha_eq ;
1699
1805
return 0 ;
@@ -1708,7 +1814,7 @@ short sidereal_time(double jd_high, double jd_low, double ut1_to_tt, enum novas_
1708
1814
st = eqeq + 0.014506 + ((((-0.0000000368 * t - 0.000029956 ) * t - 0.00000044 ) * t + 1.3915817 ) * t + 4612.156534 ) * t ;
1709
1815
1710
1816
// Form the Greenwich sidereal time.
1711
- * gst = fmod ((st / 3600.0 + theta ), DEG360 ) / 15.0 ;
1817
+ * gst = remainder ((st / 3600.0 + theta ), DEG360 ) / 15.0 ;
1712
1818
1713
1819
if (* gst < 0.0 ) * gst += DAY_HOURS ;
1714
1820
return 0 ;
@@ -1746,9 +1852,9 @@ double era(double jd_high, double jd_low) {
1746
1852
1747
1853
thet1 = 0.7790572732640 + 0.00273781191135448 * (jd_high - JD_J2000 );
1748
1854
thet2 = 0.00273781191135448 * jd_low ;
1749
- thet3 = fmod (jd_high , 1.0 ) + fmod (jd_low , 1.0 );
1855
+ thet3 = remainder (jd_high , 1.0 ) + remainder (jd_low , 1.0 );
1750
1856
1751
- theta = fmod (thet1 + thet2 + thet3 , 1.0 ) * DEG360 ;
1857
+ theta = remainder (thet1 + thet2 + thet3 , 1.0 ) * DEG360 ;
1752
1858
if (theta < 0.0 ) theta += DEG360 ;
1753
1859
1754
1860
return theta ;
@@ -1998,7 +2104,7 @@ int spin(double angle, const double *pos1, double *pos2) {
1998
2104
return -1 ;
1999
2105
}
2000
2106
2001
- if (fmod (fabs (angle - ang_last ), DEG360 ) >= 1.0e-12 ) {
2107
+ if (remainder (fabs (angle - ang_last ), DEG360 ) >= 1.0e-12 ) {
2002
2108
const double angr = angle * DEGREE ;
2003
2109
const double cosang = cos (angr );
2004
2110
const double sinang = sin (angr );
@@ -3535,10 +3641,21 @@ int fund_args(double t, novas_fundamental_args *a) {
3535
3641
return 0 ;
3536
3642
}
3537
3643
3644
+ /**
3645
+ * Returns the planetary longitude, for Mercury through Neptune, w.r.t. mean dynamical
3646
+ * ecliptic and equinox of J2000, with high order terms omitted
3647
+ * (Simon et al. 1994, 5.8.1-5.8.8).
3648
+ *
3649
+ * @param t [cy] Julian centuries since J2000
3650
+ * @param planet Novas planet id, e.g. NOVAS_MARS.
3651
+ * @return [rad] The approximate longitude of the planet in radians.
3652
+ *
3653
+ * @sa accum_prec()
3654
+ * @sa nutation()
3655
+ * @sa nutation_angles()
3656
+ * @sa NOVAS_JD_J2000
3657
+ */
3538
3658
double planet_lon (double t , enum novas_planet planet ) {
3539
- // Planetary longitudes, Mercury through Neptune, wrt mean dynamical
3540
- // ecliptic and equinox of J2000, with high order terms omitted
3541
- // (Simon et al. 1994, 5.8.1-5.8.8).
3542
3659
double lon ;
3543
3660
3544
3661
switch (planet ) {
@@ -3571,13 +3688,20 @@ double planet_lon(double t, enum novas_planet planet) {
3571
3688
return NAN ;
3572
3689
}
3573
3690
3574
- return fmod (lon , TWOPI );
3691
+ return remainder (lon , TWOPI );
3575
3692
}
3576
3693
3694
+ /**
3695
+ * Returns the general precession in longitude (Simon et al. 1994), equivalent
3696
+ * to 5028.8200 arcsec/cy at J2000.
3697
+ *
3698
+ * @param t [cy] Julian centuries since J2000
3699
+ * @return [rad] the approximate precession angle [-π:π].
3700
+ */
3577
3701
double accum_prec (double t ) {
3578
3702
// General precession in longitude (Simon et al. 1994), equivalent
3579
3703
// to 5028.8200 arcsec/cy at J2000.
3580
- return fmod ((0.024380407358 + 0.000005391235 * t ) * t , TWOPI );
3704
+ return remainder ((0.024380407358 + 0.000005391235 * t ) * t , TWOPI );
3581
3705
}
3582
3706
3583
3707
/**
@@ -4792,7 +4916,7 @@ int cal_date(double tjd, short *year, short *month, short *day, double *hour) {
4792
4916
djd = tjd + 0.5 ;
4793
4917
jd = (long ) floor (djd );
4794
4918
4795
- h = fmod (djd , 1.0 ) * DAY_HOURS ;
4919
+ h = remainder (djd , 1.0 ) * DAY_HOURS ;
4796
4920
4797
4921
k = jd + 68569L ;
4798
4922
n = 4L * k / 146097L ;
0 commit comments