CbmRecoQaTask.h

Go to the documentation of this file.

/* Copyright (C) 2024-2026 Hulubei National Institute of Physics and Nuclear Engineering - Horia Hulubei, Bucharest
   SPDX-License-Identifier: GPL-3.0-only
   Authors: Alexandru Bercuci [committer], Omveer Singh */
 
#ifndef CBMRECOQATASK_H
#define CBMRECOQATASK_H 1
 
#include "CbmKfTrackFitter.h"
#include "CbmTrack.h"
 
#include <FairTask.h>
 
#include <TDirectoryFile.h>
#include <TGeoNode.h>
 
#include <bitset>
#include <map>
#include <vector>
 
class FairMCPoint;
class CbmCluster;
class CbmEvent;
class CbmHit;
class CbmMCDataManager;
class CbmMCDataArray;
class CbmTimeSlice;
class TClonesArray;
class TH1;
class TVector3;
class CbmRecoQaTask : public FairTask {
 public:
  enum eRecoConfig
  {
    kRecoEvents = 0,  
    kRecoTracks,      
    kStsHits,         
    kTrdHits,         
    kTofHits,         
    kRichHits,        
    kMuchHits,        
    kUseMC,           
    kUseSts,          
    kUseTrd,          
    kUseTof,          
    kRecoQaNConfigs   
  };
  enum eSetup
  {
    kMcbm22 = 0,
    kMcbm24,
    kMcbm25,
    kDefault
  };
  enum class eViewType : int
  {
    kDetUnit = 0,  
    kTrkProj,      
    kPV            
  };
 
#define kNtrkProjections 6
  enum class eProjectionType : int
  {
    kXYa = 0,     
    kXYp,         
    kXdX,         
    kXdXMC,       
    kXdXT,        
    kYdY,         
    kYdYMC,       
    kYdYT,        
    kWdT,         
    kXpX,         
    kYpY,         
    kChdT,        
    kXYh,         
    kDmult,       
    kDmultMC,     
    kXYhMC,       
    kPullX,       
    kPullY,       
    kResidualX,   
    kResidualY,   
    kResidualTX,  
    kResidualTY,  
    kPVxy,        
    kPVxz,        
    kPVyz,        
    kPVmult,      
    kPVxsx,       
    kPVysy,       
    kPVzsz,       
    kXYt0,        
    kXYt1,        
    kXYt2,        
    kXYt3,        
    kXYt4,        
    kXYt5,        
    kDebug0,      
    kDebug1       
  };
 
  class Trajectory : public CbmKfTrackFitter<>::Trajectory {
   public:
    using CbmKfTrackFitter<>::Trajectory::Trajectory;
 
    static int GetPlaneId(const Node& node) { return node.userReferences[3]; }
    static void SetPlaneId(Node& node, int planeId) { node.userReferences[3] = planeId; }
    void SetPlaneId(size_t inode, int planeId) { SetPlaneId(fNodes[inode], planeId); }
  };
 
  // Generic view definition
  struct Detector;
  struct View {
    friend struct Detector;
    friend class CbmRecoQaTask;
 
    View() = default;
    View(const char* n, const char* p, std::vector<int> set) : name(n), path(p), fSelector(set) { ; }
    virtual ~View() = default;
    bool SetProjection(eProjectionType prj, float range, const char* unit);
    template<class Hit>
    bool HasAddress(const CbmHit* h, double& x, double& y, double& dx, double& dy) const;
    template<class Hit>
    bool Load(const CbmHit* h, const FairMCPoint* point, const CbmEvent* ev);
    bool Load(const Trajectory::Node* n, const FairMCPoint* point, const int dinfo = 0xab00);
    bool Load(TVector3* p, TVector3* pe = nullptr);
    std::string ToString() const;
    static std::string ToString(eProjectionType prj);
    void SetSetup(CbmRecoQaTask::eSetup setup) { fSetup = setup; }
    void SetType(CbmRecoQaTask::eViewType type) { fType = type; }
 
    std::string name           = "";    
    std::string path           = "";    
    double size[3]             = {0.};  
    double pos[3]              = {0.};  
    std::vector<int> fSelector = {};    
    std::map<eProjectionType, std::tuple<int, float, TH1*>> fProjection =
      {};                         

