FastJet 3.0.2
ILConeAlgorithm.hpp
00001 #ifndef  D0RunIIconeJets_ILCONEALGORITHM
00002 #define  D0RunIIconeJets_ILCONEALGORITHM
00003 // ---------------------------------------------------------------------------
00004 // ILConeAlgorithm.hpp
00005 //
00006 // Created: 28-JUL-2000 Francois Touze (+ Laurent Duflot)
00007 //
00008 // Purpose: Implements the Improved Legacy Cone Algorithm
00009 //
00010 // Modified:
00011 //   31-JUL-2000  Laurent Duflot
00012 //     + introduce support for additional informations (ConeJetInfo)
00013 //    1-AUG-2000  Laurent Duflot
00014 //     + seedET for midpoint jet is -999999., consistent with seedET ordering
00015 //       in ConeSplitMerge: final jets with seedET=-999999. will be midpoint 
00016 //       jets which were actually different from stable cones.
00017 //    4-Aug-2000  Laurent Duflot
00018 //     + remove unecessary copy in TemporaryJet::is_stable()
00019 //   11-Aug-2000  Laurent Duflot
00020 //     + remove using namespace std
00021 //     + add threshold on item. The input list to makeClusters() IS modified
00022 //   20-June-2002 John Krane
00023 //     + remove bug in midpoint calculation based on pT weight
00024 //     + started to add new midpoint calculation based on 4-vectors,
00025 //       but not enough info held by this class
00026 //   24-June-2002 John Krane
00027 //     + modify is_stable() to not iterate if desired 
00028 //       (to expand search cone w/out moving it)
00029 //     + added search cone logic for initial jets but not midpoint jets as per
00030 //       agreement with CDF
00031 //   19-July-2002 John Krane
00032 //     + _SEARCH_CONE size factor now provided by calreco/CalClusterReco.cpp 
00033 //     + (default = 1.0 = like original ILCone behavior)
00034 //   10-Oct-2002 John Krane
00035 //     + Check min Pt of cone with full size first, then iterate with search cone
00036 //   07-Dec-2002 John Krane
00037 //     + speed up the midpoint phi-wrap solution
00038 //   01-May-2007 Lars Sonnenschein
00039 //   extracted from D0 software framework and modified to remove subsequent dependencies
00040 //
00041 //
00042 // This file is distributed with FastJet under the terms of the GNU
00043 // General Public License (v2). Permission to do so has been granted
00044 // by Lars Sonnenschein and the D0 collaboration (see COPYING for
00045 // details)
00046 //
00047 // History of changes in FastJet compared tothe original version of
00048 // ProtoJet.hpp
00049 //
00050 // 2011-12-13  Gregory Soyez  <soyez@fastjet.fr>
00051 // 
00052 //        * added license information
00053 //
00054 // 2011-11-14  Gregory Soyez  <soyez@fastjet.fr>
00055 //
00056 //        * changed the name of a few parameters to avoid a gcc
00057 //          -Wshadow warning
00058 //
00059 // 2009-01-17  Gregory Soyez  <soyez@fastjet.fr>
00060 //
00061 //        * put the code in the fastjet::d0 namespace
00062 //
00063 // 2007-12-14  Gavin Salam  <salam@lpthe.jussieu.fr>
00064 // 
00065 //        * moved the 'std::vector<ProtoJet<Item> > ilcv' structure
00066 //          containing the info about the final jets from a local
00067 //          variable to a class variable (for integration in the
00068 //          FastJet plugin core)
00069 //
00070 // ---------------------------------------------------------------------------
00071 
00072 #include <vector>
00073 #include <list>
00074 #include <utility>  // defines pair<,>
00075 #include <map>
00076 #include <algorithm>
00077 #include <iostream>
00078 
00079 
00080 //#include "energycluster/EnergyClusterReco.hpp"
00081 //#include "energycluster/ProtoJet.hpp"
00082 #include "ProtoJet.hpp"
00083 //#include "energycluster/ConeSplitMerge.hpp"
00084 #include "ConeSplitMerge.hpp"
00085 //#include "energycluster/ConeJetInfoChunk.hpp"
00086 
00087 #include "inline_maths.h"
00088 
00089 ///////////////////////////////////////////////////////////////////////////////
00090 #include <fastjet/internal/base.hh>
00091 
00092 FASTJET_BEGIN_NAMESPACE
00093 
00094 namespace d0{
00095 
00096 using namespace inline_maths;
00097 
00098 /*
00099  NB: Some attempt at optimizing the code has been made by ordering the object
00100      along rapidity but this did not improve the speed. There are traces of 
00101      these atemps in the code, that will be cleaned up in the future.
00102  */
00103 
00104 // at most one of those !
00105 // order the input list and use lower_bound() and upper_bound() to restrict the
00106 // on item to those that could possibly be in the cone.
00107 //#define ILCA_USE_ORDERED_LIST
00108 
00109 // idem but use an intermediate multimap in hope that lower_bound() and 
00110 // upper_bound() are faster in this case.
00111 //#define ILCA_USE_MMAP
00112 
00113 
00114 #ifdef ILCA_USE_ORDERED_LIST
00115 // this class is used to order the item list along rapidity
00116 template <class Item>
00117 class rapidity_order
00118 {
00119 public:
00120   bool operator()(const Item* first, const Item* second)
00121   {
00122     return (first->y() < second->y());
00123   }
00124   bool operator()(float const & first, const Item* second)
00125   {
00126     return (first  < second->y());
00127   }
00128   bool operator()(const Item* first, float const& second)
00129   {
00130     return (first->y() < second);
00131   }
00132 };
00133 #endif
00134 
00135 
00136 //template <class Item,class ItemAddress,class IChunk>
00137 template <class Item>
00138 class ILConeAlgorithm 
00139 {
00140 
00141 public:
00142 
00143   // default constructor (default parameters are crazy: you should not use that
00144   // constructor !)
00145   ILConeAlgorithm():
00146     _CONE_RADIUS(0.),
00147     _MIN_JET_ET(99999.),
00148     _ET_MIN_RATIO(0.5),
00149     _FAR_RATIO(0.5),
00150     _SPLIT_RATIO(0.5),
00151     _DUPLICATE_DR(0.005),
00152     _DUPLICATE_DPT(0.01),
00153     _SEARCH_CONE(0.5),
00154     _PT_MIN_LEADING_PROTOJET(0), 
00155     _PT_MIN_SECOND_PROTOJET(0), 
00156     _MERGE_MAX(10000), 
00157     _PT_MIN_noMERGE_MAX(0)
00158   {;}
00159 
00160   // full constructor
00161   ILConeAlgorithm(float cone_radius, float min_jet_Et, float split_ratio,
00162                   float far_ratio=0.5, float Et_min_ratio=0.5,
00163                   bool kill_duplicate=true, float duplicate_dR=0.005, 
00164                   float duplicate_dPT=0.01, float search_factor=1.0, 
00165                   float pT_min_leading_protojet=0., float pT_min_second_protojet=0.,
00166                   int merge_max=10000, float pT_min_nomerge=0.) :
00167     // cone radius
00168     _CONE_RADIUS(cone_radius), 
00169     // minimum jet ET
00170     _MIN_JET_ET(min_jet_Et), 
00171     // stable cones must have ET > _ET_MIN_RATIO*_MIN_JET_ET at any iteration
00172     _ET_MIN_RATIO(Et_min_ratio),
00173     // precluster at least _FAR_RATIO*_CONE_RADIUS away from stable cones
00174     _FAR_RATIO(far_ratio), 
00175     // split or merge criterium           
00176     _SPLIT_RATIO(split_ratio),
00177     _DUPLICATE_DR(duplicate_dR),
00178     _DUPLICATE_DPT(duplicate_dPT),
00179     _SEARCH_CONE(cone_radius/search_factor),
00180     // kill stable cone if within _DUPLICATE_DR and delta(pT)<_DUPLICATE_DPT
00181     // of another stable cone.
00182     _KILL_DUPLICATE(kill_duplicate),
00183     _PT_MIN_LEADING_PROTOJET(pT_min_leading_protojet),
00184     _PT_MIN_SECOND_PROTOJET(pT_min_second_protojet),
00185     _MERGE_MAX(merge_max),
00186     _PT_MIN_noMERGE_MAX(pT_min_nomerge)
00187     {;}
00188 
00189   //destructor
00190   ~ILConeAlgorithm() {;}
00191 
00192   // make jet clusters using improved legacy cone algorithm
00193   //void makeClusters(const EnergyClusterReco* r,
00194   void makeClusters(
00195                     // the input list modified (ordered)
00196                     std::list<Item> &jets,
00197                     std::list<const Item*>& itemlist, 
00198                     //float zvertex,   
00199                     ////std::list<const Item*>& itemlist);   
00200                     //const EnergyClusterCollection<ItemAddress>& preclu,
00201                     //IChunk* chunkptr, ConeJetInfoChunk* infochunkptr,
00202                     const float Item_ET_Threshold);
00203 
00204   // this will hold the final jets + contents
00205   std::vector<ProtoJet<Item> > ilcv;
00206 
00207 private:
00208 
00209   float _CONE_RADIUS;
00210   float _MIN_JET_ET;
00211   float _ET_MIN_RATIO;
00212   float _FAR_RATIO;
00213   float _SPLIT_RATIO;
00214   float _DUPLICATE_DR;
00215   float _DUPLICATE_DPT;
00216   float _SEARCH_CONE;
00217 
00218   bool _KILL_DUPLICATE;
00219 
00220   float _PT_MIN_LEADING_PROTOJET;
00221   float _PT_MIN_SECOND_PROTOJET;
00222   int  _MERGE_MAX; 
00223   float _PT_MIN_noMERGE_MAX;
00224 
00225   // private class 
00226   // This is a ProtoJet with additional methods dist(), midpoint() and 
00227   // is_stable()
00228   class TemporaryJet : public ProtoJet<Item> 
00229   {
00230     
00231   public :
00232     
00233     TemporaryJet(float seedET) : ProtoJet<Item>(seedET) {;}
00234 
00235     TemporaryJet(float seedET,float y_in,float phi_in) : 
00236       ProtoJet<Item>(seedET,y_in,phi_in) {;}
00237     
00238     ~TemporaryJet() {;}
00239     
00240     float dist(TemporaryJet& jet) const 
00241     {
00242       return RDelta(this->_y,this->_phi,jet.y(),jet.phi()); 
00243     }
00244     
00245     void midpoint(const TemporaryJet& jet,float & y_out, float & phi_out) const 
00246     {
00247       // Midpoint should probably be computed w/4-vectors but don't 
00248       // have that info.  Preserving Pt-weighted calculation - JPK
00249       float pTsum = this->_pT + jet.pT();
00250       y_out = (this->_y*this->_pT + jet.y()*jet.pT())/pTsum;
00251 
00252       phi_out = (this->_phi*this->_pT + jet.phi()*jet.pT())/pTsum;
00253       // careful with phi-wrap area: convert from [0,2pi] to [-pi,pi]
00254       //ls: original D0 code, as of 23/Mar/2007
00255       //if ( abs(phi-this->_phi)>2.0 ) { // assumes cones R=1.14 or smaller, merge within 2R only  
00256       //ls: abs bug fixed 26/Mar/2007 
00257       if ( fabs(phi_out-this->_phi)>2.0 ) { // assumes cones R=1.14 or smaller, merge within 2R only  
00258         phi_out = fmod( this->_phi+PI, TWOPI);
00259         if (phi_out < 0.0) phi_out += TWOPI;
00260         phi_out -= PI;
00261 
00262         float temp=fmod( jet.phi()+PI, TWOPI);
00263         if (temp < 0.0) temp += TWOPI;
00264         temp -= PI;
00265 
00266         phi_out = (phi_out*this->_pT + temp*jet.pT()) /pTsum;
00267       }
00268 
00269       if ( phi_out < 0. ) phi_out += TWOPI;
00270     }
00271     
00272 
00273 ////////////////////////////////////////
00274 #ifdef ILCA_USE_MMAP
00275     bool is_stable(const std::multimap<float,const Item*>& items, 
00276                    float radius, float min_ET, int max_iterations=50) 
00277 #else
00278     bool is_stable(const std::list<const Item*>& itemlist, float radius, 
00279                  float min_ET, int max_iterations=50) 
00280 #endif
00281     // Note: max_iterations = 0 will just recompute the jet using the specified cone
00282     {
00283       float radius2 = radius*radius;
00284       float Rcut= 1.E-06;
00285       
00286       
00287       // ?? if(_Increase_Delta_R) Rcut= 1.E-04;
00288       bool stable= true;
00289       int trial= 0;
00290       float Yst;
00291       float PHIst;
00292       do {  
00293         trial++;
00294         //std::cout << "   trial " << trial << " " << _y << " " << _phi << std::endl; 
00295         Yst  = this->_y;
00296         PHIst= this->_phi;    
00297         //cout << "is_stable beginning do loop: this->_pT=" << this->_pT << " this->_y=" << this->_y << " this->_phi=" << this->_phi << endl;
00298         this->erase();
00299         
00300         this->setJet(Yst,PHIst,0.0);
00301         
00302         
00303 #ifdef ILCA_USE_ORDERED_LIST      
00304         std::list<const Item*>::const_iterator lower = 
00305           lower_bound(itemlist.begin(),itemlist.end(),Yst-radius,
00306                       rapidity_order<Item>());
00307         std::list<const Item*>::const_iterator upper = 
00308           upper_bound(itemlist.begin(),itemlist.end(),Yst+radius,
00309                       rapidity_order<Item>());
00310         for(std::list<const Item*>::const_iterator tk = lower; tk != upper; ++tk)      {
00311           //std::cout << " is_stable: item y=" << (*tk)->y() << " phi=" << (*tk)->phi() << " RD2=" << RD2((*tk)->y(),(*tk)->phi(),Yst,PHIst) << " " << Yst-radius << " " << Yst+radius << endl;
00312           if(RD2((*tk)->y(),(*tk)->phi(),Yst,PHIst) <= radius2) 
00313             {
00314               addItem(*tk);
00315             }
00316         }       
00317 #else
00318 #ifdef ILCA_USE_MMAP      
00319         // need to loop only on the subset with   Yst-R < y < Yst+R
00320         for ( std::multimap<float,const Item*>::const_iterator 
00321                 tk = items.lower_bound(Yst-radius);
00322               tk != items.upper_bound(Yst+radius); ++tk )
00323           {
00324             //std::cout << "     item " << (*tk)->y() << " " << (*tk)->phi() << " " << RD2((*tk)->y(),(*tk)->phi(),Yst,PHIst) << " " << Yst-radius << " " << Yst+radius << endl;
00325             if(RD2(((*tk).second)->y(),((*tk).second)->phi(),Yst,PHIst) <= radius2) 
00326               {
00327                 addItem((*tk).second);
00328               }
00329           }
00330         
00331 #else   
00332 
00333         //cout << " is_stable: itemlist.size()=" << itemlist.size() << endl;
00334         for(typename std::list<const Item*>::const_iterator tk = itemlist.begin(); tk != itemlist.end(); ++tk) 
00335           {
00336             //std::cout << "    is_stable: item (*tk)->y()=" << (*tk)->y() << " (*tk)->phi=" << (*tk)->phi() << " RD2=" << RD2((*tk)->y(),(*tk)->phi(),Yst,PHIst) << " Yst-rad=" << Yst-radius << " Yst+rad=" << Yst+radius << endl;
00337             if(RD2((*tk)->y(),(*tk)->phi(),Yst,PHIst) <= radius2) 
00338                {
00339                  //cout << "add item to *tk" << endl;
00340                 addItem(*tk);
00341               }
00342           }
00343 #endif
00344 #endif      
00345       
00346         //std::cout << "is_stable, before update: jet this->_y=" << this->_y << " _phi=" << this->_phi << " _pT=" << this->_pT << " min_ET=" << min_ET << std::endl; 
00347         this->updateJet();
00348         //std::cout << "is_stable, after update: jet this->_y=" << this->_y << " _phi=" << this->_phi << " _pT=" << this->_pT << " min_ET=" << min_ET << std::endl; 
00349         
00350         if(this->_pT < min_ET ) 
00351           {
00352             stable= false;
00353             break;
00354           } 
00355         //cout << "is_stable end while loop: this->_pT=" << this->_pT << " this->_y=" << this->_y << " this->_phi=" << this->_phi << endl;
00356       } while(RD2(this->_y,this->_phi,Yst,PHIst) >= Rcut && trial <= max_iterations);
00357       //std::cout << "   trial stable " << trial << " " << stable << std::endl; 
00358       return stable;
00359     }
00360     
00361   private :
00362     
00363   };
00364 };
00365 ///////////////////////////////////////////////////////////////////////////////
00366 //template <class Item,class ItemAddress,class IChunk>
00367 //void ILConeAlgorithm <Item,ItemAddress,IChunk>::
00368 template <class Item>
00369 void ILConeAlgorithm <Item>::
00370 //makeClusters(const EnergyClusterReco* r,
00371 makeClusters(
00372              std::list<Item> &jets,
00373              std::list<const Item*>& ilist, 
00374              //float Zvertex, 
00375              ////std::list<const Item*>& ilist) 
00376              //const EnergyClusterCollection<ItemAddress>& preclu,
00377              //IChunk* chunkptr, ConeJetInfoChunk* infochunkptr,
00378              const float Item_ET_Threshold) 
00379 {
00380   // remove items below threshold
00381   for ( typename std::list<const Item*>::iterator it = ilist.begin(); 
00382 
00383         it != ilist.end(); )
00384   {
00385     if ( (*it)->pT() < Item_ET_Threshold ) 
00386     {
00387       it = ilist.erase(it);
00388     }
00389       else ++it;
00390   }
00391 
00392   // create an energy cluster collection for jets 
00393   //EnergyClusterCollection<ItemAddress>* ptrcol;
00394   //Item* ptrcol;
00395   //r->createClusterCollection(chunkptr,ptrcol);
00396   ////std::vector<const EnergyCluster<ItemAddress>*> ecv;
00397   std::vector<const Item*> ecv;
00398   for ( typename std::list<const Item*>::iterator it = ilist.begin(); 
00399         it != ilist.end(); it++) {
00400     ecv.push_back(*it);
00401   }
00402 
00403 
00404   //preclu.getClusters(ecv);
00405   //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% need to fill here vector ecv
00406 
00407   //std::cout << " number of seed clusters: " << ecv.size() << std::endl;
00408 
00409   // skip precluster close to jets
00410   float far_def = _FAR_RATIO*_CONE_RADIUS * _FAR_RATIO*_CONE_RADIUS;
00411 
00412   // skip if jet Et is below some value
00413   float ratio= _MIN_JET_ET*_ET_MIN_RATIO;
00414 
00415 
00416 #ifdef ILCA_USE_ORDERED_LIST
00417   // sort the list in rapidity order
00418   ilist.sort(rapidity_order<Item>());
00419 #else
00420 #ifdef ILCA_USE_MMAP  
00421   // create a y ordered list of items 
00422   std::multimap<float,const Item*> items;
00423   std::list<const Item*>::const_iterator it;
00424   for(it = ilist.begin(); it != ilist.end(); ++it) 
00425   {
00426     pair<float,const Item*> p((*it)->y(),*it);
00427     items.insert(p);
00428   }
00429 #endif
00430 #endif
00431 
00432   std::vector<ProtoJet<Item> > mcoll;
00433   std::vector<TemporaryJet> scoll; 
00434 
00435 
00436   // find stable jets around seeds 
00437   //typename std::vector<const EnergyCluster<ItemAddress>* >::iterator jclu;
00438   typename std::vector<const Item*>::iterator jclu;
00439   for(jclu = ecv.begin(); jclu != ecv.end(); ++jclu) 
00440   {
00441     //const EnergyCluster<ItemAddress>* ptr= *jclu;
00442     const Item* ptr= *jclu;
00443     float p[4];
00444     ptr->p4vec(p);
00445     float Yst  = P2y(p);
00446     float PHIst= P2phi(p);
00447 
00448     // don't keep preclusters close to jet
00449     bool is_far= true;
00450     // ?? if(_Kill_Far_Clusters) {
00451     for(unsigned int i = 0; i < scoll.size(); ++i) 
00452     {
00453       float y  = scoll[i].y();
00454       float phi= scoll[i].phi();
00455       if(RD2(Yst,PHIst,y,phi) < far_def) 
00456       {
00457         is_far= false;
00458         break;
00459       }
00460     }
00461     // ?? }
00462 
00463     if(is_far) 
00464     {
00465       TemporaryJet jet(ptr->pT(),Yst,PHIst);
00466       //cout << "temporary jet: pt=" << ptr->pT() << " y=" << Yst << " phi=" << PHIst << endl;
00467 
00468       // Search cones are smaller, so contain less jet Et 
00469       // Don't throw out too many little jets during search phase!
00470       // Strategy: check Et first with full cone, then search with low-GeV min_et thresh
00471 #ifdef ILCA_USE_MMAP
00472       if(jet.is_stable(items,_CONE_RADIUS,ratio,0) && jet.is_stable(items,_SEARCH_CONE,3.0)) 
00473 #else
00474       if(jet.is_stable(ilist,_CONE_RADIUS,ratio,0) && jet.is_stable(ilist,_SEARCH_CONE,3.0)) 
00475 #endif
00476       {
00477 
00478         //cout << "temporary jet is stable" << endl;
00479 
00480 // jpk  Resize the found jets 
00481 #ifdef ILCA_USE_MMAP
00482         jet.is_stable(items,_CONE_RADIUS,ratio,0) ;
00483 #else
00484         jet.is_stable(ilist,_CONE_RADIUS,ratio,0) ;
00485 #endif
00486         //cout << "found jet has been resized if any" << endl;
00487 
00488         if ( _KILL_DUPLICATE ) 
00489         {
00490           // check if we are not finding the same jet again
00491           float distmax = 999.; int imax = -1;
00492           for(unsigned int i = 0; i < scoll.size(); ++i) 
00493           {
00494             float dist = jet.dist(scoll[i]);
00495             if ( dist < distmax )
00496             {
00497               distmax = dist;
00498               imax = i;
00499             }
00500           }
00501           if ( distmax > _DUPLICATE_DR ||
00502                fabs((jet.pT()-scoll[imax].pT())/scoll[imax].pT())>_DUPLICATE_DPT )
00503           {
00504             scoll.push_back(jet);
00505             mcoll.push_back(jet);
00506             //std::cout << " stable cone " << jet.y() << " " << jet.phi() << " " << jet.pT() << std::endl;
00507           }
00508           /*
00509             else
00510             {
00511             std::cout << " jet too close to a found jet " << distmax << " " << 
00512             fabs((jet.pT()-scoll[imax].pT())/scoll[imax].pT()) << std::endl;
00513             }
00514           */
00515         }
00516         else
00517         {
00518           scoll.push_back(jet);
00519           mcoll.push_back(jet);
00520           //std::cout << " stable cone " << jet.y() << " " << jet.phi() << " " << jet.pT() << std::endl;
00521         }
00522 
00523       }
00524     }
00525   }
00526 
00527   // find stable jets around midpoints
00528   for(unsigned int i = 0; i < scoll.size(); ++i) 
00529   {
00530     for(unsigned int k = i+1; k < scoll.size(); ++k) 
00531     {
00532       float djet= scoll[i].dist(scoll[k]);
00533       if(djet > _CONE_RADIUS && djet < 2.*_CONE_RADIUS) 
00534       {
00535         float y_mid,phi_mid;
00536         scoll[i].midpoint(scoll[k],y_mid,phi_mid);
00537         TemporaryJet jet(-999999.,y_mid,phi_mid);
00538         //std::cout << " midpoint: " << scoll[i].pT() << " " << scoll[i].info().seedET() << " " << scoll[k].pT() << " " << scoll[k].info().seedET() << " " << y_mid << " " << phi_mid << std::endl; 
00539 
00540 // midpoint jets are full size
00541 #ifdef ILCA_USE_MMAP
00542       if(jet.is_stable(items,_CONE_RADIUS,ratio)) 
00543 #else
00544       if(jet.is_stable(ilist,_CONE_RADIUS,ratio)) 
00545 #endif
00546         {
00547           mcoll.push_back(jet);
00548           //std::cout << " stable midpoint cone " << jet.y() << " " << jet.phi() << " " << jet.pT() << std::endl;
00549         }
00550       }
00551     }
00552   }
00553 
00554 
00555   // do a pT ordered splitting/merging
00556   ConeSplitMerge<Item> pjets(mcoll);
00557   ilcv.clear();
00558   pjets.split_merge(ilcv,_SPLIT_RATIO, _PT_MIN_LEADING_PROTOJET, _PT_MIN_SECOND_PROTOJET, _MERGE_MAX, _PT_MIN_noMERGE_MAX);
00559 
00560 
00561   for(unsigned int i = 0; i < ilcv.size(); ++i) 
00562   {
00563     if ( ilcv[i].pT() > _MIN_JET_ET )
00564     {
00565       //EnergyCluster<ItemAddress>* ptrclu;
00566       Item ptrclu;
00567       //r->createCluster(ptrcol,ptrclu);
00568       
00569       std::list<const Item*> tlist=ilcv[i].LItems();
00570       typename std::list<const Item*>::iterator tk;
00571       for(tk = tlist.begin(); tk != tlist.end(); ++tk) 
00572       {
00573         //ItemAddress addr= (*tk)->address();
00574         float pk[4];
00575         (*tk)->p4vec(pk);
00576         //std::cout << (*tk)->index <<" " <<  (*tk) << std::endl;
00577         //float emE= (*tk)->emE();
00578         //r->addClusterItem(ptrclu,addr,pk,emE);
00579         //ptrclu->Add(*tk);
00580         ptrclu.Add(**tk);
00581       }
00582       // print out
00583       //ptrclu->print(cout);
00584       
00585       //infochunkptr->addInfo(ilcv[i].info());
00586       jets.push_back(ptrclu);
00587     }
00588   }
00589 }
00590 
00591 }  // namespace d0
00592 
00593 FASTJET_END_NAMESPACE
00594 
00595 #endif
00596 
00597 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends