#pragma once

/*                          MOUNT CONTROL COMPONENTS LIBRARY                            */


/*                          CELESTIAL COORDINATES TRANSFORMATION ENGINE IMPLEMENTATION BASED ON ERFA LIBRARY   */


// NOTE: according to definition of astronomical azimuth it is counted from the South through the West, but
// in the ERFA the azimuth is counted from the North through the East!!!
//
// The implementation expects input azimuth as the "South"-defined one and transforms ERFA-outputs
// to the "South"-defined one respectively



#include <erfa.h>
#include <erfam.h>
#include <mutex>


#include "mcc_angle.h"
#include "mcc_ccte_iers.h"
#include "mcc_defaults.h"

namespace mcc::ccte::erfa
{


enum class MccCCTE_ERFAErrorCode : int {
    ERROR_OK = 0,
    ERROR_INVALID_INPUT_ARG,
    ERROR_julday_INVALID_YEAR,
    ERROR_julday_INVALID_MONTH,
    ERROR_julday_INVALID_DAY,
    ERROR_UNSUPPORTED_COORD_PAIR,
    ERROR_BULLETINA_OUT_OF_RANGE,
    ERROR_LEAPSECONDS_OUT_OF_RANGE,
    ERROR_DUBIOUS_YEAR,
    ERROR_UNACCEPTABLE_DATE,
    ERROR_UPDATE_LEAPSECONDS,
    ERROR_UPDATE_BULLETINA,
    ERROR_UNEXPECTED
};

}  // namespace mcc::ccte::erfa


namespace std
{

template <>
class is_error_code_enum<mcc::ccte::erfa::MccCCTE_ERFAErrorCode> : public true_type
{
};

}  // namespace std



namespace mcc::ccte::erfa
{


/*  error category definition  */

// error category
struct MccCCTE_ERFACategory : public std::error_category {
    MccCCTE_ERFACategory() : std::error_category() {}

    const char* name() const noexcept
    {
        return "CCTE-ERFA";
    }

    std::string message(int ec) const
    {
        MccCCTE_ERFAErrorCode err = static_cast<MccCCTE_ERFAErrorCode>(ec);

        switch (err) {
            case MccCCTE_ERFAErrorCode::ERROR_OK:
                return "OK";
            case MccCCTE_ERFAErrorCode::ERROR_INVALID_INPUT_ARG:
                return "invalid argument";
            case MccCCTE_ERFAErrorCode::ERROR_julday_INVALID_YEAR:
                return "invalid year number";
            case MccCCTE_ERFAErrorCode::ERROR_julday_INVALID_MONTH:
                return "invalid month number";
            case MccCCTE_ERFAErrorCode::ERROR_julday_INVALID_DAY:
                return "invalid day number";
            case MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR:
                return "unsupported coordinate pair";
            case MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE:
                return "time point is out of range";
            case MccCCTE_ERFAErrorCode::ERROR_LEAPSECONDS_OUT_OF_RANGE:
                return "time point is out of range";
            case MccCCTE_ERFAErrorCode::ERROR_DUBIOUS_YEAR:
                return "dubious year";
            case MccCCTE_ERFAErrorCode::ERROR_UNACCEPTABLE_DATE:
                return "unacceptable date";
            case MccCCTE_ERFAErrorCode::ERROR_UPDATE_LEAPSECONDS:
                return "leap seconds update error";
            case MccCCTE_ERFAErrorCode::ERROR_UPDATE_BULLETINA:
                return "bulletin A update error";
            case MccCCTE_ERFAErrorCode::ERROR_UNEXPECTED:
                return "unexpected error value";
            default:
                return "UNKNOWN";
        }
    }

    static const MccCCTE_ERFACategory& get()
    {
        static const MccCCTE_ERFACategory constInst;
        return constInst;
    }
};


inline std::error_code make_error_code(MccCCTE_ERFAErrorCode ec)
{
    return std::error_code(static_cast<int>(ec), MccCCTE_ERFACategory::get());
}


class MccCCTE_ERFA : public mcc_CCTE_interface_t<std::error_code>
{
    static constexpr double PI_2 = std::numbers::pi / 2.0;

public:
    static constexpr double DEFAULT_WAVELENGTH = 0.55;  // default observed wavelength in mkm

