FastJet  3.3.4
ClusterSequenceArea.hh
1 //FJSTARTHEADER
2 // $Id: ClusterSequenceArea.hh 4442 2020-05-05 07:50:11Z soyez $
3 //
4 // Copyright (c) 2006-2020, 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 __FASTJET_CLUSTERSEQUENCEAREA_HH__
32 #define __FASTJET_CLUSTERSEQUENCEAREA_HH__
33 
34 #include "fastjet/ClusterSequenceAreaBase.hh"
35 #include "fastjet/ClusterSequenceActiveArea.hh"
36 #include "fastjet/ClusterSequenceActiveAreaExplicitGhosts.hh"
37 #include "fastjet/ClusterSequencePassiveArea.hh"
38 #include "fastjet/ClusterSequenceVoronoiArea.hh"
39 #include "fastjet/AreaDefinition.hh"
40 
41 FASTJET_BEGIN_NAMESPACE
42 
43 /// @ingroup area_classes
44 /// \class ClusterSequenceArea
45 /// General class for user to obtain ClusterSequence with additional
46 /// area information.
47 ///
48 /// Based on the area_def, it automatically dispatches the work to the
49 /// appropriate actual ClusterSequenceAreaBase-derived-class to do the
50 /// real work.
52 public:
53  /// main constructor
54  template<class L> ClusterSequenceArea
55  (const std::vector<L> & pseudojets,
56  const JetDefinition & jet_def_in,
57  const AreaDefinition & area_def_in) : _area_def(area_def_in) {
58  initialize_and_run_cswa(pseudojets, jet_def_in);
59  }
60 
61  /// constructor with a GhostedAreaSpec
62  template<class L> ClusterSequenceArea
63  (const std::vector<L> & pseudojets,
64  const JetDefinition & jet_def_in,
65  const GhostedAreaSpec & ghost_spec) : _area_def(ghost_spec){
66  initialize_and_run_cswa(pseudojets, jet_def_in);
67  }
68 
69  /// constructor with a VoronoiAreaSpec
70  template<class L> ClusterSequenceArea
71  (const std::vector<L> & pseudojets,
72  const JetDefinition & jet_def_in,
73  const VoronoiAreaSpec & voronoi_spec) : _area_def(voronoi_spec){
74  initialize_and_run_cswa(pseudojets, jet_def_in);
75  }
76 
77  /// return a reference to the area definition
78  const AreaDefinition & area_def() const {return _area_def;}
79 
80 
81  /// return the area associated with the given jet
82  virtual double area (const PseudoJet & jet) const FASTJET_OVERRIDE {
83  return _area_base->area(jet);}
84 
85  /// return the error (uncertainty) associated with the determination
86  /// of the area of this jet
87  virtual double area_error (const PseudoJet & jet) const FASTJET_OVERRIDE {
88  return _area_base->area_error(jet);}
89 
90  /// return the 4-vector area
91  virtual PseudoJet area_4vector(const PseudoJet & jet) const FASTJET_OVERRIDE {
92  return _area_base->area_4vector(jet);}
93 
94  // /// return the total area, up to |y|<maxrap, that is free of jets
95  // virtual double empty_area(double maxrap) const {
96  // return _area_base->empty_area(maxrap);}
97  //
98  // /// return something similar to the number of pure ghost jets
99  // /// in the given rapidity range in an active area case.
100  // /// For the local implementation we return empty_area/(0.55 pi R^2),
101  // /// based on measured properties of ghost jets with kt and cam. Note
102  // /// that the number returned is a double.
103  // virtual double n_empty_jets(double maxrap) const {
104  // return _area_base->n_empty_jets(maxrap);
105 
106  /// return the total area, corresponding to the given selector, that
107  /// is free of jets
108  ///
109  /// The selector needs to have a finite area and be applicable jet by
110  /// jet (see the BackgroundEstimator and Subtractor tools for more
111  /// advanced usage)
112  virtual double empty_area(const Selector & selector) const FASTJET_OVERRIDE {
113  return _area_base->empty_area(selector);}
114 
115  /// return something similar to the number of pure ghost jets
116  /// in the given rap-phi range in an active area case.
117  /// For the local implementation we return empty_area/(0.55 pi R^2),
118  /// based on measured properties of ghost jets with kt and cam. Note
119  /// that the number returned is a double.
120  ///
121  /// The selector needs to have a finite area and be applicable jet by
122  /// jet (see the BackgroundEstimator and Subtractor tools for more
123  /// advanced usage)
124  virtual double n_empty_jets(const Selector & selector) const FASTJET_OVERRIDE {
125  return _area_base->n_empty_jets(selector);
126  }
127 
128  /// true if a jet is made exclusively of ghosts
129  virtual bool is_pure_ghost(const PseudoJet & jet) const FASTJET_OVERRIDE {
130  return _area_base->is_pure_ghost(jet);
131  }
132 
133  /// true if this ClusterSequence has explicit ghosts
134  virtual bool has_explicit_ghosts() const FASTJET_OVERRIDE {
135  return _area_base->has_explicit_ghosts();
136  }
137 
138 
139  /// overload version of what's in the ClusterSequenceAreaBase class, which
140  /// additionally checks compatibility between "selector" and region in which
141  /// ghosts are thrown.
142  ///
143  /// The selector needs to have a finite area and be applicable jet by
144  /// jet (see the BackgroundEstimator and Subtractor tools for more
145  /// advanced usage)
146  //FASTJET_DEPRECATED_MSG("ClusterSequenceArea::get_median_rho_and_sigma(...) is deprecated since FastJet 3.0. Use the BackgroundEstimator series of tools instead")
147  virtual void get_median_rho_and_sigma(const std::vector<PseudoJet> & all_jets,
148  const Selector & selector,
149  bool use_area_4vector,
150  double & median, double & sigma,
151  double & mean_area,
152  bool all_are_incl = false) const FASTJET_OVERRIDE {
153  _warn_if_range_unsuitable(selector);
154  ClusterSequenceAreaBase::_get_median_rho_and_sigma(
155  all_jets, selector, use_area_4vector,
156  median, sigma, mean_area, all_are_incl);
157  }
158 
159  /// overload version of what's in the ClusterSequenceAreaBase class,
160  /// which actually just does the same thing as the base version (but
161  /// since we've overridden the 5-argument version above, we have to
162  /// override the 4-argument version too.
163  //FASTJET_DEPRECATED_MSG("ClusterSequenceArea::get_median_rho_and_sigma(...) is deprecated since FastJet 3.0. Use the BackgroundEstimator series of tools instead")
164  virtual void get_median_rho_and_sigma(const Selector & selector,
165  bool use_area_4vector,
166  double & median, double & sigma) const FASTJET_OVERRIDE {
167  ClusterSequenceAreaBase::_get_median_rho_and_sigma(selector,use_area_4vector,
168  median,sigma);
169  }
170 
171  /// overload version of what's in the ClusterSequenceAreaBase class,
172  /// which actually just does the same thing as the base version (but
173  /// since we've overridden the multi-argument version above, we have to
174  /// override the 5-argument version too.
175  //FASTJET_DEPRECATED_MSG("ClusterSequenceArea::get_median_rho_and_sigma(...) is deprecated since FastJet 3.0. Use the BackgroundEstimator series of tools instead")
176  virtual void get_median_rho_and_sigma(const Selector & selector,
177  bool use_area_4vector,
178  double & median, double & sigma,
179  double & mean_area) const FASTJET_OVERRIDE {
180  ClusterSequenceAreaBase::_get_median_rho_and_sigma(selector,use_area_4vector,
181  median,sigma, mean_area);
182  }
183 
184 
185  /// overload version of what's in the ClusterSequenceAreaBase class, which
186  /// additionally checks compatibility between "range" and region in which
187  /// ghosts are thrown.
188  //FASTJET_DEPRECATED_MSG("ClusterSequenceArea::parabolic_pt_per_unit_area(...) is deprecated since FastJet 3.0. Use the BackgroundEstimator series of tools instead")
189  virtual void parabolic_pt_per_unit_area(double & a, double & b,
190  const Selector & selector,
191  double exclude_above=-1.0,
192  bool use_area_4vector=false) const FASTJET_OVERRIDE {
193  return _parabolic_pt_per_unit_area(a,b,selector,exclude_above,use_area_4vector);
194  }
195 
196 
197 private:
198 
199  /// print a warning if the range is unsuitable for the current
200  /// calculation of the area (e.g. because ghosts do not extend
201  /// far enough).
202  void _warn_if_range_unsuitable(const Selector & selector) const;
203 
204  template<class L> void initialize_and_run_cswa (
205  const std::vector<L> & pseudojets,
206  const JetDefinition & jet_def);
207 
209  AreaDefinition _area_def;
210  static LimitedWarning _range_warnings;
211  static LimitedWarning _explicit_ghosts_repeats_warnings;
212 
213  // the following set of private methods are all deprecated. Their
214  // role is simply to hide the corresponding methods (without the
215  // first underscore) from the public interface so that they can be
216  // used internally until all the deprecated methods are removed.
217  // DO NOT USE ANY OF THESE METHODS: THEY ARE DEPRECATED AND WILL BE
218  // REMOVED.
219  virtual void _parabolic_pt_per_unit_area(double & a, double & b,
220  const Selector & selector,
221  double exclude_above=-1.0,
222  bool use_area_4vector=false) const FASTJET_OVERRIDE {
223  _warn_if_range_unsuitable(selector);
224  ClusterSequenceAreaBase::_parabolic_pt_per_unit_area(
225  a,b,selector, exclude_above, use_area_4vector);
226  }
227 
228 };
229 
230 //----------------------------------------------------------------------
231 template<class L> void ClusterSequenceArea::initialize_and_run_cswa(
232  const std::vector<L> & pseudojets,
233  const JetDefinition & jet_def_in)
234  {
235 
236  ClusterSequenceAreaBase * _area_base_ptr;
237  switch(_area_def.area_type()) {
238  case active_area:
239  _area_base_ptr = new ClusterSequenceActiveArea(pseudojets,
240  jet_def_in,
241  _area_def.ghost_spec());
242  break;
243  case active_area_explicit_ghosts:
244  if (_area_def.ghost_spec().repeat() != 1)
245  _explicit_ghosts_repeats_warnings.warn("Requested active area with explicit ghosts with repeat != 1; only 1 set of ghosts will be used");
246  _area_base_ptr = new ClusterSequenceActiveAreaExplicitGhosts(pseudojets,
247  jet_def_in,
248  _area_def.ghost_spec());
249  break;
250  case voronoi_area:
251  _area_base_ptr = new ClusterSequenceVoronoiArea(pseudojets,
252  jet_def_in,
253  _area_def.voronoi_spec());
254  break;
255  case one_ghost_passive_area:
256  _area_base_ptr = new ClusterSequence1GhostPassiveArea(pseudojets,
257  jet_def_in,
258  _area_def.ghost_spec());
259  break;
260  case passive_area:
261  _area_base_ptr = new ClusterSequencePassiveArea(pseudojets,
262  jet_def_in,
263  _area_def.ghost_spec());
264  break;
265  default:
266  std::ostringstream err;
267  err << "Error: unrecognized area_type in ClusterSequenceArea:"
268  << _area_def.area_type();
269  throw Error(err.str());
270  //exit(-1);
271  }
272  // now copy across the information from the area base class
273  _area_base = SharedPtr<ClusterSequenceAreaBase>(_area_base_ptr);
274  transfer_from_sequence(*_area_base);
275 }
276 
277 FASTJET_END_NAMESPACE
278 
279 #endif // __FASTJET_CLUSTERSEQUENCEAREA_HH__
280 
281 
fastjet::ClusterSequenceArea::area_def
const AreaDefinition & area_def() const
return a reference to the area definition
Definition: ClusterSequenceArea.hh:78
fastjet::ClusterSequenceArea::get_median_rho_and_sigma
virtual void get_median_rho_and_sigma(const Selector &selector, bool use_area_4vector, double &median, double &sigma) const override
overload version of what's in the ClusterSequenceAreaBase class, which actually just does the same th...
Definition: ClusterSequenceArea.hh:164
fastjet::ClusterSequenceArea::area_error
virtual double area_error(const PseudoJet &jet) const override
return the error (uncertainty) associated with the determination of the area of this jet
Definition: ClusterSequenceArea.hh:87
fastjet::ClusterSequenceArea::n_empty_jets
virtual double n_empty_jets(const Selector &selector) const override
return something similar to the number of pure ghost jets in the given rap-phi range in an active are...
Definition: ClusterSequenceArea.hh:124
fastjet::JetDefinition
Definition: JetDefinition.hh:250
fastjet::VoronoiAreaSpec
Definition: AreaDefinition.hh:50
fastjet::ClusterSequenceArea::area
virtual double area(const PseudoJet &jet) const override
return the area associated with the given jet
Definition: ClusterSequenceArea.hh:82
fastjet::ClusterSequenceArea::get_median_rho_and_sigma
virtual void get_median_rho_and_sigma(const Selector &selector, bool use_area_4vector, double &median, double &sigma, double &mean_area) const override
overload version of what's in the ClusterSequenceAreaBase class, which actually just does the same th...
Definition: ClusterSequenceArea.hh:176
fastjet::ClusterSequenceArea::empty_area
virtual double empty_area(const Selector &selector) const override
return the total area, corresponding to the given selector, that is free of jets
Definition: ClusterSequenceArea.hh:112
fastjet::ClusterSequenceArea::has_explicit_ghosts
virtual bool has_explicit_ghosts() const override
true if this ClusterSequence has explicit ghosts
Definition: ClusterSequenceArea.hh:134
fastjet::SharedPtr
Definition: SharedPtr.hh:121
fastjet::AreaDefinition
Definition: AreaDefinition.hh:82
fastjet::ClusterSequenceArea
Definition: ClusterSequenceArea.hh:51
fastjet::ClusterSequenceArea::parabolic_pt_per_unit_area
virtual void parabolic_pt_per_unit_area(double &a, double &b, const Selector &selector, double exclude_above=-1.0, bool use_area_4vector=false) const override
overload version of what's in the ClusterSequenceAreaBase class, which additionally checks compatibil...
Definition: ClusterSequenceArea.hh:189
fastjet::ClusterSequenceArea::get_median_rho_and_sigma
virtual void get_median_rho_and_sigma(const std::vector< PseudoJet > &all_jets, const Selector &selector, bool use_area_4vector, double &median, double &sigma, double &mean_area, bool all_are_incl=false) const override
overload version of what's in the ClusterSequenceAreaBase class, which additionally checks compatibil...
Definition: ClusterSequenceArea.hh:147
fastjet::GhostedAreaSpec
Definition: GhostedAreaSpec.hh:65
fastjet::PseudoJet
Definition: PseudoJet.hh:67
fastjet::ClusterSequenceAreaBase
Definition: ClusterSequenceAreaBase.hh:48
fastjet::Selector
Definition: Selector.hh:149
fastjet::LimitedWarning
Definition: LimitedWarning.hh:47
fastjet::ClusterSequenceArea::is_pure_ghost
virtual bool is_pure_ghost(const PseudoJet &jet) const override
true if a jet is made exclusively of ghosts
Definition: ClusterSequenceArea.hh:129
fastjet::ClusterSequenceArea::area_4vector
virtual PseudoJet area_4vector(const PseudoJet &jet) const override
return the 4-vector area
Definition: ClusterSequenceArea.hh:91