...

mcc/mcc_ccte_erfa.h

@@ -299,12 +299,21 @@ public:
 
         using real_days_t = std::chrono::duration<double, std::ratio<86400>>;
 
-        double ut1 = jd.MJD();
+        using jd_t = decltype(jd);
-        double tt = jd.MJD();
+
+        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(jd.MJD());
+        auto dut1 = _currentState._bulletinA.DUT1(mjd);
 
         if (dut1.has_value()) {
             ut1 += std::chrono::duration_cast<real_days_t>(dut1.value()).count();
@@ -312,7 +321,7 @@ public:
             return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
         }
 
-        auto tai_utc = _currentState._leapSeconds[jd.MJD()];
+        auto tai_utc = _currentState._leapSeconds[mjd];
         if (tai_utc.has_value()) {
             tt += std::chrono::duration_cast<real_days_t>(tai_utc.value()).count();
         } else {
@@ -471,7 +480,8 @@ public:
             double eo, lst, ha, dec, az, alt;
 
             auto lst_eo = [&]() {
-                ret = eqOrigins(from_pt.time_point, &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);
                 }
@@ -566,7 +576,7 @@ public:
             return icrs2obs(from_pt, to_pt);
         }
 
-        // from apparent: copy corresponded coordinates and compute other ones
+        // 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;
@@ -647,6 +657,58 @@ public:
         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
 

mcc/mcc_generics.h

@@ -85,11 +85,19 @@ concept mcc_fp_type_like_c =
 
 template <typename T>
 concept mcc_angle_c = mcc_fp_type_like_c<T> && requires(T v, double vd) {
+    // mandatory arithmetic operations
     { v + v } -> std::same_as<T>;
     { v - v } -> std::same_as<T>;
     { v += v } -> std::same_as<T&>;
     { v -= v } -> std::same_as<T&>;
 
+    { vd + v } -> std::same_as<T>;
+    { vd - v } -> std::same_as<T>;
+    { v + vd } -> std::same_as<T>;
+    { v - vd } -> std::same_as<T>;
+    { v += vd } -> std::same_as<T&>;
+    { v -= vd } -> std::same_as<T&>;
+
     { v * vd } -> std::same_as<T>;
     { v / vd } -> std::same_as<T>;
 };
