BeamSpotProducer/plugins/BeamSpotCompatibilityChecker.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    RecoVertex/BeamSpotProducer
0004 // Class:      BeamSpotCompatibilityChecker
0005 //
0006 /**\class BeamSpotCompatibilityChecker BeamSpotCompatibilityChecker.cc RecoVertex/BeamSpotProducer/plugins/BeamSpotCompatibilityChecker.cc
0007
0008  Description: Class to check the compatibilty between the BeamSpot payload in the database and the one in the event
0009
0010  Implementation:
0011      Makes use of the Significance struct to establish how compatible are the data members of the two BeamSpots in input
0012 */
0013 //
0014 // Original Author:  Marco Musich
0015 //         Created:  Thu, 23 Apr 2020 09:00:45 GMT
0016 //
0017 //
0018
0019 // system include files
0020 #include <memory>
0021
0022 // user include files
0023 #include "CondFormats/BeamSpotObjects/interface/BeamSpotObjects.h"
0024 #include "CondFormats/DataRecord/interface/BeamSpotObjectsRcd.h"
0025 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0026 #include "DataFormats/Common/interface/Handle.h"
0027 #include "DataFormats/Scalers/interface/BeamSpotOnline.h"
0028 #include "FWCore/Framework/interface/ESHandle.h"
0029 #include "FWCore/Framework/interface/ESWatcher.h"
0030 #include "FWCore/Framework/interface/Event.h"
0031 #include "FWCore/Framework/interface/EventSetup.h"
0032 #include "FWCore/Framework/interface/Frameworkfwd.h"
0033 #include "FWCore/Framework/interface/MakerMacros.h"
0034 #include "FWCore/Framework/interface/global/EDAnalyzer.h"
0035 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0036 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0037 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0038 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0039 #include "FWCore/Utilities/interface/ESGetToken.h"
0040 #include "FWCore/Utilities/interface/InputTag.h"
0041
0042 //
0043 // class declaration
0044 //
0045
0046 // ancillary struct to check compatibility
0047 namespace {
0048
0049   struct Significance {
0050     Significance(const double& a, const double& b, const double& errA, const double& errB)
0051         : m_A(a), m_B(b), m_ErrA(errA), m_ErrB(errB) {
0052       if (m_ErrA == 0 && m_ErrB == 0) {
0053         edm::LogError("LogicalError") << "Can't calculate significance when both errors are zero!" << std::endl;
0054       }
0055       m_combinedError = std::sqrt((m_ErrA * m_ErrA) + (m_ErrB * m_ErrB));
0056     }
0057
0058     float getSig(bool verbose) {
0059       if (verbose) {
0060         edm::LogInfo("BeamSpotCompatibilityChecker")
0061             << "A= " << m_A << "+/-" << m_ErrA << "    "
0062             << "B= " << m_B << "+/-" << m_ErrB << " | delta=" << std::abs(m_A - m_B) << "+/-" << m_combinedError
0063             << std::endl;
0064       }
0065       return std::abs(m_A - m_B) / m_combinedError;
0066     }
0067
0068   private:
0069     double m_A;
0070     double m_B;
0071     double m_ErrA;
0072     double m_ErrB;
0073     double m_ErrAB;
0074     double m_combinedError;
0075   };
0076 }  // namespace
0077
0078 class BeamSpotCompatibilityChecker : public edm::global::EDAnalyzer<> {
0079 public:
0080   explicit BeamSpotCompatibilityChecker(const edm::ParameterSet&);
0081   ~BeamSpotCompatibilityChecker() override = default;
0082
0083   void analyze(edm::StreamID, edm::Event const&, edm::EventSetup const&) const final;
0084   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0085   static std::array<float, 6> compareBS(const reco::BeamSpot& BSA, const reco::BeamSpot& BSB);
0086
0087 private:
0088   // ----------member data ---------------------------
0089   static constexpr int cmToum = 10000;
0090   const double warningThreshold_;
0091   const double throwingThreshold_;
0092   const bool verbose_;
0093   const bool dbFromEvent_;
0094   const edm::EDGetTokenT<reco::BeamSpot> bsFromEventToken_;             // beamSpot from the event
0095   edm::ESGetToken<BeamSpotObjects, BeamSpotObjectsRcd> bsFromDBToken_;  // beamSpot from the DB
0096   edm::EDGetTokenT<reco::BeamSpot> dbBSfromEventToken_;                 // beamSpot from the DB (via event)
0097 };
0098
0099 //
0100 // constructors and destructor
0101 //
0102 BeamSpotCompatibilityChecker::BeamSpotCompatibilityChecker(const edm::ParameterSet& iConfig)
0103     : warningThreshold_(iConfig.getParameter<double>("warningThr")),
0104       throwingThreshold_(iConfig.getParameter<double>("errorThr")),
0105       verbose_(iConfig.getUntrackedParameter<bool>("verbose", false)),
0106       dbFromEvent_(iConfig.getParameter<bool>("dbFromEvent")),
0107       bsFromEventToken_(consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("bsFromEvent"))) {
0108   //now do what ever initialization is needed
0109   if (warningThreshold_ > throwingThreshold_) {
0110     throw cms::Exception("ConfigurationError")
0111         << __PRETTY_FUNCTION__ << "\n Warning threshold (" << warningThreshold_
0112         << ") cannot be smaller than the throwing threshold (" << throwingThreshold_ << ")" << std::endl;
0113   }
0114   if (dbFromEvent_) {
0115     edm::LogWarning("BeamSpotCompatibilityChecker")
0116         << "!!!! Warning !!!\nThe Database Beam Spot is going to be taken from the Event via BeamSpotProducer!";
0117     dbBSfromEventToken_ = consumes<reco::BeamSpot>(iConfig.getParameter<edm::InputTag>("bsFromDB"));
0118   } else {
0119     bsFromDBToken_ = esConsumes<BeamSpotObjects, BeamSpotObjectsRcd>();
0120   }
0121 }
0122
0123 //
0124 // member functions
0125 //
0126
0127 // ------------ method called for each event  ------------
0128 void BeamSpotCompatibilityChecker::analyze(edm::StreamID sid,
0129                                            edm::Event const& iEvent,
0130                                            edm::EventSetup const& iSetup) const {
0131   using namespace edm;
0132   reco::BeamSpot spotEvent, spotDB;
0133
0134   edm::Handle<reco::BeamSpot> beamSpotFromEventHandle;
0135   iEvent.getByToken(bsFromEventToken_, beamSpotFromEventHandle);
0136   spotEvent = *beamSpotFromEventHandle;
0137
0138   double evt_BSx0 = spotEvent.x0();
0139   double evt_BSy0 = spotEvent.y0();
0140   double evt_BSz0 = spotEvent.z0();
0141   double evt_Beamsigmaz = spotEvent.sigmaZ();
0142   double evt_BeamWidthX = spotEvent.BeamWidthX();
0143   double evt_BeamWidthY = spotEvent.BeamWidthY();
0144
0145   if (!dbFromEvent_) {
0146     edm::ESHandle<BeamSpotObjects> beamhandle = iSetup.getHandle(bsFromDBToken_);
0147     const BeamSpotObjects* aSpot = beamhandle.product();
0148
0149     // translate from BeamSpotObjects to reco::BeamSpot
0150     reco::BeamSpot::Point apoint(aSpot->x(), aSpot->y(), aSpot->z());
0151
0152     reco::BeamSpot::CovarianceMatrix matrix;
0153     for (int i = 0; i < reco::BeamSpot::dimension; ++i) {
0154       for (int j = 0; j < reco::BeamSpot::dimension; ++j) {
0155         matrix(i, j) = aSpot->covariance(i, j);
0156       }
0157     }
0158
0159     // this assume beam width same in x and y
0160     spotDB = reco::BeamSpot(apoint, aSpot->sigmaZ(), aSpot->dxdz(), aSpot->dydz(), aSpot->beamWidthX(), matrix);
0161   } else {
0162     // take the DB beamspot from the event (different label)
0163     edm::Handle<reco::BeamSpot> beamSpotFromDBHandle;
0164     iEvent.getByToken(dbBSfromEventToken_, beamSpotFromDBHandle);
0165     spotDB = *beamSpotFromDBHandle;
0166   }
0167
0168   double db_BSx0 = spotDB.x0();
0169   double db_BSy0 = spotDB.y0();
0170   double db_BSz0 = spotDB.z0();
0171   double db_Beamsigmaz = spotDB.sigmaZ();
0172   double db_BeamWidthX = spotDB.BeamWidthX();
0173   double db_BeamWidthY = spotDB.BeamWidthY();
0174
0175   double deltaX0 = evt_BSx0 - db_BSx0;
0176   double deltaY0 = evt_BSy0 - db_BSy0;
0177   double deltaZ0 = evt_BSz0 - db_BSz0;
0178   double deltaSigmaX = evt_BeamWidthX - db_BeamWidthX;
0179   double deltaSigmaY = evt_BeamWidthY - db_BeamWidthY;
0180   double deltaSigmaZ = evt_Beamsigmaz - db_Beamsigmaz;
0181
0182   if (verbose_) {
0183     edm::LogPrint("BeamSpotCompatibilityChecker") << "BS from Event: \n" << spotEvent << std::endl;
0184     edm::LogPrint("BeamSpotCompatibilityChecker") << "BS from DB: \n" << spotDB << std::endl;
0185   }
0186
0187   auto significances = compareBS(spotDB, spotEvent);
0188   std::vector<std::string> labels = {"x0", "y0", "z0", "sigmaX", "sigmaY", "sigmaZ"};
0189
0190   std::string msg = " |delta X0|=" + std::to_string(std::abs(deltaX0) * cmToum) +
0191                     " um |delta Y0|=" + std::to_string(std::abs(deltaY0) * cmToum) +
0192                     " um |delta Z0|=" + std::to_string(std::abs(deltaZ0) * cmToum) +
0193                     " um |delta sigmaX|=" + std::to_string(std::abs(deltaSigmaX) * cmToum) +
0194                     " um |delta sigmaY|=" + std::to_string(std::abs(deltaSigmaY) * cmToum) +
0195                     " um |delta sigmaZ|=" + std::to_string(std::abs(deltaSigmaZ)) + "cm";
0196   if (verbose_) {
0197     edm::LogPrint("BeamSpotCompatibilityChecker") << msg.c_str() << std::endl;
0198   }
0199
0200   for (unsigned int i = 0; i < 3; i++) {
0201     auto sig = significances.at(i);
0202     if (sig > throwingThreshold_) {
0203       edm::LogError("BeamSpotCompatibilityChecker") << msg.c_str() << std::endl;
0204       throw cms::Exception("BeamSpotCompatibilityChecker")
0205           << "[" << __PRETTY_FUNCTION__ << "] \n DB-Event BeamSpot " << labels.at(i) << " distance sigificance " << sig
0206           << ", exceeds the threshold of " << throwingThreshold_ << "!" << std::endl;
0207     } else if (sig > warningThreshold_) {
0208       edm::LogWarning("BeamSpotCompatibilityChecker") << msg.c_str() << std::endl;
0209       edm::LogWarning("BeamSpotCompatibilityChecker")
0210           << "[" << __PRETTY_FUNCTION__ << "]  \n  DB-Event BeamSpot " << labels.at(i) << " distance significance "
0211           << sig << ", exceeds the threshold of " << warningThreshold_ << "!" << std::endl;
0212     }
0213   }
0214 }
0215
0216 std::array<float, 6> BeamSpotCompatibilityChecker::compareBS(const reco::BeamSpot& beamSpotA,
0217                                                              const reco::BeamSpot& beamSpotB) {
0218   auto xS = Significance(beamSpotA.x0(), beamSpotB.x0(), beamSpotA.x0Error(), beamSpotB.x0Error());
0219   auto yS = Significance(beamSpotA.y0(), beamSpotB.y0(), beamSpotA.y0Error(), beamSpotB.y0Error());
0220   auto zS = Significance(beamSpotA.z0(), beamSpotB.z0(), beamSpotA.z0Error(), beamSpotB.z0Error());
0221   auto xWS = Significance(
0222       beamSpotA.BeamWidthX(), beamSpotB.BeamWidthX(), beamSpotA.BeamWidthXError(), beamSpotB.BeamWidthXError());
0223   auto yWS = Significance(
0224       beamSpotA.BeamWidthY(), beamSpotB.BeamWidthY(), beamSpotA.BeamWidthYError(), beamSpotB.BeamWidthYError());
0225   auto zWS = Significance(beamSpotA.sigmaZ(), beamSpotB.sigmaZ(), beamSpotA.sigmaZ0Error(), beamSpotB.sigmaZ0Error());
0226
0227   std::array<float, 6> ret = {
0228       {xS.getSig(false), yS.getSig(false), zS.getSig(false), xWS.getSig(false), yWS.getSig(false), zWS.getSig(false)}};
0229   return ret;
0230 }
0231
0232 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0233 void BeamSpotCompatibilityChecker::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0234   edm::ParameterSetDescription desc;
0235   desc.add<double>("warningThr", 1.)->setComment("Threshold on the signficances to emit a warning");
0236   desc.add<double>("errorThr", 3.)->setComment("Threshold on the signficances to abort the job");
0237   desc.addUntracked<bool>("verbose", false)->setComment("verbose output");
0238   desc.add<edm::InputTag>("bsFromEvent", edm::InputTag(""))
0239       ->setComment("edm::InputTag on the BeamSpot from the Event (Reference)");
0240   desc.add<bool>("dbFromEvent", false)
0241       ->setComment("Switch to take the (target) DB beamspot from the event instead of the EventSetup");
0242   desc.add<edm::InputTag>("bsFromDB", edm::InputTag(""))
0243       ->setComment("edm::InputTag on the BeamSpot from the Event (Target)\n To be used only if dbFromEvent is True");
0244   descriptions.addWithDefaultLabel(desc);
0245 }
0246
0247 //define this as a plug-in
0248 DEFINE_FWK_MODULE(BeamSpotCompatibilityChecker);