60 #include "fastjet/ClusterSequence.hh"
61 #include "fastjet/TrackJetPlugin.hh"
70 FASTJET_BEGIN_NAMESPACE
77 class TrackJetParticlePtr{
79 TrackJetParticlePtr(
int i_index,
double i_perp2)
80 : index(i_index), perp2(i_perp2){}
85 bool operator <(
const TrackJetParticlePtr &other)
const {
86 return perp2>other.perp2;
94 bool TrackJetPlugin::_first_time =
true;
96 string TrackJetPlugin::description ()
const {
98 desc <<
"TrackJet algorithm with R = " << R();
107 vector<TrackJetParticlePtr> particle_list;
109 const vector<PseudoJet> & jets = clust_seq.
jets();
111 for (vector<PseudoJet>::const_iterator mom_it = jets.begin(); mom_it != jets.end(); mom_it++){
112 particle_list.push_back(TrackJetParticlePtr(index, mom_it->perp2()));
117 stable_sort(particle_list.begin(), particle_list.end());
123 vector<PseudoJet> tuned_particles = clust_seq.
jets();
124 vector<PseudoJet> tuned_tracks = clust_seq.
jets();
125 for (vector<PseudoJet>::iterator pit = tuned_particles.begin();
126 pit != tuned_particles.end(); pit++)
127 _jet_recombiner.preprocess(*pit);
128 for (vector<PseudoJet>::iterator pit = tuned_tracks.begin();
129 pit != tuned_tracks.end(); pit++)
130 _track_recombiner.preprocess(*pit);
134 list<int> sorted_pt_index;
135 for (vector<TrackJetParticlePtr>::iterator mom_it = particle_list.begin();
136 mom_it != particle_list.end(); mom_it++)
137 sorted_pt_index.push_back(mom_it->index);
140 while (sorted_pt_index.size()){
144 int current_jet_index = sorted_pt_index.front();
145 PseudoJet current_jet = tuned_particles[current_jet_index];
146 PseudoJet current_track = tuned_tracks[current_jet_index];
149 list<int>::iterator index_it = sorted_pt_index.begin();
150 sorted_pt_index.erase(index_it);
153 index_it = sorted_pt_index.begin();
154 while (index_it != sorted_pt_index.end()){
155 const PseudoJet & current_particle = tuned_particles[*index_it];
156 const PseudoJet & current_particle_track = tuned_tracks[*index_it];
159 double distance2 = current_track.
plain_distance(current_particle_track);
160 if (distance2 <= _radius2){
164 _jet_recombiner.recombine(current_jet, current_particle, new_jet);
165 _track_recombiner.recombine(current_track, current_particle_track, new_track);
170 current_jet = new_jet;
171 current_track = new_track;
172 current_jet_index = new_jet_index;
175 sorted_pt_index.erase(index_it);
179 index_it = sorted_pt_index.begin();
192 void TrackJetPlugin::_print_banner(ostream *ostr)
const{
193 if (! _first_time)
return;
199 (*ostr) <<
"#-------------------------------------------------------------------------" << endl;
200 (*ostr) <<
"# You are running the TrackJet plugin for FastJet. It is based on " << endl;
201 (*ostr) <<
"# the implementation by Andy Buckley and Manuel Bahr that is to be " << endl;
202 (*ostr) <<
"# found in Rivet 1.1.2. See http://www.hepforge.org/downloads/rivet. " << endl;
203 (*ostr) <<
"#-------------------------------------------------------------------------" << endl;
209 FASTJET_END_NAMESPACE
double plain_distance(const PseudoJet &other) const
returns squared cylinder (rap-phi) distance between this jet and another
const std::vector< PseudoJet > & jets() const
allow the user to access the internally stored _jets() array, which contains both the initial particl...
void plugin_record_ij_recombination(int jet_i, int jet_j, double dij, int &newjet_k)
record the fact that there has been a recombination between jets()[jet_i] and jets()[jet_k], with the specified dij, and return the index (newjet_k) allocated to the new jet, whose momentum is assumed to be the 4-vector sum of that of jet_i and jet_j
void plugin_record_iB_recombination(int jet_i, double diB)
record the fact that there has been a recombination between jets()[jet_i] and the beam...
static std::ostream * fastjet_banner_stream()
returns a pointer to the stream to be used to print banners (cout by default).
Class to contain pseudojets, including minimal information of use to jet-clustering routines...