![]() |
HepMC Reference DocumentationHepMC |
00001 00002 // Matt.Dobbs@Cern.CH, Feb 2000 00003 // This example shows low to use the particle and vertex iterators 00005 // To Compile: go to the HepMC directory and type: 00006 // gmake examples/example_UsingIterators.exe 00007 // 00008 00009 #include "HepMC/IO_GenEvent.h" 00010 #include "HepMC/GenEvent.h" 00011 #include <math.h> 00012 #include <algorithm> 00013 #include <list> 00014 00016 00020 class IsPhoton { 00021 public: 00023 bool operator()( const HepMC::GenParticle* p ) { 00024 if ( p->pdg_id() == 22 00025 && p->momentum().perp() > 10. ) return 1; 00026 return 0; 00027 } 00028 }; 00029 00031 00034 class IsW_Boson { 00035 public: 00037 bool operator()( const HepMC::GenParticle* p ) { 00038 if ( abs(p->pdg_id()) == 24 ) return 1; 00039 return 0; 00040 } 00041 }; 00042 00044 00047 class IsStateFinal { 00048 public: 00050 bool operator()( const HepMC::GenParticle* p ) { 00051 if ( !p->end_vertex() && p->status()==1 ) return 1; 00052 return 0; 00053 } 00054 }; 00055 00056 int main() { 00057 { // begin scope of ascii_in 00058 // an event has been prepared in advance for this example, read it 00059 // into memory using the IO_GenEvent input strategy 00060 HepMC::IO_GenEvent ascii_in("example_UsingIterators.txt",std::ios::in); 00061 if ( ascii_in.rdstate() == std::ios::failbit ) { 00062 std::cerr << "ERROR input file example_UsingIterators.txt is needed " 00063 << "and does not exist. " 00064 << "\n Look for it in HepMC/examples, Exit." << std::endl; 00065 return 1; 00066 } 00067 00068 HepMC::GenEvent* evt = ascii_in.read_next_event(); 00069 00070 // if you wish to have a look at the event, then use evt->print(); 00071 00072 // use GenEvent::vertex_iterator to fill a list of all 00073 // vertices in the event 00074 std::list<HepMC::GenVertex*> allvertices; 00075 for ( HepMC::GenEvent::vertex_iterator v = evt->vertices_begin(); 00076 v != evt->vertices_end(); ++v ) { 00077 allvertices.push_back(*v); 00078 } 00079 00080 // we could do the same thing with the STL algorithm copy 00081 std::list<HepMC::GenVertex*> allvertices2; 00082 copy( evt->vertices_begin(), evt->vertices_end(), 00083 back_inserter(allvertices2) ); 00084 00085 // fill a list of all final state particles in the event, by requiring 00086 // that each particle satisfyies the IsStateFinal predicate 00087 IsStateFinal isfinal; 00088 std::list<HepMC::GenParticle*> finalstateparticles; 00089 for ( HepMC::GenEvent::particle_iterator p = evt->particles_begin(); 00090 p != evt->particles_end(); ++p ) { 00091 if ( isfinal(*p) ) finalstateparticles.push_back(*p); 00092 } 00093 00094 // an STL-like algorithm called HepMC::copy_if is provided in the 00095 // GenEvent.h header to do this sort of operation more easily, 00096 // you could get the identical results as above by using: 00097 std::list<HepMC::GenParticle*> finalstateparticles2; 00098 HepMC::copy_if( evt->particles_begin(), evt->particles_end(), 00099 back_inserter(finalstateparticles2), IsStateFinal() ); 00100 00101 // lets print all photons in the event that satisfy the IsPhoton criteria 00102 IsPhoton isphoton; 00103 for ( HepMC::GenEvent::particle_iterator p = evt->particles_begin(); 00104 p != evt->particles_end(); ++p ) { 00105 if ( isphoton(*p) ) (*p)->print(); 00106 } 00107 00108 // the GenVertex::particle_iterator and GenVertex::vertex_iterator 00109 // are slightly different from the GenEvent:: versions, in that 00110 // the iterator starts at the given vertex, and walks through the attached 00111 // vertex returning particles/vertices. 00112 // Thus only particles/vertices which are in the same graph as the given 00113 // vertex will be returned. A range is specified with these iterators, 00114 // the choices are: 00115 // parents, children, family, ancestors, descendants, relatives 00116 // here are some examples. 00117 00118 // use GenEvent::particle_iterator to find all W's in the event, 00119 // then 00120 // (1) for each W user the GenVertex::particle_iterator with a range of 00121 // parents to return and print the immediate mothers of these W's. 00122 // (2) for each W user the GenVertex::particle_iterator with a range of 00123 // descendants to return and print all descendants of these W's. 00124 IsW_Boson isw; 00125 for ( HepMC::GenEvent::particle_iterator p = evt->particles_begin(); 00126 p != evt->particles_end(); ++p ) { 00127 if ( isw(*p) ) { 00128 std::cout << "A W boson has been found: " << std::endl; 00129 (*p)->print(); 00130 // return all parents 00131 // we do this by pointing to the production vertex of the W 00132 // particle and asking for all particle parents of that vertex 00133 std::cout << "\t Its parents are: " << std::endl; 00134 if ( (*p)->production_vertex() ) { 00135 for ( HepMC::GenVertex::particle_iterator mother 00136 = (*p)->production_vertex()-> 00137 particles_begin(HepMC::parents); 00138 mother != (*p)->production_vertex()-> 00139 particles_end(HepMC::parents); 00140 ++mother ) { 00141 std::cout << "\t"; 00142 (*mother)->print(); 00143 } 00144 } 00145 // return all descendants 00146 // we do this by pointing to the end vertex of the W 00147 // particle and asking for all particle descendants of that vertex 00148 std::cout << "\t\t Its descendants are: " << std::endl; 00149 if ( (*p)->end_vertex() ) { 00150 for ( HepMC::GenVertex::particle_iterator des 00151 =(*p)->end_vertex()-> 00152 particles_begin(HepMC::descendants); 00153 des != (*p)->end_vertex()-> 00154 particles_end(HepMC::descendants); 00155 ++des ) { 00156 std::cout << "\t\t"; 00157 (*des)->print(); 00158 } 00159 } 00160 } 00161 } 00162 // cleanup 00163 delete evt; 00164 // in analogy to the above, similar use can be made of the 00165 // HepMC::GenVertex::vertex_iterator, which also accepts a range. 00166 } // end scope of ascii_in 00167 00168 return 0; 00169 }