FastJet  3.3.0
Recluster.hh
1 #ifndef __FASTJET_TOOLS_RECLUSTER_HH__
2 #define __FASTJET_TOOLS_RECLUSTER_HH__
3 
4 // $Id: Recluster.hh 3753 2014-12-17 15:19:55Z salam $
5 //
6 // Copyright (c) 2014, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
7 //
8 //----------------------------------------------------------------------
9 // This file is part of FastJet
10 //
11 // It is free software; you can redistribute it and/or modify it under
12 // the terms of the GNU General Public License as published by the
13 // Free Software Foundation; either version 2 of the License, or (at
14 // your option) any later version.
15 //
16 // It is distributed in the hope that it will be useful, but WITHOUT
17 // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
18 // or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
19 // License for more details.
20 //
21 // You should have received a copy of the GNU General Public License
22 // along with this code. If not, see <http://www.gnu.org/licenses/>.
23 //----------------------------------------------------------------------
24 
25 #include <fastjet/JetDefinition.hh>
26 #include <fastjet/FunctionOfPseudoJet.hh> // to derive Recluster from FOfPJ<PJ>
27 #include <iostream>
28 #include <string>
29 
30 // TODO:
31 //
32 // - maintain Voronoi areas? Requires CSAB:has_voronoi_area() {->UNASSIGNED}
33 // - make sure the description of the class is OK
34 
35 
36 
37 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
38 
39 //----------------------------------------------------------------------
40 /// @ingroup tools_generic
41 /// \class Recluster
42 /// Recluster a jet's constituents with a new jet definition.
43 ///
44 /// When Recluster is constructed from a JetDefinition, it is that
45 /// definition that will be used to obtain the new jets. The user may
46 /// then decide if the recombiner should be the one from that jet
47 /// definition or if it should be acquired from the jet being
48 /// processed (the default).
49 ///
50 /// Alternatively, Recluster can be constructed from a jet algorithm
51 /// and an optional radius. In that case the recombiner is
52 /// systematically obtained from the jet being processed (unless you
53 /// call set_acquire_recombiner(false)). If only the jet algorithm is
54 /// specified, a default radius of max_allowable_R will be assumed if
55 /// needed.
56 ///
57 /// Recluster has two possible behaviours:
58 ///
59 /// - if it is constructed with keep=keep_only_hardest the hardest
60 /// inclusive jet is returned as a "standard" jet with an
61 /// associated cluster sequence (unless there were no inclusive
62 /// jets, in which case a zero jet is returned, with no associated
63 /// cluster sequence)
64 ///
65 /// - if it is constructed with keep=keep_all
66 /// all the inclusive jets are joined into a composite jet
67 ///
68 /// [Note that since the structure of the resulting PseudoJet depends
69 /// on its usage, this class inherits from
70 /// FunctionOfPseudoJet<PseudoJet> (including a description) rather
71 /// than being a full-fledged Transformer]
72 ///
73 class Recluster : public FunctionOfPseudoJet<PseudoJet> {
74 public:
75  /// the various options for the output of Recluster
76  enum Keep{
77  /// keep only the hardest inclusive jet and return a "standard" jet with
78  /// an associated ClusterSequence [this will be the default]
80  /// keep all the inclusive jets. result() will join them into a composite
81  /// jet
82  keep_all
83  };
84 
85  /// default constructor (uses an undefined JetDefinition, and so cannot
86  /// be used directly).
87  Recluster() : _new_jet_def(), _acquire_recombiner(true),
88  _keep(keep_only_hardest), _cambridge_optimisation_enabled(true){}
89 
90  /// Constructs a Recluster object that reclusters a jet into a new jet
91  /// using a generic JetDefinition
92  ///
93  /// \param new_jet_def the jet definition applied to do the reclustering
94  /// \param acquire_recombiner
95  /// when true, the reclustering will guess the
96  /// recombiner from the input jet instead of
97  /// the one in new_jet_def. An error is then
98  /// thrown if no consistent recombiner is found
99  /// \param keep_in Recluster::keep_only_hardest: the result is
100  /// the hardest inclusive jet after reclustering,
101  /// returned as a "standard" jet.
102  /// Recluster::keep_all: the result is a
103  /// composite jet with the inclusive jets as pieces.
104  Recluster(const JetDefinition & new_jet_def,
105  bool acquire_recombiner_in = false,
106  Keep keep_in = keep_only_hardest)
107  : _new_jet_def(new_jet_def), _acquire_recombiner(acquire_recombiner_in),
108  _keep(keep_in), _cambridge_optimisation_enabled(true) {}
109 
110  /// Constructs a Recluster object that reclusters a jet into a new jet
111  /// using a JetAlgorithm and its parameters
112  ///
113  /// \param new_jet_alg the jet algorithm applied to obtain the new clustering
114  /// \param new_jet_radius the jet radius
115  /// \param keep_in Recluster::keep_only_hardest: the result is
116  /// the hardest inclusive jet after reclustering,
117  /// returned as a "standard" jet.
118  /// Recluster::keep_all: the result is a
119  /// composite jet with the inclusive jets as pieces.
120  ///
121  /// This ctor will always acquire the recombiner from the jet being
122  /// reclustered (it will throw if none can be found). If you wish
123  /// to use Recluster with an algorithm that requires an extra
124  /// parameter (like the genkt algorithm), please specify the jet
125  /// definition fully using the constructor above.
126  Recluster(JetAlgorithm new_jet_alg, double new_jet_radius, Keep keep_in = keep_only_hardest);
127 
128  /// constructor with just a jet algorithm, but no jet radius. If the
129  /// algorithm requires a jet radius, JetDefinition::max_allowable_R will be used.
130  ///
131  Recluster(JetAlgorithm new_jet_alg, Keep keep_in = keep_only_hardest);
132 
133  /// default dtor
134  virtual ~Recluster(){}
135 
136  //----------------------------------------------------------------------
137  // tweaking the behaviour and corresponding enquiry functions
138 
139  /// set whether the reclustering should attempt to acquire a
140  /// recombiner from the input jet
141  void set_acquire_recombiner(bool acquire) {_acquire_recombiner = acquire;}
142 
143  /// returns true if this reclusterer is set to acquire the
144  /// recombiner from the input jet
145  bool acquire_recombiner() const{ return _acquire_recombiner;}
146 
147 
148  /// sets whether to try to optimise reclustering with
149  /// Cambridge/Aachen algorithms (by not reclustering if the
150  /// requested C/A reclustering can be obtained by using subjets of
151  /// an input C/A jet or one composed of multiple C/A pieces from the
152  /// same clustering sequence). By default this is enabled, and
153  /// _should_ always be correct; disable it to test this statement!
154  void set_cambridge_optimisation(bool enabled){ _cambridge_optimisation_enabled = enabled;}
155  /// sets whether to try to optimise reclustering with Cambridge/Aachen algorithms (US spelling!)
156  void set_cambridge_optimization(bool enabled){ _cambridge_optimisation_enabled = enabled;}
157 
158  /// returns true if the reclusterer tries to optimise reclustering
159  /// with Cambridge/Aachen algorithms
160  bool cambridge_optimization(){return _cambridge_optimisation_enabled;}
161  bool cambridge_optimisation(){return _cambridge_optimisation_enabled;}
162 
163  /// set the behaviour with regards to keeping all resulting jets or
164  /// just the hardest.
165  void set_keep(Keep keep_in) {_keep = keep_in;}
166 
167  /// returns the current "keep" mode i.e. whether only the hardest
168  /// inclusive jet is returned or all of them (see the Keep enum above)
169  Keep keep() const{ return _keep;}
170 
171 
172  //----------------------------------------------------------------------
173  // retrieving info about the behaviour
174 
175  /// class description
176  virtual std::string description() const;
177 
178 
179  //----------------------------------------------------------------------
180  // core action of ths class
181 
182  /// runs the reclustering and sets kept and rejected to be the jets
183  /// of interest (with non-zero rho, they will have been
184  /// subtracted). Normally this will be accessed through the base
185  /// class's operator().
186  ///
187  /// \param jet the jet that gets reclustered
188  /// \return the reclustered jet
189  virtual PseudoJet result(const PseudoJet & jet) const;
190 
191  /// A lower-level method that does the actual work of reclustering
192  /// the input jet. The resulting jets are stored in output_jets.
193  /// The jet definition that has been used can be accessed from the
194  /// output_jets' ClusterSequence.
195  ///
196  /// \param input_jet the (input) jet that one wants to recluster
197  /// \param output_jets inclusive jets resulting from the new clustering
198  ///
199  /// Returns true if the C/A optimisation has been used (this means
200  /// that generate_output_jet then has to watch out for non-explicit-ghost
201  /// areas that might be leftover)
202  bool get_new_jets_and_def(const PseudoJet & input_jet,
203  std::vector<PseudoJet> & output_jets) const;
204 
205  /// given a set of inclusive jets and a jet definition used, create the
206  /// resulting PseudoJet;
207  ///
208  /// If ca_optimisation_used then special care will be taken in
209  /// deciding whether the final jet can legitimately have an area.
210  PseudoJet generate_output_jet(std::vector<PseudoJet> & incljets,
211  bool ca_optimisation_used) const;
212 
213 
214 private:
215  /// set the reclustered elements in the simple case of C/A+C/A
216  void _recluster_ca(const std::vector<PseudoJet> & all_pieces,
217  std::vector<PseudoJet> & incljets,
218  double Rfilt) const;
219 
220  /// set the reclustered elements in the generic re-clustering case
221  void _recluster_generic(const PseudoJet & jet,
222  std::vector<PseudoJet> & incljets,
223  const JetDefinition & new_jet_def,
224  bool do_areas) const;
225 
226  // a series of checks
227  //--------------------------------------------------------------------
228  /// get the pieces down to the fundamental pieces
229  bool _get_all_pieces(const PseudoJet &jet, std::vector<PseudoJet> &all_pieces) const;
230 
231  /// associate the proper recombiner taken from the underlying pieces
232  /// (an error is thrown if the pieces do no share a common
233  /// recombiner)
234  void _acquire_recombiner_from_pieces(const std::vector<PseudoJet> &all_pieces,
235  JetDefinition &new_jet_def) const;
236 
237  /// check if one can apply the simplified trick for C/A subjets
238  bool _check_ca(const std::vector<PseudoJet> &all_pieces,
239  const JetDefinition &new_jet_def) const;
240 
241  /// check if the jet (or all its pieces) have explicit ghosts
242  /// (assuming the jet has area support
243  ///
244  /// Note that if the jet has an associated cluster sequence that is no
245  /// longer valid, an error will be thrown
246  bool _check_explicit_ghosts(const std::vector<PseudoJet> &all_pieces) const;
247 
248  JetDefinition _new_jet_def; ///< the jet definition to use to extract the jets
249  bool _acquire_recombiner; ///< get the recombiner from the input
250  ///< jet rather than from _new_jet_def
251  Keep _keep; ///< dictates which inclusive jets are kept and
252  ///< which are returned (see Keep above)
253 
254  bool _cambridge_optimisation_enabled; ///<enable the checks to
255  ///< perform optimisation when
256  ///< C/A reclustering is asked
257 
258  static LimitedWarning _explicit_ghost_warning;
259 };
260 
261 FASTJET_END_NAMESPACE // defined in fastjet/internal/base.hh
262 
263 #endif // __FASTJET_TOOLS_RECLUSTER_HH__
void set_cambridge_optimisation(bool enabled)
sets whether to try to optimise reclustering with Cambridge/Aachen algorithms (by not reclustering if...
Definition: Recluster.hh:154
void set_acquire_recombiner(bool acquire)
set whether the reclustering should attempt to acquire a recombiner from the input jet ...
Definition: Recluster.hh:141
class to provide facilities for giving warnings up to some maximum number of times and to provide glo...
keep only the hardest inclusive jet and return a "standard" jet with an associated ClusterSequence [t...
Definition: Recluster.hh:79
Keep
the various options for the output of Recluster
Definition: Recluster.hh:76
virtual ~Recluster()
default dtor
Definition: Recluster.hh:134
Keep keep() const
returns the current "keep" mode i.e.
Definition: Recluster.hh:169
bool cambridge_optimization()
returns true if the reclusterer tries to optimise reclustering with Cambridge/Aachen algorithms ...
Definition: Recluster.hh:160
Recluster a jet&#39;s constituents with a new jet definition.
Definition: Recluster.hh:73
base class providing interface for a generic function of a PseudoJet
bool acquire_recombiner() const
returns true if this reclusterer is set to acquire the recombiner from the input jet ...
Definition: Recluster.hh:145
Recluster(const JetDefinition &new_jet_def, bool acquire_recombiner_in=false, Keep keep_in=keep_only_hardest)
Constructs a Recluster object that reclusters a jet into a new jet using a generic JetDefinition...
Definition: Recluster.hh:104
JetAlgorithm
the various families of jet-clustering algorithm
Class to contain pseudojets, including minimal information of use to jet-clustering routines...
Definition: PseudoJet.hh:67
Recluster()
default constructor (uses an undefined JetDefinition, and so cannot be used directly).
Definition: Recluster.hh:87
void set_cambridge_optimization(bool enabled)
sets whether to try to optimise reclustering with Cambridge/Aachen algorithms (US spelling!) ...
Definition: Recluster.hh:156
void set_keep(Keep keep_in)
set the behaviour with regards to keeping all resulting jets or just the hardest. ...
Definition: Recluster.hh:165
class that is intended to hold a full definition of the jet clusterer