CMS 3D CMS Logo

AlignmentPayloadInspectorHelper.h
Go to the documentation of this file.
1 #ifndef CONDCORE_ALIGNMENTPLUGINS_ALIGNMENTPAYLOADINSPECTORHELPER_H
2 #define CONDCORE_ALIGNMENTPLUGINS_ALIGNMENTPAYLOADINSPECTORHELPER_H
3 
4 #include <iostream>
5 #include <algorithm>
6 #include <vector>
7 #include <numeric>
8 #include <string>
9 #include "TH1.h"
10 #include "TCanvas.h"
11 #include "TPaveStats.h"
12 #include "TStyle.h"
13 #include "TList.h"
17 #include "DataFormats/Math/interface/deltaPhi.h" // for deltaPhi
18 #include "DataFormats/Math/interface/Rounding.h" // for rounding
21 
22 //#define MMDEBUG // uncomment for debugging at compile time
23 #ifdef MMDEBUG
24 #include <iostream>
25 #define COUT std::cout << "MM "
26 #else
27 #define COUT edm::LogVerbatim("")
28 #endif
29 
30 namespace AlignmentPI {
31 
32  // size of the phase-I Tracker APE payload (including both SS + DS modules)
33  static const unsigned int phase0size = 19876;
34  static const float cmToUm = 10000.f;
35  static const float tomRad = 1000.f;
36 
37  // method to zero all elements whose difference from 2Pi
38  // is less than the tolerance (2*10e-7)
39  inline double returnZeroIfNear2PI(const double phi) {
40  const double tol = 2.e-7; // default tolerance 1.e-7 doesn't account for possible variations
41  if (cms_rounding::roundIfNear0(std::abs(phi) - 2 * M_PI, tol) == 0.f) {
42  return 0.f;
43  } else {
44  return phi;
45  }
46  }
47 
48  // method to bring back around 0 all elements whose difference
49  // frm 2Pi is is less than the tolerance (in micro-rad)
50  inline double trim2PIs(const double phi, const double tolerance = 1.f) {
51  if (std::abs((std::abs(phi) - 2 * M_PI) * tomRad) < tolerance) {
52  return (std::abs(phi) - 2 * M_PI);
53  } else {
54  return phi;
55  }
56  }
57 
58  enum coordinate {
59  t_x = 1,
60  t_y = 2,
61  t_z = 3,
62  rot_alpha = 4,
63  rot_beta = 5,
64  rot_gamma = 6,
65  };
66 
67  // M.M. 2017/09/12
68  // As the matrix is symmetric, we map only 6/9 terms
69  // More terms for the extended APE can be added to the following methods
70 
71  enum index { XX = 1, XY = 2, XZ = 3, YZ = 4, YY = 5, ZZ = 6 };
72 
73  enum partitions { INVALID = 0, BPix = 1, FPix = 2, TIB = 3, TID = 4, TOB = 5, TEC = 6 };
74 
75  enum class PARTITION {
76  BPIX, // 0 Barrel Pixel
77  FPIXp, // 1 Forward Pixel Plus
78  FPIXm, // 2 Forward Pixel Minus
79  TIB, // 3 Tracker Inner Barrel
80  TIDp, // 4 Tracker Inner Disks Plus
81  TIDm, // 5 Tracker Inner Disks Minus
82  TOB, // 6 Tracker Outer Barrel
83  TECp, // 7 Tracker Endcaps Plus
84  TECm, // 8 Tracker Endcaps Minus
85  LAST = TECm
86  };
87 
88  extern const PARTITION PARTITIONS[(int)PARTITION::LAST + 1];
98 
99  inline std::ostream& operator<<(std::ostream& o, PARTITION x) {
101  }
102 
103  enum regions {
104  BPixL1o, //0 Barrel Pixel Layer 1 outer
105  BPixL1i, //1 Barrel Pixel Layer 1 inner
106  BPixL2o, //2 Barrel Pixel Layer 2 outer
107  BPixL2i, //3 Barrel Pixel Layer 2 inner
108  BPixL3o, //4 Barrel Pixel Layer 3 outer
109  BPixL3i, //5 Barrel Pixel Layer 3 inner
110  BPixL4o, //6 Barrel Pixel Layer 4 outer
111  BPixL4i, //7 Barrel Pixel Layer 4 inner
112  FPixmL1, //8 Forward Pixel Minus side Disk 1
113  FPixmL2, //9 Forward Pixel Minus side Disk 2
114  FPixmL3, //10 Forward Pixel Minus side Disk 3
115  FPixpL1, //11 Forward Pixel Plus side Disk 1
116  FPixpL2, //12 Forward Pixel Plus side Disk 2
117  FPixpL3, //13 Forward Pixel Plus side Disk 3
118  TIBL1Ro, //14 Inner Barrel Layer 1 Rphi outer
119  TIBL1Ri, //15 Inner Barrel Layer 1 Rphi inner
120  TIBL1So, //16 Inner Barrel Layer 1 Stereo outer
121  TIBL1Si, //17 Inner Barrel Layer 1 Stereo inner
122  TIBL2Ro, //18 Inner Barrel Layer 2 Rphi outer
123  TIBL2Ri, //19 Inner Barrel Layer 2 Rphi inner
124  TIBL2So, //20 Inner Barrel Layer 2 Stereo outer
125  TIBL2Si, //21 Inner Barrel Layer 2 Stereo inner
126  TIBL3o, //22 Inner Barrel Layer 3 outer
127  TIBL3i, //23 Inner Barrel Layer 3 inner
128  TIBL4o, //24 Inner Barrel Layer 4 outer
129  TIBL4i, //25 Inner Barrel Layer 4 inner
130  TOBL1Ro, //26 Outer Barrel Layer 1 Rphi outer
131  TOBL1Ri, //27 Outer Barrel Layer 1 Rphi inner
132  TOBL1So, //28 Outer Barrel Layer 1 Stereo outer
133  TOBL1Si, //29 Outer Barrel Layer 1 Stereo inner
134  TOBL2Ro, //30 Outer Barrel Layer 2 Rphi outer
135  TOBL2Ri, //31 Outer Barrel Layer 2 Rphi inner
136  TOBL2So, //32 Outer Barrel Layer 2 Stereo outer
137  TOBL2Si, //33 Outer Barrel Layer 2 Stereo inner
138  TOBL3o, //34 Outer Barrel Layer 3 outer
139  TOBL3i, //35 Outer Barrel Layer 3 inner
140  TOBL4o, //36 Outer Barrel Layer 4 outer
141  TOBL4i, //37 Outer Barrel Layer 4 inner
142  TOBL5o, //38 Outer Barrel Layer 5 outer
143  TOBL5i, //39 Outer Barrel Layer 5 inner
144  TOBL6o, //40 Outer Barrel Layer 6 outer
145  TOBL6i, //41 Outer Barrel Layer 6 inner
146  TIDmR1R, //42 Inner Disk Minus side Ring 1 Rphi
147  TIDmR1S, //43 Inner Disk Minus side Ring 1 Stereo
148  TIDmR2R, //44 Inner Disk Minus side Ring 2 Rphi
149  TIDmR2S, //45 Inner Disk Minus side Ring 2 Stereo
150  TIDmR3, //46 Inner Disk Minus side Ring 3
151  TIDpR1R, //47 Inner Disk Plus side Ring 1 Rphi
152  TIDpR1S, //48 Inner Disk Plus side Ring 1 Stereo
153  TIDpR2R, //49 Inner Disk Plus side Ring 2 Rphi
154  TIDpR2S, //50 Inner Disk Plus side Ring 2 Stereo
155  TIDpR3, //51 Inner Disk Plus side Ring 3
156  TECmR1R, //52 Endcaps Minus side Ring 1 Rphi
157  TECmR1S, //53 Endcaps Minus side Ring 1 Stereo
158  TECmR2R, //54 Encdaps Minus side Ring 2 Rphi
159  TECmR2S, //55 Endcaps Minus side Ring 2 Stereo
160  TECmR3, //56 Endcaps Minus side Ring 3
161  TECmR4, //57 Endcaps Minus side Ring 4
162  TECmR5, //58 Endcaps Minus side Ring 5
163  TECmR6, //59 Endcaps Minus side Ring 6
164  TECmR7, //60 Endcaps Minus side Ring 7
165  TECpR1R, //61 Endcaps Plus side Ring 1 Rphi
166  TECpR1S, //62 Endcaps Plus side Ring 1 Stereo
167  TECpR2R, //63 Encdaps Plus side Ring 2 Rphi
168  TECpR2S, //64 Endcaps Plus side Ring 2 Stereo
169  TECpR3, //65 Endcaps Plus side Ring 3
170  TECpR4, //66 Endcaps Plus side Ring 4
171  TECpR5, //67 Endcaps Plus side Ring 5
172  TECpR6, //68 Endcaps Plus side Ring 6
173  TECpR7, //67 Endcaps Plus side Ring 7
174  StripDoubleSide, // 70 -- not to be considered
175  NUM_OF_REGIONS // 71 -- default
176  };
177 
178  /*--------------------------------------------------------------------*/
180  /*--------------------------------------------------------------------*/
181  {
182  switch (e) {
184  return "BPixL1o";
186  return "BPixL1i";
188  return "BPixL2o";
190  return "BPixL2i";
192  return "BPixL3o";
194  return "BPixL3i";
196  return "BPixL4o";
198  return "BPixL4i";
200  return "FPixmL1";
202  return "FPixmL2";
204  return "FPixmL3";
206  return "FPixpL1";
208  return "FPixpL2";
210  return "FPixpL3";
212  return "TIBL1Ro";
214  return "TIBL1Ri";
216  return "TIBL1So";
218  return "TIBL1Si";
220  return "TIBL2Ro";
222  return "TIBL2Ri";
224  return "TIBL2So";
226  return "TIBL2Si";
227  case AlignmentPI::TIBL3o:
228  return "TIBL3o";
229  case AlignmentPI::TIBL3i:
230  return "TIBL3i";
231  case AlignmentPI::TIBL4o:
232  return "TIBL4o";
233  case AlignmentPI::TIBL4i:
234  return "TIBL4i";
236  return "TOBL1Ro";
238  return "TOBL1Ri";
240  return "TOBL1So";
242  return "TOBL1Si";
244  return "TOBL2Ro";
246  return "TOBL2Ri";
248  return "TOBL2So";
250  return "TOBL2Si";
251  case AlignmentPI::TOBL3o:
252  return "TOBL3o";
253  case AlignmentPI::TOBL3i:
254  return "TOBL3i";
255  case AlignmentPI::TOBL4o:
256  return "TOBL4o";
257  case AlignmentPI::TOBL4i:
258  return "TOBL4i";
259  case AlignmentPI::TOBL5o:
260  return "TOBL5o";
261  case AlignmentPI::TOBL5i:
262  return "TOBL5i";
263  case AlignmentPI::TOBL6o:
264  return "TOBL6o";
265  case AlignmentPI::TOBL6i:
266  return "TOBL6i";
268  return "TIDmR1R";
270  return "TIDmR1S";
272  return "TIDmR2R";
274  return "TIDmR2S";
275  case AlignmentPI::TIDmR3:
276  return "TIDmR3";
278  return "TIDpR1R";
280  return "TIDpR1S";
282  return "TIDpR2R";
284  return "TIDpR2S";
285  case AlignmentPI::TIDpR3:
286  return "TIDpR3";
288  return "TECmR1R";
290  return "TECmR1S";
292  return "TECmR2R";
294  return "TECmR2S";
295  case AlignmentPI::TECmR3:
296  return "TECmR3";
297  case AlignmentPI::TECmR4:
298  return "TECmR4";
299  case AlignmentPI::TECmR5:
300  return "TECmR5";
301  case AlignmentPI::TECmR6:
302  return "TECmR6";
303  case AlignmentPI::TECmR7:
304  return "TECmR7";
306  return "TECpR1R";
308  return "TECpR1S";
310  return "TECpR2R";
312  return "TECpR2S";
313  case AlignmentPI::TECpR3:
314  return "TECpR3";
315  case AlignmentPI::TECpR4:
316  return "TECpR4";
317  case AlignmentPI::TECpR5:
318  return "TECpR5";
319  case AlignmentPI::TECpR6:
320  return "TECpR6";
321  case AlignmentPI::TECpR7:
322  return "TECpR7";
323  default:
324  edm::LogWarning("LogicError") << "Unknown partition: " << e;
325  return "";
326  }
327  }
328 
329  /*--------------------------------------------------------------------*/
330  inline bool isBPixOuterLadder(const DetId& detid, const TrackerTopology& tTopo, bool isPhase0)
331  /*--------------------------------------------------------------------*/
332  {
333  // Using TrackerTopology
334  // Ladders have a staggered structure
335  // Non-flipped ladders are on the outer radius
336  // Phase 0: Outer ladders are odd for layer 1,3 and even for layer 2
337  // Phase 1: Outer ladders are odd for layer 1,2,3 and even for layer 4
338  bool isOuter = false;
339  int layer = tTopo.pxbLayer(detid.rawId());
340  bool odd_ladder = tTopo.pxbLadder(detid.rawId()) % 2;
341  if (isPhase0) {
342  if (layer == 2)
343  isOuter = !odd_ladder;
344  else
345  isOuter = odd_ladder;
346  } else {
347  if (layer == 4)
348  isOuter = !odd_ladder;
349  else
350  isOuter = odd_ladder;
351  }
352  return isOuter;
353  }
354 
355  // ancillary struct to manage the topology
356  // info in a more compact way
357 
358  struct topolInfo {
359  private:
360  uint32_t m_rawid;
362  int m_layer;
363  int m_side;
364  int m_ring;
365  bool m_isRphi;
368 
369  public:
370  void init();
371  void fillGeometryInfo(const DetId& detId, const TrackerTopology& tTopo, bool isPhase0);
373  bool sanityCheck();
374  void printAll();
375  virtual ~topolInfo() {}
376  };
377 
378  /*--------------------------------------------------------------------*/
379  inline void topolInfo::printAll()
380  /*--------------------------------------------------------------------*/
381  {
382  std::cout << " detId:" << m_rawid << " subdetid: " << m_subdetid << " layer: " << m_layer << " side: " << m_side
383  << " ring: " << m_ring << " isRphi:" << m_isRphi << " isDoubleSide:" << m_isDoubleSide
384  << " isInternal:" << m_isInternal << std::endl;
385  }
386 
387  /*--------------------------------------------------------------------*/
388  inline void topolInfo::init()
389  /*--------------------------------------------------------------------*/
390  {
391  m_rawid = 0;
392  m_subdetid = -1;
393  m_layer = -1;
394  m_side = -1;
395  m_ring = -1;
396  m_isRphi = false;
397  m_isDoubleSide = false;
398  m_isInternal = false;
399  };
400 
401  /*--------------------------------------------------------------------*/
403  /*--------------------------------------------------------------------*/
404  {
405  if (m_layer == 0 || (m_subdetid == 1 && m_layer > 4) || (m_subdetid == 2 && m_layer > 3)) {
406  return false;
407  } else {
408  return true;
409  }
410  }
411  /*--------------------------------------------------------------------*/
412  inline void topolInfo::fillGeometryInfo(const DetId& detId, const TrackerTopology& tTopo, bool isPhase0)
413  /*--------------------------------------------------------------------*/
414  {
415  unsigned int subdetId = static_cast<unsigned int>(detId.subdetId());
416 
417  m_rawid = detId.rawId();
418  m_subdetid = subdetId;
419 
420  if (subdetId == StripSubdetector::TIB) {
421  m_layer = tTopo.tibLayer(detId.rawId());
422  m_side = tTopo.tibSide(detId.rawId());
423  m_isRphi = tTopo.isRPhi(detId.rawId());
424  m_isDoubleSide = tTopo.tibIsDoubleSide(detId.rawId());
425  m_isInternal = tTopo.tibIsInternalString(detId.rawId());
426  } else if (subdetId == StripSubdetector::TOB) {
427  m_layer = tTopo.tobLayer(detId.rawId());
428  m_side = tTopo.tobSide(detId.rawId());
429  m_isRphi = tTopo.isRPhi(detId.rawId());
430  m_isDoubleSide = tTopo.tobIsDoubleSide(detId.rawId());
431  m_isInternal = tTopo.tobModule(detId.rawId()) % 2;
432  } else if (subdetId == StripSubdetector::TID) {
433  m_layer = tTopo.tidWheel(detId.rawId());
434  m_side = tTopo.tidSide(detId.rawId());
435  m_isRphi = tTopo.isRPhi(detId.rawId());
436  m_ring = tTopo.tidRing(detId.rawId());
437  m_isDoubleSide = tTopo.tidIsDoubleSide(detId.rawId());
438  m_isInternal = tTopo.tidModuleInfo(detId.rawId())[0];
439  } else if (subdetId == StripSubdetector::TEC) {
440  m_layer = tTopo.tecWheel(detId.rawId());
441  m_side = tTopo.tecSide(detId.rawId());
442  m_isRphi = tTopo.isRPhi(detId.rawId());
443  m_ring = tTopo.tecRing(detId.rawId());
444  m_isDoubleSide = tTopo.tecIsDoubleSide(detId.rawId());
445  m_isInternal = tTopo.tecPetalInfo(detId.rawId())[0];
446  } else if (subdetId == PixelSubdetector::PixelBarrel) {
447  m_layer = tTopo.pxbLayer(detId.rawId());
448  m_isInternal = !AlignmentPI::isBPixOuterLadder(detId, tTopo, isPhase0);
449  } else if (subdetId == PixelSubdetector::PixelEndcap) {
450  m_layer = tTopo.pxfDisk(detId.rawId());
451  m_side = tTopo.pxfSide(detId.rawId());
452  } else
453  edm::LogWarning("LogicError") << "Unknown subdetid: " << subdetId;
454  }
455 
456  // ------------ method to assign a partition based on the topology struct info ---------------
457 
458  /*--------------------------------------------------------------------*/
460  /*--------------------------------------------------------------------*/
461  {
463 
464  if (m_isDoubleSide) {
466  }
467 
468  // BPix
469  if (m_subdetid == 1) {
470  switch (m_layer) {
471  case 1:
473  break;
474  case 2:
476  break;
477  case 3:
479  break;
480  case 4:
482  break;
483  default:
484  edm::LogWarning("LogicError") << "Unknow BPix layer: " << m_layer;
485  break;
486  }
487  // FPix
488  } else if (m_subdetid == 2) {
489  switch (m_layer) {
490  case 1:
492  break;
493  case 2:
495  break;
496  case 3:
498  break;
499  default:
500  edm::LogWarning("LogicError") << "Unknow FPix disk: " << m_layer;
501  break;
502  }
503  // TIB
504  } else if (m_subdetid == 3) {
505  switch (m_layer) {
506  case 1:
507  if (m_isRphi) {
509  } else {
511  }
512  break;
513  case 2:
514  if (m_isRphi) {
516  } else {
518  }
519  break;
520  case 3:
522  break;
523  case 4:
525  break;
526  default:
527  edm::LogWarning("LogicError") << "Unknow TIB layer: " << m_layer;
528  break;
529  }
530  // TID
531  } else if (m_subdetid == 4) {
532  switch (m_ring) {
533  case 1:
534  if (m_isRphi) {
536  } else {
538  }
539  break;
540  case 2:
541  if (m_isRphi) {
543  } else {
545  }
546  break;
547  case 3:
549  break;
550  default:
551  edm::LogWarning("LogicError") << "Unknow TID wheel: " << m_layer;
552  break;
553  }
554  // TOB
555  } else if (m_subdetid == 5) {
556  switch (m_layer) {
557  case 1:
558  if (m_isRphi) {
560  } else {
562  }
563  break;
564  case 2:
565  if (m_isRphi) {
567  } else {
569  }
570  break;
571  case 3:
573  break;
574  case 4:
576  break;
577  case 5:
579  break;
580  case 6:
582  break;
583  default:
584  edm::LogWarning("LogicError") << "Unknow TOB layer: " << m_layer;
585  break;
586  }
587  // TEC
588  } else if (m_subdetid == 6) {
589  switch (m_ring) {
590  case 1:
591  if (m_isRphi) {
593  } else {
595  }
596  break;
597  case 2:
598  if (m_isRphi) {
600  } else {
602  }
603  break;
604  case 3:
606  break;
607  case 4:
609  break;
610  case 5:
612  break;
613  case 6:
615  break;
616  case 7:
618  break;
619  default:
620  edm::LogWarning("LogicError") << "Unknow TEC ring: " << m_ring;
621  break;
622  }
623  }
624 
625  return ret;
626  }
627 
628  /*--------------------------------------------------------------------*/
630  /*--------------------------------------------------------------------*/
631  {
632  switch (coord) {
633  case t_x:
634  return "x-translation";
635  case t_y:
636  return "y-translation";
637  case t_z:
638  return "z-translation";
639  case rot_alpha:
640  return "#alpha angle rotation";
641  case rot_beta:
642  return "#beta angle rotation";
643  case rot_gamma:
644  return "#gamma angle rotation";
645  default:
646  return "should never be here!";
647  }
648  }
649 
650  /*--------------------------------------------------------------------*/
652  /*--------------------------------------------------------------------*/
653  {
654  switch (i) {
655  case XX:
656  return "XX";
657  case XY:
658  return "XY";
659  case XZ:
660  return "XZ";
661  case YZ:
662  return "YX";
663  case YY:
664  return "YY";
665  case ZZ:
666  return "ZZ";
667  default:
668  return "should never be here!";
669  }
670  }
671 
672  /*--------------------------------------------------------------------*/
674  /*--------------------------------------------------------------------*/
675  {
676  switch (i) {
677  case BPix:
678  return "BPix";
679  case FPix:
680  return "FPix";
681  case TIB:
682  return "TIB";
683  case TID:
684  return "TID";
685  case TOB:
686  return "TOB";
687  case TEC:
688  return "TEC";
689  default:
690  return "should never be here!";
691  }
692  }
693 
694  /*--------------------------------------------------------------------*/
695  inline std::pair<int, int> getIndices(AlignmentPI::index i)
696  /*--------------------------------------------------------------------*/
697  {
698  switch (i) {
699  case XX:
700  return std::make_pair(0, 0);
701  case XY:
702  return std::make_pair(0, 1);
703  case XZ:
704  return std::make_pair(0, 2);
705  case YZ:
706  return std::make_pair(1, 0);
707  case YY:
708  return std::make_pair(1, 1);
709  case ZZ:
710  return std::make_pair(2, 2);
711  default:
712  return std::make_pair(-1, -1);
713  }
714  }
715 
716  /*--------------------------------------------------------------------*/
717  inline void makeNicePlotStyle(TH1* hist, int color)
718  /*--------------------------------------------------------------------*/
719  {
720  hist->SetStats(kFALSE);
721 
722  hist->GetXaxis()->SetTitleColor(color);
723  hist->SetLineColor(color);
724  hist->SetTitleSize(0.08);
725  hist->SetLineWidth(2);
726  hist->GetXaxis()->CenterTitle(true);
727  hist->GetYaxis()->CenterTitle(true);
728  hist->GetXaxis()->SetTitleFont(42);
729  hist->GetYaxis()->SetTitleFont(42);
730  hist->GetXaxis()->SetNdivisions(505);
731  hist->GetXaxis()->SetTitleSize(0.06);
732  hist->GetYaxis()->SetTitleSize(0.06);
733  hist->GetXaxis()->SetTitleOffset(1.0);
734  hist->GetYaxis()->SetTitleOffset(1.3);
735  hist->GetXaxis()->SetLabelFont(42);
736  hist->GetYaxis()->SetLabelFont(42);
737  hist->GetYaxis()->SetLabelSize(.05);
738  hist->GetXaxis()->SetLabelSize(.05);
739  }
740 
741  /*--------------------------------------------------------------------*/
742  inline void makeNiceStats(TH1F* hist, AlignmentPI::partitions part, int color)
743  /*--------------------------------------------------------------------*/
744  {
745  char buffer[255];
746  TPaveText* stat = new TPaveText(0.60, 0.75, 0.95, 0.95, "NDC");
747  sprintf(buffer, "%s \n", AlignmentPI::getStringFromPart(part).c_str());
748  stat->AddText(buffer);
749 
750  sprintf(buffer, "Entries : %i\n", (int)hist->GetEntries());
751  stat->AddText(buffer);
752 
753  if (std::abs(hist->GetMean()) > 0.01) {
754  sprintf(buffer, "Mean : %6.2f\n", hist->GetMean());
755  } else {
756  sprintf(buffer, "Mean : %6.2f e-2\n", 100 * hist->GetMean());
757  }
758  stat->AddText(buffer);
759 
760  if (std::abs(hist->GetRMS()) > 0.01) {
761  sprintf(buffer, "RMS : %6.2f\n", hist->GetRMS());
762  } else {
763  sprintf(buffer, "RMS : %6.2f e-2\n", 100 * hist->GetRMS());
764  }
765  stat->AddText(buffer);
766 
767  stat->SetLineColor(color);
768  stat->SetTextColor(color);
769  stat->SetFillColor(10);
770  stat->SetShadowColor(10);
771  stat->Draw();
772  }
773 
774  /*--------------------------------------------------------------------*/
775  inline std::pair<float, float> getTheRange(std::map<uint32_t, float> values, const float nsigma)
776  /*--------------------------------------------------------------------*/
777  {
778  float sum = std::accumulate(
779  std::begin(values), std::end(values), 0.0, [](float value, const std::map<uint32_t, float>::value_type& p) {
780  return value + p.second;
781  });
782 
783  float m = sum / values.size();
784 
785  float accum = 0.0;
786  std::for_each(std::begin(values), std::end(values), [&](const std::map<uint32_t, float>::value_type& p) {
787  accum += (p.second - m) * (p.second - m);
788  });
789 
790  float stdev = sqrt(accum / (values.size() - 1));
791 
792  if (stdev != 0.) {
793  return std::make_pair(m - nsigma * stdev, m + nsigma * stdev);
794  } else {
795  return std::make_pair(m > 0. ? 0.95 * m : 1.05 * m, m > 0 ? 1.05 * m : 0.95 * m);
796  }
797  }
798 
799  /*--------------------------------------------------------------------*/
800  inline std::pair<double, double> calculatePosition(TVirtualPad* myPad, int boundary)
801  /*--------------------------------------------------------------------*/
802  {
803  int ix1;
804  int ix2;
805  int iw = myPad->GetWw();
806  int ih = myPad->GetWh();
807  double x1p, y1p, x2p, y2p;
808  myPad->GetPadPar(x1p, y1p, x2p, y2p);
809  ix1 = (int)(iw * x1p);
810  ix2 = (int)(iw * x2p);
811  double wndc = std::min(1., (double)iw / (double)ih);
812  double rw = wndc / (double)iw;
813  double x1ndc = (double)ix1 * rw;
814  double x2ndc = (double)ix2 * rw;
815  double rx1, ry1, rx2, ry2;
816  myPad->GetRange(rx1, ry1, rx2, ry2);
817  double rx = (x2ndc - x1ndc) / (rx2 - rx1);
818  double _sx;
819  _sx = rx * (boundary - rx1) + x1ndc;
820  double _dx = _sx + 0.05;
821 
822  return std::make_pair(_sx, _dx);
823  }
824 
825  // ancillary struct to manage the barycenters
826  // info in a more compact way
827 
829  std::map<AlignmentPI::PARTITION, double> Xbarycenters;
830  std::map<AlignmentPI::PARTITION, double> Ybarycenters;
831  std::map<AlignmentPI::PARTITION, double> Zbarycenters;
832  std::map<AlignmentPI::PARTITION, double> nmodules;
833 
834  public:
835  void init();
837  void computeBarycenters(const std::vector<AlignTransform>& input,
838  const TrackerTopology& tTopo,
839  const std::map<AlignmentPI::coordinate, float>& GPR);
840  const double getNModules(AlignmentPI::PARTITION p) { return nmodules[p]; };
841 
842  // M.M. 2020/01/09
843  // introduce methods for entire partitions, summing up the two sides of the
844  // endcap detectors
845 
846  /*--------------------------------------------------------------------*/
847  const std::array<double, 6> getX()
848  /*--------------------------------------------------------------------*/
849  {
850  return {{Xbarycenters[PARTITION::BPIX],
856  };
857 
858  /*--------------------------------------------------------------------*/
859  const std::array<double, 6> getY()
860  /*--------------------------------------------------------------------*/
861  {
862  return {{Ybarycenters[PARTITION::BPIX],
868  };
869 
870  /*--------------------------------------------------------------------*/
871  const std::array<double, 6> getZ()
872  /*--------------------------------------------------------------------*/
873  {
874  return {{Zbarycenters[PARTITION::BPIX],
880  };
881  virtual ~TkAlBarycenters() {}
882  };
883 
884  /*--------------------------------------------------------------------*/
886  /*--------------------------------------------------------------------*/
887  {
889  }
890 
891  /*--------------------------------------------------------------------*/
892  inline void TkAlBarycenters::computeBarycenters(const std::vector<AlignTransform>& input,
893  const TrackerTopology& tTopo,
894  const std::map<AlignmentPI::coordinate, float>& GPR)
895  /*--------------------------------------------------------------------*/
896  {
897  // zero in the n. modules per partition...
898  for (const auto& p : PARTITIONS) {
899  nmodules[p] = 0.;
900  }
901 
902  for (const auto& ali : input) {
903  if (DetId(ali.rawId()).det() != DetId::Tracker) {
904  edm::LogWarning("TkAlBarycenters::computeBarycenters")
905  << "Encountered invalid Tracker DetId:" << ali.rawId() << " " << DetId(ali.rawId()).det()
906  << " is different from " << DetId::Tracker << " - terminating ";
907  assert(DetId(ali.rawId()).det() != DetId::Tracker);
908  }
909 
910  int subid = DetId(ali.rawId()).subdetId();
911  switch (subid) {
913  Xbarycenters[PARTITION::BPIX] += (ali.translation().x());
914  Ybarycenters[PARTITION::BPIX] += (ali.translation().y());
915  Zbarycenters[PARTITION::BPIX] += (ali.translation().z());
917  break;
919 
920  // minus side
921  if (tTopo.pxfSide(DetId(ali.rawId())) == 1) {
922  Xbarycenters[PARTITION::FPIXm] += (ali.translation().x());
923  Ybarycenters[PARTITION::FPIXm] += (ali.translation().y());
924  Zbarycenters[PARTITION::FPIXm] += (ali.translation().z());
926  } // plus side
927  else {
928  Xbarycenters[PARTITION::FPIXp] += (ali.translation().x());
929  Ybarycenters[PARTITION::FPIXp] += (ali.translation().y());
930  Zbarycenters[PARTITION::FPIXp] += (ali.translation().z());
932  }
933  break;
935  Xbarycenters[PARTITION::TIB] += (ali.translation().x());
936  Ybarycenters[PARTITION::TIB] += (ali.translation().y());
937  Zbarycenters[PARTITION::TIB] += (ali.translation().z());
939  break;
941  // minus side
942  if (tTopo.tidSide(DetId(ali.rawId())) == 1) {
943  Xbarycenters[PARTITION::TIDm] += (ali.translation().x());
944  Ybarycenters[PARTITION::TIDm] += (ali.translation().y());
945  Zbarycenters[PARTITION::TIDm] += (ali.translation().z());
947  } // plus side
948  else {
949  Xbarycenters[PARTITION::TIDp] += (ali.translation().x());
950  Ybarycenters[PARTITION::TIDp] += (ali.translation().y());
951  Zbarycenters[PARTITION::TIDp] += (ali.translation().z());
953  }
954  break;
956  Xbarycenters[PARTITION::TOB] += (ali.translation().x());
957  Ybarycenters[PARTITION::TOB] += (ali.translation().y());
958  Zbarycenters[PARTITION::TOB] += (ali.translation().z());
960  break;
962  // minus side
963  if (tTopo.tecSide(DetId(ali.rawId())) == 1) {
964  Xbarycenters[PARTITION::TECm] += (ali.translation().x());
965  Ybarycenters[PARTITION::TECm] += (ali.translation().y());
966  Zbarycenters[PARTITION::TECm] += (ali.translation().z());
968  } // plus side
969  else {
970  Xbarycenters[PARTITION::TECp] += (ali.translation().x());
971  Ybarycenters[PARTITION::TECp] += (ali.translation().y());
972  Zbarycenters[PARTITION::TECp] += (ali.translation().z());
974  }
975  break;
976  default:
977  edm::LogError("TrackerAlignment_PayloadInspector") << "Unrecognized partition " << subid << std::endl;
978  break;
979  }
980  }
981 
982  for (const auto& p : PARTITIONS) {
983  // take the arithmetic mean
984  Xbarycenters[p] /= nmodules[p];
985  Ybarycenters[p] /= nmodules[p];
986  Zbarycenters[p] /= nmodules[p];
987 
988  // add the Tracker Global Position Record
989  Xbarycenters[p] += GPR.at(AlignmentPI::t_x);
990  Ybarycenters[p] += GPR.at(AlignmentPI::t_y);
991  Zbarycenters[p] += GPR.at(AlignmentPI::t_z);
992 
993  COUT << "Partition: " << p << " n. modules: " << nmodules[p] << "|"
994  << " X: " << std::right << std::setw(12) << Xbarycenters[p] << " Y: " << std::right << std::setw(12)
995  << Ybarycenters[p] << " Z: " << std::right << std::setw(12) << Zbarycenters[p] << std::endl;
996  }
997  }
998 
999  /*--------------------------------------------------------------------*/
1001  std::vector<int>& boundaries,
1002  const std::vector<AlignTransform>& ref_ali,
1003  const std::vector<AlignTransform>& target_ali,
1004  std::unique_ptr<TH1F>& compare)
1005  /*--------------------------------------------------------------------*/
1006  {
1007  int counter = 0; /* start the counter */
1009  for (unsigned int i = 0; i < ref_ali.size(); i++) {
1010  if (ref_ali[i].rawId() == target_ali[i].rawId()) {
1011  counter++;
1012  int subid = DetId(ref_ali[i].rawId()).subdetId();
1013 
1014  auto thePart = static_cast<AlignmentPI::partitions>(subid);
1015  if (thePart != currentPart) {
1016  currentPart = thePart;
1017  boundaries.push_back(counter);
1018  }
1019 
1020  CLHEP::HepRotation target_rot(target_ali[i].rotation());
1021  CLHEP::HepRotation ref_rot(ref_ali[i].rotation());
1022 
1023  align::RotationType target_ROT(target_rot.xx(),
1024  target_rot.xy(),
1025  target_rot.xz(),
1026  target_rot.yx(),
1027  target_rot.yy(),
1028  target_rot.yz(),
1029  target_rot.zx(),
1030  target_rot.zy(),
1031  target_rot.zz());
1032 
1033  align::RotationType ref_ROT(ref_rot.xx(),
1034  ref_rot.xy(),
1035  ref_rot.xz(),
1036  ref_rot.yx(),
1037  ref_rot.yy(),
1038  ref_rot.yz(),
1039  ref_rot.zx(),
1040  ref_rot.zy(),
1041  ref_rot.zz());
1042 
1043  const std::vector<double> deltaRot = {::deltaPhi(align::toAngles(target_ROT)[0], align::toAngles(ref_ROT)[0]),
1044  ::deltaPhi(align::toAngles(target_ROT)[1], align::toAngles(ref_ROT)[1]),
1045  ::deltaPhi(align::toAngles(target_ROT)[2], align::toAngles(ref_ROT)[2])};
1046 
1047  const auto& deltaTrans = target_ali[i].translation() - ref_ali[i].translation();
1048 
1049  switch (coord) {
1050  case AlignmentPI::t_x:
1051  compare->SetBinContent(i + 1, deltaTrans.x() * AlignmentPI::cmToUm);
1052  break;
1053  case AlignmentPI::t_y:
1054  compare->SetBinContent(i + 1, deltaTrans.y() * AlignmentPI::cmToUm);
1055  break;
1056  case AlignmentPI::t_z:
1057  compare->SetBinContent(i + 1, deltaTrans.z() * AlignmentPI::cmToUm);
1058  break;
1060  compare->SetBinContent(i + 1, deltaRot[0] * AlignmentPI::tomRad);
1061  break;
1062  case AlignmentPI::rot_beta:
1063  compare->SetBinContent(i + 1, deltaRot[1] * AlignmentPI::tomRad);
1064  break;
1066  compare->SetBinContent(i + 1, deltaRot[2] * AlignmentPI::tomRad);
1067  break;
1068  default:
1069  edm::LogError("TrackerAlignment_PayloadInspector") << "Unrecognized coordinate " << coord << std::endl;
1070  break;
1071  } // switch on the coordinate
1072  } // check on the same detID
1073  } // loop on the components
1074  }
1075 
1076  /*--------------------------------------------------------------------*/
1077  inline void fillComparisonHistograms(std::vector<int>& boundaries,
1078  const std::vector<AlignTransform>& ref_ali,
1079  const std::vector<AlignTransform>& target_ali,
1080  std::unordered_map<AlignmentPI::coordinate, std::unique_ptr<TH1F> >& compare,
1081  bool diff = false,
1083  /*--------------------------------------------------------------------*/
1084  {
1085  int counter = 0; /* start the counter */
1087  for (unsigned int i = 0; i < ref_ali.size(); i++) {
1088  if (ref_ali[i].rawId() == target_ali[i].rawId()) {
1089  counter++;
1090  int subid = DetId(ref_ali[i].rawId()).subdetId();
1091 
1092  auto thePart = static_cast<AlignmentPI::partitions>(subid);
1093 
1094  // in case it has to be filtered
1095  if (checkPart > 0 && thePart != checkPart) {
1096  continue;
1097  }
1098 
1099  if (thePart != currentPart) {
1100  currentPart = thePart;
1101  boundaries.push_back(counter);
1102  }
1103 
1104  CLHEP::HepRotation target_rot(target_ali[i].rotation());
1105  CLHEP::HepRotation ref_rot(ref_ali[i].rotation());
1106 
1107  align::RotationType target_ROT(target_rot.xx(),
1108  target_rot.xy(),
1109  target_rot.xz(),
1110  target_rot.yx(),
1111  target_rot.yy(),
1112  target_rot.yz(),
1113  target_rot.zx(),
1114  target_rot.zy(),
1115  target_rot.zz());
1116 
1117  align::RotationType ref_ROT(ref_rot.xx(),
1118  ref_rot.xy(),
1119  ref_rot.xz(),
1120  ref_rot.yx(),
1121  ref_rot.yy(),
1122  ref_rot.yz(),
1123  ref_rot.zx(),
1124  ref_rot.zy(),
1125  ref_rot.zz());
1126 
1127  const std::vector<double> deltaRot = {::deltaPhi(align::toAngles(target_ROT)[0], align::toAngles(ref_ROT)[0]),
1128  ::deltaPhi(align::toAngles(target_ROT)[1], align::toAngles(ref_ROT)[1]),
1129  ::deltaPhi(align::toAngles(target_ROT)[2], align::toAngles(ref_ROT)[2])};
1130 
1131  const auto& deltaTrans = target_ali[i].translation() - ref_ali[i].translation();
1132 
1133  // fill the histograms
1134  if (diff) {
1135  compare[AlignmentPI::t_x]->Fill(deltaTrans.x() * AlignmentPI::cmToUm);
1136  compare[AlignmentPI::t_y]->Fill(deltaTrans.y() * AlignmentPI::cmToUm);
1137  compare[AlignmentPI::t_z]->Fill(deltaTrans.z() * AlignmentPI::cmToUm);
1138 
1139  compare[AlignmentPI::rot_alpha]->Fill(deltaRot[0] * AlignmentPI::tomRad);
1140  compare[AlignmentPI::rot_beta]->Fill(deltaRot[1] * AlignmentPI::tomRad);
1141  compare[AlignmentPI::rot_gamma]->Fill(deltaRot[2] * AlignmentPI::tomRad);
1142  } else {
1143  compare[AlignmentPI::t_x]->SetBinContent(i + 1, deltaTrans.x() * AlignmentPI::cmToUm);
1144  compare[AlignmentPI::t_y]->SetBinContent(i + 1, deltaTrans.y() * AlignmentPI::cmToUm);
1145  compare[AlignmentPI::t_z]->SetBinContent(i + 1, deltaTrans.z() * AlignmentPI::cmToUm);
1146 
1147  compare[AlignmentPI::rot_alpha]->SetBinContent(i + 1, deltaRot[0] * AlignmentPI::tomRad);
1148  compare[AlignmentPI::rot_beta]->SetBinContent(i + 1, deltaRot[1] * AlignmentPI::tomRad);
1149  compare[AlignmentPI::rot_gamma]->SetBinContent(i + 1, deltaRot[2] * AlignmentPI::tomRad);
1150  }
1151 
1152  } // if it's the same detid
1153  } // loop on detids
1154  }
1155 
1156 } // namespace AlignmentPI
1157 
1158 #endif
const std::array< double, 6 > getY()
const double getNModules(AlignmentPI::PARTITION p)
static constexpr auto TEC
bool tibIsDoubleSide(const DetId &id) const
bool tecIsDoubleSide(const DetId &id) const
bool tidIsDoubleSide(const DetId &id) const
unsigned int tobLayer(const DetId &id) const
unsigned int pxbLayer(const DetId &id) const
unsigned int tibSide(const DetId &id) const
void fillComparisonHistogram(const AlignmentPI::coordinate &coord, std::vector< int > &boundaries, const std::vector< AlignTransform > &ref_ali, const std::vector< AlignTransform > &target_ali, std::unique_ptr< TH1F > &compare)
unsigned int tobSide(const DetId &id) const
double trim2PIs(const double phi, const double tolerance=1.f)
std::ostream & operator<<(std::ostream &o, PARTITION x)
std::vector< unsigned int > tidModuleInfo(const DetId &id) const
void makeNicePlotStyle(TH1 *hist, int color)
void computeBarycenters(const std::vector< AlignTransform > &input, const TrackerTopology &tTopo, const std::map< AlignmentPI::coordinate, float > &GPR)
ret
prodAgent to be discontinued
unsigned int tidSide(const DetId &id) const
const double tolerance
std::vector< unsigned int > tecPetalInfo(const DetId &id) const
unsigned int tidWheel(const DetId &id) const
unsigned int tecWheel(const DetId &id) const
const PARTITION PARTITIONS[(int) PARTITION::LAST+1]
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
std::string getStringFromPart(AlignmentPI::partitions i)
bool isBPixOuterLadder(const DetId &detid, const TrackerTopology &tTopo, bool isPhase0)
bool isRPhi(const DetId &id) const
unsigned int pxbLadder(const DetId &id) const
Log< level::Error, false > LogError
assert(be >=bs)
unsigned int tecRing(const DetId &id) const
ring id
bool tobIsDoubleSide(const DetId &id) const
static const float cmToUm
static std::string const input
Definition: EdmProvDump.cc:50
void makeNiceStats(TH1F *hist, AlignmentPI::partitions part, int color)
void fillComparisonHistograms(std::vector< int > &boundaries, const std::vector< AlignTransform > &ref_ali, const std::vector< AlignTransform > &target_ali, std::unordered_map< AlignmentPI::coordinate, std::unique_ptr< TH1F > > &compare, bool diff=false, AlignmentPI::partitions checkPart=AlignmentPI::INVALID)
const std::array< double, 6 > getX()
static const unsigned int phase0size
unsigned int tecSide(const DetId &id) const
static const float tomRad
T sqrt(T t)
Definition: SSEVec.h:19
void fillGeometryInfo(const DetId &detId, const TrackerTopology &tTopo, bool isPhase0)
unsigned int pxfDisk(const DetId &id) const
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double returnZeroIfNear2PI(const double phi)
double f[11][100]
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
Definition: value.py:1
std::map< AlignmentPI::PARTITION, double > Zbarycenters
EulerAngles toAngles(const RotationType &)
Convert rotation matrix to angles about x-, y-, z-axes (frame rotation).
Definition: Utilities.cc:8
static constexpr auto TOB
std::pair< int, int > getIndices(AlignmentPI::index i)
constexpr valType roundIfNear0(valType value, double tolerance=1.e-7)
Definition: Rounding.h:11
#define M_PI
std::pair< float, float > getTheRange(std::map< uint32_t, float > values, const float nsigma)
GlobalPoint getPartitionAvg(AlignmentPI::PARTITION p)
Definition: DetId.h:17
unsigned int pxfSide(const DetId &id) const
static constexpr auto TIB
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
std::string getStringFromCoordinate(AlignmentPI::coordinate coord)
part
Definition: HCALResponse.h:20
AlignmentPI::regions filterThePartition()
bool tibIsInternalString(const DetId &id) const
std::string getStringFromIndex(AlignmentPI::index i)
std::map< AlignmentPI::PARTITION, double > Ybarycenters
std::pair< double, double > calculatePosition(TVirtualPad *myPad, int boundary)
std::map< AlignmentPI::PARTITION, double > nmodules
float x
unsigned int tidRing(const DetId &id) const
def stdev(xlist)
Definition: plotscripts.py:69
std::map< AlignmentPI::PARTITION, double > Xbarycenters
unsigned int tibLayer(const DetId &id) const
std::string getStringFromRegionEnum(AlignmentPI::regions e)
unsigned int tobModule(const DetId &id) const
Log< level::Warning, false > LogWarning
const std::array< double, 6 > getZ()
static constexpr auto TID