FastJet 3.0.2
|
00001 //---------------------------------------------------------------------- 00002 // This file distributed with FastJet has been obtained from SpartyJet 00003 // v2.20.0 by Pierre-Antoine Delsart, Kurtis L. Geerlings, Joey 00004 // Huston, Brian T. Martin and Chris Vermilion 00005 // For details, see http://www.pa.msu.edu/~huston/SpartyJet/ 00006 // http://projects.hepforge.org/spartyjet/ 00007 // 00008 // Changes from the original file are listed below. 00009 //---------------------------------------------------------------------- 00010 00011 //******************************************************************************* 00012 // Filename : JetSplitMergeTool.cxx 00013 // Author : Ambreesh Gupta 00014 // Created : Nov, 2001 00015 // 00016 // File taken from SpartyJet v2.20.0 00017 // Modifications: 00018 // removed the string name in the ctor 00019 // removed the Message m_log 00020 // replaced px() -> px, ... in the LorentzVector calls 00021 // cleaned the comments 00022 //******************************************************************************* 00023 00024 // History of changes from the original JetSplitMergeTool.cc file in 00025 // SpartyJet v2.20 00026 // 00027 // 2009-01-15 Gregory Soyez <soyez@fastjet.fr> 00028 // 00029 // * put the code in the fastjet::atlas namespace 00030 // 00031 // 2009-02-14 Gregory Soyez <soyez@fastjet.fr> 00032 // 00033 // * imported into FastJet 00034 // * removed the string name in the ctor 00035 // * removed the Message m_log 00036 // * replaced px() -> px, ... in the LorentzVector calls 00037 // * cleaned the comments 00038 00039 #include <iostream> 00040 00041 #include "JetSplitMergeTool.hh" 00042 #include "Jet.hh" 00043 #include "JetDistances.hh" 00044 #include "CommonUtils.hh" 00045 00046 #include <fastjet/internal/base.hh> 00047 00048 FASTJET_BEGIN_NAMESPACE 00049 00050 namespace atlas { 00051 00052 JetSplitMergeTool::JetSplitMergeTool() 00053 : m_f( 0.5 ) 00054 {} 00055 00056 JetSplitMergeTool::~JetSplitMergeTool() 00057 {} 00058 00059 ///////////////////////////////////////////////////////////////////////////////// 00060 //Execution / 00061 ///////////////////////////////////////////////////////////////////////////////// 00062 int JetSplitMergeTool::execute( jetcollection_t* theJets ) 00063 { 00064 m_ctr = 0; 00065 m_dctr = 0; 00066 00067 //////////////////////////////////////////////////// 00068 // From the input, collection create a list of Jet// 00069 //////////////////////////////////////////////////// 00070 m_preJet.clear(); 00071 m_jet.clear(); 00072 00073 jetcollection_t::iterator itrB = theJets->begin(); 00074 jetcollection_t::iterator itrE = theJets->end(); 00075 00076 double etot =0.; 00077 for (;itrB!=itrE;itrB++) { 00078 Jet* j = new Jet(); j->addJet(*itrB); 00079 m_ctr +=1; 00080 m_preJet.push_back(j); 00081 00082 etot += j->e(); 00083 } 00084 00085 ///////////////////// 00086 // Split Merge Jets// 00087 ///////////////////// 00088 this->split_merge(); 00089 00090 ///////////////////////////////////////////// 00091 // Empty and re-fill input jetcollection_t // 00092 ///////////////////////////////////////////// 00093 clear_list(*theJets); 00094 jetcollection_t::iterator it = m_jet.begin(); 00095 jetcollection_t::iterator itE = m_jet.end(); 00096 for(; it!=itE; ++it){ 00097 theJets->push_back(*it); 00098 } 00099 00100 return 1; 00101 } 00102 00103 /////////////////////////////////////////////////////////////////////////////// 00104 // Reconstruction algorithm specific methods / 00105 /////////////////////////////////////////////////////////////////////////////// 00106 00107 void JetSplitMergeTool::split_merge() 00108 { 00109 if ( m_preJet.size() >= 2 ) { 00110 do { 00111 sort_list_et(m_preJet); 00112 00113 jetcollection_t::iterator itr; 00114 jetcollection_t::iterator first = m_preJet.begin(); 00115 jetcollection_t::iterator last = m_preJet.end(); 00116 00117 itr=first; 00118 ++itr; 00119 bool overlap = false; 00120 00121 for (;itr != last;++itr) { 00122 double etaF = (*first)->eta(); 00123 double phiF = (*first)->phi(); 00124 double etaS = (*itr)->eta(); 00125 double phiS = (*itr)->phi(); 00126 00127 Jet* oJet = jet_from_overlap( (*first),*itr); 00128 m_ctr +=1; 00129 00130 Jet::constit_vect_t::iterator itro = oJet->firstConstituent(); 00131 Jet::constit_vect_t::iterator itroE = oJet->lastConstituent(); 00132 00133 if ( oJet->getConstituentNum() != 0 ) { 00134 overlap = true; 00135 00136 // fraction 00137 double f = sqrt(pow(oJet->px,2)+pow(oJet->py,2))/ 00138 sqrt(pow((*itr)->px,2)+pow((*itr)->py,2)); 00139 00140 // merge 00141 if ( f > m_f) { 00142 // we need to remove constituents ! 00143 Jet *j = (*first); 00144 for ( ;itro != itroE; ++itro ) j->removeConstituent(*itro); 00145 (*first)->addJet(*itr); 00146 //m_preJet.remove(*itr); 00147 delete *itr; 00148 m_preJet.erase(itr); 00149 m_dctr +=1; 00150 } 00151 00152 // split 00153 if ( f <= m_f) { 00154 for ( ;itro != itroE; ++itro ) { 00155 // Distance of first jet from ProtoJet 00156 double deta1 = etaF - (*itro)->eta(); 00157 double dphi1 = fabs(JetDistances::deltaPhi(phiF,(*itro)->phi())); 00158 double dist1 = pow( deta1 , 2 ) + pow( dphi1 , 2 ); 00159 00160 // Distance of second jet from ProtoJet 00161 double deta2 = etaS - (*itro)->eta(); 00162 double dphi2 = fabs(JetDistances::deltaPhi(phiS,(*itro)->phi())); 00163 double dist2 = pow( deta2 , 2 ) + pow( dphi2 , 2 ); 00164 00165 // Remove protojet from farther Jet 00166 if ( dist1 > dist2 ) (*first)->removeConstituent(*itro); 00167 if ( dist1 <= dist2 ) (*itr)->removeConstituent(*itro); 00168 } 00169 } 00170 // Delete overlap jet 00171 delete oJet; 00172 m_dctr +=1; 00173 break; 00174 } 00175 else { 00176 // Delete overlap jet 00177 delete oJet; 00178 m_dctr +=1; 00179 } 00180 } 00181 00182 if ( overlap == false ) { 00183 m_jet.push_back(*first); 00184 //m_preJet.remove(*first); 00185 m_preJet.erase(first); 00186 } 00187 00188 } while ( m_preJet.size() != 0 ); 00189 } 00190 else if ( m_preJet.size() == 1) { 00191 m_jet.push_back( *(m_preJet.begin()) ); 00192 } 00193 00194 } 00195 00196 ////////////////////////////////////////////////////////////////////// 00197 00198 // "True" eta and phi ASSUMING the 4-vector is filled as 00199 // ex -> e * sin(theta) * cos(phi) 00200 // ey -> e * sin(theta) * sin(phi) 00201 // ez -> e * cos(theta) 00202 // e -> e 00203 // Jet phi range is (-pi,pi]. 00204 00205 00206 double JetSplitMergeTool::etaTrue(Jet::constit_vect_t::iterator pj) 00207 { 00208 double s = ((*pj)->e() > 0) ? +1.0 : -1.0; 00209 double px = (*pj)->px; 00210 double py = (*pj)->py; 00211 double pz = (*pj)->pz; 00212 double theta = acos(pz*s/sqrt(px*px+py*py+pz*pz)); 00213 return -log(tan(theta/2.0)); 00214 } 00215 00216 double JetSplitMergeTool::phiTrue(Jet::constit_vect_t::iterator pj) 00217 { 00218 double s = ((*pj)->e() > 0) ? +1.0 : -1.0; 00219 double px = (*pj)->px; 00220 double py = (*pj)->py; 00221 return atan2(py*s,px*s); 00222 } 00223 00224 } // namespace atlas 00225 00226 FASTJET_END_NAMESPACE