@@ -101,7 +109,8 @@ concept mcc_angle_c = mcc_fp_type_like_c<T> && requires(T v, double vd) {
  *  USE OF STL std::chrono::time_point
  */
 template <typename T>
-concept mcc_time_point_c = requires(T t) { []<typename CT, typename DT>(std::chrono::time_point<CT, DT>) {}(t); };
+concept mcc_time_point_c = traits::mcc_systime_c<T>;
+// concept mcc_time_point_c = requires(T t) { []<typename CT, typename DT>(std::chrono::time_point<CT, DT>) {}(t); };
 
 
 /*                          JULIAN DAY CLASS CONCEPT                            */
@@ -123,6 +132,16 @@ concept mcc_error_c = std::convertible_to<T, bool> || requires(const T t) {
 };
 
 
+template <mcc_error_c ErrT, typename DErrT>
+static constexpr ErrT mcc_deduce_error(const DErrT& err, const ErrT& default_err)
+{
+    if constexpr (std::same_as<ErrT, DErrT>) {
+        return err;
+    } else {
+        return default_err;
+    }
+}
+
 /*                          ATMOSPHERIC REFRACTION MODEL CLASS CONCEPT                            */
 
 template <typename T>
@@ -144,6 +163,27 @@ concept mcc_celestial_point_c = requires(T t) {
 };
 
 
+static constexpr void mcc_copy_celestial_point(mcc_celestial_point_c auto const& from_pt,
+                                               mcc_celestial_point_c auto* to_pt)
+{
+    if (to_pt == nullptr) {
+        return;
+    }
+
+    using from_pt_t = std::remove_cvref_t<decltype(from_pt)>;
+    using to_pt_t = std::remove_cvref_t<decltype(*to_pt)>;
+
+    if constexpr (std::derived_from<to_pt_t, from_pt_t> && std::copyable<to_pt_t>) {
+        *to_pt = from_pt;
+        return;
+    }
+
+    to_pt->pair_kind = from_pt.pair_kind;
+    to_pt->time_point = std::chrono::time_point_cast<typename decltype(to_pt->time_point)::duration>(from_pt);
+    to_pt->X = (double)from_pt.X;
+    to_pt->Y = (double)from_pt.Y;
+}
+
 
 /*                          CELESTIAL POINT WITH APPARENT EQUATORIAL AND HORIZONTAL CLASS CONCEPT */
 
@@ -159,6 +199,34 @@ concept mcc_eqt_hrz_coord_c = mcc_celestial_point_c<T> && requires(T t) {
 };
 
 
+static constexpr void mcc_copy_eqt_hrz_coord(mcc_eqt_hrz_coord_c auto const& from_pt, mcc_eqt_hrz_coord_c auto* to_pt)
+{
+    if (to_pt == nullptr) {
+        return;
+    }
+
+    using from_pt_t = std::remove_cvref_t<decltype(from_pt)>;
+    using to_pt_t = std::remove_cvref_t<decltype(*to_pt)>;
+
+    if constexpr (std::derived_from<to_pt_t, from_pt_t> && std::copyable<to_pt_t>) {
+        *to_pt = from_pt;
+        return;
+    }
+
+    to_pt->pair_kind = from_pt.pair_kind;
+    to_pt->time_point = std::chrono::time_point_cast<typename decltype(to_pt->time_point)::duration>(from_pt);
+    to_pt->X = (double)from_pt.X;
+    to_pt->Y = (double)from_pt.Y;
+
+    to_pt->RA_APP = (double)from_pt.RA_APP;
+    to_pt->DEC_APP = (double)from_pt.DEC_APP;
+
+    to_pt->HA = (double)from_pt.HA;
+    to_pt->AZ = (double)from_pt.AZ;
+    to_pt->ZD = (double)from_pt.ZD;
+    to_pt->ALT = (double)from_pt.ALT;
+}
+
 
 /*                          CELESTIAL COORDINATES TRANSFORMATION ENGINE     */
 
@@ -180,6 +248,12 @@ struct mcc_CCTE_interface_t {
         return std::forward<SelfT>(self).timepointToAppSideral(std::move(tp), st, islocal);
     }
 
+    template <std::derived_from<mcc_CCTE_interface_t> SelfT>
+    RetT juldayToAppSideral(this SelfT&& self, mcc_julday_c auto jd, mcc_angle_c auto* st, bool islocal = false)
+    {
+        return std::forward<SelfT>(self).timepointToAppSideral(std::move(jd), st, islocal);
+    }
+
     // NOTE: ASSUMING THE AZINUTH IS COUNTED FROM THE SOUTH THROUGH THE WEST!!!
     template <std::derived_from<mcc_CCTE_interface_t> SelfT>
     RetT transformCoordinates(this SelfT&& self, mcc_celestial_point_c auto from_pt, mcc_celestial_point_c auto* to_pt)
@@ -230,8 +304,43 @@ concept mcc_pointing_target_coord_c = mcc_eqt_hrz_coord_c<T> && requires(T t) {
     requires mcc_angle_c<decltype(t.DEC_ICRS)>;  // ICRS declination
 };
 
+
+static constexpr void mcc_copy_pointing_target_coord(mcc_pointing_target_coord_c auto const& from_pt,
+                                                     mcc_pointing_target_coord_c auto* to_pt)
+{
+    if (to_pt == nullptr) {
+        return;
+    }
+
+    using from_pt_t = std::remove_cvref_t<decltype(from_pt)>;
+    using to_pt_t = std::remove_cvref_t<decltype(*to_pt)>;
+
+    if constexpr (std::derived_from<to_pt_t, from_pt_t> && std::copyable<to_pt_t>) {
+        *to_pt = from_pt;
+        return;
+    }
+
+    to_pt->pair_kind = from_pt.pair_kind;
+    to_pt->time_point = std::chrono::time_point_cast<typename decltype(to_pt->time_point)::duration>(from_pt);
+    to_pt->X = (double)from_pt.X;
+    to_pt->Y = (double)from_pt.Y;
+
+    to_pt->RA_ICRS = (double)from_pt.RA_ICRS;
+    to_pt->DEC_ICRS = (double)from_pt.DEC_ICRS;
+
+    to_pt->RA_APP = (double)from_pt.RA_APP;
+    to_pt->DEC_APP = (double)from_pt.DEC_APP;
+
+    to_pt->HA = (double)from_pt.HA;
+    to_pt->AZ = (double)from_pt.AZ;
+    to_pt->ZD = (double)from_pt.ZD;
+    to_pt->ALT = (double)from_pt.ALT;
+}
+
+
+
 template <typename T>
-concept mcc_telemetry_data_c = mcc_eqt_hrz_coord_c<T> && requires(T t) {
+concept mcc_telemetry_data_c = mcc_eqt_hrz_coord_c<T> && std::default_initializable<T> && requires(T t) {
     // target target coordinates
     requires mcc_pointing_target_coord_c<decltype(t.target)>;
 
@@ -249,6 +358,46 @@ concept mcc_telemetry_data_c = mcc_eqt_hrz_coord_c<T> && requires(T t) {
 };
 
 
+static constexpr void mcc_copy_telemetry_data(mcc_telemetry_data_c auto const& from_pt,
+                                              mcc_telemetry_data_c auto* to_pt)
+{
+    if (to_pt == nullptr) {
+        return;
+    }
+
+    using from_pt_t = std::remove_cvref_t<decltype(from_pt)>;
+    using to_pt_t = std::remove_cvref_t<decltype(*to_pt)>;
+
+    if constexpr (std::derived_from<to_pt_t, from_pt_t> && std::copyable<to_pt_t>) {
+        *to_pt = from_pt;
+        return;
+    }
+
+    to_pt->pair_kind = from_pt.pair_kind;
+    to_pt->time_point = std::chrono::time_point_cast<typename decltype(to_pt->time_point)::duration>(from_pt);
+
+    to_pt->X = (double)from_pt.X;
+    to_pt->Y = (double)from_pt.Y;
+    to_pt->speedX = (double)from_pt.speedX;
+    to_pt->speedY = (double)from_pt.speedY;
+
+    to_pt->RA_APP = (double)from_pt.RA_APP;
+    to_pt->DEC_APP = (double)from_pt.DEC_APP;
+
+    to_pt->HA = (double)from_pt.HA;
+    to_pt->AZ = (double)from_pt.AZ;
+    to_pt->ZD = (double)from_pt.ZD;
+    to_pt->ALT = (double)from_pt.ALT;
+
+    to_pt->pcmX = (double)from_pt.pcmX;
+    to_pt->pcmY = (double)from_pt.pcmY;
+
+    to_pt->refCorr = (double)from_pt.refCorr;
+
+    mcc_copy_pointing_target_coord(from_pt.target, &to_pt->target);
+}
+
+
 
 /*                          MOUNT TELEMETRY MANAGER CLASS CONCEPT                   */
 

mcc/mcc_telemetry.h

@@ -14,7 +14,13 @@ namespace mcc
 static constexpr double mcc_sideral_to_UT1_ratio = 1.002737909350795;  // sideral/UT1
 
 
-enum MccTelemetryErrorCode : int { ERROR_OK, ERROR_COORD_TRANSFROM, ERROR_PCM_COMP };
+enum MccTelemetryErrorCode : int {
+    ERROR_OK,
+    ERROR_NULLPTR,
+    ERROR_COORD_TRANSFORM,
+    ERROR_PCM_COMP,
+    ERROR_HARDWARE_GETPOS
+};
 
 }  // namespace mcc
 
@@ -50,10 +56,14 @@ struct MccTelemetryCategory : public std::error_category {
         switch (err) {
             case MccTelemetryErrorCode::ERROR_OK:
                 return "OK";
-            case MccTelemetryErrorCode::ERROR_COORD_TRANSFROM:
+            case MccTelemetryErrorCode::ERROR_NULLPTR:
+                return "nullptr input argument";
+            case MccTelemetryErrorCode::ERROR_COORD_TRANSFORM:
                 return "coordinate transformation error";
             case MccTelemetryErrorCode::ERROR_PCM_COMP:
                 return "PCM computation error";
+            case MccTelemetryErrorCode::ERROR_HARDWARE_GETPOS:
+                return "cannot get hardware position";
             default:
                 return "UNKNOWN";
         }
@@ -87,24 +97,54 @@ public:
     typedef std::error_code error_t;
 
 
-    MccTelemetry(mcc_ccte_c auto* ccte, mcc_PCM_c auto* pcm, mcc_hardware_c auto* hardware)
+    MccTelemetry(mcc_ccte_c auto* ccte, mcc_PCM_c auto* pcm, mcc_hardware_c auto* hardware) : _data()
     {
+        _data.target.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_ICRS;
+
         using ccte_t = std::remove_cvref_t<decltype(*ccte)>;
         using pcm_t = std::remove_cvref_t<decltype(*pcm)>;
         using hardware_t = std::remove_cvref_t<decltype(*hardware)>;
 
+
+        _updateTargetFunc = [ccte, this]() {
+            //
+            // ICRS coordinates of the taget must be already set
+            //
+            _data.target.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_ICRS;
+            _data.target.X = _data.target.RA_ICRS;
+            _data.target.Y = _data.target.DEC_ICRS;
+
+            // update apparent cordinates
+            auto ret = ccte->transformCoordinates(_data.target, &_data.target);
+
+            return mcc_deduce_error<error_t>(ret, MccTelemetryErrorCode::ERROR_COORD_TRANSFORM);
+        };
+
         _updateFunc = [ccte, pcm, hardware, this](MccTelemetryData* data) {
+            // first, update mount quantities
             typename hardware_t::axes_pos_t hw_pos;
             auto hw_err = hardware->getPos(&hw_pos);
+            if (hw_err) {
+                return mcc_deduce_error(hw_err, MccTelemetryErrorCode::ERROR_HARDWARE_GETPOS);
+            }
+
+            double eo;
 
-            if (!hw_err) {
             data->time_point =
                 std::chrono::time_point_cast<typename decltype(data->time_point)::duration>(hw_pos.time_point);
 
             auto ccte_err = ccte->timepointToJulday(data->time_point, &data->JD);
             if (!ccte_err) {
-                    ccte_err = ccte->timepointToAppSideral(data->time_point, &data->LST, true);
+                ccte_err = ccte->juldayToAppSideral(data->JD, &data->LST, true);
                 if (!ccte_err) {
+                    ccte_err = ccte->equationOrigins(data->JD, &eo);
+                }
+            }
+
+            if (ccte_err) {
+                return mcc_deduce_error(ccte_err, MccTelemetryErrorCode::ERROR_COORD_TRANSFORM);
+            }
+
             data->speedX = (double)hw_pos.speedX;
             data->speedY = (double)hw_pos.speedY;
 
@@ -114,19 +154,18 @@ public:
 
             auto pcm_err = pcm->computePCM(data, &pcm_res);
             if (pcm_err) {
-                            if constexpr (std::same_as<decltype(pcm_err), error_t>) {
+                return mcc_deduce_error(pcm_err, MccTelemetryErrorCode::ERROR_PCM_COMP);
-                                return pcm_err;
-                            } else {
-                                return MccTelemetryErrorCode::ERROR_PCM_COMP;
-                            }
             }
 
-                        MccCelestialPoint pt{.pair_kind = MccCoordPairKind::COORDS_KIND_AZALT,
+            data->pcmX = pcm_res.dx;
-                                             .time_point = data->time_point};
+            data->pcmY = pcm_res.dy;
+
+            MccCelestialPoint pt{.pair_kind = MccCoordPairKind::COORDS_KIND_AZALT, .time_point = data->time_point};
 
             if constexpr (mccIsEquatorialMount(pcm_t::mountType)) {
                 data->HA = (double)hw_pos.X + pcm_res.dx;
                 data->DEC_APP = (double)hw_pos.Y + pcm_res.dy;
+                data->RA_APP = (double)data->LST - (double)data->HA + eo;
 
                 data->X = data->HA;
                 data->Y = data->DEC_APP;
@@ -154,34 +193,116 @@ public:
                 if (!ccte) {
                     data->HA = pt.X;
                     data->DEC_APP = pt.Y;
+                    data->RA_APP = (double)data->LST - (double)data->HA + eo;
                 }
 
             } else {
                 static_assert(false, "UNKNOWN MOUNT TYPE!");
             }
+
+            if (!ccte_err) {
+                data->pair_kind = MccCoordPairKind::COORDS_KIND_AZZD;
+                data->X = data->AZ;
+                data->Y = data->ZD;
+
+                ccte_err = ccte->refractionCorrection(data, &data->refCorr);
+                if (!ccte_err) {
+                    // hardware encoders coordinates
+                    data->X = (double)hw_pos.X;
+                    data->Y = (double)hw_pos.Y;
+
+                    // update target (assuming target ICRS coordinates are already set)
+
+                    data->target.time_point =
+                        std::chrono::time_point_cast<typename decltype(data->target.time_point)::duration>(
+                            data->time_point);
+
+                    data->target.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_ICRS;
+                    data->target.X = data->target.RA_ICRS;
+                    data->target.Y = data->target.DEC_ICRS;
+
+                    ccte_err = ccte->transformCoordinates(data->target, &data->target);
                 }
             }
+
+            if (ccte_err) {
+                return mcc_deduce_error(ccte_err, MccTelemetryErrorCode::ERROR_COORD_TRANSFORM);
             }
+
+            if constexpr (mccIsEquatorialMount(pcm_t::mountType)) {
+                data->pair_kind = MccCoordPairKind::COORDS_KIND_HADEC_APP;
+            } else if constexpr (mccIsAltAzMount(pcm_t::mountType)) {
+                data->pair_kind = MccCoordPairKind::COORDS_KIND_AZALT;
+            } else {
+                static_assert(false, "UNKNOWN MOUNT TYPE!");
+            }
+
+
+            return MccTelemetryErrorCode::ERROR_OK;
+        };
+
+
+        _setTargetFunc = [ccte, this](MccCelestialPoint const& pt) {
+            // in the case of apparent input coordinates
+            // one must ensure the same time points
+
+            _data.target.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_ICRS;
+            _data.target.time_point =
+                std::chrono::time_point_cast<typename decltype(_data.target.time_point)::duration>(pt.time_point);
+
+            auto ret = ccte->transformCoordinates(pt, &_data.target);
+
+            if (!ret) {
+                if (pt.pair_kind == MccCoordPairKind::COORDS_KIND_RADEC_ICRS) {
+                    _data.target.RA_ICRS = _data.target.X;
+                    _data.target.DEC_ICRS = _data.target.Y;
+
+                    // update apparent coordinates
+                    ret = _updateTargetFunc();
+                } else {  // apparent coordinates were computed above
+                    // compute ICRS coordinates
+                    MccCelestialPoint cpt{.pair_kind = MccCoordPairKind::COORDS_KIND_RADEC_ICRS};
+                    ret = ccte->transformCoordinates(pt, &cpt);
+
+                    _data.target.RA_ICRS = cpt.X;
+                    _data.target.DEC_ICRS = cpt.Y;
+                }
+            }
+
+            return mcc_deduce_error<error_t>(ret, MccTelemetryErrorCode::ERROR_COORD_TRANSFORM);
         };
     }
 
 
-    error_t telemetryData(mcc_telemetry_data_c auto* data)
+    virtual ~MccTelemetry() = default;
+
+    error_t telemetryData(mcc_telemetry_data_c auto* tdata)
     {
-        error_t ret = MccTelemetryErrorCode::ERROR_OK;
+        if (tdata == nullptr) {
-
+            return MccTelemetryErrorCode::ERROR_NULLPTR;
-        // first, update mount coordinates
+        }
 
+        error_t ret = _updateFunc(&_data);
 
+        if (!ret) {
+            mcc_copy_telemetry_data(_data, tdata);
+        }
 
         return ret;
     }
 
 
-    error_t setPointingTarget(mcc_celestial_point_c auto pt) {}
+    error_t setPointingTarget(mcc_celestial_point_c auto pt)
+    {
+        // return
+    }
 
 protected:
-    std::function<error_t(MccTelemetryData*)> _updateFunc;
+    MccTelemetryData _data;
+
+    std ::function<error_t()> _updateTargetFunc{};
+    std::function<error_t(MccTelemetryData*)> _updateFunc{};
+    std::function<error_t(MccCelestialPoint const&)> _setTargetFunc{};
 };
 
 }  // namespace mcc

mcc/tests/ccte_test.cpp

@@ -47,9 +47,11 @@ int main()
 
     std::cout << "JD = " << std::setprecision(12) << (double)jd << "\n";
 
-    mcc::MccAngle lst;
+    mcc::MccAngle lst, eo;
     ret = erfa.timepointToAppSideral(now, &lst, true);
     std::cout << "LST = " << lst.sexagesimal(true) << "\n";
+    ret = erfa.equationOrigins(now, &eo);
+    std::cout << "EO = " << eo.sexagesimal(true) << "\n";
 
 
     std::cout << "\n\n";