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
|
//----------------------------------------------------------------------
/// \file
/// \page Example03 03 - using plugins
///
/// fastjet plugins example program:
/// we illustrate the plugin usage
/// here, we use the SISCone plugin though different choices are possible
/// see the output of 'fastjet-config --list-plugins' for more details
///
/// Note that when using plugins, the code needs to be linked against
/// the libfastjetplugins library (with the default monolithic
/// build. For non-monolithic build, individual libraries have to be
/// used for each plugin).
/// This is ensured in practice by calling
/// fastjet-config --libs --plugins
///
/// run it with : ./03-plugin < data/single-event.dat
///
/// Source code: 03-plugin.cc
//----------------------------------------------------------------------
//STARTHEADER
// $Id: 03-plugin.cc 2684 2011-11-14 07:41:44Z 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
#include "fastjet/ClusterSequence.hh"
#include <iostream> // needed for io
#include <cstdio> // needed for io
// include the SISCone plugin header if enabled
#include "fastjet/config.h"
#ifdef FASTJET_ENABLE_PLUGIN_SISCONE
#include "fastjet/SISConePlugin.hh"
#else
#warning "SISCone plugin not enabled. Skipping the example"
#endif // FASTJET_ENABLE_PLUGIN_SISCONE
using namespace std;
int main(){
#ifdef FASTJET_ENABLE_PLUGIN_SISCONE
// read in input particles
//----------------------------------------------------------
vector<fastjet::PseudoJet> 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 a jet definition fron a plugin.
// It basically requires declaring a JetDefinition from a pointer to
// the plugin
//
// we will use the SISCone plugin here. Its (mandatory) parameters
// are a cone radius and an overlap threshold, plus other optional
// parameters
//
// for other plugin, see individual documentations for a description
// of their parameters
//
// the list of available plugins for a given build of FastJet can be
// obtained using
// fastjet-config --list-plugins
// from the command line.
//----------------------------------------------------------
double cone_radius = 0.7;
double overlap_threshold = 0.75;
fastjet::SISConePlugin siscone(cone_radius, overlap_threshold);
fastjet::JetDefinition jet_def(& siscone);
// run the jet clustering with the above jet definition
//----------------------------------------------------------
fastjet::ClusterSequence clust_seq(input_particles, jet_def);
// get the resulting jets ordered in pt
//----------------------------------------------------------
double ptmin = 5.0;
vector<fastjet::PseudoJet> inclusive_jets = sorted_by_pt(clust_seq.inclusive_jets(ptmin));
// tell the user what was done
// - the description of the algorithm used
// - extract the inclusive jets with pt > 5 GeV
// show the output as
// {index, rap, phi, pt}
//----------------------------------------------------------
cout << "Ran " << jet_def.description() << endl;
// label the columns
printf("%5s %15s %15s %15s\n","jet #", "rapidity", "phi", "pt");
// print out the details for each jet
for (unsigned int i = 0; i < inclusive_jets.size(); i++) {
printf("%5u %15.8f %15.8f %15.8f\n",
i, inclusive_jets[i].rap(), inclusive_jets[i].phi(),
inclusive_jets[i].perp());
}
#endif
return 0;
}
|