File indexing completed on 2021-12-06 04:01:05
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
0014
0015 class SimTrack;
0016 class EncodedEventId;
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029 class SimCluster {
0030 friend std::ostream &operator<<(std::ostream &s, SimCluster const &tp);
0031
0032 public:
0033 typedef int Charge;
0034 typedef math::XYZTLorentzVectorD LorentzVector;
0035 typedef math::PtEtaPhiMLorentzVector PolarLorentzVector;
0036 typedef math::XYZPointD Point;
0037 typedef math::XYZVectorD Vector;
0038
0039
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);
0047
0048
0049 ~SimCluster();
0050
0051
0052
0053
0054
0055 int pdgId() const {
0056 if (genParticles_.empty())
0057 return g4Tracks_[0].type();
0058 else
0059 return (*genParticles_.begin())->pdgId();
0060 }
0061
0062
0063
0064
0065
0066 EncodedEventId eventId() const { return event_; }
0067
0068 uint64_t particleId() const { return particleId_; }
0069
0070
0071 void addGenParticle(const reco::GenParticleRef &ref) { genParticles_.push_back(ref); }
0072 void addG4Track(const SimTrack &t) { g4Tracks_.push_back(t); }
0073
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
0080 const reco::GenParticleRefVector &genParticles() const { return genParticles_; }
0081
0082 const std::vector<SimTrack> &g4Tracks() const { return g4Tracks_; }
0083
0084
0085 float charge() const { return g4Tracks_[0].charge(); }
0086
0087 int threeCharge() const { return lrintf(3.f * charge()); }
0088
0089
0090
0091 const math::XYZTLorentzVectorF &p4() const { return theMomentum_; }
0092
0093
0094 math::XYZVectorF momentum() const { return p4().Vect(); }
0095
0096
0097 math::XYZVectorF boostToCM() const { return p4().BoostToCM(); }
0098
0099
0100
0101 float p() const { return p4().P(); }
0102
0103
0104 float energy() const { return p4().E(); }
0105
0106
0107 float et() const { return p4().Et(); }
0108
0109
0110 float mass() const { return p4().M(); }
0111
0112
0113 float massSqr() const { return pow(mass(), 2); }
0114
0115
0116 float mt() const { return p4().Mt(); }
0117
0118
0119
0120 float mtSqr() const { return p4().Mt2(); }
0121
0122
0123
0124 float px() const { return p4().Px(); }
0125
0126
0127
0128 float py() const { return p4().Py(); }
0129
0130
0131
0132 float pz() const { return p4().Pz(); }
0133
0134
0135
0136 float pt() const { return p4().Pt(); }
0137
0138
0139
0140 float phi() const { return p4().Phi(); }
0141
0142
0143
0144 float theta() const { return p4().Theta(); }
0145
0146
0147
0148 float eta() const { return p4().Eta(); }
0149
0150
0151
0152 float rapidity() const { return p4().Rapidity(); }
0153
0154
0155 float y() const { return rapidity(); }
0156
0157
0158
0159
0160
0161 int status() const { return genParticles_.empty() ? -99 : (*genParticles_[0]).status(); }
0162
0163 static const unsigned int longLivedTag;
0164
0165
0166 bool longLived() const { return status() & longLivedTag; }
0167
0168
0169 int numberOfSimHits() const { return nsimhits_; }
0170
0171
0172 int numberOfRecHits() const { return hits_.size(); }
0173
0174
0175 void addRecHitAndFraction(uint32_t hit, float fraction) {
0176 hits_.emplace_back(hit);
0177 fractions_.emplace_back(fraction);
0178 }
0179
0180
0181 void addHitEnergy(float energy) { energies_.emplace_back(energy); }
0182
0183
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
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
0204 void clearHitsAndFractions() {
0205 std::vector<uint32_t>().swap(hits_);
0206 std::vector<float>().swap(fractions_);
0207 }
0208
0209
0210 void clearHitsEnergy() { std::vector<float>().swap(energies_); }
0211
0212
0213 float simEnergy() const { return simhit_energy_; }
0214
0215
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
0234 std::vector<SimTrack> g4Tracks_;
0235 reco::GenParticleRefVector genParticles_;
0236 };
0237
0238 #endif