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