1 #include <ecl/calibration/eclCosmicEAlgorithm.h>
2 #include <ecl/dbobjects/ECLCrystalCalib.h>
15 double eclCosmicNovoConst(
double* x,
double* par)
20 double width = par[2];
21 double sln4 = sqrt(log(4));
22 double y = par[3] * sln4;
23 double tail = -log(y + sqrt(1 + y * y)) / sln4;
25 if (TMath::Abs(tail) < 1.e-7) {
26 qc = 0.5 * TMath::Power(((x[0] - peak) / width), 2);
28 const double qa = tail * sqrt(log(4.));
29 const double qb = sinh(qa) / qa;
30 const double qx = (x[0] - peak) / width * qb;
31 const double qy = 1. + tail * qx;
34 qc = 0.5 * (TMath::Power((log(qy) / tail), 2) + tail * tail);
38 return par[0] * exp(-qc) + par[4];
41 eclCosmicEAlgorithm::eclCosmicEAlgorithm():
CalibrationAlgorithm(
"eclCosmicECollector"), cellIDLo(1), cellIDHi(8736),
42 minEntries(150), maxIterations(10), tRatioMin(0.2), tRatioMax(0.25), performFits(true), findExpValues(false), storeConst(0)
45 "Perform energy calibration of ecl crystals by fitting a Novosibirsk function to energy deposited by cosmic rays"
54 double limitTol = 0.0005;
55 double minFitLimit = 1e-25;
56 double minFitProbIter = 1e-8;
57 double constRatio = 0.5;
58 double peakMin(0.5), peakMax(1.75);
59 double peakTol = limitTol * (peakMax - peakMin);
60 double effSigMin(0.08), effSigMax(0.5);
61 double effSigTol = limitTol * (effSigMax - effSigMin);
62 double etaMin(-3.), etaMax(1.);
64 double etaTol = limitTol * (etaMax - etaMin);
65 double constTol = 0.1;
73 dummy = (TH1F*)gROOT->FindObject(
"EnvsCrysSameRing");
74 if (dummy) {
delete dummy;}
75 dummy = (TH1F*)gROOT->FindObject(
"EnvsCrysDifferentRing");
76 if (dummy) {
delete dummy;}
77 dummy = (TH1F*)gROOT->FindObject(
"IntegralVsCrysSame");
78 if (dummy) {
delete dummy;}
79 dummy = (TH1F*)gROOT->FindObject(
"IntegralVsCrysDifferent");
80 if (dummy) {
delete dummy;}
81 dummy = (TH1F*)gROOT->FindObject(
"AverageExpECrysSame");
82 if (dummy) {
delete dummy;}
83 dummy = (TH1F*)gROOT->FindObject(
"AverageExpECrysDifferent");
84 if (dummy) {
delete dummy;}
85 dummy = (TH1F*)gROOT->FindObject(
"AverageElecCalibSame");
86 if (dummy) {
delete dummy;}
87 dummy = (TH1F*)gROOT->FindObject(
"AverageElecCalibDifferent");
88 if (dummy) {
delete dummy;}
89 dummy = (TH1F*)gROOT->FindObject(
"AverageInitialCalibSame");
90 if (dummy) {
delete dummy;}
91 dummy = (TH1F*)gROOT->FindObject(
"AverageInitialCalibDifferent");
92 if (dummy) {
delete dummy;}
96 std::vector<std::shared_ptr<TH2F>> EnvsCrys;
97 EnvsCrys.push_back(getObjectPtr<TH2F>(
"EnvsCrysSameRing"));
98 EnvsCrys.push_back(getObjectPtr<TH2F>(
"EnvsCrysDifferentRing"));
100 std::vector<std::shared_ptr<TH1F>> ExpEvsCrys;
101 ExpEvsCrys.push_back(getObjectPtr<TH1F>(
"ExpEvsCrysSameRing"));
102 ExpEvsCrys.push_back(getObjectPtr<TH1F>(
"ExpEvsCrysDifferentRing"));
104 std::vector<std::shared_ptr<TH1F>> ElecCalibvsCrys;
105 ElecCalibvsCrys.push_back(getObjectPtr<TH1F>(
"ElecCalibvsCrysSameRing"));
106 ElecCalibvsCrys.push_back(getObjectPtr<TH1F>(
"ElecCalibvsCrysDifferentRing"));
108 std::vector<std::shared_ptr<TH1F>> InitialCalibvsCrys;
109 InitialCalibvsCrys.push_back(getObjectPtr<TH1F>(
"InitialCalibvsCrysSameRing"));
110 InitialCalibvsCrys.push_back(getObjectPtr<TH1F>(
"InitialCalibvsCrysDifferentRing"));
112 std::vector<std::shared_ptr<TH1F>> CalibEntriesvsCrys;
113 CalibEntriesvsCrys.push_back(getObjectPtr<TH1F>(
"CalibEntriesvsCrysSameRing"));
114 CalibEntriesvsCrys.push_back(getObjectPtr<TH1F>(
"CalibEntriesvsCrysDifferentRing"));
116 auto RawDigitAmpvsCrys = getObjectPtr<TH2F>(
"RawDigitAmpvsCrys");
121 TH1F* IntegralVsCrys[2];
122 IntegralVsCrys[0] =
new TH1F(
"IntegralVsCrysSame",
"Integral of EnVsCrys for each crystal, same theta ring;Crystal ID", 8736, 0,
124 IntegralVsCrys[1] =
new TH1F(
"IntegralVsCrysDifferent",
"Integral of EnVsCrys for each crystal, different theta rings;Crystal ID",
127 TH1F* AverageExpECrys[2];
128 AverageExpECrys[0] =
new TH1F(
"AverageExpECrysSame",
129 "Average expected E per crys from collector, same theta ring;Crystal ID;Energy (GeV)", 8736, 0, 8736);
130 AverageExpECrys[1] =
new TH1F(
"AverageExpECrysDifferent",
131 "Average expected E per crys from collector, different theta ring;Crystal ID;Energy (GeV)", 8736, 0, 8736);
133 TH1F* AverageElecCalib[2];
134 AverageElecCalib[0] =
new TH1F(
"AverageElecCalibSame",
135 "Average electronics calib const vs crys, same theta ring;Crystal ID;Calibration constant", 8736, 0, 8736);
136 AverageElecCalib[1] =
new TH1F(
"AverageElecCalibDifferent",
137 "Average electronics calib const vs crys, different theta rings;Crystal ID;Calibration constant", 8736, 0, 8736);
139 TH1F* AverageInitialCalib[2];
140 AverageInitialCalib[0] =
new TH1F(
"AverageInitialCalibSame",
141 "Average initial cosmic calib const vs crys, same theta ring;Crystal ID;Calibration constant", 8736, 0, 8736);
142 AverageInitialCalib[1] =
new TH1F(
"AverageInitialCalibDifferent",
143 "Average initial cosmic calib const vs crys, different theta rings;Crystal ID;Calibration constant", 8736, 0, 8736);
145 for (
int crysID = 0; crysID < 8736; crysID++) {
146 int histbin = crysID + 1;
147 for (
int idir = 0; idir < 2; idir++) {
148 TH1D* hEnergy = EnvsCrys[idir]->ProjectionY(
"hEnergy", histbin, histbin);
149 int Integral = hEnergy->Integral();
150 IntegralVsCrys[idir]->SetBinContent(histbin, Integral);
152 double TotEntries = CalibEntriesvsCrys[idir]->GetBinContent(histbin);
154 double expectedE = 0.;
155 if (TotEntries > 0.) {expectedE = ExpEvsCrys[idir]->GetBinContent(histbin) / TotEntries;}
156 AverageExpECrys[idir]->SetBinContent(histbin, expectedE);
157 AverageExpECrys[idir]->SetBinError(histbin, 0.);
159 double calibconst = 0.;
160 if (TotEntries > 0.) {calibconst = ElecCalibvsCrys[idir]->GetBinContent(histbin) / TotEntries;}
161 AverageElecCalib[idir]->SetBinContent(histbin, calibconst);
162 AverageElecCalib[idir]->SetBinError(histbin, 0);
165 if (TotEntries > 0.) {calibconst = InitialCalibvsCrys[idir]->GetBinContent(histbin) / TotEntries;}
166 AverageInitialCalib[idir]->SetBinContent(histbin, calibconst);
167 AverageInitialCalib[idir]->SetBinError(histbin, 0);
173 TFile* histfile =
new TFile(
"eclCosmicEAlgorithm.root",
"recreate");
174 for (
int idir = 0; idir < 2; idir++) {
175 EnvsCrys[idir]->Write();
176 IntegralVsCrys[idir]->Write();
177 AverageExpECrys[idir]->Write();
178 AverageElecCalib[idir]->Write();
179 AverageInitialCalib[idir]->Write();
181 RawDigitAmpvsCrys->Write();
186 B2INFO(
"eclCosmicEAlgorithm has not been asked to perform fits; copying input histograms and quitting");
195 bool sufficientData =
true;
197 int histbin = crysID + 1;
198 bool SameLow = IntegralVsCrys[0]->GetBinContent(histbin) <
minEntries;
199 bool DifferentLow = IntegralVsCrys[1]->GetBinContent(histbin) <
minEntries;
200 if ((SameLow && DifferentLow) || (
findExpValues && (SameLow || DifferentLow))) {
201 if (
storeConst == 1) {B2INFO(
"eclCosmicEAlgorithm: cellID " << histbin <<
" has insufficient statistics: " << IntegralVsCrys[0]->GetBinContent(histbin) <<
" and " << IntegralVsCrys[1]->GetBinContent(histbin) <<
". Requirement is " <<
minEntries);}
202 sufficientData =
false;
216 TString preName[2] = {
"SameRing",
"DifferentRing"};
218 TH1F* PeakperCrys[2];
219 PeakperCrys[0] =
new TH1F(
"PeakperCrysSame",
"Fit peak per crystal, same theta ring;Crystal ID;Peak normalized energy", 8736, 0,
221 PeakperCrys[1] =
new TH1F(
"PeakperCrysDifferent",
"Fit peak per crystal, different theta ring;Crystal ID;Peak normalized energy",
224 TH1F* SigmaperCrys[2];
225 SigmaperCrys[0] =
new TH1F(
"SigmaperCrysSame",
"Fit sigma per crysal, same theta ring;Crystal ID;Sigma (ADC)", 8736, 0, 8736);
226 SigmaperCrys[1] =
new TH1F(
"SigmaperCrysDifferent",
"Fit sigma per crysal, different theta ring;Crystal ID;Sigma (ADC)", 8736, 0,
230 EtaperCrys[0] =
new TH1F(
"EtaperCrysSame",
"Fit eta per crysal, same theta ring;Crystal ID;Eta", 8736, 0, 8736);
231 EtaperCrys[1] =
new TH1F(
"EtaperCrysDifferent",
"Fit eta per crysal, different theta ring;Crystal ID;Eta", 8736, 0, 8736);
233 TH1F* ConstperCrys[2];
234 ConstperCrys[0] =
new TH1F(
"ConstperCrysSame",
"Fit constant per crystal, same theta ring;Crystal ID;Constant", 8736, 0, 8736);
235 ConstperCrys[1] =
new TH1F(
"ConstperCrysDifferent",
"Fit constant per crystal, different theta ring;Crystal ID;Constant", 8736, 0,
238 TH1F* StatusperCrys[2];
239 StatusperCrys[0] =
new TH1F(
"StatusperCrysSame",
"Fit status, same theta ring;Crystal ID;Status", 8736, 0, 8736);
240 StatusperCrys[1] =
new TH1F(
"StatusperCrysDifferent",
"Fit status, different theta ring;Crystal ID;Status", 8736, 0, 8736);
244 hStatus[0] =
new TH1F(
"StatusSame",
"Fit status, same theta ring", 25, -5, 20);
245 hStatus[1] =
new TH1F(
"StatusDifferent",
"Fit status, different theta ring", 25, -5, 20);
247 TH1F* fracPeakUnc[2];
248 fracPeakUnc[0] =
new TH1F(
"fracPeakUncSame",
"Fractional uncertainty on peak location, same theta ring", 100, 0, 0.1);
249 fracPeakUnc[1] =
new TH1F(
"fracPeakUncDifferent",
"Fractional uncertainty on peak location, different theta ring", 100, 0, 0.1);
252 hfitProb[0] =
new TH1F(
"fitProbSame",
"Probability of fit, same theta ring", 200, 0, 0.02);
253 hfitProb[1] =
new TH1F(
"fitProbDifferent",
"Probability of fit, different theta ring", 200, 0, 0.02);
256 TH1F* ExpEnergyperCrys[2];
257 ExpEnergyperCrys[0] =
new TH1F(
"ExpEnergyperCrysSame",
"Expected energy per crystal, same theta ring;Crystal ID;Peak energy (GeV)",
259 ExpEnergyperCrys[1] =
new TH1F(
"ExpEnergyperCrysDifferent",
260 "Expected energy per crystal, different theta ring;Crystal ID;Peak energy (GeV)", 8736, 0, 8736);
262 TH1F* CalibvsCrys =
new TH1F(
"CalibvsCrys",
"Energy calibration constant from cosmics;Crystal ID;Calibration constant", 8736, 0,
267 bool allFitsOK =
true;
269 int histbin = crysID + 1;
270 for (
int idir = 0; idir < 2; idir++) {
273 TString hname = preName[idir];
274 hname +=
"Enormalized";
276 TH1D* hEnergy = EnvsCrys[idir]->ProjectionY(hname, histbin, histbin);
279 double histMin = hEnergy->GetXaxis()->GetXmin();
280 double histMax = hEnergy->GetXaxis()->GetXmax();
281 TF1* func =
new TF1(
"eclCosmicNovoConst", eclCosmicNovoConst, histMin, histMax, 5);
282 func->SetParNames(
"normalization",
"peak",
"effSigma",
"eta",
"const");
283 func->SetParLimits(1, peakMin, peakMax);
284 func->SetParLimits(2, effSigMin, effSigMax);
285 func->SetParLimits(3, etaMin, etaMax);
288 hEnergy->GetXaxis()->SetRangeUser(peakMin, peakMax);
289 int maxBin = hEnergy->GetMaximumBin();
290 double peakE = hEnergy->GetBinLowEdge(maxBin);
291 double peakEUnc = 0.;
292 double normalization = hEnergy->GetMaximum();
293 double effSigma = hEnergy->GetRMS();
294 double sigmaUnc = 0.;
295 hEnergy->GetXaxis()->SetRangeUser(histMin, histMax);
298 double fitlow = 0.25;
299 double fithigh = peakE + 2.5 * effSigma;
302 int il0 = hEnergy->GetXaxis()->FindBin(fitlow);
303 int il1 = hEnergy->GetXaxis()->FindBin(fitlow + 0.1);
304 double constant = hEnergy->Integral(il0, il1) / (1 + il1 - il0);
305 double constUnc = 0.;
312 double dIter = 0.1 * (histMax - histMin) / hEnergy->GetNbinsX();
313 double highold(0.), higholdold(0.);
315 double fitProbDefault(0.);
316 bool fitHist = IntegralVsCrys[idir]->GetBinContent(histbin) >=
minEntries;
317 bool fixConst =
false;
326 func->SetParameters(normalization, peakE, effSigma, eta, constant);
327 if (fixConst) { func->FixParameter(4, 0); }
330 hEnergy->Fit(func,
"LIQ",
"", fitlow, fithigh);
331 normalization = func->GetParameter(0);
332 peakE = func->GetParameter(1);
333 peakEUnc = func->GetParError(1);
334 effSigma = func->GetParameter(2);
335 sigmaUnc = func->GetParError(2);
336 eta = func->GetParameter(3);
337 etaUnc = func->GetParError(3);
338 constant = func->GetParameter(4);
339 constUnc = func->GetParError(4);
340 fitProbDefault = func->GetProb();
344 double peak = func->Eval(peakE) - constant;
345 double tRatio = (func->Eval(fithigh) - constant) / peak;
346 if (tRatio < tRatioMin || tRatio >
tRatioMax) {
348 higholdold = highold;
350 fithigh = func->GetX(targetY, peakE, histMax);
354 if (abs(fithigh - higholdold) < dIter) {fithigh = 0.5 * (highold + higholdold); }
357 if (nIter >
maxIterations - 3) {fithigh = 0.33333 * (fithigh + highold + higholdold); }
361 if (constant < constTol && !fixConst) {
369 B2DEBUG(200,
"cellID = " << histbin <<
" " << nIter <<
" " << preName[idir] <<
" " << peakE <<
" " << constant <<
" " << tRatio <<
378 int lowbin = hEnergy->GetXaxis()->FindBin(fitlow);
379 int highbin = hEnergy->GetXaxis()->FindBin(fithigh);
381 if (fixConst) {npar = 4;}
382 int ndeg = (highbin - lowbin) + 1 - npar;
384 double halfbinwidth = 0.5 * hEnergy->GetBinWidth(1);
385 for (
int ib = lowbin; ib <= highbin; ib++) {
386 double yexp = func->Eval(hEnergy->GetBinLowEdge(ib) + halfbinwidth);
387 double yobs = hEnergy->GetBinContent(ib);
388 double dchi2 = (yexp - yobs) * (yexp - yobs) / yexp;
391 fitProb = 0.5 * (TMath::Prob(chisq, ndeg) + fitProbDefault);
400 if (normalization < constRatio * constant) {iStatus =
noPeak;}
403 if (fitProb <= minFitLimit || (fitProb < minFitProbIter && iStatus ==
iterations)) {iStatus =
poorFit;}
406 if ((peakE < peakMin + peakTol) || (peakE > peakMax - peakTol)) {iStatus =
atLimit;}
407 if ((effSigma < effSigMin + effSigTol) || (effSigma > effSigMax - effSigTol)) {iStatus =
atLimit;}
408 if ((eta < etaMin + etaTol) || (eta > etaMax - etaTol)) {iStatus =
atLimit;}
411 if (nIter == 0) {iStatus =
notFit;}
414 PeakperCrys[idir]->SetBinContent(histbin, peakE);
415 PeakperCrys[idir]->SetBinError(histbin, peakEUnc);
416 SigmaperCrys[idir]->SetBinContent(histbin, effSigma);
417 SigmaperCrys[idir]->SetBinError(histbin, sigmaUnc);
418 EtaperCrys[idir]->SetBinContent(histbin, eta);
419 EtaperCrys[idir]->SetBinError(histbin, etaUnc);
420 ConstperCrys[idir]->SetBinContent(histbin, constant);
421 ConstperCrys[idir]->SetBinError(histbin, constUnc);
422 StatusperCrys[idir]->SetBinContent(histbin, iStatus);
423 hStatus[idir]->Fill(iStatus);
424 fracPeakUnc[idir]->Fill(peakEUnc / peakE);
425 hfitProb[idir]->Fill(fitProb);
428 B2INFO(
"cellID " << histbin <<
" " << preName[idir] <<
" status = " << iStatus <<
" fit probability = " << fitProb);
439 for (
int crysID = 0; crysID < 8736; crysID++) {
440 int histbin = crysID + 1;
441 bool atLeastOneOK =
false;
442 for (
int idir = 0; idir < 2; idir++) {
443 double fitstatus = StatusperCrys[idir]->GetBinContent(histbin);
444 double peakE = PeakperCrys[idir]->GetBinContent(histbin);
445 double peakEUnc = PeakperCrys[idir]->GetBinError(histbin);
450 B2INFO(
"eclCosmicEAlgorithm: crystal " << crysID <<
" " << preName[idir] <<
" is not a successful fit. Status = " << fitstatus);
457 double inputExpE = abs(AverageExpECrys[idir]->GetBinContent(histbin));
458 ExpEnergyperCrys[idir]->SetBinContent(histbin, inputExpE * peakE);
459 ExpEnergyperCrys[idir]->SetBinError(histbin, inputExpE * peakEUnc / peakE);
461 if (!atLeastOneOK) {allFitsOK =
false;}
469 for (
int crysID = 0; crysID < 8736; crysID++) {
470 int histbin = crysID + 1;
471 double calibConst[2] = {};
472 double calibConstUnc[2] = {999999., 999999.};
473 double weight[2] = {};
474 bool bothFitsBad =
true;
475 for (
int idir = 0; idir < 2; idir++) {
478 double peakE = PeakperCrys[idir]->GetBinContent(histbin);
479 double fracPeakEUnc = PeakperCrys[idir]->GetBinError(histbin) / peakE;
480 double inputConst = AverageInitialCalib[idir]->GetBinContent(histbin);
481 double fitstatus = StatusperCrys[idir]->GetBinContent(histbin);
482 double inputExpE = AverageExpECrys[idir]->GetBinContent(histbin);
483 if (fitstatus >=
iterations && inputConst == 0) {B2FATAL(
"eclCosmicEAlgorithm: input calibration = 0 for idir = " << idir <<
" and crysID = " << crysID);}
486 if (fitstatus >=
iterations && inputExpE > 0.) {
487 calibConst[idir] = abs(inputConst) / peakE;
488 calibConstUnc[idir] = calibConst[idir] * fracPeakEUnc / peakE;
489 weight[idir] = 1. / (calibConstUnc[idir] * calibConstUnc[idir]);
493 B2INFO(
"eclCosmicEAlgorithm: cellID " << histbin <<
" " << preName[idir] <<
" is not a successful fit. Status = " << fitstatus);
495 B2WARNING(
"eclCosmicEAlgorithm: cellID " << histbin <<
" " << preName[idir] <<
" has no expected energy. Status = " << fitstatus);
502 double averageConstUnc;
506 if (histbin >=
cellIDLo && histbin <=
cellIDHi) {B2INFO(
"eclCosmicEAlgorithm: no constant found for cellID = " << histbin <<
" status = " << StatusperCrys[0]->GetBinContent(histbin) <<
" and " << StatusperCrys[1]->GetBinContent(histbin));}
507 averageConst = -1.*abs(AverageInitialCalib[0]->GetBinContent(histbin));
508 averageConstUnc = 0.;
510 averageConst = (calibConst[0] * weight[0] + calibConst[1] * weight[1]) / (weight[0] + weight[1]);
511 averageConstUnc = 1. / sqrt(weight[0] + weight[1]);
513 CalibvsCrys->SetBinContent(histbin, averageConst);
514 CalibvsCrys->SetBinError(histbin, averageConstUnc);
520 bool DBsuccess =
false;
526 std::vector<std::string> DBname = {
"ECLExpCosmicESame",
"ECLExpCosmicEDifferent"};
527 for (
int idir = 0; idir < 2; idir++) {
528 std::vector<float> tempE;
529 std::vector<float> tempUnc;
530 for (
int crysID = 0; crysID < 8736; crysID++) {
531 int histbin = crysID + 1;
532 tempE.push_back(ExpEnergyperCrys[idir]->GetBinContent(histbin));
533 tempUnc.push_back(ExpEnergyperCrys[idir]->GetBinError(histbin));
538 B2INFO(
"eclCosmicEAlgorithm: successfully stored expected values for " << DBname[idir]);
543 std::vector<float> tempCalib;
544 std::vector<float> tempCalibUnc;
545 for (
int crysID = 0; crysID < 8736; crysID++) {
546 int histbin = crysID + 1;
547 tempCalib.push_back(CalibvsCrys->GetBinContent(histbin));
548 tempCalibUnc.push_back(CalibvsCrys->GetBinError(histbin));
553 B2INFO(
"eclCosmicEAlgorithm: successfully stored calibration constants");
558 for (
int idir = 0; idir < 2; idir++) {
559 PeakperCrys[idir]->Write();
560 SigmaperCrys[idir]->Write();
561 EtaperCrys[idir]->Write();
562 ConstperCrys[idir]->Write();
563 StatusperCrys[idir]->Write();
564 hStatus[idir]->Write();
565 fracPeakUnc[idir]->Write();
566 hfitProb[idir]->Write();
571 ExpEnergyperCrys[0]->Write();
572 ExpEnergyperCrys[1]->Write();
574 CalibvsCrys->Write();
580 dummy = (TH1F*)gROOT->FindObject(
"PeakperCrysSame");
delete dummy;
581 dummy = (TH1F*)gROOT->FindObject(
"SigmaperCrysSame");
delete dummy;
582 dummy = (TH1F*)gROOT->FindObject(
"EtaperCrysSame");
delete dummy;
583 dummy = (TH1F*)gROOT->FindObject(
"ConstperCrysSame");
delete dummy;
584 dummy = (TH1F*)gROOT->FindObject(
"StatusperCrysSame");
delete dummy;
585 dummy = (TH1F*)gROOT->FindObject(
"StatusSame");
delete dummy;
586 dummy = (TH1F*)gROOT->FindObject(
"fracPeakUncSame");
delete dummy;
587 dummy = (TH1F*)gROOT->FindObject(
"fitProbSame");
delete dummy;
588 dummy = (TH1F*)gROOT->FindObject(
"ExpEnergyperCrysSame");
delete dummy;
589 dummy = (TH1F*)gROOT->FindObject(
"PeakperCrysDifferent");
delete dummy;
590 dummy = (TH1F*)gROOT->FindObject(
"SigmaperCrysDifferent");
delete dummy;
591 dummy = (TH1F*)gROOT->FindObject(
"EtaperCrysDifferent");
delete dummy;
592 dummy = (TH1F*)gROOT->FindObject(
"ConstperCrysDifferent");
delete dummy;
593 dummy = (TH1F*)gROOT->FindObject(
"StatusperCrysDifferent");
delete dummy;
594 dummy = (TH1F*)gROOT->FindObject(
"StatusDifferent");
delete dummy;
595 dummy = (TH1F*)gROOT->FindObject(
"fracPeakUncDifferent");
delete dummy;
596 dummy = (TH1F*)gROOT->FindObject(
"fitProbDifferent");
delete dummy;
597 dummy = (TH1F*)gROOT->FindObject(
"ExpEnergyperCrysDifferent");
delete dummy;
598 dummy = (TH1F*)gROOT->FindObject(
"CalibvsCrys");
delete dummy;
603 B2INFO(
"eclCosmicEAlgorithm performed fits but was not asked to store contants");
605 }
else if (!DBsuccess) {
606 if (
findExpValues) { B2INFO(
"eclCosmicEAlgorithm: failed to store expected values"); }
607 else { B2INFO(
"eclCosmicEAlgorithm: failed to store calibration constants"); }