FastJet  3.3.0
PxConePlugin.hh
1 //FJSTARTHEADER
2 // $Id: PxConePlugin.hh 3433 2014-07-23 08:17:03Z salam $
3 //
4 // Copyright (c) 2005-2014, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
5 //
6 //----------------------------------------------------------------------
7 // This file is part of FastJet.
8 //
9 // FastJet is free software; you can redistribute it and/or modify
10 // it under the terms of the GNU General Public License as published by
11 // the Free Software Foundation; either version 2 of the License, or
12 // (at your option) any later version.
13 //
14 // The algorithms that underlie FastJet have required considerable
15 // development. They are described in the original FastJet paper,
16 // hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
17 // FastJet as part of work towards a scientific publication, please
18 // quote the version you use and include a citation to the manual and
19 // optionally also to hep-ph/0512210.
20 //
21 // FastJet is distributed in the hope that it will be useful,
22 // but WITHOUT ANY WARRANTY; without even the implied warranty of
23 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 // GNU General Public License for more details.
25 //
26 // You should have received a copy of the GNU General Public License
27 // along with FastJet. If not, see <http://www.gnu.org/licenses/>.
28 //----------------------------------------------------------------------
29 //FJENDHEADER
30 
31 #ifndef __PXCONEPLUGIN_HH__
32 #define __PXCONEPLUGIN_HH__
33 
34 #include "fastjet/JetDefinition.hh"
35 
36 // questionable whether this should be in fastjet namespace or not...
37 
38 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
39 
40 //----------------------------------------------------------------------
41 //
42 /// @ingroup plugins
43 /// \class PxConePlugin
44 /// Implementation of the PxCone algorithm (plugin for fastjet v2.1 upwards)
45 ///
46 /// PxConePlugin is a plugin for fastjet (v2.1 upwards) that provides
47 /// an interface to the fortran pxcone iterative cone algorithm with
48 /// midpoint seeds.
49 ///
50 /// Pxcone was written by Luis del Pozo and Michael H. Seymour. It is
51 /// not a "supported" program, so if you encounter problems, you are
52 /// on your own...
53 ///
54 /// Note that pxcone sometimes encounters non-stable iterations; in
55 /// such cases it returns an error -- the plugin propagates this by
56 /// throwing a fastjet::Error exception; if the user wishes to have
57 /// robust code, they should catch this exception.
58 ///
59 /// Pxcone has a hard-coded limit (by default 4000) on the maximum
60 /// number of particles and protojets; if the number of particles or
61 /// protojets exceeds this, again a fastjet::Error exception will be
62 /// thrown.
63 ///
64 /// The functionality of pxcone is described at
65 /// http://www.hep.man.ac.uk/u/wplano/ConeJet.ps
66 //
67 //----------------------------------------------------------------------
69 public:
70 
71  /// constructor for the PxConePlugin, whose arguments have the
72  /// following meaning:
73  ///
74  /// - the cone_radius is as usual in cone algorithms
75  ///
76  /// - stables cones (protojets) below min_jet_energy are discarded
77  /// before calling the splitting procedure to resolve overlaps
78  /// (called epslon in pxcone).
79  ///
80  /// - when two protojets overlap, if
81  /// (overlapping_Et)/(Et_of_softer_protojet) < overlap_threshold
82  /// the overlapping energy is split between the two protojets;
83  /// otherwise the less energetic protojet is discarded. Called
84  /// ovlim in pxcone.
85  ///
86  /// - pxcone carries out p-scheme recombination, and the resulting
87  /// jets are massless; setting E_scheme_jets = true (default
88  /// false) doesn't change the jet composition, but the final
89  /// momentum sum for the jets is carried out by direct
90  /// four-vector addition instead of p-scheme recombination.
91  ///
92  PxConePlugin (double cone_radius_in ,
93  double min_jet_energy_in = 5.0 ,
94  double overlap_threshold_in = 0.5,
95  bool E_scheme_jets_in = false) :
96  _cone_radius (cone_radius_in ),
97  _min_jet_energy (min_jet_energy_in ),
98  _overlap_threshold (overlap_threshold_in),
99  _E_scheme_jets (E_scheme_jets_in ) {}
100 
101 
102  // some functions to return info about parameters ----------------
103 
104  /// the cone radius
105  double cone_radius () const {return _cone_radius ;}
106 
107  /// minimum jet energy (protojets below this are thrown own before
108  /// merging/splitting) -- called epslon in pxcone
109  double min_jet_energy () const {return _min_jet_energy ;}
110 
111  /// Maximum fraction of overlap energy in a jet -- called ovlim in pxcone.
112  double overlap_threshold () const {return _overlap_threshold ;}
113 
114  /// if true then the final jets are returned as the E-scheme recombination
115  /// of the particle momenta (by default, pxcone returns massless jets with
116  /// a mean phi,eta type of recombination); regardless of what is
117  /// returned, the internal pxcone jet-finding procedure is
118  /// unaffected.
119  bool E_scheme_jets() const {return _E_scheme_jets ;}
120 
121 
122  // the things that are required by base class
123  virtual std::string description () const;
124  virtual void run_clustering(ClusterSequence &) const;
125  /// the plugin mechanism's standard way of accessing the jet radius
126  virtual double R() const {return cone_radius();}
127 
128 private:
129 
130  double _cone_radius ;
131  double _min_jet_energy ;
132  double _overlap_threshold ;
133 
134  bool _E_scheme_jets;
135 
136  static bool _first_time;
137 
138  /// print a banner for reference to the 3rd-party code
139  void _print_banner(std::ostream *ostr) const;
140 };
141 
142 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh
143 
144 #endif // __PXCONEPLUGIN_HH__
double cone_radius() const
the cone radius
virtual double R() const
the plugin mechanism&#39;s standard way of accessing the jet radius
deals with clustering
double min_jet_energy() const
minimum jet energy (protojets below this are thrown own before merging/splitting) – called epslon in...
double overlap_threshold() const
Maximum fraction of overlap energy in a jet – called ovlim in pxcone.
Implementation of the PxCone algorithm (plugin for fastjet v2.1 upwards)
Definition: PxConePlugin.hh:68
bool E_scheme_jets() const
if true then the final jets are returned as the E-scheme recombination of the particle momenta (by de...
PxConePlugin(double cone_radius_in, double min_jet_energy_in=5.0, double overlap_threshold_in=0.5, bool E_scheme_jets_in=false)
constructor for the PxConePlugin, whose arguments have the following meaning:
Definition: PxConePlugin.hh:92
a class that allows a user to introduce their own "plugin" jet finder