CMS 3D CMS Logo

OMTFProcessor.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <algorithm>
3 #include <strstream>
4 
9 
11 
16 
18 
23  for (auto it : theGPs)
24  delete it.second;
25 }
29  myResults.clear();
30  for (auto it : theGPs)
31  delete it.second;
32  theGPs.clear();
33 }
36 bool OMTFProcessor::configure(const OMTFConfiguration *omtfConfig, const L1TMuonOverlapParams *omtfPatterns) {
38 
39  myOmtfConfig = omtfConfig;
40 
42 
43  const l1t::LUT *chargeLUT = omtfPatterns->chargeLUT();
44  const l1t::LUT *etaLUT = omtfPatterns->etaLUT();
45  const l1t::LUT *ptLUT = omtfPatterns->ptLUT();
46  const l1t::LUT *pdfLUT = omtfPatterns->pdfLUT();
47  const l1t::LUT *meanDistPhiLUT = omtfPatterns->meanDistPhiLUT();
48 
49  unsigned int nGPs = myOmtfConfig->nGoldenPatterns();
50  unsigned int address = 0;
51  unsigned int iEta, iPt;
52  int iCharge;
53  for (unsigned int iGP = 0; iGP < nGPs; ++iGP) {
54  address = iGP;
55  iEta = etaLUT->data(address);
56  iCharge = chargeLUT->data(address) == 0 ? -1 : 1;
57  iPt = ptLUT->data(address);
58 
64  for (unsigned int iLayer = 0; iLayer < myOmtfConfig->nLayers(); ++iLayer) {
66  for (unsigned int iRefLayer = 0; iRefLayer < myOmtfConfig->nRefLayers(); ++iRefLayer) {
67  address = iRefLayer + iLayer * myOmtfConfig->nRefLayers() +
69  meanDistPhi1D[iRefLayer] = meanDistPhiLUT->data(address) - (1 << (meanDistPhiLUT->nrBitsData() - 1));
70  }
71  meanDistPhi2D[iLayer] = meanDistPhi1D;
73  for (unsigned int iRefLayer = 0; iRefLayer < myOmtfConfig->nRefLayers(); ++iRefLayer) {
74  pdf1D.assign(1 << myOmtfConfig->nPdfAddrBits(), 0);
75  for (unsigned int iPdf = 0; iPdf < (unsigned int)(1 << myOmtfConfig->nPdfAddrBits()); ++iPdf) {
76  address = iPdf + iRefLayer * (1 << myOmtfConfig->nPdfAddrBits()) +
77  iLayer * myOmtfConfig->nRefLayers() * (1 << myOmtfConfig->nPdfAddrBits()) +
79  pdf1D[iPdf] = pdfLUT->data(address);
80  }
81  pdf2D[iRefLayer] = pdf1D;
82  }
83  pdf3D[iLayer] = pdf2D;
84  }
85  Key aKey(iEta, iPt, iCharge, iGP);
86 
87  GoldenPattern *aGP = new GoldenPattern(aKey, myOmtfConfig);
88  aGP->setMeanDistPhi(meanDistPhi2D);
89  aGP->setPdf(pdf3D);
90  addGP(aGP);
91  }
92  return true;
93 }
97  if (theGPs.find(aGP->key()) != theGPs.end()) {
98  throw cms::Exception("Corrupted Golden Patterns data")
99  << "OMTFProcessor::addGP(...) "
100  << " Reading two Golden Patterns with the same key: " << aGP->key() << std::endl;
101  } else
102  theGPs[aGP->key()] = aGP;
103 
104  for (auto &itRegion : myResults) {
105  OMTFResult aResult;
106  aResult.configure(myOmtfConfig);
107  itRegion[aGP->key()] = aResult;
108  }
109 
110  return true;
111 }
115  Key aKey(0, 9, charge);
116 
117  while (theGPs.find(aKey) != theGPs.end()) {
118  GoldenPattern *aGP1 = theGPs.find(aKey)->second;
119  GoldenPattern *aGP2 = aGP1;
120  GoldenPattern *aGP3 = aGP1;
121  GoldenPattern *aGP4 = aGP1;
122 
123  ++aKey.thePtCode;
124  while (theGPs.find(aKey) == theGPs.end() && aKey.thePtCode <= 401)
125  ++aKey.thePtCode;
126  if (aKey.thePtCode <= 401 && theGPs.find(aKey) != theGPs.end())
127  aGP2 = theGPs.find(aKey)->second;
128 
129  if (aKey.thePtCode > 71) {
130  ++aKey.thePtCode;
131  while (theGPs.find(aKey) == theGPs.end() && aKey.thePtCode <= 401)
132  ++aKey.thePtCode;
133  if (aKey.thePtCode <= 401 && theGPs.find(aKey) != theGPs.end())
134  aGP3 = theGPs.find(aKey)->second;
135 
136  ++aKey.thePtCode;
137  while (theGPs.find(aKey) == theGPs.end() && aKey.thePtCode <= 401)
138  ++aKey.thePtCode;
139  if (aKey.thePtCode <= 401 && theGPs.find(aKey) != theGPs.end())
140  aGP4 = theGPs.find(aKey)->second;
141  } else {
142  aGP3 = aGP1;
143  aGP4 = aGP2;
144  }
145  //HACK. Have to clean this up.
147  ++aKey.thePtCode;
148  while (theGPs.find(aKey) == theGPs.end() && aKey.thePtCode <= 401)
149  ++aKey.thePtCode;
151 
152  GoldenPattern::vector2D meanDistPhi = aGP1->getMeanDistPhi();
153 
154  GoldenPattern::vector2D meanDistPhi1 = aGP1->getMeanDistPhi();
155  GoldenPattern::vector2D meanDistPhi2 = aGP2->getMeanDistPhi();
156  GoldenPattern::vector2D meanDistPhi3 = aGP3->getMeanDistPhi();
157  GoldenPattern::vector2D meanDistPhi4 = aGP4->getMeanDistPhi();
158 
159  for (unsigned int iLayer = 0; iLayer < myOmtfConfig->nLayers(); ++iLayer) {
160  for (unsigned int iRefLayer = 0; iRefLayer < myOmtfConfig->nRefLayers(); ++iRefLayer) {
161  meanDistPhi[iLayer][iRefLayer] += meanDistPhi2[iLayer][iRefLayer];
162  meanDistPhi[iLayer][iRefLayer] += meanDistPhi3[iLayer][iRefLayer];
163  meanDistPhi[iLayer][iRefLayer] += meanDistPhi4[iLayer][iRefLayer];
164  meanDistPhi[iLayer][iRefLayer] /= 4;
165  }
166  }
167 
168  aGP1->setMeanDistPhi(meanDistPhi);
169  aGP2->setMeanDistPhi(meanDistPhi);
170 
171  shiftGP(aGP1, meanDistPhi, meanDistPhi1);
172  shiftGP(aGP2, meanDistPhi, meanDistPhi2);
173  if (aGP3 != aGP1 && aGP4 != aGP2) {
174  aGP3->setMeanDistPhi(meanDistPhi);
175  aGP4->setMeanDistPhi(meanDistPhi);
176  shiftGP(aGP3, meanDistPhi, meanDistPhi3);
177  shiftGP(aGP4, meanDistPhi, meanDistPhi4);
178  }
179  }
180 }
184  const GoldenPattern::vector2D &meanDistPhiNew,
185  const GoldenPattern::vector2D &meanDistPhiOld) {
188  unsigned int nPdfBins = exp2(myOmtfConfig->nPdfAddrBits());
189  GoldenPattern::vector3D pdfAllRef = aGP->getPdf();
190 
191  int indexShift = 0;
192  for (unsigned int iLayer = 0; iLayer < myOmtfConfig->nLayers(); ++iLayer) {
193  for (unsigned int iRefLayer = 0; iRefLayer < myOmtfConfig->nRefLayers(); ++iRefLayer) {
194  indexShift = meanDistPhiOld[iLayer][iRefLayer] - meanDistPhiNew[iLayer][iRefLayer];
195  for (unsigned int iPdfBin = 0; iPdfBin < nPdfBins; ++iPdfBin)
196  pdfAllRef[iLayer][iRefLayer][iPdfBin] = 0;
197  for (unsigned int iPdfBin = 0; iPdfBin < nPdfBins; ++iPdfBin) {
198  if ((int)(iPdfBin) + indexShift >= 0 && iPdfBin + indexShift < nPdfBins)
199  pdfAllRef[iLayer][iRefLayer][iPdfBin + indexShift] = aGP->pdfValue(iLayer, iRefLayer, iPdfBin);
200  }
201  }
202  }
203  aGP->setPdf(pdfAllRef);
204 }
207 const std::map<Key, GoldenPattern *> &OMTFProcessor::getPatterns() const { return theGPs; }
210 const std::vector<OMTFProcessor::resultsMap> &OMTFProcessor::processInput(unsigned int iProcessor,
211  const OMTFinput &aInput) {
212  for (auto &itRegion : myResults)
213  for (auto &itKey : itRegion)
214  itKey.second.clear();
215 
218  std::bitset<128> refHitsBits = aInput.getRefHits(iProcessor);
219  if (refHitsBits.none())
220  return myResults;
221 
222  for (unsigned int iLayer = 0; iLayer < myOmtfConfig->nLayers(); ++iLayer) {
223  const OMTFinput::vector1D &layerHits = aInput.getLayerData(iLayer);
224  if (layerHits.empty())
225  continue;
227  unsigned int nTestedRefHits = myOmtfConfig->nTestRefHits();
228  for (unsigned int iRefHit = 0; iRefHit < myOmtfConfig->nRefHits(); ++iRefHit) {
229  if (!refHitsBits[iRefHit])
230  continue;
231  if (nTestedRefHits-- == 0)
232  break;
233  const RefHitDef &aRefHitDef = myOmtfConfig->getRefHitsDefs()[iProcessor][iRefHit];
234 
235  int phiRef = aInput.getLayerData(myOmtfConfig->getRefToLogicNumber()[aRefHitDef.iRefLayer])[aRefHitDef.iInput];
236  int etaRef =
237  aInput.getLayerData(myOmtfConfig->getRefToLogicNumber()[aRefHitDef.iRefLayer], true)[aRefHitDef.iInput];
238  unsigned int iRegion = aRefHitDef.iRegion;
239 
240  if (myOmtfConfig->getBendingLayers().count(iLayer))
241  phiRef = 0;
242  const OMTFinput::vector1D restrictedLayerHits = restrictInput(iProcessor, iRegion, iLayer, layerHits);
243  for (auto itGP : theGPs) {
244  GoldenPattern::layerResult aLayerResult =
245  itGP.second->process1Layer1RefLayer(aRefHitDef.iRefLayer, iLayer, phiRef, restrictedLayerHits);
246  int phiRefSt2 = itGP.second->propagateRefPhi(phiRef, etaRef, aRefHitDef.iRefLayer);
247  myResults[myOmtfConfig->nTestRefHits() - nTestedRefHits - 1][itGP.second->key()].setRefPhiRHits(
248  aRefHitDef.iRefLayer, phiRef);
249  myResults[myOmtfConfig->nTestRefHits() - nTestedRefHits - 1][itGP.second->key()].addResult(
250  aRefHitDef.iRefLayer, iLayer, aLayerResult.first, phiRefSt2, etaRef);
251  }
252  }
253  }
256  for (auto &itRefHit : myResults)
257  for (auto &itKey : itRefHit)
258  itKey.second.finalise();
259 
260  std::ostringstream myStr;
261  myStr << "iProcessor: " << iProcessor << std::endl;
262  myStr << "Input: ------------" << std::endl;
263  myStr << aInput << std::endl;
264  edm::LogInfo("OMTF processor") << myStr.str();
265 
266  return myResults;
267 }
271  unsigned int iRegion,
272  unsigned int iLayer,
273  const OMTFinput::vector1D &layerHits) {
274  OMTFinput::vector1D myHits = layerHits;
275 
276  unsigned int iStart = myOmtfConfig->getConnections()[iProcessor][iRegion][iLayer].first;
277  unsigned int iEnd = iStart + myOmtfConfig->getConnections()[iProcessor][iRegion][iLayer].second - 1;
278 
279  for (unsigned int iInput = 0; iInput < 14; ++iInput) {
280  if (iInput < iStart || iInput > iEnd)
281  myHits[iInput] = myOmtfConfig->nPhiBins();
282  }
283  return myHits;
284 }
287 void OMTFProcessor::fillCounts(unsigned int iProcessor, const OMTFinput &aInput, const SimTrack *aSimMuon) {
288  int theCharge = (abs(aSimMuon->type()) == 13) ? aSimMuon->type() / -13 : 0;
289  unsigned int iPt = RPCConst::iptFromPt(aSimMuon->momentum().pt());
292  iPt += 1;
293  if (iPt > 31)
294  iPt = 200 * 2 + 1;
295  else
296  iPt = RPCConst::ptFromIpt(iPt) * 2.0 +
297  1; //MicroGMT has 0.5 GeV step size, with lower bin edge (uGMT_pt_code - 1)*step_size
299 
301  std::bitset<128> refHitsBits = aInput.getRefHits(iProcessor);
302  if (refHitsBits.none())
303  return;
304 
305  std::ostringstream myStr;
306  myStr << "iProcessor: " << iProcessor << std::endl;
307  myStr << "Input: ------------" << std::endl;
308  myStr << aInput << std::endl;
309  edm::LogInfo("OMTF processor") << myStr.str();
310 
311  for (unsigned int iLayer = 0; iLayer < myOmtfConfig->nLayers(); ++iLayer) {
312  const OMTFinput::vector1D &layerHits = aInput.getLayerData(iLayer);
313  if (layerHits.empty())
314  continue;
316  for (unsigned int iRefHit = 0; iRefHit < myOmtfConfig->nRefHits(); ++iRefHit) {
317  if (!refHitsBits[iRefHit])
318  continue;
319  const RefHitDef &aRefHitDef = myOmtfConfig->getRefHitsDefs()[iProcessor][iRefHit];
320  int phiRef = aInput.getLayerData(myOmtfConfig->getRefToLogicNumber()[aRefHitDef.iRefLayer])[aRefHitDef.iInput];
321  unsigned int iRegion = aRefHitDef.iRegion;
322  if (myOmtfConfig->getBendingLayers().count(iLayer))
323  phiRef = 0;
324  const OMTFinput::vector1D restrictedLayerHits = restrictInput(iProcessor, iRegion, iLayer, layerHits);
325  for (auto itGP : theGPs) {
326  if (itGP.first.theCharge != theCharge)
327  continue;
328  if (itGP.first.thePtCode != iPt)
329  continue;
330  itGP.second->addCount(aRefHitDef.iRefLayer, iLayer, phiRef, restrictedLayerHits);
331  }
332  }
333  }
334 }
CoreSimTrack::momentum
const math::XYZTLorentzVectorD & momentum() const
Definition: CoreSimTrack.h:19
OMTFResult::configure
void configure(const OMTFConfiguration *omtfConfig)
Definition: OMTFResult.cc:10
OMTFProcessor::fillCounts
void fillCounts(unsigned int iProcessor, const OMTFinput &aInput, const SimTrack *aSimMuon)
Definition: OMTFProcessor.cc:287
MessageLogger.h
OMTFProcessor.h
OMTFinput::getRefHits
std::bitset< 128 > getRefHits(unsigned int iProcessor) const
Definition: OMTFinput.cc:26
GoldenPattern::setMeanDistPhi
void setMeanDistPhi(const vector2D &aMeanDistPhi)
Definition: GoldenPattern.h:58
RefHitDef
Definition: OMTFConfiguration.h:16
OMTFProcessor::resetConfiguration
void resetConfiguration()
Reset all configuration parameters.
Definition: OMTFProcessor.cc:28
OMTFProcessor::theGPs
std::map< Key, GoldenPattern * > theGPs
Map holding Golden Patterns.
Definition: OMTFProcessor.h:77
GoldenPattern::setPdf
void setPdf(const vector3D &aPdf)
Definition: GoldenPattern.h:64
OMTFConfiguration::getRefHitsDefs
const std::vector< std::vector< RefHitDef > > & getRefHitsDefs() const
Definition: OMTFConfiguration.h:118
Key::thePtCode
unsigned int thePtCode
Definition: GoldenPattern.h:35
L1TMuonOverlapParams::ptLUT
const l1t::LUT * ptLUT() const
Definition: L1TMuonOverlapParams.h:204
OMTFConfiguration::nRefHits
unsigned int nRefHits() const
Definition: OMTFConfiguration.h:94
if
if(0==first)
Definition: CAHitNtupletGeneratorKernelsImpl.h:58
OMTFProcessor::getPatterns
const std::map< Key, GoldenPattern * > & getPatterns() const
Return map of GoldenPatterns.
Definition: OMTFProcessor.cc:207
OMTFConfiguration
Definition: OMTFConfiguration.h:44
GoldenPattern::vector1D
std::vector< int > vector1D
Definition: GoldenPattern.h:46
edm::LogInfo
Log< level::Info, false > LogInfo
Definition: MessageLogger.h:125
OMTFConfiguration::nPhiBins
unsigned int nPhiBins() const
Definition: OMTFConfiguration.h:93
OMTFConfiguration::nGoldenPatterns
unsigned int nGoldenPatterns() const
Definition: OMTFConfiguration.h:99
GoldenPattern::getMeanDistPhi
const vector2D & getMeanDistPhi() const
Definition: GoldenPattern.h:60
RPCConst.h
RefHitDef::iInput
unsigned int iInput
Hit input number within a cone.
Definition: OMTFConfiguration.h:29
OMTFProcessor::addGP
bool addGP(GoldenPattern *aGP)
Definition: OMTFProcessor.cc:96
OMTFConfiguration::getConnections
const vector3D_pair & getConnections() const
Definition: OMTFConfiguration.h:120
Key
Definition: GoldenPattern.h:15
RPCConst::iptFromPt
static int iptFromPt(const double pt)
Definition: RPCConst.cc:10
GoldenPattern::vector3D
std::vector< vector2D > vector3D
Definition: GoldenPattern.h:48
OMTFProcessor::shiftGP
void shiftGP(GoldenPattern *aGP, const GoldenPattern::vector2D &meanDistPhiNew, const GoldenPattern::vector2D &meanDistPhiOld)
Definition: OMTFProcessor.cc:183
L1TMuonOverlapParams::meanDistPhiLUT
const l1t::LUT * meanDistPhiLUT() const
Definition: L1TMuonOverlapParams.h:206
OMTFProcessor::~OMTFProcessor
~OMTFProcessor()
Definition: OMTFProcessor.cc:22
OMTFProcessor::myResults
std::vector< OMTFProcessor::resultsMap > myResults
Definition: OMTFProcessor.h:82
OMTFConfiguration::nPdfAddrBits
unsigned int nPdfAddrBits() const
Definition: OMTFConfiguration.h:91
L1TMuonOverlapParams
Definition: L1TMuonOverlapParams.h:14
OMTFProcessor::resultsMap
std::map< Key, OMTFResult > resultsMap
Definition: OMTFProcessor.h:21
L1TMuonOverlapParams.h
GoldenPattern::pdfValue
int pdfValue(unsigned int iLayer, unsigned int iRefLayer, unsigned int iBin) const
Definition: GoldenPattern.h:68
OMTFResult.h
FileInPath.h
ALCARECOTkAlJpsiMuMu_cff.charge
charge
Definition: ALCARECOTkAlJpsiMuMu_cff.py:47
L1TMuonOverlapParams::pdfLUT
const l1t::LUT * pdfLUT() const
Definition: L1TMuonOverlapParams.h:205
OMTFConfiguration::nLayers
unsigned int nLayers() const
Definition: OMTFConfiguration.h:87
GoldenPattern::key
Key key() const
Definition: GoldenPattern.h:56
OMTFConfiguration::getRefToLogicNumber
const std::vector< int > & getRefToLogicNumber() const
Definition: OMTFConfiguration.h:105
OMTFinput::vector1D
std::vector< int > vector1D
Definition: OMTFinput.h:13
l1t::LUT::data
int data(unsigned int address) const
Definition: LUT.h:46
createfilelist.int
int
Definition: createfilelist.py:10
RefHitDef::iRefLayer
unsigned int iRefLayer
Reference layer logic number (0-7)
Definition: OMTFConfiguration.h:35
CoreSimTrack::type
int type() const
particle type (HEP PDT convension)
Definition: CoreSimTrack.h:22
OMTFProcessor::myOmtfConfig
const OMTFConfiguration * myOmtfConfig
Definition: OMTFProcessor.h:86
OMTFConfiguration::nTestRefHits
unsigned int nTestRefHits() const
Definition: OMTFConfiguration.h:95
AlignmentTrackSelector_cfi.theCharge
theCharge
Definition: AlignmentTrackSelector_cfi.py:20
GoldenPattern
Definition: GoldenPattern.h:44
l1t::LUT
Definition: LUT.h:29
OMTFProcessor::configure
bool configure(const OMTFConfiguration *omtfParams, const L1TMuonOverlapParams *omtfPatterns)
Fill GP map with patterns from CondFormats object.
Definition: OMTFProcessor.cc:36
SimTrack
Definition: SimTrack.h:9
OMTFConfiguration::nRefLayers
unsigned int nRefLayers() const
Definition: OMTFConfiguration.h:89
Exception
Definition: hltDiff.cc:245
OMTFProcessor::processInput
const std::vector< OMTFProcessor::resultsMap > & processInput(unsigned int iProcessor, const OMTFinput &aInput)
Definition: OMTFProcessor.cc:210
GoldenPattern::layerResult
std::pair< int, bool > layerResult
Definition: GoldenPattern.h:49
Exception.h
OMTFConfiguration::getBendingLayers
const std::set< int > & getBendingLayers() const
Definition: OMTFConfiguration.h:104
GoldenPattern::getPdf
const vector3D & getPdf() const
Definition: GoldenPattern.h:62
GoldenPattern.h
L1TowerCalibrationProducer_cfi.iEta
iEta
Definition: L1TowerCalibrationProducer_cfi.py:60
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
SimTrack.h
ParameterSet.h
OMTFinput::getLayerData
const OMTFinput::vector1D & getLayerData(unsigned int iLayer, bool giveEta=false) const
Definition: OMTFinput.cc:17
OMTFinput.h
L1TMuonOverlapParams::chargeLUT
const l1t::LUT * chargeLUT() const
Golden Patterns definitions.
Definition: L1TMuonOverlapParams.h:202
OMTFinput
Definition: OMTFinput.h:11
GoldenPattern::vector2D
std::vector< vector1D > vector2D
Definition: GoldenPattern.h:47
RefHitDef::iRegion
unsigned int iRegion
Region number assigned to this referecne hit.
Definition: OMTFConfiguration.h:32
l1t::LUT::nrBitsData
unsigned int nrBitsData() const
Definition: LUT.h:53
RPCConst::ptFromIpt
static double ptFromIpt(const int ipt)
Definition: RPCConst.cc:29
L1TMuonOverlapParams::etaLUT
const l1t::LUT * etaLUT() const
Definition: L1TMuonOverlapParams.h:203
OMTFProcessor::averagePatterns
void averagePatterns(int charge)
Definition: OMTFProcessor.cc:114
OMTFProcessor::restrictInput
OMTFinput::vector1D restrictInput(unsigned int iProcessor, unsigned int iCone, unsigned int iLayer, const OMTFinput::vector1D &layerHits)
Definition: OMTFProcessor.cc:270
OMTFResult
Definition: OMTFResult.h:9