File indexing completed on 2022-10-14 01:44:06
0001 #include "L1Trigger/TrackFindingTracklet/interface/MatchProcessor.h"
0002 #include "L1Trigger/TrackFindingTracklet/interface/Globals.h"
0003 #include "L1Trigger/TrackFindingTracklet/interface/Util.h"
0004 #include "L1Trigger/TrackFindingTracklet/interface/HistBase.h"
0005
0006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0007 #include "FWCore/Utilities/interface/Exception.h"
0008 #include "DataFormats/Math/interface/deltaPhi.h"
0009 #include "L1Trigger/TrackFindingTracklet/interface/IMATH_TrackletCalculator.h"
0010
0011 #include <filesystem>
0012
0013 using namespace std;
0014 using namespace trklet;
0015
0016 MatchProcessor::MatchProcessor(string name, Settings const& settings, Globals* global)
0017 : ProcessBase(name, settings, global),
0018 phimatchcuttable_(settings),
0019 zmatchcuttable_(settings),
0020 rphicutPStable_(settings),
0021 rphicut2Stable_(settings),
0022 rcutPStable_(settings),
0023 rcut2Stable_(settings),
0024 fullmatches_(12),
0025 rinvbendlut_(settings),
0026 luttable_(settings),
0027 inputProjBuffer_(3) {
0028 phiregion_ = name[8] - 'A';
0029
0030 layerdisk_ = initLayerDisk(3);
0031
0032 barrel_ = layerdisk_ < N_LAYER;
0033
0034 phishift_ = settings_.nphibitsstub(N_LAYER - 1) - settings_.nphibitsstub(layerdisk_);
0035 dzshift_ = settings_.nzbitsstub(0) - settings_.nzbitsstub(layerdisk_);
0036
0037 if (barrel_) {
0038 icorrshift_ = ilog2(settings_.kphi(layerdisk_) / (settings_.krbarrel() * settings_.kphider()));
0039 icorzshift_ = ilog2(settings_.kz(layerdisk_) / (settings_.krbarrel() * settings_.kzder()));
0040 } else {
0041 icorrshift_ = ilog2(settings_.kphi(layerdisk_) / (settings_.kz() * settings_.kphiderdisk()));
0042 icorzshift_ = ilog2(settings_.krprojshiftdisk() / (settings_.kz() * settings_.krder()));
0043 }
0044
0045 luttable_.initBendMatch(layerdisk_);
0046
0047 nrbits_ = 5;
0048 nphiderbits_ = 6;
0049
0050 if (!barrel_) {
0051 rinvbendlut_.initProjectionBend(
0052 global->ITC_L1L2()->der_phiD_final.K(), layerdisk_ - N_LAYER, nrbits_, nphiderbits_);
0053 }
0054
0055 nrinv_ = NRINVBITS;
0056
0057 unsigned int region = getName()[8] - 'A';
0058 assert(region < settings_.nallstubs(layerdisk_));
0059
0060 if (barrel_) {
0061 phimatchcuttable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::barrelphi, region);
0062 zmatchcuttable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::barrelz, region);
0063 } else {
0064 rphicutPStable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::diskPSphi, region);
0065 rphicut2Stable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::disk2Sphi, region);
0066 rcutPStable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::diskPSr, region);
0067 rcut2Stable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::disk2Sr, region);
0068 }
0069
0070 for (unsigned int i = 0; i < N_DSS_MOD * 2; i++) {
0071 ialphafactinner_[i] = (1 << settings_.alphashift()) * settings_.krprojshiftdisk() * settings_.half2SmoduleWidth() /
0072 (1 << (settings_.nbitsalpha() - 1)) / (settings_.rDSSinner(i) * settings_.rDSSinner(i)) /
0073 settings_.kphi();
0074 ialphafactouter_[i] = (1 << settings_.alphashift()) * settings_.krprojshiftdisk() * settings_.half2SmoduleWidth() /
0075 (1 << (settings_.nbitsalpha() - 1)) / (settings_.rDSSouter(i) * settings_.rDSSouter(i)) /
0076 settings_.kphi();
0077 }
0078
0079 nvm_ = settings_.nvmme(layerdisk_) * settings_.nallstubs(layerdisk_);
0080 nvmbins_ = settings_.nvmme(layerdisk_);
0081 nvmbits_ = settings_.nbitsvmme(layerdisk_) + settings_.nbitsallstubs(layerdisk_);
0082
0083 nMatchEngines_ = 4;
0084 for (unsigned int iME = 0; iME < nMatchEngines_; iME++) {
0085 MatchEngineUnit tmpME(settings_, barrel_, layerdisk_, luttable_);
0086 tmpME.setimeu(iME);
0087 matchengines_.push_back(tmpME);
0088 }
0089 }
0090
0091 void MatchProcessor::addOutput(MemoryBase* memory, string output) {
0092 if (settings_.writetrace()) {
0093 edm::LogVerbatim("Tracklet") << "In " << name_ << " adding output to " << memory->getName() << " to output "
0094 << output;
0095 }
0096 if (output.find("matchout") != std::string::npos) {
0097 auto* tmp = dynamic_cast<FullMatchMemory*>(memory);
0098 assert(tmp != nullptr);
0099 unsigned int iSeed = getISeed(tmp->getName());
0100 assert(iSeed < fullmatches_.size());
0101 assert(fullmatches_[iSeed] == nullptr);
0102 fullmatches_[iSeed] = tmp;
0103 return;
0104 }
0105 throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " could not find output: " << output;
0106 }
0107
0108 void MatchProcessor::addInput(MemoryBase* memory, string input) {
0109 if (settings_.writetrace()) {
0110 edm::LogVerbatim("Tracklet") << "In " << name_ << " adding input from " << memory->getName() << " to input "
0111 << input;
0112 }
0113 if (input == "allstubin") {
0114 auto* tmp = dynamic_cast<AllStubsMemory*>(memory);
0115 assert(tmp != nullptr);
0116 allstubs_ = tmp;
0117 return;
0118 }
0119 if (input == "vmstubin") {
0120 auto* tmp = dynamic_cast<VMStubsMEMemory*>(memory);
0121 assert(tmp != nullptr);
0122 vmstubs_.push_back(tmp);
0123 return;
0124 }
0125 if (input == "projin") {
0126 auto* tmp = dynamic_cast<TrackletProjectionsMemory*>(memory);
0127 assert(tmp != nullptr);
0128 inputprojs_.push_back(tmp);
0129 return;
0130 }
0131 throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " could not find input: " << input;
0132 }
0133
0134 void MatchProcessor::execute(unsigned int iSector, double phimin) {
0135 assert(vmstubs_.size() == 1);
0136
0137
0138
0139
0140
0141
0142
0143
0144
0145
0146
0147
0148
0149
0150
0151
0152 bool print = getName() == "MP_L3PHIC" && iSector == 3;
0153 print = false;
0154
0155 phimin_ = phimin;
0156
0157 Tracklet* oldTracklet = nullptr;
0158
0159 unsigned int countme = 0;
0160 unsigned int countall = 0;
0161 unsigned int countsel = 0;
0162 unsigned int countinputproj = 0;
0163
0164 unsigned int iprojmem = 0;
0165 while (iprojmem < inputprojs_.size() && inputprojs_[iprojmem]->nTracklets() == 0) {
0166 iprojmem++;
0167 }
0168
0169 unsigned int iproj = 0;
0170
0171 inputProjBuffer_.reset();
0172
0173 for (const auto& inputproj : inputprojs_) {
0174 countinputproj += inputproj->nTracklets();
0175 }
0176
0177 for (auto& matchengine : matchengines_) {
0178 matchengine.reset();
0179 }
0180
0181 ProjectionTemp tmpProj_, tmpProj__;
0182 bool good_ = false;
0183 bool good__ = false;
0184
0185 for (unsigned int istep = 0; istep < settings_.maxStep("MP"); istep++) {
0186
0187
0188
0189
0190
0191
0192
0193
0194
0195
0196
0197
0198
0199
0200
0201
0202
0203
0204 bool projdone = false;
0205
0206 bool projBuffNearFull = inputProjBuffer_.nearfull();
0207
0208 for (unsigned int iME = 0; iME < nMatchEngines_; iME++) {
0209 matchengines_[iME].setAlmostFull();
0210 }
0211
0212
0213
0214
0215 unsigned int iMEbest = 0;
0216 int bestTCID = matchengines_[0].TCID();
0217 bool meactive = matchengines_[0].active();
0218 for (unsigned int iME = 1; iME < nMatchEngines_; iME++) {
0219 meactive = meactive || matchengines_[iME].active();
0220 int tcid = matchengines_[iME].TCID();
0221 if (tcid < bestTCID) {
0222 bestTCID = tcid;
0223 iMEbest = iME;
0224 }
0225 }
0226
0227
0228
0229 if (!matchengines_[iMEbest].empty()) {
0230 std::pair<Tracklet*, const Stub*> candmatch = matchengines_[iMEbest].read();
0231
0232 const Stub* fpgastub = candmatch.second;
0233 Tracklet* tracklet = candmatch.first;
0234
0235
0236 if (oldTracklet != nullptr) {
0237
0238
0239 assert(oldTracklet->TCID() <= tracklet->TCID());
0240 }
0241 oldTracklet = tracklet;
0242
0243 bool match = matchCalculator(tracklet, fpgastub, print, istep);
0244
0245 if (settings_.debugTracklet() && match) {
0246 edm::LogVerbatim("Tracklet") << getName() << " have match";
0247 }
0248
0249 countall++;
0250 if (match)
0251 countsel++;
0252 }
0253
0254
0255
0256
0257 bool addedProjection = false;
0258 for (unsigned int iME = 0; iME < nMatchEngines_; iME++) {
0259 if (!matchengines_[iME].idle())
0260 countme++;
0261
0262 if ((!addedProjection) && matchengines_[iME].idle() && (!inputProjBuffer_.empty())) {
0263 ProjectionTemp tmpProj = inputProjBuffer_.read();
0264 VMStubsMEMemory* stubmem = vmstubs_[0];
0265
0266 if (settings_.debugTracklet()) {
0267 edm::LogVerbatim("Tracklet") << getName() << " adding projection to match engine";
0268 }
0269
0270 int nbins = (1 << N_RZBITS);
0271 if (layerdisk_ >= N_LAYER) {
0272 nbins *= 2;
0273 }
0274
0275 matchengines_[iME].init(stubmem,
0276 nbins,
0277 tmpProj.slot(),
0278 tmpProj.iphi(),
0279 tmpProj.shift(),
0280 tmpProj.projrinv(),
0281 tmpProj.projfinerz(),
0282 tmpProj.projfinephi(),
0283 tmpProj.use(0, 0),
0284 tmpProj.use(0, 1),
0285 tmpProj.use(1, 0),
0286 tmpProj.use(1, 1),
0287 tmpProj.isPSseed(),
0288 tmpProj.proj());
0289 meactive = true;
0290 addedProjection = true;
0291 } else {
0292 matchengines_[iME].step();
0293 }
0294 matchengines_[iME].processPipeline();
0295 }
0296
0297
0298
0299
0300
0301 if (good__) {
0302 inputProjBuffer_.store(tmpProj__);
0303 }
0304
0305 good__ = good_;
0306 tmpProj__ = tmpProj_;
0307
0308 good_ = false;
0309
0310 if (iprojmem < inputprojs_.size()) {
0311 TrackletProjectionsMemory* projMem = inputprojs_[iprojmem];
0312 if (!projBuffNearFull) {
0313 if (settings_.debugTracklet()) {
0314 edm::LogVerbatim("Tracklet") << getName() << " have projection in memory : " << projMem->getName();
0315 }
0316
0317 Tracklet* proj = projMem->getTracklet(iproj);
0318
0319 FPGAWord fpgaphi = proj->proj(layerdisk_).fpgaphiproj();
0320
0321 unsigned int iphi = (fpgaphi.value() >> (fpgaphi.nbits() - nvmbits_)) & (nvmbins_ - 1);
0322
0323 int nextrabits = 2;
0324 int overlapbits = nvmbits_ + nextrabits;
0325
0326 unsigned int extrabits = fpgaphi.bits(fpgaphi.nbits() - overlapbits - nextrabits, nextrabits);
0327
0328 unsigned int ivmPlus = iphi;
0329
0330 int shift = 0;
0331
0332 if (extrabits == ((1U << nextrabits) - 1) && iphi != ((1U << settings_.nbitsvmme(layerdisk_)) - 1)) {
0333 shift = 1;
0334 ivmPlus++;
0335 }
0336 unsigned int ivmMinus = iphi;
0337 if (extrabits == 0 && iphi != 0) {
0338 shift = -1;
0339 ivmMinus--;
0340 }
0341
0342 int projrinv = -1;
0343 if (barrel_) {
0344 FPGAWord phider = proj->proj(layerdisk_).fpgaphiprojder();
0345 projrinv = (1 << (nrinv_ - 1)) - 1 - (phider.value() >> (phider.nbits() - nrinv_));
0346 } else {
0347
0348
0349
0350
0351
0352 int rindex =
0353 (proj->proj(layerdisk_).fpgarzproj().value() >> (proj->proj(layerdisk_).fpgarzproj().nbits() - nrbits_)) &
0354 ((1 << nrbits_) - 1);
0355
0356 int phiprojder = proj->proj(layerdisk_).fpgaphiprojder().value();
0357
0358 int phiderindex = (phiprojder >> (proj->proj(layerdisk_).fpgaphiprojder().nbits() - nphiderbits_)) &
0359 ((1 << nphiderbits_) - 1);
0360
0361 int signindex = proj->proj(layerdisk_).fpgarzprojder().value() < 0;
0362
0363 int bendindex = (signindex << (nphiderbits_ + nrbits_)) + (rindex << (nphiderbits_)) + phiderindex;
0364
0365 projrinv = rinvbendlut_.lookup(bendindex);
0366
0367 proj->proj(layerdisk_).setBendIndex(projrinv);
0368 }
0369 assert(projrinv >= 0);
0370
0371 unsigned int slot = proj->proj(layerdisk_).fpgarzbin1projvm().value();
0372 bool second = proj->proj(layerdisk_).fpgarzbin2projvm().value();
0373
0374 unsigned int projfinephi =
0375 (fpgaphi.value() >> (fpgaphi.nbits() - (nvmbits_ + NFINEPHIBITS))) & ((1 << NFINEPHIBITS) - 1);
0376 int projfinerz = proj->proj(layerdisk_).fpgafinerzvm().value();
0377
0378 bool isPSseed = proj->PSseed();
0379
0380 int nbins = (1 << N_RZBITS);
0381 if (layerdisk_ >= N_LAYER) {
0382 nbins *= 2;
0383 }
0384
0385 VMStubsMEMemory* stubmem = vmstubs_[0];
0386 bool usefirstMinus = stubmem->nStubsBin(ivmMinus * nbins + slot) != 0;
0387 bool usesecondMinus = (second && (stubmem->nStubsBin(ivmMinus * nbins + slot + 1) != 0));
0388 bool usefirstPlus = ivmPlus != ivmMinus && stubmem->nStubsBin(ivmPlus * nbins + slot) != 0;
0389 bool usesecondPlus = ivmPlus != ivmMinus && (second && (stubmem->nStubsBin(ivmPlus * nbins + slot + 1) != 0));
0390
0391 good_ = usefirstPlus || usesecondPlus || usefirstMinus || usesecondMinus;
0392
0393 if (good_) {
0394 ProjectionTemp tmpProj(proj,
0395 slot,
0396 projrinv,
0397 projfinerz,
0398 projfinephi,
0399 ivmMinus,
0400 shift,
0401 usefirstMinus,
0402 usefirstPlus,
0403 usesecondMinus,
0404 usesecondPlus,
0405 isPSseed);
0406 tmpProj_ = tmpProj;
0407 }
0408
0409 iproj++;
0410 if (iproj == projMem->nTracklets()) {
0411 iproj = 0;
0412 do {
0413 iprojmem++;
0414 } while (iprojmem < inputprojs_.size() && inputprojs_[iprojmem]->nTracklets() == 0);
0415 }
0416
0417 } else {
0418 projdone = true && !good_ && !good__;
0419 }
0420 }
0421
0422
0423
0424
0425
0426
0427
0428 if ((projdone && !meactive && inputProjBuffer_.rptr() == inputProjBuffer_.wptr()) ||
0429 (istep == settings_.maxStep("MP") - 1)) {
0430 if (settings_.writeMonitorData("MP")) {
0431 globals_->ofstream("matchprocessor.txt") << getName() << " " << istep << " " << countall << " " << countsel
0432 << " " << countme << " " << countinputproj << endl;
0433 }
0434 break;
0435 }
0436 }
0437
0438 if (settings_.writeMonitorData("MC")) {
0439 globals_->ofstream("matchcalculator.txt") << getName() << " " << countall << " " << countsel << endl;
0440 }
0441 }
0442
0443 bool MatchProcessor::matchCalculator(Tracklet* tracklet, const Stub* fpgastub, bool, unsigned int) {
0444 const L1TStub* stub = fpgastub->l1tstub();
0445
0446 if (layerdisk_ < N_LAYER) {
0447 const Projection& proj = tracklet->proj(layerdisk_);
0448 int ir = fpgastub->r().value();
0449 int iphi = proj.fpgaphiproj().value();
0450 int icorr = (ir * proj.fpgaphiprojder().value()) >> icorrshift_;
0451 iphi += icorr;
0452
0453 int iz = proj.fpgarzproj().value();
0454 int izcor = (ir * proj.fpgarzprojder().value() + (1 << (icorzshift_ - 1))) >> icorzshift_;
0455 iz += izcor;
0456
0457 int ideltaz = fpgastub->z().value() - iz;
0458 int ideltaphi = (fpgastub->phi().value() - iphi) << phishift_;
0459
0460
0461
0462 double phi = stub->phi();
0463 double r = stub->r();
0464 double z = stub->z();
0465
0466 if (settings_.useapprox()) {
0467 double dphi = reco::reduceRange(phi - fpgastub->phiapprox(phimin_, 0.0));
0468 assert(std::abs(dphi) < 0.001);
0469 phi = fpgastub->phiapprox(phimin_, 0.0);
0470 z = fpgastub->zapprox();
0471 r = fpgastub->rapprox();
0472 }
0473
0474 if (phi < 0)
0475 phi += 2 * M_PI;
0476 phi -= phimin_;
0477
0478 double dr = r - settings_.rmean(layerdisk_);
0479 assert(std::abs(dr) < settings_.drmax());
0480
0481 double dphi = reco::reduceRange(phi - (proj.phiproj() + dr * proj.phiprojder()));
0482
0483 double dz = z - (proj.rzproj() + dr * proj.rzprojder());
0484
0485 double dphiapprox = reco::reduceRange(phi - (proj.phiprojapprox() + dr * proj.phiprojderapprox()));
0486
0487 double dzapprox = z - (proj.rzprojapprox() + dr * proj.rzprojderapprox());
0488
0489 int seedindex = tracklet->getISeed();
0490
0491 assert(phimatchcuttable_.lookup(seedindex) > 0);
0492 assert(zmatchcuttable_.lookup(seedindex) > 0);
0493
0494 if (settings_.bookHistos()) {
0495 bool truthmatch = tracklet->stubtruthmatch(stub);
0496
0497 HistBase* hists = globals_->histograms();
0498 hists->FillLayerResidual(layerdisk_ + 1,
0499 seedindex,
0500 dphiapprox * settings_.rmean(layerdisk_),
0501 ideltaphi * settings_.kphi1() * settings_.rmean(layerdisk_),
0502 (ideltaz << dzshift_) * settings_.kz(),
0503 dz,
0504 truthmatch);
0505 }
0506
0507 if (settings_.writeMonitorData("Residuals")) {
0508 double pt = 0.01 * settings_.c() * settings_.bfield() / std::abs(tracklet->rinv());
0509
0510 globals_->ofstream("layerresiduals.txt")
0511 << layerdisk_ + 1 << " " << seedindex << " " << pt << " "
0512 << ideltaphi * settings_.kphi1() * settings_.rmean(layerdisk_) << " "
0513 << dphiapprox * settings_.rmean(layerdisk_) << " "
0514 << phimatchcuttable_.lookup(seedindex) * settings_.kphi1() * settings_.rmean(layerdisk_) << " "
0515 << (ideltaz << dzshift_) * settings_.kz() << " " << dz << " "
0516 << zmatchcuttable_.lookup(seedindex) * settings_.kz() << endl;
0517 }
0518
0519 bool imatch = (std::abs(ideltaphi) <= phimatchcuttable_.lookup(seedindex)) &&
0520 (ideltaz << dzshift_ < zmatchcuttable_.lookup(seedindex)) &&
0521 (ideltaz << dzshift_ >= -zmatchcuttable_.lookup(seedindex));
0522
0523 if (settings_.debugTracklet()) {
0524 edm::LogVerbatim("Tracklet") << getName() << " imatch = " << imatch << " ideltaphi cut " << ideltaphi << " "
0525 << phimatchcuttable_.lookup(seedindex) << " ideltaz<<dzshift cut "
0526 << (ideltaz << dzshift_) << " " << zmatchcuttable_.lookup(seedindex);
0527 }
0528
0529
0530 if (std::abs(dphi) > 0.5 * settings_.dphisectorHG() || std::abs(dphiapprox) > 0.5 * settings_.dphisectorHG()) {
0531 throw cms::Exception("LogicError") << "WARNING dphi and/or dphiapprox too large : " << dphi << " " << dphiapprox
0532 << endl;
0533 }
0534
0535 bool keep = true;
0536 if (!settings_.doKF() || !settings_.doMultipleMatches()) {
0537
0538 if (imatch && tracklet->match(layerdisk_)) {
0539
0540 auto res = tracklet->resid(layerdisk_);
0541 keep = abs(ideltaphi) < abs(res.fpgaphiresid().value());
0542 imatch = keep;
0543 }
0544 }
0545
0546 if (imatch) {
0547 tracklet->addMatch(layerdisk_,
0548 ideltaphi,
0549 ideltaz,
0550 dphi,
0551 dz,
0552 dphiapprox,
0553 dzapprox,
0554 (phiregion_ << N_BITSMEMADDRESS) + fpgastub->stubindex().value(),
0555 fpgastub);
0556
0557 if (settings_.debugTracklet()) {
0558 edm::LogVerbatim("Tracklet") << "Accepted full match in layer " << getName() << " " << tracklet;
0559 }
0560
0561 int iSeed = tracklet->getISeed();
0562 assert(fullmatches_[iSeed] != nullptr);
0563 fullmatches_[iSeed]->addMatch(tracklet, fpgastub);
0564
0565 return true;
0566 } else {
0567 return false;
0568 }
0569 } else {
0570
0571
0572 assert(stub->z() * tracklet->t() > 0.0);
0573
0574 int sign = (tracklet->t() > 0.0) ? 1 : -1;
0575 int disk = sign * (layerdisk_ - N_LAYER + 1);
0576 assert(disk != 0);
0577
0578
0579
0580 int iz = fpgastub->z().value();
0581
0582 const Projection& proj = tracklet->proj(layerdisk_);
0583
0584 int iphi = proj.fpgaphiproj().value();
0585 int iphicorr = (iz * proj.fpgaphiprojder().value()) >> icorrshift_;
0586
0587 iphi += iphicorr;
0588
0589 int ir = proj.fpgarzproj().value();
0590 int ircorr = (iz * proj.fpgarzprojder().value()) >> icorzshift_;
0591 ir += ircorr;
0592
0593 int ideltaphi = fpgastub->phi().value() - iphi;
0594
0595 int irstub = fpgastub->r().value();
0596 int ialphafact = 0;
0597 if (!stub->isPSmodule()) {
0598 assert(irstub < (int)N_DSS_MOD * 2);
0599 if (layerdisk_ - N_LAYER <= 1) {
0600 ialphafact = ialphafactinner_[irstub];
0601 irstub = settings_.rDSSinner(irstub) / settings_.kr();
0602 } else {
0603 ialphafact = ialphafactouter_[irstub];
0604 irstub = settings_.rDSSouter(irstub) / settings_.kr();
0605 }
0606 }
0607
0608 int ideltar = (irstub * settings_.kr()) / settings_.krprojshiftdisk() - ir;
0609
0610 if (!stub->isPSmodule()) {
0611 int ialpha = fpgastub->alpha().value();
0612 int iphialphacor = ((ideltar * ialpha * ialphafact) >> settings_.alphashift());
0613 ideltaphi += iphialphacor;
0614 }
0615
0616
0617
0618 double phi = stub->phi();
0619 double z = stub->z();
0620 double r = stub->r();
0621
0622 if (settings_.useapprox()) {
0623 double dphi = reco::reduceRange(phi - fpgastub->phiapprox(phimin_, 0.0));
0624 assert(std::abs(dphi) < 0.001);
0625 phi = fpgastub->phiapprox(phimin_, 0.0);
0626 z = fpgastub->zapprox();
0627 r = fpgastub->rapprox();
0628 }
0629
0630 if (phi < 0)
0631 phi += 2 * M_PI;
0632 phi -= phimin_;
0633
0634 double dz = z - sign * settings_.zmean(layerdisk_ - N_LAYER);
0635
0636 if (std::abs(dz) > settings_.dzmax()) {
0637 edm::LogProblem("Tracklet") << __FILE__ << ":" << __LINE__ << " " << name_ << " " << tracklet->getISeed();
0638 edm::LogProblem("Tracklet") << "stub " << stub->z() << " disk " << disk << " " << dz;
0639 assert(std::abs(dz) < settings_.dzmax());
0640 }
0641
0642 double phiproj = proj.phiproj() + dz * proj.phiprojder();
0643 double rproj = proj.rzproj() + dz * proj.rzprojder();
0644 double deltar = r - rproj;
0645
0646 double dr = stub->r() - rproj;
0647 double drapprox = stub->r() - (proj.rzprojapprox() + dz * proj.rzprojderapprox());
0648
0649 double dphi = reco::reduceRange(phi - phiproj);
0650
0651 double dphiapprox = reco::reduceRange(phi - (proj.phiprojapprox() + dz * proj.phiprojderapprox()));
0652
0653 double drphi = dphi * stub->r();
0654 double drphiapprox = dphiapprox * stub->r();
0655
0656 if (!stub->isPSmodule()) {
0657 double alphanorm = stub->alphanorm();
0658 dphi += dr * alphanorm * settings_.half2SmoduleWidth() / stub->r2();
0659 ;
0660 dphiapprox += drapprox * alphanorm * settings_.half2SmoduleWidth() / stub->r2();
0661
0662 drphi += dr * alphanorm * settings_.half2SmoduleWidth() / stub->r();
0663 drphiapprox += dr * alphanorm * settings_.half2SmoduleWidth() / stub->r();
0664 }
0665
0666 int seedindex = tracklet->getISeed();
0667
0668 int idrphicut = rphicutPStable_.lookup(seedindex);
0669 int idrcut = rcutPStable_.lookup(seedindex);
0670 if (!stub->isPSmodule()) {
0671 idrphicut = rphicut2Stable_.lookup(seedindex);
0672 idrcut = rcut2Stable_.lookup(seedindex);
0673 }
0674
0675 double drphicut = idrphicut * settings_.kphi() * settings_.kr();
0676 double drcut = idrcut * settings_.krprojshiftdisk();
0677
0678 if (settings_.writeMonitorData("Residuals")) {
0679 double pt = 0.01 * settings_.c() * settings_.bfield() / std::abs(tracklet->rinv());
0680
0681 globals_->ofstream("diskresiduals.txt")
0682 << layerdisk_ - N_LAYER + 1 << " " << stub->isPSmodule() << " " << tracklet->layer() << " "
0683 << abs(tracklet->disk()) << " " << pt << " " << ideltaphi * settings_.kphi() * stub->r() << " " << drphiapprox
0684 << " " << drphicut << " " << ideltar * settings_.krprojshiftdisk() << " " << deltar << " " << drcut << " "
0685 << endl;
0686 }
0687
0688 bool match = (std::abs(drphi) < drphicut) && (std::abs(deltar) < drcut);
0689 bool imatch = (std::abs(ideltaphi * irstub) < idrphicut) && (std::abs(ideltar) < idrcut);
0690
0691 if (settings_.debugTracklet()) {
0692 edm::LogVerbatim("Tracklet") << "imatch match disk: " << imatch << " " << match << " " << std::abs(ideltaphi)
0693 << " " << drphicut / (settings_.kphi() * stub->r()) << " " << std::abs(ideltar)
0694 << " " << drcut / settings_.krprojshiftdisk() << " r = " << stub->r();
0695 }
0696
0697 bool keep = true;
0698 if (!settings_.doKF() || !settings_.doMultipleMatches()) {
0699
0700 if (imatch && tracklet->match(layerdisk_)) {
0701
0702 auto res = tracklet->resid(layerdisk_);
0703 keep = abs(ideltaphi) < abs(res.fpgaphiresid().value());
0704 imatch = keep;
0705 }
0706 }
0707
0708 if (imatch) {
0709 if (settings_.debugTracklet()) {
0710 edm::LogVerbatim("Tracklet") << "MatchCalculator found match in disk " << getName();
0711 }
0712
0713 if (std::abs(dphi) >= third * settings_.dphisectorHG()) {
0714 edm::LogPrint("Tracklet") << "dphi " << dphi << " ISeed " << tracklet->getISeed();
0715 }
0716 assert(std::abs(dphi) < third * settings_.dphisectorHG());
0717 assert(std::abs(dphiapprox) < third * settings_.dphisectorHG());
0718
0719 tracklet->addMatch(layerdisk_,
0720 ideltaphi,
0721 ideltar,
0722 drphi / stub->r(),
0723 dr,
0724 drphiapprox / stub->r(),
0725 drapprox,
0726 (phiregion_ << N_BITSMEMADDRESS) + fpgastub->stubindex().value(),
0727 fpgastub);
0728
0729 if (settings_.debugTracklet()) {
0730 edm::LogVerbatim("Tracklet") << "Accepted full match in disk " << getName() << " " << tracklet;
0731 }
0732
0733 int iSeed = tracklet->getISeed();
0734 assert(fullmatches_[iSeed] != nullptr);
0735 fullmatches_[iSeed]->addMatch(tracklet, fpgastub);
0736
0737 return true;
0738 } else {
0739 return false;
0740 }
0741 }
0742 }