Skip to content

Commit 47dc1fa

Browse files
jokonigdavidrohr
jokonig
authored and
davidrohr
committed
[EMCAL-565, EMCAL-566] Enable Multithreading for EMCal online calib
- To increase the speed of the online calibration, multithreading is introduced to fill the histograms - Each thread is assigned its own histogram (stored in a vector with the size = number of threads) - The data for each event is then splitted into nThread chunks and filled into the corresponding histograms - The number of threads can be set via the EMCal calib params - Merging of the individual histograms happens if the histogram is retrieved from the Data class - Local tests show a significant increase in speed of more threads are being used. The overhead is negligible
1 parent c1219d2 commit 47dc1fa

7 files changed

Lines changed: 214 additions & 67 deletions

File tree

‎Common/Utils/include/CommonUtils/BoostHistogramUtils.h‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -618,7 +618,7 @@ auto ProjectBoostHistoXFast(const boost::histogram::histogram<axes...>& hist2d,
618618
/// \return result
619619
/// 1d boost histogram from projection of the input 2d boost histogram
620620
template <typename... axes>
621-
auto ReduceBoostHistoFastSlice(boost::histogram::histogram<axes...>& hist2d, int binXLow, int binXHigh, int binYLow, int binYHigh, bool includeOverflowUnderflow)
621+
auto ReduceBoostHistoFastSlice(const boost::histogram::histogram<axes...>& hist2d, int binXLow, int binXHigh, int binYLow, int binYHigh, bool includeOverflowUnderflow)
622622
{
623623
int nXbins = binXHigh - binXLow + 1;
624624
int nYbins = binYHigh - binYLow + 1;

‎Detectors/EMCAL/calibration/include/EMCALCalibration/EMCALCalibExtractor.h‎

Lines changed: 5 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -90,11 +90,12 @@ class EMCALCalibExtractor
9090
/// \param hist histogram cell energy vs. cell ID. Main histogram for the bad channel calibration
9191
/// \param histTime histogram cell time vs. cell ID. If default argument is taken, no calibration based on the timing signal will be performed
9292
template <typename... axes>
93-
o2::emcal::BadChannelMap calibrateBadChannels(boost::histogram::histogram<axes...>& hist, const boost::histogram::histogram<axes...>& histTime = boost::histogram::make_histogram(boost::histogram::axis::variable<>{0., 1.}, boost::histogram::axis::variable<>{0., 1.}))
93+
o2::emcal::BadChannelMap calibrateBadChannels(const boost::histogram::histogram<axes...>& hist, const boost::histogram::histogram<axes...>& histTime = boost::histogram::make_histogram(boost::histogram::axis::variable<>{0., 1.}, boost::histogram::axis::variable<>{0., 1.}))
9494
{
9595
double time1 = std::chrono::duration_cast<std::chrono::nanoseconds>(std::chrono::high_resolution_clock::now().time_since_epoch()).count();
9696
std::map<int, std::pair<double, double>> slices = {{0, {0.1, 0.3}}, {1, {0.3, 0.5}}, {2, {0.5, 1.0}}, {3, {1.0, 4.0}}, {4, {4.0, 39.0}}};
9797

98+
auto histScaled = hist;
9899
if (mBCMScaleFactors) {
99100
LOG(info) << "Rescaling BCM histo";
100101
// rescale the histogram
@@ -103,13 +104,13 @@ class EMCALCalibExtractor
103104
double lowerE = hist.axis(0).bin(ebin).lower();
104105
double upperE = hist.axis(0).bin(ebin).upper();
105106
double midE = (lowerE + upperE) / 2.;
106-
hist.at(ebin, icell) = hist.at(ebin, icell) / mBCMScaleFactors->getScaleVal(icell, midE);
107+
histScaled.at(ebin, icell) = hist.at(ebin, icell) / mBCMScaleFactors->getScaleVal(icell, midE);
107108
}
108109
}
109110
}
110111

111112
// get all ofthe calibration information that we need in a struct
112-
BadChannelCalibInfo calibrationInformation = buildHitAndEnergyMean(slices, hist);
113+
BadChannelCalibInfo calibrationInformation = buildHitAndEnergyMean(slices, histScaled);
113114

114115
// only initialize this if the histo is not the default one
115116
const bool doIncludeTime = (histTime.axis(0).size() > 1 && EMCALCalibParams::Instance().useTimeInfoForCalib_bc) ? true : false;
@@ -338,7 +339,7 @@ class EMCALCalibExtractor
338339
/// \param maxTime -- max. time considered for fit
339340
/// \param restrictFitRangeToMax -- restrict the fit range to the maximum entry in the histogram in the range +-restrictFitRangeToMax (default: 25ns)
340341
template <typename... axes>
341-
o2::emcal::TimeCalibrationParams calibrateTime(boost::histogram::histogram<axes...>& hist, double minTime = 0, double maxTime = 1000, double restrictFitRangeToMax = 25)
342+
o2::emcal::TimeCalibrationParams calibrateTime(const boost::histogram::histogram<axes...>& hist, double minTime = 0, double maxTime = 1000, double restrictFitRangeToMax = 25)
342343
{
343344

344345
auto histReduced = boost::histogram::algorithm::reduce(hist, boost::histogram::algorithm::shrink(minTime, maxTime), boost::histogram::algorithm::shrink(0, mNcells));

‎Detectors/EMCAL/calibration/include/EMCALCalibration/EMCALCalibParams.h‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -49,6 +49,7 @@ struct EMCALCalibParams : public o2::conf::ConfigurableParamHelper<EMCALCalibPar
4949
float minNHitsForNHitCut = 1; ///< mean number of hits per cell that is needed to cut on the mean number of hits. Needed for high energy intervals as outliers can distort the distribution
5050
float minCellEnergy_bc = 0.1; ///< minimum cell energy considered for filling the histograms for bad channel calib. Should speedup the filling of the histogram to suppress noise
5151
float fractionEvents_bc = 1.; ///< fraction of events used in bad channel calibration
52+
size_t nThreads_bc = 4; ///< number of threads used for the bad channel calinration for filling the histograms
5253

5354
// parameters for time calibration
5455
unsigned int minNEvents_tc = 1e7; ///< minimum number of events to trigger the calibration
@@ -66,6 +67,7 @@ struct EMCALCalibParams : public o2::conf::ConfigurableParamHelper<EMCALCalibPar
6667
bool UpdateAtEndOfRunOnly_tc = false; ///< switch to enable trigger of calibration only at end of run
6768
float maxAllowedDeviationFromMax = 10; ///< maximum deviation allowed between the estimated maximum of the fit and the true maximum from the distribution. If deviation is larger then value, the fit likely failed. In this case, the true value is taken
6869
float fractionEvents_tc = 1.; ///< fraction of events used in time calibration
70+
size_t nThreads_tc = 2; ///< number of threads used for the time calinration for filling the histograms
6971

7072
// common parameters
7173
std::string calibType = "time"; ///< type of calibration to run

‎Detectors/EMCAL/calibration/include/EMCALCalibration/EMCALChannelData.h‎

Lines changed: 47 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -39,6 +39,8 @@
3939
#include <boost/histogram/ostream.hpp>
4040
// #include <boost/format.hpp>
4141

42+
#include <thread>
43+
4244
// #include <array>
4345

4446
namespace o2
@@ -59,15 +61,21 @@ class EMCALChannelData
5961
o2::emcal::Geometry* mGeometry = o2::emcal::Geometry::GetInstanceFromRunNumber(300000);
6062
int NCELLS = mGeometry->GetNCells();
6163

62-
EMCALChannelData() : mNBins(EMCALCalibParams::Instance().nBinsEnergyAxis_bc), mRange(EMCALCalibParams::Instance().maxValueEnergyAxis_bc), mNBinsTime(EMCALCalibParams::Instance().nBinsTimeAxis_bc), mRangeTimeLow(EMCALCalibParams::Instance().rangeTimeAxisLow_bc), mRangeTimeHigh(EMCALCalibParams::Instance().rangeTimeAxisHigh_bc)
64+
EMCALChannelData() : mNBins(EMCALCalibParams::Instance().nBinsEnergyAxis_bc), mRange(EMCALCalibParams::Instance().maxValueEnergyAxis_bc), mNBinsTime(EMCALCalibParams::Instance().nBinsTimeAxis_bc), mRangeTimeLow(EMCALCalibParams::Instance().rangeTimeAxisLow_bc), mRangeTimeHigh(EMCALCalibParams::Instance().rangeTimeAxisHigh_bc), mNThreads(EMCALCalibParams::Instance().nThreads_bc)
6365
{
6466

6567
// NCELLS includes DCal, treat as one calibration
6668
o2::emcal::Geometry* mGeometry = o2::emcal::Geometry::GetInstanceFromRunNumber(300000);
6769
int NCELLS = mGeometry->GetNCells();
6870

69-
mHisto = boost::histogram::make_histogram(boost::histogram::axis::regular<>(mNBins, 0., mRange), boost::histogram::axis::regular<>(NCELLS, -0.5, NCELLS - 0.5));
70-
mHistoTime = boost::histogram::make_histogram(boost::histogram::axis::regular<>(mNBinsTime, mRangeTimeHigh, mRangeTimeLow), boost::histogram::axis::regular<>(NCELLS, -0.5, NCELLS - 0.5));
71+
mVecNEntriesInHisto.resize(mNThreads);
72+
mHisto.resize(mNThreads);
73+
mHistoTime.resize(mNThreads);
74+
for (size_t i = 0; i < mNThreads; ++i) {
75+
mHisto[i] = boost::histogram::make_histogram(boost::histogram::axis::regular<>(mNBins, 0., mRange), boost::histogram::axis::regular<>(NCELLS, -0.5, NCELLS - 0.5));
76+
mHistoTime[i] = boost::histogram::make_histogram(boost::histogram::axis::regular<>(mNBinsTime, mRangeTimeHigh, mRangeTimeLow), boost::histogram::axis::regular<>(NCELLS, -0.5, NCELLS - 0.5));
77+
mVecNEntriesInHisto[i] = 0;
78+
}
7179
}
7280

7381
~EMCALChannelData() = default;
@@ -83,23 +91,40 @@ class EMCALChannelData
8391
/// \brief Fill the container with the cell ID and amplitude.
8492
void fill(const gsl::span<const o2::emcal::Cell> data);
8593
/// \brief Merge the data of two slots.
86-
void merge(const EMCALChannelData* prev);
94+
void merge(EMCALChannelData* prev);
8795
// int findBin(float v) const;
8896
/// \brief Check if enough stataistics was accumulated to perform calibration
8997
bool hasEnoughData() const;
9098
/// \brief Get current calibration histogram
91-
boostHisto& getHisto() { return mHisto; }
92-
const boostHisto& getHisto() const { return mHisto; }
99+
const boostHisto& getHisto()
100+
{
101+
// set the summed histogram to one of the existing histograms
102+
mHistoSummed = mHisto[0];
103+
// reset the histogram
104+
mHistoSummed.reset();
105+
// Sum up all entries
106+
for (const auto& h : mHisto) {
107+
mHistoSummed += h;
108+
}
109+
return mHistoSummed;
110+
}
93111

94112
/// \brief Set new calibration histogram
95-
void setHisto(boostHisto hist) { mHisto = hist; }
113+
void setHisto(boostHisto hist, int nthr = 0) { mHisto[nthr] = hist; }
96114

97-
/// \brief Get current calibration histogram with timing information
98-
boostHisto& getHistoTime() { return mHistoTime; }
99-
const boostHisto& getHistoTime() const { return mHistoTime; }
115+
/// \brief Get current calibration histogram with time information
116+
const boostHisto& getHistoTime()
117+
{
118+
mHistoTimeSummed = mHistoTime[0];
119+
mHistoTimeSummed.reset();
120+
for (const auto& h : mHistoTime) {
121+
mHistoTimeSummed += h;
122+
}
123+
return mHistoTimeSummed;
124+
}
100125

101126
/// \brief Set new calibration histogram with timing info
102-
void setHistoTime(boostHisto hist) { mHistoTime = hist; }
127+
void setHistoTime(boostHisto hist, int nthr = 0) { mHistoTime[nthr] = hist; }
103128

104129
/// \brief Peform the calibration and flag the bad channel map
105130
/// Average energy per hit histogram is fitted with a gaussian
@@ -118,30 +143,39 @@ class EMCALChannelData
118143

119144
long unsigned int getNEntriesInHisto() const { return mNEntriesInHisto; }
120145
void setNEntriesInHisto(long unsigned int n) { mNEntriesInHisto = n; }
146+
void addEntriesInHisto(long unsigned int n) { mNEntriesInHisto += n; }
121147

122148
void setGainCalibFactors(o2::emcal::GainCalibrationFactors* calibFactors)
123149
{
124150
mGainCalibFactors = calibFactors;
151+
for (unsigned int i = 0; i < mArrGainCalibFactors.size(); ++i) {
152+
mArrGainCalibFactors[i] = mGainCalibFactors->getGainCalibFactors(i);
153+
}
125154
mApplyGainCalib = true;
126155
}
127156

128157
private:
129158
float mRange = 10; ///< Maximum energy range of boost histogram (will be overwritten by values in the EMCALCalibParams)
130159
int mNBins = 1000; ///< Number of bins in the boost histogram (will be overwritten by values in the EMCALCalibParams)
131-
boostHisto mHisto; ///< 2d boost histogram with cellID vs cell energy
160+
size_t mNThreads = 1; ///< Number of threads used for filling the boost histograms
161+
std::vector<boostHisto> mHisto; ///< vector of 2d boost histogram with cellID vs cell energy
162+
boostHisto mHistoSummed; ///< summed 2d boost histogram (sum of mHisto)
132163
int mNBinsTime = 1000; ///< Number of time bins in boost histogram (cell time vs. cell ID)
133164
float mRangeTimeLow = -500; ///< lower bound of time axis of mHistoTime
134165
float mRangeTimeHigh = 500; ///< upper bound of time axis of mHistoTime
135-
boostHisto mHistoTime; ///< 2d boost histogram with cellID vs cell time
166+
std::vector<boostHisto> mHistoTime; ///< vector of 2d boost histogram with cellID vs cell time
167+
boostHisto mHistoTimeSummed; ///< Summed 2d boost histogram with cellID vs cell time
136168
int mEvents = 0; ///< event counter
137169
long unsigned int mNEntriesInHisto = 0; ///< Number of entries in the histogram
170+
std::vector<long unsigned int> mVecNEntriesInHisto; ///< Number of entries in the histogram for each thread per event
138171
boostHisto mEsumHisto; ///< contains the average energy per hit for each cell
139172
boostHisto mEsumHistoScaled; ///< contains the average energy (scaled) per hit for each cell
140173
boostHisto mCellAmplitude; ///< is the input for the calibration, hist of cell E vs. ID
141174
bool mTest = false; ///< flag to be used when running in test mode: it simplify the processing
142175
BadChannelMap mOutputBCM; ///< output bad channel map for the calibration
143176
bool mApplyGainCalib = false; ///< Switch if gain calibration is applied or not
144177
o2::emcal::GainCalibrationFactors* mGainCalibFactors; ///< Gain calibration factors applied to the data before filling the histograms
178+
std::array<double, 17664> mArrGainCalibFactors; ///< array of gain calibration factors
145179
std::shared_ptr<EMCALCalibExtractor> mCalibExtractor; ///< calib extractor
146180

147181
ClassDefNV(EMCALChannelData, 1);

‎Detectors/EMCAL/calibration/include/EMCALCalibration/EMCALTimeCalibData.h‎

Lines changed: 31 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -54,8 +54,13 @@ class EMCALTimeCalibData
5454
EMCALTimeCalibData()
5555
{
5656

57-
mTimeHisto = boost::histogram::make_histogram(boost::histogram::axis::regular<>(EMCALCalibParams::Instance().nBinsTimeAxis_tc, EMCALCalibParams::Instance().minValueTimeAxis_tc, EMCALCalibParams::Instance().maxValueTimeAxis_tc), boost::histogram::axis::regular<>(NCELLS, -0.5, NCELLS - 0.5));
58-
57+
mNThreads = EMCALCalibParams::Instance().nThreads_tc;
58+
mTimeHisto.resize(mNThreads);
59+
mVecNEntriesInHisto.resize(mNThreads);
60+
for (size_t i = 0; i < mNThreads; ++i) {
61+
mTimeHisto[i] = boost::histogram::make_histogram(boost::histogram::axis::regular<>(EMCALCalibParams::Instance().nBinsTimeAxis_tc, EMCALCalibParams::Instance().minValueTimeAxis_tc, EMCALCalibParams::Instance().maxValueTimeAxis_tc), boost::histogram::axis::regular<>(NCELLS, -0.5, NCELLS - 0.5));
62+
mVecNEntriesInHisto[i] = 0;
63+
}
5964
LOG(debug) << "initialize time histogram with " << NCELLS << " cells";
6065
}
6166

@@ -65,7 +70,7 @@ class EMCALTimeCalibData
6570
void fill(const gsl::span<const o2::emcal::Cell> data);
6671

6772
/// \brief Merge the data of two slots.
68-
void merge(const EMCALTimeCalibData* prev);
73+
void merge(EMCALTimeCalibData* prev);
6974

7075
/// \brief Check if enough data for calibration has been accumulated
7176
bool hasEnoughData() const;
@@ -84,29 +89,47 @@ class EMCALTimeCalibData
8489
long unsigned int getNEntriesInHisto() const { return mNEntriesInHisto; }
8590
/// \brief Set the number of entries in histogram
8691
void setNEntriesInHisto(long unsigned int n) { mNEntriesInHisto = n; }
92+
/// \brief Add the number of entries in histogram
93+
void addNEntriesInHisto(long unsigned int n) { mNEntriesInHisto += n; }
8794

8895
/// \brief Get current histogram
89-
boostHisto& getHisto() { return mTimeHisto; }
90-
const boostHisto& getHisto() const { return mTimeHisto; }
96+
const boostHisto& getHisto()
97+
{
98+
// set the summed histogram to one of the existing histograms
99+
mHistoSummed = mTimeHisto[0];
100+
// reset the histogram
101+
mHistoSummed.reset();
102+
// Sum up all entries
103+
for (const auto& h : mTimeHisto) {
104+
mHistoSummed += h;
105+
}
106+
return mHistoSummed;
107+
}
91108

109+
/// \brief Set gain calibration factors applied to the cell energy before filling the histograms
110+
/// \param calibFactors gain calibration object
92111
void setGainCalibFactors(o2::emcal::GainCalibrationFactors* calibFactors)
93112
{
94113
mGainCalibFactors = calibFactors;
95114
mApplyGainCalib = true;
96115
}
97116

98117
/// \brief Set new calibration histogram
99-
void setHisto(boostHisto hist) { mTimeHisto = hist; }
118+
void setHisto(boostHisto hist) { mTimeHisto[0] = hist; }
100119

120+
/// \brief print stream
101121
void PrintStream(std::ostream& stream) const;
102122

103123
/// \brief Actual function where calibration is done. Has to be called in has enough data when enough data is there
104124
o2::emcal::TimeCalibrationParams process();
105125

106126
private:
107-
boostHisto mTimeHisto; ///< histogram with cell time vs. cell ID
127+
unsigned int mNThreads = 1;
128+
std::vector<boostHisto> mTimeHisto; ///< vector of histogram with cell time vs. cell ID (size = number of threads)
129+
boostHisto mHistoSummed; ///< summed histogram (sum of all histograms in mTimeHisto)
108130
int mEvents = 0; ///< current number of events
109-
long unsigned int mNEntriesInHisto = 0; ///< number of entries in histogram
131+
long unsigned int mNEntriesInHisto = 0; ///< Number of entries in the histogram
132+
std::vector<long unsigned int> mVecNEntriesInHisto; ///< Number of entries in the histogram for each thread per event
110133
bool mApplyGainCalib = false; ///< Switch if gain calibration is applied or not
111134
o2::emcal::GainCalibrationFactors* mGainCalibFactors; ///< Gain calibration factors applied to the data before filling the histograms
112135

0 commit comments

Comments
 (0)