Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-12-14 11:44:53

0001 /*
0002  * \file AlcaBeamMonitor.cc
0003  * \author Lorenzo Uplegger/FNAL
0004  * modified by Simone Gennai INFN/Bicocca
0005  */
0006 
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/Framework/interface/ESHandle.h"
0010 #include "FWCore/ServiceRegistry/interface/Service.h"
0011 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0012 //#include "DataFormats/Scalers/interface/BeamSpotOnline.h"
0013 #include "DQM/BeamMonitor/plugins/AlcaBeamMonitor.h"
0014 #include "DataFormats/Common/interface/Handle.h"
0015 #include "DataFormats/Common/interface/View.h"
0016 #include "FWCore/Framework/interface/ConsumesCollector.h"
0017 #include "FWCore/Framework/interface/EventSetup.h"
0018 #include "FWCore/Framework/interface/LuminosityBlock.h"
0019 #include "FWCore/Framework/interface/Run.h"
0020 #include "RecoVertex/BeamSpotProducer/interface/BeamFitter.h"
0021 #include "RecoVertex/BeamSpotProducer/interface/PVFitter.h"
0022 #include <memory>
0023 
0024 #include <numeric>
0025 
0026 using namespace std;
0027 using namespace edm;
0028 using namespace reco;
0029 
0030 //----------------------------------------------------------------------------------------------------------------------
0031 AlcaBeamMonitor::AlcaBeamMonitor(const ParameterSet& ps)
0032     : monitorName_(ps.getUntrackedParameter<string>("MonitorName")),
0033       primaryVertexLabel_(consumes<VertexCollection>(ps.getUntrackedParameter<InputTag>("PrimaryVertexLabel"))),
0034       trackLabel_(consumes<reco::TrackCollection>(ps.getUntrackedParameter<InputTag>("TrackLabel"))),
0035       scalerLabel_(consumes<BeamSpot>(ps.getUntrackedParameter<InputTag>("ScalerLabel"))),
0036       beamSpotToken_(esConsumes<edm::Transition::BeginLuminosityBlock>()),
0037       perLSsaving_(ps.getUntrackedParameter<bool>("perLSsaving", false)),
0038       numberOfValuesToSave_(0) {
0039   if (!monitorName_.empty())
0040     monitorName_ = monitorName_ + "/";
0041 
0042   theBeamFitter_ = std::make_unique<BeamFitter>(ps, consumesCollector());
0043   theBeamFitter_->resetTrkVector();
0044   theBeamFitter_->resetLSRange();
0045   theBeamFitter_->resetRefTime();
0046   theBeamFitter_->resetPVFitter();
0047 
0048   thePVFitter_ = std::make_unique<PVFitter>(ps, consumesCollector());
0049 
0050   processedLumis_.clear();
0051 
0052   varNamesV_.push_back("x");
0053   varNamesV_.push_back("y");
0054   varNamesV_.push_back("z");
0055   varNamesV_.push_back("sigmaX");
0056   varNamesV_.push_back("sigmaY");
0057   varNamesV_.push_back("sigmaZ");
0058 
0059   if (!perLSsaving_) {
0060     histoByCategoryNames_.insert(pair<string, string>("run", "Coordinate"));
0061     histoByCategoryNames_.insert(pair<string, string>("run", "PrimaryVertex fit-DataBase"));
0062     histoByCategoryNames_.insert(pair<string, string>("run", "PrimaryVertex fit-BeamFit"));
0063     histoByCategoryNames_.insert(pair<string, string>("run", "PrimaryVertex fit-Scalers"));
0064     histoByCategoryNames_.insert(pair<string, string>("run", "PrimaryVertex-DataBase"));
0065     histoByCategoryNames_.insert(pair<string, string>("run", "PrimaryVertex-BeamFit"));
0066     histoByCategoryNames_.insert(pair<string, string>("run", "PrimaryVertex-Scalers"));
0067 
0068     histoByCategoryNames_.insert(pair<string, string>("lumi", "Lumibased BeamSpotFit"));
0069     histoByCategoryNames_.insert(pair<string, string>("lumi", "Lumibased PrimaryVertex"));
0070     histoByCategoryNames_.insert(pair<string, string>("lumi", "Lumibased DataBase"));
0071     histoByCategoryNames_.insert(pair<string, string>("lumi", "Lumibased Scalers"));
0072     histoByCategoryNames_.insert(pair<string, string>("lumi", "Lumibased PrimaryVertex-DataBase fit"));
0073     histoByCategoryNames_.insert(pair<string, string>("lumi", "Lumibased PrimaryVertex-Scalers fit"));
0074     histoByCategoryNames_.insert(pair<string, string>("validation", "Lumibased Scalers-DataBase fit"));
0075     histoByCategoryNames_.insert(pair<string, string>("validation", "Lumibased PrimaryVertex-DataBase"));
0076     histoByCategoryNames_.insert(pair<string, string>("validation", "Lumibased PrimaryVertex-Scalers"));
0077   }
0078 
0079   for (vector<string>::iterator itV = varNamesV_.begin(); itV != varNamesV_.end(); itV++) {
0080     for (multimap<string, string>::iterator itM = histoByCategoryNames_.begin(); itM != histoByCategoryNames_.end();
0081          itM++) {
0082       histosMap_[*itV][itM->first][itM->second] = nullptr;
0083     }
0084   }
0085 }
0086 
0087 void AlcaBeamMonitor::fillDescriptions(edm::ConfigurationDescriptions& iDesc) {
0088   edm::ParameterSetDescription ps;
0089 
0090   ps.addUntracked<std::string>("MonitorName", "YourSubsystemName");
0091   ps.addUntracked<edm::InputTag>("PrimaryVertexLabel");
0092   ps.addUntracked<edm::InputTag>("TrackLabel");
0093   ps.addUntracked<edm::InputTag>("ScalerLabel");
0094   ps.addUntracked<bool>("perLSsaving");
0095 
0096   BeamFitter::fillDescription(ps);
0097   PVFitter::fillDescription(ps);
0098 
0099   iDesc.addDefault(ps);
0100 }
0101 
0102 //----------------------------------------------------------------------------------------------------------------------
0103 void AlcaBeamMonitor::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& iSetup) {
0104   string name;
0105   string title;
0106   int firstLumi = 1;
0107   int lastLumi = 3000;
0108   for (HistosContainer::iterator itM = histosMap_.begin(); itM != histosMap_.end(); itM++) {
0109     ibooker.setCurrentFolder(monitorName_ + "Debug");
0110     for (map<string, MonitorElement*>::iterator itMM = itM->second["run"].begin(); itMM != itM->second["run"].end();
0111          itMM++) {
0112       name = string("h") + itM->first + itMM->first;
0113       title = itM->first + "_{0} " + itMM->first;
0114       if (itM->first == "x" || itM->first == "y") {
0115         if (itMM->first == "Coordinate") {
0116           itMM->second = ibooker.book1D(name, title, 1001, -0.2525, 0.2525);
0117         } else if (itMM->first == "PrimaryVertex fit-DataBase" || itMM->first == "PrimaryVertex fit-BeamFit" ||
0118                    itMM->first == "PrimaryVertex fit-Scalers" || itMM->first == "PrimaryVertex-DataBase" ||
0119                    itMM->first == "PrimaryVertex-BeamFit" || itMM->first == "PrimaryVertex-Scalers") {
0120           itMM->second = ibooker.book1D(name, title, 1001, -0.02525, 0.02525);
0121         } else {
0122           //assert(0);
0123         }
0124       } else if (itM->first == "z") {
0125         if (itMM->first == "Coordinate") {
0126           itMM->second = ibooker.book1D(name, title, 101, -5.05, 5.05);
0127         } else if (itMM->first == "PrimaryVertex fit-DataBase" || itMM->first == "PrimaryVertex fit-BeamFit" ||
0128                    itMM->first == "PrimaryVertex fit-Scalers") {
0129           itMM->second = ibooker.book1D(name, title, 101, -0.505, 0.505);
0130         } else if (itMM->first == "PrimaryVertex-DataBase" || itMM->first == "PrimaryVertex-BeamFit" ||
0131                    itMM->first == "PrimaryVertex-Scalers") {
0132           itMM->second = ibooker.book1D(name, title, 1001, -5.005, 5.005);
0133         } else {
0134           //assert(0);
0135         }
0136       } else if (itM->first == "sigmaX" || itM->first == "sigmaY") {
0137         if (itMM->first == "Coordinate") {
0138           itMM->second = ibooker.book1D(name, title, 100, 0, 0.015);
0139         } else if (itMM->first == "PrimaryVertex fit-DataBase" || itMM->first == "PrimaryVertex fit-BeamFit" ||
0140                    itMM->first == "PrimaryVertex fit-Scalers" || itMM->first == "PrimaryVertex-DataBase" ||
0141                    itMM->first == "PrimaryVertex-BeamFit" || itMM->first == "PrimaryVertex-Scalers") {
0142           itMM->second = nullptr;
0143         } else {
0144           //assert(0);
0145         }
0146       } else if (itM->first == "sigmaZ") {
0147         if (itMM->first == "Coordinate") {
0148           itMM->second = ibooker.book1D(name, title, 110, 0, 11);
0149         } else if (itMM->first == "PrimaryVertex fit-DataBase" || itMM->first == "PrimaryVertex fit-BeamFit" ||
0150                    itMM->first == "PrimaryVertex fit-Scalers" || itMM->first == "PrimaryVertex-DataBase" ||
0151                    itMM->first == "PrimaryVertex-BeamFit" || itMM->first == "PrimaryVertex-Scalers") {
0152           itMM->second = ibooker.book1D(name, title, 101, -5.05, 5.05);
0153         } else {
0154           //assert(0);
0155         }
0156       } else {
0157         //assert(0);
0158       }
0159       if (itMM->second != nullptr) {
0160         if (itMM->first == "Coordinate") {
0161           itMM->second->setAxisTitle(itM->first + "_{0} (cm)", 1);
0162         } else if (itMM->first == "PrimaryVertex fit-DataBase" || itMM->first == "PrimaryVertex fit-BeamFit" ||
0163                    itMM->first == "PrimaryVertex fit-Scalers" || itMM->first == "PrimaryVertex-DataBase" ||
0164                    itMM->first == "PrimaryVertex-BeamFit" || itMM->first == "PrimaryVertex-Scalers") {
0165           itMM->second->setAxisTitle(itMM->first + " " + itM->first + "_{0} (cm)", 1);
0166         }
0167         itMM->second->setAxisTitle("Entries", 2);
0168       }
0169     }
0170 
0171     //Making histos per Lumi
0172     // x,y,z,sigmaX,sigmaY,sigmaZ
0173     for (map<string, map<string, MonitorElement*> >::iterator itMM = itM->second.begin(); itMM != itM->second.end();
0174          itMM++) {
0175       if (itMM->first != "run") {
0176         for (map<string, MonitorElement*>::iterator itMMM = itMM->second.begin(); itMMM != itMM->second.end();
0177              itMMM++) {
0178           name = string("h") + itM->first + itMMM->first;
0179           title = itM->first + "_{0} " + itMMM->first;
0180           if (itMM->first == "lumi") {
0181             ibooker.setCurrentFolder(monitorName_ + "Debug");
0182             itMMM->second = ibooker.book1D(name, title, lastLumi - firstLumi + 1, firstLumi - 0.5, lastLumi + 0.5);
0183             itMMM->second->setEfficiencyFlag();
0184           } else if (itMM->first == "validation" && itMMM->first == "Lumibased Scalers-DataBase fit") {
0185             ibooker.setCurrentFolder(monitorName_ + "Validation");
0186             itMMM->second = ibooker.book1D(name, title, lastLumi - firstLumi + 1, firstLumi - 0.5, lastLumi + 0.5);
0187             itMMM->second->setEfficiencyFlag();
0188           } else if (itMM->first == "validation" && itMMM->first != "Lumibased Scalers-DataBase fit" &&
0189                      (itM->first == "x" || itM->first == "y" || itM->first == "z")) {
0190             ibooker.setCurrentFolder(monitorName_ + "Validation");
0191             itMMM->second = ibooker.book1D(name, title, lastLumi - firstLumi + 1, firstLumi - 0.5, lastLumi + 0.5);
0192             itMMM->second->setEfficiencyFlag();
0193           } else if (itMM->first == "validation" &&
0194                      (itM->first == "sigmaX" || itM->first == "sigmaY" || itM->first == "sigmaZ")) {
0195             ibooker.setCurrentFolder(monitorName_ + "Validation");
0196             itMMM->second = nullptr;
0197           } else {
0198             LogInfo("AlcaBeamMonitorClient") << "Unrecognized category " << itMM->first;
0199             // assert(0);
0200           }
0201           if (itMMM->second != nullptr) {
0202             if (itMMM->first.find('-') != string::npos) {
0203               itMMM->second->setAxisTitle(string("#Delta ") + itM->first + "_{0} (cm)", 2);
0204             } else {
0205               itMMM->second->setAxisTitle(itM->first + "_{0} (cm)", 2);
0206             }
0207             itMMM->second->setAxisTitle("Lumisection", 1);
0208           }
0209         }
0210       }
0211     }
0212   }
0213 
0214   // create and cd into new folder
0215   ibooker.setCurrentFolder(monitorName_ + "Validation");
0216   //Book histograms
0217   hD0Phi0_ = ibooker.bookProfile("hD0Phi0", "d_{0} vs. #phi_{0} (All Tracks)", 63, -3.15, 3.15, 100, -0.5, 0.5, "");
0218   hD0Phi0_->setAxisTitle("#phi_{0} (rad)", 1);
0219   hD0Phi0_->setAxisTitle("d_{0} (cm)", 2);
0220 
0221   ibooker.setCurrentFolder(monitorName_ + "Debug");
0222   hDxyBS_ = ibooker.book1D("hDxyBS", "dxy_{0} w.r.t. Beam spot (All Tracks)", 100, -0.1, 0.1);
0223   hDxyBS_->setAxisTitle("dxy_{0} w.r.t. Beam spot (cm)", 1);
0224 }
0225 
0226 //----------------------------------------------------------------------------------------------------------------------
0227 std::shared_ptr<alcabeammonitor::NoCache> AlcaBeamMonitor::globalBeginLuminosityBlock(const LuminosityBlock& iLumi,
0228                                                                                       const EventSetup& iSetup) const {
0229   // Always create a beamspot group for each lumi weather we have results or not! Each Beamspot will be of unknown type!
0230 
0231   vertices_.clear();
0232   beamSpotsMap_.clear();
0233   processedLumis_.push_back(iLumi.id().luminosityBlock());
0234   //Read BeamSpot from DB
0235   ESHandle<BeamSpotObjects> bsDBHandle;
0236   try {
0237     bsDBHandle = iSetup.getHandle(beamSpotToken_);
0238   } catch (cms::Exception& exception) {
0239     LogError("AlcaBeamMonitor") << exception.what();
0240     return nullptr;
0241   }
0242   if (bsDBHandle.isValid()) {  // check the product
0243     const BeamSpotObjects* spotDB = bsDBHandle.product();
0244 
0245     // translate from BeamSpotObjects to reco::BeamSpot
0246     BeamSpot::Point apoint(spotDB->x(), spotDB->y(), spotDB->z());
0247 
0248     BeamSpot::CovarianceMatrix matrix;
0249     for (int i = 0; i < 7; ++i) {
0250       for (int j = 0; j < 7; ++j) {
0251         matrix(i, j) = spotDB->covariance(i, j);
0252       }
0253     }
0254 
0255     beamSpotsMap_["DB"] =
0256         BeamSpot(apoint, spotDB->sigmaZ(), spotDB->dxdz(), spotDB->dydz(), spotDB->beamWidthX(), matrix);
0257 
0258     BeamSpot* aSpot = &(beamSpotsMap_["DB"]);
0259 
0260     aSpot->setBeamWidthY(spotDB->beamWidthY());
0261     aSpot->setEmittanceX(spotDB->emittanceX());
0262     aSpot->setEmittanceY(spotDB->emittanceY());
0263     aSpot->setbetaStar(spotDB->betaStar());
0264 
0265     if (spotDB->beamType() == 2) {
0266       aSpot->setType(reco::BeamSpot::Tracker);
0267     } else {
0268       aSpot->setType(reco::BeamSpot::Fake);
0269     }
0270     //LogInfo("AlcaBeamMonitor")
0271     //  << *aSpot << std::endl;
0272   } else {
0273     LogInfo("AlcaBeamMonitor") << "Database BeamSpot is not valid at lumi: " << iLumi.id().luminosityBlock();
0274   }
0275   return nullptr;
0276 }
0277 
0278 //----------------------------------------------------------------------------------------------------------------------
0279 void AlcaBeamMonitor::analyze(const Event& iEvent, const EventSetup& iSetup) {
0280   //------ BeamFitter
0281   theBeamFitter_->readEvent(iEvent);
0282   //------ PVFitter
0283   thePVFitter_->readEvent(iEvent);
0284 
0285   if (beamSpotsMap_.find("DB") != beamSpotsMap_.end()) {
0286     //------ Tracks
0287     Handle<reco::TrackCollection> TrackCollection;
0288     iEvent.getByToken(trackLabel_, TrackCollection);
0289     const reco::TrackCollection* tracks = TrackCollection.product();
0290     for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++track) {
0291       hD0Phi0_->Fill(track->phi(), -1 * track->dxy());
0292       hDxyBS_->Fill(-1 * track->dxy(beamSpotsMap_["DB"].position()));
0293     }
0294   }
0295 
0296   //------ Primary Vertices
0297   Handle<VertexCollection> PVCollection;
0298   if (iEvent.getByToken(primaryVertexLabel_, PVCollection)) {
0299     vertices_.push_back(*PVCollection.product());
0300   }
0301 
0302   if (beamSpotsMap_.find("SC") == beamSpotsMap_.end()) {
0303     //BeamSpot from file for this stream is = to the scalar BeamSpot
0304     Handle<BeamSpot> recoBeamSpotHandle;
0305     try {
0306       iEvent.getByToken(scalerLabel_, recoBeamSpotHandle);
0307     } catch (cms::Exception& exception) {
0308       LogInfo("AlcaBeamMonitor") << exception.what();
0309       return;
0310     }
0311     beamSpotsMap_["SC"] = *recoBeamSpotHandle;
0312     if (beamSpotsMap_["SC"].BeamWidthX() != 0) {
0313       beamSpotsMap_["SC"].setType(reco::BeamSpot::Tracker);
0314     } else {
0315       beamSpotsMap_["SC"].setType(reco::BeamSpot::Fake);
0316     }
0317   }
0318 }
0319 
0320 //----------------------------------------------------------------------------------------------------------------------
0321 void AlcaBeamMonitor::globalEndLuminosityBlock(const LuminosityBlock& iLumi, const EventSetup& iSetup) {
0322   if (theBeamFitter_->runPVandTrkFitter()) {
0323     beamSpotsMap_["BF"] = theBeamFitter_->getBeamSpot();
0324   }
0325   theBeamFitter_->resetTrkVector();
0326   theBeamFitter_->resetLSRange();
0327   theBeamFitter_->resetRefTime();
0328   theBeamFitter_->resetPVFitter();
0329 
0330   if (thePVFitter_->runFitter()) {
0331     beamSpotsMap_["PV"] = thePVFitter_->getBeamSpot();
0332   }
0333   thePVFitter_->resetAll();
0334 
0335   //    "PV,BF..."      Value,Error
0336   map<std::string, pair<double, double> > resultsMap;
0337   vector<pair<double, double> > vertexResults;
0338   MonitorElement* histo = nullptr;
0339   for (vector<string>::iterator itV = varNamesV_.begin(); itV != varNamesV_.end(); itV++) {
0340     resultsMap.clear();
0341     for (BeamSpotContainer::iterator itBS = beamSpotsMap_.begin(); itBS != beamSpotsMap_.end(); itBS++) {
0342       if (itBS->second.type() == BeamSpot::Tracker) {
0343         if (*itV == "x") {
0344           resultsMap[itBS->first] = pair<double, double>(itBS->second.x0(), itBS->second.x0Error());
0345         } else if (*itV == "y") {
0346           resultsMap[itBS->first] = pair<double, double>(itBS->second.y0(), itBS->second.y0Error());
0347         } else if (*itV == "z") {
0348           resultsMap[itBS->first] = pair<double, double>(itBS->second.z0(), itBS->second.z0Error());
0349         } else if (*itV == "sigmaX") {
0350           resultsMap[itBS->first] = pair<double, double>(itBS->second.BeamWidthX(), itBS->second.BeamWidthXError());
0351         } else if (*itV == "sigmaY") {
0352           resultsMap[itBS->first] = pair<double, double>(itBS->second.BeamWidthY(), itBS->second.BeamWidthYError());
0353         } else if (*itV == "sigmaZ") {
0354           resultsMap[itBS->first] = pair<double, double>(itBS->second.sigmaZ(), itBS->second.sigmaZ0Error());
0355         } else {
0356           LogInfo("AlcaBeamMonitor") << "The histosMap_ has been built with the name " << *itV
0357                                      << " that I can't recognize!";
0358           //assert(0);
0359         }
0360       }
0361     }
0362     vertexResults.clear();
0363     for (vector<VertexCollection>::iterator itPV = vertices_.begin(); itPV != vertices_.end(); itPV++) {
0364       if (!itPV->empty()) {
0365         for (VertexCollection::const_iterator pv = itPV->begin(); pv != itPV->end(); pv++) {
0366           if (pv->isFake() || pv->tracksSize() < 10)
0367             continue;
0368           if (*itV == "x") {
0369             vertexResults.push_back(pair<double, double>(pv->x(), pv->xError()));
0370           } else if (*itV == "y") {
0371             vertexResults.push_back(pair<double, double>(pv->y(), pv->yError()));
0372           } else if (*itV == "z") {
0373             vertexResults.push_back(pair<double, double>(pv->z(), pv->zError()));
0374           } else if (*itV != "sigmaX" && *itV != "sigmaY" && *itV != "sigmaZ") {
0375             LogInfo("AlcaBeamMonitor") << "The histosMap_ has been built with the name " << *itV
0376                                        << " that I can't recognize!";
0377             //assert(0);
0378           }
0379         }
0380       }
0381     }
0382 
0383     for (multimap<string, string>::iterator itM = histoByCategoryNames_.begin(); itM != histoByCategoryNames_.end();
0384          itM++) {
0385       if ((histo = histosMap_[*itV][itM->first][itM->second]) == nullptr)
0386         continue;
0387       if (itM->second == "Coordinate") {
0388         if (beamSpotsMap_.find("DB") != beamSpotsMap_.end()) {
0389           histo->Fill(resultsMap["DB"].first);
0390         }
0391       } else if (itM->second == "PrimaryVertex fit-DataBase") {
0392         if (resultsMap.find("PV") != resultsMap.end() && resultsMap.find("DB") != resultsMap.end()) {
0393           histo->Fill(resultsMap["PV"].first - resultsMap["DB"].first);
0394         }
0395       } else if (itM->second == "PrimaryVertex fit-BeamFit") {
0396         if (resultsMap.find("PV") != resultsMap.end() && resultsMap.find("BF") != resultsMap.end()) {
0397           histo->Fill(resultsMap["PV"].first - resultsMap["BF"].first);
0398         }
0399       } else if (itM->second == "PrimaryVertex fit-Scalers") {
0400         if (resultsMap.find("PV") != resultsMap.end() && resultsMap.find("SC") != resultsMap.end()) {
0401           histo->Fill(resultsMap["PV"].first - resultsMap["SC"].first);
0402         }
0403       } else if (itM->second == "PrimaryVertex-DataBase") {
0404         if (resultsMap.find("PV") != resultsMap.end() && resultsMap.find("DB") != resultsMap.end()) {
0405           for (vector<pair<double, double> >::iterator itPV = vertexResults.begin(); itPV != vertexResults.end();
0406                itPV++) {
0407             histo->Fill(itPV->first - resultsMap["DB"].first);
0408           }
0409         }
0410       } else if (itM->second == "PrimaryVertex-BeamFit") {
0411         if (resultsMap.find("PV") != resultsMap.end() && resultsMap.find("BF") != resultsMap.end()) {
0412           for (vector<pair<double, double> >::iterator itPV = vertexResults.begin(); itPV != vertexResults.end();
0413                itPV++) {
0414             histo->Fill(itPV->first - resultsMap["BF"].first);
0415           }
0416         }
0417       } else if (itM->second == "PrimaryVertex-Scalers") {
0418         if (resultsMap.find("PV") != resultsMap.end() && resultsMap.find("SC") != resultsMap.end()) {
0419           for (vector<pair<double, double> >::iterator itPV = vertexResults.begin(); itPV != vertexResults.end();
0420                itPV++) {
0421             histo->Fill(itPV->first - resultsMap["SC"].first);
0422           }
0423         }
0424       } else if (itM->second == "Lumibased BeamSpotFit") {
0425         if (resultsMap.find("BF") != resultsMap.end()) {
0426           histo->setBinContent(iLumi.id().luminosityBlock(), resultsMap["BF"].first);
0427           histo->setBinError(iLumi.id().luminosityBlock(), resultsMap["BF"].second);
0428         }
0429       } else if (itM->second == "Lumibased PrimaryVertex") {
0430         if (resultsMap.find("PV") != resultsMap.end()) {
0431           histo->setBinContent(iLumi.id().luminosityBlock(), resultsMap["PV"].first);
0432           histo->setBinError(iLumi.id().luminosityBlock(), resultsMap["PV"].second);
0433         }
0434       } else if (itM->second == "Lumibased DataBase") {
0435         if (resultsMap.find("DB") != resultsMap.end()) {
0436           histo->setBinContent(iLumi.id().luminosityBlock(), resultsMap["DB"].first);
0437           histo->setBinError(iLumi.id().luminosityBlock(), resultsMap["DB"].second);
0438         }
0439       } else if (itM->second == "Lumibased Scalers") {
0440         if (resultsMap.find("SC") != resultsMap.end()) {
0441           histo->setBinContent(iLumi.id().luminosityBlock(), resultsMap["SC"].first);
0442           histo->setBinError(iLumi.id().luminosityBlock(), resultsMap["SC"].second);
0443         }
0444       } else if (itM->second == "Lumibased PrimaryVertex-DataBase fit") {
0445         if (resultsMap.find("PV") != resultsMap.end() && resultsMap.find("DB") != resultsMap.end()) {
0446           histo->setBinContent(iLumi.id().luminosityBlock(), resultsMap["PV"].first - resultsMap["DB"].first);
0447           histo->setBinError(iLumi.id().luminosityBlock(),
0448                              std::sqrt(std::pow(resultsMap["PV"].second, 2) + std::pow(resultsMap["DB"].second, 2)));
0449         }
0450       } else if (itM->second == "Lumibased PrimaryVertex-Scalers fit") {
0451         if (resultsMap.find("PV") != resultsMap.end() && resultsMap.find("SC") != resultsMap.end()) {
0452           histo->setBinContent(iLumi.id().luminosityBlock(), resultsMap["PV"].first - resultsMap["SC"].first);
0453           histo->setBinError(iLumi.id().luminosityBlock(),
0454                              std::sqrt(std::pow(resultsMap["PV"].second, 2) + std::pow(resultsMap["SC"].second, 2)));
0455         }
0456       } else if (itM->second == "Lumibased Scalers-DataBase fit") {
0457         if (resultsMap.find("SC") != resultsMap.end() && resultsMap.find("DB") != resultsMap.end()) {
0458           histo->setBinContent(iLumi.id().luminosityBlock(), resultsMap["SC"].first - resultsMap["DB"].first);
0459           histo->setBinError(iLumi.id().luminosityBlock(),
0460                              std::sqrt(std::pow(resultsMap["SC"].second, 2) + std::pow(resultsMap["DB"].second, 2)));
0461         }
0462       } else if (itM->second == "Lumibased PrimaryVertex-DataBase") {
0463         if (resultsMap.find("DB") != resultsMap.end() && !vertexResults.empty()) {
0464           for (vector<pair<double, double> >::iterator itPV = vertexResults.begin(); itPV != vertexResults.end();
0465                itPV++) {
0466             histo->setBinContent(iLumi.id().luminosityBlock(), (*itPV).first - resultsMap["DB"].first);
0467             histo->setBinError(iLumi.id().luminosityBlock(),
0468                                std::sqrt(std::pow((*itPV).second, 2) + std::pow(resultsMap["DB"].second, 2)));
0469           }
0470         }
0471       } else if (itM->second == "Lumibased PrimaryVertex-Scalers") {
0472         if (resultsMap.find("SC") != resultsMap.end() && !vertexResults.empty()) {
0473           for (vector<pair<double, double> >::iterator itPV = vertexResults.begin(); itPV != vertexResults.end();
0474                itPV++) {
0475             histo->setBinContent(iLumi.id().luminosityBlock(), (*itPV).first - resultsMap["SC"].first);
0476             histo->setBinError(iLumi.id().luminosityBlock(),
0477                                std::sqrt(std::pow((*itPV).second, 2) + std::pow(resultsMap["SC"].second, 2)));
0478           }
0479         }
0480       }
0481       //      else if(itM->second == "Lumibased Scalers-DataBase"){
0482       //      if(resultsMap.find("SC") != resultsMap.end() && resultsMap.find("DB") != resultsMap.end()){
0483       //        itHHH->second->Fill(bin,resultsMap["SC"].first-resultsMap["DB"].first);
0484       //      }
0485       //    }
0486       else {
0487         LogInfo("AlcaBeamMonitor") << "The histosMap_ have a histogram named " << itM->second
0488                                    << " that I can't recognize in this loop!";
0489         //assert(0);
0490       }
0491     }
0492   }
0493 }
0494 
0495 void AlcaBeamMonitor::dqmEndRun(edm::Run const&, edm::EventSetup const&) {
0496   if (processedLumis_.empty()) {
0497     return;
0498   }
0499 
0500   const double bigNumber = 1000000.;
0501   std::sort(processedLumis_.begin(), processedLumis_.end());
0502   int firstLumi = *processedLumis_.begin();
0503   int lastLumi = *(--processedLumis_.end());
0504 
0505   for (HistosContainer::iterator itH = histosMap_.begin(); itH != histosMap_.end(); itH++) {
0506     for (map<string, map<string, MonitorElement*> >::iterator itHH = itH->second.begin(); itHH != itH->second.end();
0507          itHH++) {
0508       double min = bigNumber;
0509       double max = -bigNumber;
0510       double minDelta = bigNumber;
0511       double maxDelta = -bigNumber;
0512       //      double minDeltaProf = bigNumber;
0513       //      double maxDeltaProf = -bigNumber;
0514       if (itHH->first != "run") {
0515         for (map<string, MonitorElement*>::iterator itHHH = itHH->second.begin(); itHHH != itHH->second.end();
0516              itHHH++) {
0517           if (itHHH->second != nullptr) {
0518             for (int bin = 1; bin <= itHHH->second->getTH1()->GetNbinsX(); bin++) {
0519               if (itHHH->second->getTH1()->GetBinError(bin) != 0 || itHHH->second->getTH1()->GetBinContent(bin) != 0) {
0520                 if (itHHH->first == "Lumibased BeamSpotFit" || itHHH->first == "Lumibased PrimaryVertex" ||
0521                     itHHH->first == "Lumibased DataBase" || itHHH->first == "Lumibased Scalers") {
0522                   if (min > itHHH->second->getTH1()->GetBinContent(bin)) {
0523                     min = itHHH->second->getTH1()->GetBinContent(bin);
0524                   }
0525                   if (max < itHHH->second->getTH1()->GetBinContent(bin)) {
0526                     max = itHHH->second->getTH1()->GetBinContent(bin);
0527                   }
0528                 } else if (itHHH->first == "Lumibased PrimaryVertex-DataBase fit" ||
0529                            itHHH->first == "Lumibased PrimaryVertex-Scalers fit" ||
0530                            itHHH->first == "Lumibased Scalers-DataBase fit" ||
0531                            itHHH->first == "Lumibased PrimaryVertex-DataBase" ||
0532                            itHHH->first == "Lumibased PrimaryVertex-Scalers") {
0533                   if (minDelta > itHHH->second->getTH1()->GetBinContent(bin)) {
0534                     minDelta = itHHH->second->getTH1()->GetBinContent(bin);
0535                   }
0536                   if (maxDelta < itHHH->second->getTH1()->GetBinContent(bin)) {
0537                     maxDelta = itHHH->second->getTH1()->GetBinContent(bin);
0538                   }
0539                 } else {
0540                   LogInfo("AlcaBeamMonitorClient") << "The histosMap_ have a histogram named " << itHHH->first
0541                                                    << " that I can't recognize in this loop!";
0542                   // assert(0);
0543                 }
0544               }
0545             }
0546           }
0547         }
0548         for (map<string, MonitorElement*>::iterator itHHH = itHH->second.begin(); itHHH != itHH->second.end();
0549              itHHH++) {
0550           if (itHHH->second != nullptr) {
0551             if (itHHH->first == "Lumibased BeamSpotFit" || itHHH->first == "Lumibased PrimaryVertex" ||
0552                 itHHH->first == "Lumibased DataBase" || itHHH->first == "Lumibased Scalers") {
0553               if ((max == -bigNumber && min == bigNumber) || max - min == 0) {
0554                 itHHH->second->getTH1()->SetMinimum(itHHH->second->getTH1()->GetMinimum() - 0.01);
0555                 itHHH->second->getTH1()->SetMaximum(itHHH->second->getTH1()->GetMaximum() + 0.01);
0556               } else {
0557                 itHHH->second->getTH1()->SetMinimum(min - 0.1 * (max - min));
0558                 itHHH->second->getTH1()->SetMaximum(max + 0.1 * (max - min));
0559               }
0560             } else if (itHHH->first == "Lumibased PrimaryVertex-DataBase fit" ||
0561                        itHHH->first == "Lumibased PrimaryVertex-Scalers fit" ||
0562                        itHHH->first == "Lumibased Scalers-DataBase fit" ||
0563                        itHHH->first == "Lumibased PrimaryVertex-DataBase" ||
0564                        itHHH->first == "Lumibased PrimaryVertex-Scalers") {
0565               if ((maxDelta == -bigNumber && minDelta == bigNumber) || maxDelta - minDelta == 0) {
0566                 itHHH->second->getTH1()->SetMinimum(itHHH->second->getTH1()->GetMinimum() - 0.01);
0567                 itHHH->second->getTH1()->SetMaximum(itHHH->second->getTH1()->GetMaximum() + 0.01);
0568               } else {
0569                 itHHH->second->getTH1()->SetMinimum(minDelta - 2 * (maxDelta - minDelta));
0570                 itHHH->second->getTH1()->SetMaximum(maxDelta + 2 * (maxDelta - minDelta));
0571               }
0572             } else {
0573               LogInfo("AlcaBeamMonitorClient") << "The histosMap_ have a histogram named " << itHHH->first
0574                                                << " that I can't recognize in this loop!";
0575             }
0576             itHHH->second->getTH1()->GetXaxis()->SetRangeUser(firstLumi - 0.5, lastLumi + 0.5);
0577           }
0578         }
0579       }
0580     }
0581   }
0582 }
0583 DEFINE_FWK_MODULE(AlcaBeamMonitor);