From 1c268bb551b57464e6a76a3b82a54d964bd7260d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?R=C3=A9my=20Hannequin?= <remy.hannequin@gmail.com>
Date: Wed, 10 Apr 2024 11:55:53 +0200
Subject: [PATCH] Add interpolation method (#52)

Interpolation is very efficient and useful to compute intermediate
values for instants or quantities that are equidistant.

This adds an interpolation method for three of five terms, given an
interpolation factor.
---
 lib/astronoby.rb                  |  1 +
 lib/astronoby/util/maths.rb       | 72 +++++++++++++++++++++++++++++++
 spec/astronoby/util/maths_spec.rb | 64 +++++++++++++++++++++++++++
 3 files changed, 137 insertions(+)
 create mode 100644 lib/astronoby/util/maths.rb
 create mode 100644 spec/astronoby/util/maths_spec.rb

diff --git a/lib/astronoby.rb b/lib/astronoby.rb
index 58863bf..88f1383 100644
--- a/lib/astronoby.rb
+++ b/lib/astronoby.rb
@@ -21,6 +21,7 @@
 require "astronoby/time/greenwich_sidereal_time"
 require "astronoby/time/local_sidereal_time"
 require "astronoby/util/astrodynamics"
+require "astronoby/util/maths"
 require "astronoby/util/time"
 require "astronoby/util/trigonometry"
 require "astronoby/true_obliquity"
diff --git a/lib/astronoby/util/maths.rb b/lib/astronoby/util/maths.rb
new file mode 100644
index 0000000..31df238
--- /dev/null
+++ b/lib/astronoby/util/maths.rb
@@ -0,0 +1,72 @@
+# frozen_string_literal: true
+
+module Astronoby
+  module Util
+    module Maths
+      class << self
+        # Source:
+        #  Title: Astronomical Algorithms
+        #  Author: Jean Meeus
+        #  Edition: 2nd edition
+        #  Chapter: 3 - Interpolation
+
+        # @param values [Array<Numeric>] First term
+        # @param factor [Numeric] Interpolation factor
+        # @return [Float] Interpolated value
+        def interpolate(values, factor)
+          unless factor.between?(0, 1)
+            raise IncompatibleArgumentsError,
+              "Interpolation factor must be between 0 and 1"
+          end
+
+          if values.length == 3
+            return interpolate_3_terms(values, factor)
+          elsif values.length == 5
+            return interpolate_5_terms(values, factor)
+          end
+
+          raise IncompatibleArgumentsError,
+            "Only 3 or 5 terms are supported for interpolation"
+        end
+
+        private
+
+        # @return [Float] Interpolated value
+        def interpolate_3_terms(terms, factor)
+          y1, y2, y3 = terms
+
+          a = y2 - y1
+          b = y3 - y2
+          c = b - a
+
+          y2 + (factor / 2.0) * (a + b + factor * c)
+        end
+
+        # @return [Float] Interpolated value
+        def interpolate_5_terms(terms, factor)
+          y1, y2, y3, y4, y5 = terms
+
+          a = y2 - y1
+          b = y3 - y2
+          c = y4 - y3
+          d = y5 - y4
+
+          e = b - a
+          f = c - b
+          g = d - c
+
+          h = f - e
+          j = g - f
+
+          k = j - h
+
+          y3 +
+            factor * ((b + c) / 2.0 - (h + j) / 12.0) +
+            factor**2 * (f / 2.0 - k / 24.0) +
+            factor**3 * (h + j) / 12.0 +
+            factor**4 * k / 24.0
+        end
+      end
+    end
+  end
+end
diff --git a/spec/astronoby/util/maths_spec.rb b/spec/astronoby/util/maths_spec.rb
new file mode 100644
index 0000000..8d6b869
--- /dev/null
+++ b/spec/astronoby/util/maths_spec.rb
@@ -0,0 +1,64 @@
+# frozen_string_literal: true
+
+RSpec.describe Astronoby::Util::Maths do
+  describe "::interpolate" do
+    # Source:
+    #  Title: Astronomical Algorithms
+    #  Author: Jean Meeus
+    #  Edition: 2nd edition
+    #  Chapter: 3 - Interpolation, p.25
+    it "calculates the interpolation for 3 terms" do
+      interpolation = described_class.interpolate(
+        [0.884226, 0.877366, 0.870531],
+        0.18125
+      )
+
+      expect(interpolation.round(6)).to eq 0.876125
+    end
+
+    # Source:
+    #  Title: Astronomical Algorithms
+    #  Author: Jean Meeus
+    #  Edition: 2nd edition
+    #  Chapter: 3 - Interpolation, p.25
+    it "calculates the interpolation for 3 terms" do
+      interpolation = described_class.interpolate(
+        [1.3814294, 1.3812213, 1.3812453],
+        0.3966
+      )
+
+      expect(interpolation.round(6)).to eq 1.381203
+    end
+
+    # Source:
+    #  Title: Astronomical Algorithms
+    #  Author: Jean Meeus
+    #  Edition: 2nd edition
+    #  Chapter: 3 - Interpolation, p.25
+    it "calculates the interpolation for 5 terms" do
+      interpolation = described_class.interpolate(
+        [
+          Astronoby::Angle.from_hms(0, 54, 36.125).hours * 3600,
+          Astronoby::Angle.from_hms(0, 54, 24.606).hours * 3600,
+          Astronoby::Angle.from_hms(0, 54, 15.486).hours * 3600,
+          Astronoby::Angle.from_hms(0, 54, 8.694).hours * 3600,
+          Astronoby::Angle.from_hms(0, 54, 4.133).hours * 3600
+        ],
+        0.27777778
+      )
+
+      expect(interpolation.round(3))
+        .to eq(Astronoby::Angle.from_hms(0, 54, 13.369).hours * 3600)
+    end
+
+    context "when the interpolation factor is out of range" do
+      it "raises an error" do
+        expect { described_class.interpolate([1, 2, 3], 4) }
+          .to raise_error(
+            Astronoby::IncompatibleArgumentsError,
+            "Interpolation factor must be between 0 and 1"
+          )
+      end
+    end
+  end
+end