Merge pull request #14 from JuliaSpaceMissionDesign/dev

Update to v.1.1.2

Author: MicheleCeresoli

Project.toml

@@ -1,7 +1,7 @@
 name = "IERSConventions"
 uuid = "4e86e20e-879b-40dc-9e12-cee74f4cd199"
 authors = ["JSMD Team"]
-version = "1.2.0"
+version = "1.1.2"
 
 [deps]
 DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab"

src/delaunay.jl

@@ -49,20 +49,24 @@ centuries since J2000.
 """
 function delaunay_anomaly_moon(::IERSAModels, tdb_c::Number)
     # This expression is valid for both IERS2003 and IERS2010
-    mod2pi(
-        arcsec2rad(
-            @evalpoly(tdb_c, 485868.249036, 1717915923.2178, 31.8792, 0.051635, -0.00024470)
+    arcsec2rad(
+        jrem(
+            @evalpoly(tdb_c, 485868.249036, 1717915923.2178, 31.8792, 0.051635, -0.00024470), 
+            1296000
         )
     )
 end
 
 function delaunay_anomaly_moon(::IERSModel, tdb_c::Number)
-    arcsec2rad(mod(485868.249036 + 1717915923.2178tdb_c, 1296000))
+    arcsec2rad(jrem(485868.249036 + 1717915923.2178tdb_c, 1296000))
 end
 
 function delaunay_anomaly_moon(::IERS1996, tdb_c::Number)
-    r = mod(1325tdb_c, 1)
-    mod2pi(arcsec2rad(@evalpoly(tdb_c, 485866.733, 715922.633, 31.310, 0.064)) + 2π*r)
+    r = jrem(1325tdb_c, 1)
+    rem2pi(
+        arcsec2rad(@evalpoly(tdb_c, 485866.733, 715922.633, 31.310, 0.064)) + 2π*r, 
+        RoundNearest
+    )
 end
 
 
@@ -73,20 +77,24 @@ Compute the mean anomaly of the Sun, in radians, given time `tdb_c` expressed in
 centuries since J2000.
 """
 function delaunay_anomaly_sun(::IERSAModels, tdb_c::Number)
