1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
|
//////////////////////////////////////////////////////////////////////////
// Matt.Dobbs@Cern.CH, Feb 2000
// This example shows low to use the particle and vertex iterators
//////////////////////////////////////////////////////////////////////////
// To Compile: go to the HepMC directory and type:
// gmake examples/example_UsingIterators.exe
//
#include "HepMC/IO_GenEvent.h"
#include "HepMC/GenEvent.h"
#include <math.h>
#include <algorithm>
#include <list>
//! example class
/// \class IsPhoton
/// this predicate returns true if the input particle is a photon
/// in the central region (eta < 2.5) with pT > 10 GeV
class IsPhoton {
public:
/// returns true if the GenParticle is a photon with more than 10 GeV transverse momentum
bool operator()( const HepMC::GenParticle* p ) {
if ( p->pdg_id() == 22
&& p->momentum().perp() > 10. ) return 1;
return 0;
}
};
//! example class
/// \class IsW_Boson
/// this predicate returns true if the input particle is a W+/W-
class IsW_Boson {
public:
/// returns true if the GenParticle is a W
bool operator()( const HepMC::GenParticle* p ) {
if ( abs(p->pdg_id()) == 24 ) return 1;
return 0;
}
};
//! example class
/// \class IsStateFinal
/// this predicate returns true if the input has no decay vertex
class IsStateFinal {
public:
/// returns true if the GenParticle does not decay
bool operator()( const HepMC::GenParticle* p ) {
if ( !p->end_vertex() && p->status()==1 ) return 1;
return 0;
}
};
int main() {
{ // begin scope of ascii_in
// an event has been prepared in advance for this example, read it
// into memory using the IO_GenEvent input strategy
HepMC::IO_GenEvent ascii_in("example_UsingIterators.txt",std::ios::in);
if ( ascii_in.rdstate() == std::ios::failbit ) {
std::cerr << "ERROR input file example_UsingIterators.txt is needed "
<< "and does not exist. "
<< "\n Look for it in HepMC/examples, Exit." << std::endl;
return 1;
}
HepMC::GenEvent* evt = ascii_in.read_next_event();
// if you wish to have a look at the event, then use evt->print();
// use GenEvent::vertex_iterator to fill a list of all
// vertices in the event
std::list<HepMC::GenVertex*> allvertices;
for ( HepMC::GenEvent::vertex_iterator v = evt->vertices_begin();
v != evt->vertices_end(); ++v ) {
allvertices.push_back(*v);
}
// we could do the same thing with the STL algorithm copy
std::list<HepMC::GenVertex*> allvertices2;
copy( evt->vertices_begin(), evt->vertices_end(),
back_inserter(allvertices2) );
// fill a list of all final state particles in the event, by requiring
// that each particle satisfyies the IsStateFinal predicate
IsStateFinal isfinal;
std::list<HepMC::GenParticle*> finalstateparticles;
for ( HepMC::GenEvent::particle_iterator p = evt->particles_begin();
p != evt->particles_end(); ++p ) {
if ( isfinal(*p) ) finalstateparticles.push_back(*p);
}
// an STL-like algorithm called HepMC::copy_if is provided in the
// GenEvent.h header to do this sort of operation more easily,
// you could get the identical results as above by using:
std::list<HepMC::GenParticle*> finalstateparticles2;
HepMC::copy_if( evt->particles_begin(), evt->particles_end(),
back_inserter(finalstateparticles2), IsStateFinal() );
// lets print all photons in the event that satisfy the IsPhoton criteria
IsPhoton isphoton;
for ( HepMC::GenEvent::particle_iterator p = evt->particles_begin();
p != evt->particles_end(); ++p ) {
if ( isphoton(*p) ) (*p)->print();
}
// the GenVertex::particle_iterator and GenVertex::vertex_iterator
// are slightly different from the GenEvent:: versions, in that
// the iterator starts at the given vertex, and walks through the attached
// vertex returning particles/vertices.
// Thus only particles/vertices which are in the same graph as the given
// vertex will be returned. A range is specified with these iterators,
// the choices are:
// parents, children, family, ancestors, descendants, relatives
// here are some examples.
// use GenEvent::particle_iterator to find all W's in the event,
// then
// (1) for each W user the GenVertex::particle_iterator with a range of
// parents to return and print the immediate mothers of these W's.
// (2) for each W user the GenVertex::particle_iterator with a range of
// descendants to return and print all descendants of these W's.
IsW_Boson isw;
for ( HepMC::GenEvent::particle_iterator p = evt->particles_begin();
p != evt->particles_end(); ++p ) {
if ( isw(*p) ) {
std::cout << "A W boson has been found: " << std::endl;
(*p)->print();
// return all parents
// we do this by pointing to the production vertex of the W
// particle and asking for all particle parents of that vertex
std::cout << "\t Its parents are: " << std::endl;
if ( (*p)->production_vertex() ) {
for ( HepMC::GenVertex::particle_iterator mother
= (*p)->production_vertex()->
particles_begin(HepMC::parents);
mother != (*p)->production_vertex()->
particles_end(HepMC::parents);
++mother ) {
std::cout << "\t";
(*mother)->print();
}
}
// return all descendants
// we do this by pointing to the end vertex of the W
// particle and asking for all particle descendants of that vertex
std::cout << "\t\t Its descendants are: " << std::endl;
if ( (*p)->end_vertex() ) {
for ( HepMC::GenVertex::particle_iterator des
=(*p)->end_vertex()->
particles_begin(HepMC::descendants);
des != (*p)->end_vertex()->
particles_end(HepMC::descendants);
++des ) {
std::cout << "\t\t";
(*des)->print();
}
}
}
}
// cleanup
delete evt;
// in analogy to the above, similar use can be made of the
// HepMC::GenVertex::vertex_iterator, which also accepts a range.
} // end scope of ascii_in
return 0;
}
|
