diff --git a/lib/astronoby.rb b/lib/astronoby.rb index 9dfa7c8..c514a48 100644 --- a/lib/astronoby.rb +++ b/lib/astronoby.rb @@ -5,6 +5,8 @@ require "astronoby/angles/dms" require "astronoby/angles/hms" require "astronoby/epoch" +require "astronoby/astronomical_models/ephemeride_lunaire_parisienne" +require "astronoby/bodies/moon" require "astronoby/bodies/sun" require "astronoby/coordinates/ecliptic" require "astronoby/coordinates/equatorial" diff --git a/lib/astronoby/astronomical_models/ephemeride_lunaire_parisienne.rb b/lib/astronoby/astronomical_models/ephemeride_lunaire_parisienne.rb new file mode 100644 index 0000000..4efafd3 --- /dev/null +++ b/lib/astronoby/astronomical_models/ephemeride_lunaire_parisienne.rb @@ -0,0 +1,143 @@ +# frozen_string_literal: true + +module Astronoby + class EphemerideLunaireParisienne + # Lunar theory: Éphéméride Lunaire Parisienne + # By Jean Chapront, Michelle Chapront-Touzé + # https://en.wikipedia.org/wiki/Ephemeride_Lunaire_Parisienne + + DEGREES_UNIT = 10**-6 + + # @return [Array] Periodic terms for the Moon's longitude and distance + def self.periodic_terms_for_moon_longitude_and_distance + [ + [0, 0, 1, 0, 6288774, -20905355], + [2, 0, -1, 0, 1274027, -3699111], + [2, 0, 0, 0, 658314, -2955968], + [0, 0, 2, 0, 213618, -569925], + [0, 1, 0, 0, -185116, 48888], + [0, 0, 0, 2, -114332, -3149], + [2, 0, -2, 0, 58793, 246158], + [2, -1, -1, 0, 57066, -152138], + [2, 0, 1, 0, 53322, -170733], + [2, -1, 0, 0, 45758, -204586], + [0, 1, -1, 0, -40923, -129620], + [1, 0, 0, 0, -34720, 108743], + [0, 1, 1, 0, -30383, 104755], + [2, 0, 0, -2, 15327, 10321], + [0, 0, 1, 2, -12528, 0], + [0, 0, 1, -2, 10980, 79661], + [4, 0, -1, 0, 10675, -34782], + [0, 0, 3, 0, 10034, -23210], + [4, 0, -2, 0, 8548, -21636], + [2, 1, -1, 0, -7888, 24208], + [2, 1, 0, 0, -6766, 30824], + [1, 0, -1, 0, -5163, -8379], + [1, 1, 0, 0, 4987, -16675], + [2, -1, 1, 0, 4036, -12831], + [2, 0, 2, 0, 3994, -10445], + [4, 0, 0, 0, 3861, -11650], + [2, 0, -3, 0, 3665, 14403], + [0, 1, -2, 0, -2689, -7003], + [2, 0, -1, 2, -2602, 0], + [2, -1, -2, 0, 2390, 10056], + [1, 0, 1, 0, -2348, 6322], + [2, -2, 0, 0, 2236, -9884], + [0, 1, 2, 0, -2120, 5751], + [0, 2, 0, 0, -2069, 0], + [2, -2, -1, 0, 2048, -4950], + [2, 0, 1, -2, -1773, 4130], + [2, 0, 0, 2, -1595, 0], + [4, -1, -1, 0, 1215, -3958], + [0, 0, 2, 2, -1110, 0], + [3, 0, -1, 0, -892, 3258], + [2, 1, 1, 0, -810, 2616], + [4, -1, -2, 0, 759, -1897], + [0, 2, -1, 0, -713, -2117], + [2, 2, -1, 0, -700, 2354], + [2, 1, -2, 0, 691, 0], + [2, -1, 0, -2, 596, 0], + [4, 0, 1, 0, 549, -1423], + [0, 0, 4, 0, 537, -1117], + [4, -1, 0, 0, 520, -1571], + [1, 0, -2, 0, -487, -1739], + [2, 1, 0, -2, -399, 0], + [0, 0, 2, -2, -381, -4421], + [1, 1, 1, 0, 351, 0], + [3, 0, -2, 0, -340, 0], + [4, 0, -3, 0, 330, 0], + [2, -1, 2, 0, 327, 0], + [0, 2, 1, 0, -323, 1165], + [1, 1, -1, 0, 299, 0], + [2, 0, 3, 0, 294, 0], + [2, 0, -1, -2, 0, 8752] + ] + end + + # @return [Array] Periodic terms for the Moon's latitude + def self.periodic_terms_for_moon_latitude + [ + [0, 0, 0, 1, 5128122], + [0, 0, 1, 1, 280602], + [0, 0, 1, -1, 277693], + [2, 0, 0, -1, 173237], + [2, 0, -1, 1, 55413], + [2, 0, -1, -1, 46271], + [2, 0, 0, 1, 32573], + [0, 0, 2, 1, 17198], + [2, 0, 1, -1, 9266], + [0, 0, 2, -1, 8822], + [2, -1, 0, -1, 8216], + [2, 0, -2, -1, 4324], + [2, 0, 1, 1, 4200], + [2, 1, 0, -1, -3359], + [2, -1, -1, 1, 2463], + [2, -1, 0, 1, 2211], + [2, -1, -1, -1, 2065], + [0, 1, -1, -1, -1870], + [4, 0, -1, -1, 1828], + [0, 1, 0, 1, -1794], + [0, 0, 0, 3, -1749], + [0, 1, -1, 1, -1565], + [1, 0, 0, 1, -1491], + [0, 1, 1, 1, -1475], + [0, 1, 1, -1, -1410], + [0, 1, 0, -1, -1344], + [1, 0, 0, -1, -1335], + [0, 0, 3, 1, 1107], + [4, 0, 0, -1, 1021], + [4, 0, -1, 1, 833], + [0, 0, 1, -3, 777], + [4, 0, -2, 1, 671], + [2, 0, 0, -3, 607], + [2, 0, 2, -1, 596], + [2, -1, 1, -1, 491], + [2, 0, -2, 1, -451], + [0, 0, 3, -1, 439], + [2, 0, 2, 1, 422], + [2, 0, -3, -1, 421], + [2, 1, -1, 1, -366], + [2, 1, 0, 1, -351], + [4, 0, 0, 1, 331], + [2, -1, 1, 1, 315], + [2, -2, 0, -1, 302], + [0, 0, 1, 3, -283], + [2, 1, 1, -1, -229], + [1, 1, 0, -1, 223], + [1, 1, 0, 1, 223], + [0, 1, -2, -1, -220], + [2, 1, -1, -1, -220], + [1, 0, 1, 1, -185], + [2, -1, -2, -1, 181], + [0, 1, 2, 1, -177], + [4, 0, -2, -1, 176], + [4, -1, -1, -1, 166], + [1, 0, 1, -1, -164], + [4, 0, 1, -1, 132], + [1, 0, -1, -1, -119], + [4, -1, 0, -1, 115], + [2, -2, 0, 1, 107] + ] + end + end +end diff --git a/lib/astronoby/bodies/moon.rb b/lib/astronoby/bodies/moon.rb new file mode 100644 index 0000000..90684c0 --- /dev/null +++ b/lib/astronoby/bodies/moon.rb @@ -0,0 +1,226 @@ +# frozen_string_literal: true + +module Astronoby + class Moon + MEAN_GEOCENTRIC_DISTANCE = 385_000_560 + + def initialize(time:) + @time = time + end + + # Source: + # Title: Astronomical Algorithms + # Author: Jean Meeus + # Edition: 2nd edition + # Chapter: 47 - Position of the Moon + + # @return [Coordinates::Ecliptic] Ecliptic coordinates of the Moon + def apparent_ecliptic_coordinates + @ecliptic_coordinates ||= begin + latitude = Astronoby::Angle.from_degrees( + (latitude_terms + additive_latitude_terms) * + EphemerideLunaireParisienne::DEGREES_UNIT + ) + + longitude = mean_longitude + Astronoby::Angle.from_degrees( + (longitude_terms + additive_longitude_terms) * + EphemerideLunaireParisienne::DEGREES_UNIT + ) + nutation + + Coordinates::Ecliptic.new( + latitude: latitude, + longitude: longitude + ) + end + end + + # @return [Integer] Distance between the Earth and the Moon centers in meters + def distance + @distance ||= (MEAN_GEOCENTRIC_DISTANCE + distance_terms).round + end + + # @return [Angle] Moon's mean longitude + def mean_longitude + @mean_longitude ||= begin + degrees = + 218.3164477 + + 481267.88123421 * elapsed_centuries - + 0.0015786 * elapsed_centuries**2 + + elapsed_centuries**3 / 538841 - + elapsed_centuries**4 / 65194000 + degrees += Constants::DEGREES_PER_CIRCLE while degrees.negative? + Angle.from_degrees(degrees) + end + end + + # @return [Angle] Moon's mean elongation + def mean_elongation + @mean_elongation ||= begin + degrees = + 297.8501921 + + 445267.1114034 * elapsed_centuries - + 0.0018819 * elapsed_centuries**2 + + elapsed_centuries**3 / 545868 - + elapsed_centuries**4 / 113065000 + degrees += Constants::DEGREES_PER_CIRCLE while degrees.negative? + Angle.from_degrees(degrees) + end + end + + # @return [Angle] Moon's mean anomaly + def mean_anomaly + @mean_anomaly ||= begin + degrees = + 134.9633964 + + 477198.8675055 * elapsed_centuries + + 0.0087414 * elapsed_centuries**2 + + elapsed_centuries**3 / 69699 - + elapsed_centuries**4 / 14712000 + degrees += Constants::DEGREES_PER_CIRCLE while degrees.negative? + Angle.from_degrees(degrees) + end + end + + # @return [Angle] Moon's argument of latitude + def argument_of_latitude + @argument_of_latitude ||= begin + degrees = + 93.2720950 + + 483202.0175233 * elapsed_centuries - + 0.0036539 * elapsed_centuries**2 - + elapsed_centuries**3 / 3526000 + degrees += Constants::DEGREES_PER_CIRCLE while degrees.negative? + Angle.from_degrees(degrees) + end + end + + private + + def terrestrial_time + @terrestrial_time ||= begin + delta = Util::Time.terrestrial_universal_time_delta(@time) + @time + delta + end + end + + def julian_date + Epoch.from_time(terrestrial_time) + end + + def elapsed_centuries + (julian_date - Epoch::DEFAULT_EPOCH) / Constants::DAYS_PER_JULIAN_CENTURY + end + + def sun_mean_anomaly + @sun_mean_anomaly ||= Sun.new(time: @time).mean_anomaly + end + + def a1 + @a1 ||= begin + degrees = 119.75 + 131.849 * elapsed_centuries + degrees += Constants::DEGREES_PER_CIRCLE while degrees.negative? + Angle.from_degrees(degrees) + end + end + + def a2 + @a2 ||= begin + degrees = 53.09 + 479264.290 * elapsed_centuries + degrees += Constants::DEGREES_PER_CIRCLE while degrees.negative? + Angle.from_degrees(degrees) + end + end + + def a3 + @a3 ||= begin + degrees = 313.45 + 481266.484 * elapsed_centuries + degrees += Constants::DEGREES_PER_CIRCLE while degrees.negative? + Angle.from_degrees(degrees) + end + end + + def eccentricity_correction + @eccentricity_correction ||= + 1 - 0.002516 * elapsed_centuries - 0.0000074 * elapsed_centuries**2 + end + + def latitude_terms + @latitude_terms ||= + Astronoby::EphemerideLunaireParisienne + .periodic_terms_for_moon_latitude + .inject(0) do |sum, terms| + value = terms[4] * Math.sin( + terms[0] * mean_elongation.radians + + terms[1] * sun_mean_anomaly.radians + + terms[2] * mean_anomaly.radians + + terms[3] * argument_of_latitude.radians + ) + + value *= eccentricity_correction if terms[1].abs == 1 + value *= eccentricity_correction**2 if terms[1].abs == 2 + + sum + value + end + end + + def additive_latitude_terms + @additive_latitude_terms ||= + -2235 * mean_longitude.sin + + 382 * a3.sin + + 175 * (a1 - argument_of_latitude).sin + + 175 * (a1 + argument_of_latitude).sin + + 127 * (mean_longitude - mean_anomaly).sin - + 115 * (mean_longitude + mean_anomaly).sin + end + + def longitude_terms + @longitude_terms ||= + Astronoby::EphemerideLunaireParisienne + .periodic_terms_for_moon_longitude_and_distance + .inject(0) do |sum, terms| + value = terms[4] * Math.sin( + terms[0] * mean_elongation.radians + + terms[1] * sun_mean_anomaly.radians + + terms[2] * mean_anomaly.radians + + terms[3] * argument_of_latitude.radians + ) + + value *= eccentricity_correction if terms[1].abs == 1 + value *= eccentricity_correction**2 if terms[1].abs == 2 + + sum + value + end + end + + def additive_longitude_terms + @additive_longitude_terms ||= + 3958 * a1.sin + + 1962 * (mean_longitude - argument_of_latitude).sin + + 318 * a2.sin + end + + def distance_terms + @distance_terms ||= + EphemerideLunaireParisienne + .periodic_terms_for_moon_longitude_and_distance + .inject(0) do |sum, terms| + value = terms[5] * Math.cos( + terms[0] * mean_elongation.radians + + terms[1] * sun_mean_anomaly.radians + + terms[2] * mean_anomaly.radians + + terms[3] * argument_of_latitude.radians + ) + + value *= eccentricity_correction if terms[1].abs == 1 + value *= eccentricity_correction**2 if terms[1].abs == 2 + + sum + value + end + end + + def nutation + @nutation ||= + Astronoby::Nutation.for_ecliptic_longitude(epoch: julian_date) + end + end +end diff --git a/spec/astronoby/bodies/moon_spec.rb b/spec/astronoby/bodies/moon_spec.rb new file mode 100644 index 0000000..af32489 --- /dev/null +++ b/spec/astronoby/bodies/moon_spec.rb @@ -0,0 +1,123 @@ +# frozen_string_literal: true + +RSpec.describe Astronoby::Moon do + describe "#apparent_ecliptic_coordinates" do + # Source: + # Title: Astronomical Algorithms + # Author: Jean Meeus + # Edition: 2nd edition + # Chapter: 47 - Position of the Moon, p.342 + it "returns the apparent ecliptic coordinates for 1992-04-12" do + # Example gives 1992-04-12 00:00:00 TT (and not UTC) + time = Time.utc(1992, 4, 12, 0, 0, 0) + time -= Astronoby::Util::Time.terrestrial_universal_time_delta(time) + moon = described_class.new(time: time) + + coordinates = moon.apparent_ecliptic_coordinates + + expect(coordinates.latitude.str(:dms)).to eq "-3° 13′ 44.8543″" + # Result from the book: -3° 13′ 45″ + # Result from IMCCE: -3° 13′ 44.184″ + + expect(coordinates.longitude.str(:dms)).to eq "+133° 10′ 1.834″" + # Result from the book: 133° 10′ 02″ + # Result from IMCCE: 133° 10′ 0.157″ + end + + # Source: + # Title: Celestial Calculations + # Author: J. L. Lawrence + # Edition: MIT Press + # Chapter: 7 - The Moon, p.166 + it "returns the apparent ecliptic coordinates for 2015-01-01" do + moon = described_class.new(time: Time.new(2015, 1, 1, 22, 0, 0, "-05:00")) + + coordinates = moon.apparent_ecliptic_coordinates + + expect(coordinates.latitude.str(:dms)).to eq "-3° 50′ 46.0007″" + # Result from the book: -3° 57′ 22.5179″ (-3.956255°) + # Result from IMCCE: -3° 50′ 55.417″ + + expect(coordinates.longitude.str(:dms)).to eq "+65° 21′ 5.9223″" + # Result from the book: +65° 3′ 35.3304″ (65.059814°) + # Result from IMCCE: +65° 21′ 3.629″ + end + + # Source: + # Title: Celestial Calculations + # Author: J. L. Lawrence + # Edition: MIT Press + # Chapter: 7 - The Moon, p.185 + it "returns the apparent ecliptic coordinates for 2000-08-09" do + moon = described_class.new(time: Time.new(2000, 8, 9, 12, 0, 0, "-05:00")) + + coordinates = moon.apparent_ecliptic_coordinates + + expect(coordinates.latitude.str(:dms)).to eq "+3° 11′ 22.0683″" + # Result from the book: +3° 2′ 40.2″ (3.0445°) + # Result from IMCCE: +3° 11′ 25.819″ + + expect(coordinates.longitude.str(:dms)).to eq "+257° 17′ 32.7597″" + # Result from the book: +257° 13′ 11.784″ (257.21994°) + # Result from IMCCE: +257° 17′ 32.387″ + end + + # Source: + # Title: Celestial Calculations + # Author: J. L. Lawrence + # Edition: MIT Press + # Chapter: 7 - The Moon, p.185 + it "returns the apparent ecliptic coordinates for 2010-05-15" do + moon = described_class.new(time: Time.new(2010, 5, 15, 14, 30, 0, "-04:00")) + + coordinates = moon.apparent_ecliptic_coordinates + + expect(coordinates.latitude.str(:dms)).to eq "+2° 16′ 12.4064″" + # Result from the book: +2° 25′ 1.7975″ (2.417166°) + # Result from IMCCE: +2° 16′ 15.434″ + + expect(coordinates.longitude.str(:dms)).to eq "+76° 35′ 31.0243″" + # Result from the book: +76° 24′ 58.8924″ (76.416359°) + # Result from IMCCE: +76° 35′ 32.198″ + end + + # Source: + # Title: Practical Astronomy with your Calculator or Spreadsheet + # Authors: Peter Duffett-Smith and Jonathan Zwart + # Edition: Cambridge University Press + # Chapter: 65 - Calculating the Moon’s position, p.165 + it "returns the apparent ecliptic coordinates for 2003-09-01" do + moon = described_class.new(time: Time.utc(2003, 9, 1)) + + coordinates = moon.apparent_ecliptic_coordinates + + expect(coordinates.latitude.str(:dms)).to eq "+1° 37′ 12.7837″" + # Result from the book: +1° 42′ 57.8664″ (1.716074°) + # Result from IMCCE: +1° 37′ 9.680″ + + expect(coordinates.longitude.str(:dms)).to eq "+214° 46′ 16.0171″" + # Result from the book: +214° 52′ 3.0107″ (214.867503°) + # Result from IMCCE: +214° 46′ 16.888″ + end + end + + describe "#distance" do + # Source: + # Title: Astronomical Algorithms + # Author: Jean Meeus + # Edition: 2nd edition + # Chapter: 47 - Position of the Moon, p.342 + it "returns the distance for 1992-04-12" do + # Example gives 1992-04-12 00:00:00 TT (and not UTC) + time = Time.utc(1992, 4, 12, 0, 0, 0) + time -= Astronoby::Util::Time.terrestrial_universal_time_delta(time) + moon = described_class.new(time: time) + + distance = moon.distance + + expect(distance).to eq 368409707 + # Result from the book: 368409700 (36849.7 km) + # Result from IMCCE: 368439405 (0.002462865305 AU) + end + end +end