Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:21:08

0001 /*
0002  * OMTFProcessor.cpp
0003  *
0004  *  Created on: Oct 7, 2017
0005  *      Author: kbunkow
0006  */
0007 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/OMTFProcessor.h"
0008 #include "L1Trigger/L1TMuonOverlapPhase1/interface/MuonStub.h"
0009 #include "L1Trigger/L1TMuonOverlapPhase1/interface/MuonStubsInput.h"
0010 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/GhostBuster.h"
0011 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/GhostBusterPreferRefDt.h"
0012 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/GoldenPatternWithStat.h"
0013 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/IOMTFEmulationObserver.h"
0014 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/OMTFinput.h"
0015 #include "L1Trigger/L1TMuonOverlapPhase1/interface/Omtf/OMTFSorter.h"
0016 #include "L1Trigger/L1TMuonOverlapPhase1/interface/StubResult.h"
0017 
0018 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0019 
0020 #include <bitset>
0021 #include <cmath>
0022 #include <cstdlib>
0023 #include <iostream>
0024 #include <iomanip>
0025 #include <map>
0026 #include <string>
0027 #include <vector>
0028 
0029 #include <boost/property_tree/ptree.hpp>
0030 #include <boost/property_tree/xml_parser.hpp>
0031 
0032 ///////////////////////////////////////////////
0033 ///////////////////////////////////////////////
0034 template <class GoldenPatternType>
0035 OMTFProcessor<GoldenPatternType>::OMTFProcessor(OMTFConfiguration* omtfConfig,
0036                                                 const edm::ParameterSet& edmCfg,
0037                                                 edm::EventSetup const& evSetup,
0038                                                 const L1TMuonOverlapParams* omtfPatterns)
0039     : ProcessorBase<GoldenPatternType>(omtfConfig, omtfPatterns) {
0040   init(edmCfg, evSetup);
0041 };
0042 
0043 template <class GoldenPatternType>
0044 OMTFProcessor<GoldenPatternType>::OMTFProcessor(OMTFConfiguration* omtfConfig,
0045                                                 const edm::ParameterSet& edmCfg,
0046                                                 edm::EventSetup const& evSetup,
0047                                                 GoldenPatternVec<GoldenPatternType>&& gps)
0048     : ProcessorBase<GoldenPatternType>(omtfConfig, std::forward<GoldenPatternVec<GoldenPatternType> >(gps)) {
0049   init(edmCfg, evSetup);
0050 };
0051 
0052 template <class GoldenPatternType>
0053 OMTFProcessor<GoldenPatternType>::~OMTFProcessor() {
0054   if (useFloatingPointExtrapolation)
0055     saveExtrapolFactors();
0056 }
0057 
0058 template <class GoldenPatternType>
0059 void OMTFProcessor<GoldenPatternType>::init(const edm::ParameterSet& edmCfg, edm::EventSetup const& evSetup) {
0060   setSorter(new OMTFSorter<GoldenPatternType>(this->myOmtfConfig->getSorterType()));
0061   //initialize with the default sorter
0062 
0063   if (this->myOmtfConfig->getGhostBusterType() == "GhostBusterPreferRefDt" ||
0064       this->myOmtfConfig->getGhostBusterType() == "byLLH" || this->myOmtfConfig->getGhostBusterType() == "byFPLLH" ||
0065       this->myOmtfConfig->getGhostBusterType() == "byRefLayer") {
0066     setGhostBuster(new GhostBusterPreferRefDt(this->myOmtfConfig));
0067     edm::LogVerbatim("OMTFReconstruction") << "setting " << this->myOmtfConfig->getGhostBusterType() << std::endl;
0068   } else {
0069     setGhostBuster(new GhostBuster(this->myOmtfConfig));  //initialize with the default sorter
0070     edm::LogVerbatim("OMTFReconstruction") << "setting GhostBuster" << std::endl;
0071   }
0072 
0073   edm::LogVerbatim("OMTFReconstruction") << "fwVersion 0x" << hex << this->myOmtfConfig->fwVersion() << std::endl;
0074 
0075   useStubQualInExtr = this->myOmtfConfig->useStubQualInExtr();
0076   useEndcapStubsRInExtr = this->myOmtfConfig->useEndcapStubsRInExtr();
0077 
0078   if (edmCfg.exists("useFloatingPointExtrapolation"))
0079     useFloatingPointExtrapolation = edmCfg.getParameter<bool>("useFloatingPointExtrapolation");
0080 
0081   std::string extrapolFactorsFilename;
0082   if (edmCfg.exists("extrapolFactorsFilename")) {
0083     extrapolFactorsFilename = edmCfg.getParameter<edm::FileInPath>("extrapolFactorsFilename").fullPath();
0084   }
0085 
0086   if (this->myOmtfConfig->usePhiBExtrapolationMB1() || this->myOmtfConfig->usePhiBExtrapolationMB2()) {
0087     extrapolFactors.resize(2, std::vector<std::map<int, double> >(this->myOmtfConfig->nLayers()));
0088     extrapolFactorsNorm.resize(2, std::vector<std::map<int, int> >(this->myOmtfConfig->nLayers()));
0089 
0090     //when useFloatingPointExtrapolation is true the extrapolFactors are not used,
0091     //all calculations are done in the extrapolateDtPhiBFloatPoint
0092     if (!extrapolFactorsFilename.empty() && !useFloatingPointExtrapolation)
0093       loadExtrapolFactors(extrapolFactorsFilename);
0094   }
0095 
0096   edm::LogVerbatim("OMTFReconstruction") << "useFloatingPointExtrapolation " << useFloatingPointExtrapolation
0097                                          << std::endl;
0098   edm::LogVerbatim("OMTFReconstruction") << "extrapolFactorsFilename " << extrapolFactorsFilename << std::endl;
0099 }
0100 
0101 template <class GoldenPatternType>
0102 std::vector<l1t::RegionalMuonCand> OMTFProcessor<GoldenPatternType>::getFinalcandidates(unsigned int iProcessor,
0103                                                                                         l1t::tftype mtfType,
0104                                                                                         const AlgoMuons& algoCands) {
0105   std::vector<l1t::RegionalMuonCand> result;
0106 
0107   for (auto& myCand : algoCands) {
0108     l1t::RegionalMuonCand candidate;
0109 
0110     //the charge is only for the constrained measurement. The constrained measurement is always defined for a valid candidate
0111     if (ptAssignment) {
0112       candidate.setHwPt(myCand->getPtNNConstr());
0113       candidate.setHwSign(myCand->getChargeNNConstr() < 0 ? 1 : 0);
0114     } else {
0115       candidate.setHwPt(myCand->getPtConstr());
0116       candidate.setHwSign(myCand->getChargeConstr() < 0 ? 1 : 0);
0117     }
0118 
0119     if (mtfType == l1t::omtf_pos)
0120       candidate.setHwEta(myCand->getEtaHw());
0121     else
0122       candidate.setHwEta((-1) * myCand->getEtaHw());
0123 
0124     int phiValue = myCand->getPhi();
0125     if (phiValue >= int(this->myOmtfConfig->nPhiBins()))
0126       phiValue -= this->myOmtfConfig->nPhiBins();
0127     phiValue = this->myOmtfConfig->procPhiToGmtPhi(phiValue);
0128     candidate.setHwPhi(phiValue);
0129 
0130     candidate.setHwSignValid(1);
0131 
0132     if (myCand->getPtUnconstr() >= 0) {  //empty PtUnconstrained is -1, maybe should be corrected on the source
0133       //the upt has different hardware scale than the pt, the upt unit is 1 GeV
0134       candidate.setHwPtUnconstrained(myCand->getPtUnconstr());
0135     } else
0136       candidate.setHwPtUnconstrained(0);
0137 
0138     unsigned int quality = 12;
0139     if (this->myOmtfConfig->fwVersion() <= 6)
0140       quality = checkHitPatternValidity(myCand->getFiredLayerBits()) ? 0 | (1 << 2) | (1 << 3) : 0 | (1 << 2);  //12 : 4
0141 
0142     if (abs(myCand->getEtaHw()) == 115 &&  //115 is eta 1.25                         rrrrrrrrccccdddddd
0143         (static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000001110000000").to_ulong() ||
0144          static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000001110000000").to_ulong() ||
0145          static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000000110000000").to_ulong() ||
0146          static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000001100000000").to_ulong() ||
0147          static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000001010000000").to_ulong())) {
0148       if (this->myOmtfConfig->fwVersion() <= 6)
0149         quality = 4;
0150       else
0151         quality = 1;
0152     }
0153 
0154     if (this->myOmtfConfig->fwVersion() >= 5 && this->myOmtfConfig->fwVersion() <= 6) {
0155       if (static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000010000000011").to_ulong() ||
0156           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000100000000011").to_ulong() ||
0157           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001000000000011").to_ulong() ||
0158           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010000000000011").to_ulong() ||
0159           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000100000000000011").to_ulong() ||
0160           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000000000000011").to_ulong() ||
0161           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000000000000011").to_ulong() ||
0162           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000000000000011").to_ulong() ||
0163           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000010000001100").to_ulong() ||
0164           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000100000001100").to_ulong() ||
0165           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001000000001100").to_ulong() ||
0166           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010000000001100").to_ulong() ||
0167           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000100000000001100").to_ulong() ||
0168           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000000000001100").to_ulong() ||
0169           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000000000001100").to_ulong() ||
0170           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000000000001100").to_ulong() ||
0171           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000010000110000").to_ulong() ||
0172           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000100000110000").to_ulong() ||
0173           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001000000110000").to_ulong() ||
0174           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010000000110000").to_ulong() ||
0175           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000100000000110000").to_ulong() ||
0176           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000000000110000").to_ulong() ||
0177           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000000000110000").to_ulong() ||
0178           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000000000110000").to_ulong())
0179         quality = 1;
0180     } else if (this->myOmtfConfig->fwVersion() >= 8) {  //TODO fix the fwVersion     rrrrrrrrccccdddddd
0181       if (static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000110000000011").to_ulong() ||
0182           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000100000000011").to_ulong() ||
0183           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000010000000011").to_ulong() ||
0184           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000110000000001").to_ulong() ||
0185 
0186           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001000000001100").to_ulong() ||
0187           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000011000000001100").to_ulong() ||
0188           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010000000001100").to_ulong() ||
0189           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000011000000000100").to_ulong() ||
0190 
0191           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000011000000001").to_ulong() ||
0192           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000010000000001").to_ulong())
0193         quality = 1;
0194       else if (
0195           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000010000000101").to_ulong() ||
0196           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000010001000001").to_ulong() ||
0197           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000011000000001").to_ulong() ||
0198           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000011000000011").to_ulong() ||
0199           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000011100000001").to_ulong() ||
0200           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000100000000011").to_ulong() ||
0201           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000100001000100").to_ulong() ||
0202           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000100100000001").to_ulong() ||
0203           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000110100000001").to_ulong() ||
0204           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000111000000000").to_ulong() ||
0205           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000111000000001").to_ulong() ||
0206           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000111000000011").to_ulong() ||
0207           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001000001000100").to_ulong() ||
0208           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001010000000001").to_ulong() ||
0209           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001010000000011").to_ulong() ||
0210           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001010000000100").to_ulong() ||
0211           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001100000000001").to_ulong() ||
0212           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001100000000100").to_ulong() ||
0213           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001100000000111").to_ulong() ||
0214           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001100001000000").to_ulong() ||
0215           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001110000000100").to_ulong() ||
0216           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001110000000101").to_ulong() ||
0217           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010000000000101").to_ulong() ||
0218           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010010000000001").to_ulong() ||
0219           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010010000000100").to_ulong() ||
0220           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010010000000101").to_ulong() ||
0221           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010100000000001").to_ulong() ||
0222           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010100000000101").to_ulong() ||
0223           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000011110000000100").to_ulong() ||
0224           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000011110000000101").to_ulong() ||
0225           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000101000000010101").to_ulong() ||
0226           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000010000000001").to_ulong() ||
0227           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000011000000000").to_ulong() ||
0228           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000011000000001").to_ulong() ||
0229           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000100000000001").to_ulong() ||
0230           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000110000000000").to_ulong() ||
0231           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001001000000000100").to_ulong() ||
0232           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001001100000000100").to_ulong() ||
0233           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001010000000000100").to_ulong() ||
0234           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000000010000001").to_ulong() ||
0235           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000000011000100").to_ulong() ||
0236           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000010000000001").to_ulong() ||
0237           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000100000000001").to_ulong() ||
0238           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000011000000000").to_ulong() ||
0239           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000110000000001").to_ulong() ||
0240           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000010000000011").to_ulong() ||
0241           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000110000000011").to_ulong() ||
0242           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000011000000001100").to_ulong() ||
0243           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000011000000000100").to_ulong() ||
0244           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000010010000001").to_ulong() ||
0245           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010000000001100").to_ulong() ||
0246           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001001000001000100").to_ulong() ||
0247           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001100000000101").to_ulong() ||
0248           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000100000000101").to_ulong() ||
0249           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001100000000011").to_ulong() ||
0250           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001110000000111").to_ulong() ||
0251           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000110001000001").to_ulong() ||
0252           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001110000000011").to_ulong() ||
0253           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000000001000100").to_ulong() ||
0254           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000110001000001").to_ulong() ||
0255           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001000000000101").to_ulong() ||
0256           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001010000001000000").to_ulong() ||
0257           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001100000001000000").to_ulong() ||
0258           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("100000010000000001").to_ulong() ||
0259           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000010010000000").to_ulong() ||
0260           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000010100000001100").to_ulong() ||
0261           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000110000000011").to_ulong() ||
0262           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001000000001100").to_ulong() ||
0263           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000000000000111101").to_ulong() ||
0264           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000001100000110001").to_ulong() ||
0265           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000100000000010100").to_ulong() ||
0266           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000100000000011").to_ulong() ||
0267           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("001000110000000001").to_ulong() ||
0268           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("010000100010000001").to_ulong() ||
0269           static_cast<unsigned int>(myCand->getFiredLayerBits()) == std::bitset<18>("000100000000110000").to_ulong())
0270         quality = 8;
0271     }  //  if (abs(myCand->getEta()) == 121) quality = 4;
0272     if (abs(myCand->getEtaHw()) >= 121)
0273       quality = 0;  // changed from 4 on request from HI
0274 
0275     candidate.setHwQual(quality);
0276 
0277     std::map<int, int> trackAddr;
0278     trackAddr[0] = myCand->getFiredLayerBits();
0279     //TODO in the hardware, the uPt is sent to the uGMT at the trackAddr = (uPt << 18) + trackAddr;
0280     //check if it matters if it needs to be here as well
0281     trackAddr[1] = myCand->getRefLayer();
0282     trackAddr[2] = myCand->getDisc();
0283     trackAddr[3] = myCand->getGpResultUnconstr().getPdfSumUnconstr();
0284     if (candidate.hwPt() > 0) {
0285       candidate.setTrackAddress(trackAddr);
0286       candidate.setTFIdentifiers(iProcessor, mtfType);
0287       result.push_back(candidate);
0288     }
0289   }
0290   return result;
0291 }
0292 ///////////////////////////////////////////////////////
0293 ///////////////////////////////////////////////////////
0294 
0295 ///////////////////////////////////////////////////////
0296 ///////////////////////////////////////////////////////
0297 template <class GoldenPatternType>
0298 bool OMTFProcessor<GoldenPatternType>::checkHitPatternValidity(unsigned int hits) {
0299   ///FIXME: read the list from configuration so this can be controlled at runtime.
0300   std::vector<unsigned int> badPatterns = {
0301       99840, 34304, 3075, 36928, 12300, 98816, 98944, 33408, 66688, 66176, 7171, 20528, 33856, 35840, 4156, 34880};
0302 
0303   /*
0304 99840 01100001 1000 000000      011000011000000000
0305 34304 00100001 1000 000000      001000011000000000
0306  3075 00000011 0000 000011      000000110000000011
0307 36928 00100100 0001 000000      001001000001000000
0308 12300 00001100 0000 001100      000011000000001100
0309 98816 01100000 1000 000000      011000001000000000
0310 98944 01100000 1010 000000      011000001010000000
0311 33408 00100000 1010 000000      001000001010000000
0312 66688 01000001 0010 000000      010000010010000000
0313 66176 01000000 1010 000000      010000001010000000
0314  7171 00000111 0000 000011      000001110000000011
0315 20528 00010100 0000 110000      000101000000110000
0316 33856 00100001 0001 000000      001000010001000000
0317 35840 00100011 0000 000000      001000110000000000
0318  4156 00000100 0000 111100      000001000000111100
0319 34880 00100010 0001 000000      001000100001000000
0320    */
0321   for (auto aHitPattern : badPatterns) {
0322     if (hits == aHitPattern)
0323       return false;
0324   }
0325 
0326   return true;
0327 }
0328 ///////////////////////////////////////////////////////
0329 ///////////////////////////////////////////////////////
0330 template <class GoldenPatternType>
0331 AlgoMuons OMTFProcessor<GoldenPatternType>::sortResults(unsigned int iProcessor, l1t::tftype mtfType, int charge) {
0332   unsigned int procIndx = this->myOmtfConfig->getProcIndx(iProcessor, mtfType);
0333   return sorter->sortResults(procIndx, this->getPatterns(), charge);
0334 }
0335 
0336 template <class GoldenPatternType>
0337 int OMTFProcessor<GoldenPatternType>::extrapolateDtPhiBFloatPoint(const int& refLogicLayer,
0338                                                                   const int& refPhi,
0339                                                                   const int& refPhiB,
0340                                                                   unsigned int targetLayer,
0341                                                                   const int& targetStubPhi,
0342                                                                   const int& targetStubQuality,
0343                                                                   const int& targetStubEta,
0344                                                                   const int& targetStubR,
0345                                                                   const OMTFConfiguration* omtfConfig) {
0346   LogTrace("l1tOmtfEventPrint") << "\n"
0347                                 << __FUNCTION__ << ":" << __LINE__ << " refLogicLayer " << refLogicLayer
0348                                 << " targetLayer " << targetLayer << std::endl;
0349   LogTrace("l1tOmtfEventPrint") << "refPhi " << refPhi << " refPhiB " << refPhiB << " targetStubPhi " << targetStubPhi
0350                                 << " targetStubQuality " << targetStubQuality << std::endl;
0351 
0352   int phiExtr = 0;  //delta phi extrapolated
0353 
0354   float rRefLayer = 431.133;  //[cm], MB1 i.e. refLogicLayer = 0
0355   if (refLogicLayer == 2)
0356     rRefLayer = 512.401;  //MB2
0357   else if (refLogicLayer != 0) {
0358     return 0;
0359     //throw cms::Exception("OMTFProcessor<GoldenPatternType>::extrapolateDtPhiB: wrong refStubLogicLayer " + std::to_string(refLogicLayer) );
0360   }
0361 
0362   int reflLayerIndex = refLogicLayer == 0 ? 0 : 1;
0363 
0364   if (targetLayer == 0 || targetLayer == 2 || targetLayer == 4 || (targetLayer >= 10 && targetLayer <= 14)) {
0365     //all units are cm. Values from the CMS geometry
0366     float rTargetLayer = 512.401;  //MB2
0367 
0368     if (targetLayer == 0)
0369       rTargetLayer = 431.133;  //MB1
0370     else if (targetLayer == 4)
0371       rTargetLayer = 617.946;  //MB3
0372 
0373     else if (targetLayer == 10)
0374       rTargetLayer = 413.675;  //RB1in
0375     else if (targetLayer == 11)
0376       rTargetLayer = 448.675;  //RB1out
0377     else if (targetLayer == 12)
0378       rTargetLayer = 494.975;  //RB2in
0379     else if (targetLayer == 13)
0380       rTargetLayer = 529.975;  //RB2out
0381     else if (targetLayer == 14)
0382       rTargetLayer = 602.150;  //RB3
0383 
0384     if (useStubQualInExtr) {
0385       if (targetLayer == 0 || targetLayer == 2 || targetLayer == 4) {
0386         if (targetStubQuality == 2 || targetStubQuality == 0)
0387           rTargetLayer = rTargetLayer - 23.5 / 2;  //inner superlayer
0388         else if (targetStubQuality == 3 || targetStubQuality == 1)
0389           rTargetLayer = rTargetLayer + 23.5 / 2;  //outer superlayer
0390       }
0391     }
0392 
0393     float d = rTargetLayer - rRefLayer;
0394     //formula in the form as in the slides explaining the extrapolation algorithm
0395     //float deltaPhiExtr = d/rTargetLayer * refPhiB / omtfConfig->dtPhiBUnitsRad(); //[rad]
0396     //phiExtr = round(deltaPhiExtr / omtfConfig->omtfPhiUnit()); //[halfStrip]
0397 
0398     //formula with approximation, used to calculate extrFactor
0399     float extrFactor = d / rTargetLayer / omtfConfig->dtPhiBUnitsRad() / omtfConfig->omtfPhiUnit();
0400     phiExtr = extrFactor * (float)refPhiB;  //[halfStrip]
0401 
0402     //formula without approximation
0403     float deltaPhiExtr = atan(d / rTargetLayer * tan(refPhiB / omtfConfig->dtPhiBUnitsRad()));  //[rad]
0404     phiExtr = round(deltaPhiExtr / omtfConfig->omtfPhiUnit());                                  //[halfStrip]
0405 
0406     if (useStubQualInExtr & (targetLayer == 0 || targetLayer == 2 || targetLayer == 4)) {
0407       extrapolFactors[reflLayerIndex][targetLayer][targetStubQuality] = extrFactor;
0408       extrapolFactorsNorm[reflLayerIndex][targetLayer][targetStubQuality] = 1;
0409     } else {
0410       extrapolFactors[reflLayerIndex][targetLayer][0] = extrFactor;
0411       extrapolFactorsNorm[reflLayerIndex][targetLayer][0] = 1;
0412     }
0413 
0414     //LogTrace("l1tOmtfEventPrint") <<__FUNCTION__<<":"<<__LINE__<<" deltaPhiExtr "<<deltaPhiExtr<<" phiExtr "<<phiExtr<<std::endl;
0415 
0416     LogTrace("l1tOmtfEventPrint") << "\n"
0417                                   << __FUNCTION__ << ":" << __LINE__ << " refLogicLayer " << refLogicLayer
0418                                   << " targetLayer " << std::setw(2) << targetLayer << " targetStubQuality "
0419                                   << targetStubQuality << " extrFactor " << extrFactor << std::endl;
0420 
0421     LogTrace("l1tOmtfEventPrint") << __FUNCTION__ << ":" << __LINE__ << " refPhiB " << refPhiB << " phiExtr " << phiExtr
0422                                   << std::endl;
0423 
0424   } else if (targetLayer == 1 || targetLayer == 3 || targetLayer == 5) {
0425     int deltaPhi = targetStubPhi - refPhi;  //[halfStrip]
0426 
0427     //deltaPhi is here in phi_b hw scale
0428     deltaPhi = round(deltaPhi * omtfConfig->omtfPhiUnit() * omtfConfig->dtPhiBUnitsRad());
0429 
0430     phiExtr = refPhiB - deltaPhi;  //phiExtr is also in phi_b hw scale
0431     LogTrace("l1tOmtfEventPrint") << __FUNCTION__ << ":" << __LINE__ << " deltaPhi " << deltaPhi << " phiExtr "
0432                                   << phiExtr << std::endl;
0433   } else if ((targetLayer >= 6 && targetLayer <= 9) || (targetLayer >= 15 && targetLayer <= 17)) {
0434     //if true, for the CSC and endcap RPC the R is taken from the hit coordinates
0435 
0436     float rME = targetStubR;
0437     if (!useEndcapStubsRInExtr) {
0438       //all units are cm. This are the average R values for a given chamber (more or less middle of the chamber, but taking into account the OMTF eta range)
0439       if (targetLayer == 6 || targetLayer == 15)  //ME1/3, RE1/3,
0440         rME = 600.;
0441       else if (targetLayer == 7 || targetLayer == 15) {  //ME2/2, RE2/3,
0442         if (refLogicLayer == 0)
0443           rME = 600.;
0444         else
0445           rME = 640.;
0446       } else if (targetLayer == 8 || rME == 16) {  //ME3/2, RE3/3,
0447         if (refLogicLayer == 0)
0448           rME = 620.;
0449         else
0450           rME = 680.;
0451       } else if (targetLayer == 9) {
0452         rME = 460.;  //for the refLogicLayer = 1. refLogicLayer = 2 is impossible
0453       }
0454     }
0455 
0456     float d = rME - rRefLayer;
0457     //formula in the form as in the slides explaining the extrapolation algorithm
0458     //float deltaPhiExtr = d / rME * refPhiB / omtfConfig->dtPhiBUnitsRad();  //[rad]
0459     //phiExtr = round(deltaPhiExtr / omtfConfig->omtfPhiUnit()); //[halfStrip]
0460 
0461     //formula with approximation, used to calculate extrFactor
0462     float extrFactor = d / rME / omtfConfig->dtPhiBUnitsRad() / omtfConfig->omtfPhiUnit();
0463     phiExtr = extrFactor * refPhiB;  //[halfStrip]
0464 
0465     //formula without approximation
0466     float deltaPhiExtr = atan(d / rME * tan(refPhiB / omtfConfig->dtPhiBUnitsRad()));  //[rad]
0467     phiExtr = round(deltaPhiExtr / omtfConfig->omtfPhiUnit());                         //[halfStrip]
0468 
0469     if (useEndcapStubsRInExtr) {
0470       extrapolFactors[reflLayerIndex][targetLayer][abs(targetStubEta)] += extrFactor;
0471       extrapolFactorsNorm[reflLayerIndex][targetLayer][abs(targetStubEta)]++;
0472       //extrapolFactors[reflLayerIndex][targetLayer][0] += extrFactor;
0473       //extrapolFactorsNorm[reflLayerIndex][targetLayer][0]++;
0474     } else {
0475       extrapolFactors[reflLayerIndex][targetLayer][0] = extrFactor;
0476       extrapolFactorsNorm[reflLayerIndex][targetLayer][0] = 1;
0477     }
0478     LogTrace("l1tOmtfEventPrint") << "\n"
0479                                   << __FUNCTION__ << ":" << __LINE__ << " refLogicLayer " << refLogicLayer
0480                                   << " targetLayer " << std::setw(2) << targetLayer << " targetStubR " << targetStubR
0481                                   << " targetStubEta " << targetStubEta << " extrFactor "
0482                                   << " rRefLayer " << rRefLayer << " d " << d << " deltaPhiExtr " << deltaPhiExtr
0483                                   << " phiExtr " << phiExtr << std::endl;
0484   }
0485   //TODO restrict the range of the phiExtr and refPhiB !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
0486 
0487   return phiExtr;
0488 }
0489 
0490 template <class GoldenPatternType>
0491 int OMTFProcessor<GoldenPatternType>::extrapolateDtPhiBFixedPoint(const int& refLogicLayer,
0492                                                                   const int& refPhi,
0493                                                                   const int& refPhiB,
0494                                                                   unsigned int targetLayer,
0495                                                                   const int& targetStubPhi,
0496                                                                   const int& targetStubQuality,
0497                                                                   const int& targetStubEta,
0498                                                                   const int& targetStubR,
0499                                                                   const OMTFConfiguration* omtfConfig) {
0500   int phiExtr = 0;  //delta phi extrapolated
0501 
0502   int reflLayerIndex = refLogicLayer == 0 ? 0 : 1;
0503   int extrFactor = 0;
0504 
0505   if (targetLayer == 0 || targetLayer == 2 || targetLayer == 4) {
0506     if (useStubQualInExtr)
0507       extrFactor = extrapolFactors[reflLayerIndex][targetLayer][targetStubQuality];
0508     else
0509       extrFactor = extrapolFactors[reflLayerIndex][targetLayer][0];
0510   } else if (targetLayer == 1 || targetLayer == 3 || targetLayer == 5) {
0511     int deltaPhi = targetStubPhi - refPhi;  //[halfStrip]
0512 
0513     int scaleFactor = this->myOmtfConfig->omtfPhiUnit() * this->myOmtfConfig->dtPhiBUnitsRad() * 512;
0514     //= 305 for phase-1, 512 is multiplier so that scaleFactor is non-zero integer
0515 
0516     deltaPhi = (deltaPhi * scaleFactor) / 512;  //here deltaPhi is converted to the phi_b hw scale
0517 
0518     phiExtr = refPhiB - deltaPhi;  //phiExtr is also in phi_b hw scale
0519     //LogTrace("l1tOmtfEventPrint") <<__FUNCTION__<<":"<<__LINE__<<" deltaPhi "<<deltaPhi<<" phiExtr "<<phiExtr<<std::endl;
0520 
0521   } else if (targetLayer >= 10 && targetLayer <= 14) {
0522     extrFactor = extrapolFactors[reflLayerIndex][targetLayer][0];
0523   } else if ((targetLayer >= 6 && targetLayer <= 9) || (targetLayer >= 15 && targetLayer <= 17)) {
0524     if (useEndcapStubsRInExtr) {
0525       //if given abs(targetStubEta) value is not present in the map, it is added with default value of 0
0526       //so it should be good. The only problem is that the map can grow...
0527       extrFactor = extrapolFactors[reflLayerIndex][targetLayer][abs(targetStubEta)];
0528     } else {
0529       extrFactor = extrapolFactors[reflLayerIndex][targetLayer][0];
0530     }
0531   }
0532 
0533   if (this->myOmtfConfig->isBendingLayer(targetLayer) == false) {
0534     phiExtr = extrFactor * refPhiB / extrapolMultiplier;
0535   }
0536 
0537   LogTrace("l1tOmtfEventPrint") << "\n"
0538                                 << __FUNCTION__ << ":" << __LINE__ << " refLogicLayer " << refLogicLayer
0539                                 << " targetLayer " << targetLayer << std::endl;
0540   LogTrace("l1tOmtfEventPrint") << "refPhi " << refPhi << " refPhiB " << refPhiB << " targetStubPhi " << targetStubPhi
0541                                 << " targetStubQuality " << targetStubQuality << " targetStubEta " << targetStubEta
0542                                 << " extrFactor " << extrFactor << " phiExtr " << phiExtr << std::endl;
0543 
0544   return phiExtr;
0545 }
0546 
0547 template <class GoldenPatternType>
0548 int OMTFProcessor<GoldenPatternType>::extrapolateDtPhiB(const MuonStubPtr& refStub,
0549                                                         const MuonStubPtr& targetStub,
0550                                                         unsigned int targetLayer,
0551                                                         const OMTFConfiguration* omtfConfig) {
0552   if (useFloatingPointExtrapolation)
0553     return OMTFProcessor<GoldenPatternType>::extrapolateDtPhiBFloatPoint(refStub->logicLayer,
0554                                                                          refStub->phiHw,
0555                                                                          refStub->phiBHw,
0556                                                                          targetLayer,
0557                                                                          targetStub->phiHw,
0558                                                                          targetStub->qualityHw,
0559                                                                          targetStub->etaHw,
0560                                                                          targetStub->r,
0561                                                                          omtfConfig);
0562   return OMTFProcessor<GoldenPatternType>::extrapolateDtPhiBFixedPoint(refStub->logicLayer,
0563                                                                        refStub->phiHw,
0564                                                                        refStub->phiBHw,
0565                                                                        targetLayer,
0566                                                                        targetStub->phiHw,
0567                                                                        targetStub->qualityHw,
0568                                                                        targetStub->etaHw,
0569                                                                        targetStub->r,
0570                                                                        omtfConfig);
0571 }
0572 ///////////////////////////////////////////////
0573 ///////////////////////////////////////////////
0574 //const std::vector<OMTFProcessor::resultsMap> &
0575 template <class GoldenPatternType>
0576 void OMTFProcessor<GoldenPatternType>::processInput(unsigned int iProcessor,
0577                                                     l1t::tftype mtfType,
0578                                                     const OMTFinput& aInput,
0579                                                     std::vector<std::unique_ptr<IOMTFEmulationObserver> >& observers) {
0580   unsigned int procIndx = this->myOmtfConfig->getProcIndx(iProcessor, mtfType);
0581   for (auto& itGP : this->theGPs) {
0582     for (auto& result : itGP->getResults()[procIndx]) {
0583       result.reset();
0584     }
0585   }
0586 
0587   LogTrace("l1tOmtfEventPrint") << __FUNCTION__ << "\n"
0588                                 << __LINE__ << " iProcessor " << iProcessor << " mtfType " << mtfType << " procIndx "
0589                                 << procIndx << " ----------------------" << std::endl;
0590   //////////////////////////////////////
0591   //////////////////////////////////////
0592   std::vector<const RefHitDef*> refHitDefs;
0593 
0594   {
0595     auto refHitsBits = aInput.getRefHits(iProcessor);
0596     if (refHitsBits.none())
0597       return;  // myResults;
0598 
0599     //loop over all possible refHits, e.g. 128
0600     for (unsigned int iRefHit = 0; iRefHit < this->myOmtfConfig->nRefHits(); ++iRefHit) {
0601       if (!refHitsBits[iRefHit])
0602         continue;
0603 
0604       refHitDefs.push_back(&(this->myOmtfConfig->getRefHitsDefs()[iProcessor][iRefHit]));
0605 
0606       if (refHitDefs.size() == this->myOmtfConfig->nTestRefHits())
0607         break;
0608     }
0609   }
0610 
0611   boost::property_tree::ptree procDataTree;
0612   LogTrace("l1tOmtfEventPrint") << __FUNCTION__ << " " << __LINE__;
0613   for (unsigned int iLayer = 0; iLayer < this->myOmtfConfig->nLayers(); ++iLayer) {
0614     //debug
0615     /*for(auto& h : layerHits) {
0616       if(h != 5400)
0617         LogTrace("l1tOmtfEventPrint")<<__FUNCTION__<<" "<<__LINE__<<" iLayer "<<iLayer<<" layerHit "<<h<<std::endl;
0618     }*/
0619 
0620     for (unsigned int iRefHit = 0; iRefHit < refHitDefs.size(); iRefHit++) {
0621       const RefHitDef& aRefHitDef = *(refHitDefs[iRefHit]);
0622 
0623       unsigned int refLayerLogicNum = this->myOmtfConfig->getRefToLogicNumber()[aRefHitDef.iRefLayer];
0624       const MuonStubPtr refStub = aInput.getMuonStub(refLayerLogicNum, aRefHitDef.iInput);
0625       //int etaRef = refStub->etaHw;
0626 
0627       unsigned int iRegion = aRefHitDef.iRegion;
0628 
0629       MuonStubPtrs1D restrictedLayerStubs = this->restrictInput(iProcessor, iRegion, iLayer, aInput);
0630 
0631       //LogTrace("l1tOmtfEventPrint")<<__FUNCTION__<<" "<<__LINE__<<" iLayer "<<iLayer<<" iRefLayer "<<aRefHitDef.iRefLayer<<std::endl;
0632       //LogTrace("l1tOmtfEventPrint")<<"iLayer "<<iLayer<<" iRefHit "<<iRefHit;
0633       //LogTrace("l1tOmtfEventPrint")<<" nTestedRefHits "<<nTestedRefHits<<" aRefHitDef "<<aRefHitDef<<std::endl;
0634 
0635       std::vector<int> extrapolatedPhi(restrictedLayerStubs.size(), 0);
0636 
0637       //TODO make sure the that the iRefLayer numbers used here corresponds to this in the hwToLogicLayer_0x000X.xml
0638       if ((this->myOmtfConfig->usePhiBExtrapolationMB1() && aRefHitDef.iRefLayer == 0) ||
0639           (this->myOmtfConfig->usePhiBExtrapolationMB2() && aRefHitDef.iRefLayer == 2)) {
0640         if ((iLayer != refLayerLogicNum) && (iLayer != refLayerLogicNum + 1)) {
0641           unsigned int iStub = 0;
0642           for (auto& targetStub : restrictedLayerStubs) {
0643             if (targetStub) {
0644               extrapolatedPhi[iStub] = extrapolateDtPhiB(refStub, targetStub, iLayer, this->myOmtfConfig);
0645 
0646               LogTrace("l1tOmtfEventPrint")
0647                   << "\n"
0648                   << __FUNCTION__ << ":" << __LINE__ << " extrapolating from layer " << refLayerLogicNum
0649                   << " - iRefLayer " << aRefHitDef.iRefLayer << " to layer " << iLayer << " stub " << targetStub
0650                   << " value " << extrapolatedPhi[iStub] << std::endl;
0651 
0652               if (this->myOmtfConfig->getDumpResultToXML()) {
0653                 auto& extrapolatedPhiTree = procDataTree.add_child("extrapolatedPhi", boost::property_tree::ptree());
0654                 extrapolatedPhiTree.add("<xmlattr>.refLayer", refLayerLogicNum);
0655                 extrapolatedPhiTree.add("<xmlattr>.layer", iLayer);
0656                 extrapolatedPhiTree.add("<xmlattr>.refPhiBHw", refStub->phiBHw);
0657                 extrapolatedPhiTree.add("<xmlattr>.iStub", iStub);
0658                 extrapolatedPhiTree.add("<xmlattr>.qualityHw", targetStub->qualityHw);
0659                 extrapolatedPhiTree.add("<xmlattr>.etaHw", targetStub->etaHw);
0660                 extrapolatedPhiTree.add("<xmlattr>.phiExtr", extrapolatedPhi[iStub]);
0661 
0662                 if (this->myOmtfConfig->isBendingLayer(iLayer))
0663                   extrapolatedPhiTree.add("<xmlattr>.dist_phi", targetStub->phiBHw - extrapolatedPhi[iStub]);
0664                 else
0665                   extrapolatedPhiTree.add("<xmlattr>.dist_phi", targetStub->phiHw - extrapolatedPhi[iStub]);
0666               }
0667             }
0668             iStub++;
0669           }
0670         }
0671       }
0672 
0673       for (auto& itGP : this->theGPs) {
0674         if (itGP->key().thePt == 0)  //empty pattern
0675           continue;
0676 
0677         StubResult stubResult =
0678             itGP->process1Layer1RefLayer(aRefHitDef.iRefLayer, iLayer, restrictedLayerStubs, extrapolatedPhi, refStub);
0679 
0680         /* LogTrace("l1tOmtfEventPrint")<<__FUNCTION__<<":"<<__LINE__
0681                                      <<" layerResult: valid"<<stubResult.getValid()
0682                                      <<" pdfVal "<<stubResult.getPdfVal()
0683                                      <<std::endl;*/
0684 
0685         itGP->getResults()[procIndx][iRefHit].setStubResult(iLayer, stubResult);
0686       }
0687     }
0688   }
0689 
0690   for (unsigned int iRefHit = 0; iRefHit < refHitDefs.size(); iRefHit++) {
0691     const RefHitDef& aRefHitDef = *(refHitDefs[iRefHit]);
0692 
0693     unsigned int refLayerLogicNum = this->myOmtfConfig->getRefToLogicNumber()[aRefHitDef.iRefLayer];
0694     const MuonStubPtr refStub = aInput.getMuonStub(refLayerLogicNum, aRefHitDef.iInput);
0695 
0696     int phiRef = refStub->phiHw;
0697     int etaRef = refStub->etaHw;
0698 
0699     //calculating the phiExtrp in the case the RefLayer is MB1, to include it in the  candidate phi of candidate
0700     int phiExtrp = 0;
0701     if ((this->myOmtfConfig->usePhiBExtrapolationMB1() && aRefHitDef.iRefLayer == 0)) {
0702       //||(this->myOmtfConfig->getUsePhiBExtrapolationMB2() && aRefHitDef.iRefLayer == 2) ) {  //the extrapolation from the layer 2 to the layer 2 has no sense, so phiExtrp is 0
0703       LogTrace("l1tOmtfEventPrint") << "\n"
0704                                     << __FUNCTION__ << ":" << __LINE__
0705                                     << "extrapolating ref hit to get the phi of the candidate" << std::endl;
0706       if (useFloatingPointExtrapolation)
0707         phiExtrp = extrapolateDtPhiBFloatPoint(
0708             aRefHitDef.iRefLayer, phiRef, refStub->phiBHw, 2, 0, 6, 0, 0, this->myOmtfConfig);
0709       else
0710         phiExtrp = extrapolateDtPhiBFixedPoint(
0711             aRefHitDef.iRefLayer, phiRef, refStub->phiBHw, 2, 0, 6, 0, 0, this->myOmtfConfig);
0712     }
0713 
0714     for (auto& itGP : this->theGPs) {
0715       if (itGP->key().thePt == 0)  //empty pattern
0716         continue;
0717 
0718       int phiRefSt2 = itGP->propagateRefPhi(phiRef + phiExtrp, etaRef, aRefHitDef.iRefLayer);
0719       itGP->getResults()[procIndx][iRefHit].set(aRefHitDef.iRefLayer, phiRefSt2, etaRef, phiRef);
0720     }
0721   }
0722 
0723   //////////////////////////////////////
0724   //////////////////////////////////////
0725   {
0726     for (auto& itGP : this->theGPs) {
0727       itGP->finalise(procIndx);
0728       //debug
0729       /*for(unsigned int iRefHit = 0; iRefHit < itGP->getResults()[procIndx].size(); ++iRefHit) {
0730         if(itGP->getResults()[procIndx][iRefHit].isValid()) {
0731           LogTrace("l1tOmtfEventPrint")<<__FUNCTION__<<":"<<"__LINE__"<<itGP->getResults()[procIndx][iRefHit]<<std::endl;
0732         }
0733       }*/
0734     }
0735   }
0736 
0737   for (auto& obs : observers)
0738     obs->addProcesorData("extrapolation", procDataTree);
0739 
0740   return;
0741 }
0742 ///////////////////////////////////////////////////////
0743 ///////////////////////////////////////////////////////
0744 
0745 template <class GoldenPatternType>
0746 std::vector<l1t::RegionalMuonCand> OMTFProcessor<GoldenPatternType>::run(
0747     unsigned int iProcessor,
0748     l1t::tftype mtfType,
0749     int bx,
0750     OMTFinputMaker* inputMaker,
0751     std::vector<std::unique_ptr<IOMTFEmulationObserver> >& observers) {
0752   //uncomment if you want to check execution time of each method
0753   //boost::timer::auto_cpu_timer t("%ws wall, %us user in getProcessorCandidates\n");
0754 
0755   for (auto& obs : observers)
0756     obs->observeProcesorBegin(iProcessor, mtfType);
0757 
0758   //input is shared_ptr because the observers may need them after the run() method execution is finished
0759   std::shared_ptr<OMTFinput> input = std::make_shared<OMTFinput>(this->myOmtfConfig);
0760   inputMaker->buildInputForProcessor(input->getMuonStubs(), iProcessor, mtfType, bx, bx, observers);
0761 
0762   if (this->myOmtfConfig->cleanStubs()) {
0763     //this has sense for the pattern generation from the tracks with the secondaries
0764     //if more than one stub is in a given layer, all stubs are removed from this layer
0765     for (unsigned int iLayer = 0; iLayer < input->getMuonStubs().size(); ++iLayer) {
0766       auto& layerStubs = input->getMuonStubs()[iLayer];
0767       int count = std::count_if(layerStubs.begin(), layerStubs.end(), [](auto& ptr) { return ptr != nullptr; });
0768       if (count > 1) {
0769         for (auto& ptr : layerStubs)
0770           ptr.reset();
0771 
0772         LogTrace("OMTFReconstruction") << __FUNCTION__ << ":" << __LINE__ << "cleaning stubs in the layer " << iLayer
0773                                        << " stubs count :" << count << std::endl;
0774       }
0775     }
0776   }
0777 
0778   //LogTrace("l1tOmtfEventPrint")<<"buildInputForProce "; t.report();
0779   processInput(iProcessor, mtfType, *(input.get()), observers);
0780 
0781   //LogTrace("l1tOmtfEventPrint")<<"processInput       "; t.report();
0782   AlgoMuons algoCandidates = sortResults(iProcessor, mtfType);
0783 
0784   if (ptAssignment) {
0785     for (auto& myCand : algoCandidates) {
0786       if (myCand->isValid()) {
0787         auto pts = ptAssignment->getPts(myCand, observers);
0788         /*for (unsigned int i = 0; i < pts.size(); i++) {
0789         trackAddr[10 + i] = this->myOmtfConfig->ptGevToHw(pts[i]);
0790       }*/
0791       }
0792     }
0793   }
0794 
0795   //LogTrace("l1tOmtfEventPrint")<<"sortResults        "; t.report();
0796   // perform GB
0797   //watch out: etaBits2HwEta is used in the ghostBust to convert the AlgoMuons eta, it affect algoCandidates as they are pointers
0798   AlgoMuons gbCandidates = ghostBust(algoCandidates);
0799 
0800   //LogTrace("l1tOmtfEventPrint")<<"ghostBust"; t.report();
0801   // fill RegionalMuonCand colleciton
0802   std::vector<l1t::RegionalMuonCand> candMuons = getFinalcandidates(iProcessor, mtfType, gbCandidates);
0803 
0804   //LogTrace("l1tOmtfEventPrint")<<"getFinalcandidates "; t.report();
0805   //fill outgoing collection
0806   for (auto& candMuon : candMuons) {
0807     candMuon.setHwQual(candMuon.hwQual());
0808   }
0809 
0810   for (auto& obs : observers) {
0811     obs->observeProcesorEmulation(iProcessor, mtfType, input, algoCandidates, gbCandidates, candMuons);
0812   }
0813 
0814   return candMuons;
0815 }
0816 
0817 template <class GoldenPatternType>
0818 void OMTFProcessor<GoldenPatternType>::printInfo() const {
0819   edm::LogVerbatim("OMTFReconstruction") << __PRETTY_FUNCTION__ << std::endl;
0820 
0821   ProcessorBase<GoldenPatternType>::printInfo();
0822 }
0823 
0824 template <class GoldenPatternType>
0825 void OMTFProcessor<GoldenPatternType>::saveExtrapolFactors() {
0826   //if(this->myOmtfConfig->nProcessors() == 3) //phase2
0827   extrapolMultiplier = 512;
0828 
0829   boost::property_tree::ptree tree;
0830   auto& extrFactorsTree = tree.add("ExtrapolationFactors", "");
0831   extrFactorsTree.add("<xmlattr>.multiplier", extrapolMultiplier);
0832 
0833   edm::LogVerbatim("OMTFReconstruction") << "saving extrapolFactors to ExtrapolationFactors.xml" << std::endl;
0834   for (unsigned int iRefLayer = 0; iRefLayer < extrapolFactors.size(); iRefLayer++) {
0835     for (unsigned int iLayer = 0; iLayer < extrapolFactors[iRefLayer].size(); iLayer++) {
0836       edm::LogVerbatim("OMTFReconstruction") << " iRefLayer " << iRefLayer << " iLayer " << iLayer << std::endl;
0837 
0838       auto& layerTree = extrFactorsTree.add_child("Lut", boost::property_tree::ptree());
0839       layerTree.add("<xmlattr>.RefLayer", std::to_string(iRefLayer));
0840       layerTree.add("<xmlattr>.Layer", iLayer);
0841 
0842       if (useStubQualInExtr && (iLayer == 0 || iLayer == 2 || iLayer == 4))
0843         layerTree.add("<xmlattr>.KeyType", "quality");
0844       else if (useEndcapStubsRInExtr && ((iLayer >= 6 && iLayer <= 9) || (iLayer >= 15 && iLayer <= 17)))
0845         layerTree.add("<xmlattr>.KeyType", "eta");
0846       else
0847         layerTree.add("<xmlattr>.KeyType", "none");
0848 
0849       for (auto& extrFactors : extrapolFactors[iRefLayer][iLayer]) {
0850         int norm = 1;
0851         if (!extrapolFactorsNorm[iRefLayer][iLayer].empty())
0852           norm = extrapolFactorsNorm[iRefLayer][iLayer][extrFactors.first];
0853         auto& lutVal = layerTree.add_child("LutVal", boost::property_tree::ptree());
0854         if (useEndcapStubsRInExtr && ((iLayer >= 6 && iLayer <= 9) || (iLayer >= 15 && iLayer <= 17)))
0855           lutVal.add("<xmlattr>.key", extrFactors.first);
0856         else
0857           lutVal.add("<xmlattr>.key", extrFactors.first);
0858 
0859         double value = round(extrapolMultiplier * extrFactors.second / norm);
0860         lutVal.add("<xmlattr>.value", value);
0861 
0862         edm::LogVerbatim("OMTFReconstruction")
0863             << std::setw(4) << " key = " << extrFactors.first << " extrFactors.second " << std::setw(10)
0864             << extrFactors.second << " norm " << std::setw(6) << norm << " value/norm " << std::setw(10)
0865             << extrFactors.second / norm << " value " << value << std::endl;
0866       }
0867     }
0868   }
0869 
0870   boost::property_tree::write_xml("ExtrapolationFactors.xml",
0871                                   tree,
0872                                   std::locale(),
0873                                   boost::property_tree::xml_parser::xml_writer_make_settings<std::string>(' ', 2));
0874 }
0875 
0876 template <class GoldenPatternType>
0877 void OMTFProcessor<GoldenPatternType>::loadExtrapolFactors(const std::string& filename) {
0878   boost::property_tree::ptree tree;
0879 
0880   boost::property_tree::read_xml(filename, tree);
0881 
0882   edm::LogVerbatim("OMTFReconstruction") << "loadExtrapolFactors from file " << filename << std::endl;
0883 
0884   extrapolMultiplier = tree.get<int>("ExtrapolationFactors.<xmlattr>.multiplier");
0885   edm::LogVerbatim("OMTFReconstruction") << "extrapolMultiplier " << extrapolMultiplier << std::endl;
0886 
0887   auto& lutNodes = tree.get_child("ExtrapolationFactors");
0888   for (boost::property_tree::ptree::value_type& lutNode : lutNodes) {
0889     if (lutNode.first == "Lut") {
0890       int iRefLayer = lutNode.second.get<int>("<xmlattr>.RefLayer");
0891       int iLayer = lutNode.second.get<int>("<xmlattr>.Layer");
0892       std::string keyType = lutNode.second.get<std::string>("<xmlattr>.KeyType");
0893 
0894       edm::LogVerbatim("OMTFReconstruction")
0895           << "iRefLayer " << iRefLayer << " iLayer " << iLayer << " keyType " << keyType << std::endl;
0896 
0897       auto& valueNodes = lutNode.second;
0898       for (boost::property_tree::ptree::value_type& valueNode : valueNodes) {
0899         if (valueNode.first == "LutVal") {
0900           int key = valueNode.second.get<int>("<xmlattr>.key");
0901           float value = valueNode.second.get<float>("<xmlattr>.value");
0902           extrapolFactors.at(iRefLayer).at(iLayer)[key] = value;
0903           edm::LogVerbatim("OMTFReconstruction") << "key " << key << " value " << value << std::endl;
0904         }
0905       }
0906     }
0907   }
0908 }
0909 
0910 /////////////////////////////////////////////////////////
0911 
0912 template class OMTFProcessor<GoldenPattern>;
0913 template class OMTFProcessor<GoldenPatternWithStat>;
0914 template class OMTFProcessor<GoldenPatternWithThresh>;