[EMCAL-788] EMCAL cell selection, on-the-fly calibration

- Seletion of good EMCAL cells
  - Require FEE readout (HG and LG cells)
  - Time cut around trigger peak
  - Min. cell energy
- On-the-fly recalibration
  - Optional
  - Only in synchronous workflow

Author: mfasDa

EventVisualisation/Workflow/CMakeLists.txt

@@ -40,6 +40,8 @@ if(ALIGPU_BUILD_TYPE STREQUAL "O2"
                 O2::TOFBase
                 O2::PHOSBase
                 O2::EMCALBase
+                O2::EMCALCalib
+                O2::EMCALWorkflow
                 O2::MIDBase
                 O2::TOFWorkflowIO
                 O2::TPCReconstruction
@@ -74,6 +76,7 @@ if(ALIGPU_BUILD_TYPE STREQUAL "O2"
                 O2::TOFBase
                 O2::PHOSBase
                 O2::EMCALBase
+                O2::EMCALCalib
                 O2::MIDBase
                 O2::TOFWorkflowIO
                 O2::TPCReconstruction

EventVisualisation/Workflow/include/EveWorkflow/EveWorkflowHelper.h

@@ -55,7 +55,8 @@ class Geometry;
 namespace o2::emcal
 {
 class Geometry;
-}
+class CellRecalibrator;
+} // namespace o2::emcal
 
 namespace o2::event_visualisation
 {
@@ -207,6 +208,9 @@ class EveWorkflowHelper
   const o2::globaltracking::RecoContainer* mRecoCont = nullptr;
   const o2::globaltracking::RecoContainer* getRecoContainer() const { return mRecoCont; }
   void setRecoContainer(const o2::globaltracking::RecoContainer* rc) { mRecoCont = rc; }
+  void setEMCALCellRecalibrator(o2::emcal::CellRecalibrator* calibrator) { mEMCALCalib = calibrator; }
+  void setMaxEMCALCellTime(float maxtime) { mEMCALMaxCellTime = maxtime; }
+  void setMinEMCALCellEnergy(float minenergy) { mEMCALMinCellEnergy = minenergy; }
   TracksSet mTrackSet;
   o2::event_visualisation::VisualisationEvent mEvent;
   std::unordered_map<GID, std::size_t> mTotalDataTypes;
@@ -221,11 +225,14 @@ class EveWorkflowHelper
   o2::its::GeometryTGeo* mITSGeom;
   o2::phos::Geometry* mPHOSGeom;
   o2::emcal::Geometry* mEMCALGeom;
+  o2::emcal::CellRecalibrator* mEMCALCalib = nullptr;
 
   float mMUS2TPCTimeBins = 5.0098627;
   float mITSROFrameLengthMUS = 0; ///< ITS RO frame in mus
   float mMFTROFrameLengthMUS = 0; ///< MFT RO frame in mus
   float mTPCBin2MUS = 0;
+  float mEMCALMaxCellTime = 100.;  ///< EMCAL cell time cut (in ns)
+  float mEMCALMinCellEnergy = 0.3; ///< EMCAL cell energy cut (in GeV)
   static int BCDiffErrCount;
   const o2::vertexing::PVertexerParams* mPVParams = nullptr;
 };

EventVisualisation/Workflow/include/EveWorkflow/O2DPLDisplay.h

@@ -19,6 +19,8 @@
 #include "ReconstructionDataFormats/GlobalTrackID.h"
 #include "DataFormatsGlobalTracking/RecoContainer.h"
 #include "DetectorsBase/GRPGeomHelper.h"
+#include "EMCALCalib/CellRecalibrator.h"
+#include "EMCALWorkflow/CalibLoader.h"
 #include "EveWorkflow/DetectorData.h"
 #include "Framework/Task.h"
 #include <memory>
@@ -54,12 +56,13 @@ class O2DPLDisplaySpec : public o2::framework::Task
                    o2::dataformats::GlobalTrackID::mask_t clMask,
                    std::shared_ptr<o2::globaltracking::DataRequest> dataRequest,
                    std::shared_ptr<o2::base::GRPGeomRequest> gr,
+                   std::shared_ptr<o2::emcal::CalibLoader> emcCalibLoader,
                    const std::string& jsonPath, const std::string& ext,
                    std::chrono::milliseconds timeInterval, int numberOfFiles, int numberOfTracks,
                    bool eveHostNameMatch, int minITSTracks, int minTracks, bool filterITSROF, bool filterTime,
                    const EveWorkflowHelper::Bracket& timeBracket, bool removeTPCEta,
-                   const EveWorkflowHelper::Bracket& etaBracket, bool trackSorting, int onlyNthEvent, bool primaryVertex, int maxPrimaryVertices, bool primaryVertexTriggers, float primaryVertexMinZ, float primaryVertexMaxZ, float primaryVertexMinX, float primaryVertexMaxX, float primaryVertexMinY, float primaryVertexMaxY)
-    : mDisableWrite(disableWrite), mUseMC(useMC), mTrkMask(trkMask), mClMask(clMask), mDataRequest(dataRequest), mGGCCDBRequest(gr), mJsonPath(jsonPath), mExt(ext), mTimeInterval(timeInterval), mNumberOfFiles(numberOfFiles), mNumberOfTracks(numberOfTracks), mEveHostNameMatch(eveHostNameMatch), mMinITSTracks(minITSTracks), mMinTracks(minTracks), mFilterITSROF(filterITSROF), mFilterTime(filterTime), mTimeBracket(timeBracket), mRemoveTPCEta(removeTPCEta), mEtaBracket(etaBracket), mTrackSorting(trackSorting), mOnlyNthEvent(onlyNthEvent), mPrimaryVertexMode(primaryVertex), mMaxPrimaryVertices(maxPrimaryVertices), mPrimaryVertexTriggers(primaryVertexTriggers), mPrimaryVertexMinZ(primaryVertexMinZ), mPrimaryVertexMaxZ(primaryVertexMaxZ), mPrimaryVertexMinX(primaryVertexMinX), mPrimaryVertexMaxX(primaryVertexMaxX), mPrimaryVertexMinY(primaryVertexMinY), mPrimaryVertexMaxY(primaryVertexMaxY), mRunType(o2::parameters::GRPECS::NONE)
+                   const EveWorkflowHelper::Bracket& etaBracket, bool trackSorting, int onlyNthEvent, bool primaryVertex, int maxPrimaryVertices, bool primaryVertexTriggers, float primaryVertexMinZ, float primaryVertexMaxZ, float primaryVertexMinX, float primaryVertexMaxX, float primaryVertexMinY, float primaryVertexMaxY, float maxEMCALCellTime, float minEMCALCellEnergy)
+    : mDisableWrite(disableWrite), mUseMC(useMC), mTrkMask(trkMask), mClMask(clMask), mDataRequest(dataRequest), mGGCCDBRequest(gr), mEMCALCalibLoader(emcCalibLoader), mJsonPath(jsonPath), mExt(ext), mTimeInterval(timeInterval), mNumberOfFiles(numberOfFiles), mNumberOfTracks(numberOfTracks), mEveHostNameMatch(eveHostNameMatch), mMinITSTracks(minITSTracks), mMinTracks(minTracks), mFilterITSROF(filterITSROF), mFilterTime(filterTime), mTimeBracket(timeBracket), mRemoveTPCEta(removeTPCEta), mEtaBracket(etaBracket), mTrackSorting(trackSorting), mOnlyNthEvent(onlyNthEvent), mPrimaryVertexMode(primaryVertex), mMaxPrimaryVertices(maxPrimaryVertices), mPrimaryVertexTriggers(primaryVertexTriggers), mPrimaryVertexMinZ(primaryVertexMinZ), mPrimaryVertexMaxZ(primaryVertexMaxZ), mPrimaryVertexMinX(primaryVertexMinX), mPrimaryVertexMaxX(primaryVertexMaxX), mPrimaryVertexMinY(primaryVertexMinY), mPrimaryVertexMaxY(primaryVertexMaxY), mEMCALMaxCellTime(maxEMCALCellTime), mEMCALMinCellEnergy(minEMCALCellEnergy), mRunType(o2::parameters::GRPECS::NONE)
 
   {
     this->mTimeStamp = std::chrono::high_resolution_clock::now() - timeInterval; // first run meets condition
@@ -99,6 +102,8 @@ class O2DPLDisplaySpec : public o2::framework::Task
   float mPrimaryVertexMaxX;                // maximum x position of the primary vertex
   float mPrimaryVertexMinY;                // minimum y position of the primary vertex
   float mPrimaryVertexMaxY;                // maximum y position of the primary vertex
+  float mEMCALMaxCellTime;                 // max abs EMCAL cell time (in ns)
+  float mEMCALMinCellEnergy;               // min EMCAL cell energy (in GeV)
   int mEventCounter = 0;
   std::chrono::time_point<std::chrono::high_resolution_clock> mTimeStamp;
 
@@ -108,6 +113,8 @@ class O2DPLDisplaySpec : public o2::framework::Task
   o2::parameters::GRPECS::RunType mRunType;
   std::shared_ptr<o2::globaltracking::DataRequest> mDataRequest;
   std::shared_ptr<o2::base::GRPGeomRequest> mGGCCDBRequest;
+  std::shared_ptr<o2::emcal::CalibLoader> mEMCALCalibLoader;
+  std::unique_ptr<o2::emcal::CellRecalibrator> mEMCALCalibrator;
 };
 
 } // namespace o2::event_visualisation

EventVisualisation/Workflow/src/EveWorkflowHelper.cxx

@@ -35,6 +35,7 @@
 #include "ITSBase/GeometryTGeo.h"
 #include "PHOSBase/Geometry.h"
 #include "EMCALBase/Geometry.h"
+#include "EMCALCalib/CellRecalibrator.h"
 #include <TGeoBBox.h>
 #include <tuple>
 #include <gsl/span>
@@ -644,6 +645,30 @@ void EveWorkflowHelper::drawEMC(GID gid)
   const gsl::span cellsForTrigger(cells.data() + trig.getFirstEntry(), trig.getNumberOfObjects());
 
   for (const auto& cell : cellsForTrigger) {
+    if (!(cell.getType() == o2::emcal::ChannelType_t::HIGH_GAIN || cell.getType() == o2::emcal::ChannelType_t::LOW_GAIN)) {
+      // Select FEE cells (excluding LEDMON or TRU cells)
+      continue;
+    }
+    auto cellTime = cell.getTimeStamp();
+    auto cellEnergy = cell.getEnergy();
+    if (mEMCALCalib) {
+      // try to recalibrate cell in case calibration is available (bad channel removal, energy calibration, time calibration)
+      auto calibCell = mEMCALCalib->getCalibratedCell(cell);
+      if (!calibCell) {
+        // cell was rejected by bad channel calib
+        continue;
+      }
+      cellTime = calibCell->getTimeStamp();
+      cellEnergy = calibCell->getEnergy();
+    }
+    if (std::abs(cellTime) > mEMCALMaxCellTime) {
+      // cell rejected by time cut (pileup or noise)
+      continue;
+    }
+    if (cellEnergy < mEMCALMinCellEnergy) {
+      // cell rejected by energy cut (rejection of soft cells)
+      continue;
+    }
     const auto id = cell.getTower();
     // Point3D with x,y,z coordinates of cell with absId inside SM
     const auto relPosCell = this->mEMCALGeom->RelPosCellInSModule(id);
@@ -657,7 +682,7 @@ void EveWorkflowHelper::drawEMC(GID gid)
     TVector3 vPos(gPos.data());
 
     auto vCalo = mEvent.addCalo({.time = static_cast<float>(time),
-                                 .energy = cell.getEnergy(),
+                                 .energy = cellEnergy,
                                  .phi = (float)vPos.Phi(),
                                  .eta = (float)vPos.Eta(),
                                  .PID = 0,

EventVisualisation/Workflow/src/O2DPLDisplay.cxx

@@ -26,6 +26,8 @@
 #include "TOFBase/Geo.h"
 #include "TPCFastTransform.h"
 #include "TRDBase/Geometry.h"
+#include "EMCALCalib/CellRecalibrator.h"
+#include "EMCALWorkflow/CalibLoader.h"
 #include "GlobalTrackingWorkflowHelpers/InputHelper.h"
 #include "ReconstructionDataFormats/GlobalTrackID.h"
 #include "ReconstructionDataFormats/PrimaryVertex.h"
@@ -76,6 +78,9 @@ void customize(std::vector<ConfigParamSpec>& workflowOptions)
     {"primary-vertex-mode", VariantType::Bool, false, {"produce jsons with individual primary vertices, not total time frame data"}},
     {"max-primary-vertices", VariantType::Int, 5, {"maximum number of primary vertices to draw per time frame"}},
     {"primary-vertex-triggers", VariantType::Bool, false, {"instead of drawing vertices with tracks (and maybe calorimeter triggers), draw vertices with calorimeter triggers (and maybe tracks)"}},
+    {"no-calibrate-emcal", VariantType::Bool, false, {"Do not apply on-the-fly EMCAL calibration"}},
+    {"emcal-max-celltime", VariantType::Float, 100.f, {"Max. EMCAL cell time (in ns)"}},
+    {"emcal-min-cellenergy", VariantType::Float, 0.3f, {"Min. EMCAL cell energy (in GeV)"}},
     {"primary-vertex-min-z", VariantType::Float, -o2::constants::math::VeryBig, {"minimum z position for primary vertex"}},
     {"primary-vertex-max-z", VariantType::Float, o2::constants::math::VeryBig, {"maximum z position for primary vertex"}},
     {"primary-vertex-min-x", VariantType::Float, -o2::constants::math::VeryBig, {"minimum x position for primary vertex"}},
@@ -93,6 +98,9 @@ void O2DPLDisplaySpec::init(InitContext& ic)
   LOGF(info, "------------------------    O2DPLDisplay::init version ", o2_eve_version, "    ------------------------------------");
   mData.mConfig.configProcessing.runMC = mUseMC;
   o2::base::GRPGeomHelper::instance().setRequest(mGGCCDBRequest);
+  if (mEMCALCalibLoader) {
+    mEMCALCalibrator = std::make_unique<o2::emcal::CellRecalibrator>();
+  }
 }
 
 void O2DPLDisplaySpec::run(ProcessingContext& pc)
@@ -114,6 +122,19 @@ void O2DPLDisplaySpec::run(ProcessingContext& pc)
   o2::globaltracking::RecoContainer recoCont;
   recoCont.collectData(pc, *mDataRequest);
   updateTimeDependentParams(pc); // Make sure that this is called after the RecoContainer collect data, since some condition objects are fetched there
+  if (mEMCALCalibLoader) {
+    mEMCALCalibLoader->checkUpdates(pc);
+    if (mEMCALCalibLoader->hasUpdateBadChannelMap()) {
+      mEMCALCalibrator->setBadChannelMap(mEMCALCalibLoader->getBadChannelMap());
+    }
+    if (mEMCALCalibLoader->hasUpdateTimeCalib()) {
+      mEMCALCalibrator->setTimeCalibration(mEMCALCalibLoader->getTimeCalibration());
+    }
+    if (mEMCALCalibLoader->hasUpdateGainCalib()) {
+      mEMCALCalibrator->setGainCalibration(mEMCALCalibLoader->getGainCalibration());
+    }
+  }
+
   EveWorkflowHelper::FilterSet enabledFilters;
 
   enabledFilters.set(EveWorkflowHelper::Filter::ITSROF, this->mFilterITSROF);
@@ -122,6 +143,12 @@ void O2DPLDisplaySpec::run(ProcessingContext& pc)
   enabledFilters.set(EveWorkflowHelper::Filter::TotalNTracks, this->mNumberOfTracks != -1);
   EveWorkflowHelper helper(enabledFilters, this->mNumberOfTracks, this->mTimeBracket, this->mEtaBracket, this->mPrimaryVertexMode);
   helper.setRecoContainer(&recoCont);
+  if (mEMCALCalibrator) {
+    helper.setEMCALCellRecalibrator(mEMCALCalibrator.get());
+  }
+  helper.setMaxEMCALCellTime(mEMCALMaxCellTime);
+  helper.setMinEMCALCellEnergy(mEMCALMinCellEnergy);
+
   helper.setITSROFs();
   helper.selectTracks(&(mData.mConfig.configCalib), mClMask, mTrkMask, mTrkMask);
   helper.selectTowers();
@@ -240,6 +267,9 @@ void O2DPLDisplaySpec::finaliseCCDB(ConcreteDataMatcher& matcher, void* obj)
   if (o2::base::GRPGeomHelper::instance().finaliseCCDB(matcher, obj)) {
     return;
   }
+  if (mEMCALCalibLoader && mEMCALCalibLoader->finalizeCCDB(matcher, obj)) {
+    return;
+  }
   if (matcher == ConcreteDataMatcher("ITS", "CLUSDICT", 0)) {
     LOGF(info, "ITS cluster dictionary updated");
     mData.setITSDict((const o2::itsmft::TopologyDictionary*)obj);
@@ -372,6 +402,8 @@ WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
   auto primaryVertexMaxX = cfgc.options().get<float>("primary-vertex-max-x");
   auto primaryVertexMinY = cfgc.options().get<float>("primary-vertex-min-y");
   auto primaryVertexMaxY = cfgc.options().get<float>("primary-vertex-max-y");
+  auto maxEMCALCellTime = cfgc.options().get<float>("emcal-max-celltime");
+  auto minEMCALCellEnergy = cfgc.options().get<float>("emcal-min-cellenergy");
 
   if (numberOfTracks == -1) {
     tracksSorting = false; // do not sort if all tracks are allowed
@@ -385,11 +417,20 @@ WorkflowSpec defineDataProcessing(ConfigContext const& cfgc)
                                                               dataRequest->inputs,
                                                               true); // query only once all objects except mag.field
 
+  std::shared_ptr<o2::emcal::CalibLoader> emcalCalibLoader;
+  if (!cfgc.options().get<bool>("no-calibrate-emcal")) {
+    emcalCalibLoader = std::make_shared<o2::emcal::CalibLoader>();
+    emcalCalibLoader->enableTimeCalib(true);
+    emcalCalibLoader->enableBadChannelMap(true);
+    emcalCalibLoader->enableGainCalib(true);
+    emcalCalibLoader->defineInputSpecs(dataRequest->inputs);
+  }
+
   specs.emplace_back(DataProcessorSpec{
     "o2-eve-export",
     dataRequest->inputs,
     {},
-    AlgorithmSpec{adaptFromTask<O2DPLDisplaySpec>(disableWrite, useMC, srcTrk, srcCl, dataRequest, ggRequest, jsonFolder, ext, timeInterval, numberOfFiles, numberOfTracks, eveHostNameMatch, minITSTracks, minTracks, filterITSROF, filterTime, timeBracket, removeTPCEta, etaBracket, tracksSorting, onlyNthEvent, primaryVertexMode, maxPrimaryVertices, primaryVertexTriggers, primaryVertexMinZ, primaryVertexMaxZ, primaryVertexMinX, primaryVertexMaxX, primaryVertexMinY, primaryVertexMaxY)}});
+    AlgorithmSpec{adaptFromTask<O2DPLDisplaySpec>(disableWrite, useMC, srcTrk, srcCl, dataRequest, ggRequest, emcalCalibLoader, jsonFolder, ext, timeInterval, numberOfFiles, numberOfTracks, eveHostNameMatch, minITSTracks, minTracks, filterITSROF, filterTime, timeBracket, removeTPCEta, etaBracket, tracksSorting, onlyNthEvent, primaryVertexMode, maxPrimaryVertices, primaryVertexTriggers, primaryVertexMinZ, primaryVertexMaxZ, primaryVertexMinX, primaryVertexMaxX, primaryVertexMinY, primaryVertexMaxY, maxEMCALCellTime, minEMCALCellEnergy)}});
 
   // configure dpl timer to inject correct firstTForbit: start from the 1st orbit of TF containing 1st sampled orbit
   o2::raw::HBFUtilsInitializer hbfIni(cfgc, specs);