TDSPMath.h Source File

00001 
00009 /***************************************************************************
00010  *                                                                         *
00011  *   This program is free software; you can redistribute it and/or modify  *
00012  *   it under the terms of the GNU General Public License as published by  *
00013  *   the Free Software Foundation; either version 2 of the License, or     *
00014  *   (at your option) any later version.                                   *
00015  *                                                                         *
00016  ***************************************************************************/
00017 
00018 #ifndef TDSPMATH_H
00019 #define TDSPMATH_H
00020 
00021 
00022 #include <TMath.h>
00023 
00029 class TDSPMath {
00030  public: 
00031   
00032   static  Double_t  BesselJ0(Double_t x);
00033   static  Double_t  Q(Double_t x);
00034   static  Double_t  QSqrt(Double_t x);
00035   static  Double_t  Rayleigh(Double_t x,Double_t sigma = 1);
00036   static  Double_t  dB(Double_t x);
00037   static  Double_t  FromdB(Double_t x);
00038   static  Double_t  Lorentz(Double_t x, Double_t sigma = 1);
00039 
00040   ClassDef(TDSPMath,0)
00041 
00042 };
00043 
00044 
00045 // Because BesselJ is not yet in the root-TMath-Lib 
00046 //
00047 // This i adapted from : http://www.if.ufrj.br/teaching/compute/cpp/Bessel_functions.cpp
00048 // 
00049 
00050 inline Double_t TDSPMath::BesselJ0(Double_t x) {
00051   Double_t ax,z;
00052   Double_t xx,y,ans,ans1,ans2;      /* Acumulate polynomials in Double_t */
00053   
00054   if ((ax=TMath::Abs(x)) < 8.0) {       /* Direct rational function fit. */
00055     y=x*x;
00056     ans1=57568490574.0+y*(-13362590354.0+y*(651619640.7
00057                         +y*(-11214424.18+y*(77392.33017+y*(-184.9052456)))));
00058     ans2=57568490411.0+y*(1029532985.0+y*(9494680.718
00059                       +y*(59272.64853+y*(267.8532712+y*1.0))));
00060     ans=ans1/ans2;
00061   } else {      /* Fitting function. */
00062     z=8.0/ax;
00063     y=z*z;
00064     xx=ax-0.785398164;
00065     ans1=1.0+y*(-0.1098628627e-2+y*(0.2734510407e-4
00066                     +y*(-0.2073370639e-5+y*0.2093887211e-6)));
00067     ans2 = -0.1562499995e-1+y*(0.1430488765e-3
00068                    +y*(-0.6911147651e-5+y*(0.7621095161e-6
00069                                -y*0.934935152e-7)));
00070     ans=sqrt(0.636619772/ax)*(cos(xx)*ans1-z*sin(xx)*ans2);
00071   }
00072   return ans;
00073 }                           
00074 
00075 inline Double_t TDSPMath::Q(Double_t x) {
00076   return 0.5*TMath::Erfc(x/TMath::Sqrt(2.0));
00077 }
00078 
00079 inline Double_t TDSPMath::QSqrt(Double_t x) {
00080   return 0.5*TMath::Erfc(TMath::Sqrt(x/2.0));
00081 }
00082 
00083 inline Double_t  TDSPMath::Rayleigh(Double_t x, Double_t sigma) {
00084   return x/(sigma*sigma)*TMath::Exp(-x*x/(2*sigma*sigma));
00085 }
00086 
00087 inline Double_t  TDSPMath::dB(Double_t x) {
00088   return 10*TMath::Log10(x);
00089 }
00090 
00091 inline Double_t  TDSPMath::FromdB(Double_t x) {
00092   return TMath::Power(10,x/10);
00093 }
00094 
00095 inline Double_t  TDSPMath::Lorentz(Double_t x, Double_t sigma) {
00096   return 2.0/(TMath::Pi()*sigma*(1.0+4.0/(sigma*sigma)*x*x));
00097 }
00098 
00099 inline Double_t  dB(Double_t x) {
00100   return TDSPMath::dB(x);
00101 }
00102 
00103 inline Double_t  FromdB(Double_t x) {
00104   return TDSPMath::FromdB(x);
00105 }
00106 
00107 extern const Double_t Pi;
00108 
00109 #endif