File indexing completed on 2025-06-03 00:12:19
0001 #ifndef L1Trigger_TrackFindingTMTT_L1track3D_h
0002 #define L1Trigger_TrackFindingTMTT_L1track3D_h
0003
0004 #include "DataFormats/Math/interface/deltaPhi.h"
0005 #include "L1Trigger/TrackFindingTMTT/interface/L1trackBase.h"
0006 #include "L1Trigger/TrackFindingTMTT/interface/Settings.h"
0007 #include "L1Trigger/TrackFindingTMTT/interface/Utility.h"
0008 #include "L1Trigger/TrackFindingTMTT/interface/TP.h"
0009 #include "L1Trigger/TrackFindingTMTT/interface/Stub.h"
0010 #include "L1Trigger/TrackFindingTMTT/interface/Sector.h"
0011 #include "L1Trigger/TrackFindingTMTT/interface/HTrphi.h"
0012
0013 #include <vector>
0014 #include <string>
0015 #include <unordered_set>
0016 #include <utility>
0017
0018
0019
0020
0021
0022 namespace tmtt {
0023
0024 class L1track3D : public L1trackBase {
0025 public:
0026
0027 enum TrackletSeedType { L1L2, L2L3, L3L4, L5L6, D1D2, D3D4, L1D1, L2D1, L3L4L2, L5L6L4, L2L3D1, D1D2L2, NONE };
0028
0029 public:
0030 L1track3D(const Settings* settings,
0031 const std::vector<Stub*>& stubs,
0032 std::pair<unsigned int, unsigned int> cellLocationHT,
0033 std::pair<float, float> helixRphi,
0034 std::pair<float, float> helixRz,
0035 float helixD0,
0036 unsigned int iPhiSec,
0037 unsigned int iEtaReg,
0038 unsigned int optoLinkID,
0039 bool mergedHTcell)
0040 : L1trackBase(),
0041 settings_(settings),
0042 stubs_(stubs),
0043 stubsConst_(stubs_.begin(), stubs_.end()),
0044 cellLocationHT_(cellLocationHT),
0045 helixRphi_(helixRphi),
0046 helixRz_(helixRz),
0047 helixD0_(helixD0),
0048 iPhiSec_(iPhiSec),
0049 iEtaReg_(iEtaReg),
0050 optoLinkID_(optoLinkID),
0051 mergedHTcell_(mergedHTcell),
0052 seedLayerType_(TrackletSeedType::NONE),
0053 seedPS_(999) {
0054 nLayers_ = Utility::countLayers(settings, stubs_);
0055 matchedTP_ = Utility::matchingTP(settings,
0056 stubs_,
0057 nMatchedLayers_,
0058 matchedStubs_);
0059 }
0060
0061
0062
0063 L1track3D(const Settings* settings,
0064 const std::vector<Stub*>& stubs,
0065 std::pair<unsigned int, unsigned int> cellLocationHT,
0066 std::pair<float, float> helixRphi,
0067 std::pair<float, float> helixRz,
0068 unsigned int iPhiSec,
0069 unsigned int iEtaReg,
0070 unsigned int optoLinkID,
0071 bool mergedHTcell)
0072 : L1track3D(
0073 settings, stubs, cellLocationHT, helixRphi, helixRz, 0.0, iPhiSec, iEtaReg, optoLinkID, mergedHTcell) {}
0074
0075
0076 L1track3D(Settings* settings,
0077 std::vector<Stub*> stubs,
0078 double qOverPt,
0079 double phi0,
0080 double z0,
0081 double tanLambda,
0082 double helixD0,
0083 int iPhiSec,
0084 int iEtaReg)
0085 : settings_(settings),
0086 stubs_(std::move(stubs)),
0087 nLayers_(0),
0088 cellLocationHT_(0, 0),
0089 helixRphi_(qOverPt, phi0),
0090 helixRz_(z0, tanLambda),
0091 helixD0_(helixD0),
0092 iPhiSec_(iPhiSec),
0093 iEtaReg_(iEtaReg),
0094 optoLinkID_(0),
0095 mergedHTcell_(false),
0096 seedLayerType_(TrackletSeedType()),
0097 seedPS_(0),
0098 matchedTP_(nullptr),
0099 nMatchedLayers_(0) {}
0100
0101 ~L1track3D() override = default;
0102
0103
0104
0105
0106 void setSeedLayerType(unsigned int seedLayerType) { seedLayerType_ = static_cast<TrackletSeedType>(seedLayerType); }
0107 TrackletSeedType seedLayerType() const { return seedLayerType_; }
0108
0109
0110 void setSeedPS(unsigned int seedPS) { seedPS_ = seedPS; }
0111 unsigned int seedPS() const { return seedPS_; }
0112
0113
0114 void setBestStubs(std::unordered_set<const Stub*> bestStubs) { bestStubs_ = bestStubs; }
0115 std::unordered_set<const Stub*> bestStubs() const { return bestStubs_; }
0116
0117
0118
0119
0120 const std::vector<const Stub*>& stubsConst() const override { return stubsConst_; }
0121 const std::vector<Stub*>& stubs() const override { return stubs_; }
0122
0123 unsigned int numStubs() const override { return stubs_.size(); }
0124
0125 unsigned int numLayers() const override { return nLayers_; }
0126
0127 std::pair<unsigned int, unsigned int> cellLocationHT() const override { return cellLocationHT_; }
0128
0129 std::pair<float, float> helixRphi() const { return helixRphi_; }
0130
0131 std::pair<float, float> helixRz() const { return helixRz_; }
0132
0133
0134
0135 std::vector<float> chiPhi() {
0136 std::vector<float> result;
0137 for (const Stub* s : stubs_) {
0138 float chi_phi = reco::deltaPhi(s->phi(), this->phi0() - s->r() * this->qOverPt() * settings_->invPtToDphi());
0139 result.push_back(chi_phi);
0140 }
0141 return result;
0142 }
0143
0144 std::vector<int> chiPhiDigi() {
0145 std::vector<int> result;
0146 const float phiMult = pow(2, settings_->phiSBits()) / settings_->phiSRange();
0147 for (const float& chi_phi : this->chiPhi()) {
0148 int iDigi_chi_phi = floor(chi_phi * phiMult);
0149 result.push_back(iDigi_chi_phi);
0150 }
0151 return result;
0152 }
0153
0154 std::vector<float> chiZ() {
0155 std::vector<float> result;
0156 for (const Stub* s : stubs_) {
0157 float chi_z = s->z() - (this->z0() + s->r() * this->tanLambda());
0158 result.push_back(chi_z);
0159 }
0160 return result;
0161 }
0162
0163 std::vector<int> chiZDigi() {
0164 std::vector<int> result;
0165 const float zMult = pow(2, settings_->zBits()) / settings_->zRange();
0166 for (const float& chi_z : this->chiZ()) {
0167 int iDigi_chi_z = floor(chi_z * zMult);
0168 result.push_back(iDigi_chi_z);
0169 }
0170 return result;
0171 }
0172
0173
0174
0175 float qOverPt() const override { return helixRphi_.first; }
0176 float charge() const { return (this->qOverPt() > 0 ? 1 : -1); }
0177 float invPt() const { return std::abs(this->qOverPt()); }
0178
0179 float pt() const {
0180 constexpr float small = 1.0e-6;
0181 return 1. / (small + this->invPt());
0182 }
0183 float d0() const { return helixD0_; }
0184 float phi0() const override { return helixRphi_.second; }
0185 float z0() const { return helixRz_.first; }
0186 float tanLambda() const { return helixRz_.second; }
0187 float theta() const { return atan2(1., this->tanLambda()); }
0188 float eta() const { return -log(tan(0.5 * this->theta())); }
0189
0190
0191 float phiAtChosenR() const {
0192 return reco::deltaPhi(this->phi0() - (settings_->invPtToDphi() * settings_->chosenRofPhi()) * this->qOverPt(),
0193 0.);
0194 }
0195 float zAtChosenR() const {
0196 return (this->z0() + (settings_->chosenRofZ()) * this->tanLambda());
0197 }
0198
0199
0200 unsigned int iPhiSec() const override { return iPhiSec_; }
0201 unsigned int iEtaReg() const override { return iEtaReg_; }
0202
0203
0204 unsigned int optoLinkID() const override { return optoLinkID_; }
0205
0206
0207 bool mergedHTcell() const { return mergedHTcell_; }
0208
0209
0210
0211
0212 const TP* matchedTP() const override { return matchedTP_; }
0213
0214 const std::vector<const Stub*>& matchedStubs() const override { return matchedStubs_; }
0215
0216 unsigned int numMatchedStubs() const override { return matchedStubs_.size(); }
0217
0218 unsigned int numMatchedLayers() const override { return nMatchedLayers_; }
0219
0220 float purity() const { return numMatchedStubs() / float(numStubs()); }
0221
0222
0223
0224
0225
0226
0227
0228 bool cheat() {
0229 bool keep = false;
0230
0231 std::vector<Stub*> stubsSel;
0232 if (matchedTP_ != nullptr) {
0233 for (Stub* s : stubs_) {
0234 const TP* tp = s->assocTP();
0235 if (tp != nullptr) {
0236 if (matchedTP_->index() == tp->index()) {
0237 stubsSel.push_back(s);
0238 }
0239 }
0240 }
0241 }
0242 stubs_ = stubsSel;
0243 stubsConst_ = std::vector<const Stub*>(stubs_.begin(), stubs_.end());
0244
0245 nLayers_ = Utility::countLayers(settings_, stubs_);
0246 matchedTP_ = Utility::matchingTP(settings_,
0247 stubs_,
0248 nMatchedLayers_,
0249 matchedStubs_);
0250
0251 bool genuine = (matchedTP_ != nullptr);
0252
0253 if (genuine && matchedTP_->useForAlgEff()) {
0254 Sector secTmp(settings_, iPhiSec_, iEtaReg_);
0255 HTrphi htRphiTmp(settings_, iPhiSec_, iEtaReg_, secTmp.etaMin(), secTmp.etaMax(), secTmp.phiCentre());
0256 std::pair<unsigned int, unsigned int> trueCell = htRphiTmp.trueCell(matchedTP_);
0257
0258 std::pair<unsigned int, unsigned int> htCell = this->cellLocationHT();
0259 bool consistent = (htCell == trueCell);
0260 if (mergedHTcell_) {
0261
0262 std::pair<unsigned int, unsigned int> htCell10(htCell.first + 1, htCell.second);
0263 std::pair<unsigned int, unsigned int> htCell01(htCell.first, htCell.second + 1);
0264 std::pair<unsigned int, unsigned int> htCell11(htCell.first + 1, htCell.second + 1);
0265 if (htCell10 == trueCell)
0266 consistent = true;
0267 if (htCell01 == trueCell)
0268 consistent = true;
0269 if (htCell11 == trueCell)
0270 consistent = true;
0271 }
0272 if (consistent)
0273 keep = true;
0274 }
0275
0276 return keep;
0277 }
0278
0279 private:
0280
0281 const Settings* settings_;
0282
0283
0284 std::vector<Stub*> stubs_;
0285 std::vector<const Stub*> stubsConst_;
0286 std::unordered_set<const Stub*> bestStubs_;
0287 unsigned int nLayers_;
0288 std::pair<unsigned int, unsigned int> cellLocationHT_;
0289 std::pair<float, float> helixRphi_;
0290 std::pair<float, float> helixRz_;
0291 float helixD0_;
0292 unsigned int iPhiSec_;
0293 unsigned int iEtaReg_;
0294 unsigned int optoLinkID_;
0295 bool mergedHTcell_;
0296
0297
0298 TrackletSeedType seedLayerType_;
0299 unsigned int seedPS_;
0300
0301
0302 const TP* matchedTP_;
0303 std::vector<const Stub*> matchedStubs_;
0304 unsigned int nMatchedLayers_;
0305 };
0306
0307 }
0308
0309 #endif