Merge 6178630 into sapling-pr-archive-ktf

DataFormats/simulation/include/SimulationDataFormat/DigitizationContext.h

@@ -113,6 +113,11 @@ class DigitizationContext
   /// Check collision parts for vertex consistency.
   bool checkVertexCompatibility(bool verbose = false) const;
 
+  /// retrieves collision context for a single timeframe-id (which may be needed by simulation)
+  /// (Only copies collision context without QED information. This can be added to the result with the fillQED method
+  ///  in a second step. As a pre-condition, one should have called finalizeTimeframeStructure)
+  DigitizationContext extractSingleTimeframe(int timeframeid, std::vector<int> const& sources_to_offset);
+
   /// function reading the hits from a chain (previously initialized with initSimChains
   /// The hits pointer will be initialized (what to we do about ownership??)
   template <typename T>
@@ -128,8 +133,9 @@ class DigitizationContext
   // apply collision number cuts and potential relabeling of eventID
   void applyMaxCollisionFilter(long startOrbit, long orbitsPerTF, int maxColl);
 
-  // finalize timeframe structure (fixes the indices in mTimeFrameStartIndex)
-  void finalizeTimeframeStructure(long startOrbit, long orbitsPerTF);
+  /// finalize timeframe structure (fixes the indices in mTimeFrameStartIndex)
+  // returns the number of timeframes
+  int finalizeTimeframeStructure(long startOrbit, long orbitsPerTF);
 
   // Sample and fix interaction vertices (according to some distribution). Makes sure that same event ids
   // have to have same vertex, as well as event ids associated to same collision.
@@ -173,7 +179,7 @@ class DigitizationContext
   // for each collision we may record/fix the interaction vertex (to be used in event generation)
   std::vector<math_utils::Point3D<float>> mInteractionVertices;
 
-  // the collision records _with_ QED interleaved;
+  // the collision records **with** QED interleaved;
   std::vector<o2::InteractionTimeRecord> mEventRecordsWithQED;
   std::vector<std::vector<o2::steer::EventPart>> mEventPartsWithQED;
 

DataFormats/simulation/src/DigitizationContext.cxx

@@ -19,6 +19,7 @@
 #include <numeric> // for iota
 #include <MathUtils/Cartesian.h>
 #include <DataFormatsCalibration/MeanVertexObject.h>
+#include <filesystem>
 
 using namespace o2::steer;
 
@@ -196,10 +197,52 @@ o2::parameters::GRPObject const& DigitizationContext::getGRP() const
 
 void DigitizationContext::saveToFile(std::string_view filename) const
 {
+  // checks if the path content of filename exists ... otherwise it is created before creating the ROOT file
+  auto ensure_path_exists = [](std::string_view filename) {
+    try {
+      // Extract the directory path from the filename
+      std::filesystem::path file_path(filename);
+      std::filesystem::path dir_path = file_path.parent_path();
+
+      // Check if the directory path is empty (which means filename was just a name without path)
+      if (dir_path.empty()) {
+        // nothing to do
+        return true;
+      }
+
+      // Create directories if they do not exist
+      if (!std::filesystem::exists(dir_path)) {
+        if (std::filesystem::create_directories(dir_path)) {
+          // std::cout << "Directories created successfully: " << dir_path.string() << std::endl;
+          return true;
+        } else {
+          std::cerr << "Failed to create directories: " << dir_path.string() << std::endl;
+          return false;
+        }
+      }
+      return true;
+    } catch (const std::filesystem::filesystem_error& ex) {
+      std::cerr << "Filesystem error: " << ex.what() << std::endl;
+      return false;
+    } catch (const std::exception& ex) {
+      std::cerr << "General error: " << ex.what() << std::endl;
+      return false;
+    }
+  };
+
+  if (!ensure_path_exists(filename)) {
+    LOG(error) << "Filename contains path component which could not be created";
+    return;
+  }
+
   TFile file(filename.data(), "RECREATE");
-  auto cl = TClass::GetClass(typeid(*this));
-  file.WriteObjectAny(this, cl, "DigitizationContext");
-  file.Close();
+  if (file.IsOpen()) {
+    auto cl = TClass::GetClass(typeid(*this));
+    file.WriteObjectAny(this, cl, "DigitizationContext");
+    file.Close();
+  } else {
+    LOG(error) << "Could not write to file " << filename.data();
+  }
 }
 
 DigitizationContext* DigitizationContext::loadFromFile(std::string_view filename)
@@ -391,13 +434,15 @@ void DigitizationContext::applyMaxCollisionFilter(long startOrbit, long orbitsPe
   mEventParts = newparts;
 }
 
