From bf8bcf4eb5da167e8dac5733db10e08c38f75393 Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Sun, 4 Jul 2021 00:21:00 +0200
Subject: [PATCH 01/16] Set default format for all Epoch subtypes

---
 src/Epochs/dates.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/Epochs/dates.jl b/src/Epochs/dates.jl
index 55c846a..43258d8 100644
--- a/src/Epochs/dates.jl
+++ b/src/Epochs/dates.jl
@@ -7,7 +7,7 @@ function Epoch{S}(date::Date, time::Time, args...) where S
     return Epoch{S}(sec, time.fraction)
 end
 
-Dates.default_format(::Type{Epoch}) = EPOCH_ISO_FORMAT[]
+Dates.default_format(::Type{<:Epoch}) = EPOCH_ISO_FORMAT[]
 
 """
     Epoch(str[, format])

From 2433d5d566779b2b0dccbc8f22832974830f2aaf Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Sun, 4 Jul 2021 00:21:20 +0200
Subject: [PATCH 02/16] Improve precision of Epoch subtraction

---
 src/Epochs/operations.jl | 4 +++-
 test/epochs.jl           | 6 ++++++
 2 files changed, 9 insertions(+), 1 deletion(-)

diff --git a/src/Epochs/operations.jl b/src/Epochs/operations.jl
index 191dd1f..12e1b95 100644
--- a/src/Epochs/operations.jl
+++ b/src/Epochs/operations.jl
@@ -24,6 +24,8 @@ julia> TAIEpoch(2018, 2, 6, 20, 45, 20.0) - TAIEpoch(2018, 2, 6, 20, 45, 0.0)
 ```
 """
 function Base.:-(a::Epoch{S}, b::Epoch{S}) where S<:TimeScale
-    return ((a.second - b.second) + (a.fraction - b.fraction)) * seconds
+    second = a.second - b.second
+    sum, residual = two_sum(a.fraction, -b.fraction)
+    return (residual + sum + second) * seconds
 end
 
diff --git a/test/epochs.jl b/test/epochs.jl
index 0072bc2..aaa5104 100644
--- a/test/epochs.jl
+++ b/test/epochs.jl
@@ -83,6 +83,12 @@ import ERFA
         @test_throws MethodError t1 < t0
         @test t2 - t1 == -32.0seconds
         @test t2 < t1
+
+        today = TTEpoch(2000, 1, 1)
+        age_of_the_universe = 13.772e9years
+        big_bang = today - age_of_the_universe
+        baryons = big_bang + 1e-11seconds
+        @test baryons + age_of_the_universe - today == 1e-11seconds
     end
     @testset "Parsing" begin
         @test AstroTime.TimeScales.tryparse(1.0) === nothing

From eb398590cc069062cba51983aa5c90921f96e95c Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Sun, 4 Jul 2021 09:18:08 +0200
Subject: [PATCH 03/16] Track floating point error

---
 src/Epochs/offsets.jl    | 31 ++++++++++++++++--------------
 src/Epochs/operations.jl |  4 ++--
 src/Epochs/types.jl      | 41 ++++++++++++++++++++++++++++------------
 test/epochs.jl           |  2 +-
 4 files changed, 49 insertions(+), 29 deletions(-)

diff --git a/src/Epochs/offsets.jl b/src/Epochs/offsets.jl
index a75f8fa..0403020 100644
--- a/src/Epochs/offsets.jl
+++ b/src/Epochs/offsets.jl
@@ -15,7 +15,7 @@ julia> TTEpoch(32.184, ep)
 ```
 """
 function Epoch{S2}(offset, ep::Epoch{S1}) where {S1<:TimeScale, S2<:TimeScale}
-    second, fraction = apply_offset(ep.second, ep.fraction, offset)
+    second, fraction = apply_offset(ep.second, ep.fraction, ep.error, offset)
     Epoch{S2}(second, fraction)
 end
 
@@ -36,8 +36,8 @@ julia> TAIEpoch(ep)
 ```
 """
 function Epoch{S2}(ep::Epoch{S1}) where {S1<:TimeScale, S2<:TimeScale}
-    second, fraction = apply_offset(ep.second, ep.fraction, S1(), S2())
-    Epoch{S2}(second, fraction)
+    second, fraction, error = apply_offset(ep.second, ep.fraction, ep.error, S1(), S2())
+    Epoch{S2}(second, fraction, error)
 end
 
 """
