Line data Source code
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 0 : 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 0 : 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 0 : inline double mu_i(){return _mui;} // intermediate scale
92 : inline double mu_bar(){return _mui;}
93 0 : inline double mu_h(){return _mb/sqrt(2.0);} // high scale
94 0 : inline double lambda1(){return -_mupi2;}
95 :
96 : // expansion coefficients for RGE
97 0 : static double beta0(int nf=4){return 11.-2./3.*nf;}
98 0 : static double beta1(int nf=4){return 34.*3.-38./3.*nf;}
99 0 : static double beta2(int nf=4){return 1428.5-5033./18.*nf+325./54.*nf*nf;}
100 0 : static double gamma0(){return 16./3.;}
101 0 : static double gamma1(int nf=4){return 4./3.*(49.85498-40./9.*nf);}
102 0 : static double gamma2(int nf=4){return 64./3.*(55.07242-8.58691*nf-nf*nf/27.);} /* zeta3=1.20206 */
103 0 : static double gammap0(){return -20./3.;}
104 0 : 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 0 : static double C_F(double mu){return (4.0/3.0)*alphas(mu)/4./EvtConst::pi;}
111 :
112 : // Shape Functions
113 :
114 0 : inline double lambda_SF(){ return _lambdaSF;}
115 : double lambda_bar(double omega0);
116 0 : inline double lambda2(){return 0.12;}
117 : double mu_pi2(double omega0);
118 : inline double lambda(double){ return _mB-_mb;}
119 :
120 : // specail for gaussian SF
121 0 : 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 0 : inline double S0(double a, double r){return -gamma0()/4/a/pow(beta0(),2)*(1/r-1+log(r));}
137 0 : inline double S1(double /*a*/, double r){return gamma0()/4./pow(beta0(),2)*(
138 0 : pow(log(r),2)*beta1()/2./beta0()+(gamma1()/gamma0()-beta1()/beta0())*(1.-r+log(r))
139 : );}
140 0 : inline double S2(double a, double r){return gamma0()*a/4./pow(beta0(),2)*(
141 0 : -0.5*pow((1-r),2)*(
142 0 : pow(beta1()/beta0(),2)-beta2()/beta0()-beta1()/beta0()*gamma1()/gamma0()+gamma2()/gamma0()
143 : )
144 0 : +(pow(beta1()/beta0(),2)-beta2()/beta0())*(1-r)*log(r)
145 0 : +(beta1()/beta0()*gamma1()/gamma0()-beta2()/beta0())*(1-r+r*log(r))
146 : );}
147 0 : inline double dSudakovdepsi(double mu1, double mu2){return S2(alphas(mu1)/(4*EvtConst::pi),alphas(mu2)/alphas(mu1));}
148 0 : 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 0 : inline double dGdepsi(double mu1, double mu2){return 1./8./EvtConst::pi*(alphas(mu2)-alphas(mu1))*(gamma1()/beta0()-beta1()*gamma0()/pow(beta0(),2));}
152 0 : inline double aGamma(double mu1, double mu2, double epsi=0){return gamma0()/2/beta0()*log(alphas(mu2)/alphas(mu1))+epsi*dGdepsi( mu1, mu2);}
153 0 : inline double dgpdepsi(double mu1, double mu2){return 1./8./EvtConst::pi*(alphas(mu2)-alphas(mu1))*(gammap1()/beta0()-beta1()*gammap0()/pow(beta0(),2));}
154 0 : inline double agammap(double mu1, double mu2, double epsi=0){return gammap0()/2/beta0()*log(alphas(mu2)/alphas(mu1))+epsi*dgpdepsi( mu1, mu2);}
155 0 : 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 0 : inline double U1lo(double mu1, double mu2){return U1(mu1,mu2);}
157 0 : 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 0 : inline double alo(double mu1, double mu2){return -2*aGamma(mu1,mu2);}
159 0 : inline double anlo(double mu1, double mu2){return -2*dGdepsi(mu1,mu2);}
160 :
161 : };
162 :
163 : #endif
164 :
|
