12 #include <top/modules/TOPLaserCalibrator/TOPLaserCalibratorModule.h>
13 #include <top/modules/TOPLaserCalibrator/LaserCalibratorFit.h>
16 #include <framework/datastore/StoreArray.h>
19 #include <framework/logging/Logger.h>
23 #include <top/dataobjects/TOPDigit.h>
51 setDescription(
"T0 Calibration using the laser calibration system");
53 std::vector<double> frange = {100, 0., 1.};
55 addParam(
"dataFitOutput", m_dataFitOutput,
"Output root file for data",
57 addParam(
"mcInput", m_mcInput,
"Input root file from MC",
string(
""));
58 addParam(
"channelT0constant", m_chT0C,
"Output of channel T0 constant",
string(
""));
59 addParam(
"fitChannel", m_fitChannel,
"set 0 - 511 to a specific channel; set 0 to all channels in the fit",
61 addParam(
"barID", m_barID,
"ID of TOP module to calibrate");
62 addParam(
"refCh", m_refCh,
"reference channel of T0 constant");
63 addParam(
"fitMethod", m_fitMethod,
"gauss: single gaussian; cb: single Crystal Ball; cb2: double Crystal Ball",
string(
"gauss"));
64 addParam(
"fitRange", m_fitRange,
"fit range[nbins, xmin, xmax]", frange);
66 for (
int i = 0; i < c_NumChannels; ++i) {
71 TOPLaserCalibratorModule::~TOPLaserCalibratorModule()
75 void TOPLaserCalibratorModule::initialize()
77 B2WARNING(
"You are using an old version of the laser fitter, now deprecated. This module has been superseded by the CAF collector TOPLaserCalibratorCollector and the CAF fitter TOPLocalCalFitter.");
78 m_digits.isRequired();
81 void TOPLaserCalibratorModule::beginRun()
85 void TOPLaserCalibratorModule::event()
87 for (
auto& digit : m_digits) {
88 if (m_barID != 0 && digit.getModuleID() != m_barID)
continue;
89 unsigned channel = digit.getChannel();
90 if (channel < c_NumChannels) {
91 auto histo = m_histo[channel];
94 ss <<
"chan" << channel ;
97 string title =
"Times " + name;
98 histo =
new TH1F(name.c_str(), title.c_str(), (
int)m_fitRange[0], m_fitRange[1], m_fitRange[2]);
99 m_histo[channel] = histo;
101 if (digit.getHitQuality() != TOPDigit::EHitQuality::c_Good)
continue;
102 if (digit.getTime() < m_fitRange[1] || digit.getTime() > m_fitRange[2])
continue;
103 if (!digit.isTimeBaseCalibrated())
continue;
104 histo->Fill(digit.getTime());
109 void TOPLaserCalibratorModule::endRun()
113 void TOPLaserCalibratorModule::terminate()
115 std::vector<TH1F*> v_histo;
116 for (
int i = 0; i < c_NumChannels; ++i) {
117 v_histo.push_back(m_histo[i]);
120 double dataT[c_NumChannels] = {0};
121 double dataTErr[c_NumChannels] = {0};
122 double mcT[c_NumChannels] = {0};
128 if (m_fitChannel == c_NumChannels) {
129 for (
int i = 0; i < c_NumChannels; ++i) {
136 dataT[m_fitChannel] = t0fit.
getFitT();
142 auto mcFile =
new TFile(m_mcInput.c_str());
143 auto tree =
static_cast<TTree*
>(mcFile->Get(
"t0MC"));
147 tree->SetBranchAddress(
"channel", &channel_mc);
148 tree->SetBranchAddress(
"maxpos", &maxpos);
149 for (
int i = 0; i < tree->GetEntries(); i++) {
151 mcT[channel_mc] = maxpos;
155 for (
int i = 0; i < c_NumChannels; i++) {
156 if (i == m_refCh)
continue;
157 dataT[i] = dataT[i] - dataT[m_refCh];
158 mcT[i] = mcT[i] - mcT[m_refCh];
166 auto file =
new TFile(m_chT0C.c_str(),
"RECREATE");
167 auto otree =
new TTree(
"chT0",
"Channel T0 Const");
169 unsigned channel = 0;
171 double t0ConstRaw = 0;
172 double fittedTime = 0;
173 double fittedTimeError = 0;
175 double mcCorrection = 0;
176 double t0ConstError = 0;
178 otree->Branch(
"fittedTime", &fittedTime,
"fittedTime/D");
179 otree->Branch(
"fittedTimeError", &fittedTimeError,
"fittedTimeError/D");
180 otree->Branch(
"mcTime", &mcTime,
"mcTime/D");
181 otree->Branch(
"t0ConstRaw", &t0ConstRaw,
"t0ConstRaw/D");
182 otree->Branch(
"mcCorrection", &mcCorrection,
"mcCorrection/D");
183 otree->Branch(
"t0Const", &t0_const,
"t0_const/D");
184 otree->Branch(
"t0ConstError", &t0ConstError,
"t0ConstError/D");
185 otree->Branch(
"channel", &channel,
"channel/I");
186 otree->Branch(
"slot", &m_barID,
"slot/I");
188 for (
int i = 0; i < c_NumChannels; i++) {
190 fittedTime = dataT[i];
193 fittedTime = dataT[i] + dataT[m_refCh];
194 mcTime = mcT[i] + mcT[m_refCh];
196 fittedTimeError = dataTErr[i];
197 t0ConstRaw = dataT[i];
198 mcCorrection = mcT[i];
200 t0_const = dataT[i] - mcT[i];
202 t0ConstError = TMath::Sqrt(dataTErr[i] * dataTErr[i] + dataTErr[m_refCh] * dataTErr[m_refCh]);
205 t0ConstError = dataTErr[m_refCh];