@@ -57,14 +57,14 @@ julia> Epoch(ep, TAI)
 ```
 """
 function Epoch(ep::Epoch{S1}, scale::S2) where {S1<:TimeScale, S2<:TimeScale}
-    second, fraction = apply_offset(ep.second, ep.fraction, S1(), S2())
-    Epoch{S2}(second, fraction)
+    second, fraction, error = apply_offset(ep.second, ep.fraction, ep.error, S1(), S2())
+    Epoch{S2}(second, fraction, error)
 end
 
 function Epoch{S2}(ep::Epoch{S1}, args...) where {S1<:TimeScale, S2<:TimeScale}
     offset = getoffset(S1(), S2(), ep.second, ep.fraction, args...)
-    second, fraction = apply_offset(ep.second, ep.fraction, offset)
-    Epoch{S2}(second, fraction)
+    second, fraction, error = apply_offset(ep.second, ep.fraction, ep.error, offset)
+    Epoch{S2}(second, fraction, error)
 end
 
 Epoch{S}(ep::Epoch{S}) where {S<:TimeScale} = ep
@@ -115,10 +115,11 @@ function getoffset(ep::Epoch{S}, scale::TimeScale) where S<:TimeScale
     total_offset = 0.0
     second = ep.second
     fraction = ep.fraction
+    error = ep.error
     for i in 1:length(path) - 1
         offset::Float64 = getoffset(path[i], path[i+1], second, fraction)
         total_offset += offset
-        second, fraction = apply_offset(second, fraction, offset)
+        second, fraction, error = apply_offset(second, fraction, error, offset)
     end
     return total_offset
 end
@@ -146,19 +147,21 @@ end
 
 @inline function apply_offset(second::Int64,
                               fraction::T,
+                              error::T,
                               from::S1,
-                              to::S2)::Tuple{Int64, T} where {T, S1<:TimeScale, S2<:TimeScale}
+                              to::S2)::Tuple{Int64, T, T} where {T, S1<:TimeScale, S2<:TimeScale}
     path = find_path(from, to)
     isempty(path) && throw(NoPathError(string(from), string(to)))
-    length(path) == 2 && return _apply_offset(second, fraction, from, to)
-    return _apply_offset((second, fraction), path...)
+    length(path) == 2 && return _apply_offset(second, fraction, error, from, to)
+    return _apply_offset((second, fraction, error), path...)
 end
 
 @inline function _apply_offset(second::Int64,
                                fraction::T,
+                               error::T,
                                from::S1,
-                               to::S2)::Tuple{Int64, T} where {T, S1<:TimeScale, S2<:TimeScale}
-    return apply_offset(second, fraction, getoffset(from, to, second, fraction))
+                               to::S2)::Tuple{Int64, T, T} where {T, S1<:TimeScale, S2<:TimeScale}
+    return apply_offset(second, fraction, error, getoffset(from, to, second, fraction))
 end
 
 @generated function _apply_offset(sf, path...)
@@ -554,7 +557,7 @@ julia> getoffset(TT, TDB, 0, 0.0, π, 6371.0, 0.0)
     offset = 0.0
     for _ in 1:3
         offset = -getoffset(TDB, TT, tdb1, tdb2, elong, u, v)
-        tdb1, tdb2 = apply_offset(tt1, tt2, offset)
+        tdb1, tdb2 = apply_offset(tt1, tt2, 0.0, offset)
     end
     return offset
 end
diff --git a/src/Epochs/operations.jl b/src/Epochs/operations.jl
index 12e1b95..898b790 100644
--- a/src/Epochs/operations.jl
+++ b/src/Epochs/operations.jl
@@ -25,7 +25,7 @@ julia> TAIEpoch(2018, 2, 6, 20, 45, 20.0) - TAIEpoch(2018, 2, 6, 20, 45, 0.0)
 """
 function Base.:-(a::Epoch{S}, b::Epoch{S}) where S<:TimeScale
     second = a.second - b.second
-    sum, residual = two_sum(a.fraction, -b.fraction)
-    return (residual + sum + second) * seconds
+    fraction = (a.error - b.error) + (a.fraction - b.fraction)
+    return (fraction + second) * seconds
 end
 
diff --git a/src/Epochs/types.jl b/src/Epochs/types.jl
index cd2bf15..4545755 100644
--- a/src/Epochs/types.jl
+++ b/src/Epochs/types.jl
@@ -6,33 +6,50 @@
     return hi, lo
 end
 
+@inline function two_hilo_sum(a, b)
+    hi = a + b
+    lo = b - (hi - a)
+    return hi, lo
+end
+
+function three_sum(a, b, c)
+    s, t   = two_sum(b, c)
+    hi, u  = two_sum(a, s)
+    md, lo = two_sum(u, t)
+    hi, md = two_hilo_sum(hi, md)
+    return hi, md, lo
+end
+
 struct Epoch{S<:TimeScale, T} <: Dates.AbstractDateTime
     scale::S
     second::Int64
     fraction::T
-    function Epoch{S}(second::Int64, fraction::T) where {S<:TimeScale, T}
-        return new{S, T}(S(), second, fraction)
+    error::T
+    function Epoch{S}(second::Int64, fraction::T, error::T=zero(T)) where {S<:TimeScale, T<:AbstractFloat}
+        return new{S, T}(S(), second, fraction, error)
     end
 end
 
 Epoch{S,T}(ep::Epoch{S,T}) where {S,T} = ep
 
-@inline function apply_offset(second::Int64, fraction, offset)
-    sum, residual = two_sum(fraction, offset)
-    if !isfinite(sum)
-        fraction′ = sum
-        second′ = ifelse(sum < 0, typemin(Int64), typemax(Int64))
+@inline function apply_offset(second::Int64, fraction, error, offset)
+    if !isfinite(fraction + offset)
+        fraction′ = fraction + offset
+        second′ = ifelse(fraction′ < 0, typemin(Int64), typemax(Int64))
+        residual = zero(fraction)
     else
+        offset_secs = floor(Int64, offset)
+        sum, residual, _ = three_sum(fraction, error, offset - offset_secs)
         int_secs = floor(Int64, sum)
-        second′ = second + int_secs
-        fraction′ = sum - int_secs + residual
+        second′ = second + int_secs + offset_secs
+        fraction′ = sum - int_secs
     end
-    return second′, fraction′
+    return second′, fraction′, residual
 end
 
 function Epoch{S}(ep::Epoch{S}, Δt) where {S<:TimeScale}
-    second, fraction = apply_offset(ep.second, ep.fraction, Δt)
-    return Epoch{S}(second, fraction)
+    second, fraction, error = apply_offset(ep.second, ep.fraction, ep.error, Δt)
+    return Epoch{S}(second, fraction, error)
 end
 
 Base.show(io::IO, ep::Epoch) = print(io, DateTime(ep), " ", timescale(ep))
diff --git a/test/epochs.jl b/test/epochs.jl
index aaa5104..62e1ab5 100644
--- a/test/epochs.jl
+++ b/test/epochs.jl
@@ -84,7 +84,7 @@ import ERFA
         @test t2 - t1 == -32.0seconds
         @test t2 < t1
 
-        today = TTEpoch(2000, 1, 1)
+        today = TTEpoch(2000, 1, 1, 12, 0, 13.123)
         age_of_the_universe = 13.772e9years
         big_bang = today - age_of_the_universe
         baryons = big_bang + 1e-11seconds

From 16fd26fa731333688c6dfd400a6314d17ea8e7dd Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Mon, 5 Jul 2021 14:58:36 +0200
Subject: [PATCH 04/16] Add AccurateArithmetic submodule

---
 src/AccurateArithmetic.jl | 30 ++++++++++++++++++++++++++++++
 src/AstroTime.jl          |  1 +
 src/Epochs/Epochs.jl      |  1 +
 src/Epochs/types.jl       | 22 ----------------------
 4 files changed, 32 insertions(+), 22 deletions(-)
 create mode 100644 src/AccurateArithmetic.jl

diff --git a/src/AccurateArithmetic.jl b/src/AccurateArithmetic.jl
new file mode 100644
index 0000000..6b89524
--- /dev/null
+++ b/src/AccurateArithmetic.jl
@@ -0,0 +1,30 @@
+module AccurateArithmetic
+
+# Adapted from AccurateArithmetic.jl
+
+export two_sum, two_hilo_sum, three_sum
+
+@inline function two_sum(a, b)
+    hi = a + b
+    a1 = hi - b
+    b1 = hi - a1
+    lo = (a - a1) + (b - b1)
+    return hi, lo
+end
+
+@inline function two_hilo_sum(a, b)
+    hi = a + b
+    lo = b - (hi - a)
+    return hi, lo
+end
+
+function three_sum(a, b, c)
+    s, t   = two_sum(b, c)
+    hi, u  = two_sum(a, s)
+    md, lo = two_sum(u, t)
+    hi, md = two_hilo_sum(hi, md)
+    return hi, md, lo
+end
+
+end
+
diff --git a/src/AstroTime.jl b/src/AstroTime.jl
index 6022329..7d7d2fe 100644
--- a/src/AstroTime.jl
+++ b/src/AstroTime.jl
@@ -40,6 +40,7 @@ const EPOCH_ISO_FORMAT = Ref{Dates.DateFormat{Symbol("yyyy-mm-ddTHH:MM:SS.fff tt
                                                     Dates.Delim{Char, 1},
                                                     Dates.DatePart{'t'}}}}()
 
+include("AccurateArithmetic.jl")
 include("TimeScales.jl")
 include("Periods.jl")
 include("AstroDates.jl")
diff --git a/src/Epochs/Epochs.jl b/src/Epochs/Epochs.jl
index 6cc5f2e..d95846c 100644
--- a/src/Epochs/Epochs.jl
+++ b/src/Epochs/Epochs.jl
@@ -19,6 +19,7 @@ import ..AstroDates: Date, DateTime, Time
 import ..AstroDates: calendar, fractionofday, fractionofsecond, secondinday, subsecond
 import ..AstroDates: j2000, julian, julian_twopart
 
+using ..AccurateArithmetic: three_sum, two_sum
 using ..TimeScales: find_path
 
 export CCSDS_EPOCH, FIFTIES_EPOCH, FUTURE_INFINITY, GALILEO_EPOCH, GPS_EPOCH, J2000_EPOCH
diff --git a/src/Epochs/types.jl b/src/Epochs/types.jl
index 4545755..00f7640 100644
--- a/src/Epochs/types.jl
+++ b/src/Epochs/types.jl
@@ -1,25 +1,3 @@
-@inline function two_sum(a, b)
-    hi = a + b
-    a1 = hi - b
-    b1 = hi - a1
-    lo = (a - a1) + (b - b1)
-    return hi, lo
-end
-
-@inline function two_hilo_sum(a, b)
-    hi = a + b
-    lo = b - (hi - a)
-    return hi, lo
-end
-
-function three_sum(a, b, c)
-    s, t   = two_sum(b, c)
-    hi, u  = two_sum(a, s)
-    md, lo = two_sum(u, t)
-    hi, md = two_hilo_sum(hi, md)
-    return hi, md, lo
-end
-
 struct Epoch{S<:TimeScale, T} <: Dates.AbstractDateTime
     scale::S
     second::Int64

From 9a5dc1ae8a462bf3f16dcfaca04102805aa4a461 Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Thu, 22 Jul 2021 11:16:00 +0200
Subject: [PATCH 05/16] Refactor offset calculations

---
 src/AccurateArithmetic.jl | 43 ++++++++++++++++++++++++++++++++-------
 src/Epochs/Epochs.jl      |  2 +-
 src/Epochs/types.jl       | 15 --------------
 3 files changed, 37 insertions(+), 23 deletions(-)

diff --git a/src/AccurateArithmetic.jl b/src/AccurateArithmetic.jl
index 6b89524..e295a1f 100644
--- a/src/AccurateArithmetic.jl
+++ b/src/AccurateArithmetic.jl
@@ -2,7 +2,7 @@ module AccurateArithmetic
 
 # Adapted from AccurateArithmetic.jl
 
-export two_sum, two_hilo_sum, three_sum
+export two_sum, apply_offset
 
 @inline function two_sum(a, b)
     hi = a + b
@@ -18,12 +18,41 @@ end
     return hi, lo
 end
 
-function three_sum(a, b, c)
-    s, t   = two_sum(b, c)
-    hi, u  = two_sum(a, s)
-    md, lo = two_sum(u, t)
-    hi, md = two_hilo_sum(hi, md)
-    return hi, md, lo
+function four_sum(a, b, c, d)
+    t0, t1 = two_sum(a,  b)
+    t2, t3 = two_sum(c,  d)
+    hi, t4 = two_sum(t0, t2)
+    t5, lo = two_sum(t1, t3)
+    hm, ml = two_sum(t4, t5)
+    ml, lo = two_hilo_sum(ml, lo)
+    hm, ml = two_hilo_sum(hm, ml)
+    hi, hm = two_hilo_sum(hi,hm)
+    return hi, hm, ml, lo
+end
+
+function handle_infinity(a, b)
+    fraction = a + b
+    second = ifelse(fraction < 0, typemin(Int64), typemax(Int64))
+    residual = zero(fraction)
+    return second, fraction, residual
+end
+
+function apply_offset(s1::Int64, f1, e1, s2::Int64, f2, e2)
+    isfinite(f1 + f2) || return handle_infinity(f1, f2)
+
+    sum, residual, _ = four_sum(f1, f2, e1, e2)
+    int_seconds = floor(Int64, sum)
+    second = s1 + s2 + int_seconds
+    fraction = sum - int_seconds
+    return second, fraction, residual
+end
+
+function apply_offset(s1::Int64, f1, e1, offset)
+    isfinite(f1 + offset) || return handle_infinity(f1, offset)
+
+    s2 = floor(Int64, offset)
+    f2 = offset - s2
+    return apply_offset(s1, f1, e1, s2, f2, 0.0)
 end
 
 end
diff --git a/src/Epochs/Epochs.jl b/src/Epochs/Epochs.jl
index d95846c..cebb587 100644
--- a/src/Epochs/Epochs.jl
+++ b/src/Epochs/Epochs.jl
@@ -18,8 +18,8 @@ import ..EPOCH_ISO_FORMAT
 import ..AstroDates: Date, DateTime, Time
 import ..AstroDates: calendar, fractionofday, fractionofsecond, secondinday, subsecond
 import ..AstroDates: j2000, julian, julian_twopart
+import ..AccurateArithmetic: apply_offset, two_sum
 
-using ..AccurateArithmetic: three_sum, two_sum
 using ..TimeScales: find_path
 
 export CCSDS_EPOCH, FIFTIES_EPOCH, FUTURE_INFINITY, GALILEO_EPOCH, GPS_EPOCH, J2000_EPOCH
diff --git a/src/Epochs/types.jl b/src/Epochs/types.jl
index 00f7640..f796db6 100644
--- a/src/Epochs/types.jl
+++ b/src/Epochs/types.jl
@@ -10,21 +10,6 @@ end
 
 Epoch{S,T}(ep::Epoch{S,T}) where {S,T} = ep
 
-@inline function apply_offset(second::Int64, fraction, error, offset)
-    if !isfinite(fraction + offset)
-        fraction′ = fraction + offset
-        second′ = ifelse(fraction′ < 0, typemin(Int64), typemax(Int64))
-        residual = zero(fraction)
-    else
-        offset_secs = floor(Int64, offset)
-        sum, residual, _ = three_sum(fraction, error, offset - offset_secs)
-        int_secs = floor(Int64, sum)
-        second′ = second + int_secs + offset_secs
-        fraction′ = sum - int_secs
-    end
-    return second′, fraction′, residual
-end
-
 function Epoch{S}(ep::Epoch{S}, Δt) where {S<:TimeScale}
     second, fraction, error = apply_offset(ep.second, ep.fraction, ep.error, Δt)
     return Epoch{S}(second, fraction, error)

From d53ab59e7678bf2686fc3b44d5c1c8c3c52c443b Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Thu, 22 Jul 2021 12:45:22 +0200
Subject: [PATCH 06/16] Redesign periods to increase precision

---
 docs/src/tutorial.md      |  10 +-
 src/AccurateArithmetic.jl |   7 +-
 src/Epochs/aliases.jl     |   2 +-
 src/Epochs/julian.jl      |   6 +-
 src/Epochs/operations.jl  |   7 +-
 src/Epochs/ranges.jl      |   2 +-
 src/Periods.jl            | 214 ++++++++++++++++++--------------------
 test/epochs.jl            |   9 +-
 test/periods.jl           |  70 ++++---------
 9 files changed, 146 insertions(+), 181 deletions(-)

diff --git a/docs/src/tutorial.md b/docs/src/tutorial.md
index 5f13f18..9a5b7d8 100644
--- a/docs/src/tutorial.md
+++ b/docs/src/tutorial.md
@@ -83,7 +83,7 @@ julia> AstroTime.format(ep, "dd.mm.yyyy HH:MM ttt")
 
 ## Working with Epochs and Periods
 
-You can shift an `Epoch` in time by adding or subtracting a [`Period`](@ref) to it.
+You can shift an `Epoch` in time by adding or subtracting an [`AstroPeriod`](@ref) to it.
 
 AstroTime.jl provides a convenient way to construct periods by multiplying a value
 with a time unit.
@@ -105,7 +105,7 @@ The following time units are available:
 - `years`
 - `centuries`
 
-To shift an `Epoch` forward in time add a `Period` to it.
+To shift an `Epoch` forward in time add an `AstroPeriod` to it.
 
 ```julia
 julia> ep = UTCEpoch(2000, 1, 1)
@@ -125,7 +125,7 @@ julia> ep - 1days
 1999-12-31T00:00:00.000 UTC
 ```
 
-If you subtract two epochs you will receive the time between them as a `Period`.
+If you subtract two epochs you will receive the time between them as an `AstroPeriod`.
 
 ```julia
 julia> ep1 = UTCEpoch(2000, 1, 1)
