fastjet 2.4.5
JetSplitMergeTool.cc
Go to the documentation of this file.
00001 //*******************************************************************************
00002 // Filename : JetSplitMergeTool.cxx 
00003 // Author   : Ambreesh Gupta
00004 // Created  : Nov, 2001
00005 //
00006 // File taken from SpartyJet v2.20.0
00007 // Modifications:
00008 //   removed the string name in the ctor
00009 //   removed the Message m_log
00010 //   replaced px() -> px, ... in the LorentzVector calls
00011 //   cleaned the comments
00012 //*******************************************************************************
00013 
00014 #include <iostream>
00015 
00016 #include "JetSplitMergeTool.hh"
00017 #include "Jet.hh"
00018 #include "JetDistances.hh"
00019 #include "CommonUtils.hh"
00020 
00021 #include <fastjet/internal/base.hh>
00022 
00023 FASTJET_BEGIN_NAMESPACE
00024 
00025 namespace atlas { 
00026 
00027 JetSplitMergeTool::JetSplitMergeTool()
00028   :  m_f( 0.5 )
00029 {}
00030 
00031 JetSplitMergeTool::~JetSplitMergeTool()
00032 {}
00033 
00035 //Execution                                                                     /
00037 int JetSplitMergeTool::execute( jetcollection_t* theJets )
00038 {
00039   m_ctr = 0;
00040   m_dctr = 0;
00041 
00043   // From the input, collection create a list of Jet//
00045   m_preJet.clear();
00046   m_jet.clear();
00047 
00048   jetcollection_t::iterator itrB = theJets->begin();
00049   jetcollection_t::iterator itrE = theJets->end(); 
00050 
00051   double etot =0.;
00052   for (;itrB!=itrE;itrB++) {
00053     Jet* j = new Jet(); j->addJet(*itrB);
00054     m_ctr +=1;
00055     m_preJet.push_back(j);    
00056 
00057     etot += j->e();    
00058   }
00059 
00061   // Split Merge Jets//
00063   this->split_merge();
00064  
00066   // Empty and re-fill input jetcollection_t //
00068   clear_list(*theJets);
00069   jetcollection_t::iterator it = m_jet.begin();
00070   jetcollection_t::iterator itE = m_jet.end();
00071   for(; it!=itE; ++it){    
00072     theJets->push_back(*it);
00073   }
00074 
00075   return 1;
00076 }
00077 
00079 // Reconstruction algorithm specific methods                                  /
00081 
00082 void JetSplitMergeTool::split_merge()
00083 {
00084   if ( m_preJet.size() >= 2 ) {
00085     do {
00086       sort_list_et(m_preJet);
00087       
00088       jetcollection_t::iterator itr;
00089       jetcollection_t::iterator first = m_preJet.begin();
00090       jetcollection_t::iterator last  = m_preJet.end();
00091       
00092       itr=first;
00093       ++itr;
00094       bool overlap = false;
00095   
00096       for (;itr != last;++itr) {      
00097         double etaF = (*first)->eta();
00098         double phiF = (*first)->phi();
00099         double etaS = (*itr)->eta();
00100         double phiS = (*itr)->phi();
00101         
00102         Jet* oJet = jet_from_overlap( (*first),*itr);   
00103         m_ctr +=1; 
00104 
00105         Jet::constit_vect_t::iterator itro  = oJet->firstConstituent();
00106         Jet::constit_vect_t::iterator itroE = oJet->lastConstituent();
00107         
00108         if ( oJet->getConstituentNum() != 0 ) {
00109           overlap = true;
00110           
00111           // fraction
00112           double f = sqrt(pow(oJet->px,2)+pow(oJet->py,2))/
00113             sqrt(pow((*itr)->px,2)+pow((*itr)->py,2));
00114           
00115           // merge
00116           if ( f > m_f) {
00117             // we need to remove constituents !
00118             Jet *j = (*first);
00119             for ( ;itro != itroE; ++itro ) j->removeConstituent(*itro);
00120             (*first)->addJet(*itr);
00121             //m_preJet.remove(*itr);
00122             delete *itr;
00123             m_preJet.erase(itr);
00124             m_dctr +=1;
00125           }     
00126           
00127           // split      
00128           if ( f <= m_f) {        
00129             for ( ;itro != itroE; ++itro ) {          
00130               // Distance of first jet from ProtoJet
00131               double deta1 = etaF - (*itro)->eta();
00132               double dphi1 = fabs(JetDistances::deltaPhi(phiF,(*itro)->phi()));
00133               double dist1 = pow( deta1 , 2 ) + pow( dphi1 , 2 );
00134               
00135               // Distance of second jet from ProtoJet
00136               double deta2 = etaS - (*itro)->eta();
00137               double dphi2 = fabs(JetDistances::deltaPhi(phiS,(*itro)->phi()));
00138               double dist2 = pow( deta2 , 2 ) + pow( dphi2 , 2 );
00139               
00140               // Remove protojet from farther Jet                     
00141               if ( dist1 > dist2 ) (*first)->removeConstituent(*itro);
00142               if ( dist1 <= dist2 ) (*itr)->removeConstituent(*itro);   
00143             }
00144           }
00145           // Delete overlap jet     
00146           delete oJet;     
00147           m_dctr +=1;
00148           break; 
00149         }  
00150         else {
00151           // Delete overlap jet     
00152           delete oJet;     
00153           m_dctr +=1;
00154         }    
00155       }
00156       
00157       if ( overlap == false ) {
00158         m_jet.push_back(*first);
00159         //m_preJet.remove(*first);      
00160         m_preJet.erase(first);      
00161       }
00162       
00163     } while ( m_preJet.size() != 0 );    
00164   }
00165   else if ( m_preJet.size() == 1) {
00166     m_jet.push_back( *(m_preJet.begin()) );
00167   }
00168 
00169 }
00170 
00172 
00173 // "True" eta and phi ASSUMING the 4-vector is filled as
00174 // ex -> e * sin(theta) * cos(phi)
00175 // ey -> e * sin(theta) * sin(phi)
00176 // ez -> e * cos(theta)
00177 // e  -> e 
00178 // Jet phi range is (-pi,pi].
00179 
00180 
00181 double JetSplitMergeTool::etaTrue(Jet::constit_vect_t::iterator pj)
00182 {
00183   double s = ((*pj)->e() > 0) ? +1.0 : -1.0;
00184   double px = (*pj)->px;
00185   double py = (*pj)->py;
00186   double pz = (*pj)->pz;
00187   double theta = acos(pz*s/sqrt(px*px+py*py+pz*pz));
00188   return -log(tan(theta/2.0));
00189 }
00190 
00191 double JetSplitMergeTool::phiTrue(Jet::constit_vect_t::iterator pj)
00192 {
00193   double s = ((*pj)->e() > 0) ? +1.0 : -1.0;
00194   double px = (*pj)->px;
00195   double py = (*pj)->py;
00196   return atan2(py*s,px*s);
00197 }
00198 
00199 }  // namespace atlas
00200 
00201 FASTJET_END_NAMESPACE
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines