Belle II Software  release-05-01-25
ECLLocalRunCalibUnit.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2019 - Belle II Collaboration *
4  * *
5  * ECLLocalRunCalibUnit *
6  * *
7  * This class controls feature (mean value, standard deviation and *
8  * number of accepted events) accumulators for each cell id. *
9  * *
10  * Author: The Belle II Collaboration *
11  * Contributors: Sergei Gribanov (S.S.Gribanov@inp.nsk.su) (BINP) *
12  * *
13  * This software is provided "as is" without any warranty. *
14  **************************************************************************/
15 
16 // ECL
17 #include <ecl/modules/eclLocalRunCalibration/ECLLocalRunCalibUnit.h>
18 #include <ecl/dbobjects/ECLCrystalLocalRunCalib.h>
19 #include <ecl/dbobjects/ECLLocalRunCalibRef.h>
20 
21 using namespace Belle2;
22 // Constructor.
24  const int& ncellids,
25  const float& min_value,
26  const float& max_value,
27  const int* const ndevs):
28  m_isNegAmpl(false), m_unitData(ncellids,
30  min_value,
31  max_value,
32  ndevs))
33 {
34 }
35 // Destructor.
37 {
38 }
39 // This function will be called, if the
40 // negative amplitude values are observed
41 // in the current run.
43 {
44  m_isNegAmpl = true;
45 }
46 // Checking presence of negative
47 // amplitudes.
49 {
50  return m_isNegAmpl;
51 }
52 // Add value to
53 // accamulate mean value,
54 // standard deviation and
55 // number of accepted events.
56 void ECLLocalRunCalibUnit::add(const int& index, const float& value)
57 {
58  m_unitData[index].add(value);
59 }
60 // Calculate accumulated values.
62 {
63  for (auto& cellAcc : m_unitData) {
64  cellAcc.calc();
65  }
66 }
67 // Getter of accumulated values.
68 template <typename T>
70  std::vector<T>* vec,
71  T(ECLLocalRunCalibAcc::*getter)() const)
72 {
73  for (const auto& obj : m_unitData) {
74  vec->push_back((obj.*getter)());
75  }
76 }
77 // Make current run as reference.
78 void ECLLocalRunCalibUnit::markAsRefference(const bool& isLocal,
79  const std::string& dbName,
80  const int& run,
81  const IntervalOfValidity& iov)
82 {
83  int exp = iov.getExperimentLow();
84  ECLDBTool refPayload(isLocal, dbName.c_str(),
85  "ECLCalibRef");
86  refPayload.connect();
87  ECLLocalRunCalibRef refobj(exp, run);
88  refPayload.write(&refobj, iov);
89 }
90 // Change previous validity interval.
92  const IntervalOfValidity& iov)
93 {
94  int exp = iov.getExperimentLow();
95  int run_high = iov.getRunLow() - 1;
96  EventMetaData ev(1, run_high, exp);
97  IntervalOfValidity* prevIoV;
98  payload.read(&prevIoV, ev);
99  int run_low = prevIoV->getRunLow();
100  IntervalOfValidity ciov(exp, run_low,
101  exp, run_high);
102  payload.changeIoV(ev, ciov);
103  delete prevIoV;
104 }
105 // Save results.
107  bool isLocal,
108  const std::string& dbName,
109  const std::string& payloadName,
110  const IntervalOfValidity& iov,
111  const int& run,
112  const bool& changePrev,
113  const bool& addref)
114 {
115  int exp = iov.getExperimentLow();
117  if (m_unitData.size() > 0) {
118  int nevents = 0;
119  for (const auto& cellAcc : m_unitData) {
120  nevents = cellAcc.getNOfEvents();
121  if (nevents != 0) {
122  break;
123  }
124  }
125  std::vector<int> counts;
126  counts.reserve(m_unitData.size());
127  callAccGetter<int>(&counts,
129  std::vector<float> means;
130  means.reserve(m_unitData.size());
131  callAccGetter<float>(&means,
133  std::vector<float> stddevs;
134  stddevs.reserve(m_unitData.size());
135  callAccGetter<float>(&stddevs,
137  data.setNumberOfEvents(nevents);
138  data.setNumbersOfAcceptedEvents(counts);
139  data.setCalibVector(means, stddevs);
140  data.setExpRun(exp, run);
141  ECLDBTool payload(isLocal, dbName.c_str(),
142  payloadName.c_str());
143  payload.connect();
144  if (changePrev) {
145  changePreviousIoV(payload, iov);
146  }
147  payload.write(&data, iov);
148  if (addref) {
149  markAsRefference(isLocal, dbName,
150  run, iov);
151  }
152  }
153 }
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::ECLLocalRunCalibUnit::isNegAmpl
bool isNegAmpl() const
Check presence of negative amplitudes in the current run.
Definition: ECLLocalRunCalibUnit.cc:48
Belle2::ECLLocalRunCalibAcc::getCount
int getCount() const
Get number of accepted events.
Definition: ECLLocalRunCalibAcc.cc:123
Belle2::ECLLocalRunCalibAcc::getStdDev
float getStdDev() const
Get standard deviation.
Definition: ECLLocalRunCalibAcc.cc:133
Belle2::ECLLocalRunCalibUnit::m_unitData
std::vector< ECLLocalRunCalibAcc > m_unitData
Mean value and standard deviation accumulators for each cell id.
Definition: ECLLocalRunCalibUnit.h:138
Belle2::ECLLocalRunCalibUnit::changePreviousIoV
void changePreviousIoV(const ECLDBTool &payload, const IntervalOfValidity &iov)
Change previous validity interval.
Definition: ECLLocalRunCalibUnit.cc:91
Belle2::ECLDBTool::connect
void connect() const
Connect to a database.
Definition: ECLDBTool.cc:34
Belle2::ECLDBTool
The ECLDBTool class is designed to read / write object from / to database.
Definition: ECLDBTool.h:41
Belle2::ECLLocalRunCalibUnit::markAsRefference
void markAsRefference(const bool &isLocal, const std::string &dbName, const int &run, const IntervalOfValidity &iov)
Mark current run as reference.
Definition: ECLLocalRunCalibUnit.cc:78
Belle2::ECLDBTool::write
void write(TObject *const obj, const IntervalOfValidity &iov) const
Write object and validity interval to a database.
Definition: ECLDBTool.cc:44
Belle2::ECLLocalRunCalibAcc
ECLLocalRunCalibAcc is the class designed to accumulate mean values, standard deviation and number of...
Definition: ECLLocalRunCalibAcc.h:47
Belle2::ECLLocalRunCalibUnit::writeToDB
void writeToDB(bool isLocal, const std::string &dbName, const std::string &payloadName, const IntervalOfValidity &iov, const int &run, const bool &changePrev, const bool &addref)
Write calibration results into a database.
Definition: ECLLocalRunCalibUnit.cc:106
Belle2::ECLLocalRunCalibAcc::getMean
float getMean() const
Get mean value.
Definition: ECLLocalRunCalibAcc.cc:128
Belle2::ECLLocalRunCalibUnit::add
void add(const int &cellid, const float &value)
Add value to accumulate mean value, standard deviation and number of accepted events.
Definition: ECLLocalRunCalibUnit.cc:56
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::ECLLocalRunCalibUnit::ECLLocalRunCalibUnit
ECLLocalRunCalibUnit(const int &ncellids, const float &min_value, const float &max_value, const int *const ndevs)
Constructor.
Definition: ECLLocalRunCalibUnit.cc:23
Belle2::IntervalOfValidity::getRunLow
int getRunLow() const
Getter for lowest run number.
Definition: IntervalOfValidity.h:135
Belle2::ECLLocalRunCalibUnit::callAccGetter
void callAccGetter(std::vector< T > *vec, T(ECLLocalRunCalibAcc::*getter)() const)
Getter of accumulated values.
Definition: ECLLocalRunCalibUnit.cc:69
Belle2::ECLLocalRunCalibUnit::calc
void calc()
Calculate accumulated values.
Definition: ECLLocalRunCalibUnit.cc:61
Belle2::EventMetaData
Store event, run, and experiment numbers.
Definition: EventMetaData.h:43
Belle2::ECLLocalRunCalibUnit::enableNegAmpl
void enableNegAmpl()
This function will be called only in the case, if negative amplitudes are observed in the current run...
Definition: ECLLocalRunCalibUnit.cc:42
Belle2::ECLCrystalLocalRunCalib
ECLCrystalLocalRunCalib is designed to store results of the ECL local run calibration to database.
Definition: ECLCrystalLocalRunCalib.h:44
Belle2::ECLLocalRunCalibRef
ECLLocalRunCalibRef is designed to store reference marks to database for ECL local run calibration.
Definition: ECLLocalRunCalibRef.h:42
Belle2::ECLLocalRunCalibUnit::~ECLLocalRunCalibUnit
~ECLLocalRunCalibUnit()
Destructor.
Definition: ECLLocalRunCalibUnit.cc:36
Belle2::ECLLocalRunCalibUnit::m_isNegAmpl
bool m_isNegAmpl
m_isNegAmpl is true if there are negative amplitudes in the current run and false otherwise.
Definition: ECLLocalRunCalibUnit.h:115