Skip to content

Commit ee64753

Browse files
wenyaCernwenyaCern
committed
recover reduced task and edit full task
1 parent bf9f57a commit ee64753

File tree

5 files changed

+299
-104
lines changed

5 files changed

+299
-104
lines changed

‎PWGCF/DataModel/FemtoDerived.h‎

Lines changed: 5 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -35,6 +35,7 @@ enum CollisionBinning {
3535
kMult, //! Bin collision in number of charged tracks for mixing
3636
kMultPercentile, //! Bin collision in multiplicity percentile for mixing
3737
kMultMultPercentile, //! Bin collision in number of charged tracks and multiplicity percentile for mixing
38+
kMultPercentileqn, //! Bin collision in multiplicity percentile an qn value for mixing
3839
kNCollisionBinning
3940
};
4041

@@ -51,7 +52,8 @@ DECLARE_SOA_COLUMN(BitMaskTrackThree, bitmaskTrackThree, BitMaskType); //! Bit f
5152

5253
DECLARE_SOA_COLUMN(Downsample, downsample, bool); //! Flag for downsampling
5354

54-
DECLARE_SOA_COLUMN(QnBin, qnBin, int); //! qn bins for dividing events
55+
DECLARE_SOA_COLUMN(QnVal, qnVal, int); //! qn bins for dividing events
56+
DECLARE_SOA_COLUMN(Occupancy, occupancy, int); //! TPC occupancy
5557
} // namespace femtodreamcollision
5658

