File indexing completed on 2024-12-06 02:45:14
0001 #ifndef Candidate_ParticleState_h
0002 #define Candidate_ParticleState_h
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013 #include <Rtypes.h>
0014
0015 #include "DataFormats/Math/interface/Point3D.h"
0016 #include "DataFormats/Math/interface/Vector3D.h"
0017 #include "DataFormats/Math/interface/PtEtaPhiMass.h"
0018 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0019 #include "DataFormats/Math/interface/LorentzVector.h"
0020
0021 namespace reco {
0022
0023 class ParticleState {
0024 public:
0025
0026 typedef int Charge;
0027
0028 typedef math::XYZTLorentzVector LorentzVector;
0029
0030 typedef math::PtEtaPhiMLorentzVector PolarLorentzVector;
0031
0032 typedef math::XYZPoint Point;
0033
0034 typedef math::XYZVector Vector;
0035
0036 ParticleState() : vertex_(0, 0, 0), qx3_(0), pdgId_(0), status_(0) {}
0037
0038
0039 ParticleState(Charge q,
0040 const PtEtaPhiMass& p4,
0041 const Point& vertex = Point(0, 0, 0),
0042 int pdgId = 0,
0043 int status = 0,
0044 bool integerCharge = true)
0045 : vertex_(vertex),
0046 p4Polar_(p4.pt(), p4.eta(), p4.phi(), p4.mass()),
0047 p4Cartesian_(p4Polar_),
0048 qx3_(integerCharge ? q * 3 : q),
0049 pdgId_(pdgId),
0050 status_(status) {}
0051
0052
0053 ParticleState(Charge q,
0054 const LorentzVector& p4,
0055 const Point& vertex = Point(0, 0, 0),
0056 int pdgId = 0,
0057 int status = 0,
0058 bool integerCharge = true)
0059 : vertex_(vertex),
0060 p4Polar_(p4),
0061 p4Cartesian_(p4),
0062 qx3_(integerCharge ? q * 3 : q),
0063 pdgId_(pdgId),
0064 status_(status) {}
0065
0066
0067 ParticleState(Charge q,
0068 const PolarLorentzVector& p4,
0069 const Point& vertex = Point(0, 0, 0),
0070 int pdgId = 0,
0071 int status = 0,
0072 bool integerCharge = true)
0073 : vertex_(vertex),
0074 p4Polar_(p4),
0075 p4Cartesian_(p4),
0076 qx3_(integerCharge ? q * 3 : q),
0077 pdgId_(pdgId),
0078 status_(status) {}
0079
0080 ParticleState(Charge q,
0081 const GlobalVector& p3,
0082 float iEnergy,
0083 float imass,
0084 const Point& vertex = Point(0, 0, 0),
0085 int pdgId = 0,
0086 int status = 0,
0087 bool integerCharge = true)
0088 : vertex_(vertex),
0089 p4Polar_(p3.perp(), p3.eta(), p3.phi(), imass),
0090 p4Cartesian_(p3.x(), p3.y(), p3.z(), iEnergy),
0091 qx3_(integerCharge ? q * 3 : q),
0092 pdgId_(pdgId),
0093 status_(status) {}
0094
0095
0096 inline void setCartesian() { p4Cartesian_ = p4Polar_; }
0097
0098
0099 int charge() const { return qx3_ / 3; }
0100
0101 void setCharge(Charge q) { qx3_ = q * 3; }
0102
0103 int threeCharge() const { return qx3_; }
0104
0105 void setThreeCharge(Charge qx3) { qx3_ = qx3; }
0106
0107 const LorentzVector& p4() const { return p4Cartesian_; }
0108
0109 const PolarLorentzVector& polarP4() const { return p4Polar_; }
0110
0111 Vector momentum() const { return p4Cartesian_.Vect(); }
0112
0113
0114 Vector boostToCM() const { return p4Cartesian_.BoostToCM(); }
0115
0116 double p() const { return p4Cartesian_.P(); }
0117
0118 double energy() const { return p4Cartesian_.E(); }
0119
0120 double et() const { return (pt() <= 0) ? 0 : p4Cartesian_.Et(); }
0121
0122 double et2() const { return (pt() <= 0) ? 0 : p4Cartesian_.Et2(); }
0123
0124 double mass() const { return p4Polar_.mass(); }
0125
0126 double massSqr() const { return mass() * mass(); }
0127
0128 double mt() const { return p4Polar_.Mt(); }
0129
0130 double mtSqr() const { return p4Polar_.Mt2(); }
0131
0132 double px() const { return p4Cartesian_.Px(); }
0133
0134 double py() const { return p4Cartesian_.Py(); }
0135
0136 double pz() const { return p4Cartesian_.Pz(); }
0137
0138 double pt() const { return p4Polar_.pt(); }
0139
0140 double phi() const { return p4Polar_.phi(); }
0141
0142 double theta() const { return p4Cartesian_.Theta(); }
0143
0144 double eta() const { return p4Polar_.eta(); }
0145
0146 double rapidity() const { return p4Polar_.Rapidity(); }
0147
0148 double y() const { return rapidity(); }
0149
0150 void setP4(const LorentzVector& p4) {
0151 p4Cartesian_ = p4;
0152 p4Polar_ = p4;
0153 }
0154
0155 void setP4(const PolarLorentzVector& p4) {
0156 p4Polar_ = p4;
0157 p4Cartesian_ = p4;
0158 }
0159
0160 void setMass(double m) {
0161 p4Polar_.SetM(m);
0162 setCartesian();
0163 }
0164 void setPz(double pz) {
0165 p4Cartesian_.SetPz(pz);
0166 p4Polar_ = p4Cartesian_;
0167 }
0168
0169 const Point& vertex() const { return vertex_; }
0170
0171 double vx() const { return vertex_.X(); }
0172
0173 double vy() const { return vertex_.Y(); }
0174
0175 double vz() const { return vertex_.Z(); }
0176
0177 void setVertex(const Point& vertex) { vertex_ = vertex; }
0178
0179 int pdgId() const { return pdgId_; }
0180
0181 void setPdgId(int pdgId) { pdgId_ = pdgId; }
0182
0183 int status() const { return status_; }
0184
0185 void setStatus(int status) { status_ = status; }
0186
0187 void setLongLived() { status_ |= longLivedTag; }
0188
0189 bool longLived() const { return status_ & longLivedTag; }
0190
0191 void setMassConstraint() { status_ |= massConstraintTag; }
0192
0193 bool massConstraint() const { return status_ & massConstraintTag; }
0194
0195 private:
0196 static const unsigned int longLivedTag = 65536;
0197 static const unsigned int massConstraintTag = 131072;
0198
0199 private:
0200
0201 Point vertex_;
0202
0203
0204 PolarLorentzVector p4Polar_;
0205
0206 LorentzVector p4Cartesian_;
0207
0208
0209 Charge qx3_;
0210
0211
0212 int pdgId_;
0213
0214 int status_;
0215 };
0216
0217 }
0218
0219 #endif