...

cxx/mount_astrom.h

@@ -7,9 +7,12 @@
  *********************************/
 
 #include <chrono>
+#include <fstream>
 #include <xsimd/xsimd.hpp>
-#include "utils.h"
 
+#include "mcc_traits.h"
+#include "mount_astrom_default.h"
+#include "utils.h"
 
 namespace mcc::traits
 {
@@ -393,3 +396,349 @@ static leapsecond_db_t mcc_parse_leapsecs(std::derived_from<std::basic_istream<c
 
 
 }  // namespace mcc::astro
+
+
+namespace mcc::astrom
+{
+
+
+class MccLeapSeconds
+{
+public:
+    MccLeapSeconds()
+    {
+        // create default values
+        std::istringstream ist(defaults::MCC_DEFAULT_LEAP_SECONDS_FILE);
+
+        load(ist);
+    }
+
+    ~MccLeapSeconds() = default;
+
+
+    std::chrono::system_clock::time_point expireDate() const
+    {
+        return _expireDate;
+    }
+
+    // load from stream
+    bool load(std::derived_from<std::basic_istream<char>> auto& stream, char comment_sym = '#')
+    {
+        std::istringstream is;
+        double mjd;
+        unsigned day, month;
+        int year;
+        double tai_utc;
+
+        decltype(_expireDate) edate;
+        std::vector<leapsecond_db_elem_t> db;
+
+        for (std::string line; std::getline(stream, line);) {
+            auto sv = utils::trimSpaces(line, utils::TrimType::TRIM_LEFT);
+
+            if (sv.size()) {
+                if (sv[0] == comment_sym) {  // comment string
+                    if (std::regex_match(line, expr_date_rx)) {
+                        auto pos = line.find("on");
+                        sv = utils::trimSpaces(std::string_view{line.begin() + pos + 2, line.end()},
+                                               utils::TrimType::TRIM_LEFT);
+                        is.str({sv.begin(), sv.end()});
+                        is >> std::chrono::parse("%d %B %Y", edate);
+                        is.clear();
+                    }
+                    continue;
+                }
+            } else {
+                continue;
+            }
+
+            if (std::regex_match(line, data_rx)) {
+                is.str(line);
+                is >> mjd >> day >> month >> year >> tai_utc;
+                db.emplace_back(mjd, std::chrono::year_month_day{std::chrono::year{year} / month / day}, tai_utc);
+                // db.emplace_back(mjd,
+                //                 std::chrono::year_month_day{std::chrono::year{year}, std::chrono::month{month},
+                //                                             std::chrono::day{day}},
+                //                 tai_utc);
+                is.clear();
+
+                continue;
+            }
+        }
+
+        if (db.empty()) {  // keep previous data
+            return false;
+        }
+
+        _expireDate = std::move(edate);
+        _db = std::move(db);
+
+        return true;
+    }
+
+
+    bool load(traits::mcc_input_char_range auto const& filename, char comment_sym = '#')
+    {
+        std::ifstream fst(filename);
+
+        bool ok = fst.is_open();
+        if (!ok) {
+            return false;
+        }
+
+        ok = load(fst, comment_sym);
+
+        fst.close();
+
+        return ok;
+    }
+
+
+    std::optional<double> operator[](const std::chrono::system_clock::time_point& tp) const
+    {
+        if (tp > _expireDate) {  // ???????!!!!!!!!!!!
+            return std::nullopt;
+            // return _db.back().tai_utc;
+        }
+
+        std::chrono::year_month_day ymd{std::chrono::floor<std::chrono::days>(tp)};
+
+        for (auto const& el : _db | std::views::reverse) {
+            if (ymd >= el.ymd) {
+                return el.tai_utc;
+            }
+        }
+
+        return std::nullopt;
+    }
+
+    std::optional<double> operator[](const double& mjd) const
+    {
+        double e_mjd;
+
+        astro::mcc_julday(_expireDate, e_mjd);
+        if (mjd > e_mjd) {  // ???????!!!!!!!!!!!
+            return std::nullopt;
+            // return _db.back().tai_utc;
+        }
+
+        for (auto const& el : _db | std::views::reverse) {
+            if (mjd >= el.mjd) {
+                return el.tai_utc;
+            }
+        }
+
+        return std::nullopt;
+    }
+
+private:
+    inline static const std::regex expr_date_rx{
+        "^ *# *File +expires +on +[0-8]{1,2} "
+        "+(January|February|March|April|May|June|July|August|September|October|November|December) +[0-9]{4} *$"};
+
+    inline static const std::regex data_rx{"^ *[0-9]{5,}(\\.?[0-9]+) +[0-9]{1,2} +[0-9]{1,2} +[0-9]{4} +[0-9]{1,} *$"};
+
+    std::chrono::system_clock::time_point _expireDate{};
+
+    struct leapsecond_db_elem_t {
+        double mjd;
+        std::chrono::year_month_day ymd;
+        double tai_utc;  // TAI-UTC in seconds
+    };
+
+    std::vector<leapsecond_db_elem_t> _db{};
+};
+
+
+
+class MccIersBulletinA
+{
+public:
+    struct pole_pos_t {
+        double x, y;
+    };
+
+    MccIersBulletinA()
+    {
+        std::istringstream ist(defaults::MCC_DEFAULT_IERS_BULLETIN_A_FILE);
+
+        load(ist);
+    }
+
+    std::chrono::system_clock::time_point bulletinDate() const
+    {
+        return _date;
+    }
+
+    double TT_TAI() const
+    {
+        return _tt_tai;
+    }
+
+    // DUT1 = UT1 - UTC
+    std::optional<double> DUT1(const std::chrono::system_clock::time_point& tp) const
+    {
+        std::chrono::year_month_day ymd{std::chrono::floor<std::chrono::days>(tp)};
+
+        for (auto const& el : _db | std::views::reverse) {
+            if (ymd >= el.ymd) {
+                return el.dut1;
+            }
+        }
+
+        return std::nullopt;
+    }
+
+    std::optional<double> DUT1(const double& mjd) const
+    {
+        for (auto const& el : _db | std::views::reverse) {
+            if (mjd >= el.mjd) {
+                return el.dut1;
+            }
+        }
+
+        return std::nullopt;
+    }
+
+    std::optional<pole_pos_t> polePos(const std::chrono::system_clock::time_point& tp) const
+    {
+        std::chrono::year_month_day ymd{std::chrono::floor<std::chrono::days>(tp)};
+
+        for (auto const& el : _db | std::views::reverse) {
+            if (ymd >= el.ymd) {
+                return pole_pos_t{el.x, el.y};
+            }
+        }
+
+        return std::nullopt;
+    }
+
+    std::optional<pole_pos_t> polePos(const double& mjd) const
+    {
+        for (auto const& el : _db | std::views::reverse) {
+            if (mjd >= el.mjd) {
+                return pole_pos_t{el.x, el.y};
+            }
+        }
+
+        return std::nullopt;
+    }
+
+    bool load(std::derived_from<std::basic_istream<char>> auto& stream, char comment_sym = '*')
+    {
+        std::vector<earth_orient_db_elem_t> db;
+        enum { TAB_STATE_SEEK, TAB_STATE_START };
+        int tab_state = TAB_STATE_SEEK;
+
+        int year;
+        unsigned month, day;
+        double mjd, x, y, dut1;
+        std::istringstream is;
+        decltype(_date) bdate;
+        double tt_tai;
+
+        for (std::string line; std::getline(stream, line);) {
+            if (line.empty()) {
+                continue;
+            }
+
+            auto sv = utils::trimSpaces(line, utils::TrimType::TRIM_LEFT);
+
+            if (sv.size()) {
+                if (sv[0] == comment_sym) {  // comment string
+                    continue;
+                }
+
+                if (tab_state == TAB_STATE_START) {
+                    if (std::regex_match(sv.begin(), sv.end(), bull_tab_vals_rx)) {
+                        // is.str({sv.begin(), sv.end()});
+                        is.str(line);
+                        is >> year >> month >> day >> mjd >> x >> y >> dut1;
+                        db.emplace_back(mjd, std::chrono::year_month_day{std::chrono::year{year} / month / day}, x, y,
+                                        dut1);
+                        is.clear();
+                    } else {  // end of the table - just stop parsing
+                        break;
+                    }
+
+                    continue;
+                }
+
+                if (std::regex_match(sv.begin(), sv.end(), bull_date_rx)) {
+                    is.str({sv.begin(), sv.end()});
+                    is >> std::chrono::parse("%d %B %Y", bdate);
+                    continue;
+                }
+
+                if (std::regex_match(sv.begin(), sv.end(), bull_tt_tai_rx)) {
+                    is.str({sv.begin(), sv.end()});
+                    std::string dummy;
+                    is >> dummy >> dummy >> dummy >> dummy >> tt_tai;
+                    continue;
+                }
+
+                if (std::regex_match(sv.begin(), sv.end(), bull_tab_title_rx)) {
+                    tab_state = TAB_STATE_START;
+                    continue;
+                }
+
+            } else {  // empty string (only spaces)
+                continue;
+            }
+        }
+
+        if (db.empty()) {
+            return false;
+        }
+
+        _date = std::move(bdate);
+        _tt_tai = tt_tai;
+        _db = std::move(db);
+
+        return true;
+    }
+
+    bool load(traits::mcc_input_char_range auto const& filename, char comment_sym = '*')
+    {
+        std::ifstream fst(filename);
+
+        bool ok = fst.is_open();
+        if (!ok) {
+            return false;
+        }
+
+        ok = load(fst, comment_sym);
+
+        fst.close();
+
+        return ok;
+    }
+
+private:
+    inline static const std::regex bull_date_rx{
+        "^ *[0-9]{1,2} +(January|February|March|April|May|June|July|August|September|October|November|December) "
+        "+[0-9]{4,} +Vol\\. +[XMLCDVI]+ +No\\. +[0-9]+ *$"};
+
+    inline static const std::regex bull_tt_tai_rx{"^ *TT += +TAI +\\+ +[0-9]+\\.[0-9]+ +seconds *$"};
+
+    inline static const std::regex bull_tab_title_rx{"^ *MJD +x\\(arcsec\\) +y\\(arcsec\\) +UT1-UTC\\(sec\\) *$"};
+
+    inline static const std::regex bull_tab_vals_rx{
+        "^ *[0-9]{4,} +[0-9]{1,2} +[0-9]{1,2} +[0-9]{5,} +[0-9]+\\.[0-9]+ +[0-9]+\\.[0-9]+ +[0-9]+\\.[0-9]+ *$"};
+
+
+    std::chrono::system_clock::time_point _date;
+    double _tt_tai;
+
+    struct earth_orient_db_elem_t {
+        double mjd;
+        std::chrono::year_month_day ymd;
+        double x, y;  // Polar coordinates in arcsecs
+        double dut1;  // UT1-UTC in seconds
+    };
+
+    std::vector<earth_orient_db_elem_t> _db;
+};
+
+
+}  // namespace mcc::astrom