| 0ca57c2f |
1 | ////////////////////////////////////////////////////////////////////////// |
| 2 | // Matt.Dobbs@Cern.CH, September 1999 |
| 3 | // GenVertex within an event |
| 4 | ////////////////////////////////////////////////////////////////////////// |
| 5 | |
| 6 | #include "HepMC/GenParticle.h" |
| 7 | #include "HepMC/GenVertex.h" |
| 8 | #include "HepMC/GenEvent.h" |
| 9 | #include "HepMC/SearchVector.h" |
| 10 | #include <iomanip> // needed for formatted output |
| 11 | |
| 12 | namespace HepMC { |
| 13 | |
| 14 | GenVertex::GenVertex( const FourVector& position, |
| 15 | int id, const WeightContainer& weights ) |
| 16 | : m_position(position), m_id(id), m_weights(weights), m_event(0), |
| 17 | m_barcode(0) |
| 18 | {} |
| 19 | //{ |
| 20 | //s_counter++; |
| 21 | //} |
| 22 | |
| 23 | GenVertex::GenVertex( const GenVertex& invertex ) |
| 24 | : m_position( invertex.position() ), |
| 25 | m_particles_in(), |
| 26 | m_particles_out(), |
| 27 | m_id( invertex.id() ), |
| 28 | m_weights( invertex.weights() ), |
| 29 | m_event(0), |
| 30 | m_barcode(0) |
| 31 | { |
| 32 | /// Shallow copy: does not copy the FULL list of particle pointers. |
| 33 | /// Creates a copy of - invertex |
| 34 | /// - outgoing particles of invertex, but sets the |
| 35 | /// decay vertex of these particles to NULL |
| 36 | /// - all incoming particles which do not have a |
| 37 | /// creation vertex. |
| 38 | /// (i.e. it creates copies of all particles which it owns) |
| 39 | /// (note - impossible to copy the FULL list of particle pointers |
| 40 | /// while having the vertex |
| 41 | /// and particles in/out point-back to one another -- unless you |
| 42 | /// copy the entire tree -- which we don't want to do) |
| 43 | // |
| 44 | for ( particles_in_const_iterator |
| 45 | part1 = invertex.particles_in_const_begin(); |
| 46 | part1 != invertex.particles_in_const_end(); ++part1 ) { |
| 47 | if ( !(*part1)->production_vertex() ) { |
| 48 | GenParticle* pin = new GenParticle(**part1); |
| 49 | add_particle_in(pin); |
| 50 | } |
| 51 | } |
| 52 | for ( particles_out_const_iterator |
| 53 | part2 = invertex.particles_out_const_begin(); |
| 54 | part2 != invertex.particles_out_const_end(); part2++ ) { |
| 55 | GenParticle* pin = new GenParticle(**part2); |
| 56 | add_particle_out(pin); |
| 57 | } |
| 58 | suggest_barcode( invertex.barcode() ); |
| 59 | // |
| 60 | //s_counter++; |
| 61 | } |
| 62 | |
| 63 | GenVertex::~GenVertex() { |
| 64 | // |
| 65 | // need to delete any particles previously owned by this vertex |
| 66 | if ( parent_event() ) parent_event()->remove_barcode(this); |
| 67 | delete_adopted_particles(); |
| 68 | //s_counter--; |
| 69 | } |
| 70 | |
| 71 | void GenVertex::swap( GenVertex & other) |
| 72 | { |
| 73 | m_position.swap( other.m_position ); |
| 74 | m_particles_in.swap( other.m_particles_in ); |
| 75 | m_particles_out.swap( other.m_particles_out ); |
| 76 | std::swap( m_id, other.m_id ); |
| 77 | m_weights.swap( other.m_weights ); |
| 78 | std::swap( m_event, other.m_event ); |
| 79 | std::swap( m_barcode, other.m_barcode ); |
| 80 | } |
| 81 | |
| 82 | GenVertex& GenVertex::operator=( const GenVertex& invertex ) { |
| 83 | /// Shallow: does not copy the FULL list of particle pointers. |
| 84 | /// Creates a copy of - invertex |
| 85 | /// - outgoing particles of invertex, but sets the |
| 86 | /// decay vertex of these particles to NULL |
| 87 | /// - all incoming particles which do not have a |
| 88 | /// creation vertex. |
| 89 | /// - it does not alter *this's m_event (!) |
| 90 | /// (i.e. it creates copies of all particles which it owns) |
| 91 | /// (note - impossible to copy the FULL list of particle pointers |
| 92 | /// while having the vertex |
| 93 | /// and particles in/out point-back to one another -- unless you |
| 94 | /// copy the entire tree -- which we don't want to do) |
| 95 | /// |
| 96 | |
| 97 | // best practices implementation |
| 98 | GenVertex tmp( invertex ); |
| 99 | swap( tmp ); |
| 100 | return *this; |
| 101 | } |
| 102 | |
| 103 | bool GenVertex::operator==( const GenVertex& a ) const { |
| 104 | /// Returns true if the positions and the particles in the lists of a |
| 105 | /// and this are identical. Does not compare barcodes. |
| 106 | /// Note that it is impossible for two vertices to point to the same |
| 107 | /// particle's address, so we need to do more than just compare the |
| 108 | /// particle pointers |
| 109 | // |
| 110 | if ( a.position() != this->position() ) return false; |
| 111 | // if the size of the inlist differs, return false. |
| 112 | if ( a.particles_in_size() != this->particles_in_size() ) return false; |
| 113 | // if the size of the outlist differs, return false. |
| 114 | if ( a.particles_out_size() != this->particles_out_size() ) return false; |
| 115 | // loop over the inlist and ensure particles are identical |
| 116 | // (only do this if the lists aren't empty - we already know |
| 117 | // if one isn't, then other isn't either!) |
| 118 | if ( a.particles_in_const_begin() != a.particles_in_const_end() ) { |
| 119 | for ( GenVertex::particles_in_const_iterator |
| 120 | ia = a.particles_in_const_begin(), |
| 121 | ib = this->particles_in_const_begin(); |
| 122 | ia != a.particles_in_const_end(); ia++, ib++ ){ |
| 123 | if ( **ia != **ib ) return false; |
| 124 | } |
| 125 | } |
| 126 | // loop over the outlist and ensure particles are identical |
| 127 | // (only do this if the lists aren't empty - we already know |
| 128 | // if one isn't, then other isn't either!) |
| 129 | if ( a.particles_out_const_begin() != a.particles_out_const_end() ) { |
| 130 | for ( GenVertex::particles_out_const_iterator |
| 131 | ia = a.particles_out_const_begin(), |
| 132 | ib = this->particles_out_const_begin(); |
| 133 | ia != a.particles_out_const_end(); ia++, ib++ ){ |
| 134 | if ( **ia != **ib ) return false; |
| 135 | } |
| 136 | } |
| 137 | return true; |
| 138 | } |
| 139 | |
| 140 | bool GenVertex::operator!=( const GenVertex& a ) const { |
| 141 | // Returns true if the positions and lists of a and this are not equal. |
| 142 | return !( a == *this ); |
| 143 | } |
| 144 | |
| 145 | void GenVertex::print( std::ostream& ostr ) const { |
| 146 | // find the current stream state |
| 147 | std::ios_base::fmtflags orig = ostr.flags(); |
| 148 | std::streamsize prec = ostr.precision(); |
| 149 | if ( barcode()!=0 ) { |
| 150 | if ( position() != FourVector(0,0,0,0) ) { |
| 151 | ostr << "Vertex:"; |
| 152 | ostr.width(9); |
| 153 | ostr << barcode(); |
| 154 | ostr << " ID:"; |
| 155 | ostr.width(5); |
| 156 | ostr << id(); |
| 157 | ostr << " (X,cT)="; |
| 158 | ostr.width(9); |
| 159 | ostr.precision(2); |
| 160 | ostr.setf(std::ios::scientific, std::ios::floatfield); |
| 161 | ostr.setf(std::ios_base::showpos); |
| 162 | ostr << position().x() << ","; |
| 163 | ostr.width(9); |
| 164 | ostr << position().y() << ","; |
| 165 | ostr.width(9); |
| 166 | ostr << position().z() << ","; |
| 167 | ostr.width(9); |
| 168 | ostr << position().t(); |
| 169 | ostr.setf(std::ios::fmtflags(0), std::ios::floatfield); |
| 170 | ostr.unsetf(std::ios_base::showpos); |
| 171 | ostr << std::endl; |
| 172 | } else { |
| 173 | ostr << "GenVertex:"; |
| 174 | ostr.width(9); |
| 175 | ostr << barcode(); |
| 176 | ostr << " ID:"; |
| 177 | ostr.width(5); |
| 178 | ostr << id(); |
| 179 | ostr << " (X,cT):0"; |
| 180 | ostr << std::endl; |
| 181 | } |
| 182 | } else { |
| 183 | // If the vertex doesn't have a unique barcode assigned, then |
| 184 | // we print its memory address instead... so that the |
| 185 | // print out gives us a unique tag for the particle. |
| 186 | if ( position() != FourVector(0,0,0,0) ) { |
| 187 | ostr << "Vertex:"; |
| 188 | ostr.width(9); |
| 189 | ostr << (void*)this; |
| 190 | ostr << " ID:"; |
| 191 | ostr.width(5); |
| 192 | ostr << id(); |
| 193 | ostr << " (X,cT)="; |
| 194 | ostr.width(9); |
| 195 | ostr.precision(2); |
| 196 | ostr.setf(std::ios::scientific, std::ios::floatfield); |
| 197 | ostr.setf(std::ios_base::showpos); |
| 198 | ostr << position().x(); |
| 199 | ostr.width(9); |
| 200 | ostr << position().y(); |
| 201 | ostr.width(9); |
| 202 | ostr << position().z(); |
| 203 | ostr.width(9); |
| 204 | ostr << position().t(); |
| 205 | ostr.setf(std::ios::fmtflags(0), std::ios::floatfield); |
| 206 | ostr.unsetf(std::ios_base::showpos); |
| 207 | ostr << std::endl; |
| 208 | } else { |
| 209 | ostr << "GenVertex:"; |
| 210 | ostr.width(9); |
| 211 | ostr << (void*)this; |
| 212 | ostr << " ID:"; |
| 213 | ostr.width(5); |
| 214 | ostr << id(); |
| 215 | ostr << " (X,cT):0"; |
| 216 | ostr << std::endl; |
| 217 | } |
| 218 | } |
| 219 | |
| 220 | // print the weights if there are any |
| 221 | if ( ! weights().empty() ) { |
| 222 | ostr << " Wgts(" << weights().size() << ")="; |
| 223 | for ( WeightContainer::const_iterator wgt = weights().begin(); |
| 224 | wgt != weights().end(); wgt++ ) { ostr << *wgt << " "; } |
| 225 | ostr << std::endl; |
| 226 | } |
| 227 | // print out all the incoming, then outgoing particles |
| 228 | for ( particles_in_const_iterator part1 = particles_in_const_begin(); |
| 229 | part1 != particles_in_const_end(); part1++ ) { |
| 230 | if ( part1 == particles_in_const_begin() ) { |
| 231 | ostr << " I:"; |
| 232 | ostr.width(2); |
| 233 | ostr << m_particles_in.size(); |
| 234 | } else { ostr << " "; } |
| 235 | //(*part1)->print( ostr ); //uncomment for long debugging printout |
| 236 | ostr << **part1 << std::endl; |
| 237 | } |
| 238 | for ( particles_out_const_iterator part2 = particles_out_const_begin(); |
| 239 | part2 != particles_out_const_end(); part2++ ) { |
| 240 | if ( part2 == particles_out_const_begin() ) { |
| 241 | ostr << " O:"; |
| 242 | ostr.width(2); |
| 243 | ostr << m_particles_out.size(); |
| 244 | } else { ostr << " "; } |
| 245 | //(*part2)->print( ostr ); // uncomment for long debugging printout |
| 246 | ostr << **part2 << std::endl; |
| 247 | } |
| 248 | // restore the stream state |
| 249 | ostr.flags(orig); |
| 250 | ostr.precision(prec); |
| 251 | } |
| 252 | |
| 253 | double GenVertex::check_momentum_conservation() const { |
| 254 | /// finds the difference between the total momentum out and the total |
| 255 | /// momentum in vectors, and returns the magnitude of this vector |
| 256 | /// i.e. returns | vec{p_in} - vec{p_out} | |
| 257 | double sumpx = 0, sumpy = 0, sumpz = 0; |
| 258 | for ( particles_in_const_iterator part1 = particles_in_const_begin(); |
| 259 | part1 != particles_in_const_end(); part1++ ) { |
| 260 | sumpx += (*part1)->momentum().px(); |
| 261 | sumpy += (*part1)->momentum().py(); |
| 262 | sumpz += (*part1)->momentum().pz(); |
| 263 | } |
| 264 | for ( particles_out_const_iterator part2 = particles_out_const_begin(); |
| 265 | part2 != particles_out_const_end(); part2++ ) { |
| 266 | sumpx -= (*part2)->momentum().px(); |
| 267 | sumpy -= (*part2)->momentum().py(); |
| 268 | sumpz -= (*part2)->momentum().pz(); |
| 269 | } |
| 270 | return sqrt( sumpx*sumpx + sumpy*sumpy + sumpz*sumpz ); |
| 271 | } |
| 272 | |
| 273 | void GenVertex::add_particle_in( GenParticle* inparticle ) { |
| 274 | if ( !inparticle ) return; |
| 275 | // if inparticle previously had a decay vertex, remove it from that |
| 276 | // vertex's list |
| 277 | if ( inparticle->end_vertex() ) { |
| 278 | inparticle->end_vertex()->remove_particle_in( inparticle ); |
| 279 | } |
| 280 | m_particles_in.push_back( inparticle ); |
| 281 | inparticle->set_end_vertex_( this ); |
| 282 | } |
| 283 | |
| 284 | void GenVertex::add_particle_out( GenParticle* outparticle ) { |
| 285 | if ( !outparticle ) return; |
| 286 | // if outparticle previously had a production vertex, |
| 287 | // remove it from that vertex's list |
| 288 | if ( outparticle->production_vertex() ) { |
| 289 | outparticle->production_vertex()->remove_particle_out( outparticle ); |
| 290 | } |
| 291 | m_particles_out.push_back( outparticle ); |
| 292 | outparticle->set_production_vertex_( this ); |
| 293 | } |
| 294 | |
| 295 | GenParticle* GenVertex::remove_particle( GenParticle* particle ) { |
| 296 | /// this finds *particle in the in and/or out list and removes it from |
| 297 | /// these lists ... it DOES NOT DELETE THE PARTICLE or its relations. |
| 298 | /// you could delete the particle too as follows: |
| 299 | /// delete vtx->remove_particle( particle ); |
| 300 | /// or if the particle has an end vertex, you could: |
| 301 | /// delete vtx->remove_particle( particle )->end_vertex(); |
| 302 | /// which would delete the particle's end vertex, and thus would |
| 303 | /// also delete the particle, since the particle would be |
| 304 | /// owned by the end vertex. |
| 305 | if ( !particle ) return 0; |
| 306 | if ( particle->end_vertex() == this ) { |
| 307 | particle->set_end_vertex_( 0 ); |
| 308 | remove_particle_in(particle); |
| 309 | } |
| 310 | if ( particle->production_vertex() == this ) { |
| 311 | particle->set_production_vertex_(0); |
| 312 | remove_particle_out(particle); |
| 313 | } |
| 314 | return particle; |
| 315 | } |
| 316 | |
| 317 | void GenVertex::remove_particle_in( GenParticle* particle ) { |
| 318 | /// this finds *particle in m_particles_in and removes it from that list |
| 319 | if ( !particle ) return; |
| 320 | m_particles_in.erase( already_in_vector( &m_particles_in, particle ) ); |
| 321 | } |
| 322 | |
| 323 | void GenVertex::remove_particle_out( GenParticle* particle ) { |
| 324 | /// this finds *particle in m_particles_out and removes it from that list |
| 325 | if ( !particle ) return; |
| 326 | m_particles_out.erase( already_in_vector( &m_particles_out, particle ) ); |
| 327 | } |
| 328 | |
| 329 | void GenVertex::delete_adopted_particles() { |
| 330 | /// deletes all particles which this vertex owns |
| 331 | /// to be used by the vertex destructor and operator= |
| 332 | // |
| 333 | if ( m_particles_out.empty() && m_particles_in.empty() ) return; |
| 334 | // 1. delete all outgoing particles which don't have decay vertices. |
| 335 | // those that do become the responsibility of the decay vertex |
| 336 | // and have their productionvertex pointer set to NULL |
| 337 | for ( std::vector<GenParticle*>::iterator part1 = m_particles_out.begin(); |
| 338 | part1 != m_particles_out.end(); ) { |
| 339 | if ( !(*part1)->end_vertex() ) { |
| 340 | delete *(part1++); |
| 341 | } else { |
| 342 | (*part1)->set_production_vertex_(0); |
| 343 | ++part1; |
| 344 | } |
| 345 | } |
| 346 | m_particles_out.clear(); |
| 347 | // |
| 348 | // 2. delete all incoming particles which don't have production |
| 349 | // vertices. those that do become the responsibility of the |
| 350 | // production vertex and have their decayvertex pointer set to NULL |
| 351 | for ( std::vector<GenParticle*>::iterator part2 = m_particles_in.begin(); |
| 352 | part2 != m_particles_in.end(); ) { |
| 353 | if ( !(*part2)->production_vertex() ) { |
| 354 | delete *(part2++); |
| 355 | } else { |
| 356 | (*part2)->set_end_vertex_(0); |
| 357 | ++part2; |
| 358 | } |
| 359 | } |
| 360 | m_particles_in.clear(); |
| 361 | } |
| 362 | |
| 363 | bool GenVertex::suggest_barcode( int the_bar_code ) |
| 364 | { |
| 365 | /// allows a barcode to be suggested for this vertex. |
| 366 | /// In general it is better to let the event pick the barcode for |
| 367 | /// you, which is automatic. |
| 368 | /// Returns TRUE if the suggested barcode has been accepted (i.e. the |
| 369 | /// suggested barcode has not already been used in the event, |
| 370 | /// and so it was used). |
| 371 | /// Returns FALSE if the suggested barcode was rejected, or if the |
| 372 | /// vertex is not yet part of an event, such that it is not yet |
| 373 | /// possible to know if the suggested barcode will be accepted). |
| 374 | if ( the_bar_code >0 ) { |
| 375 | std::cerr << "GenVertex::suggest_barcode WARNING, vertex bar codes" |
| 376 | << "\n MUST be negative integers. Positive integers " |
| 377 | << "\n are reserved for particles only. Your suggestion " |
| 378 | << "\n has been rejected." << std::endl; |
| 379 | return false; |
| 380 | } |
| 381 | bool success = false; |
| 382 | if ( parent_event() ) { |
| 383 | success = parent_event()->set_barcode( this, the_bar_code ); |
| 384 | } else { set_barcode_( the_bar_code ); } |
| 385 | return success; |
| 386 | } |
| 387 | |
| 388 | void GenVertex::set_parent_event_( GenEvent* new_evt ) |
| 389 | { |
| 390 | GenEvent* orig_evt = m_event; |
| 391 | m_event = new_evt; |
| 392 | // |
| 393 | // every time a vertex's parent event changes, the map of barcodes |
| 394 | // in the new and old parent event needs to be modified to |
| 395 | // reflect this |
| 396 | if ( orig_evt != new_evt ) { |
| 397 | if (new_evt) new_evt->set_barcode( this, barcode() ); |
| 398 | if (orig_evt) orig_evt->remove_barcode( this ); |
| 399 | // we also need to loop over all the particles which are owned by |
| 400 | // this vertex, and remove their barcodes from the old event. |
| 401 | for ( particles_in_const_iterator part1=particles_in_const_begin(); |
| 402 | part1 != particles_in_const_end(); part1++ ) { |
| 403 | if ( !(*part1)->production_vertex() ) { |
| 404 | if ( orig_evt ) orig_evt->remove_barcode( *part1 ); |
| 405 | if ( new_evt ) new_evt->set_barcode( *part1, |
| 406 | (*part1)->barcode() ); |
| 407 | } |
| 408 | } |
| 409 | for ( particles_out_const_iterator |
| 410 | part2 = particles_out_const_begin(); |
| 411 | part2 != particles_out_const_end(); part2++ ) { |
| 412 | if ( orig_evt ) orig_evt->remove_barcode( *part2 ); |
| 413 | if ( new_evt ) new_evt->set_barcode( *part2, |
| 414 | (*part2)->barcode() ); |
| 415 | } |
| 416 | } |
| 417 | } |
| 418 | |
| 419 | void GenVertex::change_parent_event_( GenEvent* new_evt ) |
| 420 | { |
| 421 | // |
| 422 | // this method is for use with swap |
| 423 | // particles and vertices have already been exchanged, |
| 424 | // but the backpointer needs to be fixed |
| 425 | //GenEvent* orig_evt = m_event; |
| 426 | m_event = new_evt; |
| 427 | } |
| 428 | |
| 429 | ///////////// |
| 430 | // Static // |
| 431 | ///////////// |
| 432 | //unsigned int GenVertex::counter() { return s_counter; } |
| 433 | //unsigned int GenVertex::s_counter = 0; |
| 434 | |
| 435 | ///////////// |
| 436 | // Friends // |
| 437 | ///////////// |
| 438 | |
| 439 | /// send vertex information to ostr for printing |
| 440 | std::ostream& operator<<( std::ostream& ostr, const GenVertex& vtx ) { |
| 441 | if ( vtx.barcode()!=0 ) ostr << "BarCode " << vtx.barcode(); |
| 442 | else ostr << "Address " << &vtx; |
| 443 | ostr << " (X,cT)="; |
| 444 | if ( vtx.position() != FourVector(0,0,0,0)) { |
| 445 | ostr << vtx.position().x() << "," |
| 446 | << vtx.position().y() << "," |
| 447 | << vtx.position().z() << "," |
| 448 | << vtx.position().t(); |
| 449 | } else { ostr << 0; } |
| 450 | ostr << " #in:" << vtx.particles_in_size() |
| 451 | << " #out:" << vtx.particles_out_size(); |
| 452 | return ostr; |
| 453 | } |
| 454 | |
| 455 | ///////////////////////////// |
| 456 | // edge_iterator // (protected - for internal use only) |
| 457 | ///////////////////////////// |
| 458 | // If the user wants the functionality of the edge_iterator, he should |
| 459 | // use particle_iterator with IteratorRange = family, parents, or children |
| 460 | // |
| 461 | |
| 462 | GenVertex::edge_iterator::edge_iterator() : m_vertex(0), m_range(family), |
| 463 | m_is_inparticle_iter(false), m_is_past_end(true) |
| 464 | {} |
| 465 | |
| 466 | GenVertex::edge_iterator::edge_iterator( const GenVertex& vtx, |
| 467 | IteratorRange range ) : |
| 468 | m_vertex(&vtx), m_range(family) |
| 469 | { |
| 470 | // Note: (26.1.2000) the original version of edge_iterator inheritted |
| 471 | // from set<GenParticle*>::const_iterator() rather than using |
| 472 | // composition as it does now. |
| 473 | // The inheritted version suffered from a strange bug, which |
| 474 | // I have not fully understood --- it only occurred after many |
| 475 | // events were processed and only when I called the delete |
| 476 | // function on past events. I believe it had something to do with |
| 477 | // the past the end values, which are now robustly coded in this |
| 478 | // version as boolean members. |
| 479 | // |
| 480 | // default range is family, only other choices are children/parents |
| 481 | // descendants/ancestors not allowed & recasted ot children/parents |
| 482 | if ( range == descendants || range == children ) m_range = children; |
| 483 | if ( range == ancestors || range == parents ) m_range = parents; |
| 484 | // |
| 485 | if ( m_vertex->m_particles_in.empty() && |
| 486 | m_vertex->m_particles_out.empty() ) { |
| 487 | // Case: particles_in and particles_out is empty. |
| 488 | m_is_inparticle_iter = false; |
| 489 | m_is_past_end = true; |
| 490 | } else if ( m_range == parents && m_vertex->m_particles_in.empty() ){ |
| 491 | // Case: particles in is empty and parents is requested. |
| 492 | m_is_inparticle_iter = true; |
| 493 | m_is_past_end = true; |
| 494 | } else if ( m_range == children && m_vertex->m_particles_out.empty() ){ |
| 495 | // Case: particles out is empty and children is requested. |
| 496 | m_is_inparticle_iter = false; |
| 497 | m_is_past_end = true; |
| 498 | } else if ( m_range == children ) { |
| 499 | // Case: particles out is NOT empty, and children is requested |
| 500 | m_set_iter = m_vertex->m_particles_out.begin(); |
| 501 | m_is_inparticle_iter = false; |
| 502 | m_is_past_end = false; |
| 503 | } else if ( m_range == family && m_vertex->m_particles_in.empty() ) { |
| 504 | // Case: particles in is empty, particles out is NOT empty, |
| 505 | // and family is requested. Then skip ahead to partilces out. |
| 506 | m_set_iter = m_vertex->m_particles_out.begin(); |
| 507 | m_is_inparticle_iter = false; |
| 508 | m_is_past_end = false; |
| 509 | } else { |
| 510 | // Normal scenario: start with the first incoming particle |
| 511 | m_set_iter = m_vertex->m_particles_in.begin(); |
| 512 | m_is_inparticle_iter = true; |
| 513 | m_is_past_end = false; |
| 514 | } |
| 515 | } |
| 516 | |
| 517 | GenVertex::edge_iterator::edge_iterator( const edge_iterator& p ) { |
| 518 | *this = p; |
| 519 | } |
| 520 | |
| 521 | GenVertex::edge_iterator::~edge_iterator() {} |
| 522 | |
| 523 | GenVertex::edge_iterator& GenVertex::edge_iterator::operator=( |
| 524 | const edge_iterator& p ) { |
| 525 | m_vertex = p.m_vertex; |
| 526 | m_range = p.m_range; |
| 527 | m_set_iter = p.m_set_iter; |
| 528 | m_is_inparticle_iter = p.m_is_inparticle_iter; |
| 529 | m_is_past_end = p.m_is_past_end; |
| 530 | return *this; |
| 531 | } |
| 532 | |
| 533 | GenParticle* GenVertex::edge_iterator::operator*(void) const { |
| 534 | if ( !m_vertex || m_is_past_end ) return 0; |
| 535 | return *m_set_iter; |
| 536 | } |
| 537 | |
| 538 | GenVertex::edge_iterator& GenVertex::edge_iterator::operator++(void){ |
| 539 | // Pre-fix increment |
| 540 | // |
| 541 | // increment the set iterator (unless we're past the end value) |
| 542 | if ( m_is_past_end ) return *this; |
| 543 | ++m_set_iter; |
| 544 | // handle cases where m_set_iter points past the end |
| 545 | if ( m_range == family && m_is_inparticle_iter && |
| 546 | m_set_iter == m_vertex->m_particles_in.end() ) { |
| 547 | // at the end on in particle set, and range is family, so move to |
| 548 | // out particle set |
| 549 | m_set_iter = m_vertex->m_particles_out.begin(); |
| 550 | m_is_inparticle_iter = false; |
| 551 | } else if ( m_range == parents && |
| 552 | m_set_iter == m_vertex->m_particles_in.end() ) { |
| 553 | // at the end on in particle set, and range is parents only, so |
| 554 | // move into past the end state |
| 555 | m_is_past_end = true; |
| 556 | // might as well bail out now |
| 557 | return *this; |
| 558 | } |
| 559 | // are we iterating over input or output particles? |
| 560 | if( m_is_inparticle_iter ) { |
| 561 | // the following is not else if because we might have range=family |
| 562 | // with an empty particles_in set. |
| 563 | if ( m_set_iter == m_vertex->m_particles_in.end() ) { |
| 564 | //whenever out particles end is reached, go into past the end state |
| 565 | m_is_past_end = true; |
| 566 | } |
| 567 | } else { |
| 568 | // the following is not else if because we might have range=family |
| 569 | // with an empty particles_out set. |
| 570 | if ( m_set_iter == m_vertex->m_particles_out.end() ) { |
| 571 | //whenever out particles end is reached, go into past the end state |
| 572 | m_is_past_end = true; |
| 573 | } |
| 574 | } |
| 575 | return *this; |
| 576 | } |
| 577 | |
| 578 | GenVertex::edge_iterator GenVertex::edge_iterator::operator++(int){ |
| 579 | // Post-fix increment |
| 580 | edge_iterator returnvalue = *this; |
| 581 | ++*this; |
| 582 | return returnvalue; |
| 583 | } |
| 584 | |
| 585 | bool GenVertex::edge_iterator::is_parent() const { |
| 586 | if ( **this && (**this)->end_vertex() == m_vertex ) return true; |
| 587 | return false; |
| 588 | } |
| 589 | |
| 590 | bool GenVertex::edge_iterator::is_child() const { |
| 591 | if ( **this && (**this)->production_vertex() == m_vertex ) return true; |
| 592 | return false; |
| 593 | } |
| 594 | |
| 595 | int GenVertex::edges_size( IteratorRange range ) const { |
| 596 | if ( range == children ) return m_particles_out.size(); |
| 597 | if ( range == parents ) return m_particles_in.size(); |
| 598 | if ( range == family ) return m_particles_out.size() |
| 599 | + m_particles_in.size(); |
| 600 | return 0; |
| 601 | } |
| 602 | |
| 603 | ///////////////////// |
| 604 | // vertex_iterator // |
| 605 | ///////////////////// |
| 606 | |
| 607 | GenVertex::vertex_iterator::vertex_iterator() |
| 608 | : m_vertex(0), m_range(), m_visited_vertices(0), m_it_owns_set(0), |
| 609 | m_recursive_iterator(0) |
| 610 | {} |
| 611 | |
| 612 | GenVertex::vertex_iterator::vertex_iterator( GenVertex& vtx_root, |
| 613 | IteratorRange range ) |
| 614 | : m_vertex(&vtx_root), m_range(range) |
| 615 | { |
| 616 | // standard public constructor |
| 617 | // |
| 618 | m_visited_vertices = new std::set<const GenVertex*>; |
| 619 | m_it_owns_set = 1; |
| 620 | m_visited_vertices->insert( m_vertex ); |
| 621 | m_recursive_iterator = 0; |
| 622 | m_edge = m_vertex->edges_begin( m_range ); |
| 623 | // advance to the first good return value |
| 624 | if ( !follow_edge_() && |
| 625 | m_edge != m_vertex->edges_end( m_range )) ++*this; |
| 626 | } |
| 627 | |
| 628 | GenVertex::vertex_iterator::vertex_iterator( GenVertex& vtx_root, |
| 629 | IteratorRange range, std::set<const HepMC::GenVertex*>& visited_vertices ) : |
| 630 | m_vertex(&vtx_root), m_range(range), |
| 631 | m_visited_vertices(&visited_vertices), m_it_owns_set(0), |
| 632 | m_recursive_iterator(0) |
| 633 | { |
| 634 | // This constuctor is only to be called internally by this class |
| 635 | // for use with its recursive member pointer (m_recursive_iterator). |
| 636 | // Note: we do not need to insert m_vertex_root in the vertex - that is |
| 637 | // the responsibility of this iterator's mother, which is normally |
| 638 | // done just before calling this protected constructor. |
| 639 | m_edge = m_vertex->edges_begin( m_range ); |
| 640 | // advance to the first good return value |
| 641 | if ( !follow_edge_() && |
| 642 | m_edge != m_vertex->edges_end( m_range )) ++*this; |
| 643 | } |
| 644 | |
| 645 | GenVertex::vertex_iterator::vertex_iterator( const vertex_iterator& v_iter) |
| 646 | : m_vertex(0), m_visited_vertices(0), m_it_owns_set(0), |
| 647 | m_recursive_iterator(0) |
| 648 | { |
| 649 | *this = v_iter; |
| 650 | } |
| 651 | |
| 652 | GenVertex::vertex_iterator::~vertex_iterator() { |
| 653 | if ( m_recursive_iterator ) delete m_recursive_iterator; |
| 654 | if ( m_it_owns_set ) delete m_visited_vertices; |
| 655 | } |
| 656 | |
| 657 | GenVertex::vertex_iterator& GenVertex::vertex_iterator::operator=( |
| 658 | const vertex_iterator& v_iter ) |
| 659 | { |
| 660 | // Note: when copying a vertex_iterator that is NOT the owner |
| 661 | // of its set container, the pointer to the set is copied. Beware! |
| 662 | // (see copy_with_own_set() if you want a different set pointed to) |
| 663 | // In practise the user never needs to worry |
| 664 | // since such iterators are only intended to be used internally. |
| 665 | // |
| 666 | // destruct data member pointers |
| 667 | if ( m_recursive_iterator ) delete m_recursive_iterator; |
| 668 | m_recursive_iterator = 0; |
| 669 | if ( m_it_owns_set ) delete m_visited_vertices; |
| 670 | m_visited_vertices = 0; |
| 671 | m_it_owns_set = 0; |
| 672 | // copy the target vertex_iterator to this iterator |
| 673 | m_vertex = v_iter.m_vertex; |
| 674 | m_range = v_iter.m_range; |
| 675 | if ( v_iter.m_it_owns_set ) { |
| 676 | // i.e. this vertex will own its set if v_iter points to any |
| 677 | // vertex set regardless of whether v_iter owns the set or not! |
| 678 | m_visited_vertices = |
| 679 | new std::set<const GenVertex*>(*v_iter.m_visited_vertices); |
| 680 | m_it_owns_set = 1; |
| 681 | } else { |
| 682 | m_visited_vertices = v_iter.m_visited_vertices; |
| 683 | m_it_owns_set = 0; |
| 684 | } |
| 685 | // |
| 686 | // Note: m_vertex_root is already included in the set of |
| 687 | // tv_iter.m_visited_vertices, we do not need to insert it. |
| 688 | // |
| 689 | m_edge = v_iter.m_edge; |
| 690 | copy_recursive_iterator_( v_iter.m_recursive_iterator ); |
| 691 | return *this; |
| 692 | } |
| 693 | |
| 694 | GenVertex* GenVertex::vertex_iterator::operator*(void) const { |
| 695 | // de-reference operator |
| 696 | // |
| 697 | // if this iterator has an iterator_node, then we return the iterator |
| 698 | // node. |
| 699 | if ( m_recursive_iterator ) return **m_recursive_iterator; |
| 700 | // |
| 701 | // an iterator can only return its m_vertex -- any other vertex |
| 702 | // is returned by means of a recursive iterator_node |
| 703 | // (so this is the only place in the iterator code that a vertex |
| 704 | // is returned!) |
| 705 | if ( m_vertex ) return m_vertex; |
| 706 | return 0; |
| 707 | } |
| 708 | |
| 709 | GenVertex::vertex_iterator& GenVertex::vertex_iterator::operator++(void) { |
| 710 | // Pre-fix incremental operator |
| 711 | // |
| 712 | // check for "past the end condition" denoted by m_vertex=0 |
| 713 | if ( !m_vertex ) return *this; |
| 714 | // if at the last edge, move into the "past the end condition" |
| 715 | if ( m_edge == m_vertex->edges_end( m_range ) ) { |
| 716 | m_vertex = 0; |
| 717 | return *this; |
| 718 | } |
| 719 | // check to see if we need to create a new recursive iterator by |
| 720 | // following the current edge only if a recursive iterator doesn't |
| 721 | // already exist. If a new recursive_iterator is created, return it. |
| 722 | if ( follow_edge_() ) { |
| 723 | return *this; |
| 724 | } |
| 725 | // |
| 726 | // if a recursive iterator already exists, increment it, and return its |
| 727 | // value (unless the recursive iterator has null root_vertex [its |
| 728 | // root vertex is set to null if it has already returned its root] |
| 729 | // - in which case we delete it) |
| 730 | // and return the vertex pointed to by the edge. |
| 731 | if ( m_recursive_iterator ) { |
| 732 | ++(*m_recursive_iterator); |
| 733 | if ( **m_recursive_iterator ) { |
| 734 | return *this; |
| 735 | } else { |
| 736 | delete m_recursive_iterator; |
| 737 | m_recursive_iterator = 0; |
| 738 | } |
| 739 | } |
| 740 | // |
| 741 | // increment to the next particle edge |
| 742 | ++m_edge; |
| 743 | // if m_edge is at the end, then we have incremented through all |
| 744 | // edges, and it is time to return m_vertex, which we accomplish |
| 745 | // by returning *this |
| 746 | if ( m_edge == m_vertex->edges_end( m_range ) ) return *this; |
| 747 | // otherwise we follow the current edge by recursively ++ing. |
| 748 | return ++(*this); |
| 749 | } |
| 750 | |
| 751 | GenVertex::vertex_iterator GenVertex::vertex_iterator::operator++(int) { |
| 752 | // Post-fix increment |
| 753 | vertex_iterator returnvalue(*this); |
| 754 | ++(*this); |
| 755 | return returnvalue; |
| 756 | } |
| 757 | |
| 758 | void GenVertex::vertex_iterator::copy_with_own_set( |
| 759 | const vertex_iterator& v_iter, |
| 760 | std::set<const HepMC::GenVertex*>& visited_vertices ) { |
| 761 | /// intended for internal use only. (use with care!) |
| 762 | /// this is the same as the operator= method, but it allows the |
| 763 | /// user to specify which set container m_visited_vertices points to. |
| 764 | /// in all cases, this vertex will NOT own its set. |
| 765 | // |
| 766 | // destruct data member pointers |
| 767 | if ( m_recursive_iterator ) delete m_recursive_iterator; |
| 768 | m_recursive_iterator = 0; |
| 769 | if ( m_it_owns_set ) delete m_visited_vertices; |
| 770 | m_visited_vertices = 0; |
| 771 | m_it_owns_set = false; |
| 772 | // copy the target vertex_iterator to this iterator |
| 773 | m_vertex = v_iter.m_vertex; |
| 774 | m_range = v_iter.m_range; |
| 775 | m_visited_vertices = &visited_vertices; |
| 776 | m_it_owns_set = false; |
| 777 | m_edge = v_iter.m_edge; |
| 778 | copy_recursive_iterator_( v_iter.m_recursive_iterator ); |
| 779 | } |
| 780 | |
| 781 | GenVertex* GenVertex::vertex_iterator::follow_edge_() { |
| 782 | // follows the edge pointed to by m_edge by creating a |
| 783 | // recursive iterator for it. |
| 784 | // |
| 785 | // if a m_recursive_iterator already exists, |
| 786 | // this routine has nothing to do, |
| 787 | // if there's no m_vertex, there's no point following anything, |
| 788 | // also there's no point trying to follow a null edge. |
| 789 | if ( m_recursive_iterator || !m_vertex || !*m_edge ) return 0; |
| 790 | // |
| 791 | // if the range is parents, children, or family (i.e. <= family) |
| 792 | // then only the iterator which owns the set is allowed to create |
| 793 | // recursive iterators (i.e. recursivity is only allowed to go one |
| 794 | // layer deep) |
| 795 | if ( m_range <= family && m_it_owns_set == 0 ) return 0; |
| 796 | // |
| 797 | // M.Dobbs 2001-07-16 |
| 798 | // Take care of the very special-rare case where a particle might |
| 799 | // point to the same vertex for both production and end |
| 800 | if ( (*m_edge)->production_vertex() == |
| 801 | (*m_edge)->end_vertex() ) return 0; |
| 802 | // |
| 803 | // figure out which vertex m_edge is pointing to |
| 804 | GenVertex* vtx = ( m_edge.is_parent() ? |
| 805 | (*m_edge)->production_vertex() : |
| 806 | (*m_edge)->end_vertex() ); |
| 807 | // if the pointed to vertex doesn't exist or has already been visited, |
| 808 | // then return null |
| 809 | if ( !vtx || !(m_visited_vertices->insert(vtx).second) ) return 0; |
| 810 | // follow that edge by creating a recursive iterator |
| 811 | m_recursive_iterator = new vertex_iterator( *vtx, m_range, |
| 812 | *m_visited_vertices); |
| 813 | // and return the vertex pointed to by m_recursive_iterator |
| 814 | return **m_recursive_iterator; |
| 815 | } |
| 816 | |
| 817 | void GenVertex::vertex_iterator::copy_recursive_iterator_( |
| 818 | const vertex_iterator* recursive_v_iter ) { |
| 819 | // to properly copy the recursive iterator, we need to ensure |
| 820 | // the proper set container is transfered ... then do this |
| 821 | // operation .... you guessed it .... recursively! |
| 822 | // |
| 823 | if ( !recursive_v_iter ) return; |
| 824 | m_recursive_iterator = new vertex_iterator(); |
| 825 | m_recursive_iterator->m_vertex = recursive_v_iter->m_vertex; |
| 826 | m_recursive_iterator->m_range = recursive_v_iter->m_range; |
| 827 | m_recursive_iterator->m_visited_vertices = m_visited_vertices; |
| 828 | m_recursive_iterator->m_it_owns_set = 0; |
| 829 | m_recursive_iterator->m_edge = recursive_v_iter->m_edge; |
| 830 | m_recursive_iterator->copy_recursive_iterator_( |
| 831 | recursive_v_iter->m_recursive_iterator ); |
| 832 | } |
| 833 | |
| 834 | /////////////////////////////// |
| 835 | // particle_iterator // |
| 836 | /////////////////////////////// |
| 837 | |
| 838 | GenVertex::particle_iterator::particle_iterator() {} |
| 839 | |
| 840 | GenVertex::particle_iterator::particle_iterator( GenVertex& vertex_root, |
| 841 | IteratorRange range ) { |
| 842 | // General Purpose Constructor |
| 843 | // |
| 844 | if ( range <= family ) { |
| 845 | m_edge = GenVertex::edge_iterator( vertex_root, range ); |
| 846 | } else { |
| 847 | m_vertex_iterator = GenVertex::vertex_iterator(vertex_root, range); |
| 848 | m_edge = GenVertex::edge_iterator( **m_vertex_iterator, |
| 849 | m_vertex_iterator.range() ); |
| 850 | } |
| 851 | advance_to_first_(); |
| 852 | } |
| 853 | |
| 854 | GenVertex::particle_iterator::particle_iterator( |
| 855 | const particle_iterator& p_iter ){ |
| 856 | *this = p_iter; |
| 857 | } |
| 858 | |
| 859 | GenVertex::particle_iterator::~particle_iterator() {} |
| 860 | |
| 861 | GenVertex::particle_iterator& |
| 862 | GenVertex::particle_iterator::operator=( const particle_iterator& p_iter ) |
| 863 | { |
| 864 | m_vertex_iterator = p_iter.m_vertex_iterator; |
| 865 | m_edge = p_iter.m_edge; |
| 866 | return *this; |
| 867 | } |
| 868 | |
| 869 | GenParticle* GenVertex::particle_iterator::operator*(void) const { |
| 870 | return *m_edge; |
| 871 | } |
| 872 | |
| 873 | GenVertex::particle_iterator& |
| 874 | GenVertex::particle_iterator::operator++(void) { |
| 875 | //Pre-fix increment |
| 876 | // |
| 877 | if ( *m_edge ) { |
| 878 | ++m_edge; |
| 879 | } else if ( *m_vertex_iterator ) { // !*m_edge is implicit |
| 880 | // past end of edge, but still have more vertices to visit |
| 881 | // increment the vertex, checking that the result is valid |
| 882 | if ( !*(++m_vertex_iterator) ) return *this; |
| 883 | m_edge = GenVertex::edge_iterator( **m_vertex_iterator, |
| 884 | m_vertex_iterator.range() ); |
| 885 | } else { // !*m_edge and !*m_vertex_iterator are implicit |
| 886 | // past the end condition: do nothing |
| 887 | return *this; |
| 888 | } |
| 889 | advance_to_first_(); |
| 890 | return *this; |
| 891 | } |
| 892 | |
| 893 | GenVertex::particle_iterator GenVertex::particle_iterator::operator++(int){ |
| 894 | //Post-fix increment |
| 895 | particle_iterator returnvalue(*this); |
| 896 | ++(*this); |
| 897 | return returnvalue; |
| 898 | } |
| 899 | |
| 900 | GenParticle* GenVertex::particle_iterator::advance_to_first_() { |
| 901 | /// if the current edge is not a suitable return value ( because |
| 902 | /// it is a parent of the vertex root that itself belongs to a |
| 903 | /// different vertex ) it advances to the first suitable return value |
| 904 | if ( !*m_edge ) return *(++*this); |
| 905 | // if the range is relatives, we need to uniquely assign each particle |
| 906 | // to a single vertex so as to guarantee particles are returned |
| 907 | // exactly once. |
| 908 | if ( m_vertex_iterator.range() == relatives && |
| 909 | m_edge.is_parent() && |
| 910 | (*m_edge)->production_vertex() ) return *(++*this); |
| 911 | return *m_edge; |
| 912 | } |
| 913 | |
| 914 | /// scale the position vector |
| 915 | /// this method is only for use by GenEvent |
| 916 | /// convert_position assumes that 4th component of the position vector |
| 917 | /// is ctau rather than time and has units of length-time |
| 918 | void GenVertex::convert_position( const double& f ) { |
| 919 | m_position = FourVector( f*m_position.x(), |
| 920 | f*m_position.y(), |
| 921 | f*m_position.z(), |
| 922 | f*m_position.t() ); |
| 923 | } |
| 924 | |
| 925 | } // HepMC |
git://git.uio.no / u / mrichter / AliRoot.git / blame