CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
SiStripGainsPCLHarvester.cc
Go to the documentation of this file.
1 // system include files
2 #include <memory>
3 
4 // CMSSW includes
5 #include "CondFormats/SiStripObjects/interface/SiStripApvGain.h"
6 #include "DataFormats/DetId/interface/DetId.h"
7 #include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"
8 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
9 #include "DataFormats/TrackReco/interface/Track.h"
10 #include "DataFormats/TrackReco/interface/TrackFwd.h"
11 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
12 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
13 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
14 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
15 #include "Geometry/CommonDetUnit/interface/TrackingGeometry.h"
16 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
17 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
18 #include "Geometry/TrackerNumberingBuilder/interface/GeometricDet.h"
19 
20 // user include files
21 #include "CalibTracker/SiStripChannelGain/interface/SiStripGainsPCLHarvester.h"
22 #include "CalibTracker/SiStripChannelGain/interface/APVGainHelpers.h"
23 #include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
24 #include "CommonTools/UtilAlgos/interface/TFileService.h"
25 #include "FWCore/Framework/interface/EventSetup.h"
26 #include <iostream>
27 #include <sstream>
28 
29 // ROOT includes
30 #include <TROOT.h>
31 #include <TTree.h>
32 
33 //********************************************************************************//
34 SiStripGainsPCLHarvester::SiStripGainsPCLHarvester(const edm::ParameterSet& ps)
35  : doStoreOnDB(false), GOOD(0), BAD(0), MASKED(0), NStripAPVs(0), NPixelDets(0) {
36  m_Record = ps.getUntrackedParameter<std::string>("Record", "SiStripApvGainRcd");
37  CalibrationLevel = ps.getUntrackedParameter<int>("CalibrationLevel", 0);
38  MinNrEntries = ps.getUntrackedParameter<double>("minNrEntries", 20);
39  m_DQMdir = ps.getUntrackedParameter<std::string>("DQMdir", "AlCaReco/SiStripGains");
40  m_calibrationMode = ps.getUntrackedParameter<std::string>("calibrationMode", "StdBunch");
41  tagCondition_NClusters = ps.getUntrackedParameter<double>("NClustersForTagProd", 2E8);
42  tagCondition_GoodFrac = ps.getUntrackedParameter<double>("GoodFracForTagProd", 0.95);
43  doChargeMonitorPerPlane = ps.getUntrackedParameter<bool>("doChargeMonitorPerPlane", false);
44  VChargeHisto = ps.getUntrackedParameter<std::vector<std::string>>("ChargeHisto");
45  fit_gaussianConvolution_ = ps.getUntrackedParameter<bool>("FitGaussianConvolution", false);
46  fit_gaussianConvolutionTOBL56_ = ps.getUntrackedParameter<bool>("FitGaussianConvolutionTOBL5L6", false);
47  fit_dataDrivenRange_ = ps.getUntrackedParameter<bool>("FitDataDrivenRange", false);
48  storeGainsTree_ = ps.getUntrackedParameter<bool>("StoreGainsTree", false);
49 
50  //Set the monitoring element tag and store
51  dqm_tag_.reserve(7);
52  dqm_tag_.clear();
53  dqm_tag_.push_back("StdBunch"); // statistic collection from Standard Collision Bunch @ 3.8 T
54  dqm_tag_.push_back("StdBunch0T"); // statistic collection from Standard Collision Bunch @ 0 T
55  dqm_tag_.push_back("AagBunch"); // statistic collection from First Collision After Abort Gap @ 3.8 T
56  dqm_tag_.push_back("AagBunch0T"); // statistic collection from First Collision After Abort Gap @ 0 T
57  dqm_tag_.push_back("IsoMuon"); // statistic collection from Isolated Muon @ 3.8 T
58  dqm_tag_.push_back("IsoMuon0T"); // statistic collection from Isolated Muon @ 0 T
59  dqm_tag_.push_back("Harvest"); // statistic collection: Harvest
60 
61  tTopoToken_ = esConsumes<edm::Transition::EndRun>();
62  tkGeomToken_ = esConsumes<edm::Transition::BeginRun>();
63  gainToken_ = esConsumes<edm::Transition::BeginRun>();
64  qualityToken_ = esConsumes<edm::Transition::BeginRun>();
65 }
66 
67 //********************************************************************************//
68 // ------------ method called for each event ------------
69 void SiStripGainsPCLHarvester::beginRun(edm::Run const& run, const edm::EventSetup& iSetup) {
70  using namespace edm;
71  static constexpr float defaultGainTick = 690. / 640.;
72 
73  this->checkBookAPVColls(iSetup); // check whether APV colls are booked and do so if not yet done
74 
75  const auto gainHandle = iSetup.getHandle(gainToken_);
76  if (!gainHandle.isValid()) {
77  edm::LogError("SiStripGainPCLHarvester") << "gainHandle is not valid\n";
78  exit(0);
79  }
80 
81  const auto& stripQuality = iSetup.getData(qualityToken_);
82 
83  for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
84  std::shared_ptr<stAPVGain> APV = APVsCollOrdered[a];
85 
86  if (APV->SubDet == PixelSubdetector::PixelBarrel || APV->SubDet == PixelSubdetector::PixelEndcap)
87  continue;
88 
89  APV->isMasked = stripQuality.IsApvBad(APV->DetId, APV->APVId);
90 
91  if (gainHandle->getNumberOfTags() != 2) {
92  edm::LogError("SiStripGainPCLHarvester") << "NUMBER OF GAIN TAG IS EXPECTED TO BE 2\n";
93  fflush(stdout);
94  exit(0);
95  };
96  float newPreviousGain = gainHandle->getApvGain(APV->APVId, gainHandle->getRange(APV->DetId, 1), 1);
97  if (APV->PreviousGain != 1 and newPreviousGain != APV->PreviousGain)
98  edm::LogWarning("SiStripGainPCLHarvester") << "WARNING: ParticleGain in the global tag changed\n";
99  APV->PreviousGain = newPreviousGain;
100 
101  float newPreviousGainTick =
102  APV->isMasked ? defaultGainTick : gainHandle->getApvGain(APV->APVId, gainHandle->getRange(APV->DetId, 0), 0);
103  if (APV->PreviousGainTick != 1 and newPreviousGainTick != APV->PreviousGainTick) {
104  edm::LogWarning("SiStripGainPCLHarvester")
105  << "WARNING: TickMarkGain in the global tag changed\n"
106  << std::endl
107  << " APV->SubDet: " << APV->SubDet << " APV->APVId:" << APV->APVId << std::endl
108  << " APV->PreviousGainTick: " << APV->PreviousGainTick << " newPreviousGainTick: " << newPreviousGainTick
109  << std::endl;
110  }
111  APV->PreviousGainTick = newPreviousGainTick;
112  }
113 }
114 
115 //********************************************************************************//
116 void SiStripGainsPCLHarvester::dqmEndJob(DQMStore::IBooker& ibooker_, DQMStore::IGetter& igetter_) {
117  edm::LogInfo("SiStripGainsPCLHarvester") << "Starting harvesting statistics" << std::endl;
118 
119  std::string DQM_dir = m_DQMdir;
120 
121  std::string stag = *(std::find(dqm_tag_.begin(), dqm_tag_.end(), m_calibrationMode));
122  if (!stag.empty() && stag[0] != '_')
123  stag.insert(0, 1, '_');
124 
125  std::string cvi = DQM_dir + std::string("/Charge_Vs_Index") + stag;
126 
127  MonitorElement* Charge_Vs_Index = igetter_.get(cvi);
128 
129  if (Charge_Vs_Index == nullptr) {
130  edm::LogError("SiStripGainsPCLHarvester")
131  << "Harvesting: could not retrieve " << cvi.c_str() << ", statistics will not be summed!" << std::endl;
132  } else {
133  edm::LogInfo("SiStripGainsPCLHarvester")
134  << "Harvesting " << (Charge_Vs_Index)->getTH2S()->GetEntries() << " more clusters" << std::endl;
135  }
136 
137  algoComputeMPVandGain(Charge_Vs_Index);
138  std::unique_ptr<SiStripApvGain> theAPVGains = this->getNewObject(Charge_Vs_Index);
139 
140  // write out the APVGains record
141  edm::Service<cond::service::PoolDBOutputService> poolDbService;
142 
143  if (doStoreOnDB) {
144  if (poolDbService.isAvailable())
145  poolDbService->writeOneIOV(*theAPVGains, poolDbService->currentTime(), m_Record);
146  else
147  throw std::runtime_error("PoolDBService required.");
148  } else {
149  edm::LogInfo("SiStripGainsPCLHarvester") << "Will not produce payload!" << std::endl;
150  }
151 
152  //Collect the statistics for monitoring and validation
153  gainQualityMonitor(ibooker_, Charge_Vs_Index);
154 
155  if (storeGainsTree_) {
156  if (Charge_Vs_Index != nullptr) {
157  storeGainsTree(Charge_Vs_Index->getTH2S()->GetXaxis());
158  } else {
159  edm::LogError("SiStripGainsPCLHarvester")
160  << "Harvesting: could not retrieve " << cvi.c_str() << "Tree won't be stored" << std::endl;
161  }
162  }
163 }
164 
165 //********************************************************************************//
166 void SiStripGainsPCLHarvester::gainQualityMonitor(DQMStore::IBooker& ibooker_,
167  const MonitorElement* Charge_Vs_Index) const {
168  ibooker_.setCurrentFolder("AlCaReco/SiStripGainsHarvesting/");
169 
170  std::vector<APVGain::APVmon> new_charge_histos;
171  std::vector<std::pair<std::string, std::string>> cnames =
172  APVGain::monHnames(VChargeHisto, doChargeMonitorPerPlane, "newG2");
173  for (unsigned int i = 0; i < cnames.size(); i++) {
174  MonitorElement* monitor = ibooker_.book1DD((cnames[i]).first, (cnames[i]).second.c_str(), 100, 0., 1000.);
175  int thick = APVGain::thickness((cnames[i]).first);
176  int id = APVGain::subdetectorId((cnames[i]).first);
177  int side = APVGain::subdetectorSide((cnames[i]).first);
178  int plane = APVGain::subdetectorPlane((cnames[i]).first);
179  new_charge_histos.push_back(APVGain::APVmon(thick, id, side, plane, monitor));
180  }
181 
182  int MPVbin = 300;
183  float MPVmin = 0.;
184  float MPVmax = 600.;
185 
186  MonitorElement* MPV_Vs_EtaTIB =
187  ibooker_.book2DD("MPVvsEtaTIB", "MPV vs Eta TIB", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
188  MonitorElement* MPV_Vs_EtaTID =
189  ibooker_.book2DD("MPVvsEtaTID", "MPV vs Eta TID", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
190  MonitorElement* MPV_Vs_EtaTOB =
191  ibooker_.book2DD("MPVvsEtaTOB", "MPV vs Eta TOB", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
192  MonitorElement* MPV_Vs_EtaTEC =
193  ibooker_.book2DD("MPVvsEtaTEC", "MPV vs Eta TEC", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
194  MonitorElement* MPV_Vs_EtaTECthin =
195  ibooker_.book2DD("MPVvsEtaTEC1", "MPV vs Eta TEC-thin", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
196  MonitorElement* MPV_Vs_EtaTECthick =
197  ibooker_.book2DD("MPVvsEtaTEC2", "MPV vs Eta TEC-thick", 50, -3.0, 3.0, MPVbin, MPVmin, MPVmax);
198 
199  MonitorElement* MPV_Vs_PhiTIB =
200  ibooker_.book2DD("MPVvsPhiTIB", "MPV vs Phi TIB", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
201  MonitorElement* MPV_Vs_PhiTID =
202  ibooker_.book2DD("MPVvsPhiTID", "MPV vs Phi TID", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
203  MonitorElement* MPV_Vs_PhiTOB =
204  ibooker_.book2DD("MPVvsPhiTOB", "MPV vs Phi TOB", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
205  MonitorElement* MPV_Vs_PhiTEC =
206  ibooker_.book2DD("MPVvsPhiTEC", "MPV vs Phi TEC", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
207  MonitorElement* MPV_Vs_PhiTECthin =
208  ibooker_.book2DD("MPVvsPhiTEC1", "MPV vs Phi TEC-thin ", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
209  MonitorElement* MPV_Vs_PhiTECthick =
210  ibooker_.book2DD("MPVvsPhiTEC2", "MPV vs Phi TEC-thick", 50, -3.4, 3.4, MPVbin, MPVmin, MPVmax);
211 
212  MonitorElement* NoMPVfit = ibooker_.book2DD("NoMPVfit", "Modules with bad Landau Fit", 350, -350, 350, 240, 0, 120);
213  MonitorElement* NoMPVmasked = ibooker_.book2DD("NoMPVmasked", "Masked Modules", 350, -350, 350, 240, 0, 120);
214 
215  MonitorElement* Gains = ibooker_.book1DD("Gains", "Gains", 300, 0, 2);
216  MonitorElement* MPVs = ibooker_.book1DD("MPVs", "MPVs", MPVbin, MPVmin, MPVmax);
217  MonitorElement* MPVs320 = ibooker_.book1DD("MPV_320", "MPV 320 thickness", MPVbin, MPVmin, MPVmax);
218  MonitorElement* MPVs500 = ibooker_.book1DD("MPV_500", "MPV 500 thickness", MPVbin, MPVmin, MPVmax);
219  MonitorElement* MPVsTIB = ibooker_.book1DD("MPV_TIB", "MPV TIB", MPVbin, MPVmin, MPVmax);
220  MonitorElement* MPVsTID = ibooker_.book1DD("MPV_TID", "MPV TID", MPVbin, MPVmin, MPVmax);
221  MonitorElement* MPVsTIDP = ibooker_.book1DD("MPV_TIDP", "MPV TIDP", MPVbin, MPVmin, MPVmax);
222  MonitorElement* MPVsTIDM = ibooker_.book1DD("MPV_TIDM", "MPV TIDM", MPVbin, MPVmin, MPVmax);
223  MonitorElement* MPVsTOB = ibooker_.book1DD("MPV_TOB", "MPV TOB", MPVbin, MPVmin, MPVmax);
224  MonitorElement* MPVsTEC = ibooker_.book1DD("MPV_TEC", "MPV TEC", MPVbin, MPVmin, MPVmax);
225  MonitorElement* MPVsTECP = ibooker_.book1DD("MPV_TECP", "MPV TECP", MPVbin, MPVmin, MPVmax);
226  MonitorElement* MPVsTECM = ibooker_.book1DD("MPV_TECM", "MPV TECM", MPVbin, MPVmin, MPVmax);
227  MonitorElement* MPVsTECthin = ibooker_.book1DD("MPV_TEC1", "MPV TEC thin", MPVbin, MPVmin, MPVmax);
228  MonitorElement* MPVsTECthick = ibooker_.book1DD("MPV_TEC2", "MPV TEC thick", MPVbin, MPVmin, MPVmax);
229  MonitorElement* MPVsTECP1 = ibooker_.book1DD("MPV_TECP1", "MPV TECP thin ", MPVbin, MPVmin, MPVmax);
230  MonitorElement* MPVsTECP2 = ibooker_.book1DD("MPV_TECP2", "MPV TECP thick", MPVbin, MPVmin, MPVmax);
231  MonitorElement* MPVsTECM1 = ibooker_.book1DD("MPV_TECM1", "MPV TECM thin", MPVbin, MPVmin, MPVmax);
232  MonitorElement* MPVsTECM2 = ibooker_.book1DD("MPV_TECM2", "MPV TECM thick", MPVbin, MPVmin, MPVmax);
233 
234  MonitorElement* MPVError = ibooker_.book1DD("MPVError", "MPV Error", 150, 0, 150);
235  MonitorElement* MPVErrorVsMPV = ibooker_.book2DD("MPVErrorVsMPV", "MPV Error vs MPV", 300, 0, 600, 150, 0, 150);
236  MonitorElement* MPVErrorVsEta = ibooker_.book2DD("MPVErrorVsEta", "MPV Error vs Eta", 50, -3.0, 3.0, 150, 0, 150);
237  MonitorElement* MPVErrorVsPhi = ibooker_.book2DD("MPVErrorVsPhi", "MPV Error vs Phi", 50, -3.4, 3.4, 150, 0, 150);
238  MonitorElement* MPVErrorVsN = ibooker_.book2DD("MPVErrorVsN", "MPV Error vs N", 500, 0, 1000, 150, 0, 150);
239 
240  MonitorElement* DiffWRTPrevGainTIB = ibooker_.book1DD("DiffWRTPrevGainTIB", "Diff w.r.t. PrevGain TIB", 250, 0, 2);
241  MonitorElement* DiffWRTPrevGainTID = ibooker_.book1DD("DiffWRTPrevGainTID", "Diff w.r.t. PrevGain TID", 250, 0, 2);
242  MonitorElement* DiffWRTPrevGainTOB = ibooker_.book1DD("DiffWRTPrevGainTOB", "Diff w.r.t. PrevGain TOB", 250, 0, 2);
243  MonitorElement* DiffWRTPrevGainTEC = ibooker_.book1DD("DiffWRTPrevGainTEC", "Diff w.r.t. PrevGain TEC", 250, 0, 2);
244 
245  MonitorElement* GainVsPrevGainTIB =
246  ibooker_.book2DD("GainVsPrevGainTIB", "Gain vs PrevGain TIB", 100, 0, 2, 100, 0, 2);
247  MonitorElement* GainVsPrevGainTID =
248  ibooker_.book2DD("GainVsPrevGainTID", "Gain vs PrevGain TID", 100, 0, 2, 100, 0, 2);
249  MonitorElement* GainVsPrevGainTOB =
250  ibooker_.book2DD("GainVsPrevGainTOB", "Gain vs PrevGain TOB", 100, 0, 2, 100, 0, 2);
251  MonitorElement* GainVsPrevGainTEC =
252  ibooker_.book2DD("GainVsPrevGainTEC", "Gain vs PrevGain TEC", 100, 0, 2, 100, 0, 2);
253 
254  for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
255  std::shared_ptr<stAPVGain> APV = APVsCollOrdered[a];
256  if (APV == nullptr)
257  continue;
258 
259  unsigned int Index = APV->Index;
260  unsigned int SubDet = APV->SubDet;
261  unsigned int DetId = APV->DetId;
262  float z = APV->z;
263  float Eta = APV->Eta;
264  float R = APV->R;
265  float Phi = APV->Phi;
266  float Thickness = APV->Thickness;
267  double FitMPV = APV->FitMPV;
268  double FitMPVErr = APV->FitMPVErr;
269  double Gain = APV->Gain;
270  double NEntries = APV->NEntries;
271  double PreviousGain = APV->PreviousGain;
272 
273  if (SubDet < 3)
274  continue; // avoid to loop over Pixel det id
275 
276  if (Gain != 1.) {
277  std::vector<MonitorElement*> charge_histos = APVGain::FetchMonitor(new_charge_histos, DetId, tTopo_.get());
278 
279  if (!Charge_Vs_Index)
280  continue;
281  TH2S* chvsidx = (Charge_Vs_Index)->getTH2S();
282  int bin = chvsidx->GetXaxis()->FindBin(Index);
283  TH1D* Proj = chvsidx->ProjectionY("proj", bin, bin);
284  for (int binId = 0; binId < Proj->GetXaxis()->GetNbins(); binId++) {
285  double new_charge = Proj->GetXaxis()->GetBinCenter(binId) / Gain;
286  if (Proj->GetBinContent(binId) != 0.) {
287  for (unsigned int h = 0; h < charge_histos.size(); h++) {
288  TH1D* chisto = (charge_histos[h])->getTH1D();
289  for (int e = 0; e < Proj->GetBinContent(binId); e++)
290  chisto->Fill(new_charge);
291  }
292  }
293  }
294  }
295 
296  if (FitMPV <= 0.) { // No fit of MPV
297  if (APV->isMasked)
298  NoMPVmasked->Fill(z, R);
299  else
300  NoMPVfit->Fill(z, R);
301 
302  } else { // Fit of MPV
303  if (FitMPV > 0.)
304  Gains->Fill(Gain);
305 
306  MPVs->Fill(FitMPV);
307  if (Thickness < 0.04)
308  MPVs320->Fill(FitMPV);
309  if (Thickness > 0.04)
310  MPVs500->Fill(FitMPV);
311 
312  MPVError->Fill(FitMPVErr);
313  MPVErrorVsMPV->Fill(FitMPV, FitMPVErr);
314  MPVErrorVsEta->Fill(Eta, FitMPVErr);
315  MPVErrorVsPhi->Fill(Phi, FitMPVErr);
316  MPVErrorVsN->Fill(NEntries, FitMPVErr);
317 
318  if (SubDet == 3) {
319  MPV_Vs_EtaTIB->Fill(Eta, FitMPV);
320  MPV_Vs_PhiTIB->Fill(Phi, FitMPV);
321  MPVsTIB->Fill(FitMPV);
322 
323  } else if (SubDet == 4) {
324  MPV_Vs_EtaTID->Fill(Eta, FitMPV);
325  MPV_Vs_PhiTID->Fill(Phi, FitMPV);
326  MPVsTID->Fill(FitMPV);
327  if (Eta < 0.)
328  MPVsTIDM->Fill(FitMPV);
329  if (Eta > 0.)
330  MPVsTIDP->Fill(FitMPV);
331 
332  } else if (SubDet == 5) {
333  MPV_Vs_EtaTOB->Fill(Eta, FitMPV);
334  MPV_Vs_PhiTOB->Fill(Phi, FitMPV);
335  MPVsTOB->Fill(FitMPV);
336 
337  } else if (SubDet == 6) {
338  MPV_Vs_EtaTEC->Fill(Eta, FitMPV);
339  MPV_Vs_PhiTEC->Fill(Phi, FitMPV);
340  MPVsTEC->Fill(FitMPV);
341  if (Eta < 0.)
342  MPVsTECM->Fill(FitMPV);
343  if (Eta > 0.)
344  MPVsTECP->Fill(FitMPV);
345  if (Thickness < 0.04) {
346  MPV_Vs_EtaTECthin->Fill(Eta, FitMPV);
347  MPV_Vs_PhiTECthin->Fill(Phi, FitMPV);
348  MPVsTECthin->Fill(FitMPV);
349  if (Eta > 0.)
350  MPVsTECP1->Fill(FitMPV);
351  if (Eta < 0.)
352  MPVsTECM1->Fill(FitMPV);
353  }
354  if (Thickness > 0.04) {
355  MPV_Vs_EtaTECthick->Fill(Eta, FitMPV);
356  MPV_Vs_PhiTECthick->Fill(Phi, FitMPV);
357  MPVsTECthick->Fill(FitMPV);
358  if (Eta > 0.)
359  MPVsTECP2->Fill(FitMPV);
360  if (Eta < 0.)
361  MPVsTECM2->Fill(FitMPV);
362  }
363  }
364  }
365 
366  if (SubDet == 3 && PreviousGain != 0.)
367  DiffWRTPrevGainTIB->Fill(Gain / PreviousGain);
368  else if (SubDet == 4 && PreviousGain != 0.)
369  DiffWRTPrevGainTID->Fill(Gain / PreviousGain);
370  else if (SubDet == 5 && PreviousGain != 0.)
371  DiffWRTPrevGainTOB->Fill(Gain / PreviousGain);
372  else if (SubDet == 6 && PreviousGain != 0.)
373  DiffWRTPrevGainTEC->Fill(Gain / PreviousGain);
374 
375  if (SubDet == 3)
376  GainVsPrevGainTIB->Fill(PreviousGain, Gain);
377  else if (SubDet == 4)
378  GainVsPrevGainTID->Fill(PreviousGain, Gain);
379  else if (SubDet == 5)
380  GainVsPrevGainTOB->Fill(PreviousGain, Gain);
381  else if (SubDet == 6)
382  GainVsPrevGainTEC->Fill(PreviousGain, Gain);
383  }
384 }
385 
386 namespace {
387 
388  std::pair<double, double> findFitRange(TH1* inputHisto) {
389  const auto prevErrorIgnoreLevel = gErrorIgnoreLevel;
390  gErrorIgnoreLevel = kError;
391  auto charge_clone = std::unique_ptr<TH1D>(dynamic_cast<TH1D*>(inputHisto->Rebin(10, "charge_clone")));
392  gErrorIgnoreLevel = prevErrorIgnoreLevel;
393  float max_content = -1;
394  float xMax = -1;
395  for (int i = 1; i < charge_clone->GetNbinsX() + 1; ++i) {
396  const auto bin_content = charge_clone->GetBinContent(i);
397  const auto bin_center = charge_clone->GetXaxis()->GetBinCenter(i);
398  if ((bin_content > max_content) && (bin_center > 100.)) {
399  max_content = bin_content;
400  xMax = bin_center;
401  }
402  }
403  return std::pair<double, double>(xMax - 100., xMax + 500.);
404  }
405 
406  Double_t langaufun(Double_t* x, Double_t* par) {
407  // Numeric constants
408  Double_t invsq2pi = 0.3989422804014; // (2 pi)^(-1/2)
409  Double_t mpshift = -0.22278298; // Landau maximum location
410 
411  // Control constants
412  Double_t np = 100.0; // number of convolution steps
413  Double_t sc = 5.0; // convolution extends to +-sc Gaussian sigmas
414 
415  // Variables
416  Double_t xx;
417  Double_t mpc;
418  Double_t fland;
419  Double_t sum = 0.0;
420  Double_t xlow, xupp;
421  Double_t step;
422  Double_t i;
423 
424  // MP shift correction
425  mpc = par[1] - mpshift * par[0];
426 
427  // Range of convolution integral
428  xlow = x[0] - sc * par[3];
429  xupp = x[0] + sc * par[3];
430 
431  step = (xupp - xlow) / np;
432 
433  // Convolution integral of Landau and Gaussian by sum
434  for (i = 1.0; i <= np / 2; i++) {
435  xx = xlow + (i - .5) * step;
436  fland = TMath::Landau(xx, mpc, par[0]) / par[0];
437  sum += fland * TMath::Gaus(x[0], xx, par[3]);
438 
439  xx = xupp - (i - .5) * step;
440  fland = TMath::Landau(xx, mpc, par[0]) / par[0];
441 
442  sum += fland * TMath::Gaus(x[0], xx, par[3]);
443  }
444 
445  return (par[2] * step * sum * invsq2pi / par[3]);
446  }
447 
448  std::unique_ptr<TF1> langaufit(TH1D* his,
449  Double_t* fitrange,
450  Double_t* startvalues,
451  Double_t* parlimitslo,
452  Double_t* parlimitshi,
453  Double_t* fitparams,
454  Double_t* fiterrors,
455  Double_t* ChiSqr,
456  Int_t* NDF) {
457  Int_t i;
458  Char_t FunName[100];
459 
460  sprintf(FunName, "Fitfcn_%s", his->GetName());
461 
462  TF1* ffitold = (TF1*)gROOT->GetListOfFunctions()->FindObject(FunName);
463  if (ffitold)
464  delete ffitold;
465 
466  auto ffit = std::make_unique<TF1>(FunName, langaufun, fitrange[0], fitrange[1], 4);
467  //TF1 *ffit = new TF1(FunName,langaufun,his->GetXaxis()->GetXmin(),his->GetXaxis()->GetXmax(),4);
468  ffit->SetParameters(startvalues);
469  ffit->SetParNames("Width", "MP", "Area", "GSigma");
470 
471  for (i = 0; i < 4; i++) {
472  ffit->SetParLimits(i, parlimitslo[i], parlimitshi[i]);
473  }
474 
475  his->Fit(FunName, "QRB0"); // fit within specified range, use ParLimits, do not plot
476 
477  ffit->GetParameters(fitparams); // obtain fit parameters
478  for (i = 0; i < 4; i++) {
479  fiterrors[i] = ffit->GetParError(i); // obtain fit parameter errors
480  }
481  ChiSqr[0] = ffit->GetChisquare(); // obtain chi^2
482  NDF[0] = ffit->GetNDF(); // obtain ndf
483 
484  return ffit; // return fit function
485  }
486 } // namespace
487 
488 //********************************************************************************//
489 void SiStripGainsPCLHarvester::algoComputeMPVandGain(const MonitorElement* Charge_Vs_Index) {
490  unsigned int I = 0;
491  TH1F* Proj = nullptr;
492  static constexpr double DEF_F = -9999.;
493  double FitResults[6] = {DEF_F, DEF_F, DEF_F, DEF_F, DEF_F, DEF_F};
494  double MPVmean = 300;
495 
496  if (Charge_Vs_Index == nullptr) {
497  edm::LogError("SiStripGainsPCLHarvester")
498  << "Harvesting: could not execute algoComputeMPVandGain method because " << m_calibrationMode
499  << " statistics cannot be retrieved.\n"
500  << "Please check if input contains " << m_calibrationMode << " data." << std::endl;
501  return;
502  }
503 
504  TH2S* chvsidx = (Charge_Vs_Index)->getTH2S();
505 
506  printf("Progressing Bar :0%% 20%% 40%% 60%% 80%% 100%%\n");
507  printf("Fitting Charge Distribution :");
508  int TreeStep = APVsColl.size() / 50;
509 
510  for (auto it = APVsColl.begin(); it != APVsColl.end(); it++, I++) {
511  if (I % TreeStep == 0) {
512  printf(".");
513  fflush(stdout);
514  }
515  std::shared_ptr<stAPVGain> APV = it->second;
516  if (APV->Bin < 0)
517  APV->Bin = chvsidx->GetXaxis()->FindBin(APV->Index);
518 
519  if (APV->isMasked) {
520  APV->Gain = APV->PreviousGain;
521  MASKED++;
522  continue;
523  }
524 
525  Proj = (TH1F*)(chvsidx->ProjectionY(
526  "", chvsidx->GetXaxis()->FindBin(APV->Index), chvsidx->GetXaxis()->FindBin(APV->Index), "e"));
527  if (!Proj)
528  continue;
529 
530  if (CalibrationLevel == 0) {
531  } else if (CalibrationLevel == 1) {
532  int SecondAPVId = APV->APVId;
533  if (SecondAPVId % 2 == 0) {
534  SecondAPVId = SecondAPVId + 1;
535  } else {
536  SecondAPVId = SecondAPVId - 1;
537  }
538  std::shared_ptr<stAPVGain> APV2 = APVsColl[(APV->DetId << 4) | SecondAPVId];
539  if (APV2->Bin < 0)
540  APV2->Bin = chvsidx->GetXaxis()->FindBin(APV2->Index);
541  TH1F* Proj2 = (TH1F*)(chvsidx->ProjectionY("", APV2->Bin, APV2->Bin, "e"));
542  if (Proj2) {
543  Proj->Add(Proj2, 1);
544  delete Proj2;
545  }
546  } else if (CalibrationLevel == 2) {
547  for (unsigned int i = 0; i < 16; i++) { //loop up to 6APV for Strip and up to 16 for Pixels
548  auto tmpit = APVsColl.find((APV->DetId << 4) | i);
549  if (tmpit == APVsColl.end())
550  continue;
551  std::shared_ptr<stAPVGain> APV2 = tmpit->second;
552  if (APV2->DetId != APV->DetId || APV2->APVId == APV->APVId)
553  continue;
554  if (APV2->Bin < 0)
555  APV2->Bin = chvsidx->GetXaxis()->FindBin(APV2->Index);
556  TH1F* Proj2 = (TH1F*)(chvsidx->ProjectionY("", APV2->Bin, APV2->Bin, "e"));
557  if (Proj2) {
558  Proj->Add(Proj2, 1);
559  delete Proj2;
560  }
561  }
562  } else {
563  CalibrationLevel = 0;
564  printf("Unknown Calibration Level, will assume %i\n", CalibrationLevel);
565  }
566 
567  std::pair<double, double> fitRange{50., 5400.};
568  if (fit_dataDrivenRange_) {
569  fitRange = findFitRange(Proj);
570  }
571 
572  const bool isTOBL5L6 =
573  (DetId{APV->DetId}.subdetId() == StripSubdetector::TOB) && (tTopo_->tobLayer(APV->DetId) > 4);
574  getPeakOfLandau(Proj,
575  FitResults,
576  fitRange.first,
577  fitRange.second,
578  (isTOBL5L6 ? fit_gaussianConvolutionTOBL56_ : fit_gaussianConvolution_));
579 
580  // throw if the fit results are not set
581  assert(FitResults[0] != DEF_F);
582 
583  APV->FitMPV = FitResults[0];
584  APV->FitMPVErr = FitResults[1];
585  APV->FitWidth = FitResults[2];
586  APV->FitWidthErr = FitResults[3];
587  APV->FitChi2 = FitResults[4];
588  APV->FitNorm = FitResults[5];
589  APV->NEntries = Proj->GetEntries();
590 
591  // fall back to legacy fit in case of very low chi2 probability
592  if (APV->FitChi2 <= 0.1) {
593  edm::LogInfo("SiStripGainsPCLHarvester")
594  << "fit failed on this APV:" << APV->DetId << ":" << APV->APVId << " !" << std::endl;
595 
596  std::fill(FitResults, FitResults + 6, 0);
597  fitRange = std::make_pair(50., 5400.);
598 
599  APV->FitGrade = fitgrade::B;
600 
601  getPeakOfLandau(Proj, FitResults, fitRange.first, fitRange.second, false);
602 
603  APV->FitMPV = FitResults[0];
604  APV->FitMPVErr = FitResults[1];
605  APV->FitWidth = FitResults[2];
606  APV->FitWidthErr = FitResults[3];
607  APV->FitChi2 = FitResults[4];
608  APV->FitNorm = FitResults[5];
609  APV->NEntries = Proj->GetEntries();
610  } else {
611  APV->FitGrade = fitgrade::A;
612  }
613 
614  if (IsGoodLandauFit(FitResults)) {
615  APV->Gain = APV->FitMPV / MPVmean;
616  if (APV->SubDet > 2)
617  GOOD++;
618  } else {
619  APV->Gain = APV->PreviousGain;
620  if (APV->SubDet > 2)
621  BAD++;
622  }
623  if (APV->Gain <= 0)
624  APV->Gain = 1;
625 
626  delete Proj;
627  }
628  printf("\n");
629 }
630 
631 //********************************************************************************//
632 void SiStripGainsPCLHarvester::getPeakOfLandau(
633  TH1* InputHisto, double* FitResults, double LowRange, double HighRange, bool gaussianConvolution) {
634  // undo defaults (checked for assertion)
635  FitResults[0] = -0.5; //MPV
636  FitResults[1] = 0; //MPV error
637  FitResults[2] = -0.5; //Width
638  FitResults[3] = 0; //Width error
639  FitResults[4] = -0.5; //Fit Chi2/NDF
640  FitResults[5] = 0; //Normalization
641 
642  if (InputHisto->GetEntries() < MinNrEntries)
643  return;
644 
645  if (gaussianConvolution) {
646  // perform fit with landau convoluted with a gaussian
647  Double_t fr[2] = {LowRange, HighRange};
648  Double_t sv[4] = {25., 300., InputHisto->Integral(), 40.};
649  Double_t pllo[4] = {0.5, 100., 1.0, 0.4};
650  Double_t plhi[4] = {100., 500., 1.e7, 100.};
651  Double_t fp[4], fpe[4];
652  Double_t chisqr;
653  Int_t ndf;
654  auto fitsnr = langaufit(dynamic_cast<TH1D*>(InputHisto), fr, sv, pllo, plhi, fp, fpe, &chisqr, &ndf);
655  FitResults[0] = fitsnr->GetMaximumX(); //MPV
656  FitResults[1] = fpe[1]; //MPV error // FIXME add error propagation
657  FitResults[2] = fp[0]; //Width
658  FitResults[3] = fpe[0]; //Width error
659  FitResults[4] = chisqr / ndf; //Fit Chi2/NDF
660  FitResults[5] = fp[2];
661  } else {
662  // perform fit with standard landau
663  TF1 MyLandau("MyLandau", "landau", LowRange, HighRange);
664  MyLandau.SetParameter(1, 300);
665  InputHisto->Fit(&MyLandau, "0QR WW");
666  // MPV is parameter 1 (0=constant, 1=MPV, 2=Sigma)
667  FitResults[0] = MyLandau.GetParameter(1); //MPV
668  FitResults[1] = MyLandau.GetParError(1); //MPV error
669  FitResults[2] = MyLandau.GetParameter(2); //Width
670  FitResults[3] = MyLandau.GetParError(2); //Width error
671  FitResults[4] = MyLandau.GetChisquare() / MyLandau.GetNDF(); //Fit Chi2/NDF
672  FitResults[5] = MyLandau.GetParameter(0);
673  }
674 }
675 
676 //********************************************************************************//
677 bool SiStripGainsPCLHarvester::IsGoodLandauFit(double* FitResults) {
678  if (FitResults[0] <= 0)
679  return false;
680  // if(FitResults[1] > MaxMPVError )return false;
681  // if(FitResults[4] > MaxChi2OverNDF)return false;
682  return true;
683 }
684 
685 //********************************************************************************//
686 // ------------ method called once each job just before starting event loop ------------
687 void SiStripGainsPCLHarvester::checkBookAPVColls(const edm::EventSetup& es) {
688  auto newBareTkGeomPtr = &es.getData(tkGeomToken_);
689  if (newBareTkGeomPtr == bareTkGeomPtr_)
690  return; // already filled APVColls, nothing changed
691 
692  if (!bareTkGeomPtr_) { // pointer not yet set: called the first time => fill the APVColls
693  auto const& Det = newBareTkGeomPtr->dets();
694 
695  unsigned int Index = 0;
696 
697  for (unsigned int i = 0; i < Det.size(); i++) {
698  DetId Detid = Det[i]->geographicalId();
699  int SubDet = Detid.subdetId();
700 
701  if (SubDet == StripSubdetector::TIB || SubDet == StripSubdetector::TID || SubDet == StripSubdetector::TOB ||
702  SubDet == StripSubdetector::TEC) {
703  auto DetUnit = dynamic_cast<const StripGeomDetUnit*>(Det[i]);
704  if (!DetUnit)
705  continue;
706 
707  const StripTopology& Topo = DetUnit->specificTopology();
708  unsigned int NAPV = Topo.nstrips() / 128;
709 
710  for (unsigned int j = 0; j < NAPV; j++) {
711  auto APV = std::make_shared<stAPVGain>();
712  APV->Index = Index;
713  APV->Bin = -1;
714  APV->DetId = Detid.rawId();
715  APV->APVId = j;
716  APV->SubDet = SubDet;
717  APV->FitMPV = -1;
718  APV->FitMPVErr = -1;
719  APV->FitWidth = -1;
720  APV->FitWidthErr = -1;
721  APV->FitChi2 = -1;
722  APV->FitNorm = -1;
723  APV->FitGrade = fitgrade::NONE;
724  APV->Gain = -1;
725  APV->PreviousGain = 1;
726  APV->PreviousGainTick = 1;
727  APV->x = DetUnit->position().basicVector().x();
728  APV->y = DetUnit->position().basicVector().y();
729  APV->z = DetUnit->position().basicVector().z();
730  APV->Eta = DetUnit->position().basicVector().eta();
731  APV->Phi = DetUnit->position().basicVector().phi();
732  APV->R = DetUnit->position().basicVector().transverse();
733  APV->Thickness = DetUnit->surface().bounds().thickness();
734  APV->NEntries = 0;
735  APV->isMasked = false;
736 
737  APVsCollOrdered.push_back(APV);
738  APVsColl[(APV->DetId << 4) | APV->APVId] = APV;
739  Index++;
740  NStripAPVs++;
741  } // loop on APVs
742  } // if is Strips
743  } // loop on dets
744 
745  for (unsigned int i = 0; i < Det.size();
746  i++) { //Make two loop such that the Pixel information is added at the end --> make transition simpler
747  DetId Detid = Det[i]->geographicalId();
748  int SubDet = Detid.subdetId();
749  if (SubDet == PixelSubdetector::PixelBarrel || SubDet == PixelSubdetector::PixelEndcap) {
750  auto DetUnit = dynamic_cast<const PixelGeomDetUnit*>(Det[i]);
751  if (!DetUnit)
752  continue;
753 
754  const PixelTopology& Topo = DetUnit->specificTopology();
755  unsigned int NROCRow = Topo.nrows() / (80.);
756  unsigned int NROCCol = Topo.ncolumns() / (52.);
757 
758  for (unsigned int j = 0; j < NROCRow; j++) {
759  for (unsigned int i = 0; i < NROCCol; i++) {
760  auto APV = std::make_shared<stAPVGain>();
761  APV->Index = Index;
762  APV->Bin = -1;
763  APV->DetId = Detid.rawId();
764  APV->APVId = (j << 3 | i);
765  APV->SubDet = SubDet;
766  APV->FitMPV = -1;
767  APV->FitMPVErr = -1;
768  APV->FitWidth = -1;
769  APV->FitWidthErr = -1;
770  APV->FitChi2 = -1;
771  APV->FitGrade = fitgrade::NONE;
772  APV->Gain = -1;
773  APV->PreviousGain = 1;
774  APV->PreviousGainTick = 1;
775  APV->x = DetUnit->position().basicVector().x();
776  APV->y = DetUnit->position().basicVector().y();
777  APV->z = DetUnit->position().basicVector().z();
778  APV->Eta = DetUnit->position().basicVector().eta();
779  APV->Phi = DetUnit->position().basicVector().phi();
780  APV->R = DetUnit->position().basicVector().transverse();
781  APV->Thickness = DetUnit->surface().bounds().thickness();
782  APV->isMasked = false; //SiPixelQuality_->IsModuleBad(Detid.rawId());
783  APV->NEntries = 0;
784 
785  APVsCollOrdered.push_back(APV);
786  APVsColl[(APV->DetId << 4) | APV->APVId] = APV;
787  Index++;
788  NPixelDets++;
789 
790  } // loop on ROC cols
791  } // loop on ROC rows
792  } // if Pixel
793  } // loop on Dets
794  } //if (!bareTkGeomPtr_) ...
795  bareTkGeomPtr_ = newBareTkGeomPtr;
796 }
797 
798 //********************************************************************************//
799 bool SiStripGainsPCLHarvester::produceTagFilter(const MonitorElement* Charge_Vs_Index) {
800  // The goal of this function is to check wether or not there is enough statistics
801  // to produce a meaningful tag for the DB
802 
803  if (Charge_Vs_Index == nullptr) {
804  edm::LogError("SiStripGainsPCLHarvester")
805  << "produceTagFilter -> Return false: could not retrieve the " << m_calibrationMode << " statistics.\n"
806  << "Please check if input contains " << m_calibrationMode << " data." << std::endl;
807  return false;
808  }
809 
810  float integral = (Charge_Vs_Index)->getTH2S()->Integral();
811  if ((Charge_Vs_Index)->getTH2S()->Integral(0, NStripAPVs + 1, 0, 99999) < tagCondition_NClusters) {
812  edm::LogWarning("SiStripGainsPCLHarvester")
813  << "calibrationMode -> " << m_calibrationMode << "\n"
814  << "produceTagFilter -> Return false: Statistics is too low : " << integral << std::endl;
815  return false;
816  }
817  if ((1.0 * GOOD) / (GOOD + BAD) < tagCondition_GoodFrac) {
818  edm::LogWarning("SiStripGainsPCLHarvester")
819  << "calibrationMode -> " << m_calibrationMode << "\n"
820  << "produceTagFilter -> Return false: ratio of GOOD/TOTAL is too low: " << (1.0 * GOOD) / (GOOD + BAD)
821  << std::endl;
822  return false;
823  }
824  return true;
825 }
826 
827 //********************************************************************************//
828 std::unique_ptr<SiStripApvGain> SiStripGainsPCLHarvester::getNewObject(const MonitorElement* Charge_Vs_Index) {
829  std::unique_ptr<SiStripApvGain> obj = std::make_unique<SiStripApvGain>();
830 
831  if (!produceTagFilter(Charge_Vs_Index)) {
832  edm::LogWarning("SiStripGainsPCLHarvester")
833  << "getNewObject -> will not produce a paylaod because produceTagFilter returned false " << std::endl;
834  return obj;
835  } else {
836  doStoreOnDB = true;
837  }
838 
839  std::vector<float> theSiStripVector;
840  unsigned int PreviousDetId = 0;
841  for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
842  std::shared_ptr<stAPVGain> APV = APVsCollOrdered[a];
843  if (APV == nullptr) {
844  printf("Bug\n");
845  continue;
846  }
847  if (APV->SubDet <= 2)
848  continue;
849  if (APV->DetId != PreviousDetId) {
850  if (!theSiStripVector.empty()) {
851  SiStripApvGain::Range range(theSiStripVector.begin(), theSiStripVector.end());
852  if (!obj->put(PreviousDetId, range))
853  printf("Bug to put detId = %i\n", PreviousDetId);
854  }
855  theSiStripVector.clear();
856  PreviousDetId = APV->DetId;
857  }
858  theSiStripVector.push_back(APV->Gain);
859 
860  LogDebug("SiStripGainsPCLHarvester") << " DetId: " << APV->DetId << " APV: " << APV->APVId
861  << " Gain: " << APV->Gain << std::endl;
862  }
863  if (!theSiStripVector.empty()) {
864  SiStripApvGain::Range range(theSiStripVector.begin(), theSiStripVector.end());
865  if (!obj->put(PreviousDetId, range))
866  printf("Bug to put detId = %i\n", PreviousDetId);
867  }
868 
869  return obj;
870 }
871 
872 void SiStripGainsPCLHarvester::storeGainsTree(const TAxis* chVsIdxXaxis) const {
873  unsigned int t_Index, t_Bin, t_DetId;
874  unsigned char t_APVId, t_SubDet;
875  float t_x, t_y, t_z, t_Eta, t_R, t_Phi, t_Thickness;
876  float t_FitMPV, t_FitMPVErr, t_FitWidth, t_FitWidthErr, t_FitChi2NDF, t_FitNorm, t_FitGrade;
877  double t_Gain, t_PrevGain, t_PrevGainTick, t_NEntries;
878  bool t_isMasked;
879  auto tree = edm::Service<TFileService>()->make<TTree>("APVGain", "APVGain");
880  tree->Branch("Index", &t_Index, "Index/i");
881  tree->Branch("Bin", &t_Bin, "Bin/i");
882  tree->Branch("DetId", &t_DetId, "DetId/i");
883  tree->Branch("APVId", &t_APVId, "APVId/b");
884  tree->Branch("SubDet", &t_SubDet, "SubDet/b");
885  tree->Branch("x", &t_x, "x/F");
886  tree->Branch("y", &t_y, "y/F");
887  tree->Branch("z", &t_z, "z/F");
888  tree->Branch("Eta", &t_Eta, "Eta/F");
889  tree->Branch("R", &t_R, "R/F");
890  tree->Branch("Phi", &t_Phi, "Phi/F");
891  tree->Branch("Thickness", &t_Thickness, "Thickness/F");
892  tree->Branch("FitMPV", &t_FitMPV, "FitMPV/F");
893  tree->Branch("FitMPVErr", &t_FitMPVErr, "FitMPVErr/F");
894  tree->Branch("FitWidth", &t_FitWidth, "FitWidth/F");
895  tree->Branch("FitWidthErr", &t_FitWidthErr, "FitWidthErr/F");
896  tree->Branch("FitChi2NDF", &t_FitChi2NDF, "FitChi2NDF/F");
897  tree->Branch("FitNorm", &t_FitNorm, "FitNorm/F");
898  tree->Branch("FitGrade", &t_FitGrade, "FitGrade/F");
899  tree->Branch("Gain", &t_Gain, "Gain/D");
900  tree->Branch("PrevGain", &t_PrevGain, "PrevGain/D");
901  tree->Branch("PrevGainTick", &t_PrevGainTick, "PrevGainTick/D");
902  tree->Branch("NEntries", &t_NEntries, "NEntries/D");
903  tree->Branch("isMasked", &t_isMasked, "isMasked/O");
904 
905  for (const auto& iAPV : APVsCollOrdered) {
906  if (iAPV) {
907  t_Index = iAPV->Index;
908  t_Bin = chVsIdxXaxis->FindBin(iAPV->Index);
909  t_DetId = iAPV->DetId;
910  t_APVId = iAPV->APVId;
911  t_SubDet = iAPV->SubDet;
912  t_x = iAPV->x;
913  t_y = iAPV->y;
914  t_z = iAPV->z;
915  t_Eta = iAPV->Eta;
916  t_R = iAPV->R;
917  t_Phi = iAPV->Phi;
918  t_Thickness = iAPV->Thickness;
919  t_FitMPV = iAPV->FitMPV;
920  t_FitMPVErr = iAPV->FitMPVErr;
921  t_FitWidth = iAPV->FitWidth;
922  t_FitWidthErr = iAPV->FitWidthErr;
923  t_FitChi2NDF = iAPV->FitChi2;
924  t_FitNorm = iAPV->FitNorm;
925  t_FitGrade = iAPV->FitGrade;
926  t_Gain = iAPV->Gain;
927  t_PrevGain = iAPV->PreviousGain;
928  t_PrevGainTick = iAPV->PreviousGainTick;
929  t_NEntries = iAPV->NEntries;
930  t_isMasked = iAPV->isMasked;
931 
932  tree->Fill();
933  }
934  }
935 }
936 
937 //********************************************************************************//
938 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
939 void SiStripGainsPCLHarvester::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
940  edm::ParameterSetDescription desc;
941  desc.setUnknown();
942  descriptions.addDefault(desc);
943 }
944 
945 //********************************************************************************//
946 void SiStripGainsPCLHarvester::endRun(edm::Run const& run, edm::EventSetup const& isetup) {
947  if (!tTopo_) {
948  tTopo_ = std::make_unique<TrackerTopology>(isetup.getData(tTopoToken_));
949  }
950 }
SiStripGainsPCLHarvester::getNewObject
std::unique_ptr< SiStripApvGain > getNewObject(const MonitorElement *Charge_Vs_Index)
Definition: SiStripGainsPCLHarvester.cc:828
StripSubdetector::TEC
static constexpr auto TEC
Definition: StripSubdetector.h:19
StripTopology::nstrips
virtual int nstrips() const =0
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
mps_fire.i
i
Definition: mps_fire.py:428
SiStripGainsPCLHarvester::dqmEndJob
void dqmEndJob(DQMStore::IBooker &ibooker_, DQMStore::IGetter &igetter_) override
Definition: SiStripGainsPCLHarvester.cc:116
B
Definition: APVGainStruct.h:7
ntupleEnum.SubDet
tuple SubDet
Definition: ntupleEnum.py:14
dqm::implementation::NavigatorBase::setCurrentFolder
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:32
PixelTopology::ncolumns
virtual int ncolumns() const =0
DetId::rawId
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
PixelTopology::nrows
virtual int nrows() const =0
APVGain::subdetectorPlane
int subdetectorPlane(uint32_t, const TrackerTopology *)
Definition: APVGainHelpers.cc:119
SiStripGainsPCLHarvester::dqm_tag_
std::vector< std::string > dqm_tag_
Definition: SiStripGainsPCLHarvester.h:95
SiStripGainsPCLHarvester::beginRun
void beginRun(edm::Run const &run, edm::EventSetup const &isetup) override
Definition: SiStripGainsPCLHarvester.cc:69
SiStripGainsPCLHarvester::produceTagFilter
bool produceTagFilter(const MonitorElement *Charge_Vs_Index)
Definition: SiStripGainsPCLHarvester.cc:799
edm::LogError
Log< level::Error, false > LogError
Definition: MessageLogger.h:123
SiStripGainsPCLHarvester::storeGainsTree
void storeGainsTree(const TAxis *chVsIdxXaxis) const
Definition: SiStripGainsPCLHarvester.cc:872
spr::find
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:19
cms::cuda::assert
assert(be >=bs)
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
SiStripGainsPCLHarvester::algoComputeMPVandGain
void algoComputeMPVandGain(const MonitorElement *Charge_Vs_Index)
Definition: SiStripGainsPCLHarvester.cc:489
sistrip::SpyUtilities::range
const uint16_t range(const Frame &aFrame)
Definition: SiStripSpyUtilities.cc:76
dqm::impl::MonitorElement::Fill
void Fill(long long x)
Definition: MonitorElement.h:292
edm::EventSetup::getData
bool getData(T &iHolder) const
Definition: EventSetup.h:122
edm::second
U second(std::pair< T, U > const &p)
Definition: ParameterSet.cc:222
SiStripGainsPCLHarvester::SiStripGainsPCLHarvester
SiStripGainsPCLHarvester(const edm::ParameterSet &ps)
Definition: SiStripGainsPCLHarvester.cc:34
dqm::implementation::IBooker::book1DD
MonitorElement * book1DD(TString const &name, TString const &title, int nchX, double lowX, double highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:155
cscdqm::MASKED
Data available (reporting)
Definition: CSCDQM_Summary.h:48
edm::ConfigurationDescriptions::addDefault
void addDefault(ParameterSetDescription const &psetDescription)
Definition: ConfigurationDescriptions.cc:99
SiStripGainsPCLHarvester::fillDescriptions
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
Definition: SiStripGainsPCLHarvester.cc:939
np
int np
Definition: AMPTWrapper.h:43
SiStripGainsPCLHarvester::qualityToken_
edm::ESGetToken< SiStripQuality, SiStripQualityRcd > qualityToken_
Definition: SiStripGainsPCLHarvester.h:108
gpuVertexFinder::printf
printf("params %d %f %f %f\n", minT, eps, errmax, chi2max)
dqm::implementation::IGetter::get
virtual MonitorElement * get(std::string const &fullpath) const
Definition: DQMStore.cc:673
SiStripGainsPCLHarvester::getPeakOfLandau
void getPeakOfLandau(TH1 *InputHisto, double *FitResults, double LowRange=50, double HighRange=5400, bool gaussianConvolution=false)
Definition: SiStripGainsPCLHarvester.cc:632
edm::Service::isAvailable
bool isAvailable() const
Definition: Service.h:40
DetId::subdetId
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:48
SiStripApvGain::Range
std::pair< ContainerIterator, ContainerIterator > Range
Definition: SiStripApvGain.h:28
cond::service::PoolDBOutputService::writeOneIOV
Hash writeOneIOV(const T &payload, Time_t time, const std::string &recordName)
Definition: PoolDBOutputService.h:62
Exhume::I
const std::complex< double > I
Definition: I.h:8
hcaldqm::constants::GOOD
double const GOOD
Definition: Constants.h:13
StripSubdetector::TOB
static constexpr auto TOB
Definition: StripSubdetector.h:18
SiStripGainsPCLHarvester::VChargeHisto
std::vector< std::string > VChargeHisto
Definition: SiStripGainsPCLHarvester.h:93
SiStripGainsPCLHarvester::tTopo_
std::unique_ptr< TrackerTopology > tTopo_
Definition: SiStripGainsPCLHarvester.h:100
funct::integral
Integral< F, X >::type integral(const F &f)
Definition: Integral.h:70
dqm::legacy::MonitorElement::getTH2S
virtual TH2S * getTH2S() const
Definition: MonitorElement.h:503
SiStripGainsPCLHarvester::tTopoToken_
edm::ESGetToken< TrackerTopology, TrackerTopologyRcd > tTopoToken_
Definition: SiStripGainsPCLHarvester.h:105
edm::LogInfo
Log< level::Info, false > LogInfo
Definition: MessageLogger.h:125
DetId
Definition: DetId.h:17
SiStripGainsPCLHarvester::IsGoodLandauFit
bool IsGoodLandauFit(double *FitResults)
Definition: SiStripGainsPCLHarvester.cc:677
hcaldqm::constants::BAD
double const BAD
Definition: Constants.h:15
StripSubdetector::TIB
static constexpr auto TIB
Definition: StripSubdetector.h:16
PVValHelper::fill
void fill(std::map< std::string, TH1 * > &h, const std::string &s, double x)
Definition: PVValidationHelpers.cc:20
SiStripGainsPCLHarvester::gainQualityMonitor
void gainQualityMonitor(DQMStore::IBooker &ibooker_, const MonitorElement *Charge_Vs_Index) const
Definition: SiStripGainsPCLHarvester.cc:166
APVGain::FetchMonitor
std::vector< MonitorElement * > FetchMonitor(std::vector< APVmon >, uint32_t, const TrackerTopology *topo=nullptr)
Definition: APVGainHelpers.cc:151
APVGain::subdetectorSide
int subdetectorSide(uint32_t, const TrackerTopology *)
Definition: APVGainHelpers.cc:49
SiStripGainsPCLHarvester::bareTkGeomPtr_
const TrackerGeometry * bareTkGeomPtr_
Definition: SiStripGainsPCLHarvester.h:99
SiStripGainsPCLHarvester::endRun
void endRun(edm::Run const &run, edm::EventSetup const &isetup) override
Definition: SiStripGainsPCLHarvester.cc:946
APVGain::thickness
int thickness(uint32_t)
Definition: APVGainHelpers.cc:76
dqm::implementation::IBooker::book2DD
MonitorElement * book2DD(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:338
SiStripGainsPCLHarvester::checkBookAPVColls
virtual void checkBookAPVColls(const edm::EventSetup &setup)
Definition: SiStripGainsPCLHarvester.cc:687
SiStripGainsPCLHarvester::APVsCollOrdered
std::vector< std::shared_ptr< stAPVGain > > APVsCollOrdered
Definition: SiStripGainsPCLHarvester.h:102
a
double a
Definition: hdecay.h:119
SiStripGainsPCLHarvester::gainToken_
edm::ESGetToken< SiStripGain, SiStripGainRcd > gainToken_
Definition: SiStripGainsPCLHarvester.h:107
A
Definition: APVGainStruct.h:7
step
step
Definition: StallMonitor.cc:98
edm::EventSetup::getHandle
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
Definition: EventSetup.h:151
edm::LogWarning
Log< level::Warning, false > LogWarning
Definition: MessageLogger.h:122
APVGain::monHnames
std::vector< std::pair< std::string, std::string > > monHnames(std::vector< std::string >, bool, const char *tag)
Definition: APVGainHelpers.cc:215
h
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
MergeJob_cfg.doStoreOnDB
tuple doStoreOnDB
Definition: MergeJob_cfg.py:65
StripSubdetector::TID
static constexpr auto TID
Definition: StripSubdetector.h:17
APVGain::subdetectorId
int subdetectorId(uint32_t)
Definition: APVGainHelpers.cc:16
edm::Run
Definition: Run.h:45
SiStripGainsPCLHarvester::tkGeomToken_
edm::ESGetToken< TrackerGeometry, TrackerDigiGeometryRecord > tkGeomToken_
Definition: SiStripGainsPCLHarvester.h:106
DetId::DetId
constexpr DetId()
Create an empty or null id (also for persistence)
Definition: DetId.h:38
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:233
SiStripGainsPCLHarvester::APVsColl
std::unordered_map< unsigned int, std::shared_ptr< stAPVGain > > APVsColl
Definition: SiStripGainsPCLHarvester.h:103
Generated for CMSSW Reference Manual by doxygen 1.8.5