FastJet  3.0.2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
ClusterSequenceAreaBase.hh
1 //STARTHEADER
2 // $Id: ClusterSequenceAreaBase.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_CLUSTERSEQUENCEAREABASE_HH__
30 #define __FASTJET_CLUSTERSEQUENCEAREABASE_HH__
31 
32 #include "fastjet/ClusterSequence.hh"
33 #include "fastjet/LimitedWarning.hh"
34 #include "fastjet/Selector.hh"
35 
36 FASTJET_BEGIN_NAMESPACE
37 
38 /// @ingroup area_classes
39 /// \class ClusterSequenceAreaBase
40 /// base class that sets interface for extensions of ClusterSequence
41 /// that provide information about the area of each jet
42 ///
43 /// the virtual functions here all return 0, since no area determination
44 /// is implemented.
46 public:
47 
48  /// a constructor which just carries out the construction of the
49  /// parent class
50  template<class L> ClusterSequenceAreaBase
51  (const std::vector<L> & pseudojets,
52  const JetDefinition & jet_def_in,
53  const bool & writeout_combinations = false) :
54  ClusterSequence(pseudojets, jet_def_in, writeout_combinations) {}
55 
56 
57  /// default constructor
59 
60 
61  /// destructor
63 
64 
65  /// return the area associated with the given jet; this base class
66  /// returns 0.
67  virtual double area (const PseudoJet & ) const {return 0.0;}
68 
69  /// return the error (uncertainty) associated with the determination
70  /// of the area of this jet; this base class returns 0.
71  virtual double area_error (const PseudoJet & ) const {return 0.0;}
72 
73  /// return a PseudoJet whose 4-vector is defined by the following integral
74  ///
75  /// \int drap d\phi PseudoJet("rap,phi,pt=one") *
76  /// * Theta("rap,phi inside jet boundary")
77  ///
78  /// where PseudoJet("rap,phi,pt=one") is a 4-vector with the given
79  /// rapidity (rap), azimuth (phi) and pt=1, while Theta("rap,phi
80  /// inside jet boundary") is a function that is 1 when rap,phi
81  /// define a direction inside the jet boundary and 0 otherwise.
82  ///
83  /// This base class returns a null 4-vector.
84  virtual PseudoJet area_4vector(const PseudoJet & ) const {
85  return PseudoJet(0.0,0.0,0.0,0.0);}
86 
87  /// true if a jet is made exclusively of ghosts
88  ///
89  /// NB: most area classes do not give any explicit ghost jets, but
90  /// some do, and they should replace this function with their own
91  /// version.
92  virtual bool is_pure_ghost(const PseudoJet & ) const {
93  return false;
94  }
95 
96  /// returns true if ghosts are explicitly included within
97  /// jets for this ClusterSequence;
98  ///
99  /// Derived classes that do include explicit ghosts should provide
100  /// an alternative version of this routine and set it properly.
101  virtual bool has_explicit_ghosts() const {
102  return false;
103  }
104 
105  /// return the total area, corresponding to the given Selector, that
106  /// is free of jets, in general based on the inclusive jets.
107  ///
108  /// The selector passed as an argument has to have a finite area and
109  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
110  /// tools for more generic usages)
111  virtual double empty_area(const Selector & selector) const;
112 
113  /// return the total area, corresponding to the given Selector, that
114  /// is free of jets, based on the supplied all_jets
115  ///
116  /// The selector passed as an argument has to have a finite area and
117  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
118  /// tools for more generic usages)
119  double empty_area_from_jets(const std::vector<PseudoJet> & all_jets,
120  const Selector & selector) const;
121 
122  /// return something similar to the number of pure ghost jets
123  /// in the given selector's range in an active area case.
124  /// For the local implementation we return empty_area/(0.55 pi R^2),
125  /// based on measured properties of ghost jets with kt and cam
126  /// (cf arXiv:0802.1188).
127  ///
128  /// Note that the number returned is a double.
129  ///
130  /// The selector passed as an argument has to have a finite area and
131  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
132  /// tools for more generic usages)
133  virtual double n_empty_jets(const Selector & selector) const {
134  double R = jet_def().R();
135  return empty_area(selector)/(0.55*pi*R*R);
136  }
137 
138  /// the median of (pt/area) for jets contained within the selector
139  /// range, making use also of the info on n_empty_jets
140  ///
141  /// The selector passed as an argument has to have a finite area and
142  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
143  /// tools for more generic usages)
144  double median_pt_per_unit_area(const Selector & selector) const;
145 
146  /// the median of (pt/area_4vector) for jets contained within the
147  /// selector range, making use also of the info on n_empty_jets
148  ///
149  /// The selector passed as an argument has to have a finite area and
150  /// apply jet-by-jet
151  double median_pt_per_unit_area_4vector(const Selector & selector) const;
152 
153  /// the function that does the work for median_pt_per_unit_area and
154  /// median_pt_per_unit_area_4vector:
155  /// - something_is_area_4vect = false -> use plain area
156  /// - something_is_area_4vect = true -> use 4-vector area
157  double median_pt_per_unit_something(
158  const Selector & selector, bool use_area_4vector) const;
159 
160  /// using jets withing the selector range (and with 4-vector areas if
161  /// use_area_4vector), calculate the median pt/area, as well as an
162  /// "error" (uncertainty), which is defined as the 1-sigma
163  /// half-width of the distribution of pt/A, obtained by looking for
164  /// the point below which we have (1-0.6827)/2 of the jets
165  /// (including empty jets).
166  ///
167  /// The subtraction for a jet with uncorrected pt pt^U and area A is
168  ///
169  /// pt^S = pt^U - median*A +- sigma*sqrt(A)
170  ///
171  /// where the error is only that associated with the fluctuations
172  /// in the noise and not that associated with the noise having
173  /// caused changes in the hard-particle content of the jet.
174  ///
175  /// The selector passed as an argument has to have a finite area and
176  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
177  /// tools for more generic usages)
178  ///
179  /// NB: subtraction may also be done with 4-vector area of course,
180  /// and this is recommended for jets with larger values of R, as
181  /// long as rho has also been determined with a 4-vector area;
182  /// using a scalar area causes one to neglect terms of relative
183  /// order $R^2/8$ in the jet $p_t$.
184  virtual void get_median_rho_and_sigma(const Selector & selector,
185  bool use_area_4vector,
186  double & median, double & sigma,
187  double & mean_area) const;
188 
189  /// a more advanced version of get_median_rho_and_sigma, which allows
190  /// one to use any "view" of the event containing all jets (so that,
191  /// e.g. one might use Cam on a different resolution scale without
192  /// have to rerun the algorithm).
193  ///
194  /// By default it will assume that "all" are not inclusive jets,
195  /// so that in dealing with empty area it has to calculate
196  /// the number of empty jets based on the empty area and the
197  /// the observed <area> of jets rather than a surmised area
198  ///
199  /// Note that for small effective radii, this can cause problems
200  /// because the harder jets get an area >> <ghost-jet-area>
201  /// and so the estimate comes out all wrong. In these situations
202  /// it is highly advisable to use an area with explicit ghosts, since
203  /// then the "empty" jets are actually visible.
204  ///
205  /// The selector passed as an argument has to have a finite area and
206  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
207  /// tools for more generic usages)
208  virtual void get_median_rho_and_sigma(const std::vector<PseudoJet> & all_jets,
209  const Selector & selector,
210  bool use_area_4vector,
211  double & median, double & sigma,
212  double & mean_area,
213  bool all_are_inclusive = false) const;
214 
215  /// same as the full version of get_median_rho_and_error, but without
216  /// access to the mean_area
217  ///
218  /// The selector passed as an argument has to have a finite area and
219  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
220  /// tools for more generic usages)
221  virtual void get_median_rho_and_sigma(const Selector & selector,
222  bool use_area_4vector,
223  double & median, double & sigma) const {
224  double mean_area;
225  get_median_rho_and_sigma(selector, use_area_4vector,
226  median, sigma, mean_area);
227  }
228 
229 
230  /// fits a form pt_per_unit_area(y) = a + b*y^2 in the selector range.
231  /// exclude_above allows one to exclude large values of pt/area from fit.
232  /// (if negative, the cut is discarded)
233  /// use_area_4vector = true uses the 4vector areas.
234  ///
235  /// The selector passed as an argument has to have a finite area and
236  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
237  /// tools for more generic usages)
238  virtual void parabolic_pt_per_unit_area(double & a, double & b,
239  const Selector & selector,
240  double exclude_above=-1.0,
241  bool use_area_4vector=false) const;
242 
243  /// return a vector of all subtracted jets, using area_4vector, given rho.
244  /// Only inclusive_jets above ptmin are subtracted and returned.
245  /// the ordering is the same as that of sorted_by_pt(cs.inclusive_jets()),
246  /// i.e. not necessarily ordered in pt once subtracted
247  std::vector<PseudoJet> subtracted_jets(const double rho,
248  const double ptmin=0.0) const;
249 
250  /// return a vector of subtracted jets, using area_4vector.
251  /// Only inclusive_jets above ptmin are subtracted and returned.
252  /// the ordering is the same as that of sorted_by_pt(cs.inclusive_jets()),
253  /// i.e. not necessarily ordered in pt once subtracted
254  ///
255  /// The selector passed as an argument has to have a finite area and
256  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
257  /// tools for more generic usages)
258  std::vector<PseudoJet> subtracted_jets(const Selector & selector,
259  const double ptmin=0.0) const;
260 
261  /// return a subtracted jet, using area_4vector, given rho
262  PseudoJet subtracted_jet(const PseudoJet & jet,
263  const double rho) const;
264 
265  /// return a subtracted jet, using area_4vector; note
266  /// that this is potentially inefficient if repeatedly used for many
267  /// different jets, because rho will be recalculated each time
268  /// around.
269  ///
270  /// The selector passed as an argument has to have a finite area and
271  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
272  /// tools for more generic usages)
273  PseudoJet subtracted_jet(const PseudoJet & jet,
274  const Selector & selector) const;
275 
276  /// return the subtracted pt, given rho
277  double subtracted_pt(const PseudoJet & jet,
278  const double rho,
279  bool use_area_4vector=false) const;
280 
281  /// return the subtracted pt; note that this is
282  /// potentially inefficient if repeatedly used for many different
283  /// jets, because rho will be recalculated each time around.
284  ///
285  /// The selector passed as an argument has to have a finite area and
286  /// apply jet-by-jet (see the BackgroundEstimator and Subtractor
287  /// tools for more generic usages)
288  double subtracted_pt(const PseudoJet & jet,
289  const Selector & selector,
290  bool use_area_4vector=false) const;
291 
292 protected:
293  /// check the selector is suited for the computations i.e. applies jet by jet and has a finite area
294  void _check_selector_good_for_median(const Selector &selector) const;
295 
296 
297 private:
298  /// handle warning messages
299  static LimitedWarning _warnings;
300  static LimitedWarning _warnings_zero_area;
301  static LimitedWarning _warnings_empty_area;
302 
303  /// check the jet algorithm is suitable (and if not issue a warning)
304  void _check_jet_alg_good_for_median() const;
305 
306 };
307 
308 
309 
310 FASTJET_END_NAMESPACE
311 
312 #endif // __FASTJET_CLUSTERSEQUENCEAREABASE_HH__