FastJet  3.1.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
fastjet_areas.cc
1 
2 //STARTHEADER
3 // $Id: fastjet_areas.cc 2704 2011-11-16 11:11:10Z soyez $
4 //
5 // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
6 //
7 //----------------------------------------------------------------------
8 // This file is part of FastJet.
9 //
10 // FastJet is free software; you can redistribute it and/or modify
11 // it under the terms of the GNU General Public License as published by
12 // the Free Software Foundation; either version 2 of the License, or
13 // (at your option) any later version.
14 //
15 // The algorithms that underlie FastJet have required considerable
16 // development and are described in hep-ph/0512210. If you use
17 // FastJet as part of work towards a scientific publication, please
18 // include a citation to the FastJet paper.
19 //
20 // FastJet is distributed in the hope that it will be useful,
21 // but WITHOUT ANY WARRANTY; without even the implied warranty of
22 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
23 // GNU General Public License for more details.
24 //
25 // You should have received a copy of the GNU General Public License
26 // along with FastJet. If not, see <http://www.gnu.org/licenses/>.
27 //----------------------------------------------------------------------
28 //ENDHEADER
29 
30 
31 
32 
33 
34 
35 
36 //----------------------------------------------------------------------
37 // fastjet areas example program.
38 //
39 // Compile it with: make fastjet_areas
40 // run it with : ./fastjet_areas < data/single-event.dat
41 //
42 //----------------------------------------------------------------------
43 #include "fastjet/PseudoJet.hh"
44 #include "fastjet/ClusterSequenceArea.hh"
45 #include "fastjet/ClusterSequencePassiveArea.hh"
46 
47 // get info on how fastjet was configured
48 #include "fastjet/config.h"
49 
50 #ifdef ENABLE_PLUGIN_SISCONE
51 #include "fastjet/SISConePlugin.hh"
52 #endif
53 
54 #include<iostream> // needed for io
55 #include<sstream> // needed for internal io
56 #include<vector>
57 
58 using namespace std;
59 
60 // a declaration of a function that pretty prints a list of jets
61 void print_jets (const vector<fastjet::PseudoJet> &);
62 
63 /// an example program showing how to use fastjet
64 int main () {
65 
66  vector<fastjet::PseudoJet> input_particles;
67 
68  // read in input particles
69  double px, py , pz, E;
70  while (cin >> px >> py >> pz >> E) {
71  // create a fastjet::PseudoJet with these components and put it onto
72  // back of the input_particles vector
73  input_particles.push_back(fastjet::PseudoJet(px,py,pz,E));
74  }
75 
76  // create an object that represents your choice of jet algorithm, and
77  // the associated parameters
78  double Rparam = 1.0;
80  fastjet::JetDefinition jet_def(fastjet::kt_algorithm, Rparam, strategy);
81  //fastjet::JetDefinition jet_def(fastjet::cambridge_algorithm, Rparam, strategy);
82  //fastjet::JetDefinition jet_def(fastjet::antikt_algorithm, Rparam, strategy);
83  //fastjet::JetDefinition jet_def(new fastjet::SISConePlugin(Rparam,0.75));
84 
85  // create an object that specifies how we to define the area
86  fastjet::AreaDefinition area_def;
87  bool use_voronoi = false;
88  if (!use_voronoi) {
89  double ghost_etamax = 6.0;
90  double ghost_area = 0.01;
91  int active_area_repeats = 1;
92 
93  // now create the object that holds info about ghosts, and from that
94  // get an area definition
95  fastjet::GhostedAreaSpec ghost_spec(ghost_etamax, active_area_repeats,
96  ghost_area);
97  area_def = fastjet::AreaDefinition(fastjet::active_area,ghost_spec);
98  //area_def = fastjet::AreaDefinition(fastjet::passive_area,ghost_spec);
99  } else {
100  double effective_Rfact = 1.0;
101  area_def = fastjet::VoronoiAreaSpec(effective_Rfact);
102  }
103 
104  // run the jet clustering with the above jet definition
105  fastjet::ClusterSequenceArea clust_seq(input_particles,
106  jet_def, area_def);
107  // you can also run the individual area classes directly
108  //fastjet::ClusterSequencePassiveArea clust_seq(input_particles, jet_def,
109  // area_def.ghost_spec());
110 
111  // you may want to find out how much area in a given range (|y|<range)
112  // is empty of real jets (or corresponds to pure "ghost" jets).
113  //double range = 4.0;
114  //cout << clust_seq.empty_area(range) << endl;
115  //cout << clust_seq.n_empty_jets(range) << endl;
116 
117  // tell the user what was done
118  cout << "Jet definition was: " << jet_def.description() << endl;
119  cout << "Area definition was: " << area_def.description() << endl;
120  cout << "Strategy adopted by FastJet was "<<
121  clust_seq.strategy_string()<<endl<<endl;
122 
123  // extract the inclusive jets with pt > 5 GeV, sorted by pt
124  double ptmin = 5.0;
125  vector<fastjet::PseudoJet> inclusive_jets = clust_seq.inclusive_jets(ptmin);
126 
127  // print them out
128  cout << "Printing inclusive jets with pt > "<< ptmin<<" GeV\n";
129  cout << "---------------------------------------\n";
130  print_jets(inclusive_jets);
131  cout << endl;
132 
133 
134  cout << "Number of unclustered particles: "
135  << clust_seq.unclustered_particles().size() << endl;
136 
137 
138 }
139 
140 
141 //----------------------------------------------------------------------
142 /// a function that pretty prints a list of jets
143 void print_jets (const vector<fastjet::PseudoJet> & unsorted_jets) {
144 
145  // sort jets into increasing pt
146  vector<fastjet::PseudoJet> jets = sorted_by_pt(unsorted_jets);
147 
148  printf(" ijet rap phi Pt area +- err\n");
149  for (unsigned int j = 0; j < jets.size(); j++) {
150 
151  double area = jets[j].area();
152  double area_error = jets[j].area_error();
153 
154  printf("%5u %9.5f %8.5f %10.3f %8.3f +- %6.3f\n",j,jets[j].rap(),
155  jets[j].phi(),jets[j].perp(), area, area_error);
156  }
157 
158 
159 }
160 
161 
162 
163 
164 
165 
166 
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2
Definition: PseudoJet.cc:799
General class for user to obtain ClusterSequence with additional area information.
the longitudinally invariant kt algorithm
class that holds a generic area definition
int main()
an example program showing how to use fastjet
Definition: 01-basic.cc:50
Specification for the computation of the Voronoi jet area.
std::string description() const
return a description of the current area definition
void print_jets(const vector< fastjet::PseudoJet > &)
a function that pretty prints a list of jets
automatic selection of the best (based on N), including the LazyTiled strategies that are new to FJ3...
Parameters to configure the computation of jet areas using ghosts.
Strategy
the various options for the algorithmic strategy to adopt in clustering events with kt and cambridge ...
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
Definition: PseudoJet.hh:67
class that is intended to hold a full definition of the jet clusterer