Belle II Software  release-05-02-19
PXDDataMCGainCalibrationAlgorithm.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Benjamin Schwenker, Jonas Roetter, Maiko Takahashi *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <pxd/calibration/PXDDataMCGainCalibrationAlgorithm.h>
12 #include <pxd/dbobjects/PXDClusterChargeMapPar.h>
13 #include <pxd/dbobjects/PXDGainMapPar.h>
14 
15 #include <string>
16 #include <algorithm>
17 #include <map>
18 #include <set>
19 
20 #include <sstream>
21 #include <iostream>
22 
23 #include <boost/format.hpp>
24 #include <cmath>
25 
26 
27 //ROOT
28 #include <TRandom.h>
29 #include <TF1.h>
30 #include <TH2I.h>
31 #include <TCanvas.h>
32 
33 using namespace std;
34 using boost::format;
35 using namespace Belle2;
36 
37 
38 // Anonymous namespace for data objects used by PXDDataMCGainCalibrationAlgorithm class
39 namespace {
40 
42  int m_signal;
44  int m_run;
46  int m_exp;
47 
49  void getNumberOfBins(const std::shared_ptr<TH1I>& histo_ptr, unsigned short& nBinsU, unsigned short& nBinsV)
50  {
51  set<unsigned short> uBinSet;
52  set<unsigned short> vBinSet;
53 
54  // Loop over all bins of input histo
55  for (auto histoBin = 1; histoBin <= histo_ptr->GetNbinsX(); histoBin++) {
56  // The bin label contains the vxdid, uBin and vBin
57  string label = histo_ptr->GetXaxis()->GetBinLabel(histoBin);
58 
59  // Parse label string format to read sensorID, uBin and vBin
60  istringstream stream(label);
61  string token;
62  getline(stream, token, '_');
63  getline(stream, token, '_');
64  unsigned short uBin = std::stoi(token);
65 
66  getline(stream, token, '_');
67  unsigned short vBin = std::stoi(token);
68 
69  uBinSet.insert(uBin);
70  vBinSet.insert(vBin);
71  }
72 
73  if (uBinSet.empty() || vBinSet.empty()) {
74  B2FATAL("Not able to determine the grid size. Something is wrong with collected data.");
75  } else {
76  nBinsU = *uBinSet.rbegin() + 1;
77  nBinsV = *vBinSet.rbegin() + 1;
78  }
79  }
80 
82  unsigned short getNumberOfSensors(const std::shared_ptr<TH1I>& histo_ptr)
83  {
84  set<unsigned short> sensorSet;
85 
86  // Loop over all bins of input histo
87  for (auto histoBin = 1; histoBin <= histo_ptr->GetNbinsX(); histoBin++) {
88  // The bin label contains the vxdid, uBin and vBin
89  string label = histo_ptr->GetXaxis()->GetBinLabel(histoBin);
90 
91  // Parse label string format to read sensorID, uBin and vBin
92  istringstream stream(label);
93  string token;
94  getline(stream, token, '_');
95  VxdID sensorID(token);
96  sensorSet.insert(sensorID.getID());
97  }
98 
99  return sensorSet.size();
100  }
101 
102 }
103 
104 
105 PXDDataMCGainCalibrationAlgorithm::PXDDataMCGainCalibrationAlgorithm():
106  CalibrationAlgorithm("PXDClusterChargeCollector"),
107  minClusters(5000), noiseSigma(0.6), safetyFactor(2.0), forceContinue(false), strategy(0),
108  doCalibration(false), useChargeHistogram(false), chargePayloadName("PXDClusterChargeMapPar")
109 {
111  " ------------------------------ PXDDataMCGainCalibrationAlgorithm -----------------------------------\n"
112  " \n"
113  " Algorithm for estimating pxd median/MPV cluster charges for different position on sensor in ADU, and optionally calculate gain \n"
114  " ------------------------------------------------------------------------------------------------\n"
115  );
116 
117 }
118 
119 
121 {
122 
123  // Get counter histogram
124  auto cluster_counter = getObjectPtr<TH1I>("PXDClusterCounter");
125  if (!cluster_counter) {
126  B2INFO("Not enough Data: cluster counter does not exist ");
127  return c_NotEnoughData;
128  }
129 
130  // Extract number of sensors from counter histograms
131  auto nSensors = getNumberOfSensors(cluster_counter);
132 
133  // Extract the number of grid bins from counter histograms
134  unsigned short nBinsU = 0;
135  unsigned short nBinsV = 0;
136  getNumberOfBins(cluster_counter, nBinsU, nBinsV);
137 
138  // Check that we have collected enough Data
139  if (cluster_counter->GetEntries() < int(safetyFactor * minClusters * nSensors * nBinsU * nBinsV)) {
140  if (not forceContinue) {
141  B2INFO("Not enough Data: Only " << cluster_counter->GetEntries() << " hits were collected but " << int(safetyFactor * minClusters *
142  nSensors * nBinsU * nBinsV) << " needed!");
143  return c_NotEnoughData;
144  } else {
145  B2INFO("Continue despite low statistics: Only " << cluster_counter->GetEntries() << " hits were collected but" << int(
147  nSensors * nBinsU * nBinsV) << " would be desirable!");
148  }
149  }
150 
151  B2INFO("Start calibration using a " << nBinsU << "x" << nBinsV << " grid per sensor.");
152 
153  // This is the PXD charge calibration payload for conditions DB
154  PXDClusterChargeMapPar* chargeMapPar = new PXDClusterChargeMapPar(nBinsU, nBinsV);
155  // This is the PXD gain correction payload for conditions DB
156  PXDGainMapPar* gainMapPar = new PXDGainMapPar(nBinsU, nBinsV);
157  set<VxdID> pxdSensors;
158 
159  // Read back existing DB payloads
160  PXDClusterChargeMapPar* chargeMapPtr = nullptr;
161  PXDGainMapPar* gainMapPtr = nullptr;
162  auto dbtree = getObjectPtr<TTree>("dbtree");
163  dbtree->SetBranchAddress("run", &m_run);
164  dbtree->SetBranchAddress("exp", &m_exp);
165  dbtree->SetBranchAddress("chargeMap", &chargeMapPtr);
166  dbtree->SetBranchAddress("gainMap", &gainMapPtr);
167 
168  // Loop over all bins of input histo
169  for (auto histoBin = 1; histoBin <= cluster_counter->GetNbinsX(); histoBin++) {
170 
171  // The bin label contains the vxdid, uBin and vBin
172  string label = cluster_counter->GetXaxis()->GetBinLabel(histoBin);
173 
174  // Parse label string format to read sensorID, uBin and vBin
175  istringstream stream(label);
176  string token;
177  getline(stream, token, '_');
178  VxdID sensorID(token);
179 
180  getline(stream, token, '_');
181  unsigned short uBin = std::stoi(token);
182 
183  getline(stream, token, '_');
184  unsigned short vBin = std::stoi(token);
185 
186  // Read back the counters for number of collected clusters
187  int numberOfDataHits = cluster_counter->GetBinContent(histoBin);
188 
189  // Only perform fitting, when enough data is available
190  if (numberOfDataHits >= minClusters) {
191 
192  B2INFO("start EstimateCharge");
193  // Compute the cluster charge or gain for the part of PXD
194  auto Charge = EstimateCharge(sensorID, uBin, vBin, histoBin);
195 
196  // Store the charge or gain
197  if (!doCalibration) {
198 
199  B2INFO("EstimateCharge: sensor " << sensorID.getID() << " U " << uBin << " V " << vBin << " Charge " << Charge);
200  chargeMapPar->setContent(sensorID.getID(), uBin, vBin, Charge);
201 
202  } else {
203 
204  auto Gain = 0.0;
205 
206  if (Charge <= 0.0) {
207  B2WARNING("Retrieved negative charge for data for sensor=" << sensorID << " uBin=" << uBin << " vBin=" << vBin <<
208  ". Set gain to default value (=1.0).");
209  Gain = 1.0;
210  } else {
211  dbtree->GetEntry(0);
212  double mcCharge = chargeMapPtr->getContent(sensorID.getID(), uBin, vBin);
213  //GetChargeFromDB(sensorID, uBin, vBin, dbtree);
214  if (mcCharge <= 0.0) {
215  B2WARNING("Retrieved negative charge for MC from DB for sensor=" << sensorID << " uBin=" << uBin << " vBin=" << vBin <<
216  ". Set gain to default value (=1.0).");
217  Gain = 1.0;
218  } else {
219  Gain = Charge / mcCharge;
220  B2INFO("Estimated Gain: sensor " << sensorID.getID() << " U " << uBin << " V " << vBin << " Gain " << Gain << " = " << Charge <<
221  " / " << mcCharge);
222  }
223 
224  }
225 
226  gainMapPar->setContent(sensorID.getID(), uBin, vBin, Gain);
227 
228  }
229 
230  } else {
231 
232  B2WARNING(label << ": Number of data hits too small for fitting (" << numberOfDataHits << " < " << minClusters <<
233  "). Use default value of 0 for charge, 1 for gain.");
234  if (!doCalibration) chargeMapPar->setContent(sensorID.getID(), uBin, vBin, 0.0);
235  else gainMapPar->setContent(sensorID.getID(), uBin, vBin, 1.0);
236 
237  }
238 
239  pxdSensors.insert(sensorID);
240 
241  }
242 
243  chargeMapPtr = nullptr;
244  gainMapPtr = nullptr;
245 
246  // Save the charge map to database. Note that this will set the database object name to the same as the collector but you
247  // are free to change it.
248  if (!doCalibration) {
249  saveCalibration(chargeMapPar, chargePayloadName);
250 
251  B2INFO("PXD Cluster Charge Calibration Successful");
252  return c_OK;
253 
254  } else {
255 
256  // Post processing of gain map. It is possible that the gain
257  // computation failed on some parts. Here, we replace default
258  // values (1.0) by local averages of neighboring sensor parts.
259 
260  for (const auto& sensorID : pxdSensors) {
261  for (unsigned short vBin = 0; vBin < nBinsV; ++vBin) {
262  float meanGain = 0;
263  unsigned short nGood = 0;
264  unsigned short nBad = 0;
265  for (unsigned short uBin = 0; uBin < nBinsU; ++uBin) {
266  auto gain = gainMapPar->getContent(sensorID.getID(), uBin, vBin);
267  // Filter default gains
268  if (gain != 1.0) {
269  nGood += 1;
270  meanGain += gain;
271  } else {
272  nBad += 1;
273  }
274  }
275  B2RESULT("Gain calibration on sensor=" << sensorID << " and vBin=" << vBin << " was successful on " << nGood << "/" << nBinsU <<
276  " uBins.");
277 
278  // Check if we can repair bad calibrations with a local avarage
279  if (nGood > 0 && nBad > 0) {
280  meanGain /= nGood;
281  for (unsigned short uBin = 0; uBin < nBinsU; ++uBin) {
282  auto gain = gainMapPar->getContent(sensorID.getID(), uBin, vBin);
283  if (gain == 1.0) {
284  gainMapPar->setContent(sensorID.getID(), uBin, vBin, meanGain);
285  B2RESULT("Gain calibration on sensor=" << sensorID << ", vBin=" << vBin << " uBin " << uBin << ": Replace default gain wih average "
286  << meanGain);
287  }
288  }
289  }
290  }
291  }
292 
293  // Save the gain map to database. Note that this will set the database object name to the same as the collector but you
294  // are free to change it.
295  saveCalibration(gainMapPar, "PXDGainMapPar");
296 
297  B2INFO("PXD Gain Calibration Successful");
298  return c_OK;
299 
300  }
301 
302 }
303 
304 
305 double PXDDataMCGainCalibrationAlgorithm::EstimateCharge(VxdID sensorID, unsigned short uBin, unsigned short vBin,
306  unsigned short histoBin)
307 {
308 
309  // Construct a tree name for requested part of PXD
310  auto layerNumber = sensorID.getLayerNumber();
311  auto ladderNumber = sensorID.getLadderNumber();
312  auto sensorNumber = sensorID.getSensorNumber();
313 
314  if (!useChargeHistogram) {
315  const string treename = str(format("tree_%1%_%2%_%3%_%4%_%5%") % layerNumber % ladderNumber % sensorNumber % uBin % vBin);
316  // Vector with cluster signals from collected data
317  vector<double> signals;
318  // Fill data_signal vector from input data
319  auto tree = getObjectPtr<TTree>(treename);
320  tree->SetBranchAddress("signal", &m_signal);
321 
322  // Loop over tree
323  const auto nEntries = tree->GetEntries();
324  int incr(1);
325  if (int(nEntries / minClusters) > 2) incr = int(2 * nEntries / minClusters);
326  for (int i = 0; i < nEntries; i += incr) {
327  tree->GetEntry(i);
328  double noise = gRandom->Gaus(0.0, noiseSigma);
329  signals.push_back(m_signal + noise); //qyliu: why we introduce noise simulation here?
330  }
331 
332 
333  if (strategy == 0) {
334  double median = CalculateMedian(signals);
335  B2INFO("EstimateCharge: sensor " << sensorID.getID() << "(" << layerNumber << "," << ladderNumber << "," << sensorNumber
336  << ") U " << uBin << " V " << vBin
337  << " Charge " << median);
338  return median;
339  } else if (strategy == 1) {
340  double median = CalculateMedian(signals);
341  double landaumpv = FitLandau(signals);
342  double diff = (landaumpv - median);
343  double difff = 0.;
344  if (landaumpv > 0.) difff = (landaumpv - median) / landaumpv;
345  B2INFO("EstimateCharge: sensor " << sensorID.getID() << "(" << layerNumber << "," << ladderNumber << "," << sensorNumber
346  << ") U " << uBin << " V " << vBin
347  << " Charge " << landaumpv << " Median " << median << " diff = " << diff << "/" << difff);
348  return landaumpv; //FitLandau(signals);
349  } else if (strategy == 2) {
350  double mean = 0;
351  for (auto& each : signals)
352  mean += each;
353  if (signals.size() > 0) mean = mean / signals.size();
354  return mean;
355  } else {
356  B2FATAL("strategy unavailable, use 0 for medians, 1 for landau fit and 2 for mean!");
357  }
358 
359  } else {
360 
361  auto cluster_counter = getObjectPtr<TH2I>("PXDClusterCharge");
362  TH1D* hist_signals = cluster_counter->ProjectionY("proj", histoBin, histoBin);
363 
364  if (strategy == 0) {
365  double median = CalculateMedian(hist_signals);
366  B2INFO("EstimateCharge: sensor " << sensorID.getID() << "(" << layerNumber << "," << ladderNumber << "," << sensorNumber
367  << ") U " << uBin << " V " << vBin
368  << " Charge " << median);
369  return median;
370  }
371  if (strategy == 1) {
372  double median = CalculateMedian(hist_signals);
373  double landaumpv = FitLandau(hist_signals);
374  double diff = (landaumpv - median);
375  double difff = 0.;
376  if (landaumpv > 0.) difff = (landaumpv - median) / landaumpv;
377  B2INFO("EstimateCharge: sensor " << sensorID.getID() << "(" << layerNumber << "," << ladderNumber << "," << sensorNumber
378  << ") U " << uBin << " V " << vBin
379  << " Charge " << landaumpv << " Median " << median << " diff = " << diff << "/" << difff);
380  return landaumpv; //FitLandau(signals);
381  } else if (strategy == 2) {
382  return hist_signals->GetMean();
383  } else {
384  B2FATAL("strategy unavailable, use 0 for medians, 1 for landau fit and 2 for mean!");
385  }
386 
387  delete hist_signals;
388 
389  }
390 
391 }
392 
394 {
395  auto size = signals.size();
396 
397  if (size == 0) {
398  return 0.0; // Undefined, really.
399  } else {
400  sort(signals.begin(), signals.end());
401  if (size % 2 == 0) {
402  return (signals[size / 2 - 1] + signals[size / 2]) / 2;
403  } else {
404  return signals[size / 2];
405  }
406  }
407 }
408 
410 {
411  auto size = hist_signals->GetEntries();
412 
413  if (size == 0) return 0.0; // Undefined.
414 
415  int sum = 0;
416  for (int ibin = 0; ibin < hist_signals->GetNbinsX(); ++ibin) {
417  sum += hist_signals->GetBinContent(ibin + 1);
418  if (sum > size / 2) {
419  return hist_signals->GetBinLowEdge(ibin + 1);
420  }
421  }
422 
423  B2WARNING("Could not find median! using default value 0.0!");
424  return 0.0;
425 }
426 
427 double PXDDataMCGainCalibrationAlgorithm::FitLandau(vector<double>& signals)
428 {
429  auto size = signals.size();
430  if (size == 0) return 0.0; // Undefined, really.
431 
432  // get max and min values of vector
433  int max = *max_element(signals.begin(), signals.end());
434  int min = *min_element(signals.begin(), signals.end());
435  // make even bins
436  if ((max - min) % 2) max--;
437  int nbin = max - min;
438 
439  // create histogram to hold signals and fill it
440  TH1D* hist_signals = new TH1D("", "", nbin, min, max);
441  for (auto it = signals.begin(); it != signals.end(); ++it) {
442  hist_signals->Fill(*it);
443  }
444 
445  // create fit function
446  TF1* landau = new TF1("landau", "TMath::Landau(x,[0],[1])*[2]", min, max);
447  landau->SetParNames("MPV", "sigma", "scale");
448  landau->SetParameters(35, 8, 1000);
449  landau->SetParLimits(0, 0., 80.);
450 
451  //do fit and get results, fit range restricted to exclude low charge peak
452  float fitmin(min);
453  float fitmax(350);
454  Int_t status = hist_signals->Fit("landau", "Lq", "", fitmin, fitmax);
455  double MPV = landau->GetParameter("MPV");
456 
457  B2INFO("Fit result: " << status << " MPV " << MPV << " sigma " << landau->GetParameter("sigma")
458  << " scale " << landau->GetParameter("scale") << " chi2 " << landau->GetChisquare());
459 
460  // clean up
461  delete hist_signals;
462  delete landau;
463 
464  // check fit status
465  if (status == 0) return MPV;
466  else {
467  B2WARNING("Fit failed! using default value 0.0!");
468  return 0.0;
469  }
470 }
471 
473 {
474  auto size = hist_signals->GetEntries();
475  if (size == 0) return 0.0; // Undefined.
476 
477  int max = hist_signals->GetBinLowEdge(hist_signals->GetNbinsX() + 1);
478  int min = hist_signals->GetBinLowEdge(1);
479 
480  // create fit function
481  TF1* landau = new TF1("landau", "TMath::Landau(x,[0],[1])*[2]", min, max);
482  landau->SetParNames("MPV", "sigma", "scale");
483  landau->SetParameters(35, 8, 1000);
484  landau->SetParLimits(0, 0., 80.);
485 
486  // do fit and get results, fit range restricted to exclude low charge peak
487  float fitmin(min);
488  float fitmax(250);
489  Int_t status = hist_signals->Fit("landau", "Lq", "", fitmin, fitmax);
490  double MPV = landau->GetParameter("MPV");
491 
492  B2INFO("Fit result: " << status << " MPV " << MPV << " sigma " << landau->GetParameter("sigma")
493  << " scale " << landau->GetParameter("scale") << " chi2 " << landau->GetChisquare());
494 
495  // clean up
496  delete landau;
497 
498  // check fit status
499  if (status == 0) return MPV;
500  else {
501  B2WARNING("Fit failed! using default value 0.0!");
502  return 0.0;
503  }
504 }
505 
Belle2::PXDDataMCGainCalibrationAlgorithm::safetyFactor
float safetyFactor
Safety factor for determining whether the collected number of clusters is enough.
Definition: PXDDataMCGainCalibrationAlgorithm.h:48
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::CalibrationAlgorithm::saveCalibration
void saveCalibration(TClonesArray *data, const std::string &name)
Store DBArray payload with given name with default IOV.
Definition: CalibrationAlgorithm.cc:290
Belle2::PXDDataMCGainCalibrationAlgorithm::EstimateCharge
double EstimateCharge(VxdID sensorID, unsigned short uBin, unsigned short vBin, unsigned short histoBin)
Estimate median charge form collected clusters on part of PXD.
Definition: PXDDataMCGainCalibrationAlgorithm.cc:305
Belle2::VxdID::getLadderNumber
baseType getLadderNumber() const
Get the ladder id.
Definition: VxdID.h:108
Belle2::PXDDataMCGainCalibrationAlgorithm::CalculateMedian
double CalculateMedian(std::vector< double > &signals)
Calculate a median from unsorted signal vector. The input vector gets sorted.
Definition: PXDDataMCGainCalibrationAlgorithm.cc:393
Belle2::PXDDataMCGainCalibrationAlgorithm::calibrate
virtual EResult calibrate() override
Run algo on data.
Definition: PXDDataMCGainCalibrationAlgorithm.cc:120
Belle2::VxdID::getID
baseType getID() const
Get the unique id.
Definition: VxdID.h:104
Belle2::PXDDataMCGainCalibrationAlgorithm::chargePayloadName
std::string chargePayloadName
Payload name for Cluster Charge.
Definition: PXDDataMCGainCalibrationAlgorithm.h:65
Belle2::PXDDataMCGainCalibrationAlgorithm::forceContinue
bool forceContinue
Force continue in low statistics runs instead of returning c_NotEnoughData.
Definition: PXDDataMCGainCalibrationAlgorithm.h:51
Belle2::CalibrationAlgorithm::c_OK
@ c_OK
Finished successfuly =0 in Python.
Definition: CalibrationAlgorithm.h:51
Belle2::PXDGainMapPar
The payload class for PXD gain corrections.
Definition: PXDGainMapPar.h:53
Belle2::CalibrationAlgorithm::setDescription
void setDescription(const std::string &description)
Set algorithm description (in constructor)
Definition: CalibrationAlgorithm.h:331
Belle2::PXDClusterChargeMapPar::setContent
void setContent(unsigned short sensorID, unsigned short globalID, float value)
Set map content.
Definition: PXDClusterChargeMapPar.h:65
Belle2::PXD::getNumberOfBins
void getNumberOfBins(const std::shared_ptr< TH1I > &histo_ptr, unsigned short &nBinsU, unsigned short &nBinsV)
Helper function to extract number of bins along u side and v side from counter histogram labels.
Definition: PXDCalibrationUtilities.cc:37
Belle2::PXDClusterChargeMapPar::getContent
float getContent(unsigned short sensorID, unsigned short globalID) const
Get content.
Definition: PXDClusterChargeMapPar.h:100
Belle2::PXDDataMCGainCalibrationAlgorithm::doCalibration
bool doCalibration
flag to perform full calibration or only esitmate charge: False: only estimate charge,...
Definition: PXDDataMCGainCalibrationAlgorithm.h:59
Belle2::PXDGainMapPar::getContent
float getContent(unsigned short sensorID, unsigned short globalID) const
Get content.
Definition: PXDGainMapPar.h:113
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::PXDDataMCGainCalibrationAlgorithm::strategy
int strategy
strategy to used for gain calibration, 0 for medians, 1 for landau fit
Definition: PXDDataMCGainCalibrationAlgorithm.h:54
Belle2::VxdID::getSensorNumber
baseType getSensorNumber() const
Get the sensor id.
Definition: VxdID.h:110
Belle2::PXDGainMapPar::setContent
void setContent(unsigned short sensorID, unsigned short globalID, float value)
Set map content.
Definition: PXDGainMapPar.h:78
Belle2::CalibrationAlgorithm::EResult
EResult
The result of calibration.
Definition: CalibrationAlgorithm.h:50
Belle2::PXDClusterChargeMapPar
The payload class for PXD cluster charge calibrations.
Definition: PXDClusterChargeMapPar.h:40
Belle2::PXDDataMCGainCalibrationAlgorithm::minClusters
int minClusters
Minimum number of collected clusters for estimating median charge.
Definition: PXDDataMCGainCalibrationAlgorithm.h:42
Belle2::PXDDataMCGainCalibrationAlgorithm::noiseSigma
float noiseSigma
Artificial noise sigma for smearing cluster charge.
Definition: PXDDataMCGainCalibrationAlgorithm.h:45
Belle2::CalibrationAlgorithm::c_NotEnoughData
@ c_NotEnoughData
Needs more data =2 in Python.
Definition: CalibrationAlgorithm.h:53
Belle2::PXDDataMCGainCalibrationAlgorithm::FitLandau
double FitLandau(std::vector< double > &signals)
calculate MPV of unsorted signal vector using a Landau fit
Definition: PXDDataMCGainCalibrationAlgorithm.cc:427
Belle2::CalibrationAlgorithm
Base class for calibration algorithms.
Definition: CalibrationAlgorithm.h:47
Belle2::VxdID::getLayerNumber
baseType getLayerNumber() const
Get the layer id.
Definition: VxdID.h:106
Belle2::PXD::getNumberOfSensors
unsigned short getNumberOfSensors(const std::shared_ptr< TH1I > &histo_ptr)
Helper function to extract number of sensors from counter histogram labels.
Definition: PXDCalibrationUtilities.cc:70
Belle2::PXDDataMCGainCalibrationAlgorithm::useChargeHistogram
bool useChargeHistogram
Flag to use histogram as charge input.
Definition: PXDDataMCGainCalibrationAlgorithm.h:62