58#include "fastjet/ClusterSequence.hh"
87class FlavourRecombiner :
public DefRecomb {
109ostream & operator<<(ostream &,
PseudoJet &);
115int main (
int argc,
char ** argv) {
117 vector<PseudoJet> particles;
124 while (getline(cin,line)) {
125 if (line.substr(0,1) ==
"#") {
continue;}
126 istringstream linestream(line);
128 linestream >> px >> py >> pz >> E;
132 if (! (linestream >> btag)) btag = 0;
136 particle.set_user_index(btag);
137 particles.push_back(particle);
147 FlavourRecombiner flav_recombiner;
153 vector<PseudoJet> jets =
sorted_by_pt(cs.inclusive_jets());
155 cout <<
"Ran: " << jet_def.description() << endl << endl;
156 cout <<
"Hardest jet: " << jets[0] << endl << endl;
165 double mass_drop_threshold = 0.667;
180 while ((had_parents = this_jet.
has_parents(parent1,parent2))) {
182 if (parent1.
m() < parent2.
m())
swap(parent1,parent2);
187 if (parent1.
m() < mass_drop_threshold * this_jet.
m() &&
188 parent1.
kt_distance(parent2) > pow(rtycut,2) * this_jet.
m2()) {
198 cout <<
"Did not find suitable hard substructure in this event." << endl;
202 cout <<
"Found suitable pair of subjets: " << endl;
203 cout <<
" " << parent1 << endl;
204 cout <<
" " << parent2 << endl;
205 cout <<
"Total = " << endl;
206 cout <<
" " << this_jet << endl << endl;
219 double Rfilt = min(Rbb/2, 0.3);
221 cout <<
"Subjet separation (Rbb) = " << Rbb <<
", Rfilt = " << Rfilt << endl;
223 double dcut = pow(Rfilt/R,2);
229 cout <<
"Filtered pieces are " << endl;
230 cout <<
" " << filt_subjets[0] << endl;
231 PseudoJet filtered_total = filt_subjets[0];
232 for (
unsigned i = 1; i < nfilt && i < filt_subjets.size(); i++) {
233 cout <<
" " << filt_subjets[i] << endl;
234 flav_recombiner.plus_equal(filtered_total, filt_subjets[i]);
236 cout <<
"Filtered total is " << endl;
237 cout <<
" " << filtered_total << endl;
246 ostr <<
"pt, y, phi ="
247 <<
" " << setw(10) << jet.
perp()
248 <<
" " << setw(6) << jet.
rap()
249 <<
" " << setw(6) << jet.
phi()
250 <<
", mass = " << setw(10) << jet.
m()
int main()
an example program showing how to use fastjet
A class that will provide the recombination scheme facilities and/or allow a user to extend these fac...
virtual void recombine(const PseudoJet &pa, const PseudoJet &pb, PseudoJet &pab) const override
recombine pa and pb and put result into pab
virtual std::string description() const override
return a textual description of the recombination scheme implemented here
class that is intended to hold a full definition of the jet clusterer
Class to contain pseudojets, including minimal information of use to jet-clustering routines.
double rap() const
returns the rapidity or some large value when the rapidity is infinite
void set_user_index(const int index)
set the user_index, intended to allow the user to add simple identifying information to a particle/je...
double squared_distance(const PseudoJet &other) const
returns squared cylinder (rap-phi) distance between this jet and another
double phi() const
returns phi (in the range 0..2pi)
double perp() const
returns the scalar transverse momentum
double m2() const
returns the squared invariant mass // like CLHEP
double kt_distance(const PseudoJet &other) const
returns kt distance (R=1) between this jet and another
virtual bool has_parents(PseudoJet &parent1, PseudoJet &parent2) const
check if it is the product of a recombination, in which case return the 2 parents through the 'parent...
std::vector< PseudoJet > exclusive_subjets(const double dcut) const
return a vector of all subjets of the current jet (in the sense of the exclusive algorithm) that woul...
double m() const
returns the invariant mass (If m2() is negative then -sqrt(-m2()) is returned, as in CLHEP)
int user_index() const
return the user_index,
RecombinationScheme
The various recombination schemes.
@ E_scheme
summing the 4-momenta
@ cambridge_algorithm
the longitudinally invariant variant of the cambridge algorithm (aka Aachen algoithm).
void swap(SharedPtr< T > &a, SharedPtr< T > &b)
swapping
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2