65 evt = GenEvent(Units::GEV, Units::MM);
68 auto xsec = make_shared<HepMC3::GenCrossSection>();
75 auto heavyion = make_shared<HepMC3::GenHeavyIon>();
92 if (
GetHeader().GetEventPlaneAngle() > -999) {
96 evt.set_heavy_ion(heavyion);
110 VERBOSE(1) <<
" found " << vertices.size() <<
" vertices in the list";
118 VERBOSE(1) <<
" found no hadrons, promoting final partons to status 1";
119 for (
auto p :
evt.particles() ){
120 if (
p->children().size() == 0 ){
121 if (
p->status() !=11 ){
122 JSWARN <<
"Found a final parton with status!=11 : status=" <<
p->status() <<
". This should not happen";
151 if (!pShower->is_acyclic())
152 throw std::runtime_error(
153 "PROBLEM in JetScapeWriterHepMC: Graph is not acyclic.");
159 topsortsearch.
run(*pShower);
164 map<int, GenParticlePtr> CreatedPartons;
167 map<int, bool> vused;
168 bool foundRoot =
false;
173 if (vused.find(nIt->id()) != vused.end()){
174 throw std::runtime_error(
"PROBLEM in JetScapeWriterHepMC: Reusing a vertex.");
176 vused[nIt->id()] =
true;
190 if (nIt->indeg() == 0)
continue;
200 if (nIt->indeg() > 1) {
201 JSWARN <<
"Found more than one mother parton! Should only happen if we "
202 "added medium particles. "
203 <<
"The code should work, but proceed with caution";
206 auto inIt = nIt->in_edges_begin();
207 auto inEnd = nIt->in_edges_end();
208 for (; inIt != inEnd; ++inIt) {
209 auto phepin = CreatedPartons.find(inIt->id());
210 if (phepin != CreatedPartons.end()) {
212 v->add_particle_in(phepin->second);
222 auto in = pShower->GetParton(*inIt);
224 auto status = std::abs(hepin->status());
225 if ( status < 11 || status > 200) {
229 hepin->set_status(
status);
230 CreatedPartons[inIt->id()] = hepin;
231 v->add_particle_in(hepin);
233 if ( nIt->outdeg() == 0 ) {
240 hepout->set_status( 11 );
242 v->add_particle_out(hepout);
251 if (nIt->outdeg() == 0) {
252 if (nIt->indeg() != 1) {
254 throw std::runtime_error(
"PROBLEM in JetScapeWriterHepMC: Need exactly "
255 "one parent to clone final state partons.");
257 auto in = pShower->GetParton(*(nIt->in_edges_begin()));
262 hepout->set_status( 11 );
263 v->add_particle_out(hepout);
268 if (nIt->outdeg() > 0) {
269 auto outIt = nIt->out_edges_begin();
270 auto outEnd = nIt->out_edges_end();
271 for (; outIt != outEnd; ++outIt) {
272 if (CreatedPartons.find(outIt->id()) != CreatedPartons.end()) {
273 throw std::runtime_error(
"PROBLEM in JetScapeWriterHepMC: Trying to "
274 "recreate a preexisting GenParticle.");
276 auto out = pShower->GetParton(*outIt);
278 if ( !hepout->status()) {
280 hepout->set_status(12);
283 CreatedPartons[outIt->id()] = hepout;
284 v->add_particle_out(hepout);
293 auto hadron = h.lock();
307 HepMC3::FourVector vtxPosition(0, 0, 0, 100);
311 HepMC3::FourVector pmom(0, 0, 0, 0);
312 make_shared<GenParticle>(pmom, 0, 0);
321 if ( !hepmc->status() ) {
324 hepmc->set_status(1);
331 JSINFO <<
"JetScape HepMC Writer initialized with output file = "