Skip to content

Commit c7d9f65

Browse files
alicja-ppalicja-pp
committed
Add pair randomization, fix hardcoded pion mass
1 parent 0bba10f commit c7d9f65

5 files changed

Lines changed: 109 additions & 47 deletions

‎PWGCF/FemtoUniverse/Core/FemtoUniverseContainer.h‎

Lines changed: 27 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -28,6 +28,7 @@
2828
#include <Framework/HistogramSpec.h>
2929

3030
#include <TDatabasePDG.h>
31+
#include <TRandom2.h>
3132

3233
#include <string>
3334
#include <string_view>
@@ -229,36 +230,53 @@ class FemtoUniverseContainer
229230
/// \param part2 Particle two
230231
/// \param mult Multiplicity of the event
231232
template <bool isMC, typename T>
232-
void setPair(T const& part1, T const& part2, const int mult, bool use3dplots, float weight = 1.0f, bool isiden = false)
233+
void setPair(T const& part1, T const& part2, const int mult, bool use3dplots, float weight = 1.0f, bool isiden = false, bool randomizePair = false, double randValue = 0.5)
233234
{
234235
float femtoObs, femtoObsMC;
236+
237+
auto p1 = part1;
238+
auto p2 = part2;
239+
auto mass1 = mMassOne;
240+
auto mass2 = mMassTwo;
241+
if (randomizePair) {
242+
TRandom2* randgen = new TRandom2(0);
243+
double rand = randgen->Rndm();
244+
245+
if (rand > randValue) {
246+
p1 = part2;
247+
p2 = part1;
248+
mass1 = mMassTwo;
249+
mass2 = mMassOne;
250+
}
251+
delete randgen;
252+
}
235253
// Calculate femto observable and the mT with reconstructed information
236254
if constexpr (FemtoObs == femto_universe_container::Observable::kstar) {
237255
if (!isiden) {
238-
femtoObs = FemtoUniverseMath::getkstar(part1, mMassOne, part2, mMassTwo);
256+
femtoObs = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
239257
} else {
240-
femtoObs = 2.0 * FemtoUniverseMath::getkstar(part1, mMassOne, part2, mMassTwo);
258+
femtoObs = 2.0 * FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
241259
}
242260
}
243-
const float mT = FemtoUniverseMath::getmT(part1, mMassOne, part2, mMassTwo);
261+
const float mT = FemtoUniverseMath::getmT(p1, mass1, p2, mass2);
244262

245263
if (mHistogramRegistry) {
246-
setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kRecon>(femtoObs, mT, part1, part2, mult, use3dplots, weight);
264+
setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kRecon>(femtoObs, mT, p1, p2, mult, use3dplots, weight);
247265

248266
if constexpr (isMC) {
249267
if (part1.has_fdMCParticle() && part2.has_fdMCParticle()) {
250268
// calculate the femto observable and the mT with MC truth information
251269
if constexpr (FemtoObs == femto_universe_container::Observable::kstar) {
252270
if (!isiden) {
253-
femtoObsMC = FemtoUniverseMath::getkstar(part1.fdMCParticle(), mMassOne, part2.fdMCParticle(), mMassTwo);
271+
femtoObsMC = FemtoUniverseMath::getkstar(p1.fdMCParticle(), mass1, p2.fdMCParticle(), mass2);
254272
} else {
255-
femtoObsMC = 2.0 * FemtoUniverseMath::getkstar(part1.fdMCParticle(), mMassOne, part2.fdMCParticle(), mMassTwo);
273+
femtoObsMC = 2.0 * FemtoUniverseMath::getkstar(p1.fdMCParticle(), mass1, p2.fdMCParticle(), mass2);
256274
}
257275
}
258-
const float mTMC = FemtoUniverseMath::getmT(part1.fdMCParticle(), mMassOne, part2.fdMCParticle(), mMassTwo);
276+
const float mTMC = FemtoUniverseMath::getmT(part1.fdMCParticle(), mass1, part2.fdMCParticle(), mass2);
259277

260278
if (std::abs(part1.fdMCParticle().pdgMCTruth()) == std::abs(mPDGOne) && std::abs(part2.fdMCParticle().pdgMCTruth()) == std::abs(mPDGTwo)) { // Note: all pair-histogramms are filled with MC truth information ONLY in case of non-fake candidates
261-
setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kTruth>(femtoObsMC, mTMC, part1.fdMCParticle(), part2.fdMCParticle(), mult, use3dplots, weight);
279+
setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kTruth>(femtoObsMC, mTMC, p1.fdMCParticle(), p2.fdMCParticle(), mult, use3dplots, weight);
262280
setPairMC(femtoObsMC, femtoObs, mT, mult);
263281
} else {
264282
mHistogramRegistry->fill(HIST(FolderSuffix[EventType]) + HIST(o2::aod::femtouniverse_mc_particle::MCTypeName[o2::aod::femtouniverse_mc_particle::MCType::kTruth]) + HIST("/hFakePairsCounter"), 0);

‎PWGCF/FemtoUniverse/Core/FemtoUniverseFemtoContainer.h‎

Lines changed: 28 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -27,6 +27,7 @@
2727
#include <Framework/HistogramSpec.h>
2828

2929
#include <TDatabasePDG.h>
30+
#include <TRandom2.h>
3031

3132
#include <string>
3233
#include <string_view>
@@ -205,29 +206,47 @@ class FemtoUniverseFemtoContainer
205206
/// \param part2 Particle two
206207
/// \param mult Multiplicity of the event
207208
template <bool isMC, typename T>
208-
void setPair(T const& part1, T const& part2, const int mult, bool use3dplots)
209+
void setPair(T const& part1, T const& part2, const int mult, bool use3dplots, bool onlyPrimaryMC = false, bool randomizePair = false, double randValue = 0.5)
209210
{
210211
float femtoObs, femtoObsMC;
212+
auto p1 = part1;
213+
auto p2 = part2;
214+
215+
auto mass1 = kMassOne;
216+
auto mass2 = kMassTwo;
217+
if (randomizePair) {
218+
TRandom2* randgen = new TRandom2(0);
219+
double rand = randgen->Rndm();
220+
if (rand > randValue) {
221+
p1 = part2;
222+
p2 = part1;
223+
mass1 = kMassTwo;
224+
mass2 = kMassOne;
225+
}
226+
delete randgen;
227+
}
211228
// Calculate femto observable and the mT with reconstructed information
212229
if constexpr (kFemtoObs == femto_universe_femto_container::Observable::kstar) {
213-
femtoObs = FemtoUniverseMath::getkstar(part1, kMassOne, part2, kMassTwo);
230+
femtoObs = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
214231
}
215-
const float mT = FemtoUniverseMath::getmT(part1, kMassOne, part2, kMassTwo);
232+
const float mT = FemtoUniverseMath::getmT(p1, mass1, p2, mass2);
216233

217234
if (kHistogramRegistry) {
218-
setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kRecon>(femtoObs, mT, part1, part2, mult, use3dplots);
235+
setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kRecon>(femtoObs, mT, p1, p2, mult, use3dplots);
219236

220237
if constexpr (isMC) {
221-
if (part1.has_fdMCParticle() && part2.has_fdMCParticle()) {
238+
if (p1.has_fdMCParticle() && p1.has_fdMCParticle()) {
222239
// calculate the femto observable and the mT with MC truth information
223240
if constexpr (kFemtoObs == femto_universe_femto_container::Observable::kstar) {
224-
femtoObsMC = FemtoUniverseMath::getkstar(part1.fdMCParticle(), kMassOne, part2.fdMCParticle(), kMassTwo);
241+
femtoObsMC = FemtoUniverseMath::getkstar(p1.fdMCParticle(), mass1, p2.fdMCParticle(), mass2);
225242
}
226-
const float mTMC = FemtoUniverseMath::getmT(part1.fdMCParticle(), kMassOne, part2.fdMCParticle(), kMassTwo);
243+
const float mTMC = FemtoUniverseMath::getmT(p1.fdMCParticle(), mass1, p2.fdMCParticle(), mass2);
227244

228245
if (std::abs(part1.fdMCParticle().pdgMCTruth()) == std::abs(kPDGOne) && std::abs(part2.fdMCParticle().pdgMCTruth()) == std::abs(kPDGTwo)) { // Note: all pair-histogramms are filled with MC truth information ONLY in case of non-fake candidates
229-
setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kTruth>(femtoObsMC, mTMC, part1.fdMCParticle(), part2.fdMCParticle(), mult, use3dplots);
230-
setPairMC(femtoObsMC, femtoObs, mT, mult);
246+
if (!onlyPrimaryMC || part1.fdMCParticle().partOriginMCTruth() == o2::aod::femtouniverse_mc_particle::kPrimary && part2.fdMCParticle().partOriginMCTruth() == o2::aod::femtouniverse_mc_particle::kPrimary) {
247+
setPairBase<o2::aod::femtouniverse_mc_particle::MCType::kTruth>(femtoObsMC, mTMC, p1.fdMCParticle(), p1.fdMCParticle(), mult, use3dplots);
248+
setPairMC(femtoObsMC, femtoObs, mT, mult);
249+
}
231250
} else {
232251
kHistogramRegistry->fill(HIST(kFolderSuffix[kEventType]) + HIST(o2::aod::femtouniverse_mc_particle::MCTypeName[o2::aod::femtouniverse_mc_particle::MCType::kTruth]) + HIST("/hFakePairsCounter"), 0);
233252
}

‎PWGCF/FemtoUniverse/Core/FemtoUniversePairSHCentMultKt.h‎

Lines changed: 26 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -190,7 +190,7 @@ class PairSHCentMultKt
190190
/// \param ktval kT value
191191
template <typename T>
192192
void fillMultNumDen(T const& part1, T const& part2, uint8_t ChosenEventType,
193-
int maxl, int multval, float ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF)
193+
int maxl, int multval, float ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF, bool randomizePair = false, double randValue = 0.5)
194194
{
195195
int multbinval;
196196
int absmultval = multval;
@@ -207,7 +207,7 @@ class PairSHCentMultKt
207207
return;
208208
}
209209
// std::cout<<"multbinval "<<multbinval<<std::endl;
210-
fillkTNumDen(part1, part2, ChosenEventType, maxl, multbinval, ktval, isIdenLCMS, isqinvfill, isWeight, isIdenPRF);
210+
fillkTNumDen(part1, part2, ChosenEventType, maxl, multbinval, ktval, isIdenLCMS, isqinvfill, isWeight, isIdenPRF, randomizePair, randValue);
211211
}
212212

213213
/// Templated function to access different kT directory and call addEventPair
@@ -219,7 +219,7 @@ class PairSHCentMultKt
219219
/// \param ktval kT value
220220
template <typename T>
221221
void fillkTNumDen(T const& part1, T const& part2, uint8_t ChosenEventType,
222-
int maxl, int multval, float ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF)
222+
int maxl, int multval, float ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF, bool randomizePair = false, double randValue = 0.5)
223223
{
224224
int ktbinval = -1;
225225
if (ktval >= ktBins[0] && ktval < ktBins[1]) {
@@ -239,16 +239,16 @@ class PairSHCentMultKt
239239
} else {
240240
return;
241241
}
242-
addEventPair(part1, part2, ChosenEventType, maxl, multval, ktbinval, isIdenLCMS, isqinvfill, isWeight, isIdenPRF);
242+
addEventPair(part1, part2, ChosenEventType, maxl, multval, ktbinval, isIdenLCMS, isqinvfill, isWeight, isIdenPRF, randomizePair, randValue);
243243
}
244244

245245
/// Set the PDG codes of the two particles involved
246246
/// \param pdg1 PDG code of particle one
247247
/// \param pdg2 PDG code of particle two
248-
void setPionPairMass()
248+
void setPDGCodes(const int pdg1, const int pdg2)
249249
{
250-
mMassOne = o2::constants::physics::MassPiPlus; // FIXME: Get from the PDG service of the common header
251-
mMassTwo = o2::constants::physics::MassPiPlus; // FIXME: Get from the PDG service of the common header
250+
mMassOne = TDatabasePDG::Instance()->GetParticle(pdg1)->Mass();
251+
mMassTwo = TDatabasePDG::Instance()->GetParticle(pdg2)->Mass();
252252
}
253253

254254
/// To compute the bin value for cavariance matrix
@@ -272,16 +272,31 @@ class PairSHCentMultKt
272272
/// \param ktval kT value
273273
template <typename T>
274274
void addEventPair(T const& part1, T const& part2, uint8_t ChosenEventType,
275-
int /*maxl*/, int multval, int ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF)
275+
int /*maxl*/, int multval, int ktval, bool isIdenLCMS, bool isqinvfill, bool isWeight, bool isIdenPRF, bool randomizePair = false, double randValue = 0.5)
276276
{
277277
int fMultBin = multval;
278278
int fKtBin = ktval;
279279
std::vector<std::complex<double>> fYlmBuffer(kMaxJM);
280280
std::vector<double> f3d;
281-
setPionPairMass();
282281

283-
f3d = FemtoUniverseMath::newpairfunc(part1, mMassOne, part2, mMassTwo,
284-
isIdenLCMS, isWeight, isIdenPRF);
282+
auto p1 = part1;
283+
auto p2 = part2;
284+
auto mass1 = mMassOne;
285+
auto mass2 = mMassTwo;
286+
if (randomizePair) {
287+
TRandom2* randgen = new TRandom2(0);
288+
double rand = randgen->Rndm();
289+
290+
if (rand > randValue) {
291+
p1 = part2;
292+
p2 = part1;
293+
mass1 = mMassTwo;
294+
mass2 = mMassOne;
295+
}
296+
delete randgen;
297+
}
298+
299+
f3d = FemtoUniverseMath::newpairfunc(p1, mass1, p2, mass2, isIdenLCMS, isWeight, isIdenPRF);
285300
double varout = 0.0;
286301
double varside = 0.0;
287302
double varlong = 0.0;

‎PWGCF/FemtoUniverse/Tasks/femtoUniversePairTaskTrackTrackMultKtExtended.cxx‎

Lines changed: 13 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -88,6 +88,9 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended {
8888
Configurable<bool> confIsMC{"confIsMC", false, "Enable additional Histogramms in the case of a MonteCarlo Run"};
8989
Configurable<std::vector<float>> confTrkPIDnSigmaMax{"confTrkPIDnSigmaMax", std::vector<float>{4.f, 3.f, 2.f}, "This configurable needs to be the same as the one used in the producer task"};
9090
Configurable<bool> confUse3D{"confUse3D", false, "Enable three dimensional histogramms (to be used only for analysis with high statistics): k* vs mT vs multiplicity"};
91+
Configurable<bool> confOnlyPrimaryMCPair{"confOnlyPrimaryMCPair", false, "Fill MC pair histograms only with primary particles"};
92+
Configurable<bool> confRandomizePair{"confRandomizePair", true, "Randomize pair before filling pair histograms, like k*"};
93+
9194
} twotracksconfigs;
9295

9396
using FemtoFullParticles = soa::Join<aod::FDParticles, aod::FDExtParticles>;
@@ -536,7 +539,7 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended {
536539
float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
537540
float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2);
538541

539-
sameEventCont.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
542+
sameEventCont.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
540543
if (cfgProcessMultBins)
541544
sameEventMultCont.fill<float>(kstar, multCol, kT);
542545
}
@@ -583,14 +586,14 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended {
583586
float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass1);
584587
float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass1);
585588

586-
sameEventContPP.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
589+
sameEventContPP.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
587590
if (cfgProcessMultBins)
588591
sameEventMultContPP.fill<float>(kstar, multCol, kT);
589592
} else {
590593
float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
591594
float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2);
592595

593-
sameEventContPP.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
596+
sameEventContPP.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
594597
if (cfgProcessMultBins)
595598
sameEventMultContPP.fill<float>(kstar, multCol, kT);
596599
}
@@ -603,14 +606,14 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended {
603606
float kstar = FemtoUniverseMath::getkstar(p1, mass2, p2, mass2);
604607
float kT = FemtoUniverseMath::getkT(p1, mass2, p2, mass2);
605608

606-
sameEventContMM.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
609+
sameEventContMM.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
607610
if (cfgProcessMultBins)
608611
sameEventMultContMM.fill<float>(kstar, multCol, kT);
609612
} else {
610613
float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
611614
float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2);
612615

613-
sameEventContMM.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
616+
sameEventContMM.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
614617
if (cfgProcessMultBins)
615618
sameEventMultContMM.fill<float>(kstar, multCol, kT);
616619
}
@@ -730,7 +733,7 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended {
730733
float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
731734
float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2);
732735

733-
mixedEventCont.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
736+
mixedEventCont.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
734737
if (cfgProcessMultBins)
735738
mixedEventMultCont.fill<float>(kstar, multCol, kT);
736739

@@ -741,14 +744,14 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended {
741744
float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass1);
742745
float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass1);
743746

744-
mixedEventContPP.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
747+
mixedEventContPP.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
745748
if (cfgProcessMultBins)
746749
mixedEventMultContPP.fill<float>(kstar, multCol, kT);
747750
} else {
748751
float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
749752
float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2);
750753

751-
mixedEventContPP.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
754+
mixedEventContPP.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
752755
if (cfgProcessMultBins)
753756
mixedEventMultContPP.fill<float>(kstar, multCol, kT);
754757
}
@@ -761,14 +764,14 @@ struct FemtoUniversePairTaskTrackTrackMultKtExtended {
761764
float kstar = FemtoUniverseMath::getkstar(p1, mass2, p2, mass2);
762765
float kT = FemtoUniverseMath::getkT(p1, mass2, p2, mass2);
763766

764-
mixedEventContMM.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
767+
mixedEventContMM.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
765768
if (cfgProcessMultBins)
766769
mixedEventMultContMM.fill<float>(kstar, multCol, kT);
767770
} else {
768771
float kstar = FemtoUniverseMath::getkstar(p1, mass1, p2, mass2);
769772
float kT = FemtoUniverseMath::getkT(p1, mass1, p2, mass2);
770773

771-
mixedEventContMM.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D);
774+
mixedEventContMM.setPair<isMC>(p1, p2, multCol, twotracksconfigs.confUse3D, twotracksconfigs.confOnlyPrimaryMCPair, twotracksconfigs.confRandomizePair);
772775
if (cfgProcessMultBins)
773776
mixedEventMultContMM.fill<float>(kstar, multCol, kT);
774777
}

0 commit comments

Comments
 (0)