Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-10-14 01:44:05

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