Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-05-12 01:51:39

0001 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/OMTFinputMaker.h"
0002 #include "L1Trigger/L1TMuonOverlapPhase1/interface/MuonStub.h"
0003 #include "L1Trigger/L1TMuonOverlapPhase1/interface/ProcConfigurationBase.h"
0004 
0005 #include "DataFormats/CSCDigi/interface/CSCCorrelatedLCTDigi.h"
0006 #include "DataFormats/DetId/interface/DetId.h"
0007 #include "DataFormats/L1DTTrackFinder/interface/L1MuDTChambPhDigi.h"
0008 #include "DataFormats/L1TMuon/interface/RegionalMuonCandFwd.h"
0009 #include "DataFormats/MuonDetId/interface/CSCDetId.h"
0010 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
0011 #include "DataFormats/MuonDetId/interface/MuonSubdetId.h"
0012 #include "DataFormats/MuonDetId/interface/RPCDetId.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 
0016 #include <map>
0017 #include <sstream>
0018 #include <string>
0019 
0020 //dtThDigis is provided as argument, because in the OMTF implementation the phi and eta digis are merged (even thought it is artificial)
0021 void DtDigiToStubsConverterOmtf::addDTphiDigi(MuonStubPtrs2D& muonStubsInLayers,
0022                                               const L1MuDTChambPhDigi& digi,
0023                                               const L1MuDTChambThContainer* dtThDigis,
0024                                               unsigned int iProcessor,
0025                                               l1t::tftype procTyp) {
0026   DTChamberId detid(digi.whNum(), digi.stNum(), digi.scNum() + 1);
0027 
0028   //LogTrace("l1tOmtfEventPrint")<<__FUNCTION__<<":"<<__LINE__<<" OMTFinputMaker "<<" detid "<<detid<<endl;
0029   ///Check Trigger primitive quality
0030   ///Ts2Tag() == 0 - take only first track from DT Trigger Server
0031   ///BxCnt()  == 0 - ??
0032   ///code()>=3     - take only double layer hits, HH, HL and LL
0033   // FIXME (MK): at least Ts2Tag selection is not correct! Check it
0034   //    if (digiIt.bxNum()!= 0 || digiIt.BxCnt()!= 0 || digiIt.Ts2Tag()!= 0 || digiIt.code()<4) continue;
0035 
0036   //7 is empty digi, TODO update if the definition of the quality is changed
0037   if (digi.code() == 7 || digi.code() < config->getMinDtPhiQuality())
0038     return;
0039 
0040   unsigned int hwNumber = config->getLayerNumber(detid.rawId());
0041   if (config->getHwToLogicLayer().find(hwNumber) == config->getHwToLogicLayer().end())
0042     return;
0043 
0044   auto iter = config->getHwToLogicLayer().find(hwNumber);
0045   unsigned int iLayer = iter->second;
0046   unsigned int iInput = OMTFinputMaker::getInputNumber(config, detid.rawId(), iProcessor, procTyp);
0047   //MuonStub& stub = muonStubsInLayers[iLayer][iInput];
0048   MuonStub stub;
0049 
0050   stub.type = MuonStub::DT_PHI_ETA;
0051   stub.qualityHw = digi.code();
0052   stub.phiHw = angleConverter->getProcessorPhi(
0053       OMTFinputMaker::getProcessorPhiZero(config, iProcessor), procTyp, digi.scNum(), digi.phi());
0054   stub.etaHw = angleConverter->getGlobalEta(detid, dtThDigis, digi.bxNum());
0055 
0056   if (stub.qualityHw >= config->getMinDtPhiBQuality())
0057     stub.phiBHw = digi.phiB();
0058   else
0059     stub.phiBHw = config->nPhiBins();
0060 
0061   stub.bx = digi.bxNum();  //TODO sholdn't  it be BxCnt()?
0062   //stub.timing = digi.getTiming(); //TODO what about sub-bx timing, is is available?
0063 
0064   //stub.etaType = ?? TODO
0065   stub.logicLayer = iLayer;
0066   stub.detId = detid;
0067 
0068   OMTFinputMaker::addStub(config, muonStubsInLayers, iLayer, iInput, stub);
0069 }
0070 
0071 void DtDigiToStubsConverterOmtf::addDTetaStubs(MuonStubPtrs2D& muonStubsInLayers,
0072                                                const L1MuDTChambThDigi& thetaDigi,
0073                                                unsigned int iProcessor,
0074                                                l1t::tftype procTyp) {
0075   //in the Phase1 omtf the theta stubs are merged with the phi in the addDTphiDigi
0076 }
0077 
0078 bool DtDigiToStubsConverterOmtf::acceptDigi(const DTChamberId& dTChamberId,
0079                                             unsigned int iProcessor,
0080                                             l1t::tftype procType) {
0081   return OMTFinputMaker::acceptDtDigi(config, dTChamberId, iProcessor, procType);
0082 }
0083 
0084 void CscDigiToStubsConverterOmtf::addCSCstubs(MuonStubPtrs2D& muonStubsInLayers,
0085                                               unsigned int rawid,
0086                                               const CSCCorrelatedLCTDigi& digi,
0087                                               unsigned int iProcessor,
0088                                               l1t::tftype procTyp) {
0089   unsigned int hwNumber = config->getLayerNumber(rawid);
0090   if (config->getHwToLogicLayer().find(hwNumber) == config->getHwToLogicLayer().end())
0091     return;
0092 
0093   unsigned int iLayer = config->getHwToLogicLayer().at(hwNumber);
0094   unsigned int iInput = OMTFinputMaker::getInputNumber(config, rawid, iProcessor, procTyp);
0095 
0096   MuonStub stub;
0097   stub.type = MuonStub::CSC_PHI_ETA;
0098   stub.phiHw = angleConverter->getProcessorPhi(
0099       OMTFinputMaker::getProcessorPhiZero(config, iProcessor), procTyp, CSCDetId(rawid), digi);
0100   stub.etaHw = angleConverter->getGlobalEta(rawid, digi);
0101   stub.phiBHw = digi.getPattern();  //TODO change to phiB when implemented
0102   stub.qualityHw = digi.getQuality();
0103 
0104   stub.bx = digi.getBX() - config->cscLctCentralBx();  //TODO sholdn't  it be getBX0()?
0105   //stub.timing = digi.getTiming(); //TODO what about sub-bx timing, is is available?
0106 
0107   //stub.etaType = ?? TODO
0108   stub.logicLayer = iLayer;
0109   stub.detId = rawid;
0110 
0111   OMTFinputMaker::addStub(config, muonStubsInLayers, iLayer, iInput, stub);
0112   ///Accept CSC digis only up to eta=1.26.
0113   ///The nominal OMTF range is up to 1.24, but cutting at 1.24
0114   ///kill efficiency at the edge. 1.26 is one eta bin above nominal.
0115   //if(abs(iEta)>1.26/2.61*240) continue;
0116   //if (abs(iEta) > 115) continue;
0117 
0118   //LogTrace("l1tOmtfEventPrint")<<" ADDING CSC hit, proc: "<<iProcessor<<" iPhi : " << iPhi <<" iEta: "<< iEta << std::endl;
0119 }
0120 
0121 bool CscDigiToStubsConverterOmtf::acceptDigi(const CSCDetId& csc, unsigned int iProcessor, l1t::tftype procType) {
0122   if (procType == l1t::tftype::omtf_pos && (csc.endcap() == 2 || csc.ring() == 1 || csc.station() == 4))
0123     return false;
0124   if (procType == l1t::tftype::omtf_neg && (csc.endcap() == 1 || csc.ring() == 1 || csc.station() == 4))
0125     return false;
0126 
0127   if (procType == l1t::tftype::emtf_pos && (csc.endcap() == 2 || (csc.station() == 1 && csc.ring() == 3)))
0128     return false;
0129   if (procType == l1t::tftype::emtf_neg && (csc.endcap() == 1 || (csc.station() == 1 && csc.ring() == 3)))
0130     return false;
0131 
0132   unsigned int aSector = csc.chamber();
0133   unsigned int aMin = config->getEndcap10DegMin()[iProcessor];
0134   unsigned int aMax = config->getEndcap10DegMax()[iProcessor];
0135 
0136   if ((procType == l1t::tftype::emtf_pos || procType == l1t::tftype::emtf_neg) && csc.station() > 1 &&
0137       csc.ring() == 1) {
0138     aMin = config->getEndcap20DegMin()[iProcessor];
0139     aMax = config->getEndcap20DegMax()[iProcessor];
0140   }
0141 
0142   if (aMax > aMin && aSector >= aMin && aSector <= aMax)
0143     return true;
0144   if (aMax < aMin && (aSector >= aMin || aSector <= aMax))
0145     return true;
0146 
0147   return false;
0148 }
0149 
0150 void RpcDigiToStubsConverterOmtf::addRPCstub(MuonStubPtrs2D& muonStubsInLayers,
0151                                              const RPCDetId& roll,
0152                                              const RpcCluster& cluster,
0153                                              unsigned int iProcessor,
0154                                              l1t::tftype procTyp) {
0155   unsigned int rawid = roll.rawId();
0156 
0157   unsigned int hwNumber = config->getLayerNumber(rawid);
0158   unsigned int iLayer = config->getHwToLogicLayer().at(hwNumber);
0159   unsigned int iInput = OMTFinputMaker::getInputNumber(config, rawid, iProcessor, procTyp);
0160 
0161   //LogTrace("l1tOmtfEventPrint")<<"ADDING HIT: iLayer = " << iLayer << " iInput: " << iInput << " iPhi: " << iPhi << std::endl;
0162 
0163   //if (iLayer==17 && (iInput==0 || iInput==1)) continue;  // FIXME (MK) there is no RPC link for that input, because it is taken by DAQ link
0164 
0165   MuonStub stub;
0166   stub.type = MuonStub::RPC;
0167   stub.phiHw = angleConverter->getProcessorPhi(
0168       OMTFinputMaker::getProcessorPhiZero(config, iProcessor), procTyp, roll, cluster.firstStrip, cluster.lastStrip);
0169   stub.etaHw = angleConverter->getGlobalEtaRpc(rawid, cluster.firstStrip);
0170 
0171   stub.qualityHw = cluster.size();
0172 
0173   stub.bx = cluster.bx;
0174   stub.timing = cluster.timing;
0175 
0176   //stub.etaType = ?? TODO
0177   stub.logicLayer = iLayer;
0178   stub.detId = rawid;
0179 
0180   //This is very simple filtering of the clusters
0181   //Till Nov 2021: unfortunately performance of the firmware cannot be easily emulated from digi
0182   //(in principle would required raws, because in the firmware the clusterizaton is based on the 8-bit strip partitions)
0183   //The FW from from Nov 2021 solved this problem - option dropAllClustersIfMoreThanMax:
0184   //if any cluster is dropped in one barrel roll or endcap chamber - all are dropped for this roll/chamber.
0185   //Beside better data-to-emulator agreement it provides better eff for high pt muons
0186   if (config->getRpcDropAllClustersIfMoreThanMax()) {
0187     //two clusters were already added, so as we have the next one, we mark as dropped the one that was added before
0188     if (muonStubsInLayers[iLayer][iInput + 1]) {
0189       //if iInput+1 is not null, iInput is not null as well
0190       muonStubsInLayers[iLayer][iInput]->type = MuonStub::RPC_DROPPED;
0191       muonStubsInLayers[iLayer][iInput + 1]->type = MuonStub::RPC_DROPPED;
0192     } else if (cluster.size() > config->getRpcMaxClusterSize()) {
0193       //marking as dropped the one that was added before on the iInput
0194       if (muonStubsInLayers[iLayer][iInput]) {
0195         muonStubsInLayers[iLayer][iInput]->type = MuonStub::RPC_DROPPED;
0196 
0197         muonStubsInLayers[iLayer][iInput + 1] = std::make_shared<MuonStub>(stub);
0198         muonStubsInLayers[iLayer][iInput + 1]->type = MuonStub::RPC_DROPPED;
0199       } else {
0200         //no stub was added at this input already, so adding a stub and marking it as dropped
0201         muonStubsInLayers[iLayer].at(iInput) = std::make_shared<MuonStub>(stub);
0202         muonStubsInLayers[iLayer][iInput]->type = MuonStub::RPC_DROPPED;
0203 
0204         muonStubsInLayers[iLayer][iInput + 1] = std::make_shared<MuonStub>(stub);
0205         muonStubsInLayers[iLayer][iInput + 1]->type = MuonStub::RPC_DROPPED;
0206       }
0207     } else
0208       OMTFinputMaker::addStub(config, muonStubsInLayers, iLayer, iInput, stub);
0209   } else {
0210     if (cluster.size() <= config->getRpcMaxClusterSize())
0211       OMTFinputMaker::addStub(config, muonStubsInLayers, iLayer, iInput, stub);
0212   }
0213 
0214   std::ostringstream str;
0215   str << " RPC halfDigi "
0216       << " begin: " << cluster.firstStrip << " end: " << cluster.lastStrip << " iPhi: " << stub.phiHw
0217       << " iEta: " << stub.etaHw << " hwNumber: " << hwNumber << " iInput: " << iInput << " iLayer: " << iLayer
0218       << std::endl;
0219 
0220   LogTrace("l1tOmtfEventPrint") << str.str();
0221 }
0222 
0223 bool RpcDigiToStubsConverterOmtf::acceptDigi(const RPCDetId& rpcDetId, unsigned int iProcessor, l1t::tftype procType) {
0224   unsigned int aMin = config->getBarrelMin()[iProcessor];
0225   unsigned int aMax = config->getBarrelMax()[iProcessor];
0226 
0227   unsigned int aSector = rpcDetId.sector();
0228 
0229   ///Select RPC chambers connected to OMTF
0230   if (procType == l1t::tftype::omtf_pos &&
0231       (rpcDetId.region() < 0 || (rpcDetId.region() == 0 && rpcDetId.ring() != 2) ||
0232        (rpcDetId.region() == 0 && rpcDetId.station() == 4) ||
0233        (rpcDetId.region() == 0 && rpcDetId.station() == 2 && rpcDetId.layer() == 2 && rpcDetId.roll() == 1) ||
0234        (rpcDetId.region() == 0 && rpcDetId.station() == 3 && rpcDetId.roll() == 1) ||
0235        (rpcDetId.region() == 1 && rpcDetId.station() == 4) ||
0236        ///RPC RE1/2 temporarily not used (rpcDetId.region()==1 && rpcDetId.station()==1 && rpcDetId.ring()<2) ||
0237        (rpcDetId.region() == 1 && rpcDetId.station() > 0 && rpcDetId.ring() < 3)))
0238     return false;
0239 
0240   if (procType == l1t::tftype::omtf_neg &&
0241       (rpcDetId.region() > 0 || (rpcDetId.region() == 0 && rpcDetId.ring() != -2) ||
0242        (rpcDetId.region() == 0 && rpcDetId.station() == 4) ||
0243        (rpcDetId.region() == 0 && rpcDetId.station() == 2 && rpcDetId.layer() == 2 && rpcDetId.roll() == 1) ||
0244        (rpcDetId.region() == 0 && rpcDetId.station() == 3 && rpcDetId.roll() == 1) ||
0245        (rpcDetId.region() == -1 && rpcDetId.station() == 4) ||
0246        //RPC RE1/2 temporarily not used (rpcDetId.region()==1 && rpcDetId.station()==1 && rpcDetId.ring()<2) ||
0247        (rpcDetId.region() == -1 && rpcDetId.station() > 0 && rpcDetId.ring() < 3)))
0248     return false;
0249 
0250   if (procType == l1t::tftype::omtf_pos || procType == l1t::tftype::omtf_neg) {
0251     if (rpcDetId.region() != 0 && rpcDetId.station() == 3) {  //endcaps, layer 17
0252       unsigned int iInput = OMTFinputMaker::getInputNumber(config, rpcDetId.rawId(), iProcessor, procType);
0253       if (iInput == 0 || iInput == 1)
0254         return false;  // FIXME (MK) there is no RPC link for that input, because it is taken by DAQ link
0255     }
0256   }
0257 
0258   if (procType == l1t::tftype::bmtf && rpcDetId.region() != 0)
0259     return false;
0260 
0261   if (procType == l1t::tftype::emtf_pos &&
0262       (rpcDetId.region() <= 0 || (rpcDetId.station() == 1 && rpcDetId.ring() == 3)))
0263     return false;
0264   if (procType == l1t::tftype::emtf_neg &&
0265       (rpcDetId.region() >= 0 || (rpcDetId.station() == 1 && rpcDetId.ring() == 3)))
0266     return false;
0267   ////////////////
0268   if (rpcDetId.region() == 0)
0269     aSector = rpcDetId.sector();
0270   if (rpcDetId.region() != 0) {
0271     aSector = (rpcDetId.sector() - 1) * 6 + rpcDetId.subsector();
0272     aMin = config->getEndcap10DegMin()[iProcessor];
0273     aMax = config->getEndcap10DegMax()[iProcessor];
0274   }
0275 
0276   if (aMax > aMin && aSector >= aMin && aSector <= aMax)
0277     return true;
0278   if (aMax < aMin && (aSector >= aMin || aSector <= aMax))
0279     return true;
0280 
0281   return false;
0282 }
0283 ///////////////////////////////////////
0284 ///////////////////////////////////////
0285 bool OMTFinputMaker::acceptDtDigi(const OMTFConfiguration* config,
0286                                   const DTChamberId& dTChamberId,
0287                                   unsigned int iProcessor,
0288                                   l1t::tftype procType) {
0289   unsigned int aMin = config->getBarrelMin()[iProcessor];
0290   unsigned int aMax = config->getBarrelMax()[iProcessor];
0291 
0292   if (procType == l1t::tftype::omtf_pos && dTChamberId.wheel() != 2)
0293     return false;
0294   if (procType == l1t::tftype::omtf_neg && dTChamberId.wheel() != -2)
0295     return false;
0296   if (procType == l1t::tftype::emtf_pos || procType == l1t::tftype::emtf_neg)
0297     return false;
0298 
0299   unsigned int aSector = dTChamberId.sector();
0300 
0301   if (aMax > aMin && aSector >= aMin && aSector <= aMax)
0302     return true;
0303   if (aMax < aMin && (aSector >= aMin || aSector <= aMax))
0304     return true;
0305 
0306   return false;
0307 }
0308 
0309 ///////////////////////////////////////
0310 ///////////////////////////////////////
0311 OMTFinputMaker::OMTFinputMaker(const edm::ParameterSet& edmParameterSet,
0312                                MuStubsInputTokens& muStubsInputTokens,
0313                                const OMTFConfiguration* config,
0314                                std::unique_ptr<OmtfAngleConverter> angleConverter)
0315     : MuonStubMakerBase(config), config(config), angleConverter(std::move(angleConverter)) {
0316   if (!edmParameterSet.getParameter<bool>("dropDTPrimitives"))
0317     digiToStubsConverters.emplace_back(std::make_unique<DtDigiToStubsConverterOmtf>(
0318         config, this->angleConverter.get(), muStubsInputTokens.inputTokenDtPh, muStubsInputTokens.inputTokenDtTh));
0319 
0320   if (!edmParameterSet.getParameter<bool>("dropCSCPrimitives"))
0321     digiToStubsConverters.emplace_back(std::make_unique<CscDigiToStubsConverterOmtf>(
0322         config, this->angleConverter.get(), muStubsInputTokens.inputTokenCSC));
0323 
0324   if (!edmParameterSet.getParameter<bool>("dropRPCPrimitives"))
0325     digiToStubsConverters.emplace_back(std::make_unique<RpcDigiToStubsConverterOmtf>(
0326         config, this->angleConverter.get(), &rpcClusterization, muStubsInputTokens.inputTokenRPC));
0327 }
0328 
0329 void OMTFinputMaker::initialize(const edm::ParameterSet& edmCfg,
0330                                 const edm::EventSetup& es,
0331                                 const MuonGeometryTokens& muonGeometryTokens) {
0332   MuonStubMakerBase::initialize(edmCfg, es, muonGeometryTokens);
0333   angleConverter->checkAndUpdateGeometry(es, config, muonGeometryTokens);
0334 }
0335 
0336 ///////////////////////////////////////
0337 ///////////////////////////////////////
0338 OMTFinputMaker::~OMTFinputMaker() {}
0339 ///////////////////////////////////////
0340 ///////////////////////////////////////
0341 unsigned int OMTFinputMaker::getInputNumber(const OMTFConfiguration* config,
0342                                             unsigned int rawId,
0343                                             unsigned int iProcessor,
0344                                             l1t::tftype type) {
0345   unsigned int iInput = 99;
0346   unsigned int aSector = 99;
0347   int aMin = config->getBarrelMin()[iProcessor];
0348   int iRoll = 1;
0349   int nInputsPerSector = 2;
0350 
0351   DetId detId(rawId);
0352   if (detId.det() != DetId::Muon)
0353     edm::LogError("Critical OMTFinputMaker") << "PROBLEM: hit in unknown Det, detID: " << detId.det() << std::endl;
0354   switch (detId.subdetId()) {
0355     case MuonSubdetId::RPC: {
0356       RPCDetId rpc(rawId);
0357       if (rpc.region() == 0) {
0358         nInputsPerSector = 4;
0359         aSector = rpc.sector();
0360         ///on the 0-2pi border we need to add 1 30 deg sector
0361         ///to get the correct index
0362         if (iProcessor == 5 && aSector < 3)
0363           aMin = -1;
0364         //Use division into rolls
0365         iRoll = rpc.roll();
0366         ///Set roll number by hand to keep common input
0367         ///number shift formula for all stations
0368         if (rpc.station() == 2 && rpc.layer() == 2 && rpc.roll() == 2)
0369           iRoll = 1;
0370         ///Only one roll from station 3 is connected.
0371         if (rpc.station() == 3) {
0372           iRoll = 1;
0373           nInputsPerSector = 2;
0374         }
0375         ///At the moment do not use RPC chambers splitting into rolls for bmtf part
0376         if (type == l1t::tftype::bmtf)
0377           iRoll = 1;
0378       }
0379       if (rpc.region() != 0) {
0380         aSector = (rpc.sector() - 1) * 6 + rpc.subsector();
0381         aMin = config->getEndcap10DegMin()[iProcessor];
0382         ///on the 0-2pi border we need to add 4 10 deg sectors
0383         ///to get the correct index
0384         if (iProcessor == 5 && aSector < 5)
0385           aMin = -4;
0386       }
0387       break;
0388     }
0389     case MuonSubdetId::DT: {
0390       DTChamberId dt(rawId);
0391       aSector = dt.sector();
0392       ///on the 0-2pi border we need to add 1 30 deg sector
0393       ///to get the correct index
0394       if (iProcessor == 5 && aSector < 3)
0395         aMin = -1;
0396       break;
0397     }
0398     case MuonSubdetId::CSC: {
0399       CSCDetId csc(rawId);
0400       aSector = csc.chamber();
0401       aMin = config->getEndcap10DegMin()[iProcessor];
0402       ///on the 0-2pi border we need to add 4 10deg sectors
0403       ///to get the correct index
0404       if (iProcessor == 5 && aSector < 5)
0405         aMin = -4;
0406       ///Endcap region covers algo 10 deg sectors
0407       ///on the 0-2pi border we need to add 2 20deg sectors
0408       ///to get the correct index
0409       if ((type == l1t::tftype::emtf_pos || type == l1t::tftype::emtf_neg) && csc.station() > 1 && csc.ring() == 1) {
0410         aMin = config->getEndcap20DegMin()[iProcessor];
0411         if (iProcessor == 5 && aSector < 3)
0412           aMin = -2;
0413       }
0414       break;
0415     }
0416   }
0417 
0418   ///Assume 2 hits per chamber
0419   iInput = (aSector - aMin) * nInputsPerSector;
0420   ///Chambers divided into two rolls have rolls number 1 and 3
0421   iInput += iRoll - 1;
0422 
0423   return iInput;
0424 }
0425 ////////////////////////////////////////////
0426 ////////////////////////////////////////////
0427 
0428 //iProcessor counted from 0
0429 int OMTFinputMaker::getProcessorPhiZero(const OMTFConfiguration* config, unsigned int iProcessor) {
0430   unsigned int nPhiBins = config->nPhiBins();
0431 
0432   int phiZero = nPhiBins / 6. * (iProcessor) + nPhiBins / 24;
0433   // "0" is 15degree moved cyclically to each processor, note [0,2pi]
0434 
0435   return config->foldPhi(phiZero);
0436 }
0437 
0438 ////////////////////////////////////////////
0439 ////////////////////////////////////////////
0440 void OMTFinputMaker::addStub(const OMTFConfiguration* config,
0441                              MuonStubPtrs2D& muonStubsInLayers,
0442                              unsigned int iLayer,
0443                              unsigned int iInput,
0444                              MuonStub& stub) {
0445   //LogTrace("l1tOmtfEventPrint") << __FUNCTION__ << ":" << __LINE__ << " iInput " << iInput << " " << stub << endl;
0446   //there is a small rate of duplicated digis in the real data in the DT and CSC, the reason for this duplicates is not understood
0447   //in case of RPC the duplicated digis appear in data only for the layer 17 (RE3), where the rolls 2 and 3 has the same eta = 115 assigned, and the muon can hit both rolls
0448   //the reason of the  duplicates cannot be the fact that the same data (links) goes to neighboring boards, because this is  removed in the OMTF unpacker
0449   //the duplicates cannot be dropped, because if they appear in data, are also duplicated on the input of the algorithm, and both can be used as the reference hits
0450   if (muonStubsInLayers[iLayer][iInput] && muonStubsInLayers[iLayer][iInput]->phiHw == stub.phiHw &&
0451       muonStubsInLayers[iLayer][iInput]->phiBHw == stub.phiBHw &&
0452       muonStubsInLayers[iLayer][iInput]->etaHw == stub.etaHw) {
0453     LogTrace("OMTFReconstruction") << "addStub: the stub with exactly the same phi, phiB and eta was already added:\n"
0454                                    << "incomnig stub " << stub << "\n"
0455                                    << "existing stub " << *(muonStubsInLayers[iLayer][iInput]) << std::endl;
0456     //return;
0457   }
0458 
0459   if (muonStubsInLayers[iLayer][iInput] && muonStubsInLayers[iLayer][iInput]->phiHw != (int)config->nPhiBins())
0460     ++iInput;
0461 
0462   if (muonStubsInLayers[iLayer][iInput] && muonStubsInLayers[iLayer][iInput]->phiHw != (int)config->nPhiBins())
0463     return;
0464   //in this implementation only two first stubs are added for a given iInput
0465 
0466   muonStubsInLayers.at(iLayer).at(iInput) = std::make_shared<MuonStub>(stub);
0467 }