From 77aa39a8f52a29a04f8fddad2a240d1f6ed76569 Mon Sep 17 00:00:00 2001 From: Luisa Bergmann Date: Tue, 25 Jul 2023 14:14:46 +0200 Subject: [PATCH 01/10] Adding TRD t0 fitting procedure --- DataFormats/Detectors/TRD/CMakeLists.txt | 2 + .../TRD/include/DataFormatsTRD/CalT0.h | 23 ++ .../TRD/include/DataFormatsTRD/Constants.h | 1 + .../TRD/include/DataFormatsTRD/PHData.h | 1 + .../TRD/include/DataFormatsTRD/T0FitHistos.h | 61 ++++ .../Detectors/TRD/src/DataFormatsTRDLinkDef.h | 2 + DataFormats/Detectors/TRD/src/T0FitHistos.cxx | 75 +++++ Detectors/TRD/calibration/CMakeLists.txt | 2 + .../TRDCalibration/CalibrationParams.h | 4 + .../include/TRDCalibration/T0Fit.h | 111 +++++++ Detectors/TRD/calibration/src/T0Fit.cxx | 270 ++++++++++++++++++ .../calibration/src/TRDCalibrationLinkDef.h | 3 + Detectors/TRD/workflow/CMakeLists.txt | 1 + .../workflow/include/TRDWorkflow/T0FitSpec.h | 177 ++++++++++++ Detectors/TRD/workflow/io/CMakeLists.txt | 1 + .../include/TRDWorkflowIO/TRDPHReaderSpec.h | 53 ++++ .../TRD/workflow/io/src/TRDPHReaderSpec.cxx | 83 ++++++ .../TRD/workflow/src/trd-calib-workflow.cxx | 10 + 18 files changed, 880 insertions(+) create mode 100644 DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h create mode 100644 DataFormats/Detectors/TRD/src/T0FitHistos.cxx create mode 100644 Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h create mode 100644 Detectors/TRD/calibration/src/T0Fit.cxx create mode 100644 Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h create mode 100644 Detectors/TRD/workflow/io/include/TRDWorkflowIO/TRDPHReaderSpec.h create mode 100644 Detectors/TRD/workflow/io/src/TRDPHReaderSpec.cxx diff --git a/DataFormats/Detectors/TRD/CMakeLists.txt b/DataFormats/Detectors/TRD/CMakeLists.txt index 7697d12fec1b2..4689fa6a40fc9 100644 --- a/DataFormats/Detectors/TRD/CMakeLists.txt +++ b/DataFormats/Detectors/TRD/CMakeLists.txt @@ -14,6 +14,7 @@ o2_add_library(DataFormatsTRD src/LinkRecord.cxx src/AngularResidHistos.cxx src/GainCalibHistos.cxx + src/T0FitHistos.cxx src/DcsCcdbObjects.cxx src/Tracklet64.cxx src/RawData.cxx @@ -37,6 +38,7 @@ o2_target_root_dictionary(DataFormatsTRD include/DataFormatsTRD/DcsCcdbObjects.h include/DataFormatsTRD/AngularResidHistos.h include/DataFormatsTRD/GainCalibHistos.h + include/DataFormatsTRD/T0FitHistos.h include/DataFormatsTRD/HelperMethods.h include/DataFormatsTRD/Hit.h include/DataFormatsTRD/Digit.h diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h index 30474888a7424..91f251b14b4d5 100644 --- a/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h @@ -32,11 +32,34 @@ class CalT0 ~CalT0() = default; void setT0(int iDet, float t0) { mT0[iDet] = t0; } + void setT0av(float t0) { mT0av = t0; } float getT0(int iDet) const { return mT0[iDet]; } + float getT0av() const { return mT0av; } + float calcT0av() const + { + if (mT0.size() == 0) { + return -1; + } + float sum = 0; + int counts = 0; + for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { + if (mT0[iDet] > -5) { + sum += mT0[iDet]; + ++counts; + } + } + + if (counts > 0) { + return sum / counts; + } else { + return -2; + } + } private: std::array mT0{}; ///< calibrated T0 per TRD chamber + float mT0av{-1}; ///< calibrated average T0 ClassDefNV(CalT0, 1); }; diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/Constants.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/Constants.h index 04f05c402e667..103b830715fff 100644 --- a/DataFormats/Detectors/TRD/include/DataFormatsTRD/Constants.h +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/Constants.h @@ -77,6 +77,7 @@ constexpr double VDRIFTDEFAULT = 1.546; ///< default value for vDrift constexpr double EXBDEFAULT = 0.0; ///< default value for LorentzAngle constexpr int NBINSGAINCALIB = 320; ///< number of bins in the charge (Q0+Q1+Q2) histogram for gain calibration constexpr float MPVDEDXDEFAULT = 42.; ///< default Most Probable Value of TRD dEdx +constexpr float T0DEFAULT = 1.2; ///< default value for t0 // array size to store incoming half cru payload. constexpr int HBFBUFFERMAX = 1048576; ///< max buffer size for data read from a half cru, (all events) diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/PHData.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/PHData.h index de939ef0a0091..dbbc3161f9b02 100644 --- a/DataFormats/Detectors/TRD/include/DataFormatsTRD/PHData.h +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/PHData.h @@ -13,6 +13,7 @@ #define ALICEO2_TRD_PHDATA_H_ #include +#include "Rtypes.h" namespace o2::trd { diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h new file mode 100644 index 0000000000000..eff39bbbac438 --- /dev/null +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h @@ -0,0 +1,61 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file T0FitHistos.h +/// \brief Class to store the TRD PH values for each TRD chamber + +#ifndef ALICEO2_T0FITHISTOS_H +#define ALICEO2_T0FITHISTOS_H + +#include "DataFormatsTRD/Constants.h" +#include "DataFormatsTRD/PHData.h" +#include "Framework/InputRecord.h" +#include "Rtypes.h" +#include +#include +#include + +namespace o2 +{ +namespace trd +{ + +class T0FitHistos +{ + public: + T0FitHistos() = default; + T0FitHistos(const T0FitHistos&) = default; + ~T0FitHistos() = default; + void init(); + void reset(); + auto getDetector(int index) const { return mDet[index]; } + auto getTimeBin(int index) const { return mTB[index]; } + auto getADC(int index) const { return mADC[index]; } + auto getNEntries() const { return mNEntriesTot; } + + void fill(const std::vector data); + void merge(const T0FitHistos* prev); + void print(); + + private: + std::vector mDet{}; + std::vector mTB{}; + std::vector mADC{}; + size_t mNEntriesTot{0}; + bool mInitialized{false}; + + ClassDefNV(T0FitHistos, 1); +}; + +} // namespace trd +} // namespace o2 + +#endif // ALICEO2_T0FITHISTOS_H diff --git a/DataFormats/Detectors/TRD/src/DataFormatsTRDLinkDef.h b/DataFormats/Detectors/TRD/src/DataFormatsTRDLinkDef.h index 8f258914d5a48..3c3b93e58ce40 100644 --- a/DataFormats/Detectors/TRD/src/DataFormatsTRDLinkDef.h +++ b/DataFormats/Detectors/TRD/src/DataFormatsTRDLinkDef.h @@ -33,6 +33,7 @@ #pragma link C++ class o2::trd::NoiseStatusMCM + ; #pragma link C++ class o2::trd::AngularResidHistos + ; #pragma link C++ class o2::trd::GainCalibHistos + ; +#pragma link C++ class o2::trd::T0FitHistos + ; #pragma link C++ class o2::trd::CalVdriftExB + ; #pragma link C++ class o2::trd::CalGain + ; #pragma link C++ class o2::trd::CalT0 + ; @@ -51,6 +52,7 @@ #pragma link C++ class std::vector < o2::trd::Digit> + ; #pragma link C++ class std::vector < o2::trd::AngularResidHistos> + ; #pragma link C++ class std::vector < o2::trd::GainCalibHistos> + ; +#pragma link C++ class std::vector < o2::trd::T0FitHistos> + ; #pragma link C++ class std::vector < o2::trd::PHData> + ; #pragma link C++ class std::vector < o2::trd::KrCluster> + ; #pragma link C++ class std::vector < o2::trd::KrClusterTriggerRecord> + ; diff --git a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx new file mode 100644 index 0000000000000..3b62525393fe8 --- /dev/null +++ b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx @@ -0,0 +1,75 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file T0FitHistos.cxx +/// \brief Class to store the TRD PH values for each chamber + +#include "DataFormatsTRD/T0FitHistos.h" +#include +#include + +using namespace o2::trd; +using namespace o2::trd::constants; + +void T0FitHistos::init() +{ + mDet.resize(0); + mTB.resize(0); + mADC.resize(0); + mInitialized = true; +} + +void T0FitHistos::reset() +{ + mDet.resize(0); + mTB.resize(0); + mADC.resize(0); + mNEntriesTot = 0; +} + +void T0FitHistos::fill(const std::vector data) +{ + + if (!mInitialized) { + init(); + } + + for (auto ph : data) { + int det = ph.getDetector(); + int tb = ph.getTimebin(); + int adc = ph.getADC(); + + if (ph.getNneighbours() != 2) { + continue; + } + + mDet.push_back(det); + mTB.push_back(tb); + mADC.push_back(adc); + ++mNEntriesTot; + } +} + +void T0FitHistos::merge(const T0FitHistos* prev) +{ + auto sizePrev = (int)prev->getNEntries(); + + for (int i = 0; i < sizePrev; ++i) { + mDet.push_back(prev->getDetector(i)); + mTB.push_back(prev->getTimeBin(i)); + mADC.push_back(prev->getADC(i)); + } +} + +void T0FitHistos::print() +{ + LOG(info) << "There are " << mNEntriesTot << " entries in the container"; +} diff --git a/Detectors/TRD/calibration/CMakeLists.txt b/Detectors/TRD/calibration/CMakeLists.txt index 4a9af0739a328..52444d2855b1f 100644 --- a/Detectors/TRD/calibration/CMakeLists.txt +++ b/Detectors/TRD/calibration/CMakeLists.txt @@ -16,6 +16,7 @@ o2_add_library(TRDCalibration src/CalibratorVdExB.cxx src/CalibratorGain.cxx src/CalibratorNoise.cxx + src/T0Fit.cxx src/PadCalibCCDBBuilder.cxx src/KrClusterFinder.cxx src/DCSProcessor.cxx @@ -35,6 +36,7 @@ o2_add_library(TRDCalibration include/TRDCalibration/CalibratorVdExB.h include/TRDCalibration/CalibratorGain.h include/TRDCalibration/CalibratorNoise.h + include/TRDCalibration/T0Fit.h include/TRDCalibration/CalibrationParams.h include/TRDCalibration/PadCalibCCDBBuilder.h include/TRDCalibration/KrClusterFinder.h diff --git a/Detectors/TRD/calibration/include/TRDCalibration/CalibrationParams.h b/Detectors/TRD/calibration/include/TRDCalibration/CalibrationParams.h index 6ae7c6c727363..677673a7f85f3 100644 --- a/Detectors/TRD/calibration/include/TRDCalibration/CalibrationParams.h +++ b/Detectors/TRD/calibration/include/TRDCalibration/CalibrationParams.h @@ -39,6 +39,10 @@ struct TRDCalibParams : public o2::conf::ConfigurableParamHelper float dEdxTPCMin = 30.; float dEdxTPCMax = 70.; + // For t0 fits + size_t minEntriesTotalT0Fit = 1'000'000; ///< minimum number of entries in inclusive distribution for (meaningful) t0 fit + size_t minEntriesChamberT0Fit = 30'000; ///< minimum number of entries in one chamber for (meaningful) t0 fit + // Creation of PH plots float pileupCut = 0.7; diff --git a/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h b/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h new file mode 100644 index 0000000000000..ef08b3b9f80f6 --- /dev/null +++ b/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h @@ -0,0 +1,111 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file T0Fit.h +/// \brief Fits the TRD PH distributions to extract the t0 value +/// \author Luisa Bergmann + +#ifndef O2_TRD_T0FIT_H +#define O2_TRD_T0FIT_H + +#include "DataFormatsTRD/T0FitHistos.h" +#include "DetectorsCalibration/TimeSlotCalibration.h" +#include "DetectorsCalibration/TimeSlot.h" +#include "DataFormatsTRD/Constants.h" +#include "CCDB/CcdbObjectInfo.h" +#include "DataFormatsTRD/CalT0.h" +#include "TRDCalibration/CalibrationParams.h" + +#include "Rtypes.h" +#include "TH2.h" +#include "TH1.h" +#include "TProfile.h" +#include "TF1.h" +#include "Fit/Fitter.h" +#include "TFile.h" +#include "TTree.h" + +#include +#include +#include + +namespace o2 +{ +namespace trd +{ +//______________________________________________________________________________________________ +struct ErfLandauChi2Functor { + double operator()(const double* par) const; + std::vector x; ///< x-value (time-bin) of adc profile + std::vector y; ///< y-value (av. adc) for corresp. time-bin + float lowerBoundFit; ///< lower bound for fit + float upperBoundFit; ///< upper bound for fit +}; + +//______________________________________________________________________________________________ +class T0Fit final : public o2::calibration::TimeSlotCalibration +{ + using Slot = o2::calibration::TimeSlot; + + public: + T0Fit() = default; + ~T0Fit() final = default; + + bool hasEnoughData(const Slot& slot) const final { return slot.getContainer()->getNEntries() >= mMinEntriesTotal; } + void initOutput() final; + void finalizeSlot(Slot& slot) final; + Slot& emplaceNewSlot(bool front, TFType tStart, TFType tEnd) final; + + /// (Re-)Creates a file "trd_t0fit.root". This lets continually fill + /// a tree with the fit results. + void createOutputFile(); + + /// Close the output file. E.g. First write the tree to the file and let the + /// smart pointers take care of closing the file. + void closeOutputFile(); + + const std::vector& getCcdbObjectVector() const { return mObjectVector; } + std::vector& getCcdbObjectInfoVector() { return mInfoVector; } + + void initProcessing(); + + /// Initialize the fit values once with the previous valid ones if they are + /// available. + void retrievePrev(o2::framework::ProcessingContext& pc); + + private: + bool mInitDone{false}; ///< flag to avoid creating output etc multiple times + const TRDCalibParams& mParams{TRDCalibParams::Instance()}; ///< reference to calibration parameters + size_t mMinEntriesTotal{mParams.minEntriesTotalT0Fit}; ///< minimum total number of angular deviations + size_t mMinEntriesChamber{mParams.minEntriesChamberT0Fit}; ///< minimum number of angular deviations per chamber for accepting refitted value + bool mEnableOutput{false}; ///< enable output in a root file instead of the ccdb + std::unique_ptr mOutFile{nullptr}; ///< output file + std::unique_ptr mOutTree{nullptr}; ///< output tree + ErfLandauChi2Functor mFitFunctor; ///< used for minimization process, provides chi2 estimate + ROOT::Fit::Fitter mFitter; ///< instance of the ROOT fitter + std::array mParamsStart; ///< Starting parameters for fit + std::unique_ptr mFuncErfLandau; ///< helper function to calculate the t0 value after the fitting procedure + std::array t0_chambers; ///< t0 values of the individual chambers + float t0_average{-5}; ///< average t0 value across all chambers + + std::vector mInfoVector; ///< vector of CCDB infos; each element is filled with CCDB description of accompanying CCDB calibration object + std::vector mObjectVector; ///< vector of CCDB calibration objects; the extracted t0 per chamber and average for given slot + + std::unique_ptr adcProfIncl; ///< profile that holds inclusive PH spectrum + std::array, o2::trd::constants::MAXCHAMBER> adcProfDet; ///< array of profiles for PH spectrum of each chamber + + ClassDefNV(T0Fit, 1); +}; + +} // namespace trd +} // namespace o2 + +#endif // O2_TRD_T0FIT_H diff --git a/Detectors/TRD/calibration/src/T0Fit.cxx b/Detectors/TRD/calibration/src/T0Fit.cxx new file mode 100644 index 0000000000000..17030842bd981 --- /dev/null +++ b/Detectors/TRD/calibration/src/T0Fit.cxx @@ -0,0 +1,270 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file T0Fit.cxx +/// \brief Fits the TRD PH distributions to extract the t0 value +/// \author Luisa Bergmann + +#include "TRDCalibration/T0Fit.h" +#include "Framework/ProcessingContext.h" +#include "Framework/TimingInfo.h" +#include "Framework/InputRecord.h" +#include "DataFormatsTRD/Constants.h" +#include "TRDBase/GeometryBase.h" +#include "TStopwatch.h" +#include "CCDB/CcdbApi.h" +#include "CCDB/BasicCCDBManager.h" +#include "CommonUtils/NameConf.h" +#include "CommonUtils/MemFileHelper.h" +#include +#include + +#include "TMath.h" +#include "TCanvas.h" + +#include +#include +#include +#include + +using namespace o2::trd::constants; + +namespace o2::trd +{ +//______________________________________________________________________________________________ +double ErfLandauChi2Functor::operator()(const double* par) const +{ + // provides chi2 estimate of PH profile comparing the y-value of profiles to + // par[0]*ROOT::Math::erf(x) + par[1]*TMath::Landau(x,par[2],par[3]) + // + // par[0] : offset + // par[1] : amplitude + // par[2] : mean + // par[3] : sigma + + double sum2 = 0; + + for (int i = lowerBoundFit; i <= upperBoundFit; ++i) { + + if (TMath::Abs(y[i]) < 1e-7) { + continue; + } + + double funcVal = par[0] * TMath::Erf(x[i]) + par[1] * TMath::Landau(x[i], par[2], par[3]); + + sum2 += TMath::Power(y[i] - funcVal, 2) / y[i]; + // sum2 += TMath::Power((y[i] - funcVal)/yErr[i],2); + } + + return sum2; +} + +//______________________________________________________________________________________________ +using Slot = o2::calibration::TimeSlot; + +void T0Fit::initOutput() +{ + // reset the CCDB output vectors + mInfoVector.clear(); + mObjectVector.clear(); +} + +void T0Fit::initProcessing() +{ + if (mInitDone) { + return; + } + + adcProfIncl = std::make_unique("adcProfIncl", "adcProfIncl", 30, -0.5, 29.5); + for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { + adcProfDet[iDet] = std::make_unique(Form("adcProfDet_%i", iDet), Form("adcProfDet_%i", iDet), 30, -0.5, 29.5); + } + + mFitFunctor.lowerBoundFit = 0; + mFitFunctor.upperBoundFit = 5; + + double mParamsStart[4]; + mParamsStart[0] = 100; + mParamsStart[1] = 500; + mParamsStart[2] = 0.5; + mParamsStart[3] = 0.5; + + mFitter.SetFCN(4, mFitFunctor, mParamsStart); + + ROOT::Math::MinimizerOptions opt; + opt.SetMinimizerType("Minuit2"); + opt.SetMinimizerAlgorithm("Migrad"); + opt.SetPrintLevel(0); + opt.SetMaxFunctionCalls(1000); + opt.SetTolerance(.001); + mFitter.Config().SetMinimizerOptions(opt); + + mFuncErfLandau = std::make_unique( + "fErfLandau", [&](double* x, double* par) { return par[0] * TMath::Erf(x[0]) + par[1] * TMath::Landau(x[0], par[2], par[3]); }, mFitFunctor.lowerBoundFit, mFitFunctor.upperBoundFit, 4); + + // set tree addresses + if (mEnableOutput) { + mOutTree->Branch("t0_chambers", &t0_chambers); + mOutTree->Branch("t0_average", &t0_average); + } + + mInitDone = true; +} + +void T0Fit::retrievePrev(o2::framework::ProcessingContext& pc) +{ + // We either get a pointer to a valid object from the last ~hour or to the default object + // which is always present. The first has precedence over the latter. + auto dataCalT0 = pc.inputs().get("calt0"); + std::string msg = "Default Object"; + // We check if the object we got is the default one by comparing it to the defaults. + for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { + if (dataCalT0->getT0(iDet) != constants::T0DEFAULT) { + msg = "Previous Object"; + break; + } + } + LOG(info) << "FitInstance: From CCDB retrieved " << msg; + + // Here we set each entry regardless if it is the default or not. + for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { + t0_chambers[iDet] = dataCalT0->getT0(iDet); + } + t0_average = dataCalT0->getT0av(); +} + +void T0Fit::finalizeSlot(Slot& slot) +{ + // do actual fits for the data provided in the given slot + + TStopwatch timer; + timer.Start(); + initProcessing(); + + // get data and fill profiles + auto dataPH = slot.getContainer(); + + for (int iEntry = 0; iEntry < dataPH->getNEntries(); ++iEntry) { + int iDet = dataPH->getDetector(iEntry); + adcProfIncl->Fill(dataPH->getTimeBin(iEntry), dataPH->getADC(iEntry)); + adcProfDet[iDet]->Fill(dataPH->getTimeBin(iEntry), dataPH->getADC(iEntry)); + } + + // do fits + // inclusive distribution + mFitFunctor.x.clear(); + mFitFunctor.y.clear(); + for (int i = 1; i <= 30; ++i) { + mFitFunctor.x.push_back(i - 1); + mFitFunctor.y.push_back(adcProfIncl->GetBinContent(i)); + } + + mFitter.FitFCN(); + auto fitResult = mFitter.Result(); + + if (fitResult.IsValid()) { + mFuncErfLandau->SetParameters(fitResult.GetParams()[0], fitResult.GetParams()[1], fitResult.GetParams()[2], fitResult.GetParams()[3]); + t0_average = mFuncErfLandau->GetMaximumX(); + } else { + LOG(info) << "t0 fit for inclusive distribtion is not valid, set to -5"; + t0_average = -5; + } + + // single chambers + for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { + if (adcProfDet[iDet]->GetEntries() < mMinEntriesChamber) { + LOG(info) << "not enough entries in chamber " << iDet << " for t0 fit, set to -5"; + t0_chambers[iDet] = -5; + continue; + } + + mFitFunctor.x.clear(); + mFitFunctor.y.clear(); + for (int i = 1; i <= 30; ++i) { + mFitFunctor.x.push_back(i - 1); + mFitFunctor.y.push_back(adcProfDet[iDet]->GetBinContent(i)); + } + + mFitter.FitFCN(); + fitResult = mFitter.Result(); + + if (fitResult.IsValid()) { + mFuncErfLandau->SetParameters(fitResult.GetParams()[0], fitResult.GetParams()[1], fitResult.GetParams()[2], fitResult.GetParams()[3]); + t0_chambers[iDet] = mFuncErfLandau->GetMaximumX(); + } else { + LOG(info) << "t0 fit for chamber " << iDet << " is not valid, set to -5"; + t0_chambers[iDet] = -5; + } + } + + // fill tree + if (mEnableOutput) { + mOutTree->Fill(); + + LOGF(debug, "Fit result for inclusive distribution: t0 = ", t0_average); + for (int iDet = 0; iDet < MAXCHAMBER; ++iDet) { + LOGF(debug, "Fit result for chamber %i: t0 = ", iDet, t0_chambers[iDet]); + } + } + + // assemble CCDB object + CalT0 t0Object; + t0Object.setT0av(t0_average); + for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { + t0Object.setT0(iDet, t0_chambers[iDet]); + } + auto clName = o2::utils::MemFileHelper::getClassName(t0Object); + auto flName = o2::ccdb::CcdbApi::generateFileName(clName); + std::map metadata; // TODO: do we want to store any meta data? + long startValidity = slot.getStartTimeMS(); + mInfoVector.emplace_back("TRD/Calib/CalT0", clName, flName, metadata, startValidity, startValidity + o2::ccdb::CcdbObjectInfo::HOUR); + mObjectVector.push_back(t0Object); +} + +Slot& T0Fit::emplaceNewSlot(bool front, TFType tStart, TFType tEnd) +{ + auto& container = getSlots(); + auto& slot = front ? container.emplace_front(tStart, tEnd) : container.emplace_back(tStart, tEnd); + slot.setContainer(std::make_unique()); + return slot; +} + +void T0Fit::createOutputFile() +{ + mEnableOutput = true; + mOutFile = std::make_unique("trd_calt0.root", "RECREATE"); + if (mOutFile->IsZombie()) { + LOG(error) << "Failed to create output file!"; + mEnableOutput = false; + return; + } + mOutTree = std::make_unique("calib", "t0 values"); + LOG(info) << "Created output file trd_calt0.root"; +} + +void T0Fit::closeOutputFile() +{ + if (!mEnableOutput) { + return; + } + + try { + mOutFile->cd(); + mOutTree->Write(); + mOutTree.reset(); + mOutFile->Close(); + mOutFile.reset(); + } catch (std::exception const& e) { + LOG(error) << "Failed to write t0-value data file, reason: " << e.what(); + } + mEnableOutput = false; +} +} // namespace o2::trd diff --git a/Detectors/TRD/calibration/src/TRDCalibrationLinkDef.h b/Detectors/TRD/calibration/src/TRDCalibrationLinkDef.h index 52556ab90cebf..26a808fd15860 100644 --- a/Detectors/TRD/calibration/src/TRDCalibrationLinkDef.h +++ b/Detectors/TRD/calibration/src/TRDCalibrationLinkDef.h @@ -20,8 +20,11 @@ #pragma link C++ class o2::calibration::TimeSlotCalibration < o2::trd::AngularResidHistos> + ; #pragma link C++ class o2::calibration::TimeSlot < o2::trd::GainCalibHistos> + ; #pragma link C++ class o2::calibration::TimeSlotCalibration < o2::trd::GainCalibHistos> + ; +#pragma link C++ class o2::calibration::TimeSlot < o2::trd::T0FitHistos> + ; +#pragma link C++ class o2::calibration::TimeSlotCalibration < o2::trd::T0FitHistos> + ; #pragma link C++ class o2::trd::CalibratorNoise + ; #pragma link C++ class o2::trd::CalibratorGain + ; +#pragma link C++ class o2::trd::T0Fit + ; #pragma link C++ class o2::trd::ChannelInfoDetailed + ; #pragma link C++ class o2::trd::TrackBasedCalib + ; #pragma link C++ class o2::trd::PadCalibCCDBBuilder + ; diff --git a/Detectors/TRD/workflow/CMakeLists.txt b/Detectors/TRD/workflow/CMakeLists.txt index f71995a8fed00..069e5f11dc00e 100644 --- a/Detectors/TRD/workflow/CMakeLists.txt +++ b/Detectors/TRD/workflow/CMakeLists.txt @@ -27,6 +27,7 @@ o2_add_library(TRDWorkflow include/TRDWorkflow/TRDPulseHeightSpec.h include/TRDWorkflow/TRDGlobalTrackingQCSpec.h include/TRDWorkflow/NoiseCalibSpec.h + include/TRDWorkflow/T0FitSpec.h PUBLIC_LINK_LIBRARIES O2::Framework O2::DPLUtils O2::Steer O2::Algorithm diff --git a/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h b/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h new file mode 100644 index 0000000000000..decd9b8220f28 --- /dev/null +++ b/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h @@ -0,0 +1,177 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +#ifndef O2_TRD_T0FITSPEC_H +#define O2_TRD_T0FITSPEC_H + +/// \file T0FitSpec.h +/// \brief DPL device for steering the TRD t0 fits +/// \author Luisa Bergmann + +#include "TRDCalibration/T0Fit.h" +#include "DetectorsCalibration/Utils.h" +#include "DataFormatsTRD/PHData.h" +#include "DataFormatsTRD/CalT0.h" +#include "CommonUtils/MemFileHelper.h" +#include "Framework/Task.h" +#include "Framework/ConfigParamRegistry.h" +#include "Framework/ControlService.h" +#include "Framework/WorkflowSpec.h" +#include "CCDB/CcdbApi.h" +#include "CCDB/CcdbObjectInfo.h" +#include "Framework/CCDBParamSpec.h" +#include "DetectorsBase/GRPGeomHelper.h" +#include + +#include "TH2.h" +#include "TH1.h" + +using namespace o2::framework; + +namespace o2 +{ +namespace calibration +{ + +class T0FitDevice : public o2::framework::Task +{ + public: + T0FitDevice(std::shared_ptr req) : mCCDBRequest(req) {} + void init(o2::framework::InitContext& ic) final + { + o2::base::GRPGeomHelper::instance().setRequest(mCCDBRequest); + auto slotL = ic.options().get("sec-per-slot"); + auto delay = ic.options().get("max-delay"); + + mFitInstance = std::make_unique(); + mFitInstance->setSlotLengthInSeconds(slotL); + mFitInstance->setMaxSlotsDelay(delay); + if (ic.options().get("enable-root-output")) { + mFitInstance->createOutputFile(); + } + } + + void finaliseCCDB(o2::framework::ConcreteDataMatcher& matcher, void* obj) final + { + o2::base::GRPGeomHelper::instance().finaliseCCDB(matcher, obj); + } + + void run(o2::framework::ProcessingContext& pc) final + { + const auto& tinfo = pc.services().get(); + if (tinfo.globalRunNumberChanged) { // new run is starting + mRunStopRequested = false; + mFitInstance->retrievePrev(pc); // SOR initialization is performed here + } + if (mRunStopRequested) { + return; + } + o2::base::GRPGeomHelper::instance().checkUpdates(pc); + + auto dataT0Fit = pc.inputs().get>("input"); + o2::base::TFIDInfoHelper::fillTFIDInfo(pc, mFitInstance->getCurrentTFInfo()); + LOG(detail) << "Processing TF " << mFitInstance->getCurrentTFInfo().tfCounter << " with " << dataT0Fit.size() << " PHData entries"; + mFitInstance->process(dataT0Fit); + + if (pc.transitionState() == TransitionHandlingState::Requested) { + LOG(info) << "Run stop requested, finalizing"; + mRunStopRequested = true; + mFitInstance->checkSlotsToFinalize(o2::calibration::INFINITE_TF); + mFitInstance->closeOutputFile(); + } + sendOutput(pc.outputs()); + } + + void endOfStream(o2::framework::EndOfStreamContext& ec) final + { + if (mRunStopRequested) { + return; + } + mFitInstance->checkSlotsToFinalize(o2::calibration::INFINITE_TF); + mFitInstance->closeOutputFile(); + sendOutput(ec.outputs()); + } + + void stop() final + { + mFitInstance->closeOutputFile(); + } + + private: + std::unique_ptr mFitInstance; + std::shared_ptr mCCDBRequest; + bool mRunStopRequested = false; // flag that run was stopped (and the last output is sent) + //________________________________________________________________ + void sendOutput(DataAllocator& output) + { + // extract CCDB infos and calibration objects, convert it to TMemFile and send them to the output + + using clbUtils = o2::calibration::Utils; + const auto& payloadVec = mFitInstance->getCcdbObjectVector(); + auto& infoVec = mFitInstance->getCcdbObjectInfoVector(); // use non-const version as we update it + + assert(payloadVec.size() == infoVec.size()); + + for (uint32_t i = 0; i < payloadVec.size(); i++) { + auto& w = infoVec[i]; + auto image = o2::ccdb::CcdbApi::createObjectImage(&payloadVec[i], &w); + LOG(info) << "Sending object " << w.getPath() << "/" << w.getFileName() << " of size " << image->size() + << " bytes, valid for " << w.getStartValidityTimestamp() << " : " << w.getEndValidityTimestamp(); + + output.snapshot(Output{clbUtils::gDataOriginCDBPayload, "CALT0", i}, *image.get()); // vector + output.snapshot(Output{clbUtils::gDataOriginCDBWrapper, "CALT0", i}, w); // root-serialized + } + if (payloadVec.size()) { + mFitInstance->initOutput(); // reset the outputs once they are already sent + } + } +}; + +} // namespace calibration + +namespace framework +{ + +DataProcessorSpec getTRDT0FitSpec() +{ + using device = o2::calibration::T0FitDevice; + using clbUtils = o2::calibration::Utils; + + std::vector outputs; + outputs.emplace_back(ConcreteDataTypeMatcher{o2::calibration::Utils::gDataOriginCDBPayload, "CALT0"}, Lifetime::Sporadic); + outputs.emplace_back(ConcreteDataTypeMatcher{o2::calibration::Utils::gDataOriginCDBWrapper, "CALT0"}, Lifetime::Sporadic); + std::vector inputs; + inputs.emplace_back("input", "TRD", "PULSEHEIGHT"); + inputs.emplace_back("calt0", "TRD", "CALT0", 0, Lifetime::Condition, ccdbParamSpec("TRD/Calib/CalT0")); + auto ccdbRequest = std::make_shared(true, // orbitResetTime + true, // GRPECS=true + false, // GRPLHCIF + false, // GRPMagField + false, // askMatLUT + o2::base::GRPGeomRequest::None, // geometry + inputs); + + return DataProcessorSpec{ + "calib-t0-fit", + inputs, + outputs, + AlgorithmSpec{adaptFromTask(ccdbRequest)}, + Options{ + {"sec-per-slot", VariantType::UInt32, 900u, {"number of seconds per calibration time slot"}}, + {"max-delay", VariantType::UInt32, 2u, {"number of slots in past to consider"}}, + {"enable-root-output", VariantType::Bool, false, {"output t0 values to root file"}}, + }}; +} + +} // namespace framework +} // namespace o2 + +#endif // O2_TRD_T0FITSPEC_H diff --git a/Detectors/TRD/workflow/io/CMakeLists.txt b/Detectors/TRD/workflow/io/CMakeLists.txt index dfdae1edcc706..e27dba98b477b 100644 --- a/Detectors/TRD/workflow/io/CMakeLists.txt +++ b/Detectors/TRD/workflow/io/CMakeLists.txt @@ -20,6 +20,7 @@ o2_add_library(TRDWorkflowIO src/TRDTrackWriterSpec.cxx src/TRDTrackReaderSpec.cxx src/TRDCalibWriterSpec.cxx + src/TRDPHReaderSpec.cxx src/KrClusterWriterSpec.cxx PUBLIC_LINK_LIBRARIES O2::DataFormatsTRD O2::SimulationDataFormat O2::DPLUtils O2::GPUDataTypeHeaders O2::DataFormatsTPC) diff --git a/Detectors/TRD/workflow/io/include/TRDWorkflowIO/TRDPHReaderSpec.h b/Detectors/TRD/workflow/io/include/TRDWorkflowIO/TRDPHReaderSpec.h new file mode 100644 index 0000000000000..5e71f42f9fb19 --- /dev/null +++ b/Detectors/TRD/workflow/io/include/TRDWorkflowIO/TRDPHReaderSpec.h @@ -0,0 +1,53 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// @file TRDPHReaderSpec.h + +#ifndef O2_TRD_PHREADER +#define O2_TRD_PHREADER + +#include "TFile.h" +#include "TTree.h" + +#include "Framework/DataProcessorSpec.h" +#include "Framework/Task.h" +#include "DataFormatsTRD/PHData.h" + +namespace o2 +{ +namespace trd +{ + +class TRDPHReader : public o2::framework::Task +{ + public: + TRDPHReader() = default; + ~TRDPHReader() override = default; + void init(o2::framework::InitContext& ic) final; + void run(o2::framework::ProcessingContext& pc) final; + + private: + void connectTree(); + std::unique_ptr mFile; + std::unique_ptr mTree; + std::string mInFileName{"trd_PH.root"}; + std::string mInTreeName{"ph"}; + std::vector mPHValues, *mPHValuesPtr = &mPHValues; ///< to be used for branch address +}; + +/// create a processor spec +/// read TRD calibration data from a root file +framework::DataProcessorSpec getTRDPHReaderSpec(); + +} // namespace trd +} // namespace o2 + +#endif /* O2_TRD_PHREADER */ diff --git a/Detectors/TRD/workflow/io/src/TRDPHReaderSpec.cxx b/Detectors/TRD/workflow/io/src/TRDPHReaderSpec.cxx new file mode 100644 index 0000000000000..235c58c2ad302 --- /dev/null +++ b/Detectors/TRD/workflow/io/src/TRDPHReaderSpec.cxx @@ -0,0 +1,83 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// @file TRDPHReaderSpec.cxx + +#include "TRDWorkflowIO/TRDPHReaderSpec.h" + +#include "Framework/ControlService.h" +#include "Framework/ConfigParamRegistry.h" +#include "fairlogger/Logger.h" +#include "CommonUtils/NameConf.h" + +using namespace o2::framework; + +namespace o2 +{ +namespace trd +{ + +void TRDPHReader::init(InitContext& ic) +{ + // get the option from the init context + LOG(info) << "Init TRD PH reader!"; + mInFileName = o2::utils::Str::concat_string(o2::utils::Str::rectifyDirectory(ic.options().get("input-dir")), + ic.options().get("trd-ph-infile")); + mInTreeName = o2::utils::Str::concat_string(o2::utils::Str::rectifyDirectory(ic.options().get("input-dir")), + ic.options().get("treename")); + connectTree(); +} + +void TRDPHReader::connectTree() +{ + mTree.reset(nullptr); // in case it was already loaded + mFile.reset(TFile::Open(mInFileName.c_str())); + assert(mFile && !mFile->IsZombie()); + mTree.reset((TTree*)mFile->Get(mInTreeName.c_str())); + assert(mTree); + mTree->SetBranchAddress("values", &mPHValuesPtr); + LOG(info) << "Loaded tree from " << mInFileName << " with " << mTree->GetEntries() << " entries"; +} + +void TRDPHReader::run(ProcessingContext& pc) +{ + auto currEntry = mTree->GetReadEntry() + 1; + assert(currEntry < mTree->GetEntries()); + mTree->GetEntry(currEntry); + + LOG(info) << "Pushing vector of PH values filled with " << mPHValues.size() << " entries at tree entry " << currEntry; + pc.outputs().snapshot(Output{o2::header::gDataOriginTRD, "PULSEHEIGHT", 0, Lifetime::Timeframe}, mPHValues); + + if (mTree->GetReadEntry() + 1 >= mTree->GetEntries()) { + pc.services().get().endOfStream(); + pc.services().get().readyToQuit(QuitRequest::Me); + } +} + +DataProcessorSpec getTRDPHReaderSpec() +{ + std::vector outputs; + outputs.emplace_back(o2::header::gDataOriginTRD, "PULSEHEIGHT", 0, Lifetime::Timeframe); + + return DataProcessorSpec{ + "TRDPHReader", + Inputs{}, + outputs, + AlgorithmSpec{adaptFromTask()}, + Options{ + {"trd-ph-infile", VariantType::String, "trd_PH.root", {"Name of the input file"}}, + {"input-dir", VariantType::String, "none", {"Input directory"}}, + {"treename", VariantType::String, "ph", {"Name of top-level TTree"}}, + }}; +} + +} // namespace trd +} // namespace o2 diff --git a/Detectors/TRD/workflow/src/trd-calib-workflow.cxx b/Detectors/TRD/workflow/src/trd-calib-workflow.cxx index f3eb520873dbd..6a119d59272d7 100644 --- a/Detectors/TRD/workflow/src/trd-calib-workflow.cxx +++ b/Detectors/TRD/workflow/src/trd-calib-workflow.cxx @@ -12,9 +12,11 @@ #include "Framework/DataProcessorSpec.h" #include "TRDWorkflowIO/TRDCalibReaderSpec.h" #include "TRDWorkflowIO/TRDDigitReaderSpec.h" +#include "TRDWorkflowIO/TRDPHReaderSpec.h" #include "TRDWorkflow/VdAndExBCalibSpec.h" #include "TRDWorkflow/GainCalibSpec.h" #include "TRDWorkflow/NoiseCalibSpec.h" +#include "TRDWorkflow/T0FitSpec.h" #include "CommonUtils/ConfigurableParam.h" using namespace o2::framework; @@ -28,6 +30,7 @@ void customize(std::vector& workflowOptions) {"vDriftAndExB", o2::framework::VariantType::Bool, false, {"enable vDrift and ExB calibration"}}, {"noise", o2::framework::VariantType::Bool, false, {"enable noise and pad status calibration"}}, {"gain", o2::framework::VariantType::Bool, false, {"enable gain calibration"}}, + {"t0", o2::framework::VariantType::Bool, false, {"enable t0 fit"}}, {"calib-dds-collection-index", VariantType::Int, -1, {"allow only single collection to produce calibration objects (use -1 for no limit)"}}, {"configKeyValues", VariantType::String, "", {"Semicolon separated key=value strings"}}}; @@ -77,5 +80,12 @@ WorkflowSpec defineDataProcessing(ConfigContext const& configcontext) specs.emplace_back(getTRDGainCalibSpec()); } + if (configcontext.options().get("t0")) { + if (enableRootInp) { + specs.emplace_back(o2::trd::getTRDPHReaderSpec()); + } + specs.emplace_back(getTRDT0FitSpec()); + } + return specs; } From 4dcda76eaa3627e60465992a1fda19b10cac32fd Mon Sep 17 00:00:00 2001 From: Luisa Bergmann Date: Tue, 25 Jul 2023 18:08:22 +0200 Subject: [PATCH 02/10] cleanup/documentation --- .../TRD/include/DataFormatsTRD/CalT0.h | 4 +- .../TRD/include/DataFormatsTRD/T0FitHistos.h | 2 - DataFormats/Detectors/TRD/src/T0FitHistos.cxx | 12 ------ .../include/TRDCalibration/T0Fit.h | 10 +---- Detectors/TRD/calibration/src/T0Fit.cxx | 37 ++++--------------- .../workflow/include/TRDWorkflow/T0FitSpec.h | 6 --- 6 files changed, 12 insertions(+), 59 deletions(-) diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h index 91f251b14b4d5..f1de8a4bf9d5c 100644 --- a/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h @@ -35,7 +35,9 @@ class CalT0 void setT0av(float t0) { mT0av = t0; } float getT0(int iDet) const { return mT0[iDet]; } + //getT0av() returns the average T0 obtained by fitting the data from all chambers combined float getT0av() const { return mT0av; } + //calcT0av() returns the average T0 from all individual chambers float calcT0av() const { if (mT0.size() == 0) { @@ -59,7 +61,7 @@ class CalT0 private: std::array mT0{}; ///< calibrated T0 per TRD chamber - float mT0av{-1}; ///< calibrated average T0 + float mT0av{-1}; ///< average T0 obtained from fitting the PH data from all chambers combined ClassDefNV(CalT0, 1); }; diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h index eff39bbbac438..821a16b63c304 100644 --- a/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h @@ -34,7 +34,6 @@ class T0FitHistos T0FitHistos() = default; T0FitHistos(const T0FitHistos&) = default; ~T0FitHistos() = default; - void init(); void reset(); auto getDetector(int index) const { return mDet[index]; } auto getTimeBin(int index) const { return mTB[index]; } @@ -50,7 +49,6 @@ class T0FitHistos std::vector mTB{}; std::vector mADC{}; size_t mNEntriesTot{0}; - bool mInitialized{false}; ClassDefNV(T0FitHistos, 1); }; diff --git a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx index 3b62525393fe8..d4c3a184a5ee8 100644 --- a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx +++ b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx @@ -19,13 +19,6 @@ using namespace o2::trd; using namespace o2::trd::constants; -void T0FitHistos::init() -{ - mDet.resize(0); - mTB.resize(0); - mADC.resize(0); - mInitialized = true; -} void T0FitHistos::reset() { @@ -37,11 +30,6 @@ void T0FitHistos::reset() void T0FitHistos::fill(const std::vector data) { - - if (!mInitialized) { - init(); - } - for (auto ph : data) { int det = ph.getDetector(); int tb = ph.getTimebin(); diff --git a/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h b/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h index ef08b3b9f80f6..081e59aa462fd 100644 --- a/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h +++ b/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h @@ -25,8 +25,6 @@ #include "TRDCalibration/CalibrationParams.h" #include "Rtypes.h" -#include "TH2.h" -#include "TH1.h" #include "TProfile.h" #include "TF1.h" #include "Fit/Fitter.h" @@ -77,15 +75,10 @@ class T0Fit final : public o2::calibration::TimeSlotCalibration mOutFile{nullptr}; ///< output file std::unique_ptr mOutTree{nullptr}; ///< output tree @@ -93,6 +86,7 @@ class T0Fit final : public o2::calibration::TimeSlotCalibration mParamsStart; ///< Starting parameters for fit std::unique_ptr mFuncErfLandau; ///< helper function to calculate the t0 value after the fitting procedure + float mDummyT0{-5}; ///< dummy value for t0, to be used if fit fails or not enough statistics std::array t0_chambers; ///< t0 values of the individual chambers float t0_average{-5}; ///< average t0 value across all chambers diff --git a/Detectors/TRD/calibration/src/T0Fit.cxx b/Detectors/TRD/calibration/src/T0Fit.cxx index 17030842bd981..5967cca123e7b 100644 --- a/Detectors/TRD/calibration/src/T0Fit.cxx +++ b/Detectors/TRD/calibration/src/T0Fit.cxx @@ -28,7 +28,6 @@ #include #include "TMath.h" -#include "TCanvas.h" #include #include @@ -61,7 +60,6 @@ double ErfLandauChi2Functor::operator()(const double* par) const double funcVal = par[0] * TMath::Erf(x[i]) + par[1] * TMath::Landau(x[i], par[2], par[3]); sum2 += TMath::Power(y[i] - funcVal, 2) / y[i]; - // sum2 += TMath::Power((y[i] - funcVal)/yErr[i],2); } return sum2; @@ -119,27 +117,6 @@ void T0Fit::initProcessing() mInitDone = true; } -void T0Fit::retrievePrev(o2::framework::ProcessingContext& pc) -{ - // We either get a pointer to a valid object from the last ~hour or to the default object - // which is always present. The first has precedence over the latter. - auto dataCalT0 = pc.inputs().get("calt0"); - std::string msg = "Default Object"; - // We check if the object we got is the default one by comparing it to the defaults. - for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { - if (dataCalT0->getT0(iDet) != constants::T0DEFAULT) { - msg = "Previous Object"; - break; - } - } - LOG(info) << "FitInstance: From CCDB retrieved " << msg; - - // Here we set each entry regardless if it is the default or not. - for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { - t0_chambers[iDet] = dataCalT0->getT0(iDet); - } - t0_average = dataCalT0->getT0av(); -} void T0Fit::finalizeSlot(Slot& slot) { @@ -174,15 +151,15 @@ void T0Fit::finalizeSlot(Slot& slot) mFuncErfLandau->SetParameters(fitResult.GetParams()[0], fitResult.GetParams()[1], fitResult.GetParams()[2], fitResult.GetParams()[3]); t0_average = mFuncErfLandau->GetMaximumX(); } else { - LOG(info) << "t0 fit for inclusive distribtion is not valid, set to -5"; - t0_average = -5; + LOG(info) << "t0 fit for inclusive distribtion is not valid, set to " << mDummyT0; + t0_average = mDummyT0; } // single chambers for (int iDet = 0; iDet < constants::MAXCHAMBER; ++iDet) { - if (adcProfDet[iDet]->GetEntries() < mMinEntriesChamber) { - LOG(info) << "not enough entries in chamber " << iDet << " for t0 fit, set to -5"; - t0_chambers[iDet] = -5; + if (adcProfDet[iDet]->GetEntries() < mParams.minEntriesChamberT0Fit) { + LOG(info) << "not enough entries in chamber " << iDet << " for t0 fit, set to " << mDummyT0; + t0_chambers[iDet] = mDummyT0; continue; } @@ -200,8 +177,8 @@ void T0Fit::finalizeSlot(Slot& slot) mFuncErfLandau->SetParameters(fitResult.GetParams()[0], fitResult.GetParams()[1], fitResult.GetParams()[2], fitResult.GetParams()[3]); t0_chambers[iDet] = mFuncErfLandau->GetMaximumX(); } else { - LOG(info) << "t0 fit for chamber " << iDet << " is not valid, set to -5"; - t0_chambers[iDet] = -5; + LOG(info) << "t0 fit for chamber " << iDet << " is not valid, set to " << mDummyT0; + t0_chambers[iDet] = mDummyT0; } } diff --git a/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h b/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h index decd9b8220f28..700383936adf1 100644 --- a/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h +++ b/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h @@ -31,8 +31,6 @@ #include "DetectorsBase/GRPGeomHelper.h" #include -#include "TH2.h" -#include "TH1.h" using namespace o2::framework; @@ -67,10 +65,6 @@ class T0FitDevice : public o2::framework::Task void run(o2::framework::ProcessingContext& pc) final { const auto& tinfo = pc.services().get(); - if (tinfo.globalRunNumberChanged) { // new run is starting - mRunStopRequested = false; - mFitInstance->retrievePrev(pc); // SOR initialization is performed here - } if (mRunStopRequested) { return; } From 6ab5139d93d70169b40581f3d026fe509ee68201 Mon Sep 17 00:00:00 2001 From: ALICE Action Bot Date: Tue, 25 Jul 2023 16:10:27 +0000 Subject: [PATCH 03/10] Please consider the following formatting changes --- DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h | 4 ++-- DataFormats/Detectors/TRD/src/T0FitHistos.cxx | 1 - Detectors/TRD/calibration/src/T0Fit.cxx | 1 - Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h | 1 - 4 files changed, 2 insertions(+), 5 deletions(-) diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h index f1de8a4bf9d5c..038c2a5684d5f 100644 --- a/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/CalT0.h @@ -35,9 +35,9 @@ class CalT0 void setT0av(float t0) { mT0av = t0; } float getT0(int iDet) const { return mT0[iDet]; } - //getT0av() returns the average T0 obtained by fitting the data from all chambers combined + // getT0av() returns the average T0 obtained by fitting the data from all chambers combined float getT0av() const { return mT0av; } - //calcT0av() returns the average T0 from all individual chambers + // calcT0av() returns the average T0 from all individual chambers float calcT0av() const { if (mT0.size() == 0) { diff --git a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx index d4c3a184a5ee8..0e2bfc8880f20 100644 --- a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx +++ b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx @@ -19,7 +19,6 @@ using namespace o2::trd; using namespace o2::trd::constants; - void T0FitHistos::reset() { mDet.resize(0); diff --git a/Detectors/TRD/calibration/src/T0Fit.cxx b/Detectors/TRD/calibration/src/T0Fit.cxx index 5967cca123e7b..a155eccc045a0 100644 --- a/Detectors/TRD/calibration/src/T0Fit.cxx +++ b/Detectors/TRD/calibration/src/T0Fit.cxx @@ -117,7 +117,6 @@ void T0Fit::initProcessing() mInitDone = true; } - void T0Fit::finalizeSlot(Slot& slot) { // do actual fits for the data provided in the given slot diff --git a/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h b/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h index 700383936adf1..3f4fc7a1e69fd 100644 --- a/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h +++ b/Detectors/TRD/workflow/include/TRDWorkflow/T0FitSpec.h @@ -31,7 +31,6 @@ #include "DetectorsBase/GRPGeomHelper.h" #include - using namespace o2::framework; namespace o2 From f0d1ac9149d052fa155842a0315f602aab5ff077 Mon Sep 17 00:00:00 2001 From: Luisa Bergmann Date: Thu, 27 Jul 2023 11:04:30 +0200 Subject: [PATCH 04/10] replaced entries-variable with direct instance of calibration parameters --- Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h b/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h index 081e59aa462fd..f7daea0365446 100644 --- a/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h +++ b/Detectors/TRD/calibration/include/TRDCalibration/T0Fit.h @@ -57,7 +57,7 @@ class T0Fit final : public o2::calibration::TimeSlotCalibrationgetNEntries() >= mMinEntriesTotal; } + bool hasEnoughData(const Slot& slot) const final { return slot.getContainer()->getNEntries() >= mParams.minEntriesTotalT0Fit; } void initOutput() final; void finalizeSlot(Slot& slot) final; Slot& emplaceNewSlot(bool front, TFType tStart, TFType tEnd) final; @@ -78,7 +78,6 @@ class T0Fit final : public o2::calibration::TimeSlotCalibration mOutFile{nullptr}; ///< output file std::unique_ptr mOutTree{nullptr}; ///< output tree From 8dbf74a8a961445a406d6d6cb69c6ce71bd011b3 Mon Sep 17 00:00:00 2001 From: Luisa Bergmann Date: Thu, 27 Jul 2023 17:55:13 +0200 Subject: [PATCH 05/10] removed reset() from calibration container --- .../Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h | 1 - DataFormats/Detectors/TRD/src/T0FitHistos.cxx | 7 ------- 2 files changed, 8 deletions(-) diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h index 821a16b63c304..5418f5277362c 100644 --- a/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h @@ -34,7 +34,6 @@ class T0FitHistos T0FitHistos() = default; T0FitHistos(const T0FitHistos&) = default; ~T0FitHistos() = default; - void reset(); auto getDetector(int index) const { return mDet[index]; } auto getTimeBin(int index) const { return mTB[index]; } auto getADC(int index) const { return mADC[index]; } diff --git a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx index 0e2bfc8880f20..6d285fc257604 100644 --- a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx +++ b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx @@ -19,13 +19,6 @@ using namespace o2::trd; using namespace o2::trd::constants; -void T0FitHistos::reset() -{ - mDet.resize(0); - mTB.resize(0); - mADC.resize(0); - mNEntriesTot = 0; -} void T0FitHistos::fill(const std::vector data) { From f937ba900807f67be9e116d5cb2d6f12b081a4fd Mon Sep 17 00:00:00 2001 From: ALICE Action Bot Date: Thu, 27 Jul 2023 15:56:37 +0000 Subject: [PATCH 06/10] Please consider the following formatting changes --- DataFormats/Detectors/TRD/src/T0FitHistos.cxx | 1 - 1 file changed, 1 deletion(-) diff --git a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx index 6d285fc257604..113d22a24aaa8 100644 --- a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx +++ b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx @@ -19,7 +19,6 @@ using namespace o2::trd; using namespace o2::trd::constants; - void T0FitHistos::fill(const std::vector data) { for (auto ph : data) { From 44942ef1257d6e4e0474153a5fe7190ae84538a2 Mon Sep 17 00:00:00 2001 From: luisabergmann Date: Thu, 3 Aug 2023 12:26:47 +0200 Subject: [PATCH 07/10] Update DataFormats/Detectors/TRD/src/T0FitHistos.cxx Co-authored-by: Felix Schlepper --- DataFormats/Detectors/TRD/src/T0FitHistos.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx index 113d22a24aaa8..99b61f5aff340 100644 --- a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx +++ b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx @@ -19,7 +19,7 @@ using namespace o2::trd; using namespace o2::trd::constants; -void T0FitHistos::fill(const std::vector data) +void T0FitHistos::fill(const std::vector& data) { for (auto ph : data) { int det = ph.getDetector(); From 2382abca68ed426405596fa230ae76bd5b58514f Mon Sep 17 00:00:00 2001 From: luisabergmann Date: Thu, 3 Aug 2023 12:27:32 +0200 Subject: [PATCH 08/10] Update DataFormats/Detectors/TRD/src/T0FitHistos.cxx Co-authored-by: Felix Schlepper --- DataFormats/Detectors/TRD/src/T0FitHistos.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx index 99b61f5aff340..623f61b761188 100644 --- a/DataFormats/Detectors/TRD/src/T0FitHistos.cxx +++ b/DataFormats/Detectors/TRD/src/T0FitHistos.cxx @@ -21,7 +21,7 @@ using namespace o2::trd::constants; void T0FitHistos::fill(const std::vector& data) { - for (auto ph : data) { + for (const auto& ph : data) { int det = ph.getDetector(); int tb = ph.getTimebin(); int adc = ph.getADC(); From e563022fb0a56e41a5ee48b4012482dcad5e9109 Mon Sep 17 00:00:00 2001 From: luisabergmann Date: Thu, 3 Aug 2023 12:29:30 +0200 Subject: [PATCH 09/10] Update Detectors/TRD/calibration/src/T0Fit.cxx Co-authored-by: Felix Schlepper --- Detectors/TRD/calibration/src/T0Fit.cxx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Detectors/TRD/calibration/src/T0Fit.cxx b/Detectors/TRD/calibration/src/T0Fit.cxx index a155eccc045a0..38fc131fb30cb 100644 --- a/Detectors/TRD/calibration/src/T0Fit.cxx +++ b/Detectors/TRD/calibration/src/T0Fit.cxx @@ -150,7 +150,7 @@ void T0Fit::finalizeSlot(Slot& slot) mFuncErfLandau->SetParameters(fitResult.GetParams()[0], fitResult.GetParams()[1], fitResult.GetParams()[2], fitResult.GetParams()[3]); t0_average = mFuncErfLandau->GetMaximumX(); } else { - LOG(info) << "t0 fit for inclusive distribtion is not valid, set to " << mDummyT0; + LOG(warn) << "t0 fit for inclusive distribtion is not valid, set to " << mDummyT0; t0_average = mDummyT0; } From 2629e2dddf6ef47a8cc231f0003065968a855334 Mon Sep 17 00:00:00 2001 From: Luisa Bergmann Date: Thu, 3 Aug 2023 12:56:48 +0200 Subject: [PATCH 10/10] Updated DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h --- DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h b/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h index 5418f5277362c..889c67dd24b84 100644 --- a/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h +++ b/DataFormats/Detectors/TRD/include/DataFormatsTRD/T0FitHistos.h @@ -39,7 +39,7 @@ class T0FitHistos auto getADC(int index) const { return mADC[index]; } auto getNEntries() const { return mNEntriesTot; } - void fill(const std::vector data); + void fill(const std::vector& data); void merge(const T0FitHistos* prev); void print();