#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 fastjet::ClusterSequence &clust_seq, const vector< fastjet::PseudoJet > &jets) |
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 56 of file fastjet_example.cc.
References Best, E_scheme, kt_algorithm, and print_jets().
00056 { 00057 00058 vector<fastjet::PseudoJet> input_particles; 00059 00060 // read in input particles 00061 double px, py , pz, E; 00062 while (cin >> px >> py >> pz >> E) { 00063 // create a fastjet::PseudoJet with these components and put it onto 00064 // back of the input_particles vector 00065 input_particles.push_back(fastjet::PseudoJet(px,py,pz,E)); 00066 } 00067 00068 // create an object that represents your choice of jet algorithm and 00069 // the associated parameters 00070 double Rparam = 1.0; 00071 fastjet::Strategy strategy = fastjet::Best; 00072 fastjet::RecombinationScheme recomb_scheme = fastjet::E_scheme; 00073 fastjet::JetDefinition jet_def(fastjet::kt_algorithm, Rparam, recomb_scheme, strategy); 00074 00075 // run the jet clustering with the above jet definition 00076 fastjet::ClusterSequence clust_seq(input_particles, jet_def); 00077 00078 // tell the user what was done 00079 cout << "Ran " << jet_def.description() << endl; 00080 cout << "Strategy adopted by FastJet was "<< 00081 clust_seq.strategy_string()<<endl<<endl; 00082 00083 // extract the inclusive jets with pt > 5 GeV 00084 double ptmin = 5.0; 00085 vector<fastjet::PseudoJet> inclusive_jets = clust_seq.inclusive_jets(ptmin); 00086 00087 // print them out 00088 cout << "Printing inclusive jets with pt > "<< ptmin<<" GeV\n"; 00089 cout << "---------------------------------------\n"; 00090 print_jets(clust_seq, inclusive_jets); 00091 cout << endl; 00092 00093 // extract the exclusive jets with dcut = 25 GeV^2 00094 double dcut = 25.0; 00095 vector<fastjet::PseudoJet> exclusive_jets = clust_seq.exclusive_jets(dcut); 00096 00097 // print them out 00098 cout << "Printing exclusive jets with dcut = "<< dcut<<" GeV^2\n"; 00099 cout << "--------------------------------------------\n"; 00100 print_jets(clust_seq, exclusive_jets); 00101 00102 00103 }
void print_jets | ( | const fastjet::ClusterSequence & | clust_seq, | |
const vector< fastjet::PseudoJet > & | jets | |||
) |
a function that pretty prints a list of jets
Definition at line 108 of file fastjet_example.cc.
References fastjet::d0::inline_maths::phi(), and sorted_by_pt().
00109 { 00110 00111 // sort jets into increasing pt 00112 vector<fastjet::PseudoJet> sorted_jets = sorted_by_pt(jets); 00113 00114 // label the columns 00115 printf("%5s %15s %15s %15s %15s\n","jet #", "rapidity", 00116 "phi", "pt", "n constituents"); 00117 00118 // print out the details for each jet 00119 for (unsigned int i = 0; i < sorted_jets.size(); i++) { 00120 int n_constituents = clust_seq.constituents(sorted_jets[i]).size(); 00121 printf("%5u %15.8f %15.8f %15.8f %8u\n", 00122 i, sorted_jets[i].rap(), sorted_jets[i].phi(), 00123 sorted_jets[i].perp(), n_constituents); 00124 } 00125 00126 }