hepmc - Blame information for rev 25

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