    typedef std::error_code error_t;

    struct refract_model_t {
        static constexpr std::string_view name()
        {
            return "ERFA";
        }

        double refa, refb;
    };

    /*  use of the same type for representation of celestial and geodetic coordinates, celestial angles (e.g. P.A.),
     *  and sideral time */
    typedef MccCelestialPoint::coord_t coord_t;

    // meteo parameters (to compute refraction)
    struct meteo_t {
        typedef double temp_t;
        typedef double humid_t;
        typedef double press_t;

        temp_t temperature;  // Temperature in C
        humid_t humidity;    // humidity in % ([0.0, 1.0])
        press_t pressure;    // atmospheric presure in hPa=mB
    };

    struct engine_state_t {
        meteo_t meteo{.temperature = 0.0, .humidity = 0.5, .pressure = 1010.0};

        double wavelength = DEFAULT_WAVELENGTH;  // observed wavelength in mkm

        coord_t lat = 0.0;  // site latitude
        coord_t lon = 0.0;  // site longitude
        double elev = 0.0;  // site elevation (in meters)

        mcc::ccte::iers::MccLeapSeconds _leapSeconds{};
        mcc::ccte::iers::MccIersBulletinA _bulletinA{};
    };

    MccCCTE_ERFA() : _stateMutex(new std::mutex) {}

    MccCCTE_ERFA(engine_state_t state) : _currentState(std::move(state)), _stateMutex(new std::mutex) {}

    MccCCTE_ERFA(const MccCCTE_ERFA&) = delete;
    MccCCTE_ERFA& operator=(const MccCCTE_ERFA&) = delete;

    MccCCTE_ERFA(MccCCTE_ERFA&&) = default;
    MccCCTE_ERFA& operator=(MccCCTE_ERFA&&) = default;

    virtual ~MccCCTE_ERFA() = default;


    std::string_view nameCCTE() const
    {
        return "ERFA-CCTE-ENGINE";
    }

    // engine state related methods

    void setStateERFA(engine_state_t state)
    {
        std::lock_guard lock{*_stateMutex};

        _currentState = std::move(state);
    }

    engine_state_t getStateERFA() const
    {
        std::lock_guard lock{*_stateMutex};

        return _currentState;
    }

    void updateMeteoERFA(meteo_t meteo)
    {
        std::lock_guard lock{*_stateMutex};

        _currentState.meteo = std::move(meteo);
    }

    error_t updateLeapSeconds(std::derived_from<std::basic_istream<char>> auto& stream, char comment_sym = '#')
    {
        std::lock_guard lock{*_stateMutex};

        if (!_currentState._leapSeconds.load(stream, comment_sym)) {
            return MccCCTE_ERFAErrorCode::ERROR_UPDATE_LEAPSECONDS;
        }

        return MccCCTE_ERFAErrorCode::ERROR_OK;
    }


    error_t updateLeapSeconds(traits::mcc_input_char_range auto const& filename, char comment_sym = '#')
    {
        std::lock_guard lock{*_stateMutex};

        if (!_currentState._leapSeconds.load(filename, comment_sym)) {
            return MccCCTE_ERFAErrorCode::ERROR_UPDATE_LEAPSECONDS;
        }

        return MccCCTE_ERFAErrorCode::ERROR_OK;
    }


    error_t updateBulletinA(std::derived_from<std::basic_istream<char>> auto& stream, char comment_sym = '*')
    {
        std::lock_guard lock{*_stateMutex};

        if (!_currentState._bulletinA.load(stream, comment_sym)) {
            return MccCCTE_ERFAErrorCode::ERROR_UPDATE_BULLETINA;
        }

        return MccCCTE_ERFAErrorCode::ERROR_OK;
    }


    error_t updateBulletinA(traits::mcc_input_char_range auto const& filename, char comment_sym = '*')
    {
        std::lock_guard lock{*_stateMutex};

        if (!_currentState._bulletinA.load(filename, comment_sym)) {
            return MccCCTE_ERFAErrorCode::ERROR_UPDATE_BULLETINA;
        }

        return MccCCTE_ERFAErrorCode::ERROR_OK;
    }


    // time-related methods

