15 #include "HepMC/GenEvent.h"
16 #include "HepMC/GenCrossSection.h"
17 #include "HepMC/Version.h"
18 #include "HepMC/StreamHelpers.h"
26 const std::vector<long>& random_states,
27 Units::MomentumUnit mom,
28 Units::LengthUnit len ) :
29 m_signal_process_id(signal_process_id),
30 m_event_number(event_number),
35 m_signal_process_vertex(signal_vertex),
39 m_random_states(random_states),
41 m_particle_barcodes(),
58 const std::vector<long>& random_states,
61 Units::MomentumUnit mom,
62 Units::LengthUnit len ) :
63 m_signal_process_id(signal_process_id),
64 m_event_number(event_number),
69 m_signal_process_vertex(signal_vertex),
73 m_random_states(random_states),
75 m_particle_barcodes(),
89 Units::LengthUnit len,
90 int signal_process_id,
94 const std::vector<long>& random_states ) :
95 m_signal_process_id(signal_process_id),
96 m_event_number(event_number),
101 m_signal_process_vertex(signal_vertex),
102 m_beam_particle_1(0),
103 m_beam_particle_2(0),
105 m_random_states(random_states),
107 m_particle_barcodes(),
111 m_momentum_unit(mom),
123 Units::LengthUnit len,
124 int signal_process_id,
int event_number,
127 const std::vector<long>& random_states,
130 m_signal_process_id(signal_process_id),
131 m_event_number(event_number),
136 m_signal_process_vertex(signal_vertex),
137 m_beam_particle_1(0),
138 m_beam_particle_2(0),
140 m_random_states(random_states),
142 m_particle_barcodes(),
145 m_pdf_info( new
PdfInfo(pdf) ),
146 m_momentum_unit(mom),
157 : m_signal_process_id ( inevent.signal_process_id() ),
158 m_event_number ( inevent.event_number() ),
159 m_mpi ( inevent.mpi() ),
160 m_event_scale ( inevent.event_scale() ),
161 m_alphaQCD ( inevent.alphaQCD() ),
162 m_alphaQED ( inevent.alphaQED() ),
163 m_signal_process_vertex( ),
164 m_beam_particle_1 ( ),
165 m_beam_particle_2 ( ),
168 m_vertex_barcodes ( ),
169 m_particle_barcodes ( ),
170 m_cross_section ( inevent.cross_section() ? new
GenCrossSection(*inevent.cross_section()) : 0 ),
171 m_heavy_ion ( inevent.heavy_ion() ? new
HeavyIon(*inevent.heavy_ion()) : 0 ),
172 m_pdf_info ( inevent.pdf_info() ? new
PdfInfo(*inevent.pdf_info()) : 0 ),
173 m_momentum_unit ( inevent.momentum_unit() ),
174 m_position_unit ( inevent.length_unit() )
184 std::map<const GenVertex*,GenVertex*> map_in_to_new;
188 (*v)->position(), (*v)->id(), (*v)->weights() );
190 map_in_to_new[*v] = newvertex;
210 add_particle_in(newparticle);
214 add_particle_out(newparticle);
216 if ( oldparticle == inevent.
beam_particles().first ) beam1 = newparticle;
217 if ( oldparticle == inevent.
beam_particles().second ) beam2 = newparticle;
229 std::swap(m_signal_process_id , other.m_signal_process_id );
230 std::swap(m_event_number , other.m_event_number );
231 std::swap(m_mpi , other.m_mpi );
232 std::swap(m_event_scale , other.m_event_scale );
233 std::swap(m_alphaQCD , other.m_alphaQCD );
234 std::swap(m_alphaQED , other.m_alphaQED );
235 std::swap(m_signal_process_vertex, other.m_signal_process_vertex);
236 std::swap(m_beam_particle_1 , other.m_beam_particle_1 );
237 std::swap(m_beam_particle_2 , other.m_beam_particle_2 );
238 m_weights.
swap( other.m_weights );
239 m_random_states.swap( other.m_random_states );
240 m_vertex_barcodes.swap( other.m_vertex_barcodes );
241 m_particle_barcodes.swap( other.m_particle_barcodes );
242 std::swap(m_cross_section , other.m_cross_section );
243 std::swap(m_heavy_ion , other.m_heavy_ion );
244 std::swap(m_pdf_info , other.m_pdf_info );
245 std::swap(m_momentum_unit , other.m_momentum_unit );
246 std::swap(m_position_unit , other.m_position_unit );
250 (*vthis)->change_parent_event_(
this );
254 (*voth)->change_parent_event_( &other );
264 delete m_cross_section;
282 ostr <<
"________________________________________"
283 <<
"________________________________________\n";
286 <<
" SignalProcessGenVertex Barcode: "
292 ostr <<
" Entries this event: " <<
vertices_size() <<
" vertices, "
294 if( m_beam_particle_1 && m_beam_particle_2 ) {
295 ostr <<
" Beam Particle barcodes: " <<
beam_particles().first->barcode() <<
" "
298 ostr <<
" Beam Particles are not defined.\n";
301 ostr <<
" RndmState(" << m_random_states.size() <<
")=";
302 for ( std::vector<long>::const_iterator rs
303 = m_random_states.begin();
304 rs != m_random_states.end(); ++rs ) { ostr << *rs <<
" "; }
310 <<
" [energy] \t alphaQCD=" <<
alphaQCD()
311 <<
"\t alphaQED=" <<
alphaQED() << std::endl;
313 ostr <<
" GenParticle Legend\n";
314 ostr <<
" Barcode PDG ID "
315 <<
"( Px, Py, Pz, E )"
316 <<
" Stat DecayVtx\n";
317 ostr <<
"________________________________________"
318 <<
"________________________________________\n";
324 ostr <<
"________________________________________"
325 <<
"________________________________________" << std::endl;
329 ostr <<
"---------------------------------------------" << std::endl;
331 ostr <<
"---------------------------------------------" << std::endl;
337 if ( !vtx )
return false;
342 if ( !remove_status ) {
343 std::cerr <<
"GenEvent::add_vertex ERROR "
344 <<
"GenVertex::parent_event points to \n"
345 <<
"an event that does not point back to the "
346 <<
"GenVertex. \n This probably indicates a deeper "
347 <<
"problem. " << std::endl;
354 return ( m_vertex_barcodes.count(vtx->
barcode()) ?
true :
false );
360 if ( m_signal_process_vertex == vtx ) m_signal_process_vertex = 0;
362 return ( m_vertex_barcodes.count(vtx->
barcode()) ?
false :
true );
372 delete m_cross_section;
378 m_signal_process_id = 0;
379 m_beam_particle_1 = 0;
380 m_beam_particle_2 = 0;
386 m_weights = std::vector<double>();
387 m_random_states = std::vector<long>();
389 m_momentum_unit = Units::default_momentum_unit();
390 m_position_unit = Units::default_length_unit();
392 if ( m_vertex_barcodes.size() != 0
393 || m_particle_barcodes.size() != 0 ) {
394 std::cerr <<
"GenEvent::clear() strange result ... \n"
395 <<
"either the particle and/or the vertex map isn't empty" << std::endl;
396 std::cerr <<
"Number vtx,particle the event after deleting = "
397 << m_vertex_barcodes.size() <<
" "
398 << m_particle_barcodes.size() << std::endl;
412 GenVertex* vtx = ( m_vertex_barcodes.begin() )->second;
413 m_vertex_barcodes.erase( m_vertex_barcodes.begin() );
419 std::cerr <<
"GenEvent::delete_all_vertices strange result ... "
420 <<
"after deleting all vertices, \nthe particle and "
421 <<
"vertex maps aren't empty.\n This probably "
422 <<
"indicates deeper problems or memory leak in the "
423 <<
"code." << std::endl;
424 std::cerr <<
"Number vtx,particle the event after deleting = "
425 << m_vertex_barcodes.size() <<
" "
426 << m_particle_barcodes.size() << std::endl;
433 std::cerr <<
"GenEvent::set_barcode attempted, but the argument's"
434 <<
"\n parent_event is not this ... request rejected."
443 if ( m_particle_barcodes.count(p->
barcode()) &&
444 m_particle_barcodes[p->
barcode()] == p ) {
445 m_particle_barcodes.erase( p->
barcode() );
454 bool insert_success =
true;
455 if ( suggested_barcode > 0 ) {
456 if ( m_particle_barcodes.count(suggested_barcode) ) {
458 if ( m_particle_barcodes[suggested_barcode] == p ) {
463 insert_success =
false;
464 suggested_barcode = 0;
466 m_particle_barcodes[suggested_barcode] = p;
474 if ( suggested_barcode < 0 ) insert_success =
false;
475 if ( suggested_barcode <= 0 ) {
476 if ( !m_particle_barcodes.empty() ) {
479 suggested_barcode = m_particle_barcodes.rbegin()->first;
488 if ( suggested_barcode <= 10000 ) suggested_barcode = 10001;
491 if ( m_particle_barcodes.count(suggested_barcode) ) {
492 std::cerr <<
"GenEvent::set_barcode ERROR, this should never "
493 <<
"happen \n report bug to matt.dobbs@cern.ch"
496 m_particle_barcodes[suggested_barcode] = p;
498 return insert_success;
504 std::cerr <<
"GenEvent::set_barcode attempted, but the argument's"
505 <<
"\n parent_event is not this ... request rejected."
514 if ( m_vertex_barcodes.count(v->
barcode()) &&
515 m_vertex_barcodes[v->
barcode()] == v ) {
516 m_vertex_barcodes.erase( v->
barcode() );
526 bool insert_success =
true;
527 if ( suggested_barcode < 0 ) {
528 if ( m_vertex_barcodes.count(suggested_barcode) ) {
530 if ( m_vertex_barcodes[suggested_barcode] == v ) {
535 insert_success =
false;
536 suggested_barcode = 0;
538 m_vertex_barcodes[suggested_barcode] = v;
546 if ( suggested_barcode > 0 ) insert_success =
false;
547 if ( suggested_barcode >= 0 ) {
548 if ( !m_vertex_barcodes.empty() ) {
551 suggested_barcode = m_vertex_barcodes.rbegin()->first;
554 if ( suggested_barcode >= 0 ) suggested_barcode = -1;
557 if ( m_vertex_barcodes.count(suggested_barcode) ) {
558 std::cerr <<
"GenEvent::set_barcode ERROR, this should never "
559 <<
"happen \n report bug to matt.dobbs@cern.ch"
562 m_vertex_barcodes[suggested_barcode] = v;
564 return insert_success;
572 if(m_beam_particle_1 && m_beam_particle_2) {
576 if( m_beam_particle_1 == *p ) have1 =
true;
577 if( m_beam_particle_2 == *p ) have2 =
true;
580 if( have1 && have2 )
return true;
587 m_beam_particle_1 = bp1;
588 m_beam_particle_2 = bp2;
589 if( m_beam_particle_1 && m_beam_particle_2 )
return true;
600 os <<
" Momenutm units:" << std::setw(8) << name(
momentum_unit());
601 os <<
" Position units:" << std::setw(8) << name(
length_unit());
616 bool GenEvent::use_momentum_unit( Units::MomentumUnit newunit ) {
619 if ( m_momentum_unit != newunit ) {
620 const double factor = Units::conversion_factor( m_momentum_unit, newunit );
626 (*p)->convert_momentum(factor);
629 m_momentum_unit = newunit;
634 bool GenEvent::use_length_unit( Units::LengthUnit newunit ) {
637 if ( m_position_unit != newunit ) {
638 const double factor = Units::conversion_factor( m_position_unit, newunit );
643 (*vtx)->convert_position(factor);
646 m_position_unit = newunit;
651 bool GenEvent::use_momentum_unit( std::string& newunit ) {
652 if ( newunit ==
"MEV" )
return use_momentum_unit( Units::MEV );
653 else if( newunit ==
"GEV" )
return use_momentum_unit( Units::GEV );
654 else std::cerr <<
"GenEvent::use_momentum_unit ERROR: use either MEV or GEV\n";
658 bool GenEvent::use_length_unit( std::string& newunit ) {
659 if ( newunit ==
"MM" )
return use_length_unit( Units::MM );
660 else if( newunit ==
"CM" )
return use_length_unit( Units::CM );
661 else std::cerr <<
"GenEvent::use_length_unit ERROR: use either MM or CM\n";
667 if ( new_m ==
"MEV" ) m_momentum_unit = Units::MEV ;
668 else if( new_m ==
"GEV" ) m_momentum_unit = Units::GEV ;
669 else std::cerr <<
"GenEvent::define_units ERROR: use either MEV or GEV\n";
671 if ( new_l ==
"MM" ) m_position_unit = Units::MM ;
672 else if( new_l ==
"CM" ) m_position_unit = Units::CM ;
673 else std::cerr <<
"GenEvent::define_units ERROR: use either MM or CM\n";
685 std::ostream & GenEvent::write_beam_particles(std::ostream & os,
686 std::pair<HepMC::GenParticle *,HepMC::GenParticle *> pr )
690 detail::output( os, 0 );
692 detail::output( os, p->
barcode() );
696 detail::output( os, 0 );
698 detail::output( os, p->
barcode() );
704 std::ostream & GenEvent::write_vertex(std::ostream & os, GenVertex
const * v)
707 std::cerr <<
"GenEvent::write_vertex !v||!os, "
708 <<
"v="<< v <<
" setting badbit" << std::endl;
709 os.clear(std::ios::badbit);
714 int num_orphans_in = 0;
716 = v->particles_in_const_begin();
717 p1 != v->particles_in_const_end(); ++p1 ) {
718 if ( !(*p1)->production_vertex() ) ++num_orphans_in;
722 detail::output( os, v->barcode() );
723 detail::output( os, v->id() );
724 detail::output( os, v->position().x() );
725 detail::output( os, v->position().y() );
726 detail::output( os, v->position().z() );
727 detail::output( os, v->position().t() );
728 detail::output( os, num_orphans_in );
729 detail::output( os, (
int)v->particles_out_size() );
730 detail::output( os, (
int)v->weights().size() );
732 w != v->weights().end(); ++w ) {
733 detail::output( os, *w );
735 detail::output( os,
'\n');
738 = v->particles_in_const_begin();
739 p2 != v->particles_in_const_end(); ++p2 ) {
740 if ( !(*p2)->production_vertex() ) {
741 write_particle( os, *p2 );
746 = v->particles_out_const_begin();
747 p3 != v->particles_out_const_end(); ++p3 ) {
748 write_particle( os, *p3 );
753 std::ostream & GenEvent::write_particle( std::ostream & os, GenParticle
const * p )
756 std::cerr <<
"GenEvent::write_particle !p||!os, "
757 <<
"p="<< p <<
" setting badbit" << std::endl;
758 os.clear(std::ios::badbit);
762 detail::output( os, p->barcode() );
763 detail::output( os, p->pdg_id() );
764 detail::output( os, p->momentum().px() );
765 detail::output( os, p->momentum().py() );
766 detail::output( os, p->momentum().pz() );
767 detail::output( os, p->momentum().e() );
768 detail::output( os, p->generated_mass() );
769 detail::output( os, p->status() );
770 detail::output( os, p->polarization().theta() );
771 detail::output( os, p->polarization().phi() );
774 detail::output( os, ( p->end_vertex() ? p->end_vertex()->barcode() : 0 ) );
775 os <<
' ' << p->flow() <<
"\n";
bool remove_vertex(GenVertex *vtx)
erases vtx from evt
GenEvent * parent_event() const
pointer to the event that owns this particle
void writeVersion(std::ostream &os)
write HepMC version to os
bool particles_empty() const
return true if there are no particle barcodes
void print(std::ostream &ostr=std::cout) const
print weights
bool vertices_empty() const
return true if there are no vertex barcodes
int barcode() const
unique identifier
std::vector< HepMC::GenParticle * >::const_iterator particles_in_const_iterator
const iterator for incoming particles
The GenCrossSection class stores the generated cross section.
bool set_barcode(GenParticle *p, int suggested_barcode=false)
set the barcode - intended for use by GenParticle
void set_signal_process_vertex(GenVertex *)
set pointer to the vertex containing the signal process
int particles_size() const
how many particle barcodes exist?
int barcode() const
particle barcode
particle_const_iterator particles_begin() const
begin particle iteration
std::vector< HepMC::GenParticle * >::const_iterator particles_out_const_iterator
const iterator for outgoing particles
bool suggest_barcode(int the_bar_code)
In general there is no reason to "suggest_barcode".
non-const particle iterator
std::pair< HepMC::GenParticle *, HepMC::GenParticle * > beam_particles() const
pair of pointers to the two incoming beam particles
void delete_all_vertices()
delete all vertices owned by this event
GenCrossSection const * cross_section() const
access the GenCrossSection container if it exists
GenVertex * production_vertex() const
pointer to the production vertex
void set_barcode_(int the_bar_code)
set identifier
GenVertex contains information about decay vertices.
std::vector< double >::const_iterator const_iterator
const iterator for the weight container
void write_units(std::ostream &os=std::cout) const
The GenEvent class is the core of HepMC.
GenVertex * signal_process_vertex() const
pointer to the vertex containing the signal process
Units::MomentumUnit momentum_unit() const
Units used by the GenParticle momentum FourVector.
void swap(GenEvent &other)
swap
void set_random_states(const std::vector< long > &randomstates)
provide random state information
GenEvent(int signal_process_id=0, int event_number=0, GenVertex *signal_vertex=0, const WeightContainer &weights=std::vector< double >(), const std::vector< long > &randomstates=std::vector< long >(), Units::MomentumUnit=Units::default_momentum_unit(), Units::LengthUnit=Units::default_length_unit())
default constructor creates null pointers to HeavyIon, PdfInfo, and GenCrossSection ...
bool add_vertex(GenVertex *vtx)
adds to evt and adopts
void print_version(std::ostream &ostr=std::cout) const
dumps release version to ostr
vertex_const_iterator vertices_end() const
end vertex iteration
void set_barcode_(int the_bar_code)
for use by GenEvent only
WeightContainer & weights()
direct access to WeightContainer
void write_cross_section(std::ostream &ostr=std::cout) const
int vertices_size() const
how many vertex barcodes exist?
bool valid_beam_particles() const
test to see if we have two valid beam particles
GenEvent * parent_event() const
pointer to the event that owns this vertex
bool set_beam_particles(GenParticle *, GenParticle *)
set incoming beam particles
double cross_section_error() const
error associated with this cross section in pb
double cross_section() const
cross section in pb
void print(std::ostream &ostr=std::cout) const
dumps to ostr
void clear()
empties the entire event
particle_const_iterator particles_end() const
end particle iteration
void define_units(Units::MomentumUnit, Units::LengthUnit)
GenEvent & operator=(const GenEvent &inevent)
make a deep copy
void set_parent_event_(GenEvent *evt)
set parent event
vertex_const_iterator vertices_begin() const
begin vertex iteration
GenVertex * end_vertex() const
pointer to the decay vertex
int event_number() const
event number
int signal_process_id() const
unique signal process id
The HeavyIon class stores information about heavy ions.
Container for the Weights associated with an event or vertex.
void swap(WeightContainer &other)
swap
const std::vector< long > & random_states() const
vector of integers containing information about the random state
virtual ~GenEvent()
deletes all vertices/particles in this evt
double event_scale() const
energy scale, see hep-ph/0109068
size_type size() const
size of weight container
The GenParticle class contains information about generated particles.
double alphaQCD() const
QCD coupling, see hep-ph/0109068.
The PdfInfo class stores PDF information.
Units::LengthUnit length_unit() const
Units used by the GenVertex position FourVector.