@@ -96,6 +96,11 @@ struct HfTaskD0 {
9696 Configurable<std::string> ccdbUrl{" ccdbUrl" , " http://alice-ccdb.cern.ch" , " url of the ccdb repository" };
9797 Configurable<std::string> irSource{" irSource" , " ZNC hadronic" , " Estimator of the interaction rate (Recommended: pp --> T0VTX, Pb-Pb --> ZNC hadronic)" };
9898
99+ // UPC gap determination thresholds
100+ Configurable<float > upcFT0AThreshold{" upcFT0AThreshold" , 100 .0f , " FT0-A amplitude threshold for UPC gap determination (a.u.)" };
101+ Configurable<float > upcFT0CThreshold{" upcFT0CThreshold" , 50 .0f , " FT0-C amplitude threshold for UPC gap determination (a.u.)" };
102+ Configurable<float > upcZDCThreshold{" upcZDCThreshold" , 1 .0f , " ZDC energy threshold for UPC gap determination (a.u.)" };
103+
99104 HfEventSelection hfEvSel; // event selection and monitoring
100105
101106 ctpRateFetcher mRateFetcher ;
@@ -157,35 +162,6 @@ struct HfTaskD0 {
157162 ConfigurableAxis thnConfigAxisFT0A{" thnConfigAxisFT0A" , {1001 , -1.5 , 999.5 }, " axis for FT0-A amplitude (a.u.)" };
158163 ConfigurableAxis thnConfigAxisFT0C{" thnConfigAxisFT0C" , {1001 , -1.5 , 999.5 }, " axis for FT0-C amplitude (a.u.)" };
159164
160- // UPC gap determination thresholds
161- Configurable<float > upcFT0AThreshold{" upcFT0AThreshold" , hf_upc::defaults::FT0AThreshold, " FT0-A amplitude threshold for UPC gap determination (a.u.)" };
162- Configurable<float > upcFT0CThreshold{" upcFT0CThreshold" , hf_upc::defaults::FT0CThreshold, " FT0-C amplitude threshold for UPC gap determination (a.u.)" };
163- Configurable<float > upcZDCThreshold{" upcZDCThreshold" , hf_upc::defaults::ZDCThreshold, " ZDC energy threshold for UPC gap determination (a.u.)" };
164-
165- HfEventSelection hfEvSel; // event selection and monitoring
166-
167- HfHelper hfHelper;
168- ctpRateFetcher mRateFetcher ;
169-
170- SliceCache cache;
171- Service<o2::ccdb::BasicCCDBManager> ccdb;
172-
173- using D0Candidates = soa::Join<aod::HfCand2Prong, aod::HfSelD0>;
174- using D0CandidatesMc = soa::Join<D0Candidates, aod::HfCand2ProngMcRec>;
175- using D0CandidatesKF = soa::Join<D0Candidates, aod::HfCand2ProngKF>;
176- using D0CandidatesMcKF = soa::Join<D0CandidatesKF, aod::HfCand2ProngMcRec>;
177-
178- using D0CandidatesMl = soa::Join<D0Candidates, aod::HfMlD0>;
179- using D0CandidatesMlMc = soa::Join<D0CandidatesMl, aod::HfCand2ProngMcRec>;
180- using D0CandidatesMlKF = soa::Join<D0CandidatesMl, aod::HfCand2ProngKF>;
181- using D0CandidatesMlMcKF = soa::Join<D0CandidatesMlKF, aod::HfCand2ProngMcRec>;
182-
183- using Collisions = soa::Join<aod::Collisions, aod::EvSels>;
184- using CollisionsCent = soa::Join<aod::Collisions, aod::EvSels, aod::CentFT0Ms, aod::CentFT0Cs>;
185- using CollisionsWithMcLabels = soa::Join<aod::Collisions, aod::McCollisionLabels, aod::EvSels>;
186- using CollisionsWithMcLabelsCent = soa::Join<aod::Collisions, aod::McCollisionLabels, aod::EvSels, aod::CentFT0Ms, aod::CentFT0Cs>;
187- using TracksSelQuality = soa::Join<aod::TracksExtra, aod::TracksWMc>;
188-
189165 HistogramRegistry registry{
190166 " registry" ,
191167 {{" hPtCand" , " 2-prong candidates;candidate #it{p}_{T} (GeV/#it{c});entries" , {HistType::kTH1F , {{360 , 0 ., 36 .}}}},
@@ -262,8 +238,6 @@ struct HfTaskD0 {
262238 {" hMassReflBkgD0bar" , " 2-prong candidates (matched);#it{m}_{inv} (GeV/#it{c}^{2}); #it{p}_{T}; #it{y}" , {HistType::kTH3F , {{120 , 1.5848 , 2.1848 }, {150 , 0 ., 30 .}, {20 , -5 ., 5 .}}}},
263239 {" hMassSigBkgD0bar" , " 2-prong candidates (not checked);#it{m}_{inv} (GeV/#it{c}^{2}); #it{p}_{T}; #it{y}" , {HistType::kTH3F , {{120 , 1.5848 , 2.1848 }, {150 , 0 ., 30 .}, {20 , -5 ., 5 .}}}}}};
264240
265- HistogramRegistry qaRegistry{" QAHistos" , {}, OutputObjHandlingPolicy::AnalysisObject};
266-
267241 void init (InitContext&)
268242 {
269243 std::array<bool , 13 > doprocess{doprocessDataWithDCAFitterN, doprocessDataWithDCAFitterNCent, doprocessDataWithKFParticle, doprocessMcWithDCAFitterN, doprocessMcWithDCAFitterNCent, doprocessMcWithKFParticle, doprocessDataWithDCAFitterNMl, doprocessDataWithDCAFitterNMlCent, doprocessDataWithKFParticleMl, doprocessMcWithDCAFitterNMl, doprocessMcWithDCAFitterNMlCent, doprocessMcWithKFParticleMl, doprocessDataWithDCAFitterNMlWithUpc};
@@ -395,14 +369,14 @@ struct HfTaskD0 {
395369 registry.get <THnSparse>(HIST (" hMassVsPtVsPtBVsYVsOriginVsD0Type" ))->Sumw2 ();
396370 }
397371
398- qaRegistry .add (" Data/fitInfo/ampFT0A_vs_ampFT0C" , " FT0-A vs FT0-C amplitude;FT0-A amplitude (a.u.);FT0-C amplitude (a.u.)" , {HistType::kTH2F , {{2500 , 0 ., 250 }, {2500 , 0 ., 250 }}});
399- qaRegistry .add (" Data/zdc/energyZNA_vs_energyZNC" , " ZNA vs ZNC common energy;E_{ZNA}^{common} (a.u.);E_{ZNC}^{common} (a.u.)" , {HistType::kTH2F , {{200 , 0 ., 20 }, {200 , 0 ., 20 }}});
400- qaRegistry .add (" Data/hUpcGapAfterSelection" , " UPC gap type after selection;Gap side;Counts" , {HistType::kTH1F , {{3 , -0.5 , 2.5 }}});
401- qaRegistry .get <TH1>(HIST (" Data/hUpcGapAfterSelection" ))->GetXaxis ()->SetBinLabel (static_cast <int >(GapType::GapA) + 1 , " A" );
402- qaRegistry .get <TH1>(HIST (" Data/hUpcGapAfterSelection" ))->GetXaxis ()->SetBinLabel (static_cast <int >(GapType::GapC) + 1 , " C" );
403- qaRegistry .get <TH1>(HIST (" Data/hUpcGapAfterSelection" ))->GetXaxis ()->SetBinLabel (static_cast <int >(GapType::DoubleGap) + 1 , " Double" );
372+ registry .add (" Data/fitInfo/ampFT0A_vs_ampFT0C" , " FT0-A vs FT0-C amplitude;FT0-A amplitude (a.u.);FT0-C amplitude (a.u.)" , {HistType::kTH2F , {{2500 , 0 ., 250 }, {2500 , 0 ., 250 }}});
373+ registry .add (" Data/zdc/energyZNA_vs_energyZNC" , " ZNA vs ZNC common energy;E_{ZNA}^{common} (a.u.);E_{ZNC}^{common} (a.u.)" , {HistType::kTH2F , {{200 , 0 ., 20 }, {200 , 0 ., 20 }}});
374+ registry .add (" Data/hUpcGapAfterSelection" , " UPC gap type after selection;Gap side;Counts" , {HistType::kTH1F , {{3 , -0.5 , 2.5 }}});
375+ registry .get <TH1>(HIST (" Data/hUpcGapAfterSelection" ))->GetXaxis ()->SetBinLabel (static_cast <int >(GapType::GapA) + 1 , " A" );
376+ registry .get <TH1>(HIST (" Data/hUpcGapAfterSelection" ))->GetXaxis ()->SetBinLabel (static_cast <int >(GapType::GapC) + 1 , " C" );
377+ registry .get <TH1>(HIST (" Data/hUpcGapAfterSelection" ))->GetXaxis ()->SetBinLabel (static_cast <int >(GapType::DoubleGap) + 1 , " Double" );
404378
405- hfEvSel.addHistograms (qaRegistry );
379+ hfEvSel.addHistograms (registry );
406380
407381 ccdb->setURL (ccdbUrl);
408382 ccdb->setCaching (true );
@@ -589,7 +563,7 @@ struct HfTaskD0 {
589563 }
590564 }
591565
592- template <bool fillMl , typename CollType, typename CandType, typename BCsType>
566+ template <bool FillMl , typename CollType, typename CandType, typename BCsType>
593567 void runAnalysisPerCollisionDataWithUpc (CollType const & collisions,
594568 CandType const & candidates,
595569 BCsType const & bcs,
@@ -598,40 +572,45 @@ struct HfTaskD0 {
598572 aod::FDDs const & fdds)
599573 {
600574 for (const auto & collision : collisions) {
601- uint32_t rejectionMask{0 };
602575 float centrality{-1 .f };
603- rejectionMask = hfEvSel.getHfCollisionRejectionMaskWithUpc <true , CentralityEstimator::None, BCsType>(collision, centrality, ccdb, qaRegistry , bcs);
576+ const auto rejectionMask = hfEvSel.getHfCollisionRejectionMaskWithUpc <true , CentralityEstimator::None, BCsType>(collision, centrality, ccdb, registry , bcs);
604577 if (rejectionMask != 0 ) {
605578 continue ;
606579 }
607- auto bc = collision.template bc_as <BCsType>();
580+ const auto & bc = collision.template bc_as <BCsType>();
608581 upchelpers::FITInfo fitInfo{};
609582 udhelpers::getFITinfo (fitInfo, bc, bcs, ft0s, fv0as, fdds);
610583
611584 GapType gap = GapType::DoubleGap;
612585 if (bc.has_zdc ()) {
613- auto zdc = bc.zdc ();
614- qaRegistry .fill (HIST (" Data/fitInfo/ampFT0A_vs_ampFT0C" ), fitInfo.ampFT0A , fitInfo.ampFT0C );
615- qaRegistry .fill (HIST (" Data/zdc/energyZNA_vs_energyZNC" ), zdc.energyCommonZNA (), zdc.energyCommonZNC ());
586+ const auto & zdc = bc.zdc ();
587+ registry .fill (HIST (" Data/fitInfo/ampFT0A_vs_ampFT0C" ), fitInfo.ampFT0A , fitInfo.ampFT0C );
588+ registry .fill (HIST (" Data/zdc/energyZNA_vs_energyZNC" ), zdc.energyCommonZNA (), zdc.energyCommonZNC ());
616589 gap = hf_upc::determineGapType (fitInfo.ampFT0A , fitInfo.ampFT0C , zdc.energyCommonZNA (), zdc.energyCommonZNC (),
617590 upcFT0AThreshold, upcFT0CThreshold, upcZDCThreshold);
618- qaRegistry .fill (HIST (" Data/hUpcGapAfterSelection" ), hf_upc::gapTypeToInt (gap));
591+ registry .fill (HIST (" Data/hUpcGapAfterSelection" ), hf_upc::gapTypeToInt (gap));
619592 }
620593 if (hf_upc::isSingleSidedGap (gap)) {
621594 int const gapTypeInt = hf_upc::gapTypeToInt (gap);
622595 const auto thisCollId = collision.globalIndex ();
623596 const auto & groupedD0Candidates = candidates.sliceBy (candD0PerCollision, thisCollId);
624- float cent{centrality};
597+
598+ // Calculate occupancy and interaction rate if needed
625599 float occ{-1 .f };
600+ float ir{-1 .f };
601+ if (storeOccupancyAndIR && occEstimator != OccupancyEstimator::None) {
602+ occ = o2::hf_occupancy::getOccupancyColl (collision, occEstimator);
603+ ir = mRateFetcher .fetch (ccdb.service , bc.timestamp (), bc.runNumber (), irSource, true ) * 1 .e -3 ; // kHz
604+ }
626605
627606 for (const auto & candidate : groupedD0Candidates) {
628607 if (yCandRecoMax >= 0 . && std::abs (HfHelper::yD0 (candidate)) > yCandRecoMax) {
629608 continue ;
630609 }
631610
632- float massD0 = HfHelper::invMassD0ToPiK (candidate);
633- float massD0bar = HfHelper::invMassD0barToKPi (candidate);
634- auto ptCandidate = candidate.pt ();
611+ const float massD0 = HfHelper::invMassD0ToPiK (candidate);
612+ const float massD0bar = HfHelper::invMassD0barToKPi (candidate);
613+ const auto ptCandidate = candidate.pt ();
635614
636615 if (candidate.isSelD0 () >= selectionFlagD0) {
637616 registry.fill (HIST (" hMass" ), massD0, ptCandidate);
@@ -644,56 +623,53 @@ struct HfTaskD0 {
644623 registry.fill (HIST (" hMassVsPhi" ), massD0bar, ptCandidate, candidate.phi ());
645624 }
646625
647- // Fill THnSparse with structure matching taskDplus: [mass, pt, mlScores, cent, occ, gapType, FT0A, FT0C]
648- if constexpr (fillMl) {
649- auto fillTHnData = [&](float mass, int d0Type) {
650- std::vector<double > valuesToFill{static_cast <double >(mass), static_cast <double >(ptCandidate),
651- candidate.mlProbD0 ()[0 ], candidate.mlProbD0 ()[1 ], candidate.mlProbD0 ()[2 ],
652- static_cast <double >(HfHelper::yD0 (candidate)), static_cast <double >(d0Type)};
653- if (storeCentrality) {
654- valuesToFill.push_back (cent);
655- }
656- if (storeOccupancyAndIR) {
657- valuesToFill.push_back (occ);
658- }
659- valuesToFill.push_back (static_cast <double >(gapTypeInt));
660- valuesToFill.push_back (static_cast <double >(fitInfo.ampFT0A ));
661- valuesToFill.push_back (static_cast <double >(fitInfo.ampFT0C ));
662- registry.get <THnSparse>(HIST (" hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type" ))->Fill (valuesToFill.data ());
663- };
664-
665- if (candidate.isSelD0 () >= selectionFlagD0) {
666- fillTHnData (massD0, SigD0);
667- fillTHnData (massD0, candidate.isSelD0bar () ? ReflectedD0 : PureSigD0);
626+ // Fill THnSparse with structure matching histogram axes: [mass, pt, (mlScores if FillMl), rapidity, d0Type, (cent if storeCentrality), (occ, ir if storeOccupancyAndIR), gapType, FT0A, FT0C]
627+ auto fillTHnData = [&](float mass, int d0Type) {
628+ // Pre-calculate vector size to avoid reallocations
629+ constexpr int nAxesBase = 7 ; // mass, pt, rapidity, d0Type, gapType, FT0A, FT0C
630+ constexpr int nAxesMl = FillMl ? 3 : 0 ; // 3 ML scores if FillMl
631+ int const nAxesCent = storeCentrality ? 1 : 0 ; // centrality if storeCentrality
632+ int const nAxesOccIR = storeOccupancyAndIR ? 2 : 0 ; // occupancy and IR if storeOccupancyAndIR
633+ int const nAxesTotal = nAxesBase + nAxesMl + nAxesCent + nAxesOccIR;
634+
635+ std::vector<double > valuesToFill;
636+ valuesToFill.reserve (nAxesTotal);
637+
638+ // Fill values in order matching histogram axes
639+ valuesToFill.push_back (static_cast <double >(mass));
640+ valuesToFill.push_back (static_cast <double >(ptCandidate));
641+ if constexpr (FillMl) {
642+ valuesToFill.push_back (candidate.mlProbD0 ()[0 ]);
643+ valuesToFill.push_back (candidate.mlProbD0 ()[1 ]);
644+ valuesToFill.push_back (candidate.mlProbD0 ()[2 ]);
668645 }
669- if (candidate.isSelD0bar () >= selectionFlagD0bar) {
670- fillTHnData (massD0bar, SigD0bar);
671- fillTHnData (massD0bar, candidate.isSelD0 () ? ReflectedD0bar : PureSigD0bar);
646+ valuesToFill.push_back (static_cast <double >(HfHelper::yD0 (candidate)));
647+ valuesToFill.push_back (static_cast <double >(d0Type));
648+ if (storeCentrality) {
649+ valuesToFill.push_back (centrality);
672650 }
673- } else {
674- auto fillTHnData = [&](float mass, int d0Type) {
675- std::vector<double > valuesToFill{static_cast <double >(mass), static_cast <double >(ptCandidate),
676- static_cast <double >(HfHelper::yD0 (candidate)), static_cast <double >(d0Type)};
677- if (storeCentrality) {
678- valuesToFill.push_back (cent);
679- }
680- if (storeOccupancyAndIR) {
681- valuesToFill.push_back (occ);
682- }
683- valuesToFill.push_back (static_cast <double >(gapTypeInt));
684- valuesToFill.push_back (static_cast <double >(fitInfo.ampFT0A ));
685- valuesToFill.push_back (static_cast <double >(fitInfo.ampFT0C ));
686- registry.get <THnSparse>(HIST (" hMassVsPtVsPtBVsYVsOriginVsD0Type" ))->Fill (valuesToFill.data ());
687- };
688-
689- if (candidate.isSelD0 () >= selectionFlagD0) {
690- fillTHnData (massD0, SigD0);
691- fillTHnData (massD0, candidate.isSelD0bar () ? ReflectedD0 : PureSigD0);
651+ if (storeOccupancyAndIR) {
652+ valuesToFill.push_back (occ);
653+ valuesToFill.push_back (ir);
692654 }
693- if (candidate.isSelD0bar () >= selectionFlagD0bar) {
694- fillTHnData (massD0bar, SigD0bar);
695- fillTHnData (massD0bar, candidate.isSelD0 () ? ReflectedD0bar : PureSigD0bar);
655+ valuesToFill.push_back (static_cast <double >(gapTypeInt));
656+ valuesToFill.push_back (static_cast <double >(fitInfo.ampFT0A ));
657+ valuesToFill.push_back (static_cast <double >(fitInfo.ampFT0C ));
658+
659+ if constexpr (FillMl) {
660+ registry.get <THnSparse>(HIST (" hBdtScoreVsMassVsPtVsPtBVsYVsOriginVsD0Type" ))->Fill (valuesToFill.data ());
661+ } else {
662+ registry.get <THnSparse>(HIST (" hMassVsPtVsPtBVsYVsOriginVsD0Type" ))->Fill (valuesToFill.data ());
696663 }
664+ };
665+
666+ if (candidate.isSelD0 () >= selectionFlagD0) {
667+ fillTHnData (massD0, SigD0);
668+ fillTHnData (massD0, candidate.isSelD0bar () ? ReflectedD0 : PureSigD0);
669+ }
670+ if (candidate.isSelD0bar () >= selectionFlagD0bar) {
671+ fillTHnData (massD0bar, SigD0bar);
672+ fillTHnData (massD0bar, candidate.isSelD0 () ? ReflectedD0bar : PureSigD0bar);
697673 }
698674 }
699675 }
0 commit comments