#include "fastjet/PseudoJet.hh"
#include "fastjet/ClusterSequence.hh"
#include <iostream>
#include <sstream>
#include <vector>
#include <cstdio>
Go to the source code of this file.
Functions | |
void | print_jets (const fj::ClusterSequence &clust_seq, const vector< fj::PseudoJet > &jets) |
a function that pretty prints a list of jets | |
void | print_jet (const fj::ClusterSequence &clust_seq, const fj::PseudoJet &jet) |
print a single jet | |
void | print_jets_and_sub (const fj::ClusterSequence &clust_seq, const vector< fj::PseudoJet > &jets, double dcut) |
a function that pretty prints a list of jets | |
int | main (int argc, char **argv) |
an example program showing how to use fastjet |
int main | ( | int | argc, | |
char ** | argv | |||
) |
an example program showing how to use fastjet
Definition at line 70 of file fastjet_subjets.cc.
References cambridge_algorithm, print_jets(), and print_jets_and_sub().
00070 { 00071 00072 vector<fj::PseudoJet> input_particles; 00073 00074 // read in input particles 00075 double px, py , pz, E; 00076 while (cin >> px >> py >> pz >> E) { 00077 // create a fj::PseudoJet with these components and put it onto 00078 // back of the input_particles vector 00079 input_particles.push_back(fj::PseudoJet(px,py,pz,E)); 00080 } 00081 00082 // We'll do two subjet analyses, physically rather different 00083 double R = 1.0; 00084 //fj::JetDefinition kt_def(fj::kt_algorithm, R); 00085 fj::JetDefinition cam_def(fj::cambridge_algorithm, R); 00086 00087 // run the jet clustering with the above jet definition 00088 //fj::ClusterSequence kt_seq(input_particles, kt_def); 00089 fj::ClusterSequence cam_seq(input_particles, cam_def); 00090 00091 // tell the user what was done 00092 cout << "Ran " << cam_def.description() << endl; 00093 cout << "Strategy adopted by FastJet was "<< 00094 cam_seq.strategy_string()<<endl<<endl; 00095 00096 // extract the inclusive jets with pt > 5 GeV 00097 double ptmin = 5.0; 00098 vector<fj::PseudoJet> inclusive_jets = cam_seq.inclusive_jets(ptmin); 00099 00100 // for the subjets 00101 double smallR = 0.4; 00102 double dcut_cam = pow(smallR/R,2); 00103 00104 // print them out 00105 cout << "Printing inclusive jets (R = "<<R<<") with pt > "<< ptmin<<" GeV\n"; 00106 cout << "and their subjets with smallR = " << smallR << "\n"; 00107 cout << "---------------------------------------\n"; 00108 print_jets_and_sub(cam_seq, inclusive_jets, dcut_cam); 00109 cout << endl; 00110 00111 00112 // print them out 00113 vector<fj::PseudoJet> exclusive_jets = cam_seq.exclusive_jets(dcut_cam); 00114 cout << "Printing exclusive jets with dcut = "<< dcut_cam<<" \n"; 00115 cout << "--------------------------------------------\n"; 00116 print_jets(cam_seq, exclusive_jets); 00117 00118 00119 }
void print_jet | ( | const fj::ClusterSequence & | clust_seq, | |
const fj::PseudoJet & | jet | |||
) |
print a single jet
Definition at line 173 of file fastjet_subjets.cc.
Referenced by print_jets(), and print_jets_and_sub().
void print_jets | ( | const fj::ClusterSequence & | clust_seq, | |
const vector< fj::PseudoJet > & | jets | |||
) |
a function that pretty prints a list of jets
Definition at line 124 of file fastjet_subjets.cc.
References print_jet(), and sorted_by_pt().
00125 { 00126 00127 // sort jets into increasing pt 00128 vector<fj::PseudoJet> sorted_jets = sorted_by_pt(jets); 00129 00130 // label the columns 00131 printf("%5s %15s %15s %15s %15s\n","jet #", "rapidity", 00132 "phi", "pt", "n constituents"); 00133 00134 // print out the details for each jet 00135 for (unsigned int i = 0; i < sorted_jets.size(); i++) { 00136 printf("%5u ",i); 00137 print_jet(clust_seq, sorted_jets[i]); 00138 } 00139 00140 }
void print_jets_and_sub | ( | const fj::ClusterSequence & | clust_seq, | |
const vector< fj::PseudoJet > & | jets, | |||
double | dcut | |||
) |
a function that pretty prints a list of jets
Definition at line 145 of file fastjet_subjets.cc.
References print_jet(), and sorted_by_pt().
Referenced by main().
00147 { 00148 00149 // sort jets into increasing pt 00150 vector<fj::PseudoJet> sorted_jets = sorted_by_pt(jets); 00151 00152 // label the columns 00153 printf("%5s %15s %15s %15s %15s\n","jet #", "rapidity", 00154 "phi", "pt", "n constituents"); 00155 00156 // print out the details for each jet 00157 for (unsigned int i = 0; i < sorted_jets.size(); i++) { 00158 printf("%5u ",i); 00159 print_jet(clust_seq, sorted_jets[i]); 00160 vector<fj::PseudoJet> subjets = 00161 sorted_by_pt(clust_seq.exclusive_subjets(sorted_jets[i], dcut)); 00162 for (unsigned int j = 0; j < subjets.size(); j++) { 00163 printf(" -sub-%02u ",j); 00164 print_jet(clust_seq, subjets[j]); 00165 } 00166 } 00167 00168 }