Merge 99f8dc8 into 3799fae

src/Epochs.jl

@@ -5,10 +5,10 @@ using ERFA
 
 import Base: +, -, ==, isapprox, isless
 
-using ..TimeScales, ..Periods
+using ..TimeScales
 import ..TimeScales: acronyms
-import ..LeapSeconds: leapseconds
-
+import ..LeapSeconds: leapseconds, LSK_DATA
+using ..Periods
 export Epoch, julian, julian1, julian2, +, -, ==, isapprox, isless,
     leapseconds, jd2000, jd1950, mjd, timescale
 

src/LeapSeconds.jl

@@ -3,9 +3,31 @@ module LeapSeconds
 using ERFA
 using OptionalData
 using RemoteFiles
+include("Periods.jl")
+using .Periods:MJD, LS_1972, DRIFT
 
+export leapseconds, LSK, LSK_FILE, LSK_DATA, fractionofday
+
+
+struct changes
+    t ::Vector{Float64}
+    leapseconds ::Vector{Float64}
+    drift :: Vector{Tuple}
+end
+function changes(date, drift_dat)
+    t = Vector{Float64}()
+    leapseconds = Vector{Float64}()
+    drift = Vector{Tuple}()
+    for i in range(1,size(date)[1])
+        push!(t, Dates.datetime2julian(DateTime(date[i][1], date[i][2], 1)))
+        push!(leapseconds, date[i][3])
+        push!(drift,drift_dat[i])
+    end
+    changes(t,leapseconds,drift)
+end
+
+CHANGE = changes(LS_1972, DRIFT)
 
-export leapseconds, LSK, LSK_FILE, LSK_DATA
 
 struct LSK
     t::Vector{Float64}
@@ -41,11 +63,14 @@ function leapseconds(lsk::LSK, jd)
         return 0.0
     elseif jd < lsk.t[1]
         dt = Dates.julian2datetime(jd)
-        return ERFA.dat(Dates.year(dt), Dates.month(dt), Dates.day(dt), fractionofday(dt))
+        fd =  fractionofday(dt)
+        index = findlast(jd .>= CHANGE.t)
+        deltat = CHANGE.leapseconds[index]
+        deltat += ((jd- MJD- CHANGE.drift[index][1]) * CHANGE.drift[index][2])
+        return deltat
     else
         return lsk.leapseconds[findlast(jd .>= lsk.t)]
     end
 end
 leapseconds(jd) = leapseconds(get(LSK_DATA), jd)
-
 end

src/Periods.jl

@@ -1,7 +1,6 @@
 module Periods
 
 import Base: *, /, get, isapprox, show
-
 export TimeUnit, Second, Minute, Hour, Day, Year, Century,
     J2000, J1950, MJD,
     seconds, minutes, hours, days, years, centuries,

src/constants.jl

@@ -4,7 +4,8 @@ export MJD, J2000, J1950,
     HOURS_PER_DAY, HOURS_PER_YEAR, HOURS_PER_CENTURY,
     DAYS_PER_YEAR, DAYS_PER_CENTURY,
     YEARS_PER_CENTURY,
-    OFFSET_TT_TAI, MOD_JD_77, ELG, fairhd, DAYS_PER_MILLENNIUM, TDB0, ELB, JD_MAX, JD_MIN
+    OFFSET_TT_TAI, MOD_JD_77, ELG, fairhd, DAYS_PER_MILLENNIUM, TDB0, ELB, JD_MAX, JD_MIN,
+    DRIFT, LS_1972
 
 const MJD = 2400000.5
 const J2000 = Dates.datetime2julian(DateTime(2000, 1, 1, 12, 0, 0))
@@ -42,4 +43,40 @@ const ELB = 1.550519768e-8
 const JD_MIN = -68569.5
 const JD_MAX = 1e9
 
+export DRIFT, LS_1972
+
+const DRIFT = [
+( 37300.0, 0.0012960 ),
+( 37300.0, 0.0012960 ),
+( 37300.0, 0.0012960 ),
+( 37665.0, 0.0011232 ),
+( 37665.0, 0.0011232 ),
+( 38761.0, 0.0012960 ),
+( 38761.0, 0.0012960 ),
+( 38761.0, 0.0012960 ),
+( 38761.0, 0.0012960 ),
+( 38761.0, 0.0012960 ),
+( 38761.0, 0.0012960 ),
+( 38761.0, 0.0012960 ),
+( 39126.0, 0.0025920 ),
+( 39126.0, 0.0025920 )]
+
+LS_1972 = [
+( 1960,  1,  1.4178180 ),
+( 1961,  1,  1.4228180 ),
+( 1961,  8,  1.3728180 ),
+( 1962,  1,  1.8458580 ),
+( 1963, 11,  1.9458580 ),
+( 1964,  1,  3.2401300 ),
+( 1964,  4,  3.3401300 ),
+( 1964,  9,  3.4401300 ),
+( 1965,  1,  3.5401300 ),
+( 1965,  3,  3.6401300 ),
+( 1965,  7,  3.7401300 ),
+( 1965,  9,  3.8401300 ),
+( 1966,  1,  4.3131700 ),
+( 1968,  2,  4.2131700 )]
+
+
+
 include("fairhd.jl")

src/conversions.jl

@@ -1,10 +1,11 @@
 using ItemGraphs: ItemGraph, add_edge!, edgeitems
 
 using MuladdMacro
+using ..LeapSeconds
 using ..Periods
 export @transform
 
-using ..LeapSeconds
+
 
 abstract type Transformation end
 
@@ -561,7 +562,7 @@ end
 
 end
 
-function cal2jd(iy, im, id)
+@inline function cal2jd(iy, im, id)
     EYEAR_ALLOWED = -4799
     MON_LENGTH = [31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31]
 
@@ -591,7 +592,6 @@ function cal2jd(iy, im, id)
 end
 
 
-
 # TAI <-> UTC
 @transform UTC TAI ep begin
     jd1, jd2 = julian1(ep), julian2(ep)

test/runtests.jl

@@ -220,9 +220,15 @@ AstroTime.update()
         @test_throws ArgumentError Epochs.cal2jd(-4800, 1, 1)
         @test_throws ArgumentError Epochs.cal2jd(2000, 15, 1)
         @test_throws ArgumentError Epochs.cal2jd(2000, 1, 40)
-    
+
         @test Epochs.jd2cal(julian1(tt), julian2(tt)) == ERFA.jd2cal(julian1(tt), julian2(tt))
         @test Epochs.jd2cal(julian2(tt), julian1(tt)) == ERFA.jd2cal(julian2(tt), julian1(tt))
+
+        dt = DateTime(1963, 04, 01, 15, 30, 33)
+        jd_dt = Dates.datetime2julian(dt)
+        fraction_day = fractionofday(dt)
+        @test ERFA.dat(Dates.year(dt), Dates.month(dt), Dates.day(dt),fraction_day) ≈ leapseconds(jd_dt)
+
     end
     @testset "Leap Seconds" begin
         @test leapseconds(TTEpoch(1959,1,1)) == 0