source: HiSusy/trunk/Delphes-3.0.0/external/fastjet/JetDefinition.hh @ 1

Last change on this file since 1 was 1, checked in by zerwas, 11 years ago

first import of structure, PYTHIA8 and DELPHES

File size: 19.3 KB
Line 
1//STARTHEADER
2// $Id: JetDefinition.hh 2687 2011-11-14 11:17:51Z soyez $
3//
4// Copyright (c) 2005-2011, 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 and are described in hep-ph/0512210. If you use
16//  FastJet as part of work towards a scientific publication, please
17//  include a citation to the FastJet paper.
18//
19//  FastJet is distributed in the hope that it will be useful,
20//  but WITHOUT ANY WARRANTY; without even the implied warranty of
21//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
22//  GNU General Public License for more details.
23//
24//  You should have received a copy of the GNU General Public License
25//  along with FastJet. If not, see <http://www.gnu.org/licenses/>.
26//----------------------------------------------------------------------
27//ENDHEADER
28
29#ifndef __FASTJET_JETDEFINITION_HH__
30#define __FASTJET_JETDEFINITION_HH__
31
32#include<cassert>
33#include "fastjet/internal/numconsts.hh"
34#include "fastjet/PseudoJet.hh"
35#include<string>
36#include<memory>
37
38FASTJET_BEGIN_NAMESPACE      // defined in fastjet/internal/base.hh
39
40/// return a string containing information about the release
41//  NB: (implemented in ClusterSequence.cc but defined here because
42//  this is a visible location)
43std::string fastjet_version_string();
44
45//======================================================================
46/// the various options for the algorithmic strategy to adopt in
47/// clustering events with kt and cambridge style algorithms.
48enum Strategy {
49  /// fastest form about 500..10^4
50  N2MinHeapTiled   = -4, 
51  /// fastest from about 50..500
52  N2Tiled     = -3, 
53  /// legacy
54  N2PoorTiled = -2, 
55  /// fastest below 50
56  N2Plain     = -1, 
57  /// worse even than the usual N^3 algorithms
58  N3Dumb      =  0, 
59  /// automatic selection of the best (based on N)
60  Best        =  1, 
61  /// best of the NlnN variants -- best overall for N>10^4.
62  /// (Does not work for R>=2pi)
63  NlnN        =  2, 
64  /// legacy N ln N using 3pi coverage of cylinder.
65  /// (Does not work for R>=2pi)
66  NlnN3pi     =  3, 
67  /// legacy N ln N using 4pi coverage of cylinder
68  NlnN4pi     =  4,
69  /// Chan's closest pair method (in a variant with 4pi coverage),
70  /// for use exclusively with the Cambridge algorithm.
71  /// (Does not work for R>=2pi)
72  NlnNCam4pi   = 14,
73  /// Chan's closest pair method (in a variant with 2pi+2R coverage),
74  /// for use exclusively with the Cambridge algorithm.
75  /// (Does not work for R>=2pi)
76  NlnNCam2pi2R = 13,
77  /// Chan's closest pair method (in a variant with 2pi+minimal extra
78  /// variant), for use exclusively with the Cambridge algorithm.
79  /// (Does not work for R>=2pi)
80  NlnNCam      = 12, // 2piMultD
81  /// the plugin has been used...
82  plugin_strategy = 999
83};
84
85
86//======================================================================
87/// \enum JetAlgorithm
88/// the various families of jet-clustering algorithm
89enum JetAlgorithm {
90  /// the longitudinally invariant kt algorithm
91  kt_algorithm=0,
92  /// the longitudinally invariant variant of the cambridge algorithm
93  /// (aka Aachen algoithm).
94  cambridge_algorithm=1,
95  /// like the k_t but with distance measures
96  ///       dij = min(1/kti^2,1/ktj^2) Delta R_{ij}^2 / R^2
97  ///       diB = 1/kti^2
98  antikt_algorithm=2, 
99  /// like the k_t but with distance measures
100  ///       dij = min(kti^{2p},ktj^{2p}) Delta R_{ij}^2 / R^2
101  ///       diB = 1/kti^{2p}
102  /// where p = extra_param()
103  genkt_algorithm=3, 
104  /// a version of cambridge with a special distance measure for particles
105  /// whose pt is < extra_param()
106  cambridge_for_passive_algorithm=11,
107  /// a version of genkt with a special distance measure for particles
108  /// whose pt is < extra_param() [relevant for passive areas when p<=0]
109  genkt_for_passive_algorithm=13, 
110  //.................................................................
111  /// the e+e- kt algorithm
112  ee_kt_algorithm=50,
113  /// the e+e- genkt algorithm  (R > 2 and p=1 gives ee_kt)
114  ee_genkt_algorithm=53,
115  //.................................................................
116  /// any plugin algorithm supplied by the user
117  plugin_algorithm = 99,
118  //.................................................................
119  /// the value for the jet algorithm in a JetDefinition for which
120  /// no algorithm has yet been defined
121  undefined_jet_algorithm = 999
122};
123
124/// make standard Les Houches nomenclature JetAlgorithm (algorithm is general
125/// recipe without the parameters) backward-compatible with old JetFinder
126typedef JetAlgorithm JetFinder;
127
128/// provide other possible names for the Cambridge/Aachen algorithm
129const JetAlgorithm aachen_algorithm = cambridge_algorithm;
130const JetAlgorithm cambridge_aachen_algorithm = cambridge_algorithm;
131
132//======================================================================
133/// the various recombination schemes
134enum RecombinationScheme {
135  /// summing the 4-momenta
136  E_scheme=0,
137  /// pt weighted recombination of y,phi (and summing of pt's)
138  /// with preprocessing to make things massless by rescaling E=|\vec p|
139  pt_scheme=1,
140  /// pt^2 weighted recombination of y,phi (and summing of pt's)
141  /// with preprocessing to make things massless by rescaling E=|\vec p|
142  pt2_scheme=2,
143  /// pt weighted recombination of y,phi (and summing of pt's)
144  /// with preprocessing to make things massless by rescaling |\vec p|->=E
145  Et_scheme=3,
146  /// pt^2 weighted recombination of y,phi (and summing of pt's)
147  /// with preprocessing to make things massless by rescaling |\vec p|->=E
148  Et2_scheme=4,
149  /// pt weighted recombination of y,phi (and summing of pt's), with
150  /// no preprocessing
151  BIpt_scheme=5,
152  /// pt^2 weighted recombination of y,phi (and summing of pt's)
153  /// no preprocessing
154  BIpt2_scheme=6,
155  /// for the user's external scheme
156  external_scheme = 99
157};
158
159
160
161// forward declaration, needed in order to specify interface for the
162// plugin.
163class ClusterSequence;
164
165
166
167
168//======================================================================
169/// @ingroup basic_classes
170/// \class JetDefinition
171/// class that is intended to hold a full definition of the jet
172/// clusterer
173class JetDefinition {
174 
175public:
176
177  /// forward declaration of a class that allows the user to introduce
178  /// their own plugin
179  class Plugin;
180
181  // forward declaration of a class that will provide the
182  // recombination scheme facilities and/or allow a user to
183  // extend these facilities
184  class Recombiner;
185
186
187  /// constructor with alternative ordering or arguments -- note that
188  /// we have not provided a default jet finder, to avoid ambiguous
189  /// JetDefinition() constructor.
190  JetDefinition(JetAlgorithm jet_algorithm_in, 
191                double R_in, 
192                RecombinationScheme recomb_scheme_in = E_scheme,
193                Strategy strategy_in = Best) {
194    *this = JetDefinition(jet_algorithm_in, R_in, strategy_in, recomb_scheme_in, 1);
195  }
196
197  /// constructor for algorithms that have no free parameters
198  /// (e.g. ee_kt_algorithm)
199  JetDefinition(JetAlgorithm jet_algorithm_in, 
200                RecombinationScheme recomb_scheme_in = E_scheme,
201                Strategy strategy_in = Best) {
202    double dummyR = 0.0;
203    *this = JetDefinition(jet_algorithm_in, dummyR, strategy_in, recomb_scheme_in, 0);
204  }
205
206  /// constructor for algorithms that require R + one extra parameter to be set
207  /// (the gen-kt series for example)
208  JetDefinition(JetAlgorithm jet_algorithm_in, 
209                double R_in, 
210                double xtra_param_in,
211                RecombinationScheme recomb_scheme_in = E_scheme,
212                Strategy strategy_in = Best) {
213    *this = JetDefinition(jet_algorithm_in, R_in, strategy_in, recomb_scheme_in, 2);
214    set_extra_param(xtra_param_in);
215  }
216
217
218  /// constructor in a form that allows the user to provide a pointer
219  /// to an external recombiner class (which must remain valid for the
220  /// life of the JetDefinition object).
221  JetDefinition(JetAlgorithm jet_algorithm_in, 
222                double R_in, 
223                const Recombiner * recombiner_in,
224                Strategy strategy_in = Best) {
225    *this = JetDefinition(jet_algorithm_in, R_in, external_scheme, strategy_in);
226    _recombiner = recombiner_in;
227  }
228
229
230  /// constructor for case with 0 parameters (ee_kt_algorithm) and
231  /// and external recombiner
232  JetDefinition(JetAlgorithm jet_algorithm_in, 
233                const Recombiner * recombiner_in,
234                Strategy strategy_in = Best) {
235    *this = JetDefinition(jet_algorithm_in, external_scheme, strategy_in);
236    _recombiner = recombiner_in;
237  }
238
239  /// constructor allowing the extra parameter to be set and a pointer to
240  /// a recombiner
241  JetDefinition(JetAlgorithm jet_algorithm_in, 
242                double R_in, 
243                double xtra_param_in,
244                const Recombiner * recombiner_in,
245                Strategy strategy_in = Best) {
246    *this = JetDefinition(jet_algorithm_in, R_in, external_scheme, strategy_in);
247    _recombiner = recombiner_in;
248    set_extra_param(xtra_param_in);
249  }
250
251  /// a default constructor which creates a jet definition that is in
252  /// a well-defined internal state, but not actually usable for jet
253  /// clustering.
254  JetDefinition()  {
255    *this = JetDefinition(undefined_jet_algorithm, 1.0);
256  }
257 
258
259  // /// a default constructor
260  // JetDefinition() {
261  //   *this = JetDefinition(kt_algorithm, 1.0);
262  // }
263
264  /// constructor based on a pointer to a user's plugin; the object
265  /// pointed to must remain valid for the whole duration of existence
266  /// of the JetDefinition and any related ClusterSequences
267  JetDefinition(const Plugin * plugin_in) {
268    _plugin = plugin_in;
269    _strategy = plugin_strategy;
270    _Rparam = _plugin->R();
271    _jet_algorithm = plugin_algorithm;
272    set_recombination_scheme(E_scheme);
273  }
274
275
276  /// constructor to fully specify a jet-definition (together with
277  /// information about how algorithically to run it).
278  ///
279  /// the ordering of arguments here is old and deprecated (except
280  /// as the common constructor for internal use)
281  JetDefinition(JetAlgorithm jet_algorithm_in, 
282                double R_in, 
283                Strategy strategy_in,
284                RecombinationScheme recomb_scheme_in = E_scheme,
285                int nparameters_in = 1);
286 
287  /// R values larger than max_allowable_R are not allowed.
288  ///
289  /// We use a value of 1000, substantially smaller than
290  /// numeric_limits<double>::max(), to leave room for the convention
291  /// within PseudoJet of setting unphysical (infinite) rapidities to
292  /// +-(MaxRap + abs(pz())), where MaxRap is 10^5.
293  static const double max_allowable_R; //= 1000.0;
294
295  /// set the recombination scheme to the one provided
296  void set_recombination_scheme(RecombinationScheme);
297
298  /// set the recombiner class to the one provided
299  void set_recombiner(const Recombiner * recomb) {
300    if (_recombiner_shared()) _recombiner_shared.reset(recomb);
301    _recombiner = recomb;
302    _default_recombiner = DefaultRecombiner(external_scheme);
303  }
304
305  /// calling this tells the JetDefinition to handle the deletion of
306  /// the recombiner when it is no longer used
307  void delete_recombiner_when_unused();
308
309  /// return a pointer to the plugin
310  const Plugin * plugin() const {return _plugin;};
311
312  /// allows to let the JetDefinition handle the deletion of the
313  /// plugin when it is no longer used
314  void delete_plugin_when_unused();
315
316  /// return information about the definition...
317  JetAlgorithm jet_algorithm  () const {return _jet_algorithm  ;}
318  /// same as above for backward compatibility
319  JetAlgorithm jet_finder     () const {return _jet_algorithm  ;}
320  double    R           () const {return _Rparam      ;}
321  // a general purpose extra parameter, whose meaning depends on
322  // the algorithm, and may often be unused.
323  double    extra_param () const {return _extra_param ;}
324  Strategy  strategy    () const {return _strategy    ;}
325  RecombinationScheme recombination_scheme() const {
326    return _default_recombiner.scheme();}
327
328  /// (re)set the jet finder
329  void set_jet_algorithm(JetAlgorithm njf) {_jet_algorithm = njf;}
330  /// same as above for backward compatibility
331  void set_jet_finder(JetAlgorithm njf)    {_jet_algorithm = njf;}
332  /// (re)set the general purpose extra parameter
333  void set_extra_param(double xtra_param) {_extra_param = xtra_param;}
334
335  /// return a pointer to the currently defined recombiner.
336  ///
337  /// Warning: the pointer may be to an internal recombiner (for
338  /// default recombination schemes), in which case if the
339  /// JetDefinition becomes invalid (e.g. is deleted), the pointer
340  /// will then point to an object that no longer exists.
341  ///
342  /// Note also that if you copy a JetDefinition with a default
343  /// recombination scheme, then the two copies will have distinct
344  /// recombiners, and return different recombiner() pointers.
345  const Recombiner * recombiner() const {
346    return _recombiner == 0 ? & _default_recombiner : _recombiner;}
347
348  /// returns true if the current jet definitions shares the same
349  /// recombiner as teh one passed as an argument
350  bool has_same_recombiner(const JetDefinition &other_jd) const;
351
352  /// return a textual description of the current jet definition
353  std::string description() const;
354
355
356public:
357  //======================================================================
358  /// @ingroup advanced_usage
359  /// \class Recombiner
360  /// An abstract base class that will provide the recombination scheme
361  /// facilities and/or allow a user to extend these facilities
362  class Recombiner {
363  public:
364    /// return a textual description of the recombination scheme
365    /// implemented here
366    virtual std::string description() const = 0;
367   
368    /// recombine pa and pb and put result into pab
369    virtual void recombine(const PseudoJet & pa, const PseudoJet & pb, 
370                           PseudoJet & pab) const = 0;
371
372    /// routine called to preprocess each input jet (to make all input
373    /// jets compatible with the scheme requirements (e.g. massless).
374    virtual void preprocess(PseudoJet & ) const {};
375   
376    /// a destructor to be replaced if necessary in derived classes...
377    virtual ~Recombiner() {};
378
379    /// pa += pb in the given recombination scheme. Not virtual -- the
380    /// user should have no reason to want to redefine this!
381    inline void plus_equal(PseudoJet & pa, const PseudoJet & pb) const {
382      // put result in a temporary location in case the recombiner
383      // does something funny (ours doesn't, but who knows about the
384      // user's)
385      PseudoJet pres; 
386      recombine(pa,pb,pres);
387      pa = pres;
388    }
389
390  };
391 
392 
393  //======================================================================
394  /// @ingroup advanced_usage
395  /// \class DefaultRecombiner
396  /// A class that will provide the recombination scheme facilities and/or
397  /// allow a user to extend these facilities
398  ///
399  /// This class is derived from the (abstract) class Recombiner. It
400  /// simply "sums" PseudoJets using a specified recombination scheme
401  /// (E-scheme by default)
402  class DefaultRecombiner : public Recombiner {
403  public:
404    DefaultRecombiner(RecombinationScheme recomb_scheme = E_scheme) : 
405      _recomb_scheme(recomb_scheme) {}
406   
407    virtual std::string description() const;
408   
409    /// recombine pa and pb and put result into pab
410    virtual void recombine(const PseudoJet & pa, const PseudoJet & pb, 
411                           PseudoJet & pab) const;
412
413    virtual void preprocess(PseudoJet & p) const;
414
415    /// return the index of the recombination scheme
416    RecombinationScheme scheme() const {return _recomb_scheme;}
417   
418  private:
419    RecombinationScheme _recomb_scheme;
420  };
421
422
423  //======================================================================
424  /// @ingroup advanced_usage
425  /// \class Plugin
426  /// a class that allows a user to introduce their own "plugin" jet
427  /// finder
428  ///
429  /// Note that all the plugins provided with FastJet are derived from
430  /// this class
431  class Plugin{
432  public:
433    /// return a textual description of the jet-definition implemented
434    /// in this plugin
435    virtual std::string description() const = 0;
436   
437    /// given a ClusterSequence that has been filled up with initial
438    /// particles, the following function should fill up the rest of the
439    /// ClusterSequence, using the following member functions of
440    /// ClusterSequence:
441    ///   - plugin_do_ij_recombination(...)
442    ///   - plugin_do_iB_recombination(...)
443    virtual void run_clustering(ClusterSequence &) const = 0;
444   
445    virtual double R() const = 0;
446   
447    /// return true if there is specific support for the measurement
448    /// of passive areas, in the sense that areas determined from all
449    /// particles below the ghost separation scale will be a passive
450    /// area. [If you don't understand this, ignore it!]
451    virtual bool supports_ghosted_passive_areas() const {return false;}
452
453    /// set the ghost separation scale for passive area determinations
454    /// in future runs (strictly speaking that makes the routine
455    /// a non const, so related internal info must be stored as a mutable)
456    virtual void set_ghost_separation_scale(double scale) const;
457    virtual double ghost_separation_scale() const {return 0.0;}
458
459    /// if this returns false then a warning will be given
460    /// whenever the user requests "exclusive" jets from the
461    /// cluster sequence
462    virtual bool exclusive_sequence_meaningful() const {return false;}
463
464    /// a destructor to be replaced if necessary in derived classes...
465    virtual ~Plugin() {};
466  };
467
468private:
469
470
471  JetAlgorithm _jet_algorithm;
472  double    _Rparam;
473  double    _extra_param ; ///< parameter whose meaning varies according to context
474  Strategy  _strategy  ;
475
476  const Plugin * _plugin;
477  SharedPtr<const Plugin> _plugin_shared;
478
479  // when we use our own recombiner it's useful to point to it here
480  // so that we don't have to worry about deleting it etc...
481  DefaultRecombiner _default_recombiner;
482  const Recombiner * _recombiner;
483  SharedPtr<const Recombiner> _recombiner_shared;
484
485};
486
487
488//-------------------------------------------------------------------------------
489// helper functions to build a jet made of pieces
490//
491// These functions include an options recombiner used to compute the
492// total composite jet momentum
493// -------------------------------------------------------------------------------
494
495/// build a "CompositeJet" from the vector of its pieces
496///
497/// In this case, E-scheme recombination is assumed to compute the
498/// total momentum
499PseudoJet join(const std::vector<PseudoJet> & pieces, const JetDefinition::Recombiner & recombiner);
500
501/// build a MergedJet from a single PseudoJet
502PseudoJet join(const PseudoJet & j1, 
503               const JetDefinition::Recombiner & recombiner);
504
505/// build a MergedJet from 2 PseudoJet
506PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, 
507               const JetDefinition::Recombiner & recombiner);
508
509/// build a MergedJet from 3 PseudoJet
510PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, const PseudoJet & j3, 
511               const JetDefinition::Recombiner & recombiner);
512
513/// build a MergedJet from 4 PseudoJet
514PseudoJet join(const PseudoJet & j1, const PseudoJet & j2, const PseudoJet & j3, const PseudoJet & j4, 
515               const JetDefinition::Recombiner & recombiner);
516
517
518
519
520
521FASTJET_END_NAMESPACE
522
523#endif // __FASTJET_JETDEFINITION_HH__
Note: See TracBrowser for help on using the repository browser.