...

(start developing)

CMakeLists.txt

@@ -215,6 +215,7 @@ set(MCC_SRC
     include/mcc/mcc_pcm.h
     include/mcc/mcc_telemetry.h
     include/mcc/mcc_movement_controls.h
+    include/mcc/mcc_generic_movecontrols.h
     include/mcc/mcc_serialization_common.h
     include/mcc/mcc_deserializer.h
     include/mcc/mcc_serializer.h
@@ -265,7 +266,10 @@ endif()
 if(USE_BSPLINE_PCM)
     target_compile_definitions(${PROJECT_NAME} INTERFACE USE_BSPLINE_PCM)
     target_link_libraries(${PROJECT_NAME} INTERFACE fitpack)
-    target_link_directories(${PROJECT_NAME} INTERFACE "$<BUILD_INTERFACE:${CMAKE_BINARY_DIR}/fitpack>")
+    target_link_directories(
+        ${PROJECT_NAME}
+        INTERFACE "$<BUILD_INTERFACE:${CMAKE_BINARY_DIR}/fitpack>"
+    )
 endif()
 
 if(USE_ASIO)
@@ -275,7 +279,10 @@ endif()
 if(USE_SPDLOG)
     target_link_libraries(${PROJECT_NAME} INTERFACE spdlog::spdlog_header_only)
     # target_compile_definitions(${PROJECT_NAME} INTERFACE SPDLOG_USE_STD_FORMAT=1 SPDLOG_FMT_EXTERNAL=0)
-    target_compile_definitions(${PROJECT_NAME} INTERFACE SPDLOG_USE_STD_FORMAT=1)
+    target_compile_definitions(
+        ${PROJECT_NAME}
+        INTERFACE SPDLOG_USE_STD_FORMAT=1
+    )
 endif()
 
 if(BUILD_TESTS)

include/mcc/mcc_bsplines.h

@@ -366,4 +366,111 @@ int fitpack_eval_spl2d(const TXT& tx,
     return fitpack_eval_spl2d(tx, ty, coeffs, xv, yv, fv, kx, ky);
 }
 
-}  // namespace mcc::fitpack
+
+/*  partial derivatives  */
+
+template <std::ranges::contiguous_range TXT,
+          std::ranges::contiguous_range TYT,
+          std::ranges::contiguous_range XT,
+          std::ranges::contiguous_range YT,
+          std::ranges::contiguous_range CoeffT,
+          std::ranges::contiguous_range PderT>
+int fitpack_parder_spl2d(const TXT& tx,
+                         const TYT& ty,
+                         const CoeffT& coeffs,
+                         const XT& x,
+                         const YT& y,
+                         PderT& pder,
+                         int dx,  // partial derivatives order along X
+                         int dy,  // partial derivatives order along Y
+                         int kx = 3,
+                         int ky = 3)
+{
+    static_assert(std::same_as<std::ranges::range_value_t<TXT>, double> &&
+                      std::same_as<std::ranges::range_value_t<TYT>, double> &&
+                      std::same_as<std::ranges::range_value_t<XT>, double> &&
+                      std::same_as<std::ranges::range_value_t<YT>, double> &&
+                      std::same_as<std::ranges::range_value_t<CoeffT>, double> &&
+                      std::same_as<std::ranges::range_value_t<PderT>, double>,
+                  "Input ranges elements type must be double!");
+
+    if (kx < 0 || ky < 0 || dx < 0 || dy < 0) {
+        return 10;
+    }
+
+    int ntx = std::ranges::size(tx);
+    int nty = std::ranges::size(ty);
+
+    auto n_coeffs = (ntx - kx - 1) * (nty - ky - 1);
+    if (std::ranges::size(coeffs) < n_coeffs) {
+        return 10;
+    }
+
+    int mx = std::ranges::size(x), my = std::ranges::size(y);
+    int N = mx * my;
+
+    if (std::ranges::size(pder) < N) {
+        std::ranges::fill_n(std::back_inserter(pder), N - std::ranges::size(pder), 0.0);
+    }
+
+    // mx >=1, my >=1, 0 <= nux < kx, 0 <= nuy < ky, kwrk>=mx+my
+    // lwrk>=mx*(kx+1-nux)+my*(ky+1-nuy)+(nx-kx-1)*(ny-ky-1),
+
+    // compute sizes of working arrays according parder.f
+    int lwrk = mx * (kx + 1 - dx) + my * (ky + 1 - dy) + (ntx - kx - 1) * (nty - ky - 1);
+    std::vector<double> wrk(lwrk);
+
+    int kwrk = mx + my;
+    std::vector<int> iwrk(kwrk);
+
+
+    auto tx_ptr = const_cast<double*>(std::ranges::data(tx));
+    auto ty_ptr = const_cast<double*>(std::ranges::data(ty));
+    auto coeffs_ptr = const_cast<double*>(std::ranges::data(coeffs));
+    auto x_ptr = const_cast<double*>(std::ranges::data(x));
+    auto y_ptr = const_cast<double*>(std::ranges::data(y));
+
+    int ier = 0;
+
+    parder(tx_ptr, &ntx, ty_ptr, &nty, coeffs_ptr, &kx, &ky, &dx, &dy, x_ptr, &mx, y_ptr, &my, std::ranges::data(pder),
+           wrk.data(), &lwrk, iwrk.data(), &kwrk, &ier);
+
+    return ier;
+}
+
+// scalar version
+template <std::ranges::contiguous_range TXT,
+          std::ranges::contiguous_range TYT,
+          typename XT,
+          typename YT,
+          std::ranges::contiguous_range CoeffT,
+          typename PderT>
+int fitpack_parder_spl2d(const TXT& tx,
+                         const TYT& ty,
+                         const CoeffT& coeffs,
+                         const XT& x,
+                         const YT& y,
+                         PderT& pder,
+                         int dx,  // partial derivatives order along X
+                         int dy,  // partial derivatives order along Y
+                         int kx = 3,
+                         int ky = 3)
+{
+    static_assert(std::same_as<std::ranges::range_value_t<TXT>, double> &&
+                      std::same_as<std::ranges::range_value_t<TYT>, double> &&
+                      std::same_as<std::ranges::range_value_t<CoeffT>, double>,
+                  "Input ranges elements type must be double!");
+
+    static_assert(
+        std::convertible_to<XT, double> && std::convertible_to<YT, double> && std::convertible_to<PderT, double>,
+        "XT, YT and FuncT types must be a scalar convertible to double!");
+
+    auto xv = std::vector<double>(1, x);
+    auto yv = std::vector<double>(1, y);
+    auto pv = std::vector<double>(1, pder);
+
+    return fitpack_parder_spl2d(tx, ty, coeffs, xv, yv, pv, dx, dy, kx, ky);
+}
+
+
+}  // namespace mcc::bsplines

include/mcc/mcc_concepts.h

@@ -499,21 +499,55 @@ struct mcc_skypoint_interface_t {
     {
         return std::forward<SelfT>(self).EO(eo);
     }
+
+    // template <std::derived_from<mcc_skypoint_interface_t> SelfT>
+    // typename std::remove_cvref_t<SelfT>::dist_result_t distance(this SelfT&& self, auto const& sp)
+    // {
+    //     return std::forward<SelfT>(self).distance(sp);
+    // }
+    // template <std::derived_from<mcc_skypoint_interface_t> SelfT>
+    // SelfT& operator+=(this SelfT& self, mcc_coord_pair_c auto const& dxy)
+    // {
+    //     return self.operator+=(dxy);
+    // }
+
+    // template <std::derived_from<mcc_skypoint_interface_t> SelfT>
+    // SelfT& operator-=(this SelfT& self, mcc_coord_pair_c auto const& dxy)
+    // {
+    //     return self.operator-=(dxy);
+    // }
 };
 
 template <typename T>
-concept mcc_skypoint_c = std::derived_from<T, mcc_skypoint_interface_t> && requires(const T t_const) {
-    { t_const.epoch() } -> mcc_coord_epoch_c;
+concept mcc_skypoint_c =
+    std::derived_from<T, mcc_skypoint_interface_t> && requires(const T t_const, const T t_other_const, T t) {
+        { t_const.epoch() } -> mcc_coord_epoch_c;
 
-    // currently stored coordinates pair kind
-    { t_const.pairKind() } -> std::same_as<impl::MccCoordPairKind>;
+        // currently stored coordinates pair kind
+        { t_const.pairKind() } -> std::same_as<impl::MccCoordPairKind>;
 
-    // currently stored co-longitude coordinate
-    { t_const.co_lon() } -> std::convertible_to<double>;
+        // currently stored co-longitude coordinate
+        { t_const.co_lon() } -> std::convertible_to<double>;
 
-    // currently stored co-latitude coordinate
-    { t_const.co_lat() } -> std::convertible_to<double>;
-};
+        // currently stored co-latitude coordinate
+        { t_const.co_lat() } -> std::convertible_to<double>;
+
+        requires requires(typename T::dist_result_t res) {
+            requires mcc_angle_c<decltype(res.dist)>;  // distance on sphere
+            requires mcc_angle_c<decltype(res.dx)>;    // difference along co-longitude coordinate
+            requires mcc_angle_c<decltype(res.dy)>;    // difference along co-latitude coordinate
+            requires mcc_angle_c<decltype(res.x2)>;    // co-longitude coordinates of target sky point (in the same
+                                                       // coordinate system as 'this')
+            requires mcc_angle_c<decltype(res.y2)>;    // co-latitude coordinates of target sky point (in the same
+                                                       // coordinate system as 'this')
+        };
+
+        // distance on sphere between two sky points
+        { t_const.distance(std::declval<const T&>()) } -> std::same_as<typename T::dist_result_t>;
+
+        // { t_const - t_other_const } -> mcc_coord_pair_c;
+        // { t_const + t_other_const } -> mcc_coord_pair_c;
+    };
 
 
 
@@ -672,75 +706,12 @@ concept mcc_hardware_c = requires(T t) {
     }();  // to ensure 'mountType' can be used in compile-time context
 
 
