File indexing completed on 2024-12-02 23:44:59
0001 #ifndef SimDataFormats_CaloParticle_h
0002 #define SimDataFormats_CaloParticle_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/CaloAnalysis/interface/SimClusterFwd.h"
0009 #include "SimDataFormats/CaloHit/interface/PCaloHit.h"
0010 #include "SimDataFormats/Track/interface/SimTrack.h"
0011 #include <vector>
0012
0013 class SimTrack;
0014 class EncodedEventId;
0015
0016 class CaloParticle {
0017 friend std::ostream &operator<<(std::ostream &s, CaloParticle const &tp);
0018
0019 public:
0020 typedef int Charge;
0021 typedef math::XYZTLorentzVectorD LorentzVector;
0022 typedef math::PtEtaPhiMLorentzVector PolarLorentzVector;
0023 typedef math::XYZPointD Point;
0024 typedef math::XYZVectorD Vector;
0025
0026
0027 typedef reco::GenParticleRefVector::iterator genp_iterator;
0028 typedef std::vector<SimTrack>::const_iterator g4t_iterator;
0029 typedef SimClusterRefVector::iterator sc_iterator;
0030
0031 CaloParticle();
0032
0033 CaloParticle(const SimTrack &simtrk);
0034 CaloParticle(EncodedEventId eventID, uint32_t particleID);
0035
0036
0037 ~CaloParticle();
0038
0039
0040
0041
0042
0043 int pdgId() const {
0044 if (genParticles_.empty())
0045 return g4Tracks_[0].type();
0046 else
0047 return (*genParticles_.begin())->pdgId();
0048 }
0049
0050
0051
0052
0053
0054 EncodedEventId eventId() const { return event_; }
0055
0056 uint64_t particleId() const { return particleId_; }
0057
0058
0059 void addGenParticle(const reco::GenParticleRef &ref) { genParticles_.push_back(ref); }
0060 void addSimCluster(const SimClusterRef &ref) { simClusters_.push_back(ref); }
0061 void addG4Track(const SimTrack &t) { g4Tracks_.push_back(t); }
0062
0063 genp_iterator genParticle_begin() const { return genParticles_.begin(); }
0064 genp_iterator genParticle_end() const { return genParticles_.end(); }
0065 g4t_iterator g4Track_begin() const { return g4Tracks_.begin(); }
0066 g4t_iterator g4Track_end() const { return g4Tracks_.end(); }
0067 sc_iterator simCluster_begin() const { return simClusters_.begin(); }
0068 sc_iterator simCluster_end() const { return simClusters_.end(); }
0069
0070
0071 const reco::GenParticleRefVector &genParticles() const { return genParticles_; }
0072 const SimClusterRefVector &simClusters() const { return simClusters_; }
0073
0074 const std::vector<SimTrack> &g4Tracks() const { return g4Tracks_; }
0075
0076 void clearSimClusters() { simClusters_.clear(); }
0077
0078
0079 float charge() const { return g4Tracks_[0].charge(); }
0080
0081 int threeCharge() const { return lrintf(3.f * charge()); }
0082
0083
0084
0085 const math::XYZTLorentzVectorF &p4() const { return theMomentum_; }
0086
0087
0088 math::XYZVectorF momentum() const { return p4().Vect(); }
0089
0090
0091 math::XYZVectorF boostToCM() const { return p4().BoostToCM(); }
0092
0093
0094
0095 float p() const { return p4().P(); }
0096
0097
0098 float energy() const { return p4().E(); }
0099
0100
0101 float et() const { return p4().Et(); }
0102
0103
0104 float mass() const { return p4().M(); }
0105
0106
0107 float massSqr() const { return pow(mass(), 2); }
0108
0109
0110 float mt() const { return p4().Mt(); }
0111
0112
0113
0114 float mtSqr() const { return p4().Mt2(); }
0115
0116
0117
0118 float px() const { return p4().Px(); }
0119
0120
0121
0122 float py() const { return p4().Py(); }
0123
0124
0125
0126 float pz() const { return p4().Pz(); }
0127
0128
0129
0130 float pt() const { return p4().Pt(); }
0131
0132
0133
0134 float phi() const { return p4().Phi(); }
0135
0136
0137
0138 float theta() const { return p4().Theta(); }
0139
0140
0141
0142 float eta() const { return p4().Eta(); }
0143
0144
0145
0146 float rapidity() const { return p4().Rapidity(); }
0147
0148
0149 float y() const { return rapidity(); }
0150
0151
0152
0153
0154
0155 int status() const { return genParticles_.empty() ? -99 : (*genParticles_[0]).status(); }
0156
0157 static const unsigned int longLivedTag;
0158
0159
0160 bool longLived() const { return status() & longLivedTag; }
0161
0162
0163 int numberOfSimHits() const { return nsimhits_; }
0164
0165
0166 int numberOfRecHits() const { return hits_.size(); }
0167
0168
0169 float simTime() const { return time_; }
0170
0171
0172 void addRecHitAndFraction(uint32_t hit, float fraction) {
0173 hits_.emplace_back(hit);
0174 fractions_.emplace_back(fraction);
0175 }
0176
0177
0178 std::vector<std::pair<uint32_t, float>> hits_and_fractions() const {
0179 std::vector<std::pair<uint32_t, float>> result;
0180 for (size_t i = 0; i < hits_.size(); ++i) {
0181 result.emplace_back(hits_[i], fractions_[i]);
0182 }
0183 return result;
0184 }
0185
0186
0187 float simEnergy() const { return simhit_energy_; }
0188
0189
0190 void addSimHit(const PCaloHit &hit) {
0191 simhit_energy_ += hit.energy();
0192 ++nsimhits_;
0193 }
0194
0195
0196 void setSimTime(const float time) { time_ = time; }
0197
0198 protected:
0199 uint64_t nsimhits_{0};
0200 EncodedEventId event_;
0201
0202 uint32_t particleId_{0};
0203 float simhit_energy_{0.f};
0204 float time_{std::numeric_limits<float>::lowest()};
0205 std::vector<uint32_t> hits_;
0206 std::vector<float> fractions_;
0207
0208 math::XYZTLorentzVectorF theMomentum_;
0209
0210
0211 std::vector<SimTrack> g4Tracks_;
0212 reco::GenParticleRefVector genParticles_;
0213
0214 SimClusterRefVector simClusters_;
0215 };
0216
0217 #endif