8 #include <ecl/dbobjects/ECLDatabaseImporter.h>
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>
28 #include <framework/database/IntervalOfValidity.h>
29 #include <framework/database/Database.h>
30 #include <framework/database/DBImportArray.h>
31 #include <framework/database/DBImportObjPtr.h>
36 #include <TClonesArray.h>
44 ECLDatabaseImporter::ECLDatabaseImporter(vector<string> inputFileNames,
const std::string& name)
47 for (
auto& inputFileName : inputFileNames)
48 m_inputFileNames.push_back(inputFileName);
54 void ECLDatabaseImporter::importDigitEnergyCalibration()
57 TClonesArray digitCalibrationConstants(
"Belle2::ECLDigitEnergyConstants");
63 for (
const string& inputFile : m_inputFileNames) {
65 TFile* f = TFile::Open(inputFile.c_str(),
"READ");
67 TIter next(f->GetListOfKeys());
69 while ((key = (TKey*) next())) {
71 string histconstants = key->GetName();
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());
79 else { B2FATAL(
"Key name does not match any of the following: energy, amplitude!"); }
85 if (nFiles != 1) { B2FATAL(
"Sorry, you must only import one file at a time for now!"); }
89 for (
int bin = 1; bin <= amplitude->GetNbinsX(); ++bin) {
90 float amplitudeval = amplitude->GetBinContent(bin);
91 float energyval = energy->GetBinContent(bin);
97 Database::Instance().storeData(m_name, &digitCalibrationConstants, iov);
100 void ECLDatabaseImporter::importDigitTimeCalibration()
103 TClonesArray digitCalibrationConstants(
"Belle2::ECLDigitTimeConstants");
108 for (
const string& inputFile : m_inputFileNames) {
110 TFile* f = TFile::Open(inputFile.c_str(),
"READ");
112 TIter next(f->GetListOfKeys());
114 while ((key = (TKey*) next())) {
116 string histconstants = key->GetName();
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!"); }
126 if (nFiles != 1) { B2FATAL(
"Sorry, you must only import one file at a time for now!"); }
130 for (
int bin = 1; bin <= offset->GetNbinsX(); ++bin) {
131 float offsetval = offset->GetBinContent(bin);
137 Database::Instance().storeData(m_name, &digitCalibrationConstants, iov);
140 void ECLDatabaseImporter::importShowerCorrectorLeakageCorrections()
142 if (m_inputFileNames.size() > 1)
143 B2FATAL(
"Sorry, you must only import one file at a time for now!");
146 TFile* inputFile =
new TFile(m_inputFileNames[0].data(),
"READ");
148 if (!inputFile || inputFile->IsZombie())
149 B2FATAL(
"Could not open file " << m_inputFileNames[0]);
152 TTree* correctionTree = getRootObjectFromFile<TTree*>(inputFile,
"ParameterNtuple");
153 TTree* helperTree = getRootObjectFromFile<TTree*>(inputFile,
"ConstantNtuple");
159 int bgFractionBinNum;
164 float correctionFactor;
167 correctionTree->SetBranchAddress(m_bgFractionBinNumBranchName.c_str(), &bgFractionBinNum);
168 correctionTree->SetBranchAddress(m_regNumBranchName.c_str(), ®Num);
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);
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;
182 for (
long iEntry = 0; iEntry < correctionTree->GetEntries(); ++iEntry) {
183 correctionTree->GetEntry(iEntry);
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);
207 int numOfReg1ThetaBins;
208 int numOfReg2ThetaBins;
209 int numOfReg3ThetaBins;
216 std::vector<float> avgRecEns(m_numAvgRecEnEntries);
217 helperTree->SetBranchAddress(m_avgRecEnBranchName.c_str(), avgRecEns.data());
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);
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;
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]);
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);
288 m_numOfReg1ThetaBins,
289 m_numOfReg2ThetaBins,
290 m_numOfReg3ThetaBins,
307 void ECLDatabaseImporter::importShowerShapesSecondMomentCorrections()
309 if (m_inputFileNames.size() > 1)
310 B2FATAL(
"Sorry, you must only import one file at a time for now!");
314 TFile* inputFile =
new TFile(m_inputFileNames[0].data(),
"READ");
316 if (!inputFile || inputFile->IsZombie())
317 B2FATAL(
"Could not open file " << m_inputFileNames[0]);
320 TGraph* theta_N1_graph = getRootObjectFromFile<TGraph*>(inputFile,
"SecondMomentCorrections_theta_N1");
321 dbArray.
appendNew(ECLShower::c_nPhotons, ECLShowerShapeModule::c_thetaType, *theta_N1_graph);
324 TGraph* phi_N1_graph = getRootObjectFromFile<TGraph*>(inputFile,
"SecondMomentCorrections_phi_N1");
325 dbArray.
appendNew(ECLShower::c_nPhotons, ECLShowerShapeModule::c_phiType, *phi_N1_graph);
328 TGraph* theta_N2_graph = getRootObjectFromFile<TGraph*>(inputFile,
"SecondMomentCorrections_theta_N2");
329 dbArray.
appendNew(ECLShower::c_neutralHadron, ECLShowerShapeModule::c_thetaType, *theta_N2_graph);
332 TGraph* phi_N2_graph = getRootObjectFromFile<TGraph*>(inputFile,
"SecondMomentCorrections_phi_N2");
333 dbArray.
appendNew(ECLShower::c_neutralHadron, ECLShowerShapeModule::c_phiType, *phi_N2_graph);
349 void ECLDatabaseImporter::importShowerEnergyCorrectionTemporary()
351 if (m_inputFileNames.size() > 1)
352 B2FATAL(
"Sorry, you must only import one file at a time for now!");
355 std::filesystem::path path(m_inputFileNames[0]);
356 if (path.extension() !=
".root")
357 B2FATAL(
"Expecting a .root file. Aborting");
359 TFile* inputFile =
new TFile(m_inputFileNames[0].data(),
"READ");
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");
367 double bkgFactor = bg_histo->GetBinContent(1);
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();
374 double energyMin = theta_en_graph->GetYmin();
375 double energyMax = theta_en_graph->GetYmax();
378 B2DEBUG(28,
"Leakage DBobjects angle boundaries: thetaMin=" << thetaMin <<
" thetaMax=" << thetaMax <<
" phiMin= " << phiMin <<
379 " phiMax= " << phiMax <<
" enmin=" << energyMin <<
380 " enmax=" << energyMax);
389 if (std::abs(bkgFactor - 1.0) < 1e-9) {
392 dbPtr_theta_geo.
construct(*theta_geo_graph, thetaMin, thetaMax, energyMin, energyMax);
394 dbPtr_phi_geo.
construct(*phi_geo_graph, phiMin, phiMax, energyMin, energyMax);
396 dbPtr_theta_en.
construct(*theta_en_graph, thetaMin, thetaMax, energyMin, energyMax);
398 dbPtr_phi_en.
construct(*phi_en_graph, phiMin, phiMax, energyMin, energyMax);
401 dbPtr_theta_geo.
import(iov);
402 dbPtr_phi_geo.
import(iov);
403 dbPtr_theta_en.
import(iov);
406 if (std::abs(bkgFactor - 1.0) < 1e-9) {
409 dbPtr_theta_geo.
construct(*theta_geo_graph, thetaMin, thetaMax, energyMin, energyMax);
411 dbPtr_phi_geo.
construct(*phi_geo_graph, phiMin, phiMax, energyMin, energyMax);
413 dbPtr_theta_en.
construct(*theta_en_graph, thetaMin, thetaMax, energyMin, energyMax);
415 dbPtr_phi_en.
construct(*phi_en_graph, phiMin, phiMax, energyMin, energyMax);
418 dbPtr_theta_geo.
import(iov);
419 dbPtr_phi_geo.
import(iov);
420 dbPtr_theta_en.
import(iov);
427 void ECLDatabaseImporter::importTrackClusterMatchingThresholds()
429 if (m_inputFileNames.size() > 1)
430 B2FATAL(
"Sorry, you must only import one file at a time for now!");
433 std::filesystem::path path(m_inputFileNames[0]);
434 if (path.extension() !=
".txt")
435 B2FATAL(
"Expecting a .txt file. Aborting");
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;
445 ifstream infile(m_inputFileNames[0]);
447 while (getline(infile, line)) {
448 istringstream iss(line);
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);
464 dbPtr.
construct(m_matchingThresholdPairsFWD, m_matchingThresholdPairsBWD, m_matchingThresholdPairsBRL);
472 void ECLDatabaseImporter::importTrackClusterMatchingParameterizations()
474 if (m_inputFileNames.size() > 1)
475 B2FATAL(
"Sorry, you must only import one file at a time for now!");
478 TFile* inputFile =
new TFile(m_inputFileNames[0].data(),
"READ");
480 if (!inputFile || inputFile->IsZombie())
481 B2FATAL(
"Could not open file " << m_inputFileNames[0]);
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"};
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))));
498 dbPtr.
construct(m_parametrizationFunctions);
Class for importing array of objects to the database.
T * appendNew()
Construct a new T object at the end of the array.
bool import(const IntervalOfValidity &iov)
Import the object to database.
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.