CMS 3D CMS Logo

SiStripGainsCalibTreeWorker.cc
Go to the documentation of this file.
28 
29 #include <TChain.h>
30 
34 
35 //
39 //
40 // Original Author: L. Quertermont (calibration algorithm)
41 // Contributors: M. Verzetti (data access)
42 // A. Di Mattia (PCL multi stream processing and monitoring)
43 // M. Delcourt (monitoring)
44 // M. Musich (migration to thread-safe DQMStore access)
45 // P. David (adapted to calibtrees)
46 //
47 // Created: Wed, 12 Apr 2017 14:46:48 GMT
48 //
49 
50 class SiStripGainsCalibTreeWorker : public DQMGlobalEDAnalyzer<APVGain::APVGainHistograms> {
51 public:
53 
55  edm::Run const&,
56  edm::EventSetup const&,
57  APVGain::APVGainHistograms&) const override;
58  void dqmAnalyze(edm::Event const&, edm::EventSetup const&, APVGain::APVGainHistograms const&) const override;
59 
60  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
61 
62 private:
63  void dqmBeginRun(edm::Run const&, edm::EventSetup const&, APVGain::APVGainHistograms&) const override;
65 
66  std::vector<std::string> dqm_tag_;
67 
68  int statCollectionFromMode(const char* tag) const;
69 
72  double MinTrackEta;
73  double MaxTrackEta;
74  unsigned int MaxNrStrips;
75  unsigned int MinTrackHits;
80  bool Validation;
85  mutable bool hasProcessed_;
86 
89  std::vector<std::string> VChargeHisto;
91  // maps histograms index to topology
92  std::map<unsigned int, APVloc> theTopologyMap;
93 
96  std::string TrackPrefix_; //("track");
98  std::string CalibPrefix_; //("GainCalibration");
101  std::vector<std::string> calibTreeFileNames_;
102 
107 };
108 
110  std::vector<std::string>::const_iterator it = dqm_tag_.begin();
111  while (it != dqm_tag_.end()) {
112  if (*it == std::string(tag))
113  return it - dqm_tag_.begin();
114  it++;
115  }
116 
117  if (std::string(tag).empty())
118  return 0; // return StdBunch calibration mode for backward compatibility
119 
120  return None;
121 }
122 
125 
126 //********************************************************************************//
128  MinTrackMomentum = iConfig.getUntrackedParameter<double>("minTrackMomentum", 3.0);
129  MaxTrackMomentum = iConfig.getUntrackedParameter<double>("maxTrackMomentum", 99999.0);
130  MinTrackEta = iConfig.getUntrackedParameter<double>("minTrackEta", -5.0);
131  MaxTrackEta = iConfig.getUntrackedParameter<double>("maxTrackEta", 5.0);
132  MaxNrStrips = iConfig.getUntrackedParameter<unsigned>("maxNrStrips", 2);
133  MinTrackHits = iConfig.getUntrackedParameter<unsigned>("MinTrackHits", 8);
134  MaxTrackChiOverNdf = iConfig.getUntrackedParameter<double>("MaxTrackChiOverNdf", 3);
135  MaxTrackingIteration = iConfig.getUntrackedParameter<int>("MaxTrackingIteration", 7);
136  AllowSaturation = iConfig.getUntrackedParameter<bool>("AllowSaturation", false);
137  FirstSetOfConstants = iConfig.getUntrackedParameter<bool>("FirstSetOfConstants", true);
138  Validation = iConfig.getUntrackedParameter<bool>("Validation", false);
139  OldGainRemoving = iConfig.getUntrackedParameter<bool>("OldGainRemoving", false);
140  useCalibration = iConfig.getUntrackedParameter<bool>("UseCalibration", false);
141  doChargeMonitorPerPlane = iConfig.getUntrackedParameter<bool>("doChargeMonitorPerPlane", false);
142  m_DQMdir = iConfig.getUntrackedParameter<std::string>("DQMdir", "AlCaReco/SiStripGains");
143  m_calibrationMode = iConfig.getUntrackedParameter<std::string>("calibrationMode", "StdBunch");
144  VChargeHisto = iConfig.getUntrackedParameter<std::vector<std::string>>("ChargeHisto");
145 
146  // fill in the mapping between the histogram indices and the (id,side,plane) tuple
147  std::vector<std::pair<std::string, std::string>> hnames =
149  for (unsigned int i = 0; i < hnames.size(); i++) {
150  int id = APVGain::subdetectorId((hnames[i]).first);
151  int side = APVGain::subdetectorSide((hnames[i]).first);
152  int plane = APVGain::subdetectorPlane((hnames[i]).first);
153  int thick = APVGain::thickness((hnames[i]).first);
154  std::string s = hnames[i].first;
155 
156  auto loc = APVloc(thick, id, side, plane, s);
157  theTopologyMap.insert(std::make_pair(i, loc));
158  }
159 
160  //Set the monitoring element tag and store
161  dqm_tag_.reserve(7);
162  dqm_tag_.clear();
163  dqm_tag_.push_back("StdBunch"); // statistic collection from Standard Collision Bunch @ 3.8 T
164  dqm_tag_.push_back("StdBunch0T"); // statistic collection from Standard Collision Bunch @ 0 T
165  dqm_tag_.push_back("AagBunch"); // statistic collection from First Collision After Abort Gap @ 3.8 T
166  dqm_tag_.push_back("AagBunch0T"); // statistic collection from First Collision After Abort Gap @ 0 T
167  dqm_tag_.push_back("IsoMuon"); // statistic collection from Isolated Muon @ 3.8 T
168  dqm_tag_.push_back("IsoMuon0T"); // statistic collection from Isolated Muon @ 0 T
169  dqm_tag_.push_back("Harvest"); // statistic collection: Harvest
170 
171  hasProcessed_ = false;
172 
173  edm::ParameterSet swhallowgain_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("gain");
174  CalibPrefix_ = swhallowgain_pset.getUntrackedParameter<std::string>("prefix");
175  CalibSuffix_ = swhallowgain_pset.getUntrackedParameter<std::string>("suffix");
176  edm::ParameterSet evtinfo_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("evtinfo");
177  EventPrefix_ = evtinfo_pset.getUntrackedParameter<std::string>("prefix");
178  EventSuffix_ = evtinfo_pset.getUntrackedParameter<std::string>("suffix");
179  edm::ParameterSet track_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("tracks");
180  TrackPrefix_ = track_pset.getUntrackedParameter<std::string>("prefix");
181  TrackSuffix_ = track_pset.getUntrackedParameter<std::string>("suffix");
182 
183  calibTreeName_ = iConfig.getUntrackedParameter<std::string>("CalibTreeName");
184  calibTreeFileNames_ = iConfig.getUntrackedParameter<std::vector<std::string>>("CalibTreeFiles");
185 
187  tkGeomTokenBR_ = esConsumes<edm::Transition::BeginRun>();
188  gainToken_ = esConsumes<edm::Transition::BeginRun>();
189  qualityToken_ = esConsumes<edm::Transition::BeginRun>();
190 }
191 
192 //********************************************************************************//
194  edm::EventSetup const& iSetup,
196  using namespace edm;
197 
198  // fills the APV collections at each begin run
199  const TrackerGeometry* bareTkGeomPtr = &iSetup.getData(tkGeomTokenBR_);
200  checkBookAPVColls(bareTkGeomPtr, histograms);
201 
202  const auto gainHandle = iSetup.getHandle(gainToken_);
203  if (!gainHandle.isValid()) {
204  edm::LogError("SiStripGainPCLWorker") << "gainHandle is not valid\n";
205  exit(0);
206  }
207 
208  const auto& SiStripQuality_ = &iSetup.getData(qualityToken_);
209 
210  for (unsigned int a = 0; a < histograms.APVsCollOrdered.size(); a++) {
211  std::shared_ptr<stAPVGain> APV = histograms.APVsCollOrdered[a];
212 
214  continue;
215 
216  APV->isMasked = SiStripQuality_->IsApvBad(APV->DetId, APV->APVId);
217 
218  if (gainHandle->getNumberOfTags() != 2) {
219  edm::LogError("SiStripGainPCLWorker") << "NUMBER OF GAIN TAG IS EXPECTED TO BE 2\n";
220  fflush(stdout);
221  exit(0);
222  };
223  float newPreviousGain = gainHandle->getApvGain(APV->APVId, gainHandle->getRange(APV->DetId, 1), 1);
224  if (APV->PreviousGain != 1 and newPreviousGain != APV->PreviousGain)
225  edm::LogWarning("SiStripGainPCLWorker") << "WARNING: ParticleGain in the global tag changed\n";
226  APV->PreviousGain = newPreviousGain;
227 
228  float newPreviousGainTick = gainHandle->getApvGain(APV->APVId, gainHandle->getRange(APV->DetId, 0), 0);
229  if (APV->PreviousGainTick != 1 and newPreviousGainTick != APV->PreviousGainTick) {
230  edm::LogWarning("SiStripGainPCLWorker")
231  << "WARNING: TickMarkGain in the global tag changed\n"
232  << std::endl
233  << " APV->SubDet: " << APV->SubDet << " APV->APVId:" << APV->APVId << std::endl
234  << " APV->PreviousGainTick: " << APV->PreviousGainTick << " newPreviousGainTick: " << newPreviousGainTick
235  << std::endl;
236  }
237  APV->PreviousGainTick = newPreviousGainTick;
238  }
239 }
240 
241 //********************************************************************************//
242 // ------------ method called for each event ------------
244  edm::EventSetup const& iSetup,
245  APVGain::APVGainHistograms const& histograms) const {
246  if (!hasProcessed_) {
247  const TrackerTopology* topo = &iSetup.getData(tTopoToken_);
248 
249  for (const auto& elem : theTopologyMap) {
250  LogDebug("SiStripGainsCalibTreeWorker")
251  << elem.first << " - " << elem.second.m_string << " " << elem.second.m_subdetectorId << " "
252  << elem.second.m_subdetectorSide << " " << elem.second.m_subdetectorPlane << std::endl;
253  }
254 
255  LogDebug("SiStripGainsCalibTreeWorker") << "for mode" << m_calibrationMode << std::endl;
256 
257  // data members
258  unsigned int eventnumber = 0;
259  unsigned int runnumber = 0;
260 #ifdef ExtendedCALIBTree
261  const std::vector<bool>* TrigTech = nullptr;
262  const std::vector<double>* chargeoverpath = nullptr;
263 #endif
264  // Track data
265  const std::vector<double>* trackchi2ndof = nullptr;
266  const std::vector<float>* trackp = nullptr;
267  const std::vector<float>* trackpt = nullptr;
268  const std::vector<double>* tracketa = nullptr;
269  const std::vector<double>* trackphi = nullptr;
270  const std::vector<unsigned int>* trackhitsvalid = nullptr;
271  const std::vector<int>* trackalgo = nullptr;
272  // CalibTree data
273  const std::vector<int>* trackindex = nullptr;
274  const std::vector<unsigned int>* rawid = nullptr;
275  const std::vector<double>* localdirx = nullptr;
276  const std::vector<double>* localdiry = nullptr;
277  const std::vector<double>* localdirz = nullptr;
278  const std::vector<unsigned short>* firststrip = nullptr;
279  const std::vector<unsigned short>* nstrips = nullptr;
280  const std::vector<bool>* saturation = nullptr;
281  const std::vector<bool>* overlapping = nullptr;
282  const std::vector<bool>* farfromedge = nullptr;
283  const std::vector<unsigned int>* charge = nullptr;
284  const std::vector<double>* path = nullptr;
285  const std::vector<unsigned char>* amplitude = nullptr;
286  const std::vector<double>* gainused = nullptr;
287  const std::vector<double>* gainusedTick = nullptr;
288 
289  TChain tree{calibTreeName_.c_str()};
290  for (const auto& ctFn : calibTreeFileNames_) {
291  tree.Add(ctFn.c_str());
292  }
293 
294  tree.SetBranchAddress((EventPrefix_ + "event" + EventSuffix_).c_str(), &eventnumber, nullptr);
295  tree.SetBranchAddress((EventPrefix_ + "run" + EventSuffix_).c_str(), &runnumber, nullptr);
296 #ifdef ExtendedCALIBTree
297  tree.SetBranchAddress((EventPrefix_ + "TrigTech" + EventSuffix_).c_str(), &TrigTech, nullptr);
298  tree.SetBranchAddress((CalibPrefix_ + "chargeoverpath" + CalibSuffix_).c_str(), &chargeoverpath, nullptr);
299 #endif
300  tree.SetBranchAddress((TrackPrefix_ + "chi2ndof" + TrackSuffix_).c_str(), &trackchi2ndof, nullptr);
301  tree.SetBranchAddress((TrackPrefix_ + "momentum" + TrackSuffix_).c_str(), &trackp, nullptr);
302  tree.SetBranchAddress((TrackPrefix_ + "pt" + TrackSuffix_).c_str(), &trackpt, nullptr);
303  tree.SetBranchAddress((TrackPrefix_ + "eta" + TrackSuffix_).c_str(), &tracketa, nullptr);
304  tree.SetBranchAddress((TrackPrefix_ + "phi" + TrackSuffix_).c_str(), &trackphi, nullptr);
305  tree.SetBranchAddress((TrackPrefix_ + "hitsvalid" + TrackSuffix_).c_str(), &trackhitsvalid, nullptr);
306  tree.SetBranchAddress((TrackPrefix_ + "algo" + TrackSuffix_).c_str(), &trackalgo, nullptr);
307 
308  tree.SetBranchAddress((CalibPrefix_ + "trackindex" + CalibSuffix_).c_str(), &trackindex, nullptr);
309  tree.SetBranchAddress((CalibPrefix_ + "rawid" + CalibSuffix_).c_str(), &rawid, nullptr);
310  tree.SetBranchAddress((CalibPrefix_ + "localdirx" + CalibSuffix_).c_str(), &localdirx, nullptr);
311  tree.SetBranchAddress((CalibPrefix_ + "localdiry" + CalibSuffix_).c_str(), &localdiry, nullptr);
312  tree.SetBranchAddress((CalibPrefix_ + "localdirz" + CalibSuffix_).c_str(), &localdirz, nullptr);
313  tree.SetBranchAddress((CalibPrefix_ + "firststrip" + CalibSuffix_).c_str(), &firststrip, nullptr);
314  tree.SetBranchAddress((CalibPrefix_ + "nstrips" + CalibSuffix_).c_str(), &nstrips, nullptr);
315  tree.SetBranchAddress((CalibPrefix_ + "saturation" + CalibSuffix_).c_str(), &saturation, nullptr);
316  tree.SetBranchAddress((CalibPrefix_ + "overlapping" + CalibSuffix_).c_str(), &overlapping, nullptr);
317  tree.SetBranchAddress((CalibPrefix_ + "farfromedge" + CalibSuffix_).c_str(), &farfromedge, nullptr);
318  tree.SetBranchAddress((CalibPrefix_ + "charge" + CalibSuffix_).c_str(), &charge, nullptr);
319  tree.SetBranchAddress((CalibPrefix_ + "path" + CalibSuffix_).c_str(), &path, nullptr);
320 
321  tree.SetBranchAddress((CalibPrefix_ + "amplitude" + CalibSuffix_).c_str(), &amplitude, nullptr);
322  tree.SetBranchAddress((CalibPrefix_ + "gainused" + CalibSuffix_).c_str(), &gainused, nullptr);
323  tree.SetBranchAddress((CalibPrefix_ + "gainusedTick" + CalibSuffix_).c_str(), &gainusedTick, nullptr);
324 
325  int elepos = statCollectionFromMode(m_calibrationMode.c_str());
326  const auto nEntries = tree.GetEntries();
327  printf("Progressing Bar :0%% 20%% 40%% 60%% 80%% 100%%\n");
328  printf("Looping on the Tree :");
329  const auto treeProgStep = std::max(nEntries / 50, Long64_t(1));
330  for (Long64_t i{0}; i != nEntries; ++i) {
331  if ((i % treeProgStep) == 0) {
332  printf(".");
333  fflush(stdout);
334  }
335  tree.GetEntry(i);
336 
337  unsigned int FirstAmplitude = 0;
338  for (unsigned int i = 0; i < charge->size(); i++) {
339  FirstAmplitude += (*nstrips)[i];
340  int TI = (*trackindex)[i];
341 
342  if ((*tracketa)[TI] < MinTrackEta)
343  continue;
344  if ((*tracketa)[TI] > MaxTrackEta)
345  continue;
346  if ((*trackp)[TI] < MinTrackMomentum)
347  continue;
348  if ((*trackp)[TI] > MaxTrackMomentum)
349  continue;
350  if ((*trackhitsvalid)[TI] < MinTrackHits)
351  continue;
352  if ((*trackchi2ndof)[TI] > MaxTrackChiOverNdf)
353  continue;
354  if ((*trackalgo)[TI] > MaxTrackingIteration)
355  continue;
356 
357  std::shared_ptr<stAPVGain> APV = histograms.APVsColl.at(
358  ((*rawid)[i] << 4) |
359  ((*firststrip)[i] /
360  128)); //works for both strip and pixel thanks to firstStrip encoding for pixel in the calibTree
361 
362  if (APV->SubDet > 2 && (*farfromedge)[i] == false)
363  continue;
364  if (APV->SubDet > 2 && (*overlapping)[i] == true)
365  continue;
366  if (APV->SubDet > 2 && (*saturation)[i] && !AllowSaturation)
367  continue;
368  if (APV->SubDet > 2 && (*nstrips)[i] > MaxNrStrips)
369  continue;
370 
371  int Charge = 0;
372  if (APV->SubDet > 2 && (useCalibration || !FirstSetOfConstants)) {
373  bool Saturation = false;
374  for (unsigned int s = 0; s < (*nstrips)[i]; s++) {
375  int StripCharge = (*amplitude)[FirstAmplitude - (*nstrips)[i] + s];
377  StripCharge = (int)(StripCharge * (APV->PreviousGain / APV->CalibGain));
378  } else if (useCalibration) {
379  StripCharge = (int)(StripCharge / APV->CalibGain);
380  } else if (!FirstSetOfConstants) {
381  StripCharge = (int)(StripCharge * APV->PreviousGain);
382  }
383  if (StripCharge > 1024) {
384  StripCharge = 255;
385  Saturation = true;
386  } else if (StripCharge > 254) {
387  StripCharge = 254;
388  Saturation = true;
389  }
390  Charge += StripCharge;
391  }
392  if (Saturation && !AllowSaturation)
393  continue;
394  } else if (APV->SubDet > 2) {
395  Charge = (*charge)[i];
396  } else {
397  Charge = (*charge)[i] / 265.0; //expected scale factor between pixel and strip charge
398  }
399 
400  double ClusterChargeOverPath = ((double)Charge) / (*path)[i];
401  if (APV->SubDet > 2) {
402  if (Validation) {
403  ClusterChargeOverPath /= (*gainused)[i];
404  }
405  if (OldGainRemoving) {
406  ClusterChargeOverPath *= (*gainused)[i];
407  }
408  }
409 
410  // keep processing of pixel cluster charge until here
411  if (APV->SubDet <= 2)
412  continue;
413 
414  // real histogram for calibration
415  histograms.Charge_Vs_Index[elepos]->Fill(APV->Index, ClusterChargeOverPath);
416  LogDebug("SiStripGainsCalibTreeWorker")
417  << " for mode " << m_calibrationMode << "\n"
418  << " i " << i << " useCalibration " << useCalibration << " FirstSetOfConstants " << FirstSetOfConstants
419  << " APV->PreviousGain " << APV->PreviousGain << " APV->CalibGain " << APV->CalibGain << " APV->DetId "
420  << APV->DetId << " APV->Index " << APV->Index << " Charge " << Charge << " Path " << (*path)[i]
421  << " ClusterChargeOverPath " << ClusterChargeOverPath << std::endl;
422 
423  // Fill monitoring histograms
424  int mCharge1 = 0;
425  int mCharge2 = 0;
426  int mCharge3 = 0;
427  int mCharge4 = 0;
428  if (APV->SubDet > 2) {
429  for (unsigned int s = 0; s < (*nstrips)[i]; s++) {
430  int StripCharge = (*amplitude)[FirstAmplitude - (*nstrips)[i] + s];
431  if (StripCharge > 1024)
432  StripCharge = 255;
433  else if (StripCharge > 254)
434  StripCharge = 254;
435  mCharge1 += StripCharge;
436  mCharge2 += StripCharge;
437  mCharge3 += StripCharge;
438  mCharge4 += StripCharge;
439  }
440  // Revome gains for monitoring
441  mCharge2 *= (*gainused)[i]; // remove G2
442  mCharge3 *= (*gainusedTick)[i]; // remove G1
443  mCharge4 *= ((*gainused)[i] * (*gainusedTick)[i]); // remove G1 and G2
444  }
445 
446  LogDebug("SiStripGainsCalibTreeWorker")
447  << " full charge " << mCharge1 << " remove G2 " << mCharge2 << " remove G1 " << mCharge3 << " remove G1*G2 "
448  << mCharge4 << std::endl;
449 
450  auto indices = APVGain::FetchIndices(theTopologyMap, (*rawid)[i], topo);
451 
452  for (auto m : indices)
453  histograms.Charge_1[elepos][m]->Fill(((double)mCharge1) / (*path)[i]);
454  for (auto m : indices)
455  histograms.Charge_2[elepos][m]->Fill(((double)mCharge2) / (*path)[i]);
456  for (auto m : indices)
457  histograms.Charge_3[elepos][m]->Fill(((double)mCharge3) / (*path)[i]);
458  for (auto m : indices)
459  histograms.Charge_4[elepos][m]->Fill(((double)mCharge4) / (*path)[i]);
460 
461  if (APV->SubDet == StripSubdetector::TIB) {
462  histograms.Charge_Vs_PathlengthTIB[elepos]->Fill((*path)[i], Charge); // TIB
463 
464  } else if (APV->SubDet == StripSubdetector::TOB) {
465  histograms.Charge_Vs_PathlengthTOB[elepos]->Fill((*path)[i], Charge); // TOB
466 
467  } else if (APV->SubDet == StripSubdetector::TID) {
468  if (APV->Eta < 0) {
469  histograms.Charge_Vs_PathlengthTIDM[elepos]->Fill((*path)[i], Charge);
470  } // TID minus
471  else if (APV->Eta > 0) {
472  histograms.Charge_Vs_PathlengthTIDP[elepos]->Fill((*path)[i], Charge);
473  } // TID plus
474 
475  } else if (APV->SubDet == StripSubdetector::TEC) {
476  if (APV->Eta < 0) {
477  if (APV->Thickness < 0.04) {
478  histograms.Charge_Vs_PathlengthTECM1[elepos]->Fill((*path)[i], Charge);
479  } // TEC minus, type 1
480  else if (APV->Thickness > 0.04) {
481  histograms.Charge_Vs_PathlengthTECM2[elepos]->Fill((*path)[i], Charge);
482  } // TEC minus, type 2
483  } else if (APV->Eta > 0) {
484  if (APV->Thickness < 0.04) {
485  histograms.Charge_Vs_PathlengthTECP1[elepos]->Fill((*path)[i], Charge);
486  } // TEC plus, type 1
487  else if (APV->Thickness > 0.04) {
488  histograms.Charge_Vs_PathlengthTECP2[elepos]->Fill((*path)[i], Charge);
489  } // TEC plus, type 2
490  }
491  }
492 
493  } // END OF ON-CLUSTER LOOP
494 
495  //LogDebug("SiStripGainsCalibTreeWorker")<<" for mode"<< m_calibrationMode
496  // <<" entries in histogram:"<< histograms.Charge_Vs_Index[elepos].getEntries()
497 
498  // <<std::endl;
499  }
500 
501  hasProcessed_ = true;
502  }
503 }
504 
505 //********************************************************************************//
506 // ------------ method called once each job just before starting event loop ------------
509  if (bareTkGeomPtr) { // pointer not yet set: called the first time => fill the APVColls
510  auto const& Det = bareTkGeomPtr->dets();
511 
512  edm::LogInfo("SiStripGainsCalibTreeWorker") << " Resetting APV struct" << std::endl;
513 
514  unsigned int Index = 0;
515 
516  for (unsigned int i = 0; i < Det.size(); i++) {
517  DetId Detid = Det[i]->geographicalId();
518  int SubDet = Detid.subdetId();
519 
522  auto DetUnit = dynamic_cast<const StripGeomDetUnit*>(Det[i]);
523  if (!DetUnit)
524  continue;
525 
526  const StripTopology& Topo = DetUnit->specificTopology();
527  unsigned int NAPV = Topo.nstrips() / 128;
528 
529  for (unsigned int j = 0; j < NAPV; j++) {
530  auto APV = std::make_shared<stAPVGain>();
531  APV->Index = Index;
532  APV->Bin = -1;
533  APV->DetId = Detid.rawId();
534  APV->APVId = j;
535  APV->SubDet = SubDet;
536  APV->FitMPV = -1;
537  APV->FitMPVErr = -1;
538  APV->FitWidth = -1;
539  APV->FitWidthErr = -1;
540  APV->FitChi2 = -1;
541  APV->FitNorm = -1;
542  APV->Gain = -1;
543  APV->PreviousGain = 1;
544  APV->PreviousGainTick = 1;
545  APV->x = DetUnit->position().basicVector().x();
546  APV->y = DetUnit->position().basicVector().y();
547  APV->z = DetUnit->position().basicVector().z();
548  APV->Eta = DetUnit->position().basicVector().eta();
549  APV->Phi = DetUnit->position().basicVector().phi();
550  APV->R = DetUnit->position().basicVector().transverse();
551  APV->Thickness = DetUnit->surface().bounds().thickness();
552  APV->NEntries = 0;
553  APV->isMasked = false;
554 
555  histograms.APVsCollOrdered.push_back(APV);
556  histograms.APVsColl[(APV->DetId << 4) | APV->APVId] = APV;
557  Index++;
558  histograms.NStripAPVs++;
559  } // loop on APVs
560  } // if is Strips
561  } // loop on dets
562 
563  for (unsigned int i = 0; i < Det.size();
564  i++) { //Make two loop such that the Pixel information is added at the end --> make transition simpler
565  DetId Detid = Det[i]->geographicalId();
566  int SubDet = Detid.subdetId();
568  auto DetUnit = dynamic_cast<const PixelGeomDetUnit*>(Det[i]);
569  if (!DetUnit)
570  continue;
571 
572  const PixelTopology& Topo = DetUnit->specificTopology();
573  unsigned int NROCRow = Topo.nrows() / (80.);
574  unsigned int NROCCol = Topo.ncolumns() / (52.);
575 
576  for (unsigned int j = 0; j < NROCRow; j++) {
577  for (unsigned int i = 0; i < NROCCol; i++) {
578  auto APV = std::make_shared<stAPVGain>();
579  APV->Index = Index;
580  APV->Bin = -1;
581  APV->DetId = Detid.rawId();
582  APV->APVId = (j << 3 | i);
583  APV->SubDet = SubDet;
584  APV->FitMPV = -1;
585  APV->FitMPVErr = -1;
586  APV->FitWidth = -1;
587  APV->FitWidthErr = -1;
588  APV->FitChi2 = -1;
589  APV->Gain = -1;
590  APV->PreviousGain = 1;
591  APV->PreviousGainTick = 1;
592  APV->x = DetUnit->position().basicVector().x();
593  APV->y = DetUnit->position().basicVector().y();
594  APV->z = DetUnit->position().basicVector().z();
595  APV->Eta = DetUnit->position().basicVector().eta();
596  APV->Phi = DetUnit->position().basicVector().phi();
597  APV->R = DetUnit->position().basicVector().transverse();
598  APV->Thickness = DetUnit->surface().bounds().thickness();
599  APV->isMasked = false; //SiPixelQuality_->IsModuleBad(Detid.rawId());
600  APV->NEntries = 0;
601 
602  histograms.APVsCollOrdered.push_back(APV);
603  histograms.APVsColl[(APV->DetId << 4) | APV->APVId] = APV;
604  Index++;
605  histograms.NPixelDets++;
606 
607  } // loop on ROC cols
608  } // loop on ROC rows
609  } // if Pixel
610  } // loop on Dets
611  } //if (!bareTkGeomPtr_) ...
612 }
613 
614 //********************************************************************************//
617  desc.setUnknown();
618  descriptions.addDefault(desc);
619 }
620 
621 //********************************************************************************//
623  edm::Run const& run,
624  edm::EventSetup const& setup,
626  ibooker.cd();
627  std::string dqm_dir = m_DQMdir;
628  const char* tag = dqm_tag_[statCollectionFromMode(m_calibrationMode.c_str())].c_str();
629 
630  edm::LogInfo("SiStripGainsCalibTreeWorker")
631  << "Setting " << dqm_dir << " in DQM and booking histograms for tag " << tag << std::endl;
632 
633  ibooker.setCurrentFolder(dqm_dir);
634 
635  std::string stag(tag);
636  if (!stag.empty() && stag[0] != '_')
637  stag.insert(0, 1, '_');
638 
639  std::string cvi = std::string("Charge_Vs_Index") + stag;
640  std::string cvpTIB = std::string("Charge_Vs_PathlengthTIB") + stag;
641  std::string cvpTOB = std::string("Charge_Vs_PathlengthTOB") + stag;
642  std::string cvpTIDP = std::string("Charge_Vs_PathlengthTIDP") + stag;
643  std::string cvpTIDM = std::string("Charge_Vs_PathlengthTIDM") + stag;
644  std::string cvpTECP1 = std::string("Charge_Vs_PathlengthTECP1") + stag;
645  std::string cvpTECP2 = std::string("Charge_Vs_PathlengthTECP2") + stag;
646  std::string cvpTECM1 = std::string("Charge_Vs_PathlengthTECM1") + stag;
647  std::string cvpTECM2 = std::string("Charge_Vs_PathlengthTECM2") + stag;
648 
649  int elepos = statCollectionFromMode(tag);
650 
651  histograms.Charge_Vs_Index.reserve(dqm_tag_.size());
652  histograms.Charge_Vs_PathlengthTIB.reserve(dqm_tag_.size());
653  histograms.Charge_Vs_PathlengthTOB.reserve(dqm_tag_.size());
654  histograms.Charge_Vs_PathlengthTIDP.reserve(dqm_tag_.size());
655  histograms.Charge_Vs_PathlengthTIDM.reserve(dqm_tag_.size());
656  histograms.Charge_Vs_PathlengthTECP1.reserve(dqm_tag_.size());
657  histograms.Charge_Vs_PathlengthTECP2.reserve(dqm_tag_.size());
658  histograms.Charge_Vs_PathlengthTECM1.reserve(dqm_tag_.size());
659  histograms.Charge_Vs_PathlengthTECM2.reserve(dqm_tag_.size());
660 
661  // The cluster charge is stored by exploiting a non uniform binning in order
662  // reduce the histogram memory size. The bin width is relaxed with a falling
663  // exponential function and the bin boundaries are stored in the binYarray.
664  // The binXarray is used to provide as many bins as the APVs.
665  //
666  // More details about this implementations are here:
667  // https://indico.cern.ch/event/649344/contributions/2672267/attachments/1498323/2332518/OptimizeChHisto.pdf
668 
669  std::vector<float> binXarray;
670  binXarray.reserve(histograms.NStripAPVs + 1);
671  for (unsigned int a = 0; a <= histograms.NStripAPVs; a++) {
672  binXarray.push_back((float)a);
673  }
674 
675  std::array<float, 688> binYarray;
676  double p0 = 5.445;
677  double p1 = 0.002113;
678  double p2 = 69.01576;
679  double y = 0.;
680  for (int b = 0; b < 687; b++) {
681  binYarray[b] = y;
682  if (y <= 902.)
683  y = y + 2.;
684  else
685  y = (p0 - log(exp(p0 - p1 * y) - p2 * p1)) / p1;
686  }
687  binYarray[687] = 4000.;
688 
689  histograms.Charge_1[elepos].clear();
690  histograms.Charge_2[elepos].clear();
691  histograms.Charge_3[elepos].clear();
692  histograms.Charge_4[elepos].clear();
693 
694  auto it = histograms.Charge_Vs_Index.begin();
695  histograms.Charge_Vs_Index.insert(
696  it + elepos,
697  ibooker.book2S(cvi.c_str(), cvi.c_str(), histograms.NStripAPVs, &binXarray[0], 687, binYarray.data()));
698 
699  it = histograms.Charge_Vs_PathlengthTIB.begin();
700  histograms.Charge_Vs_PathlengthTIB.insert(it + elepos,
701  ibooker.book2S(cvpTIB.c_str(), cvpTIB.c_str(), 20, 0.3, 1.3, 250, 0, 2000));
702 
703  it = histograms.Charge_Vs_PathlengthTOB.begin();
704  histograms.Charge_Vs_PathlengthTOB.insert(it + elepos,
705  ibooker.book2S(cvpTOB.c_str(), cvpTOB.c_str(), 20, 0.3, 1.3, 250, 0, 2000));
706 
707  it = histograms.Charge_Vs_PathlengthTIDP.begin();
708  histograms.Charge_Vs_PathlengthTIDP.insert(
709  it + elepos, ibooker.book2S(cvpTIDP.c_str(), cvpTIDP.c_str(), 20, 0.3, 1.3, 250, 0, 2000));
710 
711  it = histograms.Charge_Vs_PathlengthTIDM.begin();
712  histograms.Charge_Vs_PathlengthTIDM.insert(
713  it + elepos, ibooker.book2S(cvpTIDM.c_str(), cvpTIDM.c_str(), 20, 0.3, 1.3, 250, 0, 2000));
714 
715  it = histograms.Charge_Vs_PathlengthTECP1.begin();
716  histograms.Charge_Vs_PathlengthTECP1.insert(
717  it + elepos, ibooker.book2S(cvpTECP1.c_str(), cvpTECP1.c_str(), 20, 0.3, 1.3, 250, 0, 2000));
718 
719  it = histograms.Charge_Vs_PathlengthTECP2.begin();
720  histograms.Charge_Vs_PathlengthTECP2.insert(
721  it + elepos, ibooker.book2S(cvpTECP2.c_str(), cvpTECP2.c_str(), 20, 0.3, 1.3, 250, 0, 2000));
722 
723  it = histograms.Charge_Vs_PathlengthTECM1.begin();
724  histograms.Charge_Vs_PathlengthTECM1.insert(
725  it + elepos, ibooker.book2S(cvpTECM1.c_str(), cvpTECM1.c_str(), 20, 0.3, 1.3, 250, 0, 2000));
726 
727  it = histograms.Charge_Vs_PathlengthTECM2.begin();
728  histograms.Charge_Vs_PathlengthTECM2.insert(
729  it + elepos, ibooker.book2S(cvpTECM2.c_str(), cvpTECM2.c_str(), 20, 0.3, 1.3, 250, 0, 2000));
730 
731  std::vector<std::pair<std::string, std::string>> hnames =
733  for (unsigned int i = 0; i < hnames.size(); i++) {
734  std::string htag = (hnames[i]).first + stag;
735  histograms.Charge_1[elepos].push_back(ibooker.book1DD(htag.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.));
736  }
737 
739  for (unsigned int i = 0; i < hnames.size(); i++) {
740  std::string htag = (hnames[i]).first + stag;
741  histograms.Charge_2[elepos].push_back(ibooker.book1DD(htag.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.));
742  }
743 
745  for (unsigned int i = 0; i < hnames.size(); i++) {
746  std::string htag = (hnames[i]).first + stag;
747  histograms.Charge_3[elepos].push_back(ibooker.book1DD(htag.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.));
748  }
749 
751  for (unsigned int i = 0; i < hnames.size(); i++) {
752  std::string htag = (hnames[i]).first + stag;
753  histograms.Charge_4[elepos].push_back(ibooker.book1DD(htag.c_str(), (hnames[i]).second.c_str(), 100, 0., 1000.));
754  }
755 }
756 
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &, APVGain::APVGainHistograms &) const override
static constexpr auto TEC
virtual int nstrips() const =0
void dqmAnalyze(edm::Event const &, edm::EventSetup const &, APVGain::APVGainHistograms const &) const override
std::vector< unsigned int > FetchIndices(std::map< unsigned int, APVloc >, uint32_t, const TrackerTopology *topo=nullptr)
std::vector< std::string > dqm_tag_
T const & getData(const ESGetToken< T, R > &iToken) const noexcept(false)
Definition: EventSetup.h:119
edm::ESGetToken< SiStripGain, SiStripGainRcd > gainToken_
MonitorElement * book2S(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:254
std::vector< std::string > VChargeHisto
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:36
virtual int ncolumns() const =0
edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > tkGeomTokenBR_
std::map< unsigned int, APVloc > theTopologyMap
virtual int nrows() const =0
int subdetectorPlane(uint32_t, const TrackerTopology *)
const DetContainer & dets() const override
Returm a vector of all GeomDet (including all GeomDetUnits)
Log< level::Error, false > LogError
edm::ESGetToken< SiStripQuality, SiStripQualityRcd > qualityToken_
SiStripGainsCalibTreeWorker(const edm::ParameterSet &)
T getUntrackedParameter(std::string const &, T const &) const
U second(std::pair< T, U > const &p)
MonitorElement * book1DD(TString const &name, TString const &title, int nchX, double lowX, double highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:155
int statCollectionFromMode(const char *tag) const
int iEvent
Definition: GenABIO.cc:224
void addDefault(ParameterSetDescription const &psetDescription)
edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > tkGeomToken_
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
void dqmBeginRun(edm::Run const &, edm::EventSetup const &, APVGain::APVGainHistograms &) const override
static constexpr auto TOB
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
Definition: EventSetup.h:130
edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > tTopoToken_
Log< level::Info, false > LogInfo
Definition: DetId.h:17
static constexpr auto TIB
std::vector< std::string > calibTreeFileNames_
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
int subdetectorSide(uint32_t, const TrackerTopology *)
double b
Definition: hdecay.h:118
int thickness(uint32_t)
HLT enums.
double a
Definition: hdecay.h:119
Definition: tree.py:1
Log< level::Warning, false > LogWarning
std::vector< std::pair< std::string, std::string > > monHnames(std::vector< std::string >, bool, const char *tag)
void checkBookAPVColls(const TrackerGeometry *bareTkGeomPtr, APVGain::APVGainHistograms &histograms) const
static constexpr auto TID
int subdetectorId(uint32_t)
Definition: Run.h:45
#define LogDebug(id)
def exit(msg="")