TEvtGen/EvtGenModels/EvtVubNLO.hh

   1 //--------------------------------------------------------------------------
   2 //
   3 // Environment:
   4 //      This software is part of the EvtGen package developed jointly
   5 //      for the BaBar and CLEO collaborations.  If you use all or part
   6 //      of it, please give an appropriate acknowledgement.
   7 //
   8 // Copyright Information: See EvtGen/COPYRIGHT
   9 //      Copyright (C) 1998      Caltech, UCSB
  10 //
  11 // Module: EvtGen/EvtVubNLO.hh
  12 //
  13 // Description:
  14 // Class to generate inclusive B to X_u l nu decays according to various
  15 // decay models. Implemtented are ACCM, parton-model and a QCD model.
  16 //
  17 // Modification history:
  18 //
  19 //   Sven Menke     January 17, 2001         Module created
  20 //
  21 //------------------------------------------------------------------------
  22
  23 #ifndef EVTVUBNLO_HH
  24 #define EVTVUBNLO_HH
  25
  26 #include <vector>
  27 #include "EvtGenBase/EvtDecayIncoherent.hh"
  28
  29 class EvtParticle;
  30 class RandGeneral;
  31
  32 class EvtVubNLO:public  EvtDecayIncoherent  {
  33
  34 public:
  35
  36   EvtVubNLO() {}
  37   virtual ~EvtVubNLO();
  38
  39   std::string getName();
  40
  41   EvtDecayBase* clone();
  42
  43   void initProbMax();
  44
  45   void init();
  46
  47   void decay(EvtParticle *p);
  48
  49
  50 private:
  51
  52   // cache
  53   double _lbar;
  54   double _mupi2;
  55
  56   double _mb;     // the b-quark pole mass in GeV
  57   double _mB;
  58   double _lambdaSF;
  59   double _b;      // Parameter for the Fermi Motion
  60   double _kpar;
  61   double _mui; // renormalization scale (preferred value=1.5 GeV)
  62   double _SFNorm; // SF normalization
  63   double _dGMax;  // max dGamma*p2 value;
  64   int    _nbins;
  65   int    _idSF;// which shape function?
  66   double * _masses;
  67   double * _weights;
  68
  69   double _gmax;
  70   int _ngood,_ntot;
  71
  72
  73   double tripleDiff(double pp, double pl, double pm);
  74   double SFNorm(const std::vector<double> &coeffs);
  75   static double integrand(double omega, const std::vector<double> &coeffs);
  76   double F10(const std::vector<double> &coeffs);
  77   static double F1Int(double omega,const std::vector<double> &coeffs);
  78   double F20(const std::vector<double> &coeffs);
  79   static double F2Int(double omega,const std::vector<double> &coeffs);
  80   double F30(const std::vector<double> &coeffs);
  81   static double F3Int(double omega,const std::vector<double> &coeffs);
  82   static double g1(double y, double z);
  83   static double g2(double y, double z);
  84   static double g3(double y, double z);
  85
  86   static double Gamma(double z);// Euler Gamma Function
  87   static double dgamma(double t, const std::vector<double> &c){  return pow(t,c[0]-1)*exp(-t);}
  88   static double Gamma(double z, double tmax);
  89
  90   // theory parameters
  91   inline double mu_i(){return _mui;} // intermediate scale
  92   inline double mu_bar(){return _mui;}
  93   inline double mu_h(){return _mb/sqrt(2.0);} // high scale
  94   inline double lambda1(){return -_mupi2;}
  95
  96   // expansion coefficients for RGE
  97   static double beta0(int nf=4){return 11.-2./3.*nf;}
  98   static double beta1(int nf=4){return 34.*3.-38./3.*nf;}
  99   static double beta2(int nf=4){return 1428.5-5033./18.*nf+325./54.*nf*nf;}
 100   static double gamma0(){return 16./3.;}
 101   static double gamma1(int nf=4){return 4./3.*(49.85498-40./9.*nf);}
 102   static double gamma2(int nf=4){return 64./3.*(55.07242-8.58691*nf-nf*nf/27.);} /*  zeta3=1.20206 */
 103   static double gammap0(){return -20./3.;}
 104   static double gammap1(int nf=4){return -32./3.*(6.92653-0.9899*nf);} /* ??  zeta3=1.202 */
 105
 106
 107   // running constants
 108
 109   static double alphas(double mu) ;
 110   static double C_F(double mu){return (4.0/3.0)*alphas(mu)/4./EvtConst::pi;}
 111
 112   // Shape Functions
 113
 114   inline double lambda_SF(){ return _lambdaSF;}
 115   double lambda_bar(double omega0);
 116   inline double lambda2(){return 0.12;}
 117   double mu_pi2(double omega0);
 118   inline double lambda(double mu=0){ return _mB-_mb;}
 119
 120   // specail for gaussian SF
 121   static double cGaus(double b){return pow(Gamma(1+b/2.)/Gamma((1+b)/2.),2);}
 122
 123   double M0(double mui,double omega0);
 124   static double shapeFunction(double omega, const std::vector<double> &coeffs);
 125   static double expShapeFunction(double omega, const std::vector<double> &coeffs);
 126   static double gausShapeFunction(double omega, const std::vector<double> &coeffs);
 127   // SSF (not yet implemented)
 128   double subS(const std::vector<double> &coeffs );
 129   double subT(const std::vector<double> &coeffs);
 130   double subU(const std::vector<double> &coeffs);
 131   double subV(const std::vector<double> &coeffs);
 132
 133
 134   // Sudakov
 135
 136   inline double S0(double a, double r){return -gamma0()/4/a/pow(beta0(),2)*(1/r-1+log(r));}
 137   inline double S1(double a, double r){return gamma0()/4./pow(beta0(),2)*(
 138                                                                   pow(log(r),2)*beta1()/2./beta0()+(gamma1()/gamma0()-beta1()/beta0())*(1.-r+log(r))
 139                                                                   );}
 140   inline double S2(double a, double r){return gamma0()*a/4./pow(beta0(),2)*(
 141                                                                            -0.5*pow((1-r),2)*(
 142                                                                                          pow(beta1()/beta0(),2)-beta2()/beta0()-beta1()/beta0()*gamma1()/gamma0()+gamma2()/gamma0()
 143                                                                                          )
 144                                                                            +(pow(beta1()/beta0(),2)-beta2()/beta0())*(1-r)*log(r)
 145                                                                            +(beta1()/beta0()*gamma1()/gamma0()-beta2()/beta0())*(1-r+r*log(r))
 146                                                                            );}
 147   inline double dSudakovdepsi(double mu1, double mu2){return S2(alphas(mu1)/(4*EvtConst::pi),alphas(mu2)/alphas(mu1));}
 148   inline double Sudakov(double mu1, double mu2, double epsi=0){double fp(4*EvtConst::pi);return S0(alphas(mu1)/fp,alphas(mu2)/alphas(mu1))+S1(alphas(mu1)/fp,alphas(mu2)/alphas(mu1))+epsi*dSudakovdepsi(mu1,mu2);}
 149
 150   // RG
 151   inline double dGdepsi(double mu1, double mu2){return 1./8./EvtConst::pi*(alphas(mu2)-alphas(mu1))*(gamma1()/beta0()-beta1()*gamma0()/pow(beta0(),2));}
 152   inline double aGamma(double mu1, double mu2, double epsi=0){return gamma0()/2/beta0()*log(alphas(mu2)/alphas(mu1))+epsi*dGdepsi( mu1, mu2);}
 153   inline double dgpdepsi(double mu1, double mu2){return 1./8./EvtConst::pi*(alphas(mu2)-alphas(mu1))*(gammap1()/beta0()-beta1()*gammap0()/pow(beta0(),2));}
 154   inline double agammap(double mu1, double mu2, double epsi=0){return gammap0()/2/beta0()*log(alphas(mu2)/alphas(mu1))+epsi*dgpdepsi( mu1, mu2);}
 155   inline double U1(double mu1, double mu2, double epsi=0){return exp(2*(Sudakov(mu1,mu2,epsi)-agammap(mu1,mu2,epsi)-aGamma(mu1,mu2,epsi)*log(_mb/mu1)));}
 156   inline double U1lo(double mu1, double mu2){return U1(mu1,mu2);}
 157   inline double U1nlo(double mu1, double mu2){return U1(mu1,mu2)*(1+2*(dSudakovdepsi(mu1,mu2)-dgpdepsi( mu1, mu2)-log(_mb/mu1)*dGdepsi( mu1, mu2)));}
 158   inline double alo(double mu1, double mu2){return -2*aGamma(mu1,mu2);}
 159   inline double anlo(double mu1, double mu2){return -2*dGdepsi(mu1,mu2);}
 160
 161 };
 162
 163 #endif
 164