Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 10:05:11

0001 #ifndef HistoryBase_h
0002 #define HistoryBase_h
0003 
0004 #include <set>
0005 
0006 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
0007 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticleFwd.h"
0008 #include "SimDataFormats/TrackingAnalysis/interface/TrackingVertex.h"
0009 #include "SimDataFormats/TrackingAnalysis/interface/TrackingVertexContainer.h"
0010 
0011 //! Base class to all the history types.
0012 class HistoryBase {
0013 public:
0014   //! HepMC::GenParticle trail type.
0015   typedef std::vector<const HepMC::GenParticle *> GenParticleTrail;
0016 
0017   //! reco::GenParticle trail type.
0018   typedef std::vector<const reco::GenParticle *> RecoGenParticleTrail;
0019 
0020   //! reco::GenParticle trail helper type.
0021   typedef std::set<const reco::GenParticle *> RecoGenParticleTrailHelper;
0022 
0023   //! GenVertex trail type.
0024   typedef std::vector<const HepMC::GenVertex *> GenVertexTrail;
0025 
0026   //! GenVertex trail helper type.
0027   typedef std::set<const HepMC::GenVertex *> GenVertexTrailHelper;
0028 
0029   //! SimParticle trail type.
0030   typedef std::vector<TrackingParticleRef> SimParticleTrail;
0031 
0032   //! SimVertex trail type.
0033   typedef std::vector<TrackingVertexRef> SimVertexTrail;
0034 
0035   // Default constructor
0036   HistoryBase() {
0037     // Default depth
0038     depth_ = -1;
0039   }
0040 
0041   //! Set the depth of the history.
0042   /* Set TrackHistory to given depth. Positive values
0043      constrain the number of TrackingVertex visit in the history.
0044      Negatives values set the limit of the iteration over generated
0045      information i.e. (-1 -> status 1 or -2 -> status 2 particles).
0046 
0047      /param[in] depth the history
0048   */
0049   void depth(int d) { depth_ = d; }
0050 
0051   //! Return all the simulated vertices in the history.
0052   SimVertexTrail const &simVertexTrail() const { return simVertexTrail_; }
0053 
0054   //! Return all the simulated particle in the history.
0055   SimParticleTrail const &simParticleTrail() const { return simParticleTrail_; }
0056 
0057   //! Return all generated vertex in the history.
0058   GenVertexTrail const &genVertexTrail() const { return genVertexTrail_; }
0059 
0060   //! Return all generated particle (HepMC::GenParticle) in the history.
0061   GenParticleTrail const &genParticleTrail() const { return genParticleTrail_; }
0062 
0063   //! Return all reco::GenParticle in the history.
0064   RecoGenParticleTrail const &recoGenParticleTrail() const { return recoGenParticleTrail_; }
0065 
0066   //! Return the initial tracking particle from the history.
0067   const TrackingParticleRef &simParticle() const { return simParticleTrail_[0]; }
0068 
0069   //! Return the initial tracking vertex from the history.
0070   const TrackingVertexRef &simVertex() const { return simVertexTrail_[0]; }
0071 
0072   //! Returns a pointer to most primitive status 1 or 2 particle in the
0073   //! genParticleTrail_.
0074   const HepMC::GenParticle *genParticle() const {
0075     if (genParticleTrail_.empty())
0076       return nullptr;
0077     return genParticleTrail_[genParticleTrail_.size() - 1];
0078   }
0079 
0080   //! Returns a pointer to most primitive status 1 or 2 particle in the
0081   //! recoGenParticleTrail_.
0082   const reco::GenParticle *recoGenParticle() const {
0083     if (recoGenParticleTrail_.empty())
0084       return nullptr;
0085     return recoGenParticleTrail_[recoGenParticleTrail_.size() - 1];
0086   }
0087 
0088   //! Evaluate track history using a TrackingParticleRef.
0089   /* Return false when the history cannot be determined upto a given depth.
0090      If not depth is pass to the function no restriction are apply to it.
0091 
0092      /param[in] TrackingParticleRef of a simulated track
0093      /param[in] depth of the track history
0094      /param[out] boolean that is true when history can be determined
0095   */
0096   bool evaluate(TrackingParticleRef tpr) {
0097     resetTrails(tpr);
0098     return traceSimHistory(tpr, depth_);
0099   }
0100 
0101   //! Evaluate track history using a TrackingParticleRef.
0102   /* Return false when the history cannot be determined upto a given depth.
0103      If not depth is pass to the function no restriction are apply to it.
0104 
0105      /param[in] TrackingVertexRef of a simulated vertex
0106      /param[in] depth of the track history
0107      /param[out] boolean that is true when history can be determined
0108   */
0109   bool evaluate(TrackingVertexRef tvr) {
0110     resetTrails();
0111     return traceSimHistory(tvr, depth_);
0112   }
0113 
0114 protected:
0115   // History cointainers
0116   GenVertexTrail genVertexTrail_;
0117   GenParticleTrail genParticleTrail_;
0118   RecoGenParticleTrail recoGenParticleTrail_;
0119   SimVertexTrail simVertexTrail_;
0120   SimParticleTrail simParticleTrail_;
0121 
0122   // Helper function to speedup search
0123   GenVertexTrailHelper genVertexTrailHelper_;
0124   RecoGenParticleTrailHelper recoGenParticleTrailHelper_;
0125 
0126 private:
0127   int depth_;
0128 
0129   //! Trace all the simulated information for a given reference to a
0130   //! TrackingParticle.
0131   bool traceSimHistory(TrackingParticleRef const &, int);
0132 
0133   //! Trace all the simulated information for a given reference to a
0134   //! TrackingVertex.
0135   bool traceSimHistory(TrackingVertexRef const &, int);
0136 
0137   //! Trace all the simulated information for a given pointer to a
0138   //! HepMC::GenParticle.
0139   void traceGenHistory(HepMC::GenParticle const *);
0140 
0141   //! Trace all the simulated information for a given pointer to a
0142   //! reco::GenParticle.
0143   void traceRecoGenHistory(reco::GenParticle const *);
0144 
0145   //! Trace all the simulated information for a given pointer to a GenVertex.
0146   void traceGenHistory(HepMC::GenVertex const *);
0147 
0148   //! Reset trail functions.
0149   void resetTrails() {
0150     simParticleTrail_.clear();
0151     simVertexTrail_.clear();
0152     genVertexTrail_.clear();
0153     genParticleTrail_.clear();
0154     recoGenParticleTrail_.clear();
0155     recoGenParticleTrailHelper_.clear();
0156     genVertexTrailHelper_.clear();
0157   }
0158 
0159   void resetTrails(TrackingParticleRef tpr) {
0160     resetTrails();
0161     simParticleTrail_.push_back(tpr);
0162   }
0163 };
0164 
0165 #endif