8 #include <ecl/dbobjects/ECLDatabaseImporter.h>
16 #include <boost/filesystem.hpp>
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>
30 #include <framework/database/IntervalOfValidity.h>
31 #include <framework/database/Database.h>
32 #include <framework/database/DBImportArray.h>
33 #include <framework/database/DBImportObjPtr.h>
38 #include <TClonesArray.h>
46 ECLDatabaseImporter::ECLDatabaseImporter(vector<string> inputFileNames,
const std::string& name)
49 for (
auto& inputFileName : inputFileNames)
50 m_inputFileNames.push_back(inputFileName);
56 void ECLDatabaseImporter::importDigitEnergyCalibration()
59 TClonesArray digitCalibrationConstants(
"Belle2::ECLDigitEnergyConstants");
65 for (
const string& inputFile : m_inputFileNames) {
67 TFile* f = TFile::Open(inputFile.c_str(),
"READ");
69 TIter next(f->GetListOfKeys());
71 while ((key = (TKey*) next())) {
73 string histconstants = key->GetName();
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());
81 else { B2FATAL(
"Key name does not match any of the following: energy, amplitude!"); }
87 if (nFiles != 1) { B2FATAL(
"Sorry, you must only import one file at a time for now!"); }
91 for (
int bin = 1; bin <= amplitude->GetNbinsX(); ++bin) {
92 float amplitudeval = amplitude->GetBinContent(bin);
93 float energyval = energy->GetBinContent(bin);
99 Database::Instance().storeData(m_name, &digitCalibrationConstants, iov);
102 void ECLDatabaseImporter::importDigitTimeCalibration()
105 TClonesArray digitCalibrationConstants(
"Belle2::ECLDigitTimeConstants");
110 for (
const string& inputFile : m_inputFileNames) {
112 TFile* f = TFile::Open(inputFile.c_str(),
"READ");
114 TIter next(f->GetListOfKeys());
116 while ((key = (TKey*) next())) {
118 string histconstants = key->GetName();
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!"); }
128 if (nFiles != 1) { B2FATAL(
"Sorry, you must only import one file at a time for now!"); }
132 for (
int bin = 1; bin <= offset->GetNbinsX(); ++bin) {
133 float offsetval = offset->GetBinContent(bin);
139 Database::Instance().storeData(m_name, &digitCalibrationConstants, iov);
142 void ECLDatabaseImporter::importShowerCorrectorLeakageCorrections()
144 if (m_inputFileNames.size() > 1)
145 B2FATAL(
"Sorry, you must only import one file at a time for now!");
148 TFile* inputFile =
new TFile(m_inputFileNames[0].data(),
"READ");
150 if (!inputFile || inputFile->IsZombie())
151 B2FATAL(
"Could not open file " << m_inputFileNames[0]);
154 TTree* correctionTree = getRootObjectFromFile<TTree*>(inputFile,
"ParameterNtuple");
155 TTree* helperTree = getRootObjectFromFile<TTree*>(inputFile,
"ConstantNtuple");
161 int bgFractionBinNum;
166 float correctionFactor;
169 correctionTree->SetBranchAddress(m_bgFractionBinNumBranchName.c_str(), &bgFractionBinNum);
170 correctionTree->SetBranchAddress(m_regNumBranchName.c_str(), ®Num);
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);
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;
184 for (
long iEntry = 0; iEntry < correctionTree->GetEntries(); ++iEntry) {
185 correctionTree->GetEntry(iEntry);
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);
209 int numOfReg1ThetaBins;
210 int numOfReg2ThetaBins;
211 int numOfReg3ThetaBins;
218 std::vector<float> avgRecEns(m_numAvgRecEnEntries);
219 helperTree->SetBranchAddress(m_avgRecEnBranchName.c_str(), avgRecEns.data());
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);
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;
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]);
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);
290 m_numOfReg1ThetaBins,
291 m_numOfReg2ThetaBins,
292 m_numOfReg3ThetaBins,
309 void ECLDatabaseImporter::importShowerShapesSecondMomentCorrections()
311 if (m_inputFileNames.size() > 1)
312 B2FATAL(
"Sorry, you must only import one file at a time for now!");
316 TFile* inputFile =
new TFile(m_inputFileNames[0].data(),
"READ");
318 if (!inputFile || inputFile->IsZombie())
319 B2FATAL(
"Could not open file " << m_inputFileNames[0]);
322 TGraph* theta_N1_graph = getRootObjectFromFile<TGraph*>(inputFile,
"SecondMomentCorrections_theta_N1");
323 dbArray.
appendNew(ECLShower::c_nPhotons, ECLShowerShapeModule::c_thetaType , *theta_N1_graph);
326 TGraph* phi_N1_graph = getRootObjectFromFile<TGraph*>(inputFile,
"SecondMomentCorrections_phi_N1");
327 dbArray.
appendNew(ECLShower::c_nPhotons, ECLShowerShapeModule::c_phiType , *phi_N1_graph);
330 TGraph* theta_N2_graph = getRootObjectFromFile<TGraph*>(inputFile,
"SecondMomentCorrections_theta_N2");
331 dbArray.
appendNew(ECLShower::c_neutralHadron, ECLShowerShapeModule::c_thetaType , *theta_N2_graph);
334 TGraph* phi_N2_graph = getRootObjectFromFile<TGraph*>(inputFile,
"SecondMomentCorrections_phi_N2");
335 dbArray.
appendNew(ECLShower::c_neutralHadron, ECLShowerShapeModule::c_phiType , *phi_N2_graph);
351 void ECLDatabaseImporter::importShowerEnergyCorrectionTemporary()
353 if (m_inputFileNames.size() > 1)
354 B2FATAL(
"Sorry, you must only import one file at a time for now!");
357 boost::filesystem::path path(m_inputFileNames[0]);
358 if (path.extension() !=
".root")
359 B2FATAL(
"Expecting a .root file. Aborting");
361 TFile* inputFile =
new TFile(m_inputFileNames[0].data(),
"READ");
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");
369 double bkgFactor = bg_histo->GetBinContent(1);
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();
376 double energyMin = theta_en_graph->GetYmin();
377 double energyMax = theta_en_graph->GetYmax();
380 B2DEBUG(28,
"Leakage DBobjects angle boundaries: thetaMin=" << thetaMin <<
" thetaMax=" << thetaMax <<
" phiMin= " << phiMin <<
381 " phiMax= " << phiMax <<
" enmin=" << energyMin <<
382 " enmax=" << energyMax);
391 if (std::abs(bkgFactor - 1.0) < 1e-9) {
394 dbPtr_theta_geo.
construct(*theta_geo_graph, thetaMin, thetaMax, energyMin, energyMax);
396 dbPtr_phi_geo.
construct(*phi_geo_graph, phiMin, phiMax, energyMin, energyMax);
398 dbPtr_theta_en.
construct(*theta_en_graph, thetaMin, thetaMax, energyMin, energyMax);
400 dbPtr_phi_en.
construct(*phi_en_graph, phiMin, phiMax, energyMin, energyMax);
403 dbPtr_theta_geo.
import(iov);
404 dbPtr_phi_geo.
import(iov);
405 dbPtr_theta_en.
import(iov);
408 if (std::abs(bkgFactor - 1.0) < 1e-9) {
411 dbPtr_theta_geo.
construct(*theta_geo_graph, thetaMin, thetaMax, energyMin, energyMax);
413 dbPtr_phi_geo.
construct(*phi_geo_graph, phiMin, phiMax, energyMin, energyMax);
415 dbPtr_theta_en.
construct(*theta_en_graph, thetaMin, thetaMax, energyMin, energyMax);
417 dbPtr_phi_en.
construct(*phi_en_graph, phiMin, phiMax, energyMin, energyMax);
420 dbPtr_theta_geo.
import(iov);
421 dbPtr_phi_geo.
import(iov);
422 dbPtr_theta_en.
import(iov);
429 void ECLDatabaseImporter::importTrackClusterMatchingThresholds()
431 if (m_inputFileNames.size() > 1)
432 B2FATAL(
"Sorry, you must only import one file at a time for now!");
435 boost::filesystem::path path(m_inputFileNames[0]);
436 if (path.extension() !=
".txt")
437 B2FATAL(
"Expecting a .txt file. Aborting");
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;
447 ifstream infile(m_inputFileNames[0]);
449 while (getline(infile, line)) {
450 istringstream iss(line);
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);
466 dbPtr.
construct(m_matchingThresholdPairsFWD, m_matchingThresholdPairsBWD, m_matchingThresholdPairsBRL);
474 void ECLDatabaseImporter::importTrackClusterMatchingParameterizations()
476 if (m_inputFileNames.size() > 1)
477 B2FATAL(
"Sorry, you must only import one file at a time for now!");
480 TFile* inputFile =
new TFile(m_inputFileNames[0].data(),
"READ");
482 if (!inputFile || inputFile->IsZombie())
483 B2FATAL(
"Could not open file " << m_inputFileNames[0]);
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"};
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))));
500 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.