FastJet 3.4.3
Loading...
Searching...
No Matches
03-plugin.cc
Go to the documentation of this file.
1//----------------------------------------------------------------------
2/// \file
3/// \page Example03 03 - using plugins
4///
5/// fastjet plugins example program:
6/// we illustrate the plugin usage
7/// here, we use the SISCone plugin though different choices are possible
8/// see the output of 'fastjet-config --list-plugins' for more details
9///
10/// Note that when using plugins, the code needs to be linked against
11/// the libfastjetplugins library (with the default monolithic
12/// build. For non-monolithic build, individual libraries have to be
13/// used for each plugin).
14/// This is ensured in practice by calling
15/// fastjet-config --libs --plugins
16///
17/// run it with : ./03-plugin < data/single-event.dat
18///
19/// Source code: 03-plugin.cc
20//----------------------------------------------------------------------
21
22//FJSTARTHEADER
23// $Id$
24//
25// Copyright (c) 2005-2024, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
26//
27//----------------------------------------------------------------------
28// This file is part of FastJet.
29//
30// FastJet is free software; you can redistribute it and/or modify
31// it under the terms of the GNU General Public License as published by
32// the Free Software Foundation; either version 2 of the License, or
33// (at your option) any later version.
34//
35// The algorithms that underlie FastJet have required considerable
36// development. They are described in the original FastJet paper,
37// hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
38// FastJet as part of work towards a scientific publication, please
39// quote the version you use and include a citation to the manual and
40// optionally also to hep-ph/0512210.
41//
42// FastJet is distributed in the hope that it will be useful,
43// but WITHOUT ANY WARRANTY; without even the implied warranty of
44// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
45// GNU General Public License for more details.
46//
47// You should have received a copy of the GNU General Public License
48// along with FastJet. If not, see <http://www.gnu.org/licenses/>.
49//----------------------------------------------------------------------
50//FJENDHEADER
51
52#include "fastjet/ClusterSequence.hh"
53#include <iostream> // needed for io
54#include <cstdio> // needed for io
55
56// include the SISCone plugin header if enabled
57#include "fastjet/config.h"
58#ifdef FASTJET_ENABLE_PLUGIN_SISCONE
59#include "fastjet/SISConePlugin.hh"
60#else
61#warning "SISCone plugin not enabled. Skipping the example"
62#endif // FASTJET_ENABLE_PLUGIN_SISCONE
63
64
65using namespace std;
66
67int main(){
68
69#ifdef FASTJET_ENABLE_PLUGIN_SISCONE
70
71 // read in input particles
72 //----------------------------------------------------------
73 vector<fastjet::PseudoJet> input_particles;
74
75 double px, py , pz, E;
76 while (cin >> px >> py >> pz >> E) {
77 // create a fastjet::PseudoJet with these components and put it onto
78 // back of the input_particles vector
79 input_particles.push_back(fastjet::PseudoJet(px,py,pz,E));
80 }
81
82
83 // create a jet definition fron a plugin.
84 // It basically requires declaring a JetDefinition from a pointer to
85 // the plugin
86 //
87 // we will use the SISCone plugin here. Its (mandatory) parameters
88 // are a cone radius and an overlap threshold, plus other optional
89 // parameters
90 //
91 // for other plugin, see individual documentations for a description
92 // of their parameters
93 //
94 // the list of available plugins for a given build of FastJet can be
95 // obtained using
96 // fastjet-config --list-plugins
97 // from the command line.
98 //----------------------------------------------------------
99 double cone_radius = 0.7;
100 double overlap_threshold = 0.75;
101 fastjet::SISConePlugin siscone(cone_radius, overlap_threshold);
102 fastjet::JetDefinition jet_def(& siscone);
103
104
105 // run the jet clustering with the above jet definition
106 //----------------------------------------------------------
107 fastjet::ClusterSequence clust_seq(input_particles, jet_def);
108
109
110 // get the resulting jets ordered in pt
111 //----------------------------------------------------------
112 double ptmin = 5.0;
113 vector<fastjet::PseudoJet> inclusive_jets = sorted_by_pt(clust_seq.inclusive_jets(ptmin));
114
115
116 // tell the user what was done
117 // - the description of the algorithm used
118 // - extract the inclusive jets with pt > 5 GeV
119 // show the output as
120 // {index, rap, phi, pt}
121 //----------------------------------------------------------
122 cout << "Ran " << jet_def.description() << endl;
123
124 // label the columns
125 printf("%5s %15s %15s %15s\n","jet #", "rapidity", "phi", "pt");
126
127 // print out the details for each jet
128 for (unsigned int i = 0; i < inclusive_jets.size(); i++) {
129 printf("%5u %15.8f %15.8f %15.8f\n",
130 i, inclusive_jets[i].rap(), inclusive_jets[i].phi(),
131 inclusive_jets[i].perp());
132 }
133
134#endif
135
136 return 0;
137
138}
int main()
an example program showing how to use fastjet
Definition 01-basic.cc:52
deals with clustering
class that is intended to hold a full definition of the jet clusterer
Class to contain pseudojets, including minimal information of use to jet-clustering routines.
Definition PseudoJet.hh:68
Implementation of the SISCone algorithm (plugin for fastjet v2.1 upwards)
vector< PseudoJet > sorted_by_pt(const vector< PseudoJet > &jets)
return a vector of jets sorted into decreasing kt2
Definition PseudoJet.cc:871