-
Notifications
You must be signed in to change notification settings - Fork 496
Expand file tree
/
Copy pathGeneratorTParticle.cxx
More file actions
161 lines (144 loc) · 4.72 KB
/
GeneratorTParticle.cxx
File metadata and controls
161 lines (144 loc) · 4.72 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
// Copyright 2023-2099 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.
/// @author Christian Holm Christensen <cholm@nbi.dk>
#include <Generators/GeneratorTParticle.h>
#include <Generators/GeneratorTParticleParam.h>
#include <SimulationDataFormat/MCGenProperties.h>
#include <SimConfig/SimConfig.h>
#include <fairlogger/Logger.h>
#include <TFile.h>
#include <TChain.h>
#include <TClonesArray.h>
#include <TParticle.h>
namespace o2
{
namespace eventgen
{
/*****************************************************************/
GeneratorTParticle::GeneratorTParticle()
{
setOutputSwitch("-o");
}
/*****************************************************************/
GeneratorTParticle::~GeneratorTParticle()
{
if (mChain) {
TFile* file = mChain->GetCurrentFile();
if (file) {
mChain->RecursiveRemove(file);
}
delete mChain;
}
if (mCmd.empty()) {
return;
}
removeTemp();
}
/*****************************************************************/
Bool_t GeneratorTParticle::Init()
{
mChain = new TChain(mTreeName.c_str());
mTParticles = new TClonesArray("TParticle");
mChain->SetBranchAddress(mBranchName.c_str(), &mTParticles);
if (not mCmd.empty()) {
if (mFileNames.empty()) {
// Set filename to be a temporary name
if (not makeTemp(false)) {
return false;
}
} else {
// Use the first filename as output for cmd line
if (not makeTemp(true)) {
return false;
}
}
// Build command line, Assumes command line parameter
std::string cmd = makeCmdLine();
LOG(info) << "EG command line is \"" << cmd << "\"";
// Execute the background command
if (not executeCmdLine(cmd)) {
LOG(fatal) << "Failed to spawn \"" << cmd << "\"";
return false;
}
}
for (auto filename : mFileNames) {
mChain->AddFile(filename.c_str());
}
// Clear the array of file names
mFileNames.clear();
return true;
}
/*****************************************************************/
void GeneratorTParticle::setup(const GeneratorFileOrCmdParam& param0,
const GeneratorTParticleParam& param,
const conf::SimConfig& config)
{
GeneratorFileOrCmd::setup(param0, config);
setTreeName(param.treeName);
setBranchName(param.branchName);
}
/*****************************************************************/
Bool_t GeneratorTParticle::generateEvent()
{
// If this is the first entry, and we're executing a command, then
// wait until the input file exists and actually contain some data.
if (mEntry == 0 and not mCmd.empty()) {
waitForData(mTemporary);
}
// Read in the next entry in the chain
int read = mChain->GetEntry(mEntry);
mEntry++;
// If we got an error while reading, then give error message
if (read < 0) {
LOG(error) << "Failed to read entry " << mEntry << " of chain";
}
// If we had an error or nothing was read back, then return false
if (read <= 0) {
return false;
}
return true;
}
Bool_t GeneratorTParticle::importParticles()
{
for (auto* object : *mTParticles) {
TParticle* particle = static_cast<TParticle*>(object);
auto statusCode = particle->GetStatusCode();
if (!mcgenstatus::isEncoded(statusCode)) {
statusCode = mcgenstatus::MCGenStatusEncoding(statusCode, 0)
.fullEncoding;
}
mParticles.emplace_back(particle->GetPdgCode(),
statusCode,
particle->GetFirstMother(),
particle->GetSecondMother(),
particle->GetFirstDaughter(),
particle->GetLastDaughter(),
particle->Px(),
particle->Py(),
particle->Pz(),
particle->Energy(),
particle->Vx(),
particle->Vy(),
particle->Vz(),
particle->T());
auto& tgt = mParticles[mParticles.size() - 1];
tgt.SetPolarTheta(particle->GetPolarTheta());
tgt.SetPolarPhi(particle->GetPolarPhi());
tgt.SetCalcMass(particle->GetCalcMass());
tgt.SetWeight(particle->GetWeight());
}
return true;
}
} // namespace eventgen
} // namespace o2
//
// EOF
//