TEvtGen/HepMC/GenParticle.cc

   1 //////////////////////////////////////////////////////////////////////////
   2 // Matt.Dobbs@Cern.CH, September 1999
   3 // Updated: 07.02.2000 no longer does particle point to ParticleData,
   4 //                     but rather it uses an int id which can be looked up
   5 // particle within an event coming in/out of a vertex
   6 //////////////////////////////////////////////////////////////////////////
   7 #include "HepMC/GenEvent.h"
   8 #include "HepMC/GenVertex.h"
   9 #include "HepMC/GenParticle.h"
  10 #include <iomanip>       // needed for formatted output
  11
  12 namespace HepMC {
  13
  14     GenParticle::GenParticle( void ) :
  15         m_momentum(0), m_pdg_id(0), m_status(0), m_flow(this),
  16         m_polarization(0), m_production_vertex(0), m_end_vertex(0),
  17         m_barcode(0), m_generated_mass(0.)
  18     {}
  19     //{
  20         //s_counter++;
  21     //}
  22
  23     GenParticle::GenParticle( const FourVector& momentum,
  24                         int pdg_id, int status,
  25                         const Flow& itsflow,
  26                         const Polarization& polar ) :
  27         m_momentum(momentum), m_pdg_id(pdg_id), m_status(status), m_flow(this),
  28         m_polarization(polar), m_production_vertex(0), m_end_vertex(0),
  29         m_barcode(0), m_generated_mass(momentum.m())
  30     {
  31         // Establishing *this as the owner of m_flow is done above,
  32         // then we set it equal to the other flow pattern (subtle)
  33         set_flow(itsflow);
  34         //s_counter++;
  35     }
  36
  37     GenParticle::GenParticle( const GenParticle& inparticle ) :
  38         m_momentum( inparticle.momentum() ),
  39         m_pdg_id( inparticle.pdg_id() ),
  40         m_status( inparticle.status() ),
  41         m_flow(inparticle.flow()),
  42         m_polarization( inparticle.polarization() ),
  43         m_production_vertex(0),
  44         m_end_vertex(0),
  45         m_barcode(0),
  46         m_generated_mass( inparticle.generated_mass() )
  47     {
  48         /// Shallow copy: does not copy the vertex pointers
  49         /// (note - impossible to copy vertex pointers which having the vertex
  50         ///         and particles in/out point-back to one another -- unless you
  51         ///         copy the entire tree -- which we don't want to do)
  52         set_production_vertex_( 0 );
  53         set_end_vertex_( 0 );
  54         suggest_barcode( inparticle.barcode() );
  55         //s_counter++;
  56     }
  57
  58     GenParticle::~GenParticle() {
  59         if ( parent_event() ) parent_event()->remove_barcode(this);
  60         //s_counter--;
  61     }
  62
  63     void GenParticle::swap( GenParticle & other)
  64     {
  65         // if a container has a swap method, use that for improved performance
  66         m_momentum.swap( other.m_momentum );
  67         std::swap( m_pdg_id, other.m_pdg_id );
  68         std::swap( m_status, other.m_status );
  69         m_flow.swap( other.m_flow );
  70         m_polarization.swap( other.m_polarization );
  71         std::swap( m_production_vertex, other.m_production_vertex );
  72         std::swap( m_end_vertex, other.m_end_vertex );
  73         std::swap( m_barcode, other.m_barcode );
  74         std::swap( m_generated_mass, other.m_generated_mass );
  75     }
  76
  77     GenParticle& GenParticle::operator=( const GenParticle& inparticle ) {
  78         /// Shallow: does not copy the vertex pointers
  79         /// (note - impossible to copy vertex pointers which having the vertex
  80         ///         and particles in/out point-back to one another -- unless you
  81         ///         copy the entire tree -- which we don't want to do)
  82
  83         // best practices implementation
  84         GenParticle tmp( inparticle );
  85         swap( tmp );
  86         return *this;
  87     }
  88
  89     bool GenParticle::operator==( const GenParticle& a ) const {
  90         /// consistent with the definition of the copy constructor as a shallow
  91         ///  constructor,.. this operator does not test the vertex pointers.
  92         ///  Does not compare barcodes.
  93         if ( a.momentum() != this->momentum() ) return false;
  94         if ( a.generated_mass() != this->generated_mass() ) return false;
  95         if ( a.pdg_id() != this->pdg_id() ) return false;
  96         if ( a.status() != this->status() ) return false;
  97         if ( a.m_flow != this->m_flow ) return false;
  98         if ( a.polarization() != this->polarization() ) return false;
  99         return true;
 100     }
 101
 102     bool GenParticle::operator!=( const GenParticle& a ) const {
 103         return !( a == *this );
 104     }
 105
 106     void GenParticle::print( std::ostream& ostr ) const {
 107         /// Dump this particle's full info to ostr, where by default
 108         ///  particle.print(); will dump to cout.
 109         ostr << "GenParticle: "
 110              << barcode() << " ID:" << pdg_id()
 111              << " (P,E)=" << momentum().px() << "," << momentum().py()
 112              << "," << momentum().pz() << "," << momentum().e()
 113              << " Stat:" << status();
 114         if ( production_vertex() && production_vertex()->barcode()!=0 ) {
 115             ostr << " PV:" << production_vertex()->barcode();
 116         } else ostr << " PV:" << production_vertex();
 117         if ( end_vertex() && end_vertex()->barcode()!=0 ) {
 118             ostr << " EV:" << end_vertex()->barcode();
 119         } else ostr << " EV:" << end_vertex();
 120         ostr << " Pol:" << polarization() << " F:" << m_flow << std::endl;
 121     }
 122
 123     GenEvent* GenParticle::parent_event() const {
 124         if ( production_vertex() ) return production_vertex()->parent_event();
 125         if ( end_vertex() ) return end_vertex()->parent_event();
 126         return 0;
 127     }
 128
 129     void GenParticle::set_production_vertex_( GenVertex* prodvertex )
 130     {
 131         GenEvent* its_orig_event = parent_event();
 132         m_production_vertex = prodvertex;
 133         GenEvent* its_new_event = parent_event();
 134         // Next bit of logic ensures the barcode maps are kept up to date
 135         //  in the GenEvent containers.
 136         if ( its_orig_event != its_new_event ) {
 137             if ( its_new_event ) its_new_event->set_barcode( this, barcode() );
 138             if ( its_orig_event ) its_orig_event->remove_barcode( this );
 139         }
 140     }
 141
 142     void GenParticle::set_end_vertex_( GenVertex* decayvertex )
 143     {
 144         GenEvent* its_orig_event = parent_event();
 145         m_end_vertex = decayvertex;
 146         GenEvent* its_new_event = parent_event();
 147         if ( its_orig_event != its_new_event ) {
 148             if ( its_new_event ) its_new_event->set_barcode( this, barcode() );
 149             if ( its_orig_event ) its_orig_event->remove_barcode( this );
 150         }
 151     }
 152
 153     bool GenParticle::suggest_barcode( int the_bar_code )
 154     {
 155         /// allows a barcode to be suggested for this particle.
 156         /// In general it is better to let the event pick the barcode for
 157         /// you, which is automatic.
 158         /// Returns TRUE if the suggested barcode has been accepted (i.e. the
 159         ///  suggested barcode has not already been used in the event,
 160         ///  and so it was used).
 161         /// Returns FALSE if the suggested barcode was rejected, or if the
 162         ///  particle is not yet part of an event, such that it is not yet
 163         ///  possible to know if the suggested barcode will be accepted).
 164         if ( the_bar_code <0 ) {
 165             std::cerr << "GenParticle::suggest_barcode WARNING, particle bar "
 166                       << "\n codes MUST be positive integers. Negative  "
 167                       << "\n integers are reserved for vertices only. Your "
 168                       << "\n suggestion has been rejected." << std::endl;
 169             return false;
 170         }
 171         bool success = false;
 172         if ( parent_event() ) {
 173             success = parent_event()->set_barcode( this, the_bar_code );
 174         } else { set_barcode_( the_bar_code ); }
 175         return success;
 176     }
 177
 178     /////////////
 179     // Static  //
 180     /////////////
 181     //unsigned int GenParticle::counter() { return s_counter; }
 182     //unsigned int GenParticle::s_counter = 0U;
 183
 184     /////////////
 185     // Friends //
 186     /////////////
 187
 188     /// Dump this particle's full info to ostr
 189     std::ostream& operator<<( std::ostream& ostr, const GenParticle& part ) {
 190         // find the current stream state
 191         std::ios_base::fmtflags orig = ostr.flags();
 192         std::streamsize prec = ostr.precision();
 193         ostr << " ";
 194         ostr.width(9);
 195         ostr << part.barcode();
 196         ostr.width(9);
 197         ostr << part.pdg_id() << " ";
 198         ostr.width(9);
 199         ostr.precision(2);
 200         ostr.setf(std::ios::scientific, std::ios::floatfield);
 201         ostr.setf(std::ios_base::showpos);
 202         ostr << part.momentum().px() << ",";
 203         ostr.width(9);
 204         ostr << part.momentum().py() << ",";
 205         ostr.width(9);
 206         ostr << part.momentum().pz() << ",";
 207         ostr.width(9);
 208         ostr << part.momentum().e() << " ";
 209         ostr.setf(std::ios::fmtflags(0), std::ios::floatfield);
 210         ostr.unsetf(std::ios_base::showpos);
 211         if ( part.end_vertex() && part.end_vertex()->barcode()!=0 ) {
 212             ostr.width(3);
 213             ostr << part.status() << " ";
 214             ostr.width(9);
 215             ostr << part.end_vertex()->barcode();
 216         } else if ( !part.end_vertex() ) {
 217             // There is no valid end_vertex
 218             // For consistency across different compilers, do not print anything
 219             ostr.width(3);
 220             ostr << part.status();
 221         } else {
 222             // In this case the end_vertex does not have a unique
 223             //   barcode assigned, so we choose instead to print its address
 224             ostr.width(3);
 225             ostr << part.status() << " ";
 226             ostr.width(9);
 227             ostr << (void*)part.end_vertex();
 228         }
 229         // restore the stream state
 230         ostr.flags(orig);
 231         ostr.precision(prec);
 232         return ostr;
 233     }
 234
 235
 236     double  GenParticle::generated_mass() const {
 237         return m_generated_mass;
 238     }
 239
 240     void   GenParticle::set_generated_mass( const double & m ) {
 241         m_generated_mass = m;
 242     }
 243
 244     /// scale the momentum vector and generated mass
 245     /// this method is only for use by GenEvent
 246     void GenParticle::convert_momentum( const double & f ) {
 247        m_momentum = FourVector( f*m_momentum.px(),
 248                                 f*m_momentum.py(),
 249                                 f*m_momentum.pz(),
 250                                 f*m_momentum.e() );
 251        if( m_generated_mass > 0. ) m_generated_mass = f*m_generated_mass;
 252     }
 253
 254 } // HepMC
 255