Merge 895dc91 into sapling-pr-archive-ktf

DataFormats/Detectors/GlobalTracking/src/RecoContainer.cxx

@@ -1447,7 +1447,7 @@ RecoContainer::GlobalIDSet RecoContainer::getSingleDetectorRefs(GTrackID gidx) c
     table[GTrackID::MCH] = parent0.getMCHRef();
     table[GTrackID::MID] = parent0.getMIDRef();
   }
-  return std::move(table);
+  return table;
 }
 
 //________________________________________________________

DataFormats/simulation/src/DigitizationContext.cxx

@@ -289,6 +289,11 @@ DigitizationContext* DigitizationContext::loadFromFile(std::string_view filename
 
 void DigitizationContext::fillQED(std::string_view QEDprefix, int max_events, double qedrate)
 {
+  if (mEventRecords.size() <= 1) {
+    // nothing to do
+    return;
+  }
+
   o2::steer::InteractionSampler qedInteractionSampler;
   qedInteractionSampler.setBunchFilling(mBCFilling);
 

Detectors/AOD/include/AODProducerWorkflow/AODProducerWorkflowSpec.h

@@ -482,8 +482,6 @@ class AODProducerWorkflowDPL : public Task
   // using -1 as dummies for AOD
   struct MCLabels {
     uint32_t labelID = -1;
-    uint32_t labelITS = -1;
-    uint32_t labelTPC = -1;
     uint16_t labelMask = 0;
     uint8_t fwdLabelMask = 0;
   };

Detectors/AOD/src/AODProducerWorkflowSpec.cxx

@@ -109,6 +109,7 @@ using PVertex = o2::dataformats::PrimaryVertex;
 using GIndex = o2::dataformats::VtxTrackIndex;
 using DataRequest = o2::globaltracking::DataRequest;
 using GID = o2::dataformats::GlobalTrackID;
+using DetID = o2::detectors::DetID;
 using SMatrix55Sym = ROOT::Math::SMatrix<double, 5, 5, ROOT::Math::MatRepSym<double, 5>>;
 
 namespace o2::aodproducer
@@ -1058,9 +1059,9 @@ void AODProducerWorkflowDPL::fillMCTrackLabelsTable(MCTrackLabelCursorType& mcTr
                                                     int vertexId)
 {
   // labelMask (temporary) usage:
-  //   bit 13 -- ITS/TPC or TPC/TOF labels are not equal
+  //   bit 13 -- ITS/TPC with ITS label (track of AB tracklet) different from TPC
   //   bit 14 -- isNoise() == true
-  //   bit 15 -- isFake() == true
+  //   bit 15 -- isFake() == true (defined by the fakeness of the top level global track, i.e. if TOF is present, fake means that the track of the TPC label does not contribute to TOF cluster)
   // labelID = -1 -- label is not set
 
   for (int src = GIndex::NSources; src--;) {
@@ -1084,7 +1085,7 @@ void AODProducerWorkflowDPL::fillMCTrackLabelsTable(MCTrackLabelCursorType& mcTr
 
       if (GIndex::includesSource(src, mInputSources)) {
         auto mcTruth = data.getTrackMCLabel(trackIndex);
-        MCLabels labelHolder;
+        MCLabels labelHolder{};
         if ((src == GIndex::Source::MFT) || (src == GIndex::Source::MFTMCH) || (src == GIndex::Source::MCH) || (src == GIndex::Source::MCHMID)) { // treating mft and fwd labels separately
           if (!needToStore(src == GIndex::Source::MFT ? mGIDToTableMFTID : mGIDToTableFwdID)) {
             continue;
@@ -1110,51 +1111,22 @@ void AODProducerWorkflowDPL::fillMCTrackLabelsTable(MCTrackLabelCursorType& mcTr
             continue;
           }
           if (mcTruth.isValid()) { // if not set, -1 will be stored
-            labelHolder.labelID = (mToStore[mcTruth.getSourceID()][mcTruth.getEventID()])[mcTruth.getTrackID()];
-          }
-          // treating possible mismatches and fakes for global tracks
-          auto contributorsGID = data.getSingleDetectorRefs(trackIndex);
-          bool isSetTPC = contributorsGID[GIndex::Source::TPC].isIndexSet();
-          bool isSetITS = contributorsGID[GIndex::Source::ITS].isIndexSet();
-          bool isSetTOF = contributorsGID[GIndex::Source::TOF].isIndexSet();
-          bool isTOFFake = true;
-          if (isSetTPC && (isSetITS || isSetTOF)) {
-            auto mcTruthTPC = data.getTrackMCLabel(contributorsGID[GIndex::Source::TPC]);
-            if (mcTruthTPC.isValid()) {
-              labelHolder.labelTPC = (mToStore[mcTruthTPC.getSourceID()][mcTruthTPC.getEventID()])[mcTruthTPC.getTrackID()];
-              labelHolder.labelID = labelHolder.labelTPC;
-            }
-            if (isSetITS) {
-              auto mcTruthITS = data.getTrackMCLabel(contributorsGID[GIndex::Source::ITS]);
-              if (mcTruthITS.isValid()) {
-                labelHolder.labelITS = (mToStore[mcTruthITS.getSourceID()][mcTruthITS.getEventID()])[mcTruthITS.getTrackID()];
-              }
-              if (labelHolder.labelITS != labelHolder.labelTPC) {
-                LOG(debug) << "ITS-TPC MCTruth: labelIDs do not match at " << trackIndex.getIndex() << ", src = " << src;
-                labelHolder.labelMask |= (0x1 << 13);
-              }
+            labelHolder.labelID = (mToStore[mcTruth.getSourceID()][mcTruth.getEventID()])[mcTruth.getTrackID()]; // defined by TPC if it contributes, otherwise: by ITS
+            if (mcTruth.isFake()) {
+              labelHolder.labelMask |= (0x1 << 15);
             }
-            if (isSetTOF) {
-              const auto& labelsTOF = data.getTOFClustersMCLabels()->getLabels(contributorsGID[GIndex::Source::TOF]);
-              for (auto& mcLabel : labelsTOF) {
-                if (!mcLabel.isValid()) {
-                  continue;
-                }
-                if (mcLabel == labelHolder.labelTPC) {
-                  isTOFFake = false;
-                  break;
+            if (trackIndex.includesDet(DetID::TPC) && trackIndex.getSource() != GIndex::Source::TPC) { // this is global track
+              auto contributorsGID = data.getSingleDetectorRefs(trackIndex);
+              if (contributorsGID[GIndex::Source::ITSTPC].isIndexSet()) { // there is a match to ITS tracks or ITSAB tracklet!
+                if (data.getTrackMCLabel(contributorsGID[GIndex::Source::ITSTPC]).isFake()) {
+                  labelHolder.labelMask |= (0x1 << 13);
                 }
               }
             }
-          }
-          if (mcTruth.isFake() || (isSetTOF && isTOFFake)) {
-            labelHolder.labelMask |= (0x1 << 15);
-          }
-          if (mcTruth.isNoise()) {
+          } else if (mcTruth.isNoise()) {
             labelHolder.labelMask |= (0x1 << 14);
           }
-          mcTrackLabelCursor(labelHolder.labelID,
-                             labelHolder.labelMask);
+          mcTrackLabelCursor(labelHolder.labelID, labelHolder.labelMask);
         }
       }
     }

Detectors/GlobalTracking/include/GlobalTracking/MatchTPCITS.h

@@ -50,6 +50,7 @@
 #include "ITSReconstruction/RecoGeomHelper.h"
 #include "TPCFastTransform.h"
 #include "GPUO2InterfaceRefit.h"
+#include "GPUTPCGeometry.h"
 #include "GlobalTracking/MatchTPCITSParams.h"
 #include "DataFormatsITSMFT/TopologyDictionary.h"
 #include "DataFormatsITSMFT/TrkClusRef.h"
@@ -132,6 +133,8 @@ struct TrackLocTPC : public o2::track::TrackParCov {
   int sourceID = 0;                     ///< TPC track origin in
   o2::dataformats::GlobalTrackID gid{}; // global track source ID (TPC track may be part of it)
   int matchID = MinusOne;               ///< entry (non if MinusOne) of its matchTPC struct in the mMatchesTPC
+  uint8_t lowestRow = -1;
+  uint8_t padFromEdge = -1;
   Constraint_t constraint{Constrained};
 
   float getCorrectedTime(float dt) const // return time0 corrected for extra drift (to match certain Z)
@@ -143,7 +146,7 @@ struct TrackLocTPC : public o2::track::TrackParCov {
     return constraint == Constrained ? 0.f : (constraint == ASide ? dt : -dt);
   }
 
-  ClassDefNV(TrackLocTPC, 2);
+  ClassDefNV(TrackLocTPC, 3);
 };
 
 ///< ITS track outward parameters propagated to reference X, with time bracket and index of
@@ -738,6 +741,7 @@ class MatchTPCITS
   static constexpr float MaxSnp = 0.9;                 // max snp of ITS or TPC track at xRef to be matched
   static constexpr float MaxTgp = 2.064;               // max tg corresponting to MaxSnp = MaxSnp/std::sqrt(1.-MaxSnp^2)
   static constexpr float MinTBToCleanCache = 600.;     // keep in AB ITS cluster refs cache at most this number of TPC bins
+  static const o2::gpu::GPUTPCGeometry TPCGeometry;
 
   enum TimerIDs { SWTot,
                   SWPrepITS,

Detectors/GlobalTracking/src/MatchTPCITS.cxx

@@ -9,6 +9,13 @@
 // granted to it by virtue of its status as an Intergovernmental Organization
 // or submit itself to any jurisdiction.
 
+#include "GPUO2Interface.h" // Needed for propper settings in GPUParam.h
+#include "GPUParam.h"
+#include "GPUParam.inc"
+#ifdef WITH_OPENMP
+#include <omp.h>
+#endif
+
 #include <TTree.h>
 #include <cassert>
 #include <algorithm>
@@ -50,13 +57,6 @@
 #include "ITS3Reconstruction/IOUtils.h"
 #endif
 
-#include "GPUO2Interface.h" // Needed for propper settings in GPUParam.h
-#include "GPUParam.h"
-#include "GPUParam.inc"
-#ifdef WITH_OPENMP
-#include <omp.h>
-#endif
-
 using namespace o2::globaltracking;
 
 using MatrixDSym4 = ROOT::Math::SMatrix<double, 4, 4, ROOT::Math::MatRepSym<double, 4>>;
@@ -68,6 +68,8 @@ constexpr float MatchTPCITS::Tan70, MatchTPCITS::Cos70I2, MatchTPCITS::MaxSnp, M
 
 LinksPoolMT* TPCABSeed::gLinksPool = nullptr;
 
+const o2::gpu::GPUTPCGeometry MatchTPCITS::TPCGeometry{};
+
 //______________________________________________
 MatchTPCITS::MatchTPCITS() = default;
 
@@ -428,6 +430,12 @@ int MatchTPCITS::addTPCSeed(const o2::track::TrackParCov& _tr, float t0, float t
   if (clRow > mParams->askMinTPCRow[clSect]) {
     return -9;
   }
+  const auto& clus = mTPCClusterIdxStruct->clusters[clSect][clRow][clIdx];
+  uint8_t padFromEdge = uint8_t(clus.getPad());
+  if (padFromEdge > TPCGeometry.NPads(clRow) / 2) {
+    padFromEdge = TPCGeometry.NPads(clRow) - 1 - padFromEdge;
+  }
+
   // create working copy of track param
   bool extConstrained = srcGID.getSource() != GTrackID::TPC;
   if (extConstrained) {
@@ -442,6 +450,8 @@ int MatchTPCITS::addTPCSeed(const o2::track::TrackParCov& _tr, float t0, float t
                 tpcID,
                 srcGID,
                 MinusOne,
+                clRow,
+                padFromEdge,
                 (extConstrained || tpcOrig.hasBothSidesClusters()) ? TrackLocTPC::Constrained : (tpcOrig.hasASideClustersOnly() ? TrackLocTPC::ASide : TrackLocTPC::CSide)});
   // propagate to matching Xref
   const auto& trackTune = TrackTuneParams::Instance();
@@ -2872,7 +2882,7 @@ void MatchTPCITS::dumpTPCOrig(bool acc, int tpcIndex)
   ///< fill debug tree for TPC original tracks (passing pT cut)
   mTimer[SWDBG].Start(false);
   const auto& tpcOrig = mTPCTracksArray[tpcIndex];
-  uint8_t clSect = 0, clRow = 0, prevRow = 0xff;
+  uint8_t clSect = 0, clRow = 0, prevRow = 0xff, padFromEdge = -1;
   uint32_t clIdx = 0;
   int nshared = 0;
   std::array<bool, 152> shMap{};
@@ -2888,6 +2898,11 @@ void MatchTPCITS::dumpTPCOrig(bool acc, int tpcIndex)
       prevRawShared = true;
     }
   }
+  const auto& clus = mTPCClusterIdxStruct->clusters[clSect][clRow][clIdx];
+  padFromEdge = uint8_t(clus.getPad());
+  if (padFromEdge > TPCGeometry.NPads(clRow) / 2) {
+    padFromEdge = TPCGeometry.NPads(clRow) - 1 - padFromEdge;
+  }
   int tb = tpcOrig.getTime0() * mNTPCOccBinLengthInv;
   float mltTPC = tb < 0 ? mTBinClOcc[0] : (tb >= mTBinClOcc.size() ? mTBinClOcc.back() : mTBinClOcc[tb]);
   (*mDBGOut) << "tpcOrig"
@@ -2900,6 +2915,7 @@ void MatchTPCITS::dumpTPCOrig(bool acc, int tpcIndex)
              << "time0=" << tpcOrig.getTime0()
              << "trc=" << ((o2::track::TrackParCov&)tpcOrig)
              << "minRow=" << clRow
+             << "padFromEdge=" << padFromEdge
              << "multTPC=" << mltTPC;
   if (mMCTruthON) {
     (*mDBGOut) << "tpcOrig"

Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TimeFrameGPU.h

@@ -77,7 +77,8 @@ class TimeFrameGPU : public TimeFrame
   void createCellsDevice();
   void createCellsLUTDevice();
   void createNeighboursIndexTablesDevice();
-  void createNeighboursDevice(const unsigned int& layer, std::vector<std::pair<int, int>>& neighbours);
+  void createNeighboursDevice(const unsigned int layer, const unsigned int nNeighbours);
+  void createNeighboursDevice(const unsigned int layer, std::vector<std::pair<int, int>>& neighbours);
   void createNeighboursLUTDevice(const int, const unsigned int);
   void createNeighboursDeviceArray();
   void createTrackITSExtDevice(std::vector<CellSeed>&);
@@ -151,6 +152,9 @@ class TimeFrameGPU : public TimeFrame
   gsl::span<Tracklet*> getDeviceTracklet() { return mTrackletsDevice; }
   gsl::span<CellSeed*> getDeviceCells() { return mCellsDevice; }
 
+  // Overridden getters
+  int getNumberOfCells() const;
+
  private:
   void allocMemAsync(void**, size_t, Stream*, bool); // Abstract owned and unowned memory allocations
   bool mHostRegistered = false;
@@ -252,6 +256,12 @@ inline std::vector<unsigned int> TimeFrameGPU<nLayers>::getClusterSizes()
   return sizes;
 }
 
+template <int nLayers>
+inline int TimeFrameGPU<nLayers>::getNumberOfCells() const
+{
+  return std::accumulate(mNCells.begin(), mNCells.end(), 0);
+}
+
 } // namespace gpu
 } // namespace its
 } // namespace o2

Detectors/ITSMFT/ITS/tracking/GPU/ITStrackingGPU/TrackingKernels.h

@@ -148,19 +148,19 @@ void computeCellsHandler(const Cluster** sortedClusters,
                          const int nBlocks,
                          const int nThreads);
 
-void countCellNeighboursHandler(CellSeed** cellsLayersDevice,
-                                int* neighboursLUTs,
-                                int** cellsLUTs,
-                                gpuPair<int, int>* cellNeighbours,
-                                int* neighboursIndexTable,
-                                const float maxChi2ClusterAttachment,
-                                const float bz,
-                                const int layerIndex,
-                                const unsigned int nCells,
-                                const unsigned int nCellsNext,
-                                const int maxCellNeighbours,
-                                const int nBlocks,
-                                const int nThreads);
+unsigned int countCellNeighboursHandler(CellSeed** cellsLayersDevice,
+                                        int* neighboursLUTs,
+                                        int** cellsLUTs,
+                                        gpuPair<int, int>* cellNeighbours,
+                                        int* neighboursIndexTable,
+                                        const float maxChi2ClusterAttachment,
+                                        const float bz,
+                                        const int layerIndex,
+                                        const unsigned int nCells,
+                                        const unsigned int nCellsNext,
+                                        const int maxCellNeighbours,
+                                        const int nBlocks,
+                                        const int nThreads);
 
 void computeCellNeighboursHandler(CellSeed** cellsLayersDevice,
                                   int* neighboursLUTs,

Detectors/ITSMFT/ITS/tracking/GPU/cuda/TimeFrameGPU.cu

@@ -406,7 +406,19 @@ void TimeFrameGPU<nLayers>::loadTrackSeedsDevice(std::vector<CellSeed>& seeds)
 }
 
 template <int nLayers>
-void TimeFrameGPU<nLayers>::createNeighboursDevice(const unsigned int& layer, std::vector<std::pair<int, int>>& neighbours)
+void TimeFrameGPU<nLayers>::createNeighboursDevice(const unsigned int layer, const unsigned int nNeighbours)
+{
+  START_GPU_STREAM_TIMER(mGpuStreams[0].get(), "reserving neighbours");
+  LOGP(debug, "gpu-allocation: reserving {} neighbours (pairs), for {} MB.", nNeighbours, nNeighbours * sizeof(gpuPair<int, int>) / MB);
+  allocMemAsync(reinterpret_cast<void**>(&mNeighbourPairsDevice[layer]), nNeighbours * sizeof(gpuPair<int, int>), &(mGpuStreams[0]), getExtAllocator());
+  checkGPUError(cudaMemsetAsync(mNeighbourPairsDevice[layer], -1, nNeighbours * sizeof(gpuPair<int, int>), mGpuStreams[0].get()));
+  LOGP(debug, "gpu-allocation: reserving {} neighbours, for {} MB.", nNeighbours, nNeighbours * sizeof(gpuPair<int, int>) / MB);
+  allocMemAsync(reinterpret_cast<void**>(&mNeighboursDevice[layer]), nNeighbours * sizeof(int), &(mGpuStreams[0]), getExtAllocator());
+  STOP_GPU_STREAM_TIMER(mGpuStreams[0].get());
+}
+
+template <int nLayers>
+void TimeFrameGPU<nLayers>::createNeighboursDevice(const unsigned int layer, std::vector<std::pair<int, int>>& neighbours)
 {
   START_GPU_STREAM_TIMER(mGpuStreams[0].get(), "reserving neighbours");
   mCellsNeighbours[layer].clear();