File indexing completed on 2021-07-07 22:33:41
0001 #include "L1Trigger/TrackFindingTracklet/interface/TripletEngine.h"
0002 #include "L1Trigger/TrackFindingTracklet/interface/Settings.h"
0003 #include "L1Trigger/TrackFindingTracklet/interface/Globals.h"
0004
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "FWCore/Utilities/interface/Exception.h"
0007
0008 #include <algorithm>
0009
0010 using namespace std;
0011 using namespace trklet;
0012
0013 TripletEngine::TripletEngine(string name, Settings const &settings, Globals *global)
0014 : ProcessBase(name, settings, global) {
0015 stubpairs_.clear();
0016 thirdvmstubs_.clear();
0017 layer1_ = 0;
0018 layer2_ = 0;
0019 layer3_ = 0;
0020 disk1_ = 0;
0021 disk2_ = 0;
0022 disk3_ = 0;
0023 dct1_ = 0;
0024 dct2_ = 0;
0025 dct3_ = 0;
0026 phi1_ = 0;
0027 phi2_ = 0;
0028 phi3_ = 0;
0029 z1_ = 0;
0030 z2_ = 0;
0031 z3_ = 0;
0032 r1_ = 0;
0033 r2_ = 0;
0034 r3_ = 0;
0035
0036 if (name_[4] == 'L')
0037 layer1_ = name_[5] - '0';
0038 if (name_[4] == 'D')
0039 disk1_ = name_[5] - '0';
0040 if (name_[7] == 'L')
0041 layer2_ = name_[8] - '0';
0042 if (name_[7] == 'D')
0043 disk2_ = name_[8] - '0';
0044
0045 if (layer1_ == 3 && layer2_ == 4) {
0046 layer3_ = 2;
0047 iSeed_ = 8;
0048 } else if (layer1_ == 5 && layer2_ == 6) {
0049 layer3_ = 4;
0050 iSeed_ = 9;
0051 } else if (layer1_ == 2 && layer2_ == 3) {
0052 disk3_ = 1;
0053 iSeed_ = 10;
0054 } else if (disk1_ == 1 && disk2_ == 2) {
0055 layer3_ = 2;
0056 iSeed_ = 11;
0057 } else
0058 throw cms::Exception("LogicError") << __FILE__ << " " << __LINE__ << " Invalid seeding!";
0059
0060 if ((layer2_ == 4 && layer3_ == 2) || (layer2_ == 6 && layer3_ == 4)) {
0061 secondphibits_ = settings_.nfinephi(1, iSeed_);
0062 thirdphibits_ = settings_.nfinephi(2, iSeed_);
0063 }
0064 if ((layer2_ == 3 && disk3_ == 1) || (disk2_ == 2 && layer3_ == 2)) {
0065 secondphibits_ = settings_.nfinephi(1, iSeed_);
0066 thirdphibits_ = settings_.nfinephi(2, iSeed_);
0067 }
0068 if (settings_.enableTripletTables() && !settings_.writeTripletTables())
0069 readTables();
0070 }
0071
0072 TripletEngine::~TripletEngine() {
0073 if (settings_.writeTripletTables())
0074 writeTables();
0075 }
0076
0077 void TripletEngine::addOutput(MemoryBase *memory, string output) {
0078 if (settings_.writetrace()) {
0079 edm::LogVerbatim("Tracklet") << "In " << name_ << " adding output to " << memory->getName() << " to output "
0080 << output;
0081 }
0082 if (output == "stubtripout") {
0083 auto *tmp = dynamic_cast<StubTripletsMemory *>(memory);
0084 assert(tmp != nullptr);
0085 stubtriplets_ = tmp;
0086 return;
0087 }
0088 throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " Could not find output : " << output;
0089 }
0090
0091 void TripletEngine::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 == "thirdvmstubin") {
0097 auto *tmp = dynamic_cast<VMStubsTEMemory *>(memory);
0098 assert(tmp != nullptr);
0099 thirdvmstubs_.push_back(tmp);
0100 return;
0101 }
0102 if (input.substr(0, 8) == "stubpair") {
0103 auto *tmp = dynamic_cast<StubPairsMemory *>(memory);
0104 assert(tmp != nullptr);
0105 stubpairs_.push_back(tmp);
0106 return;
0107 }
0108 throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " Could not find input : " << input;
0109 }
0110
0111 void TripletEngine::execute() {
0112 unsigned int countall = 0;
0113 unsigned int countpass = 0;
0114 unsigned int nThirdStubs = 0;
0115 count_ = 0;
0116
0117 for (unsigned int iThirdMem = 0; iThirdMem < thirdvmstubs_.size();
0118 nThirdStubs += thirdvmstubs_.at(iThirdMem)->nVMStubs(), iThirdMem++)
0119 ;
0120
0121 assert(!thirdvmstubs_.empty());
0122 assert(!stubpairs_.empty());
0123
0124 bool print = false && (getName().substr(0, 10) == "TRE_L2cL3c");
0125
0126 print = print && nThirdStubs > 0;
0127
0128 int hacksum = 0;
0129 if (print) {
0130 edm::LogVerbatim("Tracklet") << "In TripletEngine::execute : " << getName() << " " << nThirdStubs << ":";
0131 for (unsigned int i = 0; i < thirdvmstubs_.size(); ++i) {
0132 edm::LogVerbatim("Tracklet") << thirdvmstubs_.at(i)->getName() << " " << thirdvmstubs_.at(i)->nVMStubs();
0133 }
0134 int s = 0;
0135 std::string oss = "";
0136 for (unsigned int i = 0; i < stubpairs_.size(); ++i) {
0137 oss += std::to_string(stubpairs_.at(i)->nStubPairs());
0138 oss += " ";
0139 s += stubpairs_.at(i)->nStubPairs();
0140 }
0141 hacksum += nThirdStubs * s;
0142 edm::LogVerbatim("Tracklet") << oss;
0143 for (unsigned int i = 0; i < stubpairs_.size(); ++i) {
0144 edm::LogVerbatim("Tracklet") << " " << stubpairs_.at(i)->getName();
0145 }
0146 }
0147
0148 tmpSPTable_.clear();
0149
0150 for (unsigned int i = 0; i < stubpairs_.size(); ++i) {
0151 for (unsigned int j = 0; j < stubpairs_.at(i)->nStubPairs(); ++j) {
0152 if (print)
0153 edm::LogVerbatim("Tracklet") << " ***** " << stubpairs_.at(i)->getName() << " "
0154 << stubpairs_.at(i)->nStubPairs();
0155
0156 auto firstvmstub = stubpairs_.at(i)->getVMStub1(j);
0157 auto secondvmstub = stubpairs_.at(i)->getVMStub2(j);
0158
0159 if ((layer2_ == 4 && layer3_ == 2) || (layer2_ == 6 && layer3_ == 4)) {
0160 constexpr unsigned int vmbitshift = 10;
0161 int lookupbits = (int)((firstvmstub.vmbits().value() >> vmbitshift) & 1023);
0162 int newbin = (lookupbits & 127);
0163 int bin = newbin / 8;
0164
0165 int start = (bin >> 1);
0166 int last = start + (bin & 1);
0167
0168 for (int ibin = start; ibin <= last; ibin++) {
0169 for (unsigned int k = 0; k < thirdvmstubs_.size(); k++) {
0170 string vmsteSuffix = thirdvmstubs_.at(k)->getLastPartOfName();
0171 vmsteSuffix = vmsteSuffix.substr(0, vmsteSuffix.find_last_of('n'));
0172 if (stubpairs_.at(i)->getLastPartOfName() != vmsteSuffix)
0173 continue;
0174 for (unsigned int l = 0; l < thirdvmstubs_.at(k)->nVMStubsBinned(ibin); l++) {
0175 if (settings_.debugTracklet()) {
0176 edm::LogVerbatim("Tracklet") << "In " << getName() << " have third stub";
0177 }
0178
0179 if (countall >= settings_.maxStep("TRE"))
0180 break;
0181 countall++;
0182
0183 const VMStubTE &thirdvmstub = thirdvmstubs_.at(k)->getVMStubTEBinned(ibin, l);
0184
0185 assert(secondphibits_ != -1);
0186 assert(thirdphibits_ != -1);
0187
0188 unsigned int nvmsecond = settings_.nallstubs(layer2_ - 1) * settings_.nvmte(1, iSeed_);
0189 unsigned int nvmbitssecond = nbits(nvmsecond);
0190
0191 FPGAWord iphisecondbin = secondvmstub.stub()->iphivmFineBins(nvmbitssecond, secondphibits_);
0192
0193
0194
0195 FPGAWord iphithirdbin = thirdvmstub.finephi();
0196
0197 unsigned int index = (iphisecondbin.value() << thirdphibits_) + iphithirdbin.value();
0198
0199 FPGAWord secondbend = secondvmstub.bend();
0200 FPGAWord thirdbend = thirdvmstub.bend();
0201
0202 index = (index << secondbend.nbits()) + secondbend.value();
0203 index = (index << thirdbend.nbits()) + thirdbend.value();
0204
0205 if ((settings_.enableTripletTables() && !settings_.writeTripletTables()) &&
0206 (index >= table_.size() || !table_[index])) {
0207 if (settings_.debugTracklet()) {
0208 edm::LogVerbatim("Tracklet")
0209 << "Stub pair rejected because of stub pt cut bends : "
0210 << settings_.benddecode(secondvmstub.bend().value(), layer2_ - 1, secondvmstub.isPSmodule())
0211 << " " << settings_.benddecode(thirdvmstub.bend().value(), layer3_ - 1, thirdvmstub.isPSmodule());
0212 }
0213
0214
0215
0216
0217 }
0218 if (settings_.writeTripletTables()) {
0219 if (index >= table_.size())
0220 table_.resize(index + 1, false);
0221 table_[index] = true;
0222
0223 const unsigned spIndex = stubpairs_.at(i)->getIndex(j);
0224 const string &tedName = stubpairs_.at(i)->getTEDName(j);
0225 if (!tmpSPTable_.count(tedName))
0226 tmpSPTable_[tedName];
0227 if (spIndex >= tmpSPTable_.at(tedName).size())
0228 tmpSPTable_.at(tedName).resize(spIndex + 1);
0229 tmpSPTable_.at(tedName).at(spIndex).push_back(stubpairs_.at(i)->getName());
0230 }
0231
0232 if (settings_.debugTracklet())
0233 edm::LogVerbatim("Tracklet") << "Adding layer-layer pair in " << getName();
0234 if (settings_.writeMonitorData("Seeds")) {
0235 ofstream fout("seeds.txt", ofstream::app);
0236 fout << __FILE__ << ":" << __LINE__ << " " << name_ << " " << iSeed_ << endl;
0237 fout.close();
0238 }
0239 stubtriplets_->addStubs(thirdvmstub.stub(),
0240 (stubpairs_.at(i))->getVMStub1(j).stub(),
0241 (stubpairs_.at(i))->getVMStub2(j).stub());
0242
0243 countpass++;
0244 }
0245 }
0246 }
0247
0248 }
0249
0250 else if (disk2_ == 2 && layer3_ == 2) {
0251 int lookupbits = (int)((firstvmstub.vmbits().value() >> 10) & 1023);
0252 int newbin = (lookupbits & 127);
0253 int bin = newbin / 8;
0254
0255 int start = (bin >> 1);
0256 int last = start + (bin & 1);
0257
0258 if (firstvmstub.stub()->disk().value() < 0) {
0259 start = settings_.NLONGVMBINS() - last - 1;
0260 last = settings_.NLONGVMBINS() - start - 1;
0261 }
0262
0263 for (int ibin = start; ibin <= last; ibin++) {
0264 for (unsigned int k = 0; k < thirdvmstubs_.size(); k++) {
0265 string vmsteSuffix = thirdvmstubs_.at(k)->getLastPartOfName();
0266 vmsteSuffix = vmsteSuffix.substr(0, vmsteSuffix.find_last_of('n'));
0267 if (stubpairs_.at(i)->getLastPartOfName() != vmsteSuffix)
0268 continue;
0269
0270 for (unsigned int l = 0; l < thirdvmstubs_.at(k)->nVMStubsBinned(ibin); l++) {
0271 if (countall >= settings_.maxStep("TRE"))
0272 break;
0273 countall++;
0274
0275 const VMStubTE &thirdvmstub = thirdvmstubs_.at(k)->getVMStubTEBinned(ibin, l);
0276
0277 assert(secondphibits_ != -1);
0278 assert(thirdphibits_ != -1);
0279
0280 FPGAWord iphisecondbin = secondvmstub.finephi();
0281 FPGAWord iphithirdbin = thirdvmstub.finephi();
0282
0283 unsigned int index = (iphisecondbin.value() << thirdphibits_) + iphithirdbin.value();
0284
0285 FPGAWord secondbend = secondvmstub.bend();
0286 FPGAWord thirdbend = thirdvmstub.bend();
0287
0288 index = (index << secondbend.nbits()) + secondbend.value();
0289 index = (index << thirdbend.nbits()) + thirdbend.value();
0290
0291 if ((settings_.enableTripletTables() && !settings_.writeTripletTables()) &&
0292 (index >= table_.size() || !table_[index])) {
0293 if (settings_.debugTracklet()) {
0294 edm::LogVerbatim("Tracklet")
0295 << "Stub triplet rejected because of stub pt cut bends : "
0296 << settings_.benddecode(secondvmstub.bend().value(), disk2_ + 5, secondvmstub.isPSmodule()) << " "
0297 << settings_.benddecode(thirdvmstub.bend().value(), layer3_ - 1, thirdvmstub.isPSmodule());
0298 }
0299 continue;
0300 }
0301 if (settings_.writeTripletTables()) {
0302 if (index >= table_.size())
0303 table_.resize(index + 1, false);
0304 table_[index] = true;
0305
0306 const unsigned spIndex = stubpairs_.at(i)->getIndex(j);
0307 const string &tedName = stubpairs_.at(i)->getTEDName(j);
0308 if (!tmpSPTable_.count(tedName))
0309 tmpSPTable_[tedName];
0310 if (spIndex >= tmpSPTable_.at(tedName).size())
0311 tmpSPTable_.at(tedName).resize(spIndex + 1);
0312 tmpSPTable_.at(tedName).at(spIndex).push_back(stubpairs_.at(i)->getName());
0313 }
0314
0315 if (settings_.debugTracklet())
0316 edm::LogVerbatim("Tracklet") << "Adding layer-disk pair in " << getName();
0317 if (settings_.writeMonitorData("Seeds")) {
0318 ofstream fout("seeds.txt", ofstream::app);
0319 fout << __FILE__ << ":" << __LINE__ << " " << name_ << " " << iSeed_ << endl;
0320 fout.close();
0321 }
0322 stubtriplets_->addStubs(thirdvmstub.stub(),
0323 (stubpairs_.at(i))->getVMStub1(j).stub(),
0324 (stubpairs_.at(i))->getVMStub2(j).stub());
0325 countpass++;
0326 }
0327 }
0328 }
0329 }
0330
0331 else if (layer2_ == 3 && disk3_ == 1) {
0332 int lookupbits = (int)((firstvmstub.vmbits().value() >> 10) & 1023);
0333
0334 int newbin = (lookupbits & 127);
0335 int bin = newbin / 8;
0336
0337 int start = (bin >> 1);
0338 int last = start + (bin & 1);
0339
0340 for (int ibin = start; ibin <= last; ibin++) {
0341 for (unsigned int k = 0; k < thirdvmstubs_.size(); k++) {
0342 string vmsteSuffix = thirdvmstubs_.at(k)->getLastPartOfName();
0343 vmsteSuffix = vmsteSuffix.substr(0, vmsteSuffix.find_last_of('n'));
0344 if (stubpairs_.at(i)->getLastPartOfName() != vmsteSuffix)
0345 continue;
0346 for (unsigned int l = 0; l < thirdvmstubs_.at(k)->nVMStubsBinned(ibin); l++) {
0347 if (countall >= settings_.maxStep("TRE"))
0348 break;
0349 countall++;
0350
0351 const VMStubTE &thirdvmstub = thirdvmstubs_.at(k)->getVMStubTEBinned(ibin, l);
0352
0353 assert(secondphibits_ != -1);
0354 assert(thirdphibits_ != -1);
0355
0356 unsigned int nvmsecond;
0357
0358 nvmsecond = settings_.nallstubs(layer2_ - 1) * settings_.nvmte(1, iSeed_);
0359 unsigned int nvmbitssecond = nbits(nvmsecond);
0360
0361 FPGAWord iphisecondbin = secondvmstub.stub()->iphivmFineBins(nvmbitssecond, secondphibits_);
0362
0363
0364
0365 FPGAWord iphithirdbin = thirdvmstub.finephi();
0366
0367 unsigned int index = (iphisecondbin.value() << thirdphibits_) + iphithirdbin.value();
0368
0369 FPGAWord secondbend = secondvmstub.bend();
0370 FPGAWord thirdbend = thirdvmstub.bend();
0371
0372 index = (index << secondbend.nbits()) + secondbend.value();
0373 index = (index << thirdbend.nbits()) + thirdbend.value();
0374
0375 if ((settings_.enableTripletTables() && !settings_.writeTripletTables()) &&
0376 (index >= table_.size() || !table_[index])) {
0377 if (settings_.debugTracklet()) {
0378 edm::LogVerbatim("Tracklet")
0379 << "Stub pair rejected because of stub pt cut bends : "
0380 << settings_.benddecode(secondvmstub.bend().value(), layer2_ - 1, secondvmstub.isPSmodule())
0381 << " " << settings_.benddecode(thirdvmstub.bend().value(), disk3_ + 5, thirdvmstub.isPSmodule());
0382 }
0383 continue;
0384 }
0385 if (settings_.writeTripletTables()) {
0386 if (index >= table_.size())
0387 table_.resize(index + 1, false);
0388 table_[index] = true;
0389
0390 const unsigned spIndex = stubpairs_.at(i)->getIndex(j);
0391 const string &tedName = stubpairs_.at(i)->getTEDName(j);
0392 if (!tmpSPTable_.count(tedName))
0393 tmpSPTable_[tedName];
0394 if (spIndex >= tmpSPTable_.at(tedName).size())
0395 tmpSPTable_.at(tedName).resize(spIndex + 1);
0396 tmpSPTable_.at(tedName).at(spIndex).push_back(stubpairs_.at(i)->getName());
0397 }
0398
0399 if (settings_.debugTracklet())
0400 edm::LogVerbatim("Tracklet") << "Adding layer-disk pair in " << getName();
0401 if (settings_.writeMonitorData("Seeds")) {
0402 ofstream fout("seeds.txt", ofstream::app);
0403 fout << __FILE__ << ":" << __LINE__ << " " << name_ << " " << iSeed_ << endl;
0404 fout.close();
0405 }
0406 stubtriplets_->addStubs(thirdvmstub.stub(),
0407 (stubpairs_.at(i))->getVMStub1(j).stub(),
0408 (stubpairs_.at(i))->getVMStub2(j).stub());
0409 countpass++;
0410 }
0411 }
0412 }
0413 }
0414 }
0415 }
0416
0417 for (const auto &tedName : tmpSPTable_) {
0418 for (unsigned spIndex = 0; spIndex < tedName.second.size(); spIndex++) {
0419 if (tedName.second.at(spIndex).empty())
0420 continue;
0421 vector<string> entry(tedName.second.at(spIndex));
0422 sort(entry.begin(), entry.end());
0423 entry.erase(unique(entry.begin(), entry.end()), entry.end());
0424 const string &spName = entry.at(0);
0425
0426 if (!spTable_.count(tedName.first))
0427 spTable_[tedName.first];
0428 if (spIndex >= spTable_.at(tedName.first).size())
0429 spTable_.at(tedName.first).resize(spIndex + 1);
0430 if (!spTable_.at(tedName.first).at(spIndex).count(spName))
0431 spTable_.at(tedName.first).at(spIndex)[spName] = 0;
0432 spTable_.at(tedName.first).at(spIndex)[spName]++;
0433 }
0434 }
0435
0436 if (settings_.writeMonitorData("TRE")) {
0437 globals_->ofstream("tripletengine.txt") << getName() << " " << countall << " " << countpass << endl;
0438 }
0439 }
0440
0441 void TripletEngine::readTables() {
0442 ifstream fin;
0443 string tableName, word;
0444 unsigned num;
0445
0446 string tablePath = settings_.tableTREFile();
0447 unsigned int finddir = tablePath.find("table_TRE_");
0448 tableName = tablePath.substr(0, finddir) + "table_" + name_ + ".txt";
0449
0450 fin.open(tableName, ifstream::in);
0451 if (!fin) {
0452 throw cms::Exception("BadConfig") << "TripletEngine::readTables, file " << tableName << " not known";
0453 }
0454 while (!fin.eof()) {
0455 fin >> word;
0456 num = atoi(word.c_str());
0457 table_.push_back(num > 0 ? true : false);
0458 }
0459 fin.close();
0460 }
0461
0462 void TripletEngine::writeTables() {
0463 ofstream fout;
0464 stringstream tableName;
0465
0466 tableName << "table/table_" << name_ << ".txt";
0467
0468 fout.open(tableName.str(), ofstream::out);
0469 for (const auto entry : table_)
0470 fout << entry << endl;
0471 fout.close();
0472
0473 for (const auto &tedName : spTable_) {
0474 tableName.str("");
0475 tableName << "table/table_" << tedName.first << "_" << name_ << ".txt";
0476
0477 fout.open(tableName.str(), ofstream::out);
0478 for (const auto &entry : tedName.second) {
0479 for (const auto &spName : entry)
0480 fout << spName.first << ":" << spName.second << " ";
0481 fout << endl;
0482 }
0483 fout.close();
0484 }
0485 }