1515
1616// #include <TDatabasePDG.h>
1717#include " PWGLF/DataModel/LFStrangenessPIDTables.h"
18- #include " PWGLF/DataModel/LFStrangenessTables.h" //
18+ #include " PWGLF/DataModel/LFStrangenessTables.h"
19+ #include " PWGLF/DataModel/mcCentrality.h"
20+ #include " PWGLF/Utils/inelGt.h"
1921
2022#include " Common/Core/TrackSelection.h"
2123#include " Common/Core/trackUtilities.h"
@@ -64,6 +66,7 @@ using namespace o2::soa;
6466using namespace o2 ::aod::rctsel;
6567// using namespace o2::constants::physics;
6668using std::array;
69+ // FixME: Add INEL>0 selection for the derived data
6770
6871struct HigherMassResonances {
6972 SliceCache cache;
@@ -108,7 +111,7 @@ struct HigherMassResonances {
108111 // Configurable<bool> isSel8{"isSel8", false, "Event Selection 8"};
109112
110113 // Configurables for event selection
111- // Configurable<bool> isINELgt0{"isINELgt0", true, "INEL>0 selection"};
114+ Configurable<bool > isINELgt0{" isINELgt0" , true , " INEL>0 selection" };
112115 Configurable<bool > isTriggerTVX{" isTriggerTVX" , false , " TriggerTVX" };
113116 // Configurable<bool> isGoodZvtxFT0vsPV{"isGoodZvtxFT0vsPV", false, "IsGoodZvtxFT0vsPV"};
114117 // Configurable<bool> isApplyOccCut{"isApplyOccCut", true, "Apply occupancy cut"};
@@ -270,8 +273,8 @@ struct HigherMassResonances {
270273 hCutFlow->GetXaxis ()->SetBinLabel (8 , " Occupancy Cut" );
271274 hCutFlow->GetXaxis ()->SetBinLabel (9 , " rctChecker" );
272275 hCutFlow->GetXaxis ()->SetBinLabel (10 , " kIsTriggerTVX" );
273- // hCutFlow->GetXaxis()->SetBinLabel(11, "kIsGoodZvtxFT0vsPV ");
274- // hCutFlow->GetXaxis()->SetBinLabel(12, "IsINELgt0");
276+ hCutFlow->GetXaxis ()->SetBinLabel (11 , " No selection " );
277+ hCutFlow->GetXaxis ()->SetBinLabel (12 , " IsINELgt0" );
275278 // hCutFlow->GetXaxis()->SetBinLabel(13, "isVertexITSTPC");
276279 // hCutFlow->GetXaxis()->SetBinLabel(14, "isVertexTOFMatched");
277280
@@ -367,9 +370,13 @@ struct HigherMassResonances {
367370 // For MC
368371 if (doprocessGen || doprocessRec) {
369372 hMChists.add (" Genf1710" , " Gen f_{0}(1710)" , kTHnSparseF , {multiplicityAxis, ptAxis, thnAxisPOL});
373+ hMChists.add (" Genf1710Calib" , " Calibrated Gen f_{0}(1710)" , kTHnSparseF , {multiplicityAxis, ptAxis, thnAxisPOL});
370374 hMChists.add (" Genf17102" , " Gen f_{0}(1710)" , kTHnSparseF , {multiplicityAxis, ptAxis, thnAxisPOL});
375+ hMChists.add (" Genf1710Calib2" , " Calibrated Gen f_{0}(1710)" , kTHnSparseF , {multiplicityAxis, ptAxis, thnAxisPOL});
371376 hMChists.add (" Recf1710_pt1" , " Rec f_{0}(1710) p_{T}" , kTHnSparseF , {multiplicityAxis, ptAxis, glueballMassAxis, thnAxisPOL});
377+ hMChists.add (" Recf1710Calib_pt1" , " Calibrated Rec f_{0}(1710) p_{T}" , kTHnSparseF , {multiplicityAxis, ptAxis, glueballMassAxis, thnAxisPOL});
372378 hMChists.add (" Recf1710_pt2" , " Rec f_{0}(1710) p_{T}" , kTHnSparseF , {multiplicityAxis, ptAxis, glueballMassAxis, thnAxisPOL});
379+ hMChists.add (" Recf1710Calib_pt2" , " Calibrated Rec f_{0}(1710) p_{T}" , kTHnSparseF , {multiplicityAxis, ptAxis, glueballMassAxis, thnAxisPOL});
373380 hMChists.add (" h1Recsplit" , " Rec p_{T}2" , kTH1F , {ptAxis});
374381 hMChists.add (" Genf1710_mass" , " Gen f_{0}(1710) mass" , kTH1F , {glueballMassAxis});
375382 hMChists.add (" Genf1710_mass2" , " Gen f_{0}(1710) mass" , kTH1F , {glueballMassAxis});
@@ -388,6 +395,13 @@ struct HigherMassResonances {
388395 hMChists.add (" Rec_Multiplicity" , " Multiplicity in MC" , kTH1F , {multiplicityAxis});
389396 hMChists.add (" MC_mult_after_event_sel" , " Multiplicity in MC" , kTH1F , {multiplicityAxis});
390397 }
398+
399+ if (doprocessEvtLossSigLossMC) {
400+ hMChists.add (" MCcorrections/hSignalLossDenominator" , " Kstar generated before event selection" , kTH2F , {{ptAxis}, {multiplicityAxis}});
401+ hMChists.add (" MCcorrections/hSignalLossNumerator" , " Kstar generated after event selection" , kTH2F , {{ptAxis}, {multiplicityAxis}});
402+ hMChists.add (" MCcorrections/MultiplicityRec" , " Multiplicity in generated MC with at least 1 reconstruction" , kTH1F , {multiplicityAxis});
403+ hMChists.add (" MCcorrections/MultiplicityGen" , " Multiplicity in generated MC" , kTH1F , {multiplicityAxis});
404+ }
391405 }
392406
393407 template <typename Coll>
@@ -447,13 +461,13 @@ struct HigherMassResonances {
447461 // if (config.isGoodZvtxFT0vsPV && !collision.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV))
448462 // return false;
449463 // if (fillHist)
450- // rEventSelection.fill(HIST("hEventCut"), 10);
464+ rEventSelection.fill (HIST (" hEventCut" ), 10 );
451465
452- // if (config.isINELgt0 && !collision.isInelGt0()) {
453- // return false;
454- // }
455- // if (fillHist)
456- // rEventSelection.fill(HIST("hEventCut"), 11);
466+ if (config.isINELgt0 && !collision.isInelGt0 ()) {
467+ return false ;
468+ }
469+ if (fillHist)
470+ rEventSelection.fill (HIST (" hEventCut" ), 11 );
457471
458472 // if (config.isVertexITSTPC && !collision.selection_bit(o2::aod::evsel::kIsVertexITSTPC)) {
459473 // return false;
@@ -470,6 +484,63 @@ struct HigherMassResonances {
470484 return true ;
471485 }
472486
487+ template <typename Coll>
488+ bool selectionEventDerived (const Coll& collision, bool fillHist = true )
489+ {
490+ if (fillHist)
491+ rEventSelection.fill (HIST (" hEventCut" ), 0 );
492+
493+ if (std::abs (collision.posZ ()) > config.cutzvertex )
494+ return false ;
495+ if (fillHist)
496+ rEventSelection.fill (HIST (" hEventCut" ), 1 );
497+
498+ // if (config.isSel8 && !collision.sel8())
499+ // return false;
500+ if (fillHist)
501+ rEventSelection.fill (HIST (" hEventCut" ), 2 );
502+
503+ if (config.isNoTimeFrameBorder && !collision.selection_bit (aod::evsel::kNoTimeFrameBorder ))
504+ return false ;
505+ if (fillHist)
506+ rEventSelection.fill (HIST (" hEventCut" ), 3 );
507+
508+ if (config.isNoITSROFrameBorder && !collision.selection_bit (aod::evsel::kNoITSROFrameBorder ))
509+ return false ;
510+ if (fillHist)
511+ rEventSelection.fill (HIST (" hEventCut" ), 4 );
512+
513+ // if (config.isNoSameBunchPileup && (!collision.selection_bit(aod::evsel::kNoSameBunchPileup)))
514+ // return false;
515+ if (fillHist)
516+ rEventSelection.fill (HIST (" hEventCut" ), 5 );
517+
518+ if (config.isAllLayersGoodITS && !collision.selection_bit (o2::aod::evsel::kIsGoodITSLayersAll ))
519+ return false ;
520+ if (fillHist)
521+ rEventSelection.fill (HIST (" hEventCut" ), 6 );
522+
523+ // if (config.isNoCollInTimeRangeStandard && (!collision.selection_bit(o2::aod::evsel::kNoCollInTimeRangeStandard)))
524+ // return false;
525+
526+ // if (config.isApplyOccCut && (std::abs(collision.trackOccupancyInTimeRange()) > config.configOccCut))
527+ // return false;
528+ if (fillHist)
529+ rEventSelection.fill (HIST (" hEventCut" ), 7 );
530+
531+ if (rctCut.requireRCTFlagChecker && !rctChecker (collision))
532+ return false ;
533+ if (fillHist)
534+ rEventSelection.fill (HIST (" hEventCut" ), 8 );
535+
536+ if (config.isTriggerTVX && !collision.selection_bit (aod::evsel::kIsTriggerTVX ))
537+ return false ;
538+ if (fillHist)
539+ rEventSelection.fill (HIST (" hEventCut" ), 9 );
540+
541+ return true ;
542+ }
543+
473544 template <typename Collision, typename V0>
474545 bool selectionV0 (Collision const & collision, V0 const & candidate, float /* multiplicity*/ )
475546 {
@@ -714,9 +785,10 @@ struct HigherMassResonances {
714785 using TrackCandidates = soa::Filtered<soa::Join<aod::Tracks, aod::TrackSelection, aod::TracksExtra, aod::TracksDCA, aod::pidTPCFullPi, aod::pidTOFFullPi>>;
715786 using V0TrackCandidate = soa::Join<aod::V0Datas, aod::V0TOFPIDs, aod::V0TOFNSigmas>;
716787 // For Monte Carlo
717- using EventCandidatesMC = soa::Join<aod::Collisions, aod::EvSels, aod::McCollisionLabels, aod::CentFT0Cs , aod::CentFT0Ms , aod::CentFT0As , aod::CentFV0As , aod::PVMults>;
788+ using EventCandidatesMC = soa::Join<aod::Collisions, aod::EvSels, aod::McCollisionLabels, aod::CentFT0Ms , aod::CentFT0As , aod::CentFT0Cs , aod::MultZeqs , aod::FT0Mults, aod:: PVMults, aod::CentFV0As >;
718789 using TrackCandidatesMC = soa::Filtered<soa::Join<aod::Tracks, aod::TracksExtra, aod::TracksDCA, aod::TrackSelection, aod::pidTPCFullKa, aod::pidTOFFullKa, aod::pidTPCFullPi, aod::pidTOFFullPi, aod::McTrackLabels>>;
719790 using V0TrackCandidatesMC = soa::Filtered<soa::Join<aod::V0Datas, aod::V0TOFPIDs, aod::V0TOFNSigmas, aod::McV0Labels>>;
791+ using EventMCGenerated = soa::Join<aod::McCollisions, aod::McCentFT0Ms, aod::MultsExtraMC>;
720792 // zBeam direction in lab frame
721793
722794 template <typename T>
@@ -1200,14 +1272,25 @@ struct HigherMassResonances {
12001272 std::vector<bool > passKs;
12011273 ROOT::Math::PxPyPzMVector lResonanceGen1;
12021274 ROOT::Math::PxPyPzEVector lResonanceGen;
1275+ Service<o2::framework::O2DatabasePDG> pdgDB;
12031276
1204- void processGen (aod::McCollision const & mcCollision, aod::McParticles const & mcParticles, const soa::SmallGroups<EventCandidatesMC>& collisions)
1277+ // void processGen(aod::McCollision const& mcCollision, aod::McParticles const& mcParticles, const soa::SmallGroups<EventCandidatesMC>& collisions)
1278+ void processGen (EventMCGenerated::iterator const & mcCollision, aod::McParticles const & mcParticles, const soa::SmallGroups<EventCandidatesMC>& collisions)
12051279 {
12061280 // if (config.isMC == false) {
12071281 // return;
12081282 // }
12091283 hMChists.fill (HIST (" events_check" ), 0.5 );
12101284
1285+ bool isINELgt0true = false ;
1286+
1287+ if (pwglf::isINELgtNmc (mcParticles, 0 , pdgDB)) {
1288+ isINELgt0true = true ;
1289+ }
1290+ if (config.isINELgt0 && !isINELgt0true) {
1291+ return ;
1292+ }
1293+
12111294 std::vector<int64_t > selectedEvents (collisions.size ());
12121295 int nevts = 0 ;
12131296 multiplicityGen = -999.0 ;
@@ -1230,7 +1313,7 @@ struct HigherMassResonances {
12301313 continue ;
12311314 }
12321315
1233- selectedEvents[nevts++] = collision.mcCollision_as <aod::McCollisions >().globalIndex ();
1316+ selectedEvents[nevts++] = collision.mcCollision_as <EventMCGenerated >().globalIndex ();
12341317 }
12351318 selectedEvents.resize (nevts);
12361319 hMChists.fill (HIST (" events_check" ), 1.5 );
@@ -1240,6 +1323,8 @@ struct HigherMassResonances {
12401323 return ;
12411324 }
12421325 hMChists.fill (HIST (" events_check" ), 2.5 );
1326+ double genMultiplicity = mcCollision.centFT0M ();
1327+
12431328 for (const auto & mcParticle : mcParticles) {
12441329
12451330 if (std::abs (mcParticle.pdgCode ()) != config.pdgCodes [config.selectMCparticles ]) // f2(1525), f0(1710)
@@ -1290,6 +1375,7 @@ struct HigherMassResonances {
12901375 auto helicityGen1 = lResonanceGen1.Vect ().Dot (fourVecDauCM1.Vect ()) / (std::sqrt (fourVecDauCM1.Vect ().Mag2 ()) * std::sqrt (lResonanceGen1.Vect ().Mag2 ()));
12911376
12921377 hMChists.fill (HIST (" Genf1710" ), multiplicityGen, lResonanceGen.pt (), helicityGen);
1378+ hMChists.fill (HIST (" Genf1710Calib" ), genMultiplicity, lResonanceGen.pt (), helicityGen);
12931379 hMChists.fill (HIST (" Genf1710_mass" ), lResonanceGen.M ());
12941380 hMChists.fill (HIST (" GenRapidity" ), mcParticle.y ());
12951381 hMChists.fill (HIST (" GenEta" ), mcParticle.eta ());
@@ -1300,6 +1386,7 @@ struct HigherMassResonances {
13001386 }
13011387
13021388 hMChists.fill (HIST (" Genf17102" ), multiplicityGen, lResonanceGen1.pt (), helicityGen1);
1389+ hMChists.fill (HIST (" Genf1710Calib2" ), genMultiplicity, lResonanceGen1.pt (), helicityGen1);
13031390 hMChists.fill (HIST (" Genf1710_mass2" ), lResonanceGen1.M ());
13041391 hMChists.fill (HIST (" GenRapidity2" ), lResonanceGen1.Rapidity ());
13051392 hMChists.fill (HIST (" GenEta2" ), lResonanceGen1.Eta ());
@@ -1310,13 +1397,86 @@ struct HigherMassResonances {
13101397 }
13111398 PROCESS_SWITCH (HigherMassResonances, processGen, " Process Generated" , false );
13121399
1400+ void processEvtLossSigLossMC (EventMCGenerated::iterator const & mcCollision, aod::McParticles const & mcParticles, const soa::SmallGroups<EventCandidatesMC>& recCollisions)
1401+ {
1402+ // auto multiplicityRec = -1;
1403+ bool isSelectedEvent = false ;
1404+ // auto multiplicity1 = -999.;
1405+
1406+ for (const auto & RecCollision : recCollisions) {
1407+ if (!RecCollision.has_mcCollision ())
1408+ continue ;
1409+ if (!selectionEvent (RecCollision, false )) // don't fill event cut histogram
1410+ continue ;
1411+
1412+ // const auto& mcCollisionRec = RecCollision.mcCollision_as<EventMCGenerated>();
1413+ // multiplicityRec = mcCollisionRec.centFT0M();
1414+
1415+ // if (config.cSelectMultEstimator == kFT0M) {
1416+ // multiplicity1 = RecCollision.centFT0M();
1417+ // } else if (config.cSelectMultEstimator == kFT0A) {
1418+ // multiplicity1 = RecCollision.centFT0A();
1419+ // } else if (config.cSelectMultEstimator == kFT0C) {
1420+ // multiplicity1 = RecCollision.centFT0C();
1421+ // } else if (config.cSelectMultEstimator == kFV0A) {
1422+ // multiplicity1 = RecCollision.centFV0A();
1423+ // } else {
1424+ // multiplicity1 = RecCollision.centFT0M(); // default
1425+ // }
1426+
1427+ isSelectedEvent = true ;
1428+ }
1429+
1430+ bool isINELgt0true = false ;
1431+
1432+ if (pwglf::isINELgtNmc (mcParticles, 0 , pdgDB)) {
1433+ isINELgt0true = true ;
1434+ }
1435+ if (config.isINELgt0 && !isINELgt0true) {
1436+ return ;
1437+ }
1438+
1439+ if (std::abs (mcCollision.posZ ()) >= config.cutzvertex ) {
1440+ return ;
1441+ }
1442+
1443+ auto multiplicityGen = -1 ;
1444+ multiplicityGen = mcCollision.centFT0M ();
1445+ hMChists.fill (HIST (" MCcorrections/MultiplicityGen" ), multiplicityGen);
1446+
1447+ // Event loss
1448+ if (isSelectedEvent) {
1449+ hMChists.fill (HIST (" MCcorrections/MultiplicityRec" ), multiplicityGen);
1450+ }
1451+
1452+ // Generated MC
1453+ for (const auto & mcPart : mcParticles) {
1454+ if (std::abs (mcPart.y ()) >= config.rapidityMotherData || std::abs (mcPart.pdgCode ()) != config.pdgCodes [config.selectMCparticles ])
1455+ continue ;
1456+
1457+ hMChists.fill (HIST (" MCcorrections/hSignalLossDenominator" ), mcPart.pt (), multiplicityGen);
1458+ if (isSelectedEvent) {
1459+ hMChists.fill (HIST (" MCcorrections/hSignalLossNumerator" ), mcPart.pt (), multiplicityGen);
1460+ }
1461+ } // end loop on gen particles
1462+ }
1463+ PROCESS_SWITCH (HigherMassResonances, processEvtLossSigLossMC, " Process Signal Loss, Event Loss" , false );
1464+
13131465 int eventCounter = 0 ;
13141466 std::vector<int > gindex1, gindex2;
1315- void processRec (EventCandidatesMC::iterator const & collision, TrackCandidatesMC const &, V0TrackCandidatesMC const & V0s, aod::McParticles const &, aod::McCollisions const & /* mcCollisions */ )
1467+ void processRec (EventCandidatesMC::iterator const & collision, TrackCandidatesMC const &, V0TrackCandidatesMC const & V0s, aod::McParticles const &, EventMCGenerated const &)
13161468 {
13171469 // if (config.isMC == false) {
13181470 // return;
13191471 // }
1472+ if (!collision.has_mcCollision ()) {
1473+ return ;
1474+ }
1475+
1476+ double multiplicityRec = -1.0 ;
1477+ // multiplicityRec = collision.mcCollision_as<EventMCGenerated>().centFT0M();
1478+ const auto & mcCollisionRec = collision.mcCollision_as <EventMCGenerated>();
1479+ multiplicityRec = mcCollisionRec.centFT0M ();
13201480
13211481 auto multiplicity = -999.0 ;
13221482 if (config.cSelectMultEstimator == kFT0M ) {
@@ -1336,10 +1496,6 @@ struct HigherMassResonances {
13361496 }
13371497 hMChists.fill (HIST (" Rec_Multiplicity" ), multiplicity);
13381498
1339- if (!collision.has_mcCollision ()) {
1340- return ;
1341- }
1342-
13431499 hMChists.fill (HIST (" MC_mult_after_event_sel" ), multiplicity);
13441500 eventCounter++;
13451501
@@ -1472,6 +1628,7 @@ struct HigherMassResonances {
14721628 // }
14731629
14741630 hMChists.fill (HIST (" Recf1710_pt1" ), multiplicity, mothertrack1.pt (), mother1.M (), helicityRec2);
1631+ hMChists.fill (HIST (" Recf1710Calib_pt1" ), multiplicityRec, mothertrack1.pt (), mother1.M (), helicityRec2);
14751632 hMChists.fill (HIST (" RecRapidity" ), mothertrack1.y ());
14761633 hMChists.fill (HIST (" RecPhi" ), mothertrack1.phi ());
14771634 hMChists.fill (HIST (" RecEta" ), mothertrack1.eta ());
@@ -1481,6 +1638,7 @@ struct HigherMassResonances {
14811638 }
14821639
14831640 hMChists.fill (HIST (" Recf1710_pt2" ), multiplicity, mother.Pt (), mother.M (), helicityRec);
1641+ hMChists.fill (HIST (" Recf1710Calib_pt2" ), multiplicityRec, mother.Pt (), mother.M (), helicityRec);
14841642 hMChists.fill (HIST (" RecRapidity2" ), mother.Rapidity ());
14851643 hMChists.fill (HIST (" RecPhi2" ), mother.Phi ());
14861644 hMChists.fill (HIST (" RecEta2" ), mother.Eta ());
@@ -1508,7 +1666,7 @@ struct HigherMassResonances {
15081666 multiplicity = collision.centFT0M (); // default
15091667 }
15101668
1511- if (!selectionEvent (collision, true )) {
1669+ if (!selectionEventDerived (collision, true )) {
15121670 return ;
15131671 }
15141672
@@ -1617,7 +1775,7 @@ struct HigherMassResonances {
16171775 multiplicity = 0.0 ;
16181776 multiplicity = c1.centFT0M ();
16191777
1620- if (!selectionEvent (c1, false ) || !selectionEvent (c2, false )) {
1778+ if (!selectionEventDerived (c1, false ) || !selectionEventDerived (c2, false )) {
16211779 continue ;
16221780 }
16231781
0 commit comments