event finder algorithm and workflow

Author: pillot

DataFormats/Detectors/MUON/MCH/include/DataFormatsMCH/Digit.h

@@ -36,6 +36,7 @@ class Digit
   ~Digit() = default;
 
   bool operator==(const Digit&) const;
+  bool operator!=(const Digit& other) const { return !(*this == other); }
   bool operator<(const Digit&) const;
 
   // time in bunch crossing units, relative to the beginning of the TimeFrame
@@ -65,10 +66,10 @@ class Digit
   uint32_t mADC;        /// Amplitude of signal
 
   ClassDefNV(Digit, 4);
-}; //class Digit
+}; // class Digit
 
 std::ostream& operator<<(std::ostream& os, const Digit& d);
 
-} //namespace mch
-} //namespace o2
+} // namespace mch
+} // namespace o2
 #endif // ALICEO2_MCH_DIGIT_H_

DataFormats/Detectors/MUON/MCH/include/DataFormatsMCH/ROFRecord.h

@@ -64,6 +64,7 @@ class ROFRecord
            mDataRef == other.mDataRef &&
            mBCWidth == other.mBCWidth;
   }
+  bool operator!=(const ROFRecord& other) const { return !(*this == other); }
   bool operator<(const ROFRecord& other) const
   {
     if (mBCData == other.mBCData) {

DataFormats/Detectors/MUON/MCH/src/Digit.cxx

@@ -66,7 +66,8 @@ bool Digit::operator==(const Digit& other) const
          mPadID == other.mPadID &&
          mADC == other.mADC &&
          mTFtime == other.mTFtime &&
-         mNofSamples == other.mNofSamples;
+         mNofSamples == other.mNofSamples &&
+         mIsSaturated == other.mIsSaturated;
 }
 
 bool Digit::operator<(const Digit& other) const

Detectors/MUON/MCH/CMakeLists.txt

@@ -12,6 +12,7 @@
 add_subdirectory(Base)
 add_subdirectory(Contour)
 add_subdirectory(Mapping)
+add_subdirectory(Triggering)
 add_subdirectory(TimeClustering)
 add_subdirectory(PreClustering)
 add_subdirectory(Geometry)

Detectors/MUON/MCH/README.md

@@ -13,6 +13,7 @@ This is a top page for the MCH detector documentation.
 \subpage refDetectorsMUONMCHRaw
 \subpage refDetectorsMUONMCHMapping
 \subpage refDetectorsMUONMCHTracking
+\subpage refDetectorsMUONMCHTriggering
 \subpage refDetectorsMUONMCHWorkflow
 \subpage refDetectorsMUONMCHGeometry
 \subpage refDetectorsMUONMCHConditions

Detectors/MUON/MCH/Triggering/CMakeLists.txt

@@ -0,0 +1,23 @@
+# 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.
+
+o2_add_library(MCHTriggering
+        SOURCES
+        src/EventFinder.cxx
+        src/EventFinderParam.cxx
+        PUBLIC_LINK_LIBRARIES O2::CommonUtils O2::DataFormatsMCH O2::DataFormatsMID)
+
+o2_target_root_dictionary(MCHTriggering
+                          HEADERS include/MCHTriggering/EventFinderParam.h)
+
+if(BUILD_TESTING)
+  add_subdirectory(test)
+endif()

Detectors/MUON/MCH/Triggering/README.md

@@ -0,0 +1,33 @@
+<!-- doxy
+\page refDetectorsMUONMCHTriggering Triggering
+/doxy -->
+
+# EventFinder.h(cxx)
+
+Implementation of the MCH event finder algorithm.
+
+## Input / Output
+
+It takes as input the lists of MCH ROFs, MCH digits, MCH digit MC labels (if any) and MID ROFs. It returns the list of digits associated to an event, their MC labels (if any) and the list of ROFs corresponding to the events.
+
+## Short description of the algorithm
+
+The event finding is done in 4 steps:
+1) Create an empty event for each MID ROF with a trigger window opened around the MID Interaction Record (IR). The size of the window is defined via EventFinderParam (see below).
+2) Associate each MCH ROF to the first compatible event, if any. An event is compatible if the trigger window overlaps with the MCH RO window.
+3) Merge overlapping events and remove empty ones. Two events overlap if the MCH RO window associated to one event overlaps with the trigger window of another.
+4) Copy the digits (and MC labels) associated to each event in the ouput lists, merging them if the same pad is fired multiple times within the same event, and produce the output ROFs pointing to them. The time window of each ROF is set to contain only the MID IR(s) associated to the event.
+
+# EventFinderParam.h(cxx)
+
+Definition of the trigger window [min, max[, in BC unit, opened around to MID IR to gather the digits associated with this interaction.
+
+The window is configurable via the command line or an INI or JSON configuration file (cf. [workflow documentation](../Workflow/README.md)).
+
+## Example of workflow
+
+`o2-mid-tracks-reader-workflow | o2-mch-reco-workflow --triggered`
+
+This takes as input the root files with MCH digits and MID tracks and ROFs and gives as ouput a root file with the MCH tracks. The event finder device is activated in the reconstruction workflow with the option `--triggered`.
+
+See the [workflow documentation](../Workflow/README.md) for more details about the event finder workflow.

Detectors/MUON/MCH/Triggering/include/MCHTriggering/EventFinder.h

@@ -0,0 +1,84 @@
+// 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 EventFinder.h
+/// \brief Definition of a class to group MCH digits based on MID information
+///
+/// \author Philippe Pillot, Subatech
+
+#ifndef ALICEO2_MCH_EVENTFINDER_H_
+#define ALICEO2_MCH_EVENTFINDER_H_
+
+#include <map>
+#include <unordered_map>
+#include <vector>
+
+#include <gsl/span>
+
+#include "SimulationDataFormat/MCCompLabel.h"
+#include "SimulationDataFormat/MCTruthContainer.h"
+#include "DataFormatsMCH/Digit.h"
+#include "DataFormatsMCH/ROFRecord.h"
+#include "DataFormatsMID/ROFRecord.h"
+
+namespace o2
+{
+namespace mch
+{
+
+/// Class to group MCH digits based on MID information
+class EventFinder
+{
+ public:
+  EventFinder() = default;
+  ~EventFinder() = default;
+
+  EventFinder(const EventFinder&) = delete;
+  EventFinder& operator=(const EventFinder&) = delete;
+  EventFinder(EventFinder&&) = delete;
+  EventFinder& operator=(EventFinder&&) = delete;
+
+  void run(const gsl::span<const mch::ROFRecord>& mchROFs, const gsl::span<const mch::Digit>& digits,
+           const dataformats::MCLabelContainer* labels, const gsl::span<const mid::ROFRecord>& midROFs);
+
+  /// get the output MCH ROFs
+  const std::vector<mch::ROFRecord>& getOutputROFs() const { return mROFs; }
+  /// get the output MCH digits
+  const std::vector<mch::Digit>& getOutputDigits() const { return mDigits; }
+  /// get the output MCH labels
+  const dataformats::MCLabelContainer& getOutputLabels() const { return mLabels; }
+
+ private:
+  /// internal event structure
+  struct Event {
+    /// contruct an empty event with trigger window = [bcMin, bcMax[
+    Event(int64_t bcMin, int64_t bcMax) : maxRORange(bcMin)
+    {
+      trgRange[0] = bcMin;
+      trgRange[1] = bcMax;
+    }
+
+    int64_t trgRange[2]{};       ///< BC range of the MID trigger window
+    int64_t maxRORange = 0;      ///< upper limit of the MCH RO window
+    std::vector<int> iMCHROFs{}; ///< list of associated MCH ROF indices
+  };
+
+  std::map<int64_t, Event> mEvents{};       ///< sorted list of events found
+  std::vector<mch::ROFRecord> mROFs{};      ///< list of output MCH ROFs
+  std::vector<mch::Digit> mDigits{};        ///< list of output MCH digits
+  dataformats::MCLabelContainer mLabels{};  ///< container of output MCH digit labels
+  std::unordered_map<int, int> mDigitLoc{}; ///< map the digit indices of a particular event
+};
+
+} // namespace mch
+} // namespace o2
+
+#endif // ALICEO2_MCH_EVENTFINDER_H_

Detectors/MUON/MCH/Triggering/include/MCHTriggering/EventFinderParam.h

@@ -0,0 +1,38 @@
+// 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 EventFinderParam.h
+/// \brief Configurable parameters to group MCH digits based on MID information
+/// \author Philippe Pillot, Subatech
+
+#ifndef ALICEO2_MCH_EVENTFINDERPARAM_H_
+#define ALICEO2_MCH_EVENTFINDERPARAM_H_
+
+#include "CommonUtils/ConfigurableParam.h"
+#include "CommonUtils/ConfigurableParamHelper.h"
+
+namespace o2
+{
+namespace mch
+{
+
+/// Configurable parameters to group MCH digits based on MID information
+struct EventFinderParam : public o2::conf::ConfigurableParamHelper<EventFinderParam> {
+
+  int triggerRange[2] = {0, 1}; ///< BC range of the trigger window around the MID IR
+
+  O2ParamDef(EventFinderParam, "MCHTriggering");
+};
+
+} // namespace mch
+} // namespace o2
+
+#endif // ALICEO2_MCH_EVENTFINDERPARAM_H_

Detectors/MUON/MCH/Triggering/src/EventFinder.cxx

@@ -0,0 +1,120 @@
+// 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 EventFinder.cxx
+/// \brief Implementation of a class to group MCH digits based on MID information
+///
+/// \author Philippe Pillot, Subatech
+
+#include "MCHTriggering/EventFinder.h"
+
+#include <algorithm>
+#include <iterator>
+#include <utility>
+
+#include "CommonDataFormat/InteractionRecord.h"
+#include "MCHTriggering/EventFinderParam.h"
+
+namespace o2
+{
+namespace mch
+{
+
+//_________________________________________________________________________________________________
+/// run the event finder algorithm
+void EventFinder::run(const gsl::span<const mch::ROFRecord>& mchROFs,
+                      const gsl::span<const mch::Digit>& digits,
+                      const dataformats::MCLabelContainer* labels,
+                      const gsl::span<const mid::ROFRecord>& midROFs)
+{
+  mEvents.clear();
+  mROFs.clear();
+  mDigits.clear();
+  mLabels.clear();
+
+  if (mchROFs.empty() || midROFs.empty()) {
+    return;
+  }
+
+  // create empty events ordered in increasing trigger time
+  const auto& param = EventFinderParam::Instance();
+  for (const auto& midROF : midROFs) {
+    if (midROF.nEntries == 0) {
+      continue;
+    }
+    auto midBC = midROF.interactionRecord.toLong();
+    mEvents.emplace(std::make_pair(midBC, Event(midBC + param.triggerRange[0], midBC + param.triggerRange[1])));
+  }
+
+  // associate each MCH ROF to the first compatible trigger, if any
+  for (int i = 0; i < mchROFs.size(); ++i) {
+    auto mchBC = mchROFs[i].getBCData().toLong();
+    auto itEvent = mEvents.lower_bound(mchBC - param.triggerRange[1] + 1);
+    if (itEvent != mEvents.end() && itEvent->second.trgRange[0] < mchBC + mchROFs[i].getBCWidth()) {
+      itEvent->second.iMCHROFs.push_back(i);
+      itEvent->second.maxRORange = std::max(itEvent->second.maxRORange, mchBC + mchROFs[i].getBCWidth());
+    }
+  }
+
+  // merge overlapping events (when a MCH ROF is compatible with multiple trigger) and cleanup
+  for (auto itEvent = mEvents.begin(); itEvent != mEvents.end();) {
+    if (itEvent->second.iMCHROFs.empty()) {
+      itEvent = mEvents.erase(itEvent);
+    } else {
+      auto itNextEvent = std::next(itEvent);
+      if (itNextEvent != mEvents.end() && itNextEvent->second.trgRange[0] < itEvent->second.maxRORange) {
+        itEvent->second.trgRange[1] = itNextEvent->second.trgRange[1];
+        itEvent->second.maxRORange = std::max(itEvent->second.maxRORange, itNextEvent->second.maxRORange);
+        itEvent->second.iMCHROFs.insert(itEvent->second.iMCHROFs.end(),
+                                        itNextEvent->second.iMCHROFs.begin(), itNextEvent->second.iMCHROFs.end());
+        mEvents.erase(itNextEvent);
+      } else {
+        itEvent = itNextEvent;
+      }
+    }
+  }
+
+  // merge digits associated to each event and produce the output ROFs pointing to them
+  // the BC range of each ROF is set to contain only the MID IR(s) associated to the event
+  InteractionRecord ir{};
+  for (const auto& event : mEvents) {
+    mDigitLoc.clear();
+    int digitOffset = mDigits.size();
+    for (auto iMCHROF : event.second.iMCHROFs) {
+      for (int iDigit = mchROFs[iMCHROF].getFirstIdx(); iDigit <= mchROFs[iMCHROF].getLastIdx(); ++iDigit) {
+        const auto& digit = digits[iDigit];
+        auto digitLoc = mDigitLoc.emplace(((digit.getDetID() << 16) | digit.getPadID()), mDigits.size());
+        if (digitLoc.second) {
+          mDigits.emplace_back(digit);
+          if (labels != nullptr) {
+            mLabels.addElements(digitLoc.first->second, labels->getLabels(iDigit));
+          }
+        } else {
+          auto& digit0 = mDigits[digitLoc.first->second];
+          digit0.setADC(digit0.getADC() + digit.getADC());
+          digit0.setNofSamples(digit0.getNofSamples() + digit.getNofSamples());
+          digit0.setSaturated(digit0.isSaturated() || digit.isSaturated());
+          if (labels != nullptr) {
+            for (const auto& label : labels->getLabels(iDigit)) {
+              mLabels.addElementRandomAccess(digitLoc.first->second, label);
+            }
+          }
+        }
+      }
+    }
+    ir.setFromLong(event.first);
+    mROFs.emplace_back(ir, digitOffset, mDigits.size() - digitOffset,
+                       event.second.trgRange[1] - param.triggerRange[1] - event.first + 1);
+  }
+}
+
+} // namespace mch
+} // namespace o2