-    // a type that defines at least HW_MOVE_ERROR, HW_MOVE_STOPPING, HW_MOVE_STOPPED, HW_MOVE_SLEWING,
-    // HW_MOVE_ADJUSTING, HW_MOVE_TRACKING and HW_MOVE_GUIDING compile-time constants. The main purpose of this type is
-    // a possible tunning of hardware hardwareSetState-related commands and detect the stop and error states from
-    // hardware
-    //
-    // e.g. an implementations can be as follows:
-    //     enum class hardware_movement_state_t: int {HW_MOVE_ERROR = -1, HW_MOVE_STOPPED = 0, HW_MOVE_STOPPING,
-    //     HW_MOVE_SLEWING, HW_MOVE_ADJUSTING, HW_MOVE_TRACKING, HW_MOVE_GUIDING}
-    //
-    //     struct hardware_movement_state_t {
-    //         static constexpr uint16_t HW_MOVE_STOPPED = 0;
-    //         static constexpr uint16_t HW_MOVE_SLEWING = 111;
-    //         static constexpr uint16_t HW_MOVE_ADJUSTING = 222;
-    //         static constexpr uint16_t HW_MOVE_TRACKING = 333;
-    //         static constexpr uint16_t HW_MOVE_GUIDING = 444;
-    //         static constexpr uint16_t HW_MOVE_ERROR = 555;
-    //         static constexpr uint16_t HW_MOVE_STOPPING = 666;
-    //     }
-
     requires mcc_hardware_movement_state_c<typename T::hardware_movement_state_t>;
 
-    // requires requires(typename T::hardware_movement_state_t type) {
-    //     []() {
-    //         // mount axes were stopped
-    //         static constexpr auto v0 = T::hardware_movement_state_t::HW_MOVE_STOPPED;
-
-    //         // hardware was asked for slewing (move to given celestial point)
-    //         static constexpr auto v1 = T::hardware_movement_state_t::HW_MOVE_SLEWING;
-
-    //         // hardware was asked for adjusting after slewing
-    //         // (adjusting actual mount position to align with target celestial point at the end of slewing process)
-    //         static constexpr auto v2 = T::hardware_movement_state_t::HW_MOVE_ADJUSTING;
-
-    //         // hardware was asked for tracking (track target celestial point)
-    //         static constexpr auto v3 = T::hardware_movement_state_t::HW_MOVE_TRACKING;
-
-    //         // hardware was asked for guiding
-    //         // (small corrections to align actual mount position with target celestial point)
-    //         static constexpr auto v4 = T::hardware_movement_state_t::HW_MOVE_GUIDING;
-
-    //         // to detect possible hardware error
-    //         static constexpr auto v5 = T::hardware_movement_state_t::HW_MOVE_ERROR;
-    //     }();
-    // };
-
     requires mcc_hardware_state_c<typename T::hardware_state_t> && requires(typename T::hardware_state_t state) {
         requires std::same_as<decltype(state.movementState), typename T::hardware_movement_state_t>;
     };
 
-    // requires requires(typename T::hardware_state_t state) {
-    //     // encoder co-longitude and co-latiitude positions, as well as its measurement time point
-    //     // the given constrains on coordinate pair kind can be used to deduce mount type
-    //     requires mcc_coord_pair_c<decltype(state.XY)> &&
-    //                  (  // for equathorial mount:
-    //                      decltype(state.XY)::pairKind == impl::MccCoordPairKind::COORDS_KIND_HADEC_OBS ||
-    //                      // for alt-azimuthal mount:
-    //                      decltype(state.XY)::pairKind == impl::MccCoordPairKind::COORDS_KIND_AZALT ||
-    //                      decltype(state.XY)::pairKind == impl::MccCoordPairKind::COORDS_KIND_AZZD);
-
-    //     // co-longitude and co-latiitude axis angular speeds, as well as its measurement/computation time point
-    //     requires mcc_coord_pair_c<decltype(state.speedXY)> &&
-    //                  (decltype(state.XY)::pairKind == impl::MccCoordPairKind::COORDS_KIND_GENERIC ||
-    //                   decltype(state.XY)::pairKind == impl::MccCoordPairKind::COORDS_KIND_XY);
-
-
-    //     requires std::same_as<typename T::hardware_movement_state_t, decltype(state.movementState)>;
-    // };
-
-
     // set hardware state:
     { t.hardwareSetState(std::declval<typename T::hardware_state_t const&>()) } -> std::same_as<typename T::error_t>;
 
@@ -749,6 +720,8 @@ concept mcc_hardware_c = requires(T t) {
 
     // { t.hardwareStop() } -> std::same_as<typename T::error_t>;  // stop any moving
     { t.hardwareInit() } -> std::same_as<typename T::error_t>;  // initialize hardware
+
+    { t.hardwareShutdown() } -> std::same_as<typename T::error_t>;  // shutdown hardware
 };
 
 

include/mcc/mcc_coordinate.h

@@ -380,6 +380,14 @@ public:
 
     using error_t = typename CCTE_T::error_t;
 
+    struct dist_result_t {
+        MccAngle dist{};
+        MccAngle dx{};
+        MccAngle dy{};
+        MccAngle x2{};
+        MccAngle y2{};
+    };
+
     MccGenericSkyPoint() {}
 
     template <mcc_coord_pair_c PT>
@@ -685,11 +693,34 @@ public:
         return cctEngine.equationOrigins(_epoch, eo);
     }
 
+
+    dist_result_t distance(MccGenericSkyPoint const& sp) const
+    // dist_result_t distance(mcc_skypoint_c auto const& sp)
+    {
+        double x, y;
+
+        if (_pairKind == sp.pairKind() && utils::isEqual(_epoch.MJD(), sp.epoch().MJD())) {
+            x = sp.co_lon();
+            y = sp.co_lat();
+        } else {  // convert to the same coordinates kind
+            MccGenericSkyPoint p{sp};
+            p.to(_pairKind, _epoch);
+
+            x = p.co_lon();
+            y = p.co_lat();
+        }
+
+        auto d = utils::distanceOnSphere(_x, _y, x, y);
+
+        return {.dist = std::get<2>(d), .dx = std::get<0>(d), .dy = std::get<1>(d), .x2 = x, .y2 = y};
+    }
+
 protected:
     double _x{0.0}, _y{0.0};
     MccCoordPairKind _pairKind{MccCoordPairKind::COORDS_KIND_RADEC_ICRS};
     MccCelestialCoordEpoch _epoch{};  // J2000.0
 
+
     template <mcc_skypoint_c T>
     void fromOtherSkyPoint(T&& other)
     {

include/mcc/mcc_generic_mount.h

@@ -147,9 +147,9 @@ public:
                     std::tuple<MoveCntrCtorTs...> move_cntrl_ctor_ars,
                     std::tuple<LoggerCtorTs...> logger_ctor_args)
         : TELEMETRY_T(std::make_from_tuple<TELEMETRY_T>(std::move(telemetry_ctor_args))),
