Back to home page

Project CMSSW displayed by LXR

 
 

    


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

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