Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 #include "FWCore/Framework/interface/MakerMacros.h"
0002 #include "FWCore/Framework/interface/stream/EDProducer.h"
0003 #include "FWCore/Framework/interface/ESHandle.h"
0004 #include "FWCore/Utilities/interface/ESGetToken.h"
0005 #include "FWCore/Framework/interface/ModuleFactory.h"
0006 #include "FWCore/Framework/interface/ESProducer.h"
0007 #include "FWCore/Framework/interface/ESHandle.h"
0008 #include "FWCore/Framework/interface/ESProducts.h"
0009 #include "FWCore/Framework/interface/ConsumesCollector.h"
0010 
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/Frameworkfwd.h"
0013 #include "FWCore/Framework/interface/EventSetup.h"
0014 #include "FWCore/Framework/interface/Run.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 
0018 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0019 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0020 #include "Geometry/DTGeometry/interface/DTLayer.h"
0021 
0022 #include "L1Trigger/DTTriggerPhase2/interface/MuonPath.h"
0023 #include "L1Trigger/DTTriggerPhase2/interface/constants.h"
0024 
0025 #include "L1Trigger/DTTriggerPhase2/interface/MotherGrouping.h"
0026 #include "L1Trigger/DTTriggerPhase2/interface/InitialGrouping.h"
0027 #include "L1Trigger/DTTriggerPhase2/interface/HoughGrouping.h"
0028 #include "L1Trigger/DTTriggerPhase2/interface/PseudoBayesGrouping.h"
0029 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAnalyzer.h"
0030 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAnalyticAnalyzer.h"
0031 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAnalyzerInChamber.h"
0032 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAssociator.h"
0033 #include "L1Trigger/DTTriggerPhase2/interface/MPFilter.h"
0034 #include "L1Trigger/DTTriggerPhase2/interface/MPQualityEnhancerFilter.h"
0035 #include "L1Trigger/DTTriggerPhase2/interface/MPRedundantFilter.h"
0036 #include "L1Trigger/DTTriggerPhase2/interface/MPCleanHitsFilter.h"
0037 #include "L1Trigger/DTTriggerPhase2/interface/MPQualityEnhancerFilterBayes.h"
0038 #include "L1Trigger/DTTriggerPhase2/interface/GlobalCoordsObtainer.h"
0039 
0040 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
0041 #include "DataFormats/MuonDetId/interface/DTSuperLayerId.h"
0042 #include "DataFormats/MuonDetId/interface/DTLayerId.h"
0043 #include "DataFormats/MuonDetId/interface/DTWireId.h"
0044 #include "DataFormats/DTDigi/interface/DTDigiCollection.h"
0045 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTPhContainer.h"
0046 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTPhDigi.h"
0047 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTExtPhContainer.h"
0048 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTExtPhDigi.h"
0049 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTThContainer.h"
0050 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTThDigi.h"
0051 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTExtThContainer.h"
0052 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTExtThDigi.h"
0053 
0054 // DT trigger GeomUtils
0055 #include "DQM/DTMonitorModule/interface/DTTrigGeomUtils.h"
0056 
0057 //RPC TP
0058 #include "DataFormats/RPCRecHit/interface/RPCRecHitCollection.h"
0059 #include <DataFormats/MuonDetId/interface/RPCDetId.h>
0060 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"
0061 #include "L1Trigger/DTTriggerPhase2/interface/RPCIntegrator.h"
0062 
0063 #include <fstream>
0064 #include <iostream>
0065 #include <queue>
0066 #include <cmath>
0067 
0068 using namespace edm;
0069 using namespace std;
0070 using namespace cmsdt;
0071 
0072 class DTTrigPhase2Prod : public edm::stream::EDProducer<> {
0073   typedef std::map<DTChamberId, DTDigiCollection, std::less<DTChamberId>> DTDigiMap;
0074   typedef DTDigiMap::iterator DTDigiMap_iterator;
0075   typedef DTDigiMap::const_iterator DTDigiMap_const_iterator;
0076 
0077 public:
0078   //! Constructor
0079   DTTrigPhase2Prod(const edm::ParameterSet& pset);
0080 
0081   //! Destructor
0082   ~DTTrigPhase2Prod() override;
0083 
0084   //! Create Trigger Units before starting event processing
0085   void beginRun(edm::Run const& iRun, const edm::EventSetup& iEventSetup) override;
0086 
0087   //! Producer: process every event and generates trigger data
0088   void produce(edm::Event& iEvent, const edm::EventSetup& iEventSetup) override;
0089 
0090   //! endRun: finish things
0091   void endRun(edm::Run const& iRun, const edm::EventSetup& iEventSetup) override;
0092 
0093   // Methods
0094   int rango(const metaPrimitive& mp) const;
0095   bool outer(const metaPrimitive& mp) const;
0096   bool inner(const metaPrimitive& mp) const;
0097   void printmP(const std::string& ss, const metaPrimitive& mP) const;
0098   void printmPC(const std::string& ss, const metaPrimitive& mP) const;
0099   bool hasPosRF(int wh, int sec) const;
0100 
0101   // Getter-methods
0102   MP_QUALITY getMinimumQuality(void);
0103 
0104   // Setter-methods
0105   void setChiSquareThreshold(float ch2Thr);
0106   void setMinimumQuality(MP_QUALITY q);
0107 
0108   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0109 
0110   // data-members
0111   const DTGeometry* dtGeo_;
0112   edm::ESGetToken<DTGeometry, MuonGeometryRecord> dtGeomH;
0113   std::vector<std::pair<int, MuonPath>> primitives_;
0114 
0115 private:
0116   // Trigger Configuration Manager CCB validity flag
0117   bool my_CCBValid_;
0118 
0119   // BX offset used to correct DTTPG output
0120   int my_BXoffset_;
0121 
0122   // Debug Flag
0123   bool debug_;
0124   bool dump_;
0125   double dT0_correlate_TP_;
0126   int scenario_;
0127   int df_extended_;
0128   int max_index_;
0129 
0130   // ParameterSet
0131   edm::EDGetTokenT<DTDigiCollection> dtDigisToken_;
0132   edm::EDGetTokenT<RPCRecHitCollection> rpcRecHitsLabel_;
0133 
0134   // Grouping attributes and methods
0135   int algo_;  // Grouping code
0136   std::unique_ptr<MotherGrouping> grouping_obj_;
0137   std::unique_ptr<MuonPathAnalyzer> mpathanalyzer_;
0138   std::unique_ptr<MPFilter> mpathqualityenhancer_;
0139   std::unique_ptr<MPFilter> mpathqualityenhancerbayes_;
0140   std::unique_ptr<MPFilter> mpathredundantfilter_;
0141   std::unique_ptr<MPFilter> mpathhitsfilter_;
0142   std::unique_ptr<MuonPathAssociator> mpathassociator_;
0143   std::shared_ptr<GlobalCoordsObtainer> globalcoordsobtainer_;
0144 
0145   // Buffering
0146   bool activateBuffer_;
0147   int superCellhalfspacewidth_;
0148   float superCelltimewidth_;
0149   std::vector<DTDigiCollection*> distribDigis(std::queue<std::pair<DTLayerId, DTDigi>>& inQ);
0150   void processDigi(std::queue<std::pair<DTLayerId, DTDigi>>& inQ,
0151                    std::vector<std::queue<std::pair<DTLayerId, DTDigi>>*>& vec);
0152 
0153   // RPC
0154   std::unique_ptr<RPCIntegrator> rpc_integrator_;
0155   bool useRPC_;
0156 
0157   void assignIndex(std::vector<metaPrimitive>& inMPaths);
0158   void assignIndexPerBX(std::vector<metaPrimitive>& inMPaths);
0159   int assignQualityOrder(const metaPrimitive& mP) const;
0160 
0161   const std::unordered_map<int, int> qmap_;
0162 };
0163 
0164 namespace {
0165   struct {
0166     bool operator()(std::pair<DTLayerId, DTDigi> a, std::pair<DTLayerId, DTDigi> b) const {
0167       return (a.second.time() < b.second.time());
0168     }
0169   } const DigiTimeOrdering;
0170 }  // namespace
0171 
0172 DTTrigPhase2Prod::DTTrigPhase2Prod(const ParameterSet& pset)
0173     : qmap_({{8, 8}, {7, 7}, {6, 6}, {4, 4}, {3, 3}, {2, 2}, {1, 1}}) {
0174   produces<L1Phase2MuDTPhContainer>();
0175   produces<L1Phase2MuDTThContainer>();
0176   produces<L1Phase2MuDTExtPhContainer>();
0177   produces<L1Phase2MuDTExtThContainer>();
0178 
0179   debug_ = pset.getUntrackedParameter<bool>("debug");
0180   dump_ = pset.getUntrackedParameter<bool>("dump");
0181 
0182   scenario_ = pset.getParameter<int>("scenario");
0183 
0184   df_extended_ = pset.getParameter<int>("df_extended");
0185   max_index_ = pset.getParameter<int>("max_primitives") - 1;
0186 
0187   dtDigisToken_ = consumes<DTDigiCollection>(pset.getParameter<edm::InputTag>("digiTag"));
0188 
0189   rpcRecHitsLabel_ = consumes<RPCRecHitCollection>(pset.getParameter<edm::InputTag>("rpcRecHits"));
0190   useRPC_ = pset.getParameter<bool>("useRPC");
0191 
0192   // Choosing grouping scheme:
0193   algo_ = pset.getParameter<int>("algo");
0194 
0195   edm::ConsumesCollector consumesColl(consumesCollector());
0196   globalcoordsobtainer_ = std::make_shared<GlobalCoordsObtainer>(pset);
0197   globalcoordsobtainer_->generate_luts();
0198 
0199   if (algo_ == PseudoBayes) {
0200     grouping_obj_ =
0201         std::make_unique<PseudoBayesGrouping>(pset.getParameter<edm::ParameterSet>("PseudoBayesPattern"), consumesColl);
0202   } else if (algo_ == HoughTrans) {
0203     grouping_obj_ =
0204         std::make_unique<HoughGrouping>(pset.getParameter<edm::ParameterSet>("HoughGrouping"), consumesColl);
0205   } else {
0206     grouping_obj_ = std::make_unique<InitialGrouping>(pset, consumesColl);
0207   }
0208 
0209   if (algo_ == Standard) {
0210     if (debug_)
0211       LogDebug("DTTrigPhase2Prod") << "DTp2:constructor: JM analyzer";
0212     mpathanalyzer_ = std::make_unique<MuonPathAnalyticAnalyzer>(pset, consumesColl, globalcoordsobtainer_);
0213   } else {
0214     if (debug_)
0215       LogDebug("DTTrigPhase2Prod") << "DTp2:constructor: Full chamber analyzer";
0216     mpathanalyzer_ = std::make_unique<MuonPathAnalyzerInChamber>(pset, consumesColl, globalcoordsobtainer_);
0217   }
0218 
0219   // Getting buffer option
0220   activateBuffer_ = pset.getParameter<bool>("activateBuffer");
0221   superCellhalfspacewidth_ = pset.getParameter<int>("superCellspacewidth") / 2;
0222   superCelltimewidth_ = pset.getParameter<double>("superCelltimewidth");
0223 
0224   mpathqualityenhancer_ = std::make_unique<MPQualityEnhancerFilter>(pset);
0225   mpathqualityenhancerbayes_ = std::make_unique<MPQualityEnhancerFilterBayes>(pset);
0226   mpathredundantfilter_ = std::make_unique<MPRedundantFilter>(pset);
0227   mpathhitsfilter_ = std::make_unique<MPCleanHitsFilter>(pset);
0228   mpathassociator_ = std::make_unique<MuonPathAssociator>(pset, consumesColl, globalcoordsobtainer_);
0229   rpc_integrator_ = std::make_unique<RPCIntegrator>(pset, consumesColl);
0230 
0231   dtGeomH = esConsumes<DTGeometry, MuonGeometryRecord, edm::Transition::BeginRun>();
0232 }
0233 
0234 DTTrigPhase2Prod::~DTTrigPhase2Prod() {
0235   if (debug_)
0236     LogDebug("DTTrigPhase2Prod") << "DTp2: calling destructor" << std::endl;
0237 }
0238 
0239 void DTTrigPhase2Prod::beginRun(edm::Run const& iRun, const edm::EventSetup& iEventSetup) {
0240   if (debug_)
0241     LogDebug("DTTrigPhase2Prod") << "beginRun " << iRun.id().run();
0242   if (debug_)
0243     LogDebug("DTTrigPhase2Prod") << "beginRun: getting DT geometry";
0244 
0245   grouping_obj_->initialise(iEventSetup);               // Grouping object initialisation
0246   mpathanalyzer_->initialise(iEventSetup);              // Analyzer object initialisation
0247   mpathqualityenhancer_->initialise(iEventSetup);       // Filter object initialisation
0248   mpathredundantfilter_->initialise(iEventSetup);       // Filter object initialisation
0249   mpathqualityenhancerbayes_->initialise(iEventSetup);  // Filter object initialisation
0250   mpathhitsfilter_->initialise(iEventSetup);
0251   mpathassociator_->initialise(iEventSetup);  // Associator object initialisation
0252 
0253   if (auto geom = iEventSetup.getHandle(dtGeomH)) {
0254     dtGeo_ = &(*geom);
0255   }
0256 }
0257 
0258 void DTTrigPhase2Prod::produce(Event& iEvent, const EventSetup& iEventSetup) {
0259   if (debug_)
0260     LogDebug("DTTrigPhase2Prod") << "produce";
0261   edm::Handle<DTDigiCollection> dtdigis;
0262   iEvent.getByToken(dtDigisToken_, dtdigis);
0263 
0264   if (debug_)
0265     LogDebug("DTTrigPhase2Prod") << "\t Getting the RPC RecHits" << std::endl;
0266   edm::Handle<RPCRecHitCollection> rpcRecHits;
0267   iEvent.getByToken(rpcRecHitsLabel_, rpcRecHits);
0268 
0269   ////////////////////////////////
0270   // GROUPING CODE:
0271   ////////////////////////////////
0272 
0273   DTDigiMap digiMap;
0274   DTDigiCollection::DigiRangeIterator detUnitIt;
0275   for (const auto& detUnitIt : *dtdigis) {
0276     const DTLayerId& layId = detUnitIt.first;
0277     const DTChamberId chambId = layId.superlayerId().chamberId();
0278     const DTDigiCollection::Range& range = detUnitIt.second;
0279     digiMap[chambId].put(range, layId);
0280   }
0281 
0282   // generate a list muon paths for each event!!!
0283   if (debug_ && activateBuffer_)
0284     LogDebug("DTTrigPhase2Prod") << "produce - Getting and grouping digis per chamber using a buffer and super cells.";
0285   else if (debug_)
0286     LogDebug("DTTrigPhase2Prod") << "produce - Getting and grouping digis per chamber.";
0287 
0288   MuonPathPtrs muonpaths;
0289   for (const auto& ich : dtGeo_->chambers()) {
0290     // The code inside this for loop would ideally later fit inside a trigger unit (in principle, a DT station) of the future Phase 2 DT Trigger.
0291     const DTChamber* chamb = ich;
0292     DTChamberId chid = chamb->id();
0293     DTDigiMap_iterator dmit = digiMap.find(chid);
0294 
0295     if (dmit == digiMap.end())
0296       continue;
0297 
0298     if (activateBuffer_) {  // Use buffering (per chamber) or not
0299       // Import digis from the station
0300       std::vector<std::pair<DTLayerId, DTDigi>> tmpvec;
0301       tmpvec.clear();
0302 
0303       for (const auto& dtLayerIdIt : (*dmit).second) {
0304         for (DTDigiCollection::const_iterator digiIt = (dtLayerIdIt.second).first;
0305              digiIt != (dtLayerIdIt.second).second;
0306              digiIt++) {
0307           tmpvec.emplace_back(dtLayerIdIt.first, *digiIt);
0308         }
0309       }
0310 
0311       // Check to enhance CPU time usage
0312       if (tmpvec.empty())
0313         continue;
0314 
0315       // Order digis depending on TDC time and insert them into a queue (FIFO buffer). TODO: adapt for MC simulations.
0316       std::sort(tmpvec.begin(), tmpvec.end(), DigiTimeOrdering);
0317       std::queue<std::pair<DTLayerId, DTDigi>> timequeue;
0318 
0319       for (const auto& elem : tmpvec)
0320         timequeue.emplace(elem);
0321       tmpvec.clear();
0322 
0323       // Distribute the digis from the queue into supercells
0324       std::vector<DTDigiCollection*> superCells;
0325       superCells = distribDigis(timequeue);
0326 
0327       // Process each supercell & collect the resulting muonpaths (as the muonpaths std::vector is only enlarged each time
0328       // the groupings access it, it's not needed to "collect" the final products).
0329 
0330       while (!superCells.empty()) {
0331         grouping_obj_->run(iEvent, iEventSetup, *(superCells.back()), muonpaths);
0332         superCells.pop_back();
0333       }
0334     } else {
0335       grouping_obj_->run(iEvent, iEventSetup, (*dmit).second, muonpaths);
0336     }
0337   }
0338   digiMap.clear();
0339 
0340   if (dump_) {
0341     for (unsigned int i = 0; i < muonpaths.size(); i++) {
0342       stringstream ss;
0343       ss << iEvent.id().event() << "      mpath " << i << ": ";
0344       for (int lay = 0; lay < muonpaths.at(i)->nprimitives(); lay++)
0345         ss << muonpaths.at(i)->primitive(lay)->channelId() << " ";
0346       for (int lay = 0; lay < muonpaths.at(i)->nprimitives(); lay++)
0347         ss << muonpaths.at(i)->primitive(lay)->tdcTimeStamp() << " ";
0348       for (int lay = 0; lay < muonpaths.at(i)->nprimitives(); lay++)
0349         ss << muonpaths.at(i)->primitive(lay)->laterality() << " ";
0350       LogInfo("DTTrigPhase2Prod") << ss.str();
0351     }
0352   }
0353 
0354   // FILTER GROUPING
0355   MuonPathPtrs filteredmuonpaths;
0356   if (algo_ == Standard) {
0357     mpathredundantfilter_->run(iEvent, iEventSetup, muonpaths, filteredmuonpaths);
0358   } else {
0359     mpathhitsfilter_->run(iEvent, iEventSetup, muonpaths, filteredmuonpaths);
0360   }
0361 
0362   if (dump_) {
0363     for (unsigned int i = 0; i < filteredmuonpaths.size(); i++) {
0364       stringstream ss;
0365       ss << iEvent.id().event() << " filt. mpath " << i << ": ";
0366       for (int lay = 0; lay < filteredmuonpaths.at(i)->nprimitives(); lay++)
0367         ss << filteredmuonpaths.at(i)->primitive(lay)->channelId() << " ";
0368       for (int lay = 0; lay < filteredmuonpaths.at(i)->nprimitives(); lay++)
0369         ss << filteredmuonpaths.at(i)->primitive(lay)->tdcTimeStamp() << " ";
0370       LogInfo("DTTrigPhase2Prod") << ss.str();
0371     }
0372   }
0373 
0374   ///////////////////////////////////////////
0375   /// Fitting SECTION;
0376   ///////////////////////////////////////////
0377 
0378   if (debug_)
0379     LogDebug("DTTrigPhase2Prod") << "MUON PATHS found: " << muonpaths.size() << " (" << filteredmuonpaths.size()
0380                                  << ") in event " << iEvent.id().event();
0381   if (debug_)
0382     LogDebug("DTTrigPhase2Prod") << "filling NmetaPrimtives" << std::endl;
0383   std::vector<metaPrimitive> metaPrimitives;
0384   MuonPathPtrs outmpaths;
0385   if (algo_ == Standard) {
0386     if (debug_)
0387       LogDebug("DTTrigPhase2Prod") << "Fitting 1SL ";
0388     mpathanalyzer_->run(iEvent, iEventSetup, filteredmuonpaths, metaPrimitives);
0389   } else {
0390     // implementation for advanced (2SL) grouping, no filter required..
0391     if (debug_)
0392       LogDebug("DTTrigPhase2Prod") << "Fitting 2SL at once ";
0393     mpathanalyzer_->run(iEvent, iEventSetup, muonpaths, outmpaths);
0394   }
0395 
0396   if (dump_) {
0397     for (unsigned int i = 0; i < outmpaths.size(); i++) {
0398       LogInfo("DTTrigPhase2Prod") << iEvent.id().event() << " mp " << i << ": " << outmpaths.at(i)->bxTimeValue() << " "
0399                                   << outmpaths.at(i)->horizPos() << " " << outmpaths.at(i)->tanPhi() << " "
0400                                   << outmpaths.at(i)->phi() << " " << outmpaths.at(i)->phiB() << " "
0401                                   << outmpaths.at(i)->quality() << " " << outmpaths.at(i)->chiSquare();
0402     }
0403     for (unsigned int i = 0; i < metaPrimitives.size(); i++) {
0404       stringstream ss;
0405       ss << iEvent.id().event() << " mp " << i << ": ";
0406       printmP(ss.str(), metaPrimitives.at(i));
0407     }
0408   }
0409 
0410   muonpaths.clear();
0411   filteredmuonpaths.clear();
0412 
0413   /////////////////////////////////////
0414   //  FILTER SECTIONS:
0415   ////////////////////////////////////
0416 
0417   if (debug_)
0418     LogDebug("DTTrigPhase2Prod") << "declaring new vector for filtered" << std::endl;
0419 
0420   std::vector<metaPrimitive> filteredMetaPrimitives;
0421   if (algo_ == Standard)
0422     mpathqualityenhancer_->run(iEvent, iEventSetup, metaPrimitives, filteredMetaPrimitives);
0423 
0424   if (dump_) {
0425     for (unsigned int i = 0; i < filteredMetaPrimitives.size(); i++) {
0426       stringstream ss;
0427       ss << iEvent.id().event() << " filtered mp " << i << ": ";
0428       printmP(ss.str(), filteredMetaPrimitives.at(i));
0429     }
0430   }
0431 
0432   metaPrimitives.clear();
0433   metaPrimitives.erase(metaPrimitives.begin(), metaPrimitives.end());
0434 
0435   if (debug_)
0436     LogDebug("DTTrigPhase2Prod") << "DTp2 in event:" << iEvent.id().event() << " we found "
0437                                  << filteredMetaPrimitives.size() << " filteredMetaPrimitives (superlayer)"
0438                                  << std::endl;
0439   if (debug_)
0440     LogDebug("DTTrigPhase2Prod") << "filteredMetaPrimitives: starting correlations" << std::endl;
0441 
0442   /////////////////////////////////////
0443   //// CORRELATION:
0444   /////////////////////////////////////
0445 
0446   std::vector<metaPrimitive> correlatedMetaPrimitives;
0447   if (algo_ == Standard)
0448     mpathassociator_->run(iEvent, iEventSetup, dtdigis, filteredMetaPrimitives, correlatedMetaPrimitives);
0449   else {
0450     for (const auto& muonpath : outmpaths) {
0451       correlatedMetaPrimitives.emplace_back(muonpath->rawId(),
0452                                             (double)muonpath->bxTimeValue(),
0453                                             muonpath->horizPos(),
0454                                             muonpath->tanPhi(),
0455                                             muonpath->phi(),
0456                                             muonpath->phiB(),
0457                                             muonpath->phi_cmssw(),
0458                                             muonpath->phiB_cmssw(),
0459                                             muonpath->chiSquare(),
0460                                             (int)muonpath->quality(),
0461                                             muonpath->primitive(0)->channelId(),
0462                                             muonpath->primitive(0)->tdcTimeStamp(),
0463                                             muonpath->primitive(0)->laterality(),
0464                                             muonpath->primitive(1)->channelId(),
0465                                             muonpath->primitive(1)->tdcTimeStamp(),
0466                                             muonpath->primitive(1)->laterality(),
0467                                             muonpath->primitive(2)->channelId(),
0468                                             muonpath->primitive(2)->tdcTimeStamp(),
0469                                             muonpath->primitive(2)->laterality(),
0470                                             muonpath->primitive(3)->channelId(),
0471                                             muonpath->primitive(3)->tdcTimeStamp(),
0472                                             muonpath->primitive(3)->laterality(),
0473                                             muonpath->primitive(4)->channelId(),
0474                                             muonpath->primitive(4)->tdcTimeStamp(),
0475                                             muonpath->primitive(4)->laterality(),
0476                                             muonpath->primitive(5)->channelId(),
0477                                             muonpath->primitive(5)->tdcTimeStamp(),
0478                                             muonpath->primitive(5)->laterality(),
0479                                             muonpath->primitive(6)->channelId(),
0480                                             muonpath->primitive(6)->tdcTimeStamp(),
0481                                             muonpath->primitive(6)->laterality(),
0482                                             muonpath->primitive(7)->channelId(),
0483                                             muonpath->primitive(7)->tdcTimeStamp(),
0484                                             muonpath->primitive(7)->laterality());
0485     }
0486   }
0487   filteredMetaPrimitives.clear();
0488 
0489   if (debug_)
0490     LogDebug("DTTrigPhase2Prod") << "DTp2 in event:" << iEvent.id().event() << " we found "
0491                                  << correlatedMetaPrimitives.size() << " correlatedMetPrimitives (chamber)";
0492 
0493   if (dump_) {
0494     LogInfo("DTTrigPhase2Prod") << "DTp2 in event:" << iEvent.id().event() << " we found "
0495                                 << correlatedMetaPrimitives.size() << " correlatedMetPrimitives (chamber)";
0496 
0497     for (unsigned int i = 0; i < correlatedMetaPrimitives.size(); i++) {
0498       stringstream ss;
0499       ss << iEvent.id().event() << " correlated mp " << i << ": ";
0500       printmPC(ss.str(), correlatedMetaPrimitives.at(i));
0501     }
0502   }
0503 
0504   double shift_back = 0;
0505   if (scenario_ == MC)  //scope for MC
0506     shift_back = 400;
0507   else if (scenario_ == DATA)  //scope for data
0508     shift_back = 0;
0509   else if (scenario_ == SLICE_TEST)  //scope for slice test
0510     shift_back = 400;
0511 
0512   // RPC integration
0513   if (useRPC_) {
0514     rpc_integrator_->initialise(iEventSetup, shift_back);
0515     rpc_integrator_->prepareMetaPrimitives(rpcRecHits);
0516     rpc_integrator_->matchWithDTAndUseRPCTime(correlatedMetaPrimitives);
0517     rpc_integrator_->makeRPCOnlySegments();
0518     rpc_integrator_->storeRPCSingleHits();
0519     rpc_integrator_->removeRPCHitsUsed();
0520   }
0521 
0522   /// STORING RESULTs
0523   vector<L1Phase2MuDTPhDigi> outP2Ph;
0524   vector<L1Phase2MuDTExtPhDigi> outExtP2Ph;
0525   vector<L1Phase2MuDTThDigi> outP2Th;
0526   vector<L1Phase2MuDTExtThDigi> outExtP2Th;
0527 
0528   // Assigning index value
0529   assignIndex(correlatedMetaPrimitives);
0530   for (const auto& metaPrimitiveIt : correlatedMetaPrimitives) {
0531     DTChamberId chId(metaPrimitiveIt.rawId);
0532     DTSuperLayerId slId(metaPrimitiveIt.rawId);
0533     if (debug_)
0534       LogDebug("DTTrigPhase2Prod") << "looping in final vector: SuperLayerId" << chId << " x=" << metaPrimitiveIt.x
0535                                    << " quality=" << metaPrimitiveIt.quality
0536                                    << " BX=" << round(metaPrimitiveIt.t0 / 25.) << " index=" << metaPrimitiveIt.index;
0537 
0538     int sectorTP = chId.sector();
0539     //sectors 13 and 14 exist only for the outermost stations for sectors 4 and 10 respectively
0540     //due to the larger MB4 that are divided into two.
0541     if (sectorTP == 13)
0542       sectorTP = 4;
0543     if (sectorTP == 14)
0544       sectorTP = 10;
0545     sectorTP = sectorTP - 1;
0546     int sl = 0;
0547     if (metaPrimitiveIt.quality < LOWLOWQ || metaPrimitiveIt.quality == CHIGHQ) {
0548       if (inner(metaPrimitiveIt))
0549         sl = 1;
0550       else
0551         sl = 3;
0552     }
0553 
0554     if (debug_)
0555       LogDebug("DTTrigPhase2Prod") << "pushing back phase-2 dataformat carlo-federica dataformat";
0556 
0557     if (slId.superLayer() != 2) {
0558       if (df_extended_ == 1 || df_extended_ == 2) {
0559         int pathWireId[8] = {metaPrimitiveIt.wi1,
0560                              metaPrimitiveIt.wi2,
0561                              metaPrimitiveIt.wi3,
0562                              metaPrimitiveIt.wi4,
0563                              metaPrimitiveIt.wi5,
0564                              metaPrimitiveIt.wi6,
0565                              metaPrimitiveIt.wi7,
0566                              metaPrimitiveIt.wi8};
0567 
0568         int pathTDC[8] = {max((int)round(metaPrimitiveIt.tdc1 - shift_back * LHC_CLK_FREQ), -1),
0569                           max((int)round(metaPrimitiveIt.tdc2 - shift_back * LHC_CLK_FREQ), -1),
0570                           max((int)round(metaPrimitiveIt.tdc3 - shift_back * LHC_CLK_FREQ), -1),
0571                           max((int)round(metaPrimitiveIt.tdc4 - shift_back * LHC_CLK_FREQ), -1),
0572                           max((int)round(metaPrimitiveIt.tdc5 - shift_back * LHC_CLK_FREQ), -1),
0573                           max((int)round(metaPrimitiveIt.tdc6 - shift_back * LHC_CLK_FREQ), -1),
0574                           max((int)round(metaPrimitiveIt.tdc7 - shift_back * LHC_CLK_FREQ), -1),
0575                           max((int)round(metaPrimitiveIt.tdc8 - shift_back * LHC_CLK_FREQ), -1)};
0576 
0577         int pathLat[8] = {metaPrimitiveIt.lat1,
0578                           metaPrimitiveIt.lat2,
0579                           metaPrimitiveIt.lat3,
0580                           metaPrimitiveIt.lat4,
0581                           metaPrimitiveIt.lat5,
0582                           metaPrimitiveIt.lat6,
0583                           metaPrimitiveIt.lat7,
0584                           metaPrimitiveIt.lat8};
0585 
0586         // phiTP (extended DF)
0587         outExtP2Ph.emplace_back(
0588             L1Phase2MuDTExtPhDigi((int)round(metaPrimitiveIt.t0 / (float)LHC_CLK_FREQ) - shift_back,
0589                                   chId.wheel(),                                                // uwh   (m_wheel)
0590                                   sectorTP,                                                    // usc   (m_sector)
0591                                   chId.station(),                                              // ust   (m_station)
0592                                   sl,                                                          // ust   (m_station)
0593                                   (int)round(metaPrimitiveIt.phi * PHIRES_CONV),               // uphi  (m_phiAngle)
0594                                   (int)round(metaPrimitiveIt.phiB * PHIBRES_CONV),             // uphib (m_phiBending)
0595                                   metaPrimitiveIt.quality,                                     // uqua  (m_qualityCode)
0596                                   metaPrimitiveIt.index,                                       // uind  (m_segmentIndex)
0597                                   (int)round(metaPrimitiveIt.t0) - shift_back * LHC_CLK_FREQ,  // ut0   (m_t0Segment)
0598                                   (int)round(metaPrimitiveIt.chi2 * CHI2RES_CONV),             // uchi2 (m_chi2Segment)
0599                                   (int)round(metaPrimitiveIt.x * 1000),                        // ux    (m_xLocal)
0600                                   (int)round(metaPrimitiveIt.tanPhi * 1000),                   // utan  (m_tanPsi)
0601                                   (int)round(metaPrimitiveIt.phi_cmssw * PHIRES_CONV),    // uphi  (m_phiAngleCMSSW)
0602                                   (int)round(metaPrimitiveIt.phiB_cmssw * PHIBRES_CONV),  // uphib (m_phiBendingCMSSW)
0603                                   metaPrimitiveIt.rpcFlag,                                // urpc  (m_rpcFlag)
0604                                   pathWireId,
0605                                   pathTDC,
0606                                   pathLat));
0607       }
0608       if (df_extended_ == 0 || df_extended_ == 2) {
0609         // phiTP (standard DF)
0610         outP2Ph.push_back(L1Phase2MuDTPhDigi(
0611             (int)round(metaPrimitiveIt.t0 / (float)LHC_CLK_FREQ) - shift_back,
0612             chId.wheel(),                                                // uwh (m_wheel)
0613             sectorTP,                                                    // usc (m_sector)
0614             chId.station(),                                              // ust (m_station)
0615             sl,                                                          // ust (m_station)
0616             (int)round(metaPrimitiveIt.phi * PHIRES_CONV),               // uphi (_phiAngle)
0617             (int)round(metaPrimitiveIt.phiB * PHIBRES_CONV),             // uphib (m_phiBending)
0618             metaPrimitiveIt.quality,                                     // uqua (m_qualityCode)
0619             metaPrimitiveIt.index,                                       // uind (m_segmentIndex)
0620             (int)round(metaPrimitiveIt.t0) - shift_back * LHC_CLK_FREQ,  // ut0 (m_t0Segment)
0621             (int)round(metaPrimitiveIt.chi2 * CHI2RES_CONV),             // uchi2 (m_chi2Segment)
0622             metaPrimitiveIt.rpcFlag                                      // urpc (m_rpcFlag)
0623             ));
0624       }
0625     } else {
0626       if (df_extended_ == 1 || df_extended_ == 2) {
0627         int pathWireId[4] = {metaPrimitiveIt.wi1, metaPrimitiveIt.wi2, metaPrimitiveIt.wi3, metaPrimitiveIt.wi4};
0628 
0629         int pathTDC[4] = {max((int)round(metaPrimitiveIt.tdc1 - shift_back * LHC_CLK_FREQ), -1),
0630                           max((int)round(metaPrimitiveIt.tdc2 - shift_back * LHC_CLK_FREQ), -1),
0631                           max((int)round(metaPrimitiveIt.tdc3 - shift_back * LHC_CLK_FREQ), -1),
0632                           max((int)round(metaPrimitiveIt.tdc4 - shift_back * LHC_CLK_FREQ), -1)};
0633 
0634         int pathLat[4] = {metaPrimitiveIt.lat1, metaPrimitiveIt.lat2, metaPrimitiveIt.lat3, metaPrimitiveIt.lat4};
0635 
0636         // thTP (extended DF)
0637         outExtP2Th.emplace_back(
0638             L1Phase2MuDTExtThDigi((int)round(metaPrimitiveIt.t0 / (float)LHC_CLK_FREQ) - shift_back,
0639                                   chId.wheel(),                                                // uwh   (m_wheel)
0640                                   sectorTP,                                                    // usc   (m_sector)
0641                                   chId.station(),                                              // ust   (m_station)
0642                                   (int)round(metaPrimitiveIt.phi * ZRES_CONV),                 // uz    (m_zGlobal)
0643                                   (int)round(metaPrimitiveIt.phiB * KRES_CONV),                // uk    (m_kSlope)
0644                                   metaPrimitiveIt.quality,                                     // uqua  (m_qualityCode)
0645                                   metaPrimitiveIt.index,                                       // uind  (m_segmentIndex)
0646                                   (int)round(metaPrimitiveIt.t0) - shift_back * LHC_CLK_FREQ,  // ut0   (m_t0Segment)
0647                                   (int)round(metaPrimitiveIt.chi2 * CHI2RES_CONV),             // uchi2 (m_chi2Segment)
0648                                   (int)round(metaPrimitiveIt.x * 1000),                        // ux    (m_yLocal)
0649                                   (int)round(metaPrimitiveIt.phi_cmssw * PHIRES_CONV),         // uphi  (m_zCMSSW)
0650                                   (int)round(metaPrimitiveIt.phiB_cmssw * PHIBRES_CONV),       // uphib (m_kCMSSW)
0651                                   metaPrimitiveIt.rpcFlag,                                     // urpc  (m_rpcFlag)
0652                                   pathWireId,
0653                                   pathTDC,
0654                                   pathLat));
0655       }
0656       if (df_extended_ == 0 || df_extended_ == 2) {
0657         // thTP (standard DF)
0658         outP2Th.push_back(L1Phase2MuDTThDigi(
0659             (int)round(metaPrimitiveIt.t0 / (float)LHC_CLK_FREQ) - shift_back,
0660             chId.wheel(),                                                // uwh (m_wheel)
0661             sectorTP,                                                    // usc (m_sector)
0662             chId.station(),                                              // ust (m_station)
0663             (int)round(metaPrimitiveIt.phi * ZRES_CONV),                 // uz (m_zGlobal)
0664             (int)round(metaPrimitiveIt.phiB * KRES_CONV),                // uk (m_kSlope)
0665             metaPrimitiveIt.quality,                                     // uqua (m_qualityCode)
0666             metaPrimitiveIt.index,                                       // uind (m_segmentIndex)
0667             (int)round(metaPrimitiveIt.t0) - shift_back * LHC_CLK_FREQ,  // ut0 (m_t0Segment)
0668             (int)round(metaPrimitiveIt.chi2 * CHI2RES_CONV),             // uchi2 (m_chi2Segment)
0669             metaPrimitiveIt.rpcFlag                                      // urpc (m_rpcFlag)
0670             ));
0671       }
0672     }
0673   }
0674 
0675   // Storing RPC hits that were not used elsewhere
0676   if (useRPC_) {
0677     for (auto rpc_dt_digi = rpc_integrator_->rpcRecHits_translated_.begin();
0678          rpc_dt_digi != rpc_integrator_->rpcRecHits_translated_.end();
0679          rpc_dt_digi++) {
0680       outP2Ph.push_back(*rpc_dt_digi);
0681     }
0682   }
0683 
0684   // Storing Phi results
0685   if (df_extended_ == 1 || df_extended_ == 2) {
0686     std::unique_ptr<L1Phase2MuDTExtPhContainer> resultExtP2Ph(new L1Phase2MuDTExtPhContainer);
0687     resultExtP2Ph->setContainer(outExtP2Ph);
0688     iEvent.put(std::move(resultExtP2Ph));
0689   }
0690   if (df_extended_ == 0 || df_extended_ == 2) {
0691     std::unique_ptr<L1Phase2MuDTPhContainer> resultP2Ph(new L1Phase2MuDTPhContainer);
0692     resultP2Ph->setContainer(outP2Ph);
0693     iEvent.put(std::move(resultP2Ph));
0694   }
0695   outExtP2Ph.clear();
0696   outExtP2Ph.erase(outExtP2Ph.begin(), outExtP2Ph.end());
0697   outP2Ph.clear();
0698   outP2Ph.erase(outP2Ph.begin(), outP2Ph.end());
0699 
0700   // Storing Theta results
0701   if (df_extended_ == 1 || df_extended_ == 2) {
0702     std::unique_ptr<L1Phase2MuDTExtThContainer> resultExtP2Th(new L1Phase2MuDTExtThContainer);
0703     resultExtP2Th->setContainer(outExtP2Th);
0704     iEvent.put(std::move(resultExtP2Th));
0705   }
0706   if (df_extended_ == 0 || df_extended_ == 2) {
0707     std::unique_ptr<L1Phase2MuDTThContainer> resultP2Th(new L1Phase2MuDTThContainer);
0708     resultP2Th->setContainer(outP2Th);
0709     iEvent.put(std::move(resultP2Th));
0710   }
0711   outExtP2Th.clear();
0712   outExtP2Th.erase(outExtP2Th.begin(), outExtP2Th.end());
0713   outP2Th.clear();
0714   outP2Th.erase(outP2Th.begin(), outP2Th.end());
0715 }
0716 
0717 void DTTrigPhase2Prod::endRun(edm::Run const& iRun, const edm::EventSetup& iEventSetup) {
0718   grouping_obj_->finish();
0719   mpathanalyzer_->finish();
0720   mpathqualityenhancer_->finish();
0721   mpathqualityenhancerbayes_->finish();
0722   mpathredundantfilter_->finish();
0723   mpathhitsfilter_->finish();
0724   mpathassociator_->finish();
0725   rpc_integrator_->finish();
0726 };
0727 
0728 bool DTTrigPhase2Prod::outer(const metaPrimitive& mp) const {
0729   int counter = (mp.wi5 != -1) + (mp.wi6 != -1) + (mp.wi7 != -1) + (mp.wi8 != -1);
0730   return (counter > 2);
0731 }
0732 
0733 bool DTTrigPhase2Prod::inner(const metaPrimitive& mp) const {
0734   int counter = (mp.wi1 != -1) + (mp.wi2 != -1) + (mp.wi3 != -1) + (mp.wi4 != -1);
0735   return (counter > 2);
0736 }
0737 
0738 bool DTTrigPhase2Prod::hasPosRF(int wh, int sec) const { return wh > 0 || (wh == 0 && sec % 4 > 1); }
0739 
0740 void DTTrigPhase2Prod::printmP(const string& ss, const metaPrimitive& mP) const {
0741   DTSuperLayerId slId(mP.rawId);
0742   LogInfo("DTTrigPhase2Prod") << ss << (int)slId << "\t " << setw(2) << left << mP.wi1 << " " << setw(2) << left
0743                               << mP.wi2 << " " << setw(2) << left << mP.wi3 << " " << setw(2) << left << mP.wi4 << " "
0744                               << setw(5) << left << mP.tdc1 << " " << setw(5) << left << mP.tdc2 << " " << setw(5)
0745                               << left << mP.tdc3 << " " << setw(5) << left << mP.tdc4 << " " << setw(10) << right
0746                               << mP.x << " " << setw(9) << left << mP.tanPhi << " " << setw(5) << left << mP.t0 << " "
0747                               << setw(13) << left << mP.chi2 << " r:" << rango(mP);
0748 }
0749 
0750 void DTTrigPhase2Prod::printmPC(const string& ss, const metaPrimitive& mP) const {
0751   DTChamberId ChId(mP.rawId);
0752   LogInfo("DTTrigPhase2Prod") << ss << (int)ChId << "\t  " << setw(2) << left << mP.wi1 << " " << setw(2) << left
0753                               << mP.wi2 << " " << setw(2) << left << mP.wi3 << " " << setw(2) << left << mP.wi4 << " "
0754                               << setw(2) << left << mP.wi5 << " " << setw(2) << left << mP.wi6 << " " << setw(2) << left
0755                               << mP.wi7 << " " << setw(2) << left << mP.wi8 << " " << setw(5) << left << mP.tdc1 << " "
0756                               << setw(5) << left << mP.tdc2 << " " << setw(5) << left << mP.tdc3 << " " << setw(5)
0757                               << left << mP.tdc4 << " " << setw(5) << left << mP.tdc5 << " " << setw(5) << left
0758                               << mP.tdc6 << " " << setw(5) << left << mP.tdc7 << " " << setw(5) << left << mP.tdc8
0759                               << " " << setw(2) << left << mP.lat1 << " " << setw(2) << left << mP.lat2 << " "
0760                               << setw(2) << left << mP.lat3 << " " << setw(2) << left << mP.lat4 << " " << setw(2)
0761                               << left << mP.lat5 << " " << setw(2) << left << mP.lat6 << " " << setw(2) << left
0762                               << mP.lat7 << " " << setw(2) << left << mP.lat8 << " " << setw(10) << right << mP.x << " "
0763                               << setw(9) << left << mP.tanPhi << " " << setw(5) << left << mP.t0 << " " << setw(13)
0764                               << left << mP.chi2 << " r:" << rango(mP);
0765 }
0766 
0767 int DTTrigPhase2Prod::rango(const metaPrimitive& mp) const {
0768   if (mp.quality == 1 or mp.quality == 2)
0769     return 3;
0770   if (mp.quality == 3 or mp.quality == 4)
0771     return 4;
0772   return mp.quality;
0773 }
0774 
0775 void DTTrigPhase2Prod::assignIndex(std::vector<metaPrimitive>& inMPaths) {
0776   std::map<int, std::vector<metaPrimitive>> primsPerBX;
0777   for (const auto& metaPrimitive : inMPaths) {
0778     int BX = round(metaPrimitive.t0 / 25.);
0779     primsPerBX[BX].push_back(metaPrimitive);
0780   }
0781   inMPaths.clear();
0782   for (auto& prims : primsPerBX) {
0783     assignIndexPerBX(prims.second);
0784     for (const auto& primitive : prims.second)
0785       if (primitive.index <= max_index_)
0786         inMPaths.push_back(primitive);
0787   }
0788 }
0789 
0790 void DTTrigPhase2Prod::assignIndexPerBX(std::vector<metaPrimitive>& inMPaths) {
0791   // First we asociate a new index to the metaprimitive depending on quality or phiB;
0792   uint32_t rawId = -1;
0793   int numP = -1;
0794   for (auto& metaPrimitiveIt : inMPaths) {
0795     numP++;
0796     rawId = metaPrimitiveIt.rawId;
0797     int iOrder = assignQualityOrder(metaPrimitiveIt);
0798     int inf = 0;
0799     int numP2 = -1;
0800     for (auto& metaPrimitiveItN : inMPaths) {
0801       int nOrder = assignQualityOrder(metaPrimitiveItN);
0802       numP2++;
0803       if (rawId != metaPrimitiveItN.rawId)
0804         continue;
0805       if (numP2 == numP) {
0806         metaPrimitiveIt.index = inf;
0807         break;
0808       } else if (iOrder < nOrder) {
0809         inf++;
0810       } else if (iOrder > nOrder) {
0811         metaPrimitiveItN.index++;
0812       } else if (iOrder == nOrder) {
0813         if (std::abs(metaPrimitiveIt.phiB) >= std::abs(metaPrimitiveItN.phiB)) {
0814           inf++;
0815         } else if (std::abs(metaPrimitiveIt.phiB) < std::abs(metaPrimitiveItN.phiB)) {
0816           metaPrimitiveItN.index++;
0817         }
0818       }
0819     }  // ending second for
0820   }    // ending first for
0821 }
0822 
0823 int DTTrigPhase2Prod::assignQualityOrder(const metaPrimitive& mP) const {
0824   if (mP.quality > 8 || mP.quality < 1)
0825     return -1;
0826 
0827   return qmap_.find(mP.quality)->second;
0828 }
0829 
0830 std::vector<DTDigiCollection*> DTTrigPhase2Prod::distribDigis(std::queue<std::pair<DTLayerId, DTDigi>>& inQ) {
0831   std::vector<std::queue<std::pair<DTLayerId, DTDigi>>*> tmpVector;
0832   tmpVector.clear();
0833   std::vector<DTDigiCollection*> collVector;
0834   collVector.clear();
0835   while (!inQ.empty()) {
0836     processDigi(inQ, tmpVector);
0837   }
0838   for (auto& sQ : tmpVector) {
0839     DTDigiCollection tmpColl;
0840     while (!sQ->empty()) {
0841       tmpColl.insertDigi((sQ->front().first), (sQ->front().second));
0842       sQ->pop();
0843     }
0844     collVector.push_back(&tmpColl);
0845   }
0846   return collVector;
0847 }
0848 
0849 void DTTrigPhase2Prod::processDigi(std::queue<std::pair<DTLayerId, DTDigi>>& inQ,
0850                                    std::vector<std::queue<std::pair<DTLayerId, DTDigi>>*>& vec) {
0851   bool classified = false;
0852   if (!vec.empty()) {
0853     for (auto& sC : vec) {  // Conditions for entering a super cell.
0854       if ((sC->front().second.time() + superCelltimewidth_) > inQ.front().second.time()) {
0855         // Time requirement
0856         if (TMath::Abs(sC->front().second.wire() - inQ.front().second.wire()) <= superCellhalfspacewidth_) {
0857           // Spatial requirement
0858           sC->push(std::move(inQ.front()));
0859           classified = true;
0860         }
0861       }
0862     }
0863   }
0864   if (classified) {
0865     inQ.pop();
0866     return;
0867   }
0868 
0869   std::queue<std::pair<DTLayerId, DTDigi>> newQueue;
0870 
0871   std::pair<DTLayerId, DTDigi> tmpPair;
0872   tmpPair = std::move(inQ.front());
0873   newQueue.push(tmpPair);
0874   inQ.pop();
0875 
0876   vec.push_back(&newQueue);
0877 }
0878 
0879 void DTTrigPhase2Prod::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0880   // dtTriggerPhase2PrimitiveDigis
0881   edm::ParameterSetDescription desc;
0882   desc.add<edm::InputTag>("digiTag", edm::InputTag("CalibratedDigis"));
0883   desc.add<int>("trigger_with_sl", 4);
0884   desc.add<int>("timeTolerance", 999999);
0885   desc.add<double>("tanPhiTh", 1.0);
0886   desc.add<double>("tanPhiThw2max", 1.3);
0887   desc.add<double>("tanPhiThw2min", 0.5);
0888   desc.add<double>("tanPhiThw1max", 0.9);
0889   desc.add<double>("tanPhiThw1min", 0.2);
0890   desc.add<double>("tanPhiThw0", 0.5);
0891   desc.add<double>("chi2Th", 0.01);
0892   desc.add<double>("chi2corTh", 0.1);
0893   desc.add<bool>("useBX_correlation", false);
0894   desc.add<double>("dT0_correlate_TP", 25.0);
0895   desc.add<int>("dBX_correlate_TP", 0);
0896   desc.add<double>("dTanPsi_correlate_TP", 99999.0);
0897   desc.add<bool>("clean_chi2_correlation", true);
0898   desc.add<bool>("allow_confirmation", true);
0899   desc.add<double>("minx_match_2digis", 1.0);
0900   desc.add<int>("scenario", 0);
0901   desc.add<int>("df_extended", 0);
0902   desc.add<int>("max_primitives", 999);
0903   desc.add<edm::FileInPath>("ttrig_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/wire_rawId_ttrig.txt"));
0904   desc.add<edm::FileInPath>("z_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/wire_rawId_z.txt"));
0905   desc.add<edm::FileInPath>("shift_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/wire_rawId_x.txt"));
0906   desc.add<edm::FileInPath>("shift_theta_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/theta_shift.txt"));
0907   desc.add<edm::FileInPath>("global_coords_filename",
0908                             edm::FileInPath("L1Trigger/DTTriggerPhase2/data/global_coord_perp_x_phi0.txt"));
0909   desc.add<int>("algo", 0);
0910   desc.add<int>("minHits4Fit", 3);
0911   desc.add<bool>("splitPathPerSL", true);
0912   desc.addUntracked<bool>("debug", false);
0913   desc.addUntracked<bool>("dump", false);
0914   desc.add<edm::InputTag>("rpcRecHits", edm::InputTag("rpcRecHits"));
0915   desc.add<bool>("useRPC", false);
0916   desc.add<int>("bx_window", 1);
0917   desc.add<double>("phi_window", 50.0);
0918   desc.add<int>("max_quality_to_overwrite_t0", 9);
0919   desc.add<bool>("storeAllRPCHits", false);
0920   desc.add<bool>("activateBuffer", false);
0921   desc.add<double>("superCelltimewidth", 400);
0922   desc.add<int>("superCellspacewidth", 20);
0923   {
0924     edm::ParameterSetDescription psd0;
0925     psd0.addUntracked<bool>("debug", false);
0926     psd0.add<double>("angletan", 0.3);
0927     psd0.add<double>("anglebinwidth", 1.0);
0928     psd0.add<double>("posbinwidth", 2.1);
0929     psd0.add<double>("maxdeltaAngDeg", 10);
0930     psd0.add<double>("maxdeltaPos", 10);
0931     psd0.add<int>("UpperNumber", 6);
0932     psd0.add<int>("LowerNumber", 4);
0933     psd0.add<double>("MaxDistanceToWire", 0.03);
0934     psd0.add<int>("minNLayerHits", 6);
0935     psd0.add<int>("minSingleSLHitsMax", 3);
0936     psd0.add<int>("minSingleSLHitsMin", 3);
0937     psd0.add<bool>("allowUncorrelatedPatterns", true);
0938     psd0.add<int>("minUncorrelatedHits", 3);
0939     desc.add<edm::ParameterSetDescription>("HoughGrouping", psd0);
0940   }
0941   {
0942     edm::ParameterSetDescription psd0;
0943     psd0.add<edm::FileInPath>(
0944         "pattern_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/PseudoBayesPatterns_uncorrelated_v0.root"));
0945     psd0.addUntracked<bool>("debug", false);
0946     psd0.add<int>("minNLayerHits", 3);
0947     psd0.add<int>("minSingleSLHitsMax", 3);
0948     psd0.add<int>("minSingleSLHitsMin", 0);
0949     psd0.add<int>("allowedVariance", 1);
0950     psd0.add<bool>("allowDuplicates", false);
0951     psd0.add<bool>("setLateralities", true);
0952     psd0.add<bool>("allowUncorrelatedPatterns", true);
0953     psd0.add<int>("minUncorrelatedHits", 3);
0954     psd0.add<bool>("saveOnPlace", true);
0955     psd0.add<int>("maxPathsPerMatch", 256);
0956     desc.add<edm::ParameterSetDescription>("PseudoBayesPattern", psd0);
0957   }
0958   descriptions.add("dtTriggerPhase2PrimitiveDigis", desc);
0959 }
0960 
0961 DEFINE_FWK_MODULE(DTTrigPhase2Prod);