@@ -64,18 +64,14 @@ void AODToHepMC::process(Header const& collision,
6464void AODToHepMC::makeEvent (Header const & collision,
6565 Tracks const & tracks)
6666{
67- std::stringstream sb;
68- for (auto b : mBeams )
69- sb << b->id () << " ," ;
7067 LOG (debug) << " Event # " << mEvent .event_number () << " "
7168 << " (# " << mEventNo << " processed)"
7269 << " \n "
7370 << " # of particles " << mParticles .size () << " "
7471 << " (" << tracks.size () << " input tracks)"
7572 << " \n "
7673 << " # of orphans " << mOrphans .size () << " \n "
77- << " # of beams " << mBeams .size () << " ["
78- << sb.str () << " ]\n "
74+ << " # of beams " << mBeams .size () << " \n "
7975 << " # of vertexes " << mVertices .size ();
8076 if (mUseTree ) {
8177 mEvent .reserve (mParticles .size () + mBeams .size (), mVertices .size ());
@@ -105,16 +101,21 @@ void AODToHepMC::makeEvent(Header const& collision,
105101 makeBeams (collision, ip);
106102 }
107103 }
104+ // Flesh out the tracks based on daughter information.
108105 fleshOut (tracks);
109- if (mWriter )
106+ if (mWriter ) {
107+ // If we have a writer, then dump event to output file
110108 mWriter ->write_event (mEvent );
109+ }
111110}
112111// -------------------------------------------------------------------
113112void AODToHepMC::makeHeader (Header const & header)
114113{
115114 int eventNo = header.bcId ();
116- if (eventNo == mLastBC )
115+ if (eventNo == mLastBC ) {
117116 eventNo = mEventNo ;
117+ }
118+
118119 // Set the event number
119120 mEvent .set_event_number (eventNo);
120121 mEvent .weights ().push_back (header.weight ());
@@ -131,13 +132,18 @@ void AODToHepMC::makeHeader(Header const& header)
131132// -------------------------------------------------------------------
132133void AODToHepMC::makeXSection (XSections const & xsections)
133134{
134- if (not mCrossSec )
135+ if (not mCrossSec ) {
136+ // If we do not have a cross-sections object, create it
135137 mCrossSec = std::make_shared<HepMC3::GenCrossSection>();
138+ }
139+
136140 mEvent .set_cross_section (mCrossSec );
137141 mCrossSec ->set_cross_section (0 .f , 0 .f , 0 , 0 );
138142
139- if (xsections.size () <= 0 )
143+ if (xsections.size () <= 0 ) {
144+ // If we have no info, skip the rest
140145 return ;
146+ }
141147
142148 XSection xsection = xsections.iteratorAt (0 );
143149 mCrossSec ->set_cross_section (xsection.xsectGen (),
@@ -148,13 +154,18 @@ void AODToHepMC::makeXSection(XSections const& xsections)
148154// -------------------------------------------------------------------
149155void AODToHepMC::makePdfInfo (PdfInfos const & pdfs)
150156{
151- if (not mPdf )
157+ if (not mPdf ) {
158+ // If we do not have a Parton Distribution Function object, create it
152159 mPdf = std::make_shared<HepMC3::GenPdfInfo>();
160+ }
161+
153162 mEvent .set_pdf_info (mPdf );
154163 mPdf ->set (0 , 0 , 0 .f , 0 .f , 0 .f , 0 .f , 0 .f , 0 , 0 );
155164
156- if (pdfs.size () <= 0 )
165+ if (pdfs.size () <= 0 ) {
166+ // If we have no PDF info, skip the rest
157167 return ;
168+ }
158169
159170 PdfInfo pdf = pdfs.iteratorAt (0 );
160171 mPdf ->set (pdf.id1 (),
@@ -171,9 +182,11 @@ void AODToHepMC::makePdfInfo(PdfInfos const& pdfs)
171182void AODToHepMC::makeHeavyIon (HeavyIons const & heavyions,
172183 Header const & header)
173184{
174- // Generate heavy ion element if it doesn't exist
175- if ( not mIon )
185+ if ( not mIon ) {
186+ // Generate heavy ion element if it doesn't exist
176187 mIon = std::make_shared<HepMC3::GenHeavyIon>();
188+ }
189+
177190 mEvent .set_heavy_ion (mIon );
178191 mIon ->impact_parameter = header.impactParameter ();
179192 mIon ->event_plane_angle = 0 .f ;
@@ -191,26 +204,35 @@ void AODToHepMC::makeHeavyIon(HeavyIons const& heavyions,
191204 mIon ->sigma_inel_NN = 0 .f ;
192205 mIon ->centrality = 0 .f ;
193206#ifndef HEPMC3_NO_DEPRECATED
207+ // Deprecated interface with a single eccentricity
194208 mIon ->eccentricity = 0 .f ;
195209#else
210+ // Newer interface that stores multiple orders of eccentricities.
196211 mIon ->eccentricities [1 ] = 0 .f ;
197212#endif
198213
199- if (heavyions.size () <= 0 )
214+ if (heavyions.size () <= 0 ) {
215+ // If we have no heavy-ion information, skip the rest
200216 return ;
217+ }
201218
202219 HeavyIon heavyion = heavyions.iteratorAt (0 );
203220 float r = 1 ;
221+ // We need to calculate the ratio projectile to target participants
222+ // so that we may break up the number of spectators and so on. This
223+ // is because the AOD HepMC3HeavyIons table does not store the
224+ // relevant information directly.
204225 if (heavyion.npartProj () < heavyion.npartTarg () and
205- heavyion.npartTarg () > 0 )
226+ heavyion.npartTarg () > 0 ) {
206227 r = heavyion.npartProj () / heavyion.npartTarg ();
207- else if (heavyion.npartTarg () < heavyion.npartProj () and
208- heavyion.npartProj () > 0 ) {
228+ } else if (heavyion.npartTarg () < heavyion.npartProj () and
229+ heavyion.npartProj () > 0 ) {
209230 r = heavyion.npartTarg () / heavyion.npartProj ();
210231 r = (1 - r);
211232 }
212233
213- // Heavy ion parameters
234+ // Heavy ion parameters. Note that number of projectile/target
235+ // proton/neutrons are set by the ratio calculated above.
214236 mIon ->impact_parameter = heavyion.impactParameter ();
215237 mIon ->event_plane_angle = heavyion.eventPlaneAngle ();
216238 mIon ->Ncoll_hard = heavyion.ncollHard ();
@@ -235,15 +257,17 @@ void AODToHepMC::makeHeavyIon(HeavyIons const& heavyions,
235257// -------------------------------------------------------------------
236258void AODToHepMC::makeParticles (Tracks const & tracks)
237259{
238- for (auto track : tracks)
260+ for (auto track : tracks) {
239261 makeParticleRecursive (track, tracks);
262+ }
240263}
241264// -------------------------------------------------------------------
242265AODToHepMC::ParticlePtr AODToHepMC::getParticle (Track const & ref) const
243266{
244267 auto iter = mParticles .find (ref.globalIndex ());
245- if (iter == mParticles .end ())
268+ if (iter == mParticles .end ()) {
246269 return 0 ;
270+ }
247271
248272 return iter->second ;
249273}
@@ -266,14 +290,16 @@ AODToHepMC::ParticlePtr AODToHepMC::makeParticleRecursive(Track const& track,
266290 ParticlePtr particle = getParticle (track);
267291
268292 // Check if we already have the particle, and if so, return it
269- if (particle)
293+ if (particle) {
270294 return particle;
295+ }
271296
272297 // Make this particle and store it
273298 int motherStatus = 0 ;
274299 particle = makeParticle (track, motherStatus, force);
275- if (not particle)
300+ if (not particle) {
276301 return 0 ;
302+ }
277303
278304 // Store mapping from index to particle
279305 mParticles [track.globalIndex ()] = particle;
@@ -284,16 +310,20 @@ AODToHepMC::ParticlePtr AODToHepMC::makeParticleRecursive(Track const& track,
284310 for (auto mtrack : track.mothers_as <Tracks>()) {
285311 auto mother = makeParticleRecursive (mtrack, tracks, true );
286312 // If mother not found, continue
287- if (not mother)
313+ if (not mother) {
288314 continue ;
315+ }
316+
289317 // Overrride mother status based on production mechanism of daughter
290- if (motherStatus != 0 )
318+ if (motherStatus != 0 ) {
291319 mother->set_status (motherStatus);
320+ }
292321
293322 mothers.push_back (mother);
294323 // Update the production vertex if not set already
295- if (not vout)
324+ if (not vout) {
296325 vout = mother->end_vertex ();
326+ }
297327 }
298328
299329 // If we have no out vertex, and the particle isn't a beam
@@ -307,23 +337,27 @@ AODToHepMC::ParticlePtr AODToHepMC::makeParticleRecursive(Track const& track,
307337
308338 // If mothers do not have any end-vertex, add them to the found
309339 // vertex.
310- for (auto mother : mothers)
311- if (not mother->end_vertex ())
340+ for (auto mother : mothers) {
341+ if (not mother->end_vertex ()) {
312342 vout->add_particle_in (mother);
343+ }
344+ }
313345 }
314346
315347 // If we got a out-going vertex, add this particle to that
316- if (vout)
348+ if (vout) {
317349 vout->add_particle_out (particle);
350+ }
318351
319352 // If this is a beam particle, add to them
320- if (particle->status () == 4 )
353+ if (particle->status () == 4 ) {
321354 mBeams .push_back (particle);
355+ }
322356 // if if there no mothers, and this is not beam, then make
323357 // this an orphan.
324- else if (mothers.size () <= 0 )
358+ else if (mothers.size () <= 0 ) {
325359 mOrphans .push_back (particle);
326-
360+ }
327361 // return the particle
328362 return particle;
329363}
@@ -367,16 +401,18 @@ AODToHepMC::ParticlePtr AODToHepMC::makeParticle(const Track& track,
367401// -------------------------------------------------------------------
368402void AODToHepMC::fleshOut (Tracks const & tracks)
369403{
370- for (auto track : tracks)
404+ for (auto track : tracks) {
371405 fleshOutParticle (track, tracks);
406+ }
372407}
373408// -------------------------------------------------------------------
374409void AODToHepMC::fleshOutParticle (Track const & track, Tracks const & tracks)
375410{
376411 // If we are only propagating generated tracks, then we need
377412 // not process transported tracks
378- if (isIgnored (track))
413+ if (isIgnored (track)) {
379414 return ;
415+ }
380416
381417 // Check that we can find the track in our map
382418 auto particle = getParticle (track);
@@ -403,8 +439,9 @@ void AODToHepMC::fleshOutParticle(Track const& track, Tracks const& tracks)
403439 // Check that the daughther is generated by EG. If not, and we
404440 // only store generated tracks, then go on to the next daughter
405441 // (or return?).
406- if (isIgnored (dtrack))
442+ if (isIgnored (dtrack)) {
407443 continue ;
444+ }
408445
409446 auto daughter = getParticle (dtrack);
410447 if (not daughter) {
@@ -434,8 +471,9 @@ void AODToHepMC::fleshOutParticle(Track const& track, Tracks const& tracks)
434471 // such daughters.
435472 //
436473 // This check may not be needed
437- if (endVtx and endVtx->id () != prodVtx->id ())
474+ if (endVtx and endVtx->id () != prodVtx->id ()) {
438475 continue ;
476+ }
439477
440478 // If we have a current candidate end vertex, but it doesn't match
441479 // the production vertex of the daughter, then we give a warning.
@@ -456,13 +494,14 @@ void AODToHepMC::fleshOutParticle(Track const& track, Tracks const& tracks)
456494 // end vertex
457495 if (not candidate and
458496 (particle->status () == 4 or
459- particle->status () == 2 ))
497+ particle->status () == 2 )) {
460498 // Only warn for beam and decayed particles
461499 LOG (warning) << " Particle " << track.globalIndex ()
462500 << " (" << particle->id () << " )"
463501 << " w/status " << particle->status ()
464502 << " does not have an end-vertex, "
465503 << " nor was any found from its daughters" ;
504+ }
466505
467506 // If we have a candidate, set the end vertex
468507 if (candidate) {
@@ -471,18 +510,23 @@ void AODToHepMC::fleshOutParticle(Track const& track, Tracks const& tracks)
471510 }
472511 }
473512
513+ // If we have head-less daughters, add them here
474514 if (endVtx and headLess.size () > 0 ) {
475- for (auto daughter : headLess)
515+ for (auto daughter : headLess) {
476516 endVtx->add_particle_out (daughter);
517+ }
477518 }
478519}
479520// -------------------------------------------------------------------
480521void AODToHepMC::enableDump (const std::string& dump)
481522{
482- if (not dump.empty () and mWriter )
523+ if (not dump.empty () and mWriter ) {
483524 return ;
484- if (dump.empty () and not mWriter )
525+ }
526+
527+ if (dump.empty () and not mWriter ) {
485528 return ;
529+ }
486530
487531 if (not dump.empty ()) {
488532 LOG (debug) << " o2::rivet::Converter: Open output HepMC file " << dump;
@@ -495,8 +539,10 @@ void AODToHepMC::enableDump(const std::string& dump)
495539 << " Closing output HepMC file\n "
496540 << " *********************************" ;
497541 mWriter .reset ();
498- if (mOutput )
542+ if (mOutput ) {
499543 mOutput ->close ();
544+ }
545+
500546 delete mOutput ;
501547 mOutput = nullptr ;
502548 }
@@ -509,8 +555,9 @@ AODToHepMC::VertexPtr AODToHepMC::makeVertex(const Track& track)
509555 track.vz (),
510556 track.vt ());
511557 auto vtx = std::make_shared<HepMC3::GenVertex>(v);
512- if (not mUseTree )
558+ if (not mUseTree ) {
513559 mEvent .add_vertex (vtx);
560+ }
514561 return vtx;
515562}
516563// -------------------------------------------------------------------
0 commit comments