astronomy_conversions.hpp

Go to the documentation of this file.

// SPDX-License-Identifier: MIT
#pragma once
#ifndef _TIME_SHIELD_ASTRONOMY_CONVERSIONS_HPP_INCLUDED
#define _TIME_SHIELD_ASTRONOMY_CONVERSIONS_HPP_INCLUDED

 
#include "config.hpp"
#include "types.hpp"
#include "constants.hpp"
#include "MoonPhase.hpp"
#include <cmath>
 
namespace time_shield {

    inline jd_t fts_to_jd(fts_t ts) noexcept {
        // JD at Unix epoch:
        // 1970-01-01 00:00:00 UTC -> 2440587.5
        return static_cast<jd_t>(2440587.5)
             + static_cast<jd_t>(ts) / static_cast<jd_t>(SEC_PER_DAY);
    }

    inline jd_t ts_to_jd(ts_t ts) noexcept {
        return fts_to_jd(static_cast<fts_t>(ts));
    }

    inline jd_t gregorian_to_jd(double day, int64_t month, int64_t year) noexcept {
        // Algorithm source (as per your original code): krutov.org Julianday
        if (month == 1 || month == 2) {
            year -= 1;
            month += 12;
        }
        const double a = std::floor(static_cast<double>(year) / 100.0);
        const double b = 2.0 - a + std::floor(a / 4.0);
        const double jdn = std::floor(365.25 * (static_cast<double>(year) + 4716.0))
                         + std::floor(30.6000001 * (static_cast<double>(month) + 1.0))
                         + day + b - 1524.5;
        return static_cast<jd_t>(jdn);
    }
    
    inline jd_t gregorian_to_jd(
        uint32_t day,
        uint32_t month,
        uint32_t year,
        uint32_t hour,
        uint32_t minute,
        uint32_t second = 0,
        uint32_t millisecond = 0) noexcept {
        const double frac =
            (static_cast<double>(hour) / 24.0) +
            (static_cast<double>(minute) / (24.0 * 60.0)) +
            ((static_cast<double>(second) + static_cast<double>(millisecond) / 1000.0)
              / static_cast<double>(SEC_PER_DAY));
        return gregorian_to_jd(static_cast<double>(day) + frac,
                               static_cast<int64_t>(month),
                               static_cast<int64_t>(year));
    }
    
    inline mjd_t fts_to_mjd(fts_t ts) noexcept {
        return static_cast<mjd_t>(fts_to_jd(ts) - 2400000.5);
    }

    inline mjd_t ts_to_mjd(ts_t ts) noexcept {
        return static_cast<mjd_t>(fts_to_mjd(static_cast<fts_t>(ts)));
    }

    inline jdn_t gregorian_to_jdn(uint32_t day, uint32_t month, uint32_t year) noexcept {
        const uint64_t a = (14ULL - static_cast<uint64_t>(month)) / 12ULL;
        const uint64_t y = static_cast<uint64_t>(year) + 4800ULL - a;
        const uint64_t m = static_cast<uint64_t>(month) + 12ULL * a - 3ULL;
        const uint64_t jdn = static_cast<uint64_t>(day)
                           + (153ULL * m + 2ULL) / 5ULL
                           + 365ULL * y
                           + y / 4ULL
                           - y / 100ULL
                           + y / 400ULL
                           - 32045ULL;
        return static_cast<jdn_t>(jdn);
    }

    struct MoonPhaseSineCosine {
        double phase_sin = 0.0; 
        double phase_cos = 0.0; 
        double phase_angle_rad = 0.0; 
        constexpr MoonPhaseSineCosine() = default;
        constexpr MoonPhaseSineCosine(double phase_sin_value, double phase_cos_value, double phase_angle_rad_value) noexcept
            : phase_sin(phase_sin_value),
              phase_cos(phase_cos_value),
              phase_angle_rad(phase_angle_rad_value) {}
    };

    inline double moon_phase_jd_approx(fts_t ts) noexcept {
        double temp = (static_cast<double>(fts_to_jd(ts)) - 2451550.1) / 29.530588853;
        temp = temp - std::floor(temp);
        if (temp < 0.0) temp += 1.0;
        return temp;
    }

    inline double moon_phase(fts_t ts) noexcept {
        static const astronomy::MoonPhase calculator{};
        return calculator.compute(static_cast<double>(ts)).phase;
    }

    inline MoonPhaseSineCosine moon_phase_sincos(fts_t ts) noexcept {
        static const astronomy::MoonPhase calculator{};
        const auto result = calculator.compute(static_cast<double>(ts));
        return MoonPhaseSineCosine{result.phase_sin, result.phase_cos, result.phase_angle_rad};
    }

    inline double moon_illumination(fts_t ts) noexcept {
        static const astronomy::MoonPhase calculator{};
        return calculator.compute(static_cast<double>(ts)).illumination;
    }

    inline double moon_age_days_jd_approx(fts_t ts) noexcept {
        return moon_phase_jd_approx(ts) * 29.530588853;
    }

    inline double moon_age_days(fts_t ts) noexcept {
        static const astronomy::MoonPhase calculator{};
        return calculator.compute(static_cast<double>(ts)).age_days;
    }

    inline astronomy::MoonQuarterInstants moon_quarters(fts_t ts) noexcept {
        static const astronomy::MoonPhase calculator{};
        return calculator.quarter_instants_unix(static_cast<double>(ts));
    }

    inline bool is_new_moon_window(fts_t ts, double window_seconds = astronomy::MoonPhase::kDefaultQuarterWindow_s) noexcept {
        static const astronomy::MoonPhase calculator{};
        return calculator.is_new_moon_window(static_cast<double>(ts), window_seconds);
    }

    inline bool is_full_moon_window(fts_t ts, double window_seconds = astronomy::MoonPhase::kDefaultQuarterWindow_s) noexcept {
        static const astronomy::MoonPhase calculator{};
        return calculator.is_full_moon_window(static_cast<double>(ts), window_seconds);
    }

    inline bool is_first_quarter_window(fts_t ts, double window_seconds = astronomy::MoonPhase::kDefaultQuarterWindow_s) noexcept {
        static const astronomy::MoonPhase calculator{};
        return calculator.is_first_quarter_window(static_cast<double>(ts), window_seconds);
    }

    inline bool is_last_quarter_window(fts_t ts, double window_seconds = astronomy::MoonPhase::kDefaultQuarterWindow_s) noexcept {
        static const astronomy::MoonPhase calculator{};
        return calculator.is_last_quarter_window(static_cast<double>(ts), window_seconds);
    }
 
} // namespace time_shield
 
#endif // _TIME_SHIELD_ASTRONOMY_CONVERSIONS_HPP_INCLUDED