One more step

Author: Elros60

Common/Tasks/muonQa.cxx

@@ -165,6 +165,7 @@ struct muonQa {
 
   Preslice<aod::FwdTrkCl> perMuon = aod::fwdtrkcl::fwdtrackId;
   Preslice<aod::FwdTrackAssoc> fwdtrackIndicesPerCollision = aod::track_association::collisionId;
+  Preslice<aod::MFTTrackAssoc> mfttrackIndicesPerCollision = aod::track_association::collisionId;
 
   HistogramRegistry registry{"registry", {}};
   HistogramRegistry registryDCA{"registryDCA", {}};
@@ -653,20 +654,50 @@ struct muonQa {
     }
   }
 
-  template <typename TEvent, typename TFwdTracks, typename TMFTTracks, typename TFwdCls, typename AssocMuons>
-  void runMuonQA(TEvent const& collision, aod::BCsWithTimestamps const&, TFwdTracks const& muons, TMFTTracks const& mfts, TMuonCls const& clusters, AssocMuons const& fwdtrackIndices)
+  template <typename TEvent, typename TFwdTracks, typename TMFTTracks, typename TFwdCls, typename AssocMuons, typename AssocMFTs>
+  void runMuonQA(TEvent const& collision, aod::BCsWithTimestamps const&, TFwdTracks const& muons, TMFTTracks const& mfts, TMuonCls const& clusters, AssocMuons const& fwdtrackIndices, AssocMFTs const& mfttrackIndices)
   {
+    // Loop over re-associated forward tracks
+    for (auto const& muonId : fwdtrackIndices) {
+      auto muon = muonId.template fwdtrack_as<TFwdTracks>();
+      // Consider only global muons having both MFT and MCH informations
+      if (static_cast<int>(muon.trackType()) < 2) {
+        auto mfttrack = muon.template matchMFTTrack_as<TMFTTracks>(); // Get the matched MFT track
+        auto mchtrack = muon.template matchMCHTrack_as<TFwdTracks>(); // Get the matched MCH/MCH-MID trakc
+
+        // Re-do the propagation with muon tracking tools
+        // Do MFT tracks propagation
+        double chi2 = muon.chi2();
+        SMatrix5 tpars(muon.x(), muon.y(), muon.phi(), muon.tgl(), muon.signed1Pt());
+        std::vector<double> v1{muon.cXX(), muon.cXY(), muon.cYY(), muon.cPhiX(), muon.cPhiY(),
+                               muon.cPhiPhi(), muon.cTglX(), muon.cTglY(), muon.cTglPhi(), muon.cTglTgl(),
+                               muon.c1PtX(), muon.c1PtY(), muon.c1PtPhi(), muon.c1PtTgl(), muon.c1Pt21Pt2()};
+        SMatrix55 tcovs(v1.begin(), v1.end());
+        o2::track::TrackParCovFwd fwdtrack{muon.z(), tpars, tcovs, chi2};
+        o2::dataformats::GlobalFwdTrack propmuon;
+        double centerMFT[3] = {0, 0, -61.4};
+        o2::field::MagneticField* field = static_cast<o2::field::MagneticField*>(TGeoGlobalMagField::Instance()->GetField());
+        auto Bz = field->getBz(centerMFT); // Get field at centre of MFT
+        auto geoMan = o2::base::GeometryManager::meanMaterialBudget(muon.x(), muon.y(), muon.z(), collision.posX(), collision.posY(), collision.posZ());
+        auto x2x0 = static_cast<float>(geoMan.meanX2X0);
+        fwdtrack.propagateToVtxhelixWithMCS(collision.posZ(), {collision.posX(), collision.posY()}, {collision.covXX(), collision.covYY()}, Bz, x2x0);
+        propmuon.setParameters(fwdtrack.getParameters());
+        propmuon.setZ(fwdtrack.getZ());
+        propmuon.setCovariances(fwdtrack.getCovariances());
+      }
+    }
   }
 
-  void processMuonQa(MyEvents const& collisions, aod::BCsWithTimestamps const& bcs, MyMuonsWithCov const& fwdTracks, MyMFTs const& mftTracks, aod::FwdTrkCls const& muonClusters, aod::FwdTrackAssoc const& fwdtrackIndices)
+  void processMuonQa(MyEvents const& collisions, aod::BCsWithTimestamps const& bcs, MyMuonsWithCov const& fwdTracks, MyMFTs const& mftTracks, aod::FwdTrkCls const& muonClusters, aod::FwdTrackAssoc const& fwdtrackIndices, aod::MFTTrackAssoc const& mfttrackIndices)
   {
     // Update CCDB informations
     initCCDB(bcs);
 
     // Run over collisions
     for (auto& collision : collisions) {
       auto muonIdsThisCollision = fwdtrackIndices.sliceBy(fwdtrackIndicesPerCollision, collision.globalIndex());
-      runMuonQA(collision, bcs, fwdTracks, mftTracks, muonClusters, muonIdsThisCollision);
+      auto mftIdsThisCollision = mfttrackIndices.sliceBy(mfttrackIndicesPerCollision, collision.globalIndex());
+      runMuonQA(collision, bcs, fwdTracks, mftTracks, muonClusters, muonIdsThisCollision, mftIdsThisCollision);
     }
   }
   PROCESS_SWITCH(muonQa, processMuonQa, "Process to run muon QA", true);