packages/TrackingDiagnostics/TrkrNtuplizer.h

0001 #ifndef TRACKRECO_TRKRNTUPLIZER_H
0002 #define TRACKRECO_TRKRNTUPLIZER_H
0003
0004 //===============================================
0005 /// \file SvtxEvaluator.h
0006 /// \brief Compares reconstructed tracks to truth particles
0007 /// \author Michael P. McCumber (revised SVTX version)
0008 //===============================================
0009
0010 #include <fun4all/SubsysReco.h>
0011 #include <trackbase/ClusterErrorPara.h>
0012 #include <trackbase/TrkrDefs.h>
0013
0014 #include <TMatrixFfwd.h>
0015 #include <TMatrixT.h>
0016 #include <TMatrixTUtils.h>
0017
0018 #include <limits>
0019 #include <map>
0020 #include <queue>
0021 #include <set>
0022 #include <string>
0023
0024 class CDBInterface;
0025 class CDBTTree;
0026 class PHCompositeNode;
0027 class PHTimer;
0028 class TrkrCluster;
0029 class TFile;
0030 class TNtuple;
0031 class SvtxTrack;
0032 class TrackSeed;
0033 class SvtxTrackMap;
0034 class SvtxVertexMap;
0035 class TrkrClusterContainer;
0036 class ActsGeometry;
0037 class PHG4TpcGeomContainer;
0038 class GlobalVertexMap;
0039
0040 // class ClusterErrorPara;
0041
0042 /// \class SvtxEvaluator
0043 ///
0044 /// \brief Compares reconstructed tracks to truth particles
0045 ///
0046 /// Plan: This module will trace the reconstructed tracks back to
0047 /// the greatest contributor Monte Carlo particle and then
0048 /// test one against the other.
0049 ///
0050
0051 class TrkrNtuplizer : public SubsysReco
0052 {
0053  public:
0054   TrkrNtuplizer(const std::string &name = "TRKRNTUPLEIZER",
0055                 const std::string &filename = "trkrntuple.root",
0056                 const std::string &trackmapname = "SvtxTrackMap",
0057                 unsigned int nlayers_maps = 3,
0058                 unsigned int nlayers_intt = 4,
0059                 unsigned int nlayers_tpc = 48,
0060                 unsigned int nlayers_mms = 2);
0061   ~TrkrNtuplizer() override;
0062
0063   int Init(PHCompositeNode *topNode) override;
0064   int InitRun(PHCompositeNode *topNode) override;
0065   int process_event(PHCompositeNode *topNode) override;
0066   int End(PHCompositeNode *topNode) override;
0067   //  void do_primaries(bool b);
0068
0069   void do_info_eval(bool b) { _do_info_eval = b; }
0070   void do_vertex_eval(bool b) { _do_vertex_eval = b; }
0071   void do_hit_eval(bool b) { _do_hit_eval = b; }
0072   void do_cluster_eval(bool b) { _do_cluster_eval = b; }
0073   void do_clus_trk_eval(bool b) { _do_clus_trk_eval = b; }
0074   void do_track_eval(bool b) { _do_track_eval = b; }
0075   void do_dedx_calib(bool b) { _do_dedx_calib = b; }
0076   void do_tpcseed_eval(bool b) { _do_tpcseed_eval = b; }
0077   void do_siseed_eval(bool b) { _do_siseed_eval = b; }
0078   void set_first_event(int value) { _ievent = value; }
0079   void set_trkclus_seed_container(const std::string &name)
0080   {
0081     _clustrackseedcontainer = name;
0082   }
0083   SvtxTrack *best_track_from(TrkrDefs::cluskey cluster_key);
0084   std::set<SvtxTrack *> all_tracks_from(TrkrDefs::cluskey cluster_key);
0085   void create_cache_track_from_cluster();
0086   static std::vector<TrkrDefs::cluskey> get_track_ckeys(SvtxTrack *track);
0087   void segment(const int seg) { m_segment = seg; }
0088   void runnumber(const int run) { m_runnumber = run; }
0089   void job(const int job) { m_job = job; }
0090
0091  private:
0092   struct fee_info
0093   {
0094     unsigned int fee = std::numeric_limits<unsigned int>::max();
0095     unsigned int channel = std::numeric_limits<unsigned int>::max();
0096     unsigned int sampa = std::numeric_limits<unsigned int>::max();
0097   };
0098   int m_segment = 0;
0099   int m_runnumber = 0;
0100   int m_job = 0;
0101
0102   unsigned int _ievent{0};
0103   unsigned int _iseed{0};
0104   float m_fSeed{std::numeric_limits<float>::quiet_NaN()};
0105   // eval stack
0106
0107   TF1 *f_pion_plus{nullptr};
0108   TF1 *f_kaon_plus{nullptr};
0109   TF1 *f_proton_plus{nullptr};
0110   TF1 *f_pion_minus{nullptr};
0111   TF1 *f_kaon_minus{nullptr};
0112   TF1 *f_proton_minus{nullptr};
0113   float dedxcorr[2][12][3]{};
0114   float get_n1pix(TrackSeed *tpcseed);
0115
0116   static TMatrixF calculateClusterError(TrkrCluster *c, float &clusphi);
0117   void get_dca(SvtxTrack *track, SvtxVertexMap *vertexmap,
0118                float &dca3dxy, float &dca3dz,
0119                float &dca3dxysigma, float &dca3dzsigma);
0120   // TrkrClusterContainer *cluster_map{nullptr};
0121
0122   void FillCluster(Float_t fXcluster[30], TrkrDefs::cluskey cluster_key);
0123   void FillTrack(Float_t fXcluster[30], SvtxTrack *track, GlobalVertexMap *vertexmap);
0124   //----------------------------------
0125   // evaluator output ntuples
0126
0127   bool _do_info_eval{true};
0128   bool _do_vertex_eval{true};
0129   bool _do_hit_eval{true};
0130   bool _do_cluster_eval{true};
0131   bool _do_clus_trk_eval{true};
0132   bool _do_track_eval{true};
0133   bool _do_dedx_calib{false};
0134   bool _do_tpcseed_eval{false};
0135   bool _do_siseed_eval{false};
0136
0137   unsigned int _nlayers_maps{3};
0138   unsigned int _nlayers_intt{4};
0139   unsigned int _nlayers_tpc{48};
0140   unsigned int _nlayers_mms{2};
0141
0142   TNtuple *_ntp_info{nullptr};
0143   TNtuple *_ntp_vertex{nullptr};
0144   TNtuple *_ntp_hit{nullptr};
0145   TNtuple *_ntp_cluster{nullptr};
0146   TNtuple *_ntp_clus_trk{nullptr};
0147   TNtuple *_ntp_track{nullptr};
0148   TNtuple *_ntp_tpcseed{nullptr};
0149   TNtuple *_ntp_siseed{nullptr};
0150
0151   // evaluator output file
0152   std::string _filename;
0153   // Track map name
0154   std::string _trackmapname;
0155   ClusterErrorPara _ClusErrPara;
0156   TrkrClusterContainer *_cluster_map{nullptr};
0157   SvtxTrackMap *_trackmap{nullptr};
0158   ActsGeometry *_tgeometry{nullptr};
0159   PHG4TpcGeomContainer *_geom_container{nullptr};
0160   float m_ZDC_coincidence{0};
0161   float m_mbd_rate{0};
0162   float m_rawzdc{0};
0163   float m_livezdc{0};
0164   float m_scaledzdc{0};
0165   float m_rawmbd{0};
0166   float m_livembd{0};
0167   float m_scaledmbd{0};
0168   float m_rawmbdv10{0};
0169   float m_livembdv10{0};
0170   float m_scaledmbdv10{0};
0171   float m_mbd_rate1{0};
0172   float m_rawzdc1{0};
0173   float m_livezdc1{0};
0174   float m_scaledzdc1{0};
0175   float m_rawmbd1{0};
0176   float m_livembd1{0};
0177   float m_scaledmbd1{0};
0178   float m_rawmbdv101{0};
0179   float m_livembdv101{0};
0180   float m_scaledmbdv101{0};
0181   float m_rawzdclast{0};
0182   float m_rawmbdlast{0};
0183   float m_rawmbdv10last{0};
0184   uint64_t m_bcolast = std::numeric_limits<uint64_t>::quiet_NaN();
0185
0186   std::queue<int> m_rawzdc_hist;
0187   std::queue<int> m_rawmbd_hist;
0188   std::queue<int> m_rawmbdv10_hist;
0189   std::queue<uint64_t> m_bco_hist;
0190
0191   uint64_t m_bco = std::numeric_limits<uint64_t>::quiet_NaN();
0192   uint64_t m_bco1 = std::numeric_limits<uint64_t>::quiet_NaN();
0193   uint64_t m_bcotr = std::numeric_limits<uint64_t>::quiet_NaN();
0194   uint64_t m_bcotr1 = std::numeric_limits<uint64_t>::quiet_NaN();
0195   float m_totalmbd = std::numeric_limits<float>::quiet_NaN();
0196   std::vector<int> m_firedTriggers;
0197   uint64_t m_gl1BunchCrossing = std::numeric_limits<uint64_t>::quiet_NaN();
0198
0199   std::string _clustrackseedcontainer = "TpcTrackSeedContainer";
0200
0201   TFile *_tfile{nullptr};
0202   PHTimer *_timer{nullptr};
0203
0204   // output subroutines
0205   void fillOutputNtuples(PHCompositeNode *topNode);  ///< dump the evaluator information into ntuple for external analysis
0206   void printInputInfo(PHCompositeNode *topNode);     ///< print out the input object information (debugging upstream components)
0207   void printOutputInfo(PHCompositeNode *topNode);    ///< print out the ancestry information for detailed diagnosis
0208   double AdcClockPeriod = 53.0;                      // ns
0209
0210   bool _cache_track_from_cluster_exists = false;
0211   std::map<TrkrDefs::cluskey, std::set<SvtxTrack *> > _cache_all_tracks_from_cluster;
0212   std::map<TrkrDefs::cluskey, SvtxTrack *> _cache_best_track_from_cluster;
0213
0214   CDBTTree *m_cdbttree{nullptr};
0215   CDBInterface *m_cdb{nullptr};
0216
0217   int mc_sectors[12]{5, 4, 3, 2, 1, 0, 11, 10, 9, 8, 7, 6};
0218   int FEE_map[26]{4, 5, 0, 2, 1, 11, 9, 10, 8, 7, 6, 0, 1, 3, 7, 6, 5, 4, 3, 2, 0, 2, 1, 3, 5, 4};
0219   int FEE_R[26]{2, 2, 1, 1, 1, 3, 3, 3, 3, 3, 3, 2, 2, 1, 2, 2, 1, 1, 2, 2, 3, 3, 3, 3, 3, 3};
0220
0221   std::map<TrkrDefs::cluskey, fee_info> fee_map;
0222 };
0223
0224 #endif  // G4EVAL_SVTXEVALUATOR_H