5759
DECLARE_SOA_TABLE_STAGED(FDCollisions, "FDCOLLISION",
@@ -64,8 +66,8 @@ DECLARE_SOA_TABLE_STAGED(FDCollisions, "FDCOLLISION",
6466
using FDCollision = FDCollisions::iterator;
6567

6668
DECLARE_SOA_TABLE(FDExtQnCollisions, "AOD", "FDEXTCOLLISION",
67-
femtodreamcollision::QnBin);
68-
using FDExtQnCollision = FDExtQnCollisions::iterator;
69+
femtodreamcollision::QnVal,
70+
femtodreamcollision::Occupancy);
6971

7072
DECLARE_SOA_TABLE(FDColMasks, "AOD", "FDCOLMASK",
7173
femtodreamcollision::BitMaskTrackOne,

‎PWGCF/FemtoDream/Core/femtoDreamCollisionSelection.h‎

Lines changed: 95 additions & 85 deletions
Original file line numberDiff line numberDiff line change
@@ -166,24 +166,16 @@ class FemtoDreamCollisionSelection
166166
/// \return whether or not the collisions fulfills the specified selections
167167
template <typename C>
168168
bool isPileUpCollisionPbPb(C const& col,
169-
bool noSameBunchPileup, bool isGoodZvtxFT0vsPV,
170-
bool isGoodITSLayersAll, bool noCollInRofStandard,
171-
bool noHighMultCollInPrevRof, bool noCollInTimeRangeStandard,
172-
bool /*isVertexITSTPC*/,
169+
bool noSameBunchPileup, bool isGoodITSLayersAll,
173170
int tpcOccupancyMin, int tpcOccupancyMax)
174171
{
175172
const auto occupancy = col.trackOccupancyInTimeRange();
176173
if ((occupancy < tpcOccupancyMin || occupancy > tpcOccupancyMax)) {
177174
return false;
178175
}
179176
if ((noSameBunchPileup && !col.selection_bit(aod::evsel::kNoSameBunchPileup))
180-
|| (isGoodZvtxFT0vsPV && !col.selection_bit(aod::evsel::kIsGoodZvtxFT0vsPV))
181177
|| (isGoodITSLayersAll && !col.selection_bit(aod::evsel::kIsGoodITSLayersAll))
182-
|| (noCollInRofStandard && !col.selection_bit(aod::evsel::kNoCollInRofStandard))
183-
|| (noHighMultCollInPrevRof && !col.selection_bit(aod::evsel::kNoHighMultCollInPrevRof))
184-
|| (noCollInTimeRangeStandard && !col.selection_bit(aod::evsel::kNoCollInTimeRangeStandard))
185-
// || (isVertexITSTPC && !col.selection_bit(aod::evsel::kIsVertexITSTPC))
186-
) {
178+
) {
187179
return false;
188180
}
189181

@@ -209,6 +201,23 @@ class FemtoDreamCollisionSelection
209201
}
210202
}
211203

204+
/// Initializes histograms for qn bin
205+
/// \param registry Histogram registry to be passed
206+
void initQn(HistogramRegistry* registry, int mumQnBins = 10)
207+
{
208+
mHistogramQn = registry;
209+
mHistogramQn->add("Event/centFT0CBefore", "; cent", kTH1F, {{10, 0, 100}});
210+
mHistogramQn->add("Event/centFT0CAfter", "; cent", kTH1F, {{10, 0, 100}});
211+
mHistogramQn->add("Event/centVsqn", "; cent; qn", kTH2F, {{10, 0, 100}, {100, 0, 1000}});
212+
mHistogramQn->add("Event/centVsqnVsSpher", "; cent; qn; Sphericity", kTH3F, {{10, 0, 100}, {100, 0, 1000}, {100, 0, 1}});
213+
mHistogramQn->add("Event/qnBin", "; qnBin; entries", kTH1F, {{20, 0, 20}});
214+
215+
for (int iqn(0); iqn < mumQnBins; ++iqn) {
216+
mHistogramQn->add(("Qn/mult_" + std::to_string(iqn)).c_str(), "; cent; c22", kTH1F, {{100, 0, 3500}});
217+
}
218+
return;
219+
}
220+
212221
/// Initializes histograms for the flow calculation
213222
/// \param registry Histogram registry to be passed
214223
void initFlow(HistogramRegistry* registry, bool doQnSeparation, int mumQnBins = 10, int binPt = 100, int binEta = 32)
@@ -222,19 +231,13 @@ class FemtoDreamCollisionSelection
222231
mImQ2thisEvt = new TH2D("ImQ2thisEvt", "", binPt, 0., 5., binEta, -0.8, 0.8);
223232
mMQthisEvt = new TH2D("MQthisEvt", "", binPt, 0., 5., binEta, -0.8, 0.8);
224233
mMQWeightthisEvt = new TH2D("MQWeightthisEvt", "", binPt, 0., 5., binEta, -0.8, 0.8);
225-
mHistogramQn = registry;
226-
mHistogramQn->add("Event/centFT0CBefore", "; cent", kTH1F, {{10, 0, 100}});
227-
mHistogramQn->add("Event/centFT0CAfter", "; cent", kTH1F, {{10, 0, 100}});
228-
mHistogramQn->add("Event/centVsqn", "; cent; qn", kTH2F, {{10, 0, 100}, {100, 0, 1000}});
229-
mHistogramQn->add("Event/centVsqnVsSpher", "; cent; qn; Sphericity", kTH3F, {{10, 0, 100}, {100, 0, 1000}, {100, 0, 1}});
230-
mHistogramQn->add("Event/qnBin", "; qnBin; entries", kTH1F, {{20, 0, 20}});
231234

235+
mHistogramQn = registry;
232236
mHistogramQn->add<TProfile>("Event/profileC22", "; cent; c22", kTProfile, {{10, 0, 100}});
233237
mHistogramQn->add<TProfile>("Event/profileC24", "; cent; c24", kTProfile, {{10, 0, 100}});
234238
if (doQnSeparation){
235239
for (int iqn(0); iqn < mumQnBins; ++iqn) {
236240
mHistogramQn->add<TProfile>(("Qn/profileC22_" + std::to_string(iqn)).c_str(), "; cent; c22", kTProfile, {{10, 0, 100}});
237-
mHistogramQn->add(("Qn/mult_" + std::to_string(iqn)).c_str(), "; cent; c22", kTH1F, {{100, 0, 3500}});
238241
}
239242
}
240243
return;
@@ -325,9 +328,9 @@ class FemtoDreamCollisionSelection
325328

326329
/// \todo to be implemented!
327330
/// \return the 1-d qn-vector separator to 2-d
328-
std::vector<std::vector<float>> getQnBinSeparator2D(std::vector<float> flat)
331+
std::vector<std::vector<float>> getQnBinSeparator2D(std::vector<float> flat, const int numQnBins = 10)
329332
{
330-
constexpr size_t nBins = 11;
333+
size_t nBins = numQnBins+1;
331334

332335
if (flat.empty() || flat.size() % nBins != 0) {
333336
LOGP(error, "ConfQnBinSeparator size = {} is not divisible by {}",
@@ -348,23 +351,21 @@ class FemtoDreamCollisionSelection
348351

349352
/// \todo to be implemented!
350353
/// Get the bin number of qn-vector(FT0C) of an event
351-
/// \tparam T type of the collision
352-
/// \param col Collision
353354
/// \param centBinWidth centrality bin width, example: per 1%, per 10% ...
354355
/// \return bin number of qn-vector of the event
355-
template <typename T>
356-
int myqnBin(T const& col, float centMax, std::vector<float> qnBinSeparator, float fSpher, float fMult, float centBinWidth = 1.f)
356+
int myqnBin(float centrality, float centMax, std::vector<float> qnBinSeparator, bool doFillHisto, float fSpher, float qn, const int numQnBins, float mult, float centBinWidth = 1.f)
357357
{
358-
auto twoDSeparator = getQnBinSeparator2D(qnBinSeparator);
358+
auto twoDSeparator = getQnBinSeparator2D(qnBinSeparator, numQnBins);
359359
if (twoDSeparator.empty() || twoDSeparator[0][0] == -999.) {
360360
LOGP(warning, "ConfQnBinSeparator not set, using default fallback!");
361361
return -999; // safe fallback
362362
}
363363

364-
mHistogramQn->fill(HIST("Event/centFT0CBefore"), col.centFT0C());
364+
if (doFillHisto)
365+
mHistogramQn->fill(HIST("Event/centFT0CBefore"), centrality);
366+
365367
int qnBin = -999;
366-
float qn = computeqnVec(col);
367-
int mycentBin = static_cast<int>(col.centFT0C() / centBinWidth);
368+
int mycentBin = static_cast<int>(centrality / centBinWidth);
368369
if (mycentBin >= static_cast<int>(centMax / centBinWidth))
369370
return qnBin;
370371

@@ -379,47 +380,52 @@ class FemtoDreamCollisionSelection
379380
continue;
380381
}
381382
}
382-
383-
mHistogramQn->fill(HIST("Event/centFT0CAfter"), col.centFT0C());
384-
mHistogramQn->fill(HIST("Event/centVsqn"), col.centFT0C(), qn);
385-
mHistogramQn->fill(HIST("Event/centVsqnVsSpher"), col.centFT0C(), qn, fSpher);
386-
mHistogramQn->fill(HIST("Event/qnBin"), qnBin);
387-
if (qnBin >= 0 && qnBin < 10){
388-
switch (qnBin) {
389-
case 0:
390-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("0"), fMult);
391-
break;
392-
case 1:
393-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("1"), fMult);
394-
break;
395-
case 2:
396-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("2"), fMult);
397-
break;
398-
case 3:
399-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("3"), fMult);
400-
break;
401-
case 4:
402-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("4"), fMult);
403-
break;
404-
case 5:
405-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("5"), fMult);
406-
break;
407-
case 6:
408-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("6"), fMult);
409-
break;
410-
case 7:
411-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("7"), fMult);
412-
break;
413-
case 8:
414-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("8"), fMult);
415-
break;
416-
case 9:
417-
mHistogramQn->fill(HIST("Qn/mult_") + HIST("9"), fMult);
418-
break;
419-
default:
420-
return qnBin; // invalid qn bin
421-
}
383+
384+
mQnBin = qnBin;
385+
386+
if (doFillHisto){
387+
mHistogramQn->fill(HIST("Event/centFT0CAfter"), centrality);
388+
mHistogramQn->fill(HIST("Event/centVsqn"), centrality, qn);
389+
mHistogramQn->fill(HIST("Event/centVsqnVsSpher"), centrality, qn, fSpher);
390+
mHistogramQn->fill(HIST("Event/qnBin"), qnBin);
391+
if (qnBin >= 0 && qnBin < numQnBins){
392+
switch (qnBin) {
393+
case 0:
394+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("0"), mult);
395+
break;
396+
case 1:
397+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("1"), mult);
398+
break;
399+
case 2:
400+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("2"), mult);
401+
break;
402+
case 3:
403+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("3"), mult);
404+
break;
405+
case 4:
406+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("4"), mult);
407+
break;
408+
case 5:
409+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("5"), mult);
410+
break;
411+
case 6:
412+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("6"), mult);
413+
break;
414+
case 7:
415+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("7"), mult);
416+
break;
417+
case 8:
418+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("8"), mult);
419+
break;
420+
case 9:
421+
mHistogramQn->fill(HIST("Qn/mult_") + HIST("9"), mult);
422+
break;
423+
default:
424+
return qnBin; // invalid qn bin
425+
}
426+
}
422427
}
428+
423429
return qnBin;
424430
}
425431

