Skip to content

Commit 390bc61

Browse files
committed
Decay generator particles before transport in Generator::ReadEvent
This commit adds an optional decay step to Generator::ReadEvent, run after the trigger and before addTracks, and adds a unit test. - The new parameter EventDecayer.enable (default off) switches the step on. - The generator owns a lazily initialised EventDecayer and prints its report when destroyed. - The step uses the generator's position, time, momentum and energy units. - The step adds the header keys eventDecayer and eventDecayerSeed. - The step compiles out without Pythia8. Usage: o2-sim -g evtpool --configKeyValues "GeneratorEventPool.eventPoolPath=pool.root;EventDecayer.enable=true;EventDecayer.config[1]=my.cfg" https://its.cern.ch/jira/browse/O2-5903 Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
1 parent c993a26 commit 390bc61

5 files changed

Lines changed: 168 additions & 6 deletions

File tree

‎Generators/CMakeLists.txt‎

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -157,6 +157,13 @@ if(doBuildSimulation)
157157
ENVIRONMENT O2_ROOT=${CMAKE_BINARY_DIR}/stage
158158
LABELS generator
159159
PUBLIC_LINK_LIBRARIES O2::Generators)
160+
161+
o2_add_test(GeneratorEventDecayer NAME test_Generator_test_GeneratorEventDecayer
162+
SOURCES test/test_GeneratorEventDecayer.cxx
163+
COMPONENT_NAME Generator
164+
ENVIRONMENT O2_ROOT=${CMAKE_BINARY_DIR}/stage
165+
LABELS generator
166+
PUBLIC_LINK_LIBRARIES O2::Generators)
160167
endif()
161168

162169
# o2_add_test(GeneratorPythia8Param NAME test_Generator_test_GeneratorPythia8Param

‎Generators/include/Generators/EventDecayerParam.h‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -23,6 +23,7 @@ namespace eventgen
2323

2424
/// Settings of the EventDecayer, which decays generator particles with Pythia8 before transport
2525
struct EventDecayerParam : public o2::conf::ConfigurableParamHelper<EventDecayerParam> {
26+
bool enable = false; // decay step in Generator::ReadEvent
2627
std::string config[8] = {"", "", "", "", "", "", "", ""}; // Pythia8 config files, an empty slot takes DecayerPythia8.config
2728
double tau0Max = 10.; // species with a larger tau0 in mm/c stay undecayed (<= 0 switches it off)
2829
std::string skipPdg = "115 225 335 10221 10331 20223 20333 102132 102134 9010221"; // PDG codes never decayed, matched on |pdg|

‎Generators/include/Generators/Generator.h‎

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -35,6 +35,8 @@ namespace eventgen
3535
{
3636

3737
class GenTPCLoopers; // Forward declaration
38+
class EventDecayer;
39+
struct EventDecayerReport;
3840

3941
/*****************************************************************/
4042
/*****************************************************************/
@@ -181,6 +183,11 @@ class Generator : public FairGenerator
181183

182184
// Loopers generator instance
183185
o2::eventgen::GenTPCLoopers* mTPCLoopersGen = nullptr;
186+
187+
// decay step after the trigger, created on first use when EventDecayer.enable is set
188+
o2::eventgen::EventDecayer* mEventDecayer = nullptr; //!
189+
o2::eventgen::EventDecayerReport* mEventDecayerReport = nullptr; //!
190+
Bool_t decayParticles(o2::dataformats::MCEventHeader* header);
184191
#ifdef GENERATORS_WITH_TPCLOOPERS
185192
bool initTPCLoopersGen();
186193
#endif

‎Generators/src/Generator.cxx‎

Lines changed: 52 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,10 @@
3333
#include "Generators/TPCLoopers.h"
3434
#include "Generators/TPCLoopersParam.h"
3535
#endif
36+
#ifdef GENERATORS_WITH_PYTHIA8
37+
#include "Generators/EventDecayer.h"
38+
#include "Generators/EventDecayerParam.h"
39+
#endif
3640

3741
namespace o2
3842
{
@@ -117,6 +121,14 @@ Generator::~Generator()
117121
mTPCLoopersGen = nullptr;
118122
}
119123
#endif
124+
#ifdef GENERATORS_WITH_PYTHIA8
125+
if (mEventDecayerReport) {
126+
LOG(info) << "EventDecayer summary of generator " << GetName();
127+
mEventDecayerReport->print();
128+
}
129+
delete mEventDecayer;
130+
delete mEventDecayerReport;
131+
#endif
120132
}
121133

122134
/*****************************************************************/
@@ -335,19 +347,26 @@ Bool_t
335347
}
336348
}
337349

