Subversion Repositories hepmc

Rev

Rev 120 | Rev 168 | Go to most recent revision | Only display areas with differences | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 120 Rev 140
1
//////////////////////////////////////////////////////////////////////////
1
//////////////////////////////////////////////////////////////////////////
2
// Matt.Dobbs@Cern.CH, January 2000
2
// Matt.Dobbs@Cern.CH, January 2000
3
// particle's flow object
3
// particle's flow object
4
//////////////////////////////////////////////////////////////////////////
4
//////////////////////////////////////////////////////////////////////////
5
5
6
#include "HepMC/Flow.h"
6
#include "HepMC/Flow.h"
7
#include "HepMC/GenParticle.h"
7
#include "HepMC/GenParticle.h"
8
#include "HepMC/GenParticleComparison.h"
8
#include "HepMC/GenParticleComparison.h"
9
#include "HepMC/GenVertex.h"
9
#include "HepMC/GenVertex.h"
10
10
11
namespace HepMC {
11
namespace HepMC {
12
12
13
    Flow::Flow( GenParticle* particle_owner )
13
    Flow::Flow( GenParticle* particle_owner )
14
        : m_particle_owner(particle_owner)
14
        : m_particle_owner(particle_owner)
15
    {}
15
    {}
16
16
17
    Flow::Flow( const Flow& inflow ) :
17
    Flow::Flow( const Flow& inflow ) :
18
        m_particle_owner(inflow.m_particle_owner)
18
        m_particle_owner(inflow.m_particle_owner),
-
 
19
        m_icode(inflow.m_icode)
19
    {
20
    {
20
        /// copies both the m_icode AND the m_particle_owner
21
        /// copies both the m_icode AND the m_particle_owner
21
        *this = inflow;
-
 
22
    }
22
    }
23
23
24
    Flow::~Flow() {
24
    Flow::~Flow() {
25
        m_icode.clear();
25
        m_icode.clear();
-
 
26
    }
-
 
27
-
 
28
    void Flow::swap( Flow & other)
-
 
29
    {
-
 
30
        std::swap( m_particle_owner, other.m_particle_owner );
-
 
31
        m_icode.swap( other.m_icode );
26
    }
32
    }
27
33
28
    void Flow::print( std::ostream& ostr ) const {
34
    void Flow::print( std::ostream& ostr ) const {
29
        ostr << "Flow(" << m_particle_owner << "): " << *this << std::endl;
35
        ostr << "Flow(" << m_particle_owner << "): " << *this << std::endl;
30
    }
36
    }
31
   
37
   
32
    std::set<GenParticle*,GenParticleComparison> Flow::connected_partners( int code, int code_index,
38
    std::set<GenParticle*,GenParticleComparison> Flow::connected_partners( int code, int code_index,
33
                                                  int num_indices ) const {
39
                                                  int num_indices ) const {
34
        /// Returns all flow partners which have "code" in any  of the 
40
        /// Returns all flow partners which have "code" in any  of the 
35
        ///  num_indices beginning with index code_index.
41
        ///  num_indices beginning with index code_index.
36
        ///  m_particle_owner is included in the result.
42
        ///  m_particle_owner is included in the result.
37
        ///  Return is by value since the set should never be very big.
43
        ///  Return is by value since the set should never be very big.
38
        /// EXAMPLE: if you want to find all flow partners that have the same
44
        /// EXAMPLE: if you want to find all flow partners that have the same
39
        ///   code in indices 2,3,4 as particle p has in index 2, you would use:
45
        ///   code in indices 2,3,4 as particle p has in index 2, you would use:
40
        ///   set<GenParticle*> result = 
46
        ///   set<GenParticle*> result = 
41
        ///             p->flow().connected_partners(p->flow().icode(2),2,3);
47
        ///             p->flow().connected_partners(p->flow().icode(2),2,3);
42
        //
48
        //
43
        std::set<GenParticle*,GenParticleComparison> output;
49
        std::set<GenParticle*,GenParticleComparison> output;
44
        for ( int i = code_index; i!=code_index+num_indices; ++i ) {
50
        for ( int i = code_index; i!=code_index+num_indices; ++i ) {
45
            if ( icode(i)==code ) {
51
            if ( icode(i)==code ) {
46
                output.insert(m_particle_owner);
52
                output.insert(m_particle_owner);
47
                connected_partners( &output, code, code_index, num_indices );
53
                connected_partners( &output, code, code_index, num_indices );
48
                break;
54
                break;
49
            }
55
            }
50
        }
56
        }
51
        return output;
57
        return output;
52
    }
58
    }
53
59
54
    void Flow::connected_partners( std::set<GenParticle*,GenParticleComparison>* output, int code,
60
    void Flow::connected_partners( std::set<GenParticle*,GenParticleComparison>* output, int code,
55
                                   int code_index, int num_indices ) const
61
                                   int code_index, int num_indices ) const
56
    {
62
    {
57
        /// protected: for recursive use by Flow::connected_partners()
63
        /// protected: for recursive use by Flow::connected_partners()
58
        //
64
        //
59
        if ( !m_particle_owner ) return; // nothing to do
65
        if ( !m_particle_owner ) return; // nothing to do
60
        // look for connected partners joined to this m_particle_owner
66
        // look for connected partners joined to this m_particle_owner
61
        // through its end_vertex
67
        // through its end_vertex
62
        if ( m_particle_owner->end_vertex() ) {
68
        if ( m_particle_owner->end_vertex() ) {
63
            for ( GenVertex::particle_iterator p
69
            for ( GenVertex::particle_iterator p
64
                     = m_particle_owner->end_vertex()->particles_begin(family);
70
                     = m_particle_owner->end_vertex()->particles_begin(family);
65
                  p != m_particle_owner->end_vertex()->particles_end(family);
71
                  p != m_particle_owner->end_vertex()->particles_end(family);
66
                  ++p ) {
72
                  ++p ) {
67
                // if the particle has the correct flow code and is not yet in 
73
                // if the particle has the correct flow code and is not yet in 
68
                // the set, then we recursively call connected_partners
74
                // the set, then we recursively call connected_partners
69
                for ( int index = code_index; index!=code_index+num_indices;
75
                for ( int index = code_index; index!=code_index+num_indices;
70
                      ++index ){
76
                      ++index ){
71
                    if ( (*p)->flow(index)==code &&
77
                    if ( (*p)->flow(index)==code &&
72
                         output->insert(*p).second ) {
78
                         output->insert(*p).second ) {
73
                        (*p)->flow().connected_partners( output, code,
79
                        (*p)->flow().connected_partners( output, code,
74
                                                         code_index,
80
                                                         code_index,
75
                                                         num_indices );
81
                                                         num_indices );
76
                    }
82
                    }
77
                }
83
                }
78
            }
84
            }
79
        }
85
        }
80
        // same for production_vertex
86
        // same for production_vertex
81
        if ( m_particle_owner->production_vertex() ) {
87
        if ( m_particle_owner->production_vertex() ) {
82
            for ( GenVertex::particle_iterator p
88
            for ( GenVertex::particle_iterator p
83
                      = m_particle_owner->production_vertex()->
89
                      = m_particle_owner->production_vertex()->
84
                      particles_begin( family );
90
                      particles_begin( family );
85
                  p != m_particle_owner->production_vertex()->
91
                  p != m_particle_owner->production_vertex()->
86
                      particles_end( family ); ++p ) {
92
                      particles_end( family ); ++p ) {
87
                // if the particle has the correct flow code and is not yet in 
93
                // if the particle has the correct flow code and is not yet in 
88
                // the set, then we recursively call connected_partners
94
                // the set, then we recursively call connected_partners
89
                for ( int index = code_index; index!=code_index+num_indices;
95
                for ( int index = code_index; index!=code_index+num_indices;
90
                      ++index ){
96
                      ++index ){
91
                    if ( (*p)->flow(index)==code &&
97
                    if ( (*p)->flow(index)==code &&
92
                         output->insert(*p).second ) {
98
                         output->insert(*p).second ) {
93
                        (*p)->flow().connected_partners( output, code,
99
                        (*p)->flow().connected_partners( output, code,
94
                                                         code_index,
100
                                                         code_index,
95
                                                         num_indices );
101
                                                         num_indices );
96
                    }
102
                    }
97
                }
103
                }
98
            }
104
            }
99
        }
105
        }
100
    }
106
    }
101
107
102
    std::set<GenParticle*,GenParticleComparison> Flow::dangling_connected_partners( int code,
108
    std::set<GenParticle*,GenParticleComparison> Flow::dangling_connected_partners( int code,
103
                                     int code_index, int num_indices ) const {
109
                                     int code_index, int num_indices ) const {
104
        std::set<GenParticle*,GenParticleComparison> output;
110
        std::set<GenParticle*,GenParticleComparison> output;
105
        std::set<GenParticle*,GenParticleComparison> visited_particles;
111
        std::set<GenParticle*,GenParticleComparison> visited_particles;
106
        for ( int i = code_index; i!=code_index+num_indices; ++i ) {
112
        for ( int i = code_index; i!=code_index+num_indices; ++i ) {
107
            if ( icode(i)==code ) {
113
            if ( icode(i)==code ) {
108
                visited_particles.insert(m_particle_owner);
114
                visited_particles.insert(m_particle_owner);
109
                dangling_connected_partners( &output, &visited_particles, code,
115
                dangling_connected_partners( &output, &visited_particles, code,
110
                                             code_index, num_indices );
116
                                             code_index, num_indices );
111
                break;
117
                break;
112
            }
118
            }
113
        }
119
        }
114
        return output;
120
        return output;
115
    }
121
    }
116
122
117
    void Flow::dangling_connected_partners( std::set<GenParticle*,GenParticleComparison>* output,
123
    void Flow::dangling_connected_partners( std::set<GenParticle*,GenParticleComparison>* output,
118
                                            std::set<GenParticle*,GenParticleComparison>*
124
                                            std::set<GenParticle*,GenParticleComparison>*
119
                                            visited_particles,
125
                                            visited_particles,
120
                                            int code, int code_index,
126
                                            int code, int code_index,
121
                                            int num_indices ) const
127
                                            int num_indices ) const
122
    {
128
    {
123
        /// protected: for recursive use by Flow::dangling_connected_partners
129
        /// protected: for recursive use by Flow::dangling_connected_partners
124
        //
130
        //
125
        if ( !m_particle_owner ) return; // nothing to do
131
        if ( !m_particle_owner ) return; // nothing to do
126
        int count_partners = 0;
132
        int count_partners = 0;
127
        // look for connected partners joined to this m_particle_owner
133
        // look for connected partners joined to this m_particle_owner
128
        // through its end_vertex
134
        // through its end_vertex
129
        if ( m_particle_owner->end_vertex() ) {
135
        if ( m_particle_owner->end_vertex() ) {
130
            for ( GenVertex::particle_iterator p
136
            for ( GenVertex::particle_iterator p
131
                     = m_particle_owner->end_vertex()->particles_begin(family);
137
                     = m_particle_owner->end_vertex()->particles_begin(family);
132
                  p != m_particle_owner->end_vertex()->particles_end(family);
138
                  p != m_particle_owner->end_vertex()->particles_end(family);
133
                  ++p ) {
139
                  ++p ) {
134
                // if the particle has the correct flow code and is not yet in 
140
                // if the particle has the correct flow code and is not yet in 
135
                // the set, then we recursively call connected_partners
141
                // the set, then we recursively call connected_partners
136
                for ( int index = code_index; index!=code_index+num_indices;
142
                for ( int index = code_index; index!=code_index+num_indices;
137
                      ++index ){
143
                      ++index ){
138
                    if ( (*p)->flow(index)==code ) {
144
                    if ( (*p)->flow(index)==code ) {
139
                        if ( *p!=m_particle_owner ) ++count_partners;
145
                        if ( *p!=m_particle_owner ) ++count_partners;
140
                        if ( visited_particles->insert(*p).second ) {
146
                        if ( visited_particles->insert(*p).second ) {
141
                            (*p)->flow().dangling_connected_partners( output,
147
                            (*p)->flow().dangling_connected_partners( output,
142
                                                     visited_particles, code,
148
                                                     visited_particles, code,
143
                                                     code_index, num_indices );
149
                                                     code_index, num_indices );
144
150
145
                        }
151
                        }
146
                    }          
152
                    }          
147
                }
153
                }
148
            }
154
            }
149
        }
155
        }
150
        // same for production_vertex
156
        // same for production_vertex
151
        if ( m_particle_owner->production_vertex() ) {
157
        if ( m_particle_owner->production_vertex() ) {
152
            for ( GenVertex::particle_iterator p = m_particle_owner->
158
            for ( GenVertex::particle_iterator p = m_particle_owner->
153
                                                production_vertex()->
159
                                                production_vertex()->
154
                                                particles_begin( family );
160
                                                particles_begin( family );
155
                  p != m_particle_owner->production_vertex()->
161
                  p != m_particle_owner->production_vertex()->
156
                                                particles_end( family );
162
                                                particles_end( family );
157
                  ++p ) {
163
                  ++p ) {
158
                // if the particle has the correct flow code and is not yet in 
164
                // if the particle has the correct flow code and is not yet in 
159
                // the set, then we recursively call connected_partners
165
                // the set, then we recursively call connected_partners
160
                for ( int index = code_index; index!=code_index+num_indices;
166
                for ( int index = code_index; index!=code_index+num_indices;
161
                      ++index ){
167
                      ++index ){
162
                    if ( (*p)->flow(index)==code ) {
168
                    if ( (*p)->flow(index)==code ) {
163
                        if ( *p!=m_particle_owner ) ++count_partners;
169
                        if ( *p!=m_particle_owner ) ++count_partners;
164
                        if ( visited_particles->insert(*p).second ) {
170
                        if ( visited_particles->insert(*p).second ) {
165
                            (*p)->flow().dangling_connected_partners( output,
171
                            (*p)->flow().dangling_connected_partners( output,
166
                                                     visited_particles, code,
172
                                                     visited_particles, code,
167
                                                     code_index, num_indices );
173
                                                     code_index, num_indices );
168
174
169
                        }
175
                        }
170
                    }
176
                    }
171
                }
177
                }
172
            }
178
            }
173
        }
179
        }
174
        if ( count_partners <= 1 ) output->insert( m_particle_owner );
180
        if ( count_partners <= 1 ) output->insert( m_particle_owner );
175
    }
181
    }
176
       
182
       
177
    /////////////
183
    /////////////
178
    // Friends //
184
    // Friends //
179
    /////////////
185
    /////////////
180
186
181
        /// send Flow informatin to ostr for printing
187
        /// send Flow informatin to ostr for printing
182
    std::ostream& operator<<( std::ostream& ostr, const Flow& f ) {
188
    std::ostream& operator<<( std::ostream& ostr, const Flow& f ) {
183
        ostr << f.m_icode.size();
189
        ostr << f.m_icode.size();
184
        for ( std::map<int,int>::const_iterator i = f.m_icode.begin();
190
        for ( std::map<int,int>::const_iterator i = f.m_icode.begin();
185
              i != f.m_icode.end(); ++i ) {
191
              i != f.m_icode.end(); ++i ) {
186
            ostr << " " << (*i).first << " " << (*i).second;
192
            ostr << " " << (*i).first << " " << (*i).second;
187
        }
193
        }
188
        return ostr;
194
        return ostr;
189
    }
195
    }
190
196
191
} // HepMC
197
} // HepMC
192
198
193
199
194
200
195
201
196
202
197
203
198
204
199
205
200
206
201
207
202
208
203
209
204
210
205
211
206
212
207
213
208
214
209
215
210
 
216