Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-07 22:33:34

0001 #include "L1Trigger/TrackFindingTracklet/interface/Sector.h"
0002 #include "L1Trigger/TrackFindingTracklet/interface/MemoryBase.h"
0003 #include "L1Trigger/TrackFindingTracklet/interface/Settings.h"
0004 #include "L1Trigger/TrackFindingTracklet/interface/Globals.h"
0005 
0006 #include "L1Trigger/TrackFindingTracklet/interface/DTCLinkMemory.h"
0007 #include "L1Trigger/TrackFindingTracklet/interface/InputLinkMemory.h"
0008 #include "L1Trigger/TrackFindingTracklet/interface/AllStubsMemory.h"
0009 #include "L1Trigger/TrackFindingTracklet/interface/AllInnerStubsMemory.h"
0010 #include "L1Trigger/TrackFindingTracklet/interface/VMStubsTEMemory.h"
0011 #include "L1Trigger/TrackFindingTracklet/interface/VMStubsMEMemory.h"
0012 #include "L1Trigger/TrackFindingTracklet/interface/StubPairsMemory.h"
0013 #include "L1Trigger/TrackFindingTracklet/interface/StubTripletsMemory.h"
0014 #include "L1Trigger/TrackFindingTracklet/interface/TrackletParametersMemory.h"
0015 #include "L1Trigger/TrackFindingTracklet/interface/TrackletProjectionsMemory.h"
0016 #include "L1Trigger/TrackFindingTracklet/interface/AllProjectionsMemory.h"
0017 #include "L1Trigger/TrackFindingTracklet/interface/VMProjectionsMemory.h"
0018 #include "L1Trigger/TrackFindingTracklet/interface/CandidateMatchMemory.h"
0019 #include "L1Trigger/TrackFindingTracklet/interface/FullMatchMemory.h"
0020 #include "L1Trigger/TrackFindingTracklet/interface/TrackFitMemory.h"
0021 #include "L1Trigger/TrackFindingTracklet/interface/CleanTrackMemory.h"
0022 
0023 #include "L1Trigger/TrackFindingTracklet/interface/InputRouter.h"
0024 #include "L1Trigger/TrackFindingTracklet/interface/VMRouterCM.h"
0025 #include "L1Trigger/TrackFindingTracklet/interface/VMRouter.h"
0026 #include "L1Trigger/TrackFindingTracklet/interface/TrackletEngine.h"
0027 #include "L1Trigger/TrackFindingTracklet/interface/TrackletEngineDisplaced.h"
0028 #include "L1Trigger/TrackFindingTracklet/interface/TripletEngine.h"
0029 #include "L1Trigger/TrackFindingTracklet/interface/TrackletCalculator.h"
0030 #include "L1Trigger/TrackFindingTracklet/interface/TrackletProcessor.h"
0031 #include "L1Trigger/TrackFindingTracklet/interface/TrackletCalculatorDisplaced.h"
0032 #include "L1Trigger/TrackFindingTracklet/interface/ProjectionRouter.h"
0033 #include "L1Trigger/TrackFindingTracklet/interface/MatchEngine.h"
0034 #include "L1Trigger/TrackFindingTracklet/interface/MatchCalculator.h"
0035 #include "L1Trigger/TrackFindingTracklet/interface/MatchProcessor.h"
0036 #include "L1Trigger/TrackFindingTracklet/interface/FitTrack.h"
0037 #include "L1Trigger/TrackFindingTracklet/interface/PurgeDuplicate.h"
0038 #include "L1Trigger/TrackFindingTracklet/interface/Util.h"
0039 
0040 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0041 #include "FWCore/Utilities/interface/Exception.h"
0042 #include "DataFormats/Math/interface/deltaPhi.h"
0043 #include "L1Trigger/TrackFindingTracklet/interface/TrackletLUT.h"
0044 
0045 using namespace std;
0046 using namespace trklet;
0047 
0048 Sector::Sector(Settings const& settings, Globals* globals) : isector_(-1), settings_(settings), globals_(globals) {}
0049 
0050 Sector::~Sector() = default;
0051 
0052 void Sector::setSector(unsigned int isector) {
0053   assert(isector < N_SECTOR);
0054   isector_ = isector;
0055   double dphi = 2 * M_PI / N_SECTOR;
0056   double dphiHG = 0.5 * settings_.dphisectorHG() - M_PI / N_SECTOR;
0057   phimin_ = isector_ * dphi - dphiHG;
0058   phimax_ = phimin_ + dphi + 2 * dphiHG;
0059   phimin_ -= M_PI / N_SECTOR;
0060   phimax_ -= M_PI / N_SECTOR;
0061   phimin_ = reco::reduceRange(phimin_);
0062   phimax_ = reco::reduceRange(phimax_);
0063   if (phimin_ > phimax_) {
0064     phimin_ -= 2 * M_PI;
0065   }
0066 }
0067 
0068 bool Sector::addStub(L1TStub stub, string dtc) {
0069   unsigned int layerdisk = stub.layerdisk();
0070 
0071   if (layerdisk < N_LAYER && globals_->phiCorr(layerdisk) == nullptr) {
0072     globals_->phiCorr(layerdisk) = new TrackletLUT(settings_);
0073     globals_->phiCorr(layerdisk)->initPhiCorrTable(layerdisk, 3);
0074   }
0075 
0076   Stub fpgastub(stub, settings_, *globals_);
0077 
0078   if (layerdisk < N_LAYER) {
0079     FPGAWord r = fpgastub.r();
0080     int bendbin = fpgastub.bend().value();
0081     int rbin = (r.value() + (1 << (r.nbits() - 1))) >> (r.nbits() - 3);
0082     const TrackletLUT& phiCorrTable = *globals_->phiCorr(layerdisk);
0083     int iphicorr = phiCorrTable.lookup(bendbin * (1 << 3) + rbin);
0084     fpgastub.setPhiCorr(iphicorr);
0085   }
0086 
0087   int nadd = 0;
0088   for (unsigned int i = 0; i < DL_.size(); i++) {
0089     const string& name = DL_[i]->getName();
0090     if (name.find("_" + dtc) == string::npos)
0091       continue;
0092     DL_[i]->addStub(stub, fpgastub);
0093     nadd++;
0094   }
0095 
0096   assert(nadd == 1);
0097 
0098   return true;
0099 }
0100 
0101 void Sector::addMem(string memType, string memName) {
0102   if (memType == "DTCLink:") {
0103     addMemToVec(DL_, memName, settings_, phimin_, phimax_);
0104   } else if (memType == "InputLink:") {
0105     addMemToVec(IL_, memName, settings_, phimin_, phimax_);
0106   } else if (memType == "AllStubs:") {
0107     addMemToVec(AS_, memName, settings_);
0108   } else if (memType == "AllInnerStubs:") {
0109     addMemToVec(AIS_, memName, settings_);
0110   } else if (memType == "VMStubsTE:") {
0111     addMemToVec(VMSTE_, memName, settings_);
0112   } else if (memType == "VMStubsME:") {
0113     addMemToVec(VMSME_, memName, settings_);
0114   } else if (memType == "StubPairs:" || memType == "StubPairsDisplaced:") {
0115     addMemToVec(SP_, memName, settings_);
0116   } else if (memType == "StubTriplets:") {
0117     addMemToVec(ST_, memName, settings_);
0118   } else if (memType == "TrackletParameters:") {
0119     addMemToVec(TPAR_, memName, settings_);
0120   } else if (memType == "TrackletProjections:") {
0121     addMemToVec(TPROJ_, memName, settings_);
0122   } else if (memType == "AllProj:") {
0123     addMemToVec(AP_, memName, settings_);
0124   } else if (memType == "VMProjections:") {
0125     addMemToVec(VMPROJ_, memName, settings_);
0126   } else if (memType == "CandidateMatch:") {
0127     addMemToVec(CM_, memName, settings_);
0128   } else if (memType == "FullMatch:") {
0129     addMemToVec(FM_, memName, settings_);
0130   } else if (memType == "TrackFit:") {
0131     addMemToVec(TF_, memName, settings_, phimin_, phimax_);
0132   } else if (memType == "CleanTrack:") {
0133     addMemToVec(CT_, memName, settings_, phimin_, phimax_);
0134   } else {
0135     edm::LogPrint("Tracklet") << "Don't know of memory type: " << memType;
0136     exit(0);
0137   }
0138 }
0139 
0140 void Sector::addProc(string procType, string procName) {
0141   if (procType == "InputRouter:") {
0142     addProcToVec(IR_, procName, settings_, globals_);
0143   } else if (procType == "VMRouter:") {
0144     addProcToVec(VMR_, procName, settings_, globals_);
0145   } else if (procType == "VMRouterCM:") {
0146     addProcToVec(VMRCM_, procName, settings_, globals_);
0147   } else if (procType == "TrackletEngine:") {
0148     addProcToVec(TE_, procName, settings_, globals_);
0149   } else if (procType == "TrackletEngineDisplaced:") {
0150     addProcToVec(TED_, procName, settings_, globals_);
0151   } else if (procType == "TripletEngine:") {
0152     addProcToVec(TRE_, procName, settings_, globals_);
0153   } else if (procType == "TrackletCalculator:") {
0154     addProcToVec(TC_, procName, settings_, globals_);
0155   } else if (procType == "TrackletProcessor:") {
0156     addProcToVec(TP_, procName, settings_, globals_);
0157   } else if (procType == "TrackletCalculatorDisplaced:") {
0158     addProcToVec(TCD_, procName, settings_, globals_);
0159   } else if (procType == "ProjectionRouter:") {
0160     addProcToVec(PR_, procName, settings_, globals_);
0161   } else if (procType == "MatchEngine:") {
0162     addProcToVec(ME_, procName, settings_, globals_);
0163   } else if (procType == "MatchCalculator:" ||
0164              procType == "DiskMatchCalculator:") {  //TODO should not be used in configurations
0165     addProcToVec(MC_, procName, settings_, globals_);
0166   } else if (procType == "MatchProcessor:") {
0167     addProcToVec(MP_, procName, settings_, globals_);
0168   } else if (procType == "FitTrack:") {
0169     addProcToVec(FT_, procName, settings_, globals_);
0170   } else if (procType == "PurgeDuplicate:") {
0171     addProcToVec(PD_, procName, settings_, globals_);
0172   } else {
0173     edm::LogPrint("Tracklet") << "Don't know of processing type: " << procType;
0174     exit(0);
0175   }
0176 }
0177 
0178 void Sector::addWire(string mem, string procinfull, string procoutfull) {
0179   stringstream ss1(procinfull);
0180   string procin, output;
0181   getline(ss1, procin, '.');
0182   getline(ss1, output);
0183 
0184   stringstream ss2(procoutfull);
0185   string procout, input;
0186   getline(ss2, procout, '.');
0187   getline(ss2, input);
0188 
0189   MemoryBase* memory = getMem(mem);
0190 
0191   if (!procin.empty()) {
0192     ProcessBase* inProc = getProc(procin);
0193     inProc->addOutput(memory, output);
0194   }
0195 
0196   if (!procout.empty()) {
0197     ProcessBase* outProc = getProc(procout);
0198     outProc->addInput(memory, input);
0199   }
0200 }
0201 
0202 ProcessBase* Sector::getProc(string procName) {
0203   auto it = Processes_.find(procName);
0204 
0205   if (it != Processes_.end()) {
0206     return it->second;
0207   }
0208   throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " Could not find process : " << procName << endl;
0209   return nullptr;
0210 }
0211 
0212 MemoryBase* Sector::getMem(string memName) {
0213   auto it = Memories_.find(memName);
0214 
0215   if (it != Memories_.end()) {
0216     return it->second;
0217   }
0218   throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " Could not find memory : " << memName;
0219   return nullptr;
0220 }
0221 
0222 void Sector::writeDTCStubs(bool first) {
0223   for (auto& i : DL_) {
0224     i->writeStubs(first, isector_);
0225   }
0226 }
0227 
0228 void Sector::writeIRStubs(bool first) {
0229   for (auto& i : IL_) {
0230     i->writeStubs(first, isector_);
0231   }
0232 }
0233 
0234 void Sector::writeVMSTE(bool first) {
0235   for (auto& i : VMSTE_) {
0236     i->writeStubs(first, isector_);
0237   }
0238 }
0239 
0240 void Sector::writeVMSME(bool first) {
0241   for (auto& i : VMSME_) {
0242     i->writeStubs(first, isector_);
0243   }
0244 }
0245 
0246 void Sector::writeAS(bool first) {
0247   for (auto& i : AS_) {
0248     i->writeStubs(first, isector_);
0249   }
0250 }
0251 
0252 void Sector::writeAIS(bool first) {
0253   for (auto& i : AIS_) {
0254     i->writeStubs(first, isector_);
0255   }
0256 }
0257 
0258 void Sector::writeSP(bool first) {
0259   for (auto& i : SP_) {
0260     i->writeSP(first, isector_);
0261   }
0262 }
0263 
0264 void Sector::writeST(bool first) {
0265   for (auto& i : ST_) {
0266     i->writeST(first, isector_);
0267   }
0268 }
0269 
0270 void Sector::writeTPAR(bool first) {
0271   for (auto& i : TPAR_) {
0272     i->writeTPAR(first, isector_);
0273   }
0274 }
0275 
0276 void Sector::writeTPROJ(bool first) {
0277   for (auto& i : TPROJ_) {
0278     i->writeTPROJ(first, isector_);
0279   }
0280 }
0281 
0282 void Sector::writeAP(bool first) {
0283   for (auto& i : AP_) {
0284     i->writeAP(first, isector_);
0285   }
0286 }
0287 
0288 void Sector::writeVMPROJ(bool first) {
0289   for (auto& i : VMPROJ_) {
0290     i->writeVMPROJ(first, isector_);
0291   }
0292 }
0293 
0294 void Sector::writeCM(bool first) {
0295   for (auto& i : CM_) {
0296     i->writeCM(first, isector_);
0297   }
0298 }
0299 
0300 void Sector::writeMC(bool first) {
0301   for (auto& i : FM_) {
0302     i->writeMC(first, isector_);
0303   }
0304 }
0305 
0306 void Sector::writeTF(bool first) {
0307   for (auto& i : TF_) {
0308     i->writeTF(first, isector_);
0309   }
0310 }
0311 
0312 void Sector::writeCT(bool first) {
0313   for (auto& i : CT_) {
0314     i->writeCT(first, isector_);
0315   }
0316 }
0317 
0318 void Sector::clean() {
0319   for (auto& mem : MemoriesV_) {
0320     mem->clean();
0321   }
0322 }
0323 
0324 void Sector::executeIR() {
0325   for (auto& i : IR_) {
0326     i->execute();
0327   }
0328 }
0329 
0330 void Sector::executeVMR() {
0331   if (settings_.writeMonitorData("IL")) {
0332     ofstream& out = globals_->ofstream("inputlink.txt");
0333     for (auto& i : IL_) {
0334       out << i->getName() << " " << i->nStubs() << endl;
0335     }
0336   }
0337   for (auto& i : VMR_) {
0338     i->execute();
0339   }
0340   for (auto& i : VMRCM_) {
0341     i->execute();
0342   }
0343 }
0344 
0345 void Sector::executeTE() {
0346   for (auto& i : TE_) {
0347     i->execute();
0348   }
0349 }
0350 
0351 void Sector::executeTED() {
0352   for (auto& i : TED_) {
0353     i->execute();
0354   }
0355 }
0356 
0357 void Sector::executeTRE() {
0358   for (auto& i : TRE_) {
0359     i->execute();
0360   }
0361 }
0362 
0363 void Sector::executeTP() {
0364   for (auto& i : TP_) {
0365     i->execute(isector_, phimin_, phimax_);
0366   }
0367 }
0368 
0369 void Sector::executeTC() {
0370   for (auto& i : TC_) {
0371     i->execute(isector_, phimin_, phimax_);
0372   }
0373 
0374   if (settings_.writeMonitorData("TrackProjOcc")) {
0375     ofstream& out = globals_->ofstream("trackprojocc.txt");
0376     for (auto& i : TPROJ_) {
0377       out << i->getName() << " " << i->nTracklets() << endl;
0378     }
0379   }
0380 }
0381 
0382 void Sector::executeTCD() {
0383   for (auto& i : TCD_) {
0384     i->execute(isector_, phimin_, phimax_);
0385   }
0386 }
0387 
0388 void Sector::executePR() {
0389   for (auto& i : PR_) {
0390     i->execute();
0391   }
0392 }
0393 
0394 void Sector::executeME() {
0395   for (auto& i : ME_) {
0396     i->execute();
0397   }
0398 }
0399 
0400 void Sector::executeMC() {
0401   for (auto& i : MC_) {
0402     i->execute(phimin_);
0403   }
0404 }
0405 
0406 void Sector::executeMP() {
0407   for (auto& i : MP_) {
0408     i->execute(isector_, phimin_);
0409   }
0410 }
0411 
0412 void Sector::executeFT() {
0413   for (auto& i : FT_) {
0414     i->execute(isector_);
0415   }
0416 }
0417 
0418 void Sector::executePD(std::vector<Track>& tracks) {
0419   for (auto& i : PD_) {
0420     i->execute(tracks, isector_);
0421   }
0422 }
0423 
0424 std::vector<Tracklet*> Sector::getAllTracklets() const {
0425   std::vector<Tracklet*> tmp;
0426   for (auto& tpar : TPAR_) {
0427     for (unsigned int j = 0; j < tpar->nTracklets(); j++) {
0428       tmp.push_back(tpar->getTracklet(j));
0429     }
0430   }
0431   return tmp;
0432 }
0433 
0434 std::vector<const Stub*> Sector::getStubs() const {
0435   std::vector<const Stub*> tmp;
0436 
0437   for (auto& imem : IL_) {
0438     for (unsigned int istub = 0; istub < imem->nStubs(); istub++) {
0439       tmp.push_back(imem->getStub(istub));
0440     }
0441   }
0442 
0443   return tmp;
0444 }
0445 
0446 std::unordered_set<int> Sector::seedMatch(int itp) const {
0447   std::unordered_set<int> tmpSeeds;
0448   for (auto& i : TPAR_) {
0449     unsigned int nTracklet = i->nTracklets();
0450     for (unsigned int j = 0; j < nTracklet; j++) {
0451       if (i->getTracklet(j)->tpseed() == itp) {
0452         tmpSeeds.insert(i->getTracklet(j)->getISeed());
0453       }
0454     }
0455   }
0456   return tmpSeeds;
0457 }