From 17a76b290ed352c24dc2e472f46213ee03a4486f Mon Sep 17 00:00:00 2001 From: ldellost Date: Tue, 3 Nov 2026 10:15:06 +0000 Subject: [PATCH 1/2] Improve Lambda(1520) generator --- .../ini/GeneratorLF_Strangeness_ppLStar.ini | 2 +- .../PWGLF/pythia8/generator_pythia8_Lstar.C | 141 ++++++++++++------ 2 files changed, 94 insertions(+), 49 deletions(-) diff --git a/MC/config/PWGLF/ini/GeneratorLF_Strangeness_ppLStar.ini b/MC/config/PWGLF/ini/GeneratorLF_Strangeness_ppLStar.ini index 9836791ed..5f3809903 100644 --- a/MC/config/PWGLF/ini/GeneratorLF_Strangeness_ppLStar.ini +++ b/MC/config/PWGLF/ini/GeneratorLF_Strangeness_ppLStar.ini @@ -1,6 +1,6 @@ [GeneratorExternal] fileName=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator_pythia8_Lstar.C -funcName=generator_Lstar() +funcName=generator_Lstar(5) [GeneratorPythia8] config=${O2DPG_MC_CONFIG_ROOT}/MC/config/PWGLF/pythia8/generator/pythia8_inel_forceLStar.cfg \ No newline at end of file diff --git a/MC/config/PWGLF/pythia8/generator_pythia8_Lstar.C b/MC/config/PWGLF/pythia8/generator_pythia8_Lstar.C index 66cb4646b..3a732db28 100644 --- a/MC/config/PWGLF/pythia8/generator_pythia8_Lstar.C +++ b/MC/config/PWGLF/pythia8/generator_pythia8_Lstar.C @@ -1,4 +1,3 @@ - #include "Pythia8/Pythia.h" #include "Pythia8/HeavyIons.h" #include "FairGenerator.h" @@ -19,14 +18,25 @@ public: GeneratorPythia8ExtraStrangeness() = default; /// Constructor - GeneratorPythia8ExtraStrangeness(int input_pdg) + GeneratorPythia8ExtraStrangeness(int input_pdg, int inputTriggerRatio) { + + if (inputTriggerRatio <= 0) { + throw std::invalid_argument("inputTriggerRatio must be positive"); + } + + LOG(info) << "Gap trigger ratio N = " << inputTriggerRatio << " events"; + + + mGeneratedEvents = 0; + mTriggerRatio = inputTriggerRatio; + genMinPt=0.0; - genMaxPt=12.0; - genminY=-1.5; - genmaxY=1.5; - genminEta=-1.5; - genmaxEta=1.5; + genMaxPt=10.0; + genminY=-1.0; + genmaxY=1.0; + genminEta=-2.5; + genmaxEta=2.5; pdg=input_pdg; m = 0; @@ -40,14 +50,8 @@ public: xProd=0.; yProd=0.; zProd=0.; fLVHelper = std::make_unique(); - - randomizePDGsign = false; - - fSpectra = makeLevySpectrum("fSpectra", genMinPt, genMaxPt, mPythia.particleData.m0(input_pdg), 0.30, 7.); } - /// randomize the PDG code sign of core particle - void setRandomizePDGsign() { randomizePDGsign = true; } Double_t y2eta(Double_t pt, Double_t mass, Double_t y){ Double_t mt = TMath::Sqrt(mass * mass + pt * pt); @@ -108,22 +112,19 @@ TF1* makeLevySpectrum(const char* name, double ptMin, double ptMax, double mass, //_________________________________________________________________________________ /// generate uniform eta and uniform momentum void genSpectraMomentumEta(double minPt, double maxPt, double minY, double maxY){ - // random generator - std::unique_ptr ranGenerator { new TRandom3() }; - ranGenerator->SetSeed(0); - + // generate transverse momentum //const double gen_pT = ranGenerator->Uniform(minPt, maxPt); - const double gen_pT = fSpectra->GetRandom(minPt,maxPt); + const double gen_pT = fSpectra->GetRandom(minPt, maxPt, &mRandom); //Actually could be something else without loss of generality but okay - const double gen_phi = ranGenerator->Uniform(0,2*TMath::Pi()); + const double gen_phi = mRandom.Uniform(0,2*TMath::Pi()); // sample flat in rapidity, calculate eta Double_t gen_Y=10, gen_eta=10; while( gen_eta>genmaxEta || gen_etaUniform(minY,maxY); + gen_Y = mRandom.Uniform(minY,maxY); gen_eta = y2eta(gen_pT, m, gen_Y); } @@ -132,10 +133,10 @@ TF1* makeLevySpectrum(const char* name, double ptMin, double ptMax, double mass, } //__________________________________________________________________ - Pythia8::Particle createParticle(){ + Pythia8::Particle createParticle(int id){ //std::cout << "createParticle() mass " << m << " pdgCode " << pdg << std::endl; Pythia8::Particle myparticle; - myparticle.id(pdg); + myparticle.id(id); myparticle.status(11); myparticle.px(px); myparticle.py(py); @@ -149,46 +150,92 @@ TF1* makeLevySpectrum(const char* name, double ptMin, double ptMax, double mass, return myparticle; } + //__________________________________________________________________ - int randomizeSign() + // look for final-state (anti)Lambdas in the event; if there are any, remove one + // (chosen at random) from the event record and return the Lambda id, otherwise return 0 + int removeLambda() { - std::unique_ptr gen_random{new TRandom3(0)}; - const float n = gen_random->Uniform(-1, 1); + int kLambdaPdg = 3122; + auto& event = mPythia.event; - return n / abs(n); + std::vector candidates; + for (int i = 0; i < event.size(); ++i) { + bool isFromDecay = event[i].status() >= 91 && event[i].status() <= 99; + bool isAcc = std::abs(event[i].y()) < genmaxY; + + if (event[i].idAbs() == kLambdaPdg && event[i].isFinal() && isAcc && !isFromDecay) { + // not considering lambdas that are the first daughter in the event record. + bool isEdge = false; + for (int mom : event[i].motherList()) { + if (i == event[mom].daughter1() || i == event[mom].daughter2()) isEdge = true; + } + + if (!isEdge) { + candidates.push_back(i); + } + } + } + if (candidates.empty()) { + return 0; + } + + const int idx = candidates[static_cast(mRandom.Uniform(0,candidates.size()))]; + const int lambdaId = event[idx].id(); + + // remove the entry and shift the indices of the others + event.remove(idx, idx, true); + + return lambdaId; } - + //__________________________________________________________________ Bool_t generateEvent() override { // Generate PYTHIA event - mPythia.event.reset(); + if (!fSpectra) { + fSpectra = makeLevySpectrum("fSpectra", genMinPt, genMaxPt, mPythia.particleData.m0(pdg), 0.3942, 9.943); + } - /// go to next Pythia event + mPythia.event.reset(); Bool_t lPythiaOK = kFALSE; - while (!lPythiaOK){ - lPythiaOK = mPythia.next(); - } - /// reset event - //mPythia.event.reset(); + mGeneratedEvents++; + if (mGeneratedEvents % mTriggerRatio == 0) { + int lambdaId = 0; - /// create and append the desired particle - m = sampleMass(pdg); - genSpectraMomentumEta(genMinPt, genMaxPt, genminY, genmaxY); + // regenerate until the event contains a Lambda or an anti-Lambda, which is removed + while (lambdaId == 0) + { + lPythiaOK = mPythia.next(); + if (lPythiaOK) { + lambdaId = removeLambda(); + } + } - if (randomizePDGsign) - pdg *= randomizeSign(); + m = sampleMass(pdg); + genSpectraMomentumEta(genMinPt, genMaxPt, genminY, genmaxY); - Pythia8::Particle particle = createParticle(); - mPythia.event.append(particle); - mPythia.moreDecays(); + // Lambda -> Lambda*, anti-Lambda -> anti-Lambda* + const int sign = lambdaId > 0 ? 1 : -1; + Pythia8::Particle particle = createParticle(sign*pdg); + int idxL = mPythia.event.append(particle); + mPythia.moreDecays(idxL); + } else { + while (!lPythiaOK){ + lPythiaOK = mPythia.next(); + } + } return true; } + private: - + + unsigned long long mGeneratedEvents; + int mTriggerRatio; + double genMinPt; /// minimum 3-momentum for generated particles double genMaxPt; /// maximum 3-momentum for generated particles double genminY; /// minimum pseudorapidity for generated particles @@ -211,16 +258,14 @@ private: double xProd; /// x-coordinate position production vertex [cm] double yProd; /// y-coordinate position production vertex [cm] double zProd; /// z-coordinate position production vertex [cm] - - bool randomizePDGsign; /// bool to randomize the PDG code of the core particle TF1 *fSpectra = nullptr; /// TF1 to store more realistic shape of spectrum std::unique_ptr fLVHelper; + TRandom3 mRandom{0}; }; - FairGenerator *generator_Lstar() + FairGenerator *generator_Lstar(int inputTriggerRatio = 3) { - auto myGen = new GeneratorPythia8ExtraStrangeness(3124); - myGen->setRandomizePDGsign(); // randomization of PDG switched on + auto myGen = new GeneratorPythia8ExtraStrangeness(3124, inputTriggerRatio); return myGen; } \ No newline at end of file From 8bf660a0ea380863d49c4b7a8e46d6defc6d5ebd Mon Sep 17 00:00:00 2001 From: ldellost Date: Tue, 3 Nov 2026 13:14:13 +0000 Subject: [PATCH 2/2] fix the test --- .../tests/GeneratorLF_Strangeness_ppLStar.C | 28 +++++++++++++++---- 1 file changed, 23 insertions(+), 5 deletions(-) diff --git a/MC/config/PWGLF/ini/tests/GeneratorLF_Strangeness_ppLStar.C b/MC/config/PWGLF/ini/tests/GeneratorLF_Strangeness_ppLStar.C index ef04688a0..c69a84d5d 100644 --- a/MC/config/PWGLF/ini/tests/GeneratorLF_Strangeness_ppLStar.C +++ b/MC/config/PWGLF/ini/tests/GeneratorLF_Strangeness_ppLStar.C @@ -1,7 +1,8 @@ int External() { std::string path{"o2sim_Kine.root"}; - const int pdgMother = 3124; // Lambda(1520)0 (sign randomized by the generator) + const int pdgMother = 3124; // Lambda(1520)0 + const int triggerRatio = 5; // must match generator_Lstar(N) in the .ini file TFile file(path.c_str(), "READ"); if (file.IsZombie()) { @@ -20,10 +21,18 @@ int External() int nMother = 0; int nNotDecayed = 0; int nBadDecay = 0; // mothers not decaying exactly to Lambda + gamma - int nEventsWrongCount = 0; // events without exactly one injected mother + int nEventsWrongCount = 0; // events with the wrong number of injected mothers + int nGapEvents = 0; // events expected to contain the injected mother for (int i = 0; i < nEvents; i++) { tree->GetEntry(i); + + // the generator injects in every triggerRatio-th event (it counts events from 1) + const bool isGapEvent = ((i + 1) % triggerRatio == 0); + if (isGapEvent) { + nGapEvents++; + } + int nInEvent = 0; for (auto &track : *tracks) { const int pdg = track.GetPdgCode(); @@ -50,21 +59,30 @@ int External() std::cerr << "Unexpected decay of " << pdg << " (" << nDau << " daughters)\n"; } } - if (nInEvent != 1) { + + const int nExpected = isGapEvent ? 1 : 0; + if (nInEvent != nExpected) { nEventsWrongCount++; + std::cerr << "Event " << i << ": found " << nInEvent + << " Lambda(1520), expected " << nExpected << "\n"; } } std::cout << "--------------------------------\n"; std::cout << "# Events: " << nEvents << "\n"; + std::cout << "# Gap-trigger events: " << nGapEvents << " (ratio " << triggerRatio << ")\n"; std::cout << "# Lambda(1520) + anti: " << nMother << "\n"; std::cout << "# not decayed: " << nNotDecayed << "\n"; std::cout << "# unexpected decays: " << nBadDecay << "\n"; - std::cout << "# events with != 1 mother: " << nEventsWrongCount << "\n"; + std::cout << "# events wrong count: " << nEventsWrongCount << "\n"; std::cout << "--------------------------------\n"; + if (nMother != nGapEvents) { + std::cerr << "Expected one injected Lambda(1520) per gap-trigger event\n"; + return 1; + } if (nEventsWrongCount > 0) { - std::cerr << "Each event must contain exactly one injected Lambda(1520)\n"; + std::cerr << "Lambda(1520) must appear in every gap-triggered events\n"; return 1; } if (nNotDecayed > 0 || nBadDecay > 0) {