add to debug output original TPC tracks tree (#13090)

shahor02 · web-flow · commit f3f160216700 · 2024-04-29T18:23:06.000+02:00
diff --git a/Detectors/GlobalTracking/include/GlobalTracking/MatchTPCITS.h b/Detectors/GlobalTracking/include/GlobalTracking/MatchTPCITS.h
@@ -423,9 +423,10 @@ class MatchTPCITS
 
 #ifdef _ALLOW_DEBUG_TREES_
   enum DebugFlagTypes : UInt_t {
-    MatchTreeAll = 0x1 << 1,     ///< produce matching candidates tree for all candidates
-    MatchTreeAccOnly = 0x1 << 2, ///< fill the matching candidates tree only once the cut is passed
-    WinnerMatchesTree = 0x1 << 3 ///< separate debug tree for winner matches
+    MatchTreeAll = 0x1 << 1,      ///< produce matching candidates tree for all candidates
+    MatchTreeAccOnly = 0x1 << 2,  ///< fill the matching candidates tree only once the cut is passed
+    WinnerMatchesTree = 0x1 << 3, ///< separate debug tree for winner matches
+    TPCOrigTree = 0x1 << 4        ///< original TPC tracks with some aux info
   };
   ///< check if partucular flags are set
   bool isDebugFlag(UInt_t flags) const { return mDBGFlags & flags; }
@@ -450,6 +451,7 @@ class MatchTPCITS
   ///< fill matching debug tree
   void fillTPCITSmatchTree(int itsID, int tpcID, int rejFlag, float chi2 = -1., float tCorr = 0.);
   void dumpWinnerMatches();
+  void dumpTPCOrig(bool acc, int tpcIndex);
 #endif
 
  private:
@@ -464,7 +466,7 @@ class MatchTPCITS
   bool prepareFITData();
   int prepareInteractionTimes();
   int prepareTPCTracksAfterBurner();
-  void addTPCSeed(const o2::track::TrackParCov& _tr, float t0, float terr, o2::dataformats::GlobalTrackID srcGID, int tpcID);
+  int addTPCSeed(const o2::track::TrackParCov& _tr, float t0, float terr, o2::dataformats::GlobalTrackID srcGID, int tpcID);
 
   int preselectChipClusters(std::vector<int>& clVecOut, const ClusRange& clRange, const ITSChipClustersRefs& itsChipClRefs,
                             float trackY, float trackZ, float tolerY, float tolerZ) const;
diff --git a/Detectors/GlobalTracking/src/MatchTPCITS.cxx b/Detectors/GlobalTracking/src/MatchTPCITS.cxx
@@ -402,23 +402,23 @@ int MatchTPCITS::getNMatchRecordsITS(const TrackLocITS& tTPC) const
 }
 
 //______________________________________________
-void MatchTPCITS::addTPCSeed(const o2::track::TrackParCov& _tr, float t0, float terr, GTrackID srcGID, int tpcID)
+int MatchTPCITS::addTPCSeed(const o2::track::TrackParCov& _tr, float t0, float terr, GTrackID srcGID, int tpcID)
 {
   // account single TPC seed, can be from standalone TPC track or constrained track from match to TRD and/or TOF
   const float SQRT12DInv = 2. / sqrt(12.);
   if (_tr.getX() > o2::constants::geom::XTPCInnerRef + 0.1 || std::abs(_tr.getQ2Pt()) > mMinTPCTrackPtInv) {
-    return;
+    return -99;
   }
   const auto& tpcOrig = mTPCTracksArray[tpcID];
   // discard tracks w/o certain number of total or innermost pads (last cluster is innermost one)
   if (tpcOrig.getNClusterReferences() < mParams->minTPCClusters) {
-    return;
+    return -89;
   }
   uint8_t clSect = 0, clRow = 0;
   uint32_t clIdx = 0;
   tpcOrig.getClusterReference(mTPCTrackClusIdx, tpcOrig.getNClusterReferences() - 1, clSect, clRow, clIdx);
   if (clRow > mParams->askMinTPCRow[clSect]) {
-    return;
+    return -9;
   }
   // create working copy of track param
   bool extConstrained = srcGID.getSource() != GTrackID::TPC;
@@ -449,13 +449,14 @@ void MatchTPCITS::addTPCSeed(const o2::track::TrackParCov& _tr, float t0, float
   }
   if (!propagateToRefX(trc)) {
     mTPCWork.pop_back(); // discard track whose propagation to XMatchingRef failed
-    return;
+    return -1;
   }
   if (mMCTruthON) {
     mTPCLblWork.emplace_back(mTPCTrkLabels[tpcID]);
   }
   // cache work track index
   mTPCSectIndexCache[o2::math_utils::angle2Sector(trc.getAlpha())].push_back(mTPCWork.size() - 1);
+  return 0;
 }
 
 //______________________________________________
