File indexing completed on 2025-01-27 02:50:21
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include "L1Trigger/TrackFindingTracklet/interface/MatchCalculator.h"
0012 #include "L1Trigger/TrackFindingTracklet/interface/Globals.h"
0013 #include "L1Trigger/TrackFindingTracklet/interface/Util.h"
0014 #include "L1Trigger/TrackFindingTracklet/interface/CandidateMatchMemory.h"
0015 #include "L1Trigger/TrackFindingTracklet/interface/FullMatchMemory.h"
0016 #include "L1Trigger/TrackFindingTracklet/interface/AllStubsMemory.h"
0017 #include "L1Trigger/TrackFindingTracklet/interface/AllProjectionsMemory.h"
0018 #include "L1Trigger/TrackFindingTracklet/interface/Tracklet.h"
0019 #include "L1Trigger/TrackFindingTracklet/interface/Stub.h"
0020 #include "L1Trigger/TrackFindingTracklet/interface/HistBase.h"
0021
0022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0023 #include "FWCore/Utilities/interface/Exception.h"
0024 #include "DataFormats/Math/interface/deltaPhi.h"
0025
0026 #include <filesystem>
0027 #include <algorithm>
0028
0029 using namespace std;
0030 using namespace trklet;
0031
0032 MatchCalculator::MatchCalculator(string name, Settings const& settings, Globals* global)
0033 : ProcessBase(name, settings, global),
0034 phimatchcuttable_(settings),
0035 zmatchcuttable_(settings),
0036 rphicutPStable_(settings),
0037 rphicut2Stable_(settings),
0038 rcutPStable_(settings),
0039 rcut2Stable_(settings),
0040 alphainner_(settings),
0041 alphaouter_(settings),
0042 rSSinner_(settings),
0043 rSSouter_(settings) {
0044 phiregion_ = name[8] - 'A';
0045 layerdisk_ = initLayerDisk(3);
0046
0047 fullMatches_.resize(12, nullptr);
0048
0049
0050 icorrshift_ = 7;
0051
0052 if (layerdisk_ < N_PSLAYER) {
0053 icorzshift_ = -1 - settings_.PS_zderL_shift();
0054 } else {
0055 icorzshift_ = -1 - settings_.SS_zderL_shift();
0056 }
0057 phi0shift_ = 3;
0058 fact_ = 1;
0059 if (layerdisk_ >= N_PSLAYER && layerdisk_ < N_LAYER) {
0060 fact_ = (1 << (settings_.nzbitsstub(0) - settings_.nzbitsstub(5)));
0061 icorrshift_ -= (10 - settings_.nrbitsstub(layerdisk_));
0062 icorzshift_ += (settings_.nzbitsstub(0) - settings_.nzbitsstub(5) + settings_.nrbitsstub(layerdisk_) -
0063 settings_.nrbitsstub(0));
0064 phi0shift_ = 0;
0065 }
0066
0067 unsigned int region = getName()[8] - 'A';
0068 assert(region < settings_.nallstubs(layerdisk_));
0069
0070 if (layerdisk_ < N_LAYER) {
0071 phimatchcuttable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::barrelphi, region);
0072 zmatchcuttable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::barrelz, region);
0073 } else {
0074 rphicutPStable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::diskPSphi, region);
0075 rphicut2Stable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::disk2Sphi, region);
0076 rcutPStable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::diskPSr, region);
0077 rcut2Stable_.initmatchcut(layerdisk_, TrackletLUT::MatchType::disk2Sr, region);
0078 alphainner_.initmatchcut(layerdisk_, TrackletLUT::MatchType::alphainner, region);
0079 alphaouter_.initmatchcut(layerdisk_, TrackletLUT::MatchType::alphaouter, region);
0080 rSSinner_.initmatchcut(layerdisk_, TrackletLUT::MatchType::rSSinner, region);
0081 rSSouter_.initmatchcut(layerdisk_, TrackletLUT::MatchType::rSSouter, region);
0082 }
0083
0084 for (unsigned int i = 0; i < N_DSS_MOD * 2; i++) {
0085 ialphafactinner_[i] = (1 << settings_.alphashift()) * settings_.krprojshiftdisk() * settings_.half2SmoduleWidth() /
0086 (1 << (settings_.nbitsalpha() - 1)) / (settings_.rDSSinner(i) * settings_.rDSSinner(i)) /
0087 settings_.kphi();
0088 ialphafactouter_[i] = (1 << settings_.alphashift()) * settings_.krprojshiftdisk() * settings_.half2SmoduleWidth() /
0089 (1 << (settings_.nbitsalpha() - 1)) / (settings_.rDSSouter(i) * settings_.rDSSouter(i)) /
0090 settings_.kphi();
0091 }
0092 }
0093
0094 void MatchCalculator::addOutput(MemoryBase* memory, string output) {
0095 if (settings_.writetrace()) {
0096 edm::LogVerbatim("Tracklet") << "In " << name_ << " adding output to " << memory->getName() << " to output "
0097 << output;
0098 }
0099 if (output.substr(0, 8) == "matchout") {
0100 auto* tmp = dynamic_cast<FullMatchMemory*>(memory);
0101 assert(tmp != nullptr);
0102 unsigned int iSeed = getISeed(memory->getName());
0103 fullMatches_[iSeed] = tmp;
0104 return;
0105 }
0106 throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " could not find output " << output;
0107 }
0108
0109 void MatchCalculator::addInput(MemoryBase* memory, string input) {
0110 if (settings_.writetrace()) {
0111 edm::LogVerbatim("Tracklet") << "In " << name_ << " adding input from " << memory->getName() << " to input "
0112 << input;
0113 }
0114 if (input == "allstubin") {
0115 auto* tmp = dynamic_cast<AllStubsMemory*>(memory);
0116 assert(tmp != nullptr);
0117 allstubs_ = tmp;
0118 return;
0119 }
0120 if (input == "allprojin") {
0121 auto* tmp = dynamic_cast<AllProjectionsMemory*>(memory);
0122 assert(tmp != nullptr);
0123 allprojs_ = tmp;
0124 return;
0125 }
0126 if (input.substr(0, 5) == "match" && input.substr(input.size() - 2, 2) == "in") {
0127 auto* tmp = dynamic_cast<CandidateMatchMemory*>(memory);
0128 assert(tmp != nullptr);
0129 matches_.push_back(tmp);
0130 return;
0131 }
0132 throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " could not find input " << input;
0133 }
0134
0135 void MatchCalculator::execute(unsigned int iSector, double phioffset) {
0136 unsigned int countall = 0;
0137 unsigned int countsel = 0;
0138
0139
0140
0141 Tracklet* oldTracklet = nullptr;
0142
0143
0144 std::vector<std::pair<std::pair<Tracklet*, int>, const Stub*> > mergedMatches = mergeMatches(matches_);
0145
0146
0147
0148 unsigned int mergedepth = 3;
0149
0150 unsigned int maxProc = std::min(settings_.maxStep("MC") - mergedepth, (unsigned int)mergedMatches.size());
0151
0152
0153 int best_ideltaphi_barrel = 0xFFFF;
0154 int best_ideltaz_barrel = 0xFFFF;
0155 int best_ideltaphi_disk = 0xFFFF;
0156 int best_ideltar_disk = 0xFFFF;
0157 unsigned int curr_projid = -1;
0158 unsigned int next_projid = -1;
0159
0160 for (unsigned int j = 0; j < maxProc; j++) {
0161 if (settings_.debugTracklet() && j == 0) {
0162 edm::LogVerbatim("Tracklet") << getName() << " has " << mergedMatches.size() << " candidate matches";
0163 }
0164
0165 countall++;
0166
0167 const Stub* fpgastub = mergedMatches[j].second;
0168 Tracklet* tracklet = mergedMatches[j].first.first;
0169 const L1TStub* stub = fpgastub->l1tstub();
0170
0171
0172
0173 if (oldTracklet != nullptr) {
0174 assert(oldTracklet->TCID() <= tracklet->TCID());
0175 }
0176 oldTracklet = tracklet;
0177
0178 if (layerdisk_ < N_LAYER) {
0179
0180
0181 const Projection& proj = tracklet->proj(layerdisk_);
0182
0183 int ir = fpgastub->r().value();
0184 int iphi = proj.fpgaphiproj().value();
0185 int icorr = (ir * proj.fpgaphiprojder().value()) >> icorrshift_;
0186 iphi += icorr;
0187
0188 int iz = proj.fpgarzproj().value();
0189 int izcor = (ir * proj.fpgarzprojder().value() + (1 << (icorzshift_ - 1))) >> icorzshift_;
0190 iz += izcor;
0191
0192 int ideltaz = fpgastub->z().value() - iz;
0193 int ideltaphi = (fpgastub->phi().value() << phi0shift_) - (iphi << (settings_.phi0bitshift() - 1 + phi0shift_));
0194
0195
0196
0197 double phi = stub->phi() - phioffset;
0198 double r = stub->r();
0199 double z = stub->z();
0200
0201 if (settings_.useapprox()) {
0202 double dphi = reco::reducePhiRange(phi - fpgastub->phiapprox(0.0, 0.0));
0203 assert(std::abs(dphi) < 0.001);
0204 phi = fpgastub->phiapprox(0.0, 0.0);
0205 z = fpgastub->zapprox();
0206 r = fpgastub->rapprox();
0207 }
0208
0209 if (phi < 0)
0210 phi += 2 * M_PI;
0211
0212 double dr = r - settings_.rmean(layerdisk_);
0213 assert(std::abs(dr) < settings_.drmax());
0214
0215 double dphi = reco::reducePhiRange(phi - (proj.phiproj() + dr * proj.phiprojder()));
0216
0217 double dz = z - (proj.rzproj() + dr * proj.rzprojder());
0218
0219 double dphiapprox = reco::reducePhiRange(phi - (proj.phiprojapprox() + dr * proj.phiprojderapprox()));
0220
0221 double dzapprox = z - (proj.rzprojapprox() + dr * proj.rzprojderapprox());
0222
0223 int seedindex = tracklet->getISeed();
0224 unsigned int projindex = mergedMatches[j].first.second;
0225 curr_projid = next_projid;
0226 next_projid = projindex;
0227
0228
0229 bool newtracklet = (j == 0 || projindex != curr_projid);
0230 if (j == 0)
0231 best_ideltar_disk = (1 << (fpgastub->r().nbits() - 1));
0232
0233 if (newtracklet) {
0234 best_ideltaphi_barrel = (int)phimatchcuttable_.lookup(seedindex);
0235 best_ideltaz_barrel = (int)zmatchcuttable_.lookup(seedindex);
0236 }
0237
0238 assert(phimatchcuttable_.lookup(seedindex) > 0);
0239 assert(zmatchcuttable_.lookup(seedindex) > 0);
0240
0241 if (settings_.bookHistos()) {
0242 bool truthmatch = tracklet->stubtruthmatch(stub);
0243
0244 HistBase* hists = globals_->histograms();
0245 hists->FillLayerResidual(layerdisk_ + 1,
0246 seedindex,
0247 dphiapprox * settings_.rmean(layerdisk_),
0248 ideltaphi * settings_.kphi1() * settings_.rmean(layerdisk_),
0249 ideltaz * fact_ * settings_.kz(),
0250 dz,
0251 truthmatch);
0252 }
0253
0254
0255 if (std::abs(dphi) > 0.5 * settings_.dphisectorHG() || std::abs(dphiapprox) > 0.5 * settings_.dphisectorHG()) {
0256 throw cms::Exception("LogicError")
0257 << "WARNING dphi and/or dphiapprox too large : " << dphi << " " << dphiapprox << endl;
0258 }
0259
0260 if (settings_.writeMonitorData("Residuals")) {
0261 double pt = 0.01 * settings_.c() * settings_.bfield() / std::abs(tracklet->rinv());
0262
0263 globals_->ofstream("layerresiduals.txt")
0264 << layerdisk_ + 1 << " " << seedindex << " " << pt << " "
0265 << ideltaphi * settings_.kphi1() * settings_.rmean(layerdisk_) << " "
0266 << dphiapprox * settings_.rmean(layerdisk_) << " "
0267 << phimatchcuttable_.lookup(seedindex) * settings_.kphi1() * settings_.rmean(layerdisk_) << " "
0268 << ideltaz * fact_ * settings_.kz() << " " << dz << " "
0269 << zmatchcuttable_.lookup(seedindex) * settings_.kz() << endl;
0270 }
0271
0272
0273 bool imatch = (std::abs(ideltaphi) <= best_ideltaphi_barrel) && (ideltaz * fact_ < best_ideltaz_barrel) &&
0274 (ideltaz * fact_ >= -best_ideltaz_barrel);
0275
0276 if (imatch) {
0277 best_ideltaphi_barrel = std::abs(ideltaphi);
0278 best_ideltaz_barrel = std::abs(ideltaz * fact_);
0279 }
0280
0281 if (settings_.debugTracklet()) {
0282 edm::LogVerbatim("Tracklet") << getName() << " imatch = " << imatch << " ideltaphi cut " << ideltaphi << " "
0283 << phimatchcuttable_.lookup(seedindex) << " ideltaz*fact cut " << ideltaz * fact_
0284 << " " << zmatchcuttable_.lookup(seedindex);
0285 }
0286
0287 if (imatch) {
0288 countsel++;
0289
0290 tracklet->addMatch(layerdisk_,
0291 ideltaphi,
0292 ideltaz,
0293 dphi,
0294 dz,
0295 dphiapprox,
0296 dzapprox,
0297 (phiregion_ << 7) + fpgastub->stubindex().value(),
0298 mergedMatches[j].second);
0299
0300 if (settings_.debugTracklet()) {
0301 edm::LogVerbatim("Tracklet") << "Accepted full match in layer " << getName() << " " << tracklet;
0302 }
0303
0304 fullMatches_[seedindex]->addMatch(tracklet, mergedMatches[j].second);
0305 }
0306 } else {
0307
0308 assert(stub->z() * tracklet->t() > 0.0);
0309
0310 int sign = (tracklet->t() > 0.0) ? 1 : -1;
0311 int disk = sign * (layerdisk_ - (N_LAYER - 1));
0312 assert(disk != 0);
0313
0314
0315
0316 const Projection& proj = tracklet->proj(layerdisk_);
0317
0318 int iz = fpgastub->z().value();
0319 int iphi = proj.fpgaphiproj().value();
0320
0321
0322 int shifttmp = 6;
0323 int iphicorr = (iz * proj.fpgaphiprojder().value()) >> shifttmp;
0324
0325 iphi += iphicorr;
0326
0327 int ir = proj.fpgarzproj().value();
0328
0329
0330 int shifttmp2 = 7;
0331 int ircorr = (iz * proj.fpgarzprojder().value()) >> shifttmp2;
0332
0333 ir += ircorr;
0334
0335 int ideltaphi = fpgastub->phi().value() * settings_.kphi() / settings_.kphi() - iphi;
0336
0337 int irstub = fpgastub->r().value();
0338 int ialphafact = 0;
0339 if (!stub->isPSmodule()) {
0340 assert(irstub < (int)N_DSS_MOD * 2);
0341 if (abs(disk) <= 2) {
0342 ialphafact = ialphafactinner_[irstub];
0343 irstub = settings_.rDSSinner(irstub) / settings_.kr();
0344 } else {
0345 ialphafact = ialphafactouter_[irstub];
0346 irstub = settings_.rDSSouter(irstub) / settings_.kr();
0347 }
0348 }
0349
0350
0351 int ideltar = (irstub >> 1) - ir;
0352
0353 if (!stub->isPSmodule()) {
0354 int ialpha = fpgastub->alpha().value();
0355 int iphialphacor = ((ideltar * ialpha * ialphafact) >> settings_.alphashift());
0356 ideltaphi += iphialphacor;
0357 }
0358
0359
0360
0361 double phi = stub->phi() - phioffset;
0362 double z = stub->z();
0363 double r = stub->r();
0364
0365 if (settings_.useapprox()) {
0366 double dphi = reco::reducePhiRange(phi - fpgastub->phiapprox(0.0, 0.0));
0367 assert(std::abs(dphi) < 0.001);
0368 phi = fpgastub->phiapprox(0.0, 0.0);
0369 z = fpgastub->zapprox();
0370 r = fpgastub->rapprox();
0371 }
0372
0373 if (phi < 0)
0374 phi += 2 * M_PI;
0375
0376 double dz = z - sign * settings_.zmean(layerdisk_ - N_LAYER);
0377
0378 if (std::abs(dz) > settings_.dzmax()) {
0379 throw cms::Exception("LogicError") << __FILE__ << " " << __LINE__ << " " << name_ << " " << tracklet->getISeed()
0380 << "\n stub " << stub->z() << " disk " << disk << " " << dz;
0381 }
0382
0383 double phiproj = proj.phiproj() + dz * proj.phiprojder();
0384
0385 double rproj = proj.rzproj() + dz * proj.rzprojder();
0386
0387 double deltar = r - rproj;
0388
0389 double dr = stub->r() - rproj;
0390
0391 double dphi = reco::reducePhiRange(phi - phiproj);
0392
0393 double dphiapprox = reco::reducePhiRange(phi - (proj.phiprojapprox() + dz * proj.phiprojderapprox()));
0394
0395 double drapprox = stub->r() - (proj.rzprojapprox() + dz * proj.rzprojderapprox());
0396
0397 double drphi = dphi * stub->r();
0398 double drphiapprox = dphiapprox * stub->r();
0399
0400 if (!stub->isPSmodule()) {
0401 double alphanorm = stub->alphanorm();
0402 dphi += dr * alphanorm * settings_.half2SmoduleWidth() / stub->r2();
0403 dphiapprox += drapprox * alphanorm * settings_.half2SmoduleWidth() / stub->r2();
0404
0405 drphi += dr * alphanorm * settings_.half2SmoduleWidth() / stub->r();
0406 drphiapprox += dr * alphanorm * settings_.half2SmoduleWidth() / stub->r();
0407 }
0408
0409 int seedindex = tracklet->getISeed();
0410
0411 int idrphicut = rphicutPStable_.lookup(seedindex);
0412 int idrcut = rcutPStable_.lookup(seedindex);
0413 if (!stub->isPSmodule()) {
0414 idrphicut = rphicut2Stable_.lookup(seedindex);
0415 idrcut = rcut2Stable_.lookup(seedindex);
0416 }
0417
0418 unsigned int projindex = mergedMatches[j].first.second;
0419 curr_projid = next_projid;
0420 next_projid = projindex;
0421
0422 bool newtracklet = (j == 0 || projindex != curr_projid);
0423
0424 if (newtracklet) {
0425 best_ideltaphi_disk = idrphicut;
0426 best_ideltar_disk = idrcut;
0427 }
0428
0429 double drphicut = idrphicut * settings_.kphi() * settings_.kr();
0430 double drcut = idrcut * settings_.krprojshiftdisk();
0431
0432 bool match, imatch;
0433 if (std::abs(dphi) < third * settings_.dphisectorHG() &&
0434 std::abs(dphiapprox) < third * settings_.dphisectorHG()) {
0435 if (settings_.writeMonitorData("Residuals")) {
0436 double pt = 0.01 * settings_.c() * settings_.bfield() / std::abs(tracklet->rinv());
0437
0438 globals_->ofstream("diskresiduals.txt")
0439 << disk << " " << stub->isPSmodule() << " " << tracklet->layer() << " " << abs(tracklet->disk()) << " "
0440 << pt << " " << ideltaphi * settings_.kphi() * stub->r() << " " << drphiapprox << " " << drphicut << " "
0441 << ideltar * settings_.krprojshiftdisk() << " " << deltar << " " << drcut << " " << endl;
0442 }
0443
0444
0445 match = (std::abs(drphi) < drphicut) && (std::abs(deltar) < drcut);
0446
0447
0448 imatch = (std::abs(ideltaphi) * irstub < best_ideltaphi_disk) && (std::abs(ideltar) < best_ideltar_disk);
0449
0450 if (imatch) {
0451 best_ideltaphi_disk = std::abs(ideltaphi) * irstub;
0452 best_ideltar_disk = std::abs(ideltar);
0453 }
0454 } else {
0455 edm::LogProblem("Tracklet") << "WARNING dphi and/or dphiapprox too large : " << dphi << " " << dphiapprox
0456 << "dphi " << dphi << " Seed / ISeed " << tracklet->getISeed() << endl;
0457 match = false;
0458 imatch = false;
0459 }
0460 if (settings_.debugTracklet()) {
0461 edm::LogVerbatim("Tracklet") << "imatch match disk: " << imatch << " " << match << " " << std::abs(ideltaphi)
0462 << " " << drphicut / (settings_.kphi() * stub->r()) << " " << std::abs(ideltar)
0463 << " " << drcut / settings_.krprojshiftdisk() << " r = " << stub->r();
0464 }
0465
0466 if (imatch) {
0467 countsel++;
0468
0469 if (settings_.debugTracklet()) {
0470 edm::LogVerbatim("Tracklet") << "MatchCalculator found match in disk " << getName();
0471 }
0472
0473 tracklet->addMatch(layerdisk_,
0474 ideltaphi,
0475 ideltar,
0476 drphi / stub->r(),
0477 dr,
0478 drphiapprox / stub->r(),
0479 drapprox,
0480 (phiregion_ << 7) + fpgastub->stubindex().value(),
0481 fpgastub);
0482
0483 if (settings_.debugTracklet()) {
0484 edm::LogVerbatim("Tracklet") << "Accepted full match in disk " << getName() << " " << tracklet;
0485 }
0486
0487 fullMatches_[seedindex]->addMatch(tracklet, mergedMatches[j].second);
0488 }
0489 }
0490 if (countall >= settings_.maxStep("MC"))
0491 break;
0492 }
0493
0494 if (settings_.writeMonitorData("MC")) {
0495 globals_->ofstream("matchcalculator.txt") << getName() << " " << countall << " " << countsel << endl;
0496 }
0497 }
0498
0499
0500 std::vector<std::pair<std::pair<Tracklet*, int>, const Stub*> > MatchCalculator::mergeMatches(
0501 vector<CandidateMatchMemory*>& candmatch) {
0502 std::vector<std::pair<std::pair<Tracklet*, int>, const Stub*> > tmp;
0503
0504 std::vector<unsigned int> indexArray;
0505 indexArray.reserve(candmatch.size());
0506 for (unsigned int i = 0; i < candmatch.size(); i++) {
0507 indexArray.push_back(0);
0508 }
0509
0510 int bestIndex = -1;
0511 do {
0512 int bestSector = 100;
0513 int bestTCID = -1;
0514 bestIndex = -1;
0515 for (unsigned int i = 0; i < candmatch.size(); i++) {
0516 if (indexArray[i] >= candmatch[i]->nMatches()) {
0517
0518 continue;
0519 }
0520 int TCID = candmatch[i]->getMatch(indexArray[i]).first.first->TCID();
0521 int dSector = 0;
0522 if (dSector > 2)
0523 dSector -= N_SECTOR;
0524 if (dSector < -2)
0525 dSector += N_SECTOR;
0526 assert(abs(dSector) < 2);
0527 if (dSector == -1)
0528 dSector = 2;
0529 if (dSector < bestSector) {
0530 bestSector = dSector;
0531 bestTCID = TCID;
0532 bestIndex = i;
0533 }
0534 if (dSector == bestSector) {
0535 if (TCID < bestTCID || bestTCID < 0) {
0536 bestTCID = TCID;
0537 bestIndex = i;
0538 }
0539 }
0540 }
0541 if (bestIndex != -1) {
0542 tmp.push_back(candmatch[bestIndex]->getMatch(indexArray[bestIndex]));
0543 indexArray[bestIndex]++;
0544 }
0545 } while (bestIndex != -1);
0546
0547 if (layerdisk_ < N_LAYER) {
0548 int lastTCID = -1;
0549 bool error = false;
0550
0551
0552 for (unsigned int i = 1; i < tmp.size(); i++) {
0553 if (lastTCID > tmp[i].first.first->TCID()) {
0554 edm::LogProblem("Tracklet") << "Wrong TCID ordering for projections in " << getName() << " last " << lastTCID
0555 << " " << tmp[i].first.first->TCID();
0556 error = true;
0557 } else {
0558 lastTCID = tmp[i].first.first->TCID();
0559 }
0560 }
0561
0562 if (error) {
0563 for (unsigned int i = 1; i < tmp.size(); i++) {
0564 edm::LogProblem("Tracklet") << "Wrong order for in " << getName() << " " << i << " " << tmp[i].first.first
0565 << " " << tmp[i].first.first->TCID();
0566 }
0567 }
0568 }
0569
0570 for (unsigned int i = 0; i < tmp.size(); i++) {
0571 if (i > 0) {
0572
0573
0574
0575
0576 assert(tmp[i - 1].first.first->TCID() <= tmp[i].first.first->TCID());
0577 }
0578 }
0579
0580 return tmp;
0581 }