Belle II Software  release-08-01-10
ECLDatabaseImporter.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 #include <ecl/dbobjects/ECLDatabaseImporter.h>
9 
10 //STL
11 #include <sstream>
12 #include <fstream>
13 #include <string>
14 #include <filesystem>
15 
16 // ECL
17 #include <ecl/dbobjects/ECLDigitEnergyConstants.h>
18 #include <ecl/dbobjects/ECLDigitTimeConstants.h>
19 #include <ecl/modules/eclShowerShape/ECLShowerShapeModule.h>
20 #include <ecl/dbobjects/ECLShowerShapeSecondMomentCorrection.h>
21 #include <ecl/dbobjects/ECLShowerCorrectorLeakageCorrection.h>
22 #include <ecl/dbobjects/ECLShowerEnergyCorrectionTemporary.h>
23 #include <ecl/dbobjects/ECLTrackClusterMatchingParameterizations.h>
24 #include <ecl/dbobjects/ECLTrackClusterMatchingThresholds.h>
25 #include <ecl/dataobjects/ECLShower.h>
26 
27 // FRAMEWORK
28 #include <framework/database/IntervalOfValidity.h>
29 #include <framework/database/Database.h>
30 #include <framework/database/DBImportArray.h>
31 #include <framework/database/DBImportObjPtr.h>
32 
33 // ROOT
34 #include <TH1.h>
35 #include <TKey.h>
36 #include <TClonesArray.h>
37 #include <TTree.h>
38 #include <TGraph2D.h>
39 
40 // NAMESPACES
41 using namespace std;
42 using namespace Belle2;
43 
44 ECLDatabaseImporter::ECLDatabaseImporter(vector<string> inputFileNames, const std::string& name)
45 {
46  //input file names
47  for (auto& inputFileName : inputFileNames)
48  m_inputFileNames.push_back(inputFileName);
49 
50  //output file name
51  m_name = name;
52 }
53 
54 void ECLDatabaseImporter::importDigitEnergyCalibration()
55 {
56 
57  TClonesArray digitCalibrationConstants("Belle2::ECLDigitEnergyConstants");
58 
59  TH1F* energy = 0;
60  TH1F* amplitude = 0;
61  int nFiles = 0;
62 
63  for (const string& inputFile : m_inputFileNames) {
64 
65  TFile* f = TFile::Open(inputFile.c_str(), "READ");
66 
67  TIter next(f->GetListOfKeys());
68  TKey* key;
69  while ((key = (TKey*) next())) {
70 
71  string histconstants = key->GetName();
72 
73  if (histconstants.compare("energy") == 0) {
74  energy = (TH1F*)f->Get(histconstants.c_str());
75  } else if (histconstants.compare("amplitude") == 0) {
76  amplitude = (TH1F*)f->Get(histconstants.c_str());
77  }
78 
79  else { B2FATAL("Key name does not match any of the following: energy, amplitude!"); }
80  }
81 
82  nFiles++;
83  }
84 
85  if (nFiles != 1) { B2FATAL("Sorry, you must only import one file at a time for now!"); }
86 
87  // loop over the histogram to fill the TClonesArray
88  int cell = 0;
89  for (int bin = 1; bin <= amplitude->GetNbinsX(); ++bin) {
90  float amplitudeval = amplitude->GetBinContent(bin);
91  float energyval = energy->GetBinContent(bin);
92  new (digitCalibrationConstants[cell]) ECLDigitEnergyConstants(bin, amplitudeval, energyval);
93  cell++;
94  }
95 
96  IntervalOfValidity iov(0, 0, -1, -1); // IOV (0,0,-1,-1) is valid for all runs and experiments
97  Database::Instance().storeData(m_name, &digitCalibrationConstants, iov);
98 }
99 
100 void ECLDatabaseImporter::importDigitTimeCalibration()
101 {
102 
103  TClonesArray digitCalibrationConstants("Belle2::ECLDigitTimeConstants");
104 
105  TH1F* offset = 0;
106  int nFiles = 0;
107 
108  for (const string& inputFile : m_inputFileNames) {
109 
110  TFile* f = TFile::Open(inputFile.c_str(), "READ");
111 
112  TIter next(f->GetListOfKeys());
113  TKey* key;
114  while ((key = (TKey*) next())) {
115 
116  string histconstants = key->GetName();
117 
118  if (histconstants.compare("constantB") == 0) {
119  offset = (TH1F*)f->Get(histconstants.c_str());
120  } else { B2FATAL("Key name does not match any of the following: constantC!"); }
121  }
122 
123  nFiles++;
124  }
125 
126  if (nFiles != 1) { B2FATAL("Sorry, you must only import one file at a time for now!"); }
127 
128  // loop over the histogram to fill the TClonesArray
129  int cell = 0;
130  for (int bin = 1; bin <= offset->GetNbinsX(); ++bin) {
131  float offsetval = offset->GetBinContent(bin);
132  new (digitCalibrationConstants[cell]) ECLDigitTimeConstants(bin, offsetval);
133  cell++;
134  }
135 
136  IntervalOfValidity iov(0, 0, -1, -1); // IOV (0,0,-1,-1) is valid for all runs and experiments
137  Database::Instance().storeData(m_name, &digitCalibrationConstants, iov);
138 }
139 
140 void ECLDatabaseImporter::importShowerCorrectorLeakageCorrections()
141 {
142  if (m_inputFileNames.size() > 1)
143  B2FATAL("Sorry, you must only import one file at a time for now!");
144 
145  //Open file
146  TFile* inputFile = new TFile(m_inputFileNames[0].data(), "READ");
147 
148  if (!inputFile || inputFile->IsZombie())
149  B2FATAL("Could not open file " << m_inputFileNames[0]);
150 
151  //Get trees
152  TTree* correctionTree = getRootObjectFromFile<TTree*>(inputFile, "ParameterNtuple");
153  TTree* helperTree = getRootObjectFromFile<TTree*>(inputFile, "ConstantNtuple");
154 
155  //----------------------------------------------------------------------------------------------
156  //Fill ParameterNtuple vectors
157  //----------------------------------------------------------------------------------------------
158 
159  int bgFractionBinNum;
160  int regNum;
161  int phiBinNum;
162  int thetaBinNum;
163  int energyBinNum;
164  float correctionFactor;
165 
166  //Set Branch Addresses
167  correctionTree->SetBranchAddress(m_bgFractionBinNumBranchName.c_str(), &bgFractionBinNum);
168  correctionTree->SetBranchAddress(m_regNumBranchName.c_str(), &regNum);
169  correctionTree->SetBranchAddress(m_phiBinNumBranchName.c_str(), &phiBinNum);
170  correctionTree->SetBranchAddress(m_thetaBinNumBranchName.c_str(), &thetaBinNum);
171  correctionTree->SetBranchAddress(m_energyBinNumBranchName.c_str(), &energyBinNum);
172  correctionTree->SetBranchAddress(m_correctionFactorBranchName.c_str(), &correctionFactor);
173 
174  //Fill vectors
175  std::vector<int> m_bgFractionBinNum;
176  std::vector<int> m_regNum;
177  std::vector<int> m_phiBinNum;
178  std::vector<int> m_thetaBinNum;
179  std::vector<int> m_energyBinNum;
180  std::vector<float> m_correctionFactor;
181 
182  for (long iEntry = 0; iEntry < correctionTree->GetEntries(); ++iEntry) {
183  correctionTree->GetEntry(iEntry);
184 
185  m_bgFractionBinNum.push_back(bgFractionBinNum);
186  m_regNum.push_back(regNum);
187  m_phiBinNum.push_back(phiBinNum);
188  m_thetaBinNum.push_back(thetaBinNum);
189  m_energyBinNum.push_back(energyBinNum);
190  m_correctionFactor.push_back(correctionFactor);
191  }
192 
193  //----------------------------------------------------------------------------------------------
194  //Fill ConstantNtuple vectors
195  //----------------------------------------------------------------------------------------------
196 
197 
198  float lReg1Theta;
199  float hReg1Theta;
200  float lReg2Theta;
201  float hReg2Theta;
202  float lReg3Theta;
203  float hReg3Theta;
204  int numOfBfBins;
205  int numOfEnergyBins;
206  int numOfPhiBins;
207  int numOfReg1ThetaBins;
208  int numOfReg2ThetaBins;
209  int numOfReg3ThetaBins;
210  int phiPeriodicity;
211 
212  //Ugly hack to circumvent 'stack usage might be unbounded [-Wstack-usage=]' compiler warning that's caused by the use of c-type arrays.
213  //This is not for the faint of heart
214 
215  //because root GetEntry fills the whole internal array of the vector without changing it's size, we must ensure that it's the right size.
216  std::vector<float> avgRecEns(m_numAvgRecEnEntries);
217  helperTree->SetBranchAddress(m_avgRecEnBranchName.c_str(), avgRecEns.data()); //Read c-style array right into internal vector array.
218 
219  helperTree->SetBranchAddress(m_lReg1ThetaBranchName.c_str(), &lReg1Theta);
220  helperTree->SetBranchAddress(m_hReg1ThetaBranchName.c_str(), &hReg1Theta);
221  helperTree->SetBranchAddress(m_lReg2ThetaBranchName.c_str(), &lReg2Theta);
222  helperTree->SetBranchAddress(m_hReg2ThetaBranchName.c_str(), &hReg2Theta);
223  helperTree->SetBranchAddress(m_lReg3ThetaBranchName.c_str(), &lReg3Theta);
224  helperTree->SetBranchAddress(m_hReg3ThetaBranchName.c_str(), &hReg3Theta);
225  helperTree->SetBranchAddress(m_numOfBfBinsBranchName.c_str(), &numOfBfBins);
226  helperTree->SetBranchAddress(m_numOfEnergyBinsBranchName.c_str(), &numOfEnergyBins);
227  helperTree->SetBranchAddress(m_numOfPhiBinsBranchName.c_str(), &numOfPhiBins);
228  helperTree->SetBranchAddress(m_numOfReg1ThetaBinsBranchName.c_str(), &numOfReg1ThetaBins);
229  helperTree->SetBranchAddress(m_numOfReg2ThetaBinsBranchName.c_str(), &numOfReg2ThetaBins);
230  helperTree->SetBranchAddress(m_numOfReg3ThetaBinsBranchName.c_str(), &numOfReg3ThetaBins);
231  helperTree->SetBranchAddress(m_phiPeriodicityBranchName.c_str(), &phiPeriodicity);
232 
233  //Fill vectors
234  std::vector<float> m_avgRecEn;
235  std::vector<float> m_lReg1Theta;
236  std::vector<float> m_hReg1Theta;
237  std::vector<float> m_lReg2Theta;
238  std::vector<float> m_hReg2Theta;
239  std::vector<float> m_lReg3Theta;
240  std::vector<float> m_hReg3Theta;
241  std::vector<int> m_numOfBfBins;
242  std::vector<int> m_numOfEnergyBins;
243  std::vector<int> m_numOfPhiBins;
244  std::vector<int> m_numOfReg1ThetaBins;
245  std::vector<int> m_numOfReg2ThetaBins;
246  std::vector<int> m_numOfReg3ThetaBins;
247  std::vector<int> m_phiPeriodicity;
248 
249  for (long iEntry = 0; iEntry < helperTree->GetEntries(); ++iEntry) {
250  helperTree->GetEntry(iEntry);
251  for (unsigned int iIdx = 0; iIdx < avgRecEns.size(); ++iIdx) m_avgRecEn.push_back(avgRecEns[iIdx]);
252 
253  m_lReg1Theta.push_back(lReg1Theta);
254  m_hReg1Theta.push_back(hReg1Theta);
255  m_lReg2Theta.push_back(lReg2Theta);
256  m_hReg2Theta.push_back(hReg2Theta);
257  m_lReg3Theta.push_back(lReg3Theta);
258  m_hReg3Theta.push_back(hReg3Theta);
259  m_numOfBfBins.push_back(numOfBfBins);
260  m_numOfEnergyBins.push_back(numOfEnergyBins);
261  m_numOfPhiBins.push_back(numOfPhiBins);
262  m_numOfReg1ThetaBins.push_back(numOfReg1ThetaBins);
263  m_numOfReg2ThetaBins.push_back(numOfReg2ThetaBins);
264  m_numOfReg3ThetaBins.push_back(numOfReg3ThetaBins);
265  m_phiPeriodicity.push_back(phiPeriodicity);
266  }
267 
268  //----------------------------------------------------------------------------------------------
269 
270  //Construct DB object
271  DBImportObjPtr<ECLShowerCorrectorLeakageCorrection> dbPtr("ecl_shower_corrector_leakage_corrections");
272  dbPtr.construct(m_bgFractionBinNum,
273  m_regNum,
274  m_phiBinNum,
275  m_thetaBinNum,
276  m_energyBinNum,
277  m_correctionFactor,
278  m_avgRecEn,
279  m_lReg1Theta,
280  m_hReg1Theta,
281  m_lReg2Theta,
282  m_hReg2Theta,
283  m_lReg3Theta,
284  m_hReg3Theta,
285  m_numOfBfBins,
286  m_numOfEnergyBins,
287  m_numOfPhiBins,
288  m_numOfReg1ThetaBins,
289  m_numOfReg2ThetaBins,
290  m_numOfReg3ThetaBins,
291  m_phiPeriodicity);
292 
293  //Create IOV object
294  int startExp = 0;
295  int startRun = 0;
296  int endExp = -1;
297  int endRun = -1;
298  IntervalOfValidity iov(startExp, startRun, endExp, endRun);
299 
300  //Import into local db
301  dbPtr.import(iov);
302 
303  delete inputFile;
304 
305 }
306 
307 void ECLDatabaseImporter::importShowerShapesSecondMomentCorrections()
308 {
309  if (m_inputFileNames.size() > 1)
310  B2FATAL("Sorry, you must only import one file at a time for now!");
311 
312  //Extract TGraphs from file
313  DBImportArray<ECLShowerShapeSecondMomentCorrection> dbArray("ecl_shower_shape_second_moment_corrections");
314  TFile* inputFile = new TFile(m_inputFileNames[0].data(), "READ");
315 
316  if (!inputFile || inputFile->IsZombie())
317  B2FATAL("Could not open file " << m_inputFileNames[0]);
318 
319  //N1 theta
320  TGraph* theta_N1_graph = getRootObjectFromFile<TGraph*>(inputFile, "SecondMomentCorrections_theta_N1");
321  dbArray.appendNew(ECLShower::c_nPhotons, ECLShowerShapeModule::c_thetaType, *theta_N1_graph);
322 
323  //N1 phi
324  TGraph* phi_N1_graph = getRootObjectFromFile<TGraph*>(inputFile, "SecondMomentCorrections_phi_N1");
325  dbArray.appendNew(ECLShower::c_nPhotons, ECLShowerShapeModule::c_phiType, *phi_N1_graph);
326 
327  //N2 theta
328  TGraph* theta_N2_graph = getRootObjectFromFile<TGraph*>(inputFile, "SecondMomentCorrections_theta_N2");
329  dbArray.appendNew(ECLShower::c_neutralHadron, ECLShowerShapeModule::c_thetaType, *theta_N2_graph);
330 
331  //N2 phi
332  TGraph* phi_N2_graph = getRootObjectFromFile<TGraph*>(inputFile, "SecondMomentCorrections_phi_N2");
333  dbArray.appendNew(ECLShower::c_neutralHadron, ECLShowerShapeModule::c_phiType, *phi_N2_graph);
334 
335 
336  //Import to DB
337  int startExp = 0;
338  int startRun = 0;
339  int endExp = -1;
340  int endRun = -1;
341  IntervalOfValidity iov(startExp, startRun, endExp, endRun);
342 
343  //Import into local db
344  dbArray.import(iov);
345 
346  delete inputFile;
347 }
348 
349 void ECLDatabaseImporter::importShowerEnergyCorrectionTemporary()
350 {
351  if (m_inputFileNames.size() > 1)
352  B2FATAL("Sorry, you must only import one file at a time for now!");
353 
354  //Expect a root file
355  std::filesystem::path path(m_inputFileNames[0]);
356  if (path.extension() != ".root")
357  B2FATAL("Expecting a .root file. Aborting");
358 
359  TFile* inputFile = new TFile(m_inputFileNames[0].data(), "READ");
360 
361  TGraph2D* theta_geo_graph = getRootObjectFromFile<TGraph2D*>(inputFile, "LeakageCorrections_theta_geometry");
362  TGraph2D* phi_geo_graph = getRootObjectFromFile<TGraph2D*>(inputFile, "LeakageCorrections_phi_geometry");
363  TGraph2D* theta_en_graph = getRootObjectFromFile<TGraph2D*>(inputFile, "LeakageCorrections_theta_energy");
364  TGraph2D* phi_en_graph = getRootObjectFromFile<TGraph2D*>(inputFile, "LeakageCorrections_phi_energy");
365  TH1F* bg_histo = getRootObjectFromFile<TH1F*>(inputFile, "LeakageCorrections_background_fraction");
366 
367  double bkgFactor = bg_histo->GetBinContent(1);
368 
369  double thetaMin = theta_en_graph->GetXmin();
370  double thetaMax = theta_en_graph->GetXmax();
371  double phiMin = phi_en_graph->GetXmin();
372  double phiMax = phi_en_graph->GetXmax();
373 
374  double energyMin = theta_en_graph->GetYmin();
375  double energyMax = theta_en_graph->GetYmax();
376 
377 
378  B2DEBUG(28, "Leakage DBobjects angle boundaries: thetaMin=" << thetaMin << " thetaMax=" << thetaMax << " phiMin= " << phiMin <<
379  " phiMax= " << phiMax << " enmin=" << energyMin <<
380  " enmax=" << energyMax);
381 
382  // Import to DB
383  int startExp = 0;
384  int startRun = 0;
385  int endExp = -1;
386  int endRun = -1;
387  IntervalOfValidity iov(startExp, startRun, endExp, endRun);
388 
389  if (std::abs(bkgFactor - 1.0) < 1e-9) { //bkgFactor == 1 -> phase 2 backgrounds
390 
391  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_theta_geo("ECLLeakageCorrection_thetaGeometry_phase2");
392  dbPtr_theta_geo.construct(*theta_geo_graph, thetaMin, thetaMax, energyMin, energyMax);
393  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_phi_geo("ECLLeakageCorrection_phiGeometry_phase2");
394  dbPtr_phi_geo.construct(*phi_geo_graph, phiMin, phiMax, energyMin, energyMax);
395  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_theta_en("ECLLeakageCorrection_thetaEnergy_phase2");
396  dbPtr_theta_en.construct(*theta_en_graph, thetaMin, thetaMax, energyMin, energyMax);
397  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_phi_en("ECLLeakageCorrection_phiEnergy_phase2");
398  dbPtr_phi_en.construct(*phi_en_graph, phiMin, phiMax, energyMin, energyMax);
399 
400  //Import into local db
401  dbPtr_theta_geo.import(iov);
402  dbPtr_phi_geo.import(iov);
403  dbPtr_theta_en.import(iov);
404  dbPtr_phi_en.import(iov);
405  }
406  /*else (because currently phase_2 and phase_3 are same payload*/ if (std::abs(bkgFactor - 1.0) < 1e-9) {
407 
408  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_theta_geo("ECLLeakageCorrection_thetaGeometry_phase3");
409  dbPtr_theta_geo.construct(*theta_geo_graph, thetaMin, thetaMax, energyMin, energyMax);
410  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_phi_geo("ECLLeakageCorrection_phiGeometry_phase3");
411  dbPtr_phi_geo.construct(*phi_geo_graph, phiMin, phiMax, energyMin, energyMax);
412  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_theta_en("ECLLeakageCorrection_thetaEnergy_phase3");
413  dbPtr_theta_en.construct(*theta_en_graph, thetaMin, thetaMax, energyMin, energyMax);
414  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_phi_en("ECLLeakageCorrection_phiEnergy_phase3");
415  dbPtr_phi_en.construct(*phi_en_graph, phiMin, phiMax, energyMin, energyMax);
416 
417  //Import into local db
418  dbPtr_theta_geo.import(iov);
419  dbPtr_phi_geo.import(iov);
420  dbPtr_theta_en.import(iov);
421  dbPtr_phi_en.import(iov);
422  }
423 
424 
425 }
426 
427 void ECLDatabaseImporter::importTrackClusterMatchingThresholds()
428 {
429  if (m_inputFileNames.size() > 1)
430  B2FATAL("Sorry, you must only import one file at a time for now!");
431 
432  //Expect a txt file
433  std::filesystem::path path(m_inputFileNames[0]);
434  if (path.extension() != ".txt")
435  B2FATAL("Expecting a .txt file. Aborting");
436 
437  vector<pair<double, double>> m_matchingThresholdPairsFWD;
438  vector<pair<double, double>> m_matchingThresholdPairsBWD;
439  vector<pair<double, pair<double, double>>> m_matchingThresholdPairsBRL;
440  pair<double, double> m_matchingThresholdPair;
441  pair<double, pair<double, double>> m_thetaMatchingThresholdPair;
442  double pt, threshold, thetalimit;
443  string eclregion;
444 
445  ifstream infile(m_inputFileNames[0]);
446  string line;
447  while (getline(infile, line)) {
448  istringstream iss(line);
449  iss >> eclregion;
450  if (eclregion == "FWD" || eclregion == "BWD") {
451  iss >> pt >> threshold;
452  m_matchingThresholdPair = make_pair(pt, threshold);
453  if (eclregion == "FWD") m_matchingThresholdPairsFWD.push_back(m_matchingThresholdPair);
454  else m_matchingThresholdPairsBWD.push_back(m_matchingThresholdPair);
455  } else if (eclregion == "BRL") {
456  iss >> thetalimit >> pt >> threshold;
457  m_matchingThresholdPair = make_pair(pt, threshold);
458  m_thetaMatchingThresholdPair = make_pair(thetalimit, m_matchingThresholdPair);
459  m_matchingThresholdPairsBRL.push_back(m_thetaMatchingThresholdPair);
460  }
461  }
462 
463  DBImportObjPtr<ECLTrackClusterMatchingThresholds> dbPtr("ECLTrackClusterMatchingThresholds");
464  dbPtr.construct(m_matchingThresholdPairsFWD, m_matchingThresholdPairsBWD, m_matchingThresholdPairsBRL);
465 
466  IntervalOfValidity iov(0, 0, -1, -1);
467 
468  //Import into local db
469  dbPtr.import(iov);
470 }
471 
472 void ECLDatabaseImporter::importTrackClusterMatchingParameterizations()
473 {
474  if (m_inputFileNames.size() > 1)
475  B2FATAL("Sorry, you must only import one file at a time for now!");
476 
477  // Open file
478  TFile* inputFile = new TFile(m_inputFileNames[0].data(), "READ");
479 
480  if (!inputFile || inputFile->IsZombie())
481  B2FATAL("Could not open file " << m_inputFileNames[0]);
482 
483  map<string, TF1> m_parametrizationFunctions;
484  vector<string> angles = {"Theta", "Phi"};
485  vector<string> regions = {"BRL", "BWD", "FWD"};
486  vector<string> hittypes = {"CROSS", "DL", "NEAR"};
487 
488  for (const auto& angle : angles) {
489  for (const auto& region : regions) {
490  for (const auto& hittype : hittypes) {
491  m_parametrizationFunctions.insert(make_pair(angle + region + hittype, *(getRootObjectFromFile<TF1*>(inputFile,
492  "RMSParameterization" + angle + region + hittype))));
493  }
494  }
495  }
496 
497  DBImportObjPtr<ECLTrackClusterMatchingParameterizations> dbPtr("ECLTrackClusterMatchingParameterizations");
498  dbPtr.construct(m_parametrizationFunctions);
499 
500  IntervalOfValidity iov(0, 0, -1, -1);
501 
502  //Import into local db
503  dbPtr.import(iov);
504 
505  delete inputFile;
506 }
Class for importing array of objects to the database.
Definition: DBImportArray.h:25
T * appendNew()
Construct a new T object at the end of the array.
Definition: DBImportArray.h:67
bool import(const IntervalOfValidity &iov)
Import the object to database.
Definition: DBImportBase.cc:36
Class for importing a single object to the database.
void construct(Args &&... params)
Construct an object of type T in this DBImportObjPtr using the provided constructor arguments.
Energy calibration constants per digit.
Time and time resolution calibration constants per digit.
A class that describes the interval of experiments/runs for which an object in the database is valid.
Abstract base class for different kinds of events.