|
fastjet 2.4.5
|
|
fastjet 2.4.5
|
Namespaces | |
| namespace | D0RunIIconeJets_CONEJETINFO |
| namespace | inline_maths |
Classes | |
| class | ProtoJet_ET_seedET_order |
| class | ConeSplitMerge |
| class | HepEntity |
| class | ILConeAlgorithm |
| class | ProtoJet |
Functions | |
| template<class Item > | |
| void | ILConeAlgorithm< Item >makeClusters (std::list< Item > &jets, std::list< const Item * > &ilist, const float Item_ET_Threshold) |
| int | main () |
| float | RD2 (float y1, float phi1, float y2, float phi2) |
| float | RDelta (float y1, float phi1, float y2, float phi2) |
| float | P2y (float *p4vec) |
| float | P2phi (float *p4vec) |
| void fastjet::d0::ILConeAlgorithm< Item >makeClusters | ( | std::list< Item > & | jets, |
| std::list< const Item * > & | ilist, | ||
| const float | Item_ET_Threshold | ||
| ) |
| Item_ET_Threshold | /std::list<const Item*>& ilist) |
Definition at line 343 of file ILConeAlgorithm.hpp.
References P2phi(), P2y(), fastjet::d0::inline_maths::phi(), RD2(), fastjet::d0::ConeSplitMerge< Item >::split_merge(), and fastjet::d0::inline_maths::y().
{
// remove items below threshold
for ( typename std::list<const Item*>::iterator it = ilist.begin();
it != ilist.end(); )
{
if ( (*it)->pT() < Item_ET_Threshold )
{
it = ilist.erase(it);
}
else ++it;
}
// create an energy cluster collection for jets
//EnergyClusterCollection<ItemAddress>* ptrcol;
//Item* ptrcol;
//r->createClusterCollection(chunkptr,ptrcol);
std::vector<const Item*> ecv;
for ( typename std::list<const Item*>::iterator it = ilist.begin();
it != ilist.end(); it++) {
ecv.push_back(*it);
}
//preclu.getClusters(ecv);
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% need to fill here vector ecv
//std::cout << " number of seed clusters: " << ecv.size() << std::endl;
// skip precluster close to jets
float far_def = _FAR_RATIO*_CONE_RADIUS * _FAR_RATIO*_CONE_RADIUS;
// skip if jet Et is below some value
float ratio= _MIN_JET_ET*_ET_MIN_RATIO;
#ifdef ILCA_USE_ORDERED_LIST
// sort the list in rapidity order
ilist.sort(rapidity_order<Item>());
#else
#ifdef ILCA_USE_MMAP
// create a y ordered list of items
std::multimap<float,const Item*> items;
std::list<const Item*>::const_iterator it;
for(it = ilist.begin(); it != ilist.end(); ++it)
{
pair<float,const Item*> p((*it)->y(),*it);
items.insert(p);
}
#endif
#endif
std::vector<ProtoJet<Item> > mcoll;
std::vector<TemporaryJet> scoll;
// find stable jets around seeds
//typename std::vector<const EnergyCluster<ItemAddress>* >::iterator jclu;
typename std::vector<const Item*>::iterator jclu;
for(jclu = ecv.begin(); jclu != ecv.end(); ++jclu)
{
//const EnergyCluster<ItemAddress>* ptr= *jclu;
const Item* ptr= *jclu;
float p[4];
ptr->p4vec(p);
float Yst = P2y(p);
float PHIst= P2phi(p);
// don't keep preclusters close to jet
bool is_far= true;
// ?? if(_Kill_Far_Clusters) {
for(unsigned int i = 0; i < scoll.size(); ++i)
{
float y = scoll[i].y();
float phi= scoll[i].phi();
if(RD2(Yst,PHIst,y,phi) < far_def)
{
is_far= false;
break;
}
}
// ?? }
if(is_far)
{
TemporaryJet jet(ptr->pT(),Yst,PHIst);
//cout << "temporary jet: pt=" << ptr->pT() << " y=" << Yst << " phi=" << PHIst << endl;
// Search cones are smaller, so contain less jet Et
// Don't throw out too many little jets during search phase!
// Strategy: check Et first with full cone, then search with low-GeV min_et thresh
#ifdef ILCA_USE_MMAP
if(jet.is_stable(items,_CONE_RADIUS,ratio,0) && jet.is_stable(items,_SEARCH_CONE,3.0))
#else
if(jet.is_stable(ilist,_CONE_RADIUS,ratio,0) && jet.is_stable(ilist,_SEARCH_CONE,3.0))
#endif
{
//cout << "temporary jet is stable" << endl;
// jpk Resize the found jets
#ifdef ILCA_USE_MMAP
jet.is_stable(items,_CONE_RADIUS,ratio,0) ;
#else
jet.is_stable(ilist,_CONE_RADIUS,ratio,0) ;
#endif
//cout << "found jet has been resized if any" << endl;
if ( _KILL_DUPLICATE )
{
// check if we are not finding the same jet again
float distmax = 999.; int imax = -1;
for(unsigned int i = 0; i < scoll.size(); ++i)
{
float dist = jet.dist(scoll[i]);
if ( dist < distmax )
{
distmax = dist;
imax = i;
}
}
if ( distmax > _DUPLICATE_DR ||
fabs((jet.pT()-scoll[imax].pT())/scoll[imax].pT())>_DUPLICATE_DPT )
{
scoll.push_back(jet);
mcoll.push_back(jet);
//std::cout << " stable cone " << jet.y() << " " << jet.phi() << " " << jet.pT() << std::endl;
}
/*
else
{
std::cout << " jet too close to a found jet " << distmax << " " <<
fabs((jet.pT()-scoll[imax].pT())/scoll[imax].pT()) << std::endl;
}
*/
}
else
{
scoll.push_back(jet);
mcoll.push_back(jet);
//std::cout << " stable cone " << jet.y() << " " << jet.phi() << " " << jet.pT() << std::endl;
}
}
}
}
// find stable jets around midpoints
for(unsigned int i = 0; i < scoll.size(); ++i)
{
for(unsigned int k = i+1; k < scoll.size(); ++k)
{
float djet= scoll[i].dist(scoll[k]);
if(djet > _CONE_RADIUS && djet < 2.*_CONE_RADIUS)
{
float y_mid,phi_mid;
scoll[i].midpoint(scoll[k],y_mid,phi_mid);
TemporaryJet jet(-999999.,y_mid,phi_mid);
//std::cout << " midpoint: " << scoll[i].pT() << " " << scoll[i].info().seedET() << " " << scoll[k].pT() << " " << scoll[k].info().seedET() << " " << y_mid << " " << phi_mid << std::endl;
// midpoint jets are full size
#ifdef ILCA_USE_MMAP
if(jet.is_stable(items,_CONE_RADIUS,ratio))
#else
if(jet.is_stable(ilist,_CONE_RADIUS,ratio))
#endif
{
mcoll.push_back(jet);
//std::cout << " stable midpoint cone " << jet.y() << " " << jet.phi() << " " << jet.pT() << std::endl;
}
}
}
}
// do a pT ordered splitting/merging
ConeSplitMerge<Item> pjets(mcoll);
ilcv.clear();
pjets.split_merge(ilcv,_SPLIT_RATIO, _PT_MIN_LEADING_PROTOJET, _PT_MIN_SECOND_PROTOJET, _MERGE_MAX, _PT_MIN_noMERGE_MAX);
for(unsigned int i = 0; i < ilcv.size(); ++i)
{
if ( ilcv[i].pT() > _MIN_JET_ET )
{
//EnergyCluster<ItemAddress>* ptrclu;
Item ptrclu;
//r->createCluster(ptrcol,ptrclu);
std::list<const Item*> tlist=ilcv[i].LItems();
typename std::list<const Item*>::iterator tk;
for(tk = tlist.begin(); tk != tlist.end(); ++tk)
{
//ItemAddress addr= (*tk)->address();
float pk[4];
(*tk)->p4vec(pk);
//std::cout << (*tk)->index <<" " << (*tk) << std::endl;
//float emE= (*tk)->emE();
//r->addClusterItem(ptrclu,addr,pk,emE);
//ptrclu->Add(*tk);
ptrclu.Add(**tk);
}
// print out
//ptrclu->print(cout);
//infochunkptr->addInfo(ilcv[i].info());
jets.push_back(ptrclu);
}
}
}
| int fastjet::d0::main | ( | ) |
Definition at line 13 of file main.C.
References fastjet::d0::HepEntity::Fill(), and fastjet::d0::ILConeAlgorithm< Item >::makeClusters().
{
HepEntity el;
list<const HepEntity*> *ensemble = new list<const HepEntity*>;
//list<const HepEntity*> ensemble;
//fill with E, px, py, pz
el.Fill(100., 25., 25., 25., 0);
ensemble->push_back(new HepEntity(el));
el.Fill(105., 20., 30., 30., 1);
ensemble->push_back(new HepEntity(el));
el.Fill(60., 20., 20., 20., 2);
ensemble->push_back(new HepEntity(el));
el.Fill(95., 65., 10., 20., 3);
ensemble->push_back(new HepEntity(el));
el.Fill(110., 25., -25., -25., 4);
ensemble->push_back(new HepEntity(el));
el.Fill(100., 23., -25., -25., 5);
ensemble->push_back(new HepEntity(el));
el.Fill(101., 25., -20., -25., 6);
ensemble->push_back(new HepEntity(el));
el.Fill(102., 25., -25., -23., 7);
ensemble->push_back(new HepEntity(el));
cout << "list->size()=" << ensemble->size() << endl;
int i=1;
for (list<const HepEntity*>::iterator it = ensemble->begin(); it != ensemble->end(); ++it) {
cout << "4-vector " << i++ << " : E=" << (*it)->E << " pT=" << (*it)->pT() << " y=" << (*it)->y() << " phi=" << (*it)->phi() << endl;
cout << (*it) << endl;
}
float cone_radius = 0.5;
float min_jet_Et = 8.0;
float split_ratio = 0.5;
//the parameters below have been found to be set to the values given below
//in the original implementation, shouldn't be altered
float far_ratio=0.5;
float Et_min_ratio=0.5;
bool kill_duplicate=true;
float duplicate_dR=0.005;
float duplicate_dPT=0.01;
float search_factor=1.0;
float pT_min_leading_protojet=0.;
float pT_min_second_protojet=0.;
int merge_max=10000;
float pT_min_nomerge=0.;
ILConeAlgorithm<HepEntity>
ilegac(cone_radius, min_jet_Et, split_ratio,
far_ratio, Et_min_ratio, kill_duplicate, duplicate_dR,
duplicate_dPT, search_factor, pT_min_leading_protojet,
pT_min_second_protojet, merge_max, pT_min_nomerge);
float Item_ET_Threshold = 0.;
float Zvertex = 0.;
float* Item_ET_Threshold_ptr = &Item_ET_Threshold;
list<HepEntity> jets;
ilegac.makeClusters(jets, *ensemble, Item_ET_Threshold);
list<HepEntity>::iterator it;
cout << "Number of jets = " << jets.size() << endl;
for (it=jets.begin(); it!=jets.end(); ++it) {
cout << "jet: E=" << (*it).E << " pT=" << (*it).pT() << " y=" << (*it).y() << " phi=" << (*it).phi() << endl;
}
//delete elements of the ensemble particle list
//relevant to prevent memory leakage when running over many events
for (list<const HepEntity*>::iterator it = ensemble->begin(); it != ensemble->end(); ++it) {
delete *it;
}
delete ensemble;
return 0;
}
| float fastjet::d0::P2phi | ( | float * | p4vec | ) | [inline] |
Definition at line 52 of file ProtoJet.hpp.
References fastjet::d0::inline_maths::phi().
Referenced by ILConeAlgorithm< Item >makeClusters(), and fastjet::d0::ProtoJet< Item >::updateJet().
{
return phi(p4vec[0],p4vec[1]);
}
| float fastjet::d0::P2y | ( | float * | p4vec | ) | [inline] |
Definition at line 48 of file ProtoJet.hpp.
References fastjet::d0::inline_maths::y().
Referenced by ILConeAlgorithm< Item >makeClusters(), and fastjet::d0::ProtoJet< Item >::updateJet().
{
return y(p4vec[3],p4vec[2]);
}
| float fastjet::d0::RD2 | ( | float | y1, |
| float | phi1, | ||
| float | y2, | ||
| float | phi2 | ||
| ) | [inline] |
Definition at line 36 of file ProtoJet.hpp.
References fastjet::d0::inline_maths::delta_phi().
Referenced by ILConeAlgorithm< Item >makeClusters(), fastjet::d0::ILConeAlgorithm< Item >::TemporaryJet::is_stable(), and fastjet::d0::ConeSplitMerge< Item >::split_merge().
{
float dphi= delta_phi(phi1,phi2);
return (y1-y2)*(y1-y2)+dphi*dphi;
}
| float fastjet::d0::RDelta | ( | float | y1, |
| float | phi1, | ||
| float | y2, | ||
| float | phi2 | ||
| ) | [inline] |
Definition at line 42 of file ProtoJet.hpp.
References fastjet::d0::inline_maths::delta_phi().
Referenced by fastjet::d0::ILConeAlgorithm< Item >::TemporaryJet::dist().
{
float dphi= delta_phi(phi1,phi2);
return sqrt((y1-y2)*(y1-y2)+dphi*dphi);
}
1.7.4