    error_t timepointToJulday(mcc_time_point_c auto tp, mcc_julday_c auto* julday)
    {
        error_t ret = MccCCTE_ERFAErrorCode::ERROR_OK;

        if (julday == nullptr) {
            return ret;
        }

        auto days = std::chrono::floor<std::chrono::days>(tp);
        std::chrono::year_month_day ymd{days};

        double mjd0, mjd;

        int err = eraCal2jd((int)ymd.year(), (unsigned)ymd.month(), (unsigned)ymd.day(), &mjd0, &mjd);

        if (err != 0) {
            ret = err == -1   ? MccCCTE_ERFAErrorCode::ERROR_julday_INVALID_YEAR
                  : err == -2 ? MccCCTE_ERFAErrorCode::ERROR_julday_INVALID_MONTH
                  : err == -3 ? MccCCTE_ERFAErrorCode::ERROR_julday_INVALID_DAY
                              : MccCCTE_ERFAErrorCode::ERROR_UNEXPECTED;
        } else {  // partial part of day
            mjd += std::chrono::duration_cast<std::chrono::duration<double, std::ratio<86400>>>(tp - days).count();
            *julday = mjd + mjd0;
        }

        return ret;
    }

    error_t juldayToAppSideral(mcc_julday_c auto jd, mcc_angle_c auto* st, bool islocal = false)
    {
        error_t ret = MccCCTE_ERFAErrorCode::ERROR_OK;

        if (st == nullptr) {
            return ret;
        }

        using real_days_t = std::chrono::duration<double, std::ratio<86400>>;

        using jd_t = decltype(jd);

        double mjd;
        if constexpr (std::floating_point<jd_t>) {
            mjd = jd - ERFA_DJM0;
        } else {
            mjd = jd.MJD();
        }

        double ut1 = mjd;
        double tt = mjd;

        std::lock_guard lock{*_stateMutex};

        auto dut1 = _currentState._bulletinA.DUT1(mjd);

        if (dut1.has_value()) {
            ut1 += std::chrono::duration_cast<real_days_t>(dut1.value()).count();
        } else {  // out of range
            return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
        }

        auto tai_utc = _currentState._leapSeconds[mjd];
        if (tai_utc.has_value()) {
            tt += std::chrono::duration_cast<real_days_t>(tai_utc.value()).count();
        } else {
            return MccCCTE_ERFAErrorCode::ERROR_LEAPSECONDS_OUT_OF_RANGE;
        }


        auto tt_tai = _currentState._bulletinA.TT_TAI();
        tt += std::chrono::duration_cast<real_days_t>(tt_tai).count();

        *st = eraGst06a(ERFA_DJM0, ut1, ERFA_DJM0, tt);

        if (islocal) {
            // *st += _currentState.lon;
            // *st = MccAngle(*st + _currentState.lon).normalize<MccAngle::NORM_KIND_0_360>();
            *st = eraAnp(*st + _currentState.lon);
        }


        return ret;
    }

    error_t timepointToAppSideral(mcc_time_point_c auto tp, mcc_angle_c auto* st, bool islocal = false)
    {
        error_t ret = MccCCTE_ERFAErrorCode::ERROR_OK;

        if (st == nullptr) {
            return ret;
        }

        MccJulianDay julday;

        ret = timepointToJulday(tp, &julday);
        if (ret != MccCCTE_ERFAErrorCode::ERROR_OK) {
            return ret;
        }

        return juldayToAppSideral(julday, st, islocal);

        // double ut1 = julday.MJD();
        // double tt = julday.MJD();

        // std::lock_guard lock{*_stateMutex};

        // auto dut1 = _currentState._bulletinA.DUT1(julday.MJD());

        // if (dut1.has_value()) {
        //     ut1 += std::chrono::duration_cast<real_days_t>(dut1.value()).count();
        // } else {  // out of range
        //     return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
        // }

        // auto tai_utc = _currentState._leapSeconds[julday.MJD()];
        // if (tai_utc.has_value()) {
        //     tt += std::chrono::duration_cast<real_days_t>(tai_utc.value()).count();
        // } else {
        //     return MccCCTE_ERFAErrorCode::ERROR_LEAPSECONDS_OUT_OF_RANGE;
        // }


        // auto tt_tai = _currentState._bulletinA.TT_TAI();
        // tt += std::chrono::duration_cast<real_days_t>(tt_tai).count();

        // *st = eraGst06a(ERFA_DJM0, ut1, ERFA_DJM0, tt);

        // if (islocal) {
        //     *st += _currentState.lon;
        // }


        // return ret;
    }


