CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
PixelLumiDQM.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 
3 // Package: PixelLumiDQM
4 // Class: PixelLumiDQM
5 
6 // Author: Amita Raval
7 // Based on Jeroen Hegeman's code for Pixel Cluster Count Luminosity
8 
9 #include "PixelLumiDQM.h"
10 
36 
37 #include <ctime>
38 #include <fstream>
39 #include <map>
40 #include <string>
41 #include <sys/time.h>
42 #include <vector>
43 
44 const unsigned int PixelLumiDQM::lastBunchCrossing;
45 
46 // Constructors and destructor.
48  : fPixelClusterLabel(consumes<edmNew::DetSetVector<SiPixelCluster>>(
49  iConfig.getUntrackedParameter<edm::InputTag>("pixelClusterLabel", edm::InputTag("siPixelClusters")))),
50  tkGeomToken_(esConsumes()),
51  fIncludePixelClusterInfo(iConfig.getUntrackedParameter<bool>("includePixelClusterInfo", true)),
52  fIncludePixelQualCheckHistos(iConfig.getUntrackedParameter<bool>("includePixelQualCheckHistos", true)),
53  fResetIntervalInLumiSections(iConfig.getUntrackedParameter<int>("resetEveryNLumiSections", 1)),
54  fDeadModules(iConfig.getUntrackedParameter<std::vector<uint32_t>>("deadModules", std::vector<uint32_t>())),
55  fMinPixelsPerCluster(iConfig.getUntrackedParameter<int>("minNumPixelsPerCluster", 0)),
56  fMinClusterCharge(iConfig.getUntrackedParameter<double>("minChargePerCluster", 0)),
57  bunchTriggerMask(lastBunchCrossing + 1, false),
58  filledAndUnmaskedBunches(0),
59  useInnerBarrelLayer(iConfig.getUntrackedParameter<bool>("useInnerBarrelLayer", false)),
60  fLogFileName_(iConfig.getUntrackedParameter<std::string>("logFileName", "/tmp/pixel_lumi.txt")) {
61  edm::LogInfo("Configuration") << "PixelLumiDQM looking for pixel clusters in '"
62  << iConfig.getUntrackedParameter<edm::InputTag>("pixelClusterLabel",
63  edm::InputTag("siPixelClusters"))
64  << "'";
65  edm::LogInfo("Configuration") << "PixelLumiDQM storing pixel cluster info? " << fIncludePixelClusterInfo;
66  edm::LogInfo("Configuration") << "PixelLumiDQM storing pixel cluster quality check histograms? "
68 
69  if (fDeadModules.empty()) {
70  edm::LogInfo("Configuration") << "No pixel modules specified to be ignored";
71  } else {
72  edm::LogInfo("Configuration") << fDeadModules.size() << " pixel modules specified to be ignored:";
73  for (std::vector<uint32_t>::const_iterator it = fDeadModules.begin(); it != fDeadModules.end(); ++it) {
74  edm::LogInfo("Configuration") << " " << *it;
75  }
76  }
77  edm::LogInfo("Configuration") << "Ignoring pixel clusters with less than " << fMinPixelsPerCluster << " pixels";
78  edm::LogInfo("Configuration") << "Ignoring pixel clusters with charge less than " << fMinClusterCharge;
79 }
80 
82 
85  desc.setUnknown();
86  descriptions.addDefault(desc);
87 }
88 
90  // Collect all bookkeeping information.
91  fRunNo = iEvent.id().run();
92  fEvtNo = iEvent.id().event();
94  fBXNo = iEvent.bunchCrossing();
95  fTimestamp = iEvent.time().unixTime();
98  // This serves as event counter to compute luminosity from cluster counts.
99  std::map<int, PixelClusterCount>::iterator it = fNumPixelClusters.find(fBXNo);
100  if (it == fNumPixelClusters.end())
102  // Find tracker geometry.
103  const TrackerGeometry *trackerGeo = &iSetup.getData(tkGeomToken_);
105  // Find pixel clusters.
107  iEvent.getByToken(fPixelClusterLabel, pixelClusters);
108 
109  // Loop over entire tracker geometry.
110  for (TrackerGeometry::DetContainer::const_iterator i = trackerGeo->dets().begin(); i != trackerGeo->dets().end();
111  ++i) {
112  // See if this is a pixel unit(?).
113 
114  if (GeomDetEnumerators::isTrackerPixel((*i)->subDetector())) {
115  DetId detId = (*i)->geographicalId();
116  // Find all clusters on this detector module.
117  edmNew::DetSetVector<SiPixelCluster>::const_iterator iSearch = pixelClusters->find(detId);
118  if (iSearch != pixelClusters->end()) {
119  // Count the number of clusters with at least a minimum
120  // number of pixels per cluster and at least a minimum charge.
121  size_t numClusters = 0;
122  for (edmNew::DetSet<SiPixelCluster>::const_iterator itClus = iSearch->begin(); itClus != iSearch->end();
123  ++itClus) {
124  if ((itClus->size() >= fMinPixelsPerCluster) && (itClus->charge() >= fMinClusterCharge)) {
125  ++numClusters;
126  }
127  }
128  // DEBUG DEBUG DEBUG
129  assert(numClusters <= iSearch->size());
130  // DEBUG DEBUG DEBUG end
131 
132  // Add up the cluster count based on the position of this detector
133  // element.
134  if (detId.subdetId() == PixelSubdetector::PixelBarrel) {
136  int layer = detName.layerName() - kOffsetLayers;
137  fNumPixelClusters[fBXNo].numB.at(layer) += numClusters;
138  fNumPixelClusters[fBXNo].dnumB.at(layer) += sqrt(numClusters);
139  } else {
140  // DEBUG DEBUG DEBUG
142  // DEBUG DEBUG DEBUG end
143 
145  PixelEndcapNameUpgrade::HalfCylinder halfCylinder = detName.halfCylinder();
146  int disk = detName.diskName() - kOffsetDisks;
147  switch (halfCylinder) {
150  fNumPixelClusters[fBXNo].numFM.at(disk) += numClusters;
151  fNumPixelClusters[fBXNo].dnumFM.at(disk) += sqrt(numClusters);
152  break;
155  fNumPixelClusters[fBXNo].numFP.at(disk) += numClusters;
156  fNumPixelClusters[fBXNo].dnumFP.at(disk) += sqrt(numClusters);
157  break;
158  default:
159  assert(false);
160  break;
161  }
162  }
163  }
164  }
165  }
166  }
167  // ----------
168 
169  // Fill some pixel cluster quality check histograms if requested.
171  // Find pixel clusters.
173  iEvent.getByToken(fPixelClusterLabel, pixelClusters);
174 
175  bool filterDeadModules = (!fDeadModules.empty());
176  std::vector<uint32_t>::const_iterator deadModulesBegin = fDeadModules.begin();
177  std::vector<uint32_t>::const_iterator deadModulesEnd = fDeadModules.end();
178 
179  // Loop over entire tracker geometry.
180  for (TrackerGeometry::DetContainer::const_iterator i = trackerGeo->dets().begin(); i != trackerGeo->dets().end();
181  ++i) {
182  // See if this is a pixel module.
183  if (GeomDetEnumerators::isTrackerPixel((*i)->subDetector())) {
184  DetId detId = (*i)->geographicalId();
185 
186  // Skip this module if it's on the list of modules to be ignored.
187  if (filterDeadModules && find(deadModulesBegin, deadModulesEnd, detId()) != deadModulesEnd) {
188  continue;
189  }
190 
191  // Find all clusters in this module.
192  edmNew::DetSetVector<SiPixelCluster>::const_iterator iSearch = pixelClusters->find(detId);
193 
194  // Loop over all clusters in this module.
195  if (iSearch != pixelClusters->end()) {
196  for (edmNew::DetSet<SiPixelCluster>::const_iterator clus = iSearch->begin(); clus != iSearch->end(); ++clus) {
197  if ((clus->size() >= fMinPixelsPerCluster) && (clus->charge() >= fMinClusterCharge)) {
198  PixelGeomDetUnit const *theGeomDet = dynamic_cast<PixelGeomDetUnit const *>(trackerGeo->idToDet(detId));
199  PixelTopology const *topol = &(theGeomDet->specificTopology());
200  double x = clus->x();
201  double y = clus->y();
202  LocalPoint clustLP = topol->localPosition(MeasurementPoint(x, y));
203  GlobalPoint clustGP = theGeomDet->surface().toGlobal(clustLP);
204  double charge = clus->charge() / 1.e3;
205  int size = clus->size();
206 
207  if (detId.subdetId() == PixelSubdetector::PixelBarrel) {
209  unsigned int layer = detName.layerName() - kOffsetLayers;
210  if (layer < kNumLayers) {
211  std::string histName;
212  histName = "clusPosBarrel" + std::to_string(layer);
213  fHistContainerThisRun[histName]->Fill(clustGP.z(), clustGP.phi());
214  histName = "clusChargeBarrel" + std::to_string(layer);
215  fHistContainerThisRun[histName]->Fill(iEvent.bunchCrossing(), charge);
216  histName = "clusSizeBarrel" + std::to_string(layer);
217  fHistContainerThisRun[histName]->Fill(iEvent.bunchCrossing(), size);
218  } else {
219  edm::LogWarning("pixelLumi") << "higher layer number, " << layer << ", than layers";
220  }
221  } else {
222  // DEBUG DEBUG DEBUG
224  // DEBUG DEBUG DEBUG end
225 
227  unsigned int disk = detName.diskName() - kOffsetDisks;
228  if (disk < kNumDisks) {
229  std::string histName;
230  histName = "clusPosEndCap" + std::to_string(disk);
231  fHistContainerThisRun[histName]->Fill(clustGP.x(), clustGP.y());
232  histName = "clusChargeEndCap" + std::to_string(disk);
233  fHistContainerThisRun[histName]->Fill(iEvent.bunchCrossing(), charge);
234  histName = "clusSizeEndCap" + std::to_string(disk);
235  fHistContainerThisRun[histName]->Fill(iEvent.bunchCrossing(), size);
236  } else {
237  edm::LogWarning("pixelLumi")
238  << "higher disk number, " << disk << ", than disks," << kNumDisks << std::endl;
239  }
240  }
241  }
242  }
243  }
244  }
245  }
246  }
247 }
248 
250  edm::Run const &run,
251  edm::EventSetup const & /* iSetup */) {
252  edm::LogInfo("Status") << "Starting processing of run #" << run.id().run();
253 
254  // Top folder containing high-level information about pixel and HF lumi.
255  std::string folder = "PixelLumi/";
256  ibooker.setCurrentFolder(folder);
257 
258  fHistTotalRecordedLumiByLS = ibooker.book1D("totalPixelLumiByLS", "Pixel Lumi in nb vs LS", 8000, 0.5, 8000.5);
259  fHistRecordedByBxCumulative = ibooker.book1D("PXLumiByBXsum",
260  "Pixel Lumi in nb by BX Cumulative vs LS",
262  0.5,
263  float(lastBunchCrossing) + 0.5);
264 
265  std::string subfolder = folder + "lastLS/";
266  ibooker.setCurrentFolder(subfolder);
268  "PXByBXLastLumi", "Pixel By BX Last Lumi", lastBunchCrossing + 1, -0.5, float(lastBunchCrossing) + 0.5);
269 
270  subfolder = folder + "ClusterCountingDetails/";
271  ibooker.setCurrentFolder(subfolder);
272 
273  fHistnBClusVsLS[0] = ibooker.book1D("nBClusVsLS_0", "Fraction of Clusters vs LS Barrel layer 0", 8000, 0.5, 8000.5);
274  fHistnBClusVsLS[1] = ibooker.book1D("nBClusVsLS_1", "Fraction of Clusters vs LS Barrel layer 1", 8000, 0.5, 8000.5);
275  fHistnBClusVsLS[2] = ibooker.book1D("nBClusVsLS_2", "Fraction of Clusters vs LS Barrel layer 2", 8000, 0.5, 8000.5);
276  fHistnFPClusVsLS[0] = ibooker.book1D("nFPClusVsLS_0", "Fraction of Clusters vs LS Barrel layer 0", 8000, 0.5, 8000.5);
277  fHistnFPClusVsLS[1] = ibooker.book1D("nFPClusVsLS_1", "Fraction of Clusters vs LS Barrel layer 1", 8000, 0.5, 8000.5);
278  fHistnFMClusVsLS[0] = ibooker.book1D("nFMClusVsLS_0", "Fraction of Clusters vs LS Barrel layer 0", 8000, 0.5, 8000.5);
279  fHistnFMClusVsLS[1] = ibooker.book1D("nFMClusVsLS_1", "Fraction of Clusters vs LS Barrel layer 1", 8000, 0.5, 8000.5);
280  fHistBunchCrossings = ibooker.book1D(
281  "BunchCrossings", "Cumulative Bunch Crossings", lastBunchCrossing, 0.5, float(lastBunchCrossing) + 0.5);
283  "BunchCrossingsLL", "Bunch Crossings Last Lumi", lastBunchCrossing, 0.5, float(lastBunchCrossing) + 0.5);
285  "ClusterCountByBxLL", "Cluster Count by BX Last Lumi", lastBunchCrossing, 0.5, float(lastBunchCrossing) + 0.5);
287  "ClusterCountByBxSum", "Cluster Count by BX Cumulative", lastBunchCrossing, 0.5, float(lastBunchCrossing) + 0.5);
288  fHistClusByLS = ibooker.book1D("totalClusByLS", "Number of Clusters all dets vs LS", 8000, 0.5, 8000.5);
289 
290  // Add some pixel cluster quality check histograms (in a subfolder).
291  subfolder = folder + "qualityChecks/";
292  ibooker.setCurrentFolder(subfolder);
293 
295  // Create histograms for this run if not already present in our list.
296  edm::LogInfo("Status") << "Creating histograms for run #" << run.id().run();
297 
298  // Pixel cluster positions in the barrel - (z, phi).
299  for (size_t i = 0; i <= kNumLayers; ++i) {
300  std::stringstream key;
301  key << "clusPosBarrel" << i;
302  std::stringstream name;
303  name << key.str() << "_" << run.run();
304  std::stringstream title;
305  title << "Pixel cluster position - barrel layer " << i;
306  fHistContainerThisRun[key.str()] =
307  ibooker.book2D(name.str().c_str(), title.str().c_str(), 100, -30., 30., 64, -Geom::pi(), Geom::pi());
308  }
309 
310  // Pixel cluster positions in the endcaps (x, y).
311  for (size_t i = 0; i <= kNumDisks; ++i) {
312  std::stringstream key;
313  key << "clusPosEndCap" << i;
314  std::stringstream name;
315  name << key.str() << "_" << run.run();
316  std::stringstream title;
317  title << "Pixel cluster position - endcap disk " << i;
318  fHistContainerThisRun[key.str()] =
319  ibooker.book2D(name.str().c_str(), title.str().c_str(), 100, -20., 20., 100, -20., 20.);
320  }
321 
322  // Pixel cluster charge in the barrel, per bx.
323  for (size_t i = 0; i <= kNumLayers; ++i) {
324  std::stringstream key;
325  key << "clusChargeBarrel" << i;
326  std::stringstream name;
327  name << key.str() << "_" << run.run();
328  std::stringstream title;
329  title << "Pixel cluster charge - barrel layer " << i;
330  fHistContainerThisRun[key.str()] =
331  ibooker.book2D(name.str().c_str(), title.str().c_str(), 3564, .5, 3564.5, 100, 0., 100.);
332  }
333 
334  // Pixel cluster charge in the endcaps, per bx.
335  for (size_t i = 0; i <= kNumDisks; ++i) {
336  std::stringstream key;
337  key << "clusChargeEndCap" << i;
338  std::stringstream name;
339  name << key.str() << "_" << run.run();
340  std::stringstream title;
341  title << "Pixel cluster charge - endcap disk " << i;
342  fHistContainerThisRun[key.str()] =
343  ibooker.book2D(name.str().c_str(), title.str().c_str(), 3564, .5, 3564.5, 100, 0., 100.);
344  }
345 
346  // Pixel cluster size in the barrel, per bx.
347  for (size_t i = 0; i <= kNumLayers; ++i) {
348  std::stringstream key;
349  key << "clusSizeBarrel" << i;
350  std::stringstream name;
351  name << key.str() << "_" << run.run();
352  std::stringstream title;
353  title << "Pixel cluster size - barrel layer " << i;
354  fHistContainerThisRun[key.str()] =
355  ibooker.book2D(name.str().c_str(), title.str().c_str(), 3564, .5, 3564.5, 100, 0., 100.);
356  }
357 
358  // Pixel cluster size in the endcaps, per bx.
359  for (size_t i = 0; i <= kNumDisks; ++i) {
360  std::stringstream key;
361  key << "clusSizeEndCap" << i;
362  std::stringstream name;
363  name << key.str() << "_" << run.run();
364  std::stringstream title;
365  title << "Pixel cluster size - endcap disk " << i;
366  fHistContainerThisRun[key.str()] =
367  ibooker.book2D(name.str().c_str(), title.str().c_str(), 3564, .5, 3564.5, 100, 0., 100.);
368  }
369  }
370 }
371 
372 // ------------ Method called when starting to process a luminosity block.
373 // ------------
375  // Only reset and fill every fResetIntervalInLumiSections (default is 1 LS)
376  // Return unless the PREVIOUS LS was at the right modulo value
377  // (e.g. is resetinterval = 5 the rest will only be executed at LS=6
378  // NB: reset is done here so the histograms by LS are sent before resetting.
379  // NB: not being used for now since default is 1 LS. There is a bug here.
380 
381  unsigned int ls = lumiBlock.luminosityBlockAuxiliary().luminosityBlock();
382 
383  if ((ls - 1) % fResetIntervalInLumiSections == 0) {
387  }
388 }
389 
390 // ------------ Method called when ending the processing of a luminosity block.
391 // ------------
393  unsigned int ls = lumiBlock.luminosityBlockAuxiliary().luminosityBlock();
394 
395  // Only fill every fResetIntervalInLumiSections (default is 1 LS)
396  if (ls % fResetIntervalInLumiSections != 0)
397  return;
398 
399  printf("Lumi Block = %d\n", ls);
400 
401  if ((ls - 1) % fResetIntervalInLumiSections == 0) {
402  }
403 
404  unsigned int nBClus[3] = {0, 0, 0};
405  unsigned int nFPClus[2] = {0, 0};
406  unsigned int nFMClus[2] = {0, 0};
407 
408  double total_recorded = 0.;
409  double total_recorded_unc_square = 0.;
410 
411  // Obtain bunch-by-bunch cluster counts and compute totals for lumi
412  // calculation.
413  double totalcounts = 0.0;
414  double etotalcounts = 0.0;
415  double totalevents = 0.0;
416  double lumi_factor_per_bx = 0.0;
419  else
421 
422  for (std::map<int, PixelClusterCount>::iterator it = fNumPixelClusters.begin(); it != fNumPixelClusters.end(); it++) {
423  // Sum all clusters for this BX.
424  unsigned int total = (*it).second.numB.at(1) + (*it).second.numB.at(2) + (*it).second.numFP.at(0) +
425  (*it).second.numFP.at(1) + (*it).second.numFM.at(0) + (*it).second.numFM.at(1);
427  total += (*it).second.numB.at(0);
428  totalcounts += total;
429  double etotal = (*it).second.dnumB.at(1) + (*it).second.dnumB.at(2) + (*it).second.dnumFP.at(0) +
430  (*it).second.dnumFP.at(1) + (*it).second.dnumFM.at(0) + (*it).second.dnumFM.at(1);
432  etotal = (*it).second.dnumB.at(0);
433  etotalcounts += etotal;
434  etotal = sqrt(etotal);
435 
436  fHistClusterCountByBxLastLumi->setBinContent((*it).first, total);
437  fHistClusterCountByBxLastLumi->setBinError((*it).first, etotal);
440 
441  unsigned int events_per_bx = fHistBunchCrossingsLastLumi->getBinContent((*it).first);
442  totalevents += events_per_bx;
443  double average_cluster_count = events_per_bx != 0 ? double(total) / events_per_bx : 0.;
444  double average_cluster_count_unc = events_per_bx != 0 ? etotal / events_per_bx : 0.;
445  double pixel_bx_lumi_per_ls = lumi_factor_per_bx * average_cluster_count / CM2_TO_NANOBARN;
446  double pixel_bx_lumi_per_ls_unc = 0.0;
448  pixel_bx_lumi_per_ls_unc = sqrt(lumi_factor_per_bx * lumi_factor_per_bx *
449  (average_cluster_count_unc * average_cluster_count_unc +
450  (average_cluster_count * XSEC_PIXEL_CLUSTER_UNC / XSEC_PIXEL_CLUSTER) *
451  (average_cluster_count * XSEC_PIXEL_CLUSTER / XSEC_PIXEL_CLUSTER))) /
453  else
454  pixel_bx_lumi_per_ls_unc = sqrt(lumi_factor_per_bx * lumi_factor_per_bx *
455  (average_cluster_count_unc * average_cluster_count_unc +
456  (average_cluster_count * rXSEC_PIXEL_CLUSTER_UNC / rXSEC_PIXEL_CLUSTER) *
457  (average_cluster_count * rXSEC_PIXEL_CLUSTER / rXSEC_PIXEL_CLUSTER))) /
459 
460  fHistRecordedByBxLastLumi->setBinContent((*it).first, pixel_bx_lumi_per_ls);
461  fHistRecordedByBxLastLumi->setBinError((*it).first, pixel_bx_lumi_per_ls_unc);
462 
464  (*it).first, fHistRecordedByBxCumulative->getBinContent((*it).first) + pixel_bx_lumi_per_ls);
465 
466  /*
467  if(fHistRecordedByBxLastLumi->getBinContent((*it).first)!=0.)
468  fHistRecordedByBxLastLumi->getBinContent((*it).first));
469  if(fHistRecordedByBxCumulative->getBinContent((*it).first)!=0.)
470  fHistRecordedByBxCumulative->getBinContent((*it).first));
471  */
472 
473  nBClus[0] += (*it).second.numB.at(0);
474  nBClus[1] += (*it).second.numB.at(1);
475  nBClus[2] += (*it).second.numB.at(2);
476  nFPClus[0] += (*it).second.numFP.at(0);
477  nFPClus[1] += (*it).second.numFP.at(1);
478  nFMClus[0] += (*it).second.numFM.at(0);
479  nFMClus[1] += (*it).second.numFM.at(1);
480 
481  // Reset counters
482  (*it).second.Reset();
483 
484  // edm::LogWarning("pixelLumi") << "bx="<< (*it).first << " clusters=" <<
485  // (*it).second.numB.at(0));
486  }
487 
489  for (unsigned int i = 0; i <= lastBunchCrossing; i++) {
490  if (bunchTriggerMask[i]) {
492  total_recorded += fHistRecordedByBxLastLumi->getBinContent(i);
493  total_recorded_unc_square += err * err;
494  }
495  }
496 
497  // Replace the total obtained by summing over BXs with the average per BX
498  // from the total cluster count and rescale
499 
500  if (totalevents > 10) {
501  total_recorded = lumi_factor_per_bx * totalcounts / totalevents / CM2_TO_NANOBARN;
502  } else
503  total_recorded = 0.0;
504 
505  edm::LogWarning("pixelLumi") << " Total recorded " << total_recorded;
506  fHistTotalRecordedLumiByLS->setBinContent(ls, total_recorded);
507  fHistTotalRecordedLumiByLS->setBinError(ls, sqrt(total_recorded_unc_square));
508  }
509  // fill cluster counts by detector regions for sanity checks
510  unsigned int all_detectors_counts = 0;
511  for (unsigned int i = 0; i < 3; i++) {
512  all_detectors_counts += nBClus[i];
513  }
514  for (unsigned int i = 0; i < 2; i++) {
515  all_detectors_counts += nFPClus[i];
516  }
517  for (unsigned int i = 0; i < 2; i++) {
518  all_detectors_counts += nFMClus[i];
519  }
520 
521  fHistClusByLS->setBinContent(ls, all_detectors_counts);
522 
523  for (unsigned int i = 0; i < 3; i++) {
524  fHistnBClusVsLS[i]->setBinContent(ls, float(nBClus[i]) / float(all_detectors_counts));
525  }
526  for (unsigned int i = 0; i < 2; i++) {
527  fHistnFPClusVsLS[i]->setBinContent(ls, float(nFPClus[i]) / float(all_detectors_counts));
528  }
529  for (unsigned int i = 0; i < 2; i++) {
530  fHistnFMClusVsLS[i]->setBinContent(ls, float(nFMClus[i]) / float(all_detectors_counts));
531  }
532 
534 
535  timeval tv;
536  gettimeofday(&tv, nullptr);
537  tm *ts = gmtime(&tv.tv_sec);
538  char datestring[256];
539  strftime(datestring, sizeof(datestring), "%Y.%m.%d %T GMT %s", ts);
540  logFile_ << "RunNumber " << fRunNo << std::endl;
541  logFile_ << "EndTimeOfFit " << datestring << std::endl;
542  logFile_ << "LumiRange " << ls << "-" << ls << std::endl;
543  logFile_ << "Fill " << -99 << std::endl;
544  logFile_ << "ActiveBunchCrossings " << filledAndUnmaskedBunches << std::endl;
545  logFile_ << "PixelLumi " << fHistTotalRecordedLumiByLS->getBinContent(ls) * 0.98 << std::endl;
546  logFile_ << "HFLumi " << -99 << std::endl;
547  logFile_ << "Ratio " << -99 << std::endl;
548  logFile_.close();
549 }
550 
551 unsigned int PixelLumiDQM::calculateBunchMask(MonitorElement *e, std::vector<bool> &mask) {
552  unsigned int nbins = e->getNbinsX();
553  std::vector<float> ar(nbins + 1, 0.);
554  for (unsigned int i = 1; i <= nbins; i++) {
555  ar[i] = e->getBinContent(i);
556  }
557  return calculateBunchMask(ar, nbins, mask);
558 }
559 unsigned int PixelLumiDQM::calculateBunchMask(std::vector<float> &e, unsigned int nbins, std::vector<bool> &mask) {
560  // Take the cumulative cluster count histogram and find max and average of
561  // non-empty bins.
562  unsigned int active_count = 0;
563  double maxc = 0.0;
564  double ave = 0.0; // Average of non-empty bins
565  unsigned int non_empty_bins = 0;
566 
567  for (unsigned int i = 1; i <= nbins; i++) {
568  double bin = e[i];
569  if (bin != 0.0) {
570  if (maxc < bin)
571  maxc = bin;
572  ave += bin;
573  non_empty_bins++;
574  }
575  }
576 
577  ave /= non_empty_bins;
578  edm::LogWarning("pixelLumi") << "Bunch mask finder - non empty bins " << non_empty_bins
579  << " average of non empty bins " << ave << " max content of one bin " << maxc;
580  double mean = 0.;
581  double sigma = 0.;
582  if (non_empty_bins < 50) {
583  mean = maxc;
584  sigma = (maxc) / 20;
585  } else {
586  TH1F dist("dist", "dist", 500, 0., maxc + (maxc / 500.) * 20.);
587  for (unsigned int i = 1; i <= nbins; i++) {
588  double bin = e[i];
589  dist.Fill(bin);
590  }
591  TF1 fit("plgaus", "gaus");
592  dist.Fit(&fit, "", "", fmax(0., ave - (maxc - ave) / 5.), maxc);
593  mean = fit.GetParameter("Mean");
594  sigma = fit.GetParameter("Sigma");
595  }
596  edm::LogWarning("pixelLumi") << "Bunch mask will use mean" << mean << " sigma " << sigma;
597  // Active BX defined as those which have nclus within fixed standard
598  // deviations of peak.
599  for (unsigned int i = 1; i <= nbins; i++) {
600  double bin = e[i];
601  if (bin > 0. && std::abs(bin - mean) < 5. * (sigma)) {
602  mask[i] = true;
603  active_count++;
604  }
605  }
606  edm::LogWarning("pixelLumi") << "Bunch mask finds " << active_count << " active bunch crossings ";
607  // edm::LogWarning("pixelLumi") << "this is the full bx mask " ;
608  // for(unsigned int i = 1; i<= nbins; i++)
609  // edm::LogWarning("pixelLumi") << ((mask[i]) ? 1:0);
610  return active_count;
611 }
612 // Define this as a CMSSW plug-in.
RunNumber_t run() const
Definition: EventID.h:38
bool fIncludePixelQualCheckHistos
Definition: PixelLumiDQM.h:135
static constexpr double XSEC_PIXEL_CLUSTER
Definition: PixelLumiDQM.h:49
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:79
static constexpr double XSEC_PIXEL_CLUSTER_UNC
Definition: PixelLumiDQM.h:50
EventNumber_t event() const
Definition: EventID.h:40
T getUntrackedParameter(std::string const &, T const &) const
unsigned int filledAndUnmaskedBunches
Definition: PixelLumiDQM.h:168
const_iterator end(bool update=false) const
RunID const & id() const
Definition: RunBase.h:39
RunNumber_t run() const
Definition: RunBase.h:40
virtual LocalPoint localPosition(const MeasurementPoint &) const =0
static const unsigned int lastBunchCrossing
Definition: PixelLumiDQM.h:56
RunNumber_t run() const
Definition: RunID.h:36
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:32
std::vector< uint32_t > fDeadModules
Definition: PixelLumiDQM.h:142
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:539
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
def ls
Definition: eostools.py:349
virtual double getBinError(int binx) const
get uncertainty on content of bin (1-D) - See TH1::GetBinError for details
T y() const
Definition: PV3DBase.h:60
int bunchCrossing() const
Definition: EventBase.h:64
static constexpr size_t kOffsetDisks
Definition: PixelLumiDQM.h:78
std::string to_string(const V &value)
Definition: OMSAccess.h:71
edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > tkGeomToken_
Definition: PixelLumiDQM.h:116
const DetContainer & dets() const override
Returm a vector of all GeomDet (including all GeomDetUnits)
int fResetIntervalInLumiSections
Definition: PixelLumiDQM.h:136
data_type const * const_iterator
Definition: DetSetNew.h:31
MonitorElement * fHistRecordedByBxCumulative
Definition: PixelLumiDQM.h:165
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
void beginLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &) override
assert(be >=bs)
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
bool useInnerBarrelLayer
Definition: PixelLumiDQM.h:169
constexpr std::array< uint8_t, layerIndexSize > layer
LuminosityBlockNumber_t luminosityBlock() const
void Fill(long long x)
MonitorElement * fHistClusByLS
Definition: PixelLumiDQM.h:162
bool getData(T &iHolder) const
Definition: EventSetup.h:122
virtual void Reset()
Remove all data from the ME, keept the empty histogram with all its settings.
~PixelLumiDQM() override
Definition: PixelLumiDQM.cc:81
LuminosityBlockAuxiliary const & luminosityBlockAuxiliary() const override
static constexpr size_t kNumDisks
Definition: PixelLumiDQM.h:76
Measurement2DPoint MeasurementPoint
Measurement points are two-dimensional by default.
static constexpr double rXSEC_PIXEL_CLUSTER_UNC
Definition: PixelLumiDQM.h:54
int iEvent
Definition: GenABIO.cc:224
void addDefault(ParameterSetDescription const &psetDescription)
virtual int getNbinsX() const
get # of bins in X-axis
MonitorElement * fHistnFMClusVsLS[2]
Definition: PixelLumiDQM.h:157
static constexpr double SECONDS_PER_LS
Definition: PixelLumiDQM.h:46
edm::EDGetTokenT< edmNew::DetSetVector< SiPixelCluster > > fPixelClusterLabel
Definition: PixelLumiDQM.h:115
unsigned int calculateBunchMask(MonitorElement *, std::vector< bool > &)
auto const & fit
T sqrt(T t)
Definition: SSEVec.h:19
printf("params %d %f %f %f\n", minT, eps, errmax, chi2max)
T z() const
Definition: PV3DBase.h:61
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
tuple key
prepare the HTCondor submission files and eventually submit them
MonitorElement * fHistClusterCountByBxLastLumi
Definition: PixelLumiDQM.h:160
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
MonitorElement * fHistTotalRecordedLumiByLS
Definition: PixelLumiDQM.h:163
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
Definition: PixelLumiDQM.cc:83
unsigned int unixTime() const
Time in seconds since January 1, 1970.
Definition: Timestamp.h:40
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
LuminosityBlock const & getLuminosityBlock() const
Definition: Event.h:100
MonitorElement * fHistClusterCountByBxCumulative
Definition: PixelLumiDQM.h:161
void endLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &) override
MonitorElement * fHistnBClusVsLS[3]
Definition: PixelLumiDQM.h:155
virtual double getBinContent(int binx) const
get content of bin (1-D)
LuminosityBlockNumber_t luminosityBlock() const
const TrackerGeomDet * idToDet(DetId) const override
std::map< std::string, MonitorElement * > fHistContainerThisRun
Definition: PixelLumiDQM.h:138
double fMinClusterCharge
Definition: PixelLumiDQM.h:150
Log< level::Info, false > LogInfo
boost::transform_iterator< IterHelp, const_IdIter > const_iterator
std::vector< bool > bunchTriggerMask
Definition: PixelLumiDQM.h:167
Definition: DetId.h:17
bool fIncludePixelClusterInfo
Definition: PixelLumiDQM.h:134
static constexpr double rXSEC_PIXEL_CLUSTER
Definition: PixelLumiDQM.h:53
virtual void setBinContent(int binx, double content)
set content of bin (1-D)
MonitorElement * fHistBunchCrossings
Definition: PixelLumiDQM.h:158
MonitorElement * fHistBunchCrossingsLastLumi
Definition: PixelLumiDQM.h:159
MonitorElement * book2D(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:212
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
UInt_t fTimestamp
Definition: PixelLumiDQM.h:122
std::string fLogFileName_
Definition: PixelLumiDQM.h:171
PixelLumiDQM(const edm::ParameterSet &)
Definition: PixelLumiDQM.cc:47
static constexpr size_t kOffsetLayers
Definition: PixelLumiDQM.h:77
void analyze(const edm::Event &, const edm::EventSetup &) override
Definition: PixelLumiDQM.cc:89
edm::EventID id() const
Definition: EventBase.h:59
Pixel cluster – collection of neighboring pixels above threshold.
std::map< int, PixelClusterCount > fNumPixelClusters
Definition: PixelLumiDQM.h:131
MonitorElement * fHistRecordedByBxLastLumi
Definition: PixelLumiDQM.h:164
int fMinPixelsPerCluster
Definition: PixelLumiDQM.h:147
virtual void setBinError(int binx, double error)
set uncertainty on content of bin (1-D)
Log< level::Warning, false > LogWarning
constexpr double pi()
Definition: Pi.h:31
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
T x() const
Definition: PV3DBase.h:59
MonitorElement * fHistnFPClusVsLS[2]
Definition: PixelLumiDQM.h:156
ESGetTokenH3DDVariant esConsumes(std::string const &Reccord, edm::ConsumesCollector &)
Definition: DeDxTools.cc:283
int layerName() const
layer id
tuple size
Write out results.
edm::Timestamp time() const
Definition: EventBase.h:60
bool isTrackerPixel(GeomDetEnumerators::SubDetector m)
const_iterator begin(bool update=false) const
Definition: Run.h:45
static constexpr double FREQ_ORBIT
Definition: PixelLumiDQM.h:45
static constexpr double CM2_TO_NANOBARN
Definition: PixelLumiDQM.h:55
std::ofstream logFile_
Definition: PixelLumiDQM.h:173
static constexpr size_t kNumLayers
Definition: PixelLumiDQM.h:75