Skip to content

Add replacement config to generator_pythia8 to replace D*/Sigmac+ by Xic0/Omegac0 #2186

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
FDUEnrich wants to merge 3 commits into
base: master
Choose a base branch
from
Open

Add replacement config to generator_pythia8 to replace D*/Sigmac+ by Xic0/Omegac0 #2186

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
acb7c46
Fix the mess
fgrosa Jan 16, 2026
d04d7e7
Fix logic to check fraction of good decays
fgrosa Jan 17, 2026
e3f461a
Apply minor adjustements
fgrosa Jan 18, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
34 changes: 21 additions & 13 deletions MC/config/PWGHF/external/generator/generator_pythia8_embed_hf.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@

using namespace Pythia8;

#include <cmath>
namespace hf_generators
{
enum GenType : int {
Expand All @@ -35,7 +36,13 @@ public:
//}

/// Destructor
~GeneratorPythia8EmbedHF() = default;
~GeneratorPythia8EmbedHF() {
// Clean up the internally created HF generator if any
if (mGeneratorEvHF) {
delete mGeneratorEvHF;
mGeneratorEvHF = nullptr;
}
}

/// Init
bool Init() override
Expand All @@ -51,29 +58,30 @@ public:
/// \param yHadronMin minimum hadron rapidity
/// \param yHadronMax maximum hadron rapidity
/// \param hadronPdgList list of PDG codes for hadrons to be used in trigger
void setupGeneratorEvHF(int genType, bool usePtHardBins, float yQuarkMin, float yQuarkMax, float yHadronMin, float yHadronMax, std::vector<int> hadronPdgList = {}) {
/// \param quarkPdgList list of PDG codes for quarks to be enriched in the trigger
void setupGeneratorEvHF(int genType, bool usePtHardBins, float yQuarkMin, float yQuarkMax, float yHadronMin, float yHadronMax, std::vector<int> quarkPdgList = {}, std::vector<int> hadronPdgList = {}, std::vector<std::array<int,2>> partPdgToReplaceList = {}, std::vector<float> freqReplaceList = {}) {
mGeneratorEvHF = nullptr;
switch (genType)
{
case hf_generators::GapTriggeredCharm:
LOG(info) << "********** [GeneratorPythia8EmbedHF] configuring GeneratorPythia8GapTriggeredCharm **********";
LOG(info) << "********** Default number of HF signal events to be merged (updated by notifyEmbedding): " << mNumSigEvs;
mGeneratorEvHF = dynamic_cast<GeneratorPythia8GapTriggeredHF*>(GeneratorPythia8GapTriggeredCharm(/*no gap trigger*/1, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList));
mGeneratorEvHF = dynamic_cast<GeneratorPythia8GapTriggeredHF*>(GeneratorPythia8GapTriggeredCharm(/*no gap trigger*/1, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList, partPdgToReplaceList, freqReplaceList));
break;
case hf_generators::GapTriggeredBeauty:
LOG(info) << "********** [GeneratorPythia8EmbedHF] configuring GeneratorPythia8GapTriggeredBeauty **********";
LOG(info) << "********** Default number of HF signal events to be merged (updated by notifyEmbedding): " << mNumSigEvs;
mGeneratorEvHF = dynamic_cast<GeneratorPythia8GapTriggeredHF*>(GeneratorPythia8GapTriggeredBeauty(/*no gap trigger*/1, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList));
mGeneratorEvHF = dynamic_cast<GeneratorPythia8GapTriggeredHF*>(GeneratorPythia8GapTriggeredBeauty(/*no gap trigger*/1, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList, partPdgToReplaceList, freqReplaceList));
break;
case hf_generators::GapTriggeredCharmAndBeauty:
LOG(info) << "********** [GeneratorPythia8EmbedHF] configuring GeneratorPythia8GapTriggeredCharmAndBeauty **********";
LOG(info) << "********** Default number of HF signal events to be merged (updated by notifyEmbedding): " << mNumSigEvs;
mGeneratorEvHF = dynamic_cast<GeneratorPythia8GapTriggeredHF*>(GeneratorPythia8GapTriggeredCharmAndBeauty(/*no gap trigger*/1, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList));
mGeneratorEvHF = dynamic_cast<GeneratorPythia8GapTriggeredHF*>(GeneratorPythia8GapTriggeredCharmAndBeauty(/*no gap trigger*/1, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList, partPdgToReplaceList, freqReplaceList));
break;
case hf_generators::GapHF:
LOG(info) << "********** [GeneratorPythia8EmbedHF] configuring GeneratorPythia8GapHF **********";
LOG(info) << "********** Default number of HF signal events to be merged (updated by notifyEmbedding): " << mNumSigEvs;
mGeneratorEvHF = dynamic_cast<GeneratorPythia8GapTriggeredHF*>(GeneratorPythia8GapHF(/*no gap trigger*/1, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList));
mGeneratorEvHF = dynamic_cast<GeneratorPythia8GapTriggeredHF*>(GeneratorPythia8GapTriggeredBeauty(/*no gap trigger*/1, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList, partPdgToReplaceList, freqReplaceList));
break;
default:
LOG(fatal) << "********** [GeneratorPythia8EmbedHF] bad configuration, fix it! **********";
Expand Down Expand Up @@ -111,7 +119,7 @@ public:
LOG(info) << "[notifyEmbedding] ----- Collision impact parameter: " << x;

/// number of events to be embedded in a background event
mNumSigEvs = 5 + 0.886202881*std::pow(std::max(0.0f, 17.5f - x),1.7);
mNumSigEvs = static_cast<int>(std::lround(5.0 + 0.886202881 * std::pow(std::max(0.0f, 17.5f - x), 1.7)));
LOG(info) << "[notifyEmbedding] ----- generating " << mNumSigEvs << " signal events " << std::endl;
};

Expand Down Expand Up @@ -380,34 +388,34 @@ private:
};

// Charm enriched
FairGenerator * GeneratorPythia8EmbedHFCharm(bool usePtHardBins = false, float yQuarkMin = -1.5, float yQuarkMax = 1.5, float yHadronMin = -1.5, float yHadronMax = 1.5, std::vector<int> hadronPdgList = {})
FairGenerator * GeneratorPythia8EmbedHFCharm(bool usePtHardBins = false, float yQuarkMin = -1.5, float yQuarkMax = 1.5, float yHadronMin = -1.5, float yHadronMax = 1.5, std::vector<int> quarkPdgList = {}, std::vector<int> hadronPdgList = {}, std::vector<std::array<int,2>> partPdgToReplaceList = {}, std::vector<float> freqReplaceList = {})
{
auto myGen = new GeneratorPythia8EmbedHF();

/// setup the internal generator for HF events
myGen->setupGeneratorEvHF(hf_generators::GapTriggeredCharm, usePtHardBins, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList);
myGen->setupGeneratorEvHF(hf_generators::GapTriggeredCharm, usePtHardBins, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, quarkPdgList, hadronPdgList, partPdgToReplaceList, freqReplaceList);

return myGen;
}

// Beauty enriched
FairGenerator * GeneratorPythia8EmbedHFBeauty(bool usePtHardBins = false, float yQuarkMin = -1.5, float yQuarkMax = 1.5, float yHadronMin = -1.5, float yHadronMax = 1.5, std::vector<int> hadronPdgList = {})
FairGenerator * GeneratorPythia8EmbedHFBeauty(bool usePtHardBins = false, float yQuarkMin = -1.5, float yQuarkMax = 1.5, float yHadronMin = -1.5, float yHadronMax = 1.5, std::vector<int> quarkPdgList = {}, std::vector<int> hadronPdgList = {}, std::vector<std::array<int,2>> partPdgToReplaceList = {}, std::vector<float> freqReplaceList = {})
{
auto myGen = new GeneratorPythia8EmbedHF();

/// setup the internal generator for HF events
myGen->setupGeneratorEvHF(hf_generators::GapTriggeredBeauty, usePtHardBins, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList);
myGen->setupGeneratorEvHF(hf_generators::GapTriggeredBeauty, usePtHardBins, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, quarkPdgList, hadronPdgList, partPdgToReplaceList, freqReplaceList);

return myGen;
}

// Charm and beauty enriched (with same ratio)
FairGenerator * GeneratorPythia8EmbedHFCharmAndBeauty(bool usePtHardBins = false, float yQuarkMin = -1.5, float yQuarkMax = 1.5, float yHadronMin = -1.5, float yHadronMax = 1.5, std::vector<int> hadronPdgList = {})
FairGenerator * GeneratorPythia8EmbedHFCharmAndBeauty(bool usePtHardBins = false, float yQuarkMin = -1.5, float yQuarkMax = 1.5, float yHadronMin = -1.5, float yHadronMax = 1.5, std::vector<int> quarkPdgList = {}, std::vector<int> hadronPdgList = {}, std::vector<std::array<int,2>> partPdgToReplaceList = {}, std::vector<float> freqReplaceList = {})
{
auto myGen = new GeneratorPythia8EmbedHF();

/// setup the internal generator for HF events
myGen->setupGeneratorEvHF(hf_generators::GapTriggeredCharmAndBeauty, usePtHardBins, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, hadronPdgList);
myGen->setupGeneratorEvHF(hf_generators::GapTriggeredCharmAndBeauty, usePtHardBins, yQuarkMin, yQuarkMax, yHadronMin, yHadronMax, quarkPdgList, hadronPdgList, partPdgToReplaceList, freqReplaceList);

return myGen;
}
Expand Down
15 changes: 6 additions & 9 deletions MC/config/PWGHF/external/generator/generator_pythia8_gaptriggered_hf.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -258,17 +258,14 @@ protected:

bool replaceParticle(int iPartToReplace, int pdgCodeNew) {

auto mothers = mPythia.event[iPartToReplace].motherList();

std::array<int, 25> pdgDiquarks = {1103, 2101, 2103, 2203, 3101, 3103, 3201, 3203, 3303, 4101, 4103, 4201, 4203, 4301, 4303, 4403, 5101, 5103, 5201, 5203, 5301, 5303, 5401, 5403, 5503};

for (auto const& mother: mothers) {
auto pdgMother = std::abs(mPythia.event[mother].id());
if (pdgMother > 100 && std::find(pdgDiquarks.begin(), pdgDiquarks.end(), pdgMother) == pdgDiquarks.end()) {
return false;
}
// Check status code: particles with status 91-99 come from decays and should not be replaced
int statusMother = std::abs(mPythia.event[iPartToReplace].status());
if (statusMother >= 91 && statusMother <= 99) {
return false;
}

auto mothers = mPythia.event[iPartToReplace].motherList();

int charge = mPythia.event[iPartToReplace].id() / std::abs(mPythia.event[iPartToReplace].id());
float px = mPythia.event[iPartToReplace].px();
float py = mPythia.event[iPartToReplace].py();
Expand Down
8 changes: 8 additions & 0 deletions MC/config/PWGHF/ini/GeneratorHF_D2H_ccbar_and_bbbar_PbPb_replaced.ini
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,8 @@
### The external generator derives from GeneratorPythia8.
[GeneratorExternal]
fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/external/generator/generator_pythia8_embed_hf.C
funcName=GeneratorPythia8EmbedHFCharmAndBeauty(false, -1.5, 1.5, -1.5, 1.5, std::vector<int>{}, std::vector<int>{}, std::vector<std::array<int,2>>{{423,4132},{423,4232},{4212,4332}}, std::vector<float>{0.5,0.5,1})

[GeneratorPythia8]
config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGHF/pythia8/generator/pythia8_charmhadronic_decays_Mode2_hardQCD_5TeV_XicOmegaC.cfg
includePartonEvent=true
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHFTrigger_Bforced.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -90,8 +90,9 @@ int External()
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.75) // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) // we put some tolerance (e.g. due to oscillations which might change the final state)
{
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHFTrigger_LambdaBToNuclei.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -45,8 +45,9 @@ int External()
std::cout <<"# signal hadrons: " << nSignals << "\n";
std::cout <<"# signal hadrons decaying into nuclei: " << nSignalGoodDecay << "\n";

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.8) // we put some tolerance (lambdaB in MB events do not coalesce)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.2 + uncFracForcedDecays) // we put some tolerance (lambdaB in MB events do not coalesce)
{
std::cerr << "Fraction of signals decaying into nuclei: " << fracForcedDecays << ", lower than expected\n";
return 1;
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHFTrigger_bbbar.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -82,8 +82,9 @@ int External()
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.9) // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) // we put some tolerance (e.g. due to oscillations which might change the final state)
{
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHFTrigger_ccbar.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -82,8 +82,9 @@ int External()
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.85) // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) // we put some tolerance (e.g. due to oscillations which might change the final state)
{
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_bbbar_Bforced_gap5_Mode2.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -105,8 +105,9 @@ int External() {
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.85) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.3 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_bbbar_BtoDK_gap5_Mode2.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -105,8 +105,9 @@ int External() {
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.85) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_bbbar_Mode2_OmegaC.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -95,8 +95,9 @@ int External() {
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_bbbar_Mode2_XiC0_XiCplus.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -106,8 +106,9 @@ int External() {
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_bbbar_gap5.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -102,8 +102,9 @@ int External()
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.7) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_ccbar_Mode2_OmegaC.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -95,8 +95,9 @@ int External() {
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_ccbar_Mode2_XiC0_XiCplus.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -106,8 +106,9 @@ int External() {
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}
Expand Down
7 changes: 4 additions & 3 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_ccbar_and_bbbar_Mode2_ptHardBins.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -123,13 +123,14 @@ int External() {
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}

if (averagePt < 6.5) { // by testing locally it should be around 8.5 GeV/c with pthard bin 20-200 (contrary to 2-2.5 GeV/c of SoftQCD)
if (averagePt < 6) { // by testing locally it should be around 8.5 GeV/c with pthard bin 20-200 (contrary to 2-2.5 GeV/c of SoftQCD)
std::cerr << "Average pT of charmed hadrons " << averagePt << " lower than expected\n";
return 1;
}
Expand Down
5 changes: 3 additions & 2 deletions MC/config/PWGHF/ini/tests/GeneratorHF_D2H_ccbar_and_bbbar_PbPb.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -181,8 +181,9 @@ int External() {
return 1;
}

float fracForcedDecays = float(nSignalGoodDecay) / nSignals;
if (fracForcedDecays < 0.9) { // we put some tolerance (e.g. due to oscillations which might change the final state)
float fracForcedDecays = nSignals ? float(nSignalGoodDecay) / nSignals : 0.0f;
float uncFracForcedDecays = nSignals ? std::sqrt(fracForcedDecays * (1 - fracForcedDecays) / nSignals) : 1.0f;
if (1 - fracForcedDecays > 0.15 + uncFracForcedDecays) { // we put some tolerance (e.g. due to oscillations which might change the final state)
std::cerr << "Fraction of signals decaying into the correct channel " << fracForcedDecays << " lower than expected\n";
return 1;
}
Expand Down
Loading