CMS 3D CMS Logo

CorrPCCProducer.cc
Go to the documentation of this file.
1 
11 #include <memory>
12 #include <string>
13 #include <vector>
14 #include <boost/serialization/vector.hpp>
15 #include <iostream>
16 #include <map>
17 #include <utility>
42 #include "TMath.h"
43 #include "TH1.h"
44 #include "TGraph.h"
45 #include "TGraphErrors.h"
46 #include "TFile.h"
47 
51 
52 class DQMStore;
53 class MonitorElement;
54 
55 class CorrPCCProducer : public one::DQMEDAnalyzer<edm::one::WatchLuminosityBlocks> {
56 public:
57  explicit CorrPCCProducer(const edm::ParameterSet&);
58  ~CorrPCCProducer() override;
59 
60 private:
61  void beginLuminosityBlock(edm::LuminosityBlock const& lumiSeg, const edm::EventSetup& iSetup) final;
62  void endLuminosityBlock(edm::LuminosityBlock const& lumiSeg, const edm::EventSetup& iSetup) final;
63  void endRun(edm::Run const& runSeg, const edm::EventSetup& iSetup) final;
64  void endRunProduce(edm::Run& runSeg, const edm::EventSetup& iSetup) final;
65  void endJob() final;
66 
67  void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
68 
70  float getMaximum(std::vector<float>);
71  void estimateType1Frac(std::vector<float>, float&);
72  void evaluateCorrectionResiduals(std::vector<float>);
73  void calculateCorrections(std::vector<float>, std::vector<float>&, float&);
74  void resetBlock();
75 
77  std::string pccSrc_; //input file EDproducer module label
78  std::string prodInst_; //input file product instance
79 
80  std::vector<float> rawlumiBX_; //new vector containing clusters per bxid
81  std::vector<float> errOnLumiByBX_; //standard error per bx
82  std::vector<float> totalLumiByBX_; //summed lumi
83  std::vector<float> totalLumiByBXAvg_; //summed lumi
84  std::vector<float> events_; //Number of events in each BX
85  std::vector<float> correctionTemplate_;
86  std::vector<float> correctionScaleFactors_; //list of scale factors to apply.
87  float overallCorrection_; //The Overall correction to the integrated luminosity
88 
89  unsigned int iBlock = 0;
90  unsigned int minimumNumberOfEvents;
91 
92  std::map<unsigned int, LumiInfo*> lumiInfoMapPerLS;
93  std::vector<unsigned int> lumiSections;
94  std::map<std::pair<unsigned int, unsigned int>, LumiInfo*>::iterator lumiInfoMapIterator;
95  std::map<std::pair<unsigned int, unsigned int>, LumiInfo*>
96  lumiInfoMap; //map to obtain iov for lumiOb corrections to the luminosity.
97  std::map<std::pair<unsigned int, unsigned int>, unsigned int> lumiInfoCounter; // number of lumiSections in this block
98 
101  TH1F* lumiAvg_h;
105 
106  unsigned int maxLS = 3500;
110 
111  TGraphErrors* type1FracGraph;
112  TGraphErrors* type1resGraph;
113  TGraphErrors* type2resGraph;
114  TList* hlist; //list for the clusters and corrections
115  TFile* histoFile;
116 
117  float type1Frac;
118  float mean_type1_residual; //Type 1 residual
119  float mean_type2_residual; //Type 2 residual
120  float mean_type1_residual_unc; //Type 1 residual uncertainty rms
121  float mean_type2_residual_unc; //Type 2 residual uncertainty rms
122  unsigned int nTrain; //Number of bunch trains used in calc type 1 and 2 res, frac.
123  unsigned int countLumi_; //The lumisection count... the size of the lumiblock
124  unsigned int approxLumiBlockSize_; //The number of lumisections per block.
125  unsigned int thisLS; //Ending lumisection for the iov that we save with the lumiInfo object.
126 
127  double type2_a_; //amplitude for the type 2 correction
128  double type2_b_; //decay width for the type 2 correction
129 
131 
133 };
134 
135 //--------------------------------------------------------------------------------------------------
137  pccSrc_ =
138  iConfig.getParameter<edm::ParameterSet>("CorrPCCProducerParameters").getParameter<std::string>("inLumiObLabel");
139  prodInst_ =
140  iConfig.getParameter<edm::ParameterSet>("CorrPCCProducerParameters").getParameter<std::string>("ProdInst");
142  iConfig.getParameter<edm::ParameterSet>("CorrPCCProducerParameters").getParameter<int>("approxLumiBlockSize");
143  type2_a_ = iConfig.getParameter<edm::ParameterSet>("CorrPCCProducerParameters").getParameter<double>("type2_a");
144  type2_b_ = iConfig.getParameter<edm::ParameterSet>("CorrPCCProducerParameters").getParameter<double>("type2_b");
145  countLumi_ = 0;
146  minimumNumberOfEvents = 1000;
147 
153 
154  resetBlock();
155 
157 
158  edm::InputTag inputPCCTag_(pccSrc_, prodInst_);
159 
160  lumiInfoToken = consumes<LumiInfo, edm::InLumi>(inputPCCTag_);
161 
162  histoFile = new TFile("CorrectionHisto.root", "RECREATE");
163 
164  type1FracGraph = new TGraphErrors();
165  type1resGraph = new TGraphErrors();
166  type2resGraph = new TGraphErrors();
167  type1FracGraph->SetName("Type1Fraction");
168  type1resGraph->SetName("Type1Res");
169  type2resGraph->SetName("Type2Res");
170  type1FracGraph->GetYaxis()->SetTitle("Type 1 Fraction");
171  type1resGraph->GetYaxis()->SetTitle("Type 1 Residual");
172  type2resGraph->GetYaxis()->SetTitle("Type 2 Residual");
173  type1FracGraph->GetXaxis()->SetTitle("Unique LS ID");
174  type1resGraph->GetXaxis()->SetTitle("Unique LS ID");
175  type2resGraph->GetXaxis()->SetTitle("Unique LS ID");
176  type1FracGraph->SetMarkerStyle(8);
177  type1resGraph->SetMarkerStyle(8);
178  type2resGraph->SetMarkerStyle(8);
179 }
180 
181 //--------------------------------------------------------------------------------------------------
183 
184 //--------------------------------------------------------------------------------------------------
185 // This method builds the single bunch response given an exponential function (type 2 only)
187  for (unsigned int bx = 1; bx < LumiConstants::numBX; bx++) {
188  correctionTemplate_.at(bx) = type2_a_ * exp(-(float(bx) - 1) * type2_b_);
189  }
190 }
191 
192 //--------------------------------------------------------------------------------------------------
193 // Finds max lumi value
194 float CorrPCCProducer::getMaximum(std::vector<float> lumi_vector) {
195  float max_lumi = 0;
196  for (size_t i = 0; i < lumi_vector.size(); i++) {
197  if (lumi_vector.at(i) > max_lumi)
198  max_lumi = lumi_vector.at(i);
199  }
200  return max_lumi;
201 }
202 
203 //--------------------------------------------------------------------------------------------------
204 // This method takes luminosity from the last bunch in a train and makes a comparison with
205 // the follow non-active bunch crossing to estimate the spill over fraction (type 1 afterglow).
206 void CorrPCCProducer::estimateType1Frac(std::vector<float> uncorrPCCPerBX, float& type1Frac) {
207  std::vector<float> corrected_tmp_;
208  for (size_t i = 0; i < uncorrPCCPerBX.size(); i++) {
209  corrected_tmp_.push_back(uncorrPCCPerBX.at(i));
210  }
211 
212  //Apply initial type 1 correction
213  for (size_t k = 0; k < LumiConstants::numBX - 1; k++) {
214  float bin_k = corrected_tmp_.at(k);
215  corrected_tmp_.at(k + 1) = corrected_tmp_.at(k + 1) - type1Frac * bin_k;
216  }
217 
218  //Apply type 2 correction
219  for (size_t i = 0; i < LumiConstants::numBX - 1; i++) {
220  for (size_t j = i + 1; j < i + LumiConstants::numBX - 1; j++) {
221  float bin_i = corrected_tmp_.at(i);
222  if (j < LumiConstants::numBX) {
223  corrected_tmp_.at(j) = corrected_tmp_.at(j) - bin_i * correctionTemplate_.at(j - i);
224  } else {
225  corrected_tmp_.at(j - LumiConstants::numBX) =
226  corrected_tmp_.at(j - LumiConstants::numBX) - bin_i * correctionTemplate_.at(j - i);
227  }
228  }
229  }
230 
231  //Apply additional iteration for type 1 correction
232  evaluateCorrectionResiduals(corrected_tmp_);
233  type1Frac += mean_type1_residual;
234 }
235 
236 //--------------------------------------------------------------------------------------------------
237 void CorrPCCProducer::evaluateCorrectionResiduals(std::vector<float> corrected_tmp_) {
238  float lumiMax = getMaximum(corrected_tmp_);
239  float threshold = lumiMax * 0.2;
240 
243  nTrain = 0;
244  float lumi = 0;
245  std::vector<float> afterGlow;
246  TH1F type1("type1", "", 1000, -0.5, 0.5);
247  TH1F type2("type2", "", 1000, -0.5, 0.5);
248  for (size_t ibx = 2; ibx < LumiConstants::numBX - 5; ibx++) {
249  lumi = corrected_tmp_.at(ibx);
250  afterGlow.clear();
251  afterGlow.push_back(corrected_tmp_.at(ibx + 1));
252  afterGlow.push_back(corrected_tmp_.at(ibx + 2));
253 
254  //Where type 1 and type 2 residuals are computed
255  if (lumi > threshold && afterGlow[0] < threshold && afterGlow[1] < threshold) {
256  for (int index = 3; index < 6; index++) {
257  float thisAfterGlow = corrected_tmp_.at(ibx + index);
258  if (thisAfterGlow < threshold) {
259  afterGlow.push_back(thisAfterGlow);
260  } else {
261  break;
262  }
263  }
264  float thisType1 = 0;
265  float thisType2 = 0;
266  if (afterGlow.size() > 1) {
267  int nAfter = 0;
268  for (unsigned int index = 1; index < afterGlow.size(); index++) {
269  thisType2 += afterGlow[index];
270  type2.Fill(afterGlow[index] / lumi);
271  nAfter++;
272  }
273  thisType2 /= nAfter;
274  }
275  thisType1 = (afterGlow[0] - thisType2) / lumi;
276 
277  type1.Fill(thisType1);
278  nTrain += 1;
279  }
280  }
281 
282  mean_type1_residual = type1.GetMean(); //Calculate the mean value of the type 1 residual
283  mean_type2_residual = type2.GetMean(); //Calculate the mean value of the type 2 residual
284  mean_type1_residual_unc = type1.GetMeanError();
285  mean_type2_residual_unc = type2.GetMeanError();
286 
287  histoFile->cd();
288  type1.Write();
289  type2.Write();
290 }
291 
292 //--------------------------------------------------------------------------------------------------
293 void CorrPCCProducer::calculateCorrections(std::vector<float> uncorrected,
294  std::vector<float>& correctionScaleFactors_,
295  float& overallCorrection_) {
296  type1Frac = 0;
297 
298  int nTrials = 4;
299 
300  for (int trial = 0; trial < nTrials; trial++) {
301  estimateType1Frac(uncorrected, type1Frac);
302  edm::LogInfo("INFO") << "type 1 fraction after iteration " << trial << " is " << type1Frac;
303  }
304 
305  //correction should never be negative
306  type1Frac = std::max(0.0, (double)type1Frac);
307 
308  std::vector<float> corrected_tmp_;
309  for (size_t i = 0; i < uncorrected.size(); i++) {
310  corrected_tmp_.push_back(uncorrected.at(i));
311  }
312 
313  //Apply all corrections
314  for (size_t i = 0; i < LumiConstants::numBX - 1; i++) {
315  // type 1 - first (spillover from previous BXs real clusters)
316  float bin_i = corrected_tmp_.at(i);
317  corrected_tmp_.at(i + 1) = corrected_tmp_.at(i + 1) - type1Frac * bin_i;
318 
319  // type 2 - after (comes from real activation)
320  bin_i = corrected_tmp_.at(i);
321  for (size_t j = i + 1; j < i + LumiConstants::numBX - 1; j++) {
322  if (j < LumiConstants::numBX) {
323  corrected_tmp_.at(j) = corrected_tmp_.at(j) - bin_i * correctionTemplate_.at(j - i);
324  } else {
325  corrected_tmp_.at(j - LumiConstants::numBX) =
326  corrected_tmp_.at(j - LumiConstants::numBX) - bin_i * correctionTemplate_.at(j - i);
327  }
328  }
329  }
330 
331  evaluateCorrectionResiduals(corrected_tmp_);
332 
333  float integral_uncorr_clusters = 0;
334  float integral_corr_clusters = 0;
335 
336  //Calculate Per-BX correction factor and overall correction factor
337  float lumiMax = getMaximum(corrected_tmp_);
338  float threshold = lumiMax * 0.2;
339  for (size_t ibx = 0; ibx < corrected_tmp_.size(); ibx++) {
340  if (corrected_tmp_.at(ibx) > threshold) {
341  integral_uncorr_clusters += uncorrected.at(ibx);
342  integral_corr_clusters += corrected_tmp_.at(ibx);
343  }
344  if (corrected_tmp_.at(ibx) != 0.0 && uncorrected.at(ibx) != 0.0) {
345  correctionScaleFactors_.at(ibx) = corrected_tmp_.at(ibx) / uncorrected.at(ibx);
346  } else {
347  correctionScaleFactors_.at(ibx) = 0.0;
348  }
349  }
350 
351  overallCorrection_ = integral_corr_clusters / integral_uncorr_clusters;
352 }
353 
354 //--------------------------------------------------------------------------------------------------
356  countLumi_++;
357 }
358 
359 //--------------------------------------------------------------------------------------------------
361  thisLS = lumiSeg.luminosityBlock();
362 
363  edm::Handle<LumiInfo> PCCHandle;
364  lumiSeg.getByToken(lumiInfoToken, PCCHandle);
365 
366  const LumiInfo& inLumiOb = *(PCCHandle.product());
367 
368  errOnLumiByBX_ = inLumiOb.getErrorLumiAllBX();
369 
370  unsigned int totalEvents = 0;
371  for (unsigned int bx = 0; bx < LumiConstants::numBX; bx++) {
372  totalEvents += errOnLumiByBX_[bx];
373  }
374 
375  if (totalEvents < minimumNumberOfEvents) {
376  edm::LogInfo("INFO") << "number of events in this LS is too few " << totalEvents;
377  //return;
378  } else {
379  edm::LogInfo("INFO") << "Skipping Lumisection " << thisLS;
380  }
381 
383  totalLumiByBX_ = inLumiOb.getInstLumiAllBX();
384  events_ = inLumiOb.getErrorLumiAllBX();
385  lumiInfoMapPerLS[thisLS]->setInstLumiAllBX(totalLumiByBX_);
386  lumiInfoMapPerLS[thisLS]->setErrorLumiAllBX(events_);
387  lumiSections.push_back(thisLS);
388 }
389 
390 //--------------------------------------------------------------------------------------------------
391 void CorrPCCProducer::endRun(edm::Run const& runSeg, const edm::EventSetup& iSetup) {}
392 
393 //--------------------------------------------------------------------------------------------------
395  if (lumiSections.empty()) {
396  return;
397  }
398 
399  std::sort(lumiSections.begin(), lumiSections.end());
400 
401  edm::LogInfo("INFO") << "Number of Lumisections " << lumiSections.size() << " in run " << runSeg.run();
402 
403  //Determining integer number of blocks
404  float nBlocks_f = float(lumiSections.size()) / approxLumiBlockSize_;
405  unsigned int nBlocks = 1;
406  if (nBlocks_f > 1) {
407  if (nBlocks_f - lumiSections.size() / approxLumiBlockSize_ < 0.5) {
408  nBlocks = lumiSections.size() / approxLumiBlockSize_;
409  } else {
410  nBlocks = lumiSections.size() / approxLumiBlockSize_ + 1;
411  }
412  }
413 
414  float nLSPerBlock = float(lumiSections.size()) / nBlocks;
415 
416  std::vector<std::pair<unsigned int, unsigned int>> lsKeys;
417  lsKeys.clear();
418 
419  //Constructing nBlocks IOVs
420  for (unsigned iKey = 0; iKey < nBlocks; iKey++) {
421  lsKeys.push_back(std::make_pair(lumiSections[(unsigned int)(iKey * nLSPerBlock)],
422  lumiSections[(unsigned int)((iKey + 1) * nLSPerBlock) - 1]));
423  }
424 
425  lsKeys[0].first = 1;
426 
427  for (unsigned int lumiSection = 0; lumiSection < lumiSections.size(); lumiSection++) {
428  thisLS = lumiSections[lumiSection];
429 
430  std::pair<unsigned int, unsigned int> lsKey;
431 
432  bool foundKey = false;
433 
434  for (unsigned iKey = 0; iKey < nBlocks; iKey++) {
435  if ((thisLS >= lsKeys[iKey].first) && (thisLS <= lsKeys[iKey].second)) {
436  lsKey = lsKeys[iKey];
437  foundKey = true;
438  break;
439  }
440  }
441 
442  if (!foundKey) {
443  edm::LogInfo("WARNING") << "Didn't find key " << thisLS;
444  continue;
445  }
446 
447  if (lumiInfoMap.count(lsKey) == 0) {
448  lumiInfoMap[lsKey] = new LumiInfo();
449  }
450 
451  //Sum all lumi in IOV of lsKey
452  totalLumiByBX_ = lumiInfoMap[lsKey]->getInstLumiAllBX();
453  events_ = lumiInfoMap[lsKey]->getErrorLumiAllBX();
454 
455  rawlumiBX_ = lumiInfoMapPerLS[thisLS]->getInstLumiAllBX();
456  errOnLumiByBX_ = lumiInfoMapPerLS[thisLS]->getErrorLumiAllBX();
457 
458  for (unsigned int bx = 0; bx < LumiConstants::numBX; bx++) {
459  totalLumiByBX_[bx] += rawlumiBX_[bx];
460  events_[bx] += errOnLumiByBX_[bx];
461  }
462  lumiInfoMap[lsKey]->setInstLumiAllBX(totalLumiByBX_);
463  lumiInfoMap[lsKey]->setErrorLumiAllBX(events_);
464  lumiInfoCounter[lsKey]++;
465  }
466 
467  cond::Time_t thisIOV = 1;
468 
469  char* histname1 = new char[100];
470  char* histname2 = new char[100];
471  char* histname3 = new char[100];
472  char* histTitle1 = new char[100];
473  char* histTitle2 = new char[100];
474  char* histTitle3 = new char[100];
475  sprintf(histTitle1, "Type1Fraction_%d", runSeg.run());
476  sprintf(histTitle2, "Type1Res_%d", runSeg.run());
477  sprintf(histTitle3, "Type2Res_%d", runSeg.run());
478  type1FracHist = new TH1F(histTitle1, histTitle1, 1000, -0.5, 0.5);
479  type1resHist = new TH1F(histTitle2, histTitle2, 4000, -0.2, 0.2);
480  type2resHist = new TH1F(histTitle3, histTitle3, 4000, -0.2, 0.2);
481  delete[] histTitle1;
482  delete[] histTitle2;
483  delete[] histTitle3;
484 
486  totalLumiByBX_ = lumiInfoMapIterator->second->getInstLumiAllBX();
487  events_ = lumiInfoMapIterator->second->getErrorLumiAllBX();
488 
489  if (events_.empty()) {
490  continue;
491  }
492 
494  thisIOV = (cond::Time_t)(lu.value());
495 
496  sprintf(histname1,
497  "CorrectedLumiAvg_%d_%d_%d_%d",
498  runSeg.run(),
499  iBlock,
500  lumiInfoMapIterator->first.first,
501  lumiInfoMapIterator->first.second);
502  sprintf(histname2,
503  "ScaleFactorsAvg_%d_%d_%d_%d",
504  runSeg.run(),
505  iBlock,
506  lumiInfoMapIterator->first.first,
507  lumiInfoMapIterator->first.second);
508  sprintf(histname3,
509  "RawLumiAvg_%d_%d_%d_%d",
510  runSeg.run(),
511  iBlock,
512  lumiInfoMapIterator->first.first,
513  lumiInfoMapIterator->first.second);
514 
515  corrlumiAvg_h = new TH1F(histname1, "", LumiConstants::numBX, 1, LumiConstants::numBX);
516  scaleFactorAvg_h = new TH1F(histname2, "", LumiConstants::numBX, 1, LumiConstants::numBX);
517  lumiAvg_h = new TH1F(histname3, "", LumiConstants::numBX, 1, LumiConstants::numBX);
518 
519  //Averaging by the number of events
520  for (unsigned int i = 0; i < LumiConstants::numBX; i++) {
521  if (events_.at(i) != 0) {
523  } else {
524  totalLumiByBXAvg_[i] = 0.0;
525  }
526  }
527 
529 
530  for (unsigned int bx = 0; bx < LumiConstants::numBX; bx++) {
531  corrlumiAvg_h->SetBinContent(bx, totalLumiByBXAvg_[bx] * correctionScaleFactors_[bx]);
532  if (events_.at(bx) != 0) {
533  corrlumiAvg_h->SetBinError(bx,
534  totalLumiByBXAvg_[bx] * correctionScaleFactors_[bx] / TMath::Sqrt(events_.at(bx)));
535  } else {
536  corrlumiAvg_h->SetBinError(bx, 0.0);
537  }
538 
539  scaleFactorAvg_h->SetBinContent(bx, correctionScaleFactors_[bx]);
540  lumiAvg_h->SetBinContent(bx, totalLumiByBXAvg_[bx]);
541  }
542 
543  //Writing the corrections to SQL lite file for db
550 
551  if (poolDbService.isAvailable()) {
552  poolDbService->writeOne<LumiCorrections>(pccCorrections, thisIOV, "LumiCorrectionsRcd");
553  } else {
554  throw std::runtime_error("PoolDBService required.");
555  }
556 
557  delete pccCorrections;
558 
559  histoFile->cd();
560  corrlumiAvg_h->Write();
561  scaleFactorAvg_h->Write();
562  lumiAvg_h->Write();
563 
564  delete corrlumiAvg_h;
565  delete scaleFactorAvg_h;
566  delete lumiAvg_h;
567 
568  type1FracHist->Fill(type1Frac);
571 
572  for (unsigned int ils = lumiInfoMapIterator->first.first; ils < lumiInfoMapIterator->first.second + 1; ils++) {
573  if (ils > maxLS) {
574  std::cout << "ils out of maxLS range!!" << std::endl;
575  break;
576  }
583  }
584 
585  type1FracGraph->SetPoint(iBlock, thisIOV + approxLumiBlockSize_ / 2.0, type1Frac);
586  type1resGraph->SetPoint(iBlock, thisIOV + approxLumiBlockSize_ / 2.0, mean_type1_residual);
587  type2resGraph->SetPoint(iBlock, thisIOV + approxLumiBlockSize_ / 2.0, mean_type2_residual);
591 
592  edm::LogInfo("INFO")
593  << "iBlock type1Frac mean_type1_residual mean_type2_residual mean_type1_residual_unc mean_type2_residual_unc "
594  << iBlock << " " << type1Frac << " " << mean_type1_residual << " " << mean_type2_residual << " "
596 
597  type1Frac = 0.0;
598  mean_type1_residual = 0.0;
599  mean_type2_residual = 0.0;
601  mean_type2_residual_unc = 0;
602 
603  iBlock++;
604 
605  resetBlock();
606  }
607  histoFile->cd();
608  type1FracHist->Write();
609  type1resHist->Write();
610  type2resHist->Write();
611 
612  delete type1FracHist;
613  delete type1resHist;
614  delete type2resHist;
615 
616  delete[] histname1;
617  delete[] histname2;
618  delete[] histname3;
619 
621  delete lumiInfoMapIterator->second;
622  }
623  for (unsigned int lumiSection = 0; lumiSection < lumiSections.size(); lumiSection++) {
624  thisLS = lumiSections[lumiSection];
625  delete lumiInfoMapPerLS[thisLS];
626  }
627  lumiInfoMap.clear();
628  lumiInfoMapPerLS.clear();
629  lumiSections.clear();
630  lumiInfoCounter.clear();
631 }
632 
633 //--------------------------------------------------------------------------------------------------
635  for (unsigned int bx = 0; bx < LumiConstants::numBX; bx++) {
636  totalLumiByBX_[bx] = 0;
637  totalLumiByBXAvg_[bx] = 0;
638  events_[bx] = 0;
639  correctionScaleFactors_[bx] = 1.0;
640  }
641 }
642 //--------------------------------------------------------------------------------------------------
643 void CorrPCCProducer::bookHistograms(DQMStore::IBooker& ibooker, edm::Run const& iRun, edm::EventSetup const& context) {
644  ibooker.setCurrentFolder("AlCaReco/LumiPCC/");
645  Type1FracMon = ibooker.book1D("type1Fraction", "Type1Fraction;Lumisection;Type 1 Fraction", maxLS, 0, maxLS);
646  Type1ResMon = ibooker.book1D("type1Residual", "Type1Residual;Lumisection;Type 1 Residual", maxLS, 0, maxLS);
647  Type2ResMon = ibooker.book1D("type2Residual", "Type2Residual;Lumisection;Type 2 Residual", maxLS, 0, maxLS);
648 }
649 //--------------------------------------------------------------------------------------------------
651  histoFile->cd();
652  type1FracGraph->Write();
653  type1resGraph->Write();
654  type2resGraph->Write();
655  histoFile->Write();
656  histoFile->Close();
657  delete type1FracGraph;
658  delete type1resGraph;
659  delete type2resGraph;
660 }
661 
T getParameter(std::string const &) const
std::vector< float > events_
void calculateCorrections(std::vector< float >, std::vector< float > &, float &)
MonitorElement * Type1FracMon
void setBinContent(int binx, double content)
set content of bin (1-D)
RunID const & id() const
Definition: RunBase.h:39
std::map< std::pair< unsigned int, unsigned int >, LumiInfo * > lumiInfoMap
RunNumber_t run() const
Definition: RunBase.h:40
void setType1Fraction(float type1frac)
RunNumber_t run() const
Definition: RunID.h:39
bool getByToken(EDGetToken token, Handle< PROD > &result) const
unsigned int nTrain
edm::EDGetTokenT< LumiInfo > lumiInfoToken
void endRun(edm::Run const &runSeg, const edm::EventSetup &iSetup) final
void setOverallCorrection(float overallCorrection)
static const unsigned int numBX
Definition: LumiConstants.h:9
const std::vector< float > & getErrorLumiAllBX() const
Definition: LumiInfo.h:104
void endJob() final
unsigned int iBlock
unsigned int LuminosityBlockNumber_t
MonitorElement * Type1ResMon
LuminosityBlockNumber_t luminosityBlock() const
unsigned int countLumi_
void beginLuminosityBlock(edm::LuminosityBlock const &lumiSeg, const edm::EventSetup &iSetup) final
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
unsigned long long Time_t
Definition: Time.h:16
void endRunProduce(edm::Run &runSeg, const edm::EventSetup &iSetup) final
std::string pccSrc_
void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:268
unsigned int minimumNumberOfEvents
std::vector< float > totalLumiByBXAvg_
std::map< std::pair< unsigned int, unsigned int >, unsigned int > lumiInfoCounter
void evaluateCorrectionResiduals(std::vector< float >)
TGraphErrors * type1resGraph
LumiCorrections * pccCorrections
std::string prodInst_
bool isAvailable() const
Definition: Service.h:40
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:106
void writeOne(T *payload, Time_t time, const std::string &recordName, bool withlogging=false)
void setBinError(int binx, double error)
set uncertainty on content of bin (1-D)
std::vector< float > correctionScaleFactors_
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
void setType2Residual(float type2res)
const std::vector< float > & getInstLumiAllBX() const
Definition: LumiInfo.h:100
std::vector< float > correctionTemplate_
std::vector< unsigned int > lumiSections
float getMaximum(std::vector< float >)
int k[5][pyjets_maxn]
TGraphErrors * type2resGraph
std::vector< float > rawlumiBX_
T const * product() const
Definition: Handle.h:74
std::vector< float > totalLumiByBX_
edm::Service< cond::service::PoolDBOutputService > poolDbService
TGraphErrors * type1FracGraph
unsigned int approxLumiBlockSize_
void endLuminosityBlock(edm::LuminosityBlock const &lumiSeg, const edm::EventSetup &iSetup) final
void makeCorrectionTemplate()
void estimateType1Frac(std::vector< float >, float &)
unsigned int maxLS
~CorrPCCProducer() override
std::map< std::pair< unsigned int, unsigned int >, LumiInfo * >::iterator lumiInfoMapIterator
void setType1Residual(float type1res)
CorrPCCProducer(const edm::ParameterSet &)
std::vector< float > errOnLumiByBX_
MonitorElement * Type2ResMon
void setCorrectionsBX(std::vector< float > &correctBX)
Definition: Run.h:45
unsigned int thisLS
std::map< unsigned int, LumiInfo * > lumiInfoMapPerLS