D0RunIIConePlugin is a plugin for fastjet (v2.1 upwards) that provides an interface to the D0 version of Run-II iterative cone algorithm with midpoint seeds (also known as the Iterative Legacy Cone Algorithm, ILCA). More...
#include <D0RunIIConePlugin.hh>
Public Member Functions | |
D0RunIIConePlugin (double cone_radius, double min_jet_Et, double split_ratio=_DEFAULT_split_ratio) | |
A D0RunIIConePlugin constructor which sets the "free" parameters of the algorithm:. | |
double | cone_radius () const |
double | min_jet_Et () const |
double | split_ratio () const |
double | far_ratio () const |
double | Et_min_ratio () const |
bool | kill_duplicate () const |
double | duplicate_dR () const |
double | duplicate_dPT () const |
double | search_factor () const |
double | pT_min_leading_protojet () const |
double | pT_min_second_protojet () const |
int | merge_max () const |
double | pT_min_nomerge () const |
double | overlap_threshold () const |
access the split_ratio() also by the name overlap_threshold() | |
virtual std::string | description () const |
return a textual description of the jet-definition implemented in this plugin | |
virtual void | run_clustering (ClusterSequence &) const |
given a ClusterSequence that has been filled up with initial particles, the following function should fill up the rest of the ClusterSequence, using the following member functions of ClusterSequence:
| |
virtual double | R () const |
the plugin mechanism's standard way of accessing the jet radius | |
Private Attributes | |
double | _cone_radius |
double | _min_jet_Et |
double | _split_ratio |
double | _far_ratio |
double | _Et_min_ratio |
bool | _kill_duplicate |
double | _duplicate_dR |
double | _duplicate_dPT |
double | _search_factor |
double | _pT_min_leading_protojet |
double | _pT_min_second_protojet |
int | _merge_max |
double | _pT_min_nomerge |
Static Private Attributes | |
static const double | _DEFAULT_split_ratio = 0.5 |
static const double | _DEFAULT_far_ratio = 0.5 |
static const double | _DEFAULT_Et_min_ratio = 0.5 |
static const bool | _DEFAULT_kill_duplicate = true |
static const double | _DEFAULT_duplicate_dR = 0.005 |
static const double | _DEFAULT_duplicate_dPT = 0.01 |
static const double | _DEFAULT_search_factor = 1.0 |
static const double | _DEFAULT_pT_min_leading_protojet = 0. |
static const double | _DEFAULT_pT_min_second_protojet = 0. |
static const int | _DEFAULT_merge_max = 10000 |
static const double | _DEFAULT_pT_min_nomerge = 0. |
D0RunIIConePlugin is a plugin for fastjet (v2.1 upwards) that provides an interface to the D0 version of Run-II iterative cone algorithm with midpoint seeds (also known as the Iterative Legacy Cone Algorithm, ILCA).
The D0 code has been taken from Lars Sonnenschein's web-space http://www-d0.fnal.gov/~sonne/D0RunIIcone.tgz
The version of the D0 Run II code distributed here has been modified by the FastJet authors, so as to provide access to the contents of the jets (as is necessary for the plugin). This does not modify the results of the clustering.
Definition at line 56 of file D0RunIIConePlugin.hh.
fastjet::D0RunIIConePlugin::D0RunIIConePlugin | ( | double | cone_radius, | |
double | min_jet_Et, | |||
double | split_ratio = _DEFAULT_split_ratio | |||
) | [inline] |
A D0RunIIConePlugin constructor which sets the "free" parameters of the algorithm:.
The remaining parameters of the algorithm are not to be modified if the algorithm is to correspond to the one actually used by D0.
Definition at line 76 of file D0RunIIConePlugin.hh.
00078 : 00079 _cone_radius (cone_radius ), 00080 _min_jet_Et (min_jet_Et ), 00081 _split_ratio (split_ratio ), 00082 _far_ratio (_DEFAULT_far_ratio ), 00083 _Et_min_ratio (_DEFAULT_Et_min_ratio ), 00084 _kill_duplicate (_DEFAULT_kill_duplicate ), 00085 _duplicate_dR (_DEFAULT_duplicate_dR ), 00086 _duplicate_dPT (_DEFAULT_duplicate_dPT ), 00087 _search_factor (_DEFAULT_search_factor ), 00088 _pT_min_leading_protojet(_DEFAULT_pT_min_leading_protojet), 00089 _pT_min_second_protojet (_DEFAULT_pT_min_second_protojet ), 00090 _merge_max (_DEFAULT_merge_max ), 00091 _pT_min_nomerge (_DEFAULT_pT_min_nomerge ) 00092 { 00093 // nothing to be done here! 00094 }
double fastjet::D0RunIIConePlugin::cone_radius | ( | ) | const [inline] |
Definition at line 97 of file D0RunIIConePlugin.hh.
Referenced by description(), and run_clustering().
00097 { return _cone_radius ;} //= 0.5;
string fastjet::D0RunIIConePlugin::description | ( | ) | const [virtual] |
return a textual description of the jet-definition implemented in this plugin
Implements JetDefinition::Plugin.
Definition at line 59 of file D0RunIIConePlugin.cc.
References cone_radius(), min_jet_Et(), and split_ratio().
00059 { 00060 ostringstream desc; 00061 00062 desc << "D0 Run II Improved Legacy (midpoint) cone jet algorithm, with "; 00063 desc << "cone_radius = " << cone_radius () << ", " 00064 << "min_jet_Et = " << min_jet_Et () << ", " 00065 << "split_ratio = " << split_ratio (); 00066 00067 return desc.str(); 00068 }
double fastjet::D0RunIIConePlugin::duplicate_dPT | ( | ) | const [inline] |
Definition at line 104 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00104 { return _duplicate_dPT ;} // =0.01;
double fastjet::D0RunIIConePlugin::duplicate_dR | ( | ) | const [inline] |
Definition at line 103 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00103 { return _duplicate_dR ;} // =0.005;
double fastjet::D0RunIIConePlugin::Et_min_ratio | ( | ) | const [inline] |
Definition at line 101 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00101 { return _Et_min_ratio ;} // =0.5;
double fastjet::D0RunIIConePlugin::far_ratio | ( | ) | const [inline] |
Definition at line 100 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00100 { return _far_ratio ;} // =0.5;
bool fastjet::D0RunIIConePlugin::kill_duplicate | ( | ) | const [inline] |
Definition at line 102 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00102 { return _kill_duplicate ;} // =true;
int fastjet::D0RunIIConePlugin::merge_max | ( | ) | const [inline] |
Definition at line 108 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00108 { return _merge_max ;} // =10000;
double fastjet::D0RunIIConePlugin::min_jet_Et | ( | ) | const [inline] |
Definition at line 98 of file D0RunIIConePlugin.hh.
Referenced by description(), and run_clustering().
00098 { return _min_jet_Et ;} //= 8.0;
double fastjet::D0RunIIConePlugin::overlap_threshold | ( | ) | const [inline] |
access the split_ratio() also by the name overlap_threshold()
Definition at line 113 of file D0RunIIConePlugin.hh.
00113 {return split_ratio();}
double fastjet::D0RunIIConePlugin::pT_min_leading_protojet | ( | ) | const [inline] |
Definition at line 106 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00106 { return _pT_min_leading_protojet;} // =0.;
double fastjet::D0RunIIConePlugin::pT_min_nomerge | ( | ) | const [inline] |
Definition at line 109 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00109 { return _pT_min_nomerge ;} // =0.;
double fastjet::D0RunIIConePlugin::pT_min_second_protojet | ( | ) | const [inline] |
Definition at line 107 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00107 { return _pT_min_second_protojet ;} // =0.;
virtual double fastjet::D0RunIIConePlugin::R | ( | ) | const [inline, virtual] |
the plugin mechanism's standard way of accessing the jet radius
Implements JetDefinition::Plugin.
Definition at line 119 of file D0RunIIConePlugin.hh.
00119 {return cone_radius();}
void fastjet::D0RunIIConePlugin::run_clustering | ( | ClusterSequence & | ) | const [virtual] |
given a ClusterSequence that has been filled up with initial particles, the following function should fill up the rest of the ClusterSequence, using the following member functions of ClusterSequence:
..)
Implements JetDefinition::Plugin.
Definition at line 71 of file D0RunIIConePlugin.cc.
References cone_radius(), duplicate_dPT(), duplicate_dR(), Et_min_ratio(), far_ratio(), fastjet::d0::ILConeAlgorithm< Item >::ilcv, ClusterSequence::jets(), kill_duplicate(), fastjet::d0::ILConeAlgorithm< Item >::makeClusters(), merge_max(), min_jet_Et(), ClusterSequence::plugin_record_iB_recombination(), ClusterSequence::plugin_record_ij_recombination(), pT_min_leading_protojet(), pT_min_nomerge(), pT_min_second_protojet(), search_factor(), and split_ratio().
00071 { 00072 00073 // create the entities needed by the D0 code 00074 vector<HepEntity> entities(clust_seq.jets().size()); 00075 list<const HepEntity * > ensemble; 00076 for (unsigned i = 0; i < clust_seq.jets().size(); i++) { 00077 entities[i].Fill(clust_seq.jets()[i].E(), 00078 clust_seq.jets()[i].px(), 00079 clust_seq.jets()[i].py(), 00080 clust_seq.jets()[i].pz(), 00081 i); 00082 // use only the particles that do not have infinite rapidity 00083 if (abs(entities[i].pz) < entities[i].E) { 00084 ensemble.push_back(& (entities[i])); 00085 } 00086 } 00087 00088 // prepare the D0 algorithm 00089 ILConeAlgorithm<HepEntity> 00090 ilegac(cone_radius(), 00091 min_jet_Et(), 00092 split_ratio(), 00093 far_ratio(), 00094 Et_min_ratio(), 00095 kill_duplicate(), 00096 duplicate_dR(), 00097 duplicate_dPT(), 00098 search_factor(), 00099 pT_min_leading_protojet(), 00100 pT_min_second_protojet(), 00101 merge_max(), 00102 pT_min_nomerge()); 00103 00104 // run the algorithm 00105 float Item_ET_Threshold = 0.; 00106 list<HepEntity> jets; 00107 ilegac.makeClusters(jets, ensemble, Item_ET_Threshold); 00108 00109 // now transfer the information about the jets into the 00110 // FastJet structure 00111 for(int i = ilegac.ilcv.size()-1; i >= 0; i--) { 00112 00113 std::list<const HepEntity*> tlist = ilegac.ilcv[i].LItems(); 00114 std::list<const HepEntity*>::iterator tk; 00115 00116 // get first particle in list 00117 tk = tlist.begin(); 00118 00119 // if there is no particle, just discard it 00120 // Note: this unexpected behaviour has been observed when the 00121 // min_jet_Et parameter was set to 0 00122 if (tk==tlist.end()) 00123 continue; 00124 00125 int jet_k = (*tk)->index; 00126 // now merge with remaining particles in list 00127 tk++; 00128 for (; tk != tlist.end(); tk++) { 00129 int jet_i = jet_k; 00130 int jet_j = (*tk)->index; 00131 // do a fake recombination step with dij=0 00132 double dij = 0.0; 00133 clust_seq.plugin_record_ij_recombination(jet_i, jet_j, dij, jet_k); 00134 } 00135 00136 // NB: put a sensible looking d_iB just to be nice... 00137 double d_iB = clust_seq.jets()[jet_k].perp2(); 00138 clust_seq.plugin_record_iB_recombination(jet_k, d_iB); 00139 00140 } 00141 }
double fastjet::D0RunIIConePlugin::search_factor | ( | ) | const [inline] |
Definition at line 105 of file D0RunIIConePlugin.hh.
Referenced by run_clustering().
00105 { return _search_factor ;} // =1.0;
double fastjet::D0RunIIConePlugin::split_ratio | ( | ) | const [inline] |
Definition at line 99 of file D0RunIIConePlugin.hh.
Referenced by description(), and run_clustering().
00099 { return _split_ratio ;} //= 0.5;
double fastjet::D0RunIIConePlugin::_cone_radius [private] |
Definition at line 124 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_duplicate_dPT = 0.01 [static, private] |
Definition at line 148 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_duplicate_dR = 0.005 [static, private] |
Definition at line 147 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_Et_min_ratio = 0.5 [static, private] |
Definition at line 145 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_far_ratio = 0.5 [static, private] |
Definition at line 144 of file D0RunIIConePlugin.hh.
const bool fastjet::D0RunIIConePlugin::_DEFAULT_kill_duplicate = true [static, private] |
Definition at line 146 of file D0RunIIConePlugin.hh.
const int fastjet::D0RunIIConePlugin::_DEFAULT_merge_max = 10000 [static, private] |
Definition at line 152 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_pT_min_leading_protojet = 0. [static, private] |
Definition at line 150 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_pT_min_nomerge = 0. [static, private] |
Definition at line 153 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_pT_min_second_protojet = 0. [static, private] |
Definition at line 151 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_search_factor = 1.0 [static, private] |
Definition at line 149 of file D0RunIIConePlugin.hh.
const double fastjet::D0RunIIConePlugin::_DEFAULT_split_ratio = 0.5 [static, private] |
Definition at line 143 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_duplicate_dPT [private] |
Definition at line 134 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_duplicate_dR [private] |
Definition at line 133 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_Et_min_ratio [private] |
Definition at line 131 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_far_ratio [private] |
Definition at line 130 of file D0RunIIConePlugin.hh.
bool fastjet::D0RunIIConePlugin::_kill_duplicate [private] |
Definition at line 132 of file D0RunIIConePlugin.hh.
int fastjet::D0RunIIConePlugin::_merge_max [private] |
Definition at line 138 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_min_jet_Et [private] |
Definition at line 125 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_pT_min_leading_protojet [private] |
Definition at line 136 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_pT_min_nomerge [private] |
Definition at line 139 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_pT_min_second_protojet [private] |
Definition at line 137 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_search_factor [private] |
Definition at line 135 of file D0RunIIConePlugin.hh.
double fastjet::D0RunIIConePlugin::_split_ratio [private] |
Definition at line 126 of file D0RunIIConePlugin.hh.