350+
auto header = primGen->GetEvent();
351+
auto o2header = dynamic_cast<o2::dataformats::MCEventHeader*>(header);
352+
if (!header) {
353+
LOG(fatal) << "MC event header is not a 'o2::dataformats::MCEventHeader' object";
354+
return kFALSE;
355+
}
356+
357+
/** decay particles before transport **/
358+
if (!decayParticles(o2header)) {
359+
LOG(error) << "ReadEvent failed in decayParticles";
360+
return kFALSE;
361+
}
362+
338363
/** add tracks **/
339364
if (!addTracks(primGen)) {
340365
LOG(error) << "ReadEvent failed in addTracks";
341366
return kFALSE;
342367
}
343368

344369
/** update header **/
345-
auto header = primGen->GetEvent();
346-
auto o2header = dynamic_cast<o2::dataformats::MCEventHeader*>(header);
347-
if (!header) {
348-
LOG(fatal) << "MC event header is not a 'o2::dataformats::MCEventHeader' object";
349-
return kFALSE;
350-
}
351370
updateHeader(o2header);
352371
updateSubGeneratorInformation(o2header);
353372

@@ -357,6 +376,33 @@ Bool_t
357376

358377
/*****************************************************************/
359378

379+
Bool_t
380+
Generator::decayParticles(o2::dataformats::MCEventHeader* header)
381+
{
382+
#ifdef GENERATORS_WITH_PYTHIA8
383+
if (!EventDecayerParam::Instance().enable) {
384+
return kTRUE;
385+
}
386+
if (!mEventDecayer) {
387+
mEventDecayer = new EventDecayer();
388+
if (!mEventDecayer->init()) {
389+
LOG(fatal) << "Could not initialise the EventDecayer of generator " << GetName();
390+
return kFALSE;
391+
}
392+
mEventDecayerReport = new EventDecayerReport();
393+
}
394+
mEventDecayer->setUnits(mPositionUnit, mTimeUnit, mMomentumUnit, mEnergyUnit);
395+
mEventDecayerReport->merge(mEventDecayer->decay(mParticles));
396+
if (header) {
397+
header->putInfo<std::string>("eventDecayer", mEventDecayer->configSummary());
398+
header->putInfo<int>("eventDecayerSeed", mEventDecayer->seed());
399+
}
400+
#endif
401+
return kTRUE;
402+
}
403+
404+
/*****************************************************************/
405+
360406
Bool_t
361407
Generator::addTracks(FairPrimaryGenerator* primGen)
362408
{
Lines changed: 101 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,101 @@
1+
// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
2+
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
3+
// All rights not expressly granted are reserved.
4+
//
5+
// This software is distributed under the terms of the GNU General Public
6+
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
7+
//
8+
// In applying this license CERN does not waive the privileges and immunities
9+
// granted to it by virtue of its status as an Intergovernmental Organization
10+
// or submit itself to any jurisdiction.
11+
12+
#define BOOST_TEST_MODULE Test decay step of Generator
13+
#define BOOST_TEST_MAIN
14+
#define BOOST_TEST_DYN_LINK
15+
#include <boost/test/unit_test.hpp>
16+
#include <Generators/GeneratorService.h>
17+
#include <SimulationDataFormat/MCTrack.h>
18+
#include <SimulationDataFormat/MCEventHeader.h>
19+
#include <SimulationDataFormat/MCGenProperties.h>
20+
#include <CommonUtils/ConfigurableParam.h>
21+
#include <cmath>
22+
#include <cstdlib>
23+
#include <fstream>
24+
#include <string>
25+
#include <unistd.h>
26+
27+
using namespace o2::eventgen;
28+
29+
namespace
30+
{
31+
// forces Xic+ -> Xi(1530)0 pi+ with the Geant4 lifetime
32+
const std::string FORCE = "4232:mayDecay = on\n4232:onMode = off\n4232:oneChannel = 1 1.0 0 3324 211\n4232:tau0 = 0.1358\n";
33+
34+
std::string writeCfg(std::string const& text)
35+
{
36+
char name[] = "/tmp/test_GeneratorEventDecayer_XXXXXX";
37+
int fd = mkstemp(name);
38+
close(fd);
39+
std::ofstream out(name);
40+
out << text;
41+
return name;
42+
}
43+
44+
// generates one boxgen event of five Xic+ with the decay step on or off
45+
std::pair<std::vector<o2::MCTrack>, o2::dataformats::MCEventHeader> generate(bool enable)
46+
{
47+
static const std::string forcePath = writeCfg(FORCE);
48+
o2::conf::ConfigurableParam::updateFromString("BoxGun.pdg=4232;BoxGun.number=5;EventDecayer.config[1]=" + forcePath +
49+
";EventDecayer.enable=" + (enable ? "true" : "false"));
50+
GeneratorService service;
51+
service.initService("boxgen", "", o2::eventgen::NoVertexOption());
52+
return service.generateEvent();
53+
}
54+
} // namespace
55+
56+
BOOST_AUTO_TEST_CASE(decay_step_on)
57+
{
58+
auto [tracks, header] = generate(true);
59+
int nXic = 0;
60+
int nXiStar = 0;
61+
for (int i = 0; i < (int)tracks.size(); ++i) {
62+
auto const& t = tracks[i];
63+
if (t.GetPdgCode() == 4232) {
64+
++nXic;
65+
BOOST_CHECK(!t.getToBeDone());
66+
BOOST_CHECK_EQUAL(o2::mcgenstatus::getHepMCStatusCode(t.getStatusCode()), 2);
67+
}
68+
if (t.GetPdgCode() == 3324) {
69+
++nXiStar;
70+
int m = t.getMotherTrackId();
71+
BOOST_REQUIRE(m >= 0 && m < (int)tracks.size());
72+
BOOST_CHECK_EQUAL(tracks[m].GetPdgCode(), 4232);
73+
BOOST_CHECK(!t.getToBeDone());
74+
// the Xic+ decay point in cm, displaced from the origin by about its decay length
75+
double r = std::sqrt(t.Vx() * t.Vx() + t.Vy() * t.Vy() + t.Vz() * t.Vz());
76+
BOOST_CHECK(r > 0. && r < 10.);
77+
for (int d = t.getFirstDaughterTrackId(); d >= 0 && d <= t.getLastDaughterTrackId(); ++d) {
78+
BOOST_CHECK_EQUAL(tracks[d].getMotherTrackId(), i);
79+
BOOST_CHECK(tracks[d].getToBeDone());
80+
}
81+
}
82+
}
83+
BOOST_CHECK_EQUAL(nXic, 5);
84+
BOOST_CHECK_EQUAL(nXiStar, 5);
85+
bool valid = false;
86+
auto summary = header.getInfo<std::string>("eventDecayer", valid);
87+
BOOST_CHECK(valid);
88+
BOOST_CHECK(summary.find("pythia8/decays/base.cfg") != std::string::npos);
89+
}
90+
91+
BOOST_AUTO_TEST_CASE(decay_step_off)
92+
{
93+
auto [tracks, header] = generate(false);
94+
BOOST_REQUIRE_EQUAL(tracks.size(), 5);
95+
for (auto const& t : tracks) {
96+
BOOST_CHECK_EQUAL(t.GetPdgCode(), 4232);
97+
BOOST_CHECK(t.getToBeDone());
98+
BOOST_CHECK_EQUAL(t.getFirstDaughterTrackId(), -1);
99+
}
100+
BOOST_CHECK(!header.hasInfo("eventDecayer"));
101+
}

0 commit comments

Comments
 (0)