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