|
fastjet 2.4.5
|
|
fastjet 2.4.5
|
#include "fastjet/ClusterSequence.hh"#include <iostream>#include <sstream>#include <iomanip>#include <cmath>
Go to the source code of this file.
Classes | |
| class | FlavourRecombiner |
Typedefs | |
| typedef fj::JetDefinition::DefaultRecombiner | DefRecomb |
| set up a class to give standard (by default E-scheme) recombination, with additional tracking of flavour information in the user_index. | |
Functions | |
| ostream & | operator<< (ostream &, fj::PseudoJet &) |
| forward declaration for printing out info about a jet | |
| int | main (int argc, char **argv) |
| typedef fj::JetDefinition::DefaultRecombiner DefRecomb |
set up a class to give standard (by default E-scheme) recombination, with additional tracking of flavour information in the user_index.
If you use this, you must explicitly set the user index to 0 for non-flavoured particles (the default value is -1);
This will work for native algorithms, but not for all plugins
Definition at line 73 of file fastjet_boosted_higgs.cc.
| int main | ( | int | argc, |
| char ** | argv | ||
| ) |
now do the subjet decomposition;
when unpeeling a C/A jet, often only a very soft piece may break off; the mass_drop_threshold indicates how much "lighter" the heavier of the two resulting pieces must be in order for us to consider that we've really seen some form of substructure
QCD backgrounds that give larger jet masses have a component where a quite soft gluon is emitted; to eliminate part of this one can place a cut on the asymmetry of the branching;
Here the cut is expressed in terms of y, the kt-distance scaled to the squared jet mass; an easier way to see it is in terms of a requirement on the momentum fraction in the splitting: z/(1-z) and (1-z)/z > rtycut^2 [the correspondence holds only at LO]
Definition at line 98 of file fastjet_boosted_higgs.cc.
References fastjet::cambridge_algorithm, fastjet::JetDefinition::description(), fastjet::ClusterSequence::exclusive_subjets(), fastjet::ClusterSequence::has_parents(), fastjet::ClusterSequence::inclusive_jets(), fastjet::PseudoJet::m(), fastjet::PseudoJet::m2(), fastjet::d0::inline_maths::min(), fastjet::PseudoJet::set_user_index(), and fastjet::sorted_by_pt().
{
vector<fj::PseudoJet> particles;
// read in data in format px py pz E b-tag [last of these is optional]
string line;
while (getline(cin,line)) {
if (line.substr(0,1) == "#") {continue;}
istringstream linestream(line);
double px,py,pz,E;
linestream >> px >> py >> pz >> E;
// optionally read in btag information
int btag;
if (! (linestream >> btag)) btag = 0;
// construct the particle
fj::PseudoJet particle(px,py,pz,E);
particle.set_user_index(btag); // btag info goes in user index, for flavour tracking
particles.push_back(particle);
}
// set up the jet finding
double R = 1.2;
FlavourRecombiner flav_recombiner; // for tracking flavour
fj::JetDefinition jet_def(fj::cambridge_algorithm, R, &flav_recombiner);
// run the jet finding; find the hardest jet
fj::ClusterSequence cs(particles, jet_def);
vector<fj::PseudoJet> jets = sorted_by_pt(cs.inclusive_jets());
cout << "Ran: " << jet_def.description() << endl << endl;
cout << "Hardest jet: " << jets[0] << endl << endl;
//
double mass_drop_threshold = 0.667;
double rtycut = 0.3;
fj::PseudoJet this_jet = jets[0], parent1, parent2;
bool had_parents;
while ((had_parents = cs.has_parents(this_jet,parent1,parent2))) {
// make parent1 the more massive jet
if (parent1.m() < parent2.m()) swap(parent1,parent2);
//
// if we pass the conditions on the mass drop and its degree of
// asymmetry (z/(1-z) \sim kt_dist/m^2 > rtycut), then we've found
// something interesting, so exit the loop
if (parent1.m() < mass_drop_threshold * this_jet.m() &&
parent1.kt_distance(parent2) > pow(rtycut,2) * this_jet.m2()) {
break;
} else {
// otherwise try a futher decomposition on the more massive jet
this_jet = parent1;
}
}
// look to see what we found
if (had_parents) {
cout << "Found suitable pair of subjets: " << endl;
cout << " " << parent1 << endl;
cout << " " << parent2 << endl;
cout << "Total = " << endl;
cout << " " << this_jet << endl << endl;
// next we "filter" it, to remove UE & pileup contamination
//
// [there are two ways of doing this; here we directly use the
// exsiting cluster sequence and find the exclusive subjets of
// this_jet (i.e. work backwards within the cs starting from
// this_jet); alternatively one can recluster just the
// constituents of the jet]
//
// first get separation between the subjets (called Rbb -- assuming it's a Higgs!)
double Rbb = sqrt(parent1.squared_distance(parent2));
double Rfilt = min(Rbb/2, 0.3); // somewhat arbitrary choice
unsigned nfilt = 3; // number of pieces we'll take
cout << "Subjet separation (Rbb) = " << Rbb << ", Rfilt = " << Rfilt << endl;
double dcut = pow(Rfilt/R,2); // for C/A get a view at Rfilt by
// using a dcut=(Rfilt/R)^2
vector<fj::PseudoJet> filt_subjets = sorted_by_pt(cs.exclusive_subjets(this_jet, dcut));
// now print out the filtered jets and reconstruct total
// at the same time
cout << "Filtered pieces are " << endl;
cout << " " << filt_subjets[0] << endl;
fj::PseudoJet filtered_total = filt_subjets[0];
for (unsigned i = 1; i < nfilt && i < filt_subjets.size(); i++) {
cout << " " << filt_subjets[i] << endl;
flav_recombiner.plus_equal(filtered_total, filt_subjets[i]);
}
cout << "Filtered total is " << endl;
cout << " " << filtered_total << endl;
} else {
cout << "Did not find suitable hard substructure in this event." << endl;
}
}
| ostream & operator<< | ( | ostream & | ostr, |
| fj::PseudoJet & | jet | ||
| ) |
forward declaration for printing out info about a jet
does the actual work for printing out a jet
Definition at line 217 of file fastjet_boosted_higgs.cc.
References fastjet::PseudoJet::m(), fastjet::PseudoJet::perp(), fastjet::PseudoJet::phi(), fastjet::PseudoJet::rap(), and fastjet::PseudoJet::user_index().
{
ostr << "pt, y, phi ="
<< " " << setw(10) << jet.perp()
<< " " << setw(6) << jet.rap()
<< " " << setw(6) << jet.phi()
<< ", mass = " << setw(10) << jet.m()
<< ", btag = " << jet.user_index();
return ostr;
}
1.7.4