Back to home page

Project CMSSW displayed by LXR

 
 

    


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

0001 /*
0002  * \file FakeBeamMonitor.cc
0003  * \author Geng-yuan Jeng/UC Riverside
0004  *         Francisco Yumiceva/FNAL
0005  */
0006 
0007 /*
0008 The code has been modified for running average
0009 mode, and it gives results for the last NLS which is
0010 configurable.
0011 Sushil S. Chauhan /UCDavis
0012 Evan Friis        /UCDavis
0013 The last tag for working versions without this change is
0014 V00-03-25
0015 */
0016 
0017 #include "DQM/BeamMonitor/plugins/FakeBeamMonitor.h"
0018 #include "DQMServices/Core/interface/DQMStore.h"
0019 #include "FWCore/ServiceRegistry/interface/Service.h"
0020 #include "DataFormats/TrackCandidate/interface/TrackCandidate.h"
0021 #include "DataFormats/TrackCandidate/interface/TrackCandidateCollection.h"
0022 #include "DataFormats/Common/interface/View.h"
0023 #include "RecoVertex/BeamSpotProducer/interface/BSFitter.h"
0024 #include "FWCore/Framework/interface/MakerMacros.h"
0025 #include "FWCore/Framework/interface/LuminosityBlock.h"
0026 #include "FWCore/Framework/interface/Run.h"
0027 #include "FWCore/Framework/interface/ConsumesCollector.h"
0028 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0029 #include "FWCore/Common/interface/TriggerNames.h"
0030 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0031 #include "CondFormats/BeamSpotObjects/interface/BeamSpotOnlineObjects.h"
0032 #include <numeric>
0033 #include <cmath>
0034 #include <memory>
0035 #include <TMath.h>
0036 #include <TDatime.h>
0037 #include <iostream>
0038 #include <TStyle.h>
0039 #include <ctime>
0040 
0041 using namespace std;
0042 using namespace edm;
0043 
0044 void FakeBeamMonitor::formatFitTime(char* ts, const time_t& t) {
0045   //constexpr int CET(+1);
0046   constexpr int CEST(+2);
0047 
0048   //tm * ptm;
0049   //ptm = gmtime ( &t );
0050   //int year = ptm->tm_year;
0051 
0052   //get correct year from ctime
0053   time_t currentTime;
0054   struct tm* localTime;
0055   time(&currentTime);                   // Get the current time
0056   localTime = localtime(&currentTime);  // Convert the current time to the local time
0057   int year = localTime->tm_year + 1900;
0058 
0059   tm* ptm;
0060   ptm = gmtime(&t);
0061 
0062   //check if year is ok
0063   if (year <= 37)
0064     year += 2000;
0065   if (year >= 70 && year <= 137)
0066     year += 1900;
0067 
0068   if (year < 1995) {
0069     edm::LogError("BadTimeStamp") << "year reported is " << year << " !!" << std::endl;
0070     //year = 2015; //overwritten later by BeamFitter.cc for fits but needed here for TH1
0071     //edm::LogError("BadTimeStamp") << "Resetting to " <<year<<std::endl;
0072   }
0073   sprintf(ts,
0074           "%4d-%02d-%02d %02d:%02d:%02d",
0075           year,
0076           ptm->tm_mon + 1,
0077           ptm->tm_mday,
0078           (ptm->tm_hour + CEST) % 24,
0079           ptm->tm_min,
0080           ptm->tm_sec);
0081 
0082 #ifdef STRIP_TRAILING_BLANKS_IN_TIMEZONE
0083   unsigned int b = strlen(ts);
0084   while (ts[--b] == ' ') {
0085     ts[b] = 0;
0086   }
0087 #endif
0088 }
0089 
0090 static constexpr int buffTime = 23;
0091 
0092 std::string FakeBeamMonitor::getGMTstring(const time_t& timeToConvert) {
0093   char buff[32];
0094   std::strftime(buff, sizeof(buff), "%Y.%m.%d %H:%M:%S GMT", gmtime(&timeToConvert));
0095   std::string timeStr(buff);
0096   return timeStr;
0097 }
0098 
0099 //
0100 // constructors and destructor
0101 //
0102 FakeBeamMonitor::~FakeBeamMonitor() { delete rndm_; };
0103 
0104 FakeBeamMonitor::FakeBeamMonitor(const ParameterSet& ps)
0105     : dxBin_(ps.getParameter<int>("dxBin")),
0106       dxMin_(ps.getParameter<double>("dxMin")),
0107       dxMax_(ps.getParameter<double>("dxMax")),
0108 
0109       vxBin_(ps.getParameter<int>("vxBin")),
0110       vxMin_(ps.getParameter<double>("vxMin")),
0111       vxMax_(ps.getParameter<double>("vxMax")),
0112 
0113       phiBin_(ps.getParameter<int>("phiBin")),
0114       phiMin_(ps.getParameter<double>("phiMin")),
0115       phiMax_(ps.getParameter<double>("phiMax")),
0116 
0117       dzBin_(ps.getParameter<int>("dzBin")),
0118       dzMin_(ps.getParameter<double>("dzMin")),
0119       dzMax_(ps.getParameter<double>("dzMax")),
0120 
0121       countEvt_(0),
0122       countLumi_(0),
0123       nthBSTrk_(0),
0124       nFitElements_(3),
0125       resetHistos_(false),
0126       StartAverage_(false),
0127       firstAverageFit_(0),
0128       countGapLumi_(0) {
0129   monitorName_ = ps.getUntrackedParameter<string>("monitorName", "YourSubsystemName");
0130   recordName_ = ps.getUntrackedParameter<string>("recordName");
0131   intervalInSec_ = ps.getUntrackedParameter<int>("timeInterval", 920);  //40 LS X 23"
0132   fitNLumi_ = ps.getUntrackedParameter<int>("fitEveryNLumi", -1);
0133   resetFitNLumi_ = ps.getUntrackedParameter<int>("resetEveryNLumi", -1);
0134   fitPVNLumi_ = ps.getUntrackedParameter<int>("fitPVEveryNLumi", -1);
0135   resetPVNLumi_ = ps.getUntrackedParameter<int>("resetPVEveryNLumi", -1);
0136   deltaSigCut_ = ps.getUntrackedParameter<double>("deltaSignificanceCut", 15);
0137   debug_ = ps.getUntrackedParameter<bool>("Debug");
0138   onlineMode_ = ps.getUntrackedParameter<bool>("OnlineMode");
0139   min_Ntrks_ = ps.getParameter<ParameterSet>("BeamFitter").getUntrackedParameter<int>("MinimumInputTracks");
0140   maxZ_ = ps.getParameter<ParameterSet>("BeamFitter").getUntrackedParameter<double>("MaximumZ");
0141   minNrVertices_ = ps.getParameter<ParameterSet>("PVFitter").getUntrackedParameter<unsigned int>("minNrVerticesForFit");
0142   minVtxNdf_ = ps.getParameter<ParameterSet>("PVFitter").getUntrackedParameter<double>("minVertexNdf");
0143   minVtxWgt_ = ps.getParameter<ParameterSet>("PVFitter").getUntrackedParameter<double>("minVertexMeanWeight");
0144   useLockRecords_ = ps.getUntrackedParameter<bool>("useLockRecords");
0145 
0146   if (!monitorName_.empty())
0147     monitorName_ = monitorName_ + "/";
0148 
0149   if (fitNLumi_ <= 0)
0150     fitNLumi_ = 1;
0151   nFits_ = beginLumiOfBSFit_ = endLumiOfBSFit_ = beginLumiOfPVFit_ = endLumiOfPVFit_ = 0;
0152   refBStime[0] = refBStime[1] = refPVtime[0] = refPVtime[1] = 0;
0153   maxZ_ = std::fabs(maxZ_);
0154   lastlumi_ = 0;
0155   nextlumi_ = 0;
0156   processed_ = false;
0157 
0158   rndm_ = new TRandom3(0);
0159 }
0160 
0161 //--------------------------------------------------------
0162 namespace {
0163   /*The order of the enums is identical to the order in which
0164     MonitorElements are added to hs
0165    */
0166   enum Hists {
0167     k_x0_lumi,
0168     k_x0_lumi_all,
0169     k_y0_lumi,
0170     k_y0_lumi_all,
0171     k_z0_lumi,
0172     k_z0_lumi_all,
0173     k_sigmaX0_lumi,
0174     k_sigmaX0_lumi_all,
0175     k_sigmaY0_lumi,
0176     k_sigmaY0_lumi_all,
0177     k_sigmaZ0_lumi,
0178     k_sigmaZ0_lumi_all,
0179     k_x0_time,
0180     k_x0_time_all,
0181     k_y0_time,
0182     k_y0_time_all,
0183     k_z0_time,
0184     k_z0_time_all,
0185     k_sigmaX0_time,
0186     k_sigmaX0_time_all,
0187     k_sigmaY0_time,
0188     k_sigmaY0_time_all,
0189     k_sigmaZ0_time,
0190     k_sigmaZ0_time_all,
0191     k_PVx_lumi,
0192     k_PVx_lumi_all,
0193     k_PVy_lumi,
0194     k_PVy_lumi_all,
0195     k_PVz_lumi,
0196     k_PVz_lumi_all,
0197     k_PVx_time,
0198     k_PVx_time_all,
0199     k_PVy_time,
0200     k_PVy_time_all,
0201     k_PVz_time,
0202     k_PVz_time_all,
0203     kNumHists
0204   };
0205 }  // namespace
0206 
0207 void FakeBeamMonitor::dqmBeginRun(edm::Run const&, edm::EventSetup const&) {
0208   if (useLockRecords_ && onlineDbService_.isAvailable()) {
0209     onlineDbService_->lockRecords();
0210   }
0211 }
0212 
0213 void FakeBeamMonitor::bookHistograms(DQMStore::IBooker& iBooker, edm::Run const& iRun, edm::EventSetup const& iSetup) {
0214   frun = iRun.run();
0215   ftimestamp = iRun.beginTime().value();
0216   tmpTime = ftimestamp >> 32;
0217   startTime = refTime = tmpTime;
0218   char eventTime[64];
0219   formatFitTime(eventTime, tmpTime);
0220   edm::LogInfo("FakeBeamMonitor") << "TimeOffset = " << eventTime << std::endl;
0221   TDatime da(eventTime);
0222   if (debug_) {
0223     edm::LogInfo("FakeBeamMonitor") << "TimeOffset = ";
0224     da.Print();
0225   }
0226   auto daTime = da.Convert(kTRUE);
0227 
0228   // book some histograms here
0229 
0230   // create and cd into new folder
0231   iBooker.setCurrentFolder(monitorName_ + "Fit");
0232 
0233   h_nTrk_lumi = iBooker.book1D("nTrk_lumi", "Num. of selected tracks vs lumi (Fit)", 20, 0.5, 20.5);
0234   h_nTrk_lumi->setAxisTitle("Lumisection", 1);
0235   h_nTrk_lumi->setAxisTitle("Num of Tracks for Fit", 2);
0236 
0237   //store vtx vs lumi for monitoring why fits fail
0238   h_nVtx_lumi = iBooker.book1D("nVtx_lumi", "Num. of selected Vtx vs lumi (Fit)", 20, 0.5, 20.5);
0239   h_nVtx_lumi->setAxisTitle("Lumisection", 1);
0240   h_nVtx_lumi->setAxisTitle("Num of Vtx for Fit", 2);
0241 
0242   h_nVtx_lumi_all = iBooker.book1D("nVtx_lumi_all", "Num. of selected Vtx vs lumi (Fit) all", 20, 0.5, 20.5);
0243   h_nVtx_lumi_all->getTH1()->SetCanExtend(TH1::kAllAxes);
0244   h_nVtx_lumi_all->setAxisTitle("Lumisection", 1);
0245   h_nVtx_lumi_all->setAxisTitle("Num of Vtx for Fit", 2);
0246 
0247   h_d0_phi0 = iBooker.bookProfile(
0248       "d0_phi0", "d_{0} vs. #phi_{0} (Selected Tracks)", phiBin_, phiMin_, phiMax_, dxBin_, dxMin_, dxMax_, "");
0249   h_d0_phi0->setAxisTitle("#phi_{0} (rad)", 1);
0250   h_d0_phi0->setAxisTitle("d_{0} (cm)", 2);
0251 
0252   h_vx_vy = iBooker.book2D(
0253       "trk_vx_vy", "Vertex (PCA) position of selected tracks", vxBin_, vxMin_, vxMax_, vxBin_, vxMin_, vxMax_);
0254   h_vx_vy->setOption("COLZ");
0255   //   h_vx_vy->getTH1()->SetCanExtend(TH1::kAllAxes);
0256   h_vx_vy->setAxisTitle("x coordinate of input track at PCA (cm)", 1);
0257   h_vx_vy->setAxisTitle("y coordinate of input track at PCA (cm)", 2);
0258 
0259   {
0260     TDatime* da = new TDatime();
0261     gStyle->SetTimeOffset(da->Convert(kTRUE));
0262   }
0263 
0264   const int nvar_ = 6;
0265   string coord[nvar_] = {"x", "y", "z", "sigmaX", "sigmaY", "sigmaZ"};
0266   string label[nvar_] = {
0267       "x_{0} (cm)", "y_{0} (cm)", "z_{0} (cm)", "#sigma_{X_{0}} (cm)", "#sigma_{Y_{0}} (cm)", "#sigma_{Z_{0}} (cm)"};
0268 
0269   hs.clear();
0270   hs.reserve(kNumHists);
0271   for (int i = 0; i < 4; i++) {
0272     iBooker.setCurrentFolder(monitorName_ + "Fit");
0273     for (int ic = 0; ic < nvar_; ++ic) {
0274       TString histName(coord[ic]);
0275       TString histTitle(coord[ic]);
0276       string ytitle(label[ic]);
0277       string xtitle("");
0278       string options("E1");
0279       bool createHisto = true;
0280       switch (i) {
0281         case 1:  // BS vs time
0282           histName += "0_time";
0283           xtitle = "Time [UTC]";
0284           if (ic < 3)
0285             histTitle += " coordinate of beam spot vs time (Fit)";
0286           else
0287             histTitle = histTitle.Insert(5, " ") + " of beam spot vs time (Fit)";
0288           break;
0289         case 2:  // PV vs lumi
0290           if (ic < 3) {
0291             iBooker.setCurrentFolder(monitorName_ + "PrimaryVertex");
0292             histName.Insert(0, "PV");
0293             histName += "_lumi";
0294             histTitle.Insert(0, "Avg. ");
0295             histTitle += " position of primary vtx vs lumi";
0296             xtitle = "Lumisection";
0297             ytitle.insert(0, "PV");
0298             ytitle += " #pm #sigma_{PV";
0299             ytitle += coord[ic];
0300             ytitle += "} (cm)";
0301           } else
0302             createHisto = false;
0303           break;
0304         case 3:  // PV vs time
0305           if (ic < 3) {
0306             iBooker.setCurrentFolder(monitorName_ + "PrimaryVertex");
0307             histName.Insert(0, "PV");
0308             histName += "_time";
0309             histTitle.Insert(0, "Avg. ");
0310             histTitle += " position of primary vtx vs time";
0311             xtitle = "Time [UTC]";
0312             ytitle.insert(0, "PV");
0313             ytitle += " #pm #sigma_{PV";
0314             ytitle += coord[ic];
0315             ytitle += "} (cm)";
0316           } else
0317             createHisto = false;
0318           break;
0319         default:  // BS vs lumi
0320           histName += "0_lumi";
0321           xtitle = "Lumisection";
0322           if (ic < 3)
0323             histTitle += " coordinate of beam spot vs lumi (Fit)";
0324           else
0325             histTitle = histTitle.Insert(5, " ") + " of beam spot vs lumi (Fit)";
0326           break;
0327       }
0328       if (createHisto) {
0329         edm::LogInfo("FakeBeamMonitor") << "hitsName = " << histName << "; histTitle = " << histTitle << std::endl;
0330         auto tmpHs = iBooker.book1D(histName, histTitle, 40, 0.5, 40.5);
0331         hs.push_back(tmpHs);
0332         tmpHs->setAxisTitle(xtitle, 1);
0333         tmpHs->setAxisTitle(ytitle, 2);
0334         tmpHs->getTH1()->SetOption("E1");
0335         if (histName.Contains("time")) {
0336           //int nbins = (intervalInSec_/23 > 0 ? intervalInSec_/23 : 40);
0337           tmpHs->getTH1()->SetBins(intervalInSec_, 0.5, intervalInSec_ + 0.5);
0338           tmpHs->setAxisTimeDisplay(1);
0339           tmpHs->setAxisTimeFormat("%H:%M:%S", 1);
0340           tmpHs->getTH1()->GetXaxis()->SetTimeOffset(daTime);
0341         }
0342         histName += "_all";
0343         histTitle += " all";
0344         tmpHs = iBooker.book1D(histName, histTitle, 40, 0.5, 40.5);
0345         hs.push_back(tmpHs);
0346         tmpHs->getTH1()->SetCanExtend(TH1::kAllAxes);
0347         tmpHs->setAxisTitle(xtitle, 1);
0348         tmpHs->setAxisTitle(ytitle, 2);
0349         tmpHs->getTH1()->SetOption("E1");
0350         if (histName.Contains("time")) {
0351           //int nbins = (intervalInSec_/23 > 0 ? intervalInSec_/23 : 40);
0352           tmpHs->getTH1()->SetBins(intervalInSec_, 0.5, intervalInSec_ + 0.5);
0353           tmpHs->setAxisTimeDisplay(1);
0354           tmpHs->setAxisTimeFormat("%H:%M:%S", 1);
0355           tmpHs->getTH1()->GetXaxis()->SetTimeOffset(daTime);
0356         }
0357       }
0358     }
0359   }
0360   assert(hs.size() == kNumHists);
0361   assert(0 == strcmp(hs[k_sigmaY0_time]->getTH1()->GetName(), "sigmaY0_time"));
0362   assert(0 == strcmp(hs[k_PVz_lumi_all]->getTH1()->GetName(), "PVz_lumi_all"));
0363 
0364   iBooker.setCurrentFolder(monitorName_ + "Fit");
0365 
0366   h_trk_z0 = iBooker.book1D("trk_z0", "z_{0} of selected tracks", dzBin_, dzMin_, dzMax_);
0367   h_trk_z0->setAxisTitle("z_{0} of selected tracks (cm)", 1);
0368 
0369   h_vx_dz = iBooker.bookProfile(
0370       "vx_dz", "v_{x} vs. dz of selected tracks", dzBin_, dzMin_, dzMax_, dxBin_, dxMin_, dxMax_, "");
0371   h_vx_dz->setAxisTitle("dz (cm)", 1);
0372   h_vx_dz->setAxisTitle("x coordinate of input track at PCA (cm)", 2);
0373 
0374   h_vy_dz = iBooker.bookProfile(
0375       "vy_dz", "v_{y} vs. dz of selected tracks", dzBin_, dzMin_, dzMax_, dxBin_, dxMin_, dxMax_, "");
0376   h_vy_dz->setAxisTitle("dz (cm)", 1);
0377   h_vy_dz->setAxisTitle("y coordinate of input track at PCA (cm)", 2);
0378 
0379   h_x0 = iBooker.book1D("BeamMonitorFeedBack_x0", "x coordinate of beam spot (Fit)", 100, -0.01, 0.01);
0380   h_x0->setAxisTitle("x_{0} (cm)", 1);
0381   h_x0->getTH1()->SetCanExtend(TH1::kAllAxes);
0382 
0383   h_y0 = iBooker.book1D("BeamMonitorFeedBack_y0", "y coordinate of beam spot (Fit)", 100, -0.01, 0.01);
0384   h_y0->setAxisTitle("y_{0} (cm)", 1);
0385   h_y0->getTH1()->SetCanExtend(TH1::kAllAxes);
0386 
0387   h_z0 = iBooker.book1D("BeamMonitorFeedBack_z0", "z coordinate of beam spot (Fit)", dzBin_, dzMin_, dzMax_);
0388   h_z0->setAxisTitle("z_{0} (cm)", 1);
0389   h_z0->getTH1()->SetCanExtend(TH1::kAllAxes);
0390 
0391   h_sigmaX0 = iBooker.book1D("BeamMonitorFeedBack_sigmaX0", "sigma x0 of beam spot (Fit)", 100, 0, 0.05);
0392   h_sigmaX0->setAxisTitle("#sigma_{X_{0}} (cm)", 1);
0393   h_sigmaX0->getTH1()->SetCanExtend(TH1::kAllAxes);
0394 
0395   h_sigmaY0 = iBooker.book1D("BeamMonitorFeedBack_sigmaY0", "sigma y0 of beam spot (Fit)", 100, 0, 0.05);
0396   h_sigmaY0->setAxisTitle("#sigma_{Y_{0}} (cm)", 1);
0397   h_sigmaY0->getTH1()->SetCanExtend(TH1::kAllAxes);
0398 
0399   h_sigmaZ0 = iBooker.book1D("BeamMonitorFeedBack_sigmaZ0", "sigma z0 of beam spot (Fit)", 100, 0, 10);
0400   h_sigmaZ0->setAxisTitle("#sigma_{Z_{0}} (cm)", 1);
0401   h_sigmaZ0->getTH1()->SetCanExtend(TH1::kAllAxes);
0402 
0403   // Histograms of all reco tracks (without cuts):
0404   h_trkPt = iBooker.book1D("trkPt", "p_{T} of all reco'd tracks (no selection)", 200, 0., 50.);
0405   h_trkPt->setAxisTitle("p_{T} (GeV/c)", 1);
0406 
0407   h_trkVz = iBooker.book1D("trkVz", "Z coordinate of PCA of all reco'd tracks (no selection)", dzBin_, dzMin_, dzMax_);
0408   h_trkVz->setAxisTitle("V_{Z} (cm)", 1);
0409 
0410   cutFlowTable = iBooker.book1D("cutFlowTable", "Cut flow table of track selection", 9, 0, 9);
0411 
0412   // Results of previous good fit:
0413   fitResults = iBooker.book2D("fitResults", "Results of previous good beam fit", 2, 0, 2, 8, 0, 8);
0414   fitResults->setAxisTitle("Fitted Beam Spot (cm)", 1);
0415   fitResults->setBinLabel(8, "x_{0}", 2);
0416   fitResults->setBinLabel(7, "y_{0}", 2);
0417   fitResults->setBinLabel(6, "z_{0}", 2);
0418   fitResults->setBinLabel(5, "#sigma_{Z}", 2);
0419   fitResults->setBinLabel(4, "#frac{dx}{dz} (rad)", 2);
0420   fitResults->setBinLabel(3, "#frac{dy}{dz} (rad)", 2);
0421   fitResults->setBinLabel(2, "#sigma_{X}", 2);
0422   fitResults->setBinLabel(1, "#sigma_{Y}", 2);
0423   fitResults->setBinLabel(1, "Mean", 1);
0424   fitResults->setBinLabel(2, "Stat. Error", 1);
0425   fitResults->getTH1()->SetOption("text");
0426 
0427   // Histos of PrimaryVertices:
0428   iBooker.setCurrentFolder(monitorName_ + "PrimaryVertex");
0429 
0430   h_nVtx = iBooker.book1D("vtxNbr", "Reconstructed Vertices(non-fake) in all Event", 60, -0.5, 59.5);
0431   h_nVtx->setAxisTitle("Num. of reco. vertices", 1);
0432 
0433   //For one Trigger only
0434   h_nVtx_st = iBooker.book1D("vtxNbr_SelectedTriggers", "Reconstructed Vertices(non-fake) in Events", 60, -0.5, 59.5);
0435   //h_nVtx_st->setAxisTitle("Num. of reco. vertices for Un-Prescaled Jet Trigger",1);
0436 
0437   // Monitor only the PV with highest sum pt of assoc. trks:
0438   h_PVx[0] = iBooker.book1D("PVX", "x coordinate of Primary Vtx", 50, -0.01, 0.01);
0439   h_PVx[0]->setAxisTitle("PVx (cm)", 1);
0440   h_PVx[0]->getTH1()->SetCanExtend(TH1::kAllAxes);
0441 
0442   h_PVy[0] = iBooker.book1D("PVY", "y coordinate of Primary Vtx", 50, -0.01, 0.01);
0443   h_PVy[0]->setAxisTitle("PVy (cm)", 1);
0444   h_PVy[0]->getTH1()->SetCanExtend(TH1::kAllAxes);
0445 
0446   h_PVz[0] = iBooker.book1D("PVZ", "z coordinate of Primary Vtx", dzBin_, dzMin_, dzMax_);
0447   h_PVz[0]->setAxisTitle("PVz (cm)", 1);
0448 
0449   h_PVx[1] = iBooker.book1D("PVXFit", "x coordinate of Primary Vtx (Last Fit)", 50, -0.01, 0.01);
0450   h_PVx[1]->setAxisTitle("PVx (cm)", 1);
0451   h_PVx[1]->getTH1()->SetCanExtend(TH1::kAllAxes);
0452 
0453   h_PVy[1] = iBooker.book1D("PVYFit", "y coordinate of Primary Vtx (Last Fit)", 50, -0.01, 0.01);
0454   h_PVy[1]->setAxisTitle("PVy (cm)", 1);
0455   h_PVy[1]->getTH1()->SetCanExtend(TH1::kAllAxes);
0456 
0457   h_PVz[1] = iBooker.book1D("PVZFit", "z coordinate of Primary Vtx (Last Fit)", dzBin_, dzMin_, dzMax_);
0458   h_PVz[1]->setAxisTitle("PVz (cm)", 1);
0459 
0460   h_PVxz = iBooker.bookProfile("PVxz", "PVx vs. PVz", dzBin_ / 2, dzMin_, dzMax_, dxBin_ / 2, dxMin_, dxMax_, "");
0461   h_PVxz->setAxisTitle("PVz (cm)", 1);
0462   h_PVxz->setAxisTitle("PVx (cm)", 2);
0463 
0464   h_PVyz = iBooker.bookProfile("PVyz", "PVy vs. PVz", dzBin_ / 2, dzMin_, dzMax_, dxBin_ / 2, dxMin_, dxMax_, "");
0465   h_PVyz->setAxisTitle("PVz (cm)", 1);
0466   h_PVyz->setAxisTitle("PVy (cm)", 2);
0467 
0468   // Results of previous good fit:
0469   pvResults = iBooker.book2D("pvResults", "Results of fitting Primary Vertices", 2, 0, 2, 6, 0, 6);
0470   pvResults->setAxisTitle("Fitted Primary Vertex (cm)", 1);
0471   pvResults->setBinLabel(6, "PVx", 2);
0472   pvResults->setBinLabel(5, "PVy", 2);
0473   pvResults->setBinLabel(4, "PVz", 2);
0474   pvResults->setBinLabel(3, "#sigma_{X}", 2);
0475   pvResults->setBinLabel(2, "#sigma_{Y}", 2);
0476   pvResults->setBinLabel(1, "#sigma_{Z}", 2);
0477   pvResults->setBinLabel(1, "Mean", 1);
0478   pvResults->setBinLabel(2, "Stat. Error", 1);
0479   pvResults->getTH1()->SetOption("text");
0480 
0481   // Summary plots:
0482   iBooker.setCurrentFolder(monitorName_ + "EventInfo");
0483 
0484   reportSummary = iBooker.bookFloat("reportSummary");
0485   if (reportSummary)
0486     reportSummary->Fill(std::numeric_limits<double>::quiet_NaN());
0487 
0488   char histo[20];
0489   iBooker.setCurrentFolder(monitorName_ + "EventInfo/reportSummaryContents");
0490   for (int n = 0; n < nFitElements_; n++) {
0491     switch (n) {
0492       case 0:
0493         sprintf(histo, "x0_status");
0494         break;
0495       case 1:
0496         sprintf(histo, "y0_status");
0497         break;
0498       case 2:
0499         sprintf(histo, "z0_status");
0500         break;
0501     }
0502     reportSummaryContents[n] = iBooker.bookFloat(histo);
0503   }
0504 
0505   for (int i = 0; i < nFitElements_; i++) {
0506     summaryContent_[i] = 0.;
0507     reportSummaryContents[i]->Fill(std::numeric_limits<double>::quiet_NaN());
0508   }
0509 
0510   iBooker.setCurrentFolder(monitorName_ + "EventInfo");
0511 
0512   reportSummaryMap = iBooker.book2D("reportSummaryMap", "Beam Spot Summary Map", 1, 0, 1, 3, 0, 3);
0513   reportSummaryMap->setAxisTitle("", 1);
0514   reportSummaryMap->setAxisTitle("Fitted Beam Spot", 2);
0515   reportSummaryMap->setBinLabel(1, " ", 1);
0516   reportSummaryMap->setBinLabel(1, "x_{0}", 2);
0517   reportSummaryMap->setBinLabel(2, "y_{0}", 2);
0518   reportSummaryMap->setBinLabel(3, "z_{0}", 2);
0519   for (int i = 0; i < nFitElements_; i++) {
0520     reportSummaryMap->setBinContent(1, i + 1, -1.);
0521   }
0522 }
0523 
0524 //--------------------------------------------------------
0525 void FakeBeamMonitor::beginLuminosityBlock(const LuminosityBlock& lumiSeg, const EventSetup& context) {
0526   // start DB logger
0527   DBloggerReturn_ = 0;
0528   if (onlineDbService_.isAvailable()) {
0529     onlineDbService_->logger().start();
0530     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::beginLuminosityBlock - LS: " << lumiSeg.luminosityBlock()
0531                                          << " - Run: " << lumiSeg.getRun().run();
0532   }
0533 
0534   int nthlumi = lumiSeg.luminosityBlock();
0535   const edm::TimeValue_t fbegintimestamp = lumiSeg.beginTime().value();
0536   const std::time_t ftmptime = fbegintimestamp >> 32;
0537 
0538   if (countLumi_ == 0 && (!processed_)) {
0539     beginLumiOfBSFit_ = beginLumiOfPVFit_ = nthlumi;
0540     refBStime[0] = refPVtime[0] = ftmptime;
0541     mapBeginBSLS[countLumi_] = nthlumi;
0542     mapBeginPVLS[countLumi_] = nthlumi;
0543     mapBeginBSTime[countLumi_] = ftmptime;
0544     mapBeginPVTime[countLumi_] = ftmptime;
0545   }  //for the first record
0546 
0547   if (nthlumi > nextlumi_) {
0548     if (processed_) {
0549       countLumi_++;
0550       //store here them will need when we remove the first one of Last N LS
0551       mapBeginBSLS[countLumi_] = nthlumi;
0552       mapBeginPVLS[countLumi_] = nthlumi;
0553       mapBeginBSTime[countLumi_] = ftmptime;
0554       mapBeginPVTime[countLumi_] = ftmptime;
0555     }  //processed passed but not the first lumi
0556     if ((!processed_) && countLumi_ != 0) {
0557       mapBeginBSLS[countLumi_] = nthlumi;
0558       mapBeginPVLS[countLumi_] = nthlumi;
0559       mapBeginBSTime[countLumi_] = ftmptime;
0560       mapBeginPVTime[countLumi_] = ftmptime;
0561     }  //processed fails for last lumi
0562   }    //nthLumi > nextlumi
0563 
0564   if (StartAverage_) {
0565     //Just Make sure it get rest
0566     refBStime[0] = 0;
0567     refPVtime[0] = 0;
0568     beginLumiOfPVFit_ = 0;
0569     beginLumiOfBSFit_ = 0;
0570 
0571     if (debug_)
0572       edm::LogInfo("FakeBeamMonitor") << " beginLuminosityBlock:  Size of mapBeginBSLS before =  "
0573                                       << mapBeginBSLS.size() << endl;
0574     if (nthlumi >
0575         nextlumi_) {  //this make sure that it does not take into account this lumi for fitting and only look forward for new lumi
0576       //as countLumi also remains the same so map value  get overwritten once return to normal running.
0577       //even if few LS are misssing and DQM module do not sees them then it catchs up again
0578       map<int, int>::iterator itbs = mapBeginBSLS.begin();
0579       map<int, int>::iterator itpv = mapBeginPVLS.begin();
0580       map<int, std::time_t>::iterator itbstime = mapBeginBSTime.begin();
0581       map<int, std::time_t>::iterator itpvtime = mapBeginPVTime.begin();
0582 
0583       if (processed_) {  // otherwise if false then LS range of fit get messed up because we don't remove trk/pvs but we remove LS begin value . This prevent it as it happened if LS is there but no event are processed for some reason
0584         mapBeginBSLS.erase(itbs);
0585         mapBeginPVLS.erase(itpv);
0586         mapBeginBSTime.erase(itbstime);
0587         mapBeginPVTime.erase(itpvtime);
0588       }
0589       /*//not sure if want this or not ??
0590             map<int, int>::iterator itgapb=mapBeginBSLS.begin();
0591             map<int, int>::iterator itgape=mapBeginBSLS.end(); itgape--;
0592             countGapLumi_ = ( (itgape->second) - (itgapb->second) );
0593             //if we see Gap more than then 2*resetNFitLumi !!!!!!!
0594             //for example if 10-15 is fitted and if 16-25 are missing then we next fit will be for range 11-26 but BS can change in between
0595             // so better start  as fresh  and reset everything like starting in the begining!
0596             if(countGapLumi_ >= 2*resetFitNLumi_){RestartFitting(); mapBeginBSLS[countLumi_]   = nthlumi;}
0597             */
0598     }
0599 
0600     if (debug_)
0601       edm::LogInfo("FakeBeamMonitor") << " beginLuminosityBlock::  Size of mapBeginBSLS After = " << mapBeginBSLS.size()
0602                                       << endl;
0603 
0604     map<int, int>::iterator bbs = mapBeginBSLS.begin();
0605     map<int, int>::iterator bpv = mapBeginPVLS.begin();
0606     map<int, std::time_t>::iterator bbst = mapBeginBSTime.begin();
0607     map<int, std::time_t>::iterator bpvt = mapBeginPVTime.begin();
0608 
0609     if (beginLumiOfPVFit_ == 0)
0610       beginLumiOfPVFit_ = bpv->second;  //new begin time after removing the LS
0611     if (beginLumiOfBSFit_ == 0)
0612       beginLumiOfBSFit_ = bbs->second;
0613     if (refBStime[0] == 0)
0614       refBStime[0] = bbst->second;
0615     if (refPVtime[0] == 0)
0616       refPVtime[0] = bpvt->second;
0617 
0618   }  //same logic for average fit as above commented line
0619 
0620   map<int, std::time_t>::iterator nbbst = mapBeginBSTime.begin();
0621   map<int, std::time_t>::iterator nbpvt = mapBeginPVTime.begin();
0622 
0623   if (onlineMode_ && (nthlumi < nextlumi_))
0624     return;
0625 
0626   if (onlineMode_) {
0627     if (nthlumi > nextlumi_) {
0628       if (countLumi_ != 0 && processed_)
0629         FitAndFill(lumiSeg, lastlumi_, nextlumi_, nthlumi);
0630       nextlumi_ = nthlumi;
0631       edm::LogInfo("FakeBeamMonitor") << "beginLuminosityBlock:: Next Lumi to Fit: " << nextlumi_ << endl;
0632       if ((StartAverage_) && refBStime[0] == 0)
0633         refBStime[0] = nbbst->second;
0634       if ((StartAverage_) && refPVtime[0] == 0)
0635         refPVtime[0] = nbpvt->second;
0636     }
0637   } else {
0638     if (processed_)
0639       FitAndFill(lumiSeg, lastlumi_, nextlumi_, nthlumi);
0640     nextlumi_ = nthlumi;
0641     edm::LogInfo("FakeBeamMonitor") << " beginLuminosityBlock:: Next Lumi to Fit: " << nextlumi_ << endl;
0642     if ((StartAverage_) && refBStime[0] == 0)
0643       refBStime[0] = nbbst->second;
0644     if ((StartAverage_) && refPVtime[0] == 0)
0645       refPVtime[0] = nbpvt->second;
0646   }
0647 
0648   //countLumi_++;
0649   if (processed_)
0650     processed_ = false;
0651   edm::LogInfo("FakeBeamMonitor") << " beginLuminosityBlock::  Begin of Lumi: " << nthlumi << endl;
0652 }
0653 
0654 // ----------------------------------------------------------
0655 void FakeBeamMonitor::analyze(const Event& iEvent, const EventSetup& iSetup) {
0656   const int nthlumi = iEvent.luminosityBlock();
0657 
0658   if (onlineMode_ && (nthlumi < nextlumi_)) {
0659     edm::LogInfo("FakeBeamMonitor") << "analyze::  Spilt event from previous lumi section!" << std::endl;
0660     return;
0661   }
0662   if (onlineMode_ && (nthlumi > nextlumi_)) {
0663     edm::LogInfo("FakeBeamMonitor") << "analyze::  Spilt event from next lumi section!!!" << std::endl;
0664     return;
0665   }
0666 
0667   countEvt_++;
0668   //  theBeamFitter->readEvent(
0669   //      iEvent);  //Remember when track fitter read the event in the same place the PVFitter read the events !!!!!!!!!
0670 
0671   //  Handle<reco::BeamSpot> recoBeamSpotHandle;
0672   //  iEvent.getByToken(bsSrc_, recoBeamSpotHandle);
0673   //  refBS = *recoBeamSpotHandle;
0674 
0675   //  //------Cut Flow Table filled every event!--------------------------------------
0676   //  {
0677   //    // Make a copy of the cut flow table from the beam fitter.
0678   //    auto tmphisto = static_cast<TH1F*>(theBeamFitter->getCutFlow());
0679   //    cutFlowTable->getTH1()->SetBins(
0680   //        tmphisto->GetNbinsX(), tmphisto->GetXaxis()->GetXmin(), tmphisto->GetXaxis()->GetXmax());
0681   //    // Update the bin labels
0682   //    if (countEvt_ == 1)  // SetLabel just once
0683   //      for (int n = 0; n < tmphisto->GetNbinsX(); n++)
0684   //        cutFlowTable->setBinLabel(n + 1, tmphisto->GetXaxis()->GetBinLabel(n + 1), 1);
0685   //    cutFlowTable->Reset();
0686   //    cutFlowTable->getTH1()->Add(tmphisto);
0687   //  }
0688 
0689   //----Reco tracks -------------------------------------
0690   //  Handle<reco::TrackCollection> TrackCollection;
0691   //  iEvent.getByToken(tracksLabel_, TrackCollection);
0692   //  const reco::TrackCollection* tracks = TrackCollection.product();
0693   //  for (reco::TrackCollection::const_iterator track = tracks->begin(); track != tracks->end(); ++track) {
0694   //    h_trkPt->Fill(track->pt());  //no need to change  here for average bs
0695   //    h_trkVz->Fill(track->vz());
0696   //  }
0697 
0698   //-------HLT Trigger --------------------------------
0699   //  edm::Handle<TriggerResults> triggerResults;
0700   //  bool JetTrigPass = false;
0701   //  if (iEvent.getByToken(hltSrc_, triggerResults)) {
0702   //    const edm::TriggerNames& trigNames = iEvent.triggerNames(*triggerResults);
0703   //    for (unsigned int i = 0; i < triggerResults->size(); i++) {
0704   //      const std::string& trigName = trigNames.triggerName(i);
0705   //
0706   //      if (JetTrigPass)
0707   //        continue;
0708   //
0709   //      for (size_t t = 0; t < jetTrigger_.size(); ++t) {
0710   //        if (JetTrigPass)
0711   //          continue;
0712   //
0713   //        string string_search(jetTrigger_[t]);
0714   //        size_t found = trigName.find(string_search);
0715   //
0716   //        if (found != string::npos) {
0717   //          int thisTrigger_ = trigNames.triggerIndex(trigName);
0718   //          if (triggerResults->accept(thisTrigger_))
0719   //            JetTrigPass = true;
0720   //        }  //if trigger found
0721   //      }    //for(t=0;..)
0722   //    }      //for(i=0; ..)
0723   //  }        //if trigger colleciton exist)
0724 
0725   //------ Primary Vertices-------
0726   //  edm::Handle<reco::VertexCollection> PVCollection;
0727 
0728   //  if (iEvent.getByToken(pvSrc_, PVCollection)) {
0729   int nPVcount = 0;
0730   int nPVcount_ST = 0;  //For Single Trigger(hence ST)
0731 
0732   //    for (reco::VertexCollection::const_iterator pv = PVCollection->begin(); pv != PVCollection->end(); ++pv) {
0733   for (int tmp_idx = 0; tmp_idx < 10; tmp_idx++) {
0734     //--- vertex selection
0735     //    if (pv->isFake() || pv->tracksSize() == 0)
0736     //      continue;
0737     nPVcount++;  // count non fake pv:
0738 
0739     //if (JetTrigPass)
0740     nPVcount_ST++;  //non-fake pv with a specific trigger
0741 
0742     //    if (pv->ndof() < minVtxNdf_ || (pv->ndof() + 3.) / pv->tracksSize() < 2 * minVtxWgt_)
0743     //      continue;
0744 
0745     //Fill this map to store xyx for pv so that later we can remove the first one for run aver
0746     mapPVx[countLumi_].push_back(tmp_idx);
0747     mapPVy[countLumi_].push_back(tmp_idx);
0748     mapPVz[countLumi_].push_back(tmp_idx);
0749 
0750     //      if (!StartAverage_) {  //for first N LS
0751     //        h_PVx[0]->Fill(pv->x());
0752     //        h_PVy[0]->Fill(pv->y());
0753     //        h_PVz[0]->Fill(pv->z());
0754     //        h_PVxz->Fill(pv->z(), pv->x());
0755     //        h_PVyz->Fill(pv->z(), pv->y());
0756     //      }  //for first N LiS
0757     //      else {
0758     //        h_PVxz->Fill(pv->z(), pv->x());
0759     //        h_PVyz->Fill(pv->z(), pv->y());
0760     //      }
0761 
0762   }  //loop over pvs
0763 
0764   //    h_nVtx->Fill(nPVcount * 1.);  //no need to change it for average BS
0765 
0766   mapNPV[countLumi_].push_back((nPVcount_ST));
0767 
0768   //    if (!StartAverage_) {
0769   //      h_nVtx_st->Fill(nPVcount_ST * 1.);
0770   //    }
0771 
0772   //  }  //if pv collection is availaable
0773 
0774   if (StartAverage_) {
0775     map<int, std::vector<float> >::iterator itpvx = mapPVx.begin();
0776     map<int, std::vector<float> >::iterator itpvy = mapPVy.begin();
0777     map<int, std::vector<float> >::iterator itpvz = mapPVz.begin();
0778 
0779     map<int, std::vector<int> >::iterator itbspvinfo = mapNPV.begin();
0780 
0781     if ((int)mapPVx.size() > resetFitNLumi_) {  //sometimes the events is not there but LS is there!
0782       mapPVx.erase(itpvx);
0783       mapPVy.erase(itpvy);
0784       mapPVz.erase(itpvz);
0785       mapNPV.erase(itbspvinfo);
0786     }  //loop over Last N lumi collected
0787 
0788   }  //StartAverage==true
0789 
0790   processed_ = true;
0791 }
0792 
0793 //--------------------------------------------------------
0794 void FakeBeamMonitor::endLuminosityBlock(const LuminosityBlock& lumiSeg, const EventSetup& iSetup) {
0795   int nthlumi = lumiSeg.id().luminosityBlock();
0796   edm::LogInfo("FakeBeamMonitor") << "endLuminosityBlock:: Lumi of the last event before endLuminosityBlock: "
0797                                   << nthlumi << endl;
0798 
0799   if (onlineMode_ && nthlumi < nextlumi_)
0800     return;
0801   const edm::TimeValue_t fendtimestamp = lumiSeg.endTime().value();
0802   const std::time_t fendtime = fendtimestamp >> 32;
0803   tmpTime = refBStime[1] = refPVtime[1] = fendtime;
0804 
0805   // end DB logger
0806   if (onlineDbService_.isAvailable()) {
0807     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::endLuminosityBlock";
0808     onlineDbService_->logger().end(DBloggerReturn_);
0809   }
0810 }
0811 
0812 //--------------------------------------------------------
0813 void FakeBeamMonitor::FitAndFill(const LuminosityBlock& lumiSeg, int& lastlumi, int& nextlumi, int& nthlumi) {
0814   if (onlineMode_ && (nthlumi <= nextlumi))
0815     return;
0816 
0817   //set the correct run number when no event in the LS for fake output
0818   //  if ((processed_) && theBeamFitter->getRunNumber() != frun)
0819   //    theBeamFitter->setRun(frun);
0820 
0821   int currentlumi = nextlumi;
0822   edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  Lumi of the current fit: " << currentlumi << endl;
0823   lastlumi = currentlumi;
0824   endLumiOfBSFit_ = currentlumi;
0825   endLumiOfPVFit_ = currentlumi;
0826 
0827   //---------Fix for Runninv average-------------
0828   mapLSPVStoreSize[countLumi_] = 10;  //theBeamFitter->getPVvectorSize();
0829 
0830   //  if (StartAverage_) {
0831   //    std::map<int, std::size_t>::iterator rmLSPVi = mapLSPVStoreSize.begin();
0832   //    size_t SizeToRemovePV = rmLSPVi->second;
0833   //    for (std::map<int, std::size_t>::iterator rmLSPVe = mapLSPVStoreSize.end(); ++rmLSPVi != rmLSPVe;)
0834   //      rmLSPVi->second -= SizeToRemovePV;
0835   //
0836   //    theBeamFitter->resizePVvector(SizeToRemovePV);
0837   //
0838   //    map<int, std::size_t>::iterator tmpItpv = mapLSPVStoreSize.begin();
0839   //    mapLSPVStoreSize.erase(tmpItpv);
0840   //  }
0841   //  if (debug_)
0842   //    edm::LogInfo("BeamMonitor") << "FitAndFill:: Size of thePVvector After removing the PVs = "
0843   //                                << theBeamFitter->getPVvectorSize() << endl;
0844 
0845   //lets filt the PV for GUI here: It was in analyzer in preivous versiton but moved here due to absence of event in some lumis, works OK
0846   bool resetHistoFlag_ = false;
0847   if ((int)mapPVx.size() >= resetFitNLumi_ && (StartAverage_)) {
0848     h_PVx[0]->Reset();
0849     h_PVy[0]->Reset();
0850     h_PVz[0]->Reset();
0851     h_nVtx_st->Reset();
0852     resetHistoFlag_ = true;
0853   }
0854 
0855   int MaxPVs = 0;
0856   int countEvtLastNLS_ = 0;
0857   int countTotPV_ = 0;
0858 
0859   std::map<int, std::vector<int> >::iterator mnpv = mapNPV.begin();
0860   std::map<int, std::vector<float> >::iterator mpv2 = mapPVy.begin();
0861   std::map<int, std::vector<float> >::iterator mpv3 = mapPVz.begin();
0862 
0863   for (std::map<int, std::vector<float> >::iterator mpv1 = mapPVx.begin(); mpv1 != mapPVx.end();
0864        ++mpv1, ++mpv2, ++mpv3, ++mnpv) {
0865     std::vector<float>::iterator mpvs2 = (mpv2->second).begin();
0866     std::vector<float>::iterator mpvs3 = (mpv3->second).begin();
0867     for (std::vector<float>::iterator mpvs1 = (mpv1->second).begin(); mpvs1 != (mpv1->second).end();
0868          ++mpvs1, ++mpvs2, ++mpvs3) {
0869       if (resetHistoFlag_) {
0870         h_PVx[0]->Fill(*mpvs1);  //these histogram are reset after StartAverage_ flag is ON
0871         h_PVy[0]->Fill(*mpvs2);
0872         h_PVz[0]->Fill(*mpvs3);
0873       }
0874     }  //loop over second
0875 
0876     //Do the same here for nPV distr.
0877     for (std::vector<int>::iterator mnpvs = (mnpv->second).begin(); mnpvs != (mnpv->second).end(); ++mnpvs) {
0878       if ((*mnpvs > 0) && (resetHistoFlag_))
0879         h_nVtx_st->Fill((*mnpvs) * (1.0));
0880       countEvtLastNLS_++;
0881       countTotPV_ += (*mnpvs);
0882       if ((*mnpvs) > MaxPVs)
0883         MaxPVs = (*mnpvs);
0884     }  //loop over second of mapNPV
0885 
0886   }  //loop over last N lumis
0887 
0888   char tmpTitlePV[100];
0889   sprintf(tmpTitlePV, "%s %i %s %i", "Num. of reco. vertices for LS: ", beginLumiOfPVFit_, " to ", endLumiOfPVFit_);
0890   h_nVtx_st->setAxisTitle(tmpTitlePV, 1);
0891 
0892   std::vector<float> DipPVInfo_;
0893   DipPVInfo_.clear();
0894   //
0895   //  if (countTotPV_ != 0) {
0896   //    DipPVInfo_.push_back((float)countEvtLastNLS_);
0897   //    DipPVInfo_.push_back(h_nVtx_st->getMean());
0898   //    DipPVInfo_.push_back(h_nVtx_st->getMeanError());
0899   //    DipPVInfo_.push_back(h_nVtx_st->getRMS());
0900   //    DipPVInfo_.push_back(h_nVtx_st->getRMSError());
0901   //    DipPVInfo_.push_back((float)MaxPVs);
0902   //    DipPVInfo_.push_back((float)countTotPV_);
0903   //    MaxPVs = 0;
0904   //  } else {
0905   //    for (size_t i = 0; i < 7; i++) {
0906   //      if (i > 0) {
0907   //        DipPVInfo_.push_back(0.);
0908   //      } else {
0909   //        DipPVInfo_.push_back((float)countEvtLastNLS_);
0910   //      }
0911   //    }
0912   //  }
0913   //  theBeamFitter->SetPVInfo(DipPVInfo_);
0914   DipPVInfo_.push_back(rndm_->Gaus(1000., 100.));  // Events used
0915   DipPVInfo_.push_back(rndm_->Gaus(100., 10.));    // Mean PV
0916   DipPVInfo_.push_back(rndm_->Gaus(10., 5.));      // Mean PV err
0917   DipPVInfo_.push_back(rndm_->Gaus(10., 5.));      // Rms PV
0918   DipPVInfo_.push_back(rndm_->Gaus(5., 3.));       // Rms PV err
0919   DipPVInfo_.push_back(rndm_->Gaus(100., 10.));    // Max PVs
0920   countEvtLastNLS_ = 0;
0921 
0922   if (onlineMode_) {  // filling LS gap
0923     // FIXME: need to add protection for the case if the gap is at the resetting LS!
0924     const int countLS_bs = hs[k_x0_lumi]->getTH1()->GetEntries();
0925     const int countLS_pv = hs[k_PVx_lumi]->getTH1()->GetEntries();
0926     edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: countLS_bs = " << countLS_bs << " ; countLS_pv = " << countLS_pv
0927                                     << std::endl;
0928     int LSgap_bs = currentlumi / fitNLumi_ - countLS_bs;
0929     int LSgap_pv = currentlumi / fitPVNLumi_ - countLS_pv;
0930     if (currentlumi % fitNLumi_ == 0)
0931       LSgap_bs--;
0932     if (currentlumi % fitPVNLumi_ == 0)
0933       LSgap_pv--;
0934     edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  LSgap_bs = " << LSgap_bs << " ; LSgap_pv = " << LSgap_pv
0935                                     << std::endl;
0936     // filling previous fits if LS gap ever exists
0937     for (int ig = 0; ig < LSgap_bs; ig++) {
0938       hs[k_x0_lumi]->ShiftFillLast(0., 0., fitNLumi_);  //x0 , x0err, fitNLumi_;  see DQMCore....
0939       hs[k_y0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0940       hs[k_z0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0941       hs[k_sigmaX0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0942       hs[k_sigmaY0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0943       hs[k_sigmaZ0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0944       h_nVtx_lumi->ShiftFillLast(0., 0., fitNLumi_);
0945     }
0946     for (int ig = 0; ig < LSgap_pv; ig++) {
0947       hs[k_PVx_lumi]->ShiftFillLast(0., 0., fitPVNLumi_);
0948       hs[k_PVy_lumi]->ShiftFillLast(0., 0., fitPVNLumi_);
0949       hs[k_PVz_lumi]->ShiftFillLast(0., 0., fitPVNLumi_);
0950     }
0951     const int previousLS = h_nTrk_lumi->getTH1()->GetEntries();
0952     for (int i = 1; i < (currentlumi - previousLS);
0953          i++)  //if (current-previoius)= 1 then never go inside the for loop!!!!!!!!!!!
0954       h_nTrk_lumi->ShiftFillLast(nthBSTrk_);
0955   }
0956 
0957   edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: Time lapsed since last scroll = " << tmpTime - refTime << std::endl;
0958 
0959   if (testScroll(tmpTime, refTime)) {
0960     scrollTH1(hs[k_x0_time]->getTH1(), refTime);
0961     scrollTH1(hs[k_y0_time]->getTH1(), refTime);
0962     scrollTH1(hs[k_z0_time]->getTH1(), refTime);
0963     scrollTH1(hs[k_sigmaX0_time]->getTH1(), refTime);
0964     scrollTH1(hs[k_sigmaY0_time]->getTH1(), refTime);
0965     scrollTH1(hs[k_sigmaZ0_time]->getTH1(), refTime);
0966     scrollTH1(hs[k_PVx_time]->getTH1(), refTime);
0967     scrollTH1(hs[k_PVy_time]->getTH1(), refTime);
0968     scrollTH1(hs[k_PVz_time]->getTH1(), refTime);
0969   }
0970 
0971   //  bool doPVFit = false;
0972   //
0973   //  if (fitPVNLumi_ > 0) {
0974   //    if (onlineMode_) {
0975   //      if (currentlumi % fitPVNLumi_ == 0)
0976   //        doPVFit = true;
0977   //    } else if (countLumi_ % fitPVNLumi_ == 0)
0978   //      doPVFit = true;
0979   //  } else
0980   //    doPVFit = true;
0981   //
0982   //  if (doPVFit) {
0983   edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: Do PV Fitting for LS = " << beginLumiOfPVFit_ << " to "
0984                                   << endLumiOfPVFit_ << std::endl;
0985   // Primary Vertex Fit:
0986   if (h_PVx[0]->getTH1()->GetEntries() > minNrVertices_) {
0987     pvResults->Reset();
0988     char tmpTitle[50];
0989     sprintf(tmpTitle, "%s %i %s %i", "Fitted Primary Vertex (cm) of LS: ", beginLumiOfPVFit_, " to ", endLumiOfPVFit_);
0990     pvResults->setAxisTitle(tmpTitle, 1);
0991 
0992     std::unique_ptr<TF1> fgaus{new TF1("fgaus", "gaus")};
0993     double mean, width, meanErr, widthErr;
0994     fgaus->SetLineColor(4);
0995     h_PVx[0]->getTH1()->Fit(fgaus.get(), "QLM0");
0996     mean = fgaus->GetParameter(1);
0997     width = fgaus->GetParameter(2);
0998     meanErr = fgaus->GetParError(1);
0999     widthErr = fgaus->GetParError(2);
1000 
1001     hs[k_PVx_lumi]->ShiftFillLast(mean, width, fitPVNLumi_);
1002     hs[k_PVx_lumi_all]->setBinContent(currentlumi, mean);
1003     hs[k_PVx_lumi_all]->setBinError(currentlumi, width);
1004     int nthBin = tmpTime - refTime;
1005     if (nthBin < 0)
1006       edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  Event time outside current range of time histograms!"
1007                                       << std::endl;
1008     if (nthBin > 0) {
1009       hs[k_PVx_time]->setBinContent(nthBin, mean);
1010       hs[k_PVx_time]->setBinError(nthBin, width);
1011     }
1012     int jthBin = tmpTime - startTime;
1013     if (jthBin > 0) {
1014       hs[k_PVx_time_all]->setBinContent(jthBin, mean);
1015       hs[k_PVx_time_all]->setBinError(jthBin, width);
1016     }
1017     pvResults->setBinContent(1, 6, mean);
1018     pvResults->setBinContent(1, 3, width);
1019     pvResults->setBinContent(2, 6, meanErr);
1020     pvResults->setBinContent(2, 3, widthErr);
1021 
1022     {
1023       // snap shot of the fit
1024       auto tmphisto = h_PVx[0]->getTH1F();
1025       h_PVx[1]->getTH1()->SetBins(
1026           tmphisto->GetNbinsX(), tmphisto->GetXaxis()->GetXmin(), tmphisto->GetXaxis()->GetXmax());
1027       h_PVx[1]->Reset();
1028       h_PVx[1]->getTH1()->Add(tmphisto);
1029       h_PVx[1]->getTH1()->Fit(fgaus.get(), "QLM");
1030     }
1031 
1032     h_PVy[0]->getTH1()->Fit(fgaus.get(), "QLM0");
1033     mean = fgaus->GetParameter(1);
1034     width = fgaus->GetParameter(2);
1035     meanErr = fgaus->GetParError(1);
1036     widthErr = fgaus->GetParError(2);
1037     hs[k_PVy_lumi]->ShiftFillLast(mean, width, fitPVNLumi_);
1038     hs[k_PVy_lumi_all]->setBinContent(currentlumi, mean);
1039     hs[k_PVy_lumi_all]->setBinError(currentlumi, width);
1040     if (nthBin > 0) {
1041       hs[k_PVy_time]->setBinContent(nthBin, mean);
1042       hs[k_PVy_time]->setBinError(nthBin, width);
1043     }
1044     if (jthBin > 0) {
1045       hs[k_PVy_time_all]->setBinContent(jthBin, mean);
1046       hs[k_PVy_time_all]->setBinError(jthBin, width);
1047     }
1048     pvResults->setBinContent(1, 5, mean);
1049     pvResults->setBinContent(1, 2, width);
1050     pvResults->setBinContent(2, 5, meanErr);
1051     pvResults->setBinContent(2, 2, widthErr);
1052     // snap shot of the fit
1053     {
1054       auto tmphisto = h_PVy[0]->getTH1F();
1055       h_PVy[1]->getTH1()->SetBins(
1056           tmphisto->GetNbinsX(), tmphisto->GetXaxis()->GetXmin(), tmphisto->GetXaxis()->GetXmax());
1057       h_PVy[1]->Reset();
1058       h_PVy[1]->getTH1()->Add(tmphisto);
1059       h_PVy[1]->getTH1()->Fit(fgaus.get(), "QLM");
1060     }
1061 
1062     h_PVz[0]->getTH1()->Fit(fgaus.get(), "QLM0");
1063     mean = fgaus->GetParameter(1);
1064     width = fgaus->GetParameter(2);
1065     meanErr = fgaus->GetParError(1);
1066     widthErr = fgaus->GetParError(2);
1067     hs[k_PVz_lumi]->ShiftFillLast(mean, width, fitPVNLumi_);
1068     hs[k_PVz_lumi_all]->setBinContent(currentlumi, mean);
1069     hs[k_PVz_lumi_all]->setBinError(currentlumi, width);
1070     if (nthBin > 0) {
1071       hs[k_PVz_time]->setBinContent(nthBin, mean);
1072       hs[k_PVz_time]->setBinError(nthBin, width);
1073     }
1074     if (jthBin > 0) {
1075       hs[k_PVz_time_all]->setBinContent(jthBin, mean);
1076       hs[k_PVz_time_all]->setBinError(jthBin, width);
1077     }
1078     pvResults->setBinContent(1, 4, mean);
1079     pvResults->setBinContent(1, 1, width);
1080     pvResults->setBinContent(2, 4, meanErr);
1081     pvResults->setBinContent(2, 1, widthErr);
1082     {
1083       // snap shot of the fit
1084       auto tmphisto = h_PVz[0]->getTH1F();
1085       h_PVz[1]->getTH1()->SetBins(
1086           tmphisto->GetNbinsX(), tmphisto->GetXaxis()->GetXmin(), tmphisto->GetXaxis()->GetXmax());
1087       h_PVz[1]->Reset();
1088       h_PVz[1]->getTH1()->Add(tmphisto);
1089       h_PVz[1]->getTH1()->Fit(fgaus.get(), "QLM");
1090     }
1091   }  //check if found min Vertices
1092      //  }    //do PVfit
1093 
1094   if ((resetPVNLumi_ > 0 && countLumi_ == resetPVNLumi_) || StartAverage_) {
1095     beginLumiOfPVFit_ = 0;
1096     refPVtime[0] = 0;
1097   }
1098 
1099   //---------Readjustment of theBSvector, RefTime, beginLSofFit---------
1100   //  vector<BSTrkParameters> theBSvector1 = theBeamFitter->getBSvector();
1101   //  mapLSBSTrkSize[countLumi_] = (theBSvector1.size());
1102   size_t PreviousRecords = 0;  //needed to fill nth record of tracks in GUI
1103 
1104   //  if (StartAverage_) {
1105   //    size_t SizeToRemove = 0;
1106   //    std::map<int, std::size_t>::iterator rmls = mapLSBSTrkSize.begin();
1107   //    SizeToRemove = rmls->second;
1108   //    if (debug_)
1109   //      edm::LogInfo("BeamMonitor") << "  The size to remove is =  " << SizeToRemove << endl;
1110   //    int changedAfterThis = 0;
1111   //    for (std::map<int, std::size_t>::iterator rmLS = mapLSBSTrkSize.begin(); rmLS != mapLSBSTrkSize.end();
1112   //         ++rmLS, ++changedAfterThis) {
1113   //      if (changedAfterThis > 0) {
1114   //        (rmLS->second) = (rmLS->second) - SizeToRemove;
1115   //        if ((mapLSBSTrkSize.size() - (size_t)changedAfterThis) == 2)
1116   //          PreviousRecords = (rmLS->second);
1117   //      }
1118   //    }
1119   //
1120   //    theBeamFitter->resizeBSvector(SizeToRemove);
1121   //
1122   //    map<int, std::size_t>::iterator tmpIt = mapLSBSTrkSize.begin();
1123   //    mapLSBSTrkSize.erase(tmpIt);
1124   //
1125   //    std::pair<int, int> checkfitLS = theBeamFitter->getFitLSRange();
1126   //    std::pair<time_t, time_t> checkfitTime = theBeamFitter->getRefTime();
1127   //    theBeamFitter->setFitLSRange(beginLumiOfBSFit_, checkfitLS.second);
1128   //    theBeamFitter->setRefTime(refBStime[0], checkfitTime.second);
1129   //  }
1130 
1131   //Fill the track for this fit
1132   //  vector<BSTrkParameters> theBSvector = theBeamFitter->getBSvector();
1133   //  h_nTrk_lumi->ShiftFillLast(theBSvector.size());
1134   //
1135   //  if (debug_)
1136   //    edm::LogInfo("BeamMonitor") << "FitAndFill::   Size of  theBSViector.size()  After =" << theBSvector.size() << endl;
1137 
1138   //  bool countFitting = false;
1139   //  if (theBSvector.size() >= PreviousRecords && theBSvector.size() >= min_Ntrks_) {
1140   //    countFitting = true;
1141   //  }
1142 
1143   //---Fix for Cut Flow Table for Running average in a same way//the previous code  has problem for resetting!!!
1144   //  mapLSCF[countLumi_] = *theBeamFitter->getCutFlow();
1145   //  if (StartAverage_ && !mapLSCF.empty()) {
1146   //    const TH1F& cutFlowToSubtract = mapLSCF.begin()->second;
1147   //    // Subtract the last cut flow from all of the others.
1148   //    std::map<int, TH1F>::iterator cf = mapLSCF.begin();
1149   //    // Start on second entry
1150   //    for (; cf != mapLSCF.end(); ++cf) {
1151   //      cf->second.Add(&cutFlowToSubtract, -1);
1152   //    }
1153   //    theBeamFitter->subtractFromCutFlow(&cutFlowToSubtract);
1154   //    // Remove the obsolete lumi section
1155   //    mapLSCF.erase(mapLSCF.begin());
1156   //  }
1157 
1158   if (resetHistos_) {
1159     h_d0_phi0->Reset();
1160     h_vx_vy->Reset();
1161     h_vx_dz->Reset();
1162     h_vy_dz->Reset();
1163     h_trk_z0->Reset();
1164     resetHistos_ = false;
1165   }
1166 
1167   if (StartAverage_)
1168     nthBSTrk_ = PreviousRecords;  //after average proccess is ON//for 2-6 LS fit PreviousRecords is size from 2-5 LS
1169 
1170   edm::LogInfo("FakeBeamMonitor") << " The Previous Recored for this fit is  =" << nthBSTrk_ << endl;
1171 
1172   //  unsigned int itrk = 0;
1173   //  for (vector<BSTrkParameters>::const_iterator BSTrk = theBSvector.begin(); BSTrk != theBSvector.end();
1174   //       ++BSTrk, ++itrk) {
1175   //    if (itrk >= nthBSTrk_) {  //fill for this record only !!
1176   //      h_d0_phi0->Fill(BSTrk->phi0(), BSTrk->d0());
1177   //      double vx = BSTrk->vx();
1178   //      double vy = BSTrk->vy();
1179   //      double z0 = BSTrk->z0();
1180   //      h_vx_vy->Fill(vx, vy);
1181   //      h_vx_dz->Fill(z0, vx);
1182   //      h_vy_dz->Fill(z0, vy);
1183   //      h_trk_z0->Fill(z0);
1184   //    }
1185   //  }
1186 
1187   //  nthBSTrk_ = theBSvector.size();  // keep track of num of tracks filled so far
1188 
1189   edm::LogInfo("FakeBeamMonitor") << " The Current Recored for this fit is  =" << nthBSTrk_ << endl;
1190 
1191   //  if (countFitting)
1192   //    edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  Num of tracks collected = " << nthBSTrk_ << endl;
1193 
1194   if (fitNLumi_ > 0) {
1195     if (onlineMode_) {
1196       if (currentlumi % fitNLumi_ != 0) {
1197         //  for (std::map<TString,MonitorElement*>::iterator itAll = hs.begin();
1198         //       itAll != hs.end(); ++itAll) {
1199         //    if ((*itAll).first.Contains("all")) {
1200         //      (*itAll).second->setBinContent(currentlumi,0.);
1201         //      (*itAll).second->setBinError(currentlumi,0.);
1202         //    }
1203         //  }
1204         return;
1205       }
1206     } else if (countLumi_ % fitNLumi_ != 0)
1207       return;
1208   }
1209 
1210   edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: [DebugTime] refBStime[0] = " << refBStime[0]
1211                                   << "; address =  " << &refBStime[0] << std::endl;
1212   edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: [DebugTime] refBStime[1] = " << refBStime[1]
1213                                   << "; address =  " << &refBStime[1] << std::endl;
1214 
1215   //Fill for all LS even if fit fails
1216   //  h_nVtx_lumi->ShiftFillLast((theBeamFitter->getPVvectorSize()), 0., fitNLumi_);
1217   //  h_nVtx_lumi_all->setBinContent(currentlumi, (theBeamFitter->getPVvectorSize()));
1218 
1219   //  if (countFitting) {
1220   nFits_++;
1221   //    std::pair<int, int> fitLS = theBeamFitter->getFitLSRange();
1222   std::pair<int, int> fitLS(beginLumiOfBSFit_, endLumiOfBSFit_);
1223   //    edm::LogInfo("BeamMonitor") << "FitAndFill::  [BeamFitter] Do BeamSpot Fit for LS = " << fitLS.first << " to "
1224   //                                << fitLS.second << std::endl;
1225   edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  [FakeBeamMonitor] Do BeamSpot Fit for LS = " << beginLumiOfBSFit_
1226                                   << " to " << endLumiOfBSFit_ << std::endl;
1227 
1228   //Now Run the PV and Track Fitter over the collected tracks and pvs
1229   //    if (theBeamFitter->runPVandTrkFitter()) {
1230   //      reco::BeamSpot bs = theBeamFitter->getBeamSpot();
1231 
1232   //Create fake BS here
1233   reco::BeamSpot::CovarianceMatrix matrix;
1234   for (int j = 0; j < 7; ++j) {
1235     for (int k = j; k < 7; ++k) {
1236       matrix(j, k) = 0;
1237     }
1238   }
1239 
1240   // random values for fake BeamSpot
1241   float tmp_BSx = rndm_->Gaus(0.1, 0.1);            // [cm]
1242   float tmp_BSy = rndm_->Gaus(0.1, 0.1);            // [cm]
1243   float tmp_BSz = rndm_->Gaus(0.1, 0.1);            // [cm]
1244   float tmp_BSwidthX = rndm_->Gaus(0.001, 0.0005);  // [cm]
1245   float tmp_BSwidthY = rndm_->Gaus(0.001, 0.0005);  // [cm]
1246   float tmp_BSwidthZ = rndm_->Gaus(3.5, 0.5);       // [cm]
1247 
1248   reco::BeamSpot bs(reco::BeamSpot::Point(tmp_BSx, tmp_BSy, tmp_BSz),
1249                     tmp_BSwidthZ,
1250                     0,
1251                     0,
1252                     tmp_BSwidthX,
1253                     matrix,
1254                     reco::BeamSpot::Tracker);
1255   bs.setBeamWidthY(tmp_BSwidthY);
1256 
1257   if (bs.type() > 0)  // with good beamwidth fit
1258     preBS = bs;       // cache good fit results
1259 
1260   edm::LogInfo("FakeBeamMonitor") << "\n RESULTS OF DEFAULT FIT:" << endl;
1261   edm::LogInfo("FakeBeamMonitor") << bs << endl;
1262   edm::LogInfo("FakeBeamMonitor") << "[BeamFitter] fitting done \n" << endl;
1263 
1264   hs[k_x0_lumi]->ShiftFillLast(bs.x0(), bs.x0Error(), fitNLumi_);
1265   hs[k_y0_lumi]->ShiftFillLast(bs.y0(), bs.y0Error(), fitNLumi_);
1266   hs[k_z0_lumi]->ShiftFillLast(bs.z0(), bs.z0Error(), fitNLumi_);
1267   hs[k_sigmaX0_lumi]->ShiftFillLast(bs.BeamWidthX(), bs.BeamWidthXError(), fitNLumi_);
1268   hs[k_sigmaY0_lumi]->ShiftFillLast(bs.BeamWidthY(), bs.BeamWidthYError(), fitNLumi_);
1269   hs[k_sigmaZ0_lumi]->ShiftFillLast(bs.sigmaZ(), bs.sigmaZ0Error(), fitNLumi_);
1270   hs[k_x0_lumi_all]->setBinContent(currentlumi, bs.x0());
1271   hs[k_x0_lumi_all]->setBinError(currentlumi, bs.x0Error());
1272   hs[k_y0_lumi_all]->setBinContent(currentlumi, bs.y0());
1273   hs[k_y0_lumi_all]->setBinError(currentlumi, bs.y0Error());
1274   hs[k_z0_lumi_all]->setBinContent(currentlumi, bs.z0());
1275   hs[k_z0_lumi_all]->setBinError(currentlumi, bs.z0Error());
1276   hs[k_sigmaX0_lumi_all]->setBinContent(currentlumi, bs.BeamWidthX());
1277   hs[k_sigmaX0_lumi_all]->setBinError(currentlumi, bs.BeamWidthXError());
1278   hs[k_sigmaY0_lumi_all]->setBinContent(currentlumi, bs.BeamWidthY());
1279   hs[k_sigmaY0_lumi_all]->setBinError(currentlumi, bs.BeamWidthYError());
1280   hs[k_sigmaZ0_lumi_all]->setBinContent(currentlumi, bs.sigmaZ());
1281   hs[k_sigmaZ0_lumi_all]->setBinError(currentlumi, bs.sigmaZ0Error());
1282 
1283   int nthBin = tmpTime - refTime;
1284   if (nthBin > 0) {
1285     hs[k_x0_time]->setBinContent(nthBin, bs.x0());
1286     hs[k_y0_time]->setBinContent(nthBin, bs.y0());
1287     hs[k_z0_time]->setBinContent(nthBin, bs.z0());
1288     hs[k_sigmaX0_time]->setBinContent(nthBin, bs.BeamWidthX());
1289     hs[k_sigmaY0_time]->setBinContent(nthBin, bs.BeamWidthY());
1290     hs[k_sigmaZ0_time]->setBinContent(nthBin, bs.sigmaZ());
1291     hs[k_x0_time]->setBinError(nthBin, bs.x0Error());
1292     hs[k_y0_time]->setBinError(nthBin, bs.y0Error());
1293     hs[k_z0_time]->setBinError(nthBin, bs.z0Error());
1294     hs[k_sigmaX0_time]->setBinError(nthBin, bs.BeamWidthXError());
1295     hs[k_sigmaY0_time]->setBinError(nthBin, bs.BeamWidthYError());
1296     hs[k_sigmaZ0_time]->setBinError(nthBin, bs.sigmaZ0Error());
1297   }
1298 
1299   int jthBin = tmpTime - startTime;
1300   if (jthBin > 0) {
1301     hs[k_x0_time_all]->setBinContent(jthBin, bs.x0());
1302     hs[k_y0_time_all]->setBinContent(jthBin, bs.y0());
1303     hs[k_z0_time_all]->setBinContent(jthBin, bs.z0());
1304     hs[k_sigmaX0_time_all]->setBinContent(jthBin, bs.BeamWidthX());
1305     hs[k_sigmaY0_time_all]->setBinContent(jthBin, bs.BeamWidthY());
1306     hs[k_sigmaZ0_time_all]->setBinContent(jthBin, bs.sigmaZ());
1307     hs[k_x0_time_all]->setBinError(jthBin, bs.x0Error());
1308     hs[k_y0_time_all]->setBinError(jthBin, bs.y0Error());
1309     hs[k_z0_time_all]->setBinError(jthBin, bs.z0Error());
1310     hs[k_sigmaX0_time_all]->setBinError(jthBin, bs.BeamWidthXError());
1311     hs[k_sigmaY0_time_all]->setBinError(jthBin, bs.BeamWidthYError());
1312     hs[k_sigmaZ0_time_all]->setBinError(jthBin, bs.sigmaZ0Error());
1313   }
1314 
1315   h_x0->Fill(bs.x0());
1316   h_y0->Fill(bs.y0());
1317   h_z0->Fill(bs.z0());
1318   if (bs.type() > 0) {  // with good beamwidth fit
1319     h_sigmaX0->Fill(bs.BeamWidthX());
1320     h_sigmaY0->Fill(bs.BeamWidthY());
1321   }
1322   h_sigmaZ0->Fill(bs.sigmaZ());
1323 
1324   if (nthBSTrk_ >= 2 * min_Ntrks_) {
1325     double amp = std::sqrt(bs.x0() * bs.x0() + bs.y0() * bs.y0());
1326     double alpha = std::atan2(bs.y0(), bs.x0());
1327     std::unique_ptr<TF1> f1{new TF1("f1", "[0]*sin(x-[1])", -3.14, 3.14)};
1328     f1->SetParameters(amp, alpha);
1329     f1->SetParLimits(0, amp - 0.1, amp + 0.1);
1330     f1->SetParLimits(1, alpha - 0.577, alpha + 0.577);
1331     f1->SetLineColor(4);
1332     h_d0_phi0->getTProfile()->Fit(f1.get(), "QR");
1333 
1334     double mean = bs.z0();
1335     double width = bs.sigmaZ();
1336     std::unique_ptr<TF1> fgaus{new TF1("fgaus", "gaus")};
1337     fgaus->SetParameters(mean, width);
1338     fgaus->SetLineColor(4);
1339     h_trk_z0->getTH1()->Fit(fgaus.get(), "QLRM", "", mean - 3 * width, mean + 3 * width);
1340   }
1341 
1342   fitResults->Reset();
1343   std::pair<int, int> LSRange(beginLumiOfBSFit_, endLumiOfBSFit_);  //= theBeamFitter->getFitLSRange();
1344   char tmpTitle[50];
1345   sprintf(tmpTitle, "%s %i %s %i", "Fitted Beam Spot (cm) of LS: ", LSRange.first, " to ", LSRange.second);
1346   fitResults->setAxisTitle(tmpTitle, 1);
1347   fitResults->setBinContent(1, 8, bs.x0());
1348   fitResults->setBinContent(1, 7, bs.y0());
1349   fitResults->setBinContent(1, 6, bs.z0());
1350   fitResults->setBinContent(1, 5, bs.sigmaZ());
1351   fitResults->setBinContent(1, 4, bs.dxdz());
1352   fitResults->setBinContent(1, 3, bs.dydz());
1353   if (bs.type() > 0) {  // with good beamwidth fit
1354     fitResults->setBinContent(1, 2, bs.BeamWidthX());
1355     fitResults->setBinContent(1, 1, bs.BeamWidthY());
1356   } else {  // fill cached widths
1357     fitResults->setBinContent(1, 2, preBS.BeamWidthX());
1358     fitResults->setBinContent(1, 1, preBS.BeamWidthY());
1359   }
1360 
1361   fitResults->setBinContent(2, 8, bs.x0Error());
1362   fitResults->setBinContent(2, 7, bs.y0Error());
1363   fitResults->setBinContent(2, 6, bs.z0Error());
1364   fitResults->setBinContent(2, 5, bs.sigmaZ0Error());
1365   fitResults->setBinContent(2, 4, bs.dxdzError());
1366   fitResults->setBinContent(2, 3, bs.dydzError());
1367   if (bs.type() > 0) {  // with good beamwidth fit
1368     fitResults->setBinContent(2, 2, bs.BeamWidthXError());
1369     fitResults->setBinContent(2, 1, bs.BeamWidthYError());
1370   } else {  // fill cached width errors
1371     fitResults->setBinContent(2, 2, preBS.BeamWidthXError());
1372     fitResults->setBinContent(2, 1, preBS.BeamWidthYError());
1373   }
1374 
1375   // count good fit
1376   //     if (std::fabs(refBS.x0()-bs.x0())/bs.x0Error() < deltaSigCut_) { // disabled temporarily
1377   summaryContent_[0] += 1.;
1378   //     }
1379   //     if (std::fabs(refBS.y0()-bs.y0())/bs.y0Error() < deltaSigCut_) { // disabled temporarily
1380   summaryContent_[1] += 1.;
1381   //     }
1382   //     if (std::fabs(refBS.z0()-bs.z0())/bs.z0Error() < deltaSigCut_) { // disabled temporarily
1383   summaryContent_[2] += 1.;
1384   //     }
1385 
1386   // Create the BeamSpotOnlineObjects object
1387   BeamSpotOnlineObjects BSOnline;
1388   BSOnline.setLastAnalyzedLumi(LSRange.second);
1389   BSOnline.setLastAnalyzedRun(frun);
1390   BSOnline.setLastAnalyzedFill(0);  // To be updated with correct LHC Fill number
1391   BSOnline.setPosition(bs.x0(), bs.y0(), bs.z0());
1392   BSOnline.setSigmaZ(bs.sigmaZ());
1393   BSOnline.setBeamWidthX(bs.BeamWidthX());
1394   BSOnline.setBeamWidthY(bs.BeamWidthY());
1395   BSOnline.setBeamWidthXError(bs.BeamWidthXError());
1396   BSOnline.setBeamWidthYError(bs.BeamWidthYError());
1397   BSOnline.setdxdz(bs.dxdz());
1398   BSOnline.setdydz(bs.dydz());
1399   BSOnline.setType(bs.type());
1400   BSOnline.setEmittanceX(bs.emittanceX());
1401   BSOnline.setEmittanceY(bs.emittanceY());
1402   BSOnline.setBetaStar(bs.betaStar());
1403   for (int i = 0; i < 7; ++i) {
1404     for (int j = 0; j < 7; ++j) {
1405       BSOnline.setCovariance(i, j, bs.covariance(i, j));
1406     }
1407   }
1408   BSOnline.setNumTracks(50);
1409   BSOnline.setNumPVs(10);
1410   BSOnline.setUsedEvents((int)DipPVInfo_[0]);
1411   BSOnline.setMeanPV(DipPVInfo_[1]);
1412   BSOnline.setMeanErrorPV(DipPVInfo_[2]);
1413   BSOnline.setRmsPV(DipPVInfo_[3]);
1414   BSOnline.setRmsErrorPV(DipPVInfo_[4]);
1415   BSOnline.setMaxPVs((int)DipPVInfo_[5]);
1416   auto creationTime =
1417       std::chrono::duration_cast<std::chrono::microseconds>(std::chrono::system_clock::now().time_since_epoch()).count();
1418   BSOnline.setCreationTime(creationTime);
1419 
1420   // use fake timestamps from 1970.01.01 00:00:00 to 1970.01.01 00:00:01 GMT
1421   std::pair<time_t, time_t> timeForDIP = std::make_pair(0, 1);
1422   BSOnline.setStartTimeStamp(timeForDIP.first);
1423   BSOnline.setStartTime(getGMTstring(timeForDIP.first));
1424   BSOnline.setEndTimeStamp(timeForDIP.second);
1425   BSOnline.setEndTime(getGMTstring(timeForDIP.second));
1426 
1427   edm::LogInfo("FakeBeamMonitor") << "FitAndFill::[PayloadCreation] BeamSpotOnline object created: \n" << std::endl;
1428   edm::LogInfo("FakeBeamMonitor") << BSOnline << std::endl;
1429 
1430   // Create the payload for BeamSpotOnlineObjects object
1431   if (onlineDbService_.isAvailable()) {
1432     edm::LogInfo("FakeBeamMonitor") << "FitAndFill::[PayloadCreation] onlineDbService available \n" << std::endl;
1433     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - Lumi of the current fit: " << currentlumi;
1434     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - Do PV Fitting for LS = " << beginLumiOfPVFit_
1435                                          << " to " << endLumiOfPVFit_;
1436     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - [BeamFitter] Do BeamSpot Fit for LS = "
1437                                          << LSRange.first << " to " << LSRange.second;
1438     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - [FakeBeamMonitor] Do BeamSpot Fit for LS = "
1439                                          << beginLumiOfBSFit_ << " to " << endLumiOfBSFit_;
1440     onlineDbService_->logger().logInfo() << "FakeBeamMonitor - RESULTS OF DEFAULT FIT:";
1441     onlineDbService_->logger().logInfo() << "\n" << bs;
1442     onlineDbService_->logger().logInfo()
1443         << "FakeBeamMonitor::FitAndFill - [PayloadCreation] BeamSpotOnline object created:";
1444     onlineDbService_->logger().logInfo() << "\n" << BSOnline;
1445     onlineDbService_->logger().logInfo() << "FakeBeamMonitor - Additional parameters for DIP:";
1446     onlineDbService_->logger().logInfo() << "Events used in the fit: " << BSOnline.usedEvents();
1447     onlineDbService_->logger().logInfo() << "Mean PV               : " << BSOnline.meanPV();
1448     onlineDbService_->logger().logInfo() << "Mean PV Error         : " << BSOnline.meanErrorPV();
1449     onlineDbService_->logger().logInfo() << "Rms PV                : " << BSOnline.rmsPV();
1450     onlineDbService_->logger().logInfo() << "Rms PV Error          : " << BSOnline.rmsErrorPV();
1451     onlineDbService_->logger().logInfo() << "Max PVs               : " << BSOnline.maxPVs();
1452     onlineDbService_->logger().logInfo() << "StartTime             : " << BSOnline.startTime();
1453     onlineDbService_->logger().logInfo() << "StartTimeStamp        : " << BSOnline.startTimeStamp();
1454     onlineDbService_->logger().logInfo() << "EndTime               : " << BSOnline.endTime();
1455     onlineDbService_->logger().logInfo() << "EndTimeStamp          : " << BSOnline.endTimeStamp();
1456     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - [PayloadCreation] onlineDbService available";
1457     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - [PayloadCreation] SetCreationTime: "
1458                                          << creationTime << " [epoch in microseconds]";
1459     try {
1460       onlineDbService_->writeIOVForNextLumisection(BSOnline, recordName_);
1461       onlineDbService_->logger().logInfo()
1462           << "FakeBeamMonitor::FitAndFill - [PayloadCreation] writeIOVForNextLumisection executed correctly";
1463     } catch (const std::exception& e) {
1464       onlineDbService_->logger().logError() << "FakeBeamMonitor - Error writing record: " << recordName_
1465                                             << " for Run: " << frun << " - Lumi: " << LSRange.second;
1466       onlineDbService_->logger().logError() << "Error is: " << e.what();
1467       onlineDbService_->logger().logError() << "RESULTS OF DEFAULT FIT WAS:";
1468       onlineDbService_->logger().logError() << "\n" << bs;
1469       DBloggerReturn_ = 2;
1470     }
1471   }
1472   edm::LogInfo("FakeBeamMonitor") << "FitAndFill::[PayloadCreation] BeamSpotOnline payload created \n" << std::endl;
1473 
1474   //    }       //if (theBeamFitter->runPVandTrkFitter())
1475   //    else {  // beam fit fails
1476   //      reco::BeamSpot bs = theBeamFitter->getBeamSpot();
1477   //      edm::LogInfo("BeamMonitor") << "FitAndFill::   [BeamMonitor] Beam fit fails!!! \n" << endl;
1478   //      edm::LogInfo("BeamMonitor") << "FitAndFill::   [BeamMonitor] Output beam spot for DIP \n" << endl;
1479   //      edm::LogInfo("BeamMonitor") << bs << endl;
1480   //
1481   //      hs[k_sigmaX0_lumi]->ShiftFillLast(bs.BeamWidthX(), bs.BeamWidthXError(), fitNLumi_);
1482   //      hs[k_sigmaY0_lumi]->ShiftFillLast(bs.BeamWidthY(), bs.BeamWidthYError(), fitNLumi_);
1483   //      hs[k_sigmaZ0_lumi]->ShiftFillLast(bs.sigmaZ(), bs.sigmaZ0Error(), fitNLumi_);
1484   //      hs[k_x0_lumi]->ShiftFillLast(bs.x0(), bs.x0Error(), fitNLumi_);
1485   //      hs[k_y0_lumi]->ShiftFillLast(bs.y0(), bs.y0Error(), fitNLumi_);
1486   //      hs[k_z0_lumi]->ShiftFillLast(bs.z0(), bs.z0Error(), fitNLumi_);
1487   //    }  // end of beam fit fails
1488 
1489   //  }       //-------- end of countFitting------------------------------------------
1490   //  else {  // no fit
1491   //    // Overwrite Fit LS and fit time when no event processed or no track selected
1492   //    theBeamFitter->setFitLSRange(beginLumiOfBSFit_, endLumiOfBSFit_);
1493   //    theBeamFitter->setRefTime(refBStime[0], refBStime[1]);
1494   //    if (theBeamFitter->runPVandTrkFitter()) {
1495   //    }  // Dump fake beam spot for DIP
1496   //    reco::BeamSpot bs = theBeamFitter->getBeamSpot();
1497   //    edm::LogInfo("BeamMonitor") << "FitAndFill::  [BeamMonitor] No fitting \n" << endl;
1498   //    edm::LogInfo("BeamMonitor") << "FitAndFill::  [BeamMonitor] Output fake beam spot for DIP \n" << endl;
1499   //    edm::LogInfo("BeamMonitor") << bs << endl;
1500   //
1501   //    hs[k_sigmaX0_lumi]->ShiftFillLast(bs.BeamWidthX(), bs.BeamWidthXError(), fitNLumi_);
1502   //    hs[k_sigmaY0_lumi]->ShiftFillLast(bs.BeamWidthY(), bs.BeamWidthYError(), fitNLumi_);
1503   //    hs[k_sigmaZ0_lumi]->ShiftFillLast(bs.sigmaZ(), bs.sigmaZ0Error(), fitNLumi_);
1504   //    hs[k_x0_lumi]->ShiftFillLast(bs.x0(), bs.x0Error(), fitNLumi_);
1505   //    hs[k_y0_lumi]->ShiftFillLast(bs.y0(), bs.y0Error(), fitNLumi_);
1506   //    hs[k_z0_lumi]->ShiftFillLast(bs.z0(), bs.z0Error(), fitNLumi_);
1507   //  }
1508 
1509   // Fill summary report
1510   //  if (countFitting) {
1511   for (int n = 0; n < nFitElements_; n++) {
1512     reportSummaryContents[n]->Fill(summaryContent_[n] / (float)nFits_);
1513   }
1514 
1515   summarySum_ = 0;
1516   for (int ii = 0; ii < nFitElements_; ii++) {
1517     summarySum_ += summaryContent_[ii];
1518   }
1519   reportSummary_ = summarySum_ / (nFitElements_ * nFits_);
1520   if (reportSummary)
1521     reportSummary->Fill(reportSummary_);
1522 
1523   for (int bi = 0; bi < nFitElements_; bi++) {
1524     reportSummaryMap->setBinContent(1, bi + 1, summaryContent_[bi] / (float)nFits_);
1525   }
1526   //  }
1527 
1528   if ((resetFitNLumi_ > 0 &&
1529        ((onlineMode_ &&
1530          countLumi_ == resetFitNLumi_) ||  //OR it should be currentLumi_ (if in sequence then does not mattar)
1531         (!onlineMode_ && countLumi_ == resetFitNLumi_))) ||
1532       (StartAverage_)) {
1533     edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: The flag is ON for running average Beam Spot fit" << endl;
1534     StartAverage_ = true;
1535     firstAverageFit_++;
1536     resetHistos_ = true;
1537     nthBSTrk_ = 0;
1538     beginLumiOfBSFit_ = 0;
1539     refBStime[0] = 0;
1540   }
1541 }
1542 
1543 //--------------------------------------------------------
1544 void FakeBeamMonitor::RestartFitting() {  //
1545   if (debug_)
1546     edm::LogInfo("FakeBeamMonitor")
1547         << " RestartingFitting:: Restart Beami everything to a fresh start !!! because Gap is > 10 LS" << endl;
1548   //track based fit reset here
1549   resetHistos_ = true;
1550   nthBSTrk_ = 0;
1551   //  theBeamFitter->resetTrkVector();
1552   //  theBeamFitter->resetLSRange();
1553   //  theBeamFitter->resetRefTime();
1554   //  theBeamFitter->resetPVFitter();
1555   //  theBeamFitter->resetCutFlow();
1556   beginLumiOfBSFit_ = 0;
1557   refBStime[0] = 0;
1558   //pv based fit iis reset here
1559   h_PVx[0]->Reset();
1560   h_PVy[0]->Reset();
1561   h_PVz[0]->Reset();
1562   beginLumiOfPVFit_ = 0;
1563   refPVtime[0] = 0;
1564   //Clear all the Maps here
1565   mapPVx.clear();
1566   mapPVy.clear();
1567   mapPVz.clear();
1568   mapNPV.clear();
1569   mapBeginBSLS.clear();
1570   mapBeginPVLS.clear();
1571   mapBeginBSTime.clear();
1572   mapBeginPVTime.clear();
1573   mapLSBSTrkSize.clear();
1574   mapLSPVStoreSize.clear();
1575   mapLSCF.clear();
1576   countGapLumi_ = 0;
1577   countLumi_ = 0;
1578   StartAverage_ = false;
1579 }
1580 
1581 //-------------------------------------------------------
1582 void FakeBeamMonitor::dqmEndRun(const Run& r, const EventSetup& context) {
1583   if (debug_)
1584     edm::LogInfo("FakeBeamMonitor") << "dqmEndRun:: Clearing all the Maps " << endl;
1585   //Clear all the Maps here
1586   mapPVx.clear();
1587   mapPVy.clear();
1588   mapPVz.clear();
1589   mapNPV.clear();
1590   mapBeginBSLS.clear();
1591   mapBeginPVLS.clear();
1592   mapBeginBSTime.clear();
1593   mapBeginPVTime.clear();
1594   mapLSBSTrkSize.clear();
1595   mapLSPVStoreSize.clear();
1596   mapLSCF.clear();
1597   if (useLockRecords_ && onlineDbService_.isAvailable()) {
1598     onlineDbService_->releaseLocks();
1599   }
1600 }
1601 
1602 //--------------------------------------------------------
1603 void FakeBeamMonitor::scrollTH1(TH1* h, time_t ref) {
1604   char offsetTime[64];
1605   formatFitTime(offsetTime, ref);
1606   TDatime da(offsetTime);
1607   if (lastNZbin > 0) {
1608     double val = h->GetBinContent(lastNZbin);
1609     double valErr = h->GetBinError(lastNZbin);
1610     h->Reset();
1611     h->GetXaxis()->SetTimeOffset(da.Convert(kTRUE));
1612     int bin = (lastNZbin > buffTime ? buffTime : 1);
1613     h->SetBinContent(bin, val);
1614     h->SetBinError(bin, valErr);
1615   } else {
1616     h->Reset();
1617     h->GetXaxis()->SetTimeOffset(da.Convert(kTRUE));
1618   }
1619 }
1620 
1621 //--------------------------------------------------------
1622 // Method to check whether to chane histogram time offset (forward only)
1623 bool FakeBeamMonitor::testScroll(time_t& tmpTime_, time_t& refTime_) {
1624   bool scroll_ = false;
1625   if (tmpTime_ - refTime_ >= intervalInSec_) {
1626     scroll_ = true;
1627     edm::LogInfo("FakeBeamMonitor") << "testScroll::  Reset Time Offset" << std::endl;
1628     lastNZbin = intervalInSec_;
1629     for (int bin = intervalInSec_; bin >= 1; bin--) {
1630       if (hs[k_x0_time]->getBinContent(bin) > 0) {
1631         lastNZbin = bin;
1632         break;
1633       }
1634     }
1635     edm::LogInfo("FakeBeamMonitor") << "testScroll::  Last non zero bin = " << lastNZbin << std::endl;
1636     if (tmpTime_ - refTime_ >= intervalInSec_ + lastNZbin) {
1637       edm::LogInfo("FakeBeamMonitor") << "testScroll::  Time difference too large since last readout" << std::endl;
1638       lastNZbin = 0;
1639       refTime_ = tmpTime_ - buffTime;
1640     } else {
1641       edm::LogInfo("FakeBeamMonitor") << "testScroll::  Offset to last record" << std::endl;
1642       int offset = ((lastNZbin > buffTime) ? (lastNZbin - buffTime) : (lastNZbin - 1));
1643       refTime_ += offset;
1644     }
1645   }
1646   return scroll_;
1647 }
1648 
1649 DEFINE_FWK_MODULE(FakeBeamMonitor);
1650 
1651 // Local Variables:
1652 // show-trailing-whitespace: t
1653 // truncate-lines: t
1654 // End: