Belle II Software  release-05-02-19
eclMuMuEAlgorithm.cc
1 
2 #include <ecl/calibration/eclMuMuEAlgorithm.h>
3 #include <ecl/dbobjects/ECLCrystalCalib.h>
4 
5 #include "TH2F.h"
6 #include "TFile.h"
7 #include "TMath.h"
8 #include "TF1.h"
9 #include "TROOT.h"
10 
11 using namespace std;
12 using namespace Belle2;
13 using namespace ECL;
14 
16 double eclNovoConst(double* x, double* par)
17 {
18  double qc = 0;
19 
20  double peak = par[1];
21  double width = par[2];
22  double sln4 = sqrt(log(4));
23  double y = par[3] * sln4;
24  double tail = -log(y + sqrt(1 + y * y)) / sln4;
25 
26  if (TMath::Abs(tail) < 1.e-7) {
27  qc = 0.5 * TMath::Power(((x[0] - peak) / width), 2);
28  } else {
29  double qa = tail * sqrt(log(4.));
30  double qb = sinh(qa) / qa;
31  double qx = (x[0] - peak) / width * qb;
32  double qy = 1. + tail * qx;
33 
34  if (qy > 1.E-7)
35  qc = 0.5 * (TMath::Power((log(qy) / tail), 2) + tail * tail);
36  else
37  qc = 15.0;
38  }
39  return par[0] * exp(-qc) + par[4];
40 }
41 
42 
43 eclMuMuEAlgorithm::eclMuMuEAlgorithm(): CalibrationAlgorithm("eclMuMuECollector"), cellIDLo(1), cellIDHi(8736), minEntries(150),
44  maxIterations(10), tRatioMin(0.2), tRatioMax(0.25), performFits(true), findExpValues(false), storeConst(0)
45 {
47  "Perform energy calibration of ecl crystals by fitting a Novosibirsk function to energy deposited by muons"
48  );
49 }
50 
52 {
55  double limitTol = 0.0005; /*< tolerance for checking if a parameter is at the limit */
56  double minFitLimit = 1e-25; /*< cut off for labeling a fit as poor */
57  double minFitProbIter = 1e-8; /*< cut off for labeling a fit as poor if it also has many iterations */
58  double constRatio = 0.5; /*< Novosibirsk normalization must be greater than constRatio x constant term */
59  double peakMin(0.4), peakMax(2.2); /*< range for peak of normalized energy distribution */
60  double peakTol = limitTol * (peakMax - peakMin); /*< fit is at limit if it is within peakTol of min or max */
61  double effSigMin(0.02), effSigMax(0.2); /*< range for effective sigma of normalized energy distribution */
62  double effSigTol = limitTol * (effSigMax - effSigMin); /*< fit is at limit if it is within effSigTol of min or max */
63  double etaNom(-0.41); /*< Nominal tail parameter; fixed to this value for endcap and when findExValues=false */
64  double etaMin(-0.7), etaMax(-0.2); /*< Novosibirsk tail parameter range */
65  double etaTol = limitTol * (etaMax - etaMin); /*< fit is at limit if it is within etaTol of min or max */
66  double constMin(0.), constMax(10.); /*< constant term in normalized energy distribution */
67  double constTol = limitTol * constMax; /*< if constant is less than constTol, it will be fixed to 0 */
68 
70  gROOT->SetBatch();
71 
74  TH1F* dummy;
75  dummy = (TH1F*)gROOT->FindObject("IntegralVsCrysID");
76  if (dummy) {delete dummy;}
77  dummy = (TH1F*)gROOT->FindObject("AverageExpECrys");
78  if (dummy) {delete dummy;}
79  dummy = (TH1F*)gROOT->FindObject("AverageElecCalib");
80  if (dummy) {delete dummy;}
81  dummy = (TH1F*)gROOT->FindObject("AverageInitCalib");
82  if (dummy) {delete dummy;}
83 
86  auto TrkPerCrysID = getObjectPtr<TH1F>("TrkPerCrysID");
87  auto EnVsCrysID = getObjectPtr<TH2F>("EnVsCrysID");
88  auto ExpEvsCrys = getObjectPtr<TH1F>("ExpEvsCrys");
89  auto ElecCalibvsCrys = getObjectPtr<TH1F>("ElecCalibvsCrys");
90  auto InitialCalibvsCrys = getObjectPtr<TH1F>("InitialCalibvsCrys");
91  auto CalibEntriesvsCrys = getObjectPtr<TH1F>("CalibEntriesvsCrys");
92  auto RawDigitAmpvsCrys = getObjectPtr<TH2F>("RawDigitAmpvsCrys");
93  auto RawDigitTimevsCrys = getObjectPtr<TH2F>("RawDigitTimevsCrys");
94  auto hitCrysVsExtrapolatedCrys = getObjectPtr<TH2F>("hitCrysVsExtrapolatedCrys");
95 
99  TH1F* IntegralVsCrysID = new TH1F("IntegralVsCrysID", "Integral of EnVsCrysID for each crystal;crystal ID;Entries", 8736, 0, 8736);
100  TH1F* AverageExpECrys = new TH1F("AverageExpECrys", "Average expected E per crys from collector;Crystal ID;Energy (GeV)", 8736, 0,
101  8736);
102  TH1F* AverageElecCalib = new TH1F("AverageElecCalib", "Average electronics calib const vs crystal;Crystal ID;Calibration constant",
103  8736, 0, 8736);
104  TH1F* AverageInitCalib = new TH1F("AverageInitCalib", "Average initial calib const vs crystal;Crystal ID;Calibration constant",
105  8736, 0, 8736);
106 
107  for (int crysID = 0; crysID < 8736; crysID++) {
108  TH1D* hEnergy = EnVsCrysID->ProjectionY("hEnergy", crysID + 1, crysID + 1);
109  int Integral = hEnergy->Integral();
110  IntegralVsCrysID->SetBinContent(crysID + 1, Integral);
111 
112  double TotEntries = CalibEntriesvsCrys->GetBinContent(crysID + 1);
113 
114  double expectedE = 0.;
115  if (TotEntries > 0.) {expectedE = ExpEvsCrys->GetBinContent(crysID + 1) / TotEntries;}
116  AverageExpECrys->SetBinContent(crysID + 1, expectedE);
117 
118  double calibconst = 0.;
119  if (TotEntries > 0.) {calibconst = ElecCalibvsCrys->GetBinContent(crysID + 1) / TotEntries;}
120  AverageElecCalib->SetBinContent(crysID + 1, calibconst);
121 
122  calibconst = 0.;
123  if (TotEntries > 0.) {calibconst = InitialCalibvsCrys->GetBinContent(crysID + 1) / TotEntries;}
124  AverageInitCalib->SetBinContent(crysID + 1, calibconst);
125  }
126 
129  TFile* histfile = new TFile("eclMuMuEAlgorithm.root", "recreate");
130  TrkPerCrysID->Write();
131  EnVsCrysID->Write();
132  IntegralVsCrysID->Write();
133  AverageExpECrys->Write();
134  AverageElecCalib->Write();
135  AverageInitCalib->Write();
136  RawDigitAmpvsCrys->Write();
137  RawDigitTimevsCrys->Write();
138  hitCrysVsExtrapolatedCrys->Write();
139 
142  if (!performFits) {
143  B2RESULT("eclMuMuEAlgorithm has not been asked to perform fits; copying input histograms and quitting");
144  histfile->Close();
145  return c_NotEnoughData;
146  }
147 
150  bool sufficientData = true;
151  for (int crysID = cellIDLo - 1; crysID < cellIDHi; crysID++) {
152  if (IntegralVsCrysID->GetBinContent(crysID + 1) < minEntries) {
153  if (storeConst == 1) {B2RESULT("eclMuMuEAlgorithm: crystal " << crysID << " has insufficient statistics: " << IntegralVsCrysID->GetBinContent(crysID + 1) << ". Requirement is " << minEntries);}
154  sufficientData = false;
155  break;
156  }
157  }
158 
160  if (!sufficientData && storeConst == 1) {
161  histfile->Close();
162  return c_NotEnoughData;
163  }
164 
169  TH1F* CalibVsCrysID = new TH1F("CalibVsCrysID", "Calibration constant vs crystal ID;crystal ID;counts per GeV", 8736, 0, 8736);
170  TH1F* ExpEnergyperCrys = new TH1F("ExpEnergyperCrys", "Expected energy per crystal;Crystal ID;Peak energy (GeV)", 8736, 0, 8736);
171 
173  TH1F* PeakVsCrysID = new TH1F("PeakVsCrysID", "Peak of Novo fit vs crystal ID;crystal ID;Peak normalized energy", 8736, 0,
174  8736);
175  TH1F* effSigVsCrysID = new TH1F("effSigVsCrysID", "effSigma vs crystal ID;crystal ID;sigma)", 8736, 0, 8736);
176  TH1F* etaVsCrysID = new TH1F("etaVsCrysID", "eta vs crystal ID;crystal ID;Novo eta parameter", 8736, 0, 8736);
177  TH1F* constVsCrysID = new TH1F("constVsCrysID", "fit constant vs crystal ID;crystal ID;fit constant", 8736, 0, 8736);
178  TH1F* normVsCrysID = new TH1F("normVsCrysID", "Novosibirsk normalization vs crystal ID;crystal ID;normalization", 8736, 0, 8736);
179  TH1F* fitLimitVsCrysID = new TH1F("fitLimitVsCrysID", "fit range upper limit vs crystal ID;crystal ID;upper fit limit", 8736, 0,
180  8736);
181  TH1F* StatusVsCrysID = new TH1F("StatusVsCrysID", "Fit status vs crystal ID;crystal ID;Fit status", 8736, 0, 8736);
182 
184  TH1F* hStatus = new TH1F("hStatus", "Fit status", 25, -5, 20);
185  TH1F* hPeak = new TH1F("hPeak", "Peaks of normalized energy distributions, successful fits;Peak of Novosibirsk fit", 200, 0.8, 1.2);
186  TH1F* fracPeakUnc = new TH1F("fracPeakUnc", "Fractional uncertainty on peak uncertainty, successful fits;Uncertainty on peak", 100,
187  0, 0.1);
188  TH1F* nIterations = new TH1F("nIterations", "Number of times histogram was fit;Number of iterations", 20, -0.5, 19.5);
189 
190 
193  bool allFitsOK = true;
194  for (int crysID = cellIDLo - 1; crysID < cellIDHi; crysID++) {
195 
197  TString name = "Enormalized";
198  name += crysID;
199  TH1D* hEnergy = EnVsCrysID->ProjectionY(name, crysID + 1, crysID + 1);
200 
202  double histMin = hEnergy->GetXaxis()->GetXmin();
203  double histMax = hEnergy->GetXaxis()->GetXmax();
204  TF1* func = new TF1("eclNovoConst", eclNovoConst, histMin, histMax, 5);
205  func->SetParNames("normalization", "peak", "effSigma", "eta", "const");
206  func->SetParLimits(1, peakMin, peakMax);
207  func->SetParLimits(2, effSigMin, effSigMax);
208  func->SetParLimits(3, etaMin, etaMax);
209  func->SetParLimits(4, constMin, constMax);
210 
212  hEnergy->GetXaxis()->SetRangeUser(peakMin, peakMax);
213  int maxBin = hEnergy->GetMaximumBin();
214  double peakE = hEnergy->GetBinLowEdge(maxBin);
215  double peakEUnc = 0.;
216  double normalization = hEnergy->GetMaximum();
217  double normUnc = 0.;
218  double effSigma = hEnergy->GetRMS();
219  double sigmaUnc = 0.;
220  hEnergy->GetXaxis()->SetRangeUser(histMin, histMax);
221 
223  double fitlow = 0.25;
224  double fithigh = peakE + 2.5 * effSigma;
225 
227  double eta = etaNom;
228  double etaUnc = 0.;
229 
231  int il0 = hEnergy->GetXaxis()->FindBin(fitlow);
232  int il1 = hEnergy->GetXaxis()->FindBin(fitlow + 0.1);
233  double constant = hEnergy->Integral(il0, il1) / (1 + il1 - il0);
234  if (constant < 0.01 * normalization) {constant = 0.01 * normalization;}
235  double constUnc = 0.;
236 
238  double dIter = 0.1 * (histMax - histMin) / hEnergy->GetNbinsX();
239  double fitProb(0.);
240  double highold(0.), higholdold(0.);
241  bool fixConst = false;
242  int nIter = 0;
243  bool fitHist = IntegralVsCrysID->GetBinContent(crysID + 1) >= minEntries; /* fit only if enough events */
244 
247  while (fitHist) {
248 
250  func->SetParameters(normalization, peakE, effSigma, eta, constant);
251  if (fixConst) { func->FixParameter(4, 0); }
252 
254  if (crysID < 1152 || crysID > 7775 || !findExpValues) {func->FixParameter(3, etaNom);}
255 
257  hEnergy->Fit(func, "LIQ", "", fitlow, fithigh);
258  nIter++;
259  fitHist = false;
260  normalization = func->GetParameter(0);
261  normUnc = func->GetParError(0);
262  peakE = func->GetParameter(1);
263  peakEUnc = func->GetParError(1);
264  effSigma = func->GetParameter(2);
265  sigmaUnc = func->GetParError(2);
266  eta = func->GetParameter(3);
267  etaUnc = func->GetParError(3);
268  constant = func->GetParameter(4);
269  constUnc = func->GetParError(4);
270  fitProb = func->GetProb();
271 
273  double peak = func->Eval(peakE) - constant;
274  double tRatio = (func->Eval(fithigh) - constant) / peak;
275  if (tRatio < tRatioMin || tRatio > tRatioMax) {
276  double targetY = constant + 0.5 * (tRatioMin + tRatioMax) * peak;
277  higholdold = highold;
278  highold = fithigh;
279  fithigh = func->GetX(targetY, peakE, histMax);
280  fitHist = true;
281 
283  if (abs(fithigh - higholdold) < dIter) {fithigh = 0.5 * (highold + higholdold); }
284 
286  if (nIter > maxIterations - 3) {fithigh = 0.33333 * (fithigh + highold + higholdold); }
287  }
288 
290  if (constant < constTol && !fixConst) {
291  constant = 0;
292  fixConst = true;
293  fitHist = true;
294  }
295 
297  if (nIter == maxIterations) {fitHist = false;}
298  B2DEBUG(200, crysID << " " << nIter << " " << peakE << " " << constant << " " << tRatio << " " << fithigh);
299  }
300 
303  int iStatus = fitOK; // success
304  if (nIter == maxIterations) {iStatus = iterations;} // too many iterations
305 
307  if (normalization < constRatio * constant) {iStatus = noPeak;}
308 
310  if (fitProb <= minFitLimit || (fitProb < minFitProbIter && iStatus == iterations)) {iStatus = poorFit;}
311 
313  if ((peakE < peakMin + peakTol) || (peakE > peakMax - peakTol)) {iStatus = atLimit;}
314  if ((effSigma < effSigMin + effSigTol) || (effSigma > effSigMax - effSigTol)) {iStatus = atLimit;}
315  if ((eta < etaMin + etaTol) || (eta > etaMax - etaTol)) {iStatus = atLimit;}
316  if (constant > constMax - constTol) {iStatus = atLimit;}
317 
318  //** No fit
319  if (nIter == 0) {iStatus = notFit;} // not fit
320 
322  int histbin = crysID + 1;
323  PeakVsCrysID->SetBinContent(histbin, peakE);
324  PeakVsCrysID->SetBinError(histbin, peakEUnc);
325  effSigVsCrysID->SetBinContent(histbin, effSigma);
326  effSigVsCrysID->SetBinError(histbin, sigmaUnc);
327  etaVsCrysID->SetBinContent(histbin, eta);
328  etaVsCrysID->SetBinError(histbin, etaUnc);
329  constVsCrysID->SetBinContent(histbin, constant);
330  constVsCrysID->SetBinError(histbin, constUnc);
331  normVsCrysID->SetBinContent(histbin, normalization);
332  normVsCrysID->SetBinError(histbin, normUnc);
333  fitLimitVsCrysID->SetBinContent(histbin, fithigh);
334  fitLimitVsCrysID->SetBinError(histbin, 0);
335  StatusVsCrysID->SetBinContent(histbin, iStatus);
336  StatusVsCrysID->SetBinError(histbin, 0);
337 
339  hStatus->Fill(iStatus);
340  nIterations->Fill(nIter);
341  if (iStatus >= iterations) {
342  hPeak->Fill(peakE);
343  fracPeakUnc->Fill(peakEUnc / peakE);
344  }
345 
347  B2INFO("cellID " << crysID + 1 << " status = " << iStatus << " fit probability = " << fitProb);
348  histfile->cd();
349  hEnergy->Write();
350 
351  } /* end of loop over crystals */
352 
355  for (int crysID = 0; crysID < 8736; crysID++) {
356  int histbin = crysID + 1;
357  double fitstatus = StatusVsCrysID->GetBinContent(histbin);
358  double peakE = PeakVsCrysID->GetBinContent(histbin);
359  double fracpeakEUnc = PeakVsCrysID->GetBinError(histbin) / peakE;
360 
362  if (fitstatus < iterations) {
363  peakE = -1.;
364  fracpeakEUnc = 0.;
365  if (histbin >= cellIDLo && histbin <= cellIDHi) {
366  B2RESULT("eclMuMuEAlgorithm: cellID " << histbin << " is not a successful fit. Status = " << fitstatus);
367  allFitsOK = false;
368  }
369  }
370 
372  if (findExpValues) {
373  double inputExpE = abs(AverageExpECrys->GetBinContent(histbin));
374  ExpEnergyperCrys->SetBinContent(histbin, inputExpE * peakE);
375  ExpEnergyperCrys->SetBinError(histbin, fracpeakEUnc * inputExpE * peakE);
376  } else {
377 
379  double inputCalib = abs(AverageInitCalib->GetBinContent(histbin));
380  CalibVsCrysID->SetBinContent(histbin, inputCalib / peakE);
381  CalibVsCrysID->SetBinError(histbin, fracpeakEUnc * inputCalib / peakE);
382  }
383  }
384 
387  bool DBsuccess = false;
388  if (storeConst == 0 || (storeConst == 1 && allFitsOK)) {
389  DBsuccess = true;
390  if (findExpValues) {
391 
393  std::vector<float> tempE;
394  std::vector<float> tempUnc;
395  for (int crysID = 0; crysID < 8736; crysID++) {
396  tempE.push_back(ExpEnergyperCrys->GetBinContent(crysID + 1));
397  tempUnc.push_back(ExpEnergyperCrys->GetBinError(crysID + 1));
398  }
399  ECLCrystalCalib* ExpectedE = new ECLCrystalCalib();
400  ExpectedE->setCalibVector(tempE, tempUnc);
401  saveCalibration(ExpectedE, "ECLExpMuMuE");
402  B2RESULT("eclCosmicEAlgorithm: successfully stored expected energies ECLExpMuMuE");
403 
404  } else {
405 
407  std::vector<float> tempCalib;
408  std::vector<float> tempCalibUnc;
409  for (int crysID = 0; crysID < 8736; crysID++) {
410  tempCalib.push_back(CalibVsCrysID->GetBinContent(crysID + 1));
411  tempCalibUnc.push_back(CalibVsCrysID->GetBinError(crysID + 1));
412  }
413  ECLCrystalCalib* MuMuECalib = new ECLCrystalCalib();
414  MuMuECalib->setCalibVector(tempCalib, tempCalibUnc);
415  saveCalibration(MuMuECalib, "ECLCrystalEnergyMuMu");
416  B2RESULT("eclMuMuEAlgorithm: successfully stored ECLCrystalEnergyMuMu calibration constants");
417  }
418  }
419 
423  PeakVsCrysID->Write();
424  effSigVsCrysID->Write();
425  etaVsCrysID->Write();
426  constVsCrysID->Write();
427  normVsCrysID->Write();
428  fitLimitVsCrysID->Write();
429  StatusVsCrysID->Write();
430  hPeak->Write();
431  fracPeakUnc->Write();
432  nIterations->Write();
433  hStatus->Write();
434 
436  if (findExpValues) {
437  ExpEnergyperCrys->Write();
438  } else {
439  CalibVsCrysID->Write();
440  }
441  histfile->Close();
442 
445  dummy = (TH1F*)gROOT->FindObject("PeakVsCrysID"); delete dummy;
446  dummy = (TH1F*)gROOT->FindObject("effSigVsCrysID"); delete dummy;
447  dummy = (TH1F*)gROOT->FindObject("etaVsCrysID"); delete dummy;
448  dummy = (TH1F*)gROOT->FindObject("constVsCrysID"); delete dummy;
449  dummy = (TH1F*)gROOT->FindObject("normVsCrysID"); delete dummy;
450  dummy = (TH1F*)gROOT->FindObject("fitLimitVsCrysID"); delete dummy;
451  dummy = (TH1F*)gROOT->FindObject("StatusVsCrysID"); delete dummy;
452  dummy = (TH1F*)gROOT->FindObject("fitProbSame"); delete dummy;
453  dummy = (TH1F*)gROOT->FindObject("fracPeakUnc"); delete dummy;
454  dummy = (TH1F*)gROOT->FindObject("nIterations"); delete dummy;
455  dummy = (TH1F*)gROOT->FindObject("hStatus"); delete dummy;
456  dummy = (TH1F*)gROOT->FindObject("ExpEnergyperCrys"); delete dummy;
457  dummy = (TH1F*)gROOT->FindObject("CalibVsCrysID"); delete dummy;
458 
459 
462  if (storeConst == -1) {
463  B2RESULT("eclMuMuEAlgorithm performed fits but was not asked to store contants");
464  return c_Failure;
465  } else if (!DBsuccess) {
466  if (findExpValues) { B2RESULT("eclMuMuEAlgorithm: failed to store expected values"); }
467  else { B2RESULT("eclMuMuEAlgorithm: failed to store calibration constants"); }
468  return c_Failure;
469  }
470  return c_OK;
471 }
Belle2::ECL::eclMuMuEAlgorithm::notFit
int notFit
no fit performed
Definition: eclMuMuEAlgorithm.h:73
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::ECL::eclMuMuEAlgorithm::storeConst
int storeConst
controls which values are written to the database.
Definition: eclMuMuEAlgorithm.h:57
Belle2::ECL::eclMuMuEAlgorithm::tRatioMin
double tRatioMin
Fit range is adjusted so that fit at upper endpoint is between tRatioMin and tRatioMax of peak.
Definition: eclMuMuEAlgorithm.h:53
Belle2::CalibrationAlgorithm::c_OK
@ c_OK
Finished successfuly =0 in Python.
Definition: CalibrationAlgorithm.h:51
Belle2::ECLCrystalCalib
General DB object to store one calibration number per ECL crystal.
Definition: ECLCrystalCalib.h:34
Belle2::CalibrationAlgorithm::setDescription
void setDescription(const std::string &description)
Set algorithm description (in constructor)
Definition: CalibrationAlgorithm.h:331
Belle2::ECL::eclMuMuEAlgorithm::atLimit
int atLimit
a parameter is at the limit; fit not useable
Definition: eclMuMuEAlgorithm.h:70
Belle2::ECL::eclMuMuEAlgorithm::poorFit
int poorFit
low chi square; fit not useable
Definition: eclMuMuEAlgorithm.h:71
Belle2::ECL::eclMuMuEAlgorithm::findExpValues
bool findExpValues
if true, fits are used to find expected energy deposit for each crystal instead of the calibration co...
Definition: eclMuMuEAlgorithm.h:56
Belle2::ECL::eclMuMuEAlgorithm::cellIDHi
int cellIDHi
Last cellID to be fit.
Definition: eclMuMuEAlgorithm.h:50
Belle2::ECL::eclMuMuEAlgorithm::performFits
bool performFits
if false, input histograms are copied to output, but no fits are done.
Definition: eclMuMuEAlgorithm.h:55
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::ECL::eclMuMuEAlgorithm::calibrate
virtual EResult calibrate() override
..Run algorithm on events
Definition: eclMuMuEAlgorithm.cc:51
Belle2::ECL::eclMuMuEAlgorithm::cellIDLo
int cellIDLo
..Parameters to control Novosibirsk fit to energy deposited in each crystal by mu+mu- events
Definition: eclMuMuEAlgorithm.h:49
Belle2::CalibrationAlgorithm::c_Failure
@ c_Failure
Failed =3 in Python.
Definition: CalibrationAlgorithm.h:54
Belle2::ECL::eclMuMuEAlgorithm::noPeak
int noPeak
Novosibirsk component of fit is negligible; fit not useable.
Definition: eclMuMuEAlgorithm.h:72
Belle2::ECL::eclMuMuEAlgorithm::minEntries
int minEntries
All crystals to be fit must have at least minEntries events in the fit range.
Definition: eclMuMuEAlgorithm.h:51
Belle2::CalibrationAlgorithm::EResult
EResult
The result of calibration.
Definition: CalibrationAlgorithm.h:50
Belle2::CalibrationAlgorithm::c_NotEnoughData
@ c_NotEnoughData
Needs more data =2 in Python.
Definition: CalibrationAlgorithm.h:53
Belle2::ECLCrystalCalib::setCalibVector
void setCalibVector(const std::vector< float > &CalibConst, const std::vector< float > &CalibConstUnc)
Set vector of constants with uncertainties.
Definition: ECLCrystalCalib.h:48
Belle2::ECL::eclMuMuEAlgorithm::fitOK
int fitOK
fit is OK
Definition: eclMuMuEAlgorithm.h:68
Belle2::ECL::eclMuMuEAlgorithm::tRatioMax
double tRatioMax
Fit range is adjusted so that fit at upper endpoint is between tRatioMin and tRatioMax of peak.
Definition: eclMuMuEAlgorithm.h:54
Belle2::CalibrationAlgorithm
Base class for calibration algorithms.
Definition: CalibrationAlgorithm.h:47
Belle2::ECL::eclMuMuEAlgorithm::iterations
int iterations
fit reached max number of iterations, but is useable
Definition: eclMuMuEAlgorithm.h:69
Belle2::ECL::eclMuMuEAlgorithm::maxIterations
int maxIterations
no more than maxIteration iterations
Definition: eclMuMuEAlgorithm.h:52