TrackJetQAMaker/src/TrksInJetQAInJetFiller.cc

0001 // ----------------------------------------------------------------------------
0002 // 'TrksInJetQAInJetFiller.cc'
0003 // Derek Anderson
0004 // 04.03.2024
0005 //
0006 // A submodule for the TrksInJetsQAM F4A module to produce
0007 // QA histograms for tracks and more in jets
0008 // ----------------------------------------------------------------------------
0009
0010 #define TRKSINJETQAINJETFILLER_CC
0011
0012 // submodule definition
0013 #include "TrksInJetQAInJetFiller.h"
0014
0015
0016
0017 // inherited public methods ---------------------------------------------------
0018
0019 void TrksInJetQAInJetFiller::Fill(PHCompositeNode* topNode) {
0020
0021   GetNodes(topNode);
0022
0023   FillJetAndTrackQAHists(topNode);
0024   return;
0025
0026 }  // end 'Fill(PHCompositeNode* topNode)'
0027
0028
0029
0030 // private methods ------------------------------------------------------------
0031
0032 void TrksInJetQAInJetFiller::GetNode(const int node, PHCompositeNode* topNode) {
0033
0034   // jump to relevant node
0035   switch (node) {
0036
0037     case Node::Flow:
0038       m_flowStore = findNode::getClass<ParticleFlowElementContainer>(topNode, "ParticleFlowElements");
0039       if (!m_flowStore) {
0040         std::cerr << PHWHERE << ": PANIC: Couldn't grab particle flow container from node tree!" << std::endl;
0041         assert(m_flowStore);
0042       }
0043       break;
0044
0045     default:
0046       std::cerr << PHWHERE << ": WARNING: trying to grab unkown additional node..." << std::endl;
0047       break;
0048   }
0049   return;
0050
0051 }  // end 'GetNode(int, PHCompositeNode*)'
0052
0053
0054
0055 void TrksInJetQAInJetFiller::FillJetAndTrackQAHists(PHCompositeNode* topNode) {
0056
0057   // loop over jets
0058   for (
0059     uint64_t iJet = 0;
0060     iJet < m_jetMap -> size();
0061     ++iJet
0062   ) {
0063
0064     // grab jet and make sure track vector is clear
0065     Jet* jet = m_jetMap -> get_jet(iJet);
0066     m_trksInJet.clear();
0067
0068     // get all tracks "in" a jet
0069     GetCstTracks(jet, topNode);
0070     GetNonCstTracks(jet);
0071
0072     // grab jet info and fill histograms
0073     if (m_config.doJetQA) m_jetManager -> GetInfo(jet, m_trksInJet);
0074
0075     // loop over tracks in the jet
0076     for (SvtxTrack* track : m_trksInJet) {
0077
0078       // grab track info and fill histograms
0079       if (m_config.doTrackQA) m_trackManager -> GetInfo(track);
0080
0081       // fill cluster and hit histograms as needed
0082       if (m_config.doClustQA || m_config.doHitQA) {
0083         FillClustAndHitQAHists(track);
0084       }
0085     }  // end track loop
0086   }  // end jet loop
0087   return;
0088
0089 }  // end 'FillJetAndTrackQAHists(PHCompositeNode*)'
0090
0091
0092
0093 void TrksInJetQAInJetFiller::FillClustAndHitQAHists(SvtxTrack* track) {
0094
0095   // get cluster keys
0096   for (auto clustKey : ClusKeyIter(track)) {
0097
0098     // grab cluster and its info
0099     if (m_config.doClustQA) {
0100       m_clustManager -> GetInfo(
0101         m_clustMap -> findCluster(clustKey),
0102         clustKey,
0103         m_actsGeom
0104       );
0105     }
0106
0107     // get hits if needed
0108     if (m_config.doHitQA) {
0109
0110       // grab hit set and key associated with cluster key
0111       TrkrDefs::hitsetkey setKey = TrkrDefs::getHitSetKeyFromClusKey(clustKey);
0112       TrkrHitSet*         set    = m_hitMap -> findHitSet(setKey);
0113       if (!set || !(set -> size() > 0)) return;
0114
0115       // loop over all hits in hit set
0116       TrkrHitSet::ConstRange hits = set -> getHits();
0117       for (
0118         TrkrHitSet::ConstIterator itHit = hits.first;
0119         itHit != hits.second;
0120         ++itHit
0121       ) {
0122
0123         // grab hit
0124         TrkrDefs::hitkey hitKey = itHit -> first;
0125         TrkrHit*         hit    = itHit -> second;
0126
0127         // grab info and fill histograms
0128         m_hitManager -> GetInfo(hit, setKey, hitKey);
0129
0130       }  // end hit loop
0131     }
0132   }  // end cluster key loop
0133   return;
0134
0135 }  // end 'FillClustQAHists(SvtxTrack*)'
0136
0137
0138
0139 void TrksInJetQAInJetFiller::GetCstTracks(Jet* jet, PHCompositeNode* topNode) {
0140
0141   // loop over consituents
0142   auto csts = jet -> get_comp_vec();
0143   for (
0144     auto cst = csts.begin();
0145     cst != csts.end();
0146     ++cst
0147   ) {
0148
0149     // ignore cst if non-relevent type
0150     const uint32_t src = cst -> first;
0151     if ( IsCstNotRelevant(src) ) continue;
0152
0153     // if cst is track, add to list
0154     if (src == Jet::SRC::TRACK) {
0155       m_trksInJet.push_back( m_trkMap -> get(cst -> second) );
0156     }
0157
0158     // if pfo, grab track if needed
0159     if (src == Jet::SRC::PARTICLE) {
0160       PFObject*  pfo   = GetPFObject(cst -> second, topNode);
0161       SvtxTrack* track = GetTrkFromPFO(pfo);
0162       if (track) {
0163         m_trksInJet.push_back( track );
0164       }
0165     }
0166   }  // end cst loop
0167   return;
0168
0169 }  // end 'GetCstTracks(Jet* jet, PHCompositeNode* topNode)'
0170
0171
0172
0173 void TrksInJetQAInJetFiller::GetNonCstTracks(Jet* jet) {
0174
0175   // loop over tracks
0176   for (
0177     SvtxTrackMap::Iter itTrk = m_trkMap -> begin();
0178     itTrk != m_trkMap -> end();
0179     ++itTrk
0180   ) {
0181
0182     // grab track
0183     SvtxTrack* track = itTrk -> second;
0184
0185     // ignore tracks we've already added to the list
0186     if ( IsTrkInList(track -> get_id()) ) continue;
0187
0188     // FIXME this can be improved!
0189     //   - jets don't necessarily have areas of
0190     //     pi*(Rjet)^2
0191     //   - it may be better to instead check
0192     //     if a track projection falls in
0193     //     a constituent tower/cluster
0194     //   - Also track projections to a calo
0195     //     would be better to use than just
0196     //     the track
0197     const double dr = GetTrackJetDist(track, jet);
0198     if (dr < m_config.rJet) {
0199       m_trksInJet.push_back( track );
0200     }
0201   }  // end track loop
0202   return;
0203
0204 }  // end 'GetNonCstTracks(Jet* jet)'
0205
0206
0207
0208 bool TrksInJetQAInJetFiller::IsCstNotRelevant(const uint32_t type) {
0209
0210   const bool isVoid   = (type == Jet::SRC::VOID);
0211   const bool isImport = (type == Jet::SRC::HEPMC_IMPORT);
0212   const bool isProbe  = (type == Jet::SRC::JET_PROBE);
0213   return (isVoid || isImport || isProbe);
0214
0215 }  // end 'IsCstNotRelevant(uint32_t)'
0216
0217
0218
0219 bool TrksInJetQAInJetFiller::IsTrkInList(const uint32_t id) {
0220
0221   bool isAdded = false;
0222   for (SvtxTrack* trkInJet : m_trksInJet) {
0223     if (id == trkInJet -> get_id() ) {
0224       isAdded = true;
0225       break;
0226     }
0227   }
0228   return isAdded;
0229
0230 }  // end 'IsTrkInList(uint32_t)'
0231
0232
0233
0234 double TrksInJetQAInJetFiller::GetTrackJetDist(SvtxTrack* track, Jet* jet) {
0235
0236   // get delta eta
0237   const double dEta = (track -> get_eta()) - (jet -> get_eta());
0238
0239   // get delta phi
0240   double dPhi = (track -> get_phi()) - (jet -> get_phi());
0241   if (dPhi < (-1. * TMath::Pi())) dPhi += TMath::TwoPi();
0242   if (dPhi > (1. * TMath::Pi()))  dPhi -= TMath::TwoPi();
0243
0244   // return distance
0245   return std::hypot(dEta, dPhi);
0246
0247 }  // end 'GetTrackJetDist(SvtxTrack*, Jet*)'
0248
0249
0250
0251 PFObject* TrksInJetQAInJetFiller::GetPFObject(const uint32_t id, PHCompositeNode* topNode) {
0252
0253   // pointer to pfo
0254   PFObject* pfoToFind = NULL;
0255
0256   // grab pf node if needed
0257   if (!m_flowStore) GetNode(Node::Flow, topNode);
0258
0259   // loop over pfos
0260   for (
0261       ParticleFlowElementContainer::ConstIterator itFlow = m_flowStore -> getParticleFlowElements().first;
0262       itFlow != m_flowStore -> getParticleFlowElements().second;
0263       ++itFlow
0264   ) {
0265
0266     // get pfo
0267     PFObject* pfo = itFlow -> second;
0268
0269     // if has provided id, set pointer and exit
0270     if (id == pfo -> get_id()) {
0271       pfoToFind = pfo;
0272       break;
0273     }
0274   }  // end pfo loop
0275   return pfoToFind;
0276
0277 }  // end 'GetPFObject(uint32_t, PHCompositeNode*)'
0278
0279
0280
0281 SvtxTrack* TrksInJetQAInJetFiller::GetTrkFromPFO(PFObject* pfo) {
0282
0283   // pointer to track
0284   SvtxTrack* track = NULL;
0285
0286   // if pfo has track, try to grab
0287   const auto type = pfo -> get_type();
0288   if (
0289     (type == ParticleFlowElement::PFLOWTYPE::MATCHED_CHARGED_HADRON) ||
0290     (type == ParticleFlowElement::PFLOWTYPE::UNMATCHED_CHARGED_HADRON)
0291   ) {
0292     track = pfo -> get_track();
0293   }
0294   return track;
0295
0296 }  // end 'GetTrkFromPFO(PFObject*)'
0297
0298 // end ------------------------------------------------------------------------