hepmc - Blame information for rev 2
Subversion Repositories:
| Rev | Author | Line No. | Line |
|---|---|---|---|
| 2 | garren | 1 | ////////////////////////////////////////////////////////////////////////// |
| 2 | // Matt.Dobbs@Cern.CH, September 1999 | ||
| 3 | // Updated: 7.1.2000 iterators complete and working! | ||
| 4 | // Updated: 10.1.2000 GenEvent::vertex, particle iterators are made | ||
| 5 | // constant WRT this event ... note that | ||
| 6 | // GenVertex::***_iterator is not const, since it must | ||
| 7 | // be able to return a mutable pointer to itself. | ||
| 8 | // Updated: 08.2.2000 the event now holds a set of all attached vertices | ||
| 9 | // rather than just the roots of the graph | ||
| 10 | // Event record for MC generators (for use at any stage of generation) | ||
| 11 | ////////////////////////////////////////////////////////////////////////// | ||
| 12 | |||
| 13 | #include "HepMC/GenEvent.h" | ||
| 14 | #include <cstdio> // needed for formatted output using sprintf | ||
| 15 | |||
| 16 | namespace HepMC { | ||
| 17 | |||
| 18 | GenEvent::GenEvent( int signal_process_id, int event_number, | ||
| 19 | GenVertex* signal_vertex, | ||
| 20 | const WeightContainer& weights, | ||
| 21 | const std::vector<long int>& | ||
| 22 | random_states ) : | ||
| 23 | m_signal_process_id(signal_process_id), m_event_number(event_number), | ||
| 24 | m_event_scale(-1), m_alphaQCD(-1), m_alphaQED(-1), | ||
| 25 | m_signal_process_vertex(signal_vertex), m_weights(weights), | ||
| 26 | m_random_states(random_states) | ||
| 27 | { | ||
| 28 | // note: default values for m_event_scale, m_alphaQCD, m_alphaQED | ||
| 29 | // are as suggested in hep-ph/0109068, "Generic Interface..." | ||
| 30 | s_counter++; | ||
| 31 | } | ||
| 32 | |||
| 33 | GenEvent::GenEvent( const GenEvent& inevent ) | ||
| 34 | { | ||
| 35 | // deep copy | ||
| 36 | *this = inevent; | ||
| 37 | s_counter++; | ||
| 38 | } | ||
| 39 | |||
| 40 | GenEvent::~GenEvent() | ||
| 41 | { | ||
| 42 | // Deep destructor. | ||
| 43 | // deletes all vertices/particles in this evt | ||
| 44 | // | ||
| 45 | delete_all_vertices(); | ||
| 46 | s_counter--; | ||
| 47 | } | ||
| 48 | |||
| 49 | GenEvent& GenEvent::operator=( const GenEvent& inevent ) | ||
| 50 | { | ||
| 51 | // deep - makes a copy of all vertices! | ||
| 52 | // | ||
| 53 | // 1. Delete all vertices attached to this | ||
| 54 | delete_all_vertices(); | ||
| 55 | // | ||
| 56 | // 2. create a NEW copy of all vertices from inevent | ||
| 57 | // taking care to map new vertices onto the vertices being copied | ||
| 58 | // and add these new vertices to this event. | ||
| 59 | // We do not use GenVertex::operator= because that would copy | ||
| 60 | // the attached particles as well. | ||
| 61 | std::map<const GenVertex*,GenVertex*> map_in_to_new; | ||
| 62 | for ( GenEvent::vertex_const_iterator v = inevent.vertices_begin(); | ||
| 63 | v != inevent.vertices_end(); v++ ) { | ||
| 64 | GenVertex* newvertex = new GenVertex( | ||
| 65 | (*v)->position(), (*v)->id(), (*v)->weights() ); | ||
| 66 | newvertex->suggest_barcode( (*v)->barcode() ); | ||
| 67 | map_in_to_new[*v] = newvertex; | ||
| 68 | add_vertex( newvertex ); | ||
| 69 | } | ||
| 70 | // 2.b copy the signal process vertex info. | ||
| 71 | if ( inevent.signal_process_vertex() ) { | ||
| 72 | set_signal_process_vertex( | ||
| 73 | map_in_to_new[inevent.signal_process_vertex()] ); | ||
| 74 | } else set_signal_process_vertex( 0 ); | ||
| 75 | // | ||
| 76 | // 3. create a NEW copy of all particles from inevent | ||
| 77 | // taking care to attach them to the appropriate | ||
| 78 | for ( GenEvent::particle_const_iterator p = inevent.particles_begin(); | ||
| 79 | p != inevent.particles_end(); p++ ) | ||
| 80 | { | ||
| 81 | GenParticle* oldparticle = *p; | ||
| 82 | GenParticle* newparticle = new GenParticle(*oldparticle); | ||
| 83 | if ( oldparticle->end_vertex() ) { | ||
| 84 | map_in_to_new[ oldparticle->end_vertex() ]-> | ||
| 85 | add_particle_in(newparticle); | ||
| 86 | } | ||
| 87 | if ( oldparticle->production_vertex() ) { | ||
| 88 | map_in_to_new[ oldparticle->production_vertex() ]-> | ||
| 89 | add_particle_out(newparticle); | ||
| 90 | } | ||
| 91 | } | ||
| 92 | // | ||
| 93 | // 4. now that vtx/particles are copied, do everything else | ||
| 94 | set_signal_process_id( inevent.signal_process_id() ); | ||
| 95 | set_event_number( inevent.event_number() ); | ||
| 96 | set_event_scale( inevent.event_scale() ); | ||
| 97 | set_alphaQCD( inevent.alphaQCD() ); | ||
| 98 | set_alphaQED( inevent.alphaQED() ); | ||
| 99 | set_random_states( inevent.random_states() ); | ||
| 100 | weights() = inevent.weights(); | ||
| 101 | return *this; | ||
| 102 | } | ||
| 103 | |||
| 104 | void GenEvent::print( std::ostream& ostr ) const { | ||
| 105 | // dumps the content of this event to ostr | ||
| 106 | // to dump to cout use: event.print(); | ||
| 107 | // if you want to write this event to file outfile.txt you could use: | ||
| 108 | // std::ofstream outfile("outfile.txt"); event.print( outfile ); | ||
| 109 | ostr << "________________________________________" | ||
| 110 | << "________________________________________\n"; | ||
| 111 | ostr << "GenEvent: #" << event_number() | ||
| 112 | << " ID=" << signal_process_id() | ||
| 113 | << " SignalProcessGenVertex Barcode: " | ||
| 114 | << ( signal_process_vertex() ? signal_process_vertex()->barcode() | ||
| 115 | : 0 ) | ||
| 116 | << "\n"; | ||
| 117 | ostr << " Current Memory Usage: " | ||
| 118 | << GenEvent::counter() << " events, " | ||
| 119 | << GenVertex::counter() << " vertices, " | ||
| 120 | << GenParticle::counter() << " particles.\n"; | ||
| 121 | ostr << " Entries this event: " << vertices_size() << " vertices, " | ||
| 122 | << particles_size() << " particles.\n"; | ||
| 123 | // Random State | ||
| 124 | ostr << " RndmState(" << m_random_states.size() << ")="; | ||
| 125 | for ( std::vector<long int>::const_iterator rs | ||
| 126 | = m_random_states.begin(); | ||
| 127 | rs != m_random_states.end(); rs++ ) { ostr << *rs << " "; } | ||
| 128 | ostr << "\n"; | ||
| 129 | // Weights | ||
| 130 | ostr << " Wgts(" << weights().size() << ")="; | ||
| 131 | for ( WeightContainer::const_iterator wgt = weights().begin(); | ||
| 132 | wgt != weights().end(); wgt++ ) { ostr << *wgt << " "; } | ||
| 133 | ostr << "\n"; | ||
| 134 | ostr << " EventScale " << event_scale() | ||
| 135 | << " [energy] \t alphaQCD=" << alphaQCD() | ||
| 136 | << "\t alphaQED=" << alphaQED() << std::endl; | ||
| 137 | // print a legend to describe the particle info | ||
| 138 | char particle_legend[80]; | ||
| 139 | sprintf( particle_legend," %9s %8s %9s,%9s,%9s,%8s %4s %9s", | ||
| 140 | "Barcode","PDG ID","( Px","Py","Pz","E )","Stat","DecayVtx"); | ||
| 141 | ostr << " GenParticle Legend\n" | ||
| 142 | << particle_legend << "\n"; | ||
| 143 | ostr << "________________________________________" | ||
| 144 | << "________________________________________\n"; | ||
| 145 | // Print all Vertices | ||
| 146 | for ( GenEvent::vertex_const_iterator vtx = this->vertices_begin(); | ||
| 147 | vtx != this->vertices_end(); ++vtx ) { | ||
| 148 | (*vtx)->print(ostr); | ||
| 149 | } | ||
| 150 | ostr << "________________________________________" | ||
| 151 | << "________________________________________" << std::endl; | ||
| 152 | } | ||
| 153 | |||
| 154 | bool GenEvent::add_vertex( GenVertex* vtx ) { | ||
| 155 | // returns true if successful - generally will only return false | ||
| 156 | // if the inserted vertex is already included in the event. | ||
| 157 | if ( !vtx ) return 0; | ||
| 158 | // if vtx previously pointed to another GenEvent, remove it from that | ||
| 159 | // GenEvent's list | ||
| 160 | if ( vtx->parent_event() && vtx->parent_event() != this ) { | ||
| 161 | bool remove_status = vtx->parent_event()->remove_vertex( vtx ); | ||
| 162 | if ( !remove_status ) { | ||
| 163 | std::cerr << "GenEvent::add_vertex ERROR " | ||
| 164 | << "GenVertex::parent_event points to \n" | ||
| 165 | << "an event that does not point back to the " | ||
| 166 | << "GenVertex. \n This probably indicates a deeper " | ||
| 167 | << "problem. " << std::endl; | ||
| 168 | } | ||
| 169 | } | ||
| 170 | // | ||
| 171 | // setting the vertex parent also inserts the vertex into this | ||
| 172 | // event | ||
| 173 | vtx->set_parent_event_( this ); | ||
| 174 | return ( m_vertex_barcodes.count(vtx->barcode()) ? true : false ); | ||
| 175 | } | ||
| 176 | |||
| 177 | bool GenEvent::remove_vertex( GenVertex* vtx ) { | ||
| 178 | // this removes vtx from the event but does NOT delete it. | ||
| 179 | // returns True if an entry vtx existed in the table and was erased | ||
| 180 | if ( m_signal_process_vertex == vtx ) m_signal_process_vertex = 0; | ||
| 181 | if ( vtx->parent_event() == this ) vtx->set_parent_event_( 0 ); | ||
| 182 | return ( m_vertex_barcodes.count(vtx->barcode()) ? false : true ); | ||
| 183 | } | ||
| 184 | |||
| 185 | void GenEvent::delete_all_vertices() { | ||
| 186 | // deletes all vertices in the vertex container | ||
| 187 | // (i.e. all vertices owned by this event) | ||
| 188 | // The vertices are the "owners" of the particles, so as we delete | ||
| 189 | // the vertices, the vertex desctructors are automatically | ||
| 190 | // deleting their particles. | ||
| 191 | if ( vertices_empty() ) return; | ||
| 192 | // delete each vertex individually (this deletes particles as well) | ||
| 193 | while ( m_vertex_barcodes.begin() != m_vertex_barcodes.end() ) { | ||
| 194 | GenVertex* vtx = ( m_vertex_barcodes.begin() )->second; | ||
| 195 | m_vertex_barcodes.erase( m_vertex_barcodes.begin() ); | ||
| 196 | delete vtx; | ||
| 197 | } | ||
| 198 | // | ||
| 199 | // Error checking: | ||
| 200 | if ( !vertices_empty() || ! particles_empty() ) { | ||
| 201 | std::cerr << "GenEvent::delete_all_vertices strange result ... " | ||
| 202 | << "after deleting all vertices, \nthe particle and " | ||
| 203 | << "vertex maps aren't empty.\n This probably " | ||
| 204 | << "indicates deeper problems or memory leak in the " | ||
| 205 | << "code." << std::endl; | ||
| 206 | std::cerr << "Number vtx,particle the event after deleting = " | ||
| 207 | << m_vertex_barcodes.size() << " " | ||
| 208 | << m_particle_barcodes.size() << std::endl; | ||
| 209 | std::cerr << "Total Number vtx,particle in memory " | ||
| 210 | << "after method called = " | ||
| 211 | << GenVertex::counter() << "\t" | ||
| 212 | << GenParticle::counter() << std::endl; | ||
| 213 | } | ||
| 214 | } | ||
| 215 | |||
| 216 | bool GenEvent::set_barcode( GenParticle* p, int suggested_barcode ) | ||
| 217 | { | ||
| 218 | if ( p->parent_event() != this ) { | ||
| 219 | std::cerr << "GenEvent::set_barcode attempted, but the argument's" | ||
| 220 | << "\n parent_event is not this ... request rejected." | ||
| 221 | << std::endl; | ||
| 222 | return false; | ||
| 223 | } | ||
| 224 | // M.Dobbs Nov 4, 2002 | ||
| 225 | // First we must check to see if the particle already has a | ||
| 226 | // barcode which is different from the suggestion. If yes, we | ||
| 227 | // remove it from the particle map. | ||
| 228 | if ( p->barcode() != 0 && p->barcode() != suggested_barcode ) { | ||
| 229 | if ( m_particle_barcodes.count(p->barcode()) && | ||
| 230 | m_particle_barcodes[p->barcode()] == p ) { | ||
| 231 | m_particle_barcodes.erase( p->barcode() ); | ||
| 232 | } | ||
| 233 | // At this point either the particle is NOT in | ||
| 234 | // m_particle_barcodes, or else it is in the map, but | ||
| 235 | // already with the suggested barcode. | ||
| 236 | } | ||
| 237 | // | ||
| 238 | // First case --- a valid barcode has been suggested | ||
| 239 | // (valid barcodes are numbers greater than zero) | ||
| 240 | bool insert_success = true; | ||
| 241 | if ( suggested_barcode > 0 ) { | ||
| 242 | if ( m_particle_barcodes.count(suggested_barcode) ) { | ||
| 243 | // the suggested_barcode is already used. | ||
| 244 | if ( m_particle_barcodes[suggested_barcode] == p ) { | ||
| 245 | // but it was used for this particle ... so everythings ok | ||
| 246 | p->set_barcode_( suggested_barcode ); | ||
| 247 | return true; | ||
| 248 | } | ||
| 249 | insert_success = false; | ||
| 250 | suggested_barcode = 0; | ||
| 251 | } else { // suggested barcode is OK, proceed to insert | ||
| 252 | m_particle_barcodes[suggested_barcode] = p; | ||
| 253 | p->set_barcode_( suggested_barcode ); | ||
| 254 | return true; | ||
| 255 | } | ||
| 256 | } | ||
| 257 | // | ||
| 258 | // Other possibility -- a valid barcode has not been suggested, | ||
| 259 | // so one is automatically generated | ||
| 260 | if ( suggested_barcode < 0 ) insert_success = false; | ||
| 261 | if ( suggested_barcode <= 0 ) { | ||
| 262 | if ( !m_particle_barcodes.empty() ) { | ||
| 263 | // in this case we find the highest barcode that was used, | ||
| 264 | // and increment it by 1 | ||
| 265 | suggested_barcode = m_particle_barcodes.rbegin()->first; | ||
| 266 | ++suggested_barcode; | ||
| 267 | } | ||
| 268 | // For the automatically assigned barcodes, the first one | ||
| 269 | // we use is 10001 ... this is just because when we read | ||
| 270 | // events from HEPEVT, we will suggest their index as the | ||
| 271 | // barcode, and that index has maximum value 10000. | ||
| 272 | // This way we will immediately be able to recognize the hepevt | ||
| 273 | // particles from the auto-assigned ones. | ||
| 274 | if ( suggested_barcode <= 10000 ) suggested_barcode = 10001; | ||
| 275 | } | ||
| 276 | // At this point we should have a valid barcode | ||
| 277 | if ( m_particle_barcodes.count(suggested_barcode) ) { | ||
| 278 | std::cerr << "GenEvent::set_barcode ERROR, this should never " | ||
| 279 | << "happen \n report bug to matt.dobbs@cern.ch" | ||
| 280 | << std::endl; | ||
| 281 | } | ||
| 282 | m_particle_barcodes[suggested_barcode] = p; | ||
| 283 | p->set_barcode_( suggested_barcode ); | ||
| 284 | return insert_success; | ||
| 285 | } | ||
| 286 | |||
| 287 | bool GenEvent::set_barcode( GenVertex* v, int suggested_barcode ) | ||
| 288 | { | ||
| 289 | if ( v->parent_event() != this ) { | ||
| 290 | std::cerr << "GenEvent::set_barcode attempted, but the argument's" | ||
| 291 | << "\n parent_event is not this ... request rejected." | ||
| 292 | << std::endl; | ||
| 293 | return false; | ||
| 294 | } | ||
| 295 | // M.Dobbs Nov 4, 2002 | ||
| 296 | // First we must check to see if the vertex already has a | ||
| 297 | // barcode which is different from the suggestion. If yes, we | ||
| 298 | // remove it from the vertex map. | ||
| 299 | if ( v->barcode() != 0 && v->barcode() != suggested_barcode ) { | ||
| 300 | if ( m_vertex_barcodes.count(v->barcode()) && | ||
| 301 | m_vertex_barcodes[v->barcode()] == v ) { | ||
| 302 | m_vertex_barcodes.erase( v->barcode() ); | ||
| 303 | } | ||
| 304 | // At this point either the vertex is NOT in | ||
| 305 | // m_vertex_barcodes, or else it is in the map, but | ||
| 306 | // already with the suggested barcode. | ||
| 307 | } | ||
| 308 | |||
| 309 | // | ||
| 310 | // First case --- a valid barcode has been suggested | ||
| 311 | // (valid barcodes are numbers greater than zero) | ||
| 312 | bool insert_success = true; | ||
| 313 | if ( suggested_barcode < 0 ) { | ||
| 314 | if ( m_vertex_barcodes.count(suggested_barcode) ) { | ||
| 315 | // the suggested_barcode is already used. | ||
| 316 | if ( m_vertex_barcodes[suggested_barcode] == v ) { | ||
| 317 | // but it was used for this vertex ... so everythings ok | ||
| 318 | v->set_barcode_( suggested_barcode ); | ||
| 319 | return true; | ||
| 320 | } | ||
| 321 | insert_success = false; | ||
| 322 | suggested_barcode = 0; | ||
| 323 | } else { // suggested barcode is OK, proceed to insert | ||
| 324 | m_vertex_barcodes[suggested_barcode] = v; | ||
| 325 | v->set_barcode_( suggested_barcode ); | ||
| 326 | return true; | ||
| 327 | } | ||
| 328 | } | ||
| 329 | // | ||
| 330 | // Other possibility -- a valid barcode has not been suggested, | ||
| 331 | // so one is automatically generated | ||
| 332 | if ( suggested_barcode > 0 ) insert_success = false; | ||
| 333 | if ( suggested_barcode >= 0 ) { | ||
| 334 | if ( !m_vertex_barcodes.empty() ) { | ||
| 335 | // in this case we find the highest barcode that was used, | ||
| 336 | // and increment it by 1, (vertex barcodes are negative) | ||
| 337 | suggested_barcode = m_vertex_barcodes.rbegin()->first; | ||
| 338 | --suggested_barcode; | ||
| 339 | } | ||
| 340 | if ( suggested_barcode >= 0 ) suggested_barcode = -1; | ||
| 341 | } | ||
| 342 | // At this point we should have a valid barcode | ||
| 343 | if ( m_vertex_barcodes.count(suggested_barcode) ) { | ||
| 344 | std::cerr << "GenEvent::set_barcode ERROR, this should never " | ||
| 345 | << "happen \n report bug to matt.dobbs@cern.ch" | ||
| 346 | << std::endl; | ||
| 347 | } | ||
| 348 | m_vertex_barcodes[suggested_barcode] = v; | ||
| 349 | v->set_barcode_( suggested_barcode ); | ||
| 350 | return insert_success; | ||
| 351 | } | ||
| 352 | |||
| 353 | ///////////// | ||
| 354 | // Static // | ||
| 355 | ///////////// | ||
| 356 | unsigned int GenEvent::counter() { return s_counter; } | ||
| 357 | unsigned int GenEvent::s_counter = 0; | ||
| 358 | |||
| 359 | } // HepMC | ||
| 360 | |||
| 361 | |||
| 362 |