FastJet  3.1.0-beta.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
JetSplitMergeTool.cc
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