Belle II Software development
CDCDedxWireGainAlgorithm.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <cdc/calibration/CDCdEdx/CDCDedxWireGainAlgorithm.h>
10
11#include <cdc/geometry/CDCGeometryPar.h>
12
13#include <TCanvas.h>
14#include <TLine.h>
15#include <TStyle.h>
16
17using namespace Belle2;
18using namespace CDC;
19using namespace std;
20
21//-----------------------------------------------------------------
22// Implementation
23//-----------------------------------------------------------------
25 CalibrationAlgorithm("CDCDedxElectronCollector"),
26 m_isMakePlots(true),
27 m_isMerge(true),
28 m_isWireTruc(false),
29 m_dedxBins(250),
30 m_dedxMin(0.0),
31 m_dedxMax(5.0),
32 m_truncMin(0.05),
33 m_truncMax(0.75),
34 m_suffix("")
35{
36 // Set module properties
37 setDescription("A calibration algorithm for CDC dE/dx wire gains");
38}
39
40//-----------------------------------------------------------------
41// Run the calibration
42//-----------------------------------------------------------------
44{
45
46 //if not checking then copy function here;
48
49 if (!m_DBBadWires.isValid() || !m_DBWireGains.isValid())
50 B2FATAL("There is no valid payload for BadWires and/or Wirgain");
51
52 // Get data objects
53 auto ttree = getObjectPtr<TTree>("tree");
54 if (ttree->GetEntries() < 100)return c_NotEnoughData;
55
56 vector<int>* wire = 0;
57 ttree->SetBranchAddress("wire", &wire);
58
59 vector<double>* dedxhit = 0;
60 ttree->SetBranchAddress("dedxhit", &dedxhit);
61
62 // dedxhit vector to store dE/dx values for each wire
63 map<int, vector<double>> wirededx;
64
65 //IL = Inner Layer and OL = Outer Layer
66 array<TH1D*, 2> hdedxL;
67 string label[2] = {"IL", "OL"};
68 for (int il = 0; il < 2; il++)
69 hdedxL[il] = new TH1D(Form("hdedx%s_%s", label[il].data(), m_suffix.data()), "", m_dedxBins, m_dedxMin, m_dedxMax);
70
71 for (int i = 0; i < ttree->GetEntries(); ++i) {
72 ttree->GetEvent(i);
73 for (unsigned int j = 0; j < wire->size(); ++j) {
74 int jwire = wire->at(j);
75 double jhitdedx = dedxhit->at(j);
76 wirededx[jwire].push_back(jhitdedx);
77 //wire # 1279 end of inner layers
78 if (jwire < 1280) hdedxL[0]->Fill(jhitdedx);
79 else hdedxL[1]->Fill(jhitdedx);
80 }
81 }
82
83 //return if ~10% low stats or dead wires
84 int minstat = 0;
85 for (unsigned int jw = 0; jw < c_nwireCDC; ++jw)
86 if (wirededx[jw].size() <= 100) minstat++;
87
88 if (minstat > 0.10 * c_nwireCDC) return c_NotEnoughData;
89
90 //25-75 is average bin # for truncation
91 array<unsigned int, 2> minbinL, maxbinL;
92 for (int il = 0; il < 2; il++) {
93 getTruncatedBins(hdedxL[il], minbinL[il], maxbinL[il]);
94 hdedxL[il]->SetTitle(Form("%s(%s);%d;%d", label[il].data(), m_suffix.data(), minbinL[il], maxbinL[il]));
95 }
96
97 vector<double> vrel_mean;
98 vector<double> vdedx_means;
99
100 vector<TH1D*> hdedxhit(c_nwireCDC);
101
102 B2INFO("Creating CDCGeometryPar object");
104
105 vector<double> layermean(c_maxNSenseLayers);
106
107 int activelayers = 0;
108 double layeravg = 0.0;
109
110 int jwire = -1;
111 for (unsigned int il = 0; il < c_maxNSenseLayers; ++il) {
112
113 int activewires = 0;
114 layermean[il] = 0.0;
115
116 for (unsigned int iw = 0; iw < cdcgeo.nWiresInLayer(il); ++iw) {
117
118 jwire++;
119 hdedxhit[jwire] = new TH1D(Form("h%s_w%d", m_suffix.data(), jwire), "", m_dedxBins, m_dedxMin, m_dedxMax);
120
121 for (unsigned int ih = 0; ih < wirededx[jwire].size(); ++ih) {
122 hdedxhit[jwire]->Fill(wirededx[jwire][ih]);
123 }
124
125 unsigned int minbin, maxbin;
126 if (!m_isWireTruc) {
127 getTruncatedBins(hdedxhit[jwire], minbin, maxbin);
128 } else {
129 if (jwire < 1280) {
130 minbin = minbinL[0];
131 maxbin = maxbinL[0];
132 } else {
133 minbin = minbinL[1];
134 maxbin = maxbinL[1];
135 }
136 }
137 hdedxhit[jwire]->SetTitle(Form("dedxhit-dist, wire: %d (%s);%d;%d", jwire, m_suffix.data(), minbin, maxbin));
138
139 double dedxmean;
140 if (m_DBBadWires->getBadWireStatus(jwire) == kTRUE) dedxmean = 0.0;
141
142 else dedxmean = getTruncationMean(hdedxhit[jwire], minbin, maxbin);
143 vrel_mean.push_back(dedxmean);
144
145 double prewg = m_DBWireGains->getWireGain(jwire);
146 if (prewg > 0.0 && m_isMerge) {
147 vdedx_means.push_back(dedxmean * prewg);
148 B2INFO("merged-wireGain: [" << jwire << "], prewgvious = " << prewg << ", rel = " << dedxmean << ", merged = " << vdedx_means.at(
149 jwire));
150 } else vdedx_means.push_back(dedxmean);
151
152 //calculate layer average for active wires
153 if (vdedx_means.at(jwire) > 0) {
154 layermean[il] += vdedx_means.at(jwire);
155 activewires++;
156 }
157 }
158
159 if (activewires > 0) layermean[il] /= activewires;
160 else layermean[il] = 1.0;
161
162 //calculate outer layer average for active layer
163 if (il >= 8 && layermean[il] > 0) {
164 layeravg += layermean[il];
165 activelayers++;
166 }
167
168 }
169
170 //normalisation of wiregains to get outergain ~1.0
171 if (activelayers > 0) layeravg /= activelayers;
172
173 for (unsigned int iw = 0; iw < c_nwireCDC; ++iw) {
174 vrel_mean.at(iw) /= layeravg;
175 vdedx_means.at(iw) /= layeravg;
176 }
177
178 if (m_isMakePlots) {
179
180 //1. Inner and Outer layer dE/dx distributions
181 plotLayerDist(hdedxL);
182
183 //2. wiredist -> 14336 (good, bad)
184 plotWireDist(hdedxhit, vrel_mean);
185
186 //3. wiregains draw (1D, dist)
187 plotWireGain(vdedx_means, vrel_mean, layeravg);
188
189 //4. layer gain plot
190 plotLayerGain(layermean, layeravg);
191
192 //5. wiregains per layer
193 plotWGPerLayer(vdedx_means, layermean, layeravg);
194
195 //6. Statstics plot
197 }
198
199 createPayload(vdedx_means);
200 m_suffix.clear();
201 return c_OK;
202}
203
204//------------------------------------
206{
207
208 int cruns = 0;
209 for (auto expRun : getRunList()) {
210 if (cruns == 0) B2INFO("CDCDedxWireGain: start exp " << expRun.first << " and run " << expRun.second << "");
211 cruns++;
212 }
213
214 const auto erStart = getRunList()[0];
215 int estart = erStart.first;
216 int rstart = erStart.second;
217
218 const auto erEnd = getRunList()[cruns - 1];
219 int rend = erEnd.second;
220
221 updateDBObjPtrs(1, rstart, estart);
222
223 if (m_suffix.length() > 0) m_suffix = Form("%s_e%d_r%dr%d", m_suffix.data(), estart, rstart, rend);
224 else m_suffix = Form("e%d_r%dr%d", estart, rstart, rend);
225}
226
227//--------------------------
228void CDCDedxWireGainAlgorithm::createPayload(const vector<double>& vdedx_means)
229{
230
231 B2INFO("dE/dx Calibration done for " << vdedx_means.size() << " CDC wires");
232 CDCDedxWireGain* gains = new CDCDedxWireGain(vdedx_means);
233 saveCalibration(gains, "CDCDedxWireGain");
234}
235
236//--------------------------
237double CDCDedxWireGainAlgorithm::getTruncationMean(TH1D* hdedxhit, int binlow, int binhigh)
238{
239
240 //calculating truncation average
241 if (hdedxhit->Integral() < 100) return 1.0;
242
243 if (binlow <= 0 || binhigh > hdedxhit->GetNbinsX()) return 1.0;
244
245 double binweights = 0., sumofbc = 0.;
246 for (int ibin = binlow; ibin <= binhigh; ibin++) {
247 double bcdedx = hdedxhit->GetBinContent(ibin);
248 if (bcdedx > 0) {
249 binweights += (bcdedx * hdedxhit->GetBinCenter(ibin));
250 sumofbc += bcdedx;
251 }
252 }
253 if (sumofbc > 0) return binweights / sumofbc;
254 else return 1.0;
255}
256
257//--------------------------
258void CDCDedxWireGainAlgorithm::getTruncatedBins(TH1D* hdedxhit, unsigned int& binlow, unsigned int& binhigh)
259{
260
261 //calculating truncation average
262 double sum = hdedxhit->Integral();
263 if (sum <= 0 || hdedxhit->GetNbinsX() <= 0) {
264 binlow = 1; binhigh = 1;
265 return ;
266 }
267
268 binlow = 1.0; binhigh = 1.0;
269 double sumPer5 = 0.0, sumPer75 = 0.0;
270 for (int ibin = 1; ibin <= hdedxhit->GetNbinsX(); ibin++) {
271 double bcdedx = hdedxhit->GetBinContent(ibin);
272 if (sumPer5 <= m_truncMin * sum) {
273 sumPer5 += bcdedx;
274 binlow = ibin;
275 }
276 if (sumPer75 <= m_truncMax * sum) {
277 sumPer75 += bcdedx;
278 binhigh = ibin;
279 }
280 }
281 return;
282}
283
284//--------------------------
285void CDCDedxWireGainAlgorithm::plotLayerDist(array<TH1D*, 2> hdedxL)
286{
287
288 TCanvas cldedx("cldedx", "IL/OL dedxhit dist", 900, 400);
289 cldedx.Divide(2, 1);
290
291 for (int il = 0; il < 2; il++) {
292 cldedx.cd(il + 1);
293
294 int minbin = stoi(hdedxL[il]->GetXaxis()->GetTitle());
295 int maxbin = stoi(hdedxL[il]->GetYaxis()->GetTitle());
296 double lowedge = hdedxL[il]->GetXaxis()->GetBinLowEdge(minbin);
297 double upedge = hdedxL[il]->GetXaxis()->GetBinUpEdge(maxbin);
298
299 hdedxL[il]->SetFillColor(kYellow);
300 hdedxL[il]->SetTitle(Form("%s, trunc(%0.02f - %0.02f);dedxhit;entries", hdedxL[il]->GetTitle(), lowedge, upedge));
301 hdedxL[il]->Draw("histo");
302
303 TH1D* hdedxLC = (TH1D*)hdedxL[il]->Clone(Form("%s_c", hdedxL[il]->GetName()));
304 hdedxLC->GetXaxis()->SetRange(minbin, maxbin);
305 hdedxLC->SetFillColor(kAzure + 1);
306 hdedxLC->Draw("same histo");
307 }
308
309 cldedx.SaveAs(Form("cdcdedx_wgcal_layerdedx_%s.pdf", m_suffix.data()));
310}
311
312//------------------------------------
313void CDCDedxWireGainAlgorithm::plotWireDist(const vector<TH1D*>& hist, const vector<double>& vrel_mean)
314{
315
316 TCanvas ctmp(Form("cdcdedx_%s", m_suffix.data()), "", 1200, 1200);
317 ctmp.Divide(4, 4);
318 ctmp.SetBatch(kTRUE);
319
320 stringstream psname;
321 psname << Form("cdcdedx_wgcal_%s.pdf[", m_suffix.data());
322 ctmp.Print(psname.str().c_str());
323 psname.str("");
324 psname << Form("cdcdedx_wgcal_%s.pdf", m_suffix.data());
325
326 for (unsigned int iw = 0; iw < hist.size(); iw++) {
327
328 int minbin = stoi(hist[iw]->GetXaxis()->GetTitle());
329 int maxbin = stoi(hist[iw]->GetYaxis()->GetTitle());
330
331 hist[iw]->SetFillColor(kYellow - 9);
332 hist[iw]->SetTitle(Form("%s, rel. #mu_{trunc} %0.03f;dedxhit;entries", hist[iw]->GetTitle(), vrel_mean.at(iw)));
333
334 if (m_DBBadWires->getBadWireStatus(iw) == kTRUE) {
335 hist[iw]->SetLineColor(kRed);
336 hist[iw]->SetLineWidth(2);
337 }
338 ctmp.cd(iw % 16 + 1);
339 hist[iw]->DrawCopy();
340
341 TH1D* hdedxhitC = (TH1D*)hist[iw]->Clone(Form("%sC", hist[iw]->GetName()));
342 hdedxhitC->GetXaxis()->SetRange(minbin, maxbin);
343 hdedxhitC->SetFillColor(kAzure + 1);
344 hdedxhitC->DrawCopy("same histo");
345
346 if (((iw + 1) % 16 == 0) || iw == (hist.size() - 1)) {
347 ctmp.Print(psname.str().c_str());
348 ctmp.Clear("D");
349 }
350 delete hist[iw];
351 delete hdedxhitC;
352 }
353
354 psname.str("");
355 psname << Form("cdcdedx_wgcal_%s.pdf]", m_suffix.data());
356 ctmp.Print(psname.str().c_str());
357}
358
359//------------------------------------
360void CDCDedxWireGainAlgorithm::plotWireGain(const vector<double>& vdedx_means, const vector<double>& vrel_mean, double layeravg)
361{
362
363 //saving final constants in a histograms for validation
364 TCanvas cwconst("cwconst", "", 900, 500);
365 TCanvas cwconstvar("cwconstvar", "", 500, 400);
366
367 array<TH1D*, 2> hconstpw, hconstpwvar;
368
369 for (int i = 0; i < 2; i++) {
370
371 hconstpw[i] = new TH1D(Form("hconstpw_%d_%s", i, m_suffix.data()), "", c_nwireCDC, -0.5, 14335.5);
372 hconstpw[i]->SetTitle(Form("wiregain const (%s); wire numbers;<dedxhit>", m_suffix.data()));
373 if (m_isMerge
374 && i == 0) hconstpw[i]->SetTitle(Form("merged wiregain rel-const (%s), avg = %0.03f; wire numbers;<dedxhit>", m_suffix.data(),
375 layeravg));
376
377 hconstpwvar[i] = new TH1D(Form("hconstpwvar_%s", m_suffix.data()), "", 400, -0.5, 2.5);
378 hconstpwvar[i]->SetTitle(Form("wiregain const (%s); wire gains; nentries", m_suffix.data()));
379 if (m_isMerge
380 && i == 0) hconstpwvar[i]->SetTitle(Form("merged wiregain rel-const (%s), avg = %0.03f; wire gains; nentries", m_suffix.data(),
381 layeravg));
382
383 for (unsigned int iw = 0; iw < c_nwireCDC; iw++) {
384
385 double gain = vdedx_means.at(iw);
386 if (m_isMerge && i == 1) gain = vrel_mean.at(iw);
387 hconstpw[i]->SetBinContent(iw + 1, gain);
388 hconstpwvar[i]->Fill(gain);
389
390 if (iw % 500 == 0) hconstpw[i]->GetXaxis()->SetBinLabel(iw + 1, Form("w%d", iw + 1));
391 }
392
393 hconstpw[i]->SetLineColor(i * 2 + 2);
394 hconstpw[i]->LabelsOption("u", "X");
395 hconstpw[i]->GetYaxis()->SetRangeUser(-0.1, 3.5);
396 hconstpw[i]->LabelsDeflate();
397
398 hconstpwvar[i]->SetFillColor(i * 2 + 2);
399 hconstpwvar[i]->Scale(1 / hconstpwvar[i]->GetMaximum());
400
401 }
402 cwconst.cd();
403 cwconst.SetGridy(1);
404 hconstpw[0]->Draw("");
405 if (m_isMerge) hconstpw[1]->Draw("same");
406
407 cwconstvar.cd();
408 hconstpwvar[0]->Draw("hist");
409 if (m_isMerge) hconstpwvar[1]->Draw("hist same");
410
411 cwconst.SaveAs(Form("cdcdedx_wgcal_wireconst_%s.pdf", m_suffix.data()));
412 cwconstvar.SaveAs(Form("cdcdedx_wgcal_wireconstvar_%s.pdf", m_suffix.data()));
413}
414
415//------------------------------------
416void CDCDedxWireGainAlgorithm::plotLayerGain(const vector<double>& layermean, double layeravg)
417{
418
419 TH1D hlayeravg(Form("hlayeravg_%s", m_suffix.data()), "", layermean.size(), -0.5, 55.5);
420 hlayeravg.SetTitle(Form("layer gain avg (%s); layer numbers;<dedxhit>", m_suffix.data()));
421
422 for (unsigned int il = 0; il < layermean.size(); il++) {
423 hlayeravg.SetBinContent(il + 1, layermean[il]);
424 if (il % 2 == 0 || il == layermean.size() - 1) hlayeravg.GetXaxis()->SetBinLabel(il + 1, Form("L%d", il));
425 }
426
427 TCanvas clayeravg("clayeravg", "clayeravg", 800, 500);
428 clayeravg.SetGridy(1);
429 clayeravg.cd();
430 gStyle->SetOptStat("ne");
431 hlayeravg.LabelsOption("u", "X");
432 hlayeravg.SetLineColor(kBlue);
433 hlayeravg.GetYaxis()->SetRangeUser(-0.1, 3.5);
434 hlayeravg.SetTitle(Form("%s, avg = %0.03f (abs)", hlayeravg.GetTitle(), layeravg));
435 hlayeravg.LabelsDeflate();
436 hlayeravg.Draw("");
437 TLine* tl = new TLine(-0.5, layeravg, 55.5, layeravg);
438 tl->SetLineColor(kRed);
439 tl->DrawClone("same");
440 clayeravg.SaveAs(Form("cdcdedx_wgcal_layeravg_%s.pdf", m_suffix.data()));
441 delete tl;
442}
443
444//------------------------------------
445void CDCDedxWireGainAlgorithm::plotWGPerLayer(const vector<double>& vdedx_means, const vector<double>& layermean, double layeravg)
446{
447
449
450 TCanvas clconst("clconst", "", 800, 500);
451 clconst.Divide(2, 2);
452 clconst.SetBatch(kTRUE);
453
454 stringstream psnameL;
455 psnameL << Form("cdcdedx_wgcal_layerconst_%s.pdf[", m_suffix.data());
456 clconst.Print(psnameL.str().c_str());
457 psnameL.str(""); psnameL << Form("cdcdedx_wgcal_layerconst_%s.pdf", m_suffix.data());
458
459 int jwire = 0;
460
461 for (unsigned int il = 0; il < layermean.size(); ++il) {
462
463 unsigned int nwires = cdcgeo.nWiresInLayer(il);
464 TH1D hconstpl(Form("hconstpwvar_l%d_%s", il, m_suffix.data()), "", nwires, jwire, jwire + nwires);
465 hconstpl.SetTitle(Form("abs-const, layer: %d (%s); wire numbers;<dedxhit>", il, m_suffix.data()));
466
467 for (unsigned int iw = 0; iw < nwires; ++iw) {
468 hconstpl.SetBinContent(iw + 1, vdedx_means.at(jwire));
469 if (il < 32) {
470 if (iw % 10 == 0) hconstpl.GetXaxis()->SetBinLabel(iw + 1, Form("w%d", jwire));
471 } else {
472 if (iw % 15 == 0) hconstpl.GetXaxis()->SetBinLabel(iw + 1, Form("w%d", jwire));
473 }
474 jwire++;
475 }
476
477 double lmean = layermean.at(il) / layeravg;
478
479 clconst.cd(il % 4 + 1);
480 gStyle->SetOptStat("ne");
481
482 hconstpl.SetTitle(Form("%s, avg = %0.03f", hconstpl.GetTitle(), lmean));
483
484 if (il < 8) hconstpl.GetYaxis()->SetRangeUser(-0.1, 4.0);
485 else hconstpl.GetYaxis()->SetRangeUser(-0.1, 2.0);
486 hconstpl.SetFillColor(kAzure - 1);
487 hconstpl.LabelsOption("u", "X");
488 hconstpl.DrawCopy("hist");
489
490 TLine* tlc = new TLine(jwire - nwires, lmean, jwire, lmean);
491 tlc->SetLineColor(kRed);
492 tlc->DrawClone("same");
493
494 if ((il + 1) % 4 == 0 || (il + 1) == layermean.size()) {
495 clconst.Print(psnameL.str().c_str());
496 clconst.Clear("D");
497 }
498
499 delete tlc;
500 }
501
502 psnameL.str("");
503 psnameL << Form("cdcdedx_wgcal_layerconst_%s.pdf]", m_suffix.data());
504 clconst.Print(psnameL.str().c_str());
505}
506
507//--------------------------
509{
510
511 TCanvas cstats("cstats", "cstats", 800, 400);
512 cstats.SetBatch(kTRUE);
513 cstats.Divide(2, 1);
514
515 cstats.cd(1);
516 auto hestats = getObjectPtr<TH1I>("hestats");
517 if (hestats) {
518 hestats->SetName(Form("hestats_%s", m_suffix.data()));
519 hestats->SetStats(0);
520 hestats->DrawCopy("");
521 }
522
523 cstats.cd(2);
524 auto htstats = getObjectPtr<TH1I>("htstats");
525 if (htstats) {
526 htstats->SetName(Form("htstats_%s", m_suffix.data()));
527 htstats->SetStats(0);
528 htstats->DrawCopy("");
529 }
530
531 cstats.Print(Form("cdcdedx_wgcal_stats_%s.pdf", m_suffix.data()));
532}
void plotWGPerLayer(const std::vector< double > &vdedx_means, const std::vector< double > &layermean, double layeravg)
function to draw WG per layer
void plotLayerGain(const std::vector< double > &layermean, double layeravg)
function to draw layer gains
void plotLayerDist(std::array< TH1D *, 2 > hdedxL)
function to draw dE/dx for inner/outer layer
bool m_isMerge
merge payload at the time of calibration
double m_truncMax
max trunc range for mean
double m_truncMin
min trunc range for mean
void getTruncatedBins(TH1D *hdedxhit, unsigned int &binlow, unsigned int &binhigh)
function to get bins of truncation from histogram
void plotWireGain(const std::vector< double > &vdedx_means, const std::vector< double > &vrel_mean, double layeravg)
function to draw wire gains
void getExpRunInfo()
function to get extract calibration run/exp
void createPayload(const std::vector< double > &vdedx_tmeans)
function to finally store new payload after full calibration
void plotWireDist(const std::vector< TH1D * > &hist, const std::vector< double > &vrel_mean)
function to draw dE/dx histograms for each wire
bool m_isMakePlots
Save arithmetic and truncated mean for the 'dedx' values.
DBObjPtr< CDCDedxBadWires > m_DBBadWires
Bad wire DB object.
std::string m_suffix
suffix string to separate plots
double getTruncationMean(TH1D *hdedxhit, int binlow, int binhigh)
function to get mean of truncation from histogram
int m_dedxBins
number of bins for dedx histogram
DBObjPtr< CDCGeometry > m_cdcGeo
Geometry of CDC.
CDCDedxWireGainAlgorithm()
Constructor: Sets the description, the properties and the parameters of the algorithm.
DBObjPtr< CDCDedxWireGain > m_DBWireGains
Wire gain DB object.
void plotEventStats()
function to draw statstics
virtual EResult calibrate() override
Wire gain algorithm.
double m_dedxMax
max dedx range for wiregain cal
double m_dedxMin
min dedx range for wiregain cal
bool m_isWireTruc
method of trunc range for mean
dE/dx wire gain calibration constants
The Class for CDC Geometry Parameters.
unsigned nWiresInLayer(int layerId) const
Returns wire numbers in a layer.
static CDCGeometryPar & Instance(const CDCGeometry *=nullptr)
Static method to get a reference to the CDCGeometryPar instance.
void saveCalibration(TClonesArray *data, const std::string &name)
Store DBArray payload with given name with default IOV.
void updateDBObjPtrs(const unsigned int event, const int run, const int experiment)
Updates any DBObjPtrs by calling update(event) for DBStore.
void setDescription(const std::string &description)
Set algorithm description (in constructor)
const std::vector< Calibration::ExpRun > & getRunList() const
Get the list of runs for which calibration is called.
EResult
The result of calibration.
@ c_OK
Finished successfully =0 in Python.
@ c_NotEnoughData
Needs more data =2 in Python.
CalibrationAlgorithm(const std::string &collectorModuleName)
Constructor - sets the prefix for collected objects (won't be accesses until execute(....
std::shared_ptr< T > getObjectPtr(const std::string &name, const std::vector< Calibration::ExpRun > &requestedRuns)
Get calibration data object by name and list of runs, the Merge function will be called to generate t...
Abstract base class for different kinds of events.
STL namespace.