1 | //---------------------------------------------------------------------- |
---|
2 | // This file distributed with FastJet has been obtained from SpartyJet |
---|
3 | // v2.20.0 by Pierre-Antoine Delsart, Kurtis L. Geerlings, Joey |
---|
4 | // Huston, Brian T. Martin and Chris Vermilion |
---|
5 | // For details, see http://www.pa.msu.edu/~huston/SpartyJet/ |
---|
6 | // http://projects.hepforge.org/spartyjet/ |
---|
7 | // |
---|
8 | // Changes from the original file are listed below. |
---|
9 | //---------------------------------------------------------------------- |
---|
10 | |
---|
11 | //******************************************************************************* |
---|
12 | // Filename : JetSplitMergeTool.cxx |
---|
13 | // Author : Ambreesh Gupta |
---|
14 | // Created : Nov, 2001 |
---|
15 | // |
---|
16 | // File taken from SpartyJet v2.20.0 |
---|
17 | // Modifications: |
---|
18 | // removed the string name in the ctor |
---|
19 | // removed the Message m_log |
---|
20 | // replaced px() -> px, ... in the LorentzVector calls |
---|
21 | // cleaned the comments |
---|
22 | //******************************************************************************* |
---|
23 | |
---|
24 | // History of changes from the original JetSplitMergeTool.cc file in |
---|
25 | // SpartyJet v2.20 |
---|
26 | // |
---|
27 | // 2009-01-15 Gregory Soyez <soyez@fastjet.fr> |
---|
28 | // |
---|
29 | // * put the code in the fastjet::atlas namespace |
---|
30 | // |
---|
31 | // 2009-02-14 Gregory Soyez <soyez@fastjet.fr> |
---|
32 | // |
---|
33 | // * imported into FastJet |
---|
34 | // * removed the string name in the ctor |
---|
35 | // * removed the Message m_log |
---|
36 | // * replaced px() -> px, ... in the LorentzVector calls |
---|
37 | // * cleaned the comments |
---|
38 | |
---|
39 | #include <iostream> |
---|
40 | |
---|
41 | #include "JetSplitMergeTool.hh" |
---|
42 | #include "Jet.hh" |
---|
43 | #include "JetDistances.hh" |
---|
44 | #include "CommonUtils.hh" |
---|
45 | |
---|
46 | #include <fastjet/internal/base.hh> |
---|
47 | |
---|
48 | FASTJET_BEGIN_NAMESPACE |
---|
49 | |
---|
50 | namespace atlas { |
---|
51 | |
---|
52 | JetSplitMergeTool::JetSplitMergeTool() |
---|
53 | : m_f( 0.5 ) |
---|
54 | {} |
---|
55 | |
---|
56 | JetSplitMergeTool::~JetSplitMergeTool() |
---|
57 | {} |
---|
58 | |
---|
59 | ///////////////////////////////////////////////////////////////////////////////// |
---|
60 | //Execution / |
---|
61 | ///////////////////////////////////////////////////////////////////////////////// |
---|
62 | int JetSplitMergeTool::execute( jetcollection_t* theJets ) |
---|
63 | { |
---|
64 | m_ctr = 0; |
---|
65 | m_dctr = 0; |
---|
66 | |
---|
67 | //////////////////////////////////////////////////// |
---|
68 | // From the input, collection create a list of Jet// |
---|
69 | //////////////////////////////////////////////////// |
---|
70 | m_preJet.clear(); |
---|
71 | m_jet.clear(); |
---|
72 | |
---|
73 | jetcollection_t::iterator itrB = theJets->begin(); |
---|
74 | jetcollection_t::iterator itrE = theJets->end(); |
---|
75 | |
---|
76 | double etot =0.; |
---|
77 | for (;itrB!=itrE;itrB++) { |
---|
78 | Jet* j = new Jet(); j->addJet(*itrB); |
---|
79 | m_ctr +=1; |
---|
80 | m_preJet.push_back(j); |
---|
81 | |
---|
82 | etot += j->e(); |
---|
83 | } |
---|
84 | |
---|
85 | ///////////////////// |
---|
86 | // Split Merge Jets// |
---|
87 | ///////////////////// |
---|
88 | this->split_merge(); |
---|
89 | |
---|
90 | ///////////////////////////////////////////// |
---|
91 | // Empty and re-fill input jetcollection_t // |
---|
92 | ///////////////////////////////////////////// |
---|
93 | clear_list(*theJets); |
---|
94 | jetcollection_t::iterator it = m_jet.begin(); |
---|
95 | jetcollection_t::iterator itE = m_jet.end(); |
---|
96 | for(; it!=itE; ++it){ |
---|
97 | theJets->push_back(*it); |
---|
98 | } |
---|
99 | |
---|
100 | return 1; |
---|
101 | } |
---|
102 | |
---|
103 | /////////////////////////////////////////////////////////////////////////////// |
---|
104 | // Reconstruction algorithm specific methods / |
---|
105 | /////////////////////////////////////////////////////////////////////////////// |
---|
106 | |
---|
107 | void JetSplitMergeTool::split_merge() |
---|
108 | { |
---|
109 | if ( m_preJet.size() >= 2 ) { |
---|
110 | do { |
---|
111 | sort_list_et(m_preJet); |
---|
112 | |
---|
113 | jetcollection_t::iterator itr; |
---|
114 | jetcollection_t::iterator first = m_preJet.begin(); |
---|
115 | jetcollection_t::iterator last = m_preJet.end(); |
---|
116 | |
---|
117 | itr=first; |
---|
118 | ++itr; |
---|
119 | bool overlap = false; |
---|
120 | |
---|
121 | for (;itr != last;++itr) { |
---|
122 | double etaF = (*first)->eta(); |
---|
123 | double phiF = (*first)->phi(); |
---|
124 | double etaS = (*itr)->eta(); |
---|
125 | double phiS = (*itr)->phi(); |
---|
126 | |
---|
127 | Jet* oJet = jet_from_overlap( (*first),*itr); |
---|
128 | m_ctr +=1; |
---|
129 | |
---|
130 | Jet::constit_vect_t::iterator itro = oJet->firstConstituent(); |
---|
131 | Jet::constit_vect_t::iterator itroE = oJet->lastConstituent(); |
---|
132 | |
---|
133 | if ( oJet->getConstituentNum() != 0 ) { |
---|
134 | overlap = true; |
---|
135 | |
---|
136 | // fraction |
---|
137 | double f = sqrt(pow(oJet->px,2)+pow(oJet->py,2))/ |
---|
138 | sqrt(pow((*itr)->px,2)+pow((*itr)->py,2)); |
---|
139 | |
---|
140 | // merge |
---|
141 | if ( f > m_f) { |
---|
142 | // we need to remove constituents ! |
---|
143 | Jet *j = (*first); |
---|
144 | for ( ;itro != itroE; ++itro ) j->removeConstituent(*itro); |
---|
145 | (*first)->addJet(*itr); |
---|
146 | //m_preJet.remove(*itr); |
---|
147 | delete *itr; |
---|
148 | m_preJet.erase(itr); |
---|
149 | m_dctr +=1; |
---|
150 | } |
---|
151 | |
---|
152 | // split |
---|
153 | if ( f <= m_f) { |
---|
154 | for ( ;itro != itroE; ++itro ) { |
---|
155 | // Distance of first jet from ProtoJet |
---|
156 | double deta1 = etaF - (*itro)->eta(); |
---|
157 | double dphi1 = fabs(JetDistances::deltaPhi(phiF,(*itro)->phi())); |
---|
158 | double dist1 = pow( deta1 , 2 ) + pow( dphi1 , 2 ); |
---|
159 | |
---|
160 | // Distance of second jet from ProtoJet |
---|
161 | double deta2 = etaS - (*itro)->eta(); |
---|
162 | double dphi2 = fabs(JetDistances::deltaPhi(phiS,(*itro)->phi())); |
---|
163 | double dist2 = pow( deta2 , 2 ) + pow( dphi2 , 2 ); |
---|
164 | |
---|
165 | // Remove protojet from farther Jet |
---|
166 | if ( dist1 > dist2 ) (*first)->removeConstituent(*itro); |
---|
167 | if ( dist1 <= dist2 ) (*itr)->removeConstituent(*itro); |
---|
168 | } |
---|
169 | } |
---|
170 | // Delete overlap jet |
---|
171 | delete oJet; |
---|
172 | m_dctr +=1; |
---|
173 | break; |
---|
174 | } |
---|
175 | else { |
---|
176 | // Delete overlap jet |
---|
177 | delete oJet; |
---|
178 | m_dctr +=1; |
---|
179 | } |
---|
180 | } |
---|
181 | |
---|
182 | if ( overlap == false ) { |
---|
183 | m_jet.push_back(*first); |
---|
184 | //m_preJet.remove(*first); |
---|
185 | m_preJet.erase(first); |
---|
186 | } |
---|
187 | |
---|
188 | } while ( m_preJet.size() != 0 ); |
---|
189 | } |
---|
190 | else if ( m_preJet.size() == 1) { |
---|
191 | m_jet.push_back( *(m_preJet.begin()) ); |
---|
192 | } |
---|
193 | |
---|
194 | } |
---|
195 | |
---|
196 | ////////////////////////////////////////////////////////////////////// |
---|
197 | |
---|
198 | // "True" eta and phi ASSUMING the 4-vector is filled as |
---|
199 | // ex -> e * sin(theta) * cos(phi) |
---|
200 | // ey -> e * sin(theta) * sin(phi) |
---|
201 | // ez -> e * cos(theta) |
---|
202 | // e -> e |
---|
203 | // Jet phi range is (-pi,pi]. |
---|
204 | |
---|
205 | |
---|
206 | double JetSplitMergeTool::etaTrue(Jet::constit_vect_t::iterator pj) |
---|
207 | { |
---|
208 | double s = ((*pj)->e() > 0) ? +1.0 : -1.0; |
---|
209 | double px = (*pj)->px; |
---|
210 | double py = (*pj)->py; |
---|
211 | double pz = (*pj)->pz; |
---|
212 | double theta = acos(pz*s/sqrt(px*px+py*py+pz*pz)); |
---|
213 | return -log(tan(theta/2.0)); |
---|
214 | } |
---|
215 | |
---|
216 | double JetSplitMergeTool::phiTrue(Jet::constit_vect_t::iterator pj) |
---|
217 | { |
---|
218 | double s = ((*pj)->e() > 0) ? +1.0 : -1.0; |
---|
219 | double px = (*pj)->px; |
---|
220 | double py = (*pj)->py; |
---|
221 | return atan2(py*s,px*s); |
---|
222 | } |
---|
223 | |
---|
224 | } // namespace atlas |
---|
225 | |
---|
226 | FASTJET_END_NAMESPACE |
---|