Belle II Software  release-05-02-19
ECLDatabaseImporter Class Reference

ECL database importer. More...

#include <ECLDatabaseImporter.h>

Collaboration diagram for ECLDatabaseImporter:

Public Member Functions

 ECLDatabaseImporter (std::vector< std::string > inputFileNames, const std::string &m_name)
 Constructor.
 
virtual ~ECLDatabaseImporter ()
 Destructor.
 
void importDigitEnergyCalibration ()
 Import ECL energy calibration constants to the database.
 
void importDigitTimeCalibration ()
 Import ECL time calibration constants to the database.
 
void importShowerShapesSecondMomentCorrections ()
 Import ECL shower shape corrections to second moment to the database.
 
void importShowerCorrectorLeakageCorrections ()
 Import ECL leakage corrections to showers.
 
void importShowerEnergyCorrectionTemporary ()
 Import ECL corrections to showers energies. More...
 
void importTrackClusterMatchingThresholds ()
 Import threshold values for track ECL cluster matching to the database.
 
void importTrackClusterMatchingParameterizations ()
 Import parameterizations for the RMS between tracks and ECL clusters to the database.
 

Private Member Functions

template<class rootClass >
rootClass getRootObjectFromFile (TFile *file, const std::string &rootObjName) const
 Extract a root object rootObjName from a root file file. More...
 

Private Attributes

std::vector< std::string > m_inputFileNames
 Input file name.
 
std::string m_name
 Database object (output) file name.
 
std::string m_bgFractionBinNumBranchName = "bgFractionBinNum"
 Branch names for shower corrector leakage correction root file. More...
 
std::string m_regNumBranchName = "regNum"
 regNum branch name
 
std::string m_phiBinNumBranchName = "phiBinNum"
 phiBinNum branch name
 
std::string m_thetaBinNumBranchName = "thetaBinNum"
 thetaBinNum branch name
 
std::string m_energyBinNumBranchName = "energyBinNum"
 energyBinNum branch name
 
std::string m_correctionFactorBranchName = "correctionFactor"
 correctionFactor branch name
 
const int m_numAvgRecEnEntries = 15
 Number of entries in avgRecEn array. More...
 
std::string m_avgRecEnBranchName = "avgRecEn"
 avgRecEn branch name
 
std::string m_lReg1ThetaBranchName = "lReg1Theta"
 lReg1Theta branch name
 
std::string m_hReg1ThetaBranchName = "hReg1Theta"
 hReg1Theta branch name
 
std::string m_lReg2ThetaBranchName = "lReg2Theta"
 lReg2Theta branch name
 
std::string m_hReg2ThetaBranchName = "hReg2Theta"
 hReg2Theta branch name
 
std::string m_lReg3ThetaBranchName = "lReg3Theta"
 lReg3Theta branch name
 
std::string m_hReg3ThetaBranchName = "hReg3Theta"
 hReg3Theta branch name
 
std::string m_numOfBfBinsBranchName = "numOfBfBins"
 numOfBfBins branch name
 
std::string m_numOfEnergyBinsBranchName = "numOfEnergyBins"
 numOfEnergyBins branch name
 
std::string m_numOfPhiBinsBranchName = "numOfPhiBins"
 numOfPhiBins branch name
 
std::string m_numOfReg1ThetaBinsBranchName = "numOfReg1ThetaBins"
 numOfReg1ThetaBins branch name
 
std::string m_numOfReg2ThetaBinsBranchName = "numOfReg2ThetaBins"
 numOfReg2ThetaBins branch name

 
std::string m_numOfReg3ThetaBinsBranchName = "numOfReg3ThetaBins"
 numOfReg3ThetaBins branch name

 
std::string m_phiPeriodicityBranchName = "phiPeriodicity"
 phiPeriodicity branch name

 

Detailed Description

ECL database importer.

This module writes data from e.g. a ROOT histogram to database.

Definition at line 38 of file ECLDatabaseImporter.h.

Member Function Documentation

◆ getRootObjectFromFile()

rootClass getRootObjectFromFile ( TFile *  file,
const std::string &  rootObjName 
) const
inlineprivate

Extract a root object rootObjName from a root file file.

The file is assumed to be valid (pointer valid and not zombie). If rootObjName doesn't exist in file, do B2FATAL.

Definition at line 100 of file ECLDatabaseImporter.h.

◆ importShowerEnergyCorrectionTemporary()

void importShowerEnergyCorrectionTemporary ( )

Import ECL corrections to showers energies.

Temperary - there will be additional improvements done to these corrections in the future. Input file should be .txt file and have the format: generated-energy bkg-scale-factor theta-min theta-max corr-factor. The numbers should be seperated by spaces. For each line, the theta value used is the average of theta-min and theta-max

Definition at line 354 of file ECLDatabaseImporter.cc.

