[MCH] modify (pre)clustering workflows to process several ROF per TF (#5463)

Author: pillot

Detectors/MUON/MCH/Base/include/MCHBase/ClusterBlock.h

@@ -17,6 +17,7 @@
 #define ALICEO2_MCH_CLUSTERBLOCK_H_
 
 #include <iostream>
+#include <stdexcept>
 
 namespace o2
 {
@@ -51,7 +52,10 @@ struct ClusterStruct {
   /// Build the unique ID of the cluster from the chamber ID, detection element ID and cluster index
   static uint32_t buildUniqueId(int chamberId, int deId, int clusterIndex)
   {
-    return (((chamberId & 0xF) << 28) | ((deId & 0x7FF) << 17) | (clusterIndex & 0x1FFFF));
+    if ((clusterIndex & 0x1FFFF) != clusterIndex) {
+      throw std::runtime_error("invalid cluster index. Cannot build unique ID");
+    }
+    return (((chamberId & 0xF) << 28) | ((deId & 0x7FF) << 17) | clusterIndex);
   }
 };
 

Detectors/MUON/MCH/Base/include/MCHBase/PreCluster.h

@@ -29,11 +29,11 @@ namespace mch
 
 /// precluster minimal structure
 struct PreCluster {
-  uint16_t firstDigit; ///< index of first associated digit in the ordered vector of digits
-  uint16_t nDigits;    ///< number of digits attached to this precluster
+  uint32_t firstDigit; ///< index of first associated digit in the ordered vector of digits
+  uint32_t nDigits;    ///< number of digits attached to this precluster
 
   /// return the index of last associated digit in the ordered vector of digits
-  uint16_t lastDigit() const { return firstDigit + nDigits - 1; }
+  uint32_t lastDigit() const { return firstDigit + nDigits - 1; }
 
   void print(std::ostream& stream, gsl::span<const Digit> digits) const;
 };

Detectors/MUON/MCH/PreClustering/src/PreClusterFinder.cxx

@@ -178,8 +178,8 @@ void PreClusterFinder::getPreClusters(std::vector<o2::mch::PreCluster>& preClust
         }
 
         // add this precluster
-        uint16_t firstDigit = digits.size();
-        uint16_t nDigits = cluster->lastPad - cluster->firstPad + 1;
+        uint32_t firstDigit = digits.size();
+        uint32_t nDigits = cluster->lastPad - cluster->firstPad + 1;
         preClusters.push_back({firstDigit, nDigits});
 
         // add the digits of this precluster

Detectors/MUON/MCH/Workflow/CMakeLists.txt

@@ -14,7 +14,7 @@ o2_add_library(MCHWorkflow
                PUBLIC_LINK_LIBRARIES O2::Framework O2::DPLUtils O2::MCHRawDecoder Boost::program_options
                                      O2::MCHRawImplHelpers RapidJSON::RapidJSON O2::MCHMappingInterface
                                      O2::MCHPreClustering O2::MCHMappingImpl4 O2::MCHRawElecMap O2::MCHBase
-                                     O2::MCHClustering)
+                                     O2::DataFormatsMCH O2::MCHClustering)
 
 o2_add_executable(
         cru-page-reader-workflow
@@ -56,13 +56,13 @@ o2_add_executable(
         digits-reader-workflow
         SOURCES src/DigitSamplerSpec.cxx src/digits-reader-workflow.cxx
         COMPONENT_NAME mch
-        PUBLIC_LINK_LIBRARIES O2::MCHWorkflow)
+        PUBLIC_LINK_LIBRARIES O2::Framework O2::DataFormatsMCH O2::MCHBase O2::MCHMappingImpl3)
 
 o2_add_executable(
         preclusters-sink-workflow
         SOURCES src/PreClusterSinkSpec.cxx src/preclusters-sink-workflow.cxx
         COMPONENT_NAME mch
-        PUBLIC_LINK_LIBRARIES O2::MCHWorkflow)
+        PUBLIC_LINK_LIBRARIES O2::Framework O2::DataFormatsMCH O2::MCHBase O2::MCHMappingImpl3)
 
 o2_add_executable(
         preclusters-to-clusters-original-workflow
@@ -74,7 +74,7 @@ o2_add_executable(
         clusters-sink-workflow
         SOURCES src/ClusterSinkSpec.cxx src/clusters-sink-workflow.cxx
         COMPONENT_NAME mch
-        PUBLIC_LINK_LIBRARIES O2::MCHWorkflow)
+        PUBLIC_LINK_LIBRARIES O2::Framework O2::DataFormatsMCH O2::MCHBase O2::MCHMappingImpl3)
 
 o2_add_executable(
         clusters-sampler-workflow

Detectors/MUON/MCH/Workflow/README.md

@@ -19,9 +19,10 @@
 #include <fstream>
 #include <chrono>
 #include <vector>
-
 #include <stdexcept>
 
+#include <gsl/span>
+
 #include "Framework/CallbackService.h"
 #include "Framework/ConfigParamRegistry.h"
 #include "Framework/ControlService.h"
@@ -31,6 +32,7 @@
 #include "Framework/Task.h"
 #include "Framework/Logger.h"
 
+#include "DataFormatsMCH/ROFRecord.h"
 #include "MCHBase/Digit.h"
 #include "MCHBase/PreCluster.h"
 #include "MCHBase/ClusterBlock.h"
@@ -65,27 +67,59 @@ class ClusterFinderOriginalTask
   //_________________________________________________________________________________________________
   void run(framework::ProcessingContext& pc)
   {
-    /// read the preclusters and associated digits, clusterize and send the clusters
+    /// read the preclusters and associated digits, clusterize and send the clusters for all events in the TF
 
     // get the input preclusters and associated digits
+    auto preClusterROFs = pc.inputs().get<gsl::span<ROFRecord>>("preclusterrofs");
     auto preClusters = pc.inputs().get<gsl::span<PreCluster>>("preclusters");
     auto digits = pc.inputs().get<gsl::span<Digit>>("digits");
 
-    // clusterize every preclusters
-    auto tStart = std::chrono::high_resolution_clock::now();
-    mClusterFinder.reset();
-    for (const auto& preCluster : preClusters) {
-      mClusterFinder.findClusters(digits.subspan(preCluster.firstDigit, preCluster.nDigits));
-    }
-    auto tEnd = std::chrono::high_resolution_clock::now();
-    mTimeClusterFinder += tEnd - tStart;
+    //LOG(INFO) << "received time frame with " << preClusterROFs.size() << " interactions";
+
+    // create the output messages for clusters and attached digits
+    auto& clusterROFs = pc.outputs().make<std::vector<ROFRecord>>(OutputRef{"clusterrofs"});
+    auto& clusters = pc.outputs().make<std::vector<ClusterStruct>>(OutputRef{"clusters"});
+    auto& usedDigits = pc.outputs().make<std::vector<Digit>>(OutputRef{"clusterdigits"});
 
-    // send the output messages
-    pc.outputs().snapshot(Output{"MCH", "CLUSTERS", 0, Lifetime::Timeframe}, mClusterFinder.getClusters());
-    pc.outputs().snapshot(Output{"MCH", "CLUSTERDIGITS", 0, Lifetime::Timeframe}, mClusterFinder.getUsedDigits());
+    clusterROFs.reserve(preClusterROFs.size());
+    for (const auto& preClusterROF : preClusterROFs) {
+
+      //LOG(INFO) << "processing interaction: " << preClusterROF.getBCData() << "...";
+
+      // clusterize every preclusters
+      auto tStart = std::chrono::high_resolution_clock::now();
+      mClusterFinder.reset();
+      for (const auto& preCluster : preClusters.subspan(preClusterROF.getFirstIdx(), preClusterROF.getNEntries())) {
+        mClusterFinder.findClusters(digits.subspan(preCluster.firstDigit, preCluster.nDigits));
+      }
+      auto tEnd = std::chrono::high_resolution_clock::now();
+      mTimeClusterFinder += tEnd - tStart;
+
+      // fill the ouput messages
+      clusterROFs.emplace_back(preClusterROF.getBCData(), clusters.size(), mClusterFinder.getClusters().size());
+      writeClusters(clusters, usedDigits);
+    }
   }
 
  private:
+  //_________________________________________________________________________________________________
+  void writeClusters(std::vector<ClusterStruct, o2::pmr::polymorphic_allocator<ClusterStruct>>& clusters,
+                     std::vector<Digit, o2::pmr::polymorphic_allocator<Digit>>& usedDigits) const
+  {
+    /// fill the output messages with clusters and attached digits of the current event
+    /// modify the references to the attached digits according to their position in the global vector
+
+    auto clusterOffset = clusters.size();
+    clusters.insert(clusters.end(), mClusterFinder.getClusters().begin(), mClusterFinder.getClusters().end());
+
+    auto digitOffset = usedDigits.size();
+    usedDigits.insert(usedDigits.end(), mClusterFinder.getUsedDigits().begin(), mClusterFinder.getUsedDigits().end());
+
+    for (auto itCluster = clusters.begin() + clusterOffset; itCluster < clusters.end(); ++itCluster) {
+      itCluster->firstDigit += digitOffset;
+    }
+  }
+
   ClusterFinderOriginal mClusterFinder{};             ///< clusterizer
   std::chrono::duration<double> mTimeClusterFinder{}; ///< timer
 };
@@ -95,10 +129,12 @@ o2::framework::DataProcessorSpec getClusterFinderOriginalSpec()
 {
   return DataProcessorSpec{
     "ClusterFinderOriginal",
-    Inputs{InputSpec{"preclusters", "MCH", "PRECLUSTERS", 0, Lifetime::Timeframe},
+    Inputs{InputSpec{"preclusterrofs", "MCH", "PRECLUSTERROFS", 0, Lifetime::Timeframe},
+           InputSpec{"preclusters", "MCH", "PRECLUSTERS", 0, Lifetime::Timeframe},
            InputSpec{"digits", "MCH", "PRECLUSTERDIGITS", 0, Lifetime::Timeframe}},
-    Outputs{OutputSpec{"MCH", "CLUSTERS", 0, Lifetime::Timeframe},
-            OutputSpec{"MCH", "CLUSTERDIGITS", 0, Lifetime::Timeframe}},
+    Outputs{OutputSpec{{"clusterrofs"}, "MCH", "CLUSTERROFS", 0, Lifetime::Timeframe},
+            OutputSpec{{"clusters"}, "MCH", "CLUSTERS", 0, Lifetime::Timeframe},
+            OutputSpec{{"clusterdigits"}, "MCH", "CLUSTERDIGITS", 0, Lifetime::Timeframe}},
     AlgorithmSpec{adaptFromTask<ClusterFinderOriginalTask>()},
     Options{}};
 }

Detectors/MUON/MCH/Workflow/src/ClusterFinderOriginalSpec.cxx

@@ -17,11 +17,12 @@
 
 #include <iostream>
 #include <fstream>
-
 #include <array>
 #include <stdexcept>
 #include <vector>
 
+#include <gsl/span>
+
 #include "Framework/CallbackService.h"
 #include "Framework/ConfigParamRegistry.h"
 #include "Framework/ControlService.h"
@@ -30,6 +31,7 @@
 #include "Framework/Task.h"
 #include "Framework/Logger.h"
 
+#include "DataFormatsMCH/ROFRecord.h"
 #include "MCHBase/Digit.h"
 #include "MCHBase/ClusterBlock.h"
 #include "MCHMappingInterface/Segmentation.h"
@@ -72,42 +74,89 @@ class ClusterSinkTask
   //_________________________________________________________________________________________________
   void run(framework::ProcessingContext& pc)
   {
-    /// dump the clusters with associated digits of the current event
+    /// dump the clusters with associated digits of all events in the current TF
 
     // get the input clusters and associated digits
+    auto rofs = pc.inputs().get<gsl::span<ROFRecord>>("rofs");
     auto clusters = pc.inputs().get<gsl::span<ClusterStruct>>("clusters");
     auto digits = pc.inputs().get<gsl::span<Digit>>("digits");
 
-    if (mText) {
-      // write the clusters in ascii format
-      mOutputFile << clusters.size() << " clusters:" << endl;
-      for (const auto& cluster : clusters) {
-        mOutputFile << cluster << endl;
-      }
-    } else {
-      // write the number of clusters
-      int nClusters = clusters.size();
-      mOutputFile.write(reinterpret_cast<char*>(&nClusters), sizeof(int));
-
-      // write the total number of digits in these clusters
-      int nDigits = digits.size();
-      mOutputFile.write(reinterpret_cast<char*>(&nDigits), sizeof(int));
-
-      // write the clusters
-      mOutputFile.write(reinterpret_cast<const char*>(clusters.data()), clusters.size_bytes());
-
-      // write the digits (after converting the pad ID into a digit UID if requested)
-      if (mUseRun2DigitUID) {
-        std::vector<Digit> digitsCopy(digits.begin(), digits.end());
-        convertPadID2DigitUID(digitsCopy);
-        mOutputFile.write(reinterpret_cast<char*>(digitsCopy.data()), digitsCopy.size() * sizeof(Digit));
+    std::vector<ClusterStruct> eventClusters{};
+    for (const auto& rof : rofs) {
+
+      if (mText) {
+
+        // write the clusters in ascii format
+        mOutputFile << rof.getNEntries() << " clusters:" << endl;
+        for (const auto& cluster : clusters.subspan(rof.getFirstIdx(), rof.getNEntries())) {
+          mOutputFile << cluster << endl;
+        }
+
       } else {
-        mOutputFile.write(reinterpret_cast<const char*>(digits.data()), digits.size_bytes());
+
+        // get the clusters and associated digits of the current event
+        auto eventDigits = getEventClustersAndDigits(rof, clusters, digits, eventClusters);
+
+        // write the number of clusters
+        int nClusters = eventClusters.size();
+        mOutputFile.write(reinterpret_cast<char*>(&nClusters), sizeof(int));
+
+        // write the total number of digits in these clusters
+        int nDigits = eventDigits.size();
+        mOutputFile.write(reinterpret_cast<char*>(&nDigits), sizeof(int));
+
+        // write the clusters
+        mOutputFile.write(reinterpret_cast<const char*>(eventClusters.data()),
+                          eventClusters.size() * sizeof(ClusterStruct));
+
+        // write the digits (after converting the pad ID into a digit UID if requested)
+        if (mUseRun2DigitUID) {
+          std::vector<Digit> digitsCopy(eventDigits.begin(), eventDigits.end());
+          convertPadID2DigitUID(digitsCopy);
+          mOutputFile.write(reinterpret_cast<char*>(digitsCopy.data()), digitsCopy.size() * sizeof(Digit));
+        } else {
+          mOutputFile.write(reinterpret_cast<const char*>(eventDigits.data()), eventDigits.size_bytes());
+        }
       }
     }
   }
 
  private:
+  //_________________________________________________________________________________________________
+  gsl::span<const Digit> getEventClustersAndDigits(const ROFRecord& rof, gsl::span<const ClusterStruct> clusters,
+                                                   gsl::span<const Digit> digits,
+                                                   std::vector<ClusterStruct>& eventClusters) const
+  {
+    /// copy the clusters of the current event (needed to edit the clusters)
+    /// modify the references to the associated digits to start the indexing from 0
+    /// return a sub-span with the associated digits
+
+    eventClusters.clear();
+
+    if (rof.getNEntries() < 1) {
+      return {};
+    }
+
+    if (rof.getLastIdx() >= clusters.size()) {
+      throw length_error("missing clusters");
+    }
+
+    eventClusters.insert(eventClusters.end(), clusters.begin() + rof.getFirstIdx(),
+                         clusters.begin() + rof.getLastIdx() + 1);
+
+    auto digitOffset = eventClusters.front().firstDigit;
+    for (auto& cluster : eventClusters) {
+      cluster.firstDigit -= digitOffset;
+    }
+
+    auto nDigits = eventClusters.back().firstDigit + eventClusters.back().nDigits;
+    if (digitOffset + nDigits > digits.size()) {
+      throw length_error("missing digits");
+    }
+
+    return digits.subspan(digitOffset, nDigits);
+  }
+
   //_________________________________________________________________________________________________
   void convertPadID2DigitUID(std::vector<Digit>& digits)
   {
@@ -153,7 +202,8 @@ o2::framework::DataProcessorSpec getClusterSinkSpec()
 {
   return DataProcessorSpec{
     "ClusterSink",
-    Inputs{InputSpec{"clusters", "MCH", "CLUSTERS", 0, Lifetime::Timeframe},
+    Inputs{InputSpec{"rofs", "MCH", "CLUSTERROFS", 0, Lifetime::Timeframe},
+           InputSpec{"clusters", "MCH", "CLUSTERS", 0, Lifetime::Timeframe},
            InputSpec{"digits", "MCH", "CLUSTERDIGITS", 0, Lifetime::Timeframe}},
     Outputs{},
     AlgorithmSpec{adaptFromTask<ClusterSinkTask>()},