Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-05-14 00:30:02

0001 // \file AlignmentRcdChecker.cpp
0002 //
0003 // \author    : Marco Musich
0004 // Revision   : $Revision: 1.1 $
0005 // last update: $Date: 2022/05/11 14:44:00 $
0006 // by         : $Author: musich $
0007 
0008 // system includes
0009 #include <string>
0010 #include <map>
0011 #include <vector>
0012 
0013 // user includes
0014 #include "CondFormats/Alignment/interface/AlignTransform.h"
0015 #include "CondFormats/Alignment/interface/Alignments.h"
0016 #include "CondFormats/AlignmentRecord/interface/CSCAlignmentRcd.h"
0017 #include "CondFormats/AlignmentRecord/interface/DTAlignmentRcd.h"
0018 #include "CondFormats/AlignmentRecord/interface/TrackerAlignmentRcd.h"
0019 #include "DataFormats/DetId/interface/DetId.h"
0020 #include "FWCore/Framework/interface/ESHandle.h"
0021 #include "FWCore/Framework/interface/ESWatcher.h"
0022 #include "FWCore/Framework/interface/Event.h"
0023 #include "FWCore/Framework/interface/EventSetup.h"
0024 #include "FWCore/Framework/interface/MakerMacros.h"
0025 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0026 
0027 // ROOT includes
0028 #include <TMath.h>
0029 
0030 class AlignmentRcdChecker : public edm::one::EDAnalyzer<> {
0031 public:
0032   explicit AlignmentRcdChecker(const edm::ParameterSet& iConfig);
0033   ~AlignmentRcdChecker() = default;
0034 
0035   virtual void analyze(const edm::Event& evt, const edm::EventSetup& evtSetup);
0036 
0037 private:
0038   void inspectRecord(const std::string& rcdname,
0039                      const edm::Event& evt,
0040                      const Alignments* refAlignments,
0041                      const Alignments* alignments);
0042 
0043   bool verbose_;
0044   bool compareStrict_;
0045   std::string label_;
0046 
0047   const edm::ESGetToken<Alignments, TrackerAlignmentRcd> tkAliTokenRef_;
0048   const edm::ESGetToken<Alignments, TrackerAlignmentRcd> tkAliTokenNew_;
0049 };
0050 
0051 AlignmentRcdChecker::AlignmentRcdChecker(const edm::ParameterSet& iConfig)
0052     : verbose_(iConfig.getParameter<bool>("verbose")),
0053       compareStrict_(iConfig.getParameter<bool>("compareStrict")),
0054       label_(iConfig.getParameter<std::string>("label")),
0055       tkAliTokenRef_(esConsumes()),
0056       tkAliTokenNew_(esConsumes(edm::ESInputTag{"", label_})) {}
0057 
0058 void AlignmentRcdChecker::analyze(const edm::Event& evt, const edm::EventSetup& evtSetup) {
0059   const Alignments* alignmentsRef = &evtSetup.getData(tkAliTokenRef_);
0060   const Alignments* alignmentsNew = &evtSetup.getData(tkAliTokenNew_);
0061   inspectRecord("TrackerAlignmentRcd", evt, alignmentsRef, alignmentsNew);
0062 }
0063 
0064 void AlignmentRcdChecker::inspectRecord(const std::string& rcdname,
0065                                         const edm::Event& evt,
0066                                         const Alignments* refAlignments,
0067                                         const Alignments* alignments) {
0068   edm::LogPrint("inspectRecord") << rcdname << " content starting from run " << evt.run();
0069   edm::LogPrint("inspectRecord") << " with " << alignments->m_align.size() << " entries";
0070 
0071   if (refAlignments && alignments) {
0072     double meanX = 0;
0073     double rmsX = 0;
0074     double meanY = 0;
0075     double rmsY = 0;
0076     double meanZ = 0;
0077     double rmsZ = 0;
0078     double meanR = 0;
0079     double rmsR = 0;
0080     double dPhi;
0081     double meanPhi = 0;
0082     double rmsPhi = 0;
0083 
0084     std::vector<AlignTransform>::const_iterator iref = refAlignments->m_align.begin();
0085     for (std::vector<AlignTransform>::const_iterator i = alignments->m_align.begin(); i != alignments->m_align.end();
0086          ++i, ++iref) {
0087       meanX += i->translation().x() - iref->translation().x();
0088       rmsX += pow(i->translation().x() - iref->translation().x(), 2);
0089 
0090       meanY += i->translation().y() - iref->translation().y();
0091       rmsY += pow(i->translation().y() - iref->translation().y(), 2);
0092 
0093       meanZ += i->translation().z() - iref->translation().z();
0094       rmsZ += pow(i->translation().z() - iref->translation().z(), 2);
0095 
0096       meanR += i->translation().perp() - iref->translation().perp();
0097       rmsR += pow(i->translation().perp() - iref->translation().perp(), 2);
0098 
0099       dPhi = i->translation().phi() - iref->translation().phi();
0100       if (dPhi > M_PI)
0101         dPhi -= 2.0 * M_PI;
0102       if (dPhi < -M_PI)
0103         dPhi += 2.0 * M_PI;
0104 
0105       meanPhi += dPhi;
0106       rmsPhi += dPhi * dPhi;
0107     }
0108 
0109     meanX /= alignments->m_align.size();
0110     rmsX /= alignments->m_align.size();
0111     meanY /= alignments->m_align.size();
0112     rmsY /= alignments->m_align.size();
0113     meanZ /= alignments->m_align.size();
0114     rmsZ /= alignments->m_align.size();
0115     meanR /= alignments->m_align.size();
0116     rmsR /= alignments->m_align.size();
0117     meanPhi /= alignments->m_align.size();
0118     rmsPhi /= alignments->m_align.size();
0119 
0120     if (verbose_) {
0121       edm::LogPrint("inspectRecord") << "  Compared to previous record:";
0122       edm::LogPrint("inspectRecord") << "    mean X shift:   " << std::setw(12) << std::scientific
0123                                      << std::setprecision(3) << meanX << " (RMS = " << sqrt(rmsX) << ")";
0124       edm::LogPrint("inspectRecord") << "    mean Y shift:   " << std::setw(12) << std::scientific
0125                                      << std::setprecision(3) << meanY << " (RMS = " << sqrt(rmsY) << ")";
0126       edm::LogPrint("inspectRecord") << "    mean Z shift:   " << std::setw(12) << std::scientific
0127                                      << std::setprecision(3) << meanZ << " (RMS = " << sqrt(rmsZ) << ")";
0128       edm::LogPrint("inspectRecord") << "    mean R shift:   " << std::setw(12) << std::scientific
0129                                      << std::setprecision(3) << meanR << " (RMS = " << sqrt(rmsR) << ")";
0130       edm::LogPrint("inspectRecord") << "    mean Phi shift: " << std::setw(12) << std::scientific
0131                                      << std::setprecision(3) << meanPhi << " (RMS = " << sqrt(rmsPhi) << ")";
0132     }  // verbose
0133 
0134     if (compareStrict_) {
0135       // do not let any of the coordinates to fluctuate less then 1um
0136       assert(meanX < 1e-4);
0137       assert(meanY < 1e-4);
0138       assert(meanZ < 1e-4);
0139       assert(meanR < 1e-4);
0140       assert(meanPhi < 1e-4);
0141     }
0142 
0143   } else {
0144     throw cms::Exception("Missing Input Data") << "Could not retrieve the input alignments to compare \n";
0145   }
0146 }
0147 
0148 //define this as a plug-in
0149 DEFINE_FWK_MODULE(AlignmentRcdChecker);