Handle Sun's sunrise and sunset azimuth

UPGRADING.md

@@ -58,12 +58,20 @@ Returns the apparent Sun's longitude (`Angle`) at its perigee.
 
 ### `Sun#rising_time` method added (#35)
 
-Returns the UTC `Time` of the sunrise.
+Returns the UTC `Time` of the sunrise.`
+
+### `Sun#rising_azimuth` method added (#39)
+
+Returns the Sun's azimuth (`Angle`) at sunrise.
 
 ### `Sun#setting_time` method added (#35)
 
 Returns the UTC `Time` of the sunset.
 
+### `Sun#setting_azimuth` method added (#39)
+
+Returns the Sun's azimuth (`Angle`) at sunset.
+
 ### Added comparison methods to `Angle` (#21)
 
 With the inclusion of `Comparable`, comparison methods such as `#==`, `#<`,

lib/astronoby/bodies/sun.rb

@@ -56,6 +56,16 @@ def rising_time(observer:)
       ).to_gst.to_utc
     end
 
+    # @param observer [Astronoby::Observer] Observer of the event
+    # @return [Astronoby::Angle, nil] Azimuth of sunrise
+    def rising_azimuth(observer:)
+      equatorial_coordinates = ecliptic_coordinates.to_equatorial(epoch: @epoch)
+      Body.new(equatorial_coordinates).rising_azimuth(
+        latitude: observer.latitude,
+        vertical_shift: vertical_shift
+      )
+    end
+
     # @param observer [Astronoby::Observer] Observer of the event
     # @return [Time] Time of sunset
     def setting_time(observer:)
@@ -72,6 +82,16 @@ def setting_time(observer:)
       ).to_gst.to_utc
     end
 
+    # @param observer [Astronoby::Observer] Observer of the event
+    # @return [Astronoby::Angle, nil] Azimuth of sunset
+    def setting_azimuth(observer:)
+      equatorial_coordinates = ecliptic_coordinates.to_equatorial(epoch: @epoch)
+      Body.new(equatorial_coordinates).setting_azimuth(
+        latitude: observer.latitude,
+        vertical_shift: vertical_shift
+      )
+    end
+
     # @return [Numeric] Earth-Sun distance in meters
     def earth_distance
       SEMI_MAJOR_AXIS_IN_METERS / distance_angular_size_factor
@@ -172,5 +192,11 @@ def event_local_sidereal_time_for_date(date, observer, event)
         .to_lst(longitude: observer.longitude)
         .time
     end
+
+    def vertical_shift
+      Astronoby::Body::DEFAULT_REFRACTION_VERTICAL_SHIFT +
+        Astronoby::GeocentricParallax.angle(distance: earth_distance) +
+        Astronoby::Angle.as_degrees(angular_size.degrees / 2)
+    end
   end
 end

lib/astronoby/body.rb

@@ -14,78 +14,106 @@ def initialize(equatorial_coordinates)
     #  Authors: Peter Duffett-Smith and Jonathan Zwart
     #  Edition: Cambridge University Press
     #  Chapter: 33 - Rising and setting
+
+    # @param latitude [Astronoby::Angle] Latitude of the observer
+    # @param longitude [Astronoby::Angle] Longitude of the observer
+    # @param date [Date] Date of the event
+    # @param apparent [Boolean] Compute apparent or true data
+    # @param vertical_shift [Astronoby::Angle] Vertical shift correction angle
+    # @return [Time, nil] Sunrise time
     def rising_time(
       latitude:,
       longitude:,
       date:,
       apparent: true,
       vertical_shift: nil
     )
-      ratio = ratio(latitude, apparent, vertical_shift)
-      return nil unless RISING_SETTING_HOUR_ANGLE_RATIO_RANGE.cover?(ratio)
+      time_ratio = time_ratio(latitude, apparent, vertical_shift)
+      return nil unless RISING_SETTING_HOUR_ANGLE_RATIO_RANGE.cover?(time_ratio)
 
-      hour_angle = Angle.acos(ratio)
+      hour_angle = Angle.acos(time_ratio)
       local_sidereal_time = LocalSiderealTime.new(
         date: date,
-        time: @equatorial_coordinates.right_ascension.hours - hour_angle.hours,
+        time: right_ascension.hours - hour_angle.hours,
         longitude: longitude
       )
 
       local_sidereal_time.to_gst.to_utc
     end
 
     # Source:
-    #  Title: Celestial Calculations
-    #  Author: J. L. Lawrence
-    #  Edition: MIT Press
-    #  Chapter: 5 - Stars in the Nighttime Sky
-    def rising_azimuth(latitude:)
-      ar = azimuth_component(latitude)
-      return nil if ar >= 1
-
-      Angle.acos(ar)
+    #  Title: Practical Astronomy with your Calculator or Spreadsheet
+    #  Authors: Peter Duffett-Smith and Jonathan Zwart
+    #  Edition: Cambridge University Press
+    #  Chapter: 33 - Rising and setting
+
+    # @param latitude [Astronoby::Angle] Latitude of the observer
+    # @param apparent [Boolean] Compute apparent or true data
+    # @param vertical_shift [Astronoby::Angle] Vertical shift correction angle
+    # @return [Astronoby::Angle, nil] Sunrise azimuth
+    def rising_azimuth(latitude:, apparent: true, vertical_shift: nil)
+      time_ratio = time_ratio(latitude, apparent, vertical_shift)
+      return nil unless RISING_SETTING_HOUR_ANGLE_RATIO_RANGE.cover?(time_ratio)
+
+      azimuth_ratio = azimuth_ratio(latitude, apparent, vertical_shift)
+
+      Angle.acos(azimuth_ratio)
     end
 
     # Source:
     #  Title: Practical Astronomy with your Calculator or Spreadsheet
     #  Authors: Peter Duffett-Smith and Jonathan Zwart
     #  Edition: Cambridge University Press
     #  Chapter: 33 - Rising and setting
+
+    # @param latitude [Astronoby::Angle] Latitude of the observer
+    # @param longitude [Astronoby::Angle] Longitude of the observer
+    # @param date [Date] Date of the event
+    # @param apparent [Boolean] Compute apparent or true data
+    # @param vertical_shift [Astronoby::Angle] Vertical shift correction angle
+    # @return [Time, nil] Sunset time
     def setting_time(
       latitude:,
       longitude:,
       date:,
       apparent: true,
       vertical_shift: nil
     )
-      ratio = ratio(latitude, apparent, vertical_shift)
-      return nil unless RISING_SETTING_HOUR_ANGLE_RATIO_RANGE.cover?(ratio)
+      time_ratio = time_ratio(latitude, apparent, vertical_shift)
+      return nil unless RISING_SETTING_HOUR_ANGLE_RATIO_RANGE.cover?(time_ratio)
 
-      hour_angle = Angle.acos(ratio)
+      hour_angle = Angle.acos(time_ratio)
       local_sidereal_time = LocalSiderealTime.new(
         date: date,
-        time: @equatorial_coordinates.right_ascension.hours + hour_angle.hours,
+        time: right_ascension.hours + hour_angle.hours,
         longitude: longitude
       )
 
       local_sidereal_time.to_gst.to_utc
     end
 
     # Source:
-    #  Title: Celestial Calculations
-    #  Author: J. L. Lawrence
-    #  Edition: MIT Press
-    #  Chapter: 5 - Stars in the Nighttime Sky
-    def setting_azimuth(latitude:)
-      rising_az = rising_azimuth(latitude: latitude)
-      return nil if rising_az.nil?
-
-      Angle.as_degrees(360 - rising_az.degrees)
+    #  Title: Practical Astronomy with your Calculator or Spreadsheet
+    #  Authors: Peter Duffett-Smith and Jonathan Zwart
+    #  Edition: Cambridge University Press
+    #  Chapter: 33 - Rising and setting
+
+    # @param latitude [Astronoby::Angle] Latitude of the observer
+    # @param apparent [Boolean] Compute apparent or true data
+    # @param vertical_shift [Astronoby::Angle] Vertical shift correction angle
+    # @return [Astronoby::Angle, nil] Sunset azimuth
+    def setting_azimuth(latitude:, apparent: true, vertical_shift: nil)
+      time_ratio = time_ratio(latitude, apparent, vertical_shift)
+      return nil unless RISING_SETTING_HOUR_ANGLE_RATIO_RANGE.cover?(time_ratio)
+
+      azimuth_ratio = azimuth_ratio(latitude, apparent, vertical_shift)
+
+      Angle.as_degrees(360 - Angle.acos(azimuth_ratio).degrees)
     end
 
     private
 
-    def ratio(latitude, apparent, vertical_shift)
+    def time_ratio(latitude, apparent, vertical_shift)
       shift = if vertical_shift
         vertical_shift
       elsif apparent
@@ -94,13 +122,34 @@ def ratio(latitude, apparent, vertical_shift)
         Angle.zero
       end
 
-      -(shift.sin + latitude.sin * @equatorial_coordinates.declination.sin)./(
-        latitude.cos * @equatorial_coordinates.declination.cos
-      )
+      term1 = shift.sin + latitude.sin * declination.sin
+      term2 = latitude.cos * declination.cos
+
+      -term1 / term2
+    end
+
+    def azimuth_ratio(latitude, apparent, vertical_shift)
+      shift = if vertical_shift
+        vertical_shift
+      elsif apparent
+        DEFAULT_REFRACTION_VERTICAL_SHIFT
+      else
+        Angle.zero
+      end
+
+      (declination.sin + shift.sin * latitude.cos) / (shift.cos * latitude.cos)
     end
 
     def azimuth_component(latitude)
-      @equatorial_coordinates.declination.sin / latitude.cos
+      declination.sin / latitude.cos
+    end
+
+    def right_ascension
+      @equatorial_coordinates.right_ascension
+    end
+
+    def declination
+      @equatorial_coordinates.declination
     end
   end
 end

spec/astronoby/bodies/sun_spec.rb

@@ -378,6 +378,67 @@
     end
   end
 
