Line data Source code
1 : //---------------------------------------------------------------------------------
2 : //
3 : // Wilson coeficients according to A.J.Buras and M.Munz, Phys.Rev. D52, 186. (1995)
4 : // Thanks to N. Nikitine for example code for Pythia
5 : // Coefficient C8eff and C2 correction to C7eff taken from:
6 : // A.J.Buras, M.Misiak, M.Munz, S.Pokorski, Nucl.Phys. B424, 374 (1994)
7 : //
8 : // Used constants come from PDG 2004
9 : //
10 : // P. Reznicek 18.02.2005
11 : //
12 : // 04/03/2005 PR Added h-function
13 : //
14 : //---------------------------------------------------------------------------------
15 :
16 : #include "EvtGenBase/EvtPatches.hh"
17 : #include "EvtGenBase/EvtConst.hh"
18 : #include "EvtGenBase/EvtReport.hh"
19 : #include "EvtGenModels/EvtWilsonCoefficients.hh"
20 : #include <stdlib.h>
21 :
22 : // including EvtLi2Spence.F
23 : extern "C" {
24 : extern double li2spence_(double*);
25 : }
26 :
27 0 : EvtWilsonCoefficients::EvtWilsonCoefficients(){
28 : int i,j;
29 0 : double tmpa[8]={14./23.,16./23.,6./23.,-12./23.,0.4086,-0.4230,-0.8994,0.1456};
30 0 : double tmph[8]={2.2996,-1.0880,-3./7.,-1./14.,-0.6494,-0.0380,-0.0186,-0.0057};
31 0 : double tmpp[8]={0,0,-80./203.,8./33.,0.0433,0.1384,0.1648,-0.0073};
32 0 : double tmps[8]={0,0,-0.2009,-0.3579,0.0490,-0.3616,-0.3554,0.0072};
33 0 : double tmpq[8]={0,0,0,0,0.0318,0.0918,-0.2700,0.0059};
34 0 : double tmpg[8]={313063./363036.,0,0,0,-0.9135,0.0873,-0.0571,-0.0209};
35 0 : double tmpk[6][8]={{0,0,1./2.,-1./2.,0,0,0,0},
36 : {0,0,1./2.,+1./2.,0,0,0,0},
37 : {0,0,-1./14.,+1./6.,0.0510,-0.1403,-0.0113,0.0054},
38 : {0,0,-1./14.,-1./6.,0.0984,+0.1214,+0.0156,0.0026},
39 : {0,0,0,0,-0.0397,0.0117,-0.0025,+0.0304},
40 : {0,0,0,0,+0.0335,0.0239,-0.0462,-0.0112}};
41 0 : double tmpr[2][8]={{0,0,+0.8966,-0.1960,-0.2011,0.1328,-0.0292,-0.1858},
42 : {0,0,-0.1193,+0.1003,-0.0473,0.2323,-0.0133,-0.1799}};
43 0 : for(i=0;i<8;i++){
44 0 : a[i]=tmpa[i];
45 0 : h[i]=tmph[i];
46 0 : p[i]=tmpp[i];
47 0 : s[i]=tmps[i];
48 0 : q[i]=tmpq[i];
49 0 : g[i]=tmpg[i];
50 0 : for(j=0;j<6;j++) k[j][i]=tmpk[j][i];
51 0 : for(j=0;j<2;j++) r[j][i]=tmpr[j][i];
52 : }
53 0 : m_n_f=5;
54 0 : m_Lambda=0.2167;
55 0 : m_alphaMZ=0.1187;
56 0 : m_mu=4.8;
57 0 : m_M_Z=91.1876;
58 0 : m_M_t=174.3;
59 0 : m_M_W=80.425;
60 0 : m_alphaS=0;
61 0 : m_eta=0;
62 0 : m_sin2W=0.23120;
63 0 : m_ialpha=137.036;
64 0 : m_C1=m_C2=m_C3=m_C4=m_C5=m_C6=m_C7=m_C7eff0=m_C8=m_C8eff0=m_C9=m_C9tilda=m_C10=m_C10tilda=m_P0=0;
65 0 : m_A=m_B=m_C=m_D=m_E=m_F=m_Y=m_Z=m_PE=0;
66 0 : m_ksi=0;
67 0 : }
68 :
69 : void EvtWilsonCoefficients::SetRenormalizationScheme(std::string scheme){
70 0 : if(scheme=="NDR") m_ksi=0;
71 0 : else if(scheme=="HV") m_ksi=1;
72 : else{
73 0 : report(Severity::Error,"EvtGen") << "ERROR: EvtWilsonCoefficients knows only NDR and HV schemes !" << std::endl;
74 0 : ::abort();
75 : }
76 0 : }
77 :
78 : double EvtWilsonCoefficients::alphaS(double mu=4.8,int n_f=5,double Lambda=0.2167){
79 : // calculate strong coupling constant for n_f flavours and scale mu
80 0 : double beta0=11.-2./3.*n_f;
81 0 : double beta1=51.-19./3.*n_f;
82 0 : double beta2=2857.-5033./9.*n_f+325./27.*n_f*n_f;
83 0 : double lnratio=log(mu*mu/Lambda/Lambda);
84 0 : double aS=4.*EvtConst::pi/beta0/lnratio*(1.-2*beta1/beta0/beta0*log(lnratio)/lnratio+
85 0 : 4*beta1*beta1/beta0/beta0/beta0/beta0/lnratio/lnratio*((log(lnratio)-0.5)*(log(lnratio)-0.5)+beta2*beta0/8/beta1/beta1-5./4.));
86 0 : return aS;
87 : }
88 :
89 : double EvtWilsonCoefficients::Lambda(double alpha=0.1187,int n_f=5,double mu=91.1876,double epsilon=0.00005,int maxstep=1000){
90 : // calculate Lambda matching alphaS using simple iterative method
91 : int i;
92 : double difference=0;
93 0 : double Lambda=mu*0.9999999999;
94 0 : double step=-mu/20;
95 0 : for(i=0;i<maxstep && (difference=fabs(alphaS(mu,n_f,Lambda)-alpha))>=epsilon;i++){
96 0 : report(Severity::Debug,"EvtGen") << " Difference of alpha_S from " << alpha << " is " << difference << " at Lambda = " << Lambda << std::endl;
97 0 : if(alphaS(mu,n_f,Lambda)>alpha){
98 0 : if(step>0) step*=-0.4;
99 0 : if(alphaS(mu,n_f,Lambda+step-epsilon)<alphaS(mu,n_f,Lambda+step)) Lambda+=step;
100 0 : else step*=0.4;
101 : }else{
102 0 : if(step<0) step*=-0.4;
103 0 : if(Lambda+step<mu) Lambda+=step;
104 0 : else step*=0.4;
105 : }
106 : }
107 : report(Severity::Debug,"EvtGen") << " Difference of alpha_S from " << alpha << " is " << difference << " at Lambda = " << Lambda << std::endl;
108 0 : if(difference>=epsilon){
109 0 : report(Severity::Error,"EvtGen") << " ERROR: Did not converge Lambda for alpha_s = " << alpha << " , difference " << difference << " >= " << epsilon << " after " << i << " steps !" << std::endl;
110 0 : ::abort();
111 : return -1;
112 : }else{
113 0 : report(Severity::Info,"EvtGen") << " For alpha_s = " << alphaS(mu,n_f,Lambda) << " was found Lambda = " << Lambda << std::endl;
114 0 : return Lambda;
115 : }
116 : }
117 :
118 : double EvtWilsonCoefficients::eta(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
119 0 : return alphaS(M_W,n_f,Lambda)/alphaS(mu,n_f,Lambda);
120 : }
121 :
122 :
123 : EvtComplex EvtWilsonCoefficients::C1(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
124 : int i;
125 0 : EvtComplex myC1(0,0);
126 0 : for(i=0;i<8;i++) myC1+=k[0][i]*pow(eta(mu,n_f,Lambda,M_W),a[i]);
127 : return myC1;
128 0 : }
129 :
130 : EvtComplex EvtWilsonCoefficients::C2(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
131 : int i;
132 0 : EvtComplex myC2(0,0);
133 0 : for(i=0;i<8;i++) myC2+=k[1][i]*pow(eta(mu,n_f,Lambda,M_W),a[i]);
134 : return myC2;
135 0 : }
136 :
137 : EvtComplex EvtWilsonCoefficients::C3(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
138 : int i;
139 0 : EvtComplex myC3(0,0);
140 0 : for(i=0;i<8;i++) myC3+=k[2][i]*pow(eta(mu,n_f,Lambda,M_W),a[i]);
141 : return myC3;
142 0 : }
143 :
144 : EvtComplex EvtWilsonCoefficients::C4(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
145 : int i;
146 0 : EvtComplex myC4(0,0);
147 0 : for(i=0;i<8;i++) myC4+=k[3][i]*pow(eta(mu,n_f,Lambda,M_W),a[i]);
148 : return myC4;
149 0 : }
150 :
151 : EvtComplex EvtWilsonCoefficients::C5(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
152 : int i;
153 0 : EvtComplex myC5(0,0);
154 0 : for(i=0;i<8;i++) myC5+=k[4][i]*pow(eta(mu,n_f,Lambda,M_W),a[i]);
155 : return myC5;
156 0 : }
157 :
158 : EvtComplex EvtWilsonCoefficients::C6(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
159 : int i;
160 0 : EvtComplex myC6(0,0);
161 0 : for(i=0;i<8;i++) myC6+=k[5][i]*pow(eta(mu,n_f,Lambda,M_W),a[i]);
162 : return myC6;
163 0 : }
164 :
165 :
166 : EvtComplex EvtWilsonCoefficients::C7(double M_t=174.3,double M_W=80.425){
167 0 : return EvtComplex(-0.5*A(M_t*M_t/M_W/M_W),0);
168 : }
169 :
170 : EvtComplex EvtWilsonCoefficients::C8(double M_t=174.3,double M_W=80.425){
171 0 : return EvtComplex(-0.5*F(M_t*M_t/M_W/M_W),0);
172 : }
173 :
174 : EvtComplex EvtWilsonCoefficients::C7eff0(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_t=174.3,double M_W=80.425){
175 : int i;
176 0 : EvtComplex myC7eff(0,0);
177 0 : for(i=0;i<8;i++) myC7eff+=h[i]*pow(eta(mu,n_f,Lambda,M_W),a[i]);
178 0 : myC7eff*=C2(mu,n_f,Lambda,M_W);
179 0 : myC7eff+=pow(eta(mu,n_f,Lambda,M_W),16./23.)*C7(M_t,M_W);
180 0 : myC7eff+=8./3.*(pow(eta(mu,n_f,Lambda,M_W),14./23.)-pow(eta(mu,n_f,Lambda,M_W),16./23.))*C8(M_t,M_W);
181 : return myC7eff;
182 0 : }
183 :
184 : EvtComplex EvtWilsonCoefficients::C8eff0(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_t=174.3,double M_W=80.425){
185 : int i;
186 0 : EvtComplex myC8eff(0,0);
187 0 : for(i=0;i<8;i++) myC8eff+=g[i]*pow(eta(mu,n_f,Lambda,M_W),a[i]);
188 0 : myC8eff+=pow(eta(mu,n_f,Lambda,M_W),14./23.)*C8(M_t,M_W);
189 : return myC8eff;
190 0 : }
191 :
192 :
193 : EvtComplex EvtWilsonCoefficients::C10tilda(double sin2W=0.23120,double M_t=174.3,double M_W=80.425){
194 0 : return EvtComplex(-Y(M_t*M_t/M_W/M_W)/sin2W,0);
195 : }
196 :
197 : EvtComplex EvtWilsonCoefficients::C10(double sin2W=0.23120,double M_t=174.3,double M_W=80.425,double ialpha=137.036){
198 0 : return ( 1./2/EvtConst::pi/ialpha*C10tilda(sin2W,M_t,M_W) );
199 : }
200 :
201 :
202 : double EvtWilsonCoefficients::A(double x){
203 0 : return ( x*(8*x*x+5*x-7)/12/pow(x-1,3) + x*x*(2-3*x)*log(x)/2/pow(x-1,4) );
204 : }
205 :
206 : double EvtWilsonCoefficients::B(double x){
207 0 : return ( x/4/(1-x) + x/4/(x-1)/(x-1)*log(x) );
208 : }
209 :
210 : double EvtWilsonCoefficients::C(double x){
211 0 : return ( x*(x-6)/8/(x-1) + x*(3*x+2)/8/(x-1)/(x-1)*log(x) );
212 : }
213 :
214 : double EvtWilsonCoefficients::D(double x){
215 0 : return ( (-19*x*x*x+25*x*x)/36/pow(x-1,3) + x*x*(5*x*x-2*x-6)/18/pow(x-1,4)*log(x) - 4./9*log(x) );
216 : }
217 :
218 : double EvtWilsonCoefficients::E(double x){
219 0 : return ( x*(18-11*x-x*x)/12/pow(1-x,3) + x*x*(15-16*x+4*x*x)/6/pow(1-x,4)*log(x) - 2./3*log(x) );
220 : }
221 :
222 : double EvtWilsonCoefficients::F(double x){
223 0 : return ( x*(x*x-5*x-2)/4/pow(x-1,3) + 3*x*x/2/pow(x-1,4)*log(x) );
224 : }
225 :
226 : double EvtWilsonCoefficients::Y(double x){
227 0 : return (C(x)-B(x));
228 : }
229 :
230 : double EvtWilsonCoefficients::Z(double x){
231 0 : return (C(x)+1./4*D(x));
232 : }
233 :
234 :
235 : EvtComplex EvtWilsonCoefficients::C9(int ksi=0,double mu=4.8,int n_f=5,double Lambda=0.2167,double sin2W=0.23120,double M_t=174.3,double M_W=80.425,double ialpha=137.036){
236 0 : return ( 1./2/EvtConst::pi/ialpha*C9tilda(ksi,mu,n_f,Lambda,sin2W,M_t,M_W) );
237 : }
238 :
239 : EvtComplex EvtWilsonCoefficients::C9tilda(int ksi=0,double mu=4.8,int n_f=5,double Lambda=0.2167,double sin2W=0.23120,double M_t=174.3,double M_W=80.425){
240 0 : return ( P0(ksi,mu,n_f,Lambda,M_W) + Y(M_t*M_t/M_W/M_W)/sin2W - 4*Z(M_t*M_t/M_W/M_W) + PE(mu,n_f,Lambda,M_W)*E(M_t*M_t/M_W/M_W) );
241 : }
242 :
243 : EvtComplex EvtWilsonCoefficients::P0(int ksi=0,double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
244 : int i;
245 0 : EvtComplex myP0(0,0);
246 0 : for(i=0;i<8;i++) myP0+=p[i]*pow(eta(mu,n_f,Lambda,M_W),a[i]+1);
247 0 : myP0=EvtConst::pi/alphaS(M_W,n_f,Lambda)*(-0.1875+myP0);
248 0 : myP0+=1.2468-ksi*4./9.*(3*C1(mu,n_f,Lambda,M_W)+C2(mu,n_f,Lambda,M_W)-C3(mu,n_f,Lambda,M_W)-3*C4(mu,n_f,Lambda,M_W));
249 0 : for(i=0;i<8;i++) myP0+=pow(eta(mu,n_f,Lambda,M_W),a[i])*(r[ksi][i]+s[i]*eta(mu,n_f,Lambda,M_W));
250 : return myP0;
251 0 : }
252 :
253 : double EvtWilsonCoefficients::PE(double mu=4.8,int n_f=5,double Lambda=0.2167,double M_W=80.425){
254 : int i;
255 : double myPE=0.1405;
256 0 : for(i=0;i<8;i++) myPE+=q[i]*pow(eta(mu,n_f,Lambda,M_W),a[i]+1);
257 0 : return myPE;
258 : }
259 :
260 :
261 : void EvtWilsonCoefficients::CalculateAllCoefficients(){
262 0 : m_Lambda=Lambda(m_alphaMZ,m_n_f,m_M_Z);
263 0 : m_C1=C1(m_mu,m_n_f,m_Lambda,m_M_W);
264 0 : m_C2=C2(m_mu,m_n_f,m_Lambda,m_M_W);
265 0 : m_C3=C3(m_mu,m_n_f,m_Lambda,m_M_W);
266 0 : m_C4=C4(m_mu,m_n_f,m_Lambda,m_M_W);
267 0 : m_C5=C5(m_mu,m_n_f,m_Lambda,m_M_W);
268 0 : m_C6=C6(m_mu,m_n_f,m_Lambda,m_M_W);
269 0 : m_C7=C7(m_M_t,m_M_W);
270 0 : m_C8=C8(m_M_t,m_M_W);
271 0 : m_C7eff0=C7eff0(m_mu,m_n_f,m_Lambda,m_M_t,m_M_W);
272 0 : m_C8eff0=C8eff0(m_mu,m_n_f,m_Lambda,m_M_t,m_M_W);
273 0 : m_C10tilda=C10tilda(m_sin2W,m_M_t,m_M_W);
274 0 : m_C10=C10(m_sin2W,m_M_t,m_M_W,m_ialpha);
275 0 : m_A=A(m_M_t*m_M_t/m_M_W/m_M_W);
276 0 : m_B=B(m_M_t*m_M_t/m_M_W/m_M_W);
277 0 : m_C=C(m_M_t*m_M_t/m_M_W/m_M_W);
278 0 : m_D=D(m_M_t*m_M_t/m_M_W/m_M_W);
279 0 : m_E=E(m_M_t*m_M_t/m_M_W/m_M_W);
280 0 : m_F=F(m_M_t*m_M_t/m_M_W/m_M_W);
281 0 : m_Y=Y(m_M_t*m_M_t/m_M_W/m_M_W);
282 0 : m_Z=Z(m_M_t*m_M_t/m_M_W/m_M_W);
283 0 : m_C9=C9(m_ksi,m_mu,m_n_f,m_Lambda,m_sin2W,m_M_t,m_M_W,m_ialpha);
284 0 : m_C9tilda=C9tilda(m_ksi,m_mu,m_n_f,m_Lambda,m_sin2W,m_M_t,m_M_W);
285 0 : m_P0=P0(m_ksi,m_mu,m_n_f,m_Lambda,m_M_W);
286 0 : m_PE=PE(m_mu,m_n_f,m_Lambda,m_M_W);
287 0 : m_alphaS=alphaS(m_mu,m_n_f,m_Lambda);
288 0 : m_eta=eta(m_mu,m_n_f,m_Lambda,m_M_W);
289 0 : report(Severity::Info,"EvtGen") << " +---------------------------------------" << std::endl;
290 0 : report(Severity::Info,"EvtGen") << " | Table of Wilson coeficients:" << std::endl;
291 0 : report(Severity::Info,"EvtGen") << " +---------------------------------------" << std::endl;
292 0 : report(Severity::Info,"EvtGen") << " | C1 = " << m_C1 << std::endl;
293 0 : report(Severity::Info,"EvtGen") << " | C2 = " << m_C2 << std::endl;
294 0 : report(Severity::Info,"EvtGen") << " | C3 = " << m_C3 << std::endl;
295 0 : report(Severity::Info,"EvtGen") << " | C4 = " << m_C4 << std::endl;
296 0 : report(Severity::Info,"EvtGen") << " | C5 = " << m_C5 << std::endl;
297 0 : report(Severity::Info,"EvtGen") << " | C6 = " << m_C6 << std::endl;
298 0 : report(Severity::Info,"EvtGen") << " | C7 = " << m_C7 << std::endl;
299 0 : report(Severity::Info,"EvtGen") << " | C7eff0 = " << m_C7eff0 << std::endl;
300 0 : report(Severity::Info,"EvtGen") << " | C8 = " << m_C8 << std::endl;
301 0 : report(Severity::Info,"EvtGen") << " | C8eff0 = " << m_C8eff0 << std::endl;
302 0 : report(Severity::Info,"EvtGen") << " | C9 = " << m_C9 << std::endl;
303 0 : report(Severity::Info,"EvtGen") << " | C10 = " << m_C10 << std::endl;
304 0 : report(Severity::Info,"EvtGen") << " +---------------------------------------" << std::endl;
305 0 : report(Severity::Info,"EvtGen") << " | Other constants:" << std::endl;
306 0 : report(Severity::Info,"EvtGen") << " +---------------------------------------" << std::endl;
307 0 : report(Severity::Info,"EvtGen") << " | Scale = " << m_mu << " GeV" << std::endl;
308 0 : report(Severity::Info,"EvtGen") << " | Number of effective flavors = " << m_n_f << std::endl;
309 0 : report(Severity::Info,"EvtGen") << " | Corresponding to aS(M_Z)" << "=" << m_alphaMZ << " Lambda = " << m_Lambda << " GeV" << std::endl;
310 0 : report(Severity::Info,"EvtGen") << " | Strong coupling constant = " << m_alphaS << std::endl;
311 0 : report(Severity::Info,"EvtGen") << " | Electromagnetic constant = 1/" << m_ialpha << std::endl;
312 0 : report(Severity::Info,"EvtGen") << " | Top mass = " << m_M_t << " GeV" << std::endl;
313 0 : report(Severity::Info,"EvtGen") << " | W-boson mass = " << m_M_W << " GeV" << std::endl;
314 0 : report(Severity::Info,"EvtGen") << " | Z-boson mass = " << m_M_Z << " GeV" << std::endl;
315 0 : report(Severity::Info,"EvtGen") << " | Sinus squared of Weinberg angle = " << m_sin2W << std::endl;
316 0 : report(Severity::Info,"EvtGen") << " +---------------------------------------" << std::endl;
317 0 : report(Severity::Debug,"EvtGen") << " | Intermediate functions:" << std::endl;
318 0 : report(Severity::Debug,"EvtGen") << " +---------------------------------------" << std::endl;
319 0 : report(Severity::Debug,"EvtGen") << " | A = " << m_A << std::endl;
320 0 : report(Severity::Debug,"EvtGen") << " | B = " << m_B << std::endl;
321 0 : report(Severity::Debug,"EvtGen") << " | C = " << m_C << std::endl;
322 0 : report(Severity::Debug,"EvtGen") << " | D = " << m_D << std::endl;
323 0 : report(Severity::Debug,"EvtGen") << " | E = " << m_E << std::endl;
324 0 : report(Severity::Debug,"EvtGen") << " | F = " << m_F << std::endl;
325 0 : report(Severity::Debug,"EvtGen") << " | Y = " << m_Y << std::endl;
326 0 : report(Severity::Debug,"EvtGen") << " | Z = " << m_Z << std::endl;
327 0 : report(Severity::Debug,"EvtGen") << " | eta = " << m_eta << std::endl;
328 0 : report(Severity::Debug,"EvtGen") << " | C9~ = " << m_C9tilda << std::endl;
329 0 : report(Severity::Debug,"EvtGen") << " | C10~ = " << m_C10tilda << std::endl;
330 0 : report(Severity::Debug,"EvtGen") << " | P0 = " << m_P0 << std::endl;
331 0 : report(Severity::Debug,"EvtGen") << " | PE = " << m_PE << std::endl;
332 0 : report(Severity::Debug,"EvtGen") << " +--------------------------------------" << std::endl;
333 0 : }
334 :
335 : EvtComplex EvtWilsonCoefficients::hzs(double z,double shat,double mu=4.8,double M_b=4.8){
336 0 : EvtComplex i1(0,1);
337 0 : double x=4.*z*z/shat;
338 0 : if(x==0) return (8./27. - 8./9.*log(M_b/mu) - 4./9.*log(shat) + 4./9.*i1*EvtConst::pi);
339 0 : else if(x>1) return (8./27. - 8./9.*log(M_b/mu) - 8./9.*log(z) + 4./9.*x - 2./9.*(2.+x)*sqrt(x-1.) * 2*atan(1./sqrt(x-1.)));
340 0 : else return (8./27. - 8./9.*log(M_b/mu) - 8./9.*log(z) + 4./9.*x - 2./9.*(2.+x)*sqrt(1.-x) * (log(fabs(sqrt(1.-x)+1)/fabs(sqrt(1.-x)-1))-i1*EvtConst::pi));
341 0 : }
342 :
343 : double EvtWilsonCoefficients::fz(double z){
344 0 : return (1. - 8.*z*z + 8.*pow(z,6.) - pow(z,8.) - 24.*pow(z,4.)*log(z));
345 : }
346 :
347 : double EvtWilsonCoefficients::kappa(double z,double alpha_S){
348 0 : return (1. - 2.*alpha_S/3./EvtConst::pi*((EvtConst::pi*EvtConst::pi-31./4.)*(1.-z)*(1.-z) + 1.5) );
349 : }
350 :
351 : double EvtWilsonCoefficients::etatilda(double shat,double alpha_S){
352 0 : return (1. + alpha_S/EvtConst::pi*omega(shat));
353 : }
354 :
355 : double EvtWilsonCoefficients::omega(double shat){
356 : double o=0;
357 0 : o -= (2./9.)*EvtConst::pi*EvtConst::pi;
358 0 : o -= (4./3.)*li2spence_(&shat);
359 0 : o -= (2./3.)*log(shat)*log(1.-shat);
360 0 : o -= log(1.-shat)*(5.+4.*shat)/(3.+6.*shat);
361 0 : o -= log(shat)*2.*shat*(1.+shat)*(1.-2.*shat)/3./(1.-shat)/(1.-shat)/(1.+2.*shat);
362 0 : o += (5.+9.*shat-6.*shat*shat)/6./(1.-shat)/(1.+2.*shat);
363 0 : return o;
364 : }
365 :
366 : EvtComplex EvtWilsonCoefficients::C9efftilda(double z,double shat,double alpha_S,EvtComplex c1,EvtComplex c2,EvtComplex c3,EvtComplex c4,EvtComplex c5,EvtComplex c6,EvtComplex c9tilda,int ksi=0){
367 0 : EvtComplex c(0,0);
368 0 : c += (c9tilda+ksi*4./9.*(3.*c1+c2-c3-3.*c4))*etatilda(shat,alpha_S);
369 0 : c += hzs(z,shat)*(3.*c1+c2+3.*c3+c4+3.*c5+c6);
370 0 : c -= 0.5*hzs(1,shat)*(4.*c3+4.*c4+3.*c5+c6);
371 0 : c -= 0.5*hzs(0,shat)*(c3+3.*c4);
372 0 : c += 2./9.*(3.*c3+c4+3.*c5+c6);
373 0 : return c;
374 : }
375 :
376 : EvtComplex EvtWilsonCoefficients::C7b2sg(double alpha_S,double et,EvtComplex c2,double M_t=174.3,double M_W=80.425){
377 0 : EvtComplex i1(0,1);
378 0 : return (i1*alpha_S*(2./9.*pow(et,14./23.)*(0.5*F(M_t*M_t/M_W/M_W)-0.1687)-0.03*c2));
379 0 : }
380 :
381 : EvtComplex EvtWilsonCoefficients::Yld(double q2,double *ki,double *Gi,double *Mi,int ni,EvtComplex c1,EvtComplex c2,EvtComplex c3,EvtComplex c4,EvtComplex c5,EvtComplex c6,double ialpha=137.036){
382 0 : EvtComplex i1(0,1);
383 0 : EvtComplex y(0,0);
384 : int i;
385 0 : for(i=0;i<ni;i++) y+=ki[i]*Gi[i]*Mi[i]/(q2-Mi[i]*Mi[i]-i1*Mi[i]*Gi[i]);
386 0 : return (-3.*ialpha*ialpha*y*EvtConst::pi*(3.*c1+c2+3.*c3+c4+3.*c5+c6));
387 0 : }
|
