58 #include <fastjet/ClusterSequence.hh> 59 #include <fastjet/tools/MassDropTagger.hh> 60 #include <fastjet/tools/Filter.hh> 85 class FlavourRecombiner :
public DefRecomb {
90 virtual std::string description()
const {
107 ostream & operator<<(ostream &,
const PseudoJet &);
115 vector<PseudoJet> particles;
122 while (getline(cin,line)) {
123 if (line.substr(0,1) ==
"#") {
continue;}
124 istringstream linestream(line);
126 linestream >> px >> py >> pz >> E;
130 if (! (linestream >> btag)) btag = 0;
134 particle.set_user_index(btag);
135 particles.push_back(particle);
145 FlavourRecombiner flav_recombiner;
146 JetDefinition jet_def(cambridge_algorithm, R, &flav_recombiner);
151 vector<PseudoJet> jets =
sorted_by_pt(cs.inclusive_jets());
153 cout <<
"Ran: " << jet_def.description() << endl << endl;
154 cout <<
"Hardest jet: " << jets[0] << endl << endl;
169 cout <<
"No substructure found" << endl;
175 cout <<
"Found suitable pair of subjets: " << endl;
176 cout <<
" " << parent1 << endl;
177 cout <<
" " << parent2 << endl;
178 cout <<
"Total = " << endl;
179 cout <<
" " << tagged << endl;
198 double Rbb = parent1.
delta_R(parent2);
199 double Rfilt = min(Rbb/2, 0.3);
201 cout <<
"Subjet separation (Rbb) = " << Rbb <<
", Rfilt = " << Rfilt << endl;
209 const vector<PseudoJet> & filtered_pieces = filtered.
pieces();
210 cout <<
"Filtered pieces are " << endl;
211 for (
unsigned i = 0; i < nfilt && i < filtered_pieces.size(); i++) {
212 cout <<
" " << filtered_pieces[i] << endl;
214 cout <<
"Filtered total is " << endl;
215 cout <<
" " << filtered << endl;
223 ostream & operator<<(ostream & ostr,
const PseudoJet & jet) {
224 ostr <<
"pt, y, phi =" 225 <<
" " << setw(10) << jet.
perp()
226 <<
" " << setw(6) << jet.
rap()
227 <<
" " << setw(6) << jet.
phi()
228 <<
", mass = " << setw(10) << jet.
m()
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2
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...
virtual void recombine(const PseudoJet &pa, const PseudoJet &pb, PseudoJet &pab) const override
recombine pa and pb and put result into pab
Selector SelectorNHardest(unsigned int n)
select the n hardest objects
Class that helps perform filtering (Butterworth, Davison, Rubin and Salam, arXiv:0802.2470) and trimming (Krohn, Thaler and Wang, arXiv:0912.1342) on jets, optionally in conjunction with subtraction (Cacciari and Salam, arXiv:0707.1378).
Class that helps perform 2-pronged boosted tagging using the "mass-drop" technique (with asymmetry cu...
A class that will provide the recombination scheme facilities and/or allow a user to extend these fac...
double perp() const
returns the scalar transverse momentum
int main()
an example program showing how to use fastjet
virtual std::vector< PseudoJet > pieces() const
retrieve the pieces that make up the jet.
double m() const
returns the invariant mass (If m2() is negative then -sqrt(-m2()) is returned, as in CLHEP) ...
RecombinationScheme
The various recombination schemes.
double delta_R(const PseudoJet &other) const
return the cylinder (rap-phi) distance between this jet and another, .
const TransformerType::StructureType & structure_of() const
this is a helper to access any structure created by a Transformer (that is, of type Transformer::Stru...
double phi() const
returns phi (in the range 0..2pi)
double rap() const
returns the rapidity or some large value when the rapidity is infinite
virtual std::string description() const override
return a textual description of the recombination scheme implemented here
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
int user_index() const
return the user_index,
class that is intended to hold a full definition of the jet clusterer