|
||||
File indexing completed on 2025-01-18 10:10:22
0001 // @(#)root/mathmore:$Id$ 0002 // Authors: B. List 29.4.2010 0003 0004 /********************************************************************** 0005 * * 0006 * Copyright (c) 2004 ROOT Foundation, CERN/PH-SFT * 0007 * * 0008 * This library is free software; you can redistribute it and/or * 0009 * modify it under the terms of the GNU General Public License * 0010 * as published by the Free Software Foundation; either version 2 * 0011 * of the License, or (at your option) any later version. * 0012 * * 0013 * This library is distributed in the hope that it will be useful, * 0014 * but WITHOUT ANY WARRANTY; without even the implied warranty of * 0015 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU * 0016 * General Public License for more details. * 0017 * * 0018 * You should have received a copy of the GNU General Public License * 0019 * along with this library (see file COPYING); if not, write * 0020 * to the Free Software Foundation, Inc., 59 Temple Place, Suite * 0021 * 330, Boston, MA 02111-1307 USA, or contact the author. * 0022 * * 0023 **********************************************************************/ 0024 0025 // Header file for class Vavilov 0026 // 0027 // Created by: blist at Thu Apr 29 11:19:00 2010 0028 // 0029 // Last update: Thu Apr 29 11:19:00 2010 0030 // 0031 #ifndef ROOT_Math_Vavilov 0032 #define ROOT_Math_Vavilov 0033 0034 0035 0036 /** 0037 @ingroup StatFunc 0038 */ 0039 0040 0041 0042 namespace ROOT { 0043 namespace Math { 0044 0045 //____________________________________________________________________________ 0046 /** 0047 Base class describing a Vavilov distribution 0048 0049 The Vavilov distribution is defined in 0050 P.V. Vavilov: Ionization losses of high-energy heavy particles, 0051 Sov. Phys. JETP 5 (1957) 749 [Zh. Eksp. Teor. Fiz. 32 (1957) 920]. 0052 0053 The probability density function of the Vavilov distribution 0054 as function of Landau's parameter is given by: 0055 \f[ p(\lambda_L; \kappa, \beta^2) = 0056 \frac{1}{2 \pi i}\int_{c-i\infty}^{c+i\infty} \phi(s) e^{\lambda_L s} ds\f] 0057 where \f$\phi(s) = e^{C} e^{\psi(s)}\f$ 0058 with \f$ C = \kappa (1+\beta^2 \gamma )\f$ 0059 and \f$\psi(s)= s \ln \kappa + (s+\beta^2 \kappa) 0060 \cdot \left ( \int \limits_{0}^{1} 0061 \frac{1 - e^{\frac{-st}{\kappa}}}{t} \,d t- \gamma \right ) 0062 - \kappa \, e^{\frac{-s}{\kappa}}\f$. 0063 \f$ \gamma = 0.5772156649\dots\f$ is Euler's constant. 0064 0065 For the class Vavilov, 0066 Pdf returns the Vavilov distribution as function of Landau's parameter 0067 \f$\lambda_L = \lambda_V/\kappa - \ln \kappa\f$, 0068 which is the convention used in the CERNLIB routines, and in the tables 0069 by S.M. Seltzer and M.J. Berger: Energy loss stragglin of protons and mesons: 0070 Tabulation of the Vavilov distribution, pp 187-203 0071 in: National Research Council (U.S.), Committee on Nuclear Science: 0072 Studies in penetration of charged particles in matter, 0073 Nat. Akad. Sci. Publication 1133, 0074 Nucl. Sci. Series Report No. 39, 0075 Washington (Nat. Akad. Sci.) 1964, 388 pp. 0076 Available from 0077 <A HREF="http://books.google.de/books?id=kmMrAAAAYAAJ&lpg=PP9&pg=PA187#v=onepage&q&f=false">Google books</A> 0078 0079 Therefore, for small values of \f$\kappa < 0.01\f$, 0080 pdf approaches the Landau distribution. 0081 0082 For values \f$\kappa > 10\f$, the Gauss approximation should be used 0083 with \f$\mu\f$ and \f$\sigma\f$ given by Vavilov::Mean(kappa, beta2) 0084 and sqrt(Vavilov::Variance(kappa, beta2). 0085 0086 The original Vavilov pdf is obtained by 0087 v.Pdf(lambdaV/kappa-log(kappa))/kappa. 0088 0089 Two subclasses are provided: 0090 - VavilovFast uses the algorithm by 0091 A. Rotondi and P. Montagna, Fast calculation of Vavilov distribution, 0092 <A HREF="http://dx.doi.org/10.1016/0168-583X(90)90749-K">Nucl. Instr. and Meth. B47 (1990) 215-224</A>, 0093 which has been implemented in 0094 <A HREF="https://cern-tex.web.cern.ch/cern-tex/shortwrupsdir/g115/top.html"> 0095 CERNLIB (G115)</A>. 0096 0097 - VavilovAccurate uses the algorithm by 0098 B. Schorr, Programs for the Landau and the Vavilov distributions and the corresponding random numbers, 0099 <A HREF="http://dx.doi.org/10.1016/0010-4655(74)90091-5">Computer Phys. Comm. 7 (1974) 215-224</A>, 0100 which has been implemented in 0101 <A HREF="https://cern-tex.web.cern.ch/cern-tex/shortwrupsdir/g116/top.html"> 0102 CERNLIB (G116)</A>. 0103 0104 Both subclasses store coefficients needed to calculate \f$p(\lambda; \kappa, \beta^2)\f$ 0105 for fixed values of \f$\kappa\f$ and \f$\beta^2\f$. 0106 Changing these values is computationally expensive. 0107 0108 VavilovFast is about 5 times faster for the calculation of the Pdf than VavilovAccurate; 0109 initialization takes about 100 times longer than calculation of the Pdf value. 0110 For the quantile calculation, VavilovFast 0111 is 30 times faster for the initialization, and 6 times faster for 0112 subsequent calculations. Initialization for Quantile takes 0113 27 (11) times longer than subsequent calls for VavilovFast (VavilovAccurate). 0114 0115 @ingroup StatFunc 0116 0117 */ 0118 0119 0120 class Vavilov { 0121 0122 public: 0123 0124 0125 /** 0126 Default constructor 0127 */ 0128 Vavilov(); 0129 0130 /** 0131 Destructor 0132 */ 0133 virtual ~Vavilov(); 0134 0135 public: 0136 0137 /** 0138 Evaluate the Vavilov probability density function 0139 0140 @param x The Landau parameter \f$x = \lambda_L\f$ 0141 0142 */ 0143 virtual double Pdf (double x) const = 0; 0144 0145 /** 0146 Evaluate the Vavilov probability density function, 0147 and set kappa and beta2, if necessary 0148 0149 @param x The Landau parameter \f$x = \lambda_L\f$ 0150 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0151 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0152 */ 0153 virtual double Pdf (double x, double kappa, double beta2) = 0; 0154 0155 /** 0156 Evaluate the Vavilov cumulative probability density function 0157 0158 @param x The Landau parameter \f$x = \lambda_L\f$ 0159 */ 0160 virtual double Cdf (double x) const = 0; 0161 0162 /** 0163 Evaluate the Vavilov cumulative probability density function, 0164 and set kappa and beta2, if necessary 0165 0166 @param x The Landau parameter \f$x = \lambda_L\f$ 0167 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0168 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0169 */ 0170 virtual double Cdf (double x, double kappa, double beta2) = 0; 0171 0172 /** 0173 Evaluate the Vavilov complementary cumulative probability density function 0174 0175 @param x The Landau parameter \f$x = \lambda_L\f$ 0176 */ 0177 virtual double Cdf_c (double x) const = 0; 0178 0179 /** 0180 Evaluate the Vavilov complementary cumulative probability density function, 0181 and set kappa and beta2, if necessary 0182 0183 @param x The Landau parameter \f$x = \lambda_L\f$ 0184 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0185 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0186 */ 0187 virtual double Cdf_c (double x, double kappa, double beta2) = 0; 0188 0189 /** 0190 Evaluate the inverse of the Vavilov cumulative probability density function 0191 0192 @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$ 0193 */ 0194 virtual double Quantile (double z) const = 0; 0195 0196 /** 0197 Evaluate the inverse of the Vavilov cumulative probability density function, 0198 and set kappa and beta2, if necessary 0199 0200 @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$ 0201 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0202 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0203 */ 0204 virtual double Quantile (double z, double kappa, double beta2) = 0; 0205 0206 /** 0207 Evaluate the inverse of the complementary Vavilov cumulative probability density function 0208 0209 @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$ 0210 */ 0211 virtual double Quantile_c (double z) const = 0; 0212 0213 /** 0214 Evaluate the inverse of the complementary Vavilov cumulative probability density function, 0215 and set kappa and beta2, if necessary 0216 0217 @param z The argument \f$z\f$, which must be in the range \f$0 \le z \le 1\f$ 0218 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0219 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0220 */ 0221 virtual double Quantile_c (double z, double kappa, double beta2) = 0; 0222 0223 /** 0224 Change \f$\kappa\f$ and \f$\beta^2\f$ and recalculate coefficients if necessary 0225 0226 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0227 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0228 */ 0229 virtual void SetKappaBeta2 (double kappa, double beta2) = 0; 0230 0231 /** 0232 Return the minimum value of \f$\lambda\f$ for which \f$p(\lambda; \kappa, \beta^2)\f$ 0233 is nonzero in the current approximation 0234 */ 0235 virtual double GetLambdaMin() const = 0; 0236 0237 /** 0238 Return the maximum value of \f$\lambda\f$ for which \f$p(\lambda; \kappa, \beta^2)\f$ 0239 is nonzero in the current approximation 0240 */ 0241 virtual double GetLambdaMax() const = 0; 0242 0243 /** 0244 Return the current value of \f$\kappa\f$ 0245 */ 0246 virtual double GetKappa() const = 0; 0247 0248 /** 0249 Return the current value of \f$\beta^2\f$ 0250 */ 0251 virtual double GetBeta2() const = 0; 0252 0253 /** 0254 Return the value of \f$\lambda\f$ where the pdf is maximal 0255 */ 0256 virtual double Mode() const; 0257 0258 /** 0259 Return the value of \f$\lambda\f$ where the pdf is maximal function, 0260 and set kappa and beta2, if necessary 0261 0262 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0263 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0264 */ 0265 virtual double Mode(double kappa, double beta2); 0266 0267 /** 0268 Return the theoretical mean \f$\mu = \gamma-1- \ln \kappa - \beta^2\f$, 0269 where \f$\gamma = 0.5772\dots\f$ is Euler's constant 0270 */ 0271 virtual double Mean() const; 0272 0273 /** 0274 Return the theoretical variance \f$\sigma^2 = \frac{1 - \beta^2/2}{\kappa}\f$ 0275 */ 0276 virtual double Variance() const; 0277 0278 /** 0279 Return the theoretical skewness 0280 \f$\gamma_1 = \frac{1/2 - \beta^2/3}{\kappa^2 \sigma^3} \f$ 0281 */ 0282 virtual double Skewness() const; 0283 0284 /** 0285 Return the theoretical kurtosis 0286 \f$\gamma_2 = \frac{1/3 - \beta^2/4}{\kappa^3 \sigma^4}\f$ 0287 */ 0288 virtual double Kurtosis() const; 0289 0290 /** 0291 Return the theoretical Mean \f$\mu = \gamma-1- \ln \kappa - \beta^2\f$ 0292 0293 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0294 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0295 */ 0296 static double Mean(double kappa, double beta2); 0297 0298 /** 0299 Return the theoretical Variance \f$\sigma^2 = \frac{1 - \beta^2/2}{\kappa}\f$ 0300 0301 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0302 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0303 */ 0304 static double Variance(double kappa, double beta2); 0305 0306 /** 0307 Return the theoretical skewness 0308 \f$\gamma_1 = \frac{1/2 - \beta^2/3}{\kappa^2 \sigma^3} \f$ 0309 0310 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0311 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0312 */ 0313 static double Skewness(double kappa, double beta2); 0314 0315 /** 0316 Return the theoretical kurtosis 0317 \f$\gamma_2 = \frac{1/3 - \beta^2/4}{\kappa^3 \sigma^4}\f$ 0318 0319 @param kappa The parameter \f$\kappa\f$, which should be in the range \f$0.01 \le \kappa \le 10 \f$ 0320 @param beta2 The parameter \f$\beta^2\f$, which must be in the range \f$0 \le \beta^2 \le 1 \f$ 0321 */ 0322 static double Kurtosis(double kappa, double beta2); 0323 0324 0325 }; 0326 0327 } // namespace Math 0328 } // namespace ROOT 0329 0330 #endif /* ROOT_Math_Vavilov */
[ Source navigation ] | [ Diff markup ] | [ Identifier search ] | [ general search ] |
This page was automatically generated by the 2.3.7 LXR engine. The LXR team |