CaloAnalysis/interface/SimCluster.h

0001 #ifndef SimDataFormats_SimCluster_h
0002 #define SimDataFormats_SimCluster_h
0003
0004 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0005 #include "DataFormats/Math/interface/LorentzVector.h"
0006 #include "DataFormats/Math/interface/Point3D.h"
0007 #include "DataFormats/Math/interface/Vector3D.h"
0008 #include "SimDataFormats/CaloHit/interface/PCaloHit.h"
0009 #include "SimDataFormats/Track/interface/SimTrack.h"
0010 #include <vector>
0011
0012 //
0013 // Forward declarations
0014 //
0015 class SimTrack;
0016 class EncodedEventId;
0017
0018 /** @brief Monte Carlo truth information used for tracking validation.
0019  *
0020  * Object with references to the original SimTrack and parent and daughter
0021  * TrackingVertices. Simulation with high (~100) pileup was taking too much
0022  * memory so the class was slimmed down and copies of the SimHits were removed.
0023  *
0024  * @author original author unknown, re-engineering and slimming by Subir Sarkar
0025  * (subir.sarkar@cern.ch), some tweaking and documentation by Mark Grimes
0026  * (mark.grimes@bristol.ac.uk).
0027  * @date original date unknown, re-engineering Jan-May 2013
0028  */
0029 class SimCluster {
0030   friend std::ostream &operator<<(std::ostream &s, SimCluster const &tp);
0031
0032 public:
0033   typedef int Charge;                                       ///< electric charge type
0034   typedef math::XYZTLorentzVectorD LorentzVector;           ///< Lorentz vector
0035   typedef math::PtEtaPhiMLorentzVector PolarLorentzVector;  ///< Lorentz vector
0036   typedef math::XYZPointD Point;                            ///< point in the space
0037   typedef math::XYZVectorD Vector;                          ///< point in the space
0038
0039   /// reference to reco::GenParticle
0040   typedef reco::GenParticleRefVector::iterator genp_iterator;
0041   typedef std::vector<SimTrack>::const_iterator g4t_iterator;
0042
0043   SimCluster();
0044
0045   SimCluster(const SimTrack &simtrk);
0046   SimCluster(EncodedEventId eventID, uint32_t particleID);  // for PU
0047
0048   // destructor
0049   ~SimCluster();
0050
0051   /** @brief PDG ID.
0052    *
0053    * Returns the PDG ID of the first associated gen particle. If there are no
0054    * gen particles associated then it returns type() from the first SimTrack. */
0055   int pdgId() const {
0056     if (genParticles_.empty())
0057       return g4Tracks_[0].type();
0058     else
0059       return (*genParticles_.begin())->pdgId();
0060   }
0061
0062   /** @brief Signal source, crossing number.
0063    *
0064    * Note this is taken from the first SimTrack only, but there shouldn't be any
0065    * SimTracks from different crossings in the SimCluster. */
0066   EncodedEventId eventId() const { return event_; }
0067
0068   uint64_t particleId() const { return particleId_; }
0069
0070   // Setters for G4 and reco::GenParticle
0071   void addGenParticle(const reco::GenParticleRef &ref) { genParticles_.push_back(ref); }
0072   void addG4Track(const SimTrack &t) { g4Tracks_.push_back(t); }
0073   /// iterators
0074   genp_iterator genParticle_begin() const { return genParticles_.begin(); }
0075   genp_iterator genParticle_end() const { return genParticles_.end(); }
0076   g4t_iterator g4Track_begin() const { return g4Tracks_.begin(); }
0077   g4t_iterator g4Track_end() const { return g4Tracks_.end(); }
0078
0079   // Getters for Embd and Sim Tracks
0080   const reco::GenParticleRefVector &genParticles() const { return genParticles_; }
0081   // Only for clusters from the signal vertex
0082   const std::vector<SimTrack> &g4Tracks() const { return g4Tracks_; }
0083
0084   /// @brief Electric charge. Note this is taken from the first SimTrack only.
0085   float charge() const { return g4Tracks_[0].charge(); }
0086   /// Gives charge in unit of quark charge (should be 3 times "charge()")
0087   int threeCharge() const { return lrintf(3.f * charge()); }
0088
0089   /// @brief Four-momentum Lorentz vector. Note this is taken from the first
0090   /// SimTrack only.
0091   const math::XYZTLorentzVectorF &p4() const { return theMomentum_; }
0092
0093   /// @brief spatial momentum vector
0094   math::XYZVectorF momentum() const { return p4().Vect(); }
0095
0096   /// @brief Vector to boost to the particle centre of mass frame.
0097   math::XYZVectorF boostToCM() const { return p4().BoostToCM(); }
0098
0099   /// @brief Magnitude of momentum vector. Note this is taken from the first
0100   /// SimTrack only.
0101   float p() const { return p4().P(); }
0102
0103   /// @brief Energy. Note this is taken from the first SimTrack only.
0104   float energy() const { return p4().E(); }
0105
0106   /// @brief Transverse energy. Note this is taken from the first SimTrack only.
0107   float et() const { return p4().Et(); }
0108
0109   /// @brief Mass. Note this is taken from the first SimTrack only.
0110   float mass() const { return p4().M(); }
0111
0112   /// @brief Mass squared. Note this is taken from the first SimTrack only.
0113   float massSqr() const { return pow(mass(), 2); }
0114
0115   /// @brief Transverse mass. Note this is taken from the first SimTrack only.
0116   float mt() const { return p4().Mt(); }
0117
0118   /// @brief Transverse mass squared. Note this is taken from the first SimTrack
0119   /// only.
0120   float mtSqr() const { return p4().Mt2(); }
0121
0122   /// @brief x coordinate of momentum vector. Note this is taken from the first
0123   /// SimTrack only.
0124   float px() const { return p4().Px(); }
0125
0126   /// @brief y coordinate of momentum vector. Note this is taken from the first
0127   /// SimTrack only.
0128   float py() const { return p4().Py(); }
0129
0130   /// @brief z coordinate of momentum vector. Note this is taken from the first
0131   /// SimTrack only.
0132   float pz() const { return p4().Pz(); }
0133
0134   /// @brief Transverse momentum. Note this is taken from the first SimTrack
0135   /// only.
0136   float pt() const { return p4().Pt(); }
0137
0138   /// @brief Momentum azimuthal angle. Note this is taken from the first
0139   /// SimTrack only.
0140   float phi() const { return p4().Phi(); }
0141
0142   /// @brief Momentum polar angle. Note this is taken from the first SimTrack
0143   /// only.
0144   float theta() const { return p4().Theta(); }
0145
0146   /// @brief Momentum pseudorapidity. Note this is taken from the simtrack
0147   /// before the calorimeter
0148   float eta() const { return p4().Eta(); }
0149
0150   /// @brief Rapidity. Note this is taken from the simtrack before the
0151   /// calorimeter
0152   float rapidity() const { return p4().Rapidity(); }
0153
0154   /// @brief Same as rapidity().
0155   float y() const { return rapidity(); }
0156
0157   /** @brief Status word.
0158    *
0159    * Returns status() from the first gen particle, or -99 if there are no gen
0160    * particles attached. */
0161   int status() const { return genParticles_.empty() ? -99 : (*genParticles_[0]).status(); }
0162
0163   static const unsigned int longLivedTag;  ///< long lived flag
0164
0165   /// is long lived?
0166   bool longLived() const { return status() & longLivedTag; }
0167
0168   /** @brief Gives the total number of SimHits, in the cluster */
0169   int numberOfSimHits() const { return nsimhits_; }
0170
0171   /** @brief Gives the total number of SimHits, in the cluster */
0172   int numberOfRecHits() const { return hits_.size(); }
0173
0174   /** @brief add rechit with fraction */
0175   void addRecHitAndFraction(uint32_t hit, float fraction) {
0176     hits_.emplace_back(hit);
0177     fractions_.emplace_back(fraction);
0178   }
0179
0180   /** @brief add rechit energy */
0181   void addHitEnergy(float energy) { energies_.emplace_back(energy); }
0182
0183   /** @brief Returns list of rechit IDs and fractions for this SimCluster */
0184   std::vector<std::pair<uint32_t, float>> hits_and_fractions() const {
0185     std::vector<std::pair<uint32_t, float>> result;
0186     for (size_t i = 0; i < hits_.size(); ++i) {
0187       result.emplace_back(hits_[i], fractions_[i]);
0188     }
0189     return result;
0190   }
0191
0192   /** @brief Returns list of rechit IDs and energies for this SimCluster */
0193   std::vector<std::pair<uint32_t, float>> hits_and_energies() const {
0194     assert(hits_.size() == energies_.size());
0195     std::vector<std::pair<uint32_t, float>> result;
0196     result.reserve(hits_.size());
0197     for (size_t i = 0; i < hits_.size(); ++i) {
0198       result.emplace_back(hits_[i], energies_[i]);
0199     }
0200     return result;
0201   }
0202
0203   /** @brief clear the hits and fractions list */
0204   void clearHitsAndFractions() {
0205     std::vector<uint32_t>().swap(hits_);
0206     std::vector<float>().swap(fractions_);
0207   }
0208
0209   /** @brief clear the energies list */
0210   void clearHitsEnergy() { std::vector<float>().swap(energies_); }
0211
0212   /** @brief returns the accumulated sim energy in the cluster */
0213   float simEnergy() const { return simhit_energy_; }
0214
0215   /** @brief add simhit's energy to cluster */
0216   void addSimHit(const PCaloHit &hit) {
0217     simhit_energy_ += hit.energy();
0218     ++nsimhits_;
0219   }
0220
0221 private:
0222   uint64_t nsimhits_{0};
0223   EncodedEventId event_;
0224
0225   uint32_t particleId_{0};
0226   float simhit_energy_{0.f};
0227   std::vector<uint32_t> hits_;
0228   std::vector<float> fractions_;
0229   std::vector<float> energies_;
0230
0231   math::XYZTLorentzVectorF theMomentum_;
0232
0233   /// references to G4 and reco::GenParticle tracks
0234   std::vector<SimTrack> g4Tracks_;
0235   reco::GenParticleRefVector genParticles_;
0236 };
0237
0238 #endif  // SimDataFormats_SimCluster_H