File indexing completed on 2025-01-27 02:50:20
0001 #include "L1Trigger/TrackFindingTracklet/interface/HybridFit.h"
0002 #include "L1Trigger/TrackFindingTracklet/interface/Tracklet.h"
0003 #include "L1Trigger/TrackFindingTracklet/interface/Stub.h"
0004 #include "L1Trigger/TrackFindingTracklet/interface/L1TStub.h"
0005
0006 #ifdef USEHYBRID
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008 #include "DataFormats/Math/interface/deltaPhi.h"
0009
0010 using namespace std;
0011 using namespace trklet;
0012
0013 HybridFit::HybridFit(unsigned int iSector, Settings const& settings, Globals* globals) : settings_(settings) {
0014 iSector_ = iSector;
0015 globals_ = globals;
0016 }
0017
0018 void HybridFit::Fit(Tracklet* tracklet, std::vector<const Stub*>& trackstublist) {
0019 if (settings_.fakefit()) {
0020 vector<const L1TStub*> l1stubsFromFitTrack;
0021 for (unsigned int k = 0; k < trackstublist.size(); k++) {
0022 const L1TStub* L1stub = trackstublist[k]->l1tstub();
0023 l1stubsFromFitTrack.push_back(L1stub);
0024 }
0025 tracklet->setFitPars(tracklet->rinvapprox(),
0026 tracklet->phi0approx(),
0027 tracklet->d0approx(),
0028 tracklet->tapprox(),
0029 tracklet->z0approx(),
0030 0.,
0031 0.,
0032 tracklet->rinv(),
0033 tracklet->phi0(),
0034 tracklet->d0(),
0035 tracklet->t(),
0036 tracklet->z0(),
0037 0.,
0038 0.,
0039 tracklet->fpgarinv().value(),
0040 tracklet->fpgaphi0().value(),
0041 tracklet->fpgad0().value(),
0042 tracklet->fpgat().value(),
0043 tracklet->fpgaz0().value(),
0044 0,
0045 0,
0046 0,
0047 l1stubsFromFitTrack);
0048 return;
0049 }
0050
0051 std::vector<tmtt::Stub*> TMTTstubs;
0052 std::map<unsigned int, const L1TStub*> L1StubIndices;
0053 unsigned int L1stubID = 0;
0054
0055 if (globals_->tmttSettings() == nullptr) {
0056 if (settings_.printDebugKF())
0057 edm::LogVerbatim("L1track") << "Creating TMTT::Settings in HybridFit::Fit";
0058 globals_->tmttSettings() = make_unique<tmtt::Settings>();
0059 globals_->tmttSettings()->setMagneticField(settings_.bfield());
0060 }
0061
0062 const tmtt::Settings& TMTTsettings = *globals_->tmttSettings();
0063
0064 int kf_phi_sec = iSector_;
0065
0066 for (unsigned int k = 0; k < trackstublist.size(); k++) {
0067 const L1TStub* L1stubptr = trackstublist[k]->l1tstub();
0068
0069 double kfphi = L1stubptr->phi();
0070 double kfr = L1stubptr->r();
0071 double kfz = L1stubptr->z();
0072 double kfbend = L1stubptr->bend();
0073 bool psmodule = L1stubptr->isPSmodule();
0074 unsigned int iphi = L1stubptr->iphi();
0075 double alpha = L1stubptr->alpha(settings_.stripPitch(psmodule));
0076 bool isTilted = L1stubptr->isTilted();
0077
0078 bool isBarrel = trackstublist[k]->layerdisk() < N_LAYER;
0079 int kflayer;
0080
0081 if (isBarrel) {
0082 kflayer = L1stubptr->layer() + 1;
0083 if (settings_.printDebugKF())
0084 edm::LogVerbatim("L1track") << "Will create layer stub with : ";
0085 } else {
0086 kflayer = abs(L1stubptr->disk());
0087 if (kfz > 0) {
0088 kflayer += 10;
0089 } else {
0090 kflayer += 20;
0091 }
0092 if (settings_.printDebugKF())
0093 edm::LogVerbatim("L1track") << "Will create disk stub with : ";
0094 }
0095
0096 float stripPitch = settings_.stripPitch(psmodule);
0097 float stripLength = settings_.stripLength(psmodule);
0098 unsigned int nStrips = settings_.nStrips(psmodule);
0099
0100 if (settings_.printDebugKF()) {
0101 edm::LogVerbatim("L1track") << kfphi << " " << kfr << " " << kfz << " " << kfbend << " " << kflayer << " "
0102 << isBarrel << " " << psmodule << " " << isTilted << " \n"
0103 << stripPitch << " " << stripLength << " " << nStrips;
0104 }
0105
0106 unsigned int uniqueStubIndex = 1000 * L1stubID + L1stubptr->allStubIndex();
0107 tmtt::Stub* TMTTstubptr = new tmtt::Stub(&TMTTsettings,
0108 uniqueStubIndex,
0109 kfphi,
0110 kfr,
0111 kfz,
0112 kfbend,
0113 iphi,
0114 -alpha,
0115 kflayer,
0116 kf_phi_sec,
0117 psmodule,
0118 isBarrel,
0119 isTilted,
0120 stripPitch,
0121 stripLength,
0122 nStrips);
0123 TMTTstubs.push_back(TMTTstubptr);
0124 L1StubIndices[uniqueStubIndex] = L1stubptr;
0125 L1stubID++;
0126 }
0127
0128 if (settings_.printDebugKF()) {
0129 edm::LogVerbatim("L1track") << "Made TMTTstubs: trackstublist.size() = " << trackstublist.size();
0130 }
0131
0132 double kfrinv = tracklet->rinvapprox();
0133 double kfphi0 = tracklet->phi0approx();
0134 double kfz0 = tracklet->z0approx();
0135 double kft = tracklet->tapprox();
0136 double kfd0 = tracklet->d0approx();
0137
0138 if (settings_.printDebugKF()) {
0139 edm::LogVerbatim("L1track") << "tracklet phi0 = " << kfphi0 << "\n"
0140 << "iSector = " << iSector_ << "\n"
0141 << "dphisectorHG = " << settings_.dphisectorHG();
0142 }
0143
0144
0145 kfphi0 = reco::reducePhiRange(kfphi0 + iSector_ * settings_.dphisector() - 0.5 * settings_.dphisectorHG());
0146
0147 std::pair<float, float> helixrphi(kfrinv / (0.01 * settings_.c() * settings_.bfield()), kfphi0);
0148 std::pair<float, float> helixrz(kfz0, kft);
0149
0150
0151 double chargeOverPt = helixrphi.first;
0152 int mBin = std::floor(TMTTsettings.houghNbinsPt() / 2) +
0153 std::floor((TMTTsettings.houghNbinsPt() / 2) * chargeOverPt / (1. / TMTTsettings.houghMinPt()));
0154 mBin = max(min(mBin, int(TMTTsettings.houghNbinsPt() - 1)), 0);
0155 std::pair<unsigned int, unsigned int> celllocation(mBin, 1);
0156
0157
0158 const vector<double> etaRegions = TMTTsettings.etaRegions();
0159 const float chosenRofZ = TMTTsettings.chosenRofZ();
0160
0161 float kfzRef = kfz0 + chosenRofZ * kft;
0162
0163 unsigned int kf_eta_reg = 0;
0164 for (unsigned int iEtaSec = 1; iEtaSec < etaRegions.size() - 1; iEtaSec++) {
0165 const float etaMax = etaRegions[iEtaSec];
0166 const float zRefMax = chosenRofZ / tan(2. * atan(exp(-etaMax)));
0167 if (kfzRef > zRefMax)
0168 kf_eta_reg = iEtaSec;
0169 }
0170
0171 tmtt::L1track3D l1track3d(
0172 &TMTTsettings, TMTTstubs, celllocation, helixrphi, helixrz, kfd0, kf_phi_sec, kf_eta_reg, 1, false);
0173 unsigned int seedType = tracklet->getISeed();
0174 unsigned int numPS = tracklet->PSseed();
0175 l1track3d.setSeedLayerType(seedType);
0176 l1track3d.setSeedPS(numPS);
0177
0178 if (globals_->tmttKFParamsComb() == nullptr) {
0179 if (settings_.printDebugKF())
0180 edm::LogVerbatim("L1track") << "Will make KFParamsComb for " << settings_.nHelixPar() << " param fit";
0181 globals_->tmttKFParamsComb() = make_unique<tmtt::KFParamsComb>(&TMTTsettings, settings_.nHelixPar(), "KFfitter");
0182 }
0183
0184 tmtt::KFParamsComb& fitterKF = *globals_->tmttKFParamsComb();
0185
0186
0187 tmtt::L1fittedTrack fittedTrk = fitterKF.fit(l1track3d);
0188
0189 if (fittedTrk.accepted()) {
0190 if (fittedTrk.consistentSector()) {
0191 tmtt::KFTrackletTrack trk = fittedTrk.returnKFTrackletTrack();
0192
0193 if (settings_.printDebugKF())
0194 edm::LogVerbatim("L1track") << "Done with Kalman fit. Pars: pt = " << trk.pt()
0195 << ", 1/2R = " << settings_.bfield() * 3 * trk.qOverPt() / 2000
0196 << ", phi0 = " << trk.phi0() << ", eta = " << trk.eta() << ", z0 = " << trk.z0()
0197 << ", chi2 = " << trk.chi2() << ", accepted = " << trk.accepted();
0198
0199 double d0, chi2rphi, phi0, qoverpt = -999;
0200 if (trk.done_bcon()) {
0201 d0 = trk.d0_bcon();
0202 chi2rphi = trk.chi2rphi_bcon();
0203 phi0 = trk.phi0_bcon();
0204 qoverpt = trk.qOverPt_bcon();
0205 } else {
0206 d0 = trk.d0();
0207 chi2rphi = trk.chi2rphi();
0208 phi0 = trk.phi0();
0209 qoverpt = trk.qOverPt();
0210 }
0211
0212
0213 double phi0fit = reco::reducePhiRange(phi0 - iSector_ * 2 * M_PI / N_SECTOR + 0.5 * settings_.dphisectorHG());
0214 double rinvfit = 0.01 * settings_.c() * settings_.bfield() * qoverpt;
0215
0216 int irinvfit = floor(rinvfit / settings_.krinvpars());
0217 int iphi0fit = floor(phi0fit / settings_.kphi0pars());
0218 int itanlfit = floor(trk.tanLambda() / settings_.ktpars());
0219 int iz0fit = floor(trk.z0() / settings_.kz0pars());
0220 int id0fit = floor(d0 / settings_.kd0pars());
0221 int ichi2rphifit = chi2rphi / 16;
0222 int ichi2rzfit = trk.chi2rz() / 16;
0223
0224 const vector<const tmtt::Stub*>& stubsFromFit = trk.stubs();
0225 vector<const L1TStub*> l1stubsFromFit;
0226 for (const tmtt::Stub* s : stubsFromFit) {
0227 unsigned int IDf = s->index();
0228 const L1TStub* l1s = L1StubIndices.at(IDf);
0229 l1stubsFromFit.push_back(l1s);
0230 }
0231
0232 tracklet->setFitPars(rinvfit,
0233 phi0fit,
0234 d0,
0235 trk.tanLambda(),
0236 trk.z0(),
0237 chi2rphi,
0238 trk.chi2rz(),
0239 rinvfit,
0240 phi0fit,
0241 d0,
0242 trk.tanLambda(),
0243 trk.z0(),
0244 chi2rphi,
0245 trk.chi2rz(),
0246 irinvfit,
0247 iphi0fit,
0248 id0fit,
0249 itanlfit,
0250 iz0fit,
0251 ichi2rphifit,
0252 ichi2rzfit,
0253 trk.hitPattern(),
0254 l1stubsFromFit);
0255 } else {
0256 if (settings_.printDebugKF()) {
0257 edm::LogVerbatim("L1track") << "FitTrack:KF rejected track";
0258 }
0259 }
0260 } else {
0261 if (settings_.printDebugKF()) {
0262 edm::LogVerbatim("L1track") << "FitTrack:KF rejected track";
0263 }
0264 }
0265
0266 for (const tmtt::Stub* s : TMTTstubs) {
0267 delete s;
0268 }
0269 }
0270 #endif