-    mod2pi(
-        arcsec2rad(
-            @evalpoly(tdb_c, 1287104.793048, 129596581.0481, -0.5532, 0.000136, -0.00001149)
+    arcsec2rad(
+        jrem(
+            @evalpoly(tdb_c, 1287104.793048, 129596581.0481, -0.5532, 0.000136, -0.00001149), 
+            1296000
         )
     )
 end
 
 function delaunay_anomaly_sun(::IERSModel, tdb_c::Number)
-    arcsec2rad(mod(1287104.79305 + 129596581.0481tdb_c, 1296000))
+    arcsec2rad(jrem(1287104.79305 + 129596581.0481tdb_c, 1296000))
 end
 
 function delaunay_anomaly_sun(::IERS1996, tdb_c::Number)
-    r = mod(99tdb_c, 1)
-    mod2pi(arcsec2rad(@evalpoly(tdb_c, 1287099.804, 1292581.224, -0.577, -0.012)) + 2π*r)
+    r = jrem(99tdb_c, 1)
+    rem2pi(
+        arcsec2rad(@evalpoly(tdb_c, 1287099.804, 1292581.224, -0.577, -0.012)) + 2π*r, 
+        RoundNearest
+    )
 end
 
 
@@ -97,20 +105,24 @@ Compute the difference between the longitude of the Moon and the longitude of th
 node, in radians, given time `tdb_c` expressed in TDB Julian centuries since J2000.
 """
 function delaunay_longitude_diff(::IERSAModels, tdb_c::Number)
-    mod2pi(
-        arcsec2rad(
-            @evalpoly(tdb_c, 335779.526232, 1739527262.8478, -12.7512, -0.001037, +0.00000417)
+    arcsec2rad(
+        jrem(
+            @evalpoly(tdb_c, 335779.526232, 1739527262.8478, -12.7512, -0.001037, 0.00000417), 
+            1296000
         )
     )
 end
 
 function delaunay_longitude_diff(::IERSModel, tdb_c::Number)
-    arcsec2rad(mod(335779.526232 + 1739527262.8478tdb_c, 1296000))
+    arcsec2rad(jrem(335779.526232 + 1739527262.8478tdb_c, 1296000))
 end
 
 function delaunay_longitude_diff(::IERS1996, tdb_c::Number)
-    r = mod(1342tdb_c, 1)
-    mod2pi(arcsec2rad(@evalpoly(tdb_c, 335778.877, 295263.137, -13.257, 0.011)) + 2π*r)
+    r = jrem(1342tdb_c, 1)
+    rem2pi(
+        arcsec2rad(@evalpoly(tdb_c, 335778.877, 295263.137, -13.257, 0.011)) + 2π*r, 
+        RoundNearest
+    )
 end
 
 """
@@ -120,20 +132,24 @@ Compute the mean elongation of the Moon from the Sun, in radians, given time `td
 TDB Julian centuries since J2000.
 """
 function delaunay_elongation_moon(::IERSAModels, tdb_c::Number)
-    mod2pi(
-        arcsec2rad(
-            @evalpoly(tdb_c, 1072260.703692, 1602961601.2090, -6.3706, 0.006593, -0.00003169)
+    arcsec2rad(
+        jrem(
+            @evalpoly(tdb_c, 1072260.703692, 1602961601.2090, -6.3706, 0.006593, -0.00003169), 
+            1296000
         )
     )
 end
 
 function delaunay_elongation_moon(::IERSModel, tdb_c::Number)
-    arcsec2rad(mod(1072260.70369 + 1602961601.2090tdb_c, 1296000))
+    arcsec2rad(jrem(1072260.70369 + 1602961601.2090tdb_c, 1296000))
 end
 
 function delaunay_elongation_moon(::IERS1996, tdb_c::Number)
-    r = mod(1236tdb_c, 1)
-    mod2pi(arcsec2rad(@evalpoly(tdb_c, 1072261.307, 1105601.328, -6.891, 0.019)) + 2π*r)
+    r = jrem(1236tdb_c, 1)
+    rem2pi(
+        arcsec2rad(@evalpoly(tdb_c, 1072261.307, 1105601.328, -6.891, 0.019)) + 2π*r, 
+        RoundNearest
+    )
 end
 
 """
@@ -144,18 +160,26 @@ measured from the mean equinox of date, in radians, given time `tdb_c` expressed
 centuries since J2000.
 """
 function delaunay_longitude_node(::IERSAModels, tdb_c::Number)
-    mod2pi(
-        arcsec2rad(
-            @evalpoly(tdb_c, 450160.398036, -6962890.5431, 7.4722, 0.007702, -0.00005939)
+    arcsec2rad(
+        jrem(
+            @evalpoly(tdb_c, 450160.398036, -6962890.5431, 7.4722, 0.007702, -0.00005939), 
+            1296000
         )
     )
 end
 
 function delaunay_longitude_node(::IERSModel, tdb_c::Number)
-    arcsec2rad(mod(450160.398036 - 6962890.5431tdb_c, 1296000))
+    arcsec2rad(jrem(450160.398036 - 6962890.5431tdb_c, 1296000))
 end
 
 function delaunay_longitude_node(::IERS1996, tdb_c::Number)
-    r = mod(-5tdb_c, 1)
-    mod2pi(arcsec2rad(@evalpoly(tdb_c, 450160.280, -482890.539, 7.455, 0.008)) + 2π*r)
-end
+    r = jrem(-5tdb_c, 1)
+    rem2pi(
+        arcsec2rad(@evalpoly(tdb_c, 450160.280, -482890.539, 7.455, 0.008)) + 2π*r, 
+        RoundNearest
+    )
+end
+
+
+jrem(x, y) = x - y*(x ÷ y)
+jrem2pi(x) = jrem(x, 2π)

src/era.jl

@@ -18,10 +18,12 @@ since `J2000`, according to the IERS convention `m`.
 See also [`iers_era_rotm`](@ref).
 """
 function iers_era(::IERSModel, ut1_d::Number)
-
     # The function uses the fractional UT1 date to gain additional precision
     # in the computations
-    return mod2pi(2π * (mod(ut1_d, 1) + 0.7790572732640 + 0.00273781191135448ut1_d))
+
+    # Uses integer divsion because the `mod` and `rem` functions don't work with ForwardDiff
+    f = ut1_d - (ut1_d ÷ 1)
+    return rem2pi(2π * (f + 0.7790572732640 + 0.00273781191135448ut1_d), RoundDown)
 
 end
 

src/planetary.jl

@@ -68,7 +68,7 @@ Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/fame03.c) software library
 """
 function pa_mercury(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 4.402608842, 2608.7903141574))
+    return rem2pi(@evalpoly(tdb_c, 4.402608842, 2608.7903141574), RoundNearest)
 end
 
 function pa_mercury(::IERS1996, ::Number)
@@ -93,11 +93,11 @@ Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/fave03.c) software library
 """
 function pa_venus(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 3.176146697, 1021.3285546211))
+    return rem2pi(@evalpoly(tdb_c, 3.176146697, 1021.3285546211), RoundNearest)
 end
 
 function pa_venus(::IERS1996, tdb_c::Number)
-    return mod2pi(deg2rad(@evalpoly(tdb_c, 181.979800853, 58517.8156748)))
+    return rem2pi(deg2rad(@evalpoly(tdb_c, 181.979800853, 58517.8156748)), RoundNearest)
 end
 
 
@@ -113,11 +113,11 @@ Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/fae03.c) software library
 """
 function pa_earth(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 1.753470314, 628.3075849991))
+    return rem2pi(@evalpoly(tdb_c, 1.753470314, 628.3075849991), RoundNearest)
 end
 
 function pa_earth(::IERS1996, tdb_c::Number)
-    return mod2pi(deg2rad(@evalpoly(tdb_c, 100.466448494, 35999.3728521)))
+    return rem2pi(deg2rad(@evalpoly(tdb_c, 100.466448494, 35999.3728521)), RoundNearest)
 end
 
 
@@ -133,11 +133,11 @@ Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/fama03.c) software library
 """
 function pa_mars(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 6.203480913, 334.0612426700))
+    return rem2pi(@evalpoly(tdb_c, 6.203480913, 334.0612426700), RoundNearest)
 end
 
 function pa_mars(::IERS1996, tdb_c::Number)
-    return mod2pi(deg2rad(@evalpoly(tdb_c, 355.433274605, 19140.299314)))
+    return rem2pi(deg2rad(@evalpoly(tdb_c, 355.433274605, 19140.299314)), RoundNearest)
 end
 
 
@@ -153,11 +153,11 @@ Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/faju03.c) software library
 """
 function pa_jupiter(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 0.599546497, 52.9690962641))
+    return rem2pi(@evalpoly(tdb_c, 0.599546497, 52.9690962641), RoundNearest)
 end
 
 function pa_jupiter(::IERS1996, tdb_c::Number)
-    return mod2pi(deg2rad(@evalpoly(tdb_c, 34.3514839, 3034.90567464)))
+    return rem2pi(deg2rad(@evalpoly(tdb_c, 34.3514839, 3034.90567464)), RoundNearest)
 end
 
 
@@ -173,11 +173,11 @@ Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/fasa03.c) software library
 """
 function pa_saturn(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 0.874016757, 21.3299104960))
+    return rem2pi(@evalpoly(tdb_c, 0.874016757, 21.3299104960), RoundNearest)
 end
 
 function pa_saturn(::IERS1996, tdb_c::Number)
-    return mod2pi(deg2rad(@evalpoly(tdb_c, 50.0774713998, 1222.11379404)))
+    return rem2pi(deg2rad(@evalpoly(tdb_c, 50.0774713998, 1222.11379404)), RoundNearest)
 end
 
 
@@ -196,7 +196,7 @@ Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/faur03.c) software library
 """
 function pa_uranus(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 5.481293872, 7.4781598567))
+    return rem2pi(@evalpoly(tdb_c, 5.481293872, 7.4781598567), RoundNearest)
 end
 
 function pa_uranus(::IERS1996, ::Number)
@@ -224,7 +224,7 @@ Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/fane03.c) software library
 """
 function pa_neptune(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 5.311886287, 3.8133035638))
+    return rem2pi(@evalpoly(tdb_c, 5.311886287, 3.8133035638), RoundNearest)
 end
 
 function pa_neptune(::IERS1996, ::Number)
@@ -249,9 +249,9 @@ expressed in `TDB` Julian centuries since `J2000`.
 - [ERFA](https://github.com/liberfa/erfa/blob/master/src/fapa03.c) software library
 """
 function pa_precession(::IERSModel, tdb_c::Number)
-    return mod2pi(@evalpoly(tdb_c, 0, 0.024381750, 0.00000538691))
+    return jrem2pi(@evalpoly(tdb_c, 0, 0.024381750, 0.00000538691))
 end
 
 function pa_precession(::IERS1996, tdb_c::Number)
-    return mod2pi(deg2rad(@evalpoly(tdb_c, 0, 1.39697137214, 0.0003086)))
+    return jrem2pi(deg2rad(@evalpoly(tdb_c, 0, 1.39697137214, 0.0003086)))
 end

src/rate.jl

@@ -4,18 +4,11 @@ export iers_earth_rot_rate
 const ωₑ = 7.292_115_146_706_979e-5
 
 """ 
-    iers_earth_rot_rate(LOD::Number)
+    iers_earth_rot_rate(LOD::Number=0)
 
 Compute the true angular velocity of the Earth accounting for the Length of the Day, i.e., 
 the instantaneous rate of change of UT1 with respect to a uniform time scale. 
 """
-function iers_earth_rot_rate(LOD::Number)
+function iers_earth_rot_rate(LOD::Number=0)
     return ωₑ * (1 - LOD / 86400)
 end
-
-"""
-    iers_earth_rot_rate()
-    
-Compute the nominal Earth angular velocity. 
-"""
-iers_earth_rot_rate() = iers_earth_rot_rate(0.0)

src/rotations.jl

@@ -376,6 +376,9 @@ in TT seconds since `J2000`, following the IERS Conventions `m`.
 !!! note 
     Polar motion is neglected in the [`CPNd`](@ref) model. 
 
+!!! note 
+    The time derivative of the nutation and precession effects is neglected.
+
 ### References 
 - IERS Technical Note No. [36](https://www.iers.org/IERS/EN/Publications/TechnicalNotes/tn36.html) 
 """
@@ -414,6 +417,9 @@ Reference Frame (GCRF) to the International Terrestrial Reference Frame (ITRF) a
 !!! note 
     Polar motion is neglected in the [`CPNd`](@ref) model. 
 
+!!! note 
+    The time derivative of the nutation and precession effects is neglected.
+
 ### References 
 - IERS Technical Note No. [36](https://www.iers.org/IERS/EN/Publications/TechnicalNotes/tn36.html) 
 """
@@ -453,6 +459,9 @@ at time `tt_s`, expressed in TT seconds since `J2000`, following the IERS Conven
 !!! note 
     Polar motion is neglected in the [`CPNd`](@ref) model. 
 
+!!! note 
+    The time derivative of the nutation and precession effects is neglected.
+
 ### References 
 - IERS Technical Note No. [36](https://www.iers.org/IERS/EN/Publications/TechnicalNotes/tn36.html) 
 """

src/sidereal.jl

@@ -118,8 +118,9 @@ end
 # t = TT Julian centuries since J2000 
 function iers_gmst(::IERS1996, tt_c::Number)
 
-    # Transform from TT centuries to UT1 centuries
+    # Transform from TT centuries to UT1 centuries and days
     ut1_c = tt_c + offset_tt2ut1(tt_c*Tempo.CENTURY2SEC)/Tempo.CENTURY2SEC
+    ut1_d = ut1_c*Tempo.CENTURY2DAY
 
     # Coefficients for the IAU 1982 GMST-UT1 model. The first component has been adjusted 
     # of 12 hours because the input is defined with respect to noon of 01-01-2000
@@ -129,10 +130,10 @@ function iers_gmst(::IERS1996, tt_c::Number)
     D = -6.2e-6
 
     # Fractional part of UT1, in seconds
-    f = Tempo.DAY2SEC*mod(ut1_c*Tempo.CENTURY2DAY, 1)
+    f = Tempo.DAY2SEC*(ut1_d - (ut1_d ÷ 1))
 
     # Compute GMST
-    return mod2pi(2π/86400*(@evalpoly(ut1_c, A, B, C, D) + f))
+    return rem2pi(2π/86400*(@evalpoly(ut1_c, A, B, C, D) + f), RoundDown)
 
 end
 
@@ -164,7 +165,7 @@ function iers_gmst(::IERS2003, tt_c::Number, θ::Number)
     )
 
     # Compute GMST 
-    return mod2pi(θ + arcsec2rad(p))
+    return rem2pi(θ + arcsec2rad(p), RoundDown)
 
 end
 
@@ -190,6 +191,6 @@ function iers_gmst(::IERS2010, tt_c::Number, θ::Number)
     )
 
     # Compute GMST 
-    return mod2pi(θ + arcsec2rad(p))
+    return rem2pi(θ + arcsec2rad(p), RoundDown)
 
 end

test/Project.toml

@@ -1,6 +1,7 @@
 [deps]
 DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab"
 ERFA = "17511681-8477-586a-8d98-4cfd5a1f2ec3"
+ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"
 LazyArtifacts = "4af54fe1-eca0-43a8-85a7-787d91b784e3"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 ReferenceFrameRotations = "74f56ac7-18b3-5285-802d-d4bd4f104033"