-void DigitizationContext::finalizeTimeframeStructure(long startOrbit, long orbitsPerTF)
+int DigitizationContext::finalizeTimeframeStructure(long startOrbit, long orbitsPerTF)
 {
   mTimeFrameStartIndex = getTimeFrameBoundaries(mEventRecords, startOrbit, orbitsPerTF);
   LOG(info) << "Fixed " << mTimeFrameStartIndex.size() << " timeframes ";
   for (auto p : mTimeFrameStartIndex) {
     LOG(info) << p.first << " " << p.second;
   }
+
+  return mTimeFrameStartIndex.size();
 }
 
 std::unordered_map<int, int> DigitizationContext::getCollisionIndicesForSource(int source) const
@@ -483,3 +528,55 @@ void DigitizationContext::sampleInteractionVertices(o2::dataformats::MeanVertexO
     }
   }
 }
+
+DigitizationContext DigitizationContext::extractSingleTimeframe(int timeframeid, std::vector<int> const& sources_to_offset)
+{
+  DigitizationContext r; // make a return object
+  if (mTimeFrameStartIndex.size() == 0) {
+    LOG(error) << "No timeframe structure determined; Returning empty object. Please call ::finalizeTimeframeStructure before calling this function";
+    return r;
+  }
+  r.mSimPrefixes = mSimPrefixes;
+  r.mMuBC = mMuBC;
+  try {
+    auto startend = mTimeFrameStartIndex.at(timeframeid);
+
+    auto startindex = startend.first;
+    auto endindex = startend.second;
+
+    std::copy(mEventRecords.begin() + startindex, mEventRecords.begin() + endindex, std::back_inserter(r.mEventRecords));
+    std::copy(mEventParts.begin() + startindex, mEventParts.begin() + endindex, std::back_inserter(r.mEventParts));
+    if (mInteractionVertices.size() > endindex) {
+      std::copy(mInteractionVertices.begin() + startindex, mInteractionVertices.begin() + endindex, std::back_inserter(r.mInteractionVertices));
+    }
+
+    // let's assume we want to fix the ids for source = source_id
+    // Then we find the first index that has this source_id and take the corresponding number
+    // as offset. Thereafter we subtract this offset from all known event parts.
+    auto perform_offsetting = [&r](int source_id) {
+      auto indices_for_source = r.getCollisionIndicesForSource(source_id);
+      int minvalue = std::numeric_limits<int>::max();
+      for (auto& p : indices_for_source) {
+        if (p.first < minvalue) {
+          minvalue = p.first;
+        }
+      }
+      // now fix them
+      for (auto& p : indices_for_source) {
+        auto index_into_mEventParts = p.second;
+        for (auto& part : r.mEventParts[index_into_mEventParts]) {
+          if (part.sourceID == source_id) {
+            part.entryID -= minvalue;
+          }
+        }
+      }
+    };
+    for (auto source_id : sources_to_offset) {
+      perform_offsetting(source_id);
+    }
+
+  } catch (std::exception) {
+    LOG(warn) << "No such timeframe id in collision context. Returing empty object";
+  }
+  return r;
+}

Detectors/GlobalTrackingWorkflow/study/include/GlobalTrackingStudy/TrackInfoExt.h

@@ -27,6 +27,7 @@ namespace dataformats
 struct TrackInfoExt {
   o2::track::TrackParCov track;
   DCA dca{};
+  DCA dcaTPC{};
   VtxTrackIndex gid;
   MatchInfoTOF infoTOF;
   float ttime = 0;

Detectors/GlobalTrackingWorkflow/study/src/TrackingStudy.cxx

@@ -364,10 +364,30 @@ void TrackingStudySpec::process(o2::globaltracking::RecoContainer& recoData)
           continue;
         }
         {
+          o2::dataformats::DCA dcaTPC;
+          dcaTPC.set(-999.f, -999.f);
+          if (tpcTr) {
+            if (is == GTrackID::TPC) {
+              dcaTPC = dca;
+            } else {
+              o2::track::TrackParCov tmpTPC(*tpcTr);
+              if (iv < nv - 1 && is == GTrackID::TPC && tpcTr && !tpcTr->hasBothSidesClusters()) { // for unconstrained TPC tracks correct track Z
+                float corz = vdrit * (tpcTr->getTime0() * mTPCTBinMUS - pvvec[iv].getTimeStamp().getTimeStamp());
+                if (tpcTr->hasASideClustersOnly()) {
+                  corz = -corz; // A-side
+                }
+                tmpTPC.setZ(tmpTPC.getZ() + corz);
+              }
+              if (!prop->propagateToDCA(iv == nv - 1 ? vtxDummy : pvvec[iv], tmpTPC, prop->getNominalBz(), 2., o2::base::PropagatorF::MatCorrType::USEMatCorrLUT, &dcaTPC)) {
+                dcaTPC.set(-999.f, -999.f);
+              }
+            }
+          }
           auto& trcExt = trcExtVec.emplace_back();
           recoData.getTrackTime(vid, trcExt.ttime, trcExt.ttimeE);
           trcExt.track = trc;
           trcExt.dca = dca;
+          trcExt.dcaTPC = dcaTPC;
           trcExt.gid = vid;
           trcExt.xmin = xmin;
           auto gidRefs = recoData.getSingleDetectorRefs(vid);

Framework/Core/include/Framework/DataRefUtils.h

@@ -122,7 +122,7 @@ struct DataRefUtils {
         // object only depends on the state at serialization of the original object. However,
         // all objects created during deserialization are new and must be owned by the collection
         // to avoid memory leak. So we call SetOwner if it is available for the type.
-        if constexpr (has_root_setowner<T>::value) {
+        if constexpr (requires(T t) { t.SetOwner(true); }) {
           result->SetOwner(true);
         }
       });

Framework/Core/include/Framework/TypeTraits.h

@@ -147,22 +147,5 @@ class has_root_dictionary<T, typename std::enable_if<is_container<T>::value>::ty
 {
 };
 
-// Detect whether a class is a ROOT class implementing SetOwner
-// This member detector idiom is implemented using SFINAE idiom to look for
-// a 'SetOwner()' method.
-template <typename T, typename _ = void>
-struct has_root_setowner : std::false_type {
-};
-
-template <typename T>
-struct has_root_setowner<
-  T,
-  std::conditional_t<
-    false,
-    class_member_checker<
-      decltype(std::declval<T>().SetOwner(true))>,
-    void>> : public std::true_type {
-};
-
 } // namespace o2::framework
 #endif // FRAMEWORK_TYPETRAITS_H

GPU/GPUTracking/Base/GPUReconstructionCPU.cxx

@@ -216,8 +216,9 @@ int32_t GPUReconstructionCPU::RunChains()
 
   timerTotal.Start();
   if (mProcessingSettings.doublePipeline) {
-    if (EnqueuePipeline()) {
-      return 1;
+    int32_t retVal = EnqueuePipeline();
+    if (retVal) {
+      return retVal;
     }
   } else {
     if (mThreadId != GetThread()) {

GPU/GPUTracking/Definitions/GPUSettingsList.h

@@ -89,7 +89,7 @@ AddOptionRTC(extraClusterErrorSplitPadSharedSingleY2, float, 0.03f, "", 0, "Addi
 AddOptionRTC(extraClusterErrorFactorSplitPadSharedSingleY2, float, 3.0f, "", 0, "Multiplicative extra cluster error for Y2 if splitpad, shared, or single set")
 AddOptionRTC(extraClusterErrorSplitTimeSharedSingleZ2, float, 0.03f, "", 0, "Additive extra cluster error for Z2 if splittime, shared, or single set")
 AddOptionRTC(extraClusterErrorFactorSplitTimeSharedSingleZ2, float, 3.0f, "", 0, "Multiplicative extra cluster error for Z2 if splittime, shared, or single set")
-AddOptionArray(errorsCECrossing, float, 5, (0.f, 0.f, 0.f, 0.f, 0.f), "", 0, "Extra errors to add to track when crossing CE, depending on addErrorsCECrossing") // BUG: CUDA cannot yet hand AddOptionArrayRTC
+AddOptionArray(errorsCECrossing, float, 5, (0.f, 0.f, 0.f, 0.f, 0.f), "", 0, "Extra errors to add to track when crossing CE, depending on addErrorsCECrossing") // BUG: CUDA cannot yet handle AddOptionArrayRTC
 AddOptionRTC(globalTrackingYRangeUpper, float, 0.85f, "", 0, "Inner portion of y-range in slice that is not used in searching for global track candidates")
 AddOptionRTC(globalTrackingYRangeLower, float, 0.85f, "", 0, "Inner portion of y-range in slice that is not used in searching for global track candidates")
 AddOptionRTC(trackFollowingYFactor, float, 4.f, "", 0, "Weight of y residual vs z residual in tracklet constructor")

GPU/GPUTracking/Global/GPUChainTrackingDebugAndProfiling.cxx

@@ -302,7 +302,7 @@ void GPUChainTracking::SanityCheck()
 void GPUChainTracking::RunTPCClusterFilter(o2::tpc::ClusterNativeAccess* clusters, std::function<o2::tpc::ClusterNative*(size_t)> allocator, bool applyClusterCuts)
 {
   GPUTPCClusterFilter clusterFilter(*clusters);
-  o2::tpc::ClusterNative* outputBuffer;
+  o2::tpc::ClusterNative* outputBuffer = nullptr;
   for (int32_t iPhase = 0; iPhase < 2; iPhase++) {
     uint32_t countTotal = 0;
     for (uint32_t iSector = 0; iSector < GPUCA_NSLICES; iSector++) {

GPU/Workflow/src/GPUWorkflowSpec.cxx

@@ -851,7 +851,7 @@ void GPURecoWorkflowSpec::run(ProcessingContext& pc)
       }
       createEmptyOutput = !mConfParam->partialOutputForNonFatalErrors;
     } else {
-      throw std::runtime_error("tracker returned error code " + std::to_string(retVal));
+      throw std::runtime_error("GPU Reconstruction error: error code " + std::to_string(retVal));
     }
   }