CMS 3D CMS Logo
MuonPathAnalyzerPerSL.cc
Go to the documentation of this file.
1 #include "L1Trigger/DTTriggerPhase2/interface/MuonPathAnalyzerPerSL.h"
2 #include <cmath>
3 #include <memory>
4 
5 using namespace edm;
6 using namespace std;
7 using namespace cmsdt;
8 // ============================================================================
9 // Constructors and destructor
10 // ============================================================================
11 MuonPathAnalyzerPerSL::MuonPathAnalyzerPerSL(const ParameterSet &pset, edm::ConsumesCollector &iC)
12  : MuonPathAnalyzer(pset, iC),
13  bxTolerance_(30),
14  minQuality_(LOWQGHOST),
15  chiSquareThreshold_(50),
16  debug_(pset.getUntrackedParameter<bool>("debug")),
17  chi2Th_(pset.getUntrackedParameter<double>("chi2Th")),
18  tanPhiTh_(pset.getUntrackedParameter<double>("tanPhiTh")),
19  use_LSB_(pset.getUntrackedParameter<bool>("use_LSB")),
20  tanPsi_precision_(pset.getUntrackedParameter<double>("tanPsi_precision")),
21  x_precision_(pset.getUntrackedParameter<double>("x_precision")) {
22  if (debug_)
23  LogDebug("MuonPathAnalyzerPerSL") << "MuonPathAnalyzer: constructor";
24 
25  setChiSquareThreshold(chi2Th_ * 100.);
26 
27  //shift
28  int rawId;
29  shift_filename_ = pset.getParameter<edm::FileInPath>("shift_filename");
30  std::ifstream ifin3(shift_filename_.fullPath());
31  double shift;
32  if (ifin3.fail()) {
33  throw cms::Exception("Missing Input File")
34  << "MuonPathAnalyzerPerSL::MuonPathAnalyzerPerSL() - Cannot find " << shift_filename_.fullPath();
35  }
36  while (ifin3.good()) {
37  ifin3 >> rawId >> shift;
38  shiftinfo_[rawId] = shift;
39  }
40 
41  chosen_sl_ = pset.getUntrackedParameter<int>("trigger_with_sl");
42 
43  if (chosen_sl_ != 1 && chosen_sl_ != 3 && chosen_sl_ != 4) {
44  LogDebug("MuonPathAnalyzerPerSL") << "chosen sl must be 1,3 or 4(both superlayers)";
45  assert(chosen_sl_ != 1 && chosen_sl_ != 3 && chosen_sl_ != 4); //4 means run using the two superlayers
46  }
47 
48  dtGeomH = iC.esConsumes<DTGeometry, MuonGeometryRecord, edm::Transition::BeginRun>();
49 }
50 
51 MuonPathAnalyzerPerSL::~MuonPathAnalyzerPerSL() {
52  if (debug_)
53  LogDebug("MuonPathAnalyzerPerSL") << "MuonPathAnalyzer: destructor";
54 }
55 
56 // ============================================================================
57 // Main methods (initialise, run, finish)
58 // ============================================================================
59 void MuonPathAnalyzerPerSL::initialise(const edm::EventSetup &iEventSetup) {
60  if (debug_)
61  LogDebug("MuonPathAnalyzerPerSL") << "MuonPathAnalyzerPerSL::initialiase";
62 
63  const MuonGeometryRecord &geom = iEventSetup.get<MuonGeometryRecord>();
64  dtGeo_ = &geom.get(dtGeomH);
65 }
66 
67 void MuonPathAnalyzerPerSL::run(edm::Event &iEvent,
68  const edm::EventSetup &iEventSetup,
69  MuonPathPtrs &muonpaths,
70  std::vector<metaPrimitive> &metaPrimitives) {
71  if (debug_)
72  LogDebug("MuonPathAnalyzerPerSL") << "MuonPathAnalyzerPerSL: run";
73 
74  // fit per SL (need to allow for multiple outputs for a single mpath)
75  for (auto &muonpath : muonpaths) {
76  analyze(muonpath, metaPrimitives);
77  }
78 }
79 
80 void MuonPathAnalyzerPerSL::finish() {
81  if (debug_)
82  LogDebug("MuonPathAnalyzerPerSL") << "MuonPathAnalyzer: finish";
83 };
84 
85 constexpr int MuonPathAnalyzerPerSL::LAYER_ARRANGEMENTS_[NUM_LAYERS][NUM_CELL_COMB] = {
86  {0, 1, 2},
87  {1, 2, 3}, // Consecutive groups
88  {0, 1, 3},
89  {0, 2, 3} // Non-consecutive groups
90 };
91 
92 //------------------------------------------------------------------
93 //--- Métodos privados
94 //------------------------------------------------------------------
95 
96 void MuonPathAnalyzerPerSL::analyze(MuonPathPtr &inMPath, std::vector<metaPrimitive> &metaPrimitives) {
97  if (debug_)
98  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t starts";
99 
100  // LOCATE MPATH
101  int selected_Id = 0;
102  if (inMPath->primitive(0)->tdcTimeStamp() != -1)
103  selected_Id = inMPath->primitive(0)->cameraId();
104  else if (inMPath->primitive(1)->tdcTimeStamp() != -1)
105  selected_Id = inMPath->primitive(1)->cameraId();
106  else if (inMPath->primitive(2)->tdcTimeStamp() != -1)
107  selected_Id = inMPath->primitive(2)->cameraId();
108  else if (inMPath->primitive(3)->tdcTimeStamp() != -1)
109  selected_Id = inMPath->primitive(3)->cameraId();
110 
111  DTLayerId thisLId(selected_Id);
112  if (debug_)
113  LogDebug("MuonPathAnalyzerPerSL") << "Building up MuonPathSLId from rawId in the Primitive";
114  DTSuperLayerId MuonPathSLId(thisLId.wheel(), thisLId.station(), thisLId.sector(), thisLId.superLayer());
115  if (debug_)
116  LogDebug("MuonPathAnalyzerPerSL") << "The MuonPathSLId is" << MuonPathSLId;
117 
118  if (debug_)
119  LogDebug("MuonPathAnalyzerPerSL")
120  << "DTp2:analyze \t\t\t\t In analyze function checking if inMPath->isAnalyzable() " << inMPath->isAnalyzable();
121 
122  if (chosen_sl_ < 4 && thisLId.superLayer() != chosen_sl_)
123  return; // avoid running when mpath not in chosen SL (for 1SL fitting)
124 
125  auto mPath = std::make_shared<MuonPath>(inMPath);
126 
127  if (mPath->isAnalyzable()) {
128  if (debug_)
129  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t yes it is analyzable " << mPath->isAnalyzable();
130  setCellLayout(mPath->cellLayout());
131  evaluatePathQuality(mPath);
132  } else {
133  if (debug_)
134  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t no it is NOT analyzable " << mPath->isAnalyzable();
135  return;
136  }
137 
138  int wi[8], tdc[8], lat[8];
139  DTPrimitivePtr Prim0(mPath->primitive(0));
140  wi[0] = Prim0->channelId();
141  tdc[0] = Prim0->tdcTimeStamp();
142  DTPrimitivePtr Prim1(mPath->primitive(1));
143  wi[1] = Prim1->channelId();
144  tdc[1] = Prim1->tdcTimeStamp();
145  DTPrimitivePtr Prim2(mPath->primitive(2));
146  wi[2] = Prim2->channelId();
147  tdc[2] = Prim2->tdcTimeStamp();
148  DTPrimitivePtr Prim3(mPath->primitive(3));
149  wi[3] = Prim3->channelId();
150  tdc[3] = Prim3->tdcTimeStamp();
151  for (int i = 4; i < 8; i++) {
152  wi[i] = -1;
153  tdc[i] = -1;
154  lat[i] = -1;
155  }
156 
157  DTWireId wireId(MuonPathSLId, 2, 1);
158 
159  if (debug_)
160  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t checking if it passes the min quality cut "
161  << mPath->quality() << ">" << minQuality_;
162  if (mPath->quality() >= minQuality_) {
163  if (debug_)
164  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t min quality achievedCalidad: " << mPath->quality();
165  for (int i = 0; i <= 3; i++) {
166  if (debug_)
167  LogDebug("MuonPathAnalyzerPerSL")
168  << "DTp2:analyze \t\t\t\t Capa: " << mPath->primitive(i)->layerId()
169  << " Canal: " << mPath->primitive(i)->channelId() << " TDCTime: " << mPath->primitive(i)->tdcTimeStamp();
170  }
171  if (debug_)
172  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t Starting lateralities loop, totalNumValLateralities: "
173  << totalNumValLateralities_;
174 
175  double best_chi2 = 99999.;
176  double chi2_jm_tanPhi = 999;
177  double chi2_jm_x = -1;
178  double chi2_jm_t0 = -1;
179  double chi2_phi = -1;
180  double chi2_phiB = -1;
181  double chi2_chi2 = -1;
182  int chi2_quality = -1;
183  int bestLat[8];
184  for (int i = 0; i < 8; i++) {
185  bestLat[i] = -1;
186  }
187 
188  for (int i = 0; i < totalNumValLateralities_; i++) { //here
189  if (debug_)
190  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t laterality #- " << i;
191  if (debug_)
192  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t laterality #- " << i << " checking quality:";
193  if (debug_)
194  LogDebug("MuonPathAnalyzerPerSL")
195  << "DTp2:analyze \t\t\t\t\t laterality #- " << i << " checking mPath Quality=" << mPath->quality();
196  if (debug_)
197  LogDebug("MuonPathAnalyzerPerSL")
198  << "DTp2:analyze \t\t\t\t\t laterality #- " << i << " latQuality_[i].val=" << latQuality_[i].valid;
199  if (debug_)
200  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t laterality #- " << i << " before if:";
201 
202  if (latQuality_[i].valid and
203  (((mPath->quality() == HIGHQ or mPath->quality() == HIGHQGHOST) and latQuality_[i].quality == HIGHQ) or
204  ((mPath->quality() == LOWQ or mPath->quality() == LOWQGHOST) and latQuality_[i].quality == LOWQ))) {
205  if (debug_)
206  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t laterality #- " << i << " inside if";
207  mPath->setBxTimeValue(latQuality_[i].bxValue);
208  if (debug_)
209  LogDebug("MuonPathAnalyzerPerSL")
210  << "DTp2:analyze \t\t\t\t\t laterality #- " << i << " settingLateralCombination";
211  mPath->setLateralComb(lateralities_[i]);
212  if (debug_)
213  LogDebug("MuonPathAnalyzerPerSL")
214  << "DTp2:analyze \t\t\t\t\t laterality #- " << i << " done settingLateralCombination";
215 
216  // Clonamos el objeto analizado.
217  auto mpAux = std::make_shared<MuonPath>(mPath);
218  lat[0] = mpAux->lateralComb()[0];
219  lat[1] = mpAux->lateralComb()[1];
220  lat[2] = mpAux->lateralComb()[2];
221  lat[3] = mpAux->lateralComb()[3];
222 
223  int wiOk[NUM_LAYERS], tdcOk[NUM_LAYERS], latOk[NUM_LAYERS];
224  for (int lay = 0; lay < 4; lay++) {
225  if (latQuality_[i].invalidateHitIdx == lay) {
226  wiOk[lay] = -1;
227  tdcOk[lay] = -1;
228  latOk[lay] = -1;
229  } else {
230  wiOk[lay] = wi[lay];
231  tdcOk[lay] = tdc[lay];
232  latOk[lay] = lat[lay];
233  }
234  }
235 
236  int idxHitNotValid = latQuality_[i].invalidateHitIdx;
237  if (idxHitNotValid >= 0) {
238  auto dtpAux = std::make_shared<DTPrimitive>();
239  mpAux->setPrimitive(dtpAux, idxHitNotValid);
240  }
241 
242  if (debug_)
243  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t calculating parameters ";
244  calculatePathParameters(mpAux);
245  /*
246  * After calculating the parameters, and if it is a 4-hit fit,
247  * if the resultant chi2 is higher than the programmed threshold,
248  * the mpath is eliminated and we go to the next element
249  */
250  if ((mpAux->quality() == HIGHQ or mpAux->quality() == HIGHQGHOST) &&
251  mpAux->chiSquare() > chiSquareThreshold_) { //check this if!!!
252  if (debug_)
253  LogDebug("MuonPathAnalyzerPerSL")
254  << "DTp2:analyze \t\t\t\t\t HIGHQ or HIGHQGHOST but min chi2 or Q test not satisfied ";
255  } else {
256  if (debug_)
257  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t inside else, returning values: ";
258  if (debug_)
259  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t BX Time = " << mpAux->bxTimeValue();
260  if (debug_)
261  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t BX Id = " << mpAux->bxNumId();
262  if (debug_)
263  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t XCoor = " << mpAux->horizPos();
264  if (debug_)
265  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t tan(Phi)= " << mpAux->tanPhi();
266  if (debug_)
267  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t chi2= " << mpAux->chiSquare();
268  if (debug_)
269  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t lateralities = "
270  << " " << mpAux->lateralComb()[0] << " " << mpAux->lateralComb()[1] << " "
271  << mpAux->lateralComb()[2] << " " << mpAux->lateralComb()[3];
272 
273  DTChamberId ChId(MuonPathSLId.wheel(), MuonPathSLId.station(), MuonPathSLId.sector());
274 
275  double jm_tanPhi = -1. * mpAux->tanPhi(); //testing with this line
276  if (use_LSB_)
277  jm_tanPhi = floor(jm_tanPhi / tanPsi_precision_) * tanPsi_precision_;
278  double jm_x =
279  (((double)mpAux->horizPos()) / 10.) + x_precision_ * (round(shiftinfo_[wireId.rawId()] / x_precision_));
280  if (use_LSB_)
281  jm_x = ((double)round(((double)jm_x) / x_precision_)) * x_precision_;
282  //changing to chamber frame or reference:
283  double jm_t0 = mpAux->bxTimeValue();
284  int quality = mpAux->quality();
285 
286  //computing phi and phiB
287  double z = 0;
288  double z1 = Z_POS_SL;
289  double z3 = -1. * z1;
290  if (ChId.station() == 3 or ChId.station() == 4) {
291  z1 = z1 + Z_SHIFT_MB4;
292  z3 = z3 + Z_SHIFT_MB4;
293  } else if (MuonPathSLId.superLayer() == 1)
294  z = z1;
295  else if (MuonPathSLId.superLayer() == 3)
296  z = z3;
297 
298  GlobalPoint jm_x_cmssw_global = dtGeo_->chamber(ChId)->toGlobal(LocalPoint(jm_x, 0., z));
299  int thisec = MuonPathSLId.sector();
300  if (thisec == 13)
301  thisec = 4;
302  if (thisec == 14)
303  thisec = 10;
304  double phi = jm_x_cmssw_global.phi() - PHI_CONV * (thisec - 1);
305  double psi = atan(jm_tanPhi);
306  double phiB = hasPosRF(MuonPathSLId.wheel(), MuonPathSLId.sector()) ? psi - phi : -psi - phi;
307  double chi2 = mpAux->chiSquare() * 0.01; //in cmssw we need cm, 1 cm^2 = 100 mm^2
308 
309  if (debug_)
310  LogDebug("MuonPathAnalyzerPerSL")
311  << "DTp2:analyze \t\t\t\t\t\t\t\t pushing back metaPrimitive at x=" << jm_x << " tanPhi:" << jm_tanPhi
312  << " t0:" << jm_t0;
313 
314  if (mpAux->quality() == HIGHQ or
315  mpAux->quality() == HIGHQGHOST) { //keep only the values with the best chi2 among lateralities
316  if ((chi2 < best_chi2) && (std::abs(jm_tanPhi) <= tanPhiTh_)) {
317  chi2_jm_tanPhi = jm_tanPhi;
318  chi2_jm_x = (mpAux->horizPos() / 10.) + shiftinfo_[wireId.rawId()];
319  chi2_jm_t0 = mpAux->bxTimeValue();
320  chi2_phi = phi;
321  chi2_phiB = phiB;
322  chi2_chi2 = chi2;
323  best_chi2 = chi2;
324  chi2_quality = mpAux->quality();
325  for (int i = 0; i < 4; i++) {
326  bestLat[i] = lat[i];
327  }
328  }
329  } else if (std::abs(jm_tanPhi) <=
330  tanPhiTh_) { //write the metaprimitive in case no HIGHQ or HIGHQGHOST and tanPhi range
331  if (debug_)
332  LogDebug("MuonPathAnalyzerPerSL")
333  << "DTp2:analyze \t\t\t\t\t\t\t\t pushing back metaprimitive no HIGHQ or HIGHQGHOST";
334  metaPrimitives.emplace_back(metaPrimitive({MuonPathSLId.rawId(),
335  jm_t0,
336  jm_x,
337  jm_tanPhi,
338  phi,
339  phiB,
340  chi2,
341  quality,
342  wiOk[0],
343  tdcOk[0],
344  latOk[0],
345  wiOk[1],
346  tdcOk[1],
347  latOk[1],
348  wiOk[2],
349  tdcOk[2],
350  latOk[2],
351  wiOk[3],
352  tdcOk[3],
353  latOk[3],
354  wi[4],
355  tdc[4],
356  lat[4],
357  wi[5],
358  tdc[5],
359  lat[5],
360  wi[6],
361  tdc[6],
362  lat[6],
363  wi[7],
364  tdc[7],
365  lat[7],
366  -1}));
367  if (debug_)
368  LogDebug("MuonPathAnalyzerPerSL")
369  << "DTp2:analyze \t\t\t\t\t\t\t\t done pushing back metaprimitive no HIGHQ or HIGHQGHOST";
370  }
371  }
372  } else {
373  if (debug_)
374  LogDebug("MuonPathAnalyzerPerSL")
375  << "DTp2:analyze \t\t\t\t\t\t\t\t latQuality_[i].valid and (((mPath->quality()==HIGHQ or "
376  "mPath->quality()==HIGHQGHOST) and latQuality_[i].quality==HIGHQ) or ((mPath->quality() "
377  "== LOWQ or mPath->quality()==LOWQGHOST) and latQuality_[i].quality==LOWQ)) not passed";
378  }
379  }
380  if (chi2_jm_tanPhi != 999 and std::abs(chi2_jm_tanPhi) < tanPhiTh_) { //
381  if (debug_)
382  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t\t\t\t pushing back best chi2 metaPrimitive";
383  metaPrimitives.emplace_back(metaPrimitive({MuonPathSLId.rawId(),
384  chi2_jm_t0,
385  chi2_jm_x,
386  chi2_jm_tanPhi,
387  chi2_phi,
388  chi2_phiB,
389  chi2_chi2,
390  chi2_quality,
391  wi[0],
392  tdc[0],
393  bestLat[0],
394  wi[1],
395  tdc[1],
396  bestLat[1],
397  wi[2],
398  tdc[2],
399  bestLat[2],
400  wi[3],
401  tdc[3],
402  bestLat[3],
403  wi[4],
404  tdc[4],
405  bestLat[4],
406  wi[5],
407  tdc[5],
408  bestLat[5],
409  wi[6],
410  tdc[6],
411  bestLat[6],
412  wi[7],
413  tdc[7],
414  bestLat[7],
415  -1}));
416  }
417  }
418  if (debug_)
419  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t finishes";
420 }
421 
422 void MuonPathAnalyzerPerSL::setCellLayout(const int layout[NUM_LAYERS]) {
423  memcpy(cellLayout_, layout, 4 * sizeof(int));
424 
425  buildLateralities();
426 }
427 
432 void MuonPathAnalyzerPerSL::buildLateralities(void) {
433  LATERAL_CASES(*validCase)[NUM_LAYERS], sideComb[NUM_LAYERS];
434 
435  totalNumValLateralities_ = 0;
436  /* We generate all the possible lateralities combination for a given group
437  of cells */
438  for (int lowLay = LEFT; lowLay <= RIGHT; lowLay++)
439  for (int midLowLay = LEFT; midLowLay <= RIGHT; midLowLay++)
440  for (int midHigLay = LEFT; midHigLay <= RIGHT; midHigLay++)
441  for (int higLay = LEFT; higLay <= RIGHT; higLay++) {
442  sideComb[0] = static_cast<LATERAL_CASES>(lowLay);
443  sideComb[1] = static_cast<LATERAL_CASES>(midLowLay);
444  sideComb[2] = static_cast<LATERAL_CASES>(midHigLay);
445  sideComb[3] = static_cast<LATERAL_CASES>(higLay);
446 
447  /* If a laterality combination is valid, we store it */
448  if (isStraightPath(sideComb)) {
449  validCase = lateralities_ + totalNumValLateralities_;
450  memcpy(validCase, sideComb, 4 * sizeof(LATERAL_CASES));
451 
452  latQuality_[totalNumValLateralities_].valid = false;
453  latQuality_[totalNumValLateralities_].bxValue = 0;
454  latQuality_[totalNumValLateralities_].quality = NOPATH;
455  latQuality_[totalNumValLateralities_].invalidateHitIdx = -1;
456 
457  totalNumValLateralities_++;
458  }
459  }
460 }
461 
465 bool MuonPathAnalyzerPerSL::isStraightPath(LATERAL_CASES sideComb[NUM_LAYERS]) {
466  return true; //trying with all lateralities to be confirmed
467 
468  int i, ajustedLayout[NUM_LAYERS], pairDiff[3], desfase[3];
469 
470  for (i = 0; i <= 3; i++)
471  ajustedLayout[i] = cellLayout_[i] + sideComb[i];
472  for (i = 0; i <= 2; i++)
473  pairDiff[i] = ajustedLayout[i + 1] - ajustedLayout[i];
474  for (i = 0; i <= 1; i++)
475  desfase[i] = abs(pairDiff[i + 1] - pairDiff[i]);
476  desfase[2] = abs(pairDiff[2] - pairDiff[0]);
477  bool resultado = (desfase[0] > 1 or desfase[1] > 1 or desfase[2] > 1);
478 
479  return (!resultado);
480 }
481 void MuonPathAnalyzerPerSL::evaluatePathQuality(MuonPathPtr &mPath) {
482  int totalHighQ = 0, totalLowQ = 0;
483 
484  if (debug_)
485  LogDebug("MuonPathAnalyzerPerSL")
486  << "DTp2:evaluatePathQuality \t\t\t\t\t En evaluatePathQuality Evaluando PathQ. Celda base: "
487  << mPath->baseChannelId();
488  if (debug_)
489  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:evaluatePathQuality \t\t\t\t\t Total lateralidades: "
490  << totalNumValLateralities_;
491 
492  mPath->setQuality(NOPATH);
493 
494  for (int latIdx = 0; latIdx < totalNumValLateralities_; latIdx++) {
495  if (debug_)
496  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:evaluatePathQuality \t\t\t\t\t Analizando combinacion de lateralidad: "
497  << lateralities_[latIdx][0] << " " << lateralities_[latIdx][1] << " "
498  << lateralities_[latIdx][2] << " " << lateralities_[latIdx][3];
499 
500  evaluateLateralQuality(latIdx, mPath, &(latQuality_[latIdx]));
501 
502  if (latQuality_[latIdx].quality == HIGHQ) {
503  totalHighQ++;
504  if (debug_)
505  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:evaluatePathQuality \t\t\t\t\t\t Lateralidad HIGHQ";
506  }
507  if (latQuality_[latIdx].quality == LOWQ) {
508  totalLowQ++;
509  if (debug_)
510  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:evaluatePathQuality \t\t\t\t\t\t Lateralidad LOWQ";
511  }
512  }
513  /*
514  * Quality stablishment
515  */
516  if (totalHighQ == 1) {
517  mPath->setQuality(HIGHQ);
518  } else if (totalHighQ > 1) {
519  mPath->setQuality(HIGHQGHOST);
520  } else if (totalLowQ == 1) {
521  mPath->setQuality(LOWQ);
522  } else if (totalLowQ > 1) {
523  mPath->setQuality(LOWQGHOST);
524  }
525 }
526 
527 void MuonPathAnalyzerPerSL::evaluateLateralQuality(int latIdx, MuonPathPtr &mPath, LATQ_TYPE *latQuality) {
528  int layerGroup[3];
529  LATERAL_CASES sideComb[3];
530  PARTIAL_LATQ_TYPE latQResult[NUM_LAYERS] = {{false, 0}, {false, 0}, {false, 0}, {false, 0}};
531 
532  // Default values.
533  latQuality->valid = false;
534  latQuality->bxValue = 0;
535  latQuality->quality = NOPATH;
536  latQuality->invalidateHitIdx = -1;
537 
538  /* If, for a given laterality combination, the two consecutive 3-layer combinations
539  were a valid track, we will have found a right high-quality track, hence
540  it will be unnecessary to check the remaining 2 combinations.
541  In order to mimic the FPGA behaviour, we build a code that analyzes the 4 combinations
542  with an additional logic to discriminate the final quality of the track
543  */
544  for (int i = 0; i <= 3; i++) {
545  memcpy(layerGroup, LAYER_ARRANGEMENTS_[i], 3 * sizeof(int));
546 
547  // Pick the laterality combination for each cell
548  for (int j = 0; j < 3; j++)
549  sideComb[j] = lateralities_[latIdx][layerGroup[j]];
550 
551  validate(sideComb, layerGroup, mPath, &(latQResult[i]));
552  }
553  /*
554  Impose the condition for a complete laterality combination, that all combinations
555  should give the same BX vale to give a consistent track.
556  */
557  if (!sameBXValue(latQResult)) {
558  if (debug_)
559  LogDebug("MuonPathAnalyzerPerSL")
560  << "DTp2:evaluateLateralQuality \t\t\t\t\t Lateralidad DESCARTADA. Tolerancia de BX excedida";
561  return;
562  }
563 
564  // two complementary valid tracks => full muon track.
565  if ((latQResult[0].latQValid && latQResult[1].latQValid) or (latQResult[0].latQValid && latQResult[2].latQValid) or
566  (latQResult[0].latQValid && latQResult[3].latQValid) or (latQResult[1].latQValid && latQResult[2].latQValid) or
567  (latQResult[1].latQValid && latQResult[3].latQValid) or (latQResult[2].latQValid && latQResult[3].latQValid)) {
568  latQuality->valid = true;
569 
570  if (debug_)
571  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t\t Valid BXs";
572  long int sumBX = 0, numValid = 0;
573  for (int i = 0; i <= 3; i++) {
574  if (debug_)
575  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:analyze \t\t\t\t\t\t "
576  << "[" << latQResult[i].bxValue << "," << latQResult[i].latQValid << "]";
577  if (latQResult[i].latQValid) {
578  sumBX += latQResult[i].bxValue;
579  numValid++;
580  }
581  }
582 
583  // mean time of all lateralities.
584  if (numValid == 1)
585  latQuality->bxValue = sumBX;
586  else if (numValid == 2)
587  latQuality->bxValue = (sumBX * (MEANTIME_2LAT)) / std::pow(2, 15);
588  else if (numValid == 3)
589  latQuality->bxValue = (sumBX * (MEANTIME_3LAT)) / std::pow(2, 15);
590  else if (numValid == 4)
591  latQuality->bxValue = (sumBX * (MEANTIME_4LAT)) / std::pow(2, 15);
592 
593  latQuality->quality = HIGHQ;
594 
595  if (debug_)
596  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:evaluateLateralQuality \t\t\t\t\t Lateralidad ACEPTADA. HIGHQ.";
597  } else {
598  if (latQResult[0].latQValid or latQResult[1].latQValid or latQResult[2].latQValid or latQResult[3].latQValid) {
599  latQuality->valid = true;
600  latQuality->quality = LOWQ;
601  for (int i = 0; i < 4; i++)
602  if (latQResult[i].latQValid) {
603  latQuality->bxValue = latQResult[i].bxValue;
604  latQuality->invalidateHitIdx = omittedHit(i);
605  break;
606  }
607 
608  if (debug_)
609  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:evaluateLateralQuality \t\t\t\t\t Lateralidad ACEPTADA. LOWQ.";
610  } else {
611  if (debug_)
612  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:evaluateLateralQuality \t\t\t\t\t Lateralidad DESCARTADA. NOPATH.";
613  }
614  }
615 }
616 
621 void MuonPathAnalyzerPerSL::validate(LATERAL_CASES sideComb[3],
622  int layerIndex[3],
623  MuonPathPtr &mPath,
624  PARTIAL_LATQ_TYPE *latq) {
625  // Valor por defecto.
626  latq->bxValue = 0;
627  latq->latQValid = false;
628 
629  if (debug_)
630  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:validate \t\t\t\t\t\t\t In validate, checking muon path for layers: "
631  << layerIndex[0] << "/" << layerIndex[1] << "/" << layerIndex[2];
632 
633  if (debug_)
634  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:validate \t\t\t\t\t\t\t Partial lateralities: " << sideComb[0] << "/"
635  << sideComb[1] << "/" << sideComb[2];
636 
637  int validCells = 0;
638  for (int j = 0; j < 3; j++)
639  if (mPath->primitive(layerIndex[j])->isValidTime())
640  validCells++;
641 
642  if (validCells != 3) {
643  if (debug_)
644  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:validate \t\t\t\t\t\t\t There is no valid cells.";
645  return;
646  }
647 
648  if (debug_)
649  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:validate \t\t\t\t\t\t\t TDC values: "
650  << mPath->primitive(layerIndex[0])->tdcTimeStamp() << "/"
651  << mPath->primitive(layerIndex[1])->tdcTimeStamp() << "/"
652  << mPath->primitive(layerIndex[2])->tdcTimeStamp() << ".";
653 
654  if (debug_)
655  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:validate \t\t\t\t\t\t\t Valid TIMES: "
656  << mPath->primitive(layerIndex[0])->isValidTime() << "/"
657  << mPath->primitive(layerIndex[1])->isValidTime() << "/"
658  << mPath->primitive(layerIndex[2])->isValidTime() << ".";
659 
660  /* Vertical distances */
661  int dVertMI = layerIndex[1] - layerIndex[0];
662  int dVertSM = layerIndex[2] - layerIndex[1];
663 
664  /* Horizontal distances between lower/middle and middle/upper cells */
665  int dHorzMI = cellLayout_[layerIndex[1]] - cellLayout_[layerIndex[0]];
666  int dHorzSM = cellLayout_[layerIndex[2]] - cellLayout_[layerIndex[1]];
667 
668  /* Pair index of layers that we are using
669  SM => Upper + Middle
670  MI => Middle + Lower
671  We use pointers to simplify the code */
672  int *layPairSM = &layerIndex[1];
673  int *layPairMI = &layerIndex[0];
674 
675  /* Pair combination of cells to compose the equation. */
676  LATERAL_CASES smSides[2], miSides[2];
677 
678  /* Considering the index 0 of "sideComb" the laterality of the lower cells is stored,
679  we extract the laterality combiantion for SM and MI pairs */
680 
681  memcpy(smSides, &sideComb[1], 2 * sizeof(LATERAL_CASES));
682 
683  memcpy(miSides, &sideComb[0], 2 * sizeof(LATERAL_CASES));
684 
685  long int bxValue = 0;
686  int coefsAB[2] = {0, 0}, coefsCD[2] = {0, 0};
687  /* It's neccesary to be careful with that pointer's indirection. We need to
688  retrieve the lateral coeficientes (+-1) from the lower/middle and
689  middle/upper cell's lateral combinations. They are needed to evaluate the
690  existance of a possible BX value, following it's calculation equation */
691  lateralCoeficients(miSides, coefsAB);
692  lateralCoeficients(smSides, coefsCD);
693 
694  /* Each of the summs of the 'coefsCD' & 'coefsAB' give always as results 0, +-2
695  */
696 
697  int denominator = dVertMI * (coefsCD[1] + coefsCD[0]) - dVertSM * (coefsAB[1] + coefsAB[0]);
698 
699  if (denominator == 0) {
700  if (debug_)
701  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:validate \t\t\t\t\t\t\t Imposible to calculate BX.";
702  return;
703  }
704 
705  long int sumA = (long int)floor(MAXDRIFT * (dVertMI * dHorzSM - dVertSM * dHorzMI));
706  long int numerator =
707  (sumA + dVertMI * eqMainBXTerm(smSides, layPairSM, mPath) - dVertSM * eqMainBXTerm(miSides, layPairMI, mPath));
708 
709  // These magic numbers are for doing divisions in the FW.
710  // These divisions are done with a precision of 18bits.
711  if (denominator == -1 * DENOM_TYPE1)
712  bxValue = (numerator * (-1 * DIVISION_HELPER1)) / std::pow(2, NBITS);
713  else if (denominator == -1 * DENOM_TYPE2)
714  bxValue = (numerator * (-1 * DIVISION_HELPER2)) / std::pow(2, NBITS);
715  else if (denominator == -1 * DENOM_TYPE3)
716  bxValue = (numerator * (-1 * DIVISION_HELPER3)) / std::pow(2, NBITS);
717  else if (denominator == DENOM_TYPE3)
718  bxValue = (numerator * (DIVISION_HELPER3)) / std::pow(2, NBITS);
719  else if (denominator == DENOM_TYPE2)
720  bxValue = (numerator * (DIVISION_HELPER2)) / std::pow(2, NBITS);
721  else if (denominator == DENOM_TYPE1)
722  bxValue = (numerator * (DIVISION_HELPER1)) / std::pow(2, NBITS);
723  else
724  LogDebug("MuonPathAnalyzerPerSL") << "Different!";
725  if (bxValue < 0) {
726  if (debug_)
727  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:validate \t\t\t\t\t\t\t No-valid combination. Negative BX.";
728  return;
729  }
730 
731  for (int i = 0; i < 3; i++)
732  if (mPath->primitive(layerIndex[i])->isValidTime()) {
733  int diffTime = mPath->primitive(layerIndex[i])->tdcTimeStampNoOffset() - bxValue;
734 
735  if (diffTime <= 0 or diffTime > round(MAXDRIFT)) {
736  if (debug_)
737  LogDebug("MuonPathAnalyzerPerSL")
738  << "DTp2:validate \t\t\t\t\t\t\t Invalid BX value. at least one crazt TDC time";
739  return;
740  }
741  }
742  if (debug_)
743  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:validate \t\t\t\t\t\t\t BX: " << bxValue;
744 
745  /* If you reach here, the BX and partial laterality are considered are valid
746  */
747  latq->bxValue = bxValue;
748  latq->latQValid = true;
749 }
750 int MuonPathAnalyzerPerSL::eqMainBXTerm(LATERAL_CASES sideComb[2], int layerIdx[2], MuonPathPtr &mPath) {
751  int eqTerm = 0, coefs[2];
752 
753  lateralCoeficients(sideComb, coefs);
754 
755  eqTerm = coefs[0] * mPath->primitive(layerIdx[0])->tdcTimeStampNoOffset() +
756  coefs[1] * mPath->primitive(layerIdx[1])->tdcTimeStampNoOffset();
757 
758  if (debug_)
759  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:eqMainBXTerm \t\t\t\t\t In eqMainBXTerm EQTerm(BX): " << eqTerm;
760 
761  return (eqTerm);
762 }
763 int MuonPathAnalyzerPerSL::eqMainTerm(LATERAL_CASES sideComb[2], int layerIdx[2], MuonPathPtr &mPath, int bxValue) {
764  int eqTerm = 0, coefs[2];
765 
766  lateralCoeficients(sideComb, coefs);
767 
768  if (!use_LSB_)
769  eqTerm = coefs[0] * (mPath->primitive(layerIdx[0])->tdcTimeStampNoOffset() - bxValue) +
770  coefs[1] * (mPath->primitive(layerIdx[1])->tdcTimeStampNoOffset() - bxValue);
771  else
772  eqTerm = coefs[0] * floor((DRIFT_SPEED / (10 * x_precision_)) *
773  (mPath->primitive(layerIdx[0])->tdcTimeStampNoOffset() - bxValue)) +
774  coefs[1] * floor((DRIFT_SPEED / (10 * x_precision_)) *
775  (mPath->primitive(layerIdx[1])->tdcTimeStampNoOffset() - bxValue));
776 
777  if (debug_)
778  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:\t\t\t\t\t EQTerm(Main): " << eqTerm;
779 
780  return (eqTerm);
781 }
782 
783 void MuonPathAnalyzerPerSL::lateralCoeficients(LATERAL_CASES sideComb[2], int *coefs) {
784  if ((sideComb[0] == LEFT) && (sideComb[1] == LEFT)) {
785  *(coefs) = +1;
786  *(coefs + 1) = -1;
787  } else if ((sideComb[0] == LEFT) && (sideComb[1] == RIGHT)) {
788  *(coefs) = +1;
789  *(coefs + 1) = +1;
790  } else if ((sideComb[0] == RIGHT) && (sideComb[1] == LEFT)) {
791  *(coefs) = -1;
792  *(coefs + 1) = -1;
793  } else if ((sideComb[0] == RIGHT) && (sideComb[1] == RIGHT)) {
794  *(coefs) = -1;
795  *(coefs + 1) = +1;
796  }
797 }
798 
803 bool MuonPathAnalyzerPerSL::sameBXValue(PARTIAL_LATQ_TYPE *latq) {
804  bool result = true;
805 
806  if (debug_)
807  LogDebug("MuonPathAnalyzerPerSL") << "Dtp2:sameBXValue bxTolerance_: " << bxTolerance_;
808 
809  if (debug_)
810  LogDebug("MuonPathAnalyzerPerSL") << "Dtp2:sameBXValue \t\t\t\t\t\t d01:" << abs(latq[0].bxValue - latq[1].bxValue);
811  if (debug_)
812  LogDebug("MuonPathAnalyzerPerSL") << "Dtp2:sameBXValue \t\t\t\t\t\t d02:" << abs(latq[0].bxValue - latq[2].bxValue);
813  if (debug_)
814  LogDebug("MuonPathAnalyzerPerSL") << "Dtp2:sameBXValue \t\t\t\t\t\t d03:" << abs(latq[0].bxValue - latq[3].bxValue);
815  if (debug_)
816  LogDebug("MuonPathAnalyzerPerSL") << "Dtp2:sameBXValue \t\t\t\t\t\t d12:" << abs(latq[1].bxValue - latq[2].bxValue);
817  if (debug_)
818  LogDebug("MuonPathAnalyzerPerSL") << "Dtp2:sameBXValue \t\t\t\t\t\t d13:" << abs(latq[1].bxValue - latq[3].bxValue);
819  if (debug_)
820  LogDebug("MuonPathAnalyzerPerSL") << "Dtp2:sameBXValue \t\t\t\t\t\t d23:" << abs(latq[2].bxValue - latq[3].bxValue);
821 
822  bool d01, d02, d03, d12, d13, d23;
823  d01 = (abs(latq[0].bxValue - latq[1].bxValue) <= bxTolerance_) ? true : false;
824  d02 = (abs(latq[0].bxValue - latq[2].bxValue) <= bxTolerance_) ? true : false;
825  d03 = (abs(latq[0].bxValue - latq[3].bxValue) <= bxTolerance_) ? true : false;
826  d12 = (abs(latq[1].bxValue - latq[2].bxValue) <= bxTolerance_) ? true : false;
827  d13 = (abs(latq[1].bxValue - latq[3].bxValue) <= bxTolerance_) ? true : false;
828  d23 = (abs(latq[2].bxValue - latq[3].bxValue) <= bxTolerance_) ? true : false;
829 
830  /* 4 groups of partial combination of valid lateralities */
831  if ((latq[0].latQValid && latq[1].latQValid && latq[2].latQValid && latq[3].latQValid) && !(d01 && d12 && d23))
832  result = false;
833  else
834  /* 4 posible cases of 3 groups of valid partial lateralities */
835  if (((latq[0].latQValid && latq[1].latQValid && latq[2].latQValid) && !(d01 && d12)) or
836  ((latq[0].latQValid && latq[1].latQValid && latq[3].latQValid) && !(d01 && d13)) or
837  ((latq[0].latQValid && latq[2].latQValid && latq[3].latQValid) && !(d02 && d23)) or
838  ((latq[1].latQValid && latq[2].latQValid && latq[3].latQValid) && !(d12 && d23)))
839  result = false;
840  else
841  /* Lastly, the 4 possible cases of partial valid lateralities */
842 
843  if (((latq[0].latQValid && latq[1].latQValid) && !d01) or ((latq[0].latQValid && latq[2].latQValid) && !d02) or
844  ((latq[0].latQValid && latq[3].latQValid) && !d03) or ((latq[1].latQValid && latq[2].latQValid) && !d12) or
845  ((latq[1].latQValid && latq[3].latQValid) && !d13) or ((latq[2].latQValid && latq[3].latQValid) && !d23))
846  result = false;
847 
848  return result;
849 }
850 
852 void MuonPathAnalyzerPerSL::calculatePathParameters(MuonPathPtr &mPath) {
853  // The order is important.
854  if (debug_)
855  LogDebug("MuonPathAnalyzerPerSL")
856  << "DTp2:calculatePathParameters \t\t\t\t\t\t calculating calcCellDriftAndXcoor(mPath) ";
857  calcCellDriftAndXcoor(mPath);
858  if (debug_)
859  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:calculatePathParameters \t\t\t\t\t\t checking mPath->quality() "
860  << mPath->quality();
861  if (mPath->quality() == HIGHQ or mPath->quality() == HIGHQGHOST) {
862  if (debug_)
863  LogDebug("MuonPathAnalyzerPerSL")
864  << "DTp2:calculatePathParameters \t\t\t\t\t\t\t Quality test passed, now calcTanPhiXPosChamber4Hits(mPath) ";
865  calcTanPhiXPosChamber4Hits(mPath);
866  } else {
867  if (debug_)
868  LogDebug("MuonPathAnalyzerPerSL")
869  << "DTp2:calculatePathParameters \t\t\t\t\t\t\t Quality test NOT passed calcTanPhiXPosChamber3Hits(mPath) ";
870  calcTanPhiXPosChamber3Hits(mPath);
871  }
872 
873  if (debug_)
874  LogDebug("MuonPathAnalyzerPerSL") << "DTp2:calculatePathParameters \t\t\t\t\t\t calcChiSquare(mPath) ";
875  calcChiSquare(mPath);
876 }
877 
878 void MuonPathAnalyzerPerSL::calcTanPhiXPosChamber(MuonPathPtr &mPath) {
879  int layerIdx[2];
880  /*
881  To calculate path's angle are only necessary two valid primitives.
882  This method should be called only when a 'MuonPath' is determined as valid,
883  so, at least, three of its primitives must have a valid time.
884  With this two comparitions (which can be implemented easily as multiplexors
885  in the FPGA) this method ensures to catch two of those valid primitives to
886  evaluate the angle.
887 
888  The first one is below the middle line of the superlayer, while the other
889  one is above this line
890  */
891  if (mPath->primitive(0)->isValidTime())
892  layerIdx[0] = 0;
893  else
894  layerIdx[0] = 1;
895 
896  if (mPath->primitive(3)->isValidTime())
897  layerIdx[1] = 3;
898  else
899  layerIdx[1] = 2;
900 
901  /* We identify along which cells' sides the muon travels */
902  LATERAL_CASES sideComb[2];
903  sideComb[0] = (mPath->lateralComb())[layerIdx[0]];
904  sideComb[1] = (mPath->lateralComb())[layerIdx[1]];
905 
906  /* Horizontal gap between cells in cell's semi-length units */
907  int dHoriz = (mPath->cellLayout())[layerIdx[1]] - (mPath->cellLayout())[layerIdx[0]];
908 
909  /* Vertical gap between cells in cell's height units */
910  int dVert = layerIdx[1] - layerIdx[0];
911 
912  /*-----------------------------------------------------------------*/
913  /*--------------------- Phi angle calculation ---------------------*/
914  /*-----------------------------------------------------------------*/
915  float num = CELL_SEMILENGTH * dHoriz + DRIFT_SPEED * eqMainTerm(sideComb, layerIdx, mPath, mPath->bxTimeValue());
916 
917  float denom = CELL_HEIGHT * dVert;
918  float tanPhi = num / denom;
919 
920  mPath->setTanPhi(tanPhi);
921 
922  /*-----------------------------------------------------------------*/
923  /*----------------- Horizontal coord. calculation -----------------*/
924  /*-----------------------------------------------------------------*/
925 
926  /*
927  Using known coordinates, relative to superlayer axis reference, (left most
928  superlayer side, and middle line between 2nd and 3rd layers), calculating
929  horizontal coordinate implies using a basic line equation:
930  (y - y0) = (x - x0) * cotg(Phi)
931  This horizontal coordinate can be obtained setting y = 0 on last equation,
932  and also setting y0 and x0 with the values of a known muon's path cell
933  position hit.
934  It's enough to use the lower cell (layerIdx[0]) coordinates. So:
935  xC = x0 - y0 * tan(Phi)
936  */
937  float lowerXPHorizPos = mPath->xCoorCell(layerIdx[0]);
938 
939  float lowerXPVertPos = 0; // This is only the absolute value distance.
940  if (layerIdx[0] == 0)
941  lowerXPVertPos = CELL_HEIGHT + CELL_SEMIHEIGHT;
942  else
943  lowerXPVertPos = CELL_SEMIHEIGHT;
944 
945  mPath->setHorizPos(lowerXPHorizPos + lowerXPVertPos * tanPhi);
946 }
947 
951 void MuonPathAnalyzerPerSL::calcTanPhiXPosChamber4Hits(MuonPathPtr &mPath) {
952  int x_prec_inv = (int)(1. / (10. * x_precision_));
953  int numberOfBits = (int)(round(std::log(x_prec_inv) / std::log(2.)));
954  int numerator = 3 * (int)round(mPath->xCoorCell(3) / (10 * x_precision_)) +
955  (int)round(mPath->xCoorCell(2) / (10 * x_precision_)) -
956  (int)round(mPath->xCoorCell(1) / (10 * x_precision_)) -
957  3 * (int)round(mPath->xCoorCell(0) / (10 * x_precision_));
958  int CELL_HEIGHT_JM = pow(2, 15) / ((int)(10 * CELL_HEIGHT));
959  int tanPhi_x4096 = (numerator * CELL_HEIGHT_JM) >> (3 + numberOfBits);
960  mPath->setTanPhi(tanPhi_x4096 * tanPsi_precision_);
961 
962  float XPos = (mPath->xCoorCell(0) + mPath->xCoorCell(1) + mPath->xCoorCell(2) + mPath->xCoorCell(3)) / 4;
963  mPath->setHorizPos(floor(XPos / (10 * x_precision_)) * 10 * x_precision_);
964 }
965 
969 void MuonPathAnalyzerPerSL::calcTanPhiXPosChamber3Hits(MuonPathPtr &mPath) {
970  int layerIdx[2];
971  int x_prec_inv = (int)(1. / (10. * x_precision_));
972  int numberOfBits = (int)(round(std::log(x_prec_inv) / std::log(2.)));
973 
974  if (mPath->primitive(0)->isValidTime())
975  layerIdx[0] = 0;
976  else
977  layerIdx[0] = 1;
978 
979  if (mPath->primitive(3)->isValidTime())
980  layerIdx[1] = 3;
981  else
982  layerIdx[1] = 2;
983 
984  /*-----------------------------------------------------------------*/
985  /*--------------------- Phi angle calculation ---------------------*/
986  /*-----------------------------------------------------------------*/
987 
988  int tan_division_denominator_bits = 16;
989 
990  int num =
991  ((int)((int)(x_prec_inv * mPath->xCoorCell(layerIdx[1])) - (int)(x_prec_inv * mPath->xCoorCell(layerIdx[0])))
992  << (12 - numberOfBits));
993  int denominator = (layerIdx[1] - layerIdx[0]) * CELL_HEIGHT;
994  int denominator_inv = ((int)(0.5 + pow(2, tan_division_denominator_bits) / float(denominator)));
995 
996  float tanPhi = ((num * denominator_inv) >> tan_division_denominator_bits) / ((1. / tanPsi_precision_));
997 
998  mPath->setTanPhi(tanPhi);
999 
1000  /*-----------------------------------------------------------------*/
1001  /*----------------- Horizontal coord. calculation -----------------*/
1002  /*-----------------------------------------------------------------*/
1003  float XPos = 0;
1004  if (mPath->primitive(0)->isValidTime() and mPath->primitive(3)->isValidTime())
1005  XPos = (mPath->xCoorCell(0) + mPath->xCoorCell(3)) / 2;
1006  else
1007  XPos = (mPath->xCoorCell(1) + mPath->xCoorCell(2)) / 2;
1008 
1009  mPath->setHorizPos(floor(XPos / (10 * x_precision_)) * 10 * x_precision_);
1010 }
1011 
1015 void MuonPathAnalyzerPerSL::calcCellDriftAndXcoor(MuonPathPtr &mPath) {
1016  long int drift_speed_new = 889;
1017  long int drift_dist_um_x4;
1018  long int wireHorizPos_x4;
1019  long int pos_mm_x4;
1020  int x_prec_inv = (int)(1. / (10. * x_precision_));
1021 
1022  for (int i = 0; i <= 3; i++)
1023  if (mPath->primitive(i)->isValidTime()) {
1024  drift_dist_um_x4 =
1025  drift_speed_new * ((long int)mPath->primitive(i)->tdcTimeStampNoOffset() - (long int)mPath->bxTimeValue());
1026  wireHorizPos_x4 = (long)(mPath->primitive(i)->wireHorizPos() * x_prec_inv);
1027 
1028  if ((mPath->lateralComb())[i] == LEFT)
1029  pos_mm_x4 = wireHorizPos_x4 - (drift_dist_um_x4 >> 10);
1030  else
1031  pos_mm_x4 = wireHorizPos_x4 + (drift_dist_um_x4 >> 10);
1032 
1033  mPath->setXCoorCell(pos_mm_x4 * (10 * x_precision_), i);
1034  mPath->setDriftDistance(((float)(drift_dist_um_x4 >> 10)) * (10 * x_precision_), i);
1035  }
1036 }
1037 
1041 void MuonPathAnalyzerPerSL::calcChiSquare(MuonPathPtr &mPath) {
1042  int x_prec_inv = (int)(1. / (10. * x_precision_));
1043  int numberOfBits = (int)(round(std::log(x_prec_inv) / std::log(2.)));
1044  long int Z_FACTOR[NUM_LAYERS] = {-6, -2, 2, 6};
1045  for (int i = 0; i < 4; i++) {
1046  Z_FACTOR[i] = Z_FACTOR[i] * (long int)CELL_HEIGHT;
1047  }
1048  long int sum_A = 0, sum_B = 0;
1049  long int chi2_mm2_x1024 = 0;
1050  for (int i = 0; i < 4; i++) {
1051  if (mPath->primitive(i)->isValidTime()) {
1052  sum_A = (((int)(mPath->xCoorCell(i) / (10 * x_precision_))) - ((int)(mPath->horizPos() / (10 * x_precision_))))
1053  << (14 - numberOfBits);
1054  sum_B = Z_FACTOR[i] * ((int)(mPath->tanPhi() / tanPsi_precision_));
1055  chi2_mm2_x1024 += (sum_A - sum_B) * (sum_A - sum_B);
1056  }
1057  }
1058  chi2_mm2_x1024 = chi2_mm2_x1024 >> 18;
1059 
1060  mPath->setChiSquare(((double)chi2_mm2_x1024 / 1024.));
1061 }
1062 
1063 int MuonPathAnalyzerPerSL::omittedHit(int idx) {
1064  switch (idx) {
1065  case 0:
1066  return 3;
1067  case 1:
1068  return 0;
1069  case 2:
1070  return 2;
1071  case 3:
1072  return 1;
1073  }
1074 
1075  return -1;
1076 }
cmsdt::LATQ_TYPE
Definition: constants.h:157
MuonPathAnalyzerPerSL::hasPosRF
bool hasPosRF(int wh, int sec)
Definition: MuonPathAnalyzerPerSL.h:41
cmsdt::CELL_SEMIHEIGHT
constexpr float CELL_SEMIHEIGHT
Definition: constants.h:178
DTGeometry
Definition: DTGeometry.h:28
MuonPathAnalyzerPerSL::evaluateLateralQuality
void evaluateLateralQuality(int latIdx, MuonPathPtr &mPath, cmsdt::LATQ_TYPE *latQuality)
Definition: MuonPathAnalyzerPerSL.cc:527
cmsdt::metaPrimitive
Definition: constants.h:48
cmsdt::LATQ_TYPE::invalidateHitIdx
int invalidateHitIdx
Definition: constants.h:160
DTSuperLayerId
Definition: DTSuperLayerId.h:12
MuonPathAnalyzerPerSL::totalNumValLateralities_
int totalNumValLateralities_
Definition: MuonPathAnalyzerPerSL.h:88
MuonPathAnalyzerPerSL::analyze
void analyze(MuonPathPtr &inMPath, std::vector< cmsdt::metaPrimitive > &metaPrimitives)
Definition: MuonPathAnalyzerPerSL.cc:96
electrons_cff.bool
bool
Definition: electrons_cff.py:393
mps_fire.i
i
Definition: mps_fire.py:428
cmsdt::DENOM_TYPE3
constexpr int DENOM_TYPE3
Definition: constants.h:214
dqmMemoryStats.float
float
Definition: dqmMemoryStats.py:127
MuonPathAnalyzerPerSL::latQuality_
cmsdt::LATQ_TYPE latQuality_[cmsdt::NUM_LATERALITIES]
Definition: MuonPathAnalyzerPerSL.h:86
cmsdt::NOPATH
Definition: constants.h:41
makePileupJSON.denom
denom
Definition: makePileupJSON.py:147
MuonPathAnalyzerPerSL::lateralities_
cmsdt::LATERAL_CASES lateralities_[cmsdt::NUM_LATERALITIES][cmsdt::NUM_LAYERS]
Definition: MuonPathAnalyzerPerSL.h:85
cmsdt::DENOM_TYPE1
constexpr int DENOM_TYPE1
Definition: constants.h:212
MuonPathAnalyzerPerSL::MuonPathAnalyzerPerSL
MuonPathAnalyzerPerSL(const edm::ParameterSet &pset, edm::ConsumesCollector &iC)
Definition: MuonPathAnalyzerPerSL.cc:11
MuonPathAnalyzerPerSL::isStraightPath
bool isStraightPath(cmsdt::LATERAL_CASES sideComb[cmsdt::NUM_LAYERS])
Definition: MuonPathAnalyzerPerSL.cc:465
MuonPathAnalyzerPerSL::tanPsi_precision_
double tanPsi_precision_
Definition: MuonPathAnalyzerPerSL.h:99
cmsdt::DIVISION_HELPER2
constexpr int DIVISION_HELPER2
Definition: constants.h:210
edm
HLT enums.
Definition: AlignableModifier.h:19
MuonPathAnalyzerPerSL::eqMainTerm
int eqMainTerm(cmsdt::LATERAL_CASES sideComb[2], int layerIdx[2], MuonPathPtr &mPath, int bxValue)
Definition: MuonPathAnalyzerPerSL.cc:763
cmsdt::PARTIAL_LATQ_TYPE
Definition: constants.h:153
MuonPathAnalyzerPerSL::dtGeo_
const DTGeometry * dtGeo_
Definition: MuonPathAnalyzerPerSL.h:41
cmsdt::LEFT
Definition: constants.h:44
cmsdt::HIGHQ
Definition: constants.h:41
MuonPathAnalyzerPerSL::calcChiSquare
void calcChiSquare(MuonPathPtr &mPath)
Definition: MuonPathAnalyzerPerSL.cc:1041
MuonPathAnalyzerPerSL::calcTanPhiXPosChamber3Hits
void calcTanPhiXPosChamber3Hits(MuonPathPtr &mPath)
Definition: MuonPathAnalyzerPerSL.cc:969
cmsdt::PARTIAL_LATQ_TYPE::bxValue
int bxValue
Definition: constants.h:155
MuonPathAnalyzerPerSL::dtGeomH
edm::ESGetToken< DTGeometry, MuonGeometryRecord > dtGeomH
Definition: MuonPathAnalyzerPerSL.h:45
cms::cuda::assert
assert(be >=bs)
MuonPathAnalyzerPerSL::calcTanPhiXPosChamber4Hits
void calcTanPhiXPosChamber4Hits(MuonPathPtr &mPath)
Definition: MuonPathAnalyzerPerSL.cc:951
makeLayoutFileForGui.layout
layout
Definition: makeLayoutFileForGui.py:67
cmsdt::DENOM_TYPE2
constexpr int DENOM_TYPE2
Definition: constants.h:213
MuonPathAnalyzerPerSL::initialise
void initialise(const edm::EventSetup &iEventSetup) override
Definition: MuonPathAnalyzerPerSL.cc:59
DTSuperLayerId::superLayer
int superLayer() const
Return the superlayer number.
Definition: DTSuperLayerId.h:39
align::LocalPoint
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:30
MuonPathAnalyzerPerSL::calculatePathParameters
void calculatePathParameters(MuonPathPtr &mPath)
Definition: MuonPathAnalyzerPerSL.cc:852
hltPixelTracks_cff.chi2
chi2
Definition: hltPixelTracks_cff.py:25
edm::ConsumesCollector::esConsumes
auto esConsumes()
Definition: ConsumesCollector.h:97
MuonPathPtrs
std::vector< MuonPathPtr > MuonPathPtrs
Definition: MuonPath.h:122
DTGeometry::chamber
const DTChamber * chamber(const DTChamberId &id) const
Return a DTChamber given its id.
Definition: DTGeometry.cc:90
cmsdt::LOWQGHOST
Definition: constants.h:41
heavyIonCSV_trainingSettings.idx
idx
Definition: heavyIonCSV_trainingSettings.py:5
MuonPathAnalyzerPerSL::setChiSquareThreshold
void setChiSquareThreshold(float ch2Thr)
Definition: MuonPathAnalyzerPerSL.h:33
edm::FileInPath
Definition: FileInPath.h:64
MuonPathAnalyzerPerSL::chosen_sl_
int chosen_sl_
Definition: MuonPathAnalyzerPerSL.h:51
MuonPathAnalyzerPerSL::finish
void finish() override
Definition: MuonPathAnalyzerPerSL.cc:80
cmsdt
Definition: constants.h:39
cmsdt::LATERAL_CASES
LATERAL_CASES
Definition: constants.h:44
edm::EventSetup::get
T get() const
Definition: EventSetup.h:80
cmsdt::CELL_SEMILENGTH
constexpr int CELL_SEMILENGTH
Definition: constants.h:180
DTWireId
Definition: DTWireId.h:12
MuonPathAnalyzerPerSL::~MuonPathAnalyzerPerSL
~MuonPathAnalyzerPerSL() override
Definition: MuonPathAnalyzerPerSL.cc:51
DDAxes::z
MuonPathAnalyzerPerSL.h
cmsdt::Z_POS_SL
constexpr float Z_POS_SL
Definition: constants.h:244
relativeConstraints.geom
geom
Definition: relativeConstraints.py:72
cmsdt::DIVISION_HELPER3
constexpr int DIVISION_HELPER3
Definition: constants.h:211
Point3DBase< float, GlobalTag >
MuonPathAnalyzerPerSL::x_precision_
double x_precision_
Definition: MuonPathAnalyzerPerSL.h:100
DTLayerId
Definition: DTLayerId.h:12
MuonPathAnalyzerPerSL::evaluatePathQuality
void evaluatePathQuality(MuonPathPtr &mPath)
Definition: MuonPathAnalyzerPerSL.cc:481
cmsdt::PHI_CONV
constexpr double PHI_CONV
Definition: constants.h:240
MuonPathAnalyzerPerSL::use_LSB_
bool use_LSB_
Definition: MuonPathAnalyzerPerSL.h:98
GeomDet::toGlobal
GlobalPoint toGlobal(const Local2DPoint &lp) const
Conversion to the global R.F. from the R.F. of the GeomDet.
Definition: GeomDet.h:49
MuonPathAnalyzerPerSL::bxTolerance_
int bxTolerance_
Definition: MuonPathAnalyzerPerSL.h:90
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:223
cmsdt::LATQ_TYPE::quality
MP_QUALITY quality
Definition: constants.h:161
edm::ParameterSet
Definition: ParameterSet.h:47
MuonPathAnalyzerPerSL::buildLateralities
void buildLateralities(void)
Definition: MuonPathAnalyzerPerSL.cc:432
IntegrityClient_cfi.ChId
ChId
Definition: IntegrityClient_cfi.py:18
cmsdt::PARTIAL_LATQ_TYPE::latQValid
bool latQValid
Definition: constants.h:154
cmsdt::NBITS
constexpr int NBITS
Definition: constants.h:215
MuonPathAnalyzerPerSL::shiftinfo_
std::map< int, float > shiftinfo_
Definition: MuonPathAnalyzerPerSL.h:49
cmsdt::RIGHT
Definition: constants.h:44
cmsdt::HIGHQGHOST
Definition: constants.h:41
cmsdt::MEANTIME_2LAT
constexpr int MEANTIME_2LAT
Definition: constants.h:248
MuonPathAnalyzerPerSL::sameBXValue
bool sameBXValue(cmsdt::PARTIAL_LATQ_TYPE *latq)
Definition: MuonPathAnalyzerPerSL.cc:803
MuonPathAnalyzerPerSL::minQuality_
cmsdt::MP_QUALITY minQuality_
Definition: MuonPathAnalyzerPerSL.h:91
createfilelist.int
int
Definition: createfilelist.py:10
iEvent
int iEvent
Definition: GenABIO.cc:224
cmsdt::MEANTIME_3LAT
constexpr int MEANTIME_3LAT
Definition: constants.h:249
MuonPathAnalyzerPerSL::calcCellDriftAndXcoor
void calcCellDriftAndXcoor(MuonPathPtr &mPath)
Definition: MuonPathAnalyzerPerSL.cc:1015
MuonPathAnalyzerPerSL::cellLayout_
int cellLayout_[cmsdt::NUM_LAYERS]
Definition: MuonPathAnalyzerPerSL.h:97
MuonPathAnalyzerPerSL::eqMainBXTerm
int eqMainBXTerm(cmsdt::LATERAL_CASES sideComb[2], int layerIdx[2], MuonPathPtr &mPath)
Definition: MuonPathAnalyzerPerSL.cc:750
HLTTauDQMOffline_cfi.numerator
numerator
Definition: HLTTauDQMOffline_cfi.py:179
NUM_LAYERS
Definition: CSCALCTDigi.cc:14
edm::EventSetup
Definition: EventSetup.h:57
MuonPathAnalyzerPerSL::chiSquareThreshold_
float chiSquareThreshold_
Definition: MuonPathAnalyzerPerSL.h:92
MuonPathAnalyzer
Definition: MuonPathAnalyzer.h:37
cmsdt::MEANTIME_4LAT
constexpr int MEANTIME_4LAT
Definition: constants.h:250
MuonPathPtr
std::shared_ptr< MuonPath > MuonPathPtr
Definition: MuonPath.h:121
MuonPathAnalyzerPerSL::chi2Th_
double chi2Th_
Definition: MuonPathAnalyzerPerSL.h:94
cmsdt::NUM_CELL_COMB
constexpr int NUM_CELL_COMB
Definition: constants.h:202
cmsdt::CELL_HEIGHT
constexpr int CELL_HEIGHT
Definition: constants.h:177
cmsdt::Z_SHIFT_MB4
constexpr float Z_SHIFT_MB4
Definition: constants.h:243
EgammaValidation_cff.num
num
Definition: EgammaValidation_cff.py:34
MuonPathAnalyzerPerSL::lateralCoeficients
void lateralCoeficients(cmsdt::LATERAL_CASES sideComb[2], int *coefs)
Definition: MuonPathAnalyzerPerSL.cc:783
DDAxes::phi
cmsdt::DRIFT_SPEED
constexpr float DRIFT_SPEED
Definition: constants.h:182
HLTTauDQMOffline_cfi.denominator
denominator
Definition: HLTTauDQMOffline_cfi.py:180
psi
std::map< std::string, int, std::less< std::string > > psi
Definition: CountProcessesAction.h:15
MuonPathAnalyzerPerSL::setCellLayout
void setCellLayout(const int layout[cmsdt::NUM_LAYERS])
Definition: MuonPathAnalyzerPerSL.cc:422
DTChamberId::sector
int sector() const
Definition: DTChamberId.h:49
cmsdt::LOWQ
Definition: constants.h:41
std
Definition: JetResolutionObject.h:76
DetId::rawId
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
DTPrimitivePtr
std::shared_ptr< DTPrimitive > DTPrimitivePtr
Definition: DTprimitive.h:54
qcdUeDQM_cfi.quality
quality
Definition: qcdUeDQM_cfi.py:31
edm::shift
static unsigned const int shift
Definition: LuminosityBlockID.cc:7
MuonPathAnalyzerPerSL::run
void run(edm::Event &iEvent, const edm::EventSetup &iEventSetup, MuonPathPtrs &inMpath, std::vector< cmsdt::metaPrimitive > &metaPrimitives) override
Definition: MuonPathAnalyzerPerSL.cc:67
edm::Transition::BeginRun
Exception
Definition: hltDiff.cc:246
MuonPathAnalyzerPerSL::LAYER_ARRANGEMENTS_
static const int LAYER_ARRANGEMENTS_[cmsdt::NUM_LAYERS][cmsdt::NUM_CELL_COMB]
Definition: MuonPathAnalyzerPerSL.h:84
cmsdt::LATQ_TYPE::valid
bool valid
Definition: constants.h:158
cmsdt::DIVISION_HELPER1
constexpr int DIVISION_HELPER1
Definition: constants.h:209
or
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
MuonPathAnalyzerPerSL::tanPhiTh_
double tanPhiTh_
Definition: MuonPathAnalyzerPerSL.h:96
MuonPathAnalyzerPerSL::debug_
bool debug_
Definition: MuonPathAnalyzerPerSL.h:93
dqm-mbProfile.log
log
Definition: dqm-mbProfile.py:17
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29
mps_fire.result
result
Definition: mps_fire.py:311
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
RunInfoPI::valid
Definition: RunInfoPayloadInspectoHelper.h:16
DTChamberId
Definition: DTChamberId.h:14
MuonPathAnalyzerPerSL::shift_filename_
edm::FileInPath shift_filename_
Definition: MuonPathAnalyzerPerSL.h:48
dqmiolumiharvest.j
j
Definition: dqmiolumiharvest.py:66
cmsdt::LATQ_TYPE::bxValue
int bxValue
Definition: constants.h:159
edm::Event
Definition: Event.h:73
MuonGeometryRecord
Definition: MuonGeometryRecord.h:34
DTChamberId::wheel
int wheel() const
Return the wheel number.
Definition: DTChamberId.h:39
edm::ConsumesCollector
Definition: ConsumesCollector.h:45
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
MuonPathAnalyzerPerSL::omittedHit
int omittedHit(int idx)
Definition: MuonPathAnalyzerPerSL.cc:1063
edm::FileInPath::fullPath
std::string fullPath() const
Definition: FileInPath.cc:163
MuonPathAnalyzerPerSL::validate
void validate(cmsdt::LATERAL_CASES sideComb[3], int layerIndex[3], MuonPathPtr &mPath, cmsdt::PARTIAL_LATQ_TYPE *latq)
Definition: MuonPathAnalyzerPerSL.cc:621
DTChamberId::station
int station() const
Return the station number.
Definition: DTChamberId.h:42
MuonPathAnalyzerPerSL::calcTanPhiXPosChamber
void calcTanPhiXPosChamber(MuonPathPtr &mPath)
Definition: MuonPathAnalyzerPerSL.cc:878
muonDTDigis_cfi.pset
pset
Definition: muonDTDigis_cfi.py:27
cmsdt::MAXDRIFT
constexpr float MAXDRIFT
Definition: constants.h:175
Generated for CMSSW Reference Manual by doxygen 1.8.16