fastjet_example_v1_interface.cc

Go to the documentation of this file.
00001 //STARTHEADER
00002 // $Id: fastjet_example_v1_interface.cc 432 2007-01-23 12:12:59Z salam $
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 using fastjet-v1 interface
00034 //
00035 // Compile it with: make fastjet_example_v1_interface
00036 // run it with    : ./fastjet_example_v1_interface < 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 "FjPseudoJet.hh"
00043 #include "FjClusterSequence.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 FjClusterSequence &, const vector<FjPseudoJet> &);
00052 
00054 int main (int argc, char ** argv) {
00055   
00056   vector<FjPseudoJet> input_particles;
00057   
00058   // read in input particles
00059   double px, py , pz, E;
00060   while (cin >> px >> py >> pz >> E) {
00061     // create a FjPseudoJet with these components and put it onto
00062     // back of the input_particles vector
00063     input_particles.push_back(FjPseudoJet(px,py,pz,E)); 
00064   }
00065 
00066   // run the jet clustering with option R=1.0 and strategy=Best
00067   double Rparam = 1.0;
00068   FjClusterSequence clust_seq(input_particles, Rparam, Best);
00069 
00070   // tell the user what was done
00071   cout << "Strategy adopted by FastJet was "<<
00072        clust_seq.strategy_string()<<endl<<endl;
00073 
00074   // extract the inclusive jets with pt > 5 GeV, sorted by pt
00075   double ptmin = 5.0;
00076   vector<FjPseudoJet> inclusive_jets = clust_seq.inclusive_jets(ptmin);
00077 
00078   // print them out
00079   cout << "Printing inclusive jets with pt > "<< ptmin<<" GeV\n";
00080   cout << "---------------------------------------\n";
00081   print_jets(clust_seq, inclusive_jets);
00082   cout << endl;
00083 
00084   // extract the exclusive jets with dcut = 25 GeV^2 
00085   double dcut = 25.0;
00086   vector<FjPseudoJet> exclusive_jets = clust_seq.exclusive_jets(dcut);
00087 
00088   // print them out
00089   cout << "Printing exclusive jets with dcut = "<< dcut<<" GeV^2\n";
00090   cout << "--------------------------------------------\n";
00091   print_jets(clust_seq, exclusive_jets);
00092 
00093 
00094 }
00095 
00096 
00097 //----------------------------------------------------------------------
00099 void print_jets (const FjClusterSequence & clust_seq, 
00100                  const vector<FjPseudoJet> & jets) {
00101 
00102   // sort jets into increasing pt
00103   vector<FjPseudoJet> sorted_jets = sorted_by_pt(jets);  
00104 
00105   // label the columns
00106   printf("%5s %15s %15s %15s %15s\n","jet #", "rapidity", 
00107          "phi", "pt", "n constituents");
00108   
00109   // print out the details for each jet
00110   for (unsigned int i = 0; i < sorted_jets.size(); i++) {
00111     int n_constituents = clust_seq.constituents(sorted_jets[i]).size();
00112     printf("%5u %15.8f %15.8f %15.8f %8u\n",
00113            i, sorted_jets[i].rap(), sorted_jets[i].phi(),
00114            sorted_jets[i].perp(), n_constituents);
00115   }
00116 
00117 }

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