Wrap the QED event IDs at the number of QED events available

Steer/include/Steer/MCKinematicsReader.h

@@ -121,6 +121,11 @@ class MCKinematicsReader
   }
 
  private:
+  /// slow path of the track accessors: loads what is missing, or reports an event that is not there
+  void ensureTracksForSourceAndEvent(int source, int event) const;
+  [[noreturn]] static void reportMissingSource(int source, size_t available);
+  [[noreturn]] static void reportMissingEvent(const char* what, int source, int event, size_t available);
+
   void initTracksForSource(int source) const;
   void loadTracksForSourceAndEvent(int source, int eventID) const;
   void loadHeadersForSource(int source) const;
@@ -151,7 +156,9 @@ inline MCTrack const* MCKinematicsReader::getTrack(o2::MCCompLabel const& label)
 
 inline MCTrack const* MCKinematicsReader::getTrack(int source, int event, int track) const
 {
-  return &getTracks(source, event)[track];
+  auto const& tracks = getTracks(source, event);
+  // one comparison, and it covers the negative track ID of a hit whose track was not kept
+  return static_cast<size_t>(track) < tracks.size() ? &tracks[track] : nullptr;
 }
 
 inline MCTrack const* MCKinematicsReader::getTrack(int event, int track) const
@@ -161,13 +168,18 @@ inline MCTrack const* MCKinematicsReader::getTrack(int event, int track) const
 
 inline std::vector<MCTrack> const& MCKinematicsReader::getTracks(int source, int event) const
 {
-  if (mTracks[source].size() == 0) {
+  if (static_cast<size_t>(source) >= mTracks.size()) {
+    reportMissingSource(source, mTracks.size());
+  }
+  auto& perEvent = mTracks[source];
+  if (perEvent.size() == 0) {
     initTracksForSource(source);
   }
-  if (mTracks[source][event] == nullptr) {
-    loadTracksForSourceAndEvent(source, event);
+  // the event range shares the branch of the lazy load, so the fast path grows by one comparison
+  if (static_cast<size_t>(event) >= perEvent.size() || perEvent[event] == nullptr) {
+    ensureTracksForSourceAndEvent(source, event);
   }
-  return *mTracks[source][event];
+  return *perEvent[event];
 }
 
 inline std::vector<MCTrack> const& MCKinematicsReader::getTracks(int event) const
@@ -177,26 +189,41 @@ inline std::vector<MCTrack> const& MCKinematicsReader::getTracks(int event) cons
 
 inline o2::dataformats::MCEventHeader const& MCKinematicsReader::getMCEventHeader(int source, int event) const
 {
-  if (mHeaders.at(source).size() == 0) {
+  auto const& headers = mHeaders.at(source);
+  if (headers.size() == 0) {
     loadHeadersForSource(source);
   }
-  return mHeaders.at(source)[event];
+  if (static_cast<size_t>(event) >= headers.size()) {
+    reportMissingEvent("event headers", source, event, headers.size());
+  }
+  return headers[event];
 }
 
 inline gsl::span<o2::TrackReference> MCKinematicsReader::getTrackRefs(int source, int event, int track) const
 {
-  if (mIndexedTrackRefs[source].size() == 0) {
+  if (static_cast<size_t>(source) >= mIndexedTrackRefs.size()) {
+    return {};
+  }
+  auto& perEvent = mIndexedTrackRefs[source];
+  if (perEvent.size() == 0) {
     loadTrackRefsForSource(source);
   }
-  return mIndexedTrackRefs[source][event].getLabels(track);
+  if (static_cast<size_t>(event) >= perEvent.size()) {
+    return {};
+  }
+  return perEvent[event].getLabels(track);
 }
 
 inline const std::vector<o2::TrackReference>& MCKinematicsReader::getTrackRefsByEvent(int source, int event) const
 {
-  if (mIndexedTrackRefs[source].size() == 0) {
+  auto const& perEvent = mIndexedTrackRefs.at(source);
+  if (perEvent.size() == 0) {
     loadTrackRefsForSource(source);
   }
-  return mIndexedTrackRefs[source][event].getTruthArray();
+  if (static_cast<size_t>(event) >= perEvent.size()) {
+    reportMissingEvent("events of track references", source, event, perEvent.size());
+  }
+  return perEvent[event].getTruthArray();
 }
 
 inline gsl::span<o2::TrackReference> MCKinematicsReader::getTrackRefs(int event, int track) const

Steer/src/CollisionContextTool.cxx

@@ -18,6 +18,7 @@
 #include "CommonDataFormat/InteractionRecord.h"
 #include "DataFormatsCalibration/MeanVertexObject.h"
 #include "SimulationDataFormat/DigitizationContext.h"
+#include "SimulationDataFormat/MCEventLabel.h"
 #include "SimConfig/InteractionDiamondParam.h"
 #include "DataFormatsFT0/EventsPerBc.h"
 #include <cmath>
@@ -215,6 +216,25 @@ InteractionSpec parseInteractionSpec(std::string const& specifier, std::vector<I
   }
 }
 
+// The number of QED events which the collision context may cycle through.
+// The QED events are reused in a round robin over the sampled QED interactions, so the wrap has to
+// happen at the number of events that exist in the QED kinematics. Wrapping later hands out event
+// IDs which were never simulated: the hits are read modulo the file size while the MC labels keep
+// the unwrapped ID. See https://its.cern.ch/jira/browse/O2-7132
+int getQEDRoundRobinSize(InteractionSpec const& qedSpec)
+{
+  if (qedSpec.mcnumberavail <= 0) {
+    LOG(warn) << "No number of available QED events given (the MCNUMBERSTRING of --QEDinteraction); "
+              << "QED event IDs may name events which are not in the QED kinematics";
+    return qedSpec.mcnumberasked;
+  }
+  if (qedSpec.mcnumberavail > (int)o2::MCEventLabel::MaxEventID()) {
+    LOG(warn) << "The QED production has " << qedSpec.mcnumberavail << " events, more than the "
+              << o2::MCEventLabel::MaxEventID() << " an MCEventLabel can encode; QED event IDs will be truncated";
+  }
+  return qedSpec.mcnumberavail;
+}
+
 bool parseOptions(int argc, char* argv[], Options& optvalues)
 {
   namespace bpo = boost::program_options;
@@ -721,7 +741,7 @@ int main(int argc, char* argv[])
     // TODO: use bcFilling information
     auto qedSpec = parseInteractionSpec(options.qedInteraction, ispecs, options.useexistingkinematics);
     std::cout << "### IRATE " << qedSpec.interactionRate << "\n";
-    digicontext.fillQED(qedSpec.name, qedSpec.mcnumberasked, qedSpec.interactionRate);
+    digicontext.fillQED(qedSpec.name, getQEDRoundRobinSize(qedSpec), qedSpec.interactionRate);
   }
 
   if (options.printContext) {
@@ -786,7 +806,7 @@ int main(int argc, char* argv[])
         // This should probably be done inside the extraction itself
         if (digicontext.isQEDProvided()) {
           auto qedSpec = parseInteractionSpec(options.qedInteraction, ispecs, options.useexistingkinematics);
-          copy.fillQED(qedSpec.name, qedSpec.mcnumberasked, qedSpec.interactionRate);
+          copy.fillQED(qedSpec.name, getQEDRoundRobinSize(qedSpec), qedSpec.interactionRate);
         }
 
         std::stringstream str;

Steer/src/MCKinematicsReader.cxx

@@ -14,11 +14,33 @@
 #include "SimulationDataFormat/MCEventHeader.h"
 #include "SimulationDataFormat/TrackReference.h"
 #include <TChain.h>
+#include <stdexcept>
+#include <string>
 #include <vector>
 #include <fairlogger/Logger.h>
 
 using namespace o2::steer;
 
+void MCKinematicsReader::reportMissingSource(int source, size_t available)
+{
+  throw std::out_of_range("MCKinematicsReader: there are " + std::to_string(available) + " sources; source " +
+                          std::to_string(source) + " is not one of them");
+}
+
+void MCKinematicsReader::reportMissingEvent(const char* what, int source, int event, size_t available)
+{
+  throw std::out_of_range("MCKinematicsReader: source " + std::to_string(source) + " has " +
+                          std::to_string(available) + " " + what + "; there is no event " + std::to_string(event));
+}
+
+void MCKinematicsReader::ensureTracksForSourceAndEvent(int source, int event) const
+{
+  if (static_cast<size_t>(event) >= mTracks[source].size()) {
+    reportMissingEvent("events", source, event, mTracks[source].size());
+  }
+  loadTracksForSourceAndEvent(source, event);
+}
+
 MCKinematicsReader::~MCKinematicsReader()
 {
   for (auto chain : mInputChains) {