    int fMult               = 0;  
    mutable eSetup fSetup   = eSetup::kMcbm22;
    mutable eViewType fType = eViewType::kDetUnit;
 
   protected:
    bool AddProjection(eProjectionType prj, float range = -1, const char* unit = "cm");
    bool HasDebug() const;
    bool Init(const char* dname, bool mc = false);
    uint Register(TDirectoryFile* f);
    std::string makeTrange(const int scale, float& range);
    std::string makeYrange(const int scale, float& range);
    ClassDef(CbmRecoQaTask::View,
             1);  // Stand-alone detection set to which QA is applied
  };              // QA View definition
 
  // Detector unit definition
  struct Detector {
    struct Data {
      Data() = default;
      Data(ECbmDataType i, const char* n) : id(i), name(n) { ; }
      ECbmDataType id  = ECbmDataType::kUnknown;
      std::string name = "nn";
    };  // Identifier of data to be monitored
 
    Detector(ECbmModuleId did = ECbmModuleId::kNotExist);
    virtual ~Detector() = default;
    View* AddView(const char* n, const char* p, std::vector<int> set, bool debug = false);
    View* GetView(const char* n);
    View* FindView(double x, double y, double z);
    bool Init(TDirectoryFile* f, bool mc = false);
    void Print() const;
 
    ECbmModuleId id;
    Data cl, hit, point, trk;
    std::vector<View> fViews = {};
 
    ClassDef(CbmRecoQaTask::Detector,
             1);  // QA representation of a detector unit
  };              // Detection system agregate

  static constexpr size_t projSz = 4;
  static constexpr size_t stsSz  = 3;
  static constexpr size_t tofSz  = 3;
  static constexpr size_t trdSz  = 7;
  struct TrackFilter;
  struct TopoFilter {
    friend struct TrackFilter;
    TopoFilter() { ; }
    virtual ~TopoFilter() = default;
    template<ECbmModuleId d>
    bool Accept(const CbmHit* h, const CbmCluster* c0 = nullptr, const CbmCluster* c1 = nullptr) const;
 
   protected:
    std::vector<std::vector<int>> fAccept = {};
    std::vector<std::vector<int>> fReject =
      {};  
    std::vector<std::array<float, projSz>> fRef = {};  
    ClassDef(CbmRecoQaTask::TopoFilter, 1);            // Topological cuts for hits/track
  };
 
  struct TrackFilter {
    TrackFilter(ECbmModuleId d = ECbmModuleId::kNotExist) : fDet(d) { ; }
    virtual ~TrackFilter() = default;
    bool Accept(const CbmGlobalTrack* ptr, const CbmRecoQaTask* lnk) const;
    bool AddCondAccept(std::vector<int> cuts) { return AddTopoCut(cuts, 'a'); }
    bool AddCondReject(std::vector<int> cuts) { return AddTopoCut(cuts, 'r'); }
    const std::vector<std::vector<int>>& GetAcceptCuts() const { return fHitsFilter.fAccept; }
    const std::vector<std::vector<int>>& GetRejectCuts() const { return fHitsFilter.fReject; }
    const std::vector<std::array<float, projSz>>& GetProjCuts() const { return fProjFilter.fRef; }
    const TopoFilter* GetHitsFilter() const
    {
      if (!fHitsFilter.fAccept.size() && !fHitsFilter.fReject.size()) return nullptr;
      return &fHitsFilter;
    }

    bool SetProjFilter(std::vector<float> cuts);
    void SetMinHits(int n) { fNminHits = n; }
    std::string ToString() const;
 
    ECbmModuleId fDet = ECbmModuleId::kNotExist;
 
   private:
    bool AddTopoCut(std::vector<int> cuts, const char opt);
 
    int fNminHits = 0;       
    TopoFilter fHitsFilter;  
    TopoFilter fProjFilter;
    ClassDef(CbmRecoQaTask::TrackFilter, 1);  // Track cut implementation
  };                                          // TrackFilter definition
 
  struct EventFilter {
    enum class eEventCut : int
    {
      kMultTrk = 0,  
      kMultHit,      
      kTrigger,      
      kVertex,       
      kNone          
    };
    enum class eEventDef : int
    {
      kNofDetHit = 3  
    };
    EventFilter(eEventCut typ = eEventCut::kNone) : fType(typ) { ; }
    virtual ~EventFilter() = default;
    bool Accept(const CbmEvent* ptr, const CbmRecoQaTask* lnk);
    bool SetFilter(std::vector<float> cuts);
    std::string ToString() const;
 
