Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:06:52 

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_ST = 0;  //For Single Trigger(hence ST)
0730 
0731   //    for (reco::VertexCollection::const_iterator pv = PVCollection->begin(); pv != PVCollection->end(); ++pv) {
0732   for (int tmp_idx = 0; tmp_idx < 10; tmp_idx++) {
0733     //--- vertex selection
0734     //    if (pv->isFake() || pv->tracksSize() == 0)
0735     //      continue;
0736     //if (JetTrigPass)
0737     nPVcount_ST++;  //non-fake pv with a specific trigger
0738 
0739     //    if (pv->ndof() < minVtxNdf_ || (pv->ndof() + 3.) / pv->tracksSize() < 2 * minVtxWgt_)
0740     //      continue;
0741 
0742     //Fill this map to store xyx for pv so that later we can remove the first one for run aver
0743     mapPVx[countLumi_].push_back(tmp_idx);
0744     mapPVy[countLumi_].push_back(tmp_idx);
0745     mapPVz[countLumi_].push_back(tmp_idx);
0746 
0747     //      if (!StartAverage_) {  //for first N LS
0748     //        h_PVx[0]->Fill(pv->x());
0749     //        h_PVy[0]->Fill(pv->y());
0750     //        h_PVz[0]->Fill(pv->z());
0751     //        h_PVxz->Fill(pv->z(), pv->x());
0752     //        h_PVyz->Fill(pv->z(), pv->y());
0753     //      }  //for first N LiS
0754     //      else {
0755     //        h_PVxz->Fill(pv->z(), pv->x());
0756     //        h_PVyz->Fill(pv->z(), pv->y());
0757     //      }
0758 
0759   }  //loop over pvs
0760 
0761   mapNPV[countLumi_].push_back((nPVcount_ST));
0762 
0763   //    if (!StartAverage_) {
0764   //      h_nVtx_st->Fill(nPVcount_ST * 1.);
0765   //    }
0766 
0767   //  }  //if pv collection is availaable
0768 
0769   if (StartAverage_) {
0770     map<int, std::vector<float> >::iterator itpvx = mapPVx.begin();
0771     map<int, std::vector<float> >::iterator itpvy = mapPVy.begin();
0772     map<int, std::vector<float> >::iterator itpvz = mapPVz.begin();
0773 
0774     map<int, std::vector<int> >::iterator itbspvinfo = mapNPV.begin();
0775 
0776     if ((int)mapPVx.size() > resetFitNLumi_) {  //sometimes the events is not there but LS is there!
0777       mapPVx.erase(itpvx);
0778       mapPVy.erase(itpvy);
0779       mapPVz.erase(itpvz);
0780       mapNPV.erase(itbspvinfo);
0781     }  //loop over Last N lumi collected
0782 
0783   }  //StartAverage==true
0784 
0785   processed_ = true;
0786 }
0787 
0788 //--------------------------------------------------------
0789 void FakeBeamMonitor::endLuminosityBlock(const LuminosityBlock& lumiSeg, const EventSetup& iSetup) {
0790   int nthlumi = lumiSeg.id().luminosityBlock();
0791   edm::LogInfo("FakeBeamMonitor") << "endLuminosityBlock:: Lumi of the last event before endLuminosityBlock: "
0792                                   << nthlumi << endl;
0793 
0794   if (onlineMode_ && nthlumi < nextlumi_)
0795     return;
0796   const edm::TimeValue_t fendtimestamp = lumiSeg.endTime().value();
0797   const std::time_t fendtime = fendtimestamp >> 32;
0798   tmpTime = refBStime[1] = refPVtime[1] = fendtime;
0799 
0800   // end DB logger
0801   if (onlineDbService_.isAvailable()) {
0802     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::endLuminosityBlock";
0803     onlineDbService_->logger().end(DBloggerReturn_);
0804   }
0805 }
0806 
0807 //--------------------------------------------------------
0808 void FakeBeamMonitor::FitAndFill(const LuminosityBlock& lumiSeg, int& lastlumi, int& nextlumi, int& nthlumi) {
0809   if (onlineMode_ && (nthlumi <= nextlumi))
0810     return;
0811 
0812   //set the correct run number when no event in the LS for fake output
0813   //  if ((processed_) && theBeamFitter->getRunNumber() != frun)
0814   //    theBeamFitter->setRun(frun);
0815 
0816   int currentlumi = nextlumi;
0817   edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  Lumi of the current fit: " << currentlumi << endl;
0818   lastlumi = currentlumi;
0819   endLumiOfBSFit_ = currentlumi;
0820   endLumiOfPVFit_ = currentlumi;
0821 
0822   //---------Fix for Runninv average-------------
0823   mapLSPVStoreSize[countLumi_] = 10;  //theBeamFitter->getPVvectorSize();
0824 
0825   //  if (StartAverage_) {
0826   //    std::map<int, std::size_t>::iterator rmLSPVi = mapLSPVStoreSize.begin();
0827   //    size_t SizeToRemovePV = rmLSPVi->second;
0828   //    for (std::map<int, std::size_t>::iterator rmLSPVe = mapLSPVStoreSize.end(); ++rmLSPVi != rmLSPVe;)
0829   //      rmLSPVi->second -= SizeToRemovePV;
0830   //
0831   //    theBeamFitter->resizePVvector(SizeToRemovePV);
0832   //
0833   //    map<int, std::size_t>::iterator tmpItpv = mapLSPVStoreSize.begin();
0834   //    mapLSPVStoreSize.erase(tmpItpv);
0835   //  }
0836   //  if (debug_)
0837   //    edm::LogInfo("BeamMonitor") << "FitAndFill:: Size of thePVvector After removing the PVs = "
0838   //                                << theBeamFitter->getPVvectorSize() << endl;
0839 
0840   //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
0841   bool resetHistoFlag_ = false;
0842   if ((int)mapPVx.size() >= resetFitNLumi_ && (StartAverage_)) {
0843     h_PVx[0]->Reset();
0844     h_PVy[0]->Reset();
0845     h_PVz[0]->Reset();
0846     h_nVtx_st->Reset();
0847     resetHistoFlag_ = true;
0848   }
0849 
0850   int MaxPVs = 0;
0851 
0852   std::map<int, std::vector<int> >::iterator mnpv = mapNPV.begin();
0853   std::map<int, std::vector<float> >::iterator mpv2 = mapPVy.begin();
0854   std::map<int, std::vector<float> >::iterator mpv3 = mapPVz.begin();
0855 
0856   for (std::map<int, std::vector<float> >::iterator mpv1 = mapPVx.begin(); mpv1 != mapPVx.end();
0857        ++mpv1, ++mpv2, ++mpv3, ++mnpv) {
0858     std::vector<float>::iterator mpvs2 = (mpv2->second).begin();
0859     std::vector<float>::iterator mpvs3 = (mpv3->second).begin();
0860     for (std::vector<float>::iterator mpvs1 = (mpv1->second).begin(); mpvs1 != (mpv1->second).end();
0861          ++mpvs1, ++mpvs2, ++mpvs3) {
0862       if (resetHistoFlag_) {
0863         h_PVx[0]->Fill(*mpvs1);  //these histogram are reset after StartAverage_ flag is ON
0864         h_PVy[0]->Fill(*mpvs2);
0865         h_PVz[0]->Fill(*mpvs3);
0866       }
0867     }  //loop over second
0868 
0869     //Do the same here for nPV distr.
0870     for (std::vector<int>::iterator mnpvs = (mnpv->second).begin(); mnpvs != (mnpv->second).end(); ++mnpvs) {
0871       if ((*mnpvs > 0) && (resetHistoFlag_))
0872         h_nVtx_st->Fill((*mnpvs) * (1.0));
0873       if ((*mnpvs) > MaxPVs)
0874         MaxPVs = (*mnpvs);
0875     }  //loop over second of mapNPV
0876 
0877   }  //loop over last N lumis
0878 
0879   char tmpTitlePV[100];
0880   sprintf(tmpTitlePV, "%s %i %s %i", "Num. of reco. vertices for LS: ", beginLumiOfPVFit_, " to ", endLumiOfPVFit_);
0881   h_nVtx_st->setAxisTitle(tmpTitlePV, 1);
0882 
0883   std::vector<float> DipPVInfo_;
0884   DipPVInfo_.clear();
0885   DipPVInfo_.push_back(rndm_->Gaus(1000., 100.));  // Events used
0886   DipPVInfo_.push_back(rndm_->Gaus(100., 10.));    // Mean PV
0887   DipPVInfo_.push_back(rndm_->Gaus(10., 5.));      // Mean PV err
0888   DipPVInfo_.push_back(rndm_->Gaus(10., 5.));      // Rms PV
0889   DipPVInfo_.push_back(rndm_->Gaus(5., 3.));       // Rms PV err
0890   DipPVInfo_.push_back(rndm_->Gaus(100., 10.));    // Max PVs
0891 
0892   if (onlineMode_) {  // filling LS gap
0893     // FIXME: need to add protection for the case if the gap is at the resetting LS!
0894     const int countLS_bs = hs[k_x0_lumi]->getTH1()->GetEntries();
0895     const int countLS_pv = hs[k_PVx_lumi]->getTH1()->GetEntries();
0896     edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: countLS_bs = " << countLS_bs << " ; countLS_pv = " << countLS_pv
0897                                     << std::endl;
0898     int LSgap_bs = currentlumi / fitNLumi_ - countLS_bs;
0899     int LSgap_pv = currentlumi / fitPVNLumi_ - countLS_pv;
0900     if (currentlumi % fitNLumi_ == 0)
0901       LSgap_bs--;
0902     if (currentlumi % fitPVNLumi_ == 0)
0903       LSgap_pv--;
0904     edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  LSgap_bs = " << LSgap_bs << " ; LSgap_pv = " << LSgap_pv
0905                                     << std::endl;
0906     // filling previous fits if LS gap ever exists
0907     for (int ig = 0; ig < LSgap_bs; ig++) {
0908       hs[k_x0_lumi]->ShiftFillLast(0., 0., fitNLumi_);  //x0 , x0err, fitNLumi_;  see DQMCore....
0909       hs[k_y0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0910       hs[k_z0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0911       hs[k_sigmaX0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0912       hs[k_sigmaY0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0913       hs[k_sigmaZ0_lumi]->ShiftFillLast(0., 0., fitNLumi_);
0914       h_nVtx_lumi->ShiftFillLast(0., 0., fitNLumi_);
0915     }
0916     for (int ig = 0; ig < LSgap_pv; ig++) {
0917       hs[k_PVx_lumi]->ShiftFillLast(0., 0., fitPVNLumi_);
0918       hs[k_PVy_lumi]->ShiftFillLast(0., 0., fitPVNLumi_);
0919       hs[k_PVz_lumi]->ShiftFillLast(0., 0., fitPVNLumi_);
0920     }
0921     const int previousLS = h_nTrk_lumi->getTH1()->GetEntries();
0922     for (int i = 1; i < (currentlumi - previousLS);
0923          i++)  //if (current-previoius)= 1 then never go inside the for loop!!!!!!!!!!!
0924       h_nTrk_lumi->ShiftFillLast(nthBSTrk_);
0925   }
0926 
0927   edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: Time lapsed since last scroll = " << tmpTime - refTime << std::endl;
0928 
0929   if (testScroll(tmpTime, refTime)) {
0930     scrollTH1(hs[k_x0_time]->getTH1(), refTime);
0931     scrollTH1(hs[k_y0_time]->getTH1(), refTime);
0932     scrollTH1(hs[k_z0_time]->getTH1(), refTime);
0933     scrollTH1(hs[k_sigmaX0_time]->getTH1(), refTime);
0934     scrollTH1(hs[k_sigmaY0_time]->getTH1(), refTime);
0935     scrollTH1(hs[k_sigmaZ0_time]->getTH1(), refTime);
0936     scrollTH1(hs[k_PVx_time]->getTH1(), refTime);
0937     scrollTH1(hs[k_PVy_time]->getTH1(), refTime);
0938     scrollTH1(hs[k_PVz_time]->getTH1(), refTime);
0939   }
0940 
0941   //  bool doPVFit = false;
0942   //
0943   //  if (fitPVNLumi_ > 0) {
0944   //    if (onlineMode_) {
0945   //      if (currentlumi % fitPVNLumi_ == 0)
0946   //        doPVFit = true;
0947   //    } else if (countLumi_ % fitPVNLumi_ == 0)
0948   //      doPVFit = true;
0949   //  } else
0950   //    doPVFit = true;
0951   //
0952   //  if (doPVFit) {
0953   edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: Do PV Fitting for LS = " << beginLumiOfPVFit_ << " to "
0954                                   << endLumiOfPVFit_ << std::endl;
0955   // Primary Vertex Fit:
0956   if (h_PVx[0]->getTH1()->GetEntries() > minNrVertices_) {
0957     pvResults->Reset();
0958     char tmpTitle[50];
0959     sprintf(tmpTitle, "%s %i %s %i", "Fitted Primary Vertex (cm) of LS: ", beginLumiOfPVFit_, " to ", endLumiOfPVFit_);
0960     pvResults->setAxisTitle(tmpTitle, 1);
0961 
0962     std::unique_ptr<TF1> fgaus{new TF1("fgaus", "gaus")};
0963     double mean, width, meanErr, widthErr;
0964     fgaus->SetLineColor(4);
0965     h_PVx[0]->getTH1()->Fit(fgaus.get(), "QLM0");
0966     mean = fgaus->GetParameter(1);
0967     width = fgaus->GetParameter(2);
0968     meanErr = fgaus->GetParError(1);
0969     widthErr = fgaus->GetParError(2);
0970 
0971     hs[k_PVx_lumi]->ShiftFillLast(mean, width, fitPVNLumi_);
0972     hs[k_PVx_lumi_all]->setBinContent(currentlumi, mean);
0973     hs[k_PVx_lumi_all]->setBinError(currentlumi, width);
0974     int nthBin = tmpTime - refTime;
0975     if (nthBin < 0)
0976       edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  Event time outside current range of time histograms!"
0977                                       << std::endl;
0978     if (nthBin > 0) {
0979       hs[k_PVx_time]->setBinContent(nthBin, mean);
0980       hs[k_PVx_time]->setBinError(nthBin, width);
0981     }
0982     int jthBin = tmpTime - startTime;
0983     if (jthBin > 0) {
0984       hs[k_PVx_time_all]->setBinContent(jthBin, mean);
0985       hs[k_PVx_time_all]->setBinError(jthBin, width);
0986     }
0987     pvResults->setBinContent(1, 6, mean);
0988     pvResults->setBinContent(1, 3, width);
0989     pvResults->setBinContent(2, 6, meanErr);
0990     pvResults->setBinContent(2, 3, widthErr);
0991 
0992     {
0993       // snap shot of the fit
0994       auto tmphisto = h_PVx[0]->getTH1F();
0995       h_PVx[1]->getTH1()->SetBins(
0996           tmphisto->GetNbinsX(), tmphisto->GetXaxis()->GetXmin(), tmphisto->GetXaxis()->GetXmax());
0997       h_PVx[1]->Reset();
0998       h_PVx[1]->getTH1()->Add(tmphisto);
0999       h_PVx[1]->getTH1()->Fit(fgaus.get(), "QLM");
1000     }
1001 
1002     h_PVy[0]->getTH1()->Fit(fgaus.get(), "QLM0");
1003     mean = fgaus->GetParameter(1);
1004     width = fgaus->GetParameter(2);
1005     meanErr = fgaus->GetParError(1);
1006     widthErr = fgaus->GetParError(2);
1007     hs[k_PVy_lumi]->ShiftFillLast(mean, width, fitPVNLumi_);
1008     hs[k_PVy_lumi_all]->setBinContent(currentlumi, mean);
1009     hs[k_PVy_lumi_all]->setBinError(currentlumi, width);
1010     if (nthBin > 0) {
1011       hs[k_PVy_time]->setBinContent(nthBin, mean);
1012       hs[k_PVy_time]->setBinError(nthBin, width);
1013     }
1014     if (jthBin > 0) {
1015       hs[k_PVy_time_all]->setBinContent(jthBin, mean);
1016       hs[k_PVy_time_all]->setBinError(jthBin, width);
1017     }
1018     pvResults->setBinContent(1, 5, mean);
1019     pvResults->setBinContent(1, 2, width);
1020     pvResults->setBinContent(2, 5, meanErr);
1021     pvResults->setBinContent(2, 2, widthErr);
1022     // snap shot of the fit
1023     {
1024       auto tmphisto = h_PVy[0]->getTH1F();
1025       h_PVy[1]->getTH1()->SetBins(
1026           tmphisto->GetNbinsX(), tmphisto->GetXaxis()->GetXmin(), tmphisto->GetXaxis()->GetXmax());
1027       h_PVy[1]->Reset();
1028       h_PVy[1]->getTH1()->Add(tmphisto);
1029       h_PVy[1]->getTH1()->Fit(fgaus.get(), "QLM");
1030     }
1031 
1032     h_PVz[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_PVz_lumi]->ShiftFillLast(mean, width, fitPVNLumi_);
1038     hs[k_PVz_lumi_all]->setBinContent(currentlumi, mean);
1039     hs[k_PVz_lumi_all]->setBinError(currentlumi, width);
1040     if (nthBin > 0) {
1041       hs[k_PVz_time]->setBinContent(nthBin, mean);
1042       hs[k_PVz_time]->setBinError(nthBin, width);
1043     }
1044     if (jthBin > 0) {
1045       hs[k_PVz_time_all]->setBinContent(jthBin, mean);
1046       hs[k_PVz_time_all]->setBinError(jthBin, width);
1047     }
1048     pvResults->setBinContent(1, 4, mean);
1049     pvResults->setBinContent(1, 1, width);
1050     pvResults->setBinContent(2, 4, meanErr);
1051     pvResults->setBinContent(2, 1, widthErr);
1052     {
1053       // snap shot of the fit
1054       auto tmphisto = h_PVz[0]->getTH1F();
1055       h_PVz[1]->getTH1()->SetBins(
1056           tmphisto->GetNbinsX(), tmphisto->GetXaxis()->GetXmin(), tmphisto->GetXaxis()->GetXmax());
1057       h_PVz[1]->Reset();
1058       h_PVz[1]->getTH1()->Add(tmphisto);
1059       h_PVz[1]->getTH1()->Fit(fgaus.get(), "QLM");
1060     }
1061   }  //check if found min Vertices
1062      //  }    //do PVfit
1063 
1064   if ((resetPVNLumi_ > 0 && countLumi_ == resetPVNLumi_) || StartAverage_) {
1065     beginLumiOfPVFit_ = 0;
1066     refPVtime[0] = 0;
1067   }
1068 
1069   //---------Readjustment of theBSvector, RefTime, beginLSofFit---------
1070   //  vector<BSTrkParameters> theBSvector1 = theBeamFitter->getBSvector();
1071   //  mapLSBSTrkSize[countLumi_] = (theBSvector1.size());
1072   size_t PreviousRecords = 0;  //needed to fill nth record of tracks in GUI
1073 
1074   //  if (StartAverage_) {
1075   //    size_t SizeToRemove = 0;
1076   //    std::map<int, std::size_t>::iterator rmls = mapLSBSTrkSize.begin();
1077   //    SizeToRemove = rmls->second;
1078   //    if (debug_)
1079   //      edm::LogInfo("BeamMonitor") << "  The size to remove is =  " << SizeToRemove << endl;
1080   //    int changedAfterThis = 0;
1081   //    for (std::map<int, std::size_t>::iterator rmLS = mapLSBSTrkSize.begin(); rmLS != mapLSBSTrkSize.end();
1082   //         ++rmLS, ++changedAfterThis) {
1083   //      if (changedAfterThis > 0) {
1084   //        (rmLS->second) = (rmLS->second) - SizeToRemove;
1085   //        if ((mapLSBSTrkSize.size() - (size_t)changedAfterThis) == 2)
1086   //          PreviousRecords = (rmLS->second);
1087   //      }
1088   //    }
1089   //
1090   //    theBeamFitter->resizeBSvector(SizeToRemove);
1091   //
1092   //    map<int, std::size_t>::iterator tmpIt = mapLSBSTrkSize.begin();
1093   //    mapLSBSTrkSize.erase(tmpIt);
1094   //
1095   //    std::pair<int, int> checkfitLS = theBeamFitter->getFitLSRange();
1096   //    std::pair<time_t, time_t> checkfitTime = theBeamFitter->getRefTime();
1097   //    theBeamFitter->setFitLSRange(beginLumiOfBSFit_, checkfitLS.second);
1098   //    theBeamFitter->setRefTime(refBStime[0], checkfitTime.second);
1099   //  }
1100 
1101   //Fill the track for this fit
1102   //  vector<BSTrkParameters> theBSvector = theBeamFitter->getBSvector();
1103   //  h_nTrk_lumi->ShiftFillLast(theBSvector.size());
1104   //
1105   //  if (debug_)
1106   //    edm::LogInfo("BeamMonitor") << "FitAndFill::   Size of  theBSViector.size()  After =" << theBSvector.size() << endl;
1107 
1108   //  bool countFitting = false;
1109   //  if (theBSvector.size() >= PreviousRecords && theBSvector.size() >= min_Ntrks_) {
1110   //    countFitting = true;
1111   //  }
1112 
1113   //---Fix for Cut Flow Table for Running average in a same way//the previous code  has problem for resetting!!!
1114   //  mapLSCF[countLumi_] = *theBeamFitter->getCutFlow();
1115   //  if (StartAverage_ && !mapLSCF.empty()) {
1116   //    const TH1F& cutFlowToSubtract = mapLSCF.begin()->second;
1117   //    // Subtract the last cut flow from all of the others.
1118   //    std::map<int, TH1F>::iterator cf = mapLSCF.begin();
1119   //    // Start on second entry
1120   //    for (; cf != mapLSCF.end(); ++cf) {
1121   //      cf->second.Add(&cutFlowToSubtract, -1);
1122   //    }
1123   //    theBeamFitter->subtractFromCutFlow(&cutFlowToSubtract);
1124   //    // Remove the obsolete lumi section
1125   //    mapLSCF.erase(mapLSCF.begin());
1126   //  }
1127 
1128   if (resetHistos_) {
1129     h_d0_phi0->Reset();
1130     h_vx_vy->Reset();
1131     h_vx_dz->Reset();
1132     h_vy_dz->Reset();
1133     h_trk_z0->Reset();
1134     resetHistos_ = false;
1135   }
1136 
1137   if (StartAverage_)
1138     nthBSTrk_ = PreviousRecords;  //after average proccess is ON//for 2-6 LS fit PreviousRecords is size from 2-5 LS
1139 
1140   edm::LogInfo("FakeBeamMonitor") << " The Previous Recored for this fit is  =" << nthBSTrk_ << endl;
1141 
1142   //  unsigned int itrk = 0;
1143   //  for (vector<BSTrkParameters>::const_iterator BSTrk = theBSvector.begin(); BSTrk != theBSvector.end();
1144   //       ++BSTrk, ++itrk) {
1145   //    if (itrk >= nthBSTrk_) {  //fill for this record only !!
1146   //      h_d0_phi0->Fill(BSTrk->phi0(), BSTrk->d0());
1147   //      double vx = BSTrk->vx();
1148   //      double vy = BSTrk->vy();
1149   //      double z0 = BSTrk->z0();
1150   //      h_vx_vy->Fill(vx, vy);
1151   //      h_vx_dz->Fill(z0, vx);
1152   //      h_vy_dz->Fill(z0, vy);
1153   //      h_trk_z0->Fill(z0);
1154   //    }
1155   //  }
1156 
1157   //  nthBSTrk_ = theBSvector.size();  // keep track of num of tracks filled so far
1158 
1159   edm::LogInfo("FakeBeamMonitor") << " The Current Recored for this fit is  =" << nthBSTrk_ << endl;
1160 
1161   //  if (countFitting)
1162   //    edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  Num of tracks collected = " << nthBSTrk_ << endl;
1163 
1164   if (fitNLumi_ > 0) {
1165     if (onlineMode_) {
1166       if (currentlumi % fitNLumi_ != 0) {
1167         //  for (std::map<TString,MonitorElement*>::iterator itAll = hs.begin();
1168         //       itAll != hs.end(); ++itAll) {
1169         //    if ((*itAll).first.Contains("all")) {
1170         //      (*itAll).second->setBinContent(currentlumi,0.);
1171         //      (*itAll).second->setBinError(currentlumi,0.);
1172         //    }
1173         //  }
1174         return;
1175       }
1176     } else if (countLumi_ % fitNLumi_ != 0)
1177       return;
1178   }
1179 
1180   edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: [DebugTime] refBStime[0] = " << refBStime[0]
1181                                   << "; address =  " << &refBStime[0] << std::endl;
1182   edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: [DebugTime] refBStime[1] = " << refBStime[1]
1183                                   << "; address =  " << &refBStime[1] << std::endl;
1184 
1185   //Fill for all LS even if fit fails
1186   //  h_nVtx_lumi->ShiftFillLast((theBeamFitter->getPVvectorSize()), 0., fitNLumi_);
1187   //  h_nVtx_lumi_all->setBinContent(currentlumi, (theBeamFitter->getPVvectorSize()));
1188 
1189   //  if (countFitting) {
1190   nFits_++;
1191   //    std::pair<int, int> fitLS = theBeamFitter->getFitLSRange();
1192   std::pair<int, int> fitLS(beginLumiOfBSFit_, endLumiOfBSFit_);
1193   //    edm::LogInfo("BeamMonitor") << "FitAndFill::  [BeamFitter] Do BeamSpot Fit for LS = " << fitLS.first << " to "
1194   //                                << fitLS.second << std::endl;
1195   edm::LogInfo("FakeBeamMonitor") << "FitAndFill::  [FakeBeamMonitor] Do BeamSpot Fit for LS = " << beginLumiOfBSFit_
1196                                   << " to " << endLumiOfBSFit_ << std::endl;
1197 
1198   //Now Run the PV and Track Fitter over the collected tracks and pvs
1199   //    if (theBeamFitter->runPVandTrkFitter()) {
1200   //      reco::BeamSpot bs = theBeamFitter->getBeamSpot();
1201 
1202   //Create fake BS here
1203   reco::BeamSpot::CovarianceMatrix matrix;
1204   for (int j = 0; j < 7; ++j) {
1205     for (int k = j; k < 7; ++k) {
1206       matrix(j, k) = 0;
1207     }
1208   }
1209 
1210   // random values for fake BeamSpot
1211   float tmp_BSx = rndm_->Gaus(0.1, 0.1);            // [cm]
1212   float tmp_BSy = rndm_->Gaus(0.1, 0.1);            // [cm]
1213   float tmp_BSz = rndm_->Gaus(0.1, 0.1);            // [cm]
1214   float tmp_BSwidthX = rndm_->Gaus(0.001, 0.0005);  // [cm]
1215   float tmp_BSwidthY = rndm_->Gaus(0.001, 0.0005);  // [cm]
1216   float tmp_BSwidthZ = rndm_->Gaus(3.5, 0.5);       // [cm]
1217 
1218   reco::BeamSpot bs(reco::BeamSpot::Point(tmp_BSx, tmp_BSy, tmp_BSz),
1219                     tmp_BSwidthZ,
1220                     0,
1221                     0,
1222                     tmp_BSwidthX,
1223                     matrix,
1224                     reco::BeamSpot::Tracker);
1225   bs.setBeamWidthY(tmp_BSwidthY);
1226 
1227   if (bs.type() > 0)  // with good beamwidth fit
1228     preBS = bs;       // cache good fit results
1229 
1230   edm::LogInfo("FakeBeamMonitor") << "\n RESULTS OF DEFAULT FIT:" << endl;
1231   edm::LogInfo("FakeBeamMonitor") << bs << endl;
1232   edm::LogInfo("FakeBeamMonitor") << "[BeamFitter] fitting done \n" << endl;
1233 
1234   hs[k_x0_lumi]->ShiftFillLast(bs.x0(), bs.x0Error(), fitNLumi_);
1235   hs[k_y0_lumi]->ShiftFillLast(bs.y0(), bs.y0Error(), fitNLumi_);
1236   hs[k_z0_lumi]->ShiftFillLast(bs.z0(), bs.z0Error(), fitNLumi_);
1237   hs[k_sigmaX0_lumi]->ShiftFillLast(bs.BeamWidthX(), bs.BeamWidthXError(), fitNLumi_);
1238   hs[k_sigmaY0_lumi]->ShiftFillLast(bs.BeamWidthY(), bs.BeamWidthYError(), fitNLumi_);
1239   hs[k_sigmaZ0_lumi]->ShiftFillLast(bs.sigmaZ(), bs.sigmaZ0Error(), fitNLumi_);
1240   hs[k_x0_lumi_all]->setBinContent(currentlumi, bs.x0());
1241   hs[k_x0_lumi_all]->setBinError(currentlumi, bs.x0Error());
1242   hs[k_y0_lumi_all]->setBinContent(currentlumi, bs.y0());
1243   hs[k_y0_lumi_all]->setBinError(currentlumi, bs.y0Error());
1244   hs[k_z0_lumi_all]->setBinContent(currentlumi, bs.z0());
1245   hs[k_z0_lumi_all]->setBinError(currentlumi, bs.z0Error());
1246   hs[k_sigmaX0_lumi_all]->setBinContent(currentlumi, bs.BeamWidthX());
1247   hs[k_sigmaX0_lumi_all]->setBinError(currentlumi, bs.BeamWidthXError());
1248   hs[k_sigmaY0_lumi_all]->setBinContent(currentlumi, bs.BeamWidthY());
1249   hs[k_sigmaY0_lumi_all]->setBinError(currentlumi, bs.BeamWidthYError());
1250   hs[k_sigmaZ0_lumi_all]->setBinContent(currentlumi, bs.sigmaZ());
1251   hs[k_sigmaZ0_lumi_all]->setBinError(currentlumi, bs.sigmaZ0Error());
1252 
1253   int nthBin = tmpTime - refTime;
1254   if (nthBin > 0) {
1255     hs[k_x0_time]->setBinContent(nthBin, bs.x0());
1256     hs[k_y0_time]->setBinContent(nthBin, bs.y0());
1257     hs[k_z0_time]->setBinContent(nthBin, bs.z0());
1258     hs[k_sigmaX0_time]->setBinContent(nthBin, bs.BeamWidthX());
1259     hs[k_sigmaY0_time]->setBinContent(nthBin, bs.BeamWidthY());
1260     hs[k_sigmaZ0_time]->setBinContent(nthBin, bs.sigmaZ());
1261     hs[k_x0_time]->setBinError(nthBin, bs.x0Error());
1262     hs[k_y0_time]->setBinError(nthBin, bs.y0Error());
1263     hs[k_z0_time]->setBinError(nthBin, bs.z0Error());
1264     hs[k_sigmaX0_time]->setBinError(nthBin, bs.BeamWidthXError());
1265     hs[k_sigmaY0_time]->setBinError(nthBin, bs.BeamWidthYError());
1266     hs[k_sigmaZ0_time]->setBinError(nthBin, bs.sigmaZ0Error());
1267   }
1268 
1269   int jthBin = tmpTime - startTime;
1270   if (jthBin > 0) {
1271     hs[k_x0_time_all]->setBinContent(jthBin, bs.x0());
1272     hs[k_y0_time_all]->setBinContent(jthBin, bs.y0());
1273     hs[k_z0_time_all]->setBinContent(jthBin, bs.z0());
1274     hs[k_sigmaX0_time_all]->setBinContent(jthBin, bs.BeamWidthX());
1275     hs[k_sigmaY0_time_all]->setBinContent(jthBin, bs.BeamWidthY());
1276     hs[k_sigmaZ0_time_all]->setBinContent(jthBin, bs.sigmaZ());
1277     hs[k_x0_time_all]->setBinError(jthBin, bs.x0Error());
1278     hs[k_y0_time_all]->setBinError(jthBin, bs.y0Error());
1279     hs[k_z0_time_all]->setBinError(jthBin, bs.z0Error());
1280     hs[k_sigmaX0_time_all]->setBinError(jthBin, bs.BeamWidthXError());
1281     hs[k_sigmaY0_time_all]->setBinError(jthBin, bs.BeamWidthYError());
1282     hs[k_sigmaZ0_time_all]->setBinError(jthBin, bs.sigmaZ0Error());
1283   }
1284 
1285   h_x0->Fill(bs.x0());
1286   h_y0->Fill(bs.y0());
1287   h_z0->Fill(bs.z0());
1288   if (bs.type() > 0) {  // with good beamwidth fit
1289     h_sigmaX0->Fill(bs.BeamWidthX());
1290     h_sigmaY0->Fill(bs.BeamWidthY());
1291   }
1292   h_sigmaZ0->Fill(bs.sigmaZ());
1293 
1294   if (nthBSTrk_ >= 2 * min_Ntrks_) {
1295     double amp = std::sqrt(bs.x0() * bs.x0() + bs.y0() * bs.y0());
1296     double alpha = std::atan2(bs.y0(), bs.x0());
1297     std::unique_ptr<TF1> f1{new TF1("f1", "[0]*sin(x-[1])", -3.14, 3.14)};
1298     f1->SetParameters(amp, alpha);
1299     f1->SetParLimits(0, amp - 0.1, amp + 0.1);
1300     f1->SetParLimits(1, alpha - 0.577, alpha + 0.577);
1301     f1->SetLineColor(4);
1302     h_d0_phi0->getTProfile()->Fit(f1.get(), "QR");
1303 
1304     double mean = bs.z0();
1305     double width = bs.sigmaZ();
1306     std::unique_ptr<TF1> fgaus{new TF1("fgaus", "gaus")};
1307     fgaus->SetParameters(mean, width);
1308     fgaus->SetLineColor(4);
1309     h_trk_z0->getTH1()->Fit(fgaus.get(), "QLRM", "", mean - 3 * width, mean + 3 * width);
1310   }
1311 
1312   fitResults->Reset();
1313   std::pair<int, int> LSRange(beginLumiOfBSFit_, endLumiOfBSFit_);  //= theBeamFitter->getFitLSRange();
1314   char tmpTitle[50];
1315   sprintf(tmpTitle, "%s %i %s %i", "Fitted Beam Spot (cm) of LS: ", LSRange.first, " to ", LSRange.second);
1316   fitResults->setAxisTitle(tmpTitle, 1);
1317   fitResults->setBinContent(1, 8, bs.x0());
1318   fitResults->setBinContent(1, 7, bs.y0());
1319   fitResults->setBinContent(1, 6, bs.z0());
1320   fitResults->setBinContent(1, 5, bs.sigmaZ());
1321   fitResults->setBinContent(1, 4, bs.dxdz());
1322   fitResults->setBinContent(1, 3, bs.dydz());
1323   if (bs.type() > 0) {  // with good beamwidth fit
1324     fitResults->setBinContent(1, 2, bs.BeamWidthX());
1325     fitResults->setBinContent(1, 1, bs.BeamWidthY());
1326   } else {  // fill cached widths
1327     fitResults->setBinContent(1, 2, preBS.BeamWidthX());
1328     fitResults->setBinContent(1, 1, preBS.BeamWidthY());
1329   }
1330 
1331   fitResults->setBinContent(2, 8, bs.x0Error());
1332   fitResults->setBinContent(2, 7, bs.y0Error());
1333   fitResults->setBinContent(2, 6, bs.z0Error());
1334   fitResults->setBinContent(2, 5, bs.sigmaZ0Error());
1335   fitResults->setBinContent(2, 4, bs.dxdzError());
1336   fitResults->setBinContent(2, 3, bs.dydzError());
1337   if (bs.type() > 0) {  // with good beamwidth fit
1338     fitResults->setBinContent(2, 2, bs.BeamWidthXError());
1339     fitResults->setBinContent(2, 1, bs.BeamWidthYError());
1340   } else {  // fill cached width errors
1341     fitResults->setBinContent(2, 2, preBS.BeamWidthXError());
1342     fitResults->setBinContent(2, 1, preBS.BeamWidthYError());
1343   }
1344 
1345   // count good fit
1346   //     if (std::fabs(refBS.x0()-bs.x0())/bs.x0Error() < deltaSigCut_) { // disabled temporarily
1347   summaryContent_[0] += 1.;
1348   //     }
1349   //     if (std::fabs(refBS.y0()-bs.y0())/bs.y0Error() < deltaSigCut_) { // disabled temporarily
1350   summaryContent_[1] += 1.;
1351   //     }
1352   //     if (std::fabs(refBS.z0()-bs.z0())/bs.z0Error() < deltaSigCut_) { // disabled temporarily
1353   summaryContent_[2] += 1.;
1354   //     }
1355 
1356   // Create the BeamSpotOnlineObjects object
1357   BeamSpotOnlineObjects BSOnline;
1358   BSOnline.setLastAnalyzedLumi(LSRange.second);
1359   BSOnline.setLastAnalyzedRun(frun);
1360   BSOnline.setLastAnalyzedFill(0);  // To be updated with correct LHC Fill number
1361   BSOnline.setPosition(bs.x0(), bs.y0(), bs.z0());
1362   BSOnline.setSigmaZ(bs.sigmaZ());
1363   BSOnline.setBeamWidthX(bs.BeamWidthX());
1364   BSOnline.setBeamWidthY(bs.BeamWidthY());
1365   BSOnline.setBeamWidthXError(bs.BeamWidthXError());
1366   BSOnline.setBeamWidthYError(bs.BeamWidthYError());
1367   BSOnline.setdxdz(bs.dxdz());
1368   BSOnline.setdydz(bs.dydz());
1369   BSOnline.setType(bs.type());
1370   BSOnline.setEmittanceX(bs.emittanceX());
1371   BSOnline.setEmittanceY(bs.emittanceY());
1372   BSOnline.setBetaStar(bs.betaStar());
1373   for (int i = 0; i < 7; ++i) {
1374     for (int j = 0; j < 7; ++j) {
1375       BSOnline.setCovariance(i, j, bs.covariance(i, j));
1376     }
1377   }
1378   BSOnline.setNumTracks(50);
1379   BSOnline.setNumPVs(10);
1380   BSOnline.setUsedEvents((int)DipPVInfo_[0]);
1381   BSOnline.setMeanPV(DipPVInfo_[1]);
1382   BSOnline.setMeanErrorPV(DipPVInfo_[2]);
1383   BSOnline.setRmsPV(DipPVInfo_[3]);
1384   BSOnline.setRmsErrorPV(DipPVInfo_[4]);
1385   BSOnline.setMaxPVs((int)DipPVInfo_[5]);
1386   auto creationTime =
1387       std::chrono::duration_cast<std::chrono::microseconds>(std::chrono::system_clock::now().time_since_epoch()).count();
1388   BSOnline.setCreationTime(creationTime);
1389 
1390   // use fake timestamps from 1970.01.01 00:00:00 to 1970.01.01 00:00:01 GMT
1391   std::pair<time_t, time_t> timeForDIP = std::make_pair(0, 1);
1392   BSOnline.setStartTimeStamp(timeForDIP.first);
1393   BSOnline.setStartTime(getGMTstring(timeForDIP.first));
1394   BSOnline.setEndTimeStamp(timeForDIP.second);
1395   BSOnline.setEndTime(getGMTstring(timeForDIP.second));
1396 
1397   edm::LogInfo("FakeBeamMonitor") << "FitAndFill::[PayloadCreation] BeamSpotOnline object created: \n" << std::endl;
1398   edm::LogInfo("FakeBeamMonitor") << BSOnline << std::endl;
1399 
1400   // Create the payload for BeamSpotOnlineObjects object
1401   if (onlineDbService_.isAvailable()) {
1402     edm::LogInfo("FakeBeamMonitor") << "FitAndFill::[PayloadCreation] onlineDbService available \n" << std::endl;
1403     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - Lumi of the current fit: " << currentlumi;
1404     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - Do PV Fitting for LS = " << beginLumiOfPVFit_
1405                                          << " to " << endLumiOfPVFit_;
1406     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - [BeamFitter] Do BeamSpot Fit for LS = "
1407                                          << LSRange.first << " to " << LSRange.second;
1408     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - [FakeBeamMonitor] Do BeamSpot Fit for LS = "
1409                                          << beginLumiOfBSFit_ << " to " << endLumiOfBSFit_;
1410     onlineDbService_->logger().logInfo() << "FakeBeamMonitor - RESULTS OF DEFAULT FIT:";
1411     onlineDbService_->logger().logInfo() << "\n" << bs;
1412     onlineDbService_->logger().logInfo()
1413         << "FakeBeamMonitor::FitAndFill - [PayloadCreation] BeamSpotOnline object created:";
1414     onlineDbService_->logger().logInfo() << "\n" << BSOnline;
1415     onlineDbService_->logger().logInfo() << "FakeBeamMonitor - Additional parameters for DIP:";
1416     onlineDbService_->logger().logInfo() << "Events used in the fit: " << BSOnline.usedEvents();
1417     onlineDbService_->logger().logInfo() << "Mean PV               : " << BSOnline.meanPV();
1418     onlineDbService_->logger().logInfo() << "Mean PV Error         : " << BSOnline.meanErrorPV();
1419     onlineDbService_->logger().logInfo() << "Rms PV                : " << BSOnline.rmsPV();
1420     onlineDbService_->logger().logInfo() << "Rms PV Error          : " << BSOnline.rmsErrorPV();
1421     onlineDbService_->logger().logInfo() << "Max PVs               : " << BSOnline.maxPVs();
1422     onlineDbService_->logger().logInfo() << "StartTime             : " << BSOnline.startTime();
1423     onlineDbService_->logger().logInfo() << "StartTimeStamp        : " << BSOnline.startTimeStamp();
1424     onlineDbService_->logger().logInfo() << "EndTime               : " << BSOnline.endTime();
1425     onlineDbService_->logger().logInfo() << "EndTimeStamp          : " << BSOnline.endTimeStamp();
1426     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - [PayloadCreation] onlineDbService available";
1427     onlineDbService_->logger().logInfo() << "FakeBeamMonitor::FitAndFill - [PayloadCreation] SetCreationTime: "
1428                                          << creationTime << " [epoch in microseconds]";
1429     try {
1430       onlineDbService_->writeIOVForNextLumisection(BSOnline, recordName_);
1431       onlineDbService_->logger().logInfo()
1432           << "FakeBeamMonitor::FitAndFill - [PayloadCreation] writeIOVForNextLumisection executed correctly";
1433     } catch (const std::exception& e) {
1434       onlineDbService_->logger().logError() << "FakeBeamMonitor - Error writing record: " << recordName_
1435                                             << " for Run: " << frun << " - Lumi: " << LSRange.second;
1436       onlineDbService_->logger().logError() << "Error is: " << e.what();
1437       onlineDbService_->logger().logError() << "RESULTS OF DEFAULT FIT WAS:";
1438       onlineDbService_->logger().logError() << "\n" << bs;
1439       DBloggerReturn_ = 2;
1440     }
1441   }
1442   edm::LogInfo("FakeBeamMonitor") << "FitAndFill::[PayloadCreation] BeamSpotOnline payload created \n" << std::endl;
1443 
1444   //    }       //if (theBeamFitter->runPVandTrkFitter())
1445   //    else {  // beam fit fails
1446   //      reco::BeamSpot bs = theBeamFitter->getBeamSpot();
1447   //      edm::LogInfo("BeamMonitor") << "FitAndFill::   [BeamMonitor] Beam fit fails!!! \n" << endl;
1448   //      edm::LogInfo("BeamMonitor") << "FitAndFill::   [BeamMonitor] Output beam spot for DIP \n" << endl;
1449   //      edm::LogInfo("BeamMonitor") << bs << endl;
1450   //
1451   //      hs[k_sigmaX0_lumi]->ShiftFillLast(bs.BeamWidthX(), bs.BeamWidthXError(), fitNLumi_);
1452   //      hs[k_sigmaY0_lumi]->ShiftFillLast(bs.BeamWidthY(), bs.BeamWidthYError(), fitNLumi_);
1453   //      hs[k_sigmaZ0_lumi]->ShiftFillLast(bs.sigmaZ(), bs.sigmaZ0Error(), fitNLumi_);
1454   //      hs[k_x0_lumi]->ShiftFillLast(bs.x0(), bs.x0Error(), fitNLumi_);
1455   //      hs[k_y0_lumi]->ShiftFillLast(bs.y0(), bs.y0Error(), fitNLumi_);
1456   //      hs[k_z0_lumi]->ShiftFillLast(bs.z0(), bs.z0Error(), fitNLumi_);
1457   //    }  // end of beam fit fails
1458 
1459   //  }       //-------- end of countFitting------------------------------------------
1460   //  else {  // no fit
1461   //    // Overwrite Fit LS and fit time when no event processed or no track selected
1462   //    theBeamFitter->setFitLSRange(beginLumiOfBSFit_, endLumiOfBSFit_);
1463   //    theBeamFitter->setRefTime(refBStime[0], refBStime[1]);
1464   //    if (theBeamFitter->runPVandTrkFitter()) {
1465   //    }  // Dump fake beam spot for DIP
1466   //    reco::BeamSpot bs = theBeamFitter->getBeamSpot();
1467   //    edm::LogInfo("BeamMonitor") << "FitAndFill::  [BeamMonitor] No fitting \n" << endl;
1468   //    edm::LogInfo("BeamMonitor") << "FitAndFill::  [BeamMonitor] Output fake beam spot for DIP \n" << endl;
1469   //    edm::LogInfo("BeamMonitor") << bs << endl;
1470   //
1471   //    hs[k_sigmaX0_lumi]->ShiftFillLast(bs.BeamWidthX(), bs.BeamWidthXError(), fitNLumi_);
1472   //    hs[k_sigmaY0_lumi]->ShiftFillLast(bs.BeamWidthY(), bs.BeamWidthYError(), fitNLumi_);
1473   //    hs[k_sigmaZ0_lumi]->ShiftFillLast(bs.sigmaZ(), bs.sigmaZ0Error(), fitNLumi_);
1474   //    hs[k_x0_lumi]->ShiftFillLast(bs.x0(), bs.x0Error(), fitNLumi_);
1475   //    hs[k_y0_lumi]->ShiftFillLast(bs.y0(), bs.y0Error(), fitNLumi_);
1476   //    hs[k_z0_lumi]->ShiftFillLast(bs.z0(), bs.z0Error(), fitNLumi_);
1477   //  }
1478 
1479   // Fill summary report
1480   //  if (countFitting) {
1481   for (int n = 0; n < nFitElements_; n++) {
1482     reportSummaryContents[n]->Fill(summaryContent_[n] / (float)nFits_);
1483   }
1484 
1485   summarySum_ = 0;
1486   for (int ii = 0; ii < nFitElements_; ii++) {
1487     summarySum_ += summaryContent_[ii];
1488   }
1489   reportSummary_ = summarySum_ / (nFitElements_ * nFits_);
1490   if (reportSummary)
1491     reportSummary->Fill(reportSummary_);
1492 
1493   for (int bi = 0; bi < nFitElements_; bi++) {
1494     reportSummaryMap->setBinContent(1, bi + 1, summaryContent_[bi] / (float)nFits_);
1495   }
1496   //  }
1497 
1498   if ((resetFitNLumi_ > 0 &&
1499        ((onlineMode_ &&
1500          countLumi_ == resetFitNLumi_) ||  //OR it should be currentLumi_ (if in sequence then does not mattar)
1501         (!onlineMode_ && countLumi_ == resetFitNLumi_))) ||
1502       (StartAverage_)) {
1503     edm::LogInfo("FakeBeamMonitor") << "FitAndFill:: The flag is ON for running average Beam Spot fit" << endl;
1504     StartAverage_ = true;
1505     firstAverageFit_++;
1506     resetHistos_ = true;
1507     nthBSTrk_ = 0;
1508     beginLumiOfBSFit_ = 0;
1509     refBStime[0] = 0;
1510   }
1511 }
1512 
1513 //--------------------------------------------------------
1514 void FakeBeamMonitor::RestartFitting() {  //
1515   if (debug_)
1516     edm::LogInfo("FakeBeamMonitor")
1517         << " RestartingFitting:: Restart Beami everything to a fresh start !!! because Gap is > 10 LS" << endl;
1518   //track based fit reset here
1519   resetHistos_ = true;
1520   nthBSTrk_ = 0;
1521   //  theBeamFitter->resetTrkVector();
1522   //  theBeamFitter->resetLSRange();
1523   //  theBeamFitter->resetRefTime();
1524   //  theBeamFitter->resetPVFitter();
1525   //  theBeamFitter->resetCutFlow();
1526   beginLumiOfBSFit_ = 0;
1527   refBStime[0] = 0;
1528   //pv based fit iis reset here
1529   h_PVx[0]->Reset();
1530   h_PVy[0]->Reset();
1531   h_PVz[0]->Reset();
1532   beginLumiOfPVFit_ = 0;
1533   refPVtime[0] = 0;
1534   //Clear all the Maps here
1535   mapPVx.clear();
1536   mapPVy.clear();
1537   mapPVz.clear();
1538   mapNPV.clear();
1539   mapBeginBSLS.clear();
1540   mapBeginPVLS.clear();
1541   mapBeginBSTime.clear();
1542   mapBeginPVTime.clear();
1543   mapLSBSTrkSize.clear();
1544   mapLSPVStoreSize.clear();
1545   mapLSCF.clear();
1546   countGapLumi_ = 0;
1547   countLumi_ = 0;
1548   StartAverage_ = false;
1549 }
1550 
1551 //-------------------------------------------------------
1552 void FakeBeamMonitor::dqmEndRun(const Run& r, const EventSetup& context) {
1553   if (debug_)
1554     edm::LogInfo("FakeBeamMonitor") << "dqmEndRun:: Clearing all the Maps " << endl;
1555   //Clear all the Maps here
1556   mapPVx.clear();
1557   mapPVy.clear();
1558   mapPVz.clear();
1559   mapNPV.clear();
1560   mapBeginBSLS.clear();
1561   mapBeginPVLS.clear();
1562   mapBeginBSTime.clear();
1563   mapBeginPVTime.clear();
1564   mapLSBSTrkSize.clear();
1565   mapLSPVStoreSize.clear();
1566   mapLSCF.clear();
1567   if (useLockRecords_ && onlineDbService_.isAvailable()) {
1568     onlineDbService_->releaseLocks();
1569   }
1570 }
1571 
1572 //--------------------------------------------------------
1573 void FakeBeamMonitor::scrollTH1(TH1* h, time_t ref) {
1574   char offsetTime[64];
1575   formatFitTime(offsetTime, ref);
1576   TDatime da(offsetTime);
1577   if (lastNZbin > 0) {
1578     double val = h->GetBinContent(lastNZbin);
1579     double valErr = h->GetBinError(lastNZbin);
1580     h->Reset();
1581     h->GetXaxis()->SetTimeOffset(da.Convert(kTRUE));
1582     int bin = (lastNZbin > buffTime ? buffTime : 1);
1583     h->SetBinContent(bin, val);
1584     h->SetBinError(bin, valErr);
1585   } else {
1586     h->Reset();
1587     h->GetXaxis()->SetTimeOffset(da.Convert(kTRUE));
1588   }
1589 }
1590 
1591 //--------------------------------------------------------
1592 // Method to check whether to chane histogram time offset (forward only)
1593 bool FakeBeamMonitor::testScroll(time_t& tmpTime_, time_t& refTime_) {
1594   bool scroll_ = false;
1595   if (tmpTime_ - refTime_ >= intervalInSec_) {
1596     scroll_ = true;
1597     edm::LogInfo("FakeBeamMonitor") << "testScroll::  Reset Time Offset" << std::endl;
1598     lastNZbin = intervalInSec_;
1599     for (int bin = intervalInSec_; bin >= 1; bin--) {
1600       if (hs[k_x0_time]->getBinContent(bin) > 0) {
1601         lastNZbin = bin;
1602         break;
1603       }
1604     }
1605     edm::LogInfo("FakeBeamMonitor") << "testScroll::  Last non zero bin = " << lastNZbin << std::endl;
1606     if (tmpTime_ - refTime_ >= intervalInSec_ + lastNZbin) {
1607       edm::LogInfo("FakeBeamMonitor") << "testScroll::  Time difference too large since last readout" << std::endl;
1608       lastNZbin = 0;
1609       refTime_ = tmpTime_ - buffTime;
1610     } else {
1611       edm::LogInfo("FakeBeamMonitor") << "testScroll::  Offset to last record" << std::endl;
1612       int offset = ((lastNZbin > buffTime) ? (lastNZbin - buffTime) : (lastNZbin - 1));
1613       refTime_ += offset;
1614     }
1615   }
1616   return scroll_;
1617 }
1618 
1619 DEFINE_FWK_MODULE(FakeBeamMonitor);
1620 
1621 // Local Variables:
1622 // show-trailing-whitespace: t
1623 // truncate-lines: t
1624 // End:
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!