FastJet 3.0.0
fastjet_areas.cc
00001 
00002 //STARTHEADER
00003 // $Id: fastjet_areas.cc 2614 2011-09-28 15:42:11Z salam $
00004 //
00005 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
00006 //
00007 //----------------------------------------------------------------------
00008 // This file is part of FastJet.
00009 //
00010 //  FastJet is free software; you can redistribute it and/or modify
00011 //  it under the terms of the GNU General Public License as published by
00012 //  the Free Software Foundation; either version 2 of the License, or
00013 //  (at your option) any later version.
00014 //
00015 //  The algorithms that underlie FastJet have required considerable
00016 //  development and are described in hep-ph/0512210. If you use
00017 //  FastJet as part of work towards a scientific publication, please
00018 //  include a citation to the FastJet paper.
00019 //
00020 //  FastJet is distributed in the hope that it will be useful,
00021 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00022 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00023 //  GNU General Public License for more details.
00024 //
00025 //  You should have received a copy of the GNU General Public License
00026 //  along with FastJet. If not, see <http://www.gnu.org/licenses/>.
00027 //----------------------------------------------------------------------
00028 //ENDHEADER
00029 
00030 
00031 
00032 
00033 
00034 
00035 
00036 //----------------------------------------------------------------------
00037 // fastjet areas example program. 
00038 //
00039 // Compile it with: make fastjet_areas
00040 // run it with    : ./fastjet_areas < data/single-event.dat
00041 //
00042 //----------------------------------------------------------------------
00043 #include "fastjet/PseudoJet.hh"
00044 #include "fastjet/ClusterSequenceArea.hh"
00045 #include "fastjet/ClusterSequencePassiveArea.hh"
00046 
00047 // get info on how fastjet was configured
00048 #include "fastjet/config.h"
00049 
00050 #ifdef ENABLE_PLUGIN_SISCONE
00051 #include "fastjet/SISConePlugin.hh"
00052 #endif
00053 
00054 #include<iostream> // needed for io
00055 #include<sstream>  // needed for internal io
00056 #include<vector> 
00057 
00058 using namespace std;
00059 
00060 // a declaration of a function that pretty prints a list of jets
00061 void print_jets (const vector<fastjet::PseudoJet> &);
00062 
00063 /// an example program showing how to use fastjet
00064 int main (int argc, char ** argv) {
00065   
00066   vector<fastjet::PseudoJet> input_particles;
00067   
00068   // read in input particles
00069   double px, py , pz, E;
00070   while (cin >> px >> py >> pz >> E) {
00071     // create a fastjet::PseudoJet with these components and put it onto
00072     // back of the input_particles vector
00073     input_particles.push_back(fastjet::PseudoJet(px,py,pz,E)); 
00074   }
00075 
00076   // create an object that represents your choice of jet algorithm, and 
00077   // the associated parameters
00078   double Rparam = 1.0;
00079   fastjet::Strategy strategy = fastjet::Best;
00080   fastjet::JetDefinition jet_def(fastjet::kt_algorithm, Rparam, strategy);
00081   //fastjet::JetDefinition jet_def(fastjet::cambridge_algorithm, Rparam, strategy);
00082   //fastjet::JetDefinition jet_def(fastjet::antikt_algorithm, Rparam, strategy);
00083   //fastjet::JetDefinition jet_def(new fastjet::SISConePlugin(Rparam,0.75));
00084 
00085   // create an object that specifies how we to define the area
00086   fastjet::AreaDefinition area_def;
00087   bool use_voronoi = false;
00088   if (!use_voronoi) {
00089     double ghost_etamax = 6.0;
00090     double ghost_area    = 0.01;
00091     int    active_area_repeats = 1;
00092 
00093     // now create the object that holds info about ghosts, and from that
00094     // get an area definition
00095     fastjet::GhostedAreaSpec ghost_spec(ghost_etamax, active_area_repeats, 
00096                                         ghost_area);
00097     area_def = fastjet::AreaDefinition(fastjet::active_area,ghost_spec);
00098     //area_def = fastjet::AreaDefinition(fastjet::passive_area,ghost_spec);
00099   } else {
00100     double effective_Rfact = 1.0;
00101     area_def = fastjet::VoronoiAreaSpec(effective_Rfact);
00102   }
00103 
00104   // run the jet clustering with the above jet definition
00105   fastjet::ClusterSequenceArea clust_seq(input_particles, 
00106                                              jet_def, area_def);
00107   // you can also run the individual area classes directly
00108   //fastjet::ClusterSequencePassiveArea clust_seq(input_particles, jet_def, 
00109   //                                              area_def.ghost_spec());
00110 
00111   // you may want to find out how much area in a given range (|y|<range)
00112   // is empty of real jets (or corresponds to pure "ghost" jets).
00113   //double range = 4.0;
00114   //cout << clust_seq.empty_area(range) << endl;
00115   //cout << clust_seq.n_empty_jets(range) << endl;
00116 
00117   // tell the user what was done
00118   cout << "Jet definition was: " << jet_def.description() << endl;
00119   cout << "Area definition was: " << area_def.description() << endl;
00120   cout << "Strategy adopted by FastJet was "<<
00121        clust_seq.strategy_string()<<endl<<endl;
00122 
00123   // extract the inclusive jets with pt > 5 GeV, sorted by pt
00124   double ptmin = 5.0;
00125   vector<fastjet::PseudoJet> inclusive_jets = clust_seq.inclusive_jets(ptmin);
00126 
00127   // print them out
00128   cout << "Printing inclusive jets with pt > "<< ptmin<<" GeV\n";
00129   cout << "---------------------------------------\n";
00130   print_jets(inclusive_jets);
00131   cout << endl;
00132 
00133   
00134   cout << "Number of unclustered particles: " 
00135        << clust_seq.unclustered_particles().size() << endl;
00136 
00137 
00138 }
00139 
00140 
00141 //----------------------------------------------------------------------
00142 /// a function that pretty prints a list of jets
00143 void print_jets (const vector<fastjet::PseudoJet> & unsorted_jets) {
00144 
00145   // sort jets into increasing pt
00146   vector<fastjet::PseudoJet> jets = sorted_by_pt(unsorted_jets);  
00147 
00148   printf(" ijet   rap      phi        Pt         area  +-   err\n");
00149   for (size_t j = 0; j < jets.size(); j++) {
00150 
00151     double area       = jets[j].area();
00152     double area_error = jets[j].area_error();
00153 
00154     printf("%5u %9.5f %8.5f %10.3f %8.3f +- %6.3f\n",j,jets[j].rap(),
00155            jets[j].phi(),jets[j].perp(), area, area_error);
00156   }
00157 
00158 
00159 }
00160 
00161 
00162 
00163 
00164 
00165 
00166 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends