From d20291031a0cb06df88aeac0d953918e69f1d397 Mon Sep 17 00:00:00 2001
From: ericchristoffersen
 <46055145+ericchristoffersen@users.noreply.github.com>
Date: Sun, 14 Apr 2019 23:44:17 -0700
Subject: [PATCH] Fix divide by zero in slerp, also center spline length arc
 points (#3072)

---
 src/FileIO/LocationInterpolation.cpp |  8 ++++++--
 src/FileIO/LocationInterpolation.h   | 18 +++++++++++-------
 2 files changed, 17 insertions(+), 9 deletions(-)

diff --git a/src/FileIO/LocationInterpolation.cpp b/src/FileIO/LocationInterpolation.cpp
index 6ed39de34..e2982cc07 100644
--- a/src/FileIO/LocationInterpolation.cpp
+++ b/src/FileIO/LocationInterpolation.cpp
@@ -128,8 +128,12 @@ geolocation xyz::togeolocation() const
 
 xyz Slerper::Slerp(double frac)
 {
-    double scale = sin(frac * m_angle) / m_sin_angle;
-    return m_x0_norm.scale(sin((1 - frac) * m_angle) / m_sin_angle).add(m_x1_norm.scale(scale));
+    if (m_sin_angle != 0.0)
+    {
+        double scale = sin(frac * m_angle) / m_sin_angle;
+        return m_x0_norm.scale(sin((1 - frac) * m_angle) / m_sin_angle).add(m_x1_norm.scale(scale));
+    }
+    return m_x0_norm;
 }
 
 // Precompute invariant values needed to geoslerp
diff --git a/src/FileIO/LocationInterpolation.h b/src/FileIO/LocationInterpolation.h
index 80cad0fa6..ed247ecb0 100644
--- a/src/FileIO/LocationInterpolation.h
+++ b/src/FileIO/LocationInterpolation.h
@@ -633,9 +633,9 @@ public:
             d1 = workitem.d1;
 
             // Step 1
-            const double quarterspan = (d1 - d0) / 4;
-            const double inter0 = d0 + quarterspan;
-            const double inter1 = d0 + quarterspan + quarterspan;
+            const double thirdspan = (d1 - d0) / 3;
+            const double inter0 = d0 + thirdspan;
+            const double inter1 = d0 + thirdspan + thirdspan;
 
             // Step 2
             const xyz pm1 = this->Interpolate(d0);
@@ -645,18 +645,22 @@ public:
 
             // Step 3
             double linearDistance = p2.DistanceFrom(pm1);
-            if (linearDistance == 0.0)
+            if (linearDistance < 0.000001)
                 break;
 
             // Step 4
-            double quadDistance = p2.DistanceFrom(p1) + p1.DistanceFrom(p0) + p0.DistanceFrom(pm1);
+            double span0 = p0.DistanceFrom(pm1);
+            double span1 = p1.DistanceFrom(p0);
+            double span2 = p2.DistanceFrom(p1);
+
+            double arcDistance = span0 + span1 + span2;
 
             // Step 5
-            double difference = fabs(quadDistance / linearDistance) - 1.0;
+            double difference = fabs(arcDistance / linearDistance) - 1.0;
 
             // 5A: Settle for quaddistance if threshold met or no more room on worklist.
             if (difference < thresholdLimit || worklist.EmptySlots() < 3) {
-                finalLength += quadDistance;
+                finalLength += arcDistance;
             } else {
                 // 5B: otherwise push the 3 new subsegments onto worklist.
                 worklist.Push(CalcSplineLengthBracketPair(d0,     inter0));