Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:20:06

0001 #ifndef L1_EMUL_BIAS_H
0002 #define L1_EMUL_BIAS_H
0003 
0004 /*\class L1EmulBias
0005  *\description produces modified emulator data to mimmic hw problems
0006  *\usage l1 monitoring software validation
0007  *\author Nuno Leonardo (CERN)
0008  *\date 07.07
0009  */
0010 
0011 // system includes
0012 #include <memory>
0013 #include <string>
0014 #include <iostream>
0015 #include <fstream>
0016 #include <iomanip>
0017 #include <vector>
0018 #include <algorithm>
0019 
0020 // common includes
0021 #include "FWCore/Framework/interface/Frameworkfwd.h"
0022 #include "FWCore/Framework/interface/global/EDProducer.h"
0023 #include "FWCore/Framework/interface/Event.h"
0024 #include "FWCore/Framework/interface/MakerMacros.h"
0025 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0026 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0027 
0028 // l1 dataformats, d|e record includes
0029 #include "DataFormats/HcalDigi/interface/HcalTriggerPrimitiveDigi.h"
0030 #include "DataFormats/L1CSCTrackFinder/interface/L1Track.h"
0031 #include "L1Trigger/HardwareValidation/interface/DEtrait.h"
0032 
0033 // random generation
0034 #include "CLHEP/Random/RandomEngine.h"
0035 #include "CLHEP/Random/RandGaussQ.h"
0036 
0037 class L1EmulBias : public edm::global::EDProducer<> {
0038 public:
0039   explicit L1EmulBias(const edm::ParameterSet&);
0040   ~L1EmulBias() override;
0041 
0042 protected:
0043   //virtual void beginRun(edm::Run&, const edm::EventSetup&);
0044   void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0045 
0046 public:
0047   template <class T>
0048   void ModifyCollection(std::unique_ptr<T>& data, const edm::Handle<T> emul, CLHEP::HepRandomEngine*) const;
0049 
0050 private:
0051   int verbose_;
0052   int verbose() const { return verbose_; }
0053   edm::InputTag m_DEsource[dedefs::DEnsys][2];
0054   bool m_doSys[dedefs::DEnsys];
0055   std::string instName[dedefs::DEnsys][5];
0056 };
0057 
0058 /* Notes:
0059    .by default data is make identical to emul
0060    .biasing is to be implemented via specialization
0061    .bias may be defined for each data type, eg data word bit-shifting
0062    .keep template function specialization in header file
0063 */
0064 
0065 template <class T>
0066 void L1EmulBias::ModifyCollection(std::unique_ptr<T>& data, const edm::Handle<T> emul, CLHEP::HepRandomEngine*) const {
0067   data = (std::unique_ptr<T>)(const_cast<T*>(emul.product()));
0068 }
0069 
0070 template <>
0071 inline void L1EmulBias::ModifyCollection(std::unique_ptr<EcalTrigPrimDigiCollection>& data,
0072                                          const edm::Handle<EcalTrigPrimDigiCollection> emul,
0073                                          CLHEP::HepRandomEngine*) const {
0074   typedef EcalTrigPrimDigiCollection::const_iterator col_cit;
0075   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0076     EcalTriggerPrimitiveDigi col(*it);
0077     int iphi = it->id().iphi();
0078     bool reset = (iphi > 18 && iphi < 39);  //remove few supermodules
0079     for (int s = 0; s < 5; s++) {
0080       uint16_t sample = it->sample(s).raw();
0081       if (sample == 0)
0082         continue;
0083       uint16_t tmp = reset ? 0 : sample;
0084       if (reset)
0085         tmp = sample >> 1;
0086       col.setSampleValue(s, tmp);
0087       if (verbose() && sample != 0)
0088         std::cout << "[emulbias] etp " << *it << "\t sample: " << s << "  " << std::hex << sample << " -> "
0089                   << col.sample(s).raw() << std::dec << std::endl;
0090     }
0091     data->push_back(col);
0092   }
0093 }
0094 
0095 template <>
0096 inline void L1EmulBias::ModifyCollection(std::unique_ptr<HcalTrigPrimDigiCollection>& data,
0097                                          const edm::Handle<HcalTrigPrimDigiCollection> emul,
0098                                          CLHEP::HepRandomEngine*) const {
0099   typedef HcalTrigPrimDigiCollection::const_iterator col_cit;
0100   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0101     HcalTriggerPrimitiveDigi col(*it);
0102     int iphi = it->id().iphi();
0103     bool reset = (iphi > 18 && iphi < 27);  //remove few supermodules
0104     for (int s = 0; s < 5; s++) {
0105       uint16_t sample = it->sample(s).raw();
0106       if (sample == 0)
0107         continue;
0108       uint16_t tmp = reset ? 0 : sample;
0109       if (reset)
0110         tmp = sample >> 1;
0111       col.setSample(s, tmp);
0112     }
0113     data->push_back(col);
0114   }
0115 }
0116 
0117 template <>
0118 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1CaloEmCollection>& data,
0119                                          const edm::Handle<L1CaloEmCollection> emul,
0120                                          CLHEP::HepRandomEngine* engine) const {
0121   typedef L1CaloEmCollection::const_iterator col_cit;
0122   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0123     unsigned crate = it->rctCrate();
0124     unsigned raw = it->raw();
0125     bool iso = it->isolated();
0126     unsigned rdata = raw;
0127     if (crate < 4 * engine->flat())
0128       rdata = raw >> 1;
0129     L1CaloEmCand cand(rdata, crate, iso, it->index(), it->bx(), false);
0130     data->push_back(cand);
0131   }
0132   //L1CaloEmCand(uint16_t data, unsigned crate, bool iso);
0133   //L1CaloEmCand(uint16_t data, unsigned crate, bool iso, uint16_t index, int16_t bx, bool dummy);
0134 }
0135 
0136 template <>
0137 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1CaloRegionCollection>& data,
0138                                          const edm::Handle<L1CaloRegionCollection> emul,
0139                                          CLHEP::HepRandomEngine* engine) const {
0140   typedef L1CaloRegionCollection::const_iterator col_cit;
0141   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0142     unsigned crate = it->rctCrate();
0143     unsigned raw = it->et();
0144     uint16_t rdata = raw;
0145     if (crate < 4 * engine->flat())
0146       rdata = raw >> 1;
0147     L1CaloRegion cand(rdata, it->gctEta(), it->gctPhi(), it->bx());
0148     data->push_back(cand);
0149   }
0150   //L1CaloRegion(uint16_t data, unsigned ieta, unsigned iphi, int16_t bx);
0151   //Note: raw data accessor missing in dataformats!
0152 }
0153 
0154 template <>
0155 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1GctEmCandCollection>& data,
0156                                          const edm::Handle<L1GctEmCandCollection> emul,
0157                                          CLHEP::HepRandomEngine* engine) const {
0158   typedef L1GctEmCandCollection::const_iterator col_cit;
0159   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0160     unsigned raw = it->raw();
0161     uint16_t rdata = raw;
0162     if (it->phiIndex() < 4 * engine->flat())  //0-17
0163       rdata = raw >> 1;
0164     L1GctEmCand cand(rdata, it->isolated());
0165     data->push_back(cand);
0166   }
0167   //etaIndex(), etaSign() : -6 to -0, +0 to +6
0168   //L1GctEmCand(uint16_t data, bool iso);
0169 }
0170 
0171 template <>
0172 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1GctJetCandCollection>& data,
0173                                          const edm::Handle<L1GctJetCandCollection> emul,
0174                                          CLHEP::HepRandomEngine* engine) const {
0175   typedef L1GctJetCandCollection::const_iterator col_cit;
0176   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0177     unsigned raw = it->raw();
0178     uint16_t rdata = raw;
0179     if (it->phiIndex() < 4 * engine->flat())  //0-17
0180       rdata = raw >> 1;
0181     L1GctJetCand cand(rdata, it->isTau(), it->isForward());
0182     data->push_back(cand);
0183   }
0184   //L1GctJetCand(uint16_t data, bool isTau, bool isFor);
0185 }
0186 
0187 template <>
0188 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1MuRegionalCandCollection>& data,
0189                                          const edm::Handle<L1MuRegionalCandCollection> emul,
0190                                          CLHEP::HepRandomEngine* engine) const {
0191   typedef L1MuRegionalCandCollection::const_iterator col_cit;
0192   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0193     L1MuRegionalCand cand(*it);
0194     //unsigned raw = it->getDataWord();
0195     unsigned phi = it->phi_packed();
0196     if (phi > 90 && phi < 110)
0197       cand.setPtPacked((it->pt_packed()) >> 1);
0198     //raw = (raw>>2);
0199     //L1MuRegionalCand cand(raw);
0200     //cand.setType(it->type_idx());
0201     data->push_back(cand);
0202   }
0203   /* few alternatives...
0204   unsigned pt= it->pt_packed(); //0..31
0205   unsigned int qua = it->quality(); //0..7
0206   if(qua<4){cand.setPtPacked((pt>>2)&0x1f);cand.setQualityPacked((qua<<1)&0x07);}
0207   double rnd = CLHEP::RandGaussQ::shoot(engine);
0208   if(rnd>0.7) {
0209     raw_=(raw>>1);
0210     cand.setDataWord(raw_);
0211   } else if (rnd>0.3) {
0212     pt_ *= (int)(1+0.3*engine->flat());
0213     cand.setPtPacked(pt_);
0214   } else 
0215     cand.reset();
0216   unsigned raw = it->getDataWord();
0217   if(2.5<fabs(it->phiValue())<3.0)
0218     rdata = raw>>1;
0219   L1MuRegionalCand cand(rdata,it->bx());    
0220   */
0221   //L1MuRegionalCand(unsigned dataword = 0, int bx = 0);
0222   //L1MuRegionalCand(unsigned type_idx, unsigned phi, unsigned eta, unsigned pt, unsigned charge, unsigned ch_valid, unsigned finehalo, unsigned quality, int bx);
0223 }
0224 
0225 template <>
0226 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1MuDTTrackContainer>& data,
0227                                          const edm::Handle<L1MuDTTrackContainer> emul,
0228                                          CLHEP::HepRandomEngine* engine) const {
0229   typedef std::vector<L1MuDTTrackCand> TrackContainer;
0230   typedef TrackContainer::const_iterator col_cit;
0231   TrackContainer const* tracks_in = emul->getContainer();
0232   TrackContainer tracks;
0233   for (col_cit it = tracks_in->begin(); it != tracks_in->end(); it++) {
0234     L1MuDTTrackCand cand(*it);
0235     cand.setType(it->type_idx());
0236     unsigned pt = it->pt_packed();  //0..31
0237     unsigned qua = it->quality();   //0..7
0238     if (qua < 4) {
0239       cand.setPtPacked((pt >> 2) & 0x1f);
0240       cand.setQualityPacked((qua << 1) & 0x07);
0241     }
0242     tracks.push_back(cand);
0243   }
0244   data->setContainer(tracks);
0245   /*   few alternatives...
0246   unsigned phip = it->phi_packed();
0247   unsigned raw = it->getDataWord();
0248   uint16_t rdata = raw;
0249   if(2.5<fabs(it->phiValue())<3.0)
0250     rdata = raw>>1;
0251   L1MuRegionalCand cand(rdata,it->bx());    
0252   double rnd    = engine->flat();
0253   */
0254 }
0255 
0256 template <>
0257 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1MuDTChambPhContainer>& data,
0258                                          const edm::Handle<L1MuDTChambPhContainer> emul,
0259                                          CLHEP::HepRandomEngine* engine) const {
0260   typedef std::vector<L1MuDTChambPhDigi> Phi_Container;
0261   typedef Phi_Container::const_iterator col_it;
0262   Phi_Container const* tracks_in = emul->getContainer();
0263   Phi_Container tracks(tracks_in->size());
0264   int uqua;
0265   for (col_it it = tracks_in->begin(); it != tracks_in->end(); it++) {
0266     uqua = it->code();  // (int)(10*engine->flat());
0267     uqua = (uqua < 2 ? uqua + 1 : uqua);
0268     L1MuDTChambPhDigi cand(
0269         it->bxNum(), it->whNum(), it->scNum(), it->stNum(), it->phi(), it->phiB(), uqua, it->Ts2Tag(), it->BxCnt());
0270     tracks.push_back(cand);
0271   }
0272   data->setContainer(tracks);
0273 }
0274 
0275 template <>
0276 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1MuDTChambThContainer>& data,
0277                                          const edm::Handle<L1MuDTChambThContainer> emul,
0278                                          CLHEP::HepRandomEngine*) const {
0279   typedef std::vector<L1MuDTChambThDigi> Thi_Container;
0280   typedef Thi_Container::const_iterator col_cit;
0281   Thi_Container const* tracks_in = emul->getContainer();
0282   Thi_Container tracks(tracks_in->size());
0283   int uos[7], uqa[7];
0284   for (col_cit it = tracks_in->begin(); it != tracks_in->end(); it++) {
0285     for (int j = 0; j < 7; j++) {
0286       uos[j] = (it->position(j) ? 0 : 1);
0287       uqa[j] = (it->quality(j) ? 0 : 1);
0288     }
0289     int stnum = it->stNum();
0290     stnum = (stnum > 2 ? stnum - 1 : stnum);
0291     L1MuDTChambThDigi cand(it->bxNum(), it->whNum(), it->scNum(), stnum, uos, uqa);
0292     tracks.push_back(cand);
0293   }
0294   data->setContainer(tracks);
0295 }
0296 
0297 template <>
0298 inline void L1EmulBias::ModifyCollection(std::unique_ptr<LTCDigiCollection>& data,
0299                                          const edm::Handle<LTCDigiCollection> emul,
0300                                          CLHEP::HepRandomEngine*) const {
0301   typedef std::vector<LTCDigi>::const_iterator col_cit;
0302   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0303     data->push_back(*it);
0304     //note: raw data accessor missing in dataformats!
0305     //data->push_back(LTCDigi(it->data()>>1));
0306   }
0307 }
0308 
0309 template <>
0310 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1MuGMTCandCollection>& data,
0311                                          const edm::Handle<L1MuGMTCandCollection> emul,
0312                                          CLHEP::HepRandomEngine*) const {
0313   //typedef std::vector<L1MuGMTCand>          L1MuGMTCandCollection;
0314   typedef std::vector<L1MuGMTCand>::const_iterator col_cit;
0315   for (col_cit it = emul->begin(); it != emul->end(); it++) {
0316     float phiv = it->phiValue();
0317     unsigned dword = it->getDataWord();
0318     if (phiv > 2. && phiv < 4.)
0319       dword = dword >> 2;
0320     L1MuGMTCand cand(dword, it->bx());
0321     data->push_back(cand);
0322     //cand.setPtPacked(cand.ptIndex()>>1);
0323     //data->push_back(L1MuGMTCand((it->getDataWord()>>1),it->bx()));
0324   }
0325 }
0326 
0327 template <>
0328 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1MuGMTReadoutCollection>& data,
0329                                          const edm::Handle<L1MuGMTReadoutCollection> emul,
0330                                          CLHEP::HepRandomEngine*) const {
0331   typedef std::vector<L1MuGMTReadoutRecord>::const_iterator col_cit;
0332   std::vector<L1MuGMTReadoutRecord> col = emul->getRecords();
0333   for (col_cit it = col.begin(); it != col.end(); it++) {
0334     L1MuGMTReadoutRecord rec(it->getBxInEvent());
0335     rec.setBxNr(it->getBxNr());
0336     rec.setEvNr(it->getEvNr());
0337     rec.setBCERR(it->getBCERR());
0338 
0339     std::unique_ptr<L1MuRegionalCandCollection> new_dttf(new L1MuRegionalCandCollection);
0340     std::unique_ptr<L1MuRegionalCandCollection> new_rpcb(new L1MuRegionalCandCollection);
0341     std::unique_ptr<L1MuRegionalCandCollection> new_csc(new L1MuRegionalCandCollection);
0342     std::unique_ptr<L1MuRegionalCandCollection> new_rpcf(new L1MuRegionalCandCollection);
0343 
0344     L1MuRegionalCandCollection old_dttf = it->getDTBXCands();
0345     L1MuRegionalCandCollection old_rpcb = it->getBrlRPCCands();
0346     L1MuRegionalCandCollection old_csc = it->getCSCCands();
0347     L1MuRegionalCandCollection old_rpcf = it->getFwdRPCCands();
0348 
0349     typedef L1MuRegionalCandCollection::const_iterator ait;
0350     for (ait it = old_dttf.begin(); it != old_dttf.end(); it++) {
0351       L1MuRegionalCand cand(*it);
0352       if (it->quality() < 4)
0353         cand.setPtPacked((it->pt_packed() >> 2) & 0x1f);
0354       cand.setType(it->type_idx());
0355       new_dttf->push_back(cand);
0356     }
0357     for (ait it = old_rpcb.begin(); it != old_rpcb.end(); it++) {
0358       L1MuRegionalCand cand(*it);
0359       if (it->quality() < 4)
0360         cand.setPtPacked((it->pt_packed() >> 2) & 0x1f);
0361       cand.setType(it->type_idx());
0362       new_rpcb->push_back(cand);
0363     }
0364     for (ait it = old_csc.begin(); it != old_csc.end(); it++) {
0365       L1MuRegionalCand cand(*it);
0366       if (it->quality() < 4)
0367         cand.setPtPacked((it->pt_packed() >> 2) & 0x1f);
0368       cand.setType(it->type_idx());
0369       new_csc->push_back(cand);
0370     }
0371     for (ait it = old_rpcf.begin(); it != old_rpcf.end(); it++) {
0372       L1MuRegionalCand cand(*it);
0373       if (it->quality() < 4)
0374         cand.setPtPacked((it->pt_packed() >> 2) & 0x1f);
0375       cand.setType(it->type_idx());
0376       new_rpcf->push_back(cand);
0377     }
0378 
0379     for (unsigned i = 0; i < old_dttf.size(); i++)
0380       rec.setInputCand(i, new_dttf->at(i));  //dt  : 0..3
0381     for (unsigned i = 0; i < old_rpcb.size(); i++)
0382       rec.setInputCand(i + 4, new_rpcb->at(i));  //rpcb: 4..7
0383     for (unsigned i = 0; i < old_csc.size(); i++)
0384       rec.setInputCand(i + 8, new_csc->at(i));  //csc : 8..11
0385     for (unsigned i = 0; i < old_rpcf.size(); i++)
0386       rec.setInputCand(i + 12, new_rpcf->at(i));  //rpcf:12..15
0387 
0388     data->addRecord(rec);
0389   }
0390   //void addRecord(L1MuGMTReadoutRecord const& rec) {
0391 }
0392 
0393 template <>
0394 inline void L1EmulBias::ModifyCollection(std::unique_ptr<CSCCorrelatedLCTDigiCollection>& data,
0395                                          const edm::Handle<CSCCorrelatedLCTDigiCollection> emul,
0396                                          CLHEP::HepRandomEngine*) const {
0397   //typedef MuonDigiCollection<CSCDetId,CSCCorrelatedLCTDigi> CSCCorrelatedLCTDigiCollection;
0398   typedef CSCCorrelatedLCTDigiCollection::DigiRangeIterator mapIt;  //map iterator
0399   typedef CSCCorrelatedLCTDigiCollection::const_iterator vecIt;     //vec iterator
0400   //loop over data (map<idx,vec_digi>)
0401   for (mapIt mit = emul->begin(); mit != emul->end(); mit++) {
0402     //get detector index
0403     CSCDetId did = (*mit).first;
0404     //get vec_digi range(pair)  corresponding to idx of map
0405     //CSCCorrelatedLCTDigiCollection::Range ctpRange = emul->get(did)
0406     //loop over digi vector (ie between begin and end pointers in range)
0407     //for (vecIt vit = ctpRange.first; vit != ctpRange.second; vit++) {
0408     for (vecIt vit = emul->get((*mit).first).first; vit != emul->get((*mit).first).second; vit++) {
0409       ///modify digi
0410       CSCCorrelatedLCTDigi dg = *vit;
0411       //dg.clear;
0412       uint16_t tn = dg.getTrknmb();
0413       if (tn == 2)
0414         tn--;
0415       dg.setTrknmb(tn);
0416       //dg.setTrknmb   (dg.getTrknmb   ());
0417       //dg.setMPCLink  (dg.getMPCLink  ());
0418       //dg.setWireGroup(dg.getWireGroup());
0419       ///append digi
0420       data->insertDigi(did, dg);
0421     }
0422   }
0423 }
0424 
0425 template <>
0426 inline void L1EmulBias::ModifyCollection(std::unique_ptr<L1CSCTrackCollection>& data,
0427                                          const edm::Handle<L1CSCTrackCollection> emul,
0428                                          CLHEP::HepRandomEngine*) const {
0429   typedef L1CSCTrackCollection::const_iterator col_cit;
0430   //typedef std::vector<L1CSCTrack> L1CSCTrackCollection;
0431   //typedef std::pair<csc::L1Track,CSCCorrelatedLCTDigiCollection> L1CSCTrack;
0432   //typedef MuonDigiCollection<CSCDetId,CSCCorrelatedLCTDigi> CSCCorrelatedLCTDigiCollection;
0433   typedef CSCCorrelatedLCTDigiCollection::DigiRangeIterator mapIt;  //map iterator
0434   typedef CSCCorrelatedLCTDigiCollection::const_iterator vecIt;     //vec iterator
0435   CSCCorrelatedLCTDigiCollection_ ctf_trk_data_v, ctf_trk_emul_v;   //vector
0436   //loop over csc-tracks (ie pairs<l1track,digi_vec>)
0437   for (col_cit tcit = emul->begin(); tcit != emul->end(); tcit++) {
0438     csc::L1Track l1trk = tcit->first;
0439     if (l1trk.quality() < 4)
0440       l1trk.setPtPacked((l1trk.pt_packed() >> 2) & 0x1f);
0441     l1trk.setType(l1trk.type_idx());
0442     //L1MuRegionalCand reg(tcit->first.getDataWord(), tcit->first.bx());
0443     std::unique_ptr<CSCCorrelatedLCTDigiCollection> dgcoll(new CSCCorrelatedLCTDigiCollection);
0444     CSCCorrelatedLCTDigiCollection ldc = tcit->second;  //muondigicollection=map
0445     //get the lct-digi-collection (ie muon-digi-collection)
0446     //loop over data (map<idx,vec_digi>)
0447     for (mapIt mit = ldc.begin(); mit != ldc.end(); mit++) {
0448       //get vec_digi range(pair)  corresponding to idx of map
0449       //loop over digi vector (ie between begin and end pointers in range)
0450       //CSCCorrelatedLCTDigiCollection::Range ctpRange = ctp_lct_data_->get((*mit).first)
0451       CSCDetId did = (*mit).first;
0452       //for (vecIt vit = ctpRange.first; vit != ctpRange.second; vit++) {
0453       for (vecIt vit = ldc.get((*mit).first).first; vit != ldc.get((*mit).first).second; vit++) {
0454         CSCCorrelatedLCTDigi dg = *vit;
0455         uint16_t tn = dg.getTrknmb();
0456         if (tn == 2)
0457           tn--;
0458         dg.setTrknmb(tn);
0459         dgcoll->insertDigi(did, dg);
0460       }
0461     }
0462     L1CSCTrack l1csctrk = std::make_pair(l1trk, *dgcoll);
0463     data->push_back(l1csctrk);
0464   }
0465 }
0466 
0467 #endif