From 2ddbaea187c74536782778b63f8e98381bb7275b Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Thu, 14 Nov 2019 11:47:59 -0300
Subject: [PATCH 01/14] REL-3179: WIP

---
 .../gemini/ghost/ProperMotionCalculator.scala | 227 ++++++++++++++++++
 1 file changed, 227 insertions(+)
 create mode 100644 bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala

diff --git a/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala b/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala
new file mode 100644
index 0000000000..b9c4be7997
--- /dev/null
+++ b/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala
@@ -0,0 +1,227 @@
+// Copyright (c) 2016-2019 Association of Universities for Research in Astronomy, Inc. (AURA)
+// For license information see LICENSE or https://opensource.org/licenses/BSD-3-Clause
+
+package edu.gemini.spModel.gemini.ghost
+
+import java.time.Instant
+
+import scala.math.{atan2, cos, hypot, sin}
+import scalaz._
+import Scalaz._
+import edu.gemini.spModel.core.{Angle, Coordinates, Epoch, Offset}
+
+/**
+ * Time-parameterized coordinates, based on an observed position at some point in time (called
+ * the `epoch`) and measured velocities in distance (`radialVelocity`; i.e., doppler shift) and
+ * position (`properVelocity`) per year. Given this information we can compute the position at any
+ * instant in time. The references below are ''extremely'' helpful, so do check them out if you're
+  * trying to understand the implementation.
+ * @see The pretty good [[https://en.wikipedia.org/wiki/Proper_motion wikipedia]] article
+ * @see Astronomical Almanac 1984 [[https://babel.hathitrust.org/cgi/pt?id=uc1.b3754036;view=1up;seq=141 p.B39]]
+ * @see Astronomy and Astrophysics 134 (1984) [[http://articles.adsabs.harvard.edu/cgi-bin/nph-iarticle_query?bibcode=1984A%26A...134....1L&db_key=AST&page_ind=0&data_type=GIF&type=SCREEN_VIEW&classic=YES p.1-6]]
+ * @param baseCoordinates observed coordinates at `epoch`
+ * @param epoch           time of the base observation; typically `Epoch.J2000`
+ * @param properVelocity  proper velocity '''per year''' in [[RightAscension]] and [[Declination]], if any
+ * @param radialVelocity  radial velocity (km/y, positive if receding), if any
+ * @param parallax        parallax, if any
+ */
+final case class ProperMotion(
+  baseCoordinates: Coordinates,
+  epoch:           Epoch,
+  properVelocity:  Option[Offset],
+  radialVelocity:  Option[RadialVelocity],
+  parallax:        Option[Angle]
+) {
+
+  def at(i: Instant): ProperMotion =
+    plusYears(epoch.untilInstant(i))
+
+  /** Coordinates `elapsedYears` fractional epoch-years after `epoch`. */
+  def plusYears(elapsedYears: Double): ProperMotion =
+    ProperMotion(
+      ProperMotion.properMotion(
+        baseCoordinates,
+        epoch,
+        properVelocity.orEmpty,
+        radialVelocity.getOrElse(RadialVelocity.Zero).toDoubleKilometersPerSecond,
+        parallax.orEmpty,
+        elapsedYears
+      ),
+      epoch.plusYears(elapsedYears),
+      properVelocity,
+      radialVelocity,
+      parallax
+    )
+
+}
+
+
+object ProperMotion extends ProperMotionOptics {
+  import PhysicalConstants.{ AstronomicalUnit, TwoPi }
+
+  def const(cs: Coordinates): ProperMotion =
+    ProperMotion(cs, Epoch.J2000, None, None, None)
+
+  /**
+   * Proper motion correction in model units.
+   * @param baseCoordinates base coordinates
+   * @param epoch           the epoch
+   * @param properVelocity  proper velocity per epoch year
+   * @param radialVelocity  radial velocity (km/sec, positive if receding)
+   * @param parallax        parallax
+   * @param elapsedYears    elapsed time in epoch years
+   * @return Coordinates corrected for proper motion
+   */
+  def properMotion(
+    baseCoordinates: Coordinates,
+    epoch:           Epoch,
+    properVelocity:  Offset,
+    radialVelocity:  Double,
+    parallax:        Angle,
+    elapsedYears:    Double
+  ): Coordinates = {
+    val (ra, dec) = properMotionʹ(
+      baseCoordinates.toRadians,
+      epoch.scheme.lengthOfYear,
+      properVelocity.toRadians,
+      radialVelocity,
+      parallax.toMicroarcseconds.toDouble / 1000000.0,
+      elapsedYears
+    )
+    Coordinates.unsafeFromRadians(ra, dec)
+  }
+
+  // Some constants we need
+  private val secsPerDay  = 86400.0
+  private val auPerKm     = 1000.0 / AstronomicalUnit.toDouble
+  private val radsPerAsec = Angle.arcseconds.reverseGet(1).toDoubleRadians
+
+  // We need to do things with little vectors of doubles
+  private type Vec2 = (Double, Double)
+  private type Vec3 = (Double, Double, Double)
+
+  // |+| gives us addition for VecN, but we also need scalar multiplication
+  private implicit class Vec3Ops(a: Vec3) {
+    def *(d: Double): Vec3 =
+      (a._1 * d, a._2 * d, a._3 * d)
+  }
+
+  /**
+   * Proper motion correction in base units.
+   * @param baseCoordinates base (ra, dec) in radians, [0 … 2π) and (-π/2 … π/2)
+   * @param daysPerYear     length of epoch year in fractonal days
+   * @param properVelocity  proper velocity in (ra, dec) in radians per epoch year
+   * @param radialVelocity  radial velocity (km/sec, positive means away from earth)
+   * @param parallax        parallax in arcseconds (!)
+   * @param elapsedYears    elapsed time in epoch years
+   * @return (ra, dec) in radians, corrected for proper motion
+   */
+  // scalastyle:off method.length
+  private def properMotionʹ(
+    baseCoordinates: Vec2,
+    daysPerYear:     Double,
+    properVelocity:  Vec2,
+    parallax:        Double,
+    radialVelocity:  Double,
+    elapsedYears:    Double
+  ): Vec2 = {
+
+    // Break out our components
+    val (ra,   dec) = baseCoordinates
+    val (dRa, dDec) = properVelocity
+
+    // Convert to cartesian
+    val pos: Vec3 = {
+      val cd = cos(dec)
+      (cos(ra) * cd, sin(ra) * cd, sin(dec))
+    }
+
+    // Change per year due to radial velocity and parallax. The units work out to asec/y.
+    val dPos1: Vec3 =
+      pos            *
+      daysPerYear    *
+      secsPerDay     *
+      radsPerAsec    *
+      auPerKm        *
+      radialVelocity *
+      parallax
+
+    // Change per year due to proper velocity
+    val dPos2 = (
+      -dRa * pos._2 - dDec * cos(ra) * sin(dec),
+       dRa * pos._1 - dDec * sin(ra) * sin(dec),
+                      dDec *           cos(dec)
+    )
+
+    // Our new position (still in polar coordinates). `|+|` here is scalar addition provided by
+    // cats … unlike scalaz it does give you Semigroup[Double] even though it's not strictly lawful.
+    val pʹ = pos |+| ((dPos1 |+| dPos2) * elapsedYears)
+
+    // Back to spherical
+    val (x, y, z) = pʹ
+    val r    = hypot(x, y)
+    val raʹ  = if (r === 0.0) 0.0 else atan2(y, x)
+    val decʹ = if (z === 0.0) 0.0 else atan2(z, r)
+    val raʹʹ = {
+      // Normalize to [0 .. 2π)
+      val rem = raʹ % TwoPi
+      if (rem < 0.0) rem + TwoPi else rem
+    }
+    (raʹʹ, decʹ)
+
+  }
+  // scalastyle:on method.length
+
+  implicit val OrderProperMotion: Order[ProperMotion] = {
+
+    implicit val MonoidOrder: Monoid[Order[ProperMotion]] =
+      Order.whenEqualMonoid[ProperMotion]
+
+    implicit val AngleOrder: Order[Angle] =
+      Angle.SignedAngleOrder
+
+    def order[A: Order](f: ProperMotion => A): Order[ProperMotion] =
+      Order.by(f)
+
+    // This could be done with:
+    //
+    //   Order.by(p => (p.baseCoordinates, p.epoch, ...))
+    //
+    // but that would always perform comparisons for all the fields (and all
+    // their contained fields down to the leaves of the tree) all of the time.
+    // The Monoid approach on the other hand will stop at the first difference.
+    // This is premature optimization perhaps but it seems like it might make a
+    // difference when sorting a long list of targets.
+
+    order(_.baseCoordinates)  |+|
+      order(_.epoch)          |+|
+      order(_.properVelocity) |+|
+      order(_.radialVelocity) |+|
+      order(_.parallax)
+
+  }
+}
+
+trait ProperMotionOptics {
+
+  /** @group Optics */
+  val baseCoordinates: Lens[ProperMotion, Coordinates] =
+    Lens[ProperMotion](_.baseCoordinates)
+
+  /** @group Optics */
+  val epoch: Lens[ProperMotion, Epoch] =
+    Lens[ProperMotion](_.epoch)
+
+  /** @group Optics */
+  val properVelocity: Lens[ProperMotion, Option[Offset]] =
+    Lens[ProperMotion](_.properVelocity)
+
+  /** @group Optics */
+  val radialVelocity: Lens[ProperMotion, Option[RadialVelocity]] =
+    Lens[ProperMotion](_.radialVelocity)
+
+  /** @group Optics */
+  val parallax: Lens[ProperMotion, Option[Angle]] =
+    Lens[ProperMotion](_.parallax)
+
+}

From 5f9252df1a8254315e6b0fbae8d632180b2499ed Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Mon, 18 Nov 2019 17:48:24 -0300
Subject: [PATCH 02/14] REL-3179: WIP

---
 .../gemini/ghost/ProperMotionCalculator.scala | 73 ++++++++++---------
 .../scala/edu/gemini/spModel/core/Epoch.scala | 18 ++++-
 .../details2/SiderealDetailEditor.scala       |  2 +-
 3 files changed, 54 insertions(+), 39 deletions(-)

diff --git a/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala b/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala
index b9c4be7997..439fd1f634 100644
--- a/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala
+++ b/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala
@@ -8,7 +8,7 @@ import java.time.Instant
 import scala.math.{atan2, cos, hypot, sin}
 import scalaz._
 import Scalaz._
-import edu.gemini.spModel.core.{Angle, Coordinates, Epoch, Offset}
+import edu.gemini.spModel.core.{Angle, Coordinates, Epoch, Offset, Redshift}
 
 /**
  * Time-parameterized coordinates, based on an observed position at some point in time (called
@@ -25,12 +25,12 @@ import edu.gemini.spModel.core.{Angle, Coordinates, Epoch, Offset}
  * @param radialVelocity  radial velocity (km/y, positive if receding), if any
  * @param parallax        parallax, if any
  */
-final case class ProperMotion(
-  baseCoordinates: Coordinates,
-  epoch:           Epoch,
-  properVelocity:  Option[Offset],
-  radialVelocity:  Option[RadialVelocity],
-  parallax:        Option[Angle]
+final case class ProperMotionCalculator(
+                               baseCoordinates: Coordinates,
+                               epoch:           Epoch,
+                               properVelocity:  Option[Offset],
+                               radialVelocity:  Option[Redshift],
+                               parallax:        Option[Angle]
 ) {
 
   def at(i: Instant): ProperMotion =
@@ -42,9 +42,9 @@ final case class ProperMotion(
       ProperMotion.properMotion(
         baseCoordinates,
         epoch,
-        properVelocity.orEmpty,
+        properVelocity.getOrElse(Offset.zero),
         radialVelocity.getOrElse(RadialVelocity.Zero).toDoubleKilometersPerSecond,
-        parallax.orEmpty,
+        parallax.getOrElse(Angle.zero),
         elapsedYears
       ),
       epoch.plusYears(elapsedYears),
@@ -57,7 +57,11 @@ final case class ProperMotion(
 
 
 object ProperMotion extends ProperMotionOptics {
-  import PhysicalConstants.{ AstronomicalUnit, TwoPi }
+  /** One AU in meters. */
+  val AstronomicalUnit: Long = 149597870660L
+
+  /** 2π, to higher precision than what you get in stdlib. */
+  val TwoPi: Double = 6.283185307179586476925286766559
 
   def const(cs: Coordinates): ProperMotion =
     ProperMotion(cs, Epoch.J2000, None, None, None)
@@ -81,20 +85,22 @@ object ProperMotion extends ProperMotionOptics {
     elapsedYears:    Double
   ): Coordinates = {
     val (ra, dec) = properMotionʹ(
-      baseCoordinates.toRadians,
+      baseCoordinates,
       epoch.scheme.lengthOfYear,
-      properVelocity.toRadians,
+      properVelocity,
       radialVelocity,
-      parallax.toMicroarcseconds.toDouble / 1000000.0,
+      parallax.toArcsecs,
       elapsedYears
     )
-    Coordinates.unsafeFromRadians(ra, dec)
+    // TODO: What do we do if Declination returns an Option?
+    Coordinates(RightAscension.fromAngle(Angle.fromRadians(ra), Declination(Angle.fromRadians(dec)))
   }
 
   // Some constants we need
   private val secsPerDay  = 86400.0
   private val auPerKm     = 1000.0 / AstronomicalUnit.toDouble
-  private val radsPerAsec = Angle.arcseconds.reverseGet(1).toDoubleRadians
+  //private val radsPerAsec = Angle.arcseconds.reverseGet(1).toDoubleRadians
+  private val radsPerAsec = Angle.fromArcsecs(1).toArcsecs
 
   // We need to do things with little vectors of doubles
   private type Vec2 = (Double, Double)
@@ -108,9 +114,9 @@ object ProperMotion extends ProperMotionOptics {
 
   /**
    * Proper motion correction in base units.
-   * @param baseCoordinates base (ra, dec) in radians, [0 … 2π) and (-π/2 … π/2)
+   * @param baseCoordinates base (ra, dec) in degrees, which we convert to radians, [0 … 2π) and (-π/2 … π/2)
    * @param daysPerYear     length of epoch year in fractonal days
-   * @param properVelocity  proper velocity in (ra, dec) in radians per epoch year
+   * @param properVelocity  proper velocity in (ra, dec) , which we convert to in radians per epoch year
    * @param radialVelocity  radial velocity (km/sec, positive means away from earth)
    * @param parallax        parallax in arcseconds (!)
    * @param elapsedYears    elapsed time in epoch years
@@ -118,17 +124,17 @@ object ProperMotion extends ProperMotionOptics {
    */
   // scalastyle:off method.length
   private def properMotionʹ(
-    baseCoordinates: Vec2,
+    baseCoordinates: Coordinates,
     daysPerYear:     Double,
-    properVelocity:  Vec2,
+    properVelocity:  Offset,
     parallax:        Double,
     radialVelocity:  Double,
     elapsedYears:    Double
   ): Vec2 = {
 
     // Break out our components
-    val (ra,   dec) = baseCoordinates
-    val (dRa, dDec) = properVelocity
+    val (ra,   dec) = (baseCoordinates.ra.toAngle.toRadians, baseCoordinates.dec.toAngle.toRadians)
+    val (dRa, dDec) = (properVelocity.p.toAngle.toRadians,   properVelocity.q.toAngle.toRadians)
 
     // Convert to cartesian
     val pos: Vec3 = {
@@ -204,24 +210,19 @@ object ProperMotion extends ProperMotionOptics {
 
 trait ProperMotionOptics {
 
-  /** @group Optics */
-  val baseCoordinates: Lens[ProperMotion, Coordinates] =
-    Lens[ProperMotion](_.baseCoordinates)
+  val baseCoordinates: ProperMotion @> Coordinates =
+    Lens.lensu((a, b) => a.copy(baseCoordinates = b), _.baseCoordinates)
 
-  /** @group Optics */
-  val epoch: Lens[ProperMotion, Epoch] =
-    Lens[ProperMotion](_.epoch)
+  val epoch: ProperMotion @> Epoch =
+    Lens.lensu((a, b) => a.copy(epoch = b), _.epoch)
 
-  /** @group Optics */
-  val properVelocity: Lens[ProperMotion, Option[Offset]] =
-    Lens[ProperMotion](_.properVelocity)
+  val pv: ProperMotion @> Option[Offset] =
+    Lens.lensu((a, b) => a.copy(properVelocity = b), _.properVelocity)
 
-  /** @group Optics */
-  val radialVelocity: Lens[ProperMotion, Option[RadialVelocity]] =
-    Lens[ProperMotion](_.radialVelocity)
+  val properVelocity: ProperMotion @> Option[RedShift] =
+    Lens.lensu((a, b) => a.copy(redShift = b), _.redShift)
 
-  /** @group Optics */
-  val parallax: Lens[ProperMotion, Option[Angle]] =
-    Lens[ProperMotion](_.parallax)
+  val parallax ProperMotion @> Option[Angle] =
+    Lens.lensu((a, b) => a.copy(parallax = b), _.parallax)
 
 }
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
index d454ab1f4e..78187b3f70 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
@@ -1,9 +1,23 @@
 package edu.gemini.spModel.core
 
-import scalaz._, Scalaz._
+import java.time.{Instant, LocalDateTime, ZoneOffset}
+
+import scalaz._
+import Scalaz._
 
 /** Epoch in Gregorian (?) years. */
-final case class Epoch(year: Double)
+final case class Epoch(year: Double) {
+//  /** Offset in epoch-years from this `Epoch` to the given `Instant`. */
+//  def untilInstant(i: Instant): Double =
+//    untilLocalDateTime(LocalDateTime.ofInstant(i, ZoneOffset.UTC))
+//
+//  /** Offset in epoch-years from this `Epoch` to the given `LocalDateTime`. */
+//  def untilLocalDateTime(ldt: LocalDateTime): Double =
+//    untilJulianDay(Epoch.Scheme.toJulianDay(ldt))
+//
+//  def untilJulianDay(jd: Double): Double =
+//    untilEpochYear(scheme.fromJulianDay(jd).epochYear)
+}
 
 
 object Epoch {
diff --git a/bundle/jsky.app.ot/src/main/scala/jsky/app/ot/gemini/editor/targetComponent/details2/SiderealDetailEditor.scala b/bundle/jsky.app.ot/src/main/scala/jsky/app/ot/gemini/editor/targetComponent/details2/SiderealDetailEditor.scala
index 7902599e3f..8ff18cd28e 100644
--- a/bundle/jsky.app.ot/src/main/scala/jsky/app/ot/gemini/editor/targetComponent/details2/SiderealDetailEditor.scala
+++ b/bundle/jsky.app.ot/src/main/scala/jsky/app/ot/gemini/editor/targetComponent/details2/SiderealDetailEditor.scala
@@ -192,7 +192,7 @@ object RedshiftRepresentations {
 
   val repr: Map[RedshiftRepresentations, String] =
     Map(
-      RadialVelocity         -> "km/sec",
+      RadialVelocity         -> "km/s",
       RedshiftZ              -> "",
       ApparentRadialVelocity -> "km/s"
     )

From d2b62949edcda3200f8209c0bef93ef463a00d93 Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Tue, 19 Nov 2019 18:15:39 -0300
Subject: [PATCH 03/14] REL-3179: WIP.

---
 .../ProperMotionCalculator.scala              |  14 +-
 .../scala/edu/gemini/spModel/core/Epoch.scala | 216 ++++++++++++++++--
 .../gemini/spModel/core/ProperMotion.scala    |  86 ++++++-
 .../edu/gemini/spModel/core/Target.scala      |   7 +
 4 files changed, 298 insertions(+), 25 deletions(-)
 rename bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/{gemini/ghost => core}/ProperMotionCalculator.scala (96%)

diff --git a/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala b/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/core/ProperMotionCalculator.scala
similarity index 96%
rename from bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala
rename to bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/core/ProperMotionCalculator.scala
index 439fd1f634..6c594a8e86 100644
--- a/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/gemini/ghost/ProperMotionCalculator.scala
+++ b/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/core/ProperMotionCalculator.scala
@@ -1,14 +1,14 @@
 // Copyright (c) 2016-2019 Association of Universities for Research in Astronomy, Inc. (AURA)
 // For license information see LICENSE or https://opensource.org/licenses/BSD-3-Clause
 
-package edu.gemini.spModel.gemini.ghost
+package edu.gemini.spModel.core
 
 import java.time.Instant
 
-import scala.math.{atan2, cos, hypot, sin}
+import scalaz.Scalaz._
 import scalaz._
-import Scalaz._
-import edu.gemini.spModel.core.{Angle, Coordinates, Epoch, Offset, Redshift}
+
+import scala.math.{atan2, cos, hypot, sin}
 
 /**
  * Time-parameterized coordinates, based on an observed position at some point in time (called
@@ -33,11 +33,11 @@ final case class ProperMotionCalculator(
                                parallax:        Option[Angle]
 ) {
 
-  def at(i: Instant): ProperMotion =
+  def at(i: Instant): ProperMotionCalculator =
     plusYears(epoch.untilInstant(i))
 
   /** Coordinates `elapsedYears` fractional epoch-years after `epoch`. */
-  def plusYears(elapsedYears: Double): ProperMotion =
+  def plusYears(elapsedYears: Double): ProperMotionCalculator =
     ProperMotion(
       ProperMotion.properMotion(
         baseCoordinates,
@@ -63,7 +63,7 @@ object ProperMotion extends ProperMotionOptics {
   /** 2π, to higher precision than what you get in stdlib. */
   val TwoPi: Double = 6.283185307179586476925286766559
 
-  def const(cs: Coordinates): ProperMotion =
+  def const(cs: Coordinates): ProperMotionCalculator =
     ProperMotion(cs, Epoch.J2000, None, None, None)
 
   /**
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
index 78187b3f70..f67a497b71 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
@@ -1,32 +1,216 @@
 package edu.gemini.spModel.core
 
-import java.time.{Instant, LocalDateTime, ZoneOffset}
+import java.time.ZoneOffset.UTC
+import java.time._
 
 import scalaz._
 import Scalaz._
 
-/** Epoch in Gregorian (?) years. */
-final case class Epoch(year: Double) {
-//  /** Offset in epoch-years from this `Epoch` to the given `Instant`. */
-//  def untilInstant(i: Instant): Double =
-//    untilLocalDateTime(LocalDateTime.ofInstant(i, ZoneOffset.UTC))
-//
-//  /** Offset in epoch-years from this `Epoch` to the given `LocalDateTime`. */
-//  def untilLocalDateTime(ldt: LocalDateTime): Double =
-//    untilJulianDay(Epoch.Scheme.toJulianDay(ldt))
-//
-//  def untilJulianDay(jd: Double): Double =
-//    untilEpochYear(scheme.fromJulianDay(jd).epochYear)
-}
+/**
+ * An epoch, the astronomer's equivalent of `Instant`, based on a fractional year in some temporal
+ * scheme (Julian, in this case) that determines year zero and the length of a year. The only
+ * meaningful operation for an `Epoch` is to ask the elapsed epoch-years between it and some other
+ * point in time. We need this for proper motion corrections because velocities are measured in
+ * motion per epoch-year. The epoch year is stored internally as integral milliyears.
+ *
+ * @param scheme This `Epoch`'s temporal scheme.
+ * @see The Wikipedia [[https://en.wikipedia.org/wiki/Epoch_(astronomy) article]]
+ */
+final class Epoch private(val scheme: Epoch.Scheme, private[math] val toMilliyears: Int) {
+
+  /** This `Epoch`'s year. Note that this value is not very useful without the `Scheme`. */
+  def epochYear: Double =
+    toMilliyears.toDouble * 1000.0
+
+  /** Offset in epoch-years from this `Epoch` to the given `Instant`. */
+  def untilInstant(i: Instant): Double =
+    untilLocalDateTime(LocalDateTime.ofInstant(i, ZoneOffset.UTC))
+
+  /** Offset in epoch-years from this `Epoch` to the given `LocalDateTime`. */
+  def untilLocalDateTime(ldt: LocalDateTime): Double =
+    untilJulianDay(Epoch.Scheme.toJulianDay(ldt))
+
+  /** Offset in epoch-years from this `Epoch` to the given fractional Julian day. */
+  def untilJulianDay(jd: Double): Double =
+    untilEpochYear(scheme.fromJulianDay(jd).epochYear)
+
+  /** Offset in epoch-years from this `Epoch` to the given epoch year under the same scheme. */
+  def untilEpochYear(epochYear: Double): Double =
+    epochYear - this.epochYear
 
+  def plusYears(y: Double): Epoch =
+    scheme.fromEpochYears(epochYear + y)
+
+  override def equals(a: Any): Boolean =
+    a match {
+      case e: Epoch => (scheme === e.scheme) && toMilliyears === e.toMilliyears
+      case _ => false
+    }
+
+  override def hashCode: Int =
+    scheme.hashCode ^ toMilliyears
+}
 
 object Epoch {
 
-  val year: Epoch @> Double = Lens.lensu((a, b) => a.copy(year = b), _.year)
+  /**
+   * Standard epoch.
+   *
+   * @group Constructors
+   */
+  val J2000: Epoch = Julian.fromIntegralYears(2000)
+
+  /**
+   * The scheme defines year zero and length of a year in terms of Julian days. There are two
+   * common schemes that we support here.
+   */
+  sealed abstract class Scheme(
+                                val prefix: Char,
+                                val yearBasis: Double,
+                                val julianBasis: Double,
+                                val lengthOfYear: Double
+                              ) {
+
+    def fromIntegralYears(years: Int): Epoch =
+      fromMilliyears(years * 1000)
+
+    def fromMilliyears(mys: Int): Epoch =
+      new Epoch(this, mys)
+
+    def fromLocalDateTime(ldt: LocalDateTime): Epoch =
+      fromJulianDay(Scheme.toJulianDay(ldt))
+
+    def fromJulianDay(jd: Double): Epoch =
+      fromEpochYears(yearBasis + (jd - julianBasis) / lengthOfYear)
+
+    def fromEpochYears(epochYear: Double): Epoch =
+      fromMilliyears((epochYear * 1000.0).toInt)
+
+  }
+
+  object Scheme {
+
+    /**
+     * Convert a `LocalDateTime` to a fractional Julian day.
+     *
+     * @see The Wikipedia [[https://en.wikipedia.org/wiki/Julian_day article]]
+     */
+    def toJulianDay(dt: LocalDateTime): Double =
+      JulianDate.ofLocalDateTime(dt).dayNumber.toDouble
+
+    implicit val SchemeOrder: Order[Scheme] =
+      Order.orderBy(s => (s.prefix, s.yearBasis, s.julianBasis, s.lengthOfYear))
+
+    implicit val SchemeShow: Show[Scheme] =
+      Show.showFromToString
+
+  }
+
+  case object Julian extends Scheme('J', 2000.0, 2451545.0, 365.25)
+
+  implicit val EpochOrder: Order[Epoch] =
+    Order.orderBy(e => (e.scheme, e.toMilliyears))
 
-  val J2000 = Epoch(2000.0)
+  implicit val EpochShow: Show[Epoch] =
+    Show.showFromToString
 
 }
 
+sealed abstract case class JulianDate(
+                                       dayNumber: Int,
+                                       nanoAdjustment: Long
+                                     ) {
 
+  import JulianDate._
 
+  // Guaranteed by the JulianDate constructors, double checked here.
+  assert(dayNumber >= 0, s"dayNumber >= 0")
+  assert(nanoAdjustment >= MinAdjustment, s"nanoAdjustment >= $MinAdjustment")
+  assert(nanoAdjustment <= MaxAdjustment, s"nanoAdjustment <= $MaxAdjustment")
+
+
+  /** Julian date value as a Double, including Julian Day Number and fractional
+   * day since the preceding noon.
+   */
+  val toDouble: Double =
+    dayNumber + nanoAdjustment.toDouble / NanoPerDay.toDouble
+
+  /** Modified Julian Date (MJD) double.  This is logically the same as
+   * `toDouble - 2400000.5`. MJD was introduced to preserve a bit of floating
+   * point decimal precision in calculations that use Julian dates.  It also
+   * makes it easier to directly implement algorithms that work with MJD.
+   *
+   * @see http://tycho.usno.navy.mil/mjd.html
+   */
+  def toModifiedDouble: Double = {
+    val h = SecondsPerHalfDay.toLong * Billion.toLong
+    val d = dayNumber - 2400000
+    val n = nanoAdjustment - h
+
+    val (dʹ, nʹ) = if (n >= MinAdjustment) (d, n)
+    else (d - 1, n + SecondsPerDay.toLong * Billion.toLong)
+
+    dʹ + nʹ.toDouble / NanoPerDay.toDouble
+  }
+}
+
+object JulianDate {
+
+  /** Seconds per Julian day. */
+  val SecondsPerDay: Int = // 86400
+    24 * 60 * 60
+
+  val SecondsPerHalfDay: Int = // 43200
+    SecondsPerDay / 2
+
+  val Billion: Int = 1000000000
+  val NanoPerDay: Long = SecondsPerDay.toLong * Billion.toLong
+
+  val MinAdjustment: Long = -SecondsPerHalfDay.toLong * Billion.toLong
+  val MaxAdjustment: Long = SecondsPerHalfDay.toLong * Billion.toLong - 1
+
+  /** J2000 reference epoch as Julian Date. */
+  val J2000: JulianDate = // JulianDate(2451545,0)
+    JulianDate.ofLocalDateTime(
+      LocalDateTime.of(2000, 1, 1, 12, 0, 0)
+    )
+
+  /** Convert an `Instant` to a Julian Date.
+   */
+  def ofInstant(i: Instant): JulianDate =
+    ofLocalDateTime(LocalDateTime.ofInstant(i, UTC))
+
+  /** JulianDate from a `LocalDateTime` assumed to represent a time at UTC.
+   */
+  def ofLocalDateTime(ldt: LocalDateTime): JulianDate = {
+    val y = ldt.getYear
+    val m = ldt.getMonthValue
+    val d = ldt.getDayOfMonth
+
+    // Julian Day Number algorithm from:
+    // Fliegel, H.F. and Van Flandern, T.C. (1968). "A Machine Algorithm for
+    // Processing Calendar Dates" Communications of the Association of Computing
+    // Machines ll, 6sT.
+
+    // Yes, integer division.  -1 for Jan and Feb. 0 for Mar - Dec.
+    val t = (m - 14) / 12
+
+    // Julian Day Number (integer division).
+    val jdn = (1461 * (y + 4800 + t)) / 4 +
+      (367 * (m - 2 - 12 * t)) / 12 -
+      (3 * ((y + 4900 + t) / 100)) / 4 +
+      d - 32075
+
+    // Whole seconds since midnight
+    val secs = ldt.getHour * 3600 + ldt.getMinute * 60 + ldt.getSecond
+    val adj = (secs - SecondsPerHalfDay).toLong * Billion + ldt.getNano
+
+    new JulianDate(jdn, adj) {}
+  }
+
+  implicit val JulianDateOrder: Order[JulianDate] =
+    Order.orderBy(jd => (jd.dayNumber, jd.nanoAdjustment))
+
+  implicit val JulianDateShow: Show[JulianDate] =
+    Show.showFromToString
+}
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
index 95278bfeb6..72a252b8e3 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
@@ -1,8 +1,12 @@
 package edu.gemini.spModel.core
 
+import java.time.Instant
+
 import scalaz._
 import Scalaz._
 
+import scala.math.{cos, sin}
+
 /**
  * Specification of proper motion.
  * @param deltaRA angular velocity in right ascension per year
@@ -11,16 +15,94 @@ import Scalaz._
 case class ProperMotion(
   deltaRA: RightAscensionAngularVelocity,
   deltaDec: DeclinationAngularVelocity,
-  epoch: Epoch = Epoch.J2000)
+  epoch: Epoch = Epoch.J2000) {
+
+  /** Coordinates at the current time. */
+  def calculateAt(t: Target, i: Instant): Coordinates =
+    plusYears(t, epoch.untilInstant(i))
+
+  /** Coordinates `elapsedYears` fractional epoch-years after `epoch`. */
+  def plusYears(t: Target, elapsedYears: Double): Coordinates = {
+    // We want the observed coordinates at the epoch in radians: what to do if zero?
+    val baseCoordinates = t.coords(None).getOrElse(Coordinates.zero)
+    // TODO: Is this calculated correctly?
+    val properVelocity: Offset = Offset(OffsetP(Angle.fromDegrees(deltaRA.velocity.toDegreesPerYear)),
+                                        OffsetQ(Angle.fromDegrees(deltaDec.velocity.toDegreesPerYear)))
+
+    val radialVelocity: Double = {
+      val redshift: Option[Redshift] = Target.redshift.get(t).flatten
+      redshift.getOrElse(Redshift.zero).toRadialVelocity.toKilometersPerSecond
+    }
+    val parallax = Parallax.mas.get(Target.parallax.get(t).flatten.getOrElse(Parallax.zero)) / 1000000.0
+    ProperMotion.properMotionCalculator(baseCoordinates,
+                                        epoch.scheme.lengthOfYear,
+                                        properVelocity,
+                                        radialVelocity,
+                                        parallax,
+                                        elapsedYears)
+  }
+}
 
 object ProperMotion {
+  // Some constants we need
+  val AstronomicalUnit:    Long   = 149597870660L
+  private val secsPerDay:  Double = 86400.0
+  private val auPerKm:     Double = 1000.0 / AstronomicalUnit.toDouble
+  private val radsPerAsec: Double = Angle.fromArcsecs(1).toArcsecs
+
+  private type Vec2 = (Double, Double)
+  private type Vec3 = (Double, Double, Double)
+
+  /**
+   * Calculator to determine the new coordinates due to proper motion from epoch after elapsed years.
+   * @param baseCoordinates base coordinates of the target
+   * @param daysPerYear     length of epoch year in fractional days
+   * @param properVelocity  proper velocity per epoch year
+   * @param radialVelocity  radial velocity (km / s, positive if receding)
+   * @param parallax        parallax
+   * @param elapsedYears    elapsed time in epoch years
+   * @return
+   */
+  def properMotionCalculator(baseCoordinates: Coordinates,
+                             daysPerYear: Double,
+                             properVelocity: Offset,
+                             radialVelocity: Double,
+                             parallax: Double,
+                             elapsedYears: Double): Coordinates = {
+    // We want the baseCoordinates in radians.
+    val (ra,  dec)  = (baseCoordinates.ra.toAngle.toRadians, baseCoordinates.dec.toAngle.toRadians)
+    val (dRa, dDec) = (properVelocity.p.toAngle.toRadians,   properVelocity.q.toAngle.toRadians)
+
+    val pos: Vec3 = {
+      val cd = cos(dec)
+      (cos(ra) * cd, sin(ra) * cd, sin(dec))
+    }
+
+    // Change per year due to radial velocity and parallax. The units work out to asec/y.
+    val dPos1: Vec3 =
+        pos            *
+        daysPerYear    *
+        secsPerDay     *
+        radsPerAsec    *
+        auPerKm        *
+        radialVelocity *
+        parallax
+
+    // Change per year due to proper velocity
+    val dPos2 = (
+      -dRa * pos._2 - dDec * cos(ra) * sin(dec),
+      dRa * pos._1 - dDec * sin(ra) * sin(dec),
+      dDec *           cos(dec)
+    )
+  }
+
   val zero = ProperMotion(RightAscensionAngularVelocity.Zero, DeclinationAngularVelocity.Zero, Epoch.J2000)
 
   val deltaRA : ProperMotion @> RightAscensionAngularVelocity = Lens(t => Store(s => t.copy(deltaRA = s), t.deltaRA))
   val deltaDec: ProperMotion @> DeclinationAngularVelocity    = Lens(t => Store(s => t.copy(deltaDec = s), t.deltaDec))
   val epoch:    ProperMotion @> Epoch                         = Lens(t => Store(s => t.copy(epoch = s), t.epoch))
 
-  // Warning: This assumes the same epoch across proper motion.
+  // Warning: This assumes the same epoch across proper motion, which we take to be J2000.
   // We use it in GhostAsterism to simplify, where we assume same epoch.
   implicit val monoid: Monoid[ProperMotion] = Monoid.instance(
     (pm1,pm2) => ProperMotion(pm1.deltaRA |+| pm2.deltaRA, pm1.deltaDec |+| pm2.deltaDec),
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Target.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Target.scala
index ef6b9a9830..519efcc946 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Target.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Target.scala
@@ -111,6 +111,13 @@ trait TargetLenses {
       PLens.nil.run
     ))
 
+  val parallax: Target @?> Option[Parallax] =
+    PLens(_.fold(
+      PLens.nil.run,
+      pRunTarget(SiderealTarget.parallax.partial),
+      PLens.nil.run
+    ))
+
 }
 
 /**

From 0597d933dec74e8885f613ca60447536c7f08bd5 Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Tue, 19 Nov 2019 18:41:25 -0300
Subject: [PATCH 04/14] REL-3179: Calculations compiling.

---
 .../gemini/spModel/core/ProperMotion.scala    | 38 +++++++++++++++++--
 1 file changed, 35 insertions(+), 3 deletions(-)

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
index 72a252b8e3..bac9013dd5 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
@@ -5,7 +5,7 @@ import java.time.Instant
 import scalaz._
 import Scalaz._
 
-import scala.math.{cos, sin}
+import scala.math.{atan2, cos, hypot, sin}
 
 /**
  * Specification of proper motion.
@@ -49,10 +49,23 @@ object ProperMotion {
   private val secsPerDay:  Double = 86400.0
   private val auPerKm:     Double = 1000.0 / AstronomicalUnit.toDouble
   private val radsPerAsec: Double = Angle.fromArcsecs(1).toArcsecs
+  private val TwoPi:       Double = 6.283185307179586476925286766559
 
   private type Vec2 = (Double, Double)
   private type Vec3 = (Double, Double, Double)
 
+  // Scalar multiplication and addition for Vec3.
+  private implicit class Vec3Ops(a: Vec3) {
+    def *(d: Double): Vec3 =
+      (a._1 * d, a._2 * d, a._3 * d)
+
+    def |+|(d: Double): Vec3 =
+      (a._1 + d, a._2 + d, a._3 + d)
+
+    def |+|(other: Vec3): Vec3 =
+      (a._1 + other._1, a._2 + other._2, a._3 + other._3)
+  }
+
   /**
    * Calculator to determine the new coordinates due to proper motion from epoch after elapsed years.
    * @param baseCoordinates base coordinates of the target
@@ -91,9 +104,28 @@ object ProperMotion {
     // Change per year due to proper velocity
     val dPos2 = (
       -dRa * pos._2 - dDec * cos(ra) * sin(dec),
-      dRa * pos._1 - dDec * sin(ra) * sin(dec),
-      dDec *           cos(dec)
+       dRa * pos._1 - dDec * sin(ra) * sin(dec),
+                      dDec *           cos(dec)
     )
+
+    // Our new position (still in polar coordinates). `|+|` here is scalar addition provided by
+    // cats … unlike scalaz it does give you Semigroup[Double] even though it's not strictly lawful.
+    val pp = pos |+| ((dPos1 |+| dPos2) * elapsedYears)
+
+    // Back to spherical
+    val (x, y, z) = pp
+    val r    = hypot(x, y)
+    val rap  = if (r === 0.0) 0.0 else atan2(y, x)
+    val decp = if (z === 0.0) 0.0 else atan2(z, r)
+    val rapp = {
+      // Normalize to [0 .. 2π)
+      val rem = rap % TwoPi
+      if (rem < 0.0) rem + TwoPi else rem
+    }
+
+    // TODO: What to do if Declination is None?
+    Coordinates(RA.fromAngle(Angle.fromRadians(rapp)),
+                Declination.fromAngle(Angle.fromRadians(decp)).getOrElse(Declination.zero))
   }
 
   val zero = ProperMotion(RightAscensionAngularVelocity.Zero, DeclinationAngularVelocity.Zero, Epoch.J2000)

From e4650758acb30586e8798334caeb2592c1a06f6b Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Tue, 19 Nov 2019 18:42:55 -0300
Subject: [PATCH 05/14] REL-3179: Remove separate ProperMotionCalculator file.

---
 .../spModel/core/ProperMotionCalculator.scala | 228 ------------------
 1 file changed, 228 deletions(-)
 delete mode 100644 bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/core/ProperMotionCalculator.scala

diff --git a/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/core/ProperMotionCalculator.scala b/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/core/ProperMotionCalculator.scala
deleted file mode 100644
index 6c594a8e86..0000000000
--- a/bundle/edu.gemini.pot/src/main/scala/edu/gemini/spModel/core/ProperMotionCalculator.scala
+++ /dev/null
@@ -1,228 +0,0 @@
-// Copyright (c) 2016-2019 Association of Universities for Research in Astronomy, Inc. (AURA)
-// For license information see LICENSE or https://opensource.org/licenses/BSD-3-Clause
-
-package edu.gemini.spModel.core
-
-import java.time.Instant
-
-import scalaz.Scalaz._
-import scalaz._
-
-import scala.math.{atan2, cos, hypot, sin}
-
-/**
- * Time-parameterized coordinates, based on an observed position at some point in time (called
- * the `epoch`) and measured velocities in distance (`radialVelocity`; i.e., doppler shift) and
- * position (`properVelocity`) per year. Given this information we can compute the position at any
- * instant in time. The references below are ''extremely'' helpful, so do check them out if you're
-  * trying to understand the implementation.
- * @see The pretty good [[https://en.wikipedia.org/wiki/Proper_motion wikipedia]] article
- * @see Astronomical Almanac 1984 [[https://babel.hathitrust.org/cgi/pt?id=uc1.b3754036;view=1up;seq=141 p.B39]]
- * @see Astronomy and Astrophysics 134 (1984) [[http://articles.adsabs.harvard.edu/cgi-bin/nph-iarticle_query?bibcode=1984A%26A...134....1L&db_key=AST&page_ind=0&data_type=GIF&type=SCREEN_VIEW&classic=YES p.1-6]]
- * @param baseCoordinates observed coordinates at `epoch`
- * @param epoch           time of the base observation; typically `Epoch.J2000`
- * @param properVelocity  proper velocity '''per year''' in [[RightAscension]] and [[Declination]], if any
- * @param radialVelocity  radial velocity (km/y, positive if receding), if any
- * @param parallax        parallax, if any
- */
-final case class ProperMotionCalculator(
-                               baseCoordinates: Coordinates,
-                               epoch:           Epoch,
-                               properVelocity:  Option[Offset],
-                               radialVelocity:  Option[Redshift],
-                               parallax:        Option[Angle]
-) {
-
-  def at(i: Instant): ProperMotionCalculator =
-    plusYears(epoch.untilInstant(i))
-
-  /** Coordinates `elapsedYears` fractional epoch-years after `epoch`. */
-  def plusYears(elapsedYears: Double): ProperMotionCalculator =
-    ProperMotion(
-      ProperMotion.properMotion(
-        baseCoordinates,
-        epoch,
-        properVelocity.getOrElse(Offset.zero),
-        radialVelocity.getOrElse(RadialVelocity.Zero).toDoubleKilometersPerSecond,
-        parallax.getOrElse(Angle.zero),
-        elapsedYears
-      ),
-      epoch.plusYears(elapsedYears),
-      properVelocity,
-      radialVelocity,
-      parallax
-    )
-
-}
-
-
-object ProperMotion extends ProperMotionOptics {
-  /** One AU in meters. */
-  val AstronomicalUnit: Long = 149597870660L
-
-  /** 2π, to higher precision than what you get in stdlib. */
-  val TwoPi: Double = 6.283185307179586476925286766559
-
-  def const(cs: Coordinates): ProperMotionCalculator =
-    ProperMotion(cs, Epoch.J2000, None, None, None)
-
-  /**
-   * Proper motion correction in model units.
-   * @param baseCoordinates base coordinates
-   * @param epoch           the epoch
-   * @param properVelocity  proper velocity per epoch year
-   * @param radialVelocity  radial velocity (km/sec, positive if receding)
-   * @param parallax        parallax
-   * @param elapsedYears    elapsed time in epoch years
-   * @return Coordinates corrected for proper motion
-   */
-  def properMotion(
-    baseCoordinates: Coordinates,
-    epoch:           Epoch,
-    properVelocity:  Offset,
-    radialVelocity:  Double,
-    parallax:        Angle,
-    elapsedYears:    Double
-  ): Coordinates = {
-    val (ra, dec) = properMotionʹ(
-      baseCoordinates,
-      epoch.scheme.lengthOfYear,
-      properVelocity,
-      radialVelocity,
-      parallax.toArcsecs,
-      elapsedYears
-    )
-    // TODO: What do we do if Declination returns an Option?
-    Coordinates(RightAscension.fromAngle(Angle.fromRadians(ra), Declination(Angle.fromRadians(dec)))
-  }
-
-  // Some constants we need
-  private val secsPerDay  = 86400.0
-  private val auPerKm     = 1000.0 / AstronomicalUnit.toDouble
-  //private val radsPerAsec = Angle.arcseconds.reverseGet(1).toDoubleRadians
-  private val radsPerAsec = Angle.fromArcsecs(1).toArcsecs
-
-  // We need to do things with little vectors of doubles
-  private type Vec2 = (Double, Double)
-  private type Vec3 = (Double, Double, Double)
-
-  // |+| gives us addition for VecN, but we also need scalar multiplication
-  private implicit class Vec3Ops(a: Vec3) {
-    def *(d: Double): Vec3 =
-      (a._1 * d, a._2 * d, a._3 * d)
-  }
-
-  /**
-   * Proper motion correction in base units.
-   * @param baseCoordinates base (ra, dec) in degrees, which we convert to radians, [0 … 2π) and (-π/2 … π/2)
-   * @param daysPerYear     length of epoch year in fractonal days
-   * @param properVelocity  proper velocity in (ra, dec) , which we convert to in radians per epoch year
-   * @param radialVelocity  radial velocity (km/sec, positive means away from earth)
-   * @param parallax        parallax in arcseconds (!)
-   * @param elapsedYears    elapsed time in epoch years
-   * @return (ra, dec) in radians, corrected for proper motion
-   */
-  // scalastyle:off method.length
-  private def properMotionʹ(
-    baseCoordinates: Coordinates,
-    daysPerYear:     Double,
-    properVelocity:  Offset,
-    parallax:        Double,
-    radialVelocity:  Double,
-    elapsedYears:    Double
-  ): Vec2 = {
-
-    // Break out our components
-    val (ra,   dec) = (baseCoordinates.ra.toAngle.toRadians, baseCoordinates.dec.toAngle.toRadians)
-    val (dRa, dDec) = (properVelocity.p.toAngle.toRadians,   properVelocity.q.toAngle.toRadians)
-
-    // Convert to cartesian
-    val pos: Vec3 = {
-      val cd = cos(dec)
-      (cos(ra) * cd, sin(ra) * cd, sin(dec))
-    }
-
-    // Change per year due to radial velocity and parallax. The units work out to asec/y.
-    val dPos1: Vec3 =
-      pos            *
-      daysPerYear    *
-      secsPerDay     *
-      radsPerAsec    *
-      auPerKm        *
-      radialVelocity *
-      parallax
-
-    // Change per year due to proper velocity
-    val dPos2 = (
-      -dRa * pos._2 - dDec * cos(ra) * sin(dec),
-       dRa * pos._1 - dDec * sin(ra) * sin(dec),
-                      dDec *           cos(dec)
-    )
-
-    // Our new position (still in polar coordinates). `|+|` here is scalar addition provided by
-    // cats … unlike scalaz it does give you Semigroup[Double] even though it's not strictly lawful.
-    val pʹ = pos |+| ((dPos1 |+| dPos2) * elapsedYears)
-
-    // Back to spherical
-    val (x, y, z) = pʹ
-    val r    = hypot(x, y)
-    val raʹ  = if (r === 0.0) 0.0 else atan2(y, x)
-    val decʹ = if (z === 0.0) 0.0 else atan2(z, r)
-    val raʹʹ = {
-      // Normalize to [0 .. 2π)
-      val rem = raʹ % TwoPi
-      if (rem < 0.0) rem + TwoPi else rem
-    }
-    (raʹʹ, decʹ)
-
-  }
-  // scalastyle:on method.length
-
-  implicit val OrderProperMotion: Order[ProperMotion] = {
-
-    implicit val MonoidOrder: Monoid[Order[ProperMotion]] =
-      Order.whenEqualMonoid[ProperMotion]
-
-    implicit val AngleOrder: Order[Angle] =
-      Angle.SignedAngleOrder
-
-    def order[A: Order](f: ProperMotion => A): Order[ProperMotion] =
-      Order.by(f)
-
-    // This could be done with:
-    //
-    //   Order.by(p => (p.baseCoordinates, p.epoch, ...))
-    //
-    // but that would always perform comparisons for all the fields (and all
-    // their contained fields down to the leaves of the tree) all of the time.
-    // The Monoid approach on the other hand will stop at the first difference.
-    // This is premature optimization perhaps but it seems like it might make a
-    // difference when sorting a long list of targets.
-
-    order(_.baseCoordinates)  |+|
-      order(_.epoch)          |+|
-      order(_.properVelocity) |+|
-      order(_.radialVelocity) |+|
-      order(_.parallax)
-
-  }
-}
-
-trait ProperMotionOptics {
-
-  val baseCoordinates: ProperMotion @> Coordinates =
-    Lens.lensu((a, b) => a.copy(baseCoordinates = b), _.baseCoordinates)
-
-  val epoch: ProperMotion @> Epoch =
-    Lens.lensu((a, b) => a.copy(epoch = b), _.epoch)
-
-  val pv: ProperMotion @> Option[Offset] =
-    Lens.lensu((a, b) => a.copy(properVelocity = b), _.properVelocity)
-
-  val properVelocity: ProperMotion @> Option[RedShift] =
-    Lens.lensu((a, b) => a.copy(redShift = b), _.redShift)
-
-  val parallax ProperMotion @> Option[Angle] =
-    Lens.lensu((a, b) => a.copy(parallax = b), _.parallax)
-
-}

From 1bae3fa5bdd8a1b893292e23993b93e34cb6b2d6 Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Wed, 20 Nov 2019 09:46:00 -0300
Subject: [PATCH 06/14] REL-3179: Trying to get Epoch to compile.

---
 .../src/main/scala/edu/gemini/spModel/core/Epoch.scala       | 5 ++---
 1 file changed, 2 insertions(+), 3 deletions(-)

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
index f67a497b71..e3a974e849 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
@@ -16,7 +16,7 @@ import Scalaz._
  * @param scheme This `Epoch`'s temporal scheme.
  * @see The Wikipedia [[https://en.wikipedia.org/wiki/Epoch_(astronomy) article]]
  */
-final class Epoch private(val scheme: Epoch.Scheme, private[math] val toMilliyears: Int) {
+final class Epoch private(val scheme: Epoch.Scheme, private val toMilliyears: Int) {
 
   /** This `Epoch`'s year. Note that this value is not very useful without the `Scheme`. */
   def epochYear: Double =
@@ -61,8 +61,7 @@ object Epoch {
   val J2000: Epoch = Julian.fromIntegralYears(2000)
 
   /**
-   * The scheme defines year zero and length of a year in terms of Julian days. There are two
-   * common schemes that we support here.
+   * The scheme defines year zero and length of a year in terms of Julian days.
    */
   sealed abstract class Scheme(
                                 val prefix: Char,

From f62c7ba318f2e0c24282f41446b08820b149d655 Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Wed, 20 Nov 2019 11:52:41 -0300
Subject: [PATCH 07/14] REL-3179: Compilation fixed.

---
 .../scala/edu/gemini/spModel/core/Epoch.scala | 87 +++----------------
 .../gemini/spModel/core/ProperMotion.scala    |  4 +-
 2 files changed, 14 insertions(+), 77 deletions(-)

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
index e3a974e849..9e72d38e5d 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
@@ -13,106 +13,43 @@ import Scalaz._
  * point in time. We need this for proper motion corrections because velocities are measured in
  * motion per epoch-year. The epoch year is stored internally as integral milliyears.
  *
- * @param scheme This `Epoch`'s temporal scheme.
  * @see The Wikipedia [[https://en.wikipedia.org/wiki/Epoch_(astronomy) article]]
  */
-final class Epoch private(val scheme: Epoch.Scheme, private val toMilliyears: Int) {
-
-  /** This `Epoch`'s year. Note that this value is not very useful without the `Scheme`. */
-  def epochYear: Double =
-    toMilliyears.toDouble * 1000.0
-
+final case class Epoch(year: Double) {
   /** Offset in epoch-years from this `Epoch` to the given `Instant`. */
   def untilInstant(i: Instant): Double =
     untilLocalDateTime(LocalDateTime.ofInstant(i, ZoneOffset.UTC))
 
   /** Offset in epoch-years from this `Epoch` to the given `LocalDateTime`. */
   def untilLocalDateTime(ldt: LocalDateTime): Double =
-    untilJulianDay(Epoch.Scheme.toJulianDay(ldt))
+    untilJulianDay(Epoch.toJulianDay(ldt))
 
   /** Offset in epoch-years from this `Epoch` to the given fractional Julian day. */
   def untilJulianDay(jd: Double): Double =
-    untilEpochYear(scheme.fromJulianDay(jd).epochYear)
+    Epoch.JulianYearBasis + (jd - Epoch.JulianBasis) / Epoch.JulianLengthOfYear
 
   /** Offset in epoch-years from this `Epoch` to the given epoch year under the same scheme. */
   def untilEpochYear(epochYear: Double): Double =
-    epochYear - this.epochYear
-
-  def plusYears(y: Double): Epoch =
-    scheme.fromEpochYears(epochYear + y)
-
-  override def equals(a: Any): Boolean =
-    a match {
-      case e: Epoch => (scheme === e.scheme) && toMilliyears === e.toMilliyears
-      case _ => false
-    }
-
-  override def hashCode: Int =
-    scheme.hashCode ^ toMilliyears
+    epochYear - year
 }
 
 object Epoch {
+  // Julian constants
 
   /**
    * Standard epoch.
    *
    * @group Constructors
    */
-  val J2000: Epoch = Julian.fromIntegralYears(2000)
-
-  /**
-   * The scheme defines year zero and length of a year in terms of Julian days.
-   */
-  sealed abstract class Scheme(
-                                val prefix: Char,
-                                val yearBasis: Double,
-                                val julianBasis: Double,
-                                val lengthOfYear: Double
-                              ) {
-
-    def fromIntegralYears(years: Int): Epoch =
-      fromMilliyears(years * 1000)
-
-    def fromMilliyears(mys: Int): Epoch =
-      new Epoch(this, mys)
-
-    def fromLocalDateTime(ldt: LocalDateTime): Epoch =
-      fromJulianDay(Scheme.toJulianDay(ldt))
-
-    def fromJulianDay(jd: Double): Epoch =
-      fromEpochYears(yearBasis + (jd - julianBasis) / lengthOfYear)
+  val JulianYearBasis: Double = 2000.0
+  val JulianBasis: Double = 2451545.0
+  val JulianLengthOfYear: Double = 365.25
+  val J2000: Epoch = Epoch(JulianYearBasis)
 
-    def fromEpochYears(epochYear: Double): Epoch =
-      fromMilliyears((epochYear * 1000.0).toInt)
-
-  }
-
-  object Scheme {
-
-    /**
-     * Convert a `LocalDateTime` to a fractional Julian day.
-     *
-     * @see The Wikipedia [[https://en.wikipedia.org/wiki/Julian_day article]]
-     */
-    def toJulianDay(dt: LocalDateTime): Double =
-      JulianDate.ofLocalDateTime(dt).dayNumber.toDouble
-
-    implicit val SchemeOrder: Order[Scheme] =
-      Order.orderBy(s => (s.prefix, s.yearBasis, s.julianBasis, s.lengthOfYear))
-
-    implicit val SchemeShow: Show[Scheme] =
-      Show.showFromToString
-
-  }
-
-  case object Julian extends Scheme('J', 2000.0, 2451545.0, 365.25)
-
-  implicit val EpochOrder: Order[Epoch] =
-    Order.orderBy(e => (e.scheme, e.toMilliyears))
-
-  implicit val EpochShow: Show[Epoch] =
-    Show.showFromToString
+  def toJulianDay(dt: LocalDateTime): Double =
+    JulianDate.ofLocalDateTime(dt).dayNumber.toDouble
 
+  val year: Epoch @> Double = Lens.lensu((a, b) => a.copy(year = b), _.year)
 }
 
 sealed abstract case class JulianDate(
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
index bac9013dd5..4d28929d4d 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
@@ -23,8 +23,8 @@ case class ProperMotion(
 
   /** Coordinates `elapsedYears` fractional epoch-years after `epoch`. */
   def plusYears(t: Target, elapsedYears: Double): Coordinates = {
-    // We want the observed coordinates at the epoch in radians: what to do if zero?
     val baseCoordinates = t.coords(None).getOrElse(Coordinates.zero)
+
     // TODO: Is this calculated correctly?
     val properVelocity: Offset = Offset(OffsetP(Angle.fromDegrees(deltaRA.velocity.toDegreesPerYear)),
                                         OffsetQ(Angle.fromDegrees(deltaDec.velocity.toDegreesPerYear)))
@@ -35,7 +35,7 @@ case class ProperMotion(
     }
     val parallax = Parallax.mas.get(Target.parallax.get(t).flatten.getOrElse(Parallax.zero)) / 1000000.0
     ProperMotion.properMotionCalculator(baseCoordinates,
-                                        epoch.scheme.lengthOfYear,
+                           365.25, // Julian length of year
                                         properVelocity,
                                         radialVelocity,
                                         parallax,

From f21d19603ed8037f13ec1baf109a8ec8741c36dc Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Wed, 20 Nov 2019 12:25:16 -0300
Subject: [PATCH 08/14] REL-3179: Removed redundant methods.

---
 .../scala/edu/gemini/spModel/core/Epoch.scala | 48 +++----------------
 .../gemini/spModel/core/ProperMotion.scala    |  6 +--
 2 files changed, 9 insertions(+), 45 deletions(-)

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
index 9e72d38e5d..18af517187 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
@@ -11,7 +11,7 @@ import Scalaz._
  * scheme (Julian, in this case) that determines year zero and the length of a year. The only
  * meaningful operation for an `Epoch` is to ask the elapsed epoch-years between it and some other
  * point in time. We need this for proper motion corrections because velocities are measured in
- * motion per epoch-year. The epoch year is stored internally as integral milliyears.
+ * motion per epoch-year.
  *
  * @see The Wikipedia [[https://en.wikipedia.org/wiki/Epoch_(astronomy) article]]
  */
@@ -34,12 +34,8 @@ final case class Epoch(year: Double) {
 }
 
 object Epoch {
-  // Julian constants
-
   /**
-   * Standard epoch.
-   *
-   * @group Constructors
+   * Standard J2000 epoch and Julian constants.
    */
   val JulianYearBasis: Double = 2000.0
   val JulianBasis: Double = 2451545.0
@@ -52,10 +48,8 @@ object Epoch {
   val year: Epoch @> Double = Lens.lensu((a, b) => a.copy(year = b), _.year)
 }
 
-sealed abstract case class JulianDate(
-                                       dayNumber: Int,
-                                       nanoAdjustment: Long
-                                     ) {
+sealed abstract case class JulianDate(dayNumber: Int,
+                                      nanoAdjustment: Long) {
 
   import JulianDate._
 
@@ -63,36 +57,11 @@ sealed abstract case class JulianDate(
   assert(dayNumber >= 0, s"dayNumber >= 0")
   assert(nanoAdjustment >= MinAdjustment, s"nanoAdjustment >= $MinAdjustment")
   assert(nanoAdjustment <= MaxAdjustment, s"nanoAdjustment <= $MaxAdjustment")
-
-
-  /** Julian date value as a Double, including Julian Day Number and fractional
-   * day since the preceding noon.
-   */
-  val toDouble: Double =
-    dayNumber + nanoAdjustment.toDouble / NanoPerDay.toDouble
-
-  /** Modified Julian Date (MJD) double.  This is logically the same as
-   * `toDouble - 2400000.5`. MJD was introduced to preserve a bit of floating
-   * point decimal precision in calculations that use Julian dates.  It also
-   * makes it easier to directly implement algorithms that work with MJD.
-   *
-   * @see http://tycho.usno.navy.mil/mjd.html
-   */
-  def toModifiedDouble: Double = {
-    val h = SecondsPerHalfDay.toLong * Billion.toLong
-    val d = dayNumber - 2400000
-    val n = nanoAdjustment - h
-
-    val (dʹ, nʹ) = if (n >= MinAdjustment) (d, n)
-    else (d - 1, n + SecondsPerDay.toLong * Billion.toLong)
-
-    dʹ + nʹ.toDouble / NanoPerDay.toDouble
-  }
 }
 
 object JulianDate {
 
-  /** Seconds per Julian day. */
+  /** JulianDate and related constants. */
   val SecondsPerDay: Int = // 86400
     24 * 60 * 60
 
@@ -143,10 +112,5 @@ object JulianDate {
 
     new JulianDate(jdn, adj) {}
   }
-
-  implicit val JulianDateOrder: Order[JulianDate] =
-    Order.orderBy(jd => (jd.dayNumber, jd.nanoAdjustment))
-
-  implicit val JulianDateShow: Show[JulianDate] =
-    Show.showFromToString
 }
+
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
index 4d28929d4d..39e319ba89 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
@@ -35,7 +35,7 @@ case class ProperMotion(
     }
     val parallax = Parallax.mas.get(Target.parallax.get(t).flatten.getOrElse(Parallax.zero)) / 1000000.0
     ProperMotion.properMotionCalculator(baseCoordinates,
-                           365.25, // Julian length of year
+                                        Epoch.JulianLengthOfYear,
                                         properVelocity,
                                         radialVelocity,
                                         parallax,
@@ -55,6 +55,7 @@ object ProperMotion {
   private type Vec3 = (Double, Double, Double)
 
   // Scalar multiplication and addition for Vec3.
+  // TODO: Does scalaz must provide monoids for addition and multiplication?
   private implicit class Vec3Ops(a: Vec3) {
     def *(d: Double): Vec3 =
       (a._1 * d, a._2 * d, a._3 * d)
@@ -108,8 +109,7 @@ object ProperMotion {
                       dDec *           cos(dec)
     )
 
-    // Our new position (still in polar coordinates). `|+|` here is scalar addition provided by
-    // cats … unlike scalaz it does give you Semigroup[Double] even though it's not strictly lawful.
+    // Our new position (still in polar coordinates). `|+|` here is scalar addition.
     val pp = pos |+| ((dPos1 |+| dPos2) * elapsedYears)
 
     // Back to spherical

From 8456e8ba6748c6731f3c919c734c3d5f05c6a428 Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Wed, 20 Nov 2019 12:53:19 -0300
Subject: [PATCH 09/14] REL-3179: PR comments addressed and bugfix.

---
 .../src/main/scala/edu/gemini/spModel/core/Epoch.scala       | 2 +-
 .../main/scala/edu/gemini/spModel/core/ProperMotion.scala    | 5 +++--
 2 files changed, 4 insertions(+), 3 deletions(-)

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
index 18af517187..e1ee2ecf64 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
@@ -68,7 +68,7 @@ object JulianDate {
   val SecondsPerHalfDay: Int = // 43200
     SecondsPerDay / 2
 
-  val Billion: Int = 1000000000
+  val Billion: Long = 1000000000
   val NanoPerDay: Long = SecondsPerDay.toLong * Billion.toLong
 
   val MinAdjustment: Long = -SecondsPerHalfDay.toLong * Billion.toLong
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
index 39e319ba89..0a93b5cb2f 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
@@ -33,7 +33,7 @@ case class ProperMotion(
       val redshift: Option[Redshift] = Target.redshift.get(t).flatten
       redshift.getOrElse(Redshift.zero).toRadialVelocity.toKilometersPerSecond
     }
-    val parallax = Parallax.mas.get(Target.parallax.get(t).flatten.getOrElse(Parallax.zero)) / 1000000.0
+    val parallax = Parallax.mas.get(Target.parallax.get(t).flatten.getOrElse(Parallax.zero)) / ProperMotion.Million
     ProperMotion.properMotionCalculator(baseCoordinates,
                                         Epoch.JulianLengthOfYear,
                                         properVelocity,
@@ -48,8 +48,9 @@ object ProperMotion {
   val AstronomicalUnit:    Long   = 149597870660L
   private val secsPerDay:  Double = 86400.0
   private val auPerKm:     Double = 1000.0 / AstronomicalUnit.toDouble
-  private val radsPerAsec: Double = Angle.fromArcsecs(1).toArcsecs
+  private val radsPerAsec: Double = Angle.fromArcsecs(1).toRadians
   private val TwoPi:       Double = 6.283185307179586476925286766559
+  private val Million:     Double = 1000000.0
 
   private type Vec2 = (Double, Double)
   private type Vec3 = (Double, Double, Double)

From 57f17f6838583b6fdbca739172e4d94361b58134 Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Mon, 25 Nov 2019 15:57:25 -0300
Subject: [PATCH 10/14] REL-3179: Working on testing.

---
 .../scala/edu/gemini/spModel/core/Epoch.scala | 70 +----------------
 .../edu/gemini/spModel/core/JulianDate.scala  | 70 +++++++++++++++++
 .../gemini/spModel/core/ProperMotion.scala    |  1 -
 .../edu/gemini/spModel/core/Arbitraries.scala | 76 +++++++++++++++++--
 .../edu/gemini/spModel/core/EpochSpec.scala   | 44 +++++++++++
 5 files changed, 189 insertions(+), 72 deletions(-)
 create mode 100644 bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
 create mode 100644 bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
index e1ee2ecf64..29a97c0d54 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/Epoch.scala
@@ -1,6 +1,5 @@
 package edu.gemini.spModel.core
 
-import java.time.ZoneOffset.UTC
 import java.time._
 
 import scalaz._
@@ -46,71 +45,10 @@ object Epoch {
     JulianDate.ofLocalDateTime(dt).dayNumber.toDouble
 
   val year: Epoch @> Double = Lens.lensu((a, b) => a.copy(year = b), _.year)
-}
-
-sealed abstract case class JulianDate(dayNumber: Int,
-                                      nanoAdjustment: Long) {
 
-  import JulianDate._
+  implicit val EpochOrder: Order[Epoch] =
+    Order.orderBy(_.year)
 
-  // Guaranteed by the JulianDate constructors, double checked here.
-  assert(dayNumber >= 0, s"dayNumber >= 0")
-  assert(nanoAdjustment >= MinAdjustment, s"nanoAdjustment >= $MinAdjustment")
-  assert(nanoAdjustment <= MaxAdjustment, s"nanoAdjustment <= $MaxAdjustment")
+  implicit val EpochShow: Show[Epoch] =
+    Show.showFromToString
 }
-
-object JulianDate {
-
-  /** JulianDate and related constants. */
-  val SecondsPerDay: Int = // 86400
-    24 * 60 * 60
-
-  val SecondsPerHalfDay: Int = // 43200
-    SecondsPerDay / 2
-
-  val Billion: Long = 1000000000
-  val NanoPerDay: Long = SecondsPerDay.toLong * Billion.toLong
-
-  val MinAdjustment: Long = -SecondsPerHalfDay.toLong * Billion.toLong
-  val MaxAdjustment: Long = SecondsPerHalfDay.toLong * Billion.toLong - 1
-
-  /** J2000 reference epoch as Julian Date. */
-  val J2000: JulianDate = // JulianDate(2451545,0)
-    JulianDate.ofLocalDateTime(
-      LocalDateTime.of(2000, 1, 1, 12, 0, 0)
-    )
-
-  /** Convert an `Instant` to a Julian Date.
-   */
-  def ofInstant(i: Instant): JulianDate =
-    ofLocalDateTime(LocalDateTime.ofInstant(i, UTC))
-
-  /** JulianDate from a `LocalDateTime` assumed to represent a time at UTC.
-   */
-  def ofLocalDateTime(ldt: LocalDateTime): JulianDate = {
-    val y = ldt.getYear
-    val m = ldt.getMonthValue
-    val d = ldt.getDayOfMonth
-
-    // Julian Day Number algorithm from:
-    // Fliegel, H.F. and Van Flandern, T.C. (1968). "A Machine Algorithm for
-    // Processing Calendar Dates" Communications of the Association of Computing
-    // Machines ll, 6sT.
-
-    // Yes, integer division.  -1 for Jan and Feb. 0 for Mar - Dec.
-    val t = (m - 14) / 12
-
-    // Julian Day Number (integer division).
-    val jdn = (1461 * (y + 4800 + t)) / 4 +
-      (367 * (m - 2 - 12 * t)) / 12 -
-      (3 * ((y + 4900 + t) / 100)) / 4 +
-      d - 32075
-
-    // Whole seconds since midnight
-    val secs = ldt.getHour * 3600 + ldt.getMinute * 60 + ldt.getSecond
-    val adj = (secs - SecondsPerHalfDay).toLong * Billion + ldt.getNano
-
-    new JulianDate(jdn, adj) {}
-  }
-}
-
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
new file mode 100644
index 0000000000..ba4ac58bc3
--- /dev/null
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
@@ -0,0 +1,70 @@
+package edu.gemini.spModel.core
+
+import java.time.{Instant, LocalDateTime}
+import java.time.ZoneOffset.UTC
+
+sealed abstract case class JulianDate(dayNumber: Int,
+                                      nanoAdjustment: Long) {
+
+  import JulianDate._
+
+  // Guaranteed by the JulianDate constructors, double checked here.
+  assert(dayNumber >= 0, s"dayNumber >= 0")
+  assert(nanoAdjustment >= MinAdjustment, s"nanoAdjustment >= $MinAdjustment")
+  assert(nanoAdjustment <= MaxAdjustment, s"nanoAdjustment <= $MaxAdjustment")
+}
+
+object JulianDate {
+
+  /** JulianDate and related constants. */
+  val SecondsPerDay: Int = // 86400
+    24 * 60 * 60
+
+  val SecondsPerHalfDay: Int = // 43200
+    SecondsPerDay / 2
+
+  val Billion: Long = 1000000000
+  val NanoPerDay: Long = SecondsPerDay.toLong * Billion.toLong
+
+  val MinAdjustment: Long = -SecondsPerHalfDay.toLong * Billion.toLong
+  val MaxAdjustment: Long = SecondsPerHalfDay.toLong * Billion.toLong - 1
+
+  /** J2000 reference epoch as Julian Date. */
+  val J2000: JulianDate = // JulianDate(2451545,0)
+    JulianDate.ofLocalDateTime(
+      LocalDateTime.of(2000, 1, 1, 12, 0, 0)
+    )
+
+  /** Convert an `Instant` to a Julian Date.
+   */
+  def ofInstant(i: Instant): JulianDate =
+    ofLocalDateTime(LocalDateTime.ofInstant(i, UTC))
+
+  /** JulianDate from a `LocalDateTime` assumed to represent a time at UTC.
+   */
+  def ofLocalDateTime(ldt: LocalDateTime): JulianDate = {
+    val y = ldt.getYear
+    val m = ldt.getMonthValue
+    val d = ldt.getDayOfMonth
+
+    // Julian Day Number algorithm from:
+    // Fliegel, H.F. and Van Flandern, T.C. (1968). "A Machine Algorithm for
+    // Processing Calendar Dates" Communications of the Association of Computing
+    // Machines ll, 6sT.
+
+    // Yes, integer division.  -1 for Jan and Feb. 0 for Mar - Dec.
+    val t = (m - 14) / 12
+
+    // Julian Day Number (integer division).
+    val jdn = (1461 * (y + 4800 + t)) / 4 +
+      (367 * (m - 2 - 12 * t)) / 12 -
+      (3 * ((y + 4900 + t) / 100)) / 4 +
+      d - 32075
+
+    // Whole seconds since midnight
+    val secs = ldt.getHour * 3600 + ldt.getMinute * 60 + ldt.getSecond
+    val adj = (secs - SecondsPerHalfDay).toLong * Billion + ldt.getNano
+
+    new JulianDate(jdn, adj) {}
+  }
+}
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
index 0a93b5cb2f..4a4268a986 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
@@ -56,7 +56,6 @@ object ProperMotion {
   private type Vec3 = (Double, Double, Double)
 
   // Scalar multiplication and addition for Vec3.
-  // TODO: Does scalaz must provide monoids for addition and multiplication?
   private implicit class Vec3Ops(a: Vec3) {
     def *(d: Double): Vec3 =
       (a._1 * d, a._2 * d, a._3 * d)
diff --git a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
index 11a5d079cb..50ef47247e 100644
--- a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
+++ b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
@@ -1,15 +1,16 @@
 package edu.gemini.spModel.core
 
+import java.time.{LocalDate, LocalDateTime, LocalTime, Year, ZoneId}
+
 import edu.gemini.spModel.core.WavelengthConversions._
 import org.scalacheck._
 import org.scalacheck.Gen._
 import org.scalacheck.Arbitrary._
-
 import squants.motion.VelocityConversions._
 import squants.radio.IrradianceConversions._
 import squants.radio.SpectralIrradianceConversions._
-
-import scalaz.{ ==>>, Order }
+import scala.collection.JavaConverters._
+import scalaz.{==>>, Order}
 import scalaz.std.anyVal._
 
 trait Arbitraries {
@@ -101,7 +102,7 @@ trait Arbitraries {
       } yield Magnitude(value, band, error, system)
     }
 
-  implicit val arbDistribution = Arbitrary[SpectralDistribution] {
+  implicit val arbDistribution: Arbitrary[SpectralDistribution] = Arbitrary[SpectralDistribution] {
     Gen.oneOf(
       BlackBody(8000),
       BlackBody(10000),
@@ -115,7 +116,7 @@ trait Arbitraries {
       LibraryNonStar.GammaDra
     )
   }
-  implicit val arbProfile = Arbitrary[SpatialProfile] {
+  implicit val arbProfile: Arbitrary[SpatialProfile] = Arbitrary[SpatialProfile] {
     Gen.oneOf(
       PointSource,
       UniformSource,
@@ -192,3 +193,68 @@ trait Arbitraries {
       case _       => sys.error("generated empty list!?")
     })
 }
+
+// Arbitrary but reasonable dates and times.
+trait ArbTime {
+
+  implicit val arbZoneId: Arbitrary[ZoneId] =
+    Arbitrary {
+      oneOf(ZoneId.getAvailableZoneIds.asScala.toSeq).map(ZoneId.of)
+    }
+
+  implicit val arbYear: Arbitrary[Year] =
+    Arbitrary {
+      choose(2000, 2020).map(Year.of)
+    }
+
+  implicit val arbLocalDate: Arbitrary[LocalDate] =
+    Arbitrary {
+      for {
+        y <- arbitrary[Year]
+        d <- choose(1, y.length)
+      } yield LocalDate.ofYearDay(y.getValue, d)
+    }
+
+  implicit val arbLocalTime: Arbitrary[LocalTime] =
+    Arbitrary {
+      for {
+        h <- choose(0, 23)
+        m <- choose(0, 59)
+        s <- choose(0, 59)
+        n <- choose(0, 999999999)
+      } yield LocalTime.of(h, m, s, n)
+    }
+
+  implicit val arbLocalDateTime: Arbitrary[LocalDateTime] =
+    Arbitrary {
+      for {
+        d <- arbitrary[LocalDate]
+        t <- arbitrary[LocalTime]
+      } yield LocalDateTime.of(d, t)
+    }
+
+  implicit val cogLocalDate: Cogen[LocalDate] =
+    Cogen[(Int, Int)].contramap(d => (d.getYear, d.getDayOfYear))
+}
+
+object ArbTime extends ArbTime
+
+trait ArbEpoch {
+  import ArbTime._
+
+  implicit val arbEpoch: Arbitrary[Epoch] =
+    Arbitrary {
+      for {
+        ldt <- arbitrary[LocalDateTime]
+      } yield {
+        val jd: Double = JulianDate.ofLocalDateTime(ldt).dayNumber.toDouble
+        val epochYears: Double = Epoch.JulianYearBasis + jd - Epoch.JulianBasis
+        new Epoch(epochYears)
+      }
+    }
+
+  implicit val cogEpoch: Cogen[Epoch] =
+    Cogen[String].contramap(Epoch.fromString.reverseGet)
+}
+
+object ArbEpoch extends ArbEpoch
\ No newline at end of file
diff --git a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala
new file mode 100644
index 0000000000..a2ab3182ca
--- /dev/null
+++ b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala
@@ -0,0 +1,44 @@
+package edu.gemini.spModel.core
+
+import java.time.LocalDateTime
+
+import scalaz._, Scalaz._
+
+final class EpochSpec extends CatsSuite {
+  import ArbEpoch._
+  import ArbTime._
+
+  test("Epoch.eq.natural") {
+    forAll { (a: Epoch, b: Epoch) =>
+      a.equals(b) shouldEqual Eq[Epoch].eqv(a, b)
+    }
+  }
+
+  test("Epoch.show.natural") {
+    forAll { (a: Epoch) =>
+      a.toString shouldEqual Show[Epoch].show(a)
+    }
+  }
+
+  test("Epoch.until.identity") {
+    forAll { (a: Epoch) =>
+      a.untilEpochYear(a.epochYear) shouldEqual 0.0
+    }
+  }
+
+  test("Epoch.until.sanity") {
+    forAll { (a: Epoch, s: Short) =>
+      a.untilEpochYear(a.epochYear + s.toDouble) shouldEqual s.toDouble
+    }
+  }
+
+  test("Epoch.until.sanity2") {
+    forAll { (s: Epoch.Scheme, d1: LocalDateTime, d2: LocalDateTime) =>
+      val Δ1 = s.fromLocalDateTime(d1).untilLocalDateTime(d2)
+      val Δ2 = s.fromLocalDateTime(d2).epochYear - s.fromLocalDateTime(d1).epochYear
+      Δ1 shouldEqual Δ2
+    }
+  }
+
+}
+

From 76fd304caffefe73acdb3abed588e3df66ea1a96 Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Mon, 25 Nov 2019 17:05:07 -0300
Subject: [PATCH 11/14] REL-3179: EpochSpec.

---
 .../edu/gemini/spModel/core/Arbitraries.scala | 20 --------
 .../edu/gemini/spModel/core/EpochSpec.scala   | 48 +++++--------------
 2 files changed, 13 insertions(+), 55 deletions(-)

diff --git a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
index 50ef47247e..78b274002c 100644
--- a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
+++ b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
@@ -238,23 +238,3 @@ trait ArbTime {
 }
 
 object ArbTime extends ArbTime
-
-trait ArbEpoch {
-  import ArbTime._
-
-  implicit val arbEpoch: Arbitrary[Epoch] =
-    Arbitrary {
-      for {
-        ldt <- arbitrary[LocalDateTime]
-      } yield {
-        val jd: Double = JulianDate.ofLocalDateTime(ldt).dayNumber.toDouble
-        val epochYears: Double = Epoch.JulianYearBasis + jd - Epoch.JulianBasis
-        new Epoch(epochYears)
-      }
-    }
-
-  implicit val cogEpoch: Cogen[Epoch] =
-    Cogen[String].contramap(Epoch.fromString.reverseGet)
-}
-
-object ArbEpoch extends ArbEpoch
\ No newline at end of file
diff --git a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala
index a2ab3182ca..ae5a21de88 100644
--- a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala
+++ b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala
@@ -1,44 +1,22 @@
 package edu.gemini.spModel.core
 
-import java.time.LocalDateTime
-
-import scalaz._, Scalaz._
-
-final class EpochSpec extends CatsSuite {
-  import ArbEpoch._
-  import ArbTime._
-
-  test("Epoch.eq.natural") {
-    forAll { (a: Epoch, b: Epoch) =>
-      a.equals(b) shouldEqual Eq[Epoch].eqv(a, b)
+import scalaz._
+import Scalaz._
+import org.scalacheck.Prop._
+import org.specs2.ScalaCheck
+import org.specs2.mutable.Specification
+
+final class EpochSpec extends Specification with ScalaCheck with Arbitraries {
+  "Epoch" should {
+    "give an offset of 0 to its own year" ! {
+      Epoch.J2000.untilEpochYear(Epoch.J2000.year) shouldEqual 0.0
     }
   }
 
-  test("Epoch.show.natural") {
-    forAll { (a: Epoch) =>
-      a.toString shouldEqual Show[Epoch].show(a)
+  "provide sanity in offsets using short" ! {
+    forAll { (s: Short) =>
+      Epoch.J2000.untilEpochYear(Epoch.J2000.year + s.toDouble) shouldEqual s.toDouble
     }
   }
-
-  test("Epoch.until.identity") {
-    forAll { (a: Epoch) =>
-      a.untilEpochYear(a.epochYear) shouldEqual 0.0
-    }
-  }
-
-  test("Epoch.until.sanity") {
-    forAll { (a: Epoch, s: Short) =>
-      a.untilEpochYear(a.epochYear + s.toDouble) shouldEqual s.toDouble
-    }
-  }
-
-  test("Epoch.until.sanity2") {
-    forAll { (s: Epoch.Scheme, d1: LocalDateTime, d2: LocalDateTime) =>
-      val Δ1 = s.fromLocalDateTime(d1).untilLocalDateTime(d2)
-      val Δ2 = s.fromLocalDateTime(d2).epochYear - s.fromLocalDateTime(d1).epochYear
-      Δ1 shouldEqual Δ2
-    }
-  }
-
 }
 

From 676d64373063ba893b596bb221b1c17a975f765e Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Mon, 25 Nov 2019 18:30:42 -0300
Subject: [PATCH 12/14] REL-3179: JulianDateSpec.

---
 .../edu/gemini/spModel/core/JulianDate.scala  | 29 +++++++++
 .../edu/gemini/spModel/core/Arbitraries.scala | 14 +++++
 .../edu/gemini/spModel/core/EpochSpec.scala   |  2 +-
 .../gemini/spModel/core/JulianDateSpec.scala  | 62 +++++++++++++++++++
 4 files changed, 106 insertions(+), 1 deletion(-)
 create mode 100644 bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/JulianDateSpec.scala

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
index ba4ac58bc3..f69f7129e0 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
@@ -3,6 +3,9 @@ package edu.gemini.spModel.core
 import java.time.{Instant, LocalDateTime}
 import java.time.ZoneOffset.UTC
 
+import scalaz._
+import Scalaz._
+
 sealed abstract case class JulianDate(dayNumber: Int,
                                       nanoAdjustment: Long) {
 
@@ -12,6 +15,30 @@ sealed abstract case class JulianDate(dayNumber: Int,
   assert(dayNumber >= 0, s"dayNumber >= 0")
   assert(nanoAdjustment >= MinAdjustment, s"nanoAdjustment >= $MinAdjustment")
   assert(nanoAdjustment <= MaxAdjustment, s"nanoAdjustment <= $MaxAdjustment")
+
+  /** Julian date value as a Double, including Julian Day Number and fractional
+   * day since the preceding noon.
+   */
+  val toDouble: Double =
+    dayNumber + nanoAdjustment.toDouble / NanoPerDay.toDouble
+
+  /** Modified Julian Date (MJD) double.  This is logically the same as
+   * `toDouble - 2400000.5`. MJD was introduced to preserve a bit of floating
+   * point decimal precision in calculations that use Julian dates.  It also
+   * makes it easier to directly implement algorithms that work with MJD.
+   *
+   * @see http://tycho.usno.navy.mil/mjd.html
+   */
+  def toModifiedDouble: Double = {
+    val h = SecondsPerHalfDay.toLong * Billion.toLong
+    val d = dayNumber      - 2400000
+    val n = nanoAdjustment - h
+
+    val (dʹ,nʹ) = if (n >= MinAdjustment) (d, n)
+    else (d - 1, n + SecondsPerDay.toLong * Billion.toLong)
+
+    dʹ + nʹ.toDouble / NanoPerDay.toDouble
+  }
 }
 
 object JulianDate {
@@ -67,4 +94,6 @@ object JulianDate {
 
     new JulianDate(jdn, adj) {}
   }
+
+  implicit val JulianDateEquals: Equal[JulianDate] = Equal.equalA
 }
diff --git a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
index 78b274002c..53c2488e94 100644
--- a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
+++ b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/Arbitraries.scala
@@ -238,3 +238,17 @@ trait ArbTime {
 }
 
 object ArbTime extends ArbTime
+
+trait ArbJulianDate {
+  import ArbTime._
+
+  implicit val arbJulianDate: Arbitrary[JulianDate] =
+    Arbitrary {
+      arbitrary[LocalDateTime].map(JulianDate.ofLocalDateTime)
+    }
+
+  implicit val cogJulianDate: Cogen[JulianDate] =
+    Cogen[(Int, Long)].contramap(jd => (jd.dayNumber, jd.nanoAdjustment))
+}
+
+object ArbJulianDate extends ArbJulianDate
\ No newline at end of file
diff --git a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala
index ae5a21de88..a86cb4a7d3 100644
--- a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala
+++ b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/EpochSpec.scala
@@ -9,7 +9,7 @@ import org.specs2.mutable.Specification
 final class EpochSpec extends Specification with ScalaCheck with Arbitraries {
   "Epoch" should {
     "give an offset of 0 to its own year" ! {
-      Epoch.J2000.untilEpochYear(Epoch.J2000.year) shouldEqual 0.0
+      math.abs(Epoch.J2000.untilEpochYear(Epoch.J2000.year)) < 1e-9
     }
   }
 
diff --git a/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/JulianDateSpec.scala b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/JulianDateSpec.scala
new file mode 100644
index 0000000000..71abefb772
--- /dev/null
+++ b/bundle/edu.gemini.spModel.core/src/test/scala/edu/gemini/spModel/core/JulianDateSpec.scala
@@ -0,0 +1,62 @@
+package edu.gemini.spModel.core
+
+import scalaz._
+import Scalaz._
+import org.scalacheck.Properties
+import org.scalacheck.Prop.forAll
+import org.scalacheck.Prop._
+import org.scalacheck.Gen
+import org.specs2.ScalaCheck
+import org.specs2.mutable.Specification
+import AlmostEqual.AlmostEqualOps
+
+import java.time.LocalDateTime
+
+final class JulianDateSpec extends Specification with ScalaCheck with Arbitraries with Helpers {
+
+  import ArbTime._
+  import ArbJulianDate._
+
+  // Old Epoch.Scheme.toJulianDay algorithm
+  private def oldEpochSchemeJulianDate(dt: LocalDateTime): Double = {
+    import scala.math.floor
+
+    val a = floor((14.0 - dt.getMonthValue) / 12.0)
+    val y = dt.getYear + 4800.0 - a
+    val m = dt.getMonthValue + 12 * a - 3.0
+    dt.getDayOfMonth +
+      floor((153.0 * m + 2.0) / 5.0) +
+      365 * y +
+      floor(y / 4.0) -
+      floor(y / 100.0) +
+      floor(y / 400.0) -
+      32045.0
+  }
+
+  "JulianDate" should {
+        "respect Eq" !
+          forAll { (a: JulianDate, b: JulianDate) =>
+            a.equals(b) shouldEqual Equal[JulianDate].equal(a, b)
+          }
+
+    "have dayNumber match old calculation" !
+      forAll { (ldt: LocalDateTime) =>
+        val jd0 = oldEpochSchemeJulianDate(ldt)
+        val jd1 = JulianDate.ofLocalDateTime(ldt).dayNumber.toDouble
+
+        jd0 shouldEqual jd1
+      }
+
+    "Some specific dates compared to USNO calculations" ! {
+      JulianDate.ofLocalDateTime(LocalDateTime.of(1918, 11, 11, 11, 0, 0)).toDouble ~= 2421908.958333
+      JulianDate.ofLocalDateTime(LocalDateTime.of(1969, 7, 21, 2, 56, 15)).toDouble ~=  2440423.622396
+      JulianDate.ofLocalDateTime(LocalDateTime.of(2001, 9, 11, 8, 46, 0)).toDouble ~= 2452163.865278
+      JulianDate.ofLocalDateTime(LocalDateTime.of(2345, 6, 7, 12, 0, 0)).toDouble ~= 2577711.000000
+    }
+
+    "Modified JulianDate should almost equal JulianDate - 2400000.5" !
+      forAll { (j: JulianDate) =>
+        j.toModifiedDouble ~= (j.toDouble - 2400000.5)
+      }
+  }
+}
\ No newline at end of file

From 90e8504db516f2dc9e2ea5b53acc4dbfdda1b490 Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Thu, 28 Nov 2019 11:13:53 -0300
Subject: [PATCH 13/14] REL-3179: Addressing PR comments.

---
 .../scala/edu/gemini/spModel/core/JulianDate.scala     | 10 +++++-----
 .../scala/edu/gemini/spModel/core/ProperMotion.scala   |  5 ++---
 2 files changed, 7 insertions(+), 8 deletions(-)

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
index f69f7129e0..f97c0a9da9 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/JulianDate.scala
@@ -30,12 +30,12 @@ sealed abstract case class JulianDate(dayNumber: Int,
    * @see http://tycho.usno.navy.mil/mjd.html
    */
   def toModifiedDouble: Double = {
-    val h = SecondsPerHalfDay.toLong * Billion.toLong
+    val h = SecondsPerHalfDay.toLong * Billion
     val d = dayNumber      - 2400000
     val n = nanoAdjustment - h
 
     val (dʹ,nʹ) = if (n >= MinAdjustment) (d, n)
-    else (d - 1, n + SecondsPerDay.toLong * Billion.toLong)
+    else (d - 1, n + SecondsPerDay.toLong * Billion)
 
     dʹ + nʹ.toDouble / NanoPerDay.toDouble
   }
@@ -51,10 +51,10 @@ object JulianDate {
     SecondsPerDay / 2
 
   val Billion: Long = 1000000000
-  val NanoPerDay: Long = SecondsPerDay.toLong * Billion.toLong
+  val NanoPerDay: Long = SecondsPerDay.toLong * Billion
 
-  val MinAdjustment: Long = -SecondsPerHalfDay.toLong * Billion.toLong
-  val MaxAdjustment: Long = SecondsPerHalfDay.toLong * Billion.toLong - 1
+  val MinAdjustment: Long = -SecondsPerHalfDay.toLong * Billion
+  val MaxAdjustment: Long = SecondsPerHalfDay.toLong * Billion - 1
 
   /** J2000 reference epoch as Julian Date. */
   val J2000: JulianDate = // JulianDate(2451545,0)
diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
index 4a4268a986..2a477b0246 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
@@ -25,7 +25,6 @@ case class ProperMotion(
   def plusYears(t: Target, elapsedYears: Double): Coordinates = {
     val baseCoordinates = t.coords(None).getOrElse(Coordinates.zero)
 
-    // TODO: Is this calculated correctly?
     val properVelocity: Offset = Offset(OffsetP(Angle.fromDegrees(deltaRA.velocity.toDegreesPerYear)),
                                         OffsetQ(Angle.fromDegrees(deltaDec.velocity.toDegreesPerYear)))
 
@@ -123,9 +122,9 @@ object ProperMotion {
       if (rem < 0.0) rem + TwoPi else rem
     }
 
-    // TODO: What to do if Declination is None?
     Coordinates(RA.fromAngle(Angle.fromRadians(rapp)),
-                Declination.fromAngle(Angle.fromRadians(decp)).getOrElse(Declination.zero))
+                Declination.fromAngle(Angle.fromRadians(decp))
+                  .getOrElse(sys.error(s"Invalid declination: $decp radians")))
   }
 
   val zero = ProperMotion(RightAscensionAngularVelocity.Zero, DeclinationAngularVelocity.Zero, Epoch.J2000)

From cb5db85f0d4fe99c691d1838d995e66a15c7a96d Mon Sep 17 00:00:00 2001
From: sraaphorst <sraaphorst@gemini.edu>
Date: Fri, 29 Nov 2019 12:07:32 -0300
Subject: [PATCH 14/14] REL-3179: Proper motion sanity test performed.

---
 .../src/main/scala/edu/gemini/spModel/core/ProperMotion.scala | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
index 2a477b0246..87a29de68b 100644
--- a/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
+++ b/bundle/edu.gemini.spModel.core/src/main/scala/edu/gemini/spModel/core/ProperMotion.scala
@@ -83,8 +83,8 @@ object ProperMotion {
                              parallax: Double,
                              elapsedYears: Double): Coordinates = {
     // We want the baseCoordinates in radians.
-    val (ra,  dec)  = (baseCoordinates.ra.toAngle.toRadians, baseCoordinates.dec.toAngle.toRadians)
-    val (dRa, dDec) = (properVelocity.p.toAngle.toRadians,   properVelocity.q.toAngle.toRadians)
+    val (ra,  dec)  = (baseCoordinates.ra.toAngle.toSignedDegrees.toRadians, baseCoordinates.dec.toAngle.toSignedDegrees.toRadians)
+    val (dRa, dDec) = (properVelocity.p.toAngle.toSignedDegrees.toRadians,   properVelocity.q.toAngle.toSignedDegrees.toRadians)
 
     val pos: Vec3 = {
       val cd = cos(dec)