ash/system/time/astronomer_util.cc

// Copyright 2025 The Chromium Authors
// Use of this source code is governed by a BSD-style license that can be
// found in the LICENSE file.

#include "ash/system/time/astronomer_util.h"

#include <cmath>
#include <ctime>
#include <string>

#include "base/logging.h"
#include "base/time/time.h"
#include "third_party/icu/source/i18n/astro.h"

namespace ash {
namespace {

base::expected<base::Time, SunRiseSetError> GetTime(
    base::Time::Exploded& exploded,
    double hours_in_day) {
  // We only need minutes level precision.
  exploded.hour = std::floor(hours_in_day);
  exploded.minute = fmod(hours_in_day * 60, 60);
  exploded.second = 0;
  exploded.millisecond = 0;

  base::Time ret;
  if (!base::Time::FromUTCExploded(exploded, &ret)) {
    return base::unexpected(SunRiseSetError::kUnavailable);
  }
  return base::ok(ret);
}

}  // namespace

// The most of this function was generated by Gemini.
base::expected<SunRiseSetTime, SunRiseSetError>
GetSunriseSunset(const base::Time& time, double latitude, double longitude) {
  base::Time::Exploded exploded;
  time.UTCExplode(&exploded);

  std::tm t = {};
  t.tm_year = exploded.year - 1900;
  t.tm_mon = exploded.month - 1;
  t.tm_mday = exploded.day_of_month;
  std::mktime(&t);
  int day_of_year = t.tm_yday + 1;

  auto to_radians = [](double degrees) { return degrees * M_PI / 180.0; };
  auto to_degrees = [](double radians) { return radians * 180.0 / M_PI; };

  // Zenith for sunrise/sunset calculation (official is 90 degrees 50 minutes)
  const double zenith = 90.8333;

  // 2. Convert the longitude to hour value and calculate an approximate time
  double lng_hour = longitude / 15.0;
  double t_rise = day_of_year + ((6.0 - lng_hour) / 24.0);
  double t_set = day_of_year + ((18.0 - lng_hour) / 24.0);

  // 3. Calculate the Sun's mean anomaly
  double m_rise = (0.9856 * t_rise) - 3.289;
  double m_set = (0.9856 * t_set) - 3.289;

  // 4. Calculate the Sun's true longitude
  double l_rise = m_rise + (1.916 * sin(to_radians(m_rise))) +
                  (0.020 * sin(to_radians(2 * m_rise))) + 282.634;
  double l_set = m_set + (1.916 * sin(to_radians(m_set))) +
                 (0.020 * sin(to_radians(2 * m_set))) + 282.634;

  l_rise = fmod(l_rise, 360.0);
  l_set = fmod(l_set, 360.0);

  // 5a. Calculate the Sun's right ascension
  double ra_rise = to_degrees(atan(0.91764 * tan(to_radians(l_rise))));
  double ra_set = to_degrees(atan(0.91764 * tan(to_radians(l_set))));

  ra_rise = fmod(ra_rise, 360.0);
  ra_set = fmod(ra_set, 360.0);

  // 5b. Right ascension value needs to be in the same quadrant as L
  double l_quadrant_rise = floor(l_rise / 90.0) * 90.0;
  double ra_quadrant_rise = floor(ra_rise / 90.0) * 90.0;
  ra_rise = ra_rise + (l_quadrant_rise - ra_quadrant_rise);

  double l_quadrant_set = floor(l_set / 90.0) * 90.0;
  double ra_quadrant_set = floor(ra_set / 90.0) * 90.0;
  ra_set = ra_set + (l_quadrant_set - ra_quadrant_set);

  // 5c. Right ascension value needs to be converted into hours
  double ra_rise_hours = ra_rise / 15.0;
  double ra_set_hours = ra_set / 15.0;

  // 6. Calculate the Sun's declination
  double sin_dec_rise = 0.39782 * sin(to_radians(l_rise));
  double cos_dec_rise = cos(asin(sin_dec_rise));

  double sin_dec_set = 0.39782 * sin(to_radians(l_set));
  double cos_dec_set = cos(asin(sin_dec_set));

  // 7a. Calculate the Sun's local hour angle
  double cos_h_rise =
      (cos(to_radians(zenith)) - (sin_dec_rise * sin(to_radians(latitude)))) /
      (cos_dec_rise * cos(to_radians(latitude)));
  double cos_h_set =
      (cos(to_radians(zenith)) - (sin_dec_set * sin(to_radians(latitude)))) /
      (cos_dec_set * cos(to_radians(latitude)));

  if (cos_h_rise > 1.0) {
    return base::unexpected(SunRiseSetError::kNoSunRiseSet);
  }
  if (cos_h_set < -1.0) {
    return base::unexpected(SunRiseSetError::kNoSunRiseSet);
  }

  // 7b. Finish calculating H and convert into hours
  double h_rise_deg = 360.0 - to_degrees(acos(cos_h_rise));
  double h_set_deg = to_degrees(acos(cos_h_set));

  double h_rise_hours = h_rise_deg / 15.0;
  double h_set_hours = h_set_deg / 15.0;

  // 8. Calculate local mean time of rising/setting
  double t_local_rise =
      h_rise_hours + ra_rise_hours - (0.06571 * t_rise) - 6.622;
  double t_local_set = h_set_hours + ra_set_hours - (0.06571 * t_set) - 6.622;

  // 9. Adjust back to UTC
  double ut_rise = t_local_rise - lng_hour;
  double ut_set = t_local_set - lng_hour;

  // Normalize to 0-24
  double sunrise_utc_hours = fmod(ut_rise + 24.0, 24.0);
  double sunset_utc_hours = fmod(ut_set + 24.0, 24.0);

  SunRiseSetTime result;
  auto sunrise = GetTime(exploded, sunrise_utc_hours);
  if (!sunrise.has_value()) {
    return base::unexpected(sunrise.error());
  }
  result.sunrise = *sunrise;

  auto sunset = GetTime(exploded, sunset_utc_hours);
  if (!sunset.has_value()) {
    return base::unexpected(sunset.error());
  }
  result.sunset = *sunset;
  // Sunrise time may be 1 day earlier if the time was normalized to 24 hours.
  if (result.sunset < result.sunrise) {
    result.sunset += base::Days(1);
  }
  return result;
}

base::expected<SunRiseSetTime, SunRiseSetError>
GetSunriseSunsetICU(const base::Time& time, double latitude, double longitude) {
  icu::CalendarAstronomer astro(longitude, latitude);
  astro.setTime(time.InMillisecondsFSinceUnixEpoch());
  const double sun_rise_ms = astro.getSunRiseSet(/*sunrise=*/true);
  if (sun_rise_ms < 0) {
    return base::unexpected(SunRiseSetError::kNoSunRiseSet);
  }
  const double sun_set_ms = astro.getSunRiseSet(/*sunrise=*/false);
  if (sun_set_ms < 0) {
    return base::unexpected(SunRiseSetError::kNoSunRiseSet);
  }
  return base::ok(SunRiseSetTime{
      .sunrise = base::Time::FromMillisecondsSinceUnixEpoch(sun_rise_ms),
      .sunset = base::Time::FromMillisecondsSinceUnixEpoch(sun_set_ms)});
}

}  // namespace ash