R3BNeulandMilleCalDataProcessor.h

Go to the documentation of this file.

#pragma once
 
#include "R3BHuberRegression.h"
#include "R3BValueError.h"
#include <Fit/Fitter.h>
#include <R3BNeulandCalData2.h>
#include <R3BNeulandCommon.h>
#include <Rtypes.h>
#include <TF1.h>
#include <cstdlib>
#include <fmt/base.h>
#include <fmt/core.h>
#include <unordered_map>
#include <vector>
 
namespace R3B::Neuland
{
    class Cal2HitPar;
 
}
 
namespace R3B::Neuland::Calibration
{

    struct MilleCalData
    {
      public:
        MilleCalData() = default;
        explicit MilleCalData(const BarCalData& bar_cal_data)
            : module_num{ bar_cal_data.module_num }
            , left{ bar_cal_data.left.front() }
            , right{ bar_cal_data.right.front() }
        {
        }
        bool is_outlier = false;      
        int module_num = 0;           
        float fit_diff = 0.F;         
        float residual = -1.;         
        float residual_bar_pos = -1.; 
        float tsync_meas = 0.F;
        ValueErrorD position; 
        ValueErrorD z_pos;    
        CalDataSignal left;
        CalDataSignal right;
 
        ClassDefNV(MilleCalData, 1);
    };
 
    struct TrackFitPar
    {
        double offset = 0.;
        double slope = 0.;
        double p_value = 0.;
 
        static constexpr auto diff_threshold = 0.1; 
 
        auto operator!=(const TrackFitPar& other) const
        {
            return (std::abs(other.offset - offset) / offset > diff_threshold) or
                   (std::abs(other.slope - slope) / slope > diff_threshold);
        }
        void clear()
        {
            offset = 0.;
            slope = 0.;
            p_value = 0.;
        }
        ClassDefNV(TrackFitPar, 1);
    };
 
    struct TrackFitResult
    {
        TrackFitPar x_z;
        TrackFitPar y_z;
        void clear()
        {
            x_z.clear();
            y_z.clear();
        }
        ClassDefNV(TrackFitResult, 1);
    };
 
    struct NeulandTrackInfo
    {
        TrackFitResult time_data;
        TrackFitResult bar_disp_data;
        void clear()
        {
            time_data.clear();
            bar_disp_data.clear();
        }
        ClassDefNV(NeulandTrackInfo, 1);
    };
 
    struct NeulandTrackDataSet
    {
        std::vector<double> z_vals;
        std::vector<double> z_errs;
        std::vector<double> vals;
        std::vector<double> errs;
 
        void clear()
        {
            z_vals.clear();
            z_errs.clear();
            vals.clear();
            errs.clear();
        }
        [[nodiscard]] auto size() const { return z_vals.size(); }
    };
 
    class MilleDataProcessor
    {
      public:
        using FitPar = TrackFitPar;
        using FitResult = TrackFitResult;
        using TrackInfo = NeulandTrackInfo;
 
        struct TrackFitDataSet
        {
            NeulandTrackDataSet bar_x_z;
            NeulandTrackDataSet bar_y_z;
            NeulandTrackDataSet time_data_x_z;
            NeulandTrackDataSet time_data_y_z;
 
            void clear()
            {
                bar_x_z.clear();
                bar_y_z.clear();
                time_data_x_z.clear();
                time_data_y_z.clear();
            }
        };
 
        explicit MilleDataProcessor(int num_of_modules);

        auto process(const std::vector<BarCalData>& signals, const Cal2HitPar& hit_par) -> bool;
 
        auto check_is_outlier(int module_num) const -> bool;
 
        void reset();
 
        void set_p_value_cut(double val) { p_value_cut_ = val; }
 
        [[nodiscard]] auto get_data() const -> const auto& { return data_buffers_; }
        [[nodiscard]] auto get_track_info() const -> const auto& { return track_info_; }
        [[nodiscard]] auto get_track_fit_data() const -> const auto& { return track_fit_data_; }
        // [[nodiscard]] auto calculate_residual(double val, int module_num) const -> float;
 
        static auto linear_fit(const NeulandTrackDataSet& data,
                               FitPar& fit_par,
                               HuberRegressor& huber_regressor,
                               double p_value_cut) -> bool;
 
      private:
        double p_value_cut_ = DEFAULT_CALIBRATION_P_VALUE_CUT;

        using DataBufferType = std::unordered_map<int, std::vector<MilleCalData>>;
        std::unordered_map<int, std::vector<MilleCalData>> data_buffers_;
        TrackInfo track_info_;
        TrackFitDataSet track_fit_data_;
        HuberRegressor huber_regressor_;
 
        void init_data_registers(int num_of_modules);

        void remove_isolated_bar_signal();
        void fill_bar_disp_dataset();
        void fill_time_dataset();
        void reset_fit_pars();
        auto fit_planes(const Cal2HitPar& hit_par) -> bool;
        void check_fit_result();
        void set_data_buffers(DataBufferType& data_buffers,
                              const TrackInfo& track_info,
                              const Cal2HitPar& hit_par) const;
    };
} // namespace R3B::Neuland::Calibration
 
#ifndef __CLING__
template <>
class fmt::formatter<R3B::Neuland::Calibration::MilleCalData>
{
  public:
    static constexpr auto parse(format_parse_context& ctx) { return ctx.end(); }
    template <typename FmtContent>
    constexpr auto format(const R3B::Neuland::Calibration::MilleCalData& signal, FmtContent& ctn) const
    {
        return fmt::format_to(
            ctn.out(), "ModuleNum: {}, left bar: {}, right bar: {}", signal.module_num, signal.left, signal.right);
    }
};
 
template <>
class fmt::formatter<R3B::Neuland::Calibration::MilleDataProcessor::FitPar>
{
  public:
    static constexpr auto parse(format_parse_context& ctx) { return ctx.end(); }
    template <typename FmtContent>
    constexpr auto format(const R3B::Neuland::Calibration::MilleDataProcessor::FitPar& fit_par, FmtContent& ctn) const
    {
        return fmt::format_to(
            ctn.out(), "slope: {}, offset: {}, p_value: {}", fit_par.slope, fit_par.offset, fit_par.p_value);
    }
};
 
template <>
class fmt::formatter<R3B::Neuland::Calibration::MilleDataProcessor::FitResult>
{
  public:
    static constexpr auto parse(format_parse_context& ctx) { return ctx.end(); }
    template <typename FmtContent>
    constexpr auto format(const R3B::Neuland::Calibration::MilleDataProcessor::FitResult& fit_result,
                          FmtContent& ctn) const
    {
        return fmt::format_to(ctn.out(),
                              "\n"
                              "\tx_z plane: {{{}}}\n"
                              "\ty_z plane: {{{}}}",
                              fit_result.x_z,
                              fit_result.y_z);
    }
};
 
template <>
class fmt::formatter<R3B::Neuland::Calibration::NeulandTrackDataSet>
{
  public:
    static constexpr auto parse(format_parse_context& ctx) { return ctx.end(); }
    template <typename FmtContent>
    constexpr auto format(const R3B::Neuland::Calibration::NeulandTrackDataSet& fit_data, FmtContent& ctn) const
    {
        return fmt::format_to(ctn.out(),
                              "\tx = np.array({})\n"
                              "\tx_err = np.array({})\n"
                              "\ty = np.array({})\n"
                              "\ty_err = np.array({})",
                              fit_data.z_vals,
                              fit_data.z_errs,
                              fit_data.vals,
                              fit_data.errs);
    }
};
#endif