-          PZONE_CONT_T(std::make_from_tuple<PZONE_CONT_T>(pzone_cont_ctor_ars)),
-          MOVE_CNTRL_T(std::make_from_tuple<MOVE_CNTRL_T>(move_cntrl_ctor_ars)),
-          LOGGER_T(std::make_from_tuple<LOGGER_T>(logger_ctor_args))
+          PZONE_CONT_T(std::make_from_tuple<PZONE_CONT_T>(std::move(pzone_cont_ctor_ars))),
+          MOVE_CNTRL_T(std::make_from_tuple<MOVE_CNTRL_T>(std::move(move_cntrl_ctor_ars))),
+          LOGGER_T(std::make_from_tuple<LOGGER_T>(std::move(logger_ctor_args)))
     {
         // logDebug(std::format("Create MccGenericMount class instance (thread: {})", std::this_thread::get_id()));
         logDebug("Create MccGenericMount class instance (thread: {})", std::this_thread::get_id());
@@ -183,18 +183,21 @@ public:
 
     virtual ~MccGenericMount()
     {
+        // auto err = MOVE_CNTRL_T::stopMount();
+        // if (err) {
+        //     logError(formatError(err));
+        // }
+
+        // WARNING: it is assumed here that mount stopping is performed in a derived class or, it is more logicaly,
+        //          in MOVE_CNTRL_T-class
+
         // logDebug(std::format("Delete MccGenericMount class instance (thread: {})", std::this_thread::get_id()));
         logDebug("Delete MccGenericMount class instance (thread: {})", std::this_thread::get_id());
-
-        auto err = MOVE_CNTRL_T::stopMount();
-        if (err) {
-            logError(formatError(err));
-        }
     }
 
     error_t initMount()
     {
-        logInfo(std::format("Start MccGenericMount class initialization (thread: {}) ...", std::this_thread::get_id()));
+        logInfo("Start MccGenericMount class initialization (thread: {}) ...", std::this_thread::get_id());
 
         *_mountStatus = mount_status_t::MOUNT_STATUS_IDLE;
 

include/mcc/mcc_generic_movecontrols.h

@@ -0,0 +1,715 @@
+#pragma once
+
+/****************************************************************************************
+ *                                                                                      *
+ *                          MOUNT CONTROL COMPONENTS LIBRARY                            *
+ *                                                                                      *
+ *                                                                                      *
+ *                 GENERIC IMPLEMENTATION OF MOUNT MOVEMENT CONTROLS                    *
+ *                                                                                      *
+ ****************************************************************************************/
+
+
+
+#include <atomic>
+#include <fstream>
+#include <future>
+#include <thread>
+#include <type_traits>
+
+
+#include <print>
+
+#include "mcc_coordinate.h"
+#include "mcc_error.h"
+
+namespace mcc::impl
+{
+
+// mount movement-related generic errors
+enum class MccGenericMovementControlsErrorCode : int {
+    ERROR_OK,
+    ERROR_IN_PZONE,
+    ERROR_NEAR_PZONE,
+    ERROR_SLEW_TIMEOUT,
+    ERROR_STOP_TIMEOUT,
+    ERROR_HARDWARE,
+    ERROR_TELEMETRY_TIMEOUT
+};
+}  // namespace mcc::impl
+
+namespace std
+{
+template <>
+class is_error_code_enum<mcc::impl::MccGenericMovementControlsErrorCode> : public true_type
+{
+};
+}  // namespace std
+
+namespace mcc::impl
+{
+
+// error category
+struct MccGenericMovementControlsErrorCategory : std::error_category {
+    const char* name() const noexcept
+    {
+        return "MCC-GENERIC-MOVECONTRL";
+    }
+
+    std::string message(int ec) const
+    {
+        MccGenericMovementControlsErrorCode err = static_cast<MccGenericMovementControlsErrorCode>(ec);
+
+        switch (err) {
+            case MccGenericMovementControlsErrorCode::ERROR_OK:
+                return "OK";
+            case MccGenericMovementControlsErrorCode::ERROR_IN_PZONE:
+                return "target is in zone";
+            case MccGenericMovementControlsErrorCode::ERROR_NEAR_PZONE:
+                return "mount is near zone";
+            case MccGenericMovementControlsErrorCode::ERROR_SLEW_TIMEOUT:
+                return "a timeout occured while slewing";
+            case MccGenericMovementControlsErrorCode::ERROR_STOP_TIMEOUT:
+                return "a timeout occured while mount stopping";
+            case MccGenericMovementControlsErrorCode::ERROR_HARDWARE:
+                return "a hardware error occured";
+            case MccGenericMovementControlsErrorCode::ERROR_TELEMETRY_TIMEOUT:
+                return "telemetry data timeout";
+            default:
+                return "unknown";
+        }
+    }
+
+    static const MccGenericMovementControlsErrorCategory& get()
+    {
+        static const MccGenericMovementControlsErrorCategory constInst;
+        return constInst;
+    }
+};
+
+
+inline std::error_code make_error_code(MccGenericMovementControlsErrorCode ec)
+{
+    return std::error_code(static_cast<int>(ec), MccGenericMovementControlsErrorCategory::get());
+}
+
+
+struct MccGenericMovementControlsParams {
+    // timeout to telemetry updating
+    std::chrono::milliseconds telemetryTimeout{3000};
+
+    // braking acceleration after execution of mount stopping command (in rads/s^2)
+    // it must be given as non-negative value!!!
+    double brakingAccelX{0.0};
+    double brakingAccelY{0.0};
+
+    // *******    slewing mode    *******
+
+    // coordinates difference to stop slewing (in radians)
+    double slewToleranceRadius{5.0_arcsecs};
+
+    // slewing trajectory file. if empty - just skip saving
+    std::string slewingPathFilename{};
+
+
+    // *******    tracking mode    *******
+
+    // maximal target-to-mount difference for tracking process (in arcsecs)
+    // it it is greater then the current mount coordinates are used as target one
+    double trackingMaxCoordDiff{20.0};
+
+    // tracking trajectory file. if empty - just skip saving
+    std::string trackingPathFilename{};
+};
+
+
+/*                  UTILITY CLASS TO HOLD AND SAVE MOUNT MOVING TRAJECTORY                  */
+
+struct MccMovementPathFile {
+    static constexpr std::string_view commentSeq{"# "};
+    static constexpr std::string_view lineDelim{"\n"};
+
+    void setCommentSeq(traits::mcc_input_char_range auto const& s)
+    {
+        _commentSeq.clear();
+        std::ranges::copy(s, std::back_inserter(_commentSeq));
+    }
+
+    void setLineDelim(traits::mcc_input_char_range auto const& s)
+    {
+        _lineDelim.clear();
+        std::ranges::copy(s, std::back_inserter(_lineDelim));
+    }
+
+    // add comment string/strings
+    template <traits::mcc_input_char_range RT, traits::mcc_input_char_range... RTs>
+    void addComment(RT const& r, RTs const&... rs)
+    {
+        std::ranges::copy(_commentSeq, std::back_inserter(_buffer));
+        if constexpr (std::is_pointer_v<std::decay_t<RT>>) {  // const char*, char*, char[]
+            std::ranges::copy(std::string_view{r}, std::back_inserter(_buffer));
+        } else {
+            std::ranges::copy(r, std::back_inserter(_buffer));
+        }
+        std::ranges::copy(lineDelim, std::back_inserter(_buffer));
+
+        if constexpr (sizeof...(RTs)) {
+            addComment(rs...);
+        }
+    }
+
+    // comment corresponded to addToPath(mcc_telemetry_data_c auto const& tdata)
+    void addDefaultComment()
+    {
+        addComment("Format (time is in milliseconds, coordinates are in degrees, speeds are in degrees/s):");
+        addComment(
+            "    <UNIXTIME> <mount X> <mount Y> <target X> <target Y> <dX_{mount-target}> "
+            "<dY_{mount-target}> <mount-to-target-distance> <mount X-speed> <mount Y-speed> <moving state>");
+    }
+
+    // general purpose method
+    // template <std::formattable<char>... ArgTs>
+    // void addToPath(std::format_string<ArgTs...> fmt, ArgTs&&... args)
+    // {
+    //     std::format_to(std::back_inserter(_buffer), fmt, std::forward<ArgTs>(args)...);
+    //     std::ranges::copy(lineDelim, std::back_inserter(_buffer));
+    // }
+
+    // general purpose method
+    template <std::formattable<char>... ArgTs>
+    void addToPath(std::string_view fmt, ArgTs&&... args)
+    {
+        std::vformat_to(std::back_inserter(_buffer), fmt, std::make_format_args(args...));
+        std::ranges::copy(lineDelim, std::back_inserter(_buffer));
+    }
+
+    // default-implemented method
+    void addToPath(mcc_telemetry_data_c auto const& tdata)
+    {
+        // UNIX-time millisecs, mount X, mount Y, target X, target Y, dX(mount-target), dY(mount-target), dist, speedX,
+        // speedY, state
+        auto dist = tdata.mountPos.distance(tdata.targetPos);
+
+        using d_t = std::chrono::milliseconds;
+
+        auto tp = std::chrono::duration_cast<d_t>(tdata.mountPos.epoch().UTC().time_since_epoch());
+
+        const std::string_view d_fmt = "{:14.8f}";
+        const auto v = std::views::repeat(d_fmt, 9) | std::views::join_with(' ');
+        std::string fmt = "{} " + std::string(v.begin(), v.end()) + " {}";
+
+
+        int state = (int)tdata.hwState.movementState;
+        auto tp_val = tp.count();
+
+        double mnt_x = MccAngle(tdata.mountPos.co_lon()).degrees(), mnt_y = MccAngle(tdata.mountPos.co_lon()).degrees(),
+               tag_x = dist.x2.degrees(), tag_y = dist.y2.degrees(), dx = dist.dx.degrees(), dy = dist.dy.degrees(),
+               dd = dist.dist.degrees();
+
+        addToPath(std::string_view(fmt.begin(), fmt.end()), tp_val, mnt_x, mnt_y, tag_x, tag_y, dx, dy, dd,
+                  tdata.hwState.speedXY.x().degrees(), tdata.hwState.speedXY.y().degrees(), state);
+    }
+
+    void clearPath()
+    {
+        _buffer.clear();
+    }
+
+    bool saveToFile(std::string const& filename, std::ios_base::openmode mode = std::ios::out | std::ios::trunc)
+    {
+        if (filename.empty()) {
+            return true;
+        }
+
+        std::ofstream fst(filename, mode);
+
+        if (fst.is_open()) {
+            fst << _buffer;
+
+            return true;
+        } else {
+            return false;
+        }
+    }
+
+protected:
+    std::string _buffer{};
+
+    std::string _commentSeq{commentSeq};
+    std::string _lineDelim{lineDelim};
+};
+
+
+enum class MccGenericMovementControlsPolicy : int { POLICY_ASYNC, POLICY_BLOCKING };
+
+template <std::movable PARAMS_T,
+          MccGenericMovementControlsPolicy EXEC_POLICY = MccGenericMovementControlsPolicy::POLICY_ASYNC>
+class MccGenericMovementControls
+{
+public:
+    static constexpr MccGenericMovementControlsPolicy executePolicy = EXEC_POLICY;
+
+    static constexpr std::chrono::seconds defaultWaitTimeout{3};
+
+    typedef MccError error_t;
+
+    typedef PARAMS_T movement_params_t;
+
+    template <std::invocable<bool> SLEW_FUNC_T, std::invocable<> TRACK_FUNC_T, std::invocable<> STOP_FUNC_T>
+    MccGenericMovementControls(SLEW_FUNC_T&& slew_func, TRACK_FUNC_T&& track_func, STOP_FUNC_T&& stop_func)
+        : _slewFunc(std::forward<SLEW_FUNC_T>(slew_func)),
+          _trackFunc(std::forward<TRACK_FUNC_T>(track_func)),
+          _stopFunc(std::forward<STOP_FUNC_T>(stop_func))
+    {
+        if constexpr (executePolicy == MccGenericMovementControlsPolicy::POLICY_ASYNC) {
+            // *_stopMovementRequest = false;
+
+            // *_fsmState = STATE_IDLE;
+
+            // // start thread of movements
+            // _fstFuture = std::async(
+            //     [this](std::stop_token stoken) {
+            //         auto do_state = _fsmState->load();
+
+            //         while (!stoken.stop_requested()) {
+            //             std::println("\n{:*^80}\n", "  WAIT LOCK  ");
+            //             // wait here ...
+            //             _wakeupRequest->wait(false, std::memory_order_relaxed);
+
+            //             _wakeupRequest->clear();
+
+            //             std::println("\n{:*^80}\n", "  UNLOCKED  ");
+
+            //             // if (stoken.stop_requested()) {
+            //             //     break;
+            //             // }
+
+            //             do_state = _fsmState->load();
+
+            //             if (do_state & STATE_STOP) {
+            //                 // if (_fsmState->load() & STATE_STOP) {
+            //                 *_stopMovementRequest = true;
+            //                 _stopFunc();
+            //             }
+
+            //             if (do_state & STATE_SLEW) {
+            //                 // if (_fsmState->load() & STATE_SLEW) {
+            //                 *_stopMovementRequest = false;
+            //                 _slewFunc(_slewAndStop->load());
+            //             } else if (do_state & STATE_TRACK) {
+            //                 // } else if (_fsmState->load() & STATE_TRACK) {
+            //                 *_stopMovementRequest = false;
+            //                 _trackFunc();
+            //             }
+            //         }
+            //     },
+            //     _fstStopSource.get_token());
+
+            startAsyncMovementCycle();
+        }
+    }
+
+    MccGenericMovementControls(const MccGenericMovementControls&) = delete;
+    MccGenericMovementControls(MccGenericMovementControls&& other) = default;
+
+    MccGenericMovementControls& operator=(const MccGenericMovementControls&) = delete;
+    MccGenericMovementControls& operator=(MccGenericMovementControls&&) = default;
+
+    virtual ~MccGenericMovementControls()
+    {
+        if constexpr (executePolicy == MccGenericMovementControlsPolicy::POLICY_ASYNC) {
+            // *_fsmState = STATE_IDLE;
+            // _fstStopSource.request_stop();
+
+            // _wakeupRequest->test_and_set();
+            // _wakeupRequest->notify_one();
+
+            // std::this_thread::sleep_for(std::chrono::milliseconds(200));
+
+            // if (_fstFuture.valid()) {
+            //     auto status = _fstFuture.wait_for(_waitTimeout->load());
+            // }
+
+            stopAsyncMovementCycle();
+        }
+    }
+
+
+    template <typename SelfT>
+    void stopAsyncMovementCycle(this SelfT&& self)
+        requires(EXEC_POLICY == MccGenericMovementControlsPolicy::POLICY_ASYNC)
+    {
+        *self._fsmState = STATE_IDLE;
+        self._fstStopSource.request_stop();
+
+        self._wakeupRequest->test_and_set();
+        self._wakeupRequest->notify_all();
+
+        std::this_thread::sleep_for(std::chrono::milliseconds(200));
+
+        if (self._fstFuture.valid()) {
+            auto status = self._fstFuture.wait_for(self._waitTimeout->load());
+        }
+    }
+
+    template <typename SelfT>
+    void startAsyncMovementCycle(this SelfT& self)
+        requires(EXEC_POLICY == MccGenericMovementControlsPolicy::POLICY_ASYNC)
+    {
+        using self_t = std::decay_t<SelfT>;
+
+        if (self._fstFuture.valid()) {
+            self.stopAsyncMovementCycle();
+        }
+
+        *self._stopMovementRequest = false;
+
+        self._fstStopSource = std::stop_source{};
+
+        *self._fsmState = STATE_IDLE;
+
+        // start thread of movements
+        self._fstFuture = std::async(
+            std::launch::async,
+            [&self](std::stop_token stoken) mutable {
+                auto do_state = self._fsmState->load();
+
+                std::string log_str;
+                while (!stoken.stop_requested()) {
+                    log_str = std::format("\n{:*^80}\n", "  WAIT LOCK  ");
+                    if constexpr (mcc_generic_mount_c<self_t>) {
+                        self.logTrace(log_str);
+                    } else {
+                        std::println("{}", log_str);
+                    }
+
+                    // wait here ...
+                    self._wakeupRequest->wait(false, std::memory_order_relaxed);
+
+                    self._wakeupRequest->clear();
+
+                    log_str = std::format("\n{:*^80}\n", "  UNLOCKED  ");
+                    if constexpr (mcc_generic_mount_c<self_t>) {
+                        self.logTrace(log_str);
+                    } else {
+                        std::println("{}", log_str);
+                    }
+
+                    // if (stoken.stop_requested()) {
+                    //     break;
+                    // }
+
+                    do_state = self._fsmState->load();
+
+                    if (do_state & STATE_STOP) {
+                        // if (_fsmState->load() & STATE_STOP) {
+                        *self._stopMovementRequest = true;
+                        self._stopFunc();
+                    }
+
+                    if (do_state & STATE_SLEW) {
+                        // if (_fsmState->load() & STATE_SLEW) {
+                        *self._stopMovementRequest = false;
+                        self._slewFunc(self._slewAndStop->load());
+                    } else if (do_state & STATE_TRACK) {
+                        // } else if (_fsmState->load() & STATE_TRACK) {
+                        *self._stopMovementRequest = false;
+                        self._trackFunc();
+                    }
+                }
+            },
+            self._fstStopSource.get_token());
+    }
+
+    error_t slewToTarget(bool slew_and_stop)
+    {
+        // *_stopMovementRequest = false;
+
+        if constexpr (executePolicy == MccGenericMovementControlsPolicy::POLICY_ASYNC) {
+            auto prev_state = _fsmState->load();
+
+            *_fsmState = STATE_SLEW;
+
+            if (prev_state & STATE_SLEW || prev_state & STATE_TRACK) {
+                *_fsmState |= STATE_STOP;
+
+                *_stopMovementRequest = true;  // to exit from slewing or tracking
+            }
+
+            *_slewAndStop = slew_and_stop;
+
+            _wakeupRequest->test_and_set();
+            _wakeupRequest->notify_one();
+
+            return MccGenericMovementControlsErrorCode::ERROR_OK;
+        } else if constexpr (executePolicy == MccGenericMovementControlsPolicy::POLICY_BLOCKING) {
+            return _slewFunc(slew_and_stop);
+        } else {
+            static_assert(false, "UNKNOWN EXECUTION POLICY!");
+        }
+    }
+
+    error_t trackTarget()
+    {
+        // *_stopMovementRequest = false;
+
+        if constexpr (executePolicy == MccGenericMovementControlsPolicy::POLICY_ASYNC) {
+            auto prev_state = _fsmState->load();
+
+            if (!(prev_state & STATE_TRACK)) {
+                *_fsmState = STATE_TRACK;
+
+                if (prev_state & STATE_SLEW) {
+                    *_fsmState |= STATE_STOP;
+
+                    *_stopMovementRequest = true;  // to exit from slewing
+                }
+
+                _wakeupRequest->test_and_set();
+                _wakeupRequest->notify_one();
+            }  // already tracking, just ignore
+
+            return MccGenericMovementControlsErrorCode::ERROR_OK;
+        } else if constexpr (executePolicy == MccGenericMovementControlsPolicy::POLICY_BLOCKING) {
+            return _trackFunc();
+        } else {
+            static_assert(false, "UNKNOWN EXECUTION POLICY!");
+        }
+    }
+
+    error_t stopMount()
+    {
+        // *_stopMovementRequest = true;
+
+        if constexpr (executePolicy == MccGenericMovementControlsPolicy::POLICY_ASYNC) {
+            *_fsmState = STATE_STOP;
+
+            *_stopMovementRequest = true;  // to exit from slewing or tracking
+
+            _wakeupRequest->test_and_set();
+            _wakeupRequest->notify_one();
+
+            return MccGenericMovementControlsErrorCode::ERROR_OK;
+        } else if constexpr (executePolicy == MccGenericMovementControlsPolicy::POLICY_BLOCKING) {
+            return _stopFunc();
+        } else {
+            static_assert(false, "UNKNOWN EXECUTION POLICY!");
+        }
+    }
+
+    error_t setMovementParams(movement_params_t const& pars)
+    {
+        std::lock_guard lock{*_currentMovementParamsMutex};
+
+        _currentMovementParams = pars;
+
+        return MccGenericMovementControlsErrorCode::ERROR_OK;
+    }
+
+    movement_params_t getMovementParams() const
+    {
+        std::lock_guard lock{*_currentMovementParamsMutex};
+
+        return _currentMovementParams;
+    }
+
+protected:
+    std::unique_ptr<std::mutex> _currentMovementParamsMutex{new std::mutex{}};
+    PARAMS_T _currentMovementParams{};
+
+    std::unique_ptr<std::atomic_bool> _stopMovementRequest{new std::atomic_bool{false}};
+
+    std::function<error_t(bool)> _slewFunc{};
+    std::function<error_t()> _trackFunc{};
+    std::function<error_t()> _stopFunc{};
+
+    std::unique_ptr<std::atomic<std::chrono::nanoseconds>> _waitTimeout{
+        new std::atomic<std::chrono::nanoseconds>{defaultWaitTimeout}};
+
+
+
+    std::conditional_t<executePolicy == MccGenericMovementControlsPolicy::POLICY_ASYNC,
+                       std::future<void>,
+                       std::nullptr_t>
+        _fstFuture{};
+
+    std::stop_source _fstStopSource{};
+
+    std::unique_ptr<std::atomic_flag> _wakeupRequest{new std::atomic_flag};
+    std::unique_ptr<std::atomic_bool> _slewAndStop{new std::atomic_bool{false}};
+
+    enum { STATE_IDLE = 0x00, STATE_SLEW = 0x01, STATE_TRACK = 0x02, STATE_STOP = 0x04 };
+
+    std::unique_ptr<std::atomic_int> _fsmState{new std::atomic_int};
+
+
+    // template <typename SelfT>
+    // void mainCycle(this SelfT&& self, std::stop_token stoken)
+    // {
+    //     using self_t = std::decay_t<SelfT>;
+
+    //     auto do_state = _fsmState->load();
+
+    //     std::string log_str;
+
+    //     while (!stoken.stop_requested()) {
+    //         log_str = std::format("\n{:*^80}\n", "  WAIT LOCK  ");
+    //         if constexpr (mcc_generic_mount_c<self_t>) {
+    //             self.logTrace(log_str);
+    //         } else {
+    //             std::println(log_str);
+    //         }
+
+    //         // wait here ...
+    //         _wakeupRequest->wait(false, std::memory_order_relaxed);
+
+    //         _wakeupRequest->clear();
+
+    //         log_str = std::format("\n{:*^80}\n", "  UNLOCKED  ");
+    //         if constexpr (mcc_generic_mount_c<self_t>) {
+    //             self.logTrace(log_str);
+    //         } else {
+    //             std::println(log_str);
+    //         }
+
+    //         // if (stoken.stop_requested()) {
+    //         //     break;
+    //         // }
+
+    //         do_state = _fsmState->load();
+
+    //         if (do_state & STATE_STOP) {
+    //             // if (_fsmState->load() & STATE_STOP) {
+    //             *_stopMovementRequest = true;
+    //             _stopFunc();
+    //         }
+
+    //         if (do_state & STATE_SLEW) {
+    //             // if (_fsmState->load() & STATE_SLEW) {
+    //             *_stopMovementRequest = false;
+    //             _slewFunc(_slewAndStop->load());
+    //         } else if (do_state & STATE_TRACK) {
+    //             // } else if (_fsmState->load() & STATE_TRACK) {
+    //             *_stopMovementRequest = false;
+    //             _trackFunc();
+    //         }
+    //     }
+    // }
+
+    // utility methods
+
+    // the method calculates the change in coordinates of a point over a given time given the current speed and braking
+    // acceleration. a position after given 'time' interval is returned
+    auto coordsAfterTime(mcc_coord_pair_c auto const& cp,
+                         mcc_coord_pair_c auto const& speedXY,          // in radians per seconds
+                         mcc_coord_pair_c auto const& braking_accelXY,  // in radians per seconds in square
+                         traits::mcc_time_duration_c auto const& time,
+                         mcc_coord_pair_c auto* dxy = nullptr)
+    {
+        // time to stop mount with given current speed and constant braking acceleration
+        double tx_stop = std::abs(speedXY.x()) / braking_accelXY.x();
+        double ty_stop = std::abs(speedXY.y()) / braking_accelXY.y();
+
+        using secs_t = std::chrono::duration<double>;  // seconds as double
+
+        double tx = std::chrono::duration_cast<secs_t>(time).count();
+        double ty = std::chrono::duration_cast<secs_t>(time).count();
+
+        if (std::isfinite(tx_stop) && (tx > tx_stop)) {
+            tx = tx_stop;
+        }
+        if (std::isfinite(ty_stop) && (ty > ty_stop)) {
+            ty = ty_stop;
+        }
+
+        // the distance:
+        // here, one must take into account the sign of the speed!!!
+        double dx = speedXY.x() * tx - std::copysign(braking_accelXY.x(), speedXY.x()) * tx * tx / 2.0;
+        double dy = speedXY.y() * ty - std::copysign(braking_accelXY.y(), speedXY.y()) * ty * ty / 2.0;
+
+        std::remove_cvref_t<decltype(cp)> cp_res{};
+        cp_res.setEpoch(cp.epoch() + time);
+
+        cp_res.setX((double)cp.x() + dx);
+        cp_res.setY((double)cp.y() + dy);
+
+        if (dxy) {
+            dxy->setX(dx);
+            dxy->setY(dy);
+            dxy->setEpoch(cp_res.epoch());
+        }
+
+        return cp_res;
+    }
+
+    auto coordsAfterTime(mcc_skypoint_c auto const& sp,
+                         mcc_coord_pair_c auto const& speedXY,          // in radians per seconds
+                         mcc_coord_pair_c auto const& braking_accelXY,  // in radians per seconds in square
+                         traits::mcc_time_duration_c auto const& time,
+                         mcc_coord_pair_c auto* dxy = nullptr)
+    {
+        auto run_func = [&, this](auto& cp) {
+            sp.toAtSameEpoch(cp);
+
+            auto new_cp = coordsAfterTime(cp, speedXY, braking_accelXY, time, dxy);
+
+            std::remove_cvref_t<decltype(sp)> sp_res{};
+
+            sp_res.from(cp);
+
+            return sp_res;
+        };
+
+        switch (sp.pairKind()) {
+            case MccCoordPairKind::COORDS_KIND_RADEC_ICRS: {
+                MccSkyRADEC_ICRS rd;
+                return run_func(rd);
+            };
+            case MccCoordPairKind::COORDS_KIND_RADEC_OBS: {
+                MccSkyRADEC_OBS rd;
+                return run_func(rd);
+            };
+            case MccCoordPairKind::COORDS_KIND_RADEC_APP: {
+                MccSkyRADEC_APP rd;
+                return run_func(rd);
+            };
+            case MccCoordPairKind::COORDS_KIND_HADEC_OBS: {
+                MccSkyHADEC_OBS hd;
+                return run_func(hd);
+            };
+            case MccCoordPairKind::COORDS_KIND_HADEC_APP: {
+                MccSkyHADEC_APP hd;
+                return run_func(hd);
+            };
+            case MccCoordPairKind::COORDS_KIND_AZZD: {
+                MccSkyAZZD azzd;
+                return run_func(azzd);
+            };
+            case MccCoordPairKind::COORDS_KIND_AZALT: {
+                MccSkyAZALT azalt;
+                return run_func(azalt);
+            };
+            case MccCoordPairKind::COORDS_KIND_GENERIC:
+            case MccCoordPairKind::COORDS_KIND_XY: {
+                MccGenXY xy;
+                return run_func(xy);
+            };
+            default:
+                return sp;
+        }
+    }
+};
+
+
+template <std::default_initializable PARAMS_T>
+using MccGenericBlockingMovementControls =
+    MccGenericMovementControls<PARAMS_T, MccGenericMovementControlsPolicy::POLICY_BLOCKING>;
+
+template <std::default_initializable PARAMS_T>
+using MccGenericAsyncMovementControls =
+    MccGenericMovementControls<PARAMS_T, MccGenericMovementControlsPolicy::POLICY_ASYNC>;
+
+
+static_assert(mcc_movement_controls_c<MccGenericAsyncMovementControls<MccGenericMovementControlsParams>>, "!!!");
+
+}  // namespace mcc::impl

include/mcc/mcc_netserver.h

@@ -967,9 +967,6 @@ public:
 
             return output_msg.template byteRepr<typename base_t::handle_message_func_result_t>();
         };
-
-        // special functor (used in the destructor)
-        _stopMountFunc = [mount_ptr]() { mount_ptr->stopMount(); };
     }
 
     virtual ~MccGenericMountNetworkServer()
