FastJet 3.0.2
|
00001 /// \file 00002 /// \page Examples FastJet examples 00003 /// 00004 /// The FastJet examples have been organised by order of complexity, 00005 /// starting by the simplest case and introducing features one after 00006 /// another. 00007 /// - \subpage Example01 00008 /// - \subpage Example02 00009 /// - \subpage Example03 00010 /// - \subpage Example04 00011 /// - \subpage Example05 00012 /// - \subpage Example06 00013 /// - \subpage Example07 (\subpage Example07old "old version") 00014 /// - \subpage Example08 00015 /// - \subpage Example09 00016 /// - \subpage Example10 00017 /// - \subpage Example11 00018 /// - \subpage Example12 (\subpage Example12old "old version") 00019 /// - \subpage Example13 00020 /// - \subpage Example14 00021 00022 //STARTHEADER 00023 // $Id: fastjet_example.cc 2684 2011-11-14 07:41:44Z soyez $ 00024 // 00025 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez 00026 // 00027 //---------------------------------------------------------------------- 00028 // This file is part of FastJet. 00029 // 00030 // FastJet is free software; you can redistribute it and/or modify 00031 // it under the terms of the GNU General Public License as published by 00032 // the Free Software Foundation; either version 2 of the License, or 00033 // (at your option) any later version. 00034 // 00035 // The algorithms that underlie FastJet have required considerable 00036 // development and are described in hep-ph/0512210. If you use 00037 // FastJet as part of work towards a scientific publication, please 00038 // include a citation to the FastJet paper. 00039 // 00040 // FastJet is distributed in the hope that it will be useful, 00041 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00042 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00043 // GNU General Public License for more details. 00044 // 00045 // You should have received a copy of the GNU General Public License 00046 // along with FastJet. If not, see <http://www.gnu.org/licenses/>. 00047 //---------------------------------------------------------------------- 00048 //ENDHEADER 00049 00050 00051 //---------------------------------------------------------------------- 00052 // fastjet example program. 00053 // 00054 // Compile it with: make fastjet_example 00055 // run it with : ./fastjet_example < data/single-event.dat 00056 // 00057 // People who are familiar with the ktjet package are encouraged to 00058 // compare this file to the ktjet_example.cc program which does the 00059 // same thing in the ktjet framework. 00060 //---------------------------------------------------------------------- 00061 #include "fastjet/PseudoJet.hh" 00062 #include "fastjet/ClusterSequence.hh" 00063 #include<iostream> // needed for io 00064 #include<sstream> // needed for internal io 00065 #include<vector> 00066 #include <cstdio> 00067 00068 using namespace std; 00069 00070 // a declaration of a function that pretty prints a list of jets 00071 void print_jets (const vector<fastjet::PseudoJet> &); 00072 00073 /// an example program showing how to use fastjet 00074 int main () { 00075 00076 vector<fastjet::PseudoJet> input_particles; 00077 00078 // Read in input particles 00079 double px, py , pz, E; 00080 while (cin >> px >> py >> pz >> E) { 00081 // create a fastjet::PseudoJet with these components and put it onto 00082 // back of the input_particles vector 00083 input_particles.push_back(fastjet::PseudoJet(px,py,pz,E)); 00084 } 00085 00086 // create an object that represents your choice of jet algorithm and 00087 // the associated parameters 00088 double Rparam = 1.0; 00089 fastjet::Strategy strategy = fastjet::Best; 00090 fastjet::RecombinationScheme recomb_scheme = fastjet::E_scheme; 00091 fastjet::JetDefinition jet_def(fastjet::kt_algorithm, Rparam, recomb_scheme, strategy); 00092 00093 // run the jet clustering with the above jet definition 00094 fastjet::ClusterSequence clust_seq(input_particles, jet_def); 00095 00096 // tell the user what was done 00097 cout << "Ran " << jet_def.description() << endl; 00098 cout << "Strategy adopted by FastJet was "<< 00099 clust_seq.strategy_string()<<endl<<endl; 00100 00101 // extract the inclusive jets with pt > 5 GeV 00102 double ptmin = 5.0; 00103 vector<fastjet::PseudoJet> inclusive_jets = clust_seq.inclusive_jets(ptmin); 00104 00105 // print them out 00106 cout << "Printing inclusive jets with pt > "<< ptmin<<" GeV\n"; 00107 cout << "---------------------------------------\n"; 00108 print_jets(inclusive_jets); 00109 cout << endl; 00110 00111 // extract the exclusive jets with dcut = 25 GeV^2 00112 double dcut = 25.0; 00113 vector<fastjet::PseudoJet> exclusive_jets = clust_seq.exclusive_jets(dcut); 00114 00115 // print them out 00116 cout << "Printing exclusive jets with dcut = "<< dcut<<" GeV^2\n"; 00117 cout << "--------------------------------------------\n"; 00118 print_jets(exclusive_jets); 00119 00120 00121 } 00122 00123 00124 //---------------------------------------------------------------------- 00125 /// a function that pretty prints a list of jets 00126 void print_jets (const vector<fastjet::PseudoJet> & jets) { 00127 00128 // sort jets into increasing pt 00129 vector<fastjet::PseudoJet> sorted_jets = sorted_by_pt(jets); 00130 00131 // label the columns 00132 printf("%5s %15s %15s %15s %15s\n","jet #", "rapidity", 00133 "phi", "pt", "n constituents"); 00134 00135 // print out the details for each jet 00136 for (unsigned int i = 0; i < sorted_jets.size(); i++) { 00137 // the following is not super efficient since it creates an 00138 // intermediate constituents vector 00139 int n_constituents = sorted_jets[i].constituents().size(); 00140 printf("%5u %15.8f %15.8f %15.8f %8u\n", 00141 i, sorted_jets[i].rap(), sorted_jets[i].phi(), 00142 sorted_jets[i].perp(), n_constituents); 00143 } 00144 00145 }