Merge pull request #28 from JuliaSpaceMissionDesign/dev

Update to v1.3.1

.github/workflows/ci.yml

@@ -12,6 +12,7 @@ jobs:
       fail-fast: false
       matrix:
         version:
+          - '1.8'
           - '1.9'
           - '1.10'
           - '1'

Project.toml

@@ -1,7 +1,7 @@
 name = "Tempo"
 uuid = "c33777b2-e695-4fae-9135-aeae8855dd81"
 authors = ["JSMD Team"]
-version = "1.3.0"
+version = "1.3.1"
 
 [deps]
 FunctionWrappersWrappers = "77dc65aa-8811-40c2-897b-53d922fa7daf"

docs/src/tutorials/t01_epochs.md

@@ -92,12 +92,16 @@ The [`Epoch`](@ref) type supports a limited subset of basic mathematical and log
 e1 = Epoch(90.0, TT)
 e2 = Epoch(50.0, TT)
 
-e1 - e2
+Δe = e1 - e2
+
+value(Δe)
 
 e3 = Epoch(40, TAI)
 e1 - e3
 ```
-Notice that this operation can be performed only if the two epochs are defined on the same timescale.
+Notice that this operation can be performed only if the two epochs are defined on the same timescale. When computing the difference between two epochs, the result is returned in the 
+form of a [`Duration`](@ref) object. The [`value`](@ref) can then be used to retrieve the 
+actual number of seconds it represents.
 
 Epochs can also be shifted forward and backwards in time by adding or subtracting an arbitrary number of seconds: 
 ```@repl init 
@@ -163,7 +167,8 @@ eTAI = convert(TAI, e)
 
 A special remark must be made on the conversion between TAI and UTC. The offset between these two timescales is defined by a leap seconds, which are introduced to keep the UTC time scale within 0.9 seconds from UT1. Since the rotation of the Earth is irregular, it is not possible to predict when a new leap second will be introduced in the future. 
 
-The latest NAIF's leap second kernel ([LSK](https://naif.jpl.nasa.gov/pub/naif/generic_kernels/lsk)) is embedded within `Tempo` as a package artifact, which will be manually updated each time a new kernel is released, so that the user effort is minimised. Indeed, transforming an [`Epoch`](@ref) from a generic timescale to UTC is a simple as:
+A leapsecond table is embedded within `Tempo` and will be manually updated each time a new 
+leapsecond is introduced, so that the effort required from the user side is minimised. Indeed, transforming an [`Epoch`](@ref) from a generic timescale to UTC is a simple as:
 
 ```@repl init
 e = Epoch(90.0, TT)

src/Tempo.jl

@@ -38,18 +38,22 @@ include("offset.jl")
 
 export TIMESCALES, @timescale, add_timescale!, 
        TimeSystem, timescale_alias, timescale_name, timescale_id
+
 include("scales.jl")
 
 export Date, Time,
        year, month, day, find_dayinyear,
        j2000, j2000s,j2000c, hour, minute, second, DateTime
+
 include("datetime.jl")
 include("origin.jl")
 
-export Duration, value 
+export Duration, value
+
 include("duration.jl")
 
 export Epoch, j2000, j2000s, j2000c, doy, timescale, value
+
 include("epoch.jl")
 
 # Package precompilation routines

src/duration.jl

@@ -1,62 +1,83 @@
 
 """
-    Duration{T}
+    Duration{T} <: Number
 
 A `Duration` represents a period of time, split into an integer number of seconds and a 
-fractional part.
+fractional part for increased precision.
 
 ### Fields
 - `seconds`: The integer number of seconds.
 - `fraction`: The fractional part of the duration, where `T` is a subtype of `Number`.
+
+---
+
+    Duration(seconds::Number)
+
+Create a `Duration` object from a number of seconds. The type of the fractional part will 
+be inferred from the type of the input argument. 
+
+---
+
+    Duration{T}(seconds::Number)
+
+Create a `Duration` object from a number of seconds with the fractional part of type `T`. 
+
+### Examples 
+```julia-repl 
+julia> d = Duration(10.783)
+Duration{Float64}(10, 0.7829999999999995)
+
+julia> value(d) 
+10.783
+
+julia> d = Duration{BigFloat64}(10.3)
+Duration{BigFloat}(10, 0.300000000000000710542735760100185871124267578125)
+```
 """
-struct Duration{T} 
+struct Duration{T} <: Number 
     seconds::Int 
     fraction::T 
 end
 
-function Duration(seconds::T) where {T<:Number}
-    i, f = divrem(seconds, 1)
-    return Duration{T}(i, f)
+function Duration{T}(seconds::Number) where {T <: Number}
+    i,f = divrem(seconds, 1)
+    return Duration{T}(convert(Int, i), T(f))
 end
 
-function Duration(sec::Int, frac::T) where {T<:Number}
-    return Duration{T}(sec, frac)
-end
+function Duration(seconds::T) where {T <: Number}
+    Duration{T}(seconds)  
+end 
+
+ftype(::Duration{T}) where T = T
+
+"""
+    value(d::Duration)
 
+Return the duration `d`, in seconds.
+"""
 value(d::Duration{T}) where T = d.seconds + d.fraction
 
-function Base.isless(d::Duration{T}, q::Number) where T
-    return value(d) < q 
-end
+# ---
+# Type Conversions and Promotions 
 
-function Base.isless(d::Duration{T1}, d2::Duration{T2}) where {T1, T2}
-    return value(d) < value(d2)
+function Base.convert(::Type{Duration{T}}, d::Duration{S}) where {T,S}
+    return Duration(d.seconds, convert(T, d.fraction))
 end
 
-function fmasf(a, b, mul)
-    amulb = fma(mul, a, b)
-    i, f = divrem(amulb, 1)
-    return i, f 
+function Base.convert(::Type{T}, d::Duration{S}) where {T<:Number,S}
+    return Duration(d.seconds, convert(T, d.fraction))
 end
 
-function Base.:-(d1::Duration, d2::Duration) 
-    s1, f1 = d1.seconds, d1.fraction
-    s2, f2 = d2.seconds, d2.fraction
-    ds, df = divrem(f1 - f2, 1)
-    sec = s1 - s2 + ds 
-    if df < 0
-        sec -= 1
-        df += 1
-    end
-    return Duration(convert(Int, sec), df)
+function Base.promote_rule(::Type{Duration{T}}, ::Type{Duration{S}}) where {T,S}
+    return promote_rule(T, S)
 end
 
-function Base.:+(d1::Duration, d2::Duration) 
-    s1, f1 = d1.seconds, d1.fraction
-    s2, f2 = d2.seconds, d2.fraction
-    s, f = fmasf(f1, f2, 1)
-    return Duration(convert(Int, s1 + s2 + s), f)
-end
+# ----
+# Operations
+
+Base.isless(d::Duration, q::Number) = value(d) < q
+Base.isless(q::Number, d::Duration) = q < value(d)
+Base.isless(d1::Duration, d2::Duration) = value(d1) < value(d2)
 
 function Base.:+(d::Duration, x::Number)
     es, ef = d.seconds, d.fraction
@@ -65,6 +86,13 @@ function Base.:+(d::Duration, x::Number)
     return Duration(convert(Int, es + xs + s), f)
 end
 
+function Base.:+(d1::Duration, d2::Duration) 
+    s1, f1 = d1.seconds, d1.fraction
+    s2, f2 = d2.seconds, d2.fraction
+    s, f = fmasf(f1, f2, 1)
+    return Duration(convert(Int, s1 + s2 + s), f)
+end
+
 function Base.:-(d::Duration, x::Number)
     es, ef = d.seconds, d.fraction
     xs, xf = divrem(x, 1)
@@ -77,3 +105,21 @@ function Base.:-(d::Duration, x::Number)
     return Duration(convert(Int, sec), df)
 end
 
+function Base.:-(d1::Duration, d2::Duration) 
+    s1, f1 = d1.seconds, d1.fraction
+    s2, f2 = d2.seconds, d2.fraction
+    ds, df = divrem(f1 - f2, 1)
+    sec = s1 - s2 + ds 
+    if df < 0
+        sec -= 1
+        df += 1
+    end
+    return Duration(convert(Int, sec), df)
+end
+
+
+function fmasf(a, b, mul)
+    amulb = fma(mul, a, b)
+    i, f = divrem(amulb, 1)
+    return i, f 
+end

src/epoch.jl

@@ -70,8 +70,8 @@ struct Epoch{S,T}
 end
 
 function Epoch{S}(seconds::Number) where {S<:AbstractTimeScale}
-    sec, frac = divrem(seconds, 1)    
-    return Epoch{S, typeof(frac)}(S(), Duration(Int(sec), frac))
+    d = Duration(seconds)
+    return Epoch{S, ftype(d)}(S(), d)
 end
 
 Epoch(sec::Number, ::S) where {S<:AbstractTimeScale} = Epoch{S}(sec)
@@ -83,15 +83,15 @@ Epoch(dt::DateTime, ::Type{S}) where {S<:AbstractTimeScale} = Epoch{S}(j2000s(dt
 Epoch(e::Epoch) = e
 
 Epoch{S,T}(e::Epoch{S,T}) where {S,T} = e
-Epoch{S,T}(e::Epoch{S,N}) where {S, N, T} = Epoch{S}(T(j2000s(e)))
+Epoch{S,T}(e::Epoch{S,N}) where {S, N, T} = Epoch{S,T}(e.scale, convert(T, e.dur))
 
 # Construct an epoch from an ISO string and a scale
 function Epoch(s::AbstractString, scale::S) where {S <: AbstractTimeScale}
     y, m, d, H, M, sec, sf = parse_iso(s)
 
     # TODO: the precision of this could be improved
     _, jd2 = calhms2jd(y, m, d, H, M, sec + sf)
-    return Epoch(jd2 * 86400, scale)
+    return Epoch(jd2 * DAY2SEC, scale)
 end
 
 # Construct an epoch from an ISO string
@@ -201,21 +201,8 @@ function Base.:-(::Epoch{S1}, ::Epoch{S2}) where {S1, S2}
     throw(ErrorException("only epochs defined in the same timescale can be subtracted."))
 end 
 
-function Base.:+(e::Epoch{S, N}, x::Number) where {S, N}
-    return Epoch{S, N}(timescale(e), e.dur + x)
-end
-
-function Base.:-(e::Epoch{S, N}, x::Number) where {S, N}
-    return Epoch{S, N}(timescale(e), e.dur - x)
-end
-
-function Base.:+(e::Epoch{S, N}, d::Duration{<:Number}) where {S, N}
-    return Epoch{S, N}(timescale(e), e.dur + d)
-end
-
-function Base.:-(e::Epoch{S, N}, d::Duration{<:Number}) where {S, N}
-    return Epoch{S, N}(timescale(e), e.dur - d)
-end
+Base.:+(e::Epoch, x::Number) = Epoch(timescale(e), e.dur + x)
+Base.:-(e::Epoch, x::Number) = Epoch(timescale(e), e.dur - x)
 
 function (::Base.Colon)(start::Epoch, step::Number, stop::Epoch)
     step = start < stop ? step : -step
@@ -227,6 +214,9 @@ function (::Base.Colon)(start::Epoch, step::Duration, stop::Epoch)
     return (:)(start, value(step), stop)
 end
 
+(::Base.Colon)(start::Epoch, stop::Epoch) = (:)(start, 86400, stop)
+
+
 Base.isless(e1::Epoch{S}, e2::Epoch{S}) where {S} = e1.dur < e2.dur
 
 function Base.isapprox(e1::Epoch{S}, e2::Epoch{S}; kwargs...) where {S}

test/Tempo/duration.jl

@@ -1,6 +1,7 @@
 
 # TestSet for Duration
 @testset "Duration" begin
+
     # Test Duration constructor with a floating point number
     d1 = Duration(5.75)
     @test d1.seconds == 5
@@ -11,10 +12,25 @@
     @test d2.seconds == 10
     @test d2.fraction == 0.25
 
+    # Check Constructor with a Big Float
+    db1 = Duration(BigFloat(2.5422))
+
+    @test db1.seconds == 2 
+    @test db1.fraction ≈ 0.5422 atol=1e-14 rtol=1e-14
+
+    # Test duration type 
+    @test Tempo.ftype(db1) == BigFloat
+
     # Test value function
     @test value(d1) == 5.75
     @test value(d2) == 10.25
 
+    # Test duration conversion 
+    db2 = convert(BigFloat, d1)
+    @test Tempo.ftype(db2) == BigFloat 
+    @test db2.seconds == 5 
+    @test db2.fraction ≈ 0.75 atol=1e-14 rtol=1e-14
+
     # Test isless with a number
     @test d1 < 6.0 
     @test d2 < 10.5 
@@ -63,4 +79,5 @@
     @test d10.fraction == d1.fraction
     @test d11.seconds == d1.seconds
     @test d11.fraction == d1.fraction
+
 end

test/Tempo/epoch.jl

@@ -68,8 +68,10 @@
     @test_throws ErrorException e3-e1
 
     ems = e1:86400:e2 
+    ems2 = e1:e2 
     for j = 2:lastindex(ems)
         @test ems[j] == e1 + 86400*(j-1)
+        @test ems2[j] == e1 + 86400*(j-1)
     end
 
     # Based on Vallado "Fundamental of astrodynamics" page 196

test/runtests.jl

@@ -3,5 +3,5 @@ using ERFA
 using Test
 
 @testset "Tempo.jl" verbose = true begin
-    include("Tempo/Tempo.jl")
+    include(joinpath("Tempo", "Tempo.jl"))
 end;