355 {
356  if (m_inputFileNames.size() > 1)
357  B2FATAL("Sorry, you must only import one file at a time for now!");
358 
359  //Expect a root file
360  boost::filesystem::path path(m_inputFileNames[0]);
361  if (path.extension() != ".root")
362  B2FATAL("Expecting a .root file. Aborting");
363 
364  TFile* inputFile = new TFile(m_inputFileNames[0].data(), "READ");
365 
366  TGraph2D* theta_geo_graph = getRootObjectFromFile<TGraph2D*>(inputFile, "LeakageCorrections_theta_geometry");
367  TGraph2D* phi_geo_graph = getRootObjectFromFile<TGraph2D*>(inputFile, "LeakageCorrections_phi_geometry");
368  TGraph2D* theta_en_graph = getRootObjectFromFile<TGraph2D*>(inputFile, "LeakageCorrections_theta_energy");
369  TGraph2D* phi_en_graph = getRootObjectFromFile<TGraph2D*>(inputFile, "LeakageCorrections_phi_energy");
370  TH1F* bg_histo = getRootObjectFromFile<TH1F*>(inputFile, "LeakageCorrections_background_fraction");
371 
372  double bkgFactor = bg_histo->GetBinContent(1);
373 
374  double thetaMin = theta_en_graph->GetXmin();
375  double thetaMax = theta_en_graph->GetXmax();
376  double phiMin = phi_en_graph->GetXmin();
377  double phiMax = phi_en_graph->GetXmax();
378 
379  double energyMin = theta_en_graph->GetYmin();
380  double energyMax = theta_en_graph->GetYmax();
381 
382 
383  B2DEBUG(28, "Leakage DBobjects angle boundaries: thetaMin=" << thetaMin << " thetaMax=" << thetaMax << " phiMin= " << phiMin <<
384  " phiMax= " << phiMax << " enmin=" << energyMin <<
385  " enmax=" << energyMax);
386 
387  // Import to DB
388  int startExp = 0;
389  int startRun = 0;
390  int endExp = -1;
391  int endRun = -1;
392  IntervalOfValidity iov(startExp, startRun, endExp, endRun);
393 
394  if (std::abs(bkgFactor - 1.0) < 1e-9) { //bkgFactor == 1 -> phase 2 backgrounds
395 
396  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_theta_geo("ECLLeakageCorrection_thetaGeometry_phase2");
397  dbPtr_theta_geo.construct(*theta_geo_graph, thetaMin, thetaMax, energyMin, energyMax);
398  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_phi_geo("ECLLeakageCorrection_phiGeometry_phase2");
399  dbPtr_phi_geo.construct(*phi_geo_graph, phiMin, phiMax, energyMin, energyMax);
400  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_theta_en("ECLLeakageCorrection_thetaEnergy_phase2");
401  dbPtr_theta_en.construct(*theta_en_graph, thetaMin, thetaMax, energyMin, energyMax);
402  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_phi_en("ECLLeakageCorrection_phiEnergy_phase2");
403  dbPtr_phi_en.construct(*phi_en_graph, phiMin, phiMax, energyMin, energyMax);
404 
405  //Import into local db
406  dbPtr_theta_geo.import(iov);
407  dbPtr_phi_geo.import(iov);
408  dbPtr_theta_en.import(iov);
409  dbPtr_phi_en.import(iov);
410  }
411  /*else (because currently phase_2 and phase_3 are same payload*/ if (std::abs(bkgFactor - 1.0) < 1e-9) {
412 
413  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_theta_geo("ECLLeakageCorrection_thetaGeometry_phase3");
414  dbPtr_theta_geo.construct(*theta_geo_graph, thetaMin, thetaMax, energyMin, energyMax);
415  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_phi_geo("ECLLeakageCorrection_phiGeometry_phase3");
416  dbPtr_phi_geo.construct(*phi_geo_graph, phiMin, phiMax, energyMin, energyMax);
417  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_theta_en("ECLLeakageCorrection_thetaEnergy_phase3");
418  dbPtr_theta_en.construct(*theta_en_graph, thetaMin, thetaMax, energyMin, energyMax);
419  DBImportObjPtr<ECLShowerEnergyCorrectionTemporary> dbPtr_phi_en("ECLLeakageCorrection_phiEnergy_phase3");
420  dbPtr_phi_en.construct(*phi_en_graph, phiMin, phiMax, energyMin, energyMax);
421 
422  //Import into local db
423  dbPtr_theta_geo.import(iov);
424  dbPtr_phi_geo.import(iov);
425  dbPtr_theta_en.import(iov);
426  dbPtr_phi_en.import(iov);
427  }
428 
429 
430 }

Member Data Documentation

◆ m_bgFractionBinNumBranchName

std::string m_bgFractionBinNumBranchName = "bgFractionBinNum"
private

Branch names for shower corrector leakage correction root file.

bgFractionBinNum branch name

Definition at line 116 of file ECLDatabaseImporter.h.

◆ m_numAvgRecEnEntries

const int m_numAvgRecEnEntries = 15
private

Number of entries in avgRecEn array.

If this is wrong bad things will happen

Definition at line 124 of file ECLDatabaseImporter.h.


The documentation for this class was generated from the following files:
Belle2::IntervalOfValidity
A class that describes the interval of experiments/runs for which an object in the database is valid.
Definition: IntervalOfValidity.h:35
Belle2::ECLDatabaseImporter::m_inputFileNames
std::vector< std::string > m_inputFileNames
Input file name.
Definition: ECLDatabaseImporter.h:94
Belle2::DBImportObjPtr
Class for importing a single object to the database.
Definition: DBImportObjPtr.h:33