@@ -430,11 +436,11 @@ class FemtoDreamCollisionSelection
430436
/// \tparam T2 type of the tracks
431437
/// \param tracks All tracks
432438
template <typename T1, typename T2>
433-
void fillCumulants(T1 const& col, T2 const& tracks, float fHarmonic=2.f)
439+
bool fillCumulants(T1 const& col, T2 const& tracks, float fHarmonic=2.f)
434440
{
435441
int numOfTracks = col.numContrib();
436442
if (numOfTracks < 3)
437-
return ;
443+
return false;
438444

439445
mReQthisEvt->Reset();
440446
mImQthisEvt->Reset();
@@ -459,7 +465,7 @@ class FemtoDreamCollisionSelection
459465
mMQthisEvt ->Fill(pt, eta);
460466
mMQWeightthisEvt ->Fill(pt, eta, weight);
461467
}
462-
return;
468+
return true;
463469
}
464470

465471
/// \todo to be implemented!
@@ -468,9 +474,12 @@ class FemtoDreamCollisionSelection
468474
/// \param doQnSeparation to fill flow in divied qn bins
469475
/// \param qnBin should be <int> in 0-9
470476
/// \param fEtaGap eta gap for flow cumulant
471-
template <typename T>
472-
void doCumulants(T const& col, bool doQnSeparation = false, int qnBin = -999, int mumQnBinNum = 10, float fEtaGap = 0.3f, int binPt = 100, int binEta = 32)
477+
template <typename T1, typename T2>
478+
void doCumulants(T1 const& col, T2 const& tracks, float centrality, bool doQnSeparation = false, int numQnBins = 10, float fEtaGap = 0.3f, int binPt = 100, int binEta = 32)
473479
{
480+
if (!fillCumulants(col, tracks))
481+
return;
482+
474483
if (mMQthisEvt->Integral(1, binPt, 1, binEta) < 2)
475484
return;
476485

@@ -495,38 +504,38 @@ class FemtoDreamCollisionSelection
495504
TComplex negEtaQ = TComplex(negEtaRe, negEtaIm);
496505
TComplex negEtaQStar = TComplex::Conjugate(negEtaQ);
497506

498-
mHistogramQn->get<TProfile>(HIST("Event/profileC22"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
499-
if (doQnSeparation && qnBin >= 0 && qnBin < mumQnBinNum){
500-
switch (qnBin) {
507+
mHistogramQn->get<TProfile>(HIST("Event/profileC22"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
508+
if (doQnSeparation && mQnBin >= 0 && mQnBin < numQnBins){
509+
switch (mQnBin) {
501510
case 0:
502-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("0"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
511+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("0"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
503512
break;
504513
case 1:
505-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("1"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
514+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("1"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
506515
break;
507516
case 2:
508-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("2"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
517+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("2"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
509518
break;
510519
case 3:
511-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("3"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
520+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("3"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
512521
break;
513522
case 4:
514-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("4"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
523+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("4"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
515524
break;
516525
case 5:
517-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("5"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
526+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("5"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
518527
break;
519528
case 6:
520-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("6"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
529+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("6"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
521530
break;
522531
case 7:
523-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("7"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
532+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("7"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
524533
break;
525534
case 8:
526-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("8"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
535+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("8"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
527536
break;
528537
case 9:
529-
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("9"))->Fill(col.centFT0C(), (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
538+
mHistogramQn->get<TProfile>(HIST("Qn/profileC22_") + HIST("9"))->Fill(centrality, (negEtaQ*posEtaQStar).Re()/(negEtaMQ*posEtaMQ), (negEtaMQ*posEtaMQ));
530539
break;
531540
default:
532541
return; // invalid qn bin
@@ -547,6 +556,7 @@ class FemtoDreamCollisionSelection
547556
float mZvtxMax = 999.f; ///< Maximal deviation from nominal z-vertex (cm)
548557
float mMinSphericity = 0.f;
549558
float mSphericityPtmin = 0.f;
559+
int mQnBin = -999;
550560
HistogramRegistry* mHistogramQn = nullptr; ///< For flow cumulant output
551561
TH2D* mReQthisEvt = nullptr; ///< For flow cumulant in an event
552562
TH2D* mImQthisEvt = nullptr; ///< For flow cumulant in an event

‎PWGCF/FemtoDream/Core/femtoDreamContainer.h‎

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -366,7 +366,7 @@ class FemtoDreamContainer
366366
}
367367

368368
template <bool isMC, typename T1, typename T2>
369-
void setPair_qn(T1 const& part1, T2 const& part2, const float multPercentile, const int myQnBin)
369+
void setPair_qn(T1 const& part1, T2 const& part2, const float multPercentile, const int myQnBin, const int numQnBins = 10)
370370
{
371371
float femtoObs, femtoObsMC;
372372
// Calculate femto observable and the mT with reconstructed information
@@ -379,8 +379,7 @@ class FemtoDreamContainer
379379
}
380380
}
381381
const float mT = FemtoDreamMath::getmT(part1, mMassOne, part2, mMassTwo);
382-
const int numQnBins = 10;
383-
382+
384383
if (mHistogramRegistry) {
385384
setPair_qn_base<o2::aod::femtodreamMCparticle::MCType::kRecon>(femtoObs, mT, multPercentile, myQnBin, numQnBins);
386385

0 commit comments

Comments
 (0)