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