8 #include <svd/calibration/SVDCoGTimeCalibrationAlgorithm.h>
10 #include <svd/dbobjects/SVDCoGCalibrationFunction.h>
11 #include <svd/calibration/SVDCoGTimeCalibrations.h>
18 #include <framework/logging/Logger.h>
21 #include <TFitResult.h>
26 SVDCoGTimeCalibrationAlgorithm::SVDCoGTimeCalibrationAlgorithm(
const std::string& str) :
36 gROOT->SetBatch(
true);
41 std::unique_ptr<TF1> pol1(
new TF1(
"pol1",
"[0] + [1]*x", -10, 80));
42 pol1->SetParameters(-40, 0.9);
43 std::unique_ptr<TF1> pol3(
new TF1(
"pol3",
"[0] + [1]*x + [2]*x*x + [3]*x*x*x", -10, 80));
44 pol3->SetParLimits(0, -200, 0);
45 pol3->SetParLimits(1, 0, 10);
46 pol3->SetParLimits(2, -1, 0);
47 pol3->SetParLimits(3, 0, 1);
48 std::unique_ptr<TF1> pol5(
new TF1(
"pol5",
"[0] + [1]*x + [2]*x*x + [3]*x*x*x + [4]*x*x*x*x + [5]*x*x*x*x*x", -100, 100));
49 pol5->SetParameters(-50, 1.5, 0.01, 0.0001, 0.00001, 0.000001);
53 while (gSystem->GetPathInfo(Form(
"algorithm_6SampleCoG_output_rev_%d.root", cal_rev), info) == 0)
55 std::unique_ptr<TFile> f(
new TFile(Form(
"algorithm_6SampleCoG_output_rev_%d.root", cal_rev),
"RECREATE"));
57 auto m_tree =
new TTree(Form(
"rev_%d", cal_rev),
"RECREATE");
58 int layer_num, ladder_num, sensor_num, view, ndf;
59 float a, b, c, d, a_err, b_err, c_err, d_err, chi2, p;
60 m_tree->Branch(
"layer", &layer_num,
"layer/I");
61 m_tree->Branch(
"ladder", &ladder_num,
"ladder/I");
62 m_tree->Branch(
"sensor", &sensor_num,
"sensor/I");
63 m_tree->Branch(
"isU", &view,
"isU/I");
64 m_tree->Branch(
"a", &a,
"a/F");
65 m_tree->Branch(
"b", &b,
"b/F");
66 m_tree->Branch(
"c", &c,
"c/F");
67 m_tree->Branch(
"d", &d,
"d/F");
68 m_tree->Branch(
"a_err", &a_err,
"a_err/F");
69 m_tree->Branch(
"b_err", &b_err,
"b_err/F");
70 m_tree->Branch(
"c_err", &c_err,
"c_err/F");
71 m_tree->Branch(
"d_err", &d_err,
"d_err/F");
72 m_tree->Branch(
"chi2", &chi2,
"chi2/F");
73 m_tree->Branch(
"ndf", &ndf,
"ndf/I");
74 m_tree->Branch(
"p", &p,
"p/F");
76 auto __hEventT0vsCoG__ = getObjectPtr<TH3F>(
"__hEventT0vsCoG__");
77 auto __hEventT0__ = getObjectPtr<TH2F>(
"__hEventT0__");
78 auto __hEventT0NoSync__ = getObjectPtr<TH2F>(
"__hEventT0NoSync__");
79 auto __hBinToSensorMap__ = getObjectPtr<TH1F>(
"__hBinToSensorMap__");
81 for (
int ij = 0; ij < (__hBinToSensorMap__->GetNbinsX()); ij++) {
86 auto binLabel = __hBinToSensorMap__->GetXaxis()->GetBinLabel(ij + 1);
88 std::sscanf(binLabel,
"L%dL%dS%d%c", &layer_num, &ladder_num, &sensor_num, &side);
93 B2INFO(
"Projecting for Sensor: " << binLabel <<
" with Bin Number: " << ij + 1);
95 __hEventT0vsCoG__->GetZaxis()->SetRange(ij + 1, ij + 1);
96 auto hEventT0vsCoG = (TH2D*)__hEventT0vsCoG__->Project3D(
"yxe");
97 auto hEventT0 = (TH1D*)__hEventT0__->ProjectionX(
"hEventT0_tmp", ij + 1, ij + 1);
98 auto hEventT0nosync = (TH1D*)__hEventT0NoSync__->ProjectionX(
"hEventT0NoSync_tmp", ij + 1, ij + 1);
100 hEventT0vsCoG->SetName(Form(
"eventT0vsCoG__L%dL%dS%d%c", layer_num, ladder_num, sensor_num, side));
101 hEventT0->SetName(Form(
"eventT0__L%dL%dS%d%c", layer_num, ladder_num, sensor_num, side));
102 hEventT0nosync->SetName(Form(
"eventT0nosync__L%dL%dS%d%c", layer_num, ladder_num, sensor_num, side));
104 char sidePN = (side ==
'U' ?
'P' :
'N');
105 hEventT0vsCoG->SetTitle(Form(
"EventT0Sync vs rawTime in %d.%d.%d %c/%c", layer_num, ladder_num, sensor_num, side, sidePN));
106 hEventT0->SetTitle(Form(
"EventT0Sync in %d.%d.%d %c/%c", layer_num, ladder_num, sensor_num, side, sidePN));
107 hEventT0nosync->SetTitle(Form(
"EventT0NoSync in %d.%d.%d %c/%c", layer_num, ladder_num, sensor_num, side, sidePN));
109 hEventT0vsCoG->SetDirectory(0);
110 hEventT0->SetDirectory(0);
111 hEventT0nosync->SetDirectory(0);
113 B2INFO(
"Histogram: " << hEventT0vsCoG->GetName() <<
114 " Entries (n. clusters): " << hEventT0vsCoG->GetEntries());
115 if (layer_num == 3 && hEventT0vsCoG->GetEntries() <
m_minEntries) {
116 B2INFO(
"Histogram: " << hEventT0vsCoG->GetName() <<
117 " Entries (n. clusters): " << hEventT0vsCoG->GetEntries() <<
119 B2WARNING(
"Not enough data, adding one run to the collector");
121 gSystem->Unlink(Form(
"algorithm_6SampleCoG_output_rev_%d.root", cal_rev));
124 if (layer_num != 3 && hEventT0vsCoG->GetEntries() <
m_minEntries / 10) {
125 B2INFO(
"Histogram: " << hEventT0vsCoG->GetName() <<
126 " Entries (n. clusters): " << hEventT0vsCoG->GetEntries() <<
128 B2WARNING(
"Not enough data, adding one run to the collector");
130 gSystem->Unlink(Form(
"algorithm_6SampleCoG_output_rev_%d.root", cal_rev));
133 for (
int i = 1; i <= hEventT0vsCoG->GetNbinsX(); i++) {
134 for (
int j = 1; j <= hEventT0vsCoG->GetNbinsY(); j++) {
135 if (hEventT0vsCoG->GetBinContent(i, j) < max(2,
int(hEventT0vsCoG->GetEntries() * 0.001))) {
136 hEventT0vsCoG->SetBinContent(i, j, 0);
140 TProfile* pfx = hEventT0vsCoG->ProfileX();
141 std::string name =
"pfx_" + std::string(hEventT0vsCoG->GetName());
142 pfx->SetName(name.c_str());
143 TFitResultPtr tfr = pfx->Fit(
"pol3",
"RQSM");
145 pol3->GetParameters(par);
156 timeCal->set_current(1);
158 timeCal->set_pol3parameters(par[0], par[1], par[2], par[3]);
159 payload->set(layer_num, ladder_num, sensor_num,
bool(view), 1, *timeCal);
162 hEventT0vsCoG->Write();
163 hEventT0nosync->Write();
167 delete hEventT0vsCoG;
169 delete hEventT0nosync;
171 if (tfr.Get() ==
nullptr || (tfr->Status() != 0 && tfr->Status() != 4 && tfr->Status() != 4000)) {
173 B2FATAL(
"Fit to the histogram failed in SVDCoGTimeCalibrationAlgorithm. "
174 <<
"Check the 2-D histogram to clarify the reason.");
176 a = par[0]; b = par[1]; c = par[2]; d = par[3];
177 a_err = tfr->ParError(0); b_err = tfr->ParError(1); c_err = tfr->ParError(2); d_err = tfr->ParError(3);
200 float meanRawTimeL3V = 0;
202 auto rawTimeL3V = getObjectPtr<TH1F>(
"hRawTimeL3V");
209 meanRawTimeL3V = rawTimeL3V->GetMean();
216 B2INFO(
"Setting start payload boundary to be the first run ("
217 << currentRun.first <<
"," << currentRun.second <<
")");
223 <<
" to " << meanRawTimeL3V <<
". We are requesting a new payload boundary for ("
224 << currentRun.first <<
"," << currentRun.second <<
")");
Base class for calibration algorithms.
void saveCalibration(TClonesArray *data, const std::string &name)
Store DBArray payload with given name with default IOV.
void setDescription(const std::string &description)
Set algorithm description (in constructor)
EResult
The result of calibration.
@ c_OK
Finished successfuly =0 in Python.
@ c_NotEnoughData
Needs more data =2 in Python.
class to contain the CoG Time calibrations
std::string m_id
Parameter given to set the UniqueID of the payload.
std::optional< float > m_previousRawTimeMeanL3V
CoG time mean of the previous run for V side of layer 3.
float m_minEntries
Set the minimun number of entries required in the histograms of layer 3.
virtual EResult calibrate() override
Run algo on data.
virtual bool isBoundaryRequired(const Calibration::ExpRun ¤tRun) override
If the event T0 changes significantly return true.
float m_allowedTimeShift
Allowed EventT0 shift.
SVDCalibrationsBase< SVDCalibrationsScalar< SVDCoGCalibrationFunction > > t_payload
typedef for the SVDCoGCalibrationFunction payload of all SVD sensors
Abstract base class for different kinds of events.