@@ -977,8 +974,6 @@ public:
         std::stringstream st;
         st << std::this_thread::get_id();
 
-        _stopMountFunc();
-
         logInfo(std::format("Delete MccGenericMountNetworkServer class instance (thread ID = {})", st.str()));
     }
 
@@ -987,8 +982,6 @@ protected:
     MccSerializedAngleFormatPrec _coordPrec{2, 1, 7};
 
 
-    std::function<void()> _stopMountFunc{};
-
     template <typename RESULT_MSG_T, typename INPUT_MSG_T, mcc_generic_mount_c MountT>
     RESULT_MSG_T handleMessage(const INPUT_MSG_T& input_msg, MountT* mount_ptr)
         requires(

include/mcc/mcc_pcm.h

@@ -2,9 +2,15 @@
 #pragma once
 
 
-/*                          MOUNT CONTROL COMPONENTS LIBRARY                            */
+/****************************************************************************************
+ *                                                                                      *
+ *                          MOUNT CONTROL COMPONENTS LIBRARY                            *
+ *                                                                                      *
+ *                                                                                      *
+ *             A REFERENCE "POINTING-CORRECTION-MODEL" CLASS IMPLEMENTATION             *
+ *                                                                                      *
+ ****************************************************************************************/
 
-/*                          A REFERENCE "POINTING-CORRECTION-MODEL" CLASS IMPLEMENTATION                            */
 
 
 #include <mutex>
@@ -26,7 +32,8 @@ enum class MccDefaultPCMErrorCode : int {
     ERROR_INVALID_INPUTS_BISPLEV,
 #endif
     ERROR_EXCEED_MAX_ITERS,
-    ERROR_NULLPTR
+    ERROR_NULLPTR,
+    ERROR_JACINV
 };
 
 /*  error category definition  */
@@ -57,6 +64,8 @@ struct MccDefaultPCMCategory : public std::error_category {
                 return "exceed maximum of iterations number";
             case MccDefaultPCMErrorCode::ERROR_NULLPTR:
                 return "nullptr input argument";
+            case MccDefaultPCMErrorCode::ERROR_JACINV:
+                return "Jacobian is near singular";
             default:
                 return "UNKNOWN";
         }
@@ -157,6 +166,9 @@ public:
 
         std::vector<coeff_t> coeffsX{};
         std::vector<coeff_t> coeffsY{};
+
+        std::vector<coeff_t> inverseCoeffsX{};
+        std::vector<coeff_t> inverseCoeffsY{};
     };
 #endif
 
@@ -316,14 +328,65 @@ public:
 
         std::lock_guard lock(*_pcmDataMutex);
 
-        ret = _compResult(x, y, res, true);
-        if (!ret) {
+        if (_pcmData.type != MccDefaultPCMType::PCM_TYPE_BSPLINE) {  // compute using Newton-Raphson correction
+            struct {
+                double pcmX, pcmY;
+            } dfx, dfy, df;  // partial derivatives
+
+            ret = _computeFuncDeriv(x, y, res, true, &dfx, &dfy);
+            if (!ret) {
+                return ret;
+            }
+
+
+            dfx.pcmX += 1.0;  // a
+            dfx.pcmY += 1.0;  // b
+            dfy.pcmX += 1.0;  // c
+            dfy.pcmY += 1.0;  // d
+
+            // Jacobian determinant
+            auto detJ = dfx.pcmX * dfy.pcmY - dfx.pcmY * dfy.pcmX;
+            if (utils::isEqual(detJ, 0.0)) {
+                return MccDefaultPCMErrorCode::ERROR_JACINV;
+            }
+
+            //        | a  b |
+            // if A = | c  d |, then
+            //
+            //     -1            | d  -b |
+            //    A   = 1/detA * | -c  a |
+            //
+
+            df.pcmX = dfy.pcmY * res->pcmX - dfx.pcmY * res->pcmY;
+            df.pcmY = -dfy.pcmX * res->pcmX + dfx.pcmX * res->pcmY;
+
+            res->pcmX -= df.pcmX;
+            res->pcmY -= df.pcmY;
+
             if constexpr (mcc_coord_pair_c<T>) {
-                *hw_pt =
-                    MccCoordPair<typename T::x_t, typename T::y_t>{x + res->pcmX, y + res->pcmY, obs_skycoord.epoch()};
+                *hw_pt = MccCoordPair<typename T::x_t, typename T::y_t>{res->pcmX, res->pcmY, obs_skycoord.epoch()};
+            }
+        } else {  // for B-splines the result is computed directly from inverse B-spline coefficients
+            ret = _computeFuncDeriv(x, y, res);
+
+            if (!ret) {
+                if constexpr (mcc_coord_pair_c<T>) {
+                    *hw_pt = MccCoordPair<typename T::x_t, typename T::y_t>{x + res->pcmX, y + res->pcmY,
+                                                                            obs_skycoord.epoch()};
+                }
             }
         }
 
+
+        // ret = _compResult(x, y, res, true);
+        // if (!ret) {
+        //     if constexpr (mcc_coord_pair_c<T>) {
+        //         *hw_pt =
+        //             MccCoordPair<typename T::x_t, typename T::y_t>{x + res->pcmX, y + res->pcmY,
+        //             obs_skycoord.epoch()};
+        //     }
+        // }
+
         return ret;
     }
 
