 ROOT   6.08/07 Reference Guide
VavilovFast.h
Go to the documentation of this file.
1 // @(#)root/mathmore:$Id$
2 // Authors: B. List 29.4.2010
3
4  /**********************************************************************
5  * *
6  * Copyright (c) 2004 ROOT Foundation, CERN/PH-SFT *
7  * *
8  * This library is free software; you can redistribute it and/or *
9  * modify it under the terms of the GNU General Public License *
10  * as published by the Free Software Foundation; either version 2 *
11  * of the License, or (at your option) any later version. *
12  * *
13  * This library is distributed in the hope that it will be useful, *
14  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
16  * General Public License for more details. *
17  * *
18  * You should have received a copy of the GNU General Public License *
19  * along with this library (see file COPYING); if not, write *
20  * to the Free Software Foundation, Inc., 59 Temple Place, Suite *
21  * 330, Boston, MA 02111-1307 USA, or contact the author. *
22  * *
23  **********************************************************************/
24
25 // Header file for class VavilovFast
26 //
27 // Created by: blist at Thu Apr 29 11:19:00 2010
28 //
29 // Last update: Thu Apr 29 11:19:00 2010
30 //
31 #ifndef ROOT_Math_VavilovFast
32 #define ROOT_Math_VavilovFast
33
34
35 /**
36  @ingroup StatFunc
37  */
38
39
40 #include "Math/Vavilov.h"
41
42 namespace ROOT {
43 namespace Math {
44
45 //____________________________________________________________________________
46 /**
47  Class describing a Vavilov distribution.
48
49  The probability density function of the Vavilov distribution
50  as function of Landau's parameter is given by:
51  \f[ p(\lambda_L; \kappa, \beta^2) =
52  \frac{1}{2 \pi i}\int_{c-i\infty}^{c+i\infty} \phi(s) e^{\lambda_L s} ds\f]
53  where \f$\phi(s) = e^{C} e^{\psi(s)}\f$
54  with \f$C = \kappa (1+\beta^2 \gamma )\f$
55  and \f$\psi(s)= s \ln \kappa + (s+\beta^2 \kappa) 56 \cdot \left ( \int \limits_{0}^{1} 57 \frac{1 - e^{\frac{-st}{\kappa}}}{t} \,d t- \gamma \right ) 58 - \kappa \, e^{\frac{-s}{\kappa}}\f$.
59  \f$\gamma = 0.5772156649\dots\f$ is Euler's constant.
60
61  For the class VavilovFast,
62  Pdf returns the Vavilov distribution as function of Landau's parameter
63  \f$\lambda_L = \lambda_V/\kappa - \ln \kappa\f$,
64  which is the convention used in the CERNLIB routines, and in the tables
65  by S.M. Seltzer and M.J. Berger: Energy loss stragglin of protons and mesons:
66  Tabulation of the Vavilov distribution, pp 187-203
67  in: National Research Council (U.S.), Committee on Nuclear Science:
68  Studies in penetration of charged particles in matter,
69  Nat. Akad. Sci. Publication 1133,
70  Nucl. Sci. Series Report No. 39,
71  Washington (Nat. Akad. Sci.) 1964, 388 pp.
72  Available from
74
75  Therefore, for small values of \f$\kappa < 0.01\f$,
76  pdf approaches the Landau distribution.
77
78  For values \f$\kappa > 10\f$, the Gauss approximation should be used
79  with \f$\mu\f$ and \f$\sigma\f$ given by Vavilov::mean(kappa, beta2)
80  and sqrt(Vavilov::variance(kappa, beta2).
81
82  For values \f$\kappa > 10\f$, the Gauss approximation should be used
83  with \f$\mu\f$ and \f$\sigma\f$ given by Vavilov::mean(kappa, beta2)
84  and sqrt(Vavilov::variance(kappa, beta2).
85
86  The original Vavilov pdf is obtained by
87  v.Pdf(lambdaV/kappa-log(kappa))/kappa.
88
89  For detailed description see
90  A. Rotondi and P. Montagna, Fast calculation of Vavilov distribution,
91  <A HREF="http://dx.doi.org/10.1016/0168-583X(90)90749-K">Nucl. Instr. and Meth. B47 (1990) 215-224</A>,
92  which has been implemented in
93  <A HREF="https://cern-tex.web.cern.ch/cern-tex/shortwrupsdir/g115/top.html">
94  CERNLIB (G115)</A>.
95
96  The class stores coefficients needed to calculate \f$p(\lambda; \kappa, \beta^2)\f$
97  for fixed values of \f$\kappa\f$ and \f$\beta^2\f$.
98  Changing these values is computationally expensive.
99
100  The parameter \f$\kappa\f$ must be in the range \f$0.01 \le \kappa \le 12\f$.
101
102  The parameter \f$\beta^2\f$ must be in the range \f$0 \le \beta^2 \le 1\f$.
103
104  Average times on a Pentium Core2 Duo P8400 2.26GHz:
105  - 9.9us per call to SetKappaBeta2 or constructor
106  - 0.095us per call to Pdf, Cdf
107  - 3.7us per first call to Quantile after SetKappaBeta2 or constructor
108  - 0.137us per subsequent call to Quantile
109
110  Benno List, June 2010
111
112  @ingroup StatFunc
113  */
114
115
116 class VavilovFast: public Vavilov {
117
118 public:
119
120
121  /**
122  Initialize an object to calculate the Vavilov distribution
123
124  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
125  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
126  */
127
128  VavilovFast(double kappa=1, double beta2=1);
129
130
131  /**
132  Destructor
133  */
134  virtual ~VavilovFast();
135
136
137 public:
138
139  /**
140  Evaluate the Vavilov probability density function
141
142  @param x The Landau parameter \f$x = \lambda_L\f$
143  */
144  double Pdf (double x) const;
145
146  /**
147  Evaluate the Vavilov probability density function,
148  and set kappa and beta2, if necessary
149
150  @param x The Landau parameter \f$x = \lambda_L\f$
151  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
152  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
153  */
154  double Pdf (double x, double kappa, double beta2);
155
156  /**
157  Evaluate the Vavilov cumulative probability density function
158
159  @param x The Landau parameter \f$x = \lambda_L\f$
160  */
161  double Cdf (double x) const;
162
163  /**
164  Evaluate the Vavilov cumulative probability density function,
165  and set kappa and beta2, if necessary
166
167  @param x The Landau parameter \f$x = \lambda_L\f$
168  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
169  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
170  */
171  double Cdf (double x, double kappa, double beta2);
172
173  /**
174  Evaluate the Vavilov complementary cumulative probability density function
175
176  @param x The Landau parameter \f$x = \lambda_L\f$
177  */
178  double Cdf_c (double x) const;
179
180  /**
181  Evaluate the Vavilov complementary cumulative probability density function,
182  and set kappa and beta2, if necessary
183
184  @param x The Landau parameter \f$x = \lambda_L\f$
185  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
186  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
187  */
188  double Cdf_c (double x, double kappa, double beta2);
189
190  /**
191  Evaluate the inverse of the Vavilov cumulative probability density function
192
193  @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$
194  */
195  double Quantile (double z) const;
196
197  /**
198  Evaluate the inverse of the Vavilov cumulative probability density function,
199  and set kappa and beta2, if necessary
200
201  @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$
202  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
203  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
204  */
205  double Quantile (double z, double kappa, double beta2);
206
207  /**
208  Evaluate the inverse of the complementary Vavilov cumulative probability density function
209
210  @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$
211  */
212  double Quantile_c (double z) const;
213
214  /**
215  Evaluate the inverse of the complementary Vavilov cumulative probability density function,
216  and set kappa and beta2, if necessary
217
218  @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$
219  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
220  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
221  */
222  double Quantile_c (double z, double kappa, double beta2);
223
224  /**
225  Change \f$\kappa\f$ and \f$\beta^2\f$ and recalculate coefficients if necessary
226
227  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
228  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
229  */
230  virtual void SetKappaBeta2 (double kappa, double beta2);
231
232  /**
233  Return the minimum value of \f$\lambda\f$ for which \f$p(\lambda; \kappa, \beta^2)\f$
234  is nonzero in the current approximation
235  */
236  virtual double GetLambdaMin() const;
237
238  /**
239  Return the maximum value of \f$\lambda\f$ for which \f$p(\lambda; \kappa, \beta^2)\f$
240  is nonzero in the current approximation
241  */
242  virtual double GetLambdaMax() const;
243
244  /**
245  Return the current value of \f$\kappa\f$
246  */
247  virtual double GetKappa() const;
248
249  /**
250  Return the current value of \f$\beta^2\f$
251  */
252  virtual double GetBeta2() const;
253
254  /**
255  Returns a static instance of class VavilovFast
256  */
257  static VavilovFast *GetInstance();
258
259  /**
260  Returns a static instance of class VavilovFast,
261  and sets the values of kappa and beta2
262
263  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
264  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
265  */
266  static VavilovFast *GetInstance(double kappa, double beta2);
267
268
269 private:
270  double fKappa;
271  double fBeta2;
272
273  double fAC;
274  double fHC;
275  double fWCM;
276  int fItype;
277  int fNpt;
278
280
281 };
282
283  /**
284  The Vavilov probability density function
285
286  @param x The Landau parameter \f$x = \lambda_L\f$
287  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
288  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
289
290  @ingroup PdfFunc
291  */
292 double vavilov_fast_pdf (double x, double kappa, double beta2);
293
294  /**
295  The Vavilov cumulative probability density function
296
297  @param x The Landau parameter \f$x = \lambda_L\f$
298  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
299  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
300
301  @ingroup ProbFunc
302  */
303 double vavilov_fast_cdf (double x, double kappa, double beta2);
304
305  /**
306  The Vavilov complementary cumulative probability density function
307
308  @param x The Landau parameter \f$x = \lambda_L\f$
309  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
310  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
311
312  @ingroup ProbFunc
313  */
314 double vavilov_fast_cdf_c (double x, double kappa, double beta2);
315
316  /**
317  The inverse of the Vavilov cumulative probability density function
318
319  @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$
320  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
321  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
322
323  @ingroup QuantFunc
324  */
325 double vavilov_fast_quantile (double z, double kappa, double beta2);
326
327  /**
328  The inverse of the complementary Vavilov cumulative probability density function
329
330  @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$
331  @param kappa The parameter \f$\kappa\f$, which must be in the range \f$0.01 \le \kappa \le 12 \f$
332  @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$
333
334  @ingroup QuantFunc
335  */
336 double vavilov_fast_quantile_c (double z, double kappa, double beta2);
337
338 } // namespace Math
339 } // namespace ROOT
340
341 #endif /* ROOT_Math_VavilovFast */
double vavilov_fast_pdf(double x, double kappa, double beta2)
The Vavilov probability density function.
Base class describing a Vavilov distribution.
Definition: Vavilov.h:121
virtual double GetLambdaMax() const
Return the maximum value of for which is nonzero in the current approximation.
virtual ~VavilovFast()
Destructor.
Definition: VavilovFast.cxx:58
This namespace contains pre-defined functions to be used in conjuction with TExecutor::Map and TExecu...
Definition: StringConv.hxx:21
virtual double GetKappa() const
Return the current value of .
virtual void SetKappaBeta2(double kappa, double beta2)
Change and and recalculate coefficients if necessary.
Definition: VavilovFast.cxx:63
double Quantile(double z) const
Evaluate the inverse of the Vavilov cumulative probability density function.
double Cdf(double x) const
Evaluate the Vavilov cumulative probability density function.
virtual double GetLambdaMin() const
Return the minimum value of for which is nonzero in the current approximation.
VavilovFast(double kappa=1, double beta2=1)
Initialize an object to calculate the Vavilov distribution.
Definition: VavilovFast.cxx:52
Class describing a Vavilov distribution.
Definition: VavilovFast.h:116
double Cdf_c(double x) const
Evaluate the Vavilov complementary cumulative probability density function.
Double_t x[n]
Definition: legend1.C:17
double vavilov_fast_quantile_c(double z, double kappa, double beta2)
The inverse of the complementary Vavilov cumulative probability density function. ...
static VavilovFast * GetInstance()
Returns a static instance of class VavilovFast.
double vavilov_fast_cdf(double x, double kappa, double beta2)
The Vavilov cumulative probability density function.
double vavilov_fast_quantile(double z, double kappa, double beta2)
The inverse of the Vavilov cumulative probability density function.
Namespace for new Math classes and functions.
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
virtual double GetBeta2() const
Return the current value of .
static VavilovFast * fgInstance
Definition: VavilovFast.h:279
double Quantile_c(double z) const
Evaluate the inverse of the complementary Vavilov cumulative probability density function.
double Pdf(double x) const
Evaluate the Vavilov probability density function.
double vavilov_fast_cdf_c(double x, double kappa, double beta2)
The Vavilov complementary cumulative probability density function.