Skip to content

Commit dbc63b8

Browse files
committed
Add 3rd order ESE calculation and fix CI warnings
1 parent d8c5b05 commit dbc63b8

1 file changed

Lines changed: 82 additions & 46 deletions

File tree

PWGCF/Flow/Tasks/flowEseTask.cxx

Lines changed: 82 additions & 46 deletions
Original file line numberDiff line numberDiff line change
@@ -46,6 +46,8 @@
4646
#include <Math/Vector4D.h> // IWYU pragma: keep (do not replace with Math/Vector4Dfwd.h)
4747
#include <Math/Vector4Dfwd.h>
4848
#include <TF1.h>
49+
#include <TH2.h>
50+
#include <THnSparse.h>
4951
#include <TObject.h>
5052
#include <TProfile2D.h>
5153
#include <TProfile3D.h>
@@ -303,7 +305,7 @@ struct FlowEseTask {
303305
double getEseQ(TCollision const& collision)
304306
{
305307
const int harmonicIndex = cfgEseHarmonic.value - SecondHarmonic;
306-
if (collision.qvecFT0CReVec().size() <= static_cast<size_t>(harmonicIndex) || collision.qvecFT0CImVec().size() <= static_cast<size_t>(harmonicIndex)) {
308+
if (collision.qvecFT0CReVec().size() <= static_cast<std::size_t>(harmonicIndex) || collision.qvecFT0CImVec().size() <= static_cast<std::size_t>(harmonicIndex)) {
307309
LOGF(fatal, "FT0C Q-vector table does not contain harmonic %d", cfgEseHarmonic.value);
308310
}
309311
const double qx = collision.qvecFT0CReVec()[harmonicIndex];
@@ -585,52 +587,70 @@ struct FlowEseTask {
585587
template <typename TCollision, typename V0>
586588
bool selectionV0(TCollision const& collision, V0 const& candidate, int lambdaTag)
587589
{
588-
if (candidate.v0radius() < cfgv0radiusMin)
590+
if (candidate.v0radius() < cfgv0radiusMin) {
589591
return false;
592+
}
590593
if (lambdaTag) {
591-
if (std::abs(candidate.dcapostopv()) < cfgDCAPrToPVMin)
594+
if (std::abs(candidate.dcapostopv()) < cfgDCAPrToPVMin) {
592595
return false;
593-
if (std::abs(candidate.dcanegtopv()) < cfgDCAPiToPVMin)
596+
}
597+
if (std::abs(candidate.dcanegtopv()) < cfgDCAPiToPVMin) {
594598
return false;
599+
}
595600
} else {
596-
if (std::abs(candidate.dcapostopv()) < cfgDCAPiToPVMin)
601+
if (std::abs(candidate.dcapostopv()) < cfgDCAPiToPVMin) {
597602
return false;
598-
if (std::abs(candidate.dcanegtopv()) < cfgDCAPrToPVMin)
603+
}
604+
if (std::abs(candidate.dcanegtopv()) < cfgDCAPrToPVMin) {
599605
return false;
606+
}
600607
}
601-
if (candidate.v0cosPA() < cfgv0CosPA)
608+
if (candidate.v0cosPA() < cfgv0CosPA) {
602609
return false;
603-
if (std::abs(candidate.dcaV0daughters()) > cfgDCAV0Dau)
610+
}
611+
if (std::abs(candidate.dcaV0daughters()) > cfgDCAV0Dau) {
604612
return false;
605-
if (candidate.pt() < cfgV0PtMin)
613+
}
614+
if (candidate.pt() < cfgV0PtMin) {
606615
return false;
607-
if (candidate.yLambda() < cfgV0EtaMin)
616+
}
617+
if (candidate.yLambda() < cfgV0EtaMin) {
608618
return false;
609-
if (candidate.yLambda() > cfgV0EtaMax)
619+
}
620+
if (candidate.yLambda() > cfgV0EtaMax) {
610621
return false;
611-
if (candidate.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * massLambda > cfgV0LifeTime)
622+
}
623+
if (candidate.distovertotmom(collision.posX(), collision.posY(), collision.posZ()) * massLambda > cfgV0LifeTime) {
612624
return false;
625+
}
613626

614627
return true;
615628
}
616629

617630
template <typename T>
618631
bool isSelectedV0Daughter(T const& track, int pid) // pid 0: proton, pid 1: pion
619632
{
620-
if (track.tpcNClsFound() < cfgDaughTPCnclsMin)
633+
if (track.tpcNClsFound() < cfgDaughTPCnclsMin) {
621634
return false;
622-
if (pid == 0 && std::abs(track.tpcNSigmaPr()) > cfgDaughPIDCutsTPCPr)
635+
}
636+
if (pid == 0 && std::abs(track.tpcNSigmaPr()) > cfgDaughPIDCutsTPCPr) {
623637
return false;
624-
if (pid == 1 && std::abs(track.tpcNSigmaPi()) > cfgDaughPIDCutsTPCPi)
638+
}
639+
if (pid == 1 && std::abs(track.tpcNSigmaPi()) > cfgDaughPIDCutsTPCPi) {
625640
return false;
626-
if (track.eta() > cfgDaughEtaMax)
641+
}
642+
if (track.eta() > cfgDaughEtaMax) {
627643
return false;
628-
if (track.eta() < cfgDaughEtaMin)
644+
}
645+
if (track.eta() < cfgDaughEtaMin) {
629646
return false;
630-
if (pid == 0 && track.pt() < cfgDaughPrPt)
647+
}
648+
if (pid == 0 && track.pt() < cfgDaughPrPt) {
631649
return false;
632-
if (pid == 1 && track.pt() < cfgDaughPiPt)
650+
}
651+
if (pid == 1 && track.pt() < cfgDaughPiPt) {
633652
return false;
653+
}
634654

635655
return true;
636656
}
@@ -652,24 +672,33 @@ struct FlowEseTask {
652672
template <typename TrackType>
653673
bool selectionTrack(TrackType const& track)
654674
{
655-
if (track.pt() < cfgMinPt)
675+
if (track.pt() < cfgMinPt) {
656676
return false;
657-
if (std::abs(track.eta()) > cfgMaxEta)
677+
}
678+
if (std::abs(track.eta()) > cfgMaxEta) {
658679
return false;
659-
if (!track.passedITSNCls())
680+
}
681+
if (!track.passedITSNCls()) {
660682
return false;
661-
if (!track.passedITSChi2NDF())
683+
}
684+
if (!track.passedITSChi2NDF()) {
662685
return false;
663-
if (!track.passedITSHits())
686+
}
687+
if (!track.passedITSHits()) {
664688
return false;
665-
if (!track.passedTPCCrossedRowsOverNCls())
689+
}
690+
if (!track.passedTPCCrossedRowsOverNCls()) {
666691
return false;
667-
if (!track.passedTPCChi2NDF())
692+
}
693+
if (!track.passedTPCChi2NDF()) {
668694
return false;
669-
if (!track.passedDCAxy())
695+
}
696+
if (!track.passedDCAxy()) {
670697
return false;
671-
if (!track.passedDCAz())
698+
}
699+
if (!track.passedDCAz()) {
672700
return false;
701+
}
673702

674703
return true;
675704
}
@@ -720,8 +749,9 @@ struct FlowEseTask {
720749
qvecRefAInd = refAId * 4 + 3 + (nmode - 2) * cfgNQvec * 4;
721750
qvecRefBInd = refBId * 4 + 3 + (nmode - 2) * cfgNQvec * 4;
722751

723-
if (collision.qvecAmp()[detId] < MinAmplitudeThreshold || collision.qvecAmp()[refAId] < MinAmplitudeThreshold || collision.qvecAmp()[refBId] < MinAmplitudeThreshold)
752+
if (collision.qvecAmp()[detId] < MinAmplitudeThreshold || collision.qvecAmp()[refAId] < MinAmplitudeThreshold || collision.qvecAmp()[refBId] < MinAmplitudeThreshold) {
724753
return;
754+
}
725755

726756
if (nmode == CorrLevel[0]) {
727757
histos.fill(HIST("psi2/QA/EP_Det"), centrality, std::atan2(collision.qvecIm()[qvecDetInd], collision.qvecRe()[qvecDetInd]) / static_cast<float>(nmode));
@@ -864,8 +894,8 @@ struct FlowEseTask {
864894
if (!selectionTrack(trk)) {
865895
continue;
866896
}
867-
const double values[3] = {centrality, trk.pt(), std::cos(static_cast<float>(nmode) * (trk.phi() - esePlane))};
868-
histEseVn->Fill(values);
897+
const std::array<double, 3> values = {centrality, trk.pt(), std::cos(static_cast<float>(nmode) * (trk.phi() - esePlane))};
898+
histEseVn->Fill(values.data());
869899
}
870900
}
871901

@@ -920,11 +950,13 @@ struct FlowEseTask {
920950
aLambdaTag = 1;
921951
}
922952

923-
if (lambdaTag == aLambdaTag)
953+
if (lambdaTag == aLambdaTag) {
924954
continue;
955+
}
925956

926-
if (!selectionV0(collision, v0, lambdaTag))
957+
if (!selectionV0(collision, v0, lambdaTag)) {
927958
continue;
959+
}
928960

929961
if (lambdaTag) {
930962
protonVec = ROOT::Math::PxPyPzMVector(v0.pxpos(), v0.pypos(), v0.pzpos(), massPr);
@@ -996,26 +1028,27 @@ struct FlowEseTask {
9961028

9971029
if (fillEse && cfgFullCheck) {
9981030
const double mass = lambdaTag ? v0.mLambda() : v0.mAntiLambda();
999-
const double cosValues[4] = {mass, v0.pt(), angle * weight, centrality};
1000-
const double cos2Values[4] = {mass, v0.pt(), angle * angle, centrality};
1001-
const double cosSinValues[4] = {mass, v0.pt(), angle * std::sin(relphi) * weight, centrality};
1031+
const std::array<double, 4> cosValues = {mass, v0.pt(), angle * weight, centrality};
1032+
const std::array<double, 4> cos2Values = {mass, v0.pt(), angle * angle, centrality};
1033+
const std::array<double, 4> cosSinValues = {mass, v0.pt(), angle * std::sin(relphi) * weight, centrality};
10021034
if (lambdaTag) {
1003-
histEseLambdaCos->Fill(cosValues);
1004-
histEseLambdaCos2->Fill(cos2Values);
1005-
histEseLambdaCosSin->Fill(cosSinValues);
1035+
histEseLambdaCos->Fill(cosValues.data());
1036+
histEseLambdaCos2->Fill(cos2Values.data());
1037+
histEseLambdaCosSin->Fill(cosSinValues.data());
10061038
} else {
1007-
histEseAntiLambdaCos->Fill(cosValues);
1008-
histEseAntiLambdaCos2->Fill(cos2Values);
1009-
histEseAntiLambdaCosSin->Fill(cosSinValues);
1039+
histEseAntiLambdaCos->Fill(cosValues.data());
1040+
histEseAntiLambdaCos2->Fill(cos2Values.data());
1041+
histEseAntiLambdaCosSin->Fill(cosSinValues.data());
10101042
}
10111043
}
10121044
if (!fillRegular) {
10131045
continue;
10141046
}
10151047

10161048
double qvecMag = 1.0;
1017-
if (cfgUSESP)
1049+
if (cfgUSESP) {
10181050
qvecMag *= std::sqrt(std::pow(collision.qvecIm()[3 + (nmode - 2) * 28], 2) + std::pow(collision.qvecRe()[3 + (nmode - 2) * 28], 2));
1051+
}
10191052

10201053
if (nmode == CorrLevel[0] && cfgFullCheck) { ////////////
10211054
if (lambdaTag) {
@@ -1291,12 +1324,15 @@ struct FlowEseTask {
12911324
float deltaPhi = mcParticle.phi() - mcCollision.eventPlaneAngle();
12921325
// focus on bulk: e, mu, pi, k, p
12931326
int pdgCode = std::abs(mcParticle.pdgCode());
1294-
if (pdgCode != LambdaId)
1327+
if (pdgCode != LambdaId) {
12951328
continue;
1296-
if (!mcParticle.isPhysicalPrimary())
1329+
}
1330+
if (!mcParticle.isPhysicalPrimary()) {
12971331
continue;
1298-
if (std::abs(mcParticle.eta()) > EtaAcceptance) // main acceptance
1332+
}
1333+
if (std::abs(mcParticle.eta()) > EtaAcceptance) { // main acceptance
12991334
continue;
1335+
}
13001336
histos.fill(HIST("hSparseMCGenWeight"), centclass, RecoDecay::constrainAngle(deltaPhi, 0, 2), std::pow(std::cos(2.0 * RecoDecay::constrainAngle(deltaPhi, 0, 2)), 2.0), mcParticle.pt(), mcParticle.eta());
13011337
nCh++;
13021338
bool validGlobal = false;

0 commit comments

Comments
 (0)