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
0011 #include "DataFormats/Math/interface/deltaPhi.h"
0012
0013 #include <iomanip>
0014 #include <filesystem>
0015
0016 using namespace trklet;
0017 using namespace std;
0018
0019 TrackletEventProcessor::TrackletEventProcessor() : settings_(nullptr) {}
0020
0021 TrackletEventProcessor::~TrackletEventProcessor() {
0022   if (settings_ && settings_->bookHistos()) {
0023     histbase_->close();
0024   }
0025 }
0026
0027 void TrackletEventProcessor::init(Settings const& theSettings) {
0028   settings_ = &theSettings;
0029
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()) {
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_);
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   globals_->event() = &ev;
0178
0179   tracks_.clear();
0180
0181   eventnum_++;
0182   bool first = (eventnum_ == 1);
0183
0184   for (unsigned int k = 0; k < N_SECTOR; k++) {
0185     sector_->setSector(k);
0186
0187     cleanTimer_.start();
0188     sector_->clean();
0189     cleanTimer_.stop();
0190
0191     addStubTimer_.start();
0192
0193     vector<int> layerstubs(N_LAYER + N_DISK, 0);
0194     vector<int> layerstubssector(N_SECTOR * (N_LAYER + N_DISK), 0);
0195
0196     for (int j = 0; j < ev.nstubs(); j++) {
0197       const L1TStub& stub = ev.stub(j);
0198       unsigned int isector = stub.region();
0199       if (isector != k) {
0200         continue;
0201       }
0202
0203       const string& dtc = stub.DTClink();
0204
0205       layerstubs[stub.layerdisk()]++;
0206       layerstubssector[isector * (N_LAYER + N_DISK) + stub.layerdisk()]++;
0207
0208       sector_->addStub(stub, dtc);
0209     }
0210
0211     if (settings_->writeMonitorData("StubsLayerSector")) {
0212       for (unsigned int index = 0; index < layerstubssector.size(); index++) {
0213         int layerdisk = index % (N_LAYER + N_DISK);
0214         int sector = index / (N_LAYER + N_DISK);
0215         globals_->ofstream("stubslayersector.txt")
0216             << layerdisk << " " << sector << " " << layerstubssector[index] << endl;
0217       }
0218     }
0219
0220     if (settings_->writeMonitorData("StubsLayer")) {
0221       for (unsigned int layerdisk = 0; layerdisk < layerstubs.size(); layerdisk++) {
0222         globals_->ofstream("stubslayer.txt") << layerdisk << " " << layerstubs[layerdisk] << endl;
0223       }
0224     }
0225
0226     addStubTimer_.stop();
0227
0228     // ----------------------------------------------------------------------------------------
0229     // Now start the tracklet processing
0230
0231     // VM router
0232     InputRouterTimer_.start();
0233     sector_->executeIR();
0234     if (settings_->writeMem() && k == settings_->writememsect()) {
0235       sector_->writeDTCStubs(first);
0236       sector_->writeIRStubs(first);
0237     }
0238     InputRouterTimer_.stop();
0239
0240     VMRouterTimer_.start();
0241     sector_->executeVMR();
0242     if (settings_->writeMem() && k == settings_->writememsect()) {
0243       sector_->writeVMSTE(first);
0244       sector_->writeVMSME(first);
0245       sector_->writeAS(first);
0246       sector_->writeAIS(first);
0247     }
0248     VMRouterTimer_.stop();
0249
0250     // tracklet engine
0251     TETimer_.start();
0252     sector_->executeTE();
0253     TETimer_.stop();
0254
0255     // tracklet engine displaced
0256     TEDTimer_.start();
0257     sector_->executeTED();
0258     TEDTimer_.stop();
0259
0260     // triplet engine
0261     TRETimer_.start();
0262     sector_->executeTRE();
0263     if (settings_->writeMem() && k == settings_->writememsect()) {
0264       sector_->writeST(first);
0265     }
0266     TRETimer_.stop();
0267
0268     // tracklet processor (alternative implementation to TE+TC)
0269     TPTimer_.start();
0270     sector_->executeTP();
0271     TPTimer_.stop();
0272
0273     if (settings_->writeMem() && k == settings_->writememsect()) {
0274       sector_->writeSP(first);
0275     }
0276
0277     // tracklet calculator
0278     TCTimer_.start();
0279     sector_->executeTC();
0280     TCTimer_.stop();
0281
0282     int nTP = globals_->event()->nsimtracks();
0283     for (int iTP = 0; iTP < nTP; iTP++) {
0284       L1SimTrack simtrk = globals_->event()->simtrack(iTP);
0285       if (simtrk.pt() < 2.0)
0286         continue;
0287       if (std::abs(simtrk.vz()) > 15.0)
0288         continue;
0289       if (hypot(simtrk.vx(), simtrk.vy()) > 0.1)
0290         continue;
0291       bool electron = (abs(simtrk.type()) == 11);
0292       bool muon = (abs(simtrk.type()) == 13);
0293       bool pion = (abs(simtrk.type()) == 211);
0294       bool kaon = (abs(simtrk.type()) == 321);
0295       bool proton = (abs(simtrk.type()) == 2212);
0296       if (!(electron || muon || pion || kaon || proton))
0297         continue;
0298       int nlayers = 0;
0299       int ndisks = 0;
0300       int simtrackid = simtrk.trackid();
0301       unsigned int hitmask = ev.layersHit(simtrackid, nlayers, ndisks);
0302       if (nlayers + ndisks < 4)
0303         continue;
0304
0305       if (settings_->writeMonitorData("HitEff")) {
0306         static ofstream outhit("hiteff.txt");
0307         outhit << simtrk.eta() << " " << (hitmask & 1) << " " << (hitmask & 2) << " " << (hitmask & 4) << " "
0308                << (hitmask & 8) << " " << (hitmask & 16) << " " << (hitmask & 32) << " " << (hitmask & 64) << " "
0309                << (hitmask & 128) << " " << (hitmask & 256) << " " << (hitmask & 512) << " " << (hitmask & 1024)
0310                << endl;
0311       }
0312
0313       std::unordered_set<int> matchseed;
0314       std::unordered_set<int> matchseedtmp = sector_->seedMatch(iTP);
0315       matchseed.insert(matchseedtmp.begin(), matchseedtmp.end());
0316       if (settings_->bookHistos()) {
0317         for (int iseed = 0; iseed < 8; iseed++) {
0318           bool eff = matchseed.find(iseed) != matchseed.end();
0319           globals_->histograms()->fillSeedEff(iseed, simtrk.eta(), eff);
0320         }
0321       }
0322     }
0323
0324     // tracklet calculator displaced
0325     TCDTimer_.start();
0326     sector_->executeTCD();
0327     TCDTimer_.stop();
0328
0329     if (settings_->writeMem() && k == settings_->writememsect()) {
0330       sector_->writeTPAR(first);
0331       sector_->writeTPROJ(first);
0332     }
0333
0334     // projection router
0335     PRTimer_.start();
0336     sector_->executePR();
0337     if (settings_->writeMem() && k == settings_->writememsect()) {
0338       sector_->writeVMPROJ(first);
0339       sector_->writeAP(first);
0340     }
0341     PRTimer_.stop();
0342
0343     // match engine
0344     METimer_.start();
0345     sector_->executeME();
0346     if (settings_->writeMem() && k == settings_->writememsect()) {
0347       sector_->writeCM(first);
0348     }
0349     METimer_.stop();
0350
0351     // match calculator
0352     MCTimer_.start();
0353     sector_->executeMC();
0354     MCTimer_.stop();
0355
0356     // match processor (alternative to ME+MC)
0357     MPTimer_.start();
0358     sector_->executeMP();
0359     MPTimer_.stop();
0360
0361     if (settings_->writeMem() && k == settings_->writememsect()) {
0362       sector_->writeMC(first);
0363     }
0364
0365     // fit track
0366     FTTimer_.start();
0367     sector_->executeFT();
0368     if ((settings_->writeMem() || settings_->writeMonitorData("IFit")) && k == settings_->writememsect()) {
0369       sector_->writeTF(first);
0370     }
0371     FTTimer_.stop();
0372
0373     // purge duplicate
0374     PDTimer_.start();
0375     sector_->executePD(tracks_);
0376     if (((settings_->writeMem() || settings_->writeMonitorData("IFit")) && k == settings_->writememsect()) ||
0377         settings_->writeMonitorData("CT")) {
0378       sector_->writeCT(first);
0379     }
0380     PDTimer_.stop();
0381   }
0382 }
0383
0384 void TrackletEventProcessor::printSummary() {
0385   if (settings_->bookHistos()) {
0386     globals_->histograms()->close();
0387   }
0388
0389   edm::LogVerbatim("Tracklet") << "Process             Times called   Average time (ms)      Total time (s)"
0390                                << "\n"
0391                                << "Cleaning              " << setw(10) << cleanTimer_.ntimes() << setw(20)
0392                                << setprecision(3) << cleanTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0393                                << cleanTimer_.tottime() << "\n"
0394                                << "Add Stubs             " << setw(10) << addStubTimer_.ntimes() << setw(20)
0395                                << setprecision(3) << addStubTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0396                                << addStubTimer_.tottime() << "\n"
0397                                << "InputRouter           " << setw(10) << InputRouterTimer_.ntimes() << setw(20)
0398                                << setprecision(3) << InputRouterTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0399                                << InputRouterTimer_.tottime() << "\n"
0400                                << "VMRouter              " << setw(10) << VMRouterTimer_.ntimes() << setw(20)
0401                                << setprecision(3) << VMRouterTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0402                                << VMRouterTimer_.tottime();
0403   if (settings_->combined()) {
0404     edm::LogVerbatim("Tracklet") << "TrackletProcessor     " << setw(10) << TPTimer_.ntimes() << setw(20)
0405                                  << setprecision(3) << TPTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0406                                  << TPTimer_.tottime() << "\n"
0407                                  << "MatchProcessor        " << setw(10) << MPTimer_.ntimes() << setw(20)
0408                                  << setprecision(3) << MPTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0409                                  << MPTimer_.tottime();
0410   } else {
0411     edm::LogVerbatim("Tracklet") << "TrackletEngine        " << setw(10) << TETimer_.ntimes() << setw(20)
0412                                  << setprecision(3) << TETimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0413                                  << TETimer_.tottime();
0414     if (settings_->extended()) {
0415       edm::LogVerbatim("Tracklet") << "TrackletEngineDisplaced" << setw(10) << TEDTimer_.ntimes() << setw(20)
0416                                    << setprecision(3) << TEDTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0417                                    << TEDTimer_.tottime() << "\n"
0418                                    << "TripletEngine         " << setw(10) << TRETimer_.ntimes() << setw(20)
0419                                    << setprecision(3) << TRETimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0420                                    << TRETimer_.tottime() << "\n"
0421                                    << "TrackletCalculatorDisplaced" << setw(10) << TCDTimer_.ntimes() << setw(20)
0422                                    << setprecision(3) << TCDTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0423                                    << TCDTimer_.tottime();
0424     }
0425     edm::LogVerbatim("Tracklet") << "TrackletCalculator    " << setw(10) << TCTimer_.ntimes() << setw(20)
0426                                  << setprecision(3) << TCTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0427                                  << TCTimer_.tottime() << "\n"
0428                                  << "ProjectionRouter      " << setw(10) << PRTimer_.ntimes() << setw(20)
0429                                  << setprecision(3) << PRTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0430                                  << PRTimer_.tottime() << "\n"
0431                                  << "MatchEngine           " << setw(10) << METimer_.ntimes() << setw(20)
0432                                  << setprecision(3) << METimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0433                                  << METimer_.tottime() << "\n"
0434                                  << "MatchCalculator       " << setw(10) << MCTimer_.ntimes() << setw(20)
0435                                  << setprecision(3) << MCTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0436                                  << MCTimer_.tottime();
0437   }
0438   edm::LogVerbatim("Tracklet") << "FitTrack              " << setw(10) << FTTimer_.ntimes() << setw(20)
0439                                << setprecision(3) << FTTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0440                                << FTTimer_.tottime() << "\n"
0441                                << "PurgeDuplicate        " << setw(10) << PDTimer_.ntimes() << setw(20)
0442                                << setprecision(3) << PDTimer_.avgtime() * 1000.0 << setw(20) << setprecision(3)
0443                                << PDTimer_.tottime();
0444 }