40 #include "fastjet/PseudoJet.hh"
41 #include "fastjet/ClusterSequenceArea.hh"
42 #include "fastjet/Selector.hh"
43 #include "fastjet/tools/JetMedianBackgroundEstimator.hh"
44 #include "fastjet/tools/Subtractor.hh"
48 using namespace fastjet;
58 vector<PseudoJet> hard_event, full_event;
64 double particle_maxrap = 5.0;
68 while (getline(cin, line)) {
69 istringstream linestream(line);
72 if (line.substr(0,4) ==
"#END") {
break;}
73 if (line.substr(0,9) ==
"#SUBSTART") {
75 if (nsub == 1) hard_event = full_event;
78 if (line.substr(0,1) ==
"#") {
continue;}
80 linestream >> px >> py >> pz >> E;
85 if (abs(particle.rap()) <= particle_maxrap) full_event.push_back(particle);
89 if (nsub == 1) hard_event = full_event;
93 cerr <<
"Error: read empty event\n";
109 double ghost_maxrap = 6.0;
124 vector<PseudoJet> hard_jets =
sorted_by_pt(clust_seq_hard.inclusive_jets(ptmin));
125 vector<PseudoJet> full_jets =
sorted_by_pt(clust_seq_full.inclusive_jets(ptmin));
174 bkgd_estimator.set_particles(full_event);
181 cout <<
"Main clustering:" << endl;
183 cout <<
" Area: " << area_def.description() << endl;
184 cout <<
" Particles up to |y|=" << particle_maxrap << endl;
187 cout <<
"Background estimation:" << endl;
188 cout <<
" " << bkgd_estimator.description() << endl << endl;;
189 cout <<
" Giving, for the full event" << endl;
190 cout <<
" rho = " << bkgd_estimator.rho() << endl;
191 cout <<
" sigma = " << bkgd_estimator.sigma() << endl;
194 cout <<
"Jets above " << ptmin <<
" GeV in the hard event (" << hard_event.size() <<
" particles)" << endl;
195 cout <<
"---------------------------------------\n";
196 printf(
"%5s %15s %15s %15s %15s\n",
"jet #",
"rapidity",
"phi",
"pt",
"area");
197 for (
unsigned int i = 0; i < hard_jets.size(); i++) {
198 printf(
"%5u %15.8f %15.8f %15.8f %15.8f\n", i,
199 hard_jets[i].rap(), hard_jets[i].phi(), hard_jets[i].perp(),
200 hard_jets[i].area());
210 cout <<
"Jets above " << ptmin <<
" GeV in the full event (" << full_event.size() <<
" particles)" << endl;
211 cout <<
"---------------------------------------\n";
212 printf(
"%5s %15s %15s %15s %15s %15s %15s %15s\n",
"jet #",
"rapidity",
"phi",
"pt",
"area",
"rap_sub",
"phi_sub",
"pt_sub");
216 vector<PseudoJet> subtracted_jets = subtractor(full_jets);
218 for (
unsigned int i=0; i<full_jets.size(); i++){
220 if (subtracted_jets[i].perp2() >= ptmin*ptmin){
221 printf(
"%5u %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f\n", idx,
222 full_jets[i].rap(), full_jets[i].phi(), full_jets[i].perp(),
224 subtracted_jets[i].rap(), subtracted_jets[i].phi(),
225 subtracted_jets[i].perp());
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2
General class for user to obtain ClusterSequence with additional area information.
Selector SelectorNHardest(unsigned int n)
select the n hardest objects
Class that helps perform jet background subtraction.
Selector SelectorAbsRapMax(double absrapmax)
select objects with |rap| <= absrapmax
class that holds a generic area definition
int main()
an example program showing how to use fastjet
std::string description() const
return a textual description of the current jet definition
Class that encodes information about cuts and other selection criteria that can be applied to PseudoJ...
Parameters to configure the computation of jet areas using ghosts.
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
class that is intended to hold a full definition of the jet clusterer