@@ -384,6 +447,204 @@ private:
 
     std::unique_ptr<std::mutex> _pcmDataMutex{new std::mutex};
 
+
+    // compute PCM function and its partial derivatives if asked
+    template <typename DXT = std::nullptr_t, typename DYT = std::nullptr_t>
+    error_t _computeFuncDeriv(double x,
+                              double y,
+                              mcc_pcm_result_c auto* res,
+                              bool inverse = false,
+                              DXT derivX = nullptr,
+                              DYT derivY = nullptr)
+        requires((std::is_null_pointer_v<DXT> || mcc_pcm_result_c<std::remove_pointer_t<DXT>>) &&
+                 (std::is_null_pointer_v<DYT> || mcc_pcm_result_c<std::remove_pointer_t<DYT>>))
+    {
+        if constexpr (std::is_null_pointer_v<DXT> || std::is_null_pointer_v<DYT>) {
+            if (_pcmData.type != MccDefaultPCMType::PCM_TYPE_BSPLINE && inverse) {
+                return MccDefaultPCMErrorCode::ERROR_NULLPTR;
+            }
+        }
+
+        pcm_geom_coeffs_t* geom_coeffs = &_pcmData.geomCoefficients;
+
+#ifdef USE_BSPLINE_PCM
+        pcm_bspline_t* bspline = &_pcmData.bspline;
+        // pcm_bspline_t* inv_bspline = &_pcmData.inverseBspline;
+#endif
+
+#ifdef USE_BSPLINE_PCM
+        if (_pcmData.type == MccDefaultPCMType::PCM_TYPE_GEOMETRY ||
+            _pcmData.type == MccDefaultPCMType::PCM_TYPE_GEOMETRY_BSPLINE) {
+#endif
+            const auto tanY = std::tan(y);
+            const auto sinX = std::sin(x);
+            const auto cosX = std::cos(x);
+            const auto cosY = std::cos(y);
+            const auto sinY = std::sin(y);
+
+
+            if (utils::isEqual(cosY, 0.0)) {
+                res->pcmX = _pcmData.geomCoefficients.zeroPointX;
+                if constexpr (!std::is_null_pointer_v<DXT>) {
+                    derivX->pcmX = 0.0;  // dpcmX/dX
+                    derivX->pcmY = 0.0;  // dpcmX/dY
+                }
+            } else {
+                res->pcmX = geom_coeffs->zeroPointX + geom_coeffs->collimationErr / cosY +
+                            geom_coeffs->nonperpendErr * tanY - geom_coeffs->misalignErr1 * cosX * tanY +
+                            geom_coeffs->misalignErr2 * sinX * tanY + geom_coeffs->tubeFlexure * _cosPhi * sinX / cosY -
+                            geom_coeffs->DECaxisFlexure * (_cosPhi * cosX + _sinPhi * tanY);
+
+                if constexpr (!std::is_null_pointer_v<DXT>) {
+                    auto cos2Y = cosY * cosY;
+
+                    derivX->pcmX = (geom_coeffs->misalignErr1 * sinX + geom_coeffs->misalignErr2 * cosX) * tanY +
+                                   geom_coeffs->tubeFlexure * _cosPhi * cosX / cosY +
+                                   geom_coeffs->DECaxisFlexure * _cosPhi * sinX;  // dpcmX/dX
+
+                    derivX->pcmY =
+                        (geom_coeffs->collimationErr * sinY + geom_coeffs->nonperpendErr -
+                         geom_coeffs->misalignErr1 * cosX + geom_coeffs->misalignErr2 * sinX +
+                         geom_coeffs->tubeFlexure * _cosPhi * sinX * sinY - geom_coeffs->DECaxisFlexure * _sinPhi) /
+                        cos2Y;  // dpcmX/dY
+                }
+            }
+
+            res->pcmY = geom_coeffs->zeroPointY + geom_coeffs->misalignErr1 * sinX + geom_coeffs->misalignErr2 * cosX +
+                        geom_coeffs->tubeFlexure * (_cosPhi * cosX * sinY - _sinPhi * cosY);
+
+            if constexpr (!std::is_null_pointer_v<DYT>) {
+                derivY->pcmX = geom_coeffs->misalignErr1 * cosX - geom_coeffs->misalignErr2 * sinX -
+                               geom_coeffs->tubeFlexure * _cosPhi * sinX * sinY;  // dpcmY/dX
+
+                derivY->pcmY = geom_coeffs->tubeFlexure * (_cosPhi * cosX * cosY + _sinPhi * sinY);  // dpcmY/dY
+            }
+
+            if constexpr (pcmMountType == MccMountType::FORK_TYPE) {
+                if (!utils::isEqual(cosX, 0.0)) {
+                    res->pcmY += geom_coeffs->forkFlexure / cosX;
+
+                    if constexpr (!std::is_null_pointer_v<DYT>) {
+                        derivY->pcmY += geom_coeffs->forkFlexure * sinX / cosX / cosY;  // dpcmY/dY
+                    }
+                }
+            }
+#ifdef USE_BSPLINE_PCM
+        }
+#endif
+
+
+#ifdef USE_BSPLINE_PCM
+        if (_pcmData.type == MccDefaultPCMType::PCM_TYPE_BSPLINE ||
+            (_pcmData.type == MccDefaultPCMType::PCM_TYPE_GEOMETRY_BSPLINE && !inverse)) {
+            double spl_valX, spl_valY;
+
+            int ret = bsplines::fitpack_eval_spl2d(bspline->knotsX, bspline->knotsY, bspline->coeffsX, x, y, spl_valX,
+                                                   bspline->bsplDegreeX, bspline->bsplDegreeY);
+
+            if (ret) {
+                res->pcmX = std::numeric_limits<double>::quiet_NaN();
+                res->pcmY = std::numeric_limits<double>::quiet_NaN();
+                return MccDefaultPCMErrorCode::ERROR_INVALID_INPUTS_BISPLEV;
+            }
+
+
+            ret = bsplines::fitpack_eval_spl2d(bspline->knotsX, bspline->knotsY, bspline->coeffsY, x, y, spl_valY,
+                                               bspline->bsplDegreeX, bspline->bsplDegreeY);
+
+            if (ret) {
+                res->pcmX = std::numeric_limits<double>::quiet_NaN();
+                res->pcmY = std::numeric_limits<double>::quiet_NaN();
+                return MccDefaultPCMErrorCode::ERROR_INVALID_INPUTS_BISPLEV;
+            }
+
+
+            res->pcmX += spl_valX;
+            res->pcmY += spl_valY;
+        }
+
+        // compute partial derivatives of the bivariate B-spline
+        if (_pcmData.type == MccDefaultPCMType::PCM_TYPE_GEOMETRY_BSPLINE && inverse) {
+            double dspl_valX, dspl_valY;
+
+            int ret = 0;
+
+            if constexpr (!std::is_null_pointer_v<DXT>) {
+                ret = bsplines::fitpack_parder_spl2d(bspline->knotsX, bspline->knotsY, bspline->coeffsX, x, y,
+                                                     dspl_valX, 1, 0, bspline->bsplDegreeX, bspline->bsplDegreeY);
+
+                if (ret) {
+                    return MccDefaultPCMErrorCode::ERROR_INVALID_INPUTS_BISPLEV;
+                }
+
+                derivX->pcmX += dspl_valX;  // dpcmX/dX
+
+                ret = bsplines::fitpack_parder_spl2d(bspline->knotsX, bspline->knotsY, bspline->coeffsX, x, y,
+                                                     dspl_valX, 0, 1, bspline->bsplDegreeX, bspline->bsplDegreeY);
+
+                if (ret) {
+                    return MccDefaultPCMErrorCode::ERROR_INVALID_INPUTS_BISPLEV;
+                }
+
+                derivX->pcmY += dspl_valX;  // dpcmX/dY
+            }
+
+            if constexpr (!std::is_null_pointer_v<DYT>) {
+                ret = bsplines::fitpack_parder_spl2d(bspline->knotsX, bspline->knotsY, bspline->coeffsY, x, y,
+                                                     dspl_valY, 1, 0, bspline->bsplDegreeX, bspline->bsplDegreeY);
+
+                if (ret) {
+                    return MccDefaultPCMErrorCode::ERROR_INVALID_INPUTS_BISPLEV;
+                }
+
+                derivY->pcmX += dspl_valY;  // dpcmY/dX
+
+                ret = bsplines::fitpack_parder_spl2d(bspline->knotsX, bspline->knotsY, bspline->coeffsY, x, y,
+                                                     dspl_valY, 0, 1, bspline->bsplDegreeX, bspline->bsplDegreeY);
+
+                if (ret) {
+                    return MccDefaultPCMErrorCode::ERROR_INVALID_INPUTS_BISPLEV;
+                }
+
+                derivY->pcmY += dspl_valY;  // dpcmY/dY
+            }
+        }
+
+        // for inverse PCM the inverse spline coefficients are used (derivatives are not computed)!!!
+        if (_pcmData.type == MccDefaultPCMType::PCM_TYPE_BSPLINE && inverse) {
+            double spl_valX, spl_valY;
+
+            int ret = bsplines::fitpack_eval_spl2d(bspline->knotsX, bspline->knotsY, bspline->inverseCoeffsX, x, y,
+                                                   spl_valX, bspline->bsplDegreeX, bspline->bsplDegreeY);
+
+            if (ret) {
+                res->pcmX = std::numeric_limits<double>::quiet_NaN();
+                res->pcmY = std::numeric_limits<double>::quiet_NaN();
+
+                return MccDefaultPCMErrorCode::ERROR_INVALID_INPUTS_BISPLEV;
+            }
+
+
+            ret = bsplines::fitpack_eval_spl2d(bspline->knotsX, bspline->knotsY, bspline->inverseCoeffsY, x, y,
+                                               spl_valY, bspline->bsplDegreeX, bspline->bsplDegreeY);
+
+            if (ret) {
+                res->pcmX = std::numeric_limits<double>::quiet_NaN();
+                res->pcmY = std::numeric_limits<double>::quiet_NaN();
+
+                return MccDefaultPCMErrorCode::ERROR_INVALID_INPUTS_BISPLEV;
+            }
+
+
+            res->pcmX = spl_valX;
+            res->pcmY = spl_valY;
+        }
+#endif
+
+
+        return MccDefaultPCMErrorCode::ERROR_OK;
+    }
+
     error_t _compResult(double x, double y, mcc_pcm_result_c auto* res, bool inverse)
     {
         pcm_geom_coeffs_t* geom_coeffs;

include/mcc/mcc_pcm_construct.h

@@ -0,0 +1,373 @@
+#pragma once
+
+/****************************************************************************************
+ *                                                                                      *
+ *                          MOUNT CONTROL COMPONENTS LIBRARY                            *
+ *                                                                                      *
+ *                                                                                      *
+ *                        "POINTING-CORRECTION-MODEL" FITTER                            *
+ *                                                                                      *
+ ****************************************************************************************/
+
+#include <eigen3/Eigen/Dense>
+
+#include "mcc_concepts.h"
+#include "mcc_coordinate.h"
+#include "mcc_pcm.h"
+
+namespace mcc::impl
+{
+
+enum class MccDefaultPCMConstructorErrorCode : int {
+    ERROR_OK,
+    ERROR_NOT_ENOUGH_DATA,
+#ifdef USE_BSPLINE_PCM
+    ERROR_INVALID_KNOTS_NUMBER,
+    ERROR_BSPLINE_FIT
+#endif
+};
+
+// error category
+
+struct MccDefaultPCMConstructorErrorCategory : std::error_category {
+    MccDefaultPCMConstructorErrorCategory() = default;
+
+    const char* name() const noexcept
+    {
+        return "MCC-DEFAULT-PCM-CONSTRUCTOR-ERROR-CATEGORY";
+    }
+
+    std::string message(int ec) const
+    {
+        MccDefaultPCMConstructorErrorCode err = static_cast<MccDefaultPCMConstructorErrorCode>(ec);
+
+        switch (err) {
+            case MccDefaultPCMConstructorErrorCode::ERROR_OK:
+                return "OK";
+            case MccDefaultPCMConstructorErrorCode::ERROR_NOT_ENOUGH_DATA:
+                return "not enough data point";
+
+#ifdef USE_BSPLINE_PCM
+            case MccDefaultPCMConstructorErrorCode::ERROR_INVALID_KNOTS_NUMBER:
+                return "invalid number of B-spline knots";
+            case MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT:
+                return "B-spline fitting error";
+#endif
+            default:
+                return "UNKNOWN";
+        }
+    }
+
+    static const MccDefaultPCMConstructorErrorCategory& get()
+    {
+        static const MccDefaultPCMConstructorErrorCategory constInst;
+        return constInst;
+    }
+};
+
+
+inline std::error_code make_error_code(MccDefaultPCMConstructorErrorCode ec)
+{
+    return std::error_code(static_cast<int>(ec), MccDefaultPCMConstructorErrorCategory::get());
+}
+
+}  // namespace mcc::impl
+
+
+namespace std
+{
+
+template <>
+class is_error_code_enum<mcc::impl::MccDefaultPCMConstructorErrorCode> : public true_type
+{
+};
+
+
+}  // namespace std
+
+namespace mcc::impl
+{
+
+template <MccMountType MOUNT_TYPE>
+class MccPCMConstructor
+{
+public:
+    using ref_coordpair_t = std::conditional_t<mcc_is_equatorial_mount<MOUNT_TYPE>,
+                                               MccSkyHADEC_OBS,
+                                               std::conditional_t<mcc_is_altaz_mount<MOUNT_TYPE>, MccSkyAZZD, void>>;
+
+    static_assert(!std::is_void_v<ref_coordpair_t>, "UNSUPPORTED MOUNT TYPE!");
+
+    struct table_elem_t {
+        double target_colon, target_colat;
+        double hw_colon, hw_colat;
+        double colon_res, colat_res;  // target - hw
+    };
+
+
+    struct table_t {
+        std::vector<double> target_colon, target_colat;
+        std::vector<double> hw_colon, hw_colat;
+        std::vector<double> colon_res, colat_res;  // target - hw
+    };
+
+    struct compute_result_t {
+        MccDefaultPCMType pcm_type;
+        MccError error{};  // final model computation error
+
+        std::vector<double> model_colon, model_colat;  // fitted model values
+        std::vector<double> colon_res, colat_res;      // target - model
+
+#ifdef USE_BSPLINE_PCM
+        int bspline_fit_err{};  // bivariate B-spline fitting exit code (see FITPACK)
+
+        // quantities below are computed only fo pcm_type == MccDefaultPCMType::PCM_TYPE_BSPLINE
+        std::vector<double> inv_model_colon, inv_model_colat;  // fitted inverse model values
+        std::vector<double> inv_colon_res, inv_colat_res;      // encoder - model
+#endif
+    };
+
+    MccError addPoint(mcc_skypoint_c auto target_coords,
+                      mcc_coord_pair_c auto const& hw_counts,
+                      mcc_coord_epoch_c auto const& ref_epoch)
+        requires(decltype(hw_counts)::pairKind == MccCoordPairKind::COORDS_KIND_XY)
+    {
+        auto err = target_coords.to(ref_coordpair_t::pairKind, ref_epoch);
+        if (err) {
+            return mcc_deduced_err(err, MccDefaultPCMErrorCode::ERROR_CCTE);
+        }
+
+        // _table.push_back({target_coords.co_lon(), target_coords.co_lat(), hw_counts.x(), hw_counts.y(),
+        //                   target_coords.co_lon() - hw_counts.x(), target_coords.co_lat() - hw_counts.y()});
+
+        _table.target_colon.emplace_back(target_coords.co_lon());
+        _table.target_colat.emplace_back(target_coords.co_lon());
+
+        _table.hw_colon.emplace_back(hw_counts.x());
+        _table.hw_colat.emplace_back(hw_counts.y());
+
+        _table.colon_res.emplace_back(target_coords.co_lon() - hw_counts.x());
+        _table.colat_res.emplace_back(target_coords.co_lat() - hw_counts.y());
+
+        return {};
+    }
+
+
+    MccError addPoint(mcc_skypoint_c auto target_coords, mcc_coord_pair_c auto const& hw_counts)
+    {
+        auto ref_epoch = hw_counts.epoch();
+
+        return addPoint(std::move(target_coords), hw_counts, ref_epoch);
+    }
+
+
+    size_t numberOfPoints() const
+    {
+        return _table.colon_res.size();
+    }
+
+    void deletePoints()
+    {
+        _table.target_colon.clear();
+        _table.target_colat.clear();
+
+        _table.hw_colon.clear();
+        _table.hw_colat.clear();
+
+        _table.colon_res.clear();
+        _table.colat_res.clear();
+    }
+
+
+    //
+    // for B-splines interior knots along axes must be given in 'pcm_data'
+    // NOTE: the size of the interior knots array must be at least 2 as
+    //       it are interpretated as border knots and final full knots set is:
+    //       knots = [input_knots[0], input_knots[0], input_knots[0], input_knots[0], input_knots[1], input_knots[2],
+    //       ..., input_knots[N-1], input_knots[N-1], input_knots[N-1], input_knots[N-1]], where N = input_knots.size()
+    //
+    // WARNING: the input knots for inverse B-spline are ignored so the direct and inverse B-spline coefficients are
+    //          calculated on the same mesh!
+    compute_result_t computeModel(MccDefaultPCM<MOUNT_TYPE>::pcm_data_t& pcm_data)
+    {
+        compute_result_t result{.pcm_type = pcm_data.type, .error = MccDefaultPCMConstructorErrorCode::ERROR_OK};
+
+        size_t min_data_size = 2;  // 2 is for BSPLINE
+
+        if (pcm_data.type == MccDefaultPCMType::PCM_TYPE_GEOMETRY
+#ifdef USE_BSPLINE_PCM
+            || pcm_data.type == MccDefaultPCMType::PCM_TYPE_GEOMETRY_BSPLINE
+#endif
+        ) {
+            if constexpr (MOUNT_TYPE == MccMountType::FORK_TYPE) {
+                min_data_size = 9;
+            } else {
+                min_data_size = 8;
+            }
+
+            if (_table.size() < min_data_size) {
+                result.error = MccDefaultPCMConstructorErrorCode::ERROR_NOT_ENOUGH_DATA;
+                return result;
+            }
+
+            // robust linear regression with Tukey's loss function
+        }
+
+
+#ifdef USE_BSPLINE_PCM
+        if (pcm_data.type == MccDefaultPCMType::PCM_TYPE_BSPLINE ||
+            pcm_data.type == MccDefaultPCMType::PCM_TYPE_GEOMETRY_BSPLINE) {
+            if (pcm_data.bspline.knotsX.size() < 2 || pcm_data.bspline.knotsY.size() < 2) {
+                return MccDefaultPCMConstructorErrorCode::ERROR_INVALID_KNOTS_NUMBER;
+            }
+
+            double resi2x, resi2y;  // fitting residuals
+
+            std::vector<double> tx(pcm_data.bspline.knotsX.size() + 6), ty(pcm_data.bspline.knotsY.size() + 6);
+
+            size_t Ncoeffs = (tx.size() - 4) * (ty.size() - 4);
+
+            pcm_data.bspline.coeffsX.resize(Ncoeffs);
+            pcm_data.bspline.coeffsY.resize(Ncoeffs);
+
+            if (pcm_data.type == MccDefaultPCMType::PCM_TYPE_BSPLINE) {
+                if (_table.size() < min_data_size) {
+                    result.error = MccDefaultPCMConstructorErrorCode::ERROR_NOT_ENOUGH_DATA;
+                    return result;
+                }
+
+                // here both direct and inverse coefficients will be calculated
+                pcm_data.inverseBspline.coeffsX.resize(Ncoeffs);
+                pcm_data.inverseBspline.coeffsY.resize(Ncoeffs);
+
+                // direct (celestial = encoder + pcm)
+                result.bspline_fit_err = bsplines::fitpack_sphere_fit(
+                    _table.hw_colat, _table.hw_colon, _table.colon_res, 1.0, pcm_data.bspline.knotsY,
+                    pcm_data.bspline.knotsX, pcm_data.bspline.coeffsX, resi2x);
+                if (result.bspline_fit_err > 0) {
+                    result.error = MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT;
+                    return result;
+                }
+
+                result.bspline_fit_err = bsplines::fitpack_sphere_fit(
+                    _table.hw_colat, _table.hw_colon, _table.colat_res, 1.0, pcm_data.bspline.knotsY,
+                    pcm_data.bspline.knotsX, pcm_data.bspline.coeffsY, resi2y);
+                if (result.bspline_fit_err > 0) {
+                    result.error = MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT;
+                    return result;
+                }
+
+
+                // inverse (encoder = celestial + pcm)
+                std::vector<double> colon_res = _table.colon_res;
+                std::vector<double> colat_res = _table.colat_res;
+                for (size_t i = 0; i < colat_res.size(); ++i) {
+                    colon_res[i] = -colon_res[i];
+                    colat_res[i] = -colat_res[i];
+                }
+
+                result.bspline_fit_err = bsplines::fitpack_sphere_fit(
+                    _table.target_colon, _table.target_colat, colon_res, 1.0, pcm_data.bspline.knotsY,
+                    pcm_data.bspline.knotsX, pcm_data.bspline.inverseCoeffsX, resi2x);
+                if (result.bspline_fit_err > 0) {
+                    result.error = MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT;
+                    return result;
+                }
+
+                result.bspline_fit_err = bsplines::fitpack_sphere_fit(
+                    _table.target_colon, _table.target_colat, colat_res, 1.0, pcm_data.bspline.knotsY,
+                    pcm_data.bspline.knotsX, pcm_data.bspline.inverseCoeffsY, resi2y);
+                if (result.bspline_fit_err > 0) {
+                    result.error = MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT;
+                    return result;
+                }
+            } else {  // geometry + B-spline
+                // the fitting for geometrical coefficients is already done above so
+                // one must fit residuals by bivariate B-splines
+
+                std::vector<double> xres(_table.size()), yres(_table.size());
+                // for (size_t i = 0; i < _table.size(); ++i) {
+                //     xres = _table[i].target_colon;
+                //     yres = _table[i].target_colat;
+                // }
+                for (size_t i = 0; i < _table.size(); ++i) {
+                    xres = _table.target_colon[i];
+                    yres = _table.target_colat[i];
+                }
+            }
+        }
+    }
+#endif
+
+protected:
+    // std::vector<table_elem_t> _table;
+    table_t _table;
+
+    compute_result_t bsplineFitting(MccDefaultPCM<MOUNT_TYPE>::pcm_data_t& pcm_data)
+    {
+        compute_result_t result{.pcm_type = pcm_data.type, .error = MccDefaultPCMConstructorErrorCode::ERROR_OK};
+
+        if (pcm_data.bspline.knotsX.size() < 2 || pcm_data.bspline.knotsY.size() < 2) {
+            return MccDefaultPCMConstructorErrorCode::ERROR_INVALID_KNOTS_NUMBER;
+        }
+
+        double resi2x, resi2y;  // fitting residuals
+
+        std::vector<double> tx(pcm_data.bspline.knotsX.size() + 6), ty(pcm_data.bspline.knotsY.size() + 6);
+
+        size_t Ncoeffs = (tx.size() - 4) * (ty.size() - 4);
+
+        pcm_data.bspline.coeffsX.resize(Ncoeffs);
+        pcm_data.bspline.coeffsY.resize(Ncoeffs);
+
+        if (pcm_data.type == MccDefaultPCMType::PCM_TYPE_BSPLINE) {
+            // here both direct and inverse coefficients will be calculated
+            pcm_data.inverseBspline.coeffsX.resize(Ncoeffs);
+            pcm_data.inverseBspline.coeffsY.resize(Ncoeffs);
+
+            // direct (celestial = encoder + pcm)
+            result.bspline_fit_err = bsplines::fitpack_sphere_fit(_table.hw_colat, _table.hw_colon, _table.colon_res,
+                                                                  1.0, pcm_data.bspline.knotsY, pcm_data.bspline.knotsX,
+                                                                  pcm_data.bspline.coeffsX, resi2x);
+            if (result.bspline_fit_err > 0) {
+                result.error = MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT;
+                return result;
+            }
+
+            result.bspline_fit_err = bsplines::fitpack_sphere_fit(_table.hw_colat, _table.hw_colon, _table.colat_res,
+                                                                  1.0, pcm_data.bspline.knotsY, pcm_data.bspline.knotsX,
+                                                                  pcm_data.bspline.coeffsY, resi2y);
+            if (result.bspline_fit_err > 0) {
+                result.error = MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT;
+                return result;
+            }
+
+
+            // inverse (encoder = celestial + pcm)
+            std::vector<double> colon_res = _table.colon_res;
+            std::vector<double> colat_res = _table.colat_res;
+            for (size_t i = 0; i < colat_res.size(); ++i) {
+                colon_res[i] = -colon_res[i];
+                colat_res[i] = -colat_res[i];
+            }
+
+            result.bspline_fit_err = bsplines::fitpack_sphere_fit(_table.target_colon, _table.target_colat, colon_res,
+                                                                  1.0, pcm_data.bspline.knotsY, pcm_data.bspline.knotsX,
+                                                                  pcm_data.bspline.inverseCoeffsX, resi2x);
+            if (result.bspline_fit_err > 0) {
+                result.error = MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT;
+                return result;
+            }
+
+            result.bspline_fit_err = bsplines::fitpack_sphere_fit(_table.target_colon, _table.target_colat, colat_res,
+                                                                  1.0, pcm_data.bspline.knotsY, pcm_data.bspline.knotsX,
+                                                                  pcm_data.bspline.inverseCoeffsY, resi2y);
+            if (result.bspline_fit_err > 0) {
+                result.error = MccDefaultPCMConstructorErrorCode::ERROR_BSPLINE_FIT;
+                return result;
+            }
+        }
+    }
+};
+
+}  // namespace mcc::impl

include/mcc/mcc_telemetry.h

@@ -318,6 +318,8 @@ public:
 
     error_t setPointingTarget(mcc_skypoint_c auto const& sp)
     {
+        std::lock_guard lock{*_updateMutex};
+
         _enteredTargetPos = sp;
 
         return MccTelemetryErrorCode::ERROR_OK;
@@ -326,16 +328,20 @@ public:
 
     auto getPointingTarget() const
     {
+        std::lock_guard lock{*_updateMutex};
+
         return _enteredTargetPos;
     }
 
 
-    error_t getPointingTarget(mcc_skypoint_c auto* sp)
-    {
-        sp = _enteredTargetPos;
+    // error_t getPointingTarget(mcc_skypoint_c auto* sp)
+    // {
+    //     if (sp) {
+    //         *sp = _enteredTargetPos;
+    //     }
 
-        return MccTelemetryErrorCode::ERROR_OK;
-    }
+    //     return MccTelemetryErrorCode::ERROR_OK;
+    // }
 
 
     //

tests/mcc_telemetry_test.cpp

@@ -14,6 +14,9 @@ static std::uniform_real_distribution<double> y_distrib(1000, 7000);
 static std::uniform_real_distribution<double> xs_distrib(100, 700);
 static std::uniform_real_distribution<double> ys_distrib(100, 700);
 
+static std::uniform_int_distribution<uint8_t> err_distrib(0, 255);
+
+
 struct hw_t {
     static constexpr mcc::MccMountType hwMountType{mcc::MccMountType::ALTAZ_TYPE};
     static constexpr std::string_view hardwareName{"HW-TEST"};
@@ -53,13 +56,19 @@ struct hw_t {
 
         *state = hardware_state_t{.XY{x_distrib(gen), y_distrib(gen)}, .speedXY{xs_distrib(gen), ys_distrib(gen)}};
 
-        return 0;
+        // return 0;
+        return err_distrib(gen) < 5 ? 1 : 0;
     }
 
     error_t hardwareInit()
     {
         return 0;
     }
+
+    error_t hardwareShutdown()
+    {
+        return 0;
+    }
 };
 
 template <>