Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-04-15 22:41:36

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/TrapezoidalGrouping.h"
0027 #include "L1Trigger/DTTriggerPhase2/interface/HoughGrouping.h"
0028 #include "L1Trigger/DTTriggerPhase2/interface/PseudoBayesGrouping.h"
0029 #include "L1Trigger/DTTriggerPhase2/interface/LateralityProvider.h"
0030 #include "L1Trigger/DTTriggerPhase2/interface/LateralityBasicProvider.h"
0031 #include "L1Trigger/DTTriggerPhase2/interface/LateralityCoarsedProvider.h"
0032 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAnalyzer.h"
0033 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathSLFitter.h"
0034 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathCorFitter.h"
0035 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAnalyticAnalyzer.h"
0036 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAnalyzerInChamber.h"
0037 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAssociator.h"
0038 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathConfirmator.h"
0039 #include "L1Trigger/DTTriggerPhase2/interface/MPFilter.h"
0040 #include "L1Trigger/DTTriggerPhase2/interface/MPSLFilter.h"
0041 #include "L1Trigger/DTTriggerPhase2/interface/MPCorFilter.h"
0042 #include "L1Trigger/DTTriggerPhase2/interface/MPCoincidenceFilter.h"
0043 #include "L1Trigger/DTTriggerPhase2/interface/MPThetaMatching.h"
0044 #include "L1Trigger/DTTriggerPhase2/interface/MPQualityEnhancerFilter.h"
0045 #include "L1Trigger/DTTriggerPhase2/interface/MPRedundantFilter.h"
0046 #include "L1Trigger/DTTriggerPhase2/interface/MPCleanHitsFilter.h"
0047 #include "L1Trigger/DTTriggerPhase2/interface/MPQualityEnhancerFilterBayes.h"
0048 #include "L1Trigger/DTTriggerPhase2/interface/GlobalCoordsObtainer.h"
0049 
0050 #include "DataFormats/MuonDetId/interface/DTChamberId.h"
0051 #include "DataFormats/MuonDetId/interface/DTSuperLayerId.h"
0052 #include "DataFormats/MuonDetId/interface/DTLayerId.h"
0053 #include "DataFormats/MuonDetId/interface/DTWireId.h"
0054 #include "DataFormats/DTDigi/interface/DTDigiCollection.h"
0055 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTPhContainer.h"
0056 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTPhDigi.h"
0057 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTExtPhContainer.h"
0058 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTExtPhDigi.h"
0059 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTThContainer.h"
0060 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTThDigi.h"
0061 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTExtThContainer.h"
0062 #include "DataFormats/L1DTTrackFinder/interface/L1Phase2MuDTExtThDigi.h"
0063 
0064 // DT trigger GeomUtils
0065 #include "DQM/DTMonitorModule/interface/DTTrigGeomUtils.h"
0066 
0067 //RPC TP
0068 #include "DataFormats/RPCRecHit/interface/RPCRecHitCollection.h"
0069 #include <DataFormats/MuonDetId/interface/RPCDetId.h>
0070 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"
0071 #include "L1Trigger/DTTriggerPhase2/interface/RPCIntegrator.h"
0072 
0073 #include <fstream>
0074 #include <iostream>
0075 #include <queue>
0076 #include <cmath>
0077 
0078 using namespace edm;
0079 using namespace std;
0080 using namespace cmsdt;
0081 
0082 class DTTrigPhase2Prod : public edm::stream::EDProducer<> {
0083   typedef std::map<DTChamberId, DTDigiCollection, std::less<DTChamberId>> DTDigiMap;
0084   typedef DTDigiMap::iterator DTDigiMap_iterator;
0085   typedef DTDigiMap::const_iterator DTDigiMap_const_iterator;
0086 
0087 public:
0088   //! Constructor
0089   DTTrigPhase2Prod(const edm::ParameterSet& pset);
0090 
0091   //! Destructor
0092   ~DTTrigPhase2Prod() override;
0093 
0094   //! Create Trigger Units before starting event processing
0095   void beginRun(edm::Run const& iRun, const edm::EventSetup& iEventSetup) override;
0096 
0097   //! Producer: process every event and generates trigger data
0098   void produce(edm::Event& iEvent, const edm::EventSetup& iEventSetup) override;
0099 
0100   //! endRun: finish things
0101   void endRun(edm::Run const& iRun, const edm::EventSetup& iEventSetup) override;
0102 
0103   // Methods
0104   int rango(const metaPrimitive& mp) const;
0105   bool outer(const metaPrimitive& mp) const;
0106   bool inner(const metaPrimitive& mp) const;
0107   void printmP(const std::string& ss, const metaPrimitive& mP) const;
0108   void printmP(const metaPrimitive& mP) const;
0109   void printmPC(const std::string& ss, const metaPrimitive& mP) const;
0110   void printmPC(const metaPrimitive& mP) const;
0111   bool hasPosRF(int wh, int sec) const;
0112 
0113   // Getter-methods
0114   MP_QUALITY getMinimumQuality(void);
0115 
0116   // Setter-methods
0117   void setChiSquareThreshold(float ch2Thr);
0118   void setMinimumQuality(MP_QUALITY q);
0119 
0120   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0121 
0122   // data-members
0123   const DTGeometry* dtGeo_;
0124   edm::ESGetToken<DTGeometry, MuonGeometryRecord> dtGeomH;
0125   std::vector<std::pair<int, MuonPath>> primitives_;
0126 
0127 private:
0128   // Trigger Configuration Manager CCB validity flag
0129   bool my_CCBValid_;
0130 
0131   // BX offset used to correct DTTPG output
0132   int my_BXoffset_;
0133 
0134   // Debug Flag
0135   bool debug_;
0136   bool dump_;
0137   double dT0_correlate_TP_;
0138   int scenario_;
0139   int df_extended_;
0140   int co_option_;  //coincidence
0141   int co_quality_;
0142   int co_wh2option_;
0143   int th_option_;  //theta matching
0144   int th_quality_;
0145   int max_index_;
0146 
0147   bool output_mixer_;
0148   bool output_latpredictor_;
0149   bool output_slfitter_;
0150   bool output_slfilter_;
0151   bool output_confirmed_;
0152   bool output_matcher_;
0153   bool skip_processing_;
0154   bool allow_confirmation_;
0155 
0156   // ParameterSet
0157   edm::EDGetTokenT<DTDigiCollection> dtDigisToken_;
0158   edm::EDGetTokenT<RPCRecHitCollection> rpcRecHitsLabel_;
0159 
0160   // Grouping attributes and methods
0161   int algo_;  // Grouping code
0162   std::unique_ptr<MotherGrouping> grouping_obj_;
0163   std::unique_ptr<MuonPathAnalyzer> mpathanalyzer_;
0164   std::unique_ptr<LateralityProvider> latprovider_;
0165   std::unique_ptr<MPFilter> mpathqualityenhancer_;
0166   std::unique_ptr<MPFilter> mpathqualityenhancerbayes_;
0167   std::unique_ptr<MPFilter> mpathredundantfilter_;
0168   std::unique_ptr<MPFilter> mpathhitsfilter_;
0169   std::unique_ptr<MuonPathAnalyzer> mpathassociator_;
0170   std::unique_ptr<MuonPathConfirmator> mpathconfirmator_;
0171   std::unique_ptr<MPFilter> mpathcorfilter_;
0172   std::unique_ptr<MPFilter> mpathcoifilter_;
0173   std::unique_ptr<MPFilter> mpaththetamatching_;
0174   std::shared_ptr<GlobalCoordsObtainer> globalcoordsobtainer_;
0175 
0176   // Buffering
0177   bool activateBuffer_;
0178   int superCellhalfspacewidth_;
0179   float superCelltimewidth_;
0180   std::vector<DTDigiCollection*> distribDigis(std::queue<std::pair<DTLayerId, DTDigi>>& inQ);
0181   void processDigi(std::queue<std::pair<DTLayerId, DTDigi>>& inQ,
0182                    std::vector<std::queue<std::pair<DTLayerId, DTDigi>>*>& vec);
0183 
0184   // RPC
0185   std::unique_ptr<RPCIntegrator> rpc_integrator_;
0186   bool useRPC_;
0187 
0188   void assignIndex(std::vector<metaPrimitive>& inMPaths);
0189   void assignIndexPerBX(std::vector<metaPrimitive>& inMPaths);
0190   int assignQualityOrder(const metaPrimitive& mP) const;
0191 
0192   const std::unordered_map<int, int> qmap_;
0193 };
0194 
0195 namespace {
0196   struct {
0197     bool operator()(std::pair<DTLayerId, DTDigi> a, std::pair<DTLayerId, DTDigi> b) const {
0198       return (a.second.time() < b.second.time());
0199     }
0200   } const DigiTimeOrdering;
0201 }  // namespace
0202 
0203 DTTrigPhase2Prod::DTTrigPhase2Prod(const ParameterSet& pset)
0204     : qmap_({{8, 8}, {7, 7}, {6, 6}, {4, 4}, {3, 3}, {2, 2}, {1, 1}}) {
0205   produces<L1Phase2MuDTPhContainer>();
0206   produces<L1Phase2MuDTThContainer>();
0207   produces<L1Phase2MuDTExtPhContainer>();
0208   produces<L1Phase2MuDTExtThContainer>();
0209 
0210   debug_ = pset.getUntrackedParameter<bool>("debug");
0211   dump_ = pset.getUntrackedParameter<bool>("dump");
0212 
0213   scenario_ = pset.getParameter<int>("scenario");
0214 
0215   df_extended_ = pset.getParameter<int>("df_extended");
0216   co_option_ = pset.getParameter<int>("co_option");
0217   co_quality_ = pset.getParameter<int>("co_quality");
0218   co_wh2option_ = pset.getParameter<int>("co_wh2option");
0219   th_option_ = pset.getParameter<int>("th_option");
0220   th_quality_ = pset.getParameter<int>("th_quality");
0221   max_index_ = pset.getParameter<int>("max_primitives") - 1;
0222 
0223   dtDigisToken_ = consumes<DTDigiCollection>(pset.getParameter<edm::InputTag>("digiTag"));
0224 
0225   rpcRecHitsLabel_ = consumes<RPCRecHitCollection>(pset.getParameter<edm::InputTag>("rpcRecHits"));
0226   useRPC_ = pset.getParameter<bool>("useRPC");
0227 
0228   // Choosing grouping scheme:
0229   algo_ = pset.getParameter<int>("algo");
0230 
0231   // shortcuts
0232 
0233   output_mixer_ = pset.getParameter<bool>("output_mixer");
0234   output_latpredictor_ = pset.getParameter<bool>("output_latpredictor");
0235   output_slfitter_ = pset.getParameter<bool>("output_slfitter");
0236   output_slfilter_ = pset.getParameter<bool>("output_slfilter");
0237   output_confirmed_ = pset.getParameter<bool>("output_confirmed");
0238   output_matcher_ = pset.getParameter<bool>("output_matcher");
0239   allow_confirmation_ = pset.getParameter<bool>("allow_confirmation");
0240 
0241   edm::ConsumesCollector consumesColl(consumesCollector());
0242   globalcoordsobtainer_ = std::make_shared<GlobalCoordsObtainer>(pset);
0243   globalcoordsobtainer_->generate_luts();
0244 
0245   if (algo_ == PseudoBayes) {
0246     grouping_obj_ =
0247         std::make_unique<PseudoBayesGrouping>(pset.getParameter<edm::ParameterSet>("PseudoBayesPattern"), consumesColl);
0248   } else if (algo_ == HoughTrans) {
0249     grouping_obj_ =
0250         std::make_unique<HoughGrouping>(pset.getParameter<edm::ParameterSet>("HoughGrouping"), consumesColl);
0251   } else {
0252     grouping_obj_ = std::make_unique<TrapezoidalGrouping>(pset, consumesColl);
0253   }
0254 
0255   if (algo_ == Standard) {
0256     if (debug_)
0257       LogDebug("DTTrigPhase2Prod") << "DTp2:constructor: JM analyzer";
0258     mpathanalyzer_ = std::make_unique<MuonPathSLFitter>(pset, consumesColl, globalcoordsobtainer_);
0259     latprovider_ = std::make_unique<LateralityCoarsedProvider>(pset, consumesColl);
0260   } else {
0261     if (debug_)
0262       LogDebug("DTTrigPhase2Prod") << "DTp2:constructor: Full chamber analyzer";
0263     mpathanalyzer_ = std::make_unique<MuonPathAnalyzerInChamber>(pset, consumesColl, globalcoordsobtainer_);
0264   }
0265 
0266   // Getting buffer option
0267   activateBuffer_ = pset.getParameter<bool>("activateBuffer");
0268   superCellhalfspacewidth_ = pset.getParameter<int>("superCellspacewidth") / 2;
0269   superCelltimewidth_ = pset.getParameter<double>("superCelltimewidth");
0270 
0271   mpathqualityenhancer_ = std::make_unique<MPSLFilter>(pset);
0272   mpathqualityenhancerbayes_ = std::make_unique<MPQualityEnhancerFilterBayes>(pset);
0273   mpathredundantfilter_ = std::make_unique<MPRedundantFilter>(pset);
0274   mpathhitsfilter_ = std::make_unique<MPCleanHitsFilter>(pset);
0275   mpathconfirmator_ = std::make_unique<MuonPathConfirmator>(pset, consumesColl);
0276   mpathassociator_ = std::make_unique<MuonPathCorFitter>(pset, consumesColl, globalcoordsobtainer_);
0277   mpathcorfilter_ = std::make_unique<MPCorFilter>(pset);
0278   mpathcoifilter_ = std::make_unique<MPCoincidenceFilter>(pset);
0279   mpaththetamatching_ = std::make_unique<MPThetaMatching>(pset);
0280   rpc_integrator_ = std::make_unique<RPCIntegrator>(pset, consumesColl);
0281 
0282   dtGeomH = esConsumes<DTGeometry, MuonGeometryRecord, edm::Transition::BeginRun>();
0283 }
0284 
0285 DTTrigPhase2Prod::~DTTrigPhase2Prod() {
0286   if (debug_)
0287     LogDebug("DTTrigPhase2Prod") << "DTp2: calling destructor" << std::endl;
0288 }
0289 
0290 void DTTrigPhase2Prod::beginRun(edm::Run const& iRun, const edm::EventSetup& iEventSetup) {
0291   if (debug_)
0292     LogDebug("DTTrigPhase2Prod") << "beginRun " << iRun.id().run();
0293   if (debug_)
0294     LogDebug("DTTrigPhase2Prod") << "beginRun: getting DT geometry";
0295 
0296   grouping_obj_->initialise(iEventSetup);               // Grouping object initialisation
0297   mpathanalyzer_->initialise(iEventSetup);              // Analyzer object initialisation
0298   mpathqualityenhancer_->initialise(iEventSetup);       // Filter object initialisation
0299   mpathredundantfilter_->initialise(iEventSetup);       // Filter object initialisation
0300   mpathqualityenhancerbayes_->initialise(iEventSetup);  // Filter object initialisation
0301   mpathhitsfilter_->initialise(iEventSetup);
0302   mpathassociator_->initialise(iEventSetup);  // Associator object initialisation
0303   mpathcorfilter_->initialise(iEventSetup);
0304   mpathcoifilter_->initialise(iEventSetup);
0305   mpaththetamatching_->initialise(iEventSetup);
0306 
0307   if (auto geom = iEventSetup.getHandle(dtGeomH)) {
0308     dtGeo_ = &(*geom);
0309   }
0310 }
0311 
0312 void DTTrigPhase2Prod::produce(Event& iEvent, const EventSetup& iEventSetup) {
0313   if (debug_)
0314     LogDebug("DTTrigPhase2Prod") << "produce";
0315   edm::Handle<DTDigiCollection> dtdigis;
0316   iEvent.getByToken(dtDigisToken_, dtdigis);
0317 
0318   if (debug_)
0319     LogDebug("DTTrigPhase2Prod") << "\t Getting the RPC RecHits" << std::endl;
0320   edm::Handle<RPCRecHitCollection> rpcRecHits;
0321   iEvent.getByToken(rpcRecHitsLabel_, rpcRecHits);
0322 
0323   ////////////////////////////////
0324   // GROUPING CODE:
0325   ////////////////////////////////
0326 
0327   DTDigiMap digiMap;
0328   DTDigiCollection::DigiRangeIterator detUnitIt;
0329   for (const auto& detUnitIt : *dtdigis) {
0330     const DTLayerId& layId = detUnitIt.first;
0331     const DTChamberId chambId = layId.superlayerId().chamberId();
0332     const DTDigiCollection::Range& range = detUnitIt.second;
0333     digiMap[chambId].put(range, layId);
0334   }
0335 
0336   // generate a list muon paths for each event!!!
0337   if (debug_ && activateBuffer_)
0338     LogDebug("DTTrigPhase2Prod") << "produce - Getting and grouping digis per chamber using a buffer and super cells.";
0339   else if (debug_)
0340     LogDebug("DTTrigPhase2Prod") << "produce - Getting and grouping digis per chamber.";
0341 
0342   std::map<int, MuonPathPtrs> muonpaths;
0343   for (const auto& ich : dtGeo_->chambers()) {
0344     // 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.
0345     const DTChamber* chamb = ich;
0346     DTChamberId chid = chamb->id();
0347     DTDigiMap_iterator dmit = digiMap.find(chid);
0348 
0349     if (dmit == digiMap.end())
0350       continue;
0351 
0352     if (activateBuffer_) {  // Use buffering (per chamber) or not
0353       // Import digis from the station
0354       std::vector<std::pair<DTLayerId, DTDigi>> tmpvec;
0355       tmpvec.clear();
0356 
0357       for (const auto& dtLayerIdIt : (*dmit).second) {
0358         for (DTDigiCollection::const_iterator digiIt = (dtLayerIdIt.second).first;
0359              digiIt != (dtLayerIdIt.second).second;
0360              digiIt++) {
0361           tmpvec.emplace_back(dtLayerIdIt.first, *digiIt);
0362         }
0363       }
0364 
0365       // Check to enhance CPU time usage
0366       if (tmpvec.empty())
0367         continue;
0368 
0369       // Order digis depending on TDC time and insert them into a queue (FIFO buffer). TODO: adapt for MC simulations.
0370       std::sort(tmpvec.begin(), tmpvec.end(), DigiTimeOrdering);
0371       std::queue<std::pair<DTLayerId, DTDigi>> timequeue;
0372 
0373       for (const auto& elem : tmpvec)
0374         timequeue.emplace(elem);
0375       tmpvec.clear();
0376 
0377       // Distribute the digis from the queue into supercells
0378       std::vector<DTDigiCollection*> superCells;
0379       superCells = distribDigis(timequeue);
0380 
0381       // Process each supercell & collect the resulting muonpaths (as the muonpaths std::vector is only enlarged each time
0382       // the groupings access it, it's not needed to "collect" the final products).
0383 
0384       while (!superCells.empty()) {
0385         grouping_obj_->run(iEvent, iEventSetup, *(superCells.back()), muonpaths[chid.rawId()]);
0386         superCells.pop_back();
0387       }
0388     } else {
0389       grouping_obj_->run(iEvent, iEventSetup, (*dmit).second, muonpaths[chid.rawId()]);
0390     }
0391   }
0392   digiMap.clear();
0393 
0394   if (dump_) {
0395     for (auto& ch_muonpaths : muonpaths) {
0396       for (unsigned int i = 0; i < ch_muonpaths.second.size(); i++) {
0397         stringstream ss;
0398         ss << iEvent.id().event() << "      mpath " << i << ": ";
0399         for (int lay = 0; lay < ch_muonpaths.second.at(i)->nprimitives(); lay++)
0400           ss << ch_muonpaths.second.at(i)->primitive(lay)->channelId() << " ";
0401         for (int lay = 0; lay < ch_muonpaths.second.at(i)->nprimitives(); lay++)
0402           ss << ch_muonpaths.second.at(i)->primitive(lay)->tdcTimeStamp() << " ";
0403         for (int lay = 0; lay < ch_muonpaths.second.at(i)->nprimitives(); lay++)
0404           ss << ch_muonpaths.second.at(i)->primitive(lay)->laterality() << " ";
0405         LogInfo("DTTrigPhase2Prod") << ss.str();
0406       }
0407     }
0408   }
0409 
0410   std::map<int, std::vector<lat_vector>> lateralities;
0411   if (!output_mixer_) {
0412     for (auto& ch_muonpaths : muonpaths) {
0413       if (algo_ == Standard) {
0414         latprovider_->run(iEvent, iEventSetup, ch_muonpaths.second, lateralities[ch_muonpaths.first]);
0415       }
0416     }
0417   }
0418 
0419   // FILTER GROUPING
0420   std::map<int, MuonPathPtrs> filteredmuonpaths;
0421   for (auto& ch_muonpaths : muonpaths) {
0422     if (algo_ == Standard) {
0423       mpathredundantfilter_->run(iEvent, iEventSetup, ch_muonpaths.second, filteredmuonpaths[ch_muonpaths.first]);
0424     } else {
0425       mpathhitsfilter_->run(iEvent, iEventSetup, ch_muonpaths.second, filteredmuonpaths[ch_muonpaths.first]);
0426     }
0427   }
0428 
0429   if (dump_) {
0430     for (auto& ch_filteredmuonpaths : filteredmuonpaths) {
0431       for (unsigned int i = 0; i < ch_filteredmuonpaths.second.size(); i++) {
0432         stringstream ss;
0433         ss << iEvent.id().event() << " filt. mpath " << i << ": ";
0434         for (int lay = 0; lay < ch_filteredmuonpaths.second.at(i)->nprimitives(); lay++)
0435           ss << ch_filteredmuonpaths.second.at(i)->primitive(lay)->channelId() << " ";
0436         for (int lay = 0; lay < ch_filteredmuonpaths.second.at(i)->nprimitives(); lay++)
0437           ss << ch_filteredmuonpaths.second.at(i)->primitive(lay)->tdcTimeStamp() << " ";
0438         LogInfo("DTTrigPhase2Prod") << ss.str();
0439       }
0440     }
0441   }
0442 
0443   skip_processing_ = output_mixer_ || output_latpredictor_;
0444 
0445   ///////////////////////////////////////////
0446   /// Fitting SECTION;
0447   ///////////////////////////////////////////
0448 
0449   if (debug_) {
0450     for (auto& ch_muonpaths : muonpaths) {
0451       LogDebug("DTTrigPhase2Prod") << "MUON PATHS found: " << ch_muonpaths.second.size() << " ("
0452                                    << filteredmuonpaths[ch_muonpaths.first].size() << ") in event "
0453                                    << iEvent.id().event();
0454     }
0455   }
0456   if (debug_)
0457     LogDebug("DTTrigPhase2Prod") << "filling NmetaPrimtives" << std::endl;
0458   std::map<int, std::vector<metaPrimitive>> metaPrimitives;
0459   std::map<int, MuonPathPtrs> outmpaths;
0460   if (algo_ == Standard) {
0461     if (debug_)
0462       LogDebug("DTTrigPhase2Prod") << "Fitting 1SL ";
0463     for (auto& ch_muonpaths : muonpaths) {  // FIXME, do we need filtered muonpaths?
0464       if (!output_mixer_ && !output_latpredictor_)
0465         mpathanalyzer_->run(iEvent,
0466                             iEventSetup,
0467                             ch_muonpaths.second,
0468                             lateralities[ch_muonpaths.first],
0469                             metaPrimitives[ch_muonpaths.first]);
0470       else if (output_mixer_) {
0471         for (auto& inMPath : ch_muonpaths.second) {
0472           auto sl = inMPath->primitive(0)->superLayerId();  // 0, 1, 2
0473           int selected_lay = 1;
0474           if (inMPath->primitive(0)->tdcTimeStamp() != -1)
0475             selected_lay = 0;
0476           int dumLayId = inMPath->primitive(selected_lay)->cameraId();
0477           auto dtDumlayerId = DTLayerId(dumLayId);
0478           DTSuperLayerId MuonPathSLId(dtDumlayerId.wheel(), dtDumlayerId.station(), dtDumlayerId.sector(), sl + 1);
0479           if (sl == 0)
0480             metaPrimitives[ch_muonpaths.first].emplace_back(metaPrimitive({MuonPathSLId.rawId(),
0481                                                                            -1,
0482                                                                            -1,
0483                                                                            -1,
0484                                                                            -1,
0485                                                                            -1,
0486                                                                            -1,
0487                                                                            -1,
0488                                                                            -1,
0489                                                                            -1,
0490                                                                            inMPath->primitive(0)->channelId(),
0491                                                                            inMPath->primitive(0)->tdcTimeStamp(),
0492                                                                            -1,
0493                                                                            inMPath->primitive(1)->channelId(),
0494                                                                            inMPath->primitive(1)->tdcTimeStamp(),
0495                                                                            -1,
0496                                                                            inMPath->primitive(2)->channelId(),
0497                                                                            inMPath->primitive(2)->tdcTimeStamp(),
0498                                                                            -1,
0499                                                                            inMPath->primitive(3)->channelId(),
0500                                                                            inMPath->primitive(3)->tdcTimeStamp(),
0501                                                                            -1,
0502                                                                            -1,
0503                                                                            -1,
0504                                                                            -1,
0505                                                                            -1,
0506                                                                            -1,
0507                                                                            -1,
0508                                                                            -1,
0509                                                                            -1,
0510                                                                            -1,
0511                                                                            -1,
0512                                                                            -1,
0513                                                                            -1,
0514                                                                            -1}));
0515           else
0516             metaPrimitives[ch_muonpaths.first].emplace_back(metaPrimitive({MuonPathSLId.rawId(),
0517                                                                            -1,
0518                                                                            -1,
0519                                                                            -1,
0520                                                                            -1,
0521                                                                            -1,
0522                                                                            -1,
0523                                                                            -1,
0524                                                                            -1,
0525                                                                            -1,
0526                                                                            -1,
0527                                                                            -1,
0528                                                                            -1,
0529                                                                            -1,
0530                                                                            -1,
0531                                                                            -1,
0532                                                                            -1,
0533                                                                            -1,
0534                                                                            -1,
0535                                                                            -1,
0536                                                                            -1,
0537                                                                            -1,
0538                                                                            inMPath->primitive(0)->channelId(),
0539                                                                            inMPath->primitive(0)->tdcTimeStamp(),
0540                                                                            -1,
0541                                                                            inMPath->primitive(1)->channelId(),
0542                                                                            inMPath->primitive(1)->tdcTimeStamp(),
0543                                                                            -1,
0544                                                                            inMPath->primitive(2)->channelId(),
0545                                                                            inMPath->primitive(2)->tdcTimeStamp(),
0546                                                                            -1,
0547                                                                            inMPath->primitive(3)->channelId(),
0548                                                                            inMPath->primitive(3)->tdcTimeStamp(),
0549                                                                            -1,
0550                                                                            -1}));
0551         }
0552       } else if (output_latpredictor_) {
0553         int imp = -1;
0554         for (auto& inMPath : ch_muonpaths.second) {
0555           imp++;
0556           auto sl = inMPath->primitive(0)->superLayerId();  // 0, 1, 2
0557           int selected_lay = 1;
0558           if (inMPath->primitive(0)->tdcTimeStamp() != -1)
0559             selected_lay = 0;
0560           int dumLayId = inMPath->primitive(selected_lay)->cameraId();
0561           auto dtDumlayerId = DTLayerId(dumLayId);
0562           DTSuperLayerId MuonPathSLId(dtDumlayerId.wheel(), dtDumlayerId.station(), dtDumlayerId.sector(), sl + 1);
0563           for (auto& latcomb : lateralities[ch_muonpaths.first][imp]) {
0564             if (sl == 0)
0565               metaPrimitives[ch_muonpaths.first].emplace_back(metaPrimitive({MuonPathSLId.rawId(),
0566                                                                              -1,
0567                                                                              -1,
0568                                                                              -1,
0569                                                                              -1,
0570                                                                              -1,
0571                                                                              -1,
0572                                                                              -1,
0573                                                                              -1,
0574                                                                              -1,
0575                                                                              inMPath->primitive(0)->channelId(),
0576                                                                              inMPath->primitive(0)->tdcTimeStamp(),
0577                                                                              latcomb[0],
0578                                                                              inMPath->primitive(1)->channelId(),
0579                                                                              inMPath->primitive(1)->tdcTimeStamp(),
0580                                                                              latcomb[1],
0581                                                                              inMPath->primitive(2)->channelId(),
0582                                                                              inMPath->primitive(2)->tdcTimeStamp(),
0583                                                                              latcomb[2],
0584                                                                              inMPath->primitive(3)->channelId(),
0585                                                                              inMPath->primitive(3)->tdcTimeStamp(),
0586                                                                              latcomb[3],
0587                                                                              -1,
0588                                                                              -1,
0589                                                                              -1,
0590                                                                              -1,
0591                                                                              -1,
0592                                                                              -1,
0593                                                                              -1,
0594                                                                              -1,
0595                                                                              -1,
0596                                                                              -1,
0597                                                                              -1,
0598                                                                              -1,
0599                                                                              -1}));
0600             else
0601               metaPrimitives[ch_muonpaths.first].emplace_back(metaPrimitive({MuonPathSLId.rawId(),
0602                                                                              -1,
0603                                                                              -1,
0604                                                                              -1,
0605                                                                              -1,
0606                                                                              -1,
0607                                                                              -1,
0608                                                                              -1,
0609                                                                              -1,
0610                                                                              -1,
0611                                                                              -1,
0612                                                                              -1,
0613                                                                              -1,
0614                                                                              -1,
0615                                                                              -1,
0616                                                                              -1,
0617                                                                              -1,
0618                                                                              -1,
0619                                                                              -1,
0620                                                                              -1,
0621                                                                              -1,
0622                                                                              -1,
0623                                                                              inMPath->primitive(0)->channelId(),
0624                                                                              inMPath->primitive(0)->tdcTimeStamp(),
0625                                                                              latcomb[0],
0626                                                                              inMPath->primitive(1)->channelId(),
0627                                                                              inMPath->primitive(1)->tdcTimeStamp(),
0628                                                                              latcomb[1],
0629                                                                              inMPath->primitive(2)->channelId(),
0630                                                                              inMPath->primitive(2)->tdcTimeStamp(),
0631                                                                              latcomb[2],
0632                                                                              inMPath->primitive(3)->channelId(),
0633                                                                              inMPath->primitive(3)->tdcTimeStamp(),
0634                                                                              latcomb[3],
0635                                                                              -1}));
0636           }
0637         }
0638       }
0639     }
0640   } else {
0641     // implementation for advanced (2SL) grouping, no filter required..
0642     if (debug_)
0643       LogDebug("DTTrigPhase2Prod") << "Fitting 2SL at once ";
0644     for (auto& ch_muonpaths : muonpaths) {
0645       mpathanalyzer_->run(iEvent, iEventSetup, ch_muonpaths.second, outmpaths[ch_muonpaths.first]);
0646     }
0647   }
0648 
0649   skip_processing_ = skip_processing_ || output_slfitter_;
0650 
0651   if (dump_) {
0652     for (auto& ch_outmpaths : outmpaths) {
0653       for (unsigned int i = 0; i < ch_outmpaths.second.size(); i++) {
0654         LogInfo("DTTrigPhase2Prod") << iEvent.id().event() << " mp " << i << ": "
0655                                     << ch_outmpaths.second.at(i)->bxTimeValue() << " "
0656                                     << ch_outmpaths.second.at(i)->horizPos() << " "
0657                                     << ch_outmpaths.second.at(i)->tanPhi() << " " << ch_outmpaths.second.at(i)->phi()
0658                                     << " " << ch_outmpaths.second.at(i)->phiB() << " "
0659                                     << ch_outmpaths.second.at(i)->quality() << " "
0660                                     << ch_outmpaths.second.at(i)->chiSquare();
0661       }
0662     }
0663     for (auto& ch_metaPrimitives : metaPrimitives) {
0664       for (unsigned int i = 0; i < ch_metaPrimitives.second.size(); i++) {
0665         stringstream ss;
0666         ss << iEvent.id().event() << " mp " << i << ": ";
0667         printmP(ss.str(), ch_metaPrimitives.second.at(i));
0668       }
0669     }
0670   }
0671 
0672   muonpaths.clear();
0673   filteredmuonpaths.clear();
0674 
0675   /////////////////////////////////////
0676   //// CONFIRMATION:
0677   /////////////////////////////////////
0678 
0679   std::map<int, std::vector<metaPrimitive>> confirmedMetaPrimitives;
0680   for (auto& ch_metaPrimitives : metaPrimitives) {
0681     if (!skip_processing_ && allow_confirmation_)
0682       mpathconfirmator_->run(
0683           iEvent, iEventSetup, ch_metaPrimitives.second, dtdigis, confirmedMetaPrimitives[ch_metaPrimitives.first]);
0684     else
0685       for (auto& mp : ch_metaPrimitives.second) {
0686         confirmedMetaPrimitives[ch_metaPrimitives.first].push_back(mp);
0687       }
0688   }
0689 
0690   metaPrimitives.clear();
0691   skip_processing_ = skip_processing_ || output_confirmed_;
0692 
0693   /////////////////////////////////////
0694   //  FILTER SECTIONS:
0695   ////////////////////////////////////
0696 
0697   if (debug_)
0698     LogDebug("DTTrigPhase2Prod") << "declaring new vector for filtered" << std::endl;
0699 
0700   std::map<int, std::vector<metaPrimitive>> filteredMetaPrimitives;
0701   if (algo_ == Standard)
0702     for (auto& ch_confirmedMetaPrimitives : confirmedMetaPrimitives) {
0703       if (!skip_processing_)
0704         mpathqualityenhancer_->run(iEvent,
0705                                    iEventSetup,
0706                                    ch_confirmedMetaPrimitives.second,
0707                                    filteredMetaPrimitives[ch_confirmedMetaPrimitives.first]);
0708       else
0709         for (auto& mp : ch_confirmedMetaPrimitives.second) {
0710           filteredMetaPrimitives[ch_confirmedMetaPrimitives.first].push_back(mp);
0711         }
0712     }
0713   if (dump_) {
0714     for (auto& ch_filteredMetaPrimitives : filteredMetaPrimitives) {
0715       for (unsigned int i = 0; i < ch_filteredMetaPrimitives.second.size(); i++) {
0716         stringstream ss;
0717         ss << iEvent.id().event() << " filtered mp " << i << ": ";
0718         printmP(ss.str(), ch_filteredMetaPrimitives.second.at(i));
0719       }
0720     }
0721   }
0722 
0723   skip_processing_ = skip_processing_ || output_slfilter_;
0724   confirmedMetaPrimitives.clear();
0725 
0726   if (debug_)
0727     for (auto& ch_filteredMetaPrimitives : filteredMetaPrimitives) {
0728       LogDebug("DTTrigPhase2Prod") << "DTp2 in event:" << iEvent.id().event() << " we found "
0729                                    << ch_filteredMetaPrimitives.second.size() << " filteredMetaPrimitives (superlayer)"
0730                                    << std::endl;
0731     }
0732   if (debug_)
0733     LogDebug("DTTrigPhase2Prod") << "filteredMetaPrimitives: starting correlations" << std::endl;
0734 
0735   /////////////////////////////////////
0736   //// CORRELATION:
0737   /////////////////////////////////////
0738 
0739   std::map<int, std::vector<metaPrimitive>> correlatedMetaPrimitives;
0740   if (algo_ == Standard) {
0741     for (auto& ch_filteredMetaPrimitives : filteredMetaPrimitives) {
0742       if (!skip_processing_)
0743         mpathassociator_->run(iEvent,
0744                               iEventSetup,
0745                               ch_filteredMetaPrimitives.second,
0746                               correlatedMetaPrimitives[ch_filteredMetaPrimitives.first]);
0747       else
0748         for (auto& mp : ch_filteredMetaPrimitives.second) {
0749           correlatedMetaPrimitives[ch_filteredMetaPrimitives.first].push_back(mp);
0750         }
0751     }
0752   } else {
0753     for (auto& ch_outmpaths : outmpaths) {
0754       for (const auto& muonpath : ch_outmpaths.second) {
0755         correlatedMetaPrimitives[ch_outmpaths.first].emplace_back(muonpath->rawId(),
0756                                                                   (double)muonpath->bxTimeValue(),
0757                                                                   muonpath->horizPos(),
0758                                                                   muonpath->tanPhi(),
0759                                                                   muonpath->phi(),
0760                                                                   muonpath->phiB(),
0761                                                                   muonpath->phi_cmssw(),
0762                                                                   muonpath->phiB_cmssw(),
0763                                                                   muonpath->chiSquare(),
0764                                                                   (int)muonpath->quality(),
0765                                                                   muonpath->primitive(0)->channelId(),
0766                                                                   muonpath->primitive(0)->tdcTimeStamp(),
0767                                                                   muonpath->primitive(0)->laterality(),
0768                                                                   muonpath->primitive(1)->channelId(),
0769                                                                   muonpath->primitive(1)->tdcTimeStamp(),
0770                                                                   muonpath->primitive(1)->laterality(),
0771                                                                   muonpath->primitive(2)->channelId(),
0772                                                                   muonpath->primitive(2)->tdcTimeStamp(),
0773                                                                   muonpath->primitive(2)->laterality(),
0774                                                                   muonpath->primitive(3)->channelId(),
0775                                                                   muonpath->primitive(3)->tdcTimeStamp(),
0776                                                                   muonpath->primitive(3)->laterality(),
0777                                                                   muonpath->primitive(4)->channelId(),
0778                                                                   muonpath->primitive(4)->tdcTimeStamp(),
0779                                                                   muonpath->primitive(4)->laterality(),
0780                                                                   muonpath->primitive(5)->channelId(),
0781                                                                   muonpath->primitive(5)->tdcTimeStamp(),
0782                                                                   muonpath->primitive(5)->laterality(),
0783                                                                   muonpath->primitive(6)->channelId(),
0784                                                                   muonpath->primitive(6)->tdcTimeStamp(),
0785                                                                   muonpath->primitive(6)->laterality(),
0786                                                                   muonpath->primitive(7)->channelId(),
0787                                                                   muonpath->primitive(7)->tdcTimeStamp(),
0788                                                                   muonpath->primitive(7)->laterality());
0789       }
0790     }
0791   }
0792 
0793   skip_processing_ = skip_processing_ || output_matcher_;
0794 
0795   if (debug_)
0796     for (auto& ch_correlatedMetaPrimitives : correlatedMetaPrimitives) {
0797       LogDebug("DTTrigPhase2Prod") << "DTp2 in event:" << iEvent.id().event() << " we found "
0798                                    << ch_correlatedMetaPrimitives.second.size() << " correlatedMetPrimitives (chamber)";
0799     }
0800   if (dump_) {
0801     for (auto& ch_correlatedMetaPrimitives : correlatedMetaPrimitives) {
0802       LogDebug("DTTrigPhase2Prod") << "DTp2 in event:" << iEvent.id().event() << " we found "
0803                                    << ch_correlatedMetaPrimitives.second.size() << " correlatedMetPrimitives (chamber)";
0804     }
0805     for (auto& ch_correlatedMetaPrimitives : correlatedMetaPrimitives) {
0806       for (unsigned int i = 0; i < ch_correlatedMetaPrimitives.second.size(); i++) {
0807         stringstream ss;
0808         ss << iEvent.id().event() << " correlated mp " << i << ": ";
0809         printmPC(ss.str(), ch_correlatedMetaPrimitives.second.at(i));
0810       }
0811     }
0812   }
0813 
0814   // Correlated Filtering
0815   std::map<int, std::vector<metaPrimitive>> filtCorrelatedMetaPrimitives;
0816   if (algo_ == Standard) {
0817     for (auto& ch_filteredMetaPrimitives : filteredMetaPrimitives) {
0818       if (!skip_processing_)
0819         mpathcorfilter_->run(iEvent,
0820                              iEventSetup,
0821                              ch_filteredMetaPrimitives.second,
0822                              correlatedMetaPrimitives[ch_filteredMetaPrimitives.first],
0823                              filtCorrelatedMetaPrimitives[ch_filteredMetaPrimitives.first]);
0824       else {
0825         for (auto& mp : ch_filteredMetaPrimitives.second) {
0826           filtCorrelatedMetaPrimitives[ch_filteredMetaPrimitives.first].push_back(mp);
0827         }
0828         if (output_matcher_)
0829           for (auto& mp : correlatedMetaPrimitives[ch_filteredMetaPrimitives.first]) {
0830             filtCorrelatedMetaPrimitives[ch_filteredMetaPrimitives.first].push_back(mp);
0831           }
0832       }
0833     }
0834   }
0835 
0836   correlatedMetaPrimitives.clear();
0837   filteredMetaPrimitives.clear();
0838 
0839   // Coincidence (co) filter
0840 
0841   std::vector<metaPrimitive> allMetaPrimitives;
0842   for (auto& ch_filtcorrelatedMetaPrimitives : filtCorrelatedMetaPrimitives) {
0843     for (const auto& metaPrimitiveIt : ch_filtcorrelatedMetaPrimitives.second) {
0844       allMetaPrimitives.push_back(metaPrimitiveIt);
0845     }
0846   }
0847 
0848   std::map<int, std::vector<metaPrimitive>> coMetaPrimitives;
0849   if (algo_ == Standard) {
0850     for (auto& ch_filtcorrelatedMetaPrimitives : filtCorrelatedMetaPrimitives) {
0851       if (!skip_processing_)
0852         mpathcoifilter_->run(iEvent,
0853                              iEventSetup,
0854                              allMetaPrimitives,
0855                              filtCorrelatedMetaPrimitives[ch_filtcorrelatedMetaPrimitives.first],
0856                              coMetaPrimitives[ch_filtcorrelatedMetaPrimitives.first]);
0857       else {
0858         for (auto& mp : ch_filtcorrelatedMetaPrimitives.second) {
0859           coMetaPrimitives[ch_filtcorrelatedMetaPrimitives.first].push_back(mp);
0860         }
0861       }
0862     }
0863   }
0864 
0865   allMetaPrimitives.clear();
0866 
0867   // Theta (th) matching filter
0868   std::map<int, std::vector<metaPrimitive>> thMatchedMetaPrimitives;
0869   if (algo_ == Standard) {
0870     for (auto& ch_filtcoMetaPrimitives : coMetaPrimitives) {
0871       if (!skip_processing_)
0872         mpaththetamatching_->run(iEvent,
0873                                  iEventSetup,
0874                                  coMetaPrimitives[ch_filtcoMetaPrimitives.first],
0875                                  thMatchedMetaPrimitives[ch_filtcoMetaPrimitives.first]);
0876       else {
0877         for (auto& mp : ch_filtcoMetaPrimitives.second) {
0878           thMatchedMetaPrimitives[ch_filtcoMetaPrimitives.first].push_back(mp);
0879         }
0880       }
0881     }
0882   }
0883 
0884   /////////////
0885 
0886   double shift_back = 0;
0887   if (scenario_ == MC)  //scope for MC
0888     shift_back = 400;
0889   else if (scenario_ == DATA)  //scope for data
0890     shift_back = 0;
0891   else if (scenario_ == SLICE_TEST)  //scope for slice test
0892     shift_back = 400;
0893 
0894   // RPC integration
0895   if (useRPC_) {
0896     rpc_integrator_->initialise(iEventSetup, shift_back);
0897     rpc_integrator_->prepareMetaPrimitives(rpcRecHits);
0898     for (auto& ch_thMatchedMetaPrimitives : thMatchedMetaPrimitives) {
0899       rpc_integrator_->matchWithDTAndUseRPCTime(ch_thMatchedMetaPrimitives.second);  // Probably this is a FIXME
0900     }
0901     rpc_integrator_->makeRPCOnlySegments();
0902     rpc_integrator_->storeRPCSingleHits();
0903     rpc_integrator_->removeRPCHitsUsed();
0904   }
0905 
0906   /// STORING RESULTs
0907   vector<L1Phase2MuDTPhDigi> outP2Ph;
0908   vector<L1Phase2MuDTExtPhDigi> outExtP2Ph;
0909   vector<L1Phase2MuDTThDigi> outP2Th;
0910   vector<L1Phase2MuDTExtThDigi> outExtP2Th;
0911 
0912   // Assigning index value
0913   if (!skip_processing_)
0914     for (auto& ch_thMatchedMetaPrimitives : thMatchedMetaPrimitives) {
0915       assignIndex(ch_thMatchedMetaPrimitives.second);
0916     }
0917 
0918   for (auto& ch_thMatchedMetaPrimitives : thMatchedMetaPrimitives) {
0919     for (const auto& metaPrimitiveIt : ch_thMatchedMetaPrimitives.second) {
0920       DTChamberId chId(metaPrimitiveIt.rawId);
0921       DTSuperLayerId slId(metaPrimitiveIt.rawId);
0922       if (debug_)
0923         LogDebug("DTTrigPhase2Prod") << "looping in final vector: SuperLayerId" << chId << " x=" << metaPrimitiveIt.x
0924                                      << " quality=" << metaPrimitiveIt.quality
0925                                      << " BX=" << round(metaPrimitiveIt.t0 / 25.) << " index=" << metaPrimitiveIt.index;
0926 
0927       int sectorTP = chId.sector();
0928       //sectors 13 and 14 exist only for the outermost stations for sectors 4 and 10 respectively
0929       //due to the larger MB4 that are divided into two.
0930       if (sectorTP == 13)
0931         sectorTP = 4;
0932       if (sectorTP == 14)
0933         sectorTP = 10;
0934       sectorTP = sectorTP - 1;
0935       int sl = 0;
0936       if (metaPrimitiveIt.quality < LOWLOWQ || metaPrimitiveIt.quality == CHIGHQ) {
0937         if (inner(metaPrimitiveIt))
0938           sl = 1;
0939         else
0940           sl = 3;
0941       }
0942 
0943       float tp_t0 =
0944           (metaPrimitiveIt.t0 - shift_back * LHC_CLK_FREQ) * ((float)TIME_TO_TDC_COUNTS / (float)LHC_CLK_FREQ);
0945 
0946       if (debug_)
0947         LogDebug("DTTrigPhase2Prod") << "pushing back phase-2 dataformat carlo-federica dataformat";
0948 
0949       if (slId.superLayer() != 2) {
0950         if (df_extended_ == 1 || df_extended_ == 2) {
0951           int pathWireId[8] = {metaPrimitiveIt.wi1,
0952                                metaPrimitiveIt.wi2,
0953                                metaPrimitiveIt.wi3,
0954                                metaPrimitiveIt.wi4,
0955                                metaPrimitiveIt.wi5,
0956                                metaPrimitiveIt.wi6,
0957                                metaPrimitiveIt.wi7,
0958                                metaPrimitiveIt.wi8};
0959 
0960           int pathTDC[8] = {max((int)round(metaPrimitiveIt.tdc1 - shift_back * LHC_CLK_FREQ), -1),
0961                             max((int)round(metaPrimitiveIt.tdc2 - shift_back * LHC_CLK_FREQ), -1),
0962                             max((int)round(metaPrimitiveIt.tdc3 - shift_back * LHC_CLK_FREQ), -1),
0963                             max((int)round(metaPrimitiveIt.tdc4 - shift_back * LHC_CLK_FREQ), -1),
0964                             max((int)round(metaPrimitiveIt.tdc5 - shift_back * LHC_CLK_FREQ), -1),
0965                             max((int)round(metaPrimitiveIt.tdc6 - shift_back * LHC_CLK_FREQ), -1),
0966                             max((int)round(metaPrimitiveIt.tdc7 - shift_back * LHC_CLK_FREQ), -1),
0967                             max((int)round(metaPrimitiveIt.tdc8 - shift_back * LHC_CLK_FREQ), -1)};
0968 
0969           int pathLat[8] = {metaPrimitiveIt.lat1,
0970                             metaPrimitiveIt.lat2,
0971                             metaPrimitiveIt.lat3,
0972                             metaPrimitiveIt.lat4,
0973                             metaPrimitiveIt.lat5,
0974                             metaPrimitiveIt.lat6,
0975                             metaPrimitiveIt.lat7,
0976                             metaPrimitiveIt.lat8};
0977 
0978           // phiTP (extended DF)
0979           outExtP2Ph.emplace_back(
0980               L1Phase2MuDTExtPhDigi((int)round(metaPrimitiveIt.t0 / (float)LHC_CLK_FREQ) - shift_back,
0981                                     chId.wheel(),                                           // uwh   (m_wheel)
0982                                     sectorTP,                                               // usc   (m_sector)
0983                                     chId.station(),                                         // ust   (m_station)
0984                                     sl,                                                     // ust   (m_station)
0985                                     (int)round(metaPrimitiveIt.phi * PHIRES_CONV),          // uphi  (m_phiAngle)
0986                                     (int)round(metaPrimitiveIt.phiB * PHIBRES_CONV),        // uphib (m_phiBending)
0987                                     metaPrimitiveIt.quality,                                // uqua  (m_qualityCode)
0988                                     metaPrimitiveIt.index,                                  // uind  (m_segmentIndex)
0989                                     tp_t0,                                                  // ut0   (m_t0Segment)
0990                                     (int)round(metaPrimitiveIt.chi2 * CHI2RES_CONV),        // uchi2 (m_chi2Segment)
0991                                     (int)round(metaPrimitiveIt.x * 1000),                   // ux    (m_xLocal)
0992                                     (int)round(metaPrimitiveIt.tanPhi * 1000),              // utan  (m_tanPsi)
0993                                     (int)round(metaPrimitiveIt.phi_cmssw * PHIRES_CONV),    // uphi  (m_phiAngleCMSSW)
0994                                     (int)round(metaPrimitiveIt.phiB_cmssw * PHIBRES_CONV),  // uphib (m_phiBendingCMSSW)
0995                                     metaPrimitiveIt.rpcFlag,                                // urpc  (m_rpcFlag)
0996                                     pathWireId,
0997                                     pathTDC,
0998                                     pathLat));
0999         }
1000         if (df_extended_ == 0 || df_extended_ == 2) {
1001           // phiTP (standard DF)
1002           outP2Ph.push_back(L1Phase2MuDTPhDigi(
1003               (int)round(metaPrimitiveIt.t0 / (float)LHC_CLK_FREQ) - shift_back,
1004               chId.wheel(),                                     // uwh (m_wheel)
1005               sectorTP,                                         // usc (m_sector)
1006               chId.station(),                                   // ust (m_station)
1007               sl,                                               // ust (m_station)
1008               (int)round(metaPrimitiveIt.phi * PHIRES_CONV),    // uphi (_phiAngle)
1009               (int)round(metaPrimitiveIt.phiB * PHIBRES_CONV),  // uphib (m_phiBending)
1010               metaPrimitiveIt.quality,                          // uqua (m_qualityCode)
1011               metaPrimitiveIt.index,                            // uind (m_segmentIndex)
1012               tp_t0,                                            // ut0 (m_t0Segment)
1013               (int)round(metaPrimitiveIt.chi2 * CHI2RES_CONV),  // uchi2 (m_chi2Segment)
1014               metaPrimitiveIt.rpcFlag                           // urpc (m_rpcFlag)
1015               ));
1016         }
1017       } else {
1018         if (df_extended_ == 1 || df_extended_ == 2) {
1019           int pathWireId[4] = {metaPrimitiveIt.wi1, metaPrimitiveIt.wi2, metaPrimitiveIt.wi3, metaPrimitiveIt.wi4};
1020 
1021           int pathTDC[4] = {max((int)round(metaPrimitiveIt.tdc1 - shift_back * LHC_CLK_FREQ), -1),
1022                             max((int)round(metaPrimitiveIt.tdc2 - shift_back * LHC_CLK_FREQ), -1),
1023                             max((int)round(metaPrimitiveIt.tdc3 - shift_back * LHC_CLK_FREQ), -1),
1024                             max((int)round(metaPrimitiveIt.tdc4 - shift_back * LHC_CLK_FREQ), -1)};
1025 
1026           int pathLat[4] = {metaPrimitiveIt.lat1, metaPrimitiveIt.lat2, metaPrimitiveIt.lat3, metaPrimitiveIt.lat4};
1027 
1028           // thTP (extended DF)
1029           outExtP2Th.emplace_back(
1030               L1Phase2MuDTExtThDigi((int)round(metaPrimitiveIt.t0 / (float)LHC_CLK_FREQ) - shift_back,
1031                                     chId.wheel(),                                        // uwh   (m_wheel)
1032                                     sectorTP,                                            // usc   (m_sector)
1033                                     chId.station(),                                      // ust   (m_station)
1034                                     (int)round(metaPrimitiveIt.phi * ZRES_CONV),         // uz    (m_zGlobal)
1035                                     (int)round(metaPrimitiveIt.phiB * KRES_CONV),        // uk    (m_kSlope)
1036                                     metaPrimitiveIt.quality,                             // uqua  (m_qualityCode)
1037                                     metaPrimitiveIt.index,                               // uind  (m_segmentIndex)
1038                                     tp_t0,                                               // ut0   (m_t0Segment)
1039                                     (int)round(metaPrimitiveIt.chi2 * CHI2RES_CONV),     // uchi2 (m_chi2Segment)
1040                                     (int)round(metaPrimitiveIt.x * 1000),                // ux    (m_yLocal)
1041                                     (int)round(metaPrimitiveIt.phi_cmssw * ZRES_CONV),   // uphi  (m_zCMSSW)
1042                                     (int)round(metaPrimitiveIt.phiB_cmssw * KRES_CONV),  // uphib (m_kCMSSW)
1043                                     metaPrimitiveIt.rpcFlag,                             // urpc  (m_rpcFlag)
1044                                     pathWireId,
1045                                     pathTDC,
1046                                     pathLat));
1047         }
1048         if (df_extended_ == 0 || df_extended_ == 2) {
1049           // thTP (standard DF)
1050           outP2Th.push_back(L1Phase2MuDTThDigi(
1051               (int)round(metaPrimitiveIt.t0 / (float)LHC_CLK_FREQ) - shift_back,
1052               chId.wheel(),                                     // uwh (m_wheel)
1053               sectorTP,                                         // usc (m_sector)
1054               chId.station(),                                   // ust (m_station)
1055               (int)round(metaPrimitiveIt.phi * ZRES_CONV),      // uz (m_zGlobal)
1056               (int)round(metaPrimitiveIt.phiB * KRES_CONV),     // uk (m_kSlope)
1057               metaPrimitiveIt.quality,                          // uqua (m_qualityCode)
1058               metaPrimitiveIt.index,                            // uind (m_segmentIndex)
1059               tp_t0,                                            // ut0 (m_t0Segment)
1060               (int)round(metaPrimitiveIt.chi2 * CHI2RES_CONV),  // uchi2 (m_chi2Segment)
1061               metaPrimitiveIt.rpcFlag                           // urpc (m_rpcFlag)
1062               ));
1063         }
1064       }
1065     }
1066   }
1067 
1068   // Storing RPC hits that were not used elsewhere
1069   if (useRPC_) {
1070     for (auto rpc_dt_digi = rpc_integrator_->rpcRecHits_translated_.begin();
1071          rpc_dt_digi != rpc_integrator_->rpcRecHits_translated_.end();
1072          rpc_dt_digi++) {
1073       outP2Ph.push_back(*rpc_dt_digi);
1074     }
1075   }
1076 
1077   // Storing Phi results
1078   if (df_extended_ == 1 || df_extended_ == 2) {
1079     std::unique_ptr<L1Phase2MuDTExtPhContainer> resultExtP2Ph(new L1Phase2MuDTExtPhContainer);
1080     resultExtP2Ph->setContainer(outExtP2Ph);
1081     iEvent.put(std::move(resultExtP2Ph));
1082   }
1083   if (df_extended_ == 0 || df_extended_ == 2) {
1084     std::unique_ptr<L1Phase2MuDTPhContainer> resultP2Ph(new L1Phase2MuDTPhContainer);
1085     resultP2Ph->setContainer(outP2Ph);
1086     iEvent.put(std::move(resultP2Ph));
1087   }
1088   outExtP2Ph.clear();
1089   outExtP2Ph.erase(outExtP2Ph.begin(), outExtP2Ph.end());
1090   outP2Ph.clear();
1091   outP2Ph.erase(outP2Ph.begin(), outP2Ph.end());
1092 
1093   // Storing Theta results
1094   if (df_extended_ == 1 || df_extended_ == 2) {
1095     std::unique_ptr<L1Phase2MuDTExtThContainer> resultExtP2Th(new L1Phase2MuDTExtThContainer);
1096     resultExtP2Th->setContainer(outExtP2Th);
1097     iEvent.put(std::move(resultExtP2Th));
1098   }
1099   if (df_extended_ == 0 || df_extended_ == 2) {
1100     std::unique_ptr<L1Phase2MuDTThContainer> resultP2Th(new L1Phase2MuDTThContainer);
1101     resultP2Th->setContainer(outP2Th);
1102     iEvent.put(std::move(resultP2Th));
1103   }
1104   outExtP2Th.clear();
1105   outExtP2Th.erase(outExtP2Th.begin(), outExtP2Th.end());
1106   outP2Th.clear();
1107   outP2Th.erase(outP2Th.begin(), outP2Th.end());
1108 }
1109 
1110 void DTTrigPhase2Prod::endRun(edm::Run const& iRun, const edm::EventSetup& iEventSetup) {
1111   grouping_obj_->finish();
1112   mpathanalyzer_->finish();
1113   mpathqualityenhancer_->finish();
1114   mpathqualityenhancerbayes_->finish();
1115   mpathredundantfilter_->finish();
1116   mpathhitsfilter_->finish();
1117   mpathassociator_->finish();
1118   rpc_integrator_->finish();
1119 };
1120 
1121 bool DTTrigPhase2Prod::outer(const metaPrimitive& mp) const {
1122   int counter = (mp.wi5 != -1) + (mp.wi6 != -1) + (mp.wi7 != -1) + (mp.wi8 != -1);
1123   return (counter > 2);
1124 }
1125 
1126 bool DTTrigPhase2Prod::inner(const metaPrimitive& mp) const {
1127   int counter = (mp.wi1 != -1) + (mp.wi2 != -1) + (mp.wi3 != -1) + (mp.wi4 != -1);
1128   return (counter > 2);
1129 }
1130 
1131 bool DTTrigPhase2Prod::hasPosRF(int wh, int sec) const { return wh > 0 || (wh == 0 && sec % 4 > 1); }
1132 
1133 void DTTrigPhase2Prod::printmP(const string& ss, const metaPrimitive& mP) const {
1134   DTSuperLayerId slId(mP.rawId);
1135   LogInfo("DTTrigPhase2Prod") << ss << (int)slId << "\t " << setw(2) << left << mP.wi1 << " " << setw(2) << left
1136                               << mP.wi2 << " " << setw(2) << left << mP.wi3 << " " << setw(2) << left << mP.wi4 << " "
1137                               << setw(5) << left << mP.tdc1 << " " << setw(5) << left << mP.tdc2 << " " << setw(5)
1138                               << left << mP.tdc3 << " " << setw(5) << left << mP.tdc4 << " " << setw(10) << right
1139                               << mP.x << " " << setw(9) << left << mP.tanPhi << " " << setw(5) << left << mP.t0 << " "
1140                               << setw(13) << left << mP.chi2 << " r:" << rango(mP);
1141 }
1142 
1143 void DTTrigPhase2Prod::printmP(const metaPrimitive& mP) const {
1144   DTSuperLayerId slId(mP.rawId);
1145   LogInfo("DTTrigPhase2Prod") << (int)slId << "\t " << setw(2) << left << mP.wi1 << " " << setw(2) << left << mP.wi2
1146                               << " " << setw(2) << left << mP.wi3 << " " << setw(2) << left << mP.wi4 << " " << setw(5)
1147                               << left << mP.tdc1 << " " << setw(5) << left << mP.tdc2 << " " << setw(5) << left
1148                               << mP.tdc3 << " " << setw(5) << left << mP.tdc4 << " " << setw(10) << right << mP.x << " "
1149                               << setw(9) << left << mP.tanPhi << " " << setw(5) << left << mP.t0 << " " << setw(13)
1150                               << left << mP.chi2 << " r:" << rango(mP) << std::endl;
1151 }
1152 
1153 void DTTrigPhase2Prod::printmPC(const string& ss, const metaPrimitive& mP) const {
1154   DTChamberId ChId(mP.rawId);
1155   LogInfo("DTTrigPhase2Prod") << ss << (int)ChId << "\t  " << setw(2) << left << mP.wi1 << " " << setw(2) << left
1156                               << mP.wi2 << " " << setw(2) << left << mP.wi3 << " " << setw(2) << left << mP.wi4 << " "
1157                               << setw(2) << left << mP.wi5 << " " << setw(2) << left << mP.wi6 << " " << setw(2) << left
1158                               << mP.wi7 << " " << setw(2) << left << mP.wi8 << " " << setw(5) << left << mP.tdc1 << " "
1159                               << setw(5) << left << mP.tdc2 << " " << setw(5) << left << mP.tdc3 << " " << setw(5)
1160                               << left << mP.tdc4 << " " << setw(5) << left << mP.tdc5 << " " << setw(5) << left
1161                               << mP.tdc6 << " " << setw(5) << left << mP.tdc7 << " " << setw(5) << left << mP.tdc8
1162                               << " " << setw(2) << left << mP.lat1 << " " << setw(2) << left << mP.lat2 << " "
1163                               << setw(2) << left << mP.lat3 << " " << setw(2) << left << mP.lat4 << " " << setw(2)
1164                               << left << mP.lat5 << " " << setw(2) << left << mP.lat6 << " " << setw(2) << left
1165                               << mP.lat7 << " " << setw(2) << left << mP.lat8 << " " << setw(10) << right << mP.x << " "
1166                               << setw(9) << left << mP.tanPhi << " " << setw(5) << left << mP.t0 << " " << setw(13)
1167                               << left << mP.chi2 << " r:" << rango(mP);
1168 }
1169 
1170 void DTTrigPhase2Prod::printmPC(const metaPrimitive& mP) const {
1171   DTChamberId ChId(mP.rawId);
1172   LogInfo("DTTrigPhase2Prod") << (int)ChId << "\t  " << setw(2) << left << mP.wi1 << " " << setw(2) << left << mP.wi2
1173                               << " " << setw(2) << left << mP.wi3 << " " << setw(2) << left << mP.wi4 << " " << setw(2)
1174                               << left << mP.wi5 << " " << setw(2) << left << mP.wi6 << " " << setw(2) << left << mP.wi7
1175                               << " " << setw(2) << left << mP.wi8 << " " << setw(5) << left << mP.tdc1 << " " << setw(5)
1176                               << left << mP.tdc2 << " " << setw(5) << left << mP.tdc3 << " " << setw(5) << left
1177                               << mP.tdc4 << " " << setw(5) << left << mP.tdc5 << " " << setw(5) << left << mP.tdc6
1178                               << " " << setw(5) << left << mP.tdc7 << " " << setw(5) << left << mP.tdc8 << " "
1179                               << setw(2) << left << mP.lat1 << " " << setw(2) << left << mP.lat2 << " " << setw(2)
1180                               << left << mP.lat3 << " " << setw(2) << left << mP.lat4 << " " << setw(2) << left
1181                               << mP.lat5 << " " << setw(2) << left << mP.lat6 << " " << setw(2) << left << mP.lat7
1182                               << " " << setw(2) << left << mP.lat8 << " " << setw(10) << right << mP.x << " " << setw(9)
1183                               << left << mP.tanPhi << " " << setw(5) << left << mP.t0 << " " << setw(13) << left
1184                               << mP.chi2 << " r:" << rango(mP) << std::endl;
1185 }
1186 
1187 int DTTrigPhase2Prod::rango(const metaPrimitive& mp) const {
1188   if (mp.quality == 1 or mp.quality == 2)
1189     return 3;
1190   if (mp.quality == 3 or mp.quality == 4)
1191     return 4;
1192   return mp.quality;
1193 }
1194 
1195 void DTTrigPhase2Prod::assignIndex(std::vector<metaPrimitive>& inMPaths) {
1196   std::map<int, std::vector<metaPrimitive>> primsPerBX;
1197   for (const auto& metaPrimitive : inMPaths) {
1198     int BX = round(metaPrimitive.t0 / 25.);
1199     primsPerBX[BX].push_back(metaPrimitive);
1200   }
1201   inMPaths.clear();
1202   for (auto& prims : primsPerBX) {
1203     assignIndexPerBX(prims.second);
1204     for (const auto& primitive : prims.second)
1205       if (primitive.index <= max_index_)
1206         inMPaths.push_back(primitive);
1207   }
1208 }
1209 
1210 void DTTrigPhase2Prod::assignIndexPerBX(std::vector<metaPrimitive>& inMPaths) {
1211   // First we asociate a new index to the metaprimitive depending on quality or phiB;
1212   uint32_t rawId = -1;
1213   int numP = -1;
1214   for (auto& metaPrimitiveIt : inMPaths) {
1215     numP++;
1216     rawId = metaPrimitiveIt.rawId;
1217     int iOrder = assignQualityOrder(metaPrimitiveIt);
1218     int inf = 0;
1219     int numP2 = -1;
1220     for (auto& metaPrimitiveItN : inMPaths) {
1221       int nOrder = assignQualityOrder(metaPrimitiveItN);
1222       numP2++;
1223       if (rawId != metaPrimitiveItN.rawId)
1224         continue;
1225       if (numP2 == numP) {
1226         metaPrimitiveIt.index = inf;
1227         break;
1228       } else if (iOrder < nOrder) {
1229         inf++;
1230       } else if (iOrder > nOrder) {
1231         metaPrimitiveItN.index++;
1232       } else if (iOrder == nOrder) {
1233         if (std::abs(metaPrimitiveIt.phiB) >= std::abs(metaPrimitiveItN.phiB)) {
1234           inf++;
1235         } else if (std::abs(metaPrimitiveIt.phiB) < std::abs(metaPrimitiveItN.phiB)) {
1236           metaPrimitiveItN.index++;
1237         }
1238       }
1239     }  // ending second for
1240   }  // ending first for
1241 }
1242 
1243 int DTTrigPhase2Prod::assignQualityOrder(const metaPrimitive& mP) const {
1244   if (mP.quality > 8 || mP.quality < 1)
1245     return -1;
1246 
1247   return qmap_.find(mP.quality)->second;
1248 }
1249 
1250 std::vector<DTDigiCollection*> DTTrigPhase2Prod::distribDigis(std::queue<std::pair<DTLayerId, DTDigi>>& inQ) {
1251   std::vector<std::queue<std::pair<DTLayerId, DTDigi>>*> tmpVector;
1252   tmpVector.clear();
1253   std::vector<DTDigiCollection*> collVector;
1254   collVector.clear();
1255   while (!inQ.empty()) {
1256     processDigi(inQ, tmpVector);
1257   }
1258   for (auto& sQ : tmpVector) {
1259     DTDigiCollection tmpColl;
1260     while (!sQ->empty()) {
1261       tmpColl.insertDigi((sQ->front().first), (sQ->front().second));
1262       sQ->pop();
1263     }
1264     collVector.push_back(&tmpColl);
1265   }
1266   return collVector;
1267 }
1268 
1269 void DTTrigPhase2Prod::processDigi(std::queue<std::pair<DTLayerId, DTDigi>>& inQ,
1270                                    std::vector<std::queue<std::pair<DTLayerId, DTDigi>>*>& vec) {
1271   bool classified = false;
1272   if (!vec.empty()) {
1273     for (auto& sC : vec) {  // Conditions for entering a super cell.
1274       if ((sC->front().second.time() + superCelltimewidth_) > inQ.front().second.time()) {
1275         // Time requirement
1276         if (TMath::Abs(sC->front().second.wire() - inQ.front().second.wire()) <= superCellhalfspacewidth_) {
1277           // Spatial requirement
1278           sC->push(std::move(inQ.front()));
1279           classified = true;
1280         }
1281       }
1282     }
1283   }
1284   if (classified) {
1285     inQ.pop();
1286     return;
1287   }
1288 
1289   std::queue<std::pair<DTLayerId, DTDigi>> newQueue;
1290 
1291   std::pair<DTLayerId, DTDigi> tmpPair;
1292   tmpPair = std::move(inQ.front());
1293   newQueue.push(tmpPair);
1294   inQ.pop();
1295 
1296   vec.push_back(&newQueue);
1297 }
1298 
1299 void DTTrigPhase2Prod::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
1300   // dtTriggerPhase2PrimitiveDigis
1301   edm::ParameterSetDescription desc;
1302   desc.add<edm::InputTag>("digiTag", edm::InputTag("CalibratedDigis"));
1303   desc.add<int>("timeTolerance", 999999);
1304   desc.add<double>("tanPhiTh", 1.0);
1305   desc.add<double>("tanPhiThw2max", 1.3);
1306   desc.add<double>("tanPhiThw2min", 0.5);
1307   desc.add<double>("tanPhiThw1max", 0.9);
1308   desc.add<double>("tanPhiThw1min", 0.2);
1309   desc.add<double>("tanPhiThw0", 0.5);
1310   desc.add<double>("chi2Th", 0.01);
1311   desc.add<double>("chi2corTh", 0.1);
1312   desc.add<bool>("useBX_correlation", false);
1313   desc.add<double>("dT0_correlate_TP", 25.0);
1314   desc.add<int>("dBX_correlate_TP", 0);
1315   desc.add<double>("dTanPsi_correlate_TP", 99999.0);
1316   desc.add<bool>("clean_chi2_correlation", true);
1317   desc.add<bool>("allow_confirmation", true);
1318   desc.add<double>("minx_match_2digis", 1.0);
1319   desc.add<int>("scenario", 0);
1320   desc.add<int>("df_extended", 0);
1321   desc.add<int>("co_option", 0);
1322   desc.add<int>("co_quality", 0);
1323   desc.add<int>("co_wh2option", 0);
1324   desc.add<int>("th_option", 0);
1325   desc.add<int>("th_quality", 0);
1326   desc.add<int>("max_primitives", 999);
1327   desc.add<bool>("output_mixer", false);
1328   desc.add<bool>("output_latpredictor", false);
1329   desc.add<bool>("output_slfitter", false);
1330   desc.add<bool>("output_slfilter", false);
1331   desc.add<bool>("output_confirmed", false);
1332   desc.add<bool>("output_matcher", false);
1333   desc.add<edm::FileInPath>("ttrig_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/wire_rawId_ttrig.txt"));
1334   desc.add<edm::FileInPath>("z_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/wire_rawId_z.txt"));
1335   desc.add<edm::FileInPath>("lut_sl1", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/fitterlut_sl1.dat"));
1336   desc.add<edm::FileInPath>("lut_sl2", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/fitterlut_slx.dat"));
1337   desc.add<edm::FileInPath>("lut_sl3", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/fitterlut_sl3.dat"));
1338   desc.add<edm::FileInPath>("lut_2sl", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/fitterlut_2sl.dat"));
1339   desc.add<edm::FileInPath>("shift_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/wire_rawId_x.txt"));
1340   desc.add<edm::FileInPath>("maxdrift_filename",
1341                             edm::FileInPath("L1Trigger/DTTriggerPhase2/data/drift_time_per_chamber.txt"));
1342   desc.add<edm::FileInPath>("shift_theta_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/theta_shift.txt"));
1343   desc.add<edm::FileInPath>("global_coords_filename",
1344                             edm::FileInPath("L1Trigger/DTTriggerPhase2/data/global_coord_perp_x_phi0.txt"));
1345   desc.add<edm::FileInPath>("laterality_filename",
1346                             edm::FileInPath("L1Trigger/DTTriggerPhase2/data/lat_predictions.dat"));
1347   desc.add<int>("algo", 0);
1348   desc.add<int>("minHits4Fit", 3);
1349   desc.add<bool>("splitPathPerSL", true);
1350   desc.addUntracked<bool>("debug", false);
1351   desc.addUntracked<bool>("dump", false);
1352   desc.add<edm::InputTag>("rpcRecHits", edm::InputTag("rpcRecHits"));
1353   desc.add<bool>("useRPC", false);
1354   desc.add<int>("bx_window", 1);
1355   desc.add<double>("phi_window", 50.0);
1356   desc.add<int>("max_quality_to_overwrite_t0", 9);
1357   desc.add<bool>("storeAllRPCHits", false);
1358   desc.add<bool>("activateBuffer", false);
1359   desc.add<double>("superCelltimewidth", 400);
1360   desc.add<int>("superCellspacewidth", 20);
1361   {
1362     edm::ParameterSetDescription psd0;
1363     psd0.addUntracked<bool>("debug", false);
1364     psd0.add<double>("angletan", 0.3);
1365     psd0.add<double>("anglebinwidth", 1.0);
1366     psd0.add<double>("posbinwidth", 2.1);
1367     psd0.add<double>("maxdeltaAngDeg", 10);
1368     psd0.add<double>("maxdeltaPos", 10);
1369     psd0.add<int>("UpperNumber", 6);
1370     psd0.add<int>("LowerNumber", 4);
1371     psd0.add<double>("MaxDistanceToWire", 0.03);
1372     psd0.add<int>("minNLayerHits", 6);
1373     psd0.add<int>("minSingleSLHitsMax", 3);
1374     psd0.add<int>("minSingleSLHitsMin", 3);
1375     psd0.add<bool>("allowUncorrelatedPatterns", true);
1376     psd0.add<int>("minUncorrelatedHits", 3);
1377     desc.add<edm::ParameterSetDescription>("HoughGrouping", psd0);
1378   }
1379   {
1380     edm::ParameterSetDescription psd0;
1381     psd0.add<edm::FileInPath>(
1382         "pattern_filename", edm::FileInPath("L1Trigger/DTTriggerPhase2/data/PseudoBayesPatterns_uncorrelated_v0.root"));
1383     psd0.addUntracked<bool>("debug", false);
1384     psd0.add<int>("minNLayerHits", 3);
1385     psd0.add<int>("minSingleSLHitsMax", 3);
1386     psd0.add<int>("minSingleSLHitsMin", 0);
1387     psd0.add<int>("allowedVariance", 1);
1388     psd0.add<bool>("allowDuplicates", false);
1389     psd0.add<bool>("setLateralities", true);
1390     psd0.add<bool>("allowUncorrelatedPatterns", true);
1391     psd0.add<int>("minUncorrelatedHits", 3);
1392     psd0.add<bool>("saveOnPlace", true);
1393     psd0.add<int>("maxPathsPerMatch", 256);
1394     desc.add<edm::ParameterSetDescription>("PseudoBayesPattern", psd0);
1395   }
1396   descriptions.add("dtTriggerPhase2PrimitiveDigis", desc);
1397 }
1398 
1399 DEFINE_FWK_MODULE(DTTrigPhase2Prod);