Belle II Software development
CDCDedxInjectTimeAlgorithm.h
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#pragma once
10
11#include <map>
12#include <vector>
13
14#include <TH1D.h>
15
16#include <calibration/CalibrationAlgorithm.h>
17#include <framework/database/DBObjPtr.h>
18#include <cdc/dbobjects/CDCDedxInjectionTime.h>
19
20// namespace constants
21namespace numdedx {
22 static const int nrings = 2;
23 static const int nvectors = 6;
24}
25
26enum fstatus {fitOK, fitFailed};
27
28namespace Belle2 {
37
38 public:
39
44
49
53 void setMergePayload(bool value = true) {m_isMerge = value;}
54
58 void setMonitoringPlots(bool value = false) {m_ismakePlots = value;}
59
63 void setFitWidth(double value) {m_sigmaR = value;}
64
68 void setMinTracks(int value) {m_thersE = value;}
69
73 void setDedxPars(int value, double min, double max)
74 {
75 m_dedxBins = value;
76 m_dedxMin = min;
77 m_dedxMax = max;
78 }
79
83 void setChiPars(int value, double min, double max)
84 {
85 m_chiBins = value;
86 m_chiMin = min;
87 m_chiMax = max;
88 }
89
93 void fitGaussianWRange(TH1D*& temphist, fstatus& status);
94
98 void getExpRunInfo();
99
103 void defineTimeBins();
104
108 void defineHisto(std::array<std::vector<TH1D*>, numdedx::nrings>& htemp, std::string var);
109
113 void defineTimeHisto(std::array<std::array<TH1D*, 3>, numdedx::nrings>& htemp);
114
118 void checkStatistics(std::array<std::vector<TH1D*>, numdedx::nrings>& hvar);
119
123 void correctBinBias(std::map<int, std::vector<double>>& varcorr, std::map<int, std::vector<double>>& var,
124 std::map<int, std::vector<double>>& time, TH1D*& htimes);
125
129 void createPayload(std::array<double, numdedx::nrings>& scale, std::map<int, std::vector<double>>& vmeans,
130 std::map<int, std::vector<double>>& varscal, std::string svar);
131
135 void getMeanReso(std::array<std::vector<TH1D*>, numdedx::nrings>& hvar,
136 std::map<int, std::vector<double>>& vmeans, std::map<int, std::vector<double>>& vresos);
137
141 void plotEventStats();
142
146 void plotBinLevelDist(std::array<std::vector<TH1D*>, numdedx::nrings>& hvar, std::string var);
147
151 void plotRelConstants(std::map<int, std::vector<double>>& vmeans, std::map<int, std::vector<double>>& vresos,
152 std::map<int, std::vector<double>>& corr, std::string svar);
153
157 void plotTimeStat(std::array<std::vector<TH1D*>, numdedx::nrings>& htime);
158
162 void plotFinalConstants(std::map<int, std::vector<double>>& vmeans, std::map<int, std::vector<double>>& vresos,
163 std::array<double, numdedx::nrings>& scale, std::array<double, numdedx::nrings>& scale_reso);
164
168 void plotInjectionTime(std::array<std::array<TH1D*, 3>, numdedx::nrings>& hvar);
169
173 void setHistStyle(TH1D*& htemp, const int ic, const int is, const double min, const double max)
174 {
175 htemp->SetStats(0);
176 htemp->LabelsDeflate();
177 htemp->SetMarkerColor(ic);
178 htemp->SetMarkerStyle(is);
179 htemp->GetXaxis()->SetLabelOffset(-0.055);
180 htemp->GetYaxis()->SetTitleOffset(0.75);
181 htemp->SetMinimum(min);
182 htemp->SetMaximum(max);
183 }
184
188 std::string getTimeBinLabel(const double& tedges, const int& it)
189 {
190 std::string label = "";
191 if (tedges < 2e4)label = Form("%0.01f-%0.01fK", m_tedges[it] / 1e3, m_tedges[it + 1] / 1e3);
192 else if (tedges < 1e5)label = Form("%0.0f-%0.0fK", m_tedges[it] / 1e3, m_tedges[it + 1] / 1e3);
193 else label = Form("%0.01f-%0.01fM", m_tedges[it] / 1e6, m_tedges[it + 1] / 1e6);
194 return label;
195 }
196
200 void deleteHisto(std::array<std::vector<TH1D*>, numdedx::nrings>& htemp)
201 {
202 for (unsigned int ir = 0; ir < c_rings; ir++) {
203 for (unsigned int it = 0; it < m_tbins; it++) {
204 delete htemp[ir][it];
205 }
206 }
207 }
208
212 void deleteTimeHisto(std::array<std::array<TH1D*, 3>, numdedx::nrings>& htemp)
213 {
214 const int tzoom = 3;
215 for (unsigned int ir = 0; ir < c_rings; ir++) {
216 for (int wt = 0; wt < tzoom; wt++) {
217 delete htemp[ir][wt];
218 }
219 }
220 }
221
225 double getCorrection(unsigned int ring, unsigned int time, std::map<int, std::vector<double>>& vmeans);
226
227
228 protected:
229
233 virtual EResult calibrate() override;
234
235 private:
236
237 static const int c_rings = numdedx::nrings;
238 std::array<std::string, numdedx::nrings> m_sring{"ler", "her"};
239
240 std::vector<double> m_vtedges;
241 std::vector<double> m_vtlocaledges;
242
243 double* m_tedges;
244 unsigned int m_tbins;
245 double m_sigmaR;
246
248 double m_dedxMin;
249 double m_dedxMax;
250
252 double m_chiMin;
253 double m_chiMax;
254
257
261
262 std::string m_prefix;
263 std::string m_suffix;
264
265 std::vector<std::vector<double>> m_vinjPayload;
266
268 };
269
270} // namespace Belle2
void setChiPars(int value, double min, double max)
function to set chi hist parameters
void setDedxPars(int value, double min, double max)
function to set dedx hist parameters
bool m_isMerge
merge payload when rel constant
void plotBinLevelDist(std::array< std::vector< TH1D * >, numdedx::nrings > &hvar, std::string var)
function to draw dedx, chi and time dist.
std::vector< double > m_vtlocaledges
internal time vector
double m_sigmaR
fit dedx dist in sigma range
void correctBinBias(std::map< int, std::vector< double > > &varcorr, std::map< int, std::vector< double > > &var, std::map< int, std::vector< double > > &time, TH1D *&htimes)
function to correct dedx mean/reso and return corrected vector map
void defineHisto(std::array< std::vector< TH1D * >, numdedx::nrings > &htemp, std::string var)
function to define histograms for dedx and time dist.
double getCorrection(unsigned int ring, unsigned int time, std::map< int, std::vector< double > > &vmeans)
function to get the correction factor of mean
std::string m_prefix
string prefix for plot names
void getExpRunInfo()
function to get exp/run information (payload object, plotting)
void setHistStyle(TH1D *&htemp, const int ic, const int is, const double min, const double max)
function to set histogram cosmetics
void plotFinalConstants(std::map< int, std::vector< double > > &vmeans, std::map< int, std::vector< double > > &vresos, std::array< double, numdedx::nrings > &scale, std::array< double, numdedx::nrings > &scale_reso)
function to final constant from merging or abs fits
int m_countR
a hack for running functions once
double * m_tedges
internal time array (copy of vtlocaledges)
bool m_ismakePlots
produce plots for monitoring
std::string m_suffix
string suffix for object names
void setMinTracks(int value)
function to set min # of tracks in time bins (0-40ms)
static const int c_rings
injection ring constants
void createPayload(std::array< double, numdedx::nrings > &scale, std::map< int, std::vector< double > > &vmeans, std::map< int, std::vector< double > > &varscal, std::string svar)
function to store payloads after full calibration
void plotInjectionTime(std::array< std::array< TH1D *, 3 >, numdedx::nrings > &hvar)
function to injection time distributions (HER/LER in three bins)
void plotTimeStat(std::array< std::vector< TH1D * >, numdedx::nrings > &htime)
function to draw time stats
void fitGaussianWRange(TH1D *&temphist, fstatus &status)
function to perform gauss fit for input histogram
void defineTimeBins()
function to set/reset time bins
void plotRelConstants(std::map< int, std::vector< double > > &vmeans, std::map< int, std::vector< double > > &vresos, std::map< int, std::vector< double > > &corr, std::string svar)
function to relative constant from dedx fit mean and chi fit reso
void defineTimeHisto(std::array< std::array< TH1D *, 3 >, numdedx::nrings > &htemp)
function to define injection time bins histograms (monitoring only)
void deleteHisto(std::array< std::vector< TH1D * >, numdedx::nrings > &htemp)
function to delete histograms for dedx and time dist.
void setMergePayload(bool value=true)
function to decide merged vs relative calibration
void deleteTimeHisto(std::array< std::array< TH1D *, 3 >, numdedx::nrings > &htemp)
function to define injection time bins histograms (monitoring only)
void plotEventStats()
function to draw event/track statistics plots
virtual EResult calibrate() override
CDC dE/dx Injection time algorithm.
void checkStatistics(std::array< std::vector< TH1D * >, numdedx::nrings > &hvar)
check statistics for obtaining calibration const.
std::vector< double > m_vtedges
external time vector
std::vector< std::vector< double > > m_vinjPayload
vector to store payload values
bool m_isminStat
flag to merge runs for statistics thershold
CDCDedxInjectTimeAlgorithm()
Constructor: Sets the description, the properties and the parameters of the algorithm.
std::array< std::string, numdedx::nrings > m_sring
injection ring name
DBObjPtr< CDCDedxInjectionTime > m_DBInjectTime
Injection time DB object.
void setMonitoringPlots(bool value=false)
function to enable monitoring plots
int m_thersE
min tracks to start calibration
void getMeanReso(std::array< std::vector< TH1D * >, numdedx::nrings > &hvar, std::map< int, std::vector< double > > &vmeans, std::map< int, std::vector< double > > &vresos)
function to get mean and reso of histogram
std::string getTimeBinLabel(const double &tedges, const int &it)
function to return time label for histograms labeling
void setFitWidth(double value)
function to set fit range (sigma)
EResult
The result of calibration.
CalibrationAlgorithm(const std::string &collectorModuleName)
Constructor - sets the prefix for collected objects (won't be accesses until execute(....
Class for accessing objects in the database.
Definition DBObjPtr.h:21
Abstract base class for different kinds of events.