    // coordinates transformations
    template <typename ResT>
    error_t transformCoordinates(mcc_celestial_point_c auto from_pt, ResT* to_pt)
        requires(mcc_eqt_hrz_coord_c<ResT> || mcc_celestial_point_c<ResT>)
    {
        if constexpr (mcc_eqt_hrz_coord_c<ResT>) {
            return transformCoordinatesEQHR(std::move(from_pt), to_pt);
        } else if constexpr (mcc_celestial_point_c<ResT>) {
            return transformCoordinatesCP(std::move(from_pt), to_pt);
        } else {
            static_assert(false, "UNSUPPORTED TYPE!");
        }
    }

    error_t transformCoordinatesCP(mcc_celestial_point_c auto from_pt, mcc_celestial_point_c auto* to_pt)
    {
        error_t ret = MccCCTE_ERFAErrorCode::ERROR_OK;

        // MccJulianDay jd;

        if (to_pt == nullptr) {
            return ret;
        }

        // no transformations
        if (from_pt.time_point == to_pt->time_point && from_pt.pair_kind == to_pt->pair_kind) {
            to_pt->X = from_pt.X;
            to_pt->Y = from_pt.Y;

            return ret;
        }


        // special case: to ICRS from apparent
        if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_ICRS) {
            if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_AZALT ||
                from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_AZZD ||
                from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_APP ||
                from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_HADEC_APP) {
                //
            } else {
                ret = MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
            }

            ret = obs2icrs(from_pt, to_pt);

            return ret;
        }


        // special case: from ICRS to apparent
        if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_ICRS) {
            if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZALT ||
                to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZZD ||
                to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_APP ||
                to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_HADEC_APP) {
                //
                MccEqtHrzCoords eqhrz;

                ret = icrs2obs(from_pt, &eqhrz);

                if (!ret) {
                    if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZALT) {
                        to_pt->X = eqhrz.AZ;
                        to_pt->Y = eqhrz.ALT;
                    } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZZD) {
                        to_pt->X = eqhrz.AZ;
                        to_pt->Y = eqhrz.ZD;
                    } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_APP) {
                        to_pt->X = eqhrz.RA_APP;
                        to_pt->Y = eqhrz.DEC_APP;
                    } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_HADEC_APP) {
                        to_pt->X = eqhrz.HA;
                        to_pt->Y = eqhrz.DEC_APP;
                    }
                }
            } else {
                ret = MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
            }

            return ret;
        }


        // from apparent to apparent
        if (from_pt.time_point != to_pt->time_point) {  // first, convert 'from' coordinates to ICRS
            MccCelestialPoint pt{.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_ICRS};
            pt.time_point = to_pt->time_point;

            ret = obs2icrs(from_pt, &pt);

            if (!ret) {  // from ICRS to required
                return transformCoordinates(pt, to_pt);
            }
        } else {  // the same time points
            double eo, lst, ha;
            // , dec, az, alt;

            auto lst_eo = [&]() {
                // ret = eqOrigins(from_pt.time_point, &eo);
                ret = equationOrigins(from_pt.time_point, &eo);
                if (!ret) {
                    ret = timepointToAppSideral(from_pt.time_point, &lst, true);
                }
            };

            if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_APP) {
                // first, compute HA from CIO-based RA!!!
                lst_eo();
                if (!ret) {
                    // ha = MccAngle(lst - from_pt.X + eo).normalize<MccAngle::NORM_KIND_0_360>();
                    ha = MccAngle(lst - from_pt.X - eo).normalize<MccAngle::NORM_KIND_180_180>();
                    // ha = MccAngle(lst - from_pt.X - eo).normalize<MccAngle::NORM_KIND_0_360>();
                } else {
                    return ret;
                }

                if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZZD) {
                    from_pt.X = ha;
                    hadec2azalt(from_pt, to_pt);
                    to_pt->Y = PI_2 - to_pt->Y;
                } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZALT) {
                    from_pt.X = ha;
                    hadec2azalt(from_pt, to_pt);
                } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_HADEC_APP) {
                    to_pt->X = ha;
                    to_pt->Y = from_pt.Y;
                } else {
                    ret = MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
                }
            } else if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_HADEC_APP) {
                if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZZD) {
                    hadec2azalt(from_pt, to_pt);
                    to_pt->Y = PI_2 - to_pt->Y;
                } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZALT) {
                    hadec2azalt(from_pt, to_pt);
                } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_APP) {
                    lst_eo();
                    if (!ret) {
                        // to_pt->X = MccAngle(lst - from_pt.X + eo).normalize<MccAngle::NORM_KIND_0_360>();
                        to_pt->X = MccAngle(lst - from_pt.X - eo).normalize<MccAngle::NORM_KIND_0_360>();
                        to_pt->Y = from_pt.Y;
                    } else {
                        return ret;
                    }
                } else {
                    ret = MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
                }

            } else if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_AZZD) {
                from_pt.Y = PI_2 - from_pt.Y;
                from_pt.pair_kind = MccCoordPairKind::COORDS_KIND_AZALT;

                ret = transformCoordinates(std::move(from_pt), to_pt);
            } else if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_AZALT) {
                if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_HADEC_APP) {
                    azalt2hadec(from_pt, to_pt);
                } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_AZZD) {
                    to_pt->X = from_pt.X;
                    to_pt->Y = PI_2 - from_pt.Y;
                } else if (to_pt->pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_APP) {
                    azalt2hadec(from_pt, to_pt);
                    lst_eo();
                    if (!ret) {
                        // to_pt->X = MccAngle(lst - to_pt->X + eo).normalize<MccAngle::NORM_KIND_0_360>();
                        to_pt->X = MccAngle(lst - to_pt->X - eo).normalize<MccAngle::NORM_KIND_0_360>();
                    }
                } else {
                    ret = MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
                }

            } else {
                ret = MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
            }
        }

        return ret;
    }


    error_t transformCoordinatesEQHR(mcc_celestial_point_c auto from_pt, mcc_eqt_hrz_coord_c auto* to_pt)
    {
        error_t ret = MccCCTE_ERFAErrorCode::ERROR_OK;

        // MccJulianDay jd;

        if (to_pt == nullptr) {
            return ret;
        }

        using tp_pt_t = std::remove_cvref_t<decltype(*to_pt)>;
        MccGenericCelestialPoint<typename tp_pt_t::coord_t> cpt;
        cpt.time_point = from_pt.time_point;

        // the main scenario: from ICRS to apparent
        if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_ICRS) {
            return icrs2obs(from_pt, to_pt);
        }

        // from apparent: copy corresponded coordinates and compute other ones (ignore to_pt->pair_kind)
        if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_APP) {
            to_pt->RA_APP = from_pt.X;
            to_pt->DEC_APP = from_pt.Y;

            cpt.pair_kind = MccCoordPairKind::COORDS_KIND_HADEC_APP;

            ret = transformCoordinates(from_pt, &cpt);
            if (!ret) {
                to_pt->HA = cpt.X;

                auto cpt1 = cpt;
                cpt1.pair_kind = MccCoordPairKind::COORDS_KIND_AZALT;

                ret = transformCoordinates(cpt, &cpt1);
                if (!ret) {
                    to_pt->AZ = cpt1.X;
                    to_pt->ALT = cpt1.Y;
                    to_pt->ZD = PI_2 - cpt1.Y;
                }
            }
        } else if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_HADEC_APP) {
            to_pt->HA = from_pt.X;
            to_pt->DEC_APP = from_pt.Y;

            cpt.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_APP;
            ret = transformCoordinates(from_pt, &cpt);
            if (!ret) {
                to_pt->RA_APP = cpt.X;

                cpt.pair_kind = MccCoordPairKind::COORDS_KIND_AZALT;
                ret = transformCoordinates(from_pt, &cpt);
                if (!ret) {
                    to_pt->AZ = cpt.X;
                    to_pt->ALT = cpt.Y;
                    to_pt->ZD = PI_2 - cpt.Y;
                }
            }
        } else if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_AZZD) {
            to_pt->AZ = from_pt.X;
            to_pt->ZD = from_pt.Y;
            to_pt->ALT = PI_2 - from_pt.Y;

            cpt.pair_kind = MccCoordPairKind::COORDS_KIND_HADEC_APP;
            ret = transformCoordinates(from_pt, &cpt);
            if (!ret) {
                to_pt->HA = cpt.X;
                to_pt->DEC_APP = cpt.Y;

                auto cpt1 = cpt;
                cpt1.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_APP;
                ret = transformCoordinates(cpt, &cpt1);
                if (!ret) {
                    to_pt->RA_APP = cpt1.X;
                }
            }
        } else if (from_pt.pair_kind == MccCoordPairKind::COORDS_KIND_AZALT) {
            to_pt->AZ = from_pt.X;
            to_pt->ALT = from_pt.Y;
            to_pt->ZD = PI_2 - from_pt.Y;

            cpt.pair_kind = MccCoordPairKind::COORDS_KIND_HADEC_APP;
            ret = transformCoordinates(from_pt, &cpt);
            if (!ret) {
                to_pt->HA = cpt.X;
                to_pt->DEC_APP = cpt.Y;

                auto cpt1 = cpt;
                cpt1.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_APP;
                ret = transformCoordinates(cpt, &cpt1);
                if (!ret) {
                    to_pt->RA_APP = cpt1.X;
                }
            }
        } else {
            ret = MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
        }

        return ret;
    }

    // equation of the origins

    error_t equationOrigins(mcc_julday_c auto const& jd, mcc_angle_c auto* eo)
    {
        if (eo == nullptr) {
            return MccCCTE_ERFAErrorCode::ERROR_OK;
        }

        auto ret = MccCCTE_ERFAErrorCode::ERROR_OK;

        double tt, mjd;

        if constexpr (std::floating_point<decltype(jd)>) {
            mjd = jd - ERFA_DJM0;
        } else {
            mjd = jd.MJD();
        }

        std::lock_guard lock{*_stateMutex};

        using real_days_t = std::chrono::duration<double, std::ratio<86400>>;

        auto tai_utc = _currentState._leapSeconds[mjd];
        if (tai_utc.has_value()) {
            tt = jd.MJD() + std::chrono::duration_cast<real_days_t>(tai_utc.value()).count();

            auto tt_tai = _currentState._bulletinA.TT_TAI();
            tt += +std::chrono::duration_cast<real_days_t>(tt_tai).count();

            *eo = eraEo06a(ERFA_DJM0, tt);
        } else {
            ret = MccCCTE_ERFAErrorCode::ERROR_LEAPSECONDS_OUT_OF_RANGE;
        }

        return ret;
    }

    error_t equationOrigins(mcc_time_point_c auto const& tp, mcc_angle_c auto* eo)
    {
        if (eo == nullptr) {
            return MccCCTE_ERFAErrorCode::ERROR_OK;
        }

        MccJulianDay jd;
        auto ret = timepointToJulday(tp, &jd);
        if (ret) {
            return ret;
        }

        return equationOrigins(std::move(jd), eo);
    }


    // refraction

    error_t refractionModel(refract_model_t* model)
    {
        std::lock_guard lock{*_stateMutex};

        eraRefco(_currentState.meteo.pressure, _currentState.meteo.temperature, _currentState.meteo.humidity,
                 _currentState.wavelength, &model->refa, &model->refb);

        return MccCCTE_ERFAErrorCode::ERROR_OK;
    }

    error_t refractionCorrection(mcc_celestial_point_c auto pt, mcc_angle_c auto* dZ)
    {
        error_t ret = MccCCTE_ERFAErrorCode::ERROR_OK;

        if (dZ == nullptr) {
            return ret;
        }

        refract_model_t rmodel;

        refractionModel(&rmodel);

        if (pt.pair_kind == MccCoordPairKind::COORDS_KIND_AZZD) {
            if (pt.Y >= std::numbers::pi / 2.0) {
                *dZ = 35.4 / 60.0 * std::numbers::pi / 180.0;  // 35.4 arcminutes
            } else {
                auto tanZ = tan(pt.Y);
                *dZ = rmodel.refa * tanZ + rmodel.refb * tanZ * tanZ * tanZ;
                // *dZ = rmodel.refa / tanZ + rmodel.refb / (tanZ * tanZ * tanZ);
            }
        } else {
            MccCelestialPoint tr_pt{.pair_kind = MccCoordPairKind::COORDS_KIND_AZZD, .time_point = pt.time_point};
            ret = transformCoordinates(std::move(pt), &tr_pt);
            if (!ret) {
                ret = refractionCorrection(std::move(tr_pt), dZ);
            }
        }

        return ret;
    }

    /*  helper mathods  */

    auto leapSecondsExpireDate() const
    {
        return _currentState._leapSeconds.expireDate();
    }

    auto leapSecondsExpireMJD() const
    {
        return _currentState._leapSeconds.expireMJD();
    }


    auto bulletinADateRange() const
    {
        return _currentState._bulletinA.dateRange();
    }

    auto bulletinADateRangeMJD() const
    {
        return _currentState._bulletinA.dateRangeMJD();
    }


protected:
    engine_state_t _currentState{};

    std::unique_ptr<std::mutex> _stateMutex;

    error_t icrs2obs(mcc_celestial_point_c auto const& pt, mcc_eqt_hrz_coord_c auto* result)
    {
        if (result == nullptr) {
            return MccCCTE_ERFAErrorCode::ERROR_OK;
        }

        MccJulianDay jd;

        auto ret = timepointToJulday(result->time_point, &jd);
        if (ret) {
            return ret;
        }

        std::lock_guard lock{*_stateMutex};

        auto dut1 = _currentState._bulletinA.DUT1(jd.MJD());

        if (!dut1.has_value()) {
            return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
        }

        auto pol_pos = _currentState._bulletinA.polarCoords(jd.MJD());
        if (!pol_pos.has_value()) {
            return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
        }

        // const auto arcsec2rad = std::numbers::pi / 180 / 3600;
        const auto arcsec2rad = 1.0_arcsecs;
        pol_pos->x *= arcsec2rad;
        pol_pos->y *= arcsec2rad;

        double az, zd, ha, ra, dec, eo;

        int err = eraAtco13(pt.X, pt.Y, 0.0, 0.0, 0.0, 0.0, jd.MJD0, jd.MJD(), dut1->count(), _currentState.lon,
                            _currentState.lat, _currentState.elev, pol_pos->x, pol_pos->y, _currentState.meteo.pressure,
                            _currentState.meteo.temperature, _currentState.meteo.humidity, _currentState.wavelength,
                            &az, &zd, &ha, &dec, &ra, &eo);

        if (err == 1) {
            ret = MccCCTE_ERFAErrorCode::ERROR_DUBIOUS_YEAR;
        } else if (err == -1) {
            ret = MccCCTE_ERFAErrorCode::ERROR_UNACCEPTABLE_DATE;
        }


        // NOTE: according to definition of astronomical azimuth it is counted from the South through the West, but
        // in the ERFA the azimuth is counted from the North through the East!!!
        //
        result->AZ = MccAngle(az + std::numbers::pi).normalize<MccAngle::NORM_KIND_0_360>();
        result->ZD = zd;
        result->HA = ha;
        result->RA_APP = ra;
        result->DEC_APP = dec;
        result->ALT = MccCCTE_ERFA::PI_2 - zd;

        return ret;
    }


    error_t obs2icrs(mcc_celestial_point_c auto from_pt, mcc_celestial_point_c auto* to_pt)
    {
        if (to_pt == nullptr) {
            return MccCCTE_ERFAErrorCode::ERROR_OK;
        }

        MccJulianDay jd;

        auto ret = timepointToJulday(from_pt.time_point, &jd);
        if (ret) {
            return ret;
        }

        std::lock_guard lock{*_stateMutex};

        auto dut1 = _currentState._bulletinA.DUT1(jd.MJD());

        if (!dut1.has_value()) {
            return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
        }

        auto pol_pos = _currentState._bulletinA.polarCoords(jd.MJD());
        if (!pol_pos.has_value()) {
            return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
        }

        // const auto arcsec2rad = std::numbers::pi / 180 / 3600;
        const auto arcsec2rad = 1.0_arcsecs;
        pol_pos->x *= arcsec2rad;
        pol_pos->y *= arcsec2rad;

        std::string type;
        switch (from_pt.pair_kind) {
            case mcc::MccCoordPairKind::COORDS_KIND_AZZD:
                // NOTE: according to definition of astronomical azimuth it is counted from the South through the West,
                // but in the ERFA the azimuth is counted from the North through the East!!!
                //
                from_pt.X += std::numbers::pi;
                type = "A";
                break;
            case mcc::MccCoordPairKind::COORDS_KIND_AZALT:
                // NOTE: according to definition of astronomical azimuth it is counted from the South through the West,
                // but in the ERFA the azimuth is counted from the North through the East!!!
                //
                from_pt.X += std::numbers::pi;
                from_pt.Y = MccCCTE_ERFA::PI_2 - from_pt.Y;  // altitude to zenithal distance
                type = "A";
                break;
            case mcc::MccCoordPairKind::COORDS_KIND_HADEC_APP:
                type = "H";
                break;
            case mcc::MccCoordPairKind::COORDS_KIND_RADEC_APP:
                type = "R";
                break;
            default:
                return MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
        }

        int err = eraAtoc13(type.c_str(), from_pt.X, from_pt.Y, jd.MJD0, jd.MJD(), dut1->count(), _currentState.lon,
                            _currentState.lat, _currentState.elev, pol_pos->x, pol_pos->y, _currentState.meteo.pressure,
                            _currentState.meteo.temperature, _currentState.meteo.humidity, _currentState.wavelength,
                            &to_pt->X, &to_pt->Y);

        if (err == 1) {
            ret = MccCCTE_ERFAErrorCode::ERROR_DUBIOUS_YEAR;
        } else if (err == -1) {
            ret = MccCCTE_ERFAErrorCode::ERROR_UNACCEPTABLE_DATE;
        }

        return ret;
    }


    void azalt2hadec(mcc_celestial_point_c auto const& from_pt, mcc_celestial_point_c auto* to_pt)
    {
        double ha, dec;

        // NOTE: according to definition of astronomical azimuth it is counted from the South through the West, but
        // in the ERFA the azimuth is counted from the North through the East!!!
        //
        eraAe2hd(from_pt.X + std::numbers::pi, from_pt.Y, _currentState.lat, &ha, &dec);

        to_pt->X = ha;
        to_pt->Y = dec;
    }


    void hadec2azalt(mcc_celestial_point_c auto const& from_pt, mcc_celestial_point_c auto* to_pt)
    {
        double az, alt;

        eraHd2ae(from_pt.X, from_pt.Y, _currentState.lat, &az, &alt);

        // NOTE: according to definition of astronomical azimuth it is counted from the South through the West, but
        // in the ERFA the azimuth is counted from the North through the East!!!
        //
        // convert it to "from the South through the West"
        to_pt->X = MccAngle(az - std::numbers::pi).normalize<MccAngle::NORM_KIND_0_360>();
        to_pt->Y = alt;
    }

    error_t eqOrigins(mcc_time_point_c auto const& tp, double* eo)
    {
        if (eo == nullptr) {
            return MccCCTE_ERFAErrorCode::ERROR_OK;
        }

        MccJulianDay jd;
        double tt;
        auto ret = timepointToJulday(tp, &jd);
        if (ret) {
            return ret;
        }

        std::lock_guard lock{*_stateMutex};

        using real_days_t = std::chrono::duration<double, std::ratio<86400>>;

        auto tai_utc = _currentState._leapSeconds[jd.MJD()];
        if (tai_utc.has_value()) {
            tt = jd.MJD() + std::chrono::duration_cast<real_days_t>(tai_utc.value()).count();

            auto tt_tai = _currentState._bulletinA.TT_TAI();
            tt += +std::chrono::duration_cast<real_days_t>(tt_tai).count();

            *eo = eraEo06a(jd.MJD0, tt);
        } else {
            ret = MccCCTE_ERFAErrorCode::ERROR_LEAPSECONDS_OUT_OF_RANGE;
        }

        return ret;
    }
};

static_assert(mcc_ccte_c<MccCCTE_ERFA>, "");

}  // namespace mcc::ccte::erfa