    eEventCut fType = eEventCut::kNone;
 
   private:
    // track multiplicity cut definition
    int fMinTrack = 0;
    int fMaxTrack = 0;
    // hit multiplicity cut definition
    int fMultHit[(int) eEventDef::kNofDetHit] = {0};  

    void HelpMess() const;
 
    ClassDef(CbmRecoQaTask::EventFilter, 1);  // Event cut implementation
  };                                          // EventFilter definition
 
 
 public:
  CbmRecoQaTask();
  virtual ~CbmRecoQaTask() { ; }

  void Activate(ECbmModuleId det, bool set = true);
  virtual Detector* AddDetector(ECbmModuleId did);
  virtual EventFilter* AddEventFilter(EventFilter::eEventCut cut);
  virtual TrackFilter* AddTrackFilter(ECbmModuleId cut);
  virtual Detector* GetDetector(ECbmModuleId did);
  template<ECbmModuleId d>
  static size_t GetHitDetectorId(const CbmHit* h, std::vector<int>& sel);
  template<class Data>
  const Data* GetData(ECbmModuleId did, int id) const;
  virtual InitStatus Init();
  virtual void Exec(Option_t* option);
 
  virtual void Finish();
 
  // /** \brief Define the set of extra z positions where the track should be
  //  * projected in the x-y plane */
  // void SetProjections(std::vector<double> vzpj);
  void SetSetupClass(CbmRecoQaTask::eSetup setup) { fSetupClass = setup; }
  TString GetGeoTagForDetector(const TString& detector);
  std::vector<TString> GetPath(TGeoNode* node, TString, TString activeNodeName, int depth = 0,
                               const TString& path = "");
  void UseMC(bool set = true)
  {
    if (set) fuRecoConfig.set(kUseMC);
  }
 
 protected:
  static std::bitset<kRecoQaNConfigs> fuRecoConfig;
  bool IsActive(eRecoConfig feature) { return fuRecoConfig.test(feature); }
 
 private:
  CbmRecoQaTask(const CbmRecoQaTask&);
  CbmRecoQaTask& operator=(const CbmRecoQaTask&);

  template<ECbmModuleId d>
  int GetDebugInfo(const CbmHit* h);
  virtual bool FilterEvent(const CbmEvent* ptr) const;
  bool FilterHit(const CbmHit* h, ECbmModuleId d) const;
  virtual bool FilterTrack(const CbmGlobalTrack* ptr) const;
  virtual bool FilterTrack(const Trajectory* ptr) const;
  int GetNviews(eViewType type) const;
  void InitMcbm22();
  void InitMcbm24();
  void InitMcbm25();
  void InitDefault();
  template<class Data>
  std::map<ECbmModuleId, TClonesArray*> MapDataList(std::string& s) const;
 
  CbmKfTrackFitter<> fFitter;
  TClonesArray* fGTracks                          = nullptr;  
  std::map<ECbmModuleId, TClonesArray*> fTracks   = {};       
  TClonesArray* fTrackMatches                     = nullptr;  
  TClonesArray* fEvents                           = nullptr;  
  CbmTimeSlice* fTimeSlice                        = nullptr;  
  std::map<ECbmModuleId, TClonesArray*> fHits     = {};       
  std::map<ECbmModuleId, TClonesArray*> fClusters = {};       
  std::map<ECbmModuleId, CbmMCDataArray*> fPoints = {};       
  std::map<ECbmModuleId, TClonesArray*> fHitMatch = {};       
  std::vector<EventFilter> fFilterEv              = {};
  std::vector<TrackFilter> fFilterTrk             = {};
  CbmMCDataManager* cbm_mc_manager                = nullptr;
  TDirectoryFile fOutFolder                       = {"RecoQA", "CA track driven reco QA"};  
  eSetup fSetupClass                              = eSetup::kMcbm24;
  std::map<ECbmModuleId, Detector> fDetQa         = {};  
  std::map<const char*, View> fViews              = {};  
  std::vector<TVector3> fPrjPlanes = {};  
 
  ClassDef(CbmRecoQaTask, 1);  // Reconstruction QA analyzed from CA perspective
};
 
#endif  // CBMRECOQATASK_H