1 | #ifndef __FASTJET_BACKGROUND_ESTIMATOR_BASE_HH__ |
---|
2 | #define __FASTJET_BACKGROUND_ESTIMATOR_BASE_HH__ |
---|
3 | |
---|
4 | //STARTHEADER |
---|
5 | // $Id: BackgroundEstimatorBase.hh 2689 2011-11-14 14:51:06Z soyez $ |
---|
6 | // |
---|
7 | // Copyright (c) 2005-2011, Matteo Cacciari, Gavin P. Salam and Gregory Soyez |
---|
8 | // |
---|
9 | //---------------------------------------------------------------------- |
---|
10 | // This file is part of FastJet. |
---|
11 | // |
---|
12 | // FastJet is free software; you can redistribute it and/or modify |
---|
13 | // it under the terms of the GNU General Public License as published by |
---|
14 | // the Free Software Foundation; either version 2 of the License, or |
---|
15 | // (at your option) any later version. |
---|
16 | // |
---|
17 | // The algorithms that underlie FastJet have required considerable |
---|
18 | // development and are described in hep-ph/0512210. If you use |
---|
19 | // FastJet as part of work towards a scientific publication, please |
---|
20 | // include a citation to the FastJet paper. |
---|
21 | // |
---|
22 | // FastJet is distributed in the hope that it will be useful, |
---|
23 | // but WITHOUT ANY WARRANTY; without even the implied warranty of |
---|
24 | // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
---|
25 | // GNU General Public License for more details. |
---|
26 | // |
---|
27 | // You should have received a copy of the GNU General Public License |
---|
28 | // along with FastJet. If not, see <http://www.gnu.org/licenses/>. |
---|
29 | //---------------------------------------------------------------------- |
---|
30 | //ENDHEADER |
---|
31 | |
---|
32 | #include <fastjet/ClusterSequenceAreaBase.hh> |
---|
33 | #include <fastjet/FunctionOfPseudoJet.hh> |
---|
34 | #include <fastjet/Selector.hh> |
---|
35 | #include <fastjet/Error.hh> |
---|
36 | #include <iostream> |
---|
37 | |
---|
38 | FASTJET_BEGIN_NAMESPACE // defined in fastjet/internal/base.hh |
---|
39 | |
---|
40 | |
---|
41 | /// @ingroup tools_background |
---|
42 | /// \class BackgroundEstimatorBase |
---|
43 | /// |
---|
44 | /// Abstract base class that provides the basic interface for classes |
---|
45 | /// that estimate levels of background radiation in hadrion and |
---|
46 | /// heavy-ion collider events. |
---|
47 | /// |
---|
48 | /// |
---|
49 | class BackgroundEstimatorBase { |
---|
50 | public: |
---|
51 | /// @name constructors and destructors |
---|
52 | //\{ |
---|
53 | //---------------------------------------------------------------- |
---|
54 | BackgroundEstimatorBase() : _rescaling_class(0) {} |
---|
55 | //\} |
---|
56 | |
---|
57 | /// a default virtual destructor that does nothing |
---|
58 | virtual ~BackgroundEstimatorBase() {} |
---|
59 | |
---|
60 | |
---|
61 | /// @name setting a new event |
---|
62 | //\{ |
---|
63 | //---------------------------------------------------------------- |
---|
64 | |
---|
65 | /// tell the background estimator that it has a new event, composed |
---|
66 | /// of the specified particles. |
---|
67 | virtual void set_particles(const std::vector<PseudoJet> & particles) = 0; |
---|
68 | |
---|
69 | //\} |
---|
70 | |
---|
71 | /// @name retrieving fundamental information |
---|
72 | //\{ |
---|
73 | //---------------------------------------------------------------- |
---|
74 | |
---|
75 | /// get rho, the background density per unit area |
---|
76 | virtual double rho() const = 0; |
---|
77 | |
---|
78 | /// get sigma, the background fluctuations per unit area; must be |
---|
79 | /// multipled by sqrt(area) to get fluctuations for a region of a |
---|
80 | /// given area. |
---|
81 | virtual double sigma() const { |
---|
82 | throw Error("sigma() not supported for this Background Estimator"); |
---|
83 | } |
---|
84 | |
---|
85 | /// get rho, the background density per unit area, locally at the |
---|
86 | /// position of a given jet. Note that this is not const, because a |
---|
87 | /// user may then wish to query other aspects of the background that |
---|
88 | /// could depend on the position of the jet last used for a rho(jet) |
---|
89 | /// determination. |
---|
90 | virtual double rho(const PseudoJet & jet) = 0; |
---|
91 | |
---|
92 | /// get sigma, the background fluctuations per unit area, locally at |
---|
93 | /// the position of a given jet. As for rho(jet), it is non-const. |
---|
94 | virtual double sigma(const PseudoJet & /*jet*/) { |
---|
95 | throw Error("sigma(jet) not supported for this Background Estimator"); |
---|
96 | } |
---|
97 | |
---|
98 | /// returns true if this background estimator has support for |
---|
99 | /// determination of sigma |
---|
100 | virtual bool has_sigma() {return false;} |
---|
101 | //\} |
---|
102 | |
---|
103 | |
---|
104 | /// @name configuring the behaviour |
---|
105 | //\{ |
---|
106 | //---------------------------------------------------------------- |
---|
107 | |
---|
108 | /// Set a pointer to a class that calculates the rescaling factor as |
---|
109 | /// a function of the jet (position). Note that the rescaling factor |
---|
110 | /// is used both in the determination of the "global" rho (the pt/A |
---|
111 | /// of each jet is divided by this factor) and when asking for a |
---|
112 | /// local rho (the result is multiplied by this factor). |
---|
113 | /// |
---|
114 | /// The BackgroundRescalingYPolynomial class can be used to get a |
---|
115 | /// rescaling that depends just on rapidity. |
---|
116 | virtual void set_rescaling_class(const FunctionOfPseudoJet<double> * rescaling_class_in) { _rescaling_class = rescaling_class_in; } |
---|
117 | |
---|
118 | /// return the pointer to the jet density class |
---|
119 | const FunctionOfPseudoJet<double> * rescaling_class() const{ |
---|
120 | return _rescaling_class; |
---|
121 | } |
---|
122 | |
---|
123 | //\} |
---|
124 | |
---|
125 | /// @name description |
---|
126 | //\{ |
---|
127 | //---------------------------------------------------------------- |
---|
128 | |
---|
129 | /// returns a textual description of the background estimator |
---|
130 | virtual std::string description() const = 0; |
---|
131 | |
---|
132 | //\} |
---|
133 | |
---|
134 | protected: |
---|
135 | /// @name helpers for derived classes |
---|
136 | /// |
---|
137 | /// Note that these helpers are related to median-based estimation |
---|
138 | /// of the background, so there is no guarantee that they will |
---|
139 | /// remain in this base class in the long term |
---|
140 | //\{ |
---|
141 | //---------------------------------------------------------------- |
---|
142 | |
---|
143 | /// given a quantity in a vector (e.g. pt_over_area) and knowledge |
---|
144 | /// about the number of empty jets, calculate the median and |
---|
145 | /// stand_dev_if_gaussian (roughly from the 16th percentile) |
---|
146 | /// |
---|
147 | /// If do_fj2_calculation is set to true then this performs FastJet |
---|
148 | /// 2.X estimation of the standard deviation, which has a spurious |
---|
149 | /// offset in the limit of a small number of jets. |
---|
150 | void _median_and_stddev(const std::vector<double> & quantity_vector, |
---|
151 | double n_empty_jets, |
---|
152 | double & median, |
---|
153 | double & stand_dev_if_gaussian, |
---|
154 | bool do_fj2_calculation = false |
---|
155 | ) const; |
---|
156 | |
---|
157 | /// computes a percentile of a given _sorted_ vector |
---|
158 | /// \param sorted_quantity_vector the vector contains the data sample |
---|
159 | /// \param percentile the percentile (defined between 0 and 1) to compute |
---|
160 | /// \param nempty an additional number of 0's |
---|
161 | /// (considered at the beginning of |
---|
162 | /// the quantity vector) |
---|
163 | /// \param do_fj2_calculation carry out the calculation as it |
---|
164 | /// was done in fj2 (suffers from "edge effects") |
---|
165 | double _percentile(const std::vector<double> & sorted_quantity_vector, |
---|
166 | const double percentile, |
---|
167 | const double nempty=0.0, |
---|
168 | const bool do_fj2_calculation = false) const; |
---|
169 | |
---|
170 | //\} |
---|
171 | |
---|
172 | const FunctionOfPseudoJet<double> * _rescaling_class; |
---|
173 | static LimitedWarning _warnings_empty_area; |
---|
174 | }; |
---|
175 | |
---|
176 | |
---|
177 | |
---|
178 | //---------------------------------------------------------------------- |
---|
179 | /// @ingroup tools_background |
---|
180 | /// A background rescaling that is a simple polynomial in y |
---|
181 | class BackgroundRescalingYPolynomial : public FunctionOfPseudoJet<double> { |
---|
182 | public: |
---|
183 | /// construct a background rescaling polynomial of the form |
---|
184 | /// a0 + a1*y + a2*y^2 + a3*y^3 + a4*y^4 |
---|
185 | /// |
---|
186 | /// The following values give a reasonable reproduction of the |
---|
187 | /// Pythia8 tune 4C background shape for pp collisions at |
---|
188 | /// sqrt(s)=7TeV: |
---|
189 | /// |
---|
190 | /// - a0 = 1.157 |
---|
191 | /// - a1 = 0 |
---|
192 | /// - a2 = -0.0266 |
---|
193 | /// - a3 = 0 |
---|
194 | /// - a4 = 0.000048 |
---|
195 | /// |
---|
196 | BackgroundRescalingYPolynomial(double a0=1, |
---|
197 | double a1=0, |
---|
198 | double a2=0, |
---|
199 | double a3=0, |
---|
200 | double a4=0) : _a0(a0), _a1(a1), _a2(a2), _a3(a3), _a4(a4) {} |
---|
201 | |
---|
202 | /// return the rescaling factor associated with this jet |
---|
203 | virtual double result(const PseudoJet & jet) const; |
---|
204 | private: |
---|
205 | double _a0, _a1, _a2, _a3, _a4; |
---|
206 | }; |
---|
207 | |
---|
208 | |
---|
209 | |
---|
210 | |
---|
211 | |
---|
212 | FASTJET_END_NAMESPACE |
---|
213 | |
---|
214 | #endif // __BACKGROUND_ESTIMATOR_BASE_HH__ |
---|
215 | |
---|