+  describe "#rising_azimuth" do
+    it "returns an Angle" do
+      date = Date.new
+      epoch = Astronoby::Epoch.from_time(date)
+      observer = Astronoby::Observer.new(
+        latitude: Astronoby::Angle.zero,
+        longitude: Astronoby::Angle.zero
+      )
+      sun = described_class.new(epoch: epoch)
+
+      rising_azimuth = sun.rising_azimuth(observer: observer)
+
+      expect(rising_azimuth).to be_a(Astronoby::Angle)
+    end
+
+    it "returns the Sun's rising azimuth on 2015-02-05" do
+      date = Date.new(2015, 2, 5)
+      epoch = Astronoby::Epoch.from_time(date)
+      observer = Astronoby::Observer.new(
+        latitude: Astronoby::Angle.as_degrees(38),
+        longitude: Astronoby::Angle.as_degrees(-78)
+      )
+      sun = described_class.new(epoch: epoch)
+
+      rising_azimuth = sun.rising_azimuth(observer: observer)
+
+      expect(rising_azimuth&.str(:dms)).to eq "+109° 41′ 24.0917″"
+      # Time from IMCCE: +109° 53′
+    end
+
+    it "returns the Sun's rising azimuth on 1986-03-10" do
+      date = Date.new(1986, 3, 10)
+      epoch = Astronoby::Epoch.from_time(date)
+      observer = Astronoby::Observer.new(
+        latitude: Astronoby::Angle.as_degrees(42.37),
+        longitude: Astronoby::Angle.as_degrees(-71.05)
+      )
+      sun = described_class.new(epoch: epoch)
+
+      rising_azimuth = sun.rising_azimuth(observer: observer)
+
+      expect(rising_azimuth&.str(:dms)).to eq "+94° 59′ 15.7852″"
+      # Time from IMCCE: +95° 02′
+    end
+
+    it "returns the Sun's rising azimuth on 1991-03-14" do
+      date = Date.new(1991, 3, 14)
+      epoch = Astronoby::Epoch.from_time(date)
+      observer = Astronoby::Observer.new(
+        latitude: Astronoby::Angle.as_degrees(48.8566),
+        longitude: Astronoby::Angle.as_degrees(2.3522)
+      )
+      sun = described_class.new(epoch: epoch)
+
+      rising_azimuth = sun.rising_azimuth(observer: observer)
+
+      expect(rising_azimuth&.str(:dms)).to eq "+93° 26′ 26.8564″"
+      # Time from IMCCE: +93° 26′
+    end
+  end
+
   describe "#setting_time" do
     it "returns a time" do
       date = Date.new
@@ -450,4 +511,65 @@
       # Time from IMCCE: 1991-03-14T17:52:00
     end
   end
+
+  describe "#setting_azimuth" do
+    it "returns an Angle" do
+      date = Date.new
+      epoch = Astronoby::Epoch.from_time(date)
+      observer = Astronoby::Observer.new(
+        latitude: Astronoby::Angle.zero,
+        longitude: Astronoby::Angle.zero
+      )
+      sun = described_class.new(epoch: epoch)
+
+      setting_azimuth = sun.setting_azimuth(observer: observer)
+
+      expect(setting_azimuth).to be_a(Astronoby::Angle)
+    end
+
+    it "returns the Sun's setting azimuth on 2015-02-05" do
+      date = Date.new(2015, 2, 5)
+      epoch = Astronoby::Epoch.from_time(date)
+      observer = Astronoby::Observer.new(
+        latitude: Astronoby::Angle.as_degrees(38),
+        longitude: Astronoby::Angle.as_degrees(-78)
+      )
+      sun = described_class.new(epoch: epoch)
+
+      setting_azimuth = sun.setting_azimuth(observer: observer)
+
+      expect(setting_azimuth&.str(:dms)).to eq "+250° 18′ 35.9082″"
+      # Time from IMCCE: +250° 18′
+    end
+
+    it "returns the Sun's setting azimuth on 1986-03-10" do
+      date = Date.new(1986, 3, 10)
+      epoch = Astronoby::Epoch.from_time(date)
+      observer = Astronoby::Observer.new(
+        latitude: Astronoby::Angle.as_degrees(42.37),
+        longitude: Astronoby::Angle.as_degrees(-71.05)
+      )
+      sun = described_class.new(epoch: epoch)
+
+      setting_azimuth = sun.setting_azimuth(observer: observer)
+
+      expect(setting_azimuth&.str(:dms)).to eq "+265° 0′ 44.2147″"
+      # Time from IMCCE: +265° 14′
+    end
+
+    it "returns the Sun's setting azimuth on 1991-03-14" do
+      date = Date.new(1991, 3, 14)
+      epoch = Astronoby::Epoch.from_time(date)
+      observer = Astronoby::Observer.new(
+        latitude: Astronoby::Angle.as_degrees(48.8566),
+        longitude: Astronoby::Angle.as_degrees(2.3522)
+      )
+      sun = described_class.new(epoch: epoch)
+
+      setting_azimuth = sun.setting_azimuth(observer: observer)
+
+      expect(setting_azimuth&.str(:dms)).to eq "+266° 33′ 33.1435″"
+      # Time from IMCCE: +266° 52′
+    end
+  end
 end