FastJet  3.4.0
JetMedianBackgroundEstimator.hh
1 #ifndef __FASTJET_BACKGROUND_ESTIMATOR_HH__
2 #define __FASTJET_BACKGROUND_ESTIMATOR_HH__
3 
4 
5 //FJSTARTHEADER
6 // $Id$
7 //
8 // Copyright (c) 2005-2021, Matteo Cacciari, Gavin P. Salam and Gregory Soyez
9 //
10 //----------------------------------------------------------------------
11 // This file is part of FastJet.
12 //
13 // FastJet is free software; you can redistribute it and/or modify
14 // it under the terms of the GNU General Public License as published by
15 // the Free Software Foundation; either version 2 of the License, or
16 // (at your option) any later version.
17 //
18 // The algorithms that underlie FastJet have required considerable
19 // development. They are described in the original FastJet paper,
20 // hep-ph/0512210 and in the manual, arXiv:1111.6097. If you use
21 // FastJet as part of work towards a scientific publication, please
22 // quote the version you use and include a citation to the manual and
23 // optionally also to hep-ph/0512210.
24 //
25 // FastJet is distributed in the hope that it will be useful,
26 // but WITHOUT ANY WARRANTY; without even the implied warranty of
27 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
28 // GNU General Public License for more details.
29 //
30 // You should have received a copy of the GNU General Public License
31 // along with FastJet. If not, see <http://www.gnu.org/licenses/>.
32 //----------------------------------------------------------------------
33 //FJENDHEADER
34 
35 #include "fastjet/config.h"
36 #include "fastjet/ClusterSequenceAreaBase.hh"
37 #include "fastjet/AreaDefinition.hh"
38 #include "fastjet/FunctionOfPseudoJet.hh"
39 #include "fastjet/Selector.hh"
40 #include "fastjet/tools/BackgroundEstimatorBase.hh"
41 #include <iostream>
42 
43 #ifdef FASTJET_HAVE_THREAD_SAFETY
44 #include <atomic>
45 #endif
46 
47 FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh
48 
49 
50 /// @ingroup tools_background
51 /// \class JetMedianBackgroundEstimator
52 ///
53 /// Class to estimate the pt density of the background per unit area,
54 /// using the median of the distribution of pt/area from jets that
55 /// pass some selection criterion.
56 ///
57 /// Events are passed either in the form of the event particles (in
58 /// which they're clustered by the class), a ClusterSequenceArea (in
59 /// which case the jets used are those returned by "inclusive_jets()")
60 /// or directly as a set of jets.
61 ///
62 /// The selection criterion is typically a geometrical one (e.g. all
63 /// jets with |y|<2) sometimes supplemented with some kinematical
64 /// restriction (e.g. exclusion of the two hardest jets). It is passed
65 /// to the class through a Selector.
66 ///
67 /// Beware:
68 /// by default, to correctly handle partially empty events, the
69 /// class attempts to calculate an "empty area", based
70 /// (schematically) on
71 ///
72 /// range.total_area() - sum_{jets_in_range} jets.area()
73 ///
74 /// For ranges with small areas, this can be inaccurate (particularly
75 /// relevant in dense events where empty_area should be zero and ends
76 /// up not being zero).
77 ///
78 /// This calculation of empty area can be avoided if a
79 /// ClusterSequenceArea class with explicit ghosts
80 /// (ActiveAreaExplicitGhosts) is used. This is _recommended_
81 /// unless speed requirements cause you to use Voronoi areas. For
82 /// speedy background estimation you could also consider using
83 /// GridMedianBackgroundEstimator.
84 ///
85 ///
87 public:
88  /// @name constructors and destructors
89  //\{
90  //----------------------------------------------------------------
91  /// Constructor that sets the rho range as well as the jet
92  /// definition and area definition to be used to cluster the
93  /// particles. Prior to the estimation of rho, one has to provide
94  /// the particles to cluster using set_particles(...)
95  ///
96  /// \param rho_range the Selector specifying which jets will be considered
97  /// \param jet_def the jet definition to use for the clustering
98  /// \param area_def the area definition to use for the clustering
99  JetMedianBackgroundEstimator(const Selector &rho_range,
100  const JetDefinition &jet_def,
101  const AreaDefinition &area_def);
102 
103  /// ctor from a ClusterSequenceAreaBase with area
104  ///
105  /// \param rho_range the Selector specifying which jets will be considered
106  /// \param csa the ClusterSequenceArea to use
107  ///
108  /// Pre-conditions:
109  /// - one should be able to estimate the "empty area" (i.e. the area
110  /// not occupied by jets). This is feasible if at least one of the following
111  /// conditions is satisfied:
112  /// ( i) the ClusterSequence has explicit ghosts
113  /// (ii) the range has a computable area.
114  /// - the jet algorithm must be suited for median computation
115  /// (otherwise a warning will be issues)
116  ///
117  /// Note that selectors with e.g. hardest-jets exclusion do not have
118  /// a well-defined area. For this reasons, it is STRONGLY advised to
119  /// use an area with explicit ghosts.
120  JetMedianBackgroundEstimator(const Selector &rho_range,
121  const ClusterSequenceAreaBase &csa);
122 
123 
124  /// Default constructor that optionally sets the rho range. The
125  /// configuration must be done later calling
126  /// set_cluster_sequence(...) or set_jets(...).
127  ///
128  /// \param rho_range the Selector specifying which jets will be considered
129  ///
130  JetMedianBackgroundEstimator(const Selector &rho_range = SelectorIdentity())
131  : _rho_range(rho_range), _jet_def(JetDefinition()),
132  _enable_rho_m(true){ reset(); }
133 
134 
135  /// default dtor
137 
138  //\}
139 
140 
141  /// @name setting a new event
142  //\{
143  //----------------------------------------------------------------
144 
145  /// tell the background estimator that it has a new event, composed
146  /// of the specified particles.
147  virtual void set_particles(const std::vector<PseudoJet> & particles) FASTJET_OVERRIDE;
148 
149  // tell the background estimator that it has a new event, composed
150  // of the specified particles and use the supplied seed for the
151  // generation of ghosts. If the seed is empty, it is ignored.
152  virtual void set_particles_with_seed(const std::vector<PseudoJet> & particles, const std::vector<int> & seed) FASTJET_OVERRIDE;
153 
154  /// (re)set the cluster sequence (with area support) to be used by
155  /// future calls to rho() etc.
156  ///
157  /// \param csa the cluster sequence area
158  ///
159  /// Pre-conditions:
160  /// - one should be able to estimate the "empty area" (i.e. the area
161  /// not occupied by jets). This is feasible if at least one of the following
162  /// conditions is satisfied:
163  /// ( i) the ClusterSequence has explicit ghosts
164  /// (ii) the range selected has a computable area.
165  /// - the jet algorithm must be suited for median computation
166  /// (otherwise a warning will be issues)
167  ///
168  /// Note that selectors with e.g. hardest-jets exclusion do not have
169  /// a well-defined area. For this reasons, it is STRONGLY advised to
170  /// use an area with explicit ghosts.
171  void set_cluster_sequence(const ClusterSequenceAreaBase & csa);
172 
173  /// (re)set the jets (which must have area support) to be used by future
174  /// calls to rho() etc.; for the conditions that must be satisfied
175  /// by the jets, see the Constructor that takes jets.
176  void set_jets(const std::vector<PseudoJet> &jets);
177 
178  /// (re)set the selector to be used for future calls to rho() etc.
179  void set_selector(const Selector & rho_range_selector) {
180  _rho_range = rho_range_selector;
181  _set_cache_unavailable();
182  }
183 
184  /// determine whether the automatic calculation of rho_m and sigma_m
185  /// is enabled (by default true)
186  void set_compute_rho_m(bool enable){
187  _enable_rho_m = enable;
188  _set_cache_unavailable();
189  }
190 
191  //\}
192 
193 
194  /// return a pointer to a copy of this BGE; the user is responsible
195  /// for eventually deleting the resulting object.
196  BackgroundEstimatorBase * copy() const FASTJET_OVERRIDE {
197  return new JetMedianBackgroundEstimator(*this);
198  };
199 
200 
201  /// @name retrieving fundamental information
202  //\{
203  //----------------------------------------------------------------
204  /// get the full set of background properties
205  ///
206  /// For background estimators using a local ranges, this throws an
207  /// error (use estimate(jet) instead)
208  /// In the presence of a rescaling, the rescaling is not included
209  BackgroundEstimate estimate() const FASTJET_OVERRIDE;
210 
211  /// get the full set of background properties for a given reference jet
212  BackgroundEstimate estimate(const PseudoJet &jet) const FASTJET_OVERRIDE;
213 
214  /// get rho, the median background density per unit area
215  double rho() const FASTJET_OVERRIDE;
216 
217  /// get sigma, the background fluctuations per unit area
218  double sigma() const FASTJET_OVERRIDE;
219 
220  /// get rho, the median background density per unit area, locally at
221  /// the position of a given jet.
222  ///
223  /// If the Selector associated with the range takes a reference jet
224  /// (i.e. is relocatable), then for subsequent operations the
225  /// Selector has that jet set as its reference.
226  double rho(const PseudoJet & jet) FASTJET_OVERRIDE;
227 
228  /// get sigma, the background fluctuations per unit area,
229  /// locally at the position of a given jet.
230  ///
231  /// If the Selector associated with the range takes a reference jet
232  /// (i.e. is relocatable), then for subsequent operations the
233  /// Selector has that jet set as its reference.
234  double sigma(const PseudoJet &jet) FASTJET_OVERRIDE;
235 
236  /// returns true if this background estimator has support for
237  /// determination of sigma
238  virtual bool has_sigma() const FASTJET_OVERRIDE {return true;}
239 
240  //----------------------------------------------------------------
241  // now do the same thing for rho_m and sigma_m
242 
243  /// returns rho_m, the purely longitudinal, particle-mass-induced
244  /// component of the background density per unit area
245  virtual double rho_m() const FASTJET_OVERRIDE;
246 
247  /// returns sigma_m, a measure of the fluctuations in the purely
248  /// longitudinal, particle-mass-induced component of the background
249  /// density per unit area; must be multipled by sqrt(area) to get
250  /// fluctuations for a region of a given area.
251  virtual double sigma_m() const FASTJET_OVERRIDE;
252 
253  /// Returns rho_m locally at the jet position. As for rho(jet), it is non-const.
254  virtual double rho_m(const PseudoJet & /*jet*/) FASTJET_OVERRIDE;
255 
256  /// Returns sigma_m locally at the jet position. As for rho(jet), it is non-const.
257  virtual double sigma_m(const PseudoJet & /*jet*/) FASTJET_OVERRIDE;
258 
259  /// Returns true if this background estimator has support for
260  /// determination of rho_m.
261  ///
262  /// In te presence of a density class, support for rho_m is
263  /// automatically disabled
264  ///
265  /// Note that support for sigma_m is automatic is one has sigma and
266  /// rho_m support.
267  virtual bool has_rho_m() const FASTJET_OVERRIDE {return _enable_rho_m && (_jet_density_class == 0);}
268  //\}
269 
270  /// @name retrieving additional useful information
271  //\{
272  //----------------------------------------------------------------
273  /// Returns the mean area of the jets used to actually compute the
274  /// background properties in the last call of rho() or sigma()
275  /// If the configuration has changed in the meantime, throw an error.
276  double mean_area() const;
277 
278  /// returns the number of jets used to actually compute the
279  /// background properties in the last call of rho() or sigma()
280  /// If the configuration has changed in the meantime, throw an error.
281  unsigned int n_jets_used() const;
282 
283  /// returns the jets used to actually compute the background
284  /// properties
285  std::vector<PseudoJet> jets_used() const;
286 
287  /// Returns the estimate of the area (within the range defined by
288  /// the selector) that is not occupied by jets. The value is that
289  /// for the last call of rho() or sigma()
290  /// If the configuration has changed in the meantime, throw an error.
291  ///
292  /// The answer is defined to be zero if the area calculation
293  /// involved explicit ghosts; if the area calculation was an active
294  /// area, then use is made of the active area's internal list of
295  /// pure ghost jets (taking those that pass the selector); otherwise
296  /// it is based on the difference between the selector's total area
297  /// and the area of the jets that pass the selector.
298  ///
299  /// The result here is just the cached result of the corresponding
300  /// call to the ClusterSequenceAreaBase function.
301  double empty_area() const;
302 
303  /// Returns the number of empty jets used when computing the
304  /// background properties. The value is that for the last call of
305  /// rho() or sigma().
306  /// If the configuration has changed in the meantime, throw an error.
307  ///
308  /// If the area has explicit ghosts the result is zero; for active
309  /// areas it is the number of internal pure ghost jets that pass the
310  /// selector; otherwise it is deduced from the empty area, divided by
311  /// \f$ 0.55 \pi R^2 \f$ (the average pure-ghost-jet area).
312  ///
313  /// The result here is just the cached result of the corresponding
314  /// call to the ClusterSequenceAreaBase function.
315  double n_empty_jets() const;
316 
317  //}
318 
319 
320  /// @name configuring behaviour
321  //\{
322  //----------------------------------------------------------------
323 
324  /// Resets the class to its default state, including the choice to
325  /// use 4-vector areas.
326  ///
327  void reset();
328 
329  /// By default when calculating pt/Area for a jet, it is the
330  /// transverse component of the 4-vector area that is used in the ratiof \f$p_t/A\f$.
331  /// Calling this function with a "false" argument causes the scalar area to
332  /// be used instead.
333  ///
334  /// While the difference between the two choices is usually small,
335  /// for high-precision work it is usually the 4-vector area that is
336  /// to be preferred.
337  ///
338  /// \param use_it whether one uses the 4-vector area or not (true by default)
339  void set_use_area_4vector(bool use_it = true){
340  _use_area_4vector = use_it;
341  _set_cache_unavailable();
342  }
343 
344  /// check if the estimator uses the 4-vector area or the scalar area
345  bool use_area_4vector() const{ return _use_area_4vector;}
346 
347  /// The FastJet v2.X sigma calculation had a small spurious offset
348  /// in the limit of a small number of jets. This is fixed by default
349  /// in versions 3 upwards. The old behaviour can be obtained with a
350  /// call to this function.
351  void set_provide_fj2_sigma(bool provide_fj2_sigma = true) {
352  _provide_fj2_sigma = provide_fj2_sigma;
353  _set_cache_unavailable();
354  }
355 
356  /// Set a pointer to a class that calculates the quantity whose
357  /// median will be calculated; if the pointer is null then pt/area
358  /// is used (as occurs also if this function is not called).
359  ///
360  /// Note that this is still <i>preliminary</i> in FastJet 3.0 and
361  /// that backward compatibility is not guaranteed in future releases
362  /// of FastJet
363  void set_jet_density_class(const FunctionOfPseudoJet<double> * jet_density_class);
364 
365  /// return the pointer to the jet density class
367  return _jet_density_class;
368  }
369 
370  /// Set a pointer to a class that calculates the rescaling factor as
371  /// a function of the jet (position). Note that the rescaling factor
372  /// is used both in the determination of the "global" rho (the pt/A
373  /// of each jet is divided by this factor) and when asking for a
374  /// local rho (the result is multiplied by this factor).
375  ///
376  /// The BackgroundRescalingYPolynomial class can be used to get a
377  /// rescaling that depends just on rapidity.
378  virtual void set_rescaling_class(const FunctionOfPseudoJet<double> * rescaling_class_in) FASTJET_OVERRIDE {
379  BackgroundEstimatorBase::set_rescaling_class(rescaling_class_in);
380  _set_cache_unavailable();
381  }
382 
383  //\}
384 
385  /// @name description
386  //\{
387  //----------------------------------------------------------------
388 
389  /// returns a textual description of the background estimator
390  std::string description() const FASTJET_OVERRIDE;
391 
392  //\}
393 
394  /// an internal class to hold results of the calculation
395  /// that are to be assigned to the "extras" part of a BackgroundEstimate
396  class Extras : public BackgroundEstimate::Extras {
397  public:
398  Extras()
399  : _reference_jet(PseudoJet()), _n_jets_used(0),
400  _n_empty_jets(0.0), _empty_area(0.0) {}
401 
402  /// returns the current reference jet
403  PseudoJet reference_jet() const {return _reference_jet;}
404 
405  /// returns the number of jets used to estimate the background
406  unsigned int n_jets_used() const {return _n_jets_used;}
407 
408  /// returns the number of empty (pure-ghost) jets
409  double n_empty_jets() const {return _n_empty_jets;}
410 
411  /// returns the empty (pure-ghost/unclustered) area!
412  double empty_area() const {return _empty_area;}
413 
414  void set_reference_jet(const PseudoJet &reference_jet_in){
415  _reference_jet = reference_jet_in;
416  }
417  void set_n_jets_used(int n_jets_used_in){ _n_jets_used=n_jets_used_in;}
418  void set_n_empty_jets(double n_empty_jets_in){ _n_empty_jets=n_empty_jets_in;}
419  void set_empty_area(double empty_area_in){ _empty_area=empty_area_in;}
420 
421 
422  protected:
423  PseudoJet _reference_jet; ///< current reference jet
424  unsigned int _n_jets_used; ///< number of jets used to estimate the background
425  double _n_empty_jets; ///< number of empty (pure-ghost) jets
426  double _empty_area; ///< the empty (pure-ghost/unclustered) area!
427  };
428 
429 
430 private:
431 
432  /// compute the background properties for a given jet (excluding
433  /// rescaling factors) and return a corresponding BackgroundEstimate
434  ///
435  /// this leaves the cache (and the status flags) unchanged
436  BackgroundEstimate _compute(const PseudoJet &jet) const;
437 
438  //------
439  // the next calls are meant for the case where the cache can be
440  // filled once and for all, i.e. cases where the selector does NOT
441  // take a reference
442 
443  /// fill the cache with the given estimate
444  void _cache_no_overwrite(const BackgroundEstimate &estimate) const;
445 
446  /// fill the cache with a computed estimate
447  void _compute_and_cache_no_overwrite() const;
448 
449  //------
450  // the next calls are meant for the case where the selector does
451  // take a reference and the cache needs to be refilled whenever one
452  // calls this background estimate with a different reference jet
453 
454  /// fill the cache with the given estimate
455  void _cache(const BackgroundEstimate &estimate) const;
456 
457  /// update the cache if need be and return the background
458  /// estimate. This is meant to be called for cases with a local
459  /// range (selector that takesa reference)
460  BackgroundEstimate _compute_and_cache_if_needed(const PseudoJet &jet) const;
461  //-----
462 
463  /// check that the underlying structure is still alive
464  /// throw an error otherwise
465  void _check_csa_alive() const;
466 
467  /// check that the algorithm used for the clustering is adapted for
468  /// background estimation (i.e. either kt or C/A)
469  /// Issue a warning otherwise
470  void _check_jet_alg_good_for_median() const;
471 
472  // the basic parameters of this class (passed through the variou ctors)
473  Selector _rho_range; ///< range to compute the background in
474  JetDefinition _jet_def; ///< the jet def to use for teh clustering
475  AreaDefinition _area_def; ///< the area def to use for teh clustering
476  std::vector<PseudoJet> _included_jets; ///< jets to be used
477 
478  // the tunable parameters of the class
479  bool _use_area_4vector;
480  bool _provide_fj2_sigma;
481  const FunctionOfPseudoJet<double> * _jet_density_class;
482  //SharedPtr<BackgroundRescalingBase> _rescaling_class_sharedptr;
483  bool _enable_rho_m;
484 
485  // internal variables
486  SharedPtr<PseudoJetStructureBase> _csi; ///< allows to check if _csa is still valid
487 
488  /// handle warning messages
489  static LimitedWarning _warnings;
490  static LimitedWarning _warnings_zero_area;
491  static LimitedWarning _warnings_preliminary;
492 };
493 
494 
495 
496 
497 //----------------------------------------------------------------------
498 /// @ingroup tools_background
499 /// \class BackgroundJetPtDensity
500 /// Class that implements pt/area_4vector.perp() for background estimation
501 /// <i>(this is a preliminary class)</i>.
503 public:
504  virtual double result(const PseudoJet & jet) const {
505  return jet.perp() / jet.area_4vector().perp();
506  }
507  virtual std::string description() const {return "BackgroundJetPtDensity";}
508 };
509 
510 
511 //----------------------------------------------------------------------
512 /// @ingroup tools_background
513 /// \class BackgroundJetScalarPtDensity
514 /// Class that implements (scalar pt sum of jet)/(scalar area of jet)
515 /// for background estimation <i>(this is a preliminary class)</i>.
516 ///
517 /// Optionally it can return a quantity based on the sum of pt^n,
518 /// e.g. for use in subtracting fragementation function moments.
520 public:
521  /// Default constructor provides background estimation with scalar pt sum
522  BackgroundJetScalarPtDensity() : _pt_power(1) {}
523 
524  /// Constructor to provide background estimation based on
525  /// \f$ sum_{i\in jet} p_{ti}^{n} \f$
526  BackgroundJetScalarPtDensity(double n) : _pt_power(n) {}
527 
528  virtual double result(const PseudoJet & jet) const;
529 
530  virtual std::string description() const;
531 
532 private:
533  double _pt_power;
534 };
535 
536 //----------------------------------------------------------------------
537 /// @ingroup tools_background
538 /// \class BackgroundJetPtMDensity
539 /// Class that implements
540 /// \f$ \frac{1}{A} \sum_{i \in jet} (\sqrt{p_{ti}^2+m^2} - p_{ti}) \f$
541 /// for background estimation <i>(this is a preliminary class)</i>.
542 ///
543 ///
544 /// This is useful for correcting jet masses in cases where the event
545 /// involves massive particles.
547 public:
548  virtual double result(const PseudoJet & jet) const {
549  std::vector<PseudoJet> constituents = jet.constituents();
550  double scalar_ptm = 0;
551  for (unsigned i = 0; i < constituents.size(); i++) {
552  scalar_ptm += constituents[i].mperp() - constituents[i].perp();
553  }
554  return scalar_ptm / jet.area();
555  }
556 
557  virtual std::string description() const {return "BackgroundPtMDensity";}
558 };
559 
560 
561 FASTJET_END_NAMESPACE
562 
563 #endif // __BACKGROUND_ESTIMATOR_HH__
564 
class that holds a generic area definition
/// a class that holds the result of the calculation
Abstract base class that provides the basic interface for classes that estimate levels of background ...
Class that implements pt/area_4vector.perp() for background estimation (this is a preliminary class).
virtual std::string description() const
returns a description of the function (an empty string by default)
virtual double result(const PseudoJet &jet) const
the action of the function this has to be overloaded in derived classes
Class that implements for background estimation (this is a preliminary class).
virtual std::string description() const
returns a description of the function (an empty string by default)
virtual double result(const PseudoJet &jet) const
the action of the function this has to be overloaded in derived classes
Class that implements (scalar pt sum of jet)/(scalar area of jet) for background estimation (this is ...
BackgroundJetScalarPtDensity(double n)
Constructor to provide background estimation based on .
BackgroundJetScalarPtDensity()
Default constructor provides background estimation with scalar pt sum.
base class that sets interface for extensions of ClusterSequence that provide information about the a...
class that is intended to hold a full definition of the jet clusterer
an internal class to hold results of the calculation that are to be assigned to the "extras" part of ...
PseudoJet reference_jet() const
returns the current reference jet
unsigned int _n_jets_used
number of jets used to estimate the background
unsigned int n_jets_used() const
returns the number of jets used to estimate the background
double _n_empty_jets
number of empty (pure-ghost) jets
double n_empty_jets() const
returns the number of empty (pure-ghost) jets
double empty_area() const
returns the empty (pure-ghost/unclustered) area!
double _empty_area
the empty (pure-ghost/unclustered) area!
Class to estimate the pt density of the background per unit area, using the median of the distributio...
void set_provide_fj2_sigma(bool provide_fj2_sigma=true)
The FastJet v2.X sigma calculation had a small spurious offset in the limit of a small number of jets...
BackgroundEstimatorBase * copy() const override
return a pointer to a copy of this BGE; the user is responsible for eventually deleting the resulting...
const FunctionOfPseudoJet< double > * jet_density_class() const
return the pointer to the jet density class
bool use_area_4vector() const
check if the estimator uses the 4-vector area or the scalar area
void set_use_area_4vector(bool use_it=true)
By default when calculating pt/Area for a jet, it is the transverse component of the 4-vector area th...
JetMedianBackgroundEstimator(const Selector &rho_range=SelectorIdentity())
Default constructor that optionally sets the rho range.
void set_compute_rho_m(bool enable)
determine whether the automatic calculation of rho_m and sigma_m is enabled (by default true)
virtual void set_rescaling_class(const FunctionOfPseudoJet< double > *rescaling_class_in) override
Set a pointer to a class that calculates the rescaling factor as a function of the jet (position).
void set_selector(const Selector &rho_range_selector)
(re)set the selector to be used for future calls to rho() etc.
class to provide facilities for giving warnings up to some maximum number of times and to provide glo...
Class to contain pseudojets, including minimal information of use to jet-clustering routines.
Definition: PseudoJet.hh:68
virtual std::vector< PseudoJet > constituents() const
retrieve the constituents.
Definition: PseudoJet.cc:692
double perp() const
returns the scalar transverse momentum
Definition: PseudoJet.hh:158
virtual double area() const
return the jet (scalar) area.
Definition: PseudoJet.cc:829
virtual PseudoJet area_4vector() const
return the jet 4-vector area.
Definition: PseudoJet.cc:844
Class that encodes information about cuts and other selection criteria that can be applied to PseudoJ...
Definition: Selector.hh:149
An implementation of shared pointers that is broadly similar to C++11 shared_ptr (https://en....
Definition: SharedPtr.hh:341