@@ -234,7 +234,8 @@ struct JetSpectraEseTask {
234234 static constexpr EventPlaneFiller PsiFillerEse = {true , false };
235235 static constexpr EventPlaneFiller PsiFillerFalse = {false , false };
236236 TRandom3* fRndm = new TRandom3(0 );
237- static constexpr int NumSubSmpl = 10 ;
237+ static constexpr int NumSubSmpl = 5 ;
238+ static constexpr int NumSavedRhoFitEvents = 5 ;
238239 std::array<std::shared_ptr<THnSparse>, NumSubSmpl> hSameSub;
239240
240241 void applySystematicPreset ()
@@ -359,6 +360,10 @@ struct JetSpectraEseTask {
359360 registry.get <TH1 >(HIST (" trackQA/hRhoTrackCounter" ))->GetXaxis ()->SetBinLabel (3 , " Tracks for rho(phi)" );
360361
361362 registry.add (" trackQA/hPhiPtsum" , " jet sumpt;sum p_{T};entries" , {HistType::kTH1F , {{40 , 0 ., o2::constants::math::TwoPI}}});
363+ for (int i = 0 ; i < NumSavedRhoFitEvents; ++i) {
364+ const auto eventName = fmt::format (" rhoPhiFitEvents/event_{}" , i);
365+ registry.add (eventName, " event;#varphi;#Sigma #it{p}_{T} (GeV/#it{c})" , HistType::kTH1F , {{1 , 0 ., o2::constants::math::TwoPI}});
366+ }
362367 registry.add (" eventQA/hfitPar0" , " " , {HistType::kTH2F , {{centAxis}, {fitAxis0}}});
363368 registry.add (" eventQA/hfitPar1" , " " , {HistType::kTH2F , {{centAxis}, {fitAxis13}}});
364369 registry.add (" eventQA/hfitPar2" , " " , {HistType::kTH2F , {{centAxis}, {fitAxis24}}});
@@ -659,7 +664,8 @@ struct JetSpectraEseTask {
659664 for (const auto & jet : jets) {
660665 if (!jetfindingutilities::isInEtaAcceptance (jet, jetEtaMin, jetEtaMax, trackEtaMin, trackEtaMax))
661666 continue ;
662- if (!isAcceptedJet<aod::JetTracks>(jet)) {
667+ // if (!isAcceptedJet<aod::JetTracks>(jet)) {
668+ if (!isAcceptedJet<soa::Join<aod::JetTracks, aod::JTrackPIs>>(jet)) {
663669 continue ;
664670 }
665671 auto vCorr = corr (collision, jet);
@@ -793,7 +799,7 @@ struct JetSpectraEseTask {
793799 for (const auto & jet : jets1) {
794800 if (!jetfindingutilities::isInEtaAcceptance (jet, jetEtaMin, jetEtaMax, trackEtaMin, trackEtaMax))
795801 continue ;
796- if (!isAcceptedJet<aod::JetTracks>(jet)) {
802+ if (!isAcceptedJet<soa::Join< aod::JetTracks, aod::JTrackPIs> >(jet)) {
797803 continue ;
798804 }
799805 auto vCorr = corr (c1, jet);
@@ -1198,13 +1204,13 @@ struct JetSpectraEseTask {
11981204 auto vec1{qVecNoESE<DetID::FT0A , P.hist >(vec)};
11991205 epMap[" FT0A" ] = computeEP (vec1, 0 , 2.0 );
12001206 ep3Map[" FT0A" ] = computeEP (vec1, 0 , 3.0 );
1207+ auto vec2{qVecNoESE<DetID::FT0C , false >(vec)};
1208+ epMap[" FT0C" ] = computeEP (vec2, 0 , 2.0 );
1209+ ep3Map[" FT0C" ] = computeEP (vec2, 0 , 3.0 );
1210+ epMap[" FV0A" ] = computeEP (qVecNoESE<DetID::FV0A , false >(vec), 0 , 2.0 );
1211+ epMap[" TPCpos" ] = computeEP (qVecNoESE<DetID::TPCpos, false >(vec), 1 , 2.0 );
1212+ epMap[" TPCneg" ] = computeEP (qVecNoESE<DetID::TPCneg, false >(vec), 1 , 2.0 );
12011213 if constexpr (P.psi ) {
1202- auto vec2{qVecNoESE<DetID::FT0C , false >(vec)};
1203- epMap[" FT0C" ] = computeEP (vec2, 0 , 2.0 );
1204- ep3Map[" FT0C" ] = computeEP (vec2, 0 , 3.0 );
1205- epMap[" FV0A" ] = computeEP (qVecNoESE<DetID::FV0A , false >(vec), 0 , 2.0 );
1206- epMap[" TPCpos" ] = computeEP (qVecNoESE<DetID::TPCpos, false >(vec), 1 , 2.0 );
1207- epMap[" TPCneg" ] = computeEP (qVecNoESE<DetID::TPCneg, false >(vec), 1 , 2.0 );
12081214 if constexpr (P.hist )
12091215 fillEP (/* std::make_index_sequence<5>{},*/ vec, epMap);
12101216 auto cosPsi = [](float psiX, float psiY) { return (static_cast <double >(psiX) == InvalidValue || static_cast <double >(psiY) == InvalidValue) ? InvalidValue : std::cos (2.0 * (psiX - psiY)); };
@@ -1312,6 +1318,24 @@ struct JetSpectraEseTask {
13121318 return true ;
13131319 }
13141320
1321+ std::shared_ptr<TH1 > getRhoPhiFitEvent (int eventIndex)
1322+ {
1323+ switch (eventIndex) {
1324+ case 0 :
1325+ return registry.get <TH1 >(HIST (" rhoPhiFitEvents/event_0" ));
1326+ case 1 :
1327+ return registry.get <TH1 >(HIST (" rhoPhiFitEvents/event_1" ));
1328+ case 2 :
1329+ return registry.get <TH1 >(HIST (" rhoPhiFitEvents/event_2" ));
1330+ case 3 :
1331+ return registry.get <TH1 >(HIST (" rhoPhiFitEvents/event_3" ));
1332+ case 4 :
1333+ return registry.get <TH1 >(HIST (" rhoPhiFitEvents/event_4" ));
1334+ default :
1335+ return nullptr ;
1336+ }
1337+ }
1338+
13151339 template <bool fillHist, typename Col, typename TTracks, typename Jets>
13161340 std::unique_ptr<TF1 > fitRho (const Col& col, const EventPlane& ep, TTracks const & tracks, Jets const & jets)
13171341 {
@@ -1360,7 +1384,7 @@ struct JetSpectraEseTask {
13601384 modulationFit->FixParameter (2 , (ep.psi2 < 0 ) ? RecoDecay::constrainAngle (ep.psi2 ) : ep.psi2 );
13611385 modulationFit->FixParameter (4 , (ep.psi3 < 0 ) ? RecoDecay::constrainAngle (ep.psi3 ) : ep.psi3 );
13621386
1363- hPhiPt->Fit (modulationFit.get (), " Q " , " " , 0 , o2::constants::math::TwoPI);
1387+ hPhiPt->Fit (modulationFit.get (), " QN " , " " , 0 , o2::constants::math::TwoPI);
13641388
13651389 if constexpr (fillHist) {
13661390 registry.fill (HIST (" eventQA/hfitPar0" ), col.centFT0M (), modulationFit->GetParameter (0 ));
@@ -1391,6 +1415,56 @@ struct JetSpectraEseTask {
13911415 if constexpr (fillHist) {
13921416 registry.fill (HIST (" eventQA/hPValueCentCDF" ), col.centFT0M (), cDF);
13931417 registry.fill (HIST (" eventQA/hCentChi2Ndf" ), col.centFT0M (), chi2 / (static_cast <float >(nDF)));
1418+
1419+ static int numSavedRhoFitEvents{0 };
1420+ if (numSavedRhoFitEvents < NumSavedRhoFitEvents) {
1421+ auto savedEvent = getRhoPhiFitEvent (numSavedRhoFitEvents);
1422+ savedEvent->SetBins (hPhiPt->GetNbinsX (), 0 ., o2::constants::math::TwoPI);
1423+ for (int bin = 0 ; bin <= hPhiPt->GetNbinsX () + 1 ; ++bin) {
1424+ savedEvent->SetBinContent (bin, hPhiPt->GetBinContent (bin));
1425+ savedEvent->SetBinError (bin, hPhiPt->GetBinError (bin));
1426+ }
1427+ savedEvent->SetEntries (hPhiPt->GetEntries ());
1428+
1429+ const double rho0 = modulationFit->GetParameter (0 );
1430+ const double v2 = modulationFit->GetParameter (1 );
1431+ const double psi2 = modulationFit->GetParameter (2 );
1432+ const double v3 = modulationFit->GetParameter (3 );
1433+ const double psi3 = modulationFit->GetParameter (4 );
1434+
1435+ auto rho0Function = std::make_unique<TF1 >(" rho_0" , " pol0" , 0 ., o2::constants::math::TwoPI, TF1 ::EAddToList::kNo );
1436+ rho0Function->SetParameter (0 , rho0);
1437+ rho0Function->SetParError (0 , modulationFit->GetParError (0 ));
1438+ rho0Function->SetParName (0 , " rho_0" );
1439+ rho0Function->SetLineStyle (2 );
1440+
1441+ auto rhoPhiFunction = std::unique_ptr<TF1 >(static_cast <TF1 *>(modulationFit->Clone (" rho_phi" )));
1442+ rhoPhiFunction->SetParNames (" rho_0" , " v_2" , " Psi_2" , " v_3" , " Psi_3" );
1443+ rhoPhiFunction->SetLineWidth (2 );
1444+ rhoPhiFunction->ResetBit (TF1 ::kNotDraw );
1445+
1446+ auto rhoV2Function = std::make_unique<TF1 >(" rho_v2" , " [0] * (1. + 2. * [1] * std::cos(2. * (x - [2])))" , 0 ., o2::constants::math::TwoPI, TF1 ::EAddToList::kNo );
1447+ rhoV2Function->SetParameters (rho0, v2, psi2);
1448+ rhoV2Function->SetParError (0 , modulationFit->GetParError (0 ));
1449+ rhoV2Function->SetParError (1 , modulationFit->GetParError (1 ));
1450+ rhoV2Function->SetParError (2 , modulationFit->GetParError (2 ));
1451+ rhoV2Function->SetParNames (" rho_0" , " v_2" , " Psi_2" );
1452+ rhoV2Function->SetLineStyle (7 );
1453+
1454+ auto rhoV3Function = std::make_unique<TF1 >(" rho_v3" , " [0] * (1. + 2. * [1] * std::cos(3. * (x - [2])))" , 0 ., o2::constants::math::TwoPI, TF1 ::EAddToList::kNo );
1455+ rhoV3Function->SetParameters (rho0, v3, psi3);
1456+ rhoV3Function->SetParError (0 , modulationFit->GetParError (0 ));
1457+ rhoV3Function->SetParError (1 , modulationFit->GetParError (3 ));
1458+ rhoV3Function->SetParError (2 , modulationFit->GetParError (4 ));
1459+ rhoV3Function->SetParNames (" rho_0" , " v_3" , " Psi_3" );
1460+ rhoV3Function->SetLineStyle (9 );
1461+
1462+ savedEvent->GetListOfFunctions ()->Add (rho0Function.release ());
1463+ savedEvent->GetListOfFunctions ()->Add (rhoPhiFunction.release ());
1464+ savedEvent->GetListOfFunctions ()->Add (rhoV2Function.release ());
1465+ savedEvent->GetListOfFunctions ()->Add (rhoV3Function.release ());
1466+ ++numSavedRhoFitEvents;
1467+ }
13941468 }
13951469
13961470 return modulationFit;
@@ -1622,7 +1696,7 @@ struct JetSpectraEseTask {
16221696 for (const auto & jet : jets) {
16231697 if (!jetfindingutilities::isInEtaAcceptance (jet, jetEtaMin, jetEtaMax, trackEtaMin, trackEtaMax))
16241698 continue ;
1625- if (!isAcceptedJet<aod::JetTracks>(jet)) {
1699+ if (!isAcceptedJet<soa::Join< aod::JetTracks, aod::JTrackPIs> >(jet)) {
16261700 continue ;
16271701 }
16281702 auto jetPtBkg = jet.pt () - (collision.rho () * jet.area ());
0 commit comments