diff --git a/PWGCF/Femto/Core/closePairRejection.h b/PWGCF/Femto/Core/closePairRejection.h index 9fa4abf3388..c84f9beeae2 100644 --- a/PWGCF/Femto/Core/closePairRejection.h +++ b/PWGCF/Femto/Core/closePairRejection.h @@ -25,10 +25,12 @@ #include "Framework/HistogramSpec.h" #include +#include #include #include #include #include +#include #include #include @@ -68,6 +70,7 @@ struct ConfCpr : o2::framework::ConfigurableGroup { o2::framework::Configurable kinematicMax{"kinematicMax", -1.f, "Maximum kstar/Q3 of pair/triplet for plotting (Set to negative value to turn off the cut)"}; o2::framework::ConfigurableAxis binningDeta{"binningDeta", {{250, -0.5, 0.5}}, "deta"}; o2::framework::ConfigurableAxis binningDphistar{"binningDphistar", {{250, -0.5, 0.5}}, "dphi"}; + o2::framework::Configurable seed{"seed", -1, "Seed to randomize particle 1 and particle 2. Set to negative value to deactivate. Set to 0 to generate unique seed in time."}; }; constexpr const char PrefixCprTrackTrack[] = "CprTrackTrack"; @@ -174,6 +177,17 @@ class CloseTrackRejection mPlotAverage = confCpr.plotAverage.value; mPlotAllRadii = confCpr.plotAllRadii.value; + if (confCpr.seed.value >= 0) { + uint64_t randomSeed; + mRandomizeTracks = true; + if (confCpr.seed.value == 0) { + randomSeed = static_cast(std::chrono::duration_cast(std::chrono::system_clock::now().time_since_epoch()).count()); + } else { + randomSeed = static_cast(confCpr.seed.value); + } + mRng = std::mt19937(randomSeed); + } + // check if we need to apply any cut a plot is requested mIsActivated = mCutAverage || mCutAnyRadius || mPlotAverage || mPlotAllRadii; @@ -210,11 +224,19 @@ class CloseTrackRejection mDphistar.fill(0.f); mDphistarMask.fill(false); - mDeta = track1.eta() - track2.eta(); + bool swapTracks = false; + if (mRandomizeTracks) { + swapTracks = (mSwapDist(mRng) == 1); + } + + auto const& t1 = swapTracks ? track2 : track1; + auto const& t2 = swapTracks ? track1 : track2; + + mDeta = t1.eta() - t2.eta(); for (size_t i = 0; i < TpcRadii.size(); i++) { - auto phistar1 = phistar(mMagField, TpcRadii[i], mChargeAbsTrack1 * track1.signedPt(), track1.phi()); - auto phistar2 = phistar(mMagField, TpcRadii[i], mChargeAbsTrack2 * track2.signedPt(), track2.phi()); + auto phistar1 = phistar(mMagField, TpcRadii[i], mChargeAbsTrack1 * t1.signedPt(), t1.phi()); + auto phistar2 = phistar(mMagField, TpcRadii[i], mChargeAbsTrack2 * t2.signedPt(), t2.phi()); if (phistar1 && phistar2) { mDphistar.at(i) = RecoDecay::constrainAngle(phistar1.value() - phistar2.value(), -o2::constants::math::PI); // constrain angular difference between -pi and pi mDphistarMask.at(i) = true; @@ -342,6 +364,10 @@ class CloseTrackRejection float mDeta = 0.f; std::array mDphistar = {0.f}; std::array mDphistarMask = {false}; + + bool mRandomizeTracks = false; + std::mt19937 mRng; + std::uniform_int_distribution mSwapDist{0, 1}; }; template