fastjet_example.cc

Go to the documentation of this file.
00001 //STARTHEADER
00002 // $Id: fastjet_example.cc 958 2007-11-28 17:43:48Z cacciari $
00003 //
00004 // Copyright (c) 2005-2006, Matteo Cacciari and Gavin Salam
00005 //
00006 //----------------------------------------------------------------------
00007 // This file is part of FastJet.
00008 //
00009 //  FastJet is free software; you can redistribute it and/or modify
00010 //  it under the terms of the GNU General Public License as published by
00011 //  the Free Software Foundation; either version 2 of the License, or
00012 //  (at your option) any later version.
00013 //
00014 //  The algorithms that underlie FastJet have required considerable
00015 //  development and are described in hep-ph/0512210. If you use
00016 //  FastJet as part of work towards a scientific publication, please
00017 //  include a citation to the FastJet paper.
00018 //
00019 //  FastJet is distributed in the hope that it will be useful,
00020 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00021 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00022 //  GNU General Public License for more details.
00023 //
00024 //  You should have received a copy of the GNU General Public License
00025 //  along with FastJet; if not, write to the Free Software
00026 //  Foundation, Inc.:
00027 //      59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00028 //----------------------------------------------------------------------
00029 //ENDHEADER
00030 
00031 
00032 //----------------------------------------------------------------------
00033 // fastjet example program. 
00034 //
00035 // Compile it with: make fastjet_example
00036 // run it with    : ./fastjet_example < data/single-event.dat
00037 //
00038 // People who are familiar with the ktjet package are encouraged to
00039 // compare this file to the ktjet_example.cc program which does the
00040 // same thing in the ktjet framework.
00041 //----------------------------------------------------------------------
00042 #include "fastjet/PseudoJet.hh"
00043 #include "fastjet/ClusterSequence.hh"
00044 #include<iostream> // needed for io
00045 #include<sstream>  // needed for internal io
00046 #include<vector> 
00047 
00048 using namespace std;
00049 
00050 // a declaration of a function that pretty prints a list of jets
00051 void print_jets (const fastjet::ClusterSequence &, 
00052                  const vector<fastjet::PseudoJet> &);
00053 
00055 int main (int argc, char ** argv) {
00056   
00057   vector<fastjet::PseudoJet> input_particles;
00058   
00059   // read in input particles
00060   double px, py , pz, E;
00061   while (cin >> px >> py >> pz >> E) {
00062     // create a fastjet::PseudoJet with these components and put it onto
00063     // back of the input_particles vector
00064     input_particles.push_back(fastjet::PseudoJet(px,py,pz,E)); 
00065   }
00066   
00067   // create an object that represents your choice of jet algorithm and 
00068   // the associated parameters
00069   double Rparam = 1.0;
00070   fastjet::Strategy strategy = fastjet::Best;
00071   fastjet::RecombinationScheme recomb_scheme = fastjet::E_scheme;
00072   fastjet::JetDefinition jet_def(fastjet::kt_algorithm, Rparam, recomb_scheme, strategy);
00073 
00074   // run the jet clustering with the above jet definition
00075   fastjet::ClusterSequence clust_seq(input_particles, jet_def);
00076 
00077   // tell the user what was done
00078   cout << "Ran " << jet_def.description() << endl;
00079   cout << "Strategy adopted by FastJet was "<<
00080        clust_seq.strategy_string()<<endl<<endl;
00081 
00082   // extract the inclusive jets with pt > 5 GeV
00083   double ptmin = 5.0;
00084   vector<fastjet::PseudoJet> inclusive_jets = clust_seq.inclusive_jets(ptmin);
00085 
00086   // print them out
00087   cout << "Printing inclusive jets with pt > "<< ptmin<<" GeV\n";
00088   cout << "---------------------------------------\n";
00089   print_jets(clust_seq, inclusive_jets);
00090   cout << endl;
00091 
00092   // extract the exclusive jets with dcut = 25 GeV^2 
00093   double dcut = 25.0;
00094   vector<fastjet::PseudoJet> exclusive_jets = clust_seq.exclusive_jets(dcut);
00095 
00096   // print them out
00097   cout << "Printing exclusive jets with dcut = "<< dcut<<" GeV^2\n";
00098   cout << "--------------------------------------------\n";
00099   print_jets(clust_seq, exclusive_jets);
00100 
00101 
00102 }
00103 
00104 
00105 //----------------------------------------------------------------------
00107 void print_jets (const fastjet::ClusterSequence & clust_seq, 
00108                  const vector<fastjet::PseudoJet> & jets) {
00109 
00110   // sort jets into increasing pt
00111   vector<fastjet::PseudoJet> sorted_jets = sorted_by_pt(jets);  
00112 
00113   // label the columns
00114   printf("%5s %15s %15s %15s %15s\n","jet #", "rapidity", 
00115          "phi", "pt", "n constituents");
00116   
00117   // print out the details for each jet
00118   for (unsigned int i = 0; i < sorted_jets.size(); i++) {
00119     int n_constituents = clust_seq.constituents(sorted_jets[i]).size();
00120     printf("%5u %15.8f %15.8f %15.8f %8u\n",
00121            i, sorted_jets[i].rap(), sorted_jets[i].phi(),
00122            sorted_jets[i].perp(), n_constituents);
00123   }
00124 
00125 }

Generated on Tue Dec 18 17:05:02 2007 for fastjet by  doxygen 1.5.2