@@ -138,8 +138,8 @@ julia> ep2 - ep1
 86400.0 seconds
 ```
 
-You can also construct a `Period` with a different time unit from
-another `Period`.
+You can also construct an `AstroPeriod` with a different time unit from
+another `AstroPeriod`.
 
 ```julia
 julia> dt = 86400.0seconds
diff --git a/src/AccurateArithmetic.jl b/src/AccurateArithmetic.jl
index e295a1f..0355332 100644
--- a/src/AccurateArithmetic.jl
+++ b/src/AccurateArithmetic.jl
@@ -30,15 +30,14 @@ function four_sum(a, b, c, d)
     return hi, hm, ml, lo
 end
 
-function handle_infinity(a, b)
-    fraction = a + b
+function handle_infinity(fraction)
     second = ifelse(fraction < 0, typemin(Int64), typemax(Int64))
     residual = zero(fraction)
     return second, fraction, residual
 end
 
 function apply_offset(s1::Int64, f1, e1, s2::Int64, f2, e2)
-    isfinite(f1 + f2) || return handle_infinity(f1, f2)
+    isfinite(f1 + f2) || return handle_infinity(f1 + f2)
 
     sum, residual, _ = four_sum(f1, f2, e1, e2)
     int_seconds = floor(Int64, sum)
@@ -48,7 +47,7 @@ function apply_offset(s1::Int64, f1, e1, s2::Int64, f2, e2)
 end
 
 function apply_offset(s1::Int64, f1, e1, offset)
-    isfinite(f1 + offset) || return handle_infinity(f1, offset)
+    isfinite(f1 + offset) || return handle_infinity(f1 + offset)
 
     s2 = floor(Int64, offset)
     f2 = offset - s2
diff --git a/src/Epochs/aliases.jl b/src/Epochs/aliases.jl
index f5fc154..4696ebf 100644
--- a/src/Epochs/aliases.jl
+++ b/src/Epochs/aliases.jl
@@ -32,7 +32,7 @@ for (scale, acronym) in zip(TimeScales.NAMES, TimeScales.ACRONYMS)
         $epoch(::AbstractString)
 
         """
-            $($name)(jd1::T, jd2::T=zero(T); origin=:j2000) where T<:Period
+            $($name)(jd1::T, jd2::T=zero(T); origin=:j2000) where T<:AstroPeriod
 
         Construct a $($name) from a Julian date (optionally split into
         `jd1` and `jd2`). `origin` determines the variant of Julian
diff --git a/src/Epochs/julian.jl b/src/Epochs/julian.jl
index 7845df3..0e7922e 100644
--- a/src/Epochs/julian.jl
+++ b/src/Epochs/julian.jl
@@ -2,7 +2,7 @@ const J2000_TO_JULIAN = 2.451545e6days
 const J2000_TO_MJD = 51544.5days
 
 """
-    Epoch{S}(jd1::T, jd2::T=zero(T); origin=:j2000) where {S, T<:Period}
+    Epoch{S}(jd1::T, jd2::T=zero(T); origin=:j2000) where {S, T<:AstroPeriod}
 
 Construct an `Epoch` with time scale `S` from a Julian date
 (optionally split into `jd1` and `jd2`). `origin` determines the
@@ -22,7 +22,7 @@ julia> Epoch{InternationalAtomicTime}(2.451545e6days, origin=:julian)
 2000-01-01T12:00:00.000 TAI
 ```
 """
-function Epoch{S}(jd1::T, jd2::T=zero(T), args...; origin=:j2000) where {S, T<:Period}
+function Epoch{S}(jd1::T, jd2::T=zero(T), args...; origin=:j2000) where {S, T<:AstroPeriod}
     if jd2 > jd1
         jd1, jd2 = jd2, jd1
     end
@@ -52,7 +52,7 @@ end
     julian_period([T,] ep::Epoch; origin=:j2000, scale=timescale(ep), unit=days)
 
 Return the period since Julian Epoch `origin` within the time scale `scale` expressed in
-`unit` for a given epoch `ep`. The result is a [`Period`](@ref) object by default.
+`unit` for a given epoch `ep`. The result is an [`AstroPeriod`](@ref) object by default.
 If the type argument `T` is present, the result is converted to `T` instead.
 
 ### Example ###
diff --git a/src/Epochs/operations.jl b/src/Epochs/operations.jl
index 898b790..a6c0d6c 100644
--- a/src/Epochs/operations.jl
+++ b/src/Epochs/operations.jl
@@ -8,8 +8,11 @@ end
 
 Base.isless(ep1::Epoch, ep2::Epoch) = isless(value(ep1 - ep2), 0.0)
 
-Base.:+(ep::Epoch{S}, p::Period) where {S} = Epoch{S}(ep, value(seconds(p)))
-Base.:-(ep::Epoch{S}, p::Period) where {S} = Epoch{S}(ep, -value(seconds(p)))
+function Base.:+(ep::Epoch{S}, p::AstroPeriod) where {S}
+    second, fraction, error = apply_offset(ep.second, ep.fraction, ep.error, p.second, p.fraction, p.error)
+    return Epoch{S}(second, fraction, error)
+end
+Base.:-(ep::Epoch, p::AstroPeriod) = ep + (-p)
 
 """
     -(a::Epoch, b::Epoch)
diff --git a/src/Epochs/ranges.jl b/src/Epochs/ranges.jl
index f8b8d76..31bd60e 100644
--- a/src/Epochs/ranges.jl
+++ b/src/Epochs/ranges.jl
@@ -1,6 +1,6 @@
 (::Base.Colon)(start::Epoch{S}, stop::Epoch{S}) where {S} = (:)(start, 1.0days, stop)
 
-function (::Base.Colon)(start::Epoch{S}, step::Period, stop::Epoch{S}) where S
+function (::Base.Colon)(start::Epoch{S}, step::AstroPeriod, stop::Epoch{S}) where S
     step = seconds(step)
     step = start < stop ? step : -step
     StepRangeLen(start, step, floor(Int, value(stop-start)/value(step))+1)
diff --git a/src/Periods.jl b/src/Periods.jl
index 94d5b44..9c6fe3e 100644
--- a/src/Periods.jl
+++ b/src/Periods.jl
@@ -1,23 +1,16 @@
 module Periods
 
+using ..AccurateArithmetic: apply_offset, handle_infinity
+
 export TimeUnit, Second, Minute, Hour, Day, Year, Century,
     seconds, minutes, hours, days, years, centuries,
-    Period, -, *, /, +, value, unit,
+    AstroPeriod,
+    -, *, /, +, value, unit,
     SECONDS_PER_MINUTE,
     SECONDS_PER_HOUR,
     SECONDS_PER_DAY,
     SECONDS_PER_YEAR,
-    SECONDS_PER_CENTURY,
-    MINUTES_PER_HOUR,
-    MINUTES_PER_DAY,
-    MINUTES_PER_YEAR,
-    MINUTES_PER_CENTURY,
-    HOURS_PER_DAY,
-    HOURS_PER_YEAR,
-    HOURS_PER_CENTURY,
-    DAYS_PER_YEAR,
-    DAYS_PER_CENTURY,
-    YEARS_PER_CENTURY
+    SECONDS_PER_CENTURY
 
 const SECONDS_PER_MINUTE   = 60.0
 const SECONDS_PER_HOUR     = 60.0 * 60.0
@@ -25,20 +18,6 @@ const SECONDS_PER_DAY      = 60.0 * 60.0 * 24.0
 const SECONDS_PER_YEAR     = 60.0 * 60.0 * 24.0 * 365.25
 const SECONDS_PER_CENTURY  = 60.0 * 60.0 * 24.0 * 365.25 * 100.0
 
-const MINUTES_PER_HOUR     = 60.0
-const MINUTES_PER_DAY      = 60.0 * 24.0
-const MINUTES_PER_YEAR     = 60.0 * 24.0 * 365.25
-const MINUTES_PER_CENTURY  = 60.0 * 24.0 * 365.25 * 100.0
-
-const HOURS_PER_DAY        = 24.0
-const HOURS_PER_YEAR       = 24.0 * 365.25
-const HOURS_PER_CENTURY    = 24.0 * 365.25 * 100.0
-
-const DAYS_PER_YEAR        = 365.25
-const DAYS_PER_CENTURY     = 365.25 * 100.0
-
-const YEARS_PER_CENTURY    = 100.0
-
 """
 All time units are subtypes of the abstract type `TimeUnit`.
 The following time units are defined:
@@ -68,11 +47,18 @@ const centuries = Century()
 
 Base.broadcastable(u::TimeUnit) = Ref(u)
 
+factor(::Second) = 1.0
+factor(::Minute) = SECONDS_PER_MINUTE
+factor(::Hour) = SECONDS_PER_HOUR
+factor(::Day) = SECONDS_PER_DAY
+factor(::Year) = SECONDS_PER_YEAR
+factor(::Century) = SECONDS_PER_CENTURY
+
 """
-    Period{U, T}(unit, Δt) where {U<:TimeUnit, T}
+    AstroPeriod{U, T}(unit, Δt) where {U<:TimeUnit, T}
 
-A `Period` object represents a time interval of `Δt` with a [`TimeUnit`](@ref) of `unit`.
-Periods should be constructed via the shorthand syntax shown in the examples below.
+An `AstroPeriod` object represents a time interval of `Δt` with a [`TimeUnit`](@ref) of
+`unit`.  Periods should be constructed via the shorthand syntax shown in the examples below.
 
 ### Examples ###
 
@@ -84,118 +70,118 @@ julia> 1.0minutes
 1.0 minutes
 
 julia> 12hours
-12 hours
+12.0 hours
 
 julia> days_per_year = 365
 365
 julia> days_per_year * days
-365 days
+365.0 days
 
 julia> 10.0years
 10.0 years
 
 julia> 1centuries
-1 century
+1.0 centuries
 ```
 """
-struct Period{U<:TimeUnit, T}
+struct AstroPeriod{U<:TimeUnit, T}
     unit::U
-    Δt::T
+    second::Int64
+    fraction::T
+    error::T
+end
+
+function AstroPeriod(unit, dt)
+    isfinite(dt) || return AstroPeriod(unit, handle_infinity(dt)...)
+
+    seconds = dt * factor(unit)
+    int_seconds = floor(Int64, seconds)
+    fraction = seconds - int_seconds
+    return AstroPeriod(unit, int_seconds, fraction, zero(fraction))
 end
 
-value(p::Period) = p.Δt
-unit(p::Period) = p.unit
-Base.zero(p::Period) = Period(unit(p), zero(value(p)))
-Base.zero(p::Type{<:Period{U}}) where {U} = Period(U(), 0.0)
-Base.zero(p::Type{<:Period{U,T}}) where {U, T} = Period(U(), zero(T))
-Base.eltype(p::Period) = typeof(value(p))
-Base.eltype(p::Type{<:Period{U,T}}) where {U, T} = T
-
-name(::Second, val::Integer) = ifelse(val == one(val), "second", "seconds")
-name(::Second, ::Any) = "seconds"
-name(::Minute, val::Integer) = ifelse(val == one(val), "minute", "minutes")
-name(::Minute, ::Any) = "minutes"
-name(::Hour, val::Integer) = ifelse(val == one(val), "hour", "hours")
-name(::Hour, ::Any) = "hours"
-name(::Day, val::Integer) = ifelse(val == one(val), "day", "days")
-name(::Day, ::Any) = "days"
-name(::Year, val::Integer) = ifelse(val == one(val), "year", "years")
-name(::Year, ::Any) = "years"
-name(::Century, val::Integer) = ifelse(val == one(val), "century", "centuries")
-name(::Century, ::Any) = "centuries"
-
-function Base.show(io::IO, p::Period)
+
+(u::TimeUnit)(p::AstroPeriod) = AstroPeriod(u, p.second, p.fraction, p.error)
+
+"""
+    unit(p::AstroPeriod)
+
+Return the unit of the period `p`.
+
+### Examples ###
+
+```jldoctest; setup = :(using AstroTime)
+julia> unit(3.0seconds)
+Second()
+```
+"""
+unit(p::AstroPeriod) = p.unit
+
+"""
+    value(p::AstroPeriod)
+
+Return the unitless value of the period `p`.
+
+### Examples ###
+
+```jldoctest; setup = :(using AstroTime)
+julia> value(3.0seconds)
+3.0
+```
+"""
+value(p::AstroPeriod) = (p.fraction + p.second) / factor(unit(p))
+
+Base.zero(p::AstroPeriod) = AstroPeriod(unit(p), zero(value(p)))
+Base.zero(p::Type{<:AstroPeriod{U}}) where {U} = AstroPeriod(U(), 0.0)
+Base.zero(p::Type{<:AstroPeriod{U,T}}) where {U, T} = AstroPeriod(U(), zero(T))
+Base.eltype(p::AstroPeriod) = typeof(value(p))
+Base.eltype(p::Type{<:AstroPeriod{U,T}}) where {U, T} = T
+
+name(::Second) = "seconds"
+name(::Minute) = "minutes"
+name(::Hour) = "hours"
+name(::Day) = "days"
+name(::Year) = "years"
+name(::Century) = "centuries"
+
+function Base.show(io::IO, p::AstroPeriod)
     u = unit(p)
     v = value(p)
-    print(io, v, " ", name(u, v))
+    print(io, v, " ", name(u))
 end
 
-Base.:*(Δt::T, unit::TimeUnit) where {T<:Number} = Period(unit, Δt)
-Base.:*(unit::TimeUnit, Δt::T) where {T<:Number} = Period(unit, Δt)
+Base.:*(dt::Number, unit::TimeUnit) = AstroPeriod(unit, dt)
+Base.:*(unit::TimeUnit, dt::Number) = AstroPeriod(unit, dt)
 Base.:*(A::TimeUnit, B::AbstractArray) = broadcast(*, A, B)
 Base.:*(A::AbstractArray, B::TimeUnit) = broadcast(*, A, B)
 
-Base.:+(p1::Period{U}, p2::Period{U}) where {U}= Period(unit(p1), p1.Δt + p2.Δt)
-Base.:-(p1::Period{U}, p2::Period{U}) where {U}= Period(unit(p1), p1.Δt - p2.Δt)
-Base.:-(p::Period) = Period(unit(p), -p.Δt)
-Base.:*(x, p::Period) = Period(unit(p), p.Δt * x)
-Base.:*(p::Period, x) = Period(unit(p), p.Δt * x)
-Base.:/(p::Period, x) = Period(unit(p), p.Δt / x)
+Base.:-(p::AstroPeriod) = AstroPeriod(unit(p), -p.second, -p.fraction, -p.error)
+
+function Base.:+(p1::AstroPeriod{U}, p2::AstroPeriod{U}) where U
+    second, fraction, error = apply_offset(p1.second, p1.fraction, p1.error, p2.second, p2.fraction, p2.error)
+    return AstroPeriod(U(), second, fraction, error)
+end
+
+Base.:-(p1::AstroPeriod, p2::AstroPeriod) = p1 + (-p2)
+Base.:*(x, p::AstroPeriod) = AstroPeriod(unit(p), value(p) * x)
+Base.:*(p::AstroPeriod, x) = AstroPeriod(unit(p), value(p) * x)
+Base.:/(p::AstroPeriod, x) = AstroPeriod(unit(p), value(p) / x)
 
-Base.isless(p1::Period{U}, p2::Period{U}) where {U} = isless(value(p1), value(p2))
-Base.isapprox(p1::Period{U}, p2::Period{U}) where {U} = value(p1) ≈ value(p2)
+Base.isless(p1::AstroPeriod{U}, p2::AstroPeriod{U}) where {U} = isless(value(p1), value(p2))
+Base.:(==)(p1::AstroPeriod{U}, p2::AstroPeriod{U}) where {U} = value(p1) == value(p2)
+function Base.isapprox(p1::AstroPeriod{U}, p2::AstroPeriod{U}; kwargs...) where {U}
+    return isapprox(value(p1), value(p2); kwargs...)
+end
 
-(::Base.Colon)(start::Period{U,T}, stop::Period{U,T}) where {U,T} = (:)(start, one(T) * U(), stop)
+(::Base.Colon)(start::AstroPeriod{U,T}, stop::AstroPeriod{U,T}) where {U,T} = (:)(start, one(T) * U(), stop)
 
-function (::Base.Colon)(start::Period{U}, step::Period{U}, stop::Period{U}) where {U}
+function (::Base.Colon)(start::AstroPeriod{U}, step::AstroPeriod{U}, stop::AstroPeriod{U}) where {U}
     step = start < stop ? step : -step
     StepRangeLen(start, step, floor(Int, value(stop-start)/value(step))+1)
 end
 
-Period{U,T}(p::Period{U,T}) where {U,T} = p
-
-Base.step(r::StepRangeLen{<:Period}) = r.step
-
-(::Second)(p::Period{Second})   = p
-(::Second)(p::Period{Minute})   = Period(seconds, p.Δt * SECONDS_PER_MINUTE)
-(::Second)(p::Period{Hour})     = Period(seconds, p.Δt * SECONDS_PER_HOUR)
-(::Second)(p::Period{Day})      = Period(seconds, p.Δt * SECONDS_PER_DAY)
-(::Second)(p::Period{Year})     = Period(seconds, p.Δt * SECONDS_PER_YEAR)
-(::Second)(p::Period{Century})  = Period(seconds, p.Δt * SECONDS_PER_CENTURY)
-
-(::Minute)(p::Period{Second})   = Period(minutes, p.Δt / SECONDS_PER_MINUTE)
-(::Minute)(p::Period{Minute})   = p
-(::Minute)(p::Period{Hour})     = Period(minutes, p.Δt * MINUTES_PER_HOUR)
-(::Minute)(p::Period{Day})      = Period(minutes, p.Δt * MINUTES_PER_DAY)
-(::Minute)(p::Period{Year})     = Period(minutes, p.Δt * MINUTES_PER_YEAR)
-(::Minute)(p::Period{Century})  = Period(minutes, p.Δt * MINUTES_PER_CENTURY)
-
-(::Hour)(p::Period{Second})     = Period(hours, p.Δt / SECONDS_PER_HOUR)
-(::Hour)(p::Period{Minute})     = Period(hours, p.Δt / MINUTES_PER_HOUR)
-(::Hour)(p::Period{Hour})       = p
-(::Hour)(p::Period{Day})        = Period(hours, p.Δt * HOURS_PER_DAY)
-(::Hour)(p::Period{Year})       = Period(hours, p.Δt * HOURS_PER_YEAR)
-(::Hour)(p::Period{Century})    = Period(hours, p.Δt * HOURS_PER_CENTURY)
-
-(::Day)(p::Period{Second})      = Period(days, p.Δt / SECONDS_PER_DAY)
-(::Day)(p::Period{Minute})      = Period(days, p.Δt / MINUTES_PER_DAY)
-(::Day)(p::Period{Hour})        = Period(days, p.Δt / HOURS_PER_DAY)
-(::Day)(p::Period{Day})         = p
-(::Day)(p::Period{Year})        = Period(days, p.Δt * DAYS_PER_YEAR)
-(::Day)(p::Period{Century})     = Period(days, p.Δt * DAYS_PER_CENTURY)
-
-(::Year)(p::Period{Second})     = Period(years, p.Δt / SECONDS_PER_YEAR)
-(::Year)(p::Period{Minute})     = Period(years, p.Δt / MINUTES_PER_YEAR)
-(::Year)(p::Period{Hour})       = Period(years, p.Δt / HOURS_PER_YEAR)
-(::Year)(p::Period{Day})        = Period(years, p.Δt / DAYS_PER_YEAR)
-(::Year)(p::Period{Year})       = p
-(::Year)(p::Period{Century})    = Period(years, p.Δt * YEARS_PER_CENTURY)
-
-(::Century)(p::Period{Second})  = Period(centuries, p.Δt / SECONDS_PER_CENTURY)
-(::Century)(p::Period{Minute})  = Period(centuries, p.Δt / MINUTES_PER_CENTURY)
-(::Century)(p::Period{Hour})    = Period(centuries, p.Δt / HOURS_PER_CENTURY)
-(::Century)(p::Period{Day})     = Period(centuries, p.Δt / DAYS_PER_CENTURY)
-(::Century)(p::Period{Year})    = Period(centuries, p.Δt / YEARS_PER_CENTURY)
-(::Century)(p::Period{Century}) = p
+AstroPeriod{U,T}(p::AstroPeriod{U,T}) where {U,T} = p
+
+Base.step(r::StepRangeLen{<:AstroPeriod}) = r.step
 
 end
diff --git a/test/epochs.jl b/test/epochs.jl
index 62e1ab5..8d2bd91 100644
--- a/test/epochs.jl
+++ b/test/epochs.jl
@@ -89,6 +89,13 @@ import ERFA
         big_bang = today - age_of_the_universe
         baryons = big_bang + 1e-11seconds
         @test baryons + age_of_the_universe - today == 1e-11seconds
+
+        reception_time = TDBEpoch("2021-07-01T00:00:00.00")
+        rtlt_a = seconds(1.5days)
+        rtlt_b = rtlt_a + 1e-6seconds
+        transmission_time_a = reception_time + rtlt_a
+        transmission_time_b = reception_time + rtlt_b
+        @test transmission_time_b - transmission_time_a == 1e-6seconds
     end
     @testset "Parsing" begin
         @test AstroTime.TimeScales.tryparse(1.0) === nothing
@@ -239,7 +246,7 @@ import ERFA
         ep_f64 = TAIEpoch(2000, 1, 1)
         ep_err = TAIEpoch(ep_f64.second, 1.0 ± 1.1)
         Δt = (30 ± 0.1) * seconds
-        @test typeof(Δt) == Period{Second,Measurement{Float64}}
+        @test typeof(Δt) == AstroPeriod{Second,Measurement{Float64}}
         @test typeof(ep_f64) == Epoch{InternationalAtomicTime,Float64}
         @test typeof(ep_err) == Epoch{InternationalAtomicTime,Measurement{Float64}}
         @test typeof(ep_f64 + Δt) == Epoch{InternationalAtomicTime,Measurement{Float64}}
diff --git a/test/periods.jl b/test/periods.jl
index bc314a2..eb9a5ef 100644
--- a/test/periods.jl
+++ b/test/periods.jl
@@ -5,12 +5,12 @@
     d = 1.0days
     y = 1.0years
     c = 1.0centuries
-    @test s == Period(seconds, 1.0)
-    @test m == Period(minutes, 1.0)
-    @test h == Period(hours, 1.0)
-    @test d == Period(days, 1.0)
-    @test y == Period(years, 1.0)
-    @test c == Period(centuries, 1.0)
+    @test s == AstroPeriod(seconds, 1.0)
+    @test m == AstroPeriod(minutes, 1.0)
+    @test h == AstroPeriod(hours, 1.0)
+    @test d == AstroPeriod(days, 1.0)
+    @test y == AstroPeriod(years, 1.0)
+    @test c == AstroPeriod(centuries, 1.0)
 
     @test seconds(s) == 1.0seconds
     @test seconds(m) == 60.0seconds
@@ -54,7 +54,7 @@
     @test centuries(y) == (1.0 / 100.0)centuries
     @test centuries(c) == 1.0centuries
 
-    @test zero(Period{Year}) == 0.0years
+    @test zero(AstroPeriod{Year}) == 0.0years
     @test zero(1years) == 0years
     @test zero(1.0years) == 0.0years
 
@@ -83,50 +83,20 @@
     float_rng = 1.0seconds:3.0seconds
     @test step(float_rng) == 1.0seconds
     @test collect(float_rng) == [1.0seconds, 2.0seconds, 3.0seconds]
-    @test Period{Second,Float64}(1.0seconds) == 1.0seconds
+    @test AstroPeriod{Second,Float64}(1.0seconds) == 1.0seconds
 
-    @test Periods.name(seconds, 1) == "second"
-    @test Periods.name(seconds, 2) == "seconds"
-    @test Periods.name(seconds, 1.0) == "seconds"
-    @test Periods.name(minutes, 1) == "minute"
-    @test Periods.name(minutes, 2) == "minutes"
-    @test Periods.name(minutes, 1.0) == "minutes"
-    @test Periods.name(hours, 1) == "hour"
-    @test Periods.name(hours, 2) == "hours"
-    @test Periods.name(hours, 1.0) == "hours"
-    @test Periods.name(days, 1) == "day"
-    @test Periods.name(days, 2) == "days"
-    @test Periods.name(days, 1.0) == "days"
-    @test Periods.name(years, 1) == "year"
-    @test Periods.name(years, 2) == "years"
-    @test Periods.name(years, 1.0) == "years"
-    @test Periods.name(centuries, 1) == "century"
-    @test Periods.name(centuries, 2) == "centuries"
-    @test Periods.name(centuries, 1.0) == "centuries"
+    @test Periods.name(seconds) == "seconds"
+    @test Periods.name(minutes) == "minutes"
+    @test Periods.name(hours) == "hours"
+    @test Periods.name(days) == "days"
+    @test Periods.name(years) == "years"
+    @test Periods.name(centuries) == "centuries"
 
-    @test string(1seconds) == "1 second"
-    @test string(Int32(1) * seconds) == "1 second"
-    @test string(1.0seconds) == "1.0 seconds"
-    @test string(2seconds) == "2 seconds"
-    @test string(1minutes) == "1 minute"
-    @test string(Int32(1) * minutes) == "1 minute"
-    @test string(1.0minutes) == "1.0 minutes"
-    @test string(2minutes) == "2 minutes"
-    @test string(1hours) == "1 hour"
-    @test string(Int32(1) * hours) == "1 hour"
-    @test string(1.0hours) == "1.0 hours"
-    @test string(2hours) == "2 hours"
-    @test string(1days) == "1 day"
-    @test string(Int32(1) * days) == "1 day"
-    @test string(1.0days) == "1.0 days"
-    @test string(2days) == "2 days"
-    @test string(1years) == "1 year"
-    @test string(Int32(1) * years) == "1 year"
-    @test string(1.0years) == "1.0 years"
-    @test string(2years) == "2 years"
-    @test string(1centuries) == "1 century"
-    @test string(Int32(1) * centuries) == "1 century"
-    @test string(1.0centuries) == "1.0 centuries"
-    @test string(2centuries) == "2 centuries"
+    @test string(1seconds) == "1.0 seconds"
+    @test string(1minutes) == "1.0 minutes"
+    @test string(1hours) == "1.0 hours"
+    @test string(1days) == "1.0 days"
+    @test string(1years) == "1.0 years"
+    @test string(1centuries) == "1.0 centuries"
 end
 

From 4e0960f87127c0b6e303611aac3660bb779bdc26 Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Thu, 22 Jul 2021 12:57:59 +0200
Subject: [PATCH 07/16] Update CI

---
 .github/workflows/CI.yml | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml
index cb0e0b7..8280020 100644
--- a/.github/workflows/CI.yml
+++ b/.github/workflows/CI.yml
@@ -2,10 +2,10 @@ name: CI
 on:
   pull_request:
     branches:
-      - master
+      - main
   push:
     branches:
-      - master
+      - main
     tags: '*'
 jobs:
   test:

From cf067db47eab97ada3026b5d3e3625f40ac3487e Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Thu, 22 Jul 2021 13:03:30 +0200
Subject: [PATCH 08/16] Bump minimum Julia version to 1.3

---
 .github/workflows/CI.yml | 2 +-
 Project.toml             | 2 +-
 2 files changed, 2 insertions(+), 2 deletions(-)

diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml
index 8280020..775efab 100644
--- a/.github/workflows/CI.yml
+++ b/.github/workflows/CI.yml
@@ -15,7 +15,7 @@ jobs:
       fail-fast: false
       matrix:
         version:
-          - '1.0' # Replace this with the minimum Julia version that your package supports. E.g. if your package requires Julia 1.5 or higher, change this to '1.5'.
+          - '1.3' # Replace this with the minimum Julia version that your package supports. E.g. if your package requires Julia 1.5 or higher, change this to '1.5'.
           - '1' # Leave this line unchanged. '1' will automatically expand to the latest stable 1.x release of Julia.
           - 'nightly'
         os:
diff --git a/Project.toml b/Project.toml
index b39449d..8c8cdec 100644
--- a/Project.toml
+++ b/Project.toml
@@ -20,7 +20,7 @@ MacroTools = "0.5"
 Measurements = "2.2"
 MuladdMacro = "0.2"
 Reexport = "0.2, 1"
-julia = "1"
+julia = "1.3"
 
 [extras]
 ERFA = "17511681-8477-586a-8d98-4cfd5a1f2ec3"

From e6627813abfe93cd5518770b04277b3f9ebc89a1 Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Thu, 22 Jul 2021 13:22:04 +0200
Subject: [PATCH 09/16] Coverage

---
 src/AccurateArithmetic.jl | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/src/AccurateArithmetic.jl b/src/AccurateArithmetic.jl
index 0355332..6e8edc7 100644
--- a/src/AccurateArithmetic.jl
+++ b/src/AccurateArithmetic.jl
@@ -32,8 +32,7 @@ end
 
 function handle_infinity(fraction)
     second = ifelse(fraction < 0, typemin(Int64), typemax(Int64))
-    residual = zero(fraction)
-    return second, fraction, residual
+    return second, fraction, zero(fraction)
 end
 
 function apply_offset(s1::Int64, f1, e1, s2::Int64, f2, e2)

From d88b53a85659bb7613401ba4fc4a56185b7ac9b6 Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Thu, 22 Jul 2021 13:34:26 +0200
Subject: [PATCH 10/16] Cleanup

---
 src/AstroTime.jl    | 1 -
 src/Epochs/dates.jl | 2 --
 test/runtests.jl    | 2 +-
 3 files changed, 1 insertion(+), 4 deletions(-)

diff --git a/src/AstroTime.jl b/src/AstroTime.jl
index 7d7d2fe..a5d67af 100644
--- a/src/AstroTime.jl
+++ b/src/AstroTime.jl
@@ -236,7 +236,6 @@ macro timescale(scale::Symbol, parent=nothing, oneway=false)
             Dates.Second,
             Dates.Millisecond,
         )
-        Dates.default_format(::Type{$epoch_type}) = Dates.default_format(AstroDates.DateTime)
 
         $reg_expr
 
diff --git a/src/Epochs/dates.jl b/src/Epochs/dates.jl
index 43258d8..8f329e0 100644
--- a/src/Epochs/dates.jl
+++ b/src/Epochs/dates.jl
@@ -35,8 +35,6 @@ Epoch(str::AbstractString, format::Dates.DateFormat=Dates.default_format(Epoch))
 
 Epoch(str::AbstractString, format::AbstractString) = Epoch(str, Dates.DateFormat(format))
 
-Dates.default_format(::Type{Epoch{S}}) where {S} = Dates.default_format(AstroDates.DateTime)
-
 """
     Epoch{S}(str[, format]) where S
 
diff --git a/test/runtests.jl b/test/runtests.jl
index a618410..536137b 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -60,7 +60,7 @@ end
         @test TAIEpoch(scet, astronomical_unit) ≈ tai_exp
         @test SCETEpoch(tai_exp, astronomical_unit) ≈ scet
 
-        @test Dates.default_format(SCETEpoch) == AstroTime.ASTRO_ISO_FORMAT[]
+        @test Dates.default_format(SCETEpoch) == AstroTime.EPOCH_ISO_FORMAT[]
 
         @test string(Dummy) == "Dummy"
         @test typeof(Dummy) == DummyScale

From 3bb0ac3b7ad0608d69eddf5db7b894e270fe6b69 Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Thu, 22 Jul 2021 15:59:43 +0200
Subject: [PATCH 11/16] Fix ambiguities with `Dates`

---
 src/Periods.jl   | 3 +--
 test/epochs.jl   | 2 +-
 test/periods.jl  | 4 ++--
 test/runtests.jl | 3 +--
 4 files changed, 5 insertions(+), 7 deletions(-)

diff --git a/src/Periods.jl b/src/Periods.jl
index 9c6fe3e..7088490 100644
--- a/src/Periods.jl
+++ b/src/Periods.jl
@@ -2,9 +2,8 @@ module Periods
 
 using ..AccurateArithmetic: apply_offset, handle_infinity
 
-export TimeUnit, Second, Minute, Hour, Day, Year, Century,
+export AstroPeriod, TimeUnit,
     seconds, minutes, hours, days, years, centuries,
-    AstroPeriod,
     -, *, /, +, value, unit,
     SECONDS_PER_MINUTE,
     SECONDS_PER_HOUR,
diff --git a/test/epochs.jl b/test/epochs.jl
index 8d2bd91..e7e88b0 100644
--- a/test/epochs.jl
+++ b/test/epochs.jl
@@ -246,7 +246,7 @@ import ERFA
         ep_f64 = TAIEpoch(2000, 1, 1)
         ep_err = TAIEpoch(ep_f64.second, 1.0 ± 1.1)
         Δt = (30 ± 0.1) * seconds
-        @test typeof(Δt) == AstroPeriod{Second,Measurement{Float64}}
+        @test typeof(Δt) == AstroPeriod{AstroTime.Periods.Second,Measurement{Float64}}
         @test typeof(ep_f64) == Epoch{InternationalAtomicTime,Float64}
         @test typeof(ep_err) == Epoch{InternationalAtomicTime,Measurement{Float64}}
         @test typeof(ep_f64 + Δt) == Epoch{InternationalAtomicTime,Measurement{Float64}}
diff --git a/test/periods.jl b/test/periods.jl
index eb9a5ef..bd133a3 100644
--- a/test/periods.jl
+++ b/test/periods.jl
@@ -54,7 +54,7 @@
     @test centuries(y) == (1.0 / 100.0)centuries
     @test centuries(c) == 1.0centuries
 
-    @test zero(AstroPeriod{Year}) == 0.0years
+    @test zero(AstroPeriod{AstroTime.Periods.Year}) == 0.0years
     @test zero(1years) == 0years
     @test zero(1.0years) == 0.0years
 
@@ -83,7 +83,7 @@
     float_rng = 1.0seconds:3.0seconds
     @test step(float_rng) == 1.0seconds
     @test collect(float_rng) == [1.0seconds, 2.0seconds, 3.0seconds]
-    @test AstroPeriod{Second,Float64}(1.0seconds) == 1.0seconds
+    @test AstroPeriod{AstroTime.Periods.Second,Float64}(1.0seconds) == 1.0seconds
 
     @test Periods.name(seconds) == "seconds"
     @test Periods.name(minutes) == "minutes"
diff --git a/test/runtests.jl b/test/runtests.jl
index 536137b..431c040 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -1,6 +1,5 @@
 using AstroTime
-# TODO: Fix clashing `Period` exports
-import Dates
+using Dates
 using Test
 
 AstroTime.load_test_eop()

From b81ebb837a622d767285c23d07a8b891dedf48d9 Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Fri, 30 Jul 2021 17:53:30 +0200
Subject: [PATCH 12/16] Document and harmonize UTC API

---
 src/Epochs/ranges.jl |  2 +-
 src/Epochs/utc.jl    | 81 ++++++++++++++++++++++++++++++++++++++------
 src/Periods.jl       |  2 +-
 test/epochs.jl       | 11 +++---
 4 files changed, 79 insertions(+), 17 deletions(-)

diff --git a/src/Epochs/ranges.jl b/src/Epochs/ranges.jl
index 31bd60e..44a105a 100644
--- a/src/Epochs/ranges.jl
+++ b/src/Epochs/ranges.jl
@@ -1,4 +1,4 @@
-(::Base.Colon)(start::Epoch{S}, stop::Epoch{S}) where {S} = (:)(start, 1.0days, stop)
+(::Base.Colon)(start::Epoch{S}, stop::Epoch{S}) where {S} = (:)(start, 1.0seconds, stop)
 
 function (::Base.Colon)(start::Epoch{S}, step::AstroPeriod, stop::Epoch{S}) where S
     step = seconds(step)
diff --git a/src/Epochs/utc.jl b/src/Epochs/utc.jl
index 80d0381..2bdf3fe 100644
--- a/src/Epochs/utc.jl
+++ b/src/Epochs/utc.jl
@@ -2,30 +2,83 @@ const LEAP_J2000 = round.(Int, (LeapSeconds.LS_EPOCHS .- value(J2000_TO_MJD)) *
 const LEAP_TAI = LEAP_J2000 .+ round.(Int, LeapSeconds.LEAP_SECONDS) .- 1
 const LEAP_TAI_SET = Set(LEAP_TAI)
 
-function from_utc(str::AbstractString;
-        dateformat::Dates.DateFormat=Dates.default_format(AstroDates.DateTime),
+"""
+    from_utc(str::AbstractString, dateformat::Dates.DateFormat; scale=TAI)
+    from_utc(dt::Dates.DateTime; scale=TAI)
+    from_utc(year, month, day, hour=0, minute=0, second=0, fraction=0.0; scale=TAI)
+    from_utc(year, month, day, hour, minute, seconds; scale=TAI)
+
+Create an `Epoch` in `scale` based on a UTC timestamp, `Dates.DateTime` or date and
+time components.
+
+### Examples ###
+
+```jldoctest; setup = :(using AstroTime; import Dates)
+julia> from_utc(2016, 12, 31, 23, 59, 60, 0.0)
+2017-01-01T00:00:36.000 TAI
+
+julia> from_utc(2016, 12, 31, 23, 59, 60.0)
+2017-01-01T00:00:36.000 TAI
+
+julia> from_utc("2016-12-31T23:59:60.0")
+2017-01-01T00:00:36.000 TAI
+
+julia> from_utc("2016-12-31T23:59:60.0", scale=TDB)
+2017-01-01T00:01:08.183 TDB
+```
+"""
+from_utc
+
+"""
+    to_utc(ep)
+    to_utc(::Type{DateTime}, ep)
+    to_utc(::Type{Dates.DateTime}, ep)
+    to_utc(::Type{String}, ep, dateformat=Dates.default_format(DateTime))
+
+Create a UTC timestamp or `Dates.DateTime` from an `Epoch` `ep`.
+
+### Examples ###
+
+```jldoctest; setup = :(using AstroTime; import Dates)
+julia> tai = from_utc(Dates.DateTime(2018, 2, 6, 20, 45, 0, 0))
+2018-02-06T20:45:37.000 TAI
+
+julia> to_utc(tai)
+"2018-02-06T20:45:00.000"
+
+julia> to_utc(String, tai, Dates.dateformat"yyyy-mm-dd")
+"2018-02-06"
+
+julia> to_utc(Dates.DateTime, tai)
+2018-02-06T20:45:00
+```
+"""
+to_utc
+
+function from_utc(str::AbstractString,
+        dateformat::Dates.DateFormat=Dates.default_format(AstroDates.DateTime);
         scale::TimeScale=TAI)
     dt = AstroDates.DateTime(str, dateformat)
-    return from_utc(dt, scale)
+    return from_utc(dt; scale)
 end
 
 function from_utc(year::Integer, month::Integer, day::Integer,
     hour::Integer=0, minute::Integer=0, second::Integer=0, fraction=0.0;
     scale::TimeScale=TAI)
     dt = DateTime(year, month, day, hour, minute, second, fraction)
-    return from_utc(dt, scale)
+    return from_utc(dt; scale)
 end
 
 function from_utc(year::Integer, month::Integer, day::Integer,
-    hour::Integer, minute::Integer, second;
+    hour::Integer, minute::Integer, seconds;
     scale::TimeScale=TAI)
-    dt = DateTime(year, month, day, hour, minute, second)
-    return from_utc(dt, scale)
+    dt = DateTime(year, month, day, hour, minute, seconds)
+    return from_utc(dt; scale)
 end
 
-from_utc(dt::Dates.DateTime, scale::S=TAI) where {S} = from_utc(DateTime(dt), scale)
+from_utc(dt::Dates.DateTime; scale::S=TAI) where {S} = from_utc(DateTime(dt); scale)
 
-function from_utc(dt::DateTime, scale::S=TAI) where S
+function from_utc(dt::DateTime; scale::S=TAI) where S
     ep = TAIEpoch(dt)
     idx = searchsortedlast(LEAP_J2000, ep.second)
     if idx == 0
@@ -58,6 +111,11 @@ function to_utc(::Type{DateTime}, ep)
     return DateTime(d, t)
 end
 
+function to_utc(::Type{Dates.DateTime}, ep)
+    dt = to_utc(DateTime, ep)
+    return Dates.DateTime(dt)
+end
+
 function to_utc(::Type{String}, ep, dateformat=Dates.default_format(DateTime))
     dt = to_utc(DateTime, ep)
     return Dates.format(dt, dateformat)
@@ -82,9 +140,10 @@ julia> now(TDBEpoch)
 ```
 """
 function Dates.now(::Type{Epoch{S}}) where S
-    return from_utc(Dates.now(Dates.UTC), S())
+    return from_utc(Dates.now(Dates.UTC); scale=S())
 end
 
 function Dates.now(::Type{Epoch})
-    return from_utc(Dates.now(Dates.UTC), TAI)
+    return from_utc(Dates.now(Dates.UTC); scale=TAI)
 end
+
diff --git a/src/Periods.jl b/src/Periods.jl
index 7088490..caf939e 100644
--- a/src/Periods.jl
+++ b/src/Periods.jl
@@ -111,7 +111,7 @@ Return the unit of the period `p`.
 
 ```jldoctest; setup = :(using AstroTime)
 julia> unit(3.0seconds)
-Second()
+AstroTime.Periods.Second()
 ```
 """
 unit(p::AstroPeriod) = p.unit
diff --git a/test/epochs.jl b/test/epochs.jl
index e7e88b0..351d083 100644
--- a/test/epochs.jl
+++ b/test/epochs.jl
@@ -41,6 +41,9 @@ import ERFA
         @test to_utc(during_act) == during
         @test to_utc(after_act) == after
 
+        dt = Dates.DateTime(2016, 12, 31, 23, 59, 59, 300)
+        @test dt == to_utc(Dates.DateTime, from_utc(dt))
+
         sixties = AstroTime.DateTime(1961, 3, 5, 23, 4, 12.0)
         sixties_exp = (second=-1225198547, fraction=0.5057117799999999)
         sixties_act = from_utc(sixties)
@@ -233,11 +236,11 @@ import ERFA
         @test Dates.DateTime(ep) == Dates.DateTime(y, m, d, hr, mn, 59, 371)
     end
     @testset "Ranges" begin
-        rng = TAIEpoch(2018, 1, 1):TAIEpoch(2018, 2, 1)
-        @test step(rng) == 86400.0seconds
-        @test length(rng) == 32
+        rng = TAIEpoch(2018, 1, 1):TAIEpoch(2018, 1, 2)
+        @test step(rng) == 1seconds
+        @test length(rng) == 86401
         @test first(rng) == TAIEpoch(2018, 1, 1)
-        @test last(rng) == TAIEpoch(2018, 2, 1)
+        @test last(rng) == TAIEpoch(2018, 1, 2)
         rng = TAIEpoch(2018, 1, 1):13seconds:TAIEpoch(2018, 1, 1, 0, 1)
         @test step(rng) == 13seconds
         @test last(rng) == TAIEpoch(2018, 1, 1, 0, 0, 52.0)

From 1c66906942cd5ceecad45a41958ccf1c1cb3170e Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Fri, 30 Jul 2021 17:54:31 +0200
Subject: [PATCH 13/16] Update documentation

---
 README.md            |   6 +-
 docs/Project.toml    |   2 +-
 docs/src/index.md    |   6 +-
 docs/src/tutorial.md | 205 ++++++++++++++++++++++++++++++-------------
 4 files changed, 152 insertions(+), 67 deletions(-)

diff --git a/README.md b/README.md
index c4fdd1a..e163877 100644
--- a/README.md
+++ b/README.md
@@ -24,11 +24,11 @@ julia> import Pkg; Pkg.add("AstroTime")
 # Create an Epoch based on the TT (Terrestial Time) scale
 tt = TTEpoch("2018-01-01T12:00:00")
 
-# Transform to UTC (Universal Time Coordinated)
-utc = UTCEpoch(tt)
+# Transform to TAI (International Atomic Time)
+tai = TAIEpoch(tt)
 
 # Transform to TDB (Barycentric Dynamical Time)
-utc = TDBEpoch(utc)
+tdb = TDBEpoch(tai)
 
 # Shift an Epoch by one day
 another_day = tt + 1days
diff --git a/docs/Project.toml b/docs/Project.toml
index 0d8fbc7..1a6d309 100644
--- a/docs/Project.toml
+++ b/docs/Project.toml
@@ -2,4 +2,4 @@
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 
 [compat]
-Documenter = "~0.26"
+Documenter = "~0.27"
diff --git a/docs/src/index.md b/docs/src/index.md
index d5ec7d8..58559b8 100644
--- a/docs/src/index.md
+++ b/docs/src/index.md
@@ -19,11 +19,11 @@ julia> import Pkg; Pkg.add("AstroTime")
 # Create an Epoch based on the TT (Terrestial Time) scale
 tt = TTEpoch("2018-01-01T12:00:00")
 
-# Transform to UTC (Universal Time Coordinated)
-utc = UTCEpoch(tt)
+# Transform to TAI (International Atomic Time)
+tai = TAIEpoch(tt)
 
 # Transform to TDB (Barycentric Dynamical Time)
-utc = TDBEpoch(utc)
+tdb = TDBEpoch(tai)
 
 # Shift an Epoch by one day
 another_day = tt + 1days
diff --git a/docs/src/tutorial.md b/docs/src/tutorial.md
index 9a5b7d8..bf35f14 100644
--- a/docs/src/tutorial.md
+++ b/docs/src/tutorial.md
@@ -2,12 +2,17 @@
 
 This tutorial will walk you through the features and functionality of AstroTime.jl.
 Everything in this package revolves around the `Epoch` data type.
-`Epochs` are a high-precision, time-scale aware version of the [`DateTime`](https://docs.julialang.org/en/v1.0/stdlib/Dates) type from Julia's standard library.
-This means that while `DateTime` timestamps are always assumed to be based on Universal Time (UT), `Epochs` can be created in several pre-defined time scales or custom user-defined time scales.
+`Epochs` are a high-precision, time-scale aware version of the
+[`DateTime`](https://docs.julialang.org/en/v1.0/stdlib/Dates) type from Julia's standard
+library.
+This means that while `DateTime` timestamps are always assumed to be based on Universal
+Time (UT), `Epochs` can be created in several pre-defined time scales or custom user-defined
+time scales.
 
 ## Creating Epochs
 
-You construct `Epoch` instances similar to `DateTime` instance, for example by using date and time components.
+You construct `Epoch` instances similar to `DateTime` instances, for example by using date
+and time components.
 The main difference is that you need to supply the time scale to be used.
 Out of the box, the following time scales are defined:
 
@@ -18,6 +23,11 @@ Out of the box, the following time scales are defined:
 - [`TCB`](@ref): [Barycentric Coordinate Time](https://en.wikipedia.org/wiki/Barycentric_Coordinate_Time)
 - [`TDB`](@ref): [Barycentric Dynamical Time](https://en.wikipedia.org/wiki/Barycentric_Dynamical_Time)
 
+Conspicuously missing from this list is [Coordinated Universal Time (UTC)](https://en.wikipedia.org/wiki/Coordinated_Universal_Time).
+While AstroTime.jl does support UTC, it requires special treatment due to the discontinuities
+in the time scale from the introduction of leap seconds.
+See [UTC and Leap Seconds](@ref) for more details.
+
 [^1]:
     Transformations to and from UT1 depend on the measured quantity ΔUT1 which is
     published in [IERS](https://www.iers.org) tables on a weekly basis. AstroTime.jl can
@@ -27,13 +37,16 @@ Out of the box, the following time scales are defined:
 ```julia
 using AstroTime
 
-ep = Epoch{CoordinatedUniversalTime}(2018, 2, 6, 20, 45, 0.0)
+ep = Epoch{InternationalAtomicTime}(2018, 2, 6, 20, 45, 0.0)
 
 # The following shorthand syntax also works
-ep = UTCEpoch(2018, 2, 6, 20, 45, 0.0)
+ep = TAIEpoch(2018, 2, 6, 20, 45, 0.0)
 
 # Or in another time scale
-ep = TAIEpoch(2018, 2, 6, 20, 45, 0.0)
+ep = TTEpoch(2018, 2, 6, 20, 45, 0.0)
+
+# Or use UTC with leap second handling
+ep = from_utc(2018, 2, 6, 20, 45, 0.0)
 ```
 
 You can also parse an `Epoch` from a string.
@@ -41,44 +54,45 @@ AstroTime.jl uses the [`DateFormat`](https://docs.julialang.org/en/stable/stdlib
 For example:
 
 ```julia
-ep = UTCEpoch("2018-02-06T20:45:00.000", "yyyy-mm-ddTHH:MM:SS.sss")
+ep = TAIEpoch("2018-02-06T20:45:00.000", "yyyy-mm-ddTHH:MM:SS.fff")
 
-# The format string above `yyyy-mm-ddTHH:MM:SS.sss` is also the default format.
+# The format string above `yyyy-mm-ddTHH:MM:SS.fff` is also the default format.
 # Thus, this also works...
-ep = UTCEpoch("2018-02-06T20:45:00.000")
+ep = TAIEpoch("2018-02-06T20:45:00.000")
 
 import Dates
 
 # You can also reuse the format string
 df = Dates.dateformat"dd.mm.yyyy HH:MM"
 
-utc = UTCEpoch("06.02.2018 20:45", df)
+utc = from_utc("06.02.2018 20:45", df)
 tai = TAIEpoch("06.02.2018 20:45", df)
 ```
 
-There are two additional character codes supported.
+There are three additional character codes supported.
 
+- `f`: This character code is parsed as the fraction of the current second and supports an arbitrary number of decimal places.
 - `t`: This character code is parsed as the time scale.
 - `D`: This character code is parsed as the day number within a year.
 
 ```julia
 # The time scale can be omitted from the constructor because it is already
 # defined in the input string
-julia> Epoch("2018-02-06T20:45:00.000 UTC", "yyyy-mm-ddTHH:MM:SS.sss ttt")
-2018-02-06T20:45:00.000 UTC
+julia> Epoch("2018-02-06T20:45:00.000 TAI", "yyyy-mm-ddTHH:MM:SS.fff ttt")
+2018-02-06T20:45:00.000 TAI
 
 # February 6 is the 37th day of the year
-julia> UTCEpoch("2018-037T20:45:00.000", "yyyy-DDDTHH:MM:SS.sss")
-2018-02-06T20:45:00.000 UTC
+julia> TAIEpoch("2018-037T20:45:00.000", "yyyy-DDDTHH:MM:SS.fff")
+2018-02-06T20:45:00.000 TAI
 ```
 
 When printing `Epochs`, you can format the output in the same way.
 
 ```julia
-julia> ep = UTCEpoch(2018, 2, 6, 20, 45, 0.0)
-2018-02-06T20:45:00.000 UTC
+julia> ep = TAIEpoch(2018, 2, 6, 20, 45, 0.0)
+2018-02-06T20:45:00.000 TAI
 julia> AstroTime.format(ep, "dd.mm.yyyy HH:MM ttt")
-06.02.2018 20:45 UTC
+06.02.2018 20:45 TAI
 ```
 
 ## Working with Epochs and Periods
@@ -108,31 +122,31 @@ The following time units are available:
 To shift an `Epoch` forward in time add an `AstroPeriod` to it.
 
 ```julia
-julia> ep = UTCEpoch(2000, 1, 1)
-2000-01-01T00:00:00.000 UTC
+julia> ep = TAIEpoch(2000, 1, 1)
+2000-01-01T00:00:00.000 TAI
 
 julia> ep + 1days
-2000-01-02T00:00:00.000 UTC
+2000-01-02T00:00:00.000 TAI
 ```
 
 Or subtract it to shift the `Epoch` backwards.
 
 ```julia
-julia> ep = UTCEpoch(2000, 1, 1)
-2000-01-01T00:00:00.000 UTC
+julia> ep = TAIEpoch(2000, 1, 1)
+2000-01-01T00:00:00.000 TAI
 
 julia> ep - 1days
-1999-12-31T00:00:00.000 UTC
+1999-12-31T00:00:00.000 TAI
 ```
 
 If you subtract two epochs you will receive the time between them as an `AstroPeriod`.
 
 ```julia
-julia> ep1 = UTCEpoch(2000, 1, 1)
-2000-01-01T00:00:00.000 UTC
+julia> ep1 = TAIEpoch(2000, 1, 1)
+2000-01-01T00:00:00.000 TAI
 
-julia> ep2 = UTCEpoch(2000, 1, 2)
-2000-01-02T00:00:00.000 UTC
+julia> ep2 = TAIEpoch(2000, 1, 2)
+2000-01-02T00:00:00.000 TAI
 
 julia> ep2 - ep1
 86400.0 seconds
@@ -159,23 +173,94 @@ julia> value(days(dt))
 1.0
 ```
 
+## Ranges
+
+You can also construct ranges of `Epoch`s. The default step size one second.
+
+```julia
+julia> TAIEpoch(2021, 7, 30, 17, 34, 30.0):TAIEpoch(2021, 7, 30, 17, 34, 31.0)
+2021-07-30T17:34:30.000 TAI:1.0 seconds:2021-07-30T17:34:31.000 TAI
+```
+
+Or you can adjust the step size with any of the units supported.
+
+```julia
+julia> collect(TAIEpoch(2000, 1, 1):1days:TAIEpoch(2000, 1, 5))
+5-element Vector{TAIEpoch{Float64}}:
+ 2000-01-01T00:00:00.000 TAI
+ 2000-01-02T00:00:00.000 TAI
+ 2000-01-03T00:00:00.000 TAI
+ 2000-01-04T00:00:00.000 TAI
+ 2000-01-05T00:00:00.000 TAI
+```
+
 ## Converting Between Time Scales
 
 You convert an `Epoch` to another time scale by constructing a new `Epoch` with the
 target time scale from it.
 
 ```julia
-julia> utc = UTCEpoch(2018, 2, 6, 20, 45, 0.0)
-2018-02-06T20:45:00.000 UTC
+julia> tai = TAIEpoch(2018, 2, 6, 20, 45, 0.0)
+2018-02-06T20:45:00.000 TAI
+
+julia> tt = TTEpoch(tai) # Convert to TT
+2018-02-06T20:45:32.184 TT
+```
+
+### UTC and Leap Seconds
+
+UTC is the primary civil time standard and aims to provide a time scale based on TAI and
+uniform SI seconds that is at the same time aligned with UT1 which is based on solar time
+and governed by the rotation of the Earth. The problem is that Earth's rotation speed is
+much more irregular compared to atomic clocks which define the SI second.
+Over the past decades, Earth's rotation has continuously slowed and thus TAI has been
+running ahead of UT1.
+
+Leap seconds are inserted into the UTC time scale to keep it within 0.9 seconds of UT1.
+This introduces ambiguities in AstroTime.jl's data model (see [#50](@ref)).
+As a consequence, `UTCEpoch`s are not supported.
+Nevertheless, UTC is supported as an I/O format for timestamps through the
+[`from_utc`](@ref) and [`to_utc`](@ref) functions.
+
+The last leap second was introduced at the end of December 31, 2016. You can create a
+`TAIEpoch` (or other `Epoch`s) from a UTC date with proper leap second handling:
+
+```julia
+julia> from_utc(2016, 12, 31, 23, 59, 60.0)
+2017-01-01T00:00:36.000 TAI
+
+julia> from_utc("2016-12-31T23:59:60.0")
+2017-01-01T00:00:36.000 TAI
+
+julia> from_utc("2016-12-31T23:59:60.0", scale=TDB)
+2017-01-01T00:01:08.183 TDB
+```
 
-julia> tai = TAIEpoch(utc) # Convert to TAI
+You can also use `Dates.DateTime` but note that you cannot represent a leap second
+date with it.
+
+```julia
+julia> tai = from_utc(Dates.DateTime(2018, 2, 6, 20, 45, 0, 0))
 2018-02-06T20:45:37.000 TAI
 ```
 
+And go back to UTC:
+
+```julia
+julia> to_utc(tai)
+"2018-02-06T20:45:00.000"
+
+julia> to_utc(String, tai, Dates.dateformat"yyyy-mm-dd")
+"2018-02-06"
+
+julia> to_utc(Dates.DateTime, tai)
+2018-02-06T20:45:00
+```
+
 ### High-Precision Conversions and Custom Offsets
 
 Some time scale transformations depend on measured quantities which cannot be accurately
-predicted (e.g. UTC to UT1) or there are different algortihms which offer variable levels
+predicted (e.g. UT1) or there are different algorithms which offer variable levels
 of accuracy.
 For the former, AstroTime.jl can download the required data automatically from the internet.
 You need to run `AstroTime.update()` periodically (weekly) to keep this data up-to-date.
@@ -274,24 +359,24 @@ julia> import Dates; Dates.DateTime(ep)
 AstroTime.jl enables you to create your own first-class time scales via the [`@timescale`](@ref) macro.
 The macro will define the necessary structs and register the new time scale.
 
-Let's start with a simple example and assume that you want to define `GMT` as an alias for `UTC`.
+Let's start with a simple example and assume that you want to define `EphemerisTime` as an alias for `TDB`.
 You need to provide the name of the time scale and optionally a "parent" time scale to which it is linked.
 
 ```julia
-@timescale GMT UTC
+@timescale EphemerisTime TDB
 ```
 
 At this point, you can already use the new time scale to create epochs.
 
 ```julia
-julia> GMT
-GMT
+julia> EphemerisTime
+EphemerisTime
 
-julia> typeof(GMT)
-GMTScale
+julia> typeof(EphemerisTime)
+EphemerisTimeScale
 
-julia> gmt = GMTEpoch(2000, 1, 1)
-2000-01-01T00:00:00.000 GMT
+julia> et = EphemerisTimeEpoch(2000, 1, 1)
+2000-01-01T00:00:00.000 EphemerisTime
 ```
 
 Conversion to other `Epoch` types will not yet work for the newly created time
@@ -300,29 +385,28 @@ If you are unsure which methods are needed, you can try to transform the epoch
 and the resulting error message will provide a hint.
 
 ```julia
-julia> UTCEpoch(gmt)
-ERROR: No conversion 'GMT->UTC' available. If one of these is a custom time scale,
-you may need to define `AstroTime.Epochs.getoffset(::GMTScale, ::CoordinatedUniversalTime,
-second, fraction, args...)`.
+julia> TDBEpoch(et)
+ERROR: No conversion 'EphemerisTime->TDB' available. If one of these is a custom time scale,
+you may need to define `AstroTime.Epochs.getoffset(::EphemerisTimeScale, ::BarycentricDynamicalTime, second, fraction, args...)`.
 ```
 
-To enable transformations between `GMT` and `UTC` in both directions you need
+To enable transformations between `EphemerisTime` and `TDB` in both directions you need
 to define the following methods.
-Since `GMT` is the same offset as `UTC`, these can just return zero.
+Since `EphemerisTime` and `TDB` are identical, the offset between them is zero.
 
 ```julia
-AstroTime.Epochs.getoffset(::GMTType, ::CoordinatedUniversalTime, second, fraction) = 0.0
-AstroTime.Epochs.getoffset(::CoordinatedUniversalTime, ::GMTType, second, fraction) = 0.0
+AstroTime.Epochs.getoffset(::EphemerisTimeType, ::CoordinatedUniversalTime, second, fraction) = 0.0
+AstroTime.Epochs.getoffset(::CoordinatedUniversalTime, ::EphemerisTimeType, second, fraction) = 0.0
 ```
 
-You can now use `GMTEpoch` like any other epoch type, e.g.
+You can now use `EphemerisTimeEpoch` like any other epoch type, e.g.
 
 ```julia
-julia> ep = UTCEpoch(2000, 1, 1)
-2000-01-01T00:00:00.000 UTC
+julia> ep = TDBEpoch(2000, 1, 1)
+2000-01-01T00:00:00.000 TDB
 
-julia> GMTEpoch(ep)
-2000-01-01T00:00:00.000 GMT
+julia> EphemerisTimeEpoch(ep)
+2000-01-01T00:00:00.000 EphemerisTime
 ```
 
 For a more complex example, let's reimplement the Geocentric Coordinate
@@ -365,13 +449,13 @@ which is the distance of the spacecraft from Earth.
 ```julia
 const speed_of_light = 299792458.0 # m/s
 
-@timescale SCET UTC
+@timescale SCET TAI
 
-function AstroTime.Epochs.getoffset(::SCETType, ::CoordinatedUniversalTime,
+function AstroTime.Epochs.getoffset(::SCETScale, ::InternationalAtomicTime,
                                     second, fraction, distance)
     return distance / speed_of_light
 end
-function AstroTime.Epochs.getoffset(::CoordinatedUniversalTime, ::SCETType,
+function AstroTime.Epochs.getoffset(::InternationalAtomicTime, ::SCETScale,
                                     second, fraction, distance)
     return -distance / speed_of_light
 end
@@ -385,11 +469,11 @@ For example, for a spacecraft that is one astronomical unit away from Earth:
 julia> astronomical_unit = 149597870700.0 # m
 149597870700.0
 
-julia> ep = UTCEpoch(2000, 1, 1)
-2000-01-01T00:00:00.000 UTC
+julia> ep = TAIEpoch(2000, 1, 1)
+2000-01-01T00:00:00.000 TAI
 
 julia> SCETEpoch(ep, astronomical_unit)
-2000-01-01T00:08:19.005 SCET
+1999-12-31T23:51:40.995 SCET
 ```
 
 !!! note
@@ -401,7 +485,8 @@ default graph of time scales by defining a time scale without a parent.
 ```julia
 julia> @timescale Disjoint
 
-julia> typeof(Disjoint) = DisjointScale
+julia> typeof(Disjoint)
+DisjointScale
 ```
 
 By defining additional time scales connected to this scale and the appropriate

From 7638187f9bd7d056442bb5ff23761fdea95bf720 Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Fri, 30 Jul 2021 17:59:07 +0200
Subject: [PATCH 14/16] Use old kw-arg syntax

---
 src/Epochs/utc.jl | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/src/Epochs/utc.jl b/src/Epochs/utc.jl
index 2bdf3fe..4a29741 100644
--- a/src/Epochs/utc.jl
+++ b/src/Epochs/utc.jl
@@ -59,21 +59,21 @@ function from_utc(str::AbstractString,
         dateformat::Dates.DateFormat=Dates.default_format(AstroDates.DateTime);
         scale::TimeScale=TAI)
     dt = AstroDates.DateTime(str, dateformat)
-    return from_utc(dt; scale)
+    return from_utc(dt; scale=scale)
 end
 
 function from_utc(year::Integer, month::Integer, day::Integer,
     hour::Integer=0, minute::Integer=0, second::Integer=0, fraction=0.0;
     scale::TimeScale=TAI)
     dt = DateTime(year, month, day, hour, minute, second, fraction)
-    return from_utc(dt; scale)
+    return from_utc(dt; scale=scale)
 end
 
 function from_utc(year::Integer, month::Integer, day::Integer,
     hour::Integer, minute::Integer, seconds;
     scale::TimeScale=TAI)
     dt = DateTime(year, month, day, hour, minute, seconds)
-    return from_utc(dt; scale)
+    return from_utc(dt; scale=scale)
 end
 
 from_utc(dt::Dates.DateTime; scale::S=TAI) where {S} = from_utc(DateTime(dt); scale)

From ddc21e7c5d3b761975a690d5dd157cb33bfed5aa Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Fri, 30 Jul 2021 18:02:00 +0200
Subject: [PATCH 15/16] Try again

---
 src/Epochs/utc.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/Epochs/utc.jl b/src/Epochs/utc.jl
index 4a29741..a27e8de 100644
--- a/src/Epochs/utc.jl
+++ b/src/Epochs/utc.jl
@@ -76,7 +76,7 @@ function from_utc(year::Integer, month::Integer, day::Integer,
     return from_utc(dt; scale=scale)
 end
 
-from_utc(dt::Dates.DateTime; scale::S=TAI) where {S} = from_utc(DateTime(dt); scale)
+from_utc(dt::Dates.DateTime; scale::S=TAI) where {S} = from_utc(DateTime(dt); scale=scale)
 
 function from_utc(dt::DateTime; scale::S=TAI) where S
     ep = TAIEpoch(dt)

From 4fdb65a014603ace5db69ca0e2438bc4c382e7cd Mon Sep 17 00:00:00 2001
From: Helge Eichhorn <git@helgeeichhorn.de>
Date: Fri, 30 Jul 2021 18:09:42 +0200
Subject: [PATCH 16/16] Update compat and bump version

---
 Project.toml | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/Project.toml b/Project.toml
index 8c8cdec..de0b3df 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,6 +1,6 @@
 name = "AstroTime"
 uuid = "c61b5328-d09d-5e37-a9a8-0eb41c39009c"
-version = "0.6.2"
+version = "0.7.0"
 
 [deps]
 Dates = "ade2ca70-3891-5945-98fb-dc099432e06a"
@@ -12,12 +12,12 @@ MuladdMacro = "46d2c3a1-f734-5fdb-9937-b9b9aeba4221"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 
 [compat]
-ERFA = "0.5, 0.6"
+ERFA = "0.5, 0.6, 1.0"
 EarthOrientation = "0.7"
 ItemGraphs = "0.4"
 LeapSeconds = "1.1"
 MacroTools = "0.5"
-Measurements = "2.2"
+Measurements = "2"
 MuladdMacro = "0.2"
 Reexport = "0.2, 1"
 julia = "1.3"