TrackFindingTracklet/src/TrackletEventProcessor.cc

0001 #include "L1Trigger/TrackFindingTracklet/interface/TrackletEventProcessor.h"
0002 #include "L1Trigger/TrackFindingTracklet/interface/SLHCEvent.h"
0003 #include "L1Trigger/TrackFindingTracklet/interface/Globals.h"
0004 #include "L1Trigger/TrackFindingTracklet/interface/SLHCEvent.h"
0005 #include "L1Trigger/TrackFindingTracklet/interface/Sector.h"
0006 #include "L1Trigger/TrackFindingTracklet/interface/HistBase.h"
0007 #include "L1Trigger/TrackFindingTracklet/interface/Track.h"
0008 #include "L1Trigger/TrackFindingTracklet/interface/TrackletConfigBuilder.h"
0009 #include "L1Trigger/TrackFindingTracklet/interface/IMATH_TrackletCalculator.h"
0010 #include "L1Trigger/TrackFindingTracklet/interface/StubStreamData.h"
0011
0012 #include "DataFormats/Math/interface/deltaPhi.h"
0013
0014 #include <iomanip>
0015 #include <filesystem>
0016
0017 using namespace trklet;
0018 using namespace std;
0019
0020 TrackletEventProcessor::TrackletEventProcessor() : settings_(nullptr) {}
0021
0022 TrackletEventProcessor::~TrackletEventProcessor() {
0023   if (settings_ && settings_->bookHistos()) {
0024     histbase_->close();
0025   }
0026 }
0027
0028 void TrackletEventProcessor::init(Settings const& theSettings, const tt::Setup* setup) {
0029   settings_ = &theSettings;
0030   globals_ = make_unique<Globals>(*settings_);
0031
0032   //Verify consistency
0033   if (settings_->kphi0pars() != globals_->ITC_L1L2()->phi0_final.K()) {
0034     throw cms::Exception("Inconsistency") << "phi0 conversion parameter inconsistency\n";
0035   }
0036
0037   if (settings_->krinvpars() != globals_->ITC_L1L2()->rinv_final.K()) {
0038     throw cms::Exception("Inconsistency") << "ring conversion parameter inconsistency\n";
0039   }
0040
0041   if (settings_->ktpars() != globals_->ITC_L1L2()->t_final.K()) {
0042     throw cms::Exception("Inconsistency") << "t conversion parameter inconsistency\n";
0043   }
0044
0045   if (settings_->kphider() != globals_->ITC_L1L2()->der_phiL_final.K()) {
0046     throw cms::Exception("Inconsistency")
0047         << "t conversion parameter inconsistency:" << settings_->kphider() / globals_->ITC_L1L2()->der_phiL_final.K()
0048         << "\n";
0049   }
0050
0051   if (settings_->debugTracklet()) {
0052     edm::LogVerbatim("Tracklet") << "========================================================= \n"
0053                                  << "Conversion factors for global coordinates: \n"
0054                                  << "z    kz            = " << settings_->kz() << "\n"
0055                                  << "r    kr            = " << settings_->kr() << "\n"
0056                                  << "phi  kphi1         = " << settings_->kphi1() << "\n"
0057                                  << "========================================================= \n"
0058                                  << "Conversion factors for track(let) parameters: \n"
0059                                  << "rinv krinvpars     = " << settings_->krinvpars() << "\n"
0060                                  << "phi0 kphi0pars     = " << settings_->kphi0pars() << "\n"
0061                                  << "d0   kd0pars       = " << settings_->kd0pars() << "\n"
0062                                  << "t    ktpars        = " << settings_->ktpars() << "\n"
0063                                  << "z0   kz0pars       = " << settings_->kz0pars() << "\n"
0064                                  << "========================================================= \n"
0065                                  << "phi0bitshift = " << settings_->phi0bitshift() << "\n"
0066                                  << "d0bitshift   = ??? \n"
0067                                  << "=========================================================";
0068   }
0069
0070   if (settings_->bookHistos()) {
0071     histbase_ = new HistBase;
0072     histbase_->open();
0073     histbase_->bookLayerResidual();
0074     histbase_->bookDiskResidual();
0075     histbase_->bookTrackletParams();
0076     histbase_->bookSeedEff();
0077
0078     globals_->histograms() = histbase_;
0079   }
0080
0081   sector_ = make_unique<Sector>(*settings_, globals_.get());
0082
0083   if (settings_->extended() || settings_->reduced()) {
0084     ifstream inmem(settings_->memoryModulesFile().c_str());
0085     assert(inmem.good());
0086
0087     ifstream inproc(settings_->processingModulesFile().c_str());
0088     assert(inproc.good());
0089
0090     ifstream inwire(settings_->wiresFile().c_str());
0091     assert(inwire.good());
0092
0093     configure(inwire, inmem, inproc);
0094
0095   } else {
0096     TrackletConfigBuilder config(*settings_, setup);
0097
0098     //Write configurations to file.
0099     if (settings_->writeConfig()) {
0100       std::ofstream wires = openfile(settings_->tablePath(), "wires.dat", __FILE__, __LINE__);
0101       std::ofstream memorymodules = openfile(settings_->tablePath(), "memorymodules.dat", __FILE__, __LINE__);
0102       std::ofstream processingmodules = openfile(settings_->tablePath(), "processingmodules.dat", __FILE__, __LINE__);
0103
0104       config.writeAll(wires, memorymodules, processingmodules);
0105     }
0106
0107     std::stringstream wires;
0108     std::stringstream memorymodules;
0109     std::stringstream processingmodules;
0110
0111     config.writeAll(wires, memorymodules, processingmodules);
0112     configure(wires, memorymodules, processingmodules);
0113   }
0114 }
0115
0116 void TrackletEventProcessor::configure(istream& inwire, istream& inmem, istream& inproc) {
0117   // get the memory modules
0118   if (settings_->debugTracklet()) {
0119     edm::LogVerbatim("Tracklet") << "Will read memory modules";
0120   }
0121
0122   while (inmem.good()) {
0123     string memType, memName, size;
0124     inmem >> memType >> memName >> size;
0125     if (!inmem.good())
0126       continue;
0127     if (settings_->writetrace()) {
0128       edm::LogVerbatim("Tracklet") << "Read memory: " << memType << " " << memName;
0129     }
0130     sector_->addMem(memType, memName);
0131   }
0132
0133   // get the processing modules
0134   if (settings_->debugTracklet()) {
0135     edm::LogVerbatim("Tracklet") << "Will read processing modules";
0136   }
0137
0138   while (inproc.good()) {
0139     string procType, procName;
0140     inproc >> procType >> procName;
0141     if (!inproc.good())
0142       continue;
0143     if (settings_->writetrace()) {
0144       edm::LogVerbatim("Tracklet") << "Read process: " << procType << " " << procName;
0145     }
0146     sector_->addProc(procType, procName);
0147   }
0148
0149   // get the wiring information
0150   if (settings_->debugTracklet()) {
0151     edm::LogVerbatim("Tracklet") << "Will read wiring information";
0152   }
0153
0154   while (inwire.good()) {
0155     string line;
0156     getline(inwire, line);
0157     if (!inwire.good())
0158       continue;
0159     if (settings_->writetrace()) {
0160       edm::LogVerbatim("Tracklet") << "Line : " << line;
0161     }
0162     stringstream ss(line);
0163     string mem, tmp1, procin, tmp2, procout;
0164     ss >> mem >> tmp1 >> procin;
0165     if (procin == "output=>") {
0166       procin = "";
0167       ss >> procout;
0168     } else {
0169       ss >> tmp2 >> procout;
0170     }
0171
0172     sector_->addWire(mem, procin, procout);
0173   }
0174 }
0175
0176 void TrackletEventProcessor::event(SLHCEvent& ev,
0177                                    vector<vector<string>>& streamsTrackRaw,
0178                                    vector<vector<StubStreamData>>& streamsStubRaw) {
0179   globals_->event() = &ev;
0180
0181   tracks_.clear();
0182   eventnum_++;
0183   bool first = (eventnum_ == 1);
0184
0185   for (unsigned int k = 0; k < N_SECTOR; k++) {
0186     sector_->setSector(k);
0187
0188     cleanTimer_.start();
0189     sector_->clean();
0190     cleanTimer_.stop();
0191
0192     addStubTimer_.start();
0193
0194     vector<int> layerstubs(N_LAYER + N_DISK, 0);
0195     vector<int> layerstubssector(N_SECTOR * (N_LAYER + N_DISK), 0);
0196
0197     for (int j = 0; j < ev.nstubs(); j++) {
0198       const L1TStub& stub = ev.stub(j);
0199       unsigned int isector = stub.region();
0200       if (isector != k) {
0201         continue;
0202       }
0203
0204       const string& dtc = stub.DTClink();
0205
0206       layerstubs[stub.layerdisk()]++;
0207       layerstubssector[isector * (N_LAYER + N_DISK) + stub.layerdisk()]++;
0208
0209       sector_->addStub(stub, dtc);
0210     }
0211
0212     if (settings_->writeMonitorData("StubsLayerSector")) {
0213       for (unsigned int index = 0; index < layerstubssector.size(); index++) {
0214         int layerdisk = index % (N_LAYER + N_DISK);
0215         int sector = index / (N_LAYER + N_DISK);
0216         globals_->ofstream("stubslayersector.txt")
0217             << layerdisk << " " << sector << " " << layerstubssector[index] << endl;
0218       }
0219     }
0220
0221     if (settings_->writeMonitorData("StubsLayer")) {
0222       for (unsigned int layerdisk = 0; layerdisk < layerstubs.size(); layerdisk++) {
0223         globals_->ofstream("stubslayer.txt") << layerdisk << " " << layerstubs[layerdisk] << endl;
0224       }
0225     }
0226
0227     addStubTimer_.stop();
0228
0229     // ----------------------------------------------------------------------------------------
0230     // Now start the tracklet processing
0231
0232     // VM router
0233     InputRouterTimer_.start();
0234     sector_->executeIR();
0235     if (settings_->writeMem() && k == settings_->writememsect()) {
0236       sector_->writeDTCStubs(first);
0237       sector_->writeIRStubs(first);
0238     }
0239     InputRouterTimer_.stop();
0240
0241     VMRouterTimer_.start();
0242     sector_->executeVMR();
0243     if (settings_->writeMem() && k == settings_->writememsect()) {
0244       sector_->writeVMSTE(first);
0245       sector_->writeVMSME(first);
0246       sector_->writeAS(first);
0247       sector_->writeAIS(first);
0248     }
0249     VMRouterTimer_.stop();
0250
0251     // tracklet engine
0252     TETimer_.start();
0253     sector_->executeTE();
0254     TETimer_.stop();
0255
0256     // tracklet engine displaced
0257     TEDTimer_.start();
0258     sector_->executeTED();
0259     TEDTimer_.stop();
0260
0261     // triplet engine
0262     TRETimer_.start();
0263     sector_->executeTRE();
0264     if (settings_->writeMem() && k == settings_->writememsect()) {
0265       sector_->writeST(first);
0266     }
0267     TRETimer_.stop();
0268
0269     // tracklet processor (alternative implementation to TE+TC)
0270     TPTimer_.start();
0271     sector_->executeTP();
0272     TPTimer_.stop();
0273
0274     if (settings_->writeMem() && k == settings_->writememsect()) {
0275       sector_->writeSP(first);
0276     }
0277
0278     // tracklet calculator
0279     TCTimer_.start();
0280     sector_->executeTC();
0281     TCTimer_.stop();
0282
0283     if (settings_->writeMonitorData("HitEff") || settings_->bookHistos()) {
0284       int nTP = globals_->event()->nsimtracks();
0285       for (int iTP = 0; iTP < nTP; iTP++) {
0286         L1SimTrack simtrk = globals_->event()->simtrack(iTP);
0287         if (simtrk.pt() < 2.0)
0288           continue;
0289         if (std::abs(simtrk.vz()) > 15.0)
0290           continue;
0291         if (hypot(simtrk.vx(), simtrk.vy()) > 0.1)
0292           continue;
0293         bool electron = (abs(simtrk.type()) == 11);
0294         bool muon = (abs(simtrk.type()) == 13);
0295         bool pion = (abs(simtrk.type()) == 211);
0296         bool kaon = (abs(simtrk.type()) == 321);
0297         bool proton = (abs(simtrk.type()) == 2212);
0298         if (!(electron || muon || pion || kaon || proton))
0299           continue;
0300         int nlayers = 0;
0301         int ndisks = 0;
0302         int simtrackid = simtrk.trackid();
0303         unsigned int hitmask = 0;
0304         hitmask = ev.layersHit(simtrackid, nlayers, ndisks);  // FIX CPU use.
0305         if (nlayers + ndisks < 4)
0306           continue;
0307
0308         if (settings_->writeMonitorData("HitEff")) {
0309           static ofstream outhit("hiteff.txt");
0310           outhit << simtrk.eta() << " " << (hitmask & 1) << " " << (hitmask & 2) << " " << (hitmask & 4) << " "
0311                  << (hitmask & 8) << " " << (hitmask & 16) << " " << (hitmask & 32) << " " << (hitmask & 64) << " "
0312                  << (hitmask & 128) << " " << (hitmask & 256) << " " << (hitmask & 512) << " " << (hitmask & 1024)
0313                  << endl;
0314         }
0315
0316         std::unordered_set<int> matchseed;
0317         std::unordered_set<int> matchseedtmp = sector_->seedMatch(iTP);
0318         matchseed.insert(matchseedtmp.begin(), matchseedtmp.end());
0319         if (settings_->bookHistos()) {
0320           for (int iseed = 0; iseed < 8; iseed++) {
0321             bool eff = matchseed.find(iseed) != matchseed.end();
0322             globals_->histograms()->fillSeedEff(iseed, simtrk.eta(), eff);
0323           }
0324         }
0325       }
0326     }
0327
0328     // tracklet calculator displaced
0329     TCDTimer_.start();
0330     sector_->executeTCD();
0331     TCDTimer_.stop();
0332
0333     // tracklet processor displaced
0334     TPDTimer_.start();
0335     sector_->executeTPD();
0336     TPDTimer_.stop();
0337
0338     if (settings_->writeMem() && k == settings_->writememsect()) {
0339       sector_->writeTPAR(first);
0340       sector_->writeTPROJ(first);
0341     }
0342
0343     // projection router
0344     PRTimer_.start();
0345     sector_->executePR();
0346     if (settings_->writeMem() && k == settings_->writememsect()) {
0347       sector_->writeVMPROJ(first);
0348       sector_->writeAP(first);
0349     }
0350     PRTimer_.stop();
0351
0352     // match engine
0353     METimer_.start();
0354     sector_->executeME();
0355     if (settings_->writeMem() && k == settings_->writememsect()) {
0356       sector_->writeCM(first);
0357     }
0358     METimer_.stop();
0359
0360     // match calculator
0361     MCTimer_.start();
0362     sector_->executeMC();
0363     MCTimer_.stop();
0364
0365     // match processor (alternative to ME+MC)
0366     MPTimer_.start();
0367     sector_->executeMP();
0368     MPTimer_.stop();
0369
0370     if (settings_->writeMem() && k == settings_->writememsect()) {
0371       sector_->writeMC(first);
0372     }
0373
0374     // fit track
0375     FTTimer_.start();
0376     sector_->executeFT(streamsTrackRaw, streamsStubRaw);
0377     if ((settings_->writeMem() || settings_->writeMonitorData("IFit")) && k == settings_->writememsect()) {
0378       sector_->writeTF(first);
0379     }
0380     FTTimer_.stop();
0381
0382     // purge duplicate
0383     PDTimer_.start();
0384     sector_->executePD(tracks_);
0385     if (((settings_->writeMem() || settings_->writeMonitorData("IFit")) && k == settings_->writememsect()) ||
0386         settings_->writeMonitorData("CT")) {
0387       sector_->writeCT(first);
0388     }
0389     PDTimer_.stop();
0390   }
0391 }
0392
0393 void TrackletEventProcessor::printSummary() {
0394   if (settings_->bookHistos()) {
0395     globals_->histograms()->close();
0396   }
0397
0398   edm::LogVerbatim("Tracklet") << "Process             Times called   Average time (ms)      Total time (s)"
0399                                << "\n"
0400                                << "Cleaning              " << setw(10) << cleanTimer_.ntimes() << setw(20)
0401                                << setprecision(3) << cleanTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0402                                << cleanTimer_.tottime() << "\n"
0403                                << "Add Stubs             " << setw(10) << addStubTimer_.ntimes() << setw(20)
0404                                << setprecision(3) << addStubTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0405                                << addStubTimer_.tottime() << "\n"
0406                                << "InputRouter           " << setw(10) << InputRouterTimer_.ntimes() << setw(20)
0407                                << setprecision(3) << InputRouterTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0408                                << InputRouterTimer_.tottime() << "\n"
0409                                << "VMRouter              " << setw(10) << VMRouterTimer_.ntimes() << setw(20)
0410                                << setprecision(3) << VMRouterTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0411                                << VMRouterTimer_.tottime();
0412   if (settings_->combined()) {
0413     edm::LogVerbatim("Tracklet") << "TrackletProcessor     " << setw(10) << TPTimer_.ntimes() << setw(20)
0414                                  << setprecision(3) << TPTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0415                                  << TPTimer_.tottime() << "\n"
0416                                  << "MatchProcessor        " << setw(10) << MPTimer_.ntimes() << setw(20)
0417                                  << setprecision(3) << MPTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0418                                  << MPTimer_.tottime();
0419   } else {
0420     edm::LogVerbatim("Tracklet") << "TrackletEngine        " << setw(10) << TETimer_.ntimes() << setw(20)
0421                                  << setprecision(3) << TETimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0422                                  << TETimer_.tottime();
0423     if (settings_->extended()) {
0424       edm::LogVerbatim("Tracklet") << "TrackletEngineDisplaced" << setw(10) << TEDTimer_.ntimes() << setw(20)
0425                                    << setprecision(3) << TEDTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0426                                    << TEDTimer_.tottime() << "\n"
0427                                    << "TripletEngine         " << setw(10) << TRETimer_.ntimes() << setw(20)
0428                                    << setprecision(3) << TRETimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0429                                    << TRETimer_.tottime() << "\n"
0430                                    << "TrackletCalculatorDisplaced" << setw(10) << TCDTimer_.ntimes() << setw(20)
0431                                    << setprecision(3) << TCDTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0432                                    << TCDTimer_.tottime() << "\n"
0433                                    << TCDTimer_.tottime() << "\n"
0434                                    << "TrackletProcessorDisplaced" << setw(10) << TPDTimer_.ntimes() << setw(20)
0435                                    << setprecision(3) << TPDTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0436                                    << TPDTimer_.tottime();
0437     }
0438     edm::LogVerbatim("Tracklet") << "TrackletCalculator    " << setw(10) << TCTimer_.ntimes() << setw(20)
0439                                  << setprecision(3) << TCTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0440                                  << TCTimer_.tottime() << "\n"
0441                                  << "ProjectionRouter      " << setw(10) << PRTimer_.ntimes() << setw(20)
0442                                  << setprecision(3) << PRTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0443                                  << PRTimer_.tottime() << "\n"
0444                                  << "MatchEngine           " << setw(10) << METimer_.ntimes() << setw(20)
0445                                  << setprecision(3) << METimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0446                                  << METimer_.tottime() << "\n"
0447                                  << "MatchCalculator       " << setw(10) << MCTimer_.ntimes() << setw(20)
0448                                  << setprecision(3) << MCTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0449                                  << MCTimer_.tottime();
0450   }
0451   edm::LogVerbatim("Tracklet") << "FitTrack              " << setw(10) << FTTimer_.ntimes() << setw(20)
0452                                << setprecision(3) << FTTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0453                                << FTTimer_.tottime() << "\n"
0454                                << "PurgeDuplicate        " << setw(10) << PDTimer_.ntimes() << setw(20)
0455                                << setprecision(3) << PDTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0456                                << PDTimer_.tottime();
0457 }