4545#include < Framework/HistogramSpec.h>
4646#include < Framework/InitContext.h>
4747#include < Framework/runDataProcessing.h>
48+ #include < PID/PIDTOFParamService.h>
4849#include < ReconstructionDataFormats/DCA.h>
4950#include < ReconstructionDataFormats/PID.h>
5051#include < ReconstructionDataFormats/Track.h>
5152
52- // #include <Math/Vector4D.h> // IWYU pragma: keep (do not replace with Math/Vector4Dfwd.h)
53- // #include <Math/Vector4Dfwd.h>
5453#include < TH1.h>
5554
5655#include < array>
6059#include < random>
6160#include < string>
6261#include < string_view>
62+ #include < unordered_map>
6363#include < utility>
6464#include < vector>
6565
@@ -78,7 +78,7 @@ struct taggingHFE {
7878
7979 using MyTracks = soa::Join<aod::TracksIU, aod::TracksExtra, aod::TracksCovIU,
8080 aod::pidTPCFullEl, aod::pidTPCFullPi, aod::pidTPCFullKa, aod::pidTPCFullPr,
81- aod::pidTOFFullEl, aod::pidTOFFullPi, aod::pidTOFFullKa, aod::pidTOFFullPr, aod::pidTOFbeta>;
81+ aod::pidTOFFullEl, aod::pidTOFFullPi, aod::pidTOFFullKa, aod::pidTOFFullPr, aod::pidTOFbeta, aod::TOFSignal, aod::TOFEvTime >;
8282 using MyTracksWithMCLabel = soa::Join<MyTracks, aod::McTrackLabels, aod::mcTPCTuneOnData>;
8383
8484 // using MyV0s = soa::Join<aod::V0Datas, aod::V0Covs, aod::V0TOFNSigmas, aod::V0CoreMCLabels>;
@@ -101,6 +101,7 @@ struct taggingHFE {
101101 Configurable<float > d_bz_input{" d_bz_input" , -999 , " bz field in kG, -999 is automatic" };
102102 Configurable<int > cfgPdgLepton{" cfgPdgLepton" , 11 , " pdg code of desired lepton: 11 or 13" };
103103 Configurable<float > cfgDownSampling{" cfgDownSampling" , 1.1 , " down sampling for fake matches" };
104+ Configurable<bool > useTOFNSigmaDeltaBC{" useTOFNSigmaDeltaBC" , true , " Flag to shift delta BC for TOF n sigma (only with TTCA)" };
104105
105106 struct : ConfigurableGroup {
106107 std::string prefix = " dcaFitterGroup_eK" ;
@@ -122,7 +123,7 @@ struct taggingHFE {
122123
123124 struct : ConfigurableGroup {
124125 std::string prefix = " electronCut" ;
125- Configurable<float > cfg_min_pt_track{" cfg_min_pt_track" , 0.4 , " min pT for single track" };
126+ Configurable<float > cfg_min_pt_track{" cfg_min_pt_track" , 0.1 , " min pT for single track" };
126127 Configurable<float > cfg_max_pt_track{" cfg_max_pt_track" , 1e+10 , " max pT for single track" };
127128 Configurable<float > cfg_min_eta_track{" cfg_min_eta_track" , -0.8 , " min eta for single track" };
128129 Configurable<float > cfg_max_eta_track{" cfg_max_eta_track" , +0.8 , " max eta for single track" };
@@ -167,7 +168,6 @@ struct taggingHFE {
167168 Configurable<float > cfg_min_TOFNsigmaKa{" cfg_min_TOFNsigmaKa" , -3 , " min n sigma ka in TOF" };
168169 Configurable<float > cfg_max_TOFNsigmaKa{" cfg_max_TOFNsigmaKa" , +3 , " max n sigma ka in TOF" };
169170 Configurable<bool > requirePiKa{" requirePiKa" , false , " require hadron to be pion or kaon" }; // proton is not involved in semileptonic decay of HF hadrons often.
170- Configurable<bool > requireKa{" requireKa" , false , " require hadron to be kaon" }; // Mostly, kaon is involved in semileptonic decay of charm hadrons.
171171 } hadronCut;
172172
173173 struct : ConfigurableGroup {
@@ -272,14 +272,15 @@ struct taggingHFE {
272272 } lCPairCut;
273273
274274 o2::aod::rctsel::RCTFlagsChecker rctChecker;
275+ Service<o2::pid::tof::TOFResponse> mTOFResponse ;
275276
276277 HistogramRegistry fRegistry {" fRegistry" };
277278 static constexpr std::string_view hadron_names[6 ] = {" LF/" , " Jpsi/" , " D0/" , " Dpm/" , " Ds/" , " Lc/" };
278279 static constexpr std::string_view pair_names[3 ] = {" e_Kpm/" , " e_K0S/" , " e_Lambda/" };
279280 static constexpr std::string_view hTypes[4 ] = {" findable/" , " correct/" , " fake/" , " miss/" };
280281 static constexpr std::string_view promptTypes[2 ] = {" prompt/" , " nonprompt/" };
281282
282- void init (o2::framework::InitContext&)
283+ void init (o2::framework::InitContext& initContext )
283284 {
284285 // if (doprocessSA && doprocessTTCA) {
285286 // LOGF(fatal, "Cannot enable doprocessWithoutFTTCA and doprocessWithFTTCA at the same time. Please choose one.");
@@ -291,6 +292,9 @@ struct taggingHFE {
291292 ccdb->setFatalWhenNull (false );
292293 rctChecker.init (eventCut.cfgRCTLabel .value , eventCut.cfgCheckZDC .value , eventCut.cfgTreatLimitedAcceptanceAsBad .value );
293294
295+ LOGF (info, " intializing TOFResponse" );
296+ mTOFResponse ->initSetup (ccdb, initContext);
297+
294298 std::random_device seed_gen;
295299 engine = std::mt19937 (seed_gen ());
296300 dist01 = std::uniform_real_distribution<float >(0 .0f , 1 .0f );
@@ -452,21 +456,15 @@ struct taggingHFE {
452456 fRegistry .add (" Cascade/hMassOmega" , " #Omega mass;m_{#LambdaK} (GeV/c^{2})" , kTH1F , {{100 , 1.62 , 1.72 }}, false );
453457 }
454458
455- template <typename TTrack>
456- bool isKaon (TTrack const & track)
457- {
458- bool is_ka_included_TPC = hadronCut.cfg_min_TPCNsigmaKa < track.tpcNSigmaKa () && track.tpcNSigmaKa () < hadronCut.cfg_max_TPCNsigmaKa ;
459- bool is_ka_included_TOF = track.hasTOF () ? (hadronCut.cfg_min_TOFNsigmaKa < track.tofNSigmaKa () && track.tofNSigmaKa () < hadronCut.cfg_max_TOFNsigmaKa ) : true ;
460- return is_ka_included_TPC && is_ka_included_TOF;
461- }
462-
463- template <typename TTrack>
464- bool isKaon_or_isPion (TTrack const & track)
459+ template <typename TCollision, typename TTrack>
460+ bool isKaon_or_isPion (TCollision const & collision, TTrack const & track)
465461 {
462+ float tofNSigmaPi = mapTOFNsigmaPiReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())];
463+ float tofNSigmaKa = mapTOFNsigmaKaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())];
466464 bool is_ka_included_TPC = hadronCut.cfg_min_TPCNsigmaKa < track.tpcNSigmaKa () && track.tpcNSigmaKa () < hadronCut.cfg_max_TPCNsigmaKa ;
467- bool is_ka_included_TOF = track.hasTOF () ? (hadronCut.cfg_min_TOFNsigmaKa < track. tofNSigmaKa () && track. tofNSigmaKa () < hadronCut.cfg_max_TOFNsigmaKa ) : true ;
465+ bool is_ka_included_TOF = track.hasTOF () ? (hadronCut.cfg_min_TOFNsigmaKa < tofNSigmaKa && tofNSigmaKa < hadronCut.cfg_max_TOFNsigmaKa ) : true ;
468466 bool is_pi_included_TPC = hadronCut.cfg_min_TPCNsigmaPi < track.tpcNSigmaPi () && track.tpcNSigmaPi () < hadronCut.cfg_max_TPCNsigmaPi ;
469- bool is_pi_included_TOF = track.hasTOF () ? (hadronCut.cfg_min_TOFNsigmaPi < track. tofNSigmaPi () && track. tofNSigmaPi () < hadronCut.cfg_max_TOFNsigmaPi ) : true ;
467+ bool is_pi_included_TOF = track.hasTOF () ? (hadronCut.cfg_min_TOFNsigmaPi < tofNSigmaPi && tofNSigmaPi < hadronCut.cfg_max_TOFNsigmaPi ) : true ;
470468 return (is_ka_included_TPC && is_ka_included_TOF) || (is_pi_included_TPC && is_pi_included_TOF);
471469 }
472470
@@ -479,6 +477,14 @@ struct taggingHFE {
479477 return is_pi_included_TPC && is_pi_included_TOF;
480478 }
481479
480+ template <typename TTrack>
481+ bool isKaon (TTrack const & track)
482+ {
483+ bool is_ka_included_TPC = hadronCut.cfg_min_TPCNsigmaKa < track.tpcNSigmaKa () && track.tpcNSigmaKa () < hadronCut.cfg_max_TPCNsigmaKa ;
484+ bool is_ka_included_TOF = track.hasTOF () ? (hadronCut.cfg_min_TOFNsigmaKa < track.tofNSigmaKa () && track.tofNSigmaKa () < hadronCut.cfg_max_TOFNsigmaKa ) : true ;
485+ return is_ka_included_TPC && is_ka_included_TOF;
486+ }
487+
482488 template <typename TTrack>
483489 bool isProton (TTrack const & track)
484490 {
@@ -559,8 +565,8 @@ struct taggingHFE {
559565 return true ;
560566 }
561567
562- template <typename TTrack, typename TTrackParCov>
563- bool isSelectedHadron (TTrack const & track, TTrackParCov const & trackParCov, const float dcaXY, const float dcaZ)
568+ template <typename TCollision, typename TTrack, typename TTrackParCov>
569+ bool isSelectedHadron (TCollision const & collision, TTrack const & track, TTrackParCov const & trackParCov, const float dcaXY, const float dcaZ)
564570 {
565571 if (!track.hasITS () || !track.hasTPC ()) {
566572 return false ;
@@ -614,11 +620,7 @@ struct taggingHFE {
614620 return false ;
615621 }
616622
617- if (hadronCut.requirePiKa && !isKaon_or_isPion (track)) {
618- return false ;
619- }
620-
621- if (hadronCut.requireKa && !isKaon (track)) {
623+ if (hadronCut.requirePiKa && !isKaon_or_isPion (collision, track)) {
622624 return false ;
623625 }
624626
@@ -827,6 +829,87 @@ struct taggingHFE {
827829 return true ;
828830 }
829831
832+ std::map<std::pair<int , int >, float > mapTOFNsigmaPiReassociated; // map pair(collisionId, trackId) -> tof n sigma pi
833+ std::map<std::pair<int , int >, float > mapTOFNsigmaKaReassociated; // map pair(collisionId, trackId) -> tof n sigma ka
834+ std::map<std::pair<int , int >, float > mapTOFBetaReassociated; // map pair(collisionId, trackId) -> tof beta
835+ std::unordered_map<int , double > mapCollisionTime;
836+ std::unordered_map<int , double > mapCollisionTimeError;
837+
838+ template <bool withTTCA, typename TCollisions, typename TBCs, typename TTracks, typename TTrackAssoc>
839+ void calculateTOFNSigmaWithReassociation (TCollisions const & collisions, TBCs const &, TTracks const & tracks, TTrackAssoc const & trackIndices)
840+ {
841+ if (useTOFNSigmaDeltaBC) {
842+ if constexpr (withTTCA) {
843+ for (const auto & collision : collisions) {
844+ if (mapCollisionTime.find (collision.globalIndex ()) == mapCollisionTime.end ()) {
845+ continue ;
846+ }
847+ auto bcCollision = collision.template bc_as <TBCs>();
848+ auto trackIdsThisCollision = trackIndices.sliceBy (trackIndicesPerCollision, collision.globalIndex ());
849+ for (const auto & trackId : trackIdsThisCollision) {
850+ auto track = trackId.template track_as <TTracks>();
851+ if (!track.hasITS () || !track.hasTPC ()) { // apply only minimal cut
852+ continue ;
853+ }
854+
855+ if (track.hasTOF () && track.has_collision ()) { // TTCA may use orphan tracks.
856+ auto bcTrack = track.template collision_as <TCollisions>().template bc_as <TBCs>();
857+ float tofNSigmaPi = mTOFResponse ->nSigma <o2::track::PID::Pion>(track.tofSignalInAnotherBC (bcTrack.globalBC (), bcCollision.globalBC ()), track.tofExpMom (), track.length (), track.p (), track.eta (), mapCollisionTime[collision.globalIndex ()], mapCollisionTimeError[collision.globalIndex ()]);
858+ float tofNSigmaKa = mTOFResponse ->nSigma <o2::track::PID::Kaon>(track.tofSignalInAnotherBC (bcTrack.globalBC (), bcCollision.globalBC ()), track.tofExpMom (), track.length (), track.p (), track.eta (), mapCollisionTime[collision.globalIndex ()], mapCollisionTimeError[collision.globalIndex ()]);
859+ float beta = track.length () / (track.tofSignalInAnotherBC (bcTrack.globalBC (), bcCollision.globalBC ()) - mapCollisionTime[collision.globalIndex ()]) / (TMath::C () * 1e+2 * 1e-12 );
860+ mapTOFNsigmaPiReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = tofNSigmaPi;
861+ mapTOFNsigmaKaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = tofNSigmaKa;
862+ mapTOFBetaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = beta;
863+ } else {
864+ mapTOFNsigmaPiReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.tofNSigmaPi ();
865+ mapTOFNsigmaKaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.tofNSigmaKa ();
866+ mapTOFBetaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.beta ();
867+ }
868+ } // end of track loop
869+ } // end of collision loop
870+ } else {
871+ for (const auto & collision : collisions) {
872+ auto tracks_per_coll = tracks.sliceBy (perCol, collision.globalIndex ());
873+ for (const auto & track : tracks_per_coll) {
874+ if (!track.hasITS () || !track.hasTPC ()) { // apply only minimal cut
875+ continue ;
876+ }
877+ mapTOFNsigmaPiReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.tofNSigmaPi ();
878+ mapTOFNsigmaKaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.tofNSigmaKa ();
879+ mapTOFBetaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.beta ();
880+ }
881+ } // end of track loop
882+ } // end of collision loop
883+ } else {
884+ if constexpr (withTTCA) {
885+ for (const auto & collision : collisions) {
886+ auto trackIdsThisCollision = trackIndices.sliceBy (trackIndicesPerCollision, collision.globalIndex ());
887+ for (const auto & trackId : trackIdsThisCollision) {
888+ auto track = trackId.template track_as <TTracks>();
889+ if (!track.hasITS () || !track.hasTPC ()) { // apply only minimal cut
890+ continue ;
891+ }
892+ mapTOFNsigmaPiReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.tofNSigmaPi ();
893+ mapTOFNsigmaKaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.tofNSigmaKa ();
894+ mapTOFBetaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.beta ();
895+ } // end of track loop
896+ } // end of collision loop
897+ } else {
898+ for (const auto & collision : collisions) {
899+ auto tracks_per_coll = tracks.sliceBy (perCol, collision.globalIndex ());
900+ for (const auto & track : tracks_per_coll) {
901+ if (!track.hasITS () || !track.hasTPC ()) { // apply only minimal cut
902+ continue ;
903+ }
904+ mapTOFNsigmaPiReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.tofNSigmaPi ();
905+ mapTOFNsigmaKaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.tofNSigmaKa ();
906+ mapTOFBetaReassociated[std::make_pair (collision.globalIndex (), track.globalIndex ())] = track.beta ();
907+ }
908+ } // end of track loop
909+ } // end of collision loop
910+ }
911+ }
912+
830913 template <bool isMC, typename TBCs, typename TCollisions, typename TTracks, typename TTrackAssoc, typename TV0s, typename TCascades, typename TMCCollisions, typename TMCParticles>
831914 void runPairing (TBCs const &, TCollisions const & collisions, TTracks const & tracks, TTrackAssoc const & trackIndices, TV0s const & v0s, TCascades const & cascades, TMCCollisions const &, TMCParticles const & mcParticles)
832915 {
@@ -898,7 +981,7 @@ struct taggingHFE {
898981
899982 fRegistry .fill (HIST (" Electron/hs" ), trackParCov.getPt (), trackParCov.getEta (), RecoDecay::constrainAngle (trackParCov.getPhi (), 0 , 1U ));
900983 fRegistry .fill (HIST (" Electron/hTPCdEdx" ), track.tpcInnerParam (), track.mcTunedTPCSignal ());
901- fRegistry .fill (HIST (" Electron/hTOFbeta" ), track.p (), track.beta () );
984+ fRegistry .fill (HIST (" Electron/hTOFbeta" ), track.p (), mapTOFBetaReassociated[ std::make_pair (collision. globalIndex (), track.globalIndex ())] );
902985 if (track.sign () > 0 ) { // positron
903986 positronIds.emplace_back (trackId.trackId ());
904987 } else { // electron
@@ -912,10 +995,10 @@ struct taggingHFE {
912995 dcaXY = mDcaInfoCov .getY ();
913996 dcaZ = mDcaInfoCov .getZ ();
914997
915- if (isSelectedHadron (track, trackParCov, dcaXY, dcaZ)) { // electrons can be included in hadron sample.
998+ if (isSelectedHadron (collision, track, trackParCov, dcaXY, dcaZ)) { // electrons can be included in hadron sample.
916999 fRegistry .fill (HIST (" Hadron/hs" ), trackParCov.getPt (), trackParCov.getEta (), RecoDecay::constrainAngle (trackParCov.getPhi (), 0 , 1U ));
9171000 fRegistry .fill (HIST (" Hadron/hTPCdEdx" ), track.tpcInnerParam (), track.mcTunedTPCSignal ());
918- fRegistry .fill (HIST (" Hadron/hTOFbeta" ), track.p (), track.beta () );
1001+ fRegistry .fill (HIST (" Hadron/hTOFbeta" ), track.p (), mapTOFBetaReassociated[ std::make_pair (collision. globalIndex (), track.globalIndex ())] );
9191002 if (track.sign () > 0 ) { // K+
9201003 kaonPlusIds.emplace_back (trackId.trackId ());
9211004 } else { // K-
@@ -1070,9 +1153,13 @@ struct taggingHFE {
10701153 }
10711154 }
10721155
1156+ float tofNSigmaPi = mapTOFNsigmaPiReassociated[std::make_pair (collision.globalIndex (), kaon.globalIndex ())];
1157+ float tofNSigmaKa = mapTOFNsigmaKaReassociated[std::make_pair (collision.globalIndex (), kaon.globalIndex ())];
1158+
10731159 emmllhpair (leptonTable.lastIndex (),
10741160 trackParCov.getQ2Pt (), trackParCov.getEta (), dcaXY_kaon, dcaZ_kaon, trackParCov.getSigmaY2 (), trackParCov.getSigmaZY (), trackParCov.getSigmaZ2 (),
1075- kaon.tpcNSigmaKa (), kaon.tofNSigmaKa (),
1161+ kaon.tpcNSigmaPi (), tofNSigmaPi,
1162+ kaon.tpcNSigmaKa (), tofNSigmaKa,
10761163 eKpair.mass , eKpair.dca2legs , eKpair.cospa , eKpair.cospaXY ,
10771164 eKpair.lxyz , eKpair.lxyzErr ,
10781165 eKpair.lxy , eKpair.lxyErr ,
@@ -1386,9 +1473,13 @@ struct taggingHFE {
13861473 }
13871474 }
13881475
1476+ float tofNSigmaPi = mapTOFNsigmaPiReassociated[std::make_pair (collision.globalIndex (), kaon.globalIndex ())];
1477+ float tofNSigmaKa = mapTOFNsigmaKaReassociated[std::make_pair (collision.globalIndex (), kaon.globalIndex ())];
1478+
13891479 emmllhpair (leptonTable.lastIndex (),
13901480 trackParCov.getQ2Pt (), trackParCov.getEta (), dcaXY_kaon, dcaZ_kaon, trackParCov.getSigmaY2 (), trackParCov.getSigmaZY (), trackParCov.getSigmaZ2 (),
1391- kaon.tpcNSigmaKa (), kaon.tofNSigmaKa (),
1481+ kaon.tpcNSigmaPi (), tofNSigmaPi,
1482+ kaon.tpcNSigmaKa (), tofNSigmaKa,
13921483 eKpair.mass , eKpair.dca2legs , eKpair.cospa , eKpair.cospaXY ,
13931484 eKpair.lxyz , eKpair.lxyzErr ,
13941485 eKpair.lxy , eKpair.lxyErr ,
@@ -1814,8 +1905,26 @@ struct taggingHFE {
18141905
18151906 void processMC (FilteredMyCollisionsWithMCLabel const & collisions, aod::BCsWithTimestamps const & bcs, MyTracksWithMCLabel const & tracks, aod::TrackAssoc const & trackIndices, filteredV0s const & v0s, filteredMyCascades const & cascades, aod::McCollisions const & mcCollisions, aod::McParticles const & mcParticles)
18161907 {
1908+ initCCDB (bcs.iteratorAt (0 ));
1909+ mTOFResponse ->processSetup (bcs.iteratorAt (0 ));
1910+
1911+ for (const auto & track : tracks) {
1912+ if (mapCollisionTime.find (track.collisionId ()) == mapCollisionTime.end ()) {
1913+ // LOGF(info, "track.collisionId() = %d, track.tofEvTime() = %f, track.tofEvTimeErr() = %f", track.collisionId(), track.tofEvTime(), track.tofEvTimeErr());
1914+ mapCollisionTime[track.collisionId ()] = track.tofEvTime ();
1915+ mapCollisionTimeError[track.collisionId ()] = track.tofEvTimeErr ();
1916+ }
1917+ }
1918+ calculateTOFNSigmaWithReassociation<true >(collisions, bcs, tracks, trackIndices);
1919+
18171920 runPairing<true >(bcs, collisions, tracks, trackIndices, v0s, cascades, mcCollisions, mcParticles);
18181921 runGen (mcCollisions, mcParticles);
1922+
1923+ mapCollisionTime.clear ();
1924+ mapCollisionTimeError.clear ();
1925+ mapTOFNsigmaPiReassociated.clear ();
1926+ mapTOFNsigmaKaReassociated.clear ();
1927+ mapTOFBetaReassociated.clear ();
18191928 }
18201929 PROCESS_SWITCH (taggingHFE, processMC, " process with TTCA" , true );
18211930
@@ -1824,5 +1933,6 @@ struct taggingHFE {
18241933};
18251934WorkflowSpec defineDataProcessing (ConfigContext const & cfgc)
18261935{
1936+ o2::pid::tof::TOFResponseImpl::metadataInfo.initMetadata (cfgc);
18271937 return WorkflowSpec{adaptAnalysisTask<taggingHFE>(cfgc, TaskName{" tagging-hfe" })};
18281938}
0 commit comments