FastJet 3.0.2
|
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