Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-11-15 23:40:40

0001 //-------------------------------------------------
0002 //
0003 //   Class: L1MuBMEtaProcessor
0004 //
0005 //   Description: Eta Processor
0006 //
0007 //                An Eta Processor consists of :
0008 //                a receiver unit,
0009 //                one Eta Track Finder (ETF) and
0010 //                one Eta Matching Unit (EMU)
0011 //
0012 //
0013 //   Author :
0014 //   N. Neumeister            CERN EP
0015 //   J. Troconiz              UAM Madrid
0016 //Modifications:
0017 // G.Karathanasis     U.Athens
0018 //
0019 //--------------------------------------------------
0020 
0021 //-----------------------
0022 // This Class's Header --
0023 //-----------------------
0024 
0025 #include "L1Trigger/L1TMuonBarrel/src/L1MuBMEtaProcessor.h"
0026 
0027 //---------------
0028 // C++ Headers --
0029 //---------------
0030 
0031 #include <iostream>
0032 #include <iomanip>
0033 #include <bitset>
0034 
0035 //-------------------------------
0036 // Collaborating Class Headers --
0037 //-------------------------------
0038 
0039 #include "L1Trigger/L1TMuonBarrel/interface/L1MuBMTFConfig.h"
0040 #include "L1Trigger/L1TMuonBarrel/src/L1MuBMSectorProcessor.h"
0041 #include "L1Trigger/L1TMuonBarrel/interface/L1MuBMTrackFinder.h"
0042 #include "CondFormats/L1TObjects/interface/L1MuDTEtaPattern.h"
0043 #include "L1Trigger/L1TMuonBarrel/interface/L1MuBMTEtaPatternLut.h"
0044 #include "CondFormats/DataRecord/interface/L1MuDTEtaPatternLutRcd.h"
0045 #include "L1Trigger/L1TMuonBarrel/interface/L1MuBMTEtaPatternLut.h"
0046 #include "CondFormats/DataRecord/interface/L1MuDTQualPatternLutRcd.h"
0047 #include "CondFormats/L1TObjects/interface/L1MuDTTFMasks.h"
0048 #include "CondFormats/DataRecord/interface/L1MuDTTFMasksRcd.h"
0049 
0050 #include "DataFormats/L1TMuon/interface/L1MuBMTrack.h"
0051 #include "DataFormats/L1TMuon/interface/L1MuBMTrackSegEta.h"
0052 #include "DataFormats/L1TMuon/interface/BMTF/L1MuBMSecProcId.h"
0053 
0054 using namespace std;
0055 
0056 // --------------------------------
0057 //       class L1MuBMEtaProcessor
0058 //---------------------------------
0059 
0060 //----------------
0061 // Constructors --
0062 //----------------
0063 
0064 L1MuBMEtaProcessor::L1MuBMEtaProcessor(L1MuBMTrackFinder& tf, int id, edm::ConsumesCollector&& iC)
0065     : m_tf(tf),
0066       m_epid(id),
0067       m_foundPattern(0),
0068       m_tseta(15),
0069       m_bmtfParamsToken(iC.esConsumes()),
0070       m_DTDigiToken(iC.consumes<L1MuDTChambThContainer>(tf.config().getBMThetaDigiInputTag())) {
0071   m_tseta.reserve(15);
0072 }
0073 
0074 //--------------
0075 // Destructor --
0076 //--------------
0077 
0078 L1MuBMEtaProcessor::~L1MuBMEtaProcessor() {}
0079 
0080 //--------------
0081 // Operations --
0082 //--------------
0083 
0084 //
0085 // run Eta Processor
0086 //
0087 void L1MuBMEtaProcessor::run(int bx, const edm::Event& e, const edm::EventSetup& c) {
0088   auto const& params = c.getData(m_bmtfParamsToken);
0089   if (m_tf.config().getEtaTF()) {
0090     receiveData(bx, e, params);
0091     runEtaTrackFinder(params);
0092   }
0093 
0094   receiveAddresses();
0095   runEtaMatchingUnit(params);
0096 
0097   assign();
0098 }
0099 
0100 //
0101 // reset Eta Processor
0102 //
0103 void L1MuBMEtaProcessor::reset() {
0104   vector<const L1MuBMTrackSegEta*>::iterator iter = m_tseta.begin();
0105   while (iter != m_tseta.end()) {
0106     if (*iter) {
0107       delete *iter;
0108       *iter = nullptr;
0109     }
0110     iter++;
0111   }
0112 
0113   m_tseta.clear();
0114 
0115   for (int i = 0; i < 12; i++) {
0116     m_eta[i] = 99;
0117     m_fine[i] = false;
0118     m_pattern[i] = 0;
0119     m_address[i] = 0;
0120     m_TrackCand[i] = nullptr;
0121     m_TracKCand[i] = nullptr;
0122   }
0123 
0124   m_foundPattern.clear();
0125 
0126   m_mask = true;
0127 }
0128 
0129 //
0130 // print track candidates found in Eta Processor
0131 //
0132 void L1MuBMEtaProcessor::print() const {
0133   bool empty1 = true;
0134   for (int i = 0; i < 15; i++) {
0135     empty1 &= (m_tseta[i] == nullptr || m_tseta[i]->empty());
0136   }
0137 
0138   bool empty2 = true;
0139   for (int i = 0; i < 12; i++) {
0140     empty2 &= (m_address[i] == 0);
0141   }
0142 
0143   if (!empty1 || !empty2) {
0144     cout << "Eta processor " << m_epid << " : " << endl;
0145 
0146     // print local pattern
0147     if (!empty1) {
0148       cout << "Local pattern : " << endl;
0149       for (int i = 0; i < 15; i++) {
0150         if ((i + 5) % 5 == 0)
0151           cout << "station " << m_tseta[i]->station() << " : ";
0152         bitset<7> pos(m_tseta[i]->position());
0153         bitset<7> qua(m_tseta[i]->quality());
0154         for (int j = 6; j >= 0; j--) {
0155           cout << pos[j] + qua[j];
0156         }
0157         cout << " ";
0158         if ((i + 1) % 5 == 0)
0159           cout << endl;
0160       }
0161       cout << "Found patterns :" << endl;
0162       vector<int>::const_iterator iter;
0163       for (iter = m_foundPattern.begin(); iter != m_foundPattern.end(); iter++) {
0164         const L1MuDTEtaPattern p = theEtaPatternLUT.getPattern(*iter);
0165         int qualitycode = p.quality();
0166         cout << "ID = " << setw(4) << p.id() << "  "
0167              << "eta = " << setw(3) << p.eta() << "  "
0168              << "quality = " << setw(2) << qualitycode << " (" << quality(qualitycode, 1) << " "
0169              << quality(qualitycode, 2) << " " << quality(qualitycode, 3) << ")";
0170         for (int i = 0; i < 12; i++) {
0171           if (m_pattern[i] == p.id())
0172             cout << " <--";
0173         }
0174         cout << endl;
0175       }
0176     }
0177 
0178     cout << "Received addresses : " << endl;
0179     for (int i = 0; i < 12; i++)
0180       cout << setw(3) << m_address[i] << " ";
0181     cout << endl;
0182 
0183     if (!empty1) {
0184       cout << "Matched patterns : " << endl;
0185       for (int i = 0; i < 12; i++) {
0186         if (m_fine[i]) {
0187           const L1MuDTEtaPattern p = theEtaPatternLUT.getPattern(m_pattern[i]);
0188           int fineeta = p.eta();
0189           int coarseeta = theQualPatternLUT.getCoarseEta(i / 2 + 1, m_address[i]);
0190           cout << "Index = " << setw(2) << i << ", "
0191                << "address = " << setw(2) << m_address[i] << " --> "
0192                << "pattern = " << setw(4) << m_pattern[i] << " "
0193                << "eta (coarse) = " << setw(3) << coarseeta << " "
0194                << "eta (fine) = " << setw(3) << fineeta << " "
0195                << "quality = " << setw(2) << p.quality() << endl;
0196         }
0197       }
0198     }
0199 
0200     cout << "Eta values and fine bits : " << endl;
0201     for (int i = 0; i < 12; i++)
0202       cout << setw(3) << m_eta[i] << " ";
0203     cout << endl;
0204     for (int i = 0; i < 12; i++)
0205       cout << setw(3) << m_fine[i] << " ";
0206     cout << endl;
0207   }
0208 }
0209 
0210 //
0211 // receive data ( 15*3 BBMX eta trigger primitives )
0212 //
0213 void L1MuBMEtaProcessor::receiveData(int bx, const edm::Event& e, const L1TMuonBarrelParams& bmtfParams) {
0214   msks = bmtfParams.l1mudttfmasks;
0215   theEtaPatternLUT.m_lut = bmtfParams.lutparams_.eta_lut_;  //l1mudttfetaplut;  // ETF look-up table
0216   theQualPatternLUT.m_lut = bmtfParams.lutparams_.qp_lut_;  //l1mudttfqualplut; // EMU look-up tables
0217 
0218   edm::Handle<L1MuDTChambThContainer> dttrig;
0219   e.getByToken(m_DTDigiToken, dttrig);
0220 
0221   // const int bx_offset = dttrig->correctBX();
0222   int bx_offset = 0;
0223   bx = bx + bx_offset;
0224 
0225   //
0226   // get 5*3 eta track segments
0227   //
0228   int sector = m_epid;
0229   for (int stat = 1; stat <= 3; stat++) {
0230     for (int wheel = -2; wheel <= 2; wheel++) {
0231       L1MuDTChambThDigi const* tseta = dttrig->chThetaSegm(wheel, stat, sector, bx);
0232       bitset<7> pos;
0233       bitset<7> qual;
0234 
0235       int lwheel = wheel + 1;
0236       if (wheel < 0)
0237         lwheel = wheel - 1;
0238 
0239       bool masked = false;
0240       if (stat == 1)
0241         masked = msks.get_etsoc_chdis_st1(lwheel, sector);
0242       if (stat == 2)
0243         masked = msks.get_etsoc_chdis_st2(lwheel, sector);
0244       if (stat == 3)
0245         masked = msks.get_etsoc_chdis_st3(lwheel, sector);
0246 
0247       if (!masked)
0248         m_mask = false;
0249 
0250       if (tseta && !masked) {
0251         if (wheel == -2 || wheel == -1 ||
0252             (wheel == 0 && (sector == 0 || sector == 3 || sector == 4 || sector == 7 || sector == 8 || sector == 11))) {
0253           for (int i = 0; i < 7; i++) {
0254             if (tseta->position(i))
0255               pos.set(6 - i);
0256             if (tseta->quality(i))
0257               qual.set(6 - i);
0258           }
0259         } else {
0260           for (int i = 0; i < 7; i++) {
0261             if (tseta->position(i))
0262               pos.set(i);
0263             if (tseta->quality(i))
0264               qual.set(i);
0265           }
0266         }
0267       }
0268 
0269       const L1MuBMTrackSegEta* tmpts =
0270           new L1MuBMTrackSegEta(wheel, sector, stat, pos.to_ulong(), qual.to_ulong(), bx - bx_offset);
0271       m_tseta.push_back(tmpts);
0272     }
0273   }
0274 }
0275 
0276 //
0277 // receive track addresses from 6 Sector Processors
0278 //
0279 void L1MuBMEtaProcessor::receiveAddresses() {
0280   // get track address code of all track segments
0281   // 6*2 times 5 bits; valid range [1-22]
0282 
0283   int sector = m_epid;
0284 
0285   int i = 0;
0286   for (int wheel = -3; wheel <= 3; wheel++) {
0287     if (wheel == 0)
0288       continue;
0289     L1MuBMSecProcId tmpspid(wheel, sector);
0290     for (int number = 0; number < 2; number++) {
0291       L1MuBMTrack& cand = m_tf.sp(tmpspid)->track(number);
0292       L1MuBMTrack& canD = m_tf.sp(tmpspid)->tracK(number);
0293       m_address[i] = cand.address().trackAddressCode();
0294       if (!cand.empty()) {
0295         m_TrackCand[i] = &cand;
0296         m_TracKCand[i] = &canD;
0297       }
0298       i++;
0299     }
0300   }
0301 }
0302 
0303 //
0304 // run Eta Track Finder (ETF)
0305 //
0306 void L1MuBMEtaProcessor::runEtaTrackFinder(const L1TMuonBarrelParams& bmtfParams) {
0307   theEtaPatternLUT.m_lut = bmtfParams.lutparams_.eta_lut_;  //l1mudttfetaplut;  // ETF look-up table
0308 
0309   // check if there are any data
0310   bool empty = true;
0311   for (int i = 0; i < 15; i++) {
0312     empty &= m_tseta[i]->empty();
0313   }
0314   if (empty)
0315     return;
0316 
0317   // Pattern comparator:
0318   // loop over all patterns and compare with local chamber pattern
0319   // result : list of valid pattern IDs ( m_foundPattern )
0320   L1MuBMTEtaPatternLut::ETFLut_iter it = theEtaPatternLUT.begin();
0321   while (it != theEtaPatternLUT.end()) {
0322     const L1MuDTEtaPattern pattern = (*it).second;
0323     int qualitycode = pattern.quality();
0324 
0325     bool good = true;
0326 
0327     for (int station = 0; station < 3; station++) {
0328       int q = quality(qualitycode, station + 1);
0329       int wheel = pattern.wheel(station + 1);
0330       int bin = pattern.position(station + 1);
0331       if (bin == 0)
0332         continue;
0333       bitset<7> pos = m_tseta[wheel + 2 + station * 5]->position();
0334       bitset<7> qual = m_tseta[wheel + 2 + station * 5]->quality();
0335       // compare position
0336       good &= pos.test(bin - 1);
0337       // compare quality
0338       if (q == 2)
0339         good &= qual.test(bin - 1);
0340     }
0341 
0342     if (good)
0343       m_foundPattern.push_back(pattern.id());
0344 
0345     it++;
0346   }
0347 }
0348 
0349 //
0350 // run Eta Matching Unit (EMU)
0351 //
0352 void L1MuBMEtaProcessor::runEtaMatchingUnit(const L1TMuonBarrelParams& bmtfParams) {
0353   theQualPatternLUT.m_lut = bmtfParams.lutparams_.qp_lut_;  //l1mudttfqualplut; // EMU look-up tables
0354 
0355   // loop over all addresses
0356   for (int i = 0; i < 12; i++) {
0357     int adr = m_address[i];
0358     if (adr == 0)
0359       continue;
0360     int sp = i / 2 + 1;  //sector processor [1,6]
0361 
0362     // assign coarse eta value
0363     if (!m_mask)
0364       m_eta[i] = theQualPatternLUT.getCoarseEta(sp, adr);
0365 
0366     if (m_foundPattern.empty())
0367       continue;
0368 
0369     // get list of qualified patterns ordered by quality
0370     // and compare with found patterns
0371     const vector<short>& qualifiedPatterns = theQualPatternLUT.getQualifiedPatterns(sp, adr);
0372     vector<short>::const_iterator iter;
0373     vector<int>::const_iterator f_iter;
0374     for (iter = qualifiedPatterns.begin(); iter != qualifiedPatterns.end(); iter++) {
0375       f_iter = find(m_foundPattern.begin(), m_foundPattern.end(), (*iter));
0376       // found
0377       if (f_iter != m_foundPattern.end()) {
0378         const L1MuDTEtaPattern p = theEtaPatternLUT.getPattern(*f_iter);
0379         // assign fine eta value
0380         m_fine[i] = true;
0381         m_eta[i] = p.eta();  // improved eta
0382         ;
0383         m_pattern[i] = (*f_iter);
0384         break;
0385       }
0386     }
0387   }
0388 
0389   // if both tracks from one sector processor deliver the same track address
0390   // both tracks get only a coarse eta value!
0391 
0392   // loop over sector processors
0393   for (int i = 0; i < 6; i++) {
0394     int idx1 = 2 * i;
0395     int idx2 = 2 * i + 1;
0396     int adr1 = m_address[idx1];
0397     int adr2 = m_address[idx2];
0398     if (adr1 == 0 || adr2 == 0)
0399       continue;
0400     if (adr1 == adr2 && !m_mask) {
0401       // both tracks get coarse (default) eta value
0402       m_eta[idx1] = theQualPatternLUT.getCoarseEta(i + 1, adr1);
0403       m_pattern[idx1] = 0;
0404       m_fine[idx1] = false;
0405       m_eta[idx2] = theQualPatternLUT.getCoarseEta(i + 1, adr2);
0406 
0407       m_pattern[idx2] = 0;
0408       m_fine[idx2] = false;
0409     }
0410   }
0411 }
0412 
0413 //
0414 // assign values to track candidates
0415 //
0416 void L1MuBMEtaProcessor::assign() {
0417   for (int i = 0; i < 12; i++) {
0418     if (m_TrackCand[i]) {
0419       if (m_eta[i] != 99) {
0420         m_TrackCand[i]->setEta(m_eta[i]);
0421         m_TracKCand[i]->setEta(m_eta[i]);
0422       } else {
0423         //        if ( i/2 != 2 ) cerr << "L1MuBMEtaProcessor: assign invalid eta" << " " << m_address[i] << endl;
0424       }
0425       if (m_fine[i]) {
0426         m_TrackCand[i]->setFineEtaBit();
0427         m_TracKCand[i]->setFineEtaBit();
0428         // find all contributing track segments
0429         const L1MuDTEtaPattern p = theEtaPatternLUT.getPattern(m_pattern[i]);
0430         vector<const L1MuBMTrackSegEta*> TSeta;
0431         const L1MuBMTrackSegEta* ts = nullptr;
0432         for (int stat = 0; stat < 3; stat++) {
0433           int wh = p.wheel(stat + 1);
0434           int pos = p.position(stat + 1);
0435           if (pos == 0)
0436             continue;
0437           ts = m_tseta[wh + 2 + stat * 5];
0438           TSeta.push_back(ts);
0439         }
0440         m_TrackCand[i]->setTSeta(TSeta);
0441         m_TracKCand[i]->setTSeta(TSeta);
0442       }
0443     }
0444   }
0445 }
0446 
0447 //
0448 // get quality:  id [0,26], stat [1,3]
0449 //
0450 int L1MuBMEtaProcessor::quality(int id, int stat) {
0451   // quality codes as defined in CMS Note 2001/027
0452   // This QualityPatterns are used to have a defined Quality-Identifier for
0453   // all possible found tracks.
0454   // Therefore three integer numbers ( Station 1, 2, 3 from left to right )
0455   // can have numbers like 0, 1 or 2
0456   //    0 ... no hit is given
0457   //    1 ... a LTRG is given
0458   //    2 ... a HTRG is given
0459 
0460   const int qualcode[27][3] = {{0, 0, 0}, {1, 0, 0}, {0, 1, 0}, {0, 0, 1}, {2, 0, 0}, {0, 2, 0}, {0, 0, 2},
0461                                {1, 1, 0}, {1, 0, 1}, {0, 1, 1}, {2, 1, 0}, {1, 2, 0}, {2, 0, 1}, {1, 0, 2},
0462                                {0, 2, 1}, {0, 1, 2}, {2, 2, 0}, {2, 0, 2}, {0, 2, 2}, {1, 1, 1}, {2, 1, 1},
0463                                {1, 2, 1}, {1, 1, 2}, {2, 2, 1}, {2, 1, 2}, {1, 2, 2}, {2, 2, 2}};
0464 
0465   return qualcode[id][stat - 1];
0466 }