@@ -485,6 +486,32 @@ bool MatchTPCITS::prepareTPCData()
     mTPCSectIndexCache[sec].reserve(100 + 1.2 * ntrW / o2::constants::math::NSectors);
   }
 
+  mTPCRefitter = std::make_unique<o2::gpu::GPUO2InterfaceRefit>(mTPCClusterIdxStruct, mTPCCorrMapsHelper, mBz, mTPCTrackClusIdx.data(), 0, mTPCRefitterShMap.data(), mTPCRefitterOccMap.data(), mTPCRefitterOccMap.size(), nullptr, o2::base::Propagator::Instance());
+  mTPCRefitter->setTrackReferenceX(900); // disable propagation after refit by setting reference to value > 500
+  mNTPCOccBinLength = mTPCRefitter->getParam()->rec.tpc.occupancyMapTimeBins;
+  mTBinClOcc.clear();
+  if (mNTPCOccBinLength > 1 && mTPCRefitterOccMap.size()) {
+    mNTPCOccBinLengthInv = 1. / mNTPCOccBinLength;
+    int nTPCBins = mNHBPerTF * o2::constants::lhc::LHCMaxBunches / 8, ninteg = 0;
+    int nTPCOccBins = nTPCBins * mNTPCOccBinLengthInv, sumBins = std::max(1, int(o2::constants::lhc::LHCMaxBunches / 8 * mNTPCOccBinLengthInv));
+    mTBinClOcc.resize(nTPCOccBins);
+    std::vector<float> mltHistTB(nTPCOccBins);
+    float sm = 0., tb = 0.5 * mNTPCOccBinLength;
+    for (int i = 0; i < nTPCOccBins; i++) {
+      mltHistTB[i] = mTPCRefitter->getParam()->GetUnscaledMult(tb);
+      tb += mNTPCOccBinLength;
+    }
+    for (int i = nTPCOccBins; i--;) {
+      sm += mltHistTB[i];
+      if (i + sumBins < nTPCOccBins) {
+        sm -= mltHistTB[i + sumBins];
+      }
+      mTBinClOcc[i] = sm;
+    }
+  } else {
+    mTBinClOcc.resize(1);
+  }
+
   auto creator = [this](auto& trk, GTrackID gid, float time0, float terr) {
     if constexpr (isITSTrack<decltype(trk)>()) {
       // do nothing, ITS tracks will be processed in a direct loop over ROFs
@@ -494,20 +521,27 @@ bool MatchTPCITS::prepareTPCData()
     } else if (std::abs(trk.getQ2Pt()) > mMinTPCTrackPtInv) {
       return true;
     }
+    int resAdd = -100;
+    int tpcIndex = -1;
     if constexpr (isTPCTrack<decltype(trk)>()) {
       // unconstrained TPC track, with t0 = TrackTPC.getTime0+0.5*(DeltaFwd-DeltaBwd) and terr = 0.5*(DeltaFwd+DeltaBwd) in TimeBins
       if (!this->mSkipTPCOnly && trk.getNClusters() > 0) {
-        this->addTPCSeed(trk, this->tpcTimeBin2MUS(time0), this->tpcTimeBin2MUS(terr), gid, gid.getIndex());
+        resAdd = this->addTPCSeed(trk, this->tpcTimeBin2MUS(time0), this->tpcTimeBin2MUS(terr), gid, (tpcIndex = gid.getIndex()));
       }
     }
     if constexpr (isTPCTOFTrack<decltype(trk)>()) {
       // TPC track constrained by TOF time, time and its error in \mus
-      this->addTPCSeed(trk, time0, terr, gid, this->mRecoCont->getTPCContributorGID(gid));
+      resAdd = this->addTPCSeed(trk, time0, terr, gid, (tpcIndex = this->mRecoCont->getTPCContributorGID(gid)));
     }
     if constexpr (isTRDTrack<decltype(trk)>()) {
       // TPC track constrained by TRD trigger time, time and its error in \mus
-      this->addTPCSeed(trk, time0, terr, gid, this->mRecoCont->getTPCContributorGID(gid));
+      resAdd = this->addTPCSeed(trk, time0, terr, gid, (tpcIndex = this->mRecoCont->getTPCContributorGID(gid)));
     }
+#ifdef _ALLOW_DEBUG_TREES_
+    if (resAdd > -10 && mDBGOut && isDebugFlag(TPCOrigTree)) {
+      dumpTPCOrig(resAdd == 0, tpcIndex);
+    }
+#endif
     // note: TPCTRDTPF tracks are actually TRD track with extra TOF cluster
     return true;
   };
@@ -576,31 +610,6 @@ bool MatchTPCITS::prepareTPCData()
     mITSROFofTPCBin[ib] = itsROF;
   }
 */
-  mTPCRefitter = std::make_unique<o2::gpu::GPUO2InterfaceRefit>(mTPCClusterIdxStruct, mTPCCorrMapsHelper, mBz, mTPCTrackClusIdx.data(), 0, mTPCRefitterShMap.data(), mTPCRefitterOccMap.data(), mTPCRefitterOccMap.size(), nullptr, o2::base::Propagator::Instance());
-  mTPCRefitter->setTrackReferenceX(900); // disable propagation after refit by setting reference to value > 500
-  mNTPCOccBinLength = mTPCRefitter->getParam()->rec.tpc.occupancyMapTimeBins;
-  mTBinClOcc.clear();
-  if (mNTPCOccBinLength > 1 && mTPCRefitterOccMap.size()) {
-    mNTPCOccBinLengthInv = 1. / mNTPCOccBinLength;
-    int nTPCBins = mNHBPerTF * o2::constants::lhc::LHCMaxBunches / 8, ninteg = 0;
-    int nTPCOccBins = nTPCBins * mNTPCOccBinLengthInv, sumBins = std::max(1, int(o2::constants::lhc::LHCMaxBunches / 8 * mNTPCOccBinLengthInv));
-    mTBinClOcc.resize(nTPCOccBins);
-    std::vector<float> mltHistTB(nTPCOccBins);
-    float sm = 0., tb = 0.5 * mNTPCOccBinLength;
-    for (int i = 0; i < nTPCOccBins; i++) {
-      mltHistTB[i] = mTPCRefitter->getParam()->GetUnscaledMult(tb);
-      tb += mNTPCOccBinLength;
-    }
-    for (int i = nTPCOccBins; i--;) {
-      sm += mltHistTB[i];
-      if (i + sumBins < nTPCOccBins) {
-        sm -= mltHistTB[i + sumBins];
-      }
-      mTBinClOcc[i] = sm;
-    }
-  } else {
-    mTBinClOcc.resize(1);
-  }
   mInteractionMUSLUT.clear();
   mInteractionMUSLUT.resize(maxTime + 3 * o2::constants::lhc::LHCOrbitMUS, -1);
   mTimer[SWPrepTPC].Stop();
@@ -2816,6 +2825,36 @@ void MatchTPCITS::setDebugFlag(UInt_t flag, bool on)
   }
 }
 
+//_________________________________________________________
+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;
+  uint32_t clIdx = 0;
+  tpcOrig.getClusterReference(mTPCTrackClusIdx, tpcOrig.getNClusterReferences() - 1, clSect, clRow, clIdx);
+  int tb = tpcOrig.getTime0() * mNTPCOccBinLengthInv;
+  float mltTPC = tb < 0 ? mTBinClOcc[0] : (tb >= mTBinClOcc.size() ? mTBinClOcc.back() : mTBinClOcc[tb]);
+  (*mDBGOut) << "tpcOrig"
+             << "tf=" << mTFCount
+             << "index=" << tpcIndex
+             << "acc=" << acc
+             << "chi2TPC=" << tpcOrig.getChi2()
+             << "nClus=" << tpcOrig.getNClusters()
+             << "time0=" << tpcOrig.getTime0()
+             << "trc=" << ((o2::track::TrackParCov&)tpcOrig)
+             << "minRow=" << clRow
+             << "multTPC=" << mltTPC;
+  if (mMCTruthON) {
+    (*mDBGOut) << "tpcOrig"
+               << "tpcLbl=" << mTPCTrkLabels[tpcIndex];
+  }
+  (*mDBGOut) << "tpcOrig"
+             << "\n";
+  mTimer[SWDBG].Stop();
+}
+
 //_________________________________________________________
 void MatchTPCITS::fillTPCITSmatchTree(int itsID, int tpcID, int rejFlag, float chi2, float tCorr)
 {
