Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-09-18 22:46:38

0001 #ifndef CONDCORE_ALIGNMENTPLUGINS_ALIGNMENTPAYLOADINSPECTORHELPER_H
0002 #define CONDCORE_ALIGNMENTPLUGINS_ALIGNMENTPAYLOADINSPECTORHELPER_H
0003 
0004 #include <iostream>
0005 #include <algorithm>
0006 #include <vector>
0007 #include <numeric>
0008 #include <string>
0009 #include "TH1.h"
0010 #include "TCanvas.h"
0011 #include "TPaveStats.h"
0012 #include "TStyle.h"
0013 #include "TList.h"
0014 #include "Alignment/CommonAlignment/interface/Utilities.h"
0015 #include "CondFormats/Alignment/interface/Alignments.h"
0016 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0017 #include "DataFormats/Math/interface/deltaPhi.h"  // for deltaPhi
0018 #include "DataFormats/Math/interface/Rounding.h"  // for rounding
0019 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0020 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0021 
0022 //#define MMDEBUG // uncomment for debugging at compile time
0023 #ifdef MMDEBUG
0024 #include <iostream>
0025 #define COUT std::cout << "MM "
0026 #else
0027 #define COUT edm::LogVerbatim("")
0028 #endif
0029 
0030 namespace AlignmentPI {
0031 
0032   // size of the phase-I Tracker APE payload (including both SS + DS modules)
0033   static const unsigned int phase0size = 19876;
0034   static const float cmToUm = 10000.f;
0035   static const float tomRad = 1000.f;
0036 
0037   // method to zero all elements whose difference from 2Pi
0038   // is less than the tolerance (2*10e-7)
0039   inline double returnZeroIfNear2PI(const double phi) {
0040     const double tol = 2.e-7;  // default tolerance 1.e-7 doesn't account for possible variations
0041     if (cms_rounding::roundIfNear0(std::abs(phi) - 2 * M_PI, tol) == 0.f) {
0042       return 0.f;
0043     } else {
0044       return phi;
0045     }
0046   }
0047 
0048   // method to bring back around 0 all elements whose  difference
0049   // frm 2Pi is is less than the tolerance (in micro-rad)
0050   inline double trim2PIs(const double phi, const double tolerance = 1.f) {
0051     if (std::abs((std::abs(phi) - 2 * M_PI) * tomRad) < tolerance) {
0052       return (std::abs(phi) - 2 * M_PI);
0053     } else {
0054       return phi;
0055     }
0056   }
0057 
0058   enum coordinate {
0059     t_x = 1,
0060     t_y = 2,
0061     t_z = 3,
0062     rot_alpha = 4,
0063     rot_beta = 5,
0064     rot_gamma = 6,
0065   };
0066 
0067   // M.M. 2017/09/12
0068   // As the matrix is symmetric, we map only 6/9 terms
0069   // More terms for the extended APE can be added to the following methods
0070 
0071   enum index { XX = 1, XY = 2, XZ = 3, YZ = 4, YY = 5, ZZ = 6 };
0072 
0073   enum partitions { BPix = 1, FPix = 2, TIB = 3, TID = 4, TOB = 5, TEC = 6 };
0074 
0075   enum class PARTITION {
0076     BPIX,   // 0 Barrel Pixel
0077     FPIXp,  // 1 Forward Pixel Plus
0078     FPIXm,  // 2 Forward Pixel Minus
0079     TIB,    // 3 Tracker Inner Barrel
0080     TIDp,   // 4 Tracker Inner Disks Plus
0081     TIDm,   // 5 Tracker Inner Disks Minus
0082     TOB,    // 6 Tracker Outer Barrel
0083     TECp,   // 7 Tracker Endcaps Plus
0084     TECm,   // 8 Tracker Endcaps Minus
0085     LAST = TECm
0086   };
0087 
0088   extern const PARTITION PARTITIONS[(int)PARTITION::LAST + 1];
0089   const PARTITION PARTITIONS[] = {PARTITION::BPIX,
0090                                   PARTITION::FPIXp,
0091                                   PARTITION::FPIXm,
0092                                   PARTITION::TIB,
0093                                   PARTITION::TIDp,
0094                                   PARTITION::TIDm,
0095                                   PARTITION::TOB,
0096                                   PARTITION::TECp,
0097                                   PARTITION::TECm};
0098 
0099   inline std::ostream& operator<<(std::ostream& o, PARTITION x) {
0100     return o << std::underlying_type<PARTITION>::type(x);
0101   }
0102 
0103   enum regions {
0104     BPixL1o,          //0  Barrel Pixel Layer 1 outer
0105     BPixL1i,          //1  Barrel Pixel Layer 1 inner
0106     BPixL2o,          //2  Barrel Pixel Layer 2 outer
0107     BPixL2i,          //3  Barrel Pixel Layer 2 inner
0108     BPixL3o,          //4  Barrel Pixel Layer 3 outer
0109     BPixL3i,          //5  Barrel Pixel Layer 3 inner
0110     BPixL4o,          //6  Barrel Pixel Layer 4 outer
0111     BPixL4i,          //7  Barrel Pixel Layer 4 inner
0112     FPixmL1,          //8  Forward Pixel Minus side Disk 1
0113     FPixmL2,          //9 Forward Pixel Minus side Disk 2
0114     FPixmL3,          //10 Forward Pixel Minus side Disk 3
0115     FPixpL1,          //11 Forward Pixel Plus side Disk 1
0116     FPixpL2,          //12 Forward Pixel Plus side Disk 2
0117     FPixpL3,          //13 Forward Pixel Plus side Disk 3
0118     TIBL1Ro,          //14 Inner Barrel Layer 1 Rphi outer
0119     TIBL1Ri,          //15 Inner Barrel Layer 1 Rphi inner
0120     TIBL1So,          //16 Inner Barrel Layer 1 Stereo outer
0121     TIBL1Si,          //17 Inner Barrel Layer 1 Stereo inner
0122     TIBL2Ro,          //18 Inner Barrel Layer 2 Rphi outer
0123     TIBL2Ri,          //19 Inner Barrel Layer 2 Rphi inner
0124     TIBL2So,          //20 Inner Barrel Layer 2 Stereo outer
0125     TIBL2Si,          //21 Inner Barrel Layer 2 Stereo inner
0126     TIBL3o,           //22 Inner Barrel Layer 3 outer
0127     TIBL3i,           //23 Inner Barrel Layer 3 inner
0128     TIBL4o,           //24 Inner Barrel Layer 4 outer
0129     TIBL4i,           //25 Inner Barrel Layer 4 inner
0130     TOBL1Ro,          //26 Outer Barrel Layer 1 Rphi outer
0131     TOBL1Ri,          //27 Outer Barrel Layer 1 Rphi inner
0132     TOBL1So,          //28 Outer Barrel Layer 1 Stereo outer
0133     TOBL1Si,          //29 Outer Barrel Layer 1 Stereo inner
0134     TOBL2Ro,          //30 Outer Barrel Layer 2 Rphi outer
0135     TOBL2Ri,          //31 Outer Barrel Layer 2 Rphi inner
0136     TOBL2So,          //32 Outer Barrel Layer 2 Stereo outer
0137     TOBL2Si,          //33 Outer Barrel Layer 2 Stereo inner
0138     TOBL3o,           //34 Outer Barrel Layer 3 outer
0139     TOBL3i,           //35 Outer Barrel Layer 3 inner
0140     TOBL4o,           //36 Outer Barrel Layer 4 outer
0141     TOBL4i,           //37 Outer Barrel Layer 4 inner
0142     TOBL5o,           //38 Outer Barrel Layer 5 outer
0143     TOBL5i,           //39 Outer Barrel Layer 5 inner
0144     TOBL6o,           //40 Outer Barrel Layer 6 outer
0145     TOBL6i,           //41 Outer Barrel Layer 6 inner
0146     TIDmR1R,          //42 Inner Disk Minus side Ring 1 Rphi
0147     TIDmR1S,          //43 Inner Disk Minus side Ring 1 Stereo
0148     TIDmR2R,          //44 Inner Disk Minus side Ring 2 Rphi
0149     TIDmR2S,          //45 Inner Disk Minus side Ring 2 Stereo
0150     TIDmR3,           //46 Inner Disk Minus side Ring 3
0151     TIDpR1R,          //47 Inner Disk Plus side Ring 1 Rphi
0152     TIDpR1S,          //48 Inner Disk Plus side Ring 1 Stereo
0153     TIDpR2R,          //49 Inner Disk Plus side Ring 2 Rphi
0154     TIDpR2S,          //50 Inner Disk Plus side Ring 2 Stereo
0155     TIDpR3,           //51 Inner Disk Plus side Ring 3
0156     TECmR1R,          //52 Endcaps Minus side Ring 1 Rphi
0157     TECmR1S,          //53 Endcaps Minus side Ring 1 Stereo
0158     TECmR2R,          //54 Encdaps Minus side Ring 2 Rphi
0159     TECmR2S,          //55 Endcaps Minus side Ring 2 Stereo
0160     TECmR3,           //56 Endcaps Minus side Ring 3
0161     TECmR4,           //57 Endcaps Minus side Ring 4
0162     TECmR5,           //58 Endcaps Minus side Ring 5
0163     TECmR6,           //59 Endcaps Minus side Ring 6
0164     TECmR7,           //60 Endcaps Minus side Ring 7
0165     TECpR1R,          //61 Endcaps Plus side Ring 1 Rphi
0166     TECpR1S,          //62 Endcaps Plus side Ring 1 Stereo
0167     TECpR2R,          //63 Encdaps Plus side Ring 2 Rphi
0168     TECpR2S,          //64 Endcaps Plus side Ring 2 Stereo
0169     TECpR3,           //65 Endcaps Plus side Ring 3
0170     TECpR4,           //66 Endcaps Plus side Ring 4
0171     TECpR5,           //67 Endcaps Plus side Ring 5
0172     TECpR6,           //68 Endcaps Plus side Ring 6
0173     TECpR7,           //67 Endcaps Plus side Ring 7
0174     StripDoubleSide,  // 70 -- not to be considered
0175     NUM_OF_REGIONS    // 71 -- default
0176   };
0177 
0178   /*--------------------------------------------------------------------*/
0179   inline std::string getStringFromRegionEnum(AlignmentPI::regions e)
0180   /*--------------------------------------------------------------------*/
0181   {
0182     switch (e) {
0183       case AlignmentPI::BPixL1o:
0184         return "BPixL1o";
0185       case AlignmentPI::BPixL1i:
0186         return "BPixL1i";
0187       case AlignmentPI::BPixL2o:
0188         return "BPixL2o";
0189       case AlignmentPI::BPixL2i:
0190         return "BPixL2i";
0191       case AlignmentPI::BPixL3o:
0192         return "BPixL3o";
0193       case AlignmentPI::BPixL3i:
0194         return "BPixL3i";
0195       case AlignmentPI::BPixL4o:
0196         return "BPixL4o";
0197       case AlignmentPI::BPixL4i:
0198         return "BPixL4i";
0199       case AlignmentPI::FPixmL1:
0200         return "FPixmL1";
0201       case AlignmentPI::FPixmL2:
0202         return "FPixmL2";
0203       case AlignmentPI::FPixmL3:
0204         return "FPixmL3";
0205       case AlignmentPI::FPixpL1:
0206         return "FPixpL1";
0207       case AlignmentPI::FPixpL2:
0208         return "FPixpL2";
0209       case AlignmentPI::FPixpL3:
0210         return "FPixpL3";
0211       case AlignmentPI::TIBL1Ro:
0212         return "TIBL1Ro";
0213       case AlignmentPI::TIBL1Ri:
0214         return "TIBL1Ri";
0215       case AlignmentPI::TIBL1So:
0216         return "TIBL1So";
0217       case AlignmentPI::TIBL1Si:
0218         return "TIBL1Si";
0219       case AlignmentPI::TIBL2Ro:
0220         return "TIBL2Ro";
0221       case AlignmentPI::TIBL2Ri:
0222         return "TIBL2Ri";
0223       case AlignmentPI::TIBL2So:
0224         return "TIBL2So";
0225       case AlignmentPI::TIBL2Si:
0226         return "TIBL2Si";
0227       case AlignmentPI::TIBL3o:
0228         return "TIBL3o";
0229       case AlignmentPI::TIBL3i:
0230         return "TIBL3i";
0231       case AlignmentPI::TIBL4o:
0232         return "TIBL4o";
0233       case AlignmentPI::TIBL4i:
0234         return "TIBL4i";
0235       case AlignmentPI::TOBL1Ro:
0236         return "TOBL1Ro";
0237       case AlignmentPI::TOBL1Ri:
0238         return "TOBL1Ri";
0239       case AlignmentPI::TOBL1So:
0240         return "TOBL1So";
0241       case AlignmentPI::TOBL1Si:
0242         return "TOBL1Si";
0243       case AlignmentPI::TOBL2Ro:
0244         return "TOBL2Ro";
0245       case AlignmentPI::TOBL2Ri:
0246         return "TOBL2Ri";
0247       case AlignmentPI::TOBL2So:
0248         return "TOBL2So";
0249       case AlignmentPI::TOBL2Si:
0250         return "TOBL2Si";
0251       case AlignmentPI::TOBL3o:
0252         return "TOBL3o";
0253       case AlignmentPI::TOBL3i:
0254         return "TOBL3i";
0255       case AlignmentPI::TOBL4o:
0256         return "TOBL4o";
0257       case AlignmentPI::TOBL4i:
0258         return "TOBL4i";
0259       case AlignmentPI::TOBL5o:
0260         return "TOBL5o";
0261       case AlignmentPI::TOBL5i:
0262         return "TOBL5i";
0263       case AlignmentPI::TOBL6o:
0264         return "TOBL6o";
0265       case AlignmentPI::TOBL6i:
0266         return "TOBL6i";
0267       case AlignmentPI::TIDmR1R:
0268         return "TIDmR1R";
0269       case AlignmentPI::TIDmR1S:
0270         return "TIDmR1S";
0271       case AlignmentPI::TIDmR2R:
0272         return "TIDmR2R";
0273       case AlignmentPI::TIDmR2S:
0274         return "TIDmR2S";
0275       case AlignmentPI::TIDmR3:
0276         return "TIDmR3";
0277       case AlignmentPI::TIDpR1R:
0278         return "TIDpR1R";
0279       case AlignmentPI::TIDpR1S:
0280         return "TIDpR1S";
0281       case AlignmentPI::TIDpR2R:
0282         return "TIDpR2R";
0283       case AlignmentPI::TIDpR2S:
0284         return "TIDpR2S";
0285       case AlignmentPI::TIDpR3:
0286         return "TIDpR3";
0287       case AlignmentPI::TECmR1R:
0288         return "TECmR1R";
0289       case AlignmentPI::TECmR1S:
0290         return "TECmR1S";
0291       case AlignmentPI::TECmR2R:
0292         return "TECmR2R";
0293       case AlignmentPI::TECmR2S:
0294         return "TECmR2S";
0295       case AlignmentPI::TECmR3:
0296         return "TECmR3";
0297       case AlignmentPI::TECmR4:
0298         return "TECmR4";
0299       case AlignmentPI::TECmR5:
0300         return "TECmR5";
0301       case AlignmentPI::TECmR6:
0302         return "TECmR6";
0303       case AlignmentPI::TECmR7:
0304         return "TECmR7";
0305       case AlignmentPI::TECpR1R:
0306         return "TECpR1R";
0307       case AlignmentPI::TECpR1S:
0308         return "TECpR1S";
0309       case AlignmentPI::TECpR2R:
0310         return "TECpR2R";
0311       case AlignmentPI::TECpR2S:
0312         return "TECpR2S";
0313       case AlignmentPI::TECpR3:
0314         return "TECpR3";
0315       case AlignmentPI::TECpR4:
0316         return "TECpR4";
0317       case AlignmentPI::TECpR5:
0318         return "TECpR5";
0319       case AlignmentPI::TECpR6:
0320         return "TECpR6";
0321       case AlignmentPI::TECpR7:
0322         return "TECpR7";
0323       default:
0324         edm::LogWarning("LogicError") << "Unknown partition: " << e;
0325         return "";
0326     }
0327   }
0328 
0329   /*--------------------------------------------------------------------*/
0330   inline bool isBPixOuterLadder(const DetId& detid, const TrackerTopology& tTopo, bool isPhase0)
0331   /*--------------------------------------------------------------------*/
0332   {
0333     bool isOuter = false;
0334     int layer = tTopo.pxbLayer(detid.rawId());
0335     bool odd_ladder = tTopo.pxbLadder(detid.rawId()) % 2;
0336     if (isPhase0) {
0337       if (layer == 2)
0338         isOuter = !odd_ladder;
0339       else
0340         isOuter = odd_ladder;
0341     } else {
0342       if (layer == 4)
0343         isOuter = odd_ladder;
0344       else
0345         isOuter = !odd_ladder;
0346     }
0347     return isOuter;
0348   }
0349 
0350   // ancillary struct to manage the topology
0351   // info in a more compact way
0352 
0353   struct topolInfo {
0354   private:
0355     uint32_t m_rawid;
0356     int m_subdetid;
0357     int m_layer;
0358     int m_side;
0359     int m_ring;
0360     bool m_isRphi;
0361     bool m_isDoubleSide;
0362     bool m_isInternal;
0363 
0364   public:
0365     void init();
0366     void fillGeometryInfo(const DetId& detId, const TrackerTopology& tTopo, bool isPhase0);
0367     AlignmentPI::regions filterThePartition();
0368     bool sanityCheck();
0369     void printAll();
0370     virtual ~topolInfo() {}
0371   };
0372 
0373   /*--------------------------------------------------------------------*/
0374   inline void topolInfo::printAll()
0375   /*--------------------------------------------------------------------*/
0376   {
0377     std::cout << " detId:" << m_rawid << " subdetid: " << m_subdetid << " layer: " << m_layer << " side: " << m_side
0378               << " ring: " << m_ring << " isRphi:" << m_isRphi << " isDoubleSide:" << m_isDoubleSide
0379               << " isInternal:" << m_isInternal << std::endl;
0380   }
0381 
0382   /*--------------------------------------------------------------------*/
0383   inline void topolInfo::init()
0384   /*--------------------------------------------------------------------*/
0385   {
0386     m_rawid = 0;
0387     m_subdetid = -1;
0388     m_layer = -1;
0389     m_side = -1;
0390     m_ring = -1;
0391     m_isRphi = false;
0392     m_isDoubleSide = false;
0393     m_isInternal = false;
0394   };
0395 
0396   /*--------------------------------------------------------------------*/
0397   inline bool topolInfo::sanityCheck()
0398   /*--------------------------------------------------------------------*/
0399   {
0400     if (m_layer == 0 || (m_subdetid == 1 && m_layer > 4) || (m_subdetid == 2 && m_layer > 3)) {
0401       return false;
0402     } else {
0403       return true;
0404     }
0405   }
0406   /*--------------------------------------------------------------------*/
0407   inline void topolInfo::fillGeometryInfo(const DetId& detId, const TrackerTopology& tTopo, bool isPhase0)
0408   /*--------------------------------------------------------------------*/
0409   {
0410     unsigned int subdetId = static_cast<unsigned int>(detId.subdetId());
0411 
0412     m_rawid = detId.rawId();
0413     m_subdetid = subdetId;
0414 
0415     if (subdetId == StripSubdetector::TIB) {
0416       m_layer = tTopo.tibLayer(detId.rawId());
0417       m_side = tTopo.tibSide(detId.rawId());
0418       m_isRphi = tTopo.isRPhi(detId.rawId());
0419       m_isDoubleSide = tTopo.tibIsDoubleSide(detId.rawId());
0420       m_isInternal = tTopo.tibIsInternalString(detId.rawId());
0421     } else if (subdetId == StripSubdetector::TOB) {
0422       m_layer = tTopo.tobLayer(detId.rawId());
0423       m_side = tTopo.tobSide(detId.rawId());
0424       m_isRphi = tTopo.isRPhi(detId.rawId());
0425       m_isDoubleSide = tTopo.tobIsDoubleSide(detId.rawId());
0426       m_isInternal = tTopo.tobModule(detId.rawId()) % 2;
0427     } else if (subdetId == StripSubdetector::TID) {
0428       m_layer = tTopo.tidWheel(detId.rawId());
0429       m_side = tTopo.tidSide(detId.rawId());
0430       m_isRphi = tTopo.isRPhi(detId.rawId());
0431       m_ring = tTopo.tidRing(detId.rawId());
0432       m_isDoubleSide = tTopo.tidIsDoubleSide(detId.rawId());
0433       m_isInternal = tTopo.tidModuleInfo(detId.rawId())[0];
0434     } else if (subdetId == StripSubdetector::TEC) {
0435       m_layer = tTopo.tecWheel(detId.rawId());
0436       m_side = tTopo.tecSide(detId.rawId());
0437       m_isRphi = tTopo.isRPhi(detId.rawId());
0438       m_ring = tTopo.tecRing(detId.rawId());
0439       m_isDoubleSide = tTopo.tecIsDoubleSide(detId.rawId());
0440       m_isInternal = tTopo.tecPetalInfo(detId.rawId())[0];
0441     } else if (subdetId == PixelSubdetector::PixelBarrel) {
0442       m_layer = tTopo.pxbLayer(detId.rawId());
0443       m_isInternal = !AlignmentPI::isBPixOuterLadder(detId, tTopo, isPhase0);
0444     } else if (subdetId == PixelSubdetector::PixelEndcap) {
0445       m_layer = tTopo.pxfDisk(detId.rawId());
0446       m_side = tTopo.pxfSide(detId.rawId());
0447     } else
0448       edm::LogWarning("LogicError") << "Unknown subdetid: " << subdetId;
0449   }
0450 
0451   // ------------ method to assign a partition based on the topology struct info ---------------
0452 
0453   /*--------------------------------------------------------------------*/
0454   inline AlignmentPI::regions topolInfo::filterThePartition()
0455   /*--------------------------------------------------------------------*/
0456   {
0457     AlignmentPI::regions ret = AlignmentPI::NUM_OF_REGIONS;
0458 
0459     if (m_isDoubleSide) {
0460       return AlignmentPI::StripDoubleSide;
0461     }
0462 
0463     // BPix
0464     if (m_subdetid == 1) {
0465       switch (m_layer) {
0466         case 1:
0467           m_isInternal > 0 ? ret = AlignmentPI::BPixL1o : ret = AlignmentPI::BPixL1i;
0468           break;
0469         case 2:
0470           m_isInternal > 0 ? ret = AlignmentPI::BPixL2o : ret = AlignmentPI::BPixL2i;
0471           break;
0472         case 3:
0473           m_isInternal > 0 ? ret = AlignmentPI::BPixL3o : ret = AlignmentPI::BPixL3i;
0474           break;
0475         case 4:
0476           m_isInternal > 0 ? ret = AlignmentPI::BPixL4o : ret = AlignmentPI::BPixL4i;
0477           break;
0478         default:
0479           edm::LogWarning("LogicError") << "Unknow BPix layer: " << m_layer;
0480           break;
0481       }
0482       // FPix
0483     } else if (m_subdetid == 2) {
0484       switch (m_layer) {
0485         case 1:
0486           m_side > 1 ? ret = AlignmentPI::FPixpL1 : ret = AlignmentPI::FPixmL1;
0487           break;
0488         case 2:
0489           m_side > 1 ? ret = AlignmentPI::FPixpL2 : ret = AlignmentPI::FPixmL2;
0490           break;
0491         case 3:
0492           m_side > 1 ? ret = AlignmentPI::FPixpL3 : ret = AlignmentPI::FPixmL3;
0493           break;
0494         default:
0495           edm::LogWarning("LogicError") << "Unknow FPix disk: " << m_layer;
0496           break;
0497       }
0498       // TIB
0499     } else if (m_subdetid == 3) {
0500       switch (m_layer) {
0501         case 1:
0502           if (m_isRphi) {
0503             m_isInternal > 0 ? ret = AlignmentPI::TIBL1Ro : ret = AlignmentPI::TIBL1Ri;
0504           } else {
0505             m_isInternal > 0 ? ret = AlignmentPI::TIBL1So : ret = AlignmentPI::TIBL1Si;
0506           }
0507           break;
0508         case 2:
0509           if (m_isRphi) {
0510             m_isInternal > 0 ? ret = AlignmentPI::TIBL2Ro : ret = AlignmentPI::TIBL2Ri;
0511           } else {
0512             m_isInternal > 0 ? ret = AlignmentPI::TIBL2So : ret = AlignmentPI::TIBL2Si;
0513           }
0514           break;
0515         case 3:
0516           m_isInternal > 0 ? ret = AlignmentPI::TIBL3o : ret = AlignmentPI::TIBL3i;
0517           break;
0518         case 4:
0519           m_isInternal > 0 ? ret = AlignmentPI::TIBL4o : ret = AlignmentPI::TIBL4i;
0520           break;
0521         default:
0522           edm::LogWarning("LogicError") << "Unknow TIB layer: " << m_layer;
0523           break;
0524       }
0525       // TID
0526     } else if (m_subdetid == 4) {
0527       switch (m_ring) {
0528         case 1:
0529           if (m_isRphi) {
0530             m_side > 1 ? ret = AlignmentPI::TIDpR1R : ret = AlignmentPI::TIDmR1R;
0531           } else {
0532             m_side > 1 ? ret = AlignmentPI::TIDpR1S : ret = AlignmentPI::TIDmR1S;
0533           }
0534           break;
0535         case 2:
0536           if (m_isRphi) {
0537             m_side > 1 ? ret = AlignmentPI::TIDpR2R : ret = AlignmentPI::TIDmR2R;
0538           } else {
0539             m_side > 1 ? ret = AlignmentPI::TIDpR2S : ret = AlignmentPI::TIDmR2S;
0540           }
0541           break;
0542         case 3:
0543           m_side > 1 ? ret = AlignmentPI::TIDpR3 : ret = AlignmentPI::TIDmR3;
0544           break;
0545         default:
0546           edm::LogWarning("LogicError") << "Unknow TID wheel: " << m_layer;
0547           break;
0548       }
0549       // TOB
0550     } else if (m_subdetid == 5) {
0551       switch (m_layer) {
0552         case 1:
0553           if (m_isRphi) {
0554             m_isInternal > 0 ? ret = AlignmentPI::TOBL1Ro : ret = AlignmentPI::TOBL1Ri;
0555           } else {
0556             m_isInternal > 0 ? ret = AlignmentPI::TOBL1So : ret = AlignmentPI::TOBL1Si;
0557           }
0558           break;
0559         case 2:
0560           if (m_isRphi) {
0561             m_isInternal > 0 ? ret = AlignmentPI::TOBL2Ro : ret = AlignmentPI::TOBL2Ri;
0562           } else {
0563             m_isInternal > 0 ? ret = AlignmentPI::TOBL2So : ret = AlignmentPI::TOBL2Si;
0564           }
0565           break;
0566         case 3:
0567           m_isInternal > 0 ? ret = AlignmentPI::TOBL3o : ret = AlignmentPI::TOBL3i;
0568           break;
0569         case 4:
0570           m_isInternal > 0 ? ret = AlignmentPI::TOBL4o : ret = AlignmentPI::TOBL4i;
0571           break;
0572         case 5:
0573           m_isInternal > 0 ? ret = AlignmentPI::TOBL5o : ret = AlignmentPI::TOBL5i;
0574           break;
0575         case 6:
0576           m_isInternal > 0 ? ret = AlignmentPI::TOBL6o : ret = AlignmentPI::TOBL6i;
0577           break;
0578         default:
0579           edm::LogWarning("LogicError") << "Unknow TOB layer: " << m_layer;
0580           break;
0581       }
0582       // TEC
0583     } else if (m_subdetid == 6) {
0584       switch (m_ring) {
0585         case 1:
0586           if (m_isRphi) {
0587             m_side > 1 ? ret = AlignmentPI::TECpR1R : ret = AlignmentPI::TECmR1R;
0588           } else {
0589             m_side > 1 ? ret = AlignmentPI::TECpR1S : ret = AlignmentPI::TECmR1S;
0590           }
0591           break;
0592         case 2:
0593           if (m_isRphi) {
0594             m_side > 1 ? ret = AlignmentPI::TECpR2R : ret = AlignmentPI::TECmR2R;
0595           } else {
0596             m_side > 1 ? ret = AlignmentPI::TECpR2S : ret = AlignmentPI::TECmR2S;
0597           }
0598           break;
0599         case 3:
0600           m_side > 1 ? ret = AlignmentPI::TECpR3 : ret = AlignmentPI::TECmR3;
0601           break;
0602         case 4:
0603           m_side > 1 ? ret = AlignmentPI::TECpR4 : ret = AlignmentPI::TECmR4;
0604           break;
0605         case 5:
0606           m_side > 1 ? ret = AlignmentPI::TECpR5 : ret = AlignmentPI::TECmR5;
0607           break;
0608         case 6:
0609           m_side > 1 ? ret = AlignmentPI::TECpR6 : ret = AlignmentPI::TECmR6;
0610           break;
0611         case 7:
0612           m_side > 1 ? ret = AlignmentPI::TECpR7 : ret = AlignmentPI::TECmR7;
0613           break;
0614         default:
0615           edm::LogWarning("LogicError") << "Unknow TEC ring: " << m_ring;
0616           break;
0617       }
0618     }
0619 
0620     return ret;
0621   }
0622 
0623   /*--------------------------------------------------------------------*/
0624   inline std::string getStringFromCoordinate(AlignmentPI::coordinate coord)
0625   /*--------------------------------------------------------------------*/
0626   {
0627     switch (coord) {
0628       case t_x:
0629         return "x-translation";
0630       case t_y:
0631         return "y-translation";
0632       case t_z:
0633         return "z-translation";
0634       case rot_alpha:
0635         return "#alpha angle rotation";
0636       case rot_beta:
0637         return "#beta angle rotation";
0638       case rot_gamma:
0639         return "#gamma angle rotation";
0640       default:
0641         return "should never be here!";
0642     }
0643   }
0644 
0645   /*--------------------------------------------------------------------*/
0646   inline std::string getStringFromIndex(AlignmentPI::index i)
0647   /*--------------------------------------------------------------------*/
0648   {
0649     switch (i) {
0650       case XX:
0651         return "XX";
0652       case XY:
0653         return "XY";
0654       case XZ:
0655         return "XZ";
0656       case YZ:
0657         return "YX";
0658       case YY:
0659         return "YY";
0660       case ZZ:
0661         return "ZZ";
0662       default:
0663         return "should never be here!";
0664     }
0665   }
0666 
0667   /*--------------------------------------------------------------------*/
0668   inline std::string getStringFromPart(AlignmentPI::partitions i)
0669   /*--------------------------------------------------------------------*/
0670   {
0671     switch (i) {
0672       case BPix:
0673         return "BPix";
0674       case FPix:
0675         return "FPix";
0676       case TIB:
0677         return "TIB";
0678       case TID:
0679         return "TID";
0680       case TOB:
0681         return "TOB";
0682       case TEC:
0683         return "TEC";
0684       default:
0685         return "should never be here!";
0686     }
0687   }
0688 
0689   /*--------------------------------------------------------------------*/
0690   inline std::pair<int, int> getIndices(AlignmentPI::index i)
0691   /*--------------------------------------------------------------------*/
0692   {
0693     switch (i) {
0694       case XX:
0695         return std::make_pair(0, 0);
0696       case XY:
0697         return std::make_pair(0, 1);
0698       case XZ:
0699         return std::make_pair(0, 2);
0700       case YZ:
0701         return std::make_pair(1, 0);
0702       case YY:
0703         return std::make_pair(1, 1);
0704       case ZZ:
0705         return std::make_pair(2, 2);
0706       default:
0707         return std::make_pair(-1, -1);
0708     }
0709   }
0710 
0711   /*--------------------------------------------------------------------*/
0712   inline void makeNicePlotStyle(TH1* hist, int color)
0713   /*--------------------------------------------------------------------*/
0714   {
0715     hist->SetStats(kFALSE);
0716 
0717     hist->GetXaxis()->SetTitleColor(color);
0718     hist->SetLineColor(color);
0719     hist->SetTitleSize(0.08);
0720     hist->SetLineWidth(2);
0721     hist->GetXaxis()->CenterTitle(true);
0722     hist->GetYaxis()->CenterTitle(true);
0723     hist->GetXaxis()->SetTitleFont(42);
0724     hist->GetYaxis()->SetTitleFont(42);
0725     hist->GetXaxis()->SetNdivisions(505);
0726     hist->GetXaxis()->SetTitleSize(0.06);
0727     hist->GetYaxis()->SetTitleSize(0.06);
0728     hist->GetXaxis()->SetTitleOffset(1.0);
0729     hist->GetYaxis()->SetTitleOffset(1.3);
0730     hist->GetXaxis()->SetLabelFont(42);
0731     hist->GetYaxis()->SetLabelFont(42);
0732     hist->GetYaxis()->SetLabelSize(.05);
0733     hist->GetXaxis()->SetLabelSize(.05);
0734   }
0735 
0736   /*--------------------------------------------------------------------*/
0737   inline void makeNiceStats(TH1F* hist, AlignmentPI::partitions part, int color)
0738   /*--------------------------------------------------------------------*/
0739   {
0740     char buffer[255];
0741     TPaveText* stat = new TPaveText(0.60, 0.75, 0.95, 0.95, "NDC");
0742     sprintf(buffer, "%s \n", AlignmentPI::getStringFromPart(part).c_str());
0743     stat->AddText(buffer);
0744 
0745     sprintf(buffer, "Entries : %i\n", (int)hist->GetEntries());
0746     stat->AddText(buffer);
0747 
0748     if (std::abs(hist->GetMean()) > 0.01) {
0749       sprintf(buffer, "Mean    : %6.2f\n", hist->GetMean());
0750     } else {
0751       sprintf(buffer, "Mean    : %6.2f e-2\n", 100 * hist->GetMean());
0752     }
0753     stat->AddText(buffer);
0754 
0755     if (std::abs(hist->GetRMS()) > 0.01) {
0756       sprintf(buffer, "RMS     : %6.2f\n", hist->GetRMS());
0757     } else {
0758       sprintf(buffer, "RMS     : %6.2f e-2\n", 100 * hist->GetRMS());
0759     }
0760     stat->AddText(buffer);
0761 
0762     stat->SetLineColor(color);
0763     stat->SetTextColor(color);
0764     stat->SetFillColor(10);
0765     stat->SetShadowColor(10);
0766     stat->Draw();
0767   }
0768 
0769   /*--------------------------------------------------------------------*/
0770   inline std::pair<float, float> getTheRange(std::map<uint32_t, float> values, const float nsigma)
0771   /*--------------------------------------------------------------------*/
0772   {
0773     float sum = std::accumulate(
0774         std::begin(values), std::end(values), 0.0, [](float value, const std::map<uint32_t, float>::value_type& p) {
0775           return value + p.second;
0776         });
0777 
0778     float m = sum / values.size();
0779 
0780     float accum = 0.0;
0781     std::for_each(std::begin(values), std::end(values), [&](const std::map<uint32_t, float>::value_type& p) {
0782       accum += (p.second - m) * (p.second - m);
0783     });
0784 
0785     float stdev = sqrt(accum / (values.size() - 1));
0786 
0787     if (stdev != 0.) {
0788       return std::make_pair(m - nsigma * stdev, m + nsigma * stdev);
0789     } else {
0790       return std::make_pair(m > 0. ? 0.95 * m : 1.05 * m, m > 0 ? 1.05 * m : 0.95 * m);
0791     }
0792   }
0793 
0794   /*--------------------------------------------------------------------*/
0795   inline std::pair<double, double> calculatePosition(TVirtualPad* myPad, int boundary)
0796   /*--------------------------------------------------------------------*/
0797   {
0798     int ix1;
0799     int ix2;
0800     int iw = myPad->GetWw();
0801     int ih = myPad->GetWh();
0802     double x1p, y1p, x2p, y2p;
0803     myPad->GetPadPar(x1p, y1p, x2p, y2p);
0804     ix1 = (int)(iw * x1p);
0805     ix2 = (int)(iw * x2p);
0806     double wndc = std::min(1., (double)iw / (double)ih);
0807     double rw = wndc / (double)iw;
0808     double x1ndc = (double)ix1 * rw;
0809     double x2ndc = (double)ix2 * rw;
0810     double rx1, ry1, rx2, ry2;
0811     myPad->GetRange(rx1, ry1, rx2, ry2);
0812     double rx = (x2ndc - x1ndc) / (rx2 - rx1);
0813     double _sx;
0814     _sx = rx * (boundary - rx1) + x1ndc;
0815     double _dx = _sx + 0.05;
0816 
0817     return std::make_pair(_sx, _dx);
0818   }
0819 
0820   // ancillary struct to manage the barycenters
0821   // info in a more compact way
0822 
0823   struct TkAlBarycenters {
0824     std::map<AlignmentPI::PARTITION, double> Xbarycenters;
0825     std::map<AlignmentPI::PARTITION, double> Ybarycenters;
0826     std::map<AlignmentPI::PARTITION, double> Zbarycenters;
0827     std::map<AlignmentPI::PARTITION, double> nmodules;
0828 
0829   public:
0830     void init();
0831     GlobalPoint getPartitionAvg(AlignmentPI::PARTITION p);
0832     void computeBarycenters(const std::vector<AlignTransform>& input,
0833                             const TrackerTopology& tTopo,
0834                             const std::map<AlignmentPI::coordinate, float>& GPR);
0835     const double getNModules(AlignmentPI::PARTITION p) { return nmodules[p]; };
0836 
0837     // M.M. 2020/01/09
0838     // introduce methods for entire partitions, summing up the two sides of the
0839     // endcap detectors
0840 
0841     /*--------------------------------------------------------------------*/
0842     const std::array<double, 6> getX()
0843     /*--------------------------------------------------------------------*/
0844     {
0845       return {{Xbarycenters[PARTITION::BPIX],
0846                (Xbarycenters[PARTITION::FPIXm] + Xbarycenters[PARTITION::FPIXp]) / 2.,
0847                Xbarycenters[PARTITION::TIB],
0848                (Xbarycenters[PARTITION::TIDm] + Xbarycenters[PARTITION::TIDp]) / 2,
0849                Xbarycenters[PARTITION::TOB],
0850                (Xbarycenters[PARTITION::TECm] + Xbarycenters[PARTITION::TECp]) / 2}};
0851     };
0852 
0853     /*--------------------------------------------------------------------*/
0854     const std::array<double, 6> getY()
0855     /*--------------------------------------------------------------------*/
0856     {
0857       return {{Ybarycenters[PARTITION::BPIX],
0858                (Ybarycenters[PARTITION::FPIXm] + Ybarycenters[PARTITION::FPIXp]) / 2.,
0859                Ybarycenters[PARTITION::TIB],
0860                (Ybarycenters[PARTITION::TIDm] + Ybarycenters[PARTITION::TIDp]) / 2,
0861                Ybarycenters[PARTITION::TOB],
0862                (Ybarycenters[PARTITION::TECm] + Ybarycenters[PARTITION::TECp]) / 2}};
0863     };
0864 
0865     /*--------------------------------------------------------------------*/
0866     const std::array<double, 6> getZ()
0867     /*--------------------------------------------------------------------*/
0868     {
0869       return {{Zbarycenters[PARTITION::BPIX],
0870                (Zbarycenters[PARTITION::FPIXm] + Zbarycenters[PARTITION::FPIXp]) / 2.,
0871                Zbarycenters[PARTITION::TIB],
0872                (Zbarycenters[PARTITION::TIDm] + Zbarycenters[PARTITION::TIDp]) / 2,
0873                Zbarycenters[PARTITION::TOB],
0874                (Zbarycenters[PARTITION::TECm] + Zbarycenters[PARTITION::TECp]) / 2}};
0875     };
0876     virtual ~TkAlBarycenters() {}
0877   };
0878 
0879   /*--------------------------------------------------------------------*/
0880   inline GlobalPoint TkAlBarycenters::getPartitionAvg(AlignmentPI::PARTITION p)
0881   /*--------------------------------------------------------------------*/
0882   {
0883     return GlobalPoint(Xbarycenters[p], Ybarycenters[p], Zbarycenters[p]);
0884   }
0885 
0886   /*--------------------------------------------------------------------*/
0887   inline void TkAlBarycenters::computeBarycenters(const std::vector<AlignTransform>& input,
0888                                                   const TrackerTopology& tTopo,
0889                                                   const std::map<AlignmentPI::coordinate, float>& GPR)
0890   /*--------------------------------------------------------------------*/
0891   {
0892     // zero in the n. modules per partition...
0893     for (const auto& p : PARTITIONS) {
0894       nmodules[p] = 0.;
0895     }
0896 
0897     for (const auto& ali : input) {
0898       if (DetId(ali.rawId()).det() != DetId::Tracker) {
0899         edm::LogWarning("TkAlBarycenters::computeBarycenters")
0900             << "Encountered invalid Tracker DetId:" << ali.rawId() << " " << DetId(ali.rawId()).det()
0901             << " is different from " << DetId::Tracker << "  - terminating ";
0902         assert(DetId(ali.rawId()).det() != DetId::Tracker);
0903       }
0904 
0905       int subid = DetId(ali.rawId()).subdetId();
0906       switch (subid) {
0907         case PixelSubdetector::PixelBarrel:
0908           Xbarycenters[PARTITION::BPIX] += (ali.translation().x());
0909           Ybarycenters[PARTITION::BPIX] += (ali.translation().y());
0910           Zbarycenters[PARTITION::BPIX] += (ali.translation().z());
0911           nmodules[PARTITION::BPIX]++;
0912           break;
0913         case PixelSubdetector::PixelEndcap:
0914 
0915           // minus side
0916           if (tTopo.pxfSide(DetId(ali.rawId())) == 1) {
0917             Xbarycenters[PARTITION::FPIXm] += (ali.translation().x());
0918             Ybarycenters[PARTITION::FPIXm] += (ali.translation().y());
0919             Zbarycenters[PARTITION::FPIXm] += (ali.translation().z());
0920             nmodules[PARTITION::FPIXm]++;
0921           }  // plus side
0922           else {
0923             Xbarycenters[PARTITION::FPIXp] += (ali.translation().x());
0924             Ybarycenters[PARTITION::FPIXp] += (ali.translation().y());
0925             Zbarycenters[PARTITION::FPIXp] += (ali.translation().z());
0926             nmodules[PARTITION::FPIXp]++;
0927           }
0928           break;
0929         case StripSubdetector::TIB:
0930           Xbarycenters[PARTITION::TIB] += (ali.translation().x());
0931           Ybarycenters[PARTITION::TIB] += (ali.translation().y());
0932           Zbarycenters[PARTITION::TIB] += (ali.translation().z());
0933           nmodules[PARTITION::TIB]++;
0934           break;
0935         case StripSubdetector::TID:
0936           // minus side
0937           if (tTopo.tidSide(DetId(ali.rawId())) == 1) {
0938             Xbarycenters[PARTITION::TIDm] += (ali.translation().x());
0939             Ybarycenters[PARTITION::TIDm] += (ali.translation().y());
0940             Zbarycenters[PARTITION::TIDm] += (ali.translation().z());
0941             nmodules[PARTITION::TIDm]++;
0942           }  // plus side
0943           else {
0944             Xbarycenters[PARTITION::TIDp] += (ali.translation().x());
0945             Ybarycenters[PARTITION::TIDp] += (ali.translation().y());
0946             Zbarycenters[PARTITION::TIDp] += (ali.translation().z());
0947             nmodules[PARTITION::TIDp]++;
0948           }
0949           break;
0950         case StripSubdetector::TOB:
0951           Xbarycenters[PARTITION::TOB] += (ali.translation().x());
0952           Ybarycenters[PARTITION::TOB] += (ali.translation().y());
0953           Zbarycenters[PARTITION::TOB] += (ali.translation().z());
0954           nmodules[PARTITION::TOB]++;
0955           break;
0956         case StripSubdetector::TEC:
0957           // minus side
0958           if (tTopo.tecSide(DetId(ali.rawId())) == 1) {
0959             Xbarycenters[PARTITION::TECm] += (ali.translation().x());
0960             Ybarycenters[PARTITION::TECm] += (ali.translation().y());
0961             Zbarycenters[PARTITION::TECm] += (ali.translation().z());
0962             nmodules[PARTITION::TECm]++;
0963           }  // plus side
0964           else {
0965             Xbarycenters[PARTITION::TECp] += (ali.translation().x());
0966             Ybarycenters[PARTITION::TECp] += (ali.translation().y());
0967             Zbarycenters[PARTITION::TECp] += (ali.translation().z());
0968             nmodules[PARTITION::TECp]++;
0969           }
0970           break;
0971         default:
0972           edm::LogError("TrackerAlignment_PayloadInspector") << "Unrecognized partition " << subid << std::endl;
0973           break;
0974       }
0975     }
0976 
0977     for (const auto& p : PARTITIONS) {
0978       // take the arithmetic mean
0979       Xbarycenters[p] /= nmodules[p];
0980       Ybarycenters[p] /= nmodules[p];
0981       Zbarycenters[p] /= nmodules[p];
0982 
0983       // add the Tracker Global Position Record
0984       Xbarycenters[p] += GPR.at(AlignmentPI::t_x);
0985       Ybarycenters[p] += GPR.at(AlignmentPI::t_y);
0986       Zbarycenters[p] += GPR.at(AlignmentPI::t_z);
0987 
0988       COUT << "Partition: " << p << " n. modules: " << nmodules[p] << "|"
0989            << " X: " << std::right << std::setw(12) << Xbarycenters[p] << " Y: " << std::right << std::setw(12)
0990            << Ybarycenters[p] << " Z: " << std::right << std::setw(12) << Zbarycenters[p] << std::endl;
0991     }
0992   }
0993 
0994   /*--------------------------------------------------------------------*/
0995   inline void fillComparisonHistogram(const AlignmentPI::coordinate& coord,
0996                                       std::vector<int>& boundaries,
0997                                       const std::vector<AlignTransform>& ref_ali,
0998                                       const std::vector<AlignTransform>& target_ali,
0999                                       std::unique_ptr<TH1F>& compare)
1000   /*--------------------------------------------------------------------*/
1001   {
1002     int counter = 0; /* start the counter */
1003     AlignmentPI::partitions currentPart = AlignmentPI::BPix;
1004     for (unsigned int i = 0; i < ref_ali.size(); i++) {
1005       if (ref_ali[i].rawId() == target_ali[i].rawId()) {
1006         counter++;
1007         int subid = DetId(ref_ali[i].rawId()).subdetId();
1008 
1009         auto thePart = static_cast<AlignmentPI::partitions>(subid);
1010         if (thePart != currentPart) {
1011           currentPart = thePart;
1012           boundaries.push_back(counter);
1013         }
1014 
1015         CLHEP::HepRotation target_rot(target_ali[i].rotation());
1016         CLHEP::HepRotation ref_rot(ref_ali[i].rotation());
1017 
1018         align::RotationType target_ROT(target_rot.xx(),
1019                                        target_rot.xy(),
1020                                        target_rot.xz(),
1021                                        target_rot.yx(),
1022                                        target_rot.yy(),
1023                                        target_rot.yz(),
1024                                        target_rot.zx(),
1025                                        target_rot.zy(),
1026                                        target_rot.zz());
1027 
1028         align::RotationType ref_ROT(ref_rot.xx(),
1029                                     ref_rot.xy(),
1030                                     ref_rot.xz(),
1031                                     ref_rot.yx(),
1032                                     ref_rot.yy(),
1033                                     ref_rot.yz(),
1034                                     ref_rot.zx(),
1035                                     ref_rot.zy(),
1036                                     ref_rot.zz());
1037 
1038         const std::vector<double> deltaRot = {::deltaPhi(align::toAngles(target_ROT)[0], align::toAngles(ref_ROT)[0]),
1039                                               ::deltaPhi(align::toAngles(target_ROT)[1], align::toAngles(ref_ROT)[1]),
1040                                               ::deltaPhi(align::toAngles(target_ROT)[2], align::toAngles(ref_ROT)[2])};
1041 
1042         const auto& deltaTrans = target_ali[i].translation() - ref_ali[i].translation();
1043 
1044         switch (coord) {
1045           case AlignmentPI::t_x:
1046             compare->SetBinContent(i + 1, deltaTrans.x() * AlignmentPI::cmToUm);
1047             break;
1048           case AlignmentPI::t_y:
1049             compare->SetBinContent(i + 1, deltaTrans.y() * AlignmentPI::cmToUm);
1050             break;
1051           case AlignmentPI::t_z:
1052             compare->SetBinContent(i + 1, deltaTrans.z() * AlignmentPI::cmToUm);
1053             break;
1054           case AlignmentPI::rot_alpha:
1055             compare->SetBinContent(i + 1, deltaRot[0] * AlignmentPI::tomRad);
1056             break;
1057           case AlignmentPI::rot_beta:
1058             compare->SetBinContent(i + 1, deltaRot[1] * AlignmentPI::tomRad);
1059             break;
1060           case AlignmentPI::rot_gamma:
1061             compare->SetBinContent(i + 1, deltaRot[2] * AlignmentPI::tomRad);
1062             break;
1063           default:
1064             edm::LogError("TrackerAlignment_PayloadInspector") << "Unrecognized coordinate " << coord << std::endl;
1065             break;
1066         }  // switch on the coordinate
1067       }    // check on the same detID
1068     }      // loop on the components
1069   }
1070 
1071   /*--------------------------------------------------------------------*/
1072   inline void fillComparisonHistograms(std::vector<int>& boundaries,
1073                                        const std::vector<AlignTransform>& ref_ali,
1074                                        const std::vector<AlignTransform>& target_ali,
1075                                        std::unordered_map<AlignmentPI::coordinate, std::unique_ptr<TH1F> >& compare)
1076   /*--------------------------------------------------------------------*/
1077   {
1078     int counter = 0; /* start the counter */
1079     AlignmentPI::partitions currentPart = AlignmentPI::BPix;
1080     for (unsigned int i = 0; i < ref_ali.size(); i++) {
1081       if (ref_ali[i].rawId() == target_ali[i].rawId()) {
1082         counter++;
1083         int subid = DetId(ref_ali[i].rawId()).subdetId();
1084 
1085         auto thePart = static_cast<AlignmentPI::partitions>(subid);
1086         if (thePart != currentPart) {
1087           currentPart = thePart;
1088           boundaries.push_back(counter);
1089         }
1090 
1091         CLHEP::HepRotation target_rot(target_ali[i].rotation());
1092         CLHEP::HepRotation ref_rot(ref_ali[i].rotation());
1093 
1094         align::RotationType target_ROT(target_rot.xx(),
1095                                        target_rot.xy(),
1096                                        target_rot.xz(),
1097                                        target_rot.yx(),
1098                                        target_rot.yy(),
1099                                        target_rot.yz(),
1100                                        target_rot.zx(),
1101                                        target_rot.zy(),
1102                                        target_rot.zz());
1103 
1104         align::RotationType ref_ROT(ref_rot.xx(),
1105                                     ref_rot.xy(),
1106                                     ref_rot.xz(),
1107                                     ref_rot.yx(),
1108                                     ref_rot.yy(),
1109                                     ref_rot.yz(),
1110                                     ref_rot.zx(),
1111                                     ref_rot.zy(),
1112                                     ref_rot.zz());
1113 
1114         const std::vector<double> deltaRot = {::deltaPhi(align::toAngles(target_ROT)[0], align::toAngles(ref_ROT)[0]),
1115                                               ::deltaPhi(align::toAngles(target_ROT)[1], align::toAngles(ref_ROT)[1]),
1116                                               ::deltaPhi(align::toAngles(target_ROT)[2], align::toAngles(ref_ROT)[2])};
1117 
1118         const auto& deltaTrans = target_ali[i].translation() - ref_ali[i].translation();
1119 
1120         // fill the histograms
1121         compare[AlignmentPI::t_x]->SetBinContent(i + 1, deltaTrans.x() * AlignmentPI::cmToUm);
1122         compare[AlignmentPI::t_y]->SetBinContent(i + 1, deltaTrans.y() * AlignmentPI::cmToUm);
1123         compare[AlignmentPI::t_z]->SetBinContent(i + 1, deltaTrans.z() * AlignmentPI::cmToUm);
1124 
1125         compare[AlignmentPI::rot_alpha]->SetBinContent(i + 1, deltaRot[0] * AlignmentPI::tomRad);
1126         compare[AlignmentPI::rot_beta]->SetBinContent(i + 1, deltaRot[1] * AlignmentPI::tomRad);
1127         compare[AlignmentPI::rot_gamma]->SetBinContent(i + 1, deltaRot[2] * AlignmentPI::tomRad);
1128 
1129       }  // if it's the same detid
1130     }    // loop on detids
1131   }
1132 
1133 }  // namespace AlignmentPI
1134 
1135 #endif