Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:36:19

0001 /*!
0002   \file TrackerAlignmentErrorExtended_PayloadInspector
0003   \Payload Inspector Plugin for Tracker Alignment 
0004   \author M. Musich
0005   \version $Revision: 1.0 $
0006   \date $Date: 2017/07/10 10:59:24 $
0007 */
0008 
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010 
0011 #include "CondCore/Utilities/interface/PayloadInspectorModule.h"
0012 #include "CondCore/Utilities/interface/PayloadInspector.h"
0013 #include "CondCore/CondDB/interface/Time.h"
0014 
0015 // the data format of the condition to be inspected
0016 #include "CondFormats/Alignment/interface/Alignments.h"
0017 #include "DataFormats/DetId/interface/DetId.h"
0018 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0019 
0020 //#include "CondFormats/Alignment/interface/Definitions.h"
0021 #include "CLHEP/Vector/RotationInterfaces.h"
0022 #include "Alignment/CommonAlignment/interface/Utilities.h"
0023 
0024 // needed for the tracker map
0025 #include "CommonTools/TrackerMap/interface/TrackerMap.h"
0026 
0027 // needed for mapping
0028 #include "CondCore/AlignmentPlugins/interface/AlignmentPayloadInspectorHelper.h"
0029 #include "CalibTracker/StandaloneTrackerTopology/interface/StandaloneTrackerTopology.h"
0030 #include "CalibTracker/SiPixelESProducers/interface/SiPixelDetInfoFileReader.h"
0031 #include "DQM/TrackerRemapper/interface/Phase1PixelSummaryMap.h"
0032 
0033 #include <boost/range/adaptor/indexed.hpp>
0034 #include <iomanip>  // std::setprecision
0035 #include <iostream>
0036 #include <memory>
0037 #include <sstream>
0038 
0039 // include ROOT
0040 #include "TH2F.h"
0041 #include "TGaxis.h"
0042 #include "TLegend.h"
0043 #include "TCanvas.h"
0044 #include "TLine.h"
0045 #include "TStyle.h"
0046 #include "TLatex.h"
0047 #include "TPave.h"
0048 #include "TMarker.h"
0049 #include "TPaveStats.h"
0050 
0051 namespace {
0052 
0053   using namespace cond::payloadInspector;
0054 
0055   // M.M. 2017/09/29
0056   // Hardcoded Tracker Global Position Record
0057   // Without accessing the ES, it is not possible to access to the GPR with the PI technology,
0058   // so this needs to be hardcoded.
0059   // Anyway it is not likely to change until a new Tracker is installed.
0060   // Details at:
0061   // - https://indico.cern.ch/event/238026/contributions/513928/attachments/400000/556192/mm_TkAlMeeting_28_03_2013.pdf
0062   // - https://twiki.cern.ch/twiki/bin/view/CMS/TkAlignmentPixelPosition
0063 
0064   const std::map<AlignmentPI::coordinate, float> hardcodeGPR = {
0065       {AlignmentPI::t_x, -9.00e-02}, {AlignmentPI::t_y, -1.10e-01}, {AlignmentPI::t_z, -1.70e-01}};
0066 
0067   //*******************************************/
0068   // Size of the movement over all partitions,
0069   // one at a time
0070   //******************************************//
0071 
0072   template <int ntags, IOVMultiplicity nIOVs, bool doOnlyPixel>
0073   class TrackerAlignmentCompareAll : public PlotImage<Alignments, nIOVs, ntags> {
0074   public:
0075     TrackerAlignmentCompareAll()
0076         : PlotImage<Alignments, nIOVs, ntags>("comparison of all coordinates between two geometries") {}
0077 
0078     bool fill() override {
0079       TGaxis::SetExponentOffset(-0.12, 0.01, "y");  // Y offset
0080 
0081       // trick to deal with the multi-ioved tag and two tag case at the same time
0082       auto theIOVs = PlotBase::getTag<0>().iovs;
0083       auto tagname1 = PlotBase::getTag<0>().name;
0084       std::string tagname2 = "";
0085       auto firstiov = theIOVs.front();
0086       std::tuple<cond::Time_t, cond::Hash> lastiov;
0087 
0088       // we don't support (yet) comparison with more than 2 tags
0089       assert(this->m_plotAnnotations.ntags < 3);
0090 
0091       if (this->m_plotAnnotations.ntags == 2) {
0092         auto tag2iovs = PlotBase::getTag<1>().iovs;
0093         tagname2 = PlotBase::getTag<1>().name;
0094         lastiov = tag2iovs.front();
0095       } else {
0096         lastiov = theIOVs.back();
0097       }
0098 
0099       std::shared_ptr<Alignments> last_payload = this->fetchPayload(std::get<1>(lastiov));
0100       std::shared_ptr<Alignments> first_payload = this->fetchPayload(std::get<1>(firstiov));
0101 
0102       std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
0103       std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
0104 
0105       std::vector<AlignTransform> ref_ali = first_payload->m_align;
0106       std::vector<AlignTransform> target_ali = last_payload->m_align;
0107 
0108       // Use remove_if along with a lambda expression to remove elements based on the condition (subid > 2)
0109       if (doOnlyPixel) {
0110         ref_ali.erase(std::remove_if(ref_ali.begin(),
0111                                      ref_ali.end(),
0112                                      [](const AlignTransform &transform) {
0113                                        int subid = DetId(transform.rawId()).subdetId();
0114                                        return subid > 2;
0115                                      }),
0116                       ref_ali.end());
0117 
0118         target_ali.erase(std::remove_if(target_ali.begin(),
0119                                         target_ali.end(),
0120                                         [](const AlignTransform &transform) {
0121                                           int subid = DetId(transform.rawId()).subdetId();
0122                                           return subid > 2;
0123                                         }),
0124                          target_ali.end());
0125       }
0126       TCanvas canvas("Alignment Comparison", "Alignment Comparison", 2000, 1200);
0127       canvas.Divide(3, 2);
0128 
0129       if (ref_ali.size() != target_ali.size()) {
0130         edm::LogError("TrackerAlignment_PayloadInspector")
0131             << "the size of the reference alignment (" << ref_ali.size()
0132             << ") is different from the one of the target (" << target_ali.size()
0133             << ")! You are probably trying to compare different underlying geometries. Exiting";
0134         return false;
0135       }
0136 
0137       const bool ph2 = (ref_ali.size() > AlignmentPI::phase1size);
0138 
0139       // check that the geomtery is a tracker one
0140       const char *path_toTopologyXML = nullptr;
0141       if (ph2) {
0142         path_toTopologyXML = "Geometry/TrackerCommonData/data/PhaseII/trackerParameters.xml";
0143       } else {
0144         path_toTopologyXML = (ref_ali.size() == AlignmentPI::phase0size)
0145                                  ? "Geometry/TrackerCommonData/data/trackerParameters.xml"
0146                                  : "Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml";
0147       }
0148 
0149       TrackerTopology tTopo =
0150           StandaloneTrackerTopology::fromTrackerParametersXMLFile(edm::FileInPath(path_toTopologyXML).fullPath());
0151 
0152       for (const auto &ali : ref_ali) {
0153         auto mydetid = ali.rawId();
0154         if (DetId(mydetid).det() != DetId::Tracker) {
0155           edm::LogWarning("TrackerAlignment_PayloadInspector")
0156               << "Encountered invalid Tracker DetId:" << DetId(mydetid).rawId() << " (" << DetId(mydetid).det()
0157               << ") is different from " << DetId::Tracker << " (is DoubleSide: " << tTopo.tidIsDoubleSide(mydetid)
0158               << "); subdetId " << DetId(mydetid).subdetId() << " - terminating ";
0159           return false;
0160         }
0161       }
0162 
0163       const std::vector<AlignmentPI::coordinate> coords = {AlignmentPI::t_x,
0164                                                            AlignmentPI::t_y,
0165                                                            AlignmentPI::t_z,
0166                                                            AlignmentPI::rot_alpha,
0167                                                            AlignmentPI::rot_beta,
0168                                                            AlignmentPI::rot_gamma};
0169 
0170       std::unordered_map<AlignmentPI::coordinate, std::unique_ptr<TH1F> > diffs;
0171 
0172       // generate the map of histograms
0173       for (const auto &coord : coords) {
0174         auto s_coord = AlignmentPI::getStringFromCoordinate(coord);
0175         std::string unit =
0176             (coord == AlignmentPI::t_x || coord == AlignmentPI::t_y || coord == AlignmentPI::t_z) ? "[#mum]" : "[mrad]";
0177 
0178         diffs[coord] = std::make_unique<TH1F>(Form("comparison_%s", s_coord.c_str()),
0179                                               Form(";Detector Id index; #Delta%s %s", s_coord.c_str(), unit.c_str()),
0180                                               ref_ali.size(),
0181                                               -0.5,
0182                                               ref_ali.size() - 0.5);
0183       }
0184 
0185       // fill all the histograms together
0186       std::map<int, AlignmentPI::partitions> boundaries;
0187       boundaries.insert({0, AlignmentPI::BPix});  // always start with BPix, not filled in the loop
0188       AlignmentPI::fillComparisonHistograms(boundaries, ref_ali, target_ali, diffs);
0189 
0190       unsigned int subpad{1};
0191       TLegend legend = TLegend(0.17, 0.84, 0.95, 0.94);
0192       legend.SetTextSize(0.023);
0193       for (const auto &coord : coords) {
0194         canvas.cd(subpad);
0195         canvas.cd(subpad)->SetTopMargin(0.06);
0196         canvas.cd(subpad)->SetLeftMargin(0.17);
0197         canvas.cd(subpad)->SetRightMargin(0.05);
0198         canvas.cd(subpad)->SetBottomMargin(0.15);
0199         AlignmentPI::makeNicePlotStyle(diffs[coord].get(), kBlack);
0200         auto max = diffs[coord]->GetMaximum();
0201         auto min = diffs[coord]->GetMinimum();
0202         auto range = std::abs(max) > std::abs(min) ? std::abs(max) : std::abs(min);
0203         if (range == 0.f)
0204           range = 0.1;
0205         //auto newMax = (max > 0.) ? max*1.2 : max*0.8;
0206 
0207         diffs[coord]->GetYaxis()->SetRangeUser(-range * 1.5, range * 1.5);
0208         diffs[coord]->GetYaxis()->SetTitleOffset(1.5);
0209         diffs[coord]->SetMarkerStyle(20);
0210         diffs[coord]->SetMarkerSize(0.5);
0211         diffs[coord]->Draw("P");
0212 
0213         if (subpad == 1) { /* fill the legend only at the first pass */
0214           if (this->m_plotAnnotations.ntags == 2) {
0215             legend.SetHeader("#bf{Two Tags Comparison}", "C");  // option "C" allows to center the header
0216             legend.AddEntry(
0217                 diffs[coord].get(),
0218                 ("#splitline{" + tagname1 + " : " + firstIOVsince + "}{" + tagname2 + " : " + lastIOVsince + "}")
0219                     .c_str(),
0220                 "PL");
0221           } else {
0222             legend.SetHeader(("tag: #bf{" + tagname1 + "}").c_str(), "C");  // option "C" allows to center the header
0223             legend.AddEntry(diffs[coord].get(),
0224                             ("#splitline{IOV since: " + firstIOVsince + "}{IOV since: " + lastIOVsince + "}").c_str(),
0225                             "PL");
0226           }
0227         }
0228         subpad++;
0229       }
0230 
0231       canvas.Update();
0232       canvas.cd();
0233       canvas.Modified();
0234 
0235       TLine l[6][boundaries.size()];
0236       TLatex tSubdet[6];
0237       for (unsigned int i = 0; i < 6; i++) {
0238         tSubdet[i].SetTextColor(kRed);
0239         tSubdet[i].SetNDC();
0240         tSubdet[i].SetTextAlign(21);
0241         tSubdet[i].SetTextSize(doOnlyPixel ? 0.05 : 0.03);
0242         tSubdet[i].SetTextAngle(90);
0243       }
0244 
0245       subpad = 0;
0246       for (const auto &coord : coords) {
0247         auto s_coord = AlignmentPI::getStringFromCoordinate(coord);
0248         canvas.cd(subpad + 1);
0249         for (const auto &line : boundaries | boost::adaptors::indexed(0)) {
0250           const auto &index = line.index();
0251           const auto value = line.value();
0252           l[subpad][index] = TLine(diffs[coord]->GetBinLowEdge(value.first),
0253                                    canvas.cd(subpad + 1)->GetUymin(),
0254                                    diffs[coord]->GetBinLowEdge(value.first),
0255                                    canvas.cd(subpad + 1)->GetUymax() * 0.84);
0256           l[subpad][index].SetLineWidth(1);
0257           l[subpad][index].SetLineStyle(9);
0258           l[subpad][index].SetLineColor(2);
0259           l[subpad][index].Draw("same");
0260         }
0261 
0262         for (const auto &elem : boundaries | boost::adaptors::indexed(0)) {
0263           const auto &lm = canvas.cd(subpad + 1)->GetLeftMargin();
0264           const auto &rm = 1 - canvas.cd(subpad + 1)->GetRightMargin();
0265           const auto &frac = float(elem.value().first) / ref_ali.size();
0266 
0267           LogDebug("TrackerAlignmentCompareAll")
0268               << __PRETTY_FUNCTION__ << " left margin:  " << lm << " right margin: " << rm << " fraction: " << frac;
0269 
0270           float theX_ = lm + (rm - lm) * frac + ((elem.index() > 0 || doOnlyPixel) ? 0.025 : 0.01);
0271 
0272           tSubdet[subpad].DrawLatex(
0273               theX_, 0.23, Form("%s", AlignmentPI::getStringFromPart(elem.value().second, /*is phase2?*/ ph2).c_str()));
0274         }
0275 
0276         auto ltx = TLatex();
0277         ltx.SetTextFont(62);
0278         ltx.SetTextSize(0.042);
0279         ltx.SetTextAlign(11);
0280         ltx.DrawLatexNDC(canvas.cd(subpad + 1)->GetLeftMargin(),
0281                          1 - canvas.cd(subpad + 1)->GetTopMargin() + 0.01,
0282                          ("Tracker Alignment Compare : #color[4]{" + s_coord + "}").c_str());
0283         legend.Draw("same");
0284         subpad++;
0285       }  // loop on the coordinates
0286 
0287       std::string fileName(this->m_imageFileName);
0288       canvas.SaveAs(fileName.c_str());
0289       //canvas.SaveAs("out.root");
0290 
0291       return true;
0292     }
0293   };
0294 
0295   typedef TrackerAlignmentCompareAll<1, MULTI_IOV, false> TrackerAlignmentComparatorSingleTag;
0296   typedef TrackerAlignmentCompareAll<2, SINGLE_IOV, false> TrackerAlignmentComparatorTwoTags;
0297 
0298   typedef TrackerAlignmentCompareAll<1, MULTI_IOV, true> PixelAlignmentComparatorSingleTag;
0299   typedef TrackerAlignmentCompareAll<2, SINGLE_IOV, true> PixelAlignmentComparatorTwoTags;
0300 
0301   //*******************************************//
0302   // Size of the movement over all partitions,
0303   // one coordinate (x,y,z,...) at a time
0304   //******************************************//
0305 
0306   template <AlignmentPI::coordinate coord, int ntags, IOVMultiplicity nIOVs>
0307   class TrackerAlignmentComparatorBase : public PlotImage<Alignments, nIOVs, ntags> {
0308   public:
0309     TrackerAlignmentComparatorBase()
0310         : PlotImage<Alignments, nIOVs, ntags>("comparison of " + AlignmentPI::getStringFromCoordinate(coord) +
0311                                               " coordinate between two geometries") {}
0312 
0313     bool fill() override {
0314       TGaxis::SetExponentOffset(-0.12, 0.01, "y");  // Y offset
0315 
0316       // trick to deal with the multi-ioved tag and two tag case at the same time
0317       auto theIOVs = PlotBase::getTag<0>().iovs;
0318       auto tagname1 = PlotBase::getTag<0>().name;
0319       std::string tagname2 = "";
0320       auto firstiov = theIOVs.front();
0321       std::tuple<cond::Time_t, cond::Hash> lastiov;
0322 
0323       // we don't support (yet) comparison with more than 2 tags
0324       assert(this->m_plotAnnotations.ntags < 3);
0325 
0326       if (this->m_plotAnnotations.ntags == 2) {
0327         auto tag2iovs = PlotBase::getTag<1>().iovs;
0328         tagname2 = PlotBase::getTag<1>().name;
0329         lastiov = tag2iovs.front();
0330       } else {
0331         lastiov = theIOVs.back();
0332       }
0333 
0334       std::shared_ptr<Alignments> last_payload = this->fetchPayload(std::get<1>(lastiov));
0335       std::shared_ptr<Alignments> first_payload = this->fetchPayload(std::get<1>(firstiov));
0336 
0337       std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
0338       std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
0339 
0340       std::vector<AlignTransform> ref_ali = first_payload->m_align;
0341       std::vector<AlignTransform> target_ali = last_payload->m_align;
0342 
0343       TCanvas canvas("Alignment Comparison", "Alignment Comparison", 1200, 1200);
0344 
0345       if (ref_ali.size() != target_ali.size()) {
0346         edm::LogError("TrackerAlignment_PayloadInspector")
0347             << "the size of the reference alignment (" << ref_ali.size()
0348             << ") is different from the one of the target (" << target_ali.size()
0349             << ")! You are probably trying to compare different underlying geometries. Exiting";
0350         return false;
0351       }
0352 
0353       // check that the geomtery is a tracker one
0354       const char *path_toTopologyXML = (ref_ali.size() == AlignmentPI::phase0size)
0355                                            ? "Geometry/TrackerCommonData/data/trackerParameters.xml"
0356                                            : "Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml";
0357       TrackerTopology tTopo =
0358           StandaloneTrackerTopology::fromTrackerParametersXMLFile(edm::FileInPath(path_toTopologyXML).fullPath());
0359 
0360       for (const auto &ali : ref_ali) {
0361         auto mydetid = ali.rawId();
0362         if (DetId(mydetid).det() != DetId::Tracker) {
0363           edm::LogWarning("TrackerAlignment_PayloadInspector")
0364               << "Encountered invalid Tracker DetId:" << DetId(mydetid).rawId() << " (" << DetId(mydetid).det()
0365               << ") is different from " << DetId::Tracker << " (is DoubleSide: " << tTopo.tidIsDoubleSide(mydetid)
0366               << "); subdetId " << DetId(mydetid).subdetId() << " - terminating ";
0367           return false;
0368         }
0369       }
0370 
0371       auto s_coord = AlignmentPI::getStringFromCoordinate(coord);
0372       std::string unit =
0373           (coord == AlignmentPI::t_x || coord == AlignmentPI::t_y || coord == AlignmentPI::t_z) ? "[#mum]" : "[mrad]";
0374 
0375       //std::unique_ptr<TH1F> compare = std::unique_ptr<TH1F>(new TH1F("comparison",Form("Comparison of %s;DetId index; #Delta%s %s",s_coord.c_str(),s_coord.c_str(),unit.c_str()),ref_ali.size(),-0.5,ref_ali.size()-0.5));
0376       std::unique_ptr<TH1F> compare =
0377           std::make_unique<TH1F>("comparison",
0378                                  Form(";Detector Id index; #Delta%s %s", s_coord.c_str(), unit.c_str()),
0379                                  ref_ali.size(),
0380                                  -0.5,
0381                                  ref_ali.size() - 0.5);
0382 
0383       // fill the histograms
0384       std::map<int, AlignmentPI::partitions> boundaries;
0385       boundaries.insert({0, AlignmentPI::BPix});  // always start with BPix, not filled in the loop
0386       AlignmentPI::fillComparisonHistogram(coord, boundaries, ref_ali, target_ali, compare);
0387 
0388       canvas.cd();
0389 
0390       canvas.SetTopMargin(0.06);
0391       canvas.SetLeftMargin(0.17);
0392       canvas.SetRightMargin(0.05);
0393       canvas.SetBottomMargin(0.15);
0394       AlignmentPI::makeNicePlotStyle(compare.get(), kBlack);
0395       auto max = compare->GetMaximum();
0396       auto min = compare->GetMinimum();
0397       auto range = std::abs(max) > std::abs(min) ? std::abs(max) : std::abs(min);
0398       if (range == 0.f)
0399         range = 0.1;
0400       //auto newMax = (max > 0.) ? max*1.2 : max*0.8;
0401 
0402       compare->GetYaxis()->SetRangeUser(-range * 1.5, range * 1.5);
0403       compare->GetYaxis()->SetTitleOffset(1.5);
0404       compare->SetMarkerStyle(20);
0405       compare->SetMarkerSize(0.5);
0406       compare->Draw("P");
0407 
0408       canvas.Update();
0409       canvas.cd();
0410 
0411       TLine l[boundaries.size()];
0412       for (const auto &line : boundaries | boost::adaptors::indexed(0)) {
0413         const auto &index = line.index();
0414         const auto value = line.value();
0415         l[index] = TLine(compare->GetBinLowEdge(value.first),
0416                          canvas.cd()->GetUymin(),
0417                          compare->GetBinLowEdge(value.first),
0418                          canvas.cd()->GetUymax());
0419         l[index].SetLineWidth(1);
0420         l[index].SetLineStyle(9);
0421         l[index].SetLineColor(2);
0422         l[index].Draw("same");
0423       }
0424 
0425       TLatex tSubdet;
0426       tSubdet.SetNDC();
0427       tSubdet.SetTextAlign(21);
0428       tSubdet.SetTextSize(0.027);
0429       tSubdet.SetTextAngle(90);
0430 
0431       for (const auto &elem : boundaries) {
0432         tSubdet.SetTextColor(kRed);
0433         auto myPair = AlignmentPI::calculatePosition(gPad, compare->GetBinLowEdge(elem.first));
0434         float theX_ = elem.first != 0 ? myPair.first + 0.025 : myPair.first + 0.01;
0435         const bool isPhase2 = (ref_ali.size() > AlignmentPI::phase1size);
0436         tSubdet.DrawLatex(theX_, 0.20, Form("%s", AlignmentPI::getStringFromPart(elem.second, isPhase2).c_str()));
0437       }
0438 
0439       TLegend legend = TLegend(0.17, 0.86, 0.95, 0.94);
0440       if (this->m_plotAnnotations.ntags == 2) {
0441         legend.SetHeader("#bf{Two Tags Comparison}", "C");  // option "C" allows to center the header
0442         legend.AddEntry(
0443             compare.get(),
0444             ("#splitline{" + tagname1 + " : " + firstIOVsince + "}{" + tagname2 + " : " + lastIOVsince + "}").c_str(),
0445             "PL");
0446       } else {
0447         legend.SetHeader(("tag: #bf{" + tagname1 + "}").c_str(), "C");  // option "C" allows to center the header
0448         legend.AddEntry(compare.get(),
0449                         ("#splitline{IOV since: " + firstIOVsince + "}{IOV since: " + lastIOVsince + "}").c_str(),
0450                         "PL");
0451       }
0452       legend.SetTextSize(0.020);
0453       legend.Draw("same");
0454 
0455       auto ltx = TLatex();
0456       ltx.SetTextFont(62);
0457       ltx.SetTextSize(0.042);
0458       ltx.SetTextAlign(11);
0459       ltx.DrawLatexNDC(gPad->GetLeftMargin(),
0460                        1 - gPad->GetTopMargin() + 0.01,
0461                        ("Tracker Alignment Compare :#color[4]{" + s_coord + "}").c_str());
0462 
0463       std::string fileName(this->m_imageFileName);
0464       canvas.SaveAs(fileName.c_str());
0465 
0466       return true;
0467     }
0468   };
0469 
0470   template <AlignmentPI::coordinate coord>
0471   using TrackerAlignmentCompare = TrackerAlignmentComparatorBase<coord, 1, MULTI_IOV>;
0472 
0473   template <AlignmentPI::coordinate coord>
0474   using TrackerAlignmentCompareTwoTags = TrackerAlignmentComparatorBase<coord, 2, SINGLE_IOV>;
0475 
0476   typedef TrackerAlignmentCompare<AlignmentPI::t_x> TrackerAlignmentCompareX;
0477   typedef TrackerAlignmentCompare<AlignmentPI::t_y> TrackerAlignmentCompareY;
0478   typedef TrackerAlignmentCompare<AlignmentPI::t_z> TrackerAlignmentCompareZ;
0479 
0480   typedef TrackerAlignmentCompare<AlignmentPI::rot_alpha> TrackerAlignmentCompareAlpha;
0481   typedef TrackerAlignmentCompare<AlignmentPI::rot_beta> TrackerAlignmentCompareBeta;
0482   typedef TrackerAlignmentCompare<AlignmentPI::rot_gamma> TrackerAlignmentCompareGamma;
0483 
0484   typedef TrackerAlignmentCompareTwoTags<AlignmentPI::t_x> TrackerAlignmentCompareXTwoTags;
0485   typedef TrackerAlignmentCompareTwoTags<AlignmentPI::t_y> TrackerAlignmentCompareYTwoTags;
0486   typedef TrackerAlignmentCompareTwoTags<AlignmentPI::t_z> TrackerAlignmentCompareZTwoTags;
0487 
0488   typedef TrackerAlignmentCompareTwoTags<AlignmentPI::rot_alpha> TrackerAlignmentCompareAlphaTwoTags;
0489   typedef TrackerAlignmentCompareTwoTags<AlignmentPI::rot_beta> TrackerAlignmentCompareBetaTwoTags;
0490   typedef TrackerAlignmentCompareTwoTags<AlignmentPI::rot_gamma> TrackerAlignmentCompareGammaTwoTags;
0491 
0492   //*******************************************//
0493   // Summary canvas per subdetector
0494   //******************************************//
0495 
0496   template <int ntags, IOVMultiplicity nIOVs, AlignmentPI::partitions q>
0497   class TrackerAlignmentSummaryBase : public PlotImage<Alignments, nIOVs, ntags> {
0498   public:
0499     TrackerAlignmentSummaryBase()
0500         : PlotImage<Alignments, nIOVs, ntags>("Comparison of all coordinates between two geometries for " +
0501                                               getStringFromPart(q)) {}
0502 
0503     bool fill() override {
0504       // trick to deal with the multi-ioved tag and two tag case at the same time
0505       auto theIOVs = PlotBase::getTag<0>().iovs;
0506       auto tagname1 = PlotBase::getTag<0>().name;
0507       std::string tagname2 = "";
0508       auto firstiov = theIOVs.front();
0509       std::tuple<cond::Time_t, cond::Hash> lastiov;
0510 
0511       // we don't support (yet) comparison with more than 2 tags
0512       assert(this->m_plotAnnotations.ntags < 3);
0513 
0514       if (this->m_plotAnnotations.ntags == 2) {
0515         auto tag2iovs = PlotBase::getTag<1>().iovs;
0516         tagname2 = PlotBase::getTag<1>().name;
0517         lastiov = tag2iovs.front();
0518       } else {
0519         lastiov = theIOVs.back();
0520       }
0521 
0522       std::shared_ptr<Alignments> last_payload = this->fetchPayload(std::get<1>(lastiov));
0523       std::shared_ptr<Alignments> first_payload = this->fetchPayload(std::get<1>(firstiov));
0524 
0525       std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
0526       std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
0527 
0528       std::vector<AlignTransform> ref_ali = first_payload->m_align;
0529       std::vector<AlignTransform> target_ali = last_payload->m_align;
0530 
0531       if (ref_ali.size() != target_ali.size()) {
0532         edm::LogError("TrackerAlignment_PayloadInspector")
0533             << "the size of the reference alignment (" << ref_ali.size()
0534             << ") is different from the one of the target (" << target_ali.size()
0535             << ")! You are probably trying to compare different underlying geometries. Exiting";
0536         return false;
0537       }
0538 
0539       // check that the geomtery is a tracker one
0540       const char *path_toTopologyXML = (ref_ali.size() == AlignmentPI::phase0size)
0541                                            ? "Geometry/TrackerCommonData/data/trackerParameters.xml"
0542                                            : "Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml";
0543       TrackerTopology tTopo =
0544           StandaloneTrackerTopology::fromTrackerParametersXMLFile(edm::FileInPath(path_toTopologyXML).fullPath());
0545 
0546       for (const auto &ali : ref_ali) {
0547         auto mydetid = ali.rawId();
0548         if (DetId(mydetid).det() != DetId::Tracker) {
0549           edm::LogWarning("TrackerAlignment_PayloadInspector")
0550               << "Encountered invalid Tracker DetId:" << DetId(mydetid).rawId() << " (" << DetId(mydetid).det()
0551               << ") is different from " << DetId::Tracker << " (is DoubleSide: " << tTopo.tidIsDoubleSide(mydetid)
0552               << "); subdetId " << DetId(mydetid).subdetId() << " - terminating ";
0553           return false;
0554         }
0555       }
0556 
0557       TCanvas canvas("Alignment Comparison", "Alignment Comparison", 1800, 1200);
0558       canvas.Divide(3, 2);
0559 
0560       std::unordered_map<AlignmentPI::coordinate, std::unique_ptr<TH1F> > diffs;
0561       std::vector<AlignmentPI::coordinate> coords = {AlignmentPI::t_x,
0562                                                      AlignmentPI::t_y,
0563                                                      AlignmentPI::t_z,
0564                                                      AlignmentPI::rot_alpha,
0565                                                      AlignmentPI::rot_beta,
0566                                                      AlignmentPI::rot_gamma};
0567 
0568       for (const auto &coord : coords) {
0569         auto s_coord = AlignmentPI::getStringFromCoordinate(coord);
0570         std::string unit =
0571             (coord == AlignmentPI::t_x || coord == AlignmentPI::t_y || coord == AlignmentPI::t_z) ? "[#mum]" : "[mrad]";
0572 
0573         diffs[coord] = std::make_unique<TH1F>(Form("hDiff_%s", s_coord.c_str()),
0574                                               Form(";#Delta%s %s;n. of modules", s_coord.c_str(), unit.c_str()),
0575                                               1001,
0576                                               -500.5,
0577                                               500.5);
0578       }
0579 
0580       // fill the comparison histograms
0581       std::map<int, AlignmentPI::partitions> boundaries;
0582       AlignmentPI::fillComparisonHistograms(boundaries, ref_ali, target_ali, diffs, true, q);
0583 
0584       int c_index = 1;
0585 
0586       //TLegend (Double_t x1, Double_t y1, Double_t x2, Double_t y2, const char *header="", Option_t *option="brNDC")
0587       auto legend = std::make_unique<TLegend>(0.14, 0.88, 0.96, 0.99);
0588       if (this->m_plotAnnotations.ntags == 2) {
0589         legend->SetHeader("#bf{Two Tags Comparison}", "C");  // option "C" allows to center the header
0590         legend->AddEntry(
0591             diffs[AlignmentPI::t_x].get(),
0592             ("#splitline{" + tagname1 + " : " + firstIOVsince + "}{" + tagname2 + " : " + lastIOVsince + "}").c_str(),
0593             "PL");
0594       } else {
0595         legend->SetHeader(("tag: #bf{" + tagname1 + "}").c_str(), "C");  // option "C" allows to center the header
0596         legend->AddEntry(diffs[AlignmentPI::t_x].get(),
0597                          ("#splitline{IOV since: " + firstIOVsince + "}{IOV since: " + lastIOVsince + "}").c_str(),
0598                          "PL");
0599       }
0600       legend->SetTextSize(0.025);
0601 
0602       for (const auto &coord : coords) {
0603         canvas.cd(c_index)->SetLogy();
0604         canvas.cd(c_index)->SetTopMargin(0.01);
0605         canvas.cd(c_index)->SetBottomMargin(0.15);
0606         canvas.cd(c_index)->SetLeftMargin(0.14);
0607         canvas.cd(c_index)->SetRightMargin(0.04);
0608         diffs[coord]->SetLineWidth(2);
0609         AlignmentPI::makeNicePlotStyle(diffs[coord].get(), kBlack);
0610 
0611         //float x_max = diffs[coord]->GetXaxis()->GetBinCenter(diffs[coord]->FindLastBinAbove(0.));
0612         //float x_min = diffs[coord]->GetXaxis()->GetBinCenter(diffs[coord]->FindFirstBinAbove(0.));
0613         //float extremum = std::abs(x_max) > std::abs(x_min) ? std::abs(x_max) : std::abs(x_min);
0614         //diffs[coord]->GetXaxis()->SetRangeUser(-extremum*2,extremum*2);
0615 
0616         int i_max = diffs[coord]->FindLastBinAbove(0.);
0617         int i_min = diffs[coord]->FindFirstBinAbove(0.);
0618         diffs[coord]->GetXaxis()->SetRange(std::max(1, i_min - 10), std::min(i_max + 10, diffs[coord]->GetNbinsX()));
0619         diffs[coord]->SetMaximum(diffs[coord]->GetMaximum() * 5);
0620         diffs[coord]->Draw("HIST");
0621         AlignmentPI::makeNiceStats(diffs[coord].get(), q, kBlack);
0622 
0623         legend->Draw("same");
0624 
0625         c_index++;
0626       }
0627 
0628       std::string fileName(this->m_imageFileName);
0629       canvas.SaveAs(fileName.c_str());
0630 
0631       return true;
0632     }
0633   };
0634 
0635   typedef TrackerAlignmentSummaryBase<1, MULTI_IOV, AlignmentPI::BPix> TrackerAlignmentSummaryBPix;
0636   typedef TrackerAlignmentSummaryBase<1, MULTI_IOV, AlignmentPI::FPix> TrackerAlignmentSummaryFPix;
0637   typedef TrackerAlignmentSummaryBase<1, MULTI_IOV, AlignmentPI::TIB> TrackerAlignmentSummaryTIB;
0638 
0639   typedef TrackerAlignmentSummaryBase<1, MULTI_IOV, AlignmentPI::TID> TrackerAlignmentSummaryTID;
0640   typedef TrackerAlignmentSummaryBase<1, MULTI_IOV, AlignmentPI::TOB> TrackerAlignmentSummaryTOB;
0641   typedef TrackerAlignmentSummaryBase<1, MULTI_IOV, AlignmentPI::TEC> TrackerAlignmentSummaryTEC;
0642 
0643   typedef TrackerAlignmentSummaryBase<2, SINGLE_IOV, AlignmentPI::BPix> TrackerAlignmentSummaryBPixTwoTags;
0644   typedef TrackerAlignmentSummaryBase<2, SINGLE_IOV, AlignmentPI::FPix> TrackerAlignmentSummaryFPixTwoTags;
0645   typedef TrackerAlignmentSummaryBase<2, SINGLE_IOV, AlignmentPI::TIB> TrackerAlignmentSummaryTIBTwoTags;
0646 
0647   typedef TrackerAlignmentSummaryBase<2, SINGLE_IOV, AlignmentPI::TID> TrackerAlignmentSummaryTIDTwoTags;
0648   typedef TrackerAlignmentSummaryBase<2, SINGLE_IOV, AlignmentPI::TOB> TrackerAlignmentSummaryTOBTwoTags;
0649   typedef TrackerAlignmentSummaryBase<2, SINGLE_IOV, AlignmentPI::TEC> TrackerAlignmentSummaryTECTwoTags;
0650 
0651   /************************************************
0652    Full Pixel Tracker Map Comparison of coordinates
0653   *************************************************/
0654   template <AlignmentPI::coordinate coord, int ntags, IOVMultiplicity nIOVs>
0655   class PixelAlignmentComparisonMapBase : public PlotImage<Alignments, nIOVs, ntags> {
0656   public:
0657     PixelAlignmentComparisonMapBase()
0658         : PlotImage<Alignments, nIOVs, ntags>("SiPixel Comparison Map of " +
0659                                               AlignmentPI::getStringFromCoordinate(coord)) {
0660       label_ = "PixelAlignmentComparisonMap" + AlignmentPI::getStringFromCoordinate(coord);
0661       payloadString = "Tracker Alignment";
0662     }
0663 
0664     bool fill() override {
0665       gStyle->SetPalette(1);
0666 
0667       // trick to deal with the multi-ioved tag and two tag case at the same time
0668       auto theIOVs = PlotBase::getTag<0>().iovs;
0669       auto tagname1 = PlotBase::getTag<0>().name;
0670       std::string tagname2 = "";
0671       auto firstiov = theIOVs.front();
0672       std::tuple<cond::Time_t, cond::Hash> lastiov;
0673 
0674       // we don't support (yet) comparison with more than 2 tags
0675       assert(this->m_plotAnnotations.ntags < 3);
0676 
0677       if (this->m_plotAnnotations.ntags == 2) {
0678         auto tag2iovs = PlotBase::getTag<1>().iovs;
0679         tagname2 = PlotBase::getTag<1>().name;
0680         lastiov = tag2iovs.front();
0681       } else {
0682         lastiov = theIOVs.back();
0683       }
0684 
0685       std::shared_ptr<Alignments> last_payload = this->fetchPayload(std::get<1>(lastiov));
0686       std::shared_ptr<Alignments> first_payload = this->fetchPayload(std::get<1>(firstiov));
0687 
0688       std::string lastIOVsince = std::to_string(std::get<0>(lastiov));
0689       std::string firstIOVsince = std::to_string(std::get<0>(firstiov));
0690 
0691       const std::vector<AlignTransform> &ref_ali = first_payload->m_align;
0692       const std::vector<AlignTransform> &target_ali = last_payload->m_align;
0693 
0694       if (last_payload.get() && first_payload.get()) {
0695         Phase1PixelSummaryMap fullMap(
0696             "",
0697             fmt::sprintf("%s %s", payloadString, AlignmentPI::getStringFromCoordinate(coord)),
0698             fmt::sprintf(
0699                 "#Delta %s [%s]", AlignmentPI::getStringFromCoordinate(coord), (coord <= 3) ? "#mum" : "mrad"));
0700         fullMap.createTrackerBaseMap();
0701 
0702         if (this->isPhase0(ref_ali) || this->isPhase0(target_ali)) {
0703           edm::LogError(label_) << "Pixel Tracker Alignment maps are not supported for non-Phase1 Pixel geometries !";
0704           TCanvas canvas("Canv", "Canv", 1200, 1000);
0705           AlignmentPI::displayNotSupported(canvas, 0);
0706           std::string fileName(this->m_imageFileName);
0707           canvas.SaveAs(fileName.c_str());
0708           return false;
0709         }
0710 
0711         // fill the map of differences
0712         std::map<uint32_t, double> diffPerDetid;
0713         this->fillPerDetIdDiff(coord, ref_ali, target_ali, diffPerDetid);
0714 
0715         // now fill the tracker map
0716         for (const auto &elem : diffPerDetid) {
0717           // reject Strips
0718           int subid = DetId(elem.first).subdetId();
0719           if (subid > 2) {
0720             continue;
0721           }
0722           fullMap.fillTrackerMap(elem.first, elem.second);
0723         }
0724 
0725         // limit the axis range (in case of need)
0726         //fullMap.setZAxisRange(-50.f,50.f);
0727 
0728         TCanvas canvas("Canv", "Canv", 3000, 2000);
0729         fullMap.printTrackerMap(canvas);
0730 
0731         auto ltx = TLatex();
0732         ltx.SetTextFont(62);
0733         ltx.SetTextSize(0.025);
0734         ltx.SetTextAlign(11);
0735 
0736         ltx.DrawLatexNDC(gPad->GetLeftMargin() + 0.01,
0737                          gPad->GetBottomMargin() + 0.01,
0738                          ("#color[4]{" + tagname1 + "}, IOV: #color[4]{" + firstIOVsince + "} vs #color[4]{" +
0739                           tagname2 + "}, IOV: #color[4]{" + lastIOVsince + "}")
0740                              .c_str());
0741 
0742         std::string fileName(this->m_imageFileName);
0743         canvas.SaveAs(fileName.c_str());
0744       }
0745       return true;
0746     }
0747 
0748   protected:
0749     std::string payloadString;
0750     std::string label_;
0751 
0752   private:
0753     //_________________________________________________
0754     bool isPhase0(std::vector<AlignTransform> theAlis) {
0755       SiPixelDetInfoFileReader reader =
0756           SiPixelDetInfoFileReader(edm::FileInPath(SiPixelDetInfoFileReader::kPh0DefaultFile).fullPath());
0757       const auto &p0detIds = reader.getAllDetIds();
0758 
0759       std::vector<uint32_t> ownDetIds;
0760       std::transform(theAlis.begin(), theAlis.end(), std::back_inserter(ownDetIds), [](AlignTransform ali) -> uint32_t {
0761         return ali.rawId();
0762       });
0763 
0764       for (const auto &det : ownDetIds) {
0765         // if found at least one phase-0 detId early return
0766         if (std::find(p0detIds.begin(), p0detIds.end(), det) != p0detIds.end()) {
0767           return true;
0768         }
0769       }
0770       return false;
0771     }
0772 
0773     /*--------------------------------------------------------------------*/
0774     void fillPerDetIdDiff(const AlignmentPI::coordinate &myCoord,
0775                           const std::vector<AlignTransform> &ref_ali,
0776                           const std::vector<AlignTransform> &target_ali,
0777                           std::map<uint32_t, double> &diff)
0778     /*--------------------------------------------------------------------*/
0779     {
0780       for (unsigned int i = 0; i < ref_ali.size(); i++) {
0781         uint32_t detid = ref_ali[i].rawId();
0782         if (ref_ali[i].rawId() == target_ali[i].rawId()) {
0783           CLHEP::HepRotation target_rot(target_ali[i].rotation());
0784           CLHEP::HepRotation ref_rot(ref_ali[i].rotation());
0785 
0786           align::RotationType target_ROT(target_rot.xx(),
0787                                          target_rot.xy(),
0788                                          target_rot.xz(),
0789                                          target_rot.yx(),
0790                                          target_rot.yy(),
0791                                          target_rot.yz(),
0792                                          target_rot.zx(),
0793                                          target_rot.zy(),
0794                                          target_rot.zz());
0795 
0796           align::RotationType ref_ROT(ref_rot.xx(),
0797                                       ref_rot.xy(),
0798                                       ref_rot.xz(),
0799                                       ref_rot.yx(),
0800                                       ref_rot.yy(),
0801                                       ref_rot.yz(),
0802                                       ref_rot.zx(),
0803                                       ref_rot.zy(),
0804                                       ref_rot.zz());
0805 
0806           const std::vector<double> deltaRot = {
0807               ::deltaPhi(align::toAngles(target_ROT)[0], align::toAngles(ref_ROT)[0]),
0808               ::deltaPhi(align::toAngles(target_ROT)[1], align::toAngles(ref_ROT)[1]),
0809               ::deltaPhi(align::toAngles(target_ROT)[2], align::toAngles(ref_ROT)[2])};
0810 
0811           const auto &deltaTrans = target_ali[i].translation() - ref_ali[i].translation();
0812 
0813           switch (myCoord) {
0814             case AlignmentPI::t_x:
0815               diff.insert({detid, deltaTrans.x() * AlignmentPI::cmToUm});
0816               break;
0817             case AlignmentPI::t_y:
0818               diff.insert({detid, deltaTrans.y() * AlignmentPI::cmToUm});
0819               break;
0820             case AlignmentPI::t_z:
0821               diff.insert({detid, deltaTrans.z() * AlignmentPI::cmToUm});
0822               break;
0823             case AlignmentPI::rot_alpha:
0824               diff.insert({detid, deltaRot[0] * AlignmentPI::tomRad});
0825               break;
0826             case AlignmentPI::rot_beta:
0827               diff.insert({detid, deltaRot[1] * AlignmentPI::tomRad});
0828               break;
0829             case AlignmentPI::rot_gamma:
0830               diff.insert({detid, deltaRot[2] * AlignmentPI::tomRad});
0831               break;
0832             default:
0833               edm::LogError("TrackerAlignment_PayloadInspector") << "Unrecognized coordinate " << myCoord << std::endl;
0834               break;
0835           }  // switch on the coordinate
0836         }    // check on the same detID
0837       }      // loop on the components
0838     }
0839   };
0840 
0841   template <AlignmentPI::coordinate coord>
0842   using PixelAlignmentCompareMap = PixelAlignmentComparisonMapBase<coord, 1, MULTI_IOV>;
0843 
0844   template <AlignmentPI::coordinate coord>
0845   using PixelAlignmentCompareMapTwoTags = PixelAlignmentComparisonMapBase<coord, 2, SINGLE_IOV>;
0846 
0847   typedef PixelAlignmentCompareMap<AlignmentPI::t_x> PixelAlignmentCompareMapX;
0848   typedef PixelAlignmentCompareMap<AlignmentPI::t_y> PixelAlignmentCompareMapY;
0849   typedef PixelAlignmentCompareMap<AlignmentPI::t_z> PixelAlignmentCompareMapZ;
0850 
0851   typedef PixelAlignmentCompareMap<AlignmentPI::rot_alpha> PixelAlignmentCompareMapAlpha;
0852   typedef PixelAlignmentCompareMap<AlignmentPI::rot_beta> PixelAlignmentCompareMapBeta;
0853   typedef PixelAlignmentCompareMap<AlignmentPI::rot_gamma> PixelAlignmentCompareMapGamma;
0854 
0855   typedef PixelAlignmentCompareMapTwoTags<AlignmentPI::t_x> PixelAlignmentCompareMapXTwoTags;
0856   typedef PixelAlignmentCompareMapTwoTags<AlignmentPI::t_y> PixelAlignmentCompareMapYTwoTags;
0857   typedef PixelAlignmentCompareMapTwoTags<AlignmentPI::t_z> PixelAlignmentCompareMapZTwoTags;
0858 
0859   typedef PixelAlignmentCompareMapTwoTags<AlignmentPI::rot_alpha> PixelAlignmentCompareMapAlphaTwoTags;
0860   typedef PixelAlignmentCompareMapTwoTags<AlignmentPI::rot_beta> PixelAlignmentCompareMapBetaTwoTags;
0861   typedef PixelAlignmentCompareMapTwoTags<AlignmentPI::rot_gamma> PixelAlignmentCompareMapGammaTwoTags;
0862 
0863   //*******************************************//
0864   // History of the position of the BPix Barycenter
0865   //******************************************//
0866 
0867   template <AlignmentPI::coordinate coord>
0868   class BPixBarycenterHistory : public HistoryPlot<Alignments, float> {
0869   public:
0870     BPixBarycenterHistory()
0871         : HistoryPlot<Alignments, float>(
0872               " Barrel Pixel " + AlignmentPI::getStringFromCoordinate(coord) + " positions vs time",
0873               AlignmentPI::getStringFromCoordinate(coord) + " position [cm]") {}
0874     ~BPixBarycenterHistory() override = default;
0875 
0876     float getFromPayload(Alignments &payload) override {
0877       std::vector<AlignTransform> alignments = payload.m_align;
0878 
0879       float barycenter = 0.;
0880       float nmodules(0.);
0881       for (const auto &ali : alignments) {
0882         if (DetId(ali.rawId()).det() != DetId::Tracker) {
0883           edm::LogWarning("TrackerAlignment_PayloadInspector")
0884               << "Encountered invalid Tracker DetId:" << ali.rawId() << " " << DetId(ali.rawId()).det()
0885               << " is different from " << DetId::Tracker << "  - terminating ";
0886           return false;
0887         }
0888 
0889         int subid = DetId(ali.rawId()).subdetId();
0890         if (subid != PixelSubdetector::PixelBarrel)
0891           continue;
0892 
0893         nmodules++;
0894         switch (coord) {
0895           case AlignmentPI::t_x:
0896             barycenter += (ali.translation().x());
0897             break;
0898           case AlignmentPI::t_y:
0899             barycenter += (ali.translation().y());
0900             break;
0901           case AlignmentPI::t_z:
0902             barycenter += (ali.translation().z());
0903             break;
0904           default:
0905             edm::LogError("TrackerAlignment_PayloadInspector") << "Unrecognized coordinate " << coord << std::endl;
0906             break;
0907         }  // switch on the coordinate (only X,Y,Z are interesting)
0908       }    // ends loop on the alignments
0909 
0910       edm::LogInfo("TrackerAlignment_PayloadInspector") << "barycenter (" << barycenter << ")/n. modules (" << nmodules
0911                                                         << ") =  " << barycenter / nmodules << std::endl;
0912 
0913       // take the mean
0914       barycenter /= nmodules;
0915 
0916       // applied GPR correction to move barycenter to global CMS coordinates
0917       barycenter += hardcodeGPR.at(coord);
0918 
0919       return barycenter;
0920 
0921     }  // payload
0922   };
0923 
0924   typedef BPixBarycenterHistory<AlignmentPI::t_x> X_BPixBarycenterHistory;
0925   typedef BPixBarycenterHistory<AlignmentPI::t_y> Y_BPixBarycenterHistory;
0926   typedef BPixBarycenterHistory<AlignmentPI::t_z> Z_BPixBarycenterHistory;
0927 
0928   /************************************************
0929     Display of Tracker Detector barycenters
0930   *************************************************/
0931   class TrackerAlignmentBarycenters : public PlotImage<Alignments, SINGLE_IOV> {
0932   public:
0933     TrackerAlignmentBarycenters() : PlotImage<Alignments, SINGLE_IOV>("Display of Tracker Alignment Barycenters") {}
0934 
0935     bool fill() override {
0936       auto tag = PlotBase::getTag<0>();
0937       auto iov = tag.iovs.front();
0938       const auto &tagname = PlotBase::getTag<0>().name;
0939       std::shared_ptr<Alignments> payload = fetchPayload(std::get<1>(iov));
0940       unsigned int run = std::get<0>(iov);
0941 
0942       TCanvas canvas("Tracker Alignment Barycenter Summary", "Tracker Alignment Barycenter summary", 1600, 1000);
0943       canvas.cd();
0944 
0945       canvas.SetTopMargin(0.07);
0946       canvas.SetBottomMargin(0.06);
0947       canvas.SetLeftMargin(0.15);
0948       canvas.SetRightMargin(0.03);
0949       canvas.Modified();
0950       canvas.SetGrid();
0951 
0952       auto h2_BarycenterParameters =
0953           std::make_unique<TH2F>("Parameters", "SubDetector Barycenter summary", 6, 0.0, 6.0, 6, 0, 6.);
0954 
0955       auto h2_uncBarycenterParameters =
0956           std::make_unique<TH2F>("Parameters2", "SubDetector Barycenter summary", 6, 0.0, 6.0, 6, 0, 6.);
0957 
0958       h2_BarycenterParameters->SetStats(false);
0959       h2_BarycenterParameters->SetTitle(nullptr);
0960       h2_uncBarycenterParameters->SetStats(false);
0961       h2_uncBarycenterParameters->SetTitle(nullptr);
0962 
0963       std::vector<AlignTransform> alignments = payload->m_align;
0964 
0965       isPhase0 = (alignments.size() == AlignmentPI::phase0size) ? true : false;
0966 
0967       // check that the geomtery is a tracker one
0968       const char *path_toTopologyXML = isPhase0 ? "Geometry/TrackerCommonData/data/trackerParameters.xml"
0969                                                 : "Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml";
0970 
0971       TrackerTopology tTopo =
0972           StandaloneTrackerTopology::fromTrackerParametersXMLFile(edm::FileInPath(path_toTopologyXML).fullPath());
0973 
0974       AlignmentPI::TkAlBarycenters barycenters;
0975       // compute uncorrected barycenter
0976       barycenters.computeBarycenters(
0977           alignments, tTopo, {{AlignmentPI::t_x, 0.0}, {AlignmentPI::t_y, 0.0}, {AlignmentPI::t_z, 0.0}});
0978 
0979       auto Xbarycenters = barycenters.getX();
0980       auto Ybarycenters = barycenters.getY();
0981       auto Zbarycenters = barycenters.getZ();
0982 
0983       // compute barycenter corrected for the GPR
0984       barycenters.computeBarycenters(alignments, tTopo, hardcodeGPR);
0985 
0986       auto c_Xbarycenters = barycenters.getX();
0987       auto c_Ybarycenters = barycenters.getY();
0988       auto c_Zbarycenters = barycenters.getZ();
0989 
0990       h2_BarycenterParameters->GetXaxis()->SetBinLabel(1, "X [cm]");
0991       h2_BarycenterParameters->GetXaxis()->SetBinLabel(2, "Y [cm]");
0992       h2_BarycenterParameters->GetXaxis()->SetBinLabel(3, "Z [cm]");
0993       h2_BarycenterParameters->GetXaxis()->SetBinLabel(4, "X_{no GPR} [cm]");
0994       h2_BarycenterParameters->GetXaxis()->SetBinLabel(5, "Y_{no GPR} [cm]");
0995       h2_BarycenterParameters->GetXaxis()->SetBinLabel(6, "Z_{no GPR} [cm]");
0996 
0997       bool isLikelyMC(false);
0998       int checkX =
0999           std::count_if(Xbarycenters.begin(), Xbarycenters.begin() + 2, [](float a) { return (std::abs(a) >= 1.e-4); });
1000       int checkY =
1001           std::count_if(Ybarycenters.begin(), Ybarycenters.begin() + 2, [](float a) { return (std::abs(a) >= 1.e-4); });
1002       int checkZ =
1003           std::count_if(Zbarycenters.begin(), Zbarycenters.begin() + 2, [](float a) { return (std::abs(a) >= 1.e-4); });
1004 
1005       // if all the coordinate barycenters for both BPix and FPix are below 10um
1006       // this is very likely a MC payload
1007       if ((checkX + checkY + checkZ) == 0 && run == 1)
1008         isLikelyMC = true;
1009 
1010       unsigned int yBin = 6;
1011       for (unsigned int i = 0; i < 6; i++) {
1012         auto thePart = static_cast<AlignmentPI::partitions>(i + 1);
1013         std::string theLabel = getStringFromPart(thePart);
1014         h2_BarycenterParameters->GetYaxis()->SetBinLabel(yBin, theLabel.c_str());
1015         if (!isLikelyMC) {
1016           h2_BarycenterParameters->SetBinContent(1, yBin, c_Xbarycenters[i]);
1017           h2_BarycenterParameters->SetBinContent(2, yBin, c_Ybarycenters[i]);
1018           h2_BarycenterParameters->SetBinContent(3, yBin, c_Zbarycenters[i]);
1019         }
1020 
1021         h2_uncBarycenterParameters->SetBinContent(4, yBin, Xbarycenters[i]);
1022         h2_uncBarycenterParameters->SetBinContent(5, yBin, Ybarycenters[i]);
1023         h2_uncBarycenterParameters->SetBinContent(6, yBin, Zbarycenters[i]);
1024         yBin--;
1025       }
1026 
1027       h2_BarycenterParameters->GetXaxis()->LabelsOption("h");
1028       h2_BarycenterParameters->GetYaxis()->SetLabelSize(0.05);
1029       h2_BarycenterParameters->GetXaxis()->SetLabelSize(0.05);
1030       h2_BarycenterParameters->SetMarkerSize(1.5);
1031       h2_BarycenterParameters->Draw("TEXT");
1032 
1033       h2_uncBarycenterParameters->SetMarkerSize(1.5);
1034       h2_uncBarycenterParameters->SetMarkerColor(kRed);
1035       h2_uncBarycenterParameters->Draw("TEXTsame");
1036 
1037       TLatex t1;
1038       t1.SetNDC();
1039       t1.SetTextAlign(26);
1040       t1.SetTextSize(0.045);
1041       t1.DrawLatex(0.5, 0.96, Form("TkAl Barycenters, Tag: #color[4]{%s}, IOV #color[4]{%i}", tagname.c_str(), run));
1042       t1.SetTextSize(0.025);
1043 
1044       std::string fileName(m_imageFileName);
1045       canvas.SaveAs(fileName.c_str());
1046 
1047       return true;
1048     }
1049 
1050   private:
1051     bool isPhase0;
1052   };
1053 
1054   /************************************************
1055     Comparator of Tracker Detector barycenters
1056   *************************************************/
1057   template <int ntags, IOVMultiplicity nIOVs>
1058   class TrackerAlignmentBarycentersComparatorBase : public PlotImage<Alignments, nIOVs, ntags> {
1059   public:
1060     TrackerAlignmentBarycentersComparatorBase()
1061         : PlotImage<Alignments, nIOVs, ntags>("Comparison of Tracker Alignment Barycenters") {}
1062 
1063     bool fill() override {
1064       // trick to deal with the multi-ioved tag and two tag case at the same time
1065       auto theIOVs = PlotBase::getTag<0>().iovs;
1066       auto tagname1 = PlotBase::getTag<0>().name;
1067       std::string tagname2 = "";
1068       auto firstiov = theIOVs.front();
1069       std::tuple<cond::Time_t, cond::Hash> lastiov;
1070 
1071       // we don't support (yet) comparison with more than 2 tags
1072       assert(this->m_plotAnnotations.ntags < 3);
1073 
1074       if (this->m_plotAnnotations.ntags == 2) {
1075         auto tag2iovs = PlotBase::getTag<1>().iovs;
1076         tagname2 = PlotBase::getTag<1>().name;
1077         lastiov = tag2iovs.front();
1078       } else {
1079         lastiov = theIOVs.back();
1080       }
1081 
1082       unsigned int first_run = std::get<0>(firstiov);
1083       unsigned int last_run = std::get<0>(lastiov);
1084 
1085       std::shared_ptr<Alignments> last_payload = this->fetchPayload(std::get<1>(lastiov));
1086       std::vector<AlignTransform> last_alignments = last_payload->m_align;
1087 
1088       std::shared_ptr<Alignments> first_payload = this->fetchPayload(std::get<1>(firstiov));
1089       std::vector<AlignTransform> first_alignments = first_payload->m_align;
1090 
1091       isInitialPhase0 = (first_alignments.size() == AlignmentPI::phase0size) ? true : false;
1092       isFinalPhase0 = (last_alignments.size() == AlignmentPI::phase0size) ? true : false;
1093 
1094       // check that the geomtery is a tracker one
1095       const char *path_toTopologyXML = isInitialPhase0 ? "Geometry/TrackerCommonData/data/trackerParameters.xml"
1096                                                        : "Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml";
1097 
1098       TrackerTopology tTopo_f =
1099           StandaloneTrackerTopology::fromTrackerParametersXMLFile(edm::FileInPath(path_toTopologyXML).fullPath());
1100 
1101       path_toTopologyXML = isFinalPhase0 ? "Geometry/TrackerCommonData/data/trackerParameters.xml"
1102                                          : "Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml";
1103 
1104       TrackerTopology tTopo_l =
1105           StandaloneTrackerTopology::fromTrackerParametersXMLFile(edm::FileInPath(path_toTopologyXML).fullPath());
1106 
1107       TCanvas canvas("Tracker Alignment Barycenter Summary", "Tracker Alignment Barycenter summary", 1200, 800);
1108       canvas.cd();
1109 
1110       canvas.SetTopMargin(0.07);
1111       canvas.SetBottomMargin(0.06);
1112       canvas.SetLeftMargin(0.15);
1113       canvas.SetRightMargin(0.03);
1114       canvas.Modified();
1115       canvas.SetGrid();
1116 
1117       auto h2_BarycenterDiff =
1118           std::make_unique<TH2F>("Parameters diff", "SubDetector Barycenter Difference", 3, 0.0, 3.0, 6, 0, 6.);
1119 
1120       h2_BarycenterDiff->SetStats(false);
1121       h2_BarycenterDiff->SetTitle(nullptr);
1122       h2_BarycenterDiff->GetXaxis()->SetBinLabel(1, "X [#mum]");
1123       h2_BarycenterDiff->GetXaxis()->SetBinLabel(2, "Y [#mum]");
1124       h2_BarycenterDiff->GetXaxis()->SetBinLabel(3, "Z [#mum]");
1125 
1126       AlignmentPI::TkAlBarycenters l_barycenters;
1127       l_barycenters.computeBarycenters(last_alignments, tTopo_l, hardcodeGPR);
1128 
1129       AlignmentPI::TkAlBarycenters f_barycenters;
1130       f_barycenters.computeBarycenters(first_alignments, tTopo_f, hardcodeGPR);
1131 
1132       unsigned int yBin = 6;
1133       for (unsigned int i = 0; i < 6; i++) {
1134         auto thePart = static_cast<AlignmentPI::partitions>(i + 1);
1135         std::string theLabel = getStringFromPart(thePart);
1136         h2_BarycenterDiff->GetYaxis()->SetBinLabel(yBin, theLabel.c_str());
1137         h2_BarycenterDiff->SetBinContent(
1138             1, yBin, (l_barycenters.getX()[i] - f_barycenters.getX()[i]) * AlignmentPI::cmToUm);
1139         h2_BarycenterDiff->SetBinContent(
1140             2, yBin, (l_barycenters.getY()[i] - f_barycenters.getY()[i]) * AlignmentPI::cmToUm);
1141         h2_BarycenterDiff->SetBinContent(
1142             3, yBin, (l_barycenters.getZ()[i] - f_barycenters.getZ()[i]) * AlignmentPI::cmToUm);
1143         yBin--;
1144       }
1145 
1146       h2_BarycenterDiff->GetXaxis()->LabelsOption("h");
1147       h2_BarycenterDiff->GetYaxis()->SetLabelSize(0.05);
1148       h2_BarycenterDiff->GetXaxis()->SetLabelSize(0.05);
1149       h2_BarycenterDiff->SetMarkerSize(1.5);
1150       h2_BarycenterDiff->SetMarkerColor(kRed);
1151       h2_BarycenterDiff->Draw("TEXT");
1152 
1153       TLatex t1;
1154       t1.SetNDC();
1155       t1.SetTextAlign(26);
1156       t1.SetTextSize(0.05);
1157       t1.DrawLatex(0.5, 0.96, Form("Tracker Alignment Barycenters Diff, IOV %i - IOV %i", last_run, first_run));
1158       t1.SetTextSize(0.025);
1159 
1160       std::string fileName(this->m_imageFileName);
1161       canvas.SaveAs(fileName.c_str());
1162 
1163       return true;
1164     }
1165 
1166   private:
1167     bool isInitialPhase0;
1168     bool isFinalPhase0;
1169   };
1170 
1171   using TrackerAlignmentBarycentersCompare = TrackerAlignmentBarycentersComparatorBase<1, MULTI_IOV>;
1172   using TrackerAlignmentBarycentersCompareTwoTags = TrackerAlignmentBarycentersComparatorBase<2, SINGLE_IOV>;
1173 
1174   /************************************************
1175     Comparator of Pixel Tracker Detector barycenters
1176   *************************************************/
1177   template <int ntags, IOVMultiplicity nIOVs>
1178   class PixelBarycentersComparatorBase : public PlotImage<Alignments, nIOVs, ntags> {
1179   public:
1180     PixelBarycentersComparatorBase() : PlotImage<Alignments, nIOVs, ntags>("Comparison of Pixel Barycenters") {}
1181 
1182     bool fill() override {
1183       // trick to deal with the multi-ioved tag and two tag case at the same time
1184       auto theIOVs = PlotBase::getTag<0>().iovs;
1185       auto tagname1 = PlotBase::getTag<0>().name;
1186       std::string tagname2 = "";
1187       auto firstiov = theIOVs.front();
1188       std::tuple<cond::Time_t, cond::Hash> lastiov;
1189 
1190       // we don't support (yet) comparison with more than 2 tags
1191       assert(this->m_plotAnnotations.ntags < 3);
1192 
1193       if (this->m_plotAnnotations.ntags == 2) {
1194         auto tag2iovs = PlotBase::getTag<1>().iovs;
1195         tagname2 = PlotBase::getTag<1>().name;
1196         lastiov = tag2iovs.front();
1197       } else {
1198         lastiov = theIOVs.back();
1199       }
1200 
1201       unsigned int first_run = std::get<0>(firstiov);
1202       unsigned int last_run = std::get<0>(lastiov);
1203 
1204       std::shared_ptr<Alignments> last_payload = this->fetchPayload(std::get<1>(lastiov));
1205       std::vector<AlignTransform> last_alignments = last_payload->m_align;
1206 
1207       std::shared_ptr<Alignments> first_payload = this->fetchPayload(std::get<1>(firstiov));
1208       std::vector<AlignTransform> first_alignments = first_payload->m_align;
1209 
1210       TCanvas canvas("Pixel Barycenter Summary", "Pixel Barycenter summary", 1200, 1200);
1211       canvas.Divide(2, 2);
1212       canvas.cd();
1213 
1214       TLatex t1;
1215       t1.SetNDC();
1216       t1.SetTextAlign(26);
1217       t1.SetTextSize(0.03);
1218       t1.DrawLatex(0.5,
1219                    0.97,
1220                    ("Pixel Barycenters comparison, IOV: #color[2]{" + std::to_string(first_run) +
1221                     "} vs IOV: #color[4]{" + std::to_string(last_run) + "}")
1222                        .c_str());
1223       t1.SetTextSize(0.025);
1224 
1225       for (unsigned int c = 1; c <= 4; c++) {
1226         canvas.cd(c)->SetTopMargin(0.07);
1227         canvas.cd(c)->SetBottomMargin(0.12);
1228         canvas.cd(c)->SetLeftMargin(0.15);
1229         canvas.cd(c)->SetRightMargin(0.03);
1230         canvas.cd(c)->Modified();
1231         canvas.cd(c)->SetGrid();
1232       }
1233 
1234       std::array<std::string, 3> structures = {{"FPIX-", "BPIX", "FPIX+"}};
1235       std::array<std::unique_ptr<TH2F>, 3> histos;
1236 
1237       isInitialPhase0 = (first_alignments.size() == AlignmentPI::phase0size) ? true : false;
1238       isFinalPhase0 = (last_alignments.size() == AlignmentPI::phase0size) ? true : false;
1239 
1240       // check that the geomtery is a tracker one
1241       const char *path_toTopologyXML = isInitialPhase0 ? "Geometry/TrackerCommonData/data/trackerParameters.xml"
1242                                                        : "Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml";
1243 
1244       TrackerTopology tTopo_f =
1245           StandaloneTrackerTopology::fromTrackerParametersXMLFile(edm::FileInPath(path_toTopologyXML).fullPath());
1246 
1247       AlignmentPI::TkAlBarycenters myInitialBarycenters;
1248       //myInitialBarycenters.computeBarycenters(first_alignments,tTopo_f,hardcodeGPR);
1249       myInitialBarycenters.computeBarycenters(
1250           first_alignments, tTopo_f, {{AlignmentPI::t_x, 0.0}, {AlignmentPI::t_y, 0.0}, {AlignmentPI::t_z, 0.0}});
1251 
1252       path_toTopologyXML = isFinalPhase0 ? "Geometry/TrackerCommonData/data/trackerParameters.xml"
1253                                          : "Geometry/TrackerCommonData/data/PhaseI/trackerParameters.xml";
1254 
1255       TrackerTopology tTopo_l =
1256           StandaloneTrackerTopology::fromTrackerParametersXMLFile(edm::FileInPath(path_toTopologyXML).fullPath());
1257 
1258       AlignmentPI::TkAlBarycenters myFinalBarycenters;
1259       //myFinalBarycenters.computeBarycenters(last_alignments,tTopo_l,hardcodeGPR);
1260       myFinalBarycenters.computeBarycenters(
1261           last_alignments, tTopo_l, {{AlignmentPI::t_x, 0.0}, {AlignmentPI::t_y, 0.0}, {AlignmentPI::t_z, 0.0}});
1262 
1263       if (isFinalPhase0 != isInitialPhase0) {
1264         edm::LogWarning("TrackerAlignment_PayloadInspector")
1265             << "the size of the reference alignment (" << first_alignments.size()
1266             << ") is different from the one of the target (" << last_alignments.size()
1267             << ")! You are probably trying to compare different underlying geometries.";
1268       }
1269 
1270       unsigned int index(0);
1271       for (const auto &piece : structures) {
1272         const char *name = piece.c_str();
1273         histos[index] = std::make_unique<TH2F>(
1274             name,
1275             Form("%s x-y Barycenter Difference;x_{%s}-x_{TOB} [mm];y_{%s}-y_{TOB} [mm]", name, name, name),
1276             100,
1277             -3.,
1278             3.,
1279             100,
1280             -3.,
1281             3.);
1282 
1283         histos[index]->SetStats(false);
1284         histos[index]->SetTitle(nullptr);
1285         histos[index]->GetYaxis()->SetLabelSize(0.05);
1286         histos[index]->GetXaxis()->SetLabelSize(0.05);
1287         histos[index]->GetYaxis()->SetTitleSize(0.06);
1288         histos[index]->GetXaxis()->SetTitleSize(0.06);
1289         histos[index]->GetYaxis()->CenterTitle();
1290         histos[index]->GetXaxis()->CenterTitle();
1291         histos[index]->GetXaxis()->SetTitleOffset(0.9);
1292         index++;
1293       }
1294 
1295       auto h2_ZBarycenterDiff = std::make_unique<TH2F>(
1296           "Pixel_z_diff", "Pixel z-Barycenter Difference;; z_{Pixel-Ideal} -z_{TOB} [mm]", 3, -0.5, 2.5, 100, -10., 10.);
1297       h2_ZBarycenterDiff->SetStats(false);
1298       h2_ZBarycenterDiff->SetTitle(nullptr);
1299       h2_ZBarycenterDiff->GetXaxis()->SetBinLabel(1, "FPIX -");
1300       h2_ZBarycenterDiff->GetXaxis()->SetBinLabel(2, "BPIX");
1301       h2_ZBarycenterDiff->GetXaxis()->SetBinLabel(3, "FPIX +");
1302       h2_ZBarycenterDiff->GetYaxis()->SetLabelSize(0.05);
1303       h2_ZBarycenterDiff->GetXaxis()->SetLabelSize(0.07);
1304       h2_ZBarycenterDiff->GetYaxis()->SetTitleSize(0.06);
1305       h2_ZBarycenterDiff->GetXaxis()->SetTitleSize(0.06);
1306       h2_ZBarycenterDiff->GetYaxis()->CenterTitle();
1307       h2_ZBarycenterDiff->GetXaxis()->CenterTitle();
1308       h2_ZBarycenterDiff->GetYaxis()->SetTitleOffset(1.1);
1309 
1310       std::function<GlobalPoint(int)> cutFunctorInitial = [&myInitialBarycenters](int index) {
1311         switch (index) {
1312           case 1:
1313             return myInitialBarycenters.getPartitionAvg(AlignmentPI::PARTITION::FPIXm);
1314           case 2:
1315             return myInitialBarycenters.getPartitionAvg(AlignmentPI::PARTITION::BPIX);
1316           case 3:
1317             return myInitialBarycenters.getPartitionAvg(AlignmentPI::PARTITION::FPIXp);
1318           default:
1319             return GlobalPoint(0, 0, 0);
1320         }
1321       };
1322 
1323       std::function<GlobalPoint(int)> cutFunctorFinal = [&myFinalBarycenters](int index) {
1324         switch (index) {
1325           case 1:
1326             return myFinalBarycenters.getPartitionAvg(AlignmentPI::PARTITION::FPIXm);
1327           case 2:
1328             return myFinalBarycenters.getPartitionAvg(AlignmentPI::PARTITION::BPIX);
1329           case 3:
1330             return myFinalBarycenters.getPartitionAvg(AlignmentPI::PARTITION::FPIXp);
1331           default:
1332             return GlobalPoint(0, 0, 0);
1333         }
1334       };
1335 
1336       float x0i, x0f, y0i, y0f;
1337 
1338       t1.SetNDC(kFALSE);
1339       t1.SetTextSize(0.047);
1340       for (unsigned int c = 1; c <= 3; c++) {
1341         x0i = cutFunctorInitial(c).x() * 10;  // transform cm to mm (x10)
1342         x0f = cutFunctorFinal(c).x() * 10;
1343         y0i = cutFunctorInitial(c).y() * 10;
1344         y0f = cutFunctorFinal(c).y() * 10;
1345 
1346         canvas.cd(c);
1347         histos[c - 1]->Draw();
1348 
1349         COUT << "initial x,y " << std::left << std::setw(7) << structures[c - 1] << " (" << x0i << "," << y0i << ") mm"
1350              << std::endl;
1351         COUT << "final   x,y " << std::left << std::setw(7) << structures[c - 1] << " (" << x0f << "," << y0f << ") mm"
1352              << std::endl;
1353 
1354         TMarker *initial = new TMarker(x0i, y0i, 21);
1355         TMarker *final = new TMarker(x0f, y0f, 20);
1356 
1357         initial->SetMarkerColor(kRed);
1358         final->SetMarkerColor(kBlue);
1359         initial->SetMarkerSize(2.5);
1360         final->SetMarkerSize(2.5);
1361         t1.SetTextColor(kRed);
1362         initial->Draw();
1363         t1.DrawLatex(x0i, y0i + (y0i > y0f ? 0.3 : -0.5), Form("(%.2f,%.2f)", x0i, y0i));
1364         final->Draw("same");
1365         t1.SetTextColor(kBlue);
1366         t1.DrawLatex(x0f, y0f + (y0i > y0f ? -0.5 : 0.3), Form("(%.2f,%.2f)", x0f, y0f));
1367       }
1368 
1369       // fourth pad is a special case for the z coordinate
1370       canvas.cd(4);
1371       h2_ZBarycenterDiff->Draw();
1372       float z0i, z0f;
1373 
1374       // numbers do agree with https://twiki.cern.ch/twiki/bin/view/CMSPublic/TkAlignmentPerformancePhaseIStartUp17#Pixel_Barycentre_Positions
1375 
1376       std::array<double, 3> hardcodeIdealZPhase0 = {{-41.94909, 0., 41.94909}};  // units are cm
1377       std::array<double, 3> hardcodeIdealZPhase1 = {{-39.82911, 0., 39.82911}};  // units are cm
1378 
1379       for (unsigned int c = 1; c <= 3; c++) {
1380         // less than pretty but needed to remove the z position of the FPix barycenters != 0
1381 
1382         z0i =
1383             (cutFunctorInitial(c).z() - (isInitialPhase0 ? hardcodeIdealZPhase0[c - 1] : hardcodeIdealZPhase1[c - 1])) *
1384             10;  // convert to mm
1385         z0f =
1386             (cutFunctorFinal(c).z() - (isFinalPhase0 ? hardcodeIdealZPhase0[c - 1] : hardcodeIdealZPhase1[c - 1])) * 10;
1387 
1388         TMarker *initial = new TMarker(c - 1, z0i, 21);
1389         TMarker *final = new TMarker(c - 1, z0f, 20);
1390 
1391         COUT << "initial   z " << std::left << std::setw(7) << structures[c - 1] << " " << z0i << " mm" << std::endl;
1392         COUT << "final     z " << std::left << std::setw(7) << structures[c - 1] << " " << z0f << " mm" << std::endl;
1393 
1394         initial->SetMarkerColor(kRed);
1395         final->SetMarkerColor(kBlue);
1396         initial->SetMarkerSize(2.5);
1397         final->SetMarkerSize(2.5);
1398         initial->Draw();
1399         t1.SetTextColor(kRed);
1400         t1.DrawLatex(c - 1, z0i + (z0i > z0f ? 1. : -1.5), Form("(%.2f)", z0i));
1401         final->Draw("same");
1402         t1.SetTextColor(kBlue);
1403         t1.DrawLatex(c - 1, z0f + (z0i > z0f ? -1.5 : 1), Form("(%.2f)", z0f));
1404       }
1405 
1406       std::string fileName(this->m_imageFileName);
1407       canvas.SaveAs(fileName.c_str());
1408 
1409       return true;
1410     }
1411 
1412   private:
1413     bool isInitialPhase0;
1414     bool isFinalPhase0;
1415   };
1416 
1417   using PixelBarycentersCompare = PixelBarycentersComparatorBase<1, MULTI_IOV>;
1418   using PixelBarycentersCompareTwoTags = PixelBarycentersComparatorBase<2, SINGLE_IOV>;
1419 
1420 }  // namespace
1421 
1422 PAYLOAD_INSPECTOR_MODULE(TrackerAlignment) {
1423   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentComparatorSingleTag);
1424   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentComparatorTwoTags);
1425   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentComparatorSingleTag);
1426   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentComparatorTwoTags);
1427   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareX);
1428   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareY);
1429   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareZ);
1430   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareAlpha);
1431   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareBeta);
1432   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareGamma);
1433   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareXTwoTags);
1434   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareYTwoTags);
1435   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareZTwoTags);
1436   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareAlphaTwoTags);
1437   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareBetaTwoTags);
1438   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentCompareGammaTwoTags);
1439   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapX);
1440   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapY);
1441   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapZ);
1442   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapAlpha);
1443   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapBeta);
1444   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapGamma);
1445   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapXTwoTags);
1446   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapYTwoTags);
1447   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapZTwoTags);
1448   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapAlphaTwoTags);
1449   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapBetaTwoTags);
1450   PAYLOAD_INSPECTOR_CLASS(PixelAlignmentCompareMapGammaTwoTags);
1451   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryBPix);
1452   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryFPix);
1453   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryTIB);
1454   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryTID);
1455   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryTOB);
1456   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryTEC);
1457   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryBPixTwoTags);
1458   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryFPixTwoTags);
1459   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryTIBTwoTags);
1460   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryTIDTwoTags);
1461   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryTOBTwoTags);
1462   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentSummaryTECTwoTags);
1463   PAYLOAD_INSPECTOR_CLASS(X_BPixBarycenterHistory);
1464   PAYLOAD_INSPECTOR_CLASS(Y_BPixBarycenterHistory);
1465   PAYLOAD_INSPECTOR_CLASS(Z_BPixBarycenterHistory);
1466   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentBarycenters);
1467   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentBarycentersCompare);
1468   PAYLOAD_INSPECTOR_CLASS(TrackerAlignmentBarycentersCompareTwoTags);
1469   PAYLOAD_INSPECTOR_CLASS(PixelBarycentersCompare);
1470   PAYLOAD_INSPECTOR_CLASS(PixelBarycentersCompareTwoTags);
1471 }