CGEM BOSS 6.6.5.g
BESIII Offline Software System
Loading...
Searching...
No Matches
EvtbTosllAmp Class Referenceabstract

#include <EvtbTosllAmp.hh>

+ Inheritance diagram for EvtbTosllAmp:

Public Member Functions

virtual void CalcAmp (EvtParticle *parent, EvtAmp &amp, EvtbTosllFF *formFactors)=0
 
double CalcMaxProb (EvtId parent, EvtId meson, EvtId lepton, EvtId nudaug, EvtbTosllFF *formFactors, double &poleSize)
 
EvtComplex GetC7Eff (double q2, bool nnlo=true)
 
EvtComplex GetC9Eff (double q2, bool nnlo=true, bool btod=false)
 
EvtComplex GetC10Eff (double q2, bool nnlo=true)
 
double dGdsProb (double mb, double ms, double ml, double s)
 
double dGdsdupProb (double mb, double ms, double ml, double s, double u)
 

Detailed Description

Definition at line 31 of file EvtbTosllAmp.hh.

Member Function Documentation

◆ CalcAmp()

virtual void EvtbTosllAmp::CalcAmp ( EvtParticle parent,
EvtAmp amp,
EvtbTosllFF formFactors 
)
pure virtual

◆ CalcMaxProb()

double EvtbTosllAmp::CalcMaxProb ( EvtId  parent,
EvtId  meson,
EvtId  lepton,
EvtId  nudaug,
EvtbTosllFF formFactors,
double &  poleSize 
)

Definition at line 37 of file EvtbTosllAmp.cc.

40 {
41
42 //This routine takes the arguements parent, meson, and lepton
43 //number, and a form factor model, and returns a maximum
44 //probability for this semileptonic form factor model. A
45 //brute force method is used. The 2D cos theta lepton and
46 //q2 phase space is probed.
47
48 //Start by declaring a particle at rest.
49
50 //It only makes sense to have a scalar parent. For now.
51 //This should be generalized later.
52
53 EvtScalarParticle *scalar_part;
54 EvtParticle *root_part;
55
56 scalar_part=new EvtScalarParticle;
57
58 //cludge to avoid generating random numbers!
59 scalar_part->noLifeTime();
60
61 EvtVector4R p_init;
62
63 p_init.set(EvtPDL::getMass(parent),0.0,0.0,0.0);
64 scalar_part->init(parent,p_init);
65 root_part=(EvtParticle *)scalar_part;
66 root_part->setDiagonalSpinDensity();
67
68 EvtParticle *daughter, *lep1, *lep2;
69
70 EvtAmp amp;
71
72 EvtId listdaug[3];
73 listdaug[0] = meson;
74 listdaug[1] = lepton1;
75 listdaug[2] = lepton2;
76
77 amp.init(parent,3,listdaug);
78
79 root_part->makeDaughters(3,listdaug);
80 daughter=root_part->getDaug(0);
81 lep1=root_part->getDaug(1);
82 lep2=root_part->getDaug(2);
83
84 //cludge to avoid generating random numbers!
85 daughter->noLifeTime();
86 lep1->noLifeTime();
87 lep2->noLifeTime();
88
89
90 //Initial particle is unpolarized, well it is a scalar so it is
91 //trivial
93 rho.SetDiag(root_part->getSpinStates());
94
95 double mass[3];
96
97 double m = root_part->mass();
98
99 EvtVector4R p4meson, p4lepton1, p4lepton2, p4w;
100 double q2max;
101
102 double q2, elepton, plepton;
103 int i,j;
104 double erho,prho,costl;
105
106 double maxfoundprob = 0.0;
107 double prob = -10.0;
108 int massiter;
109
110 double maxpole=0;
111
112 for (massiter=0;massiter<3;massiter++){
113
114 mass[0] = EvtPDL::getMeanMass(meson);
115 mass[1] = EvtPDL::getMeanMass(lepton1);
116 mass[2] = EvtPDL::getMeanMass(lepton2);
117 if ( massiter==1 ) {
118 mass[0] = EvtPDL::getMinMass(meson);
119 }
120 if ( massiter==2 ) {
121 mass[0] = EvtPDL::getMaxMass(meson);
122 if ( (mass[0]+mass[1]+mass[2])>m) mass[0]=m-mass[1]-mass[2]-0.00001;
123 }
124
125 q2max = (m-mass[0])*(m-mass[0]);
126
127 //loop over q2
128 //cout << "m " << m << "mass[0] " << mass[0] << " q2max "<< q2max << endl;
129 for (i=0;i<25;i++) {
130 //want to avoid picking up the tail of the photon propagator
131 q2 = ((i+1.5)*q2max)/26.0;
132
133 if (i==0) q2=4*(mass[1]*mass[1]);
134
135 erho = ( m*m + mass[0]*mass[0] - q2 )/(2.0*m);
136
137 prho = sqrt(erho*erho-mass[0]*mass[0]);
138
139 p4meson.set(erho,0.0,0.0,-1.0*prho);
140 p4w.set(m-erho,0.0,0.0,prho);
141
142 //This is in the W rest frame
143 elepton = (q2+mass[1]*mass[1])/(2.0*sqrt(q2));
144 plepton = sqrt(elepton*elepton-mass[1]*mass[1]);
145
146 double probctl[3];
147
148 for (j=0;j<3;j++) {
149
150 costl = 0.99*(j - 1.0);
151
152 //These are in the W rest frame. Need to boost out into
153 //the B frame.
154 p4lepton1.set(elepton,0.0,
155 plepton*sqrt(1.0-costl*costl),plepton*costl);
156 p4lepton2.set(elepton,0.0,
157 -1.0*plepton*sqrt(1.0-costl*costl),-1.0*plepton*costl);
158
159 EvtVector4R boost((m-erho),0.0,0.0,1.0*prho);
160 p4lepton1=boostTo(p4lepton1,boost);
161 p4lepton2=boostTo(p4lepton2,boost);
162
163 //Now initialize the daughters...
164
165 daughter->init(meson,p4meson);
166 lep1->init(lepton1,p4lepton1);
167 lep2->init(lepton2,p4lepton2);
168
169 CalcAmp(root_part,amp,FormFactors);
170
171 //Now find the probability at this q2 and cos theta lepton point
172 //and compare to maxfoundprob.
173
174 //Do a little magic to get the probability!!
175
176 //cout <<"amp:"<<amp.getSpinDensity()<<endl;
177
178 prob = rho.NormalizedProb(amp.getSpinDensity());
179
180 //cout << "prob:"<<q2<<" "<<costl<<" "<<prob<<endl;
181
182 probctl[j]=prob;
183 }
184
185 //probclt contains prob at ctl=-1,0,1.
186 //prob=a+b*ctl+c*ctl^2
187
188 double a=probctl[1];
189 double b=0.5*(probctl[2]-probctl[0]);
190 double c=0.5*(probctl[2]+probctl[0])-probctl[1];
191
192 prob=probctl[0];
193 if (probctl[1]>prob) prob=probctl[1];
194 if (probctl[2]>prob) prob=probctl[2];
195
196 if (fabs(c)>1e-20){
197 double ctlx=-0.5*b/c;
198 if (fabs(ctlx)<1.0){
199 double probtmp=a+b*ctlx+c*ctlx*ctlx;
200 if (probtmp>prob) prob=probtmp;
201 }
202
203 }
204
205 //report(DEBUG,"EvtGen") << "prob,probctl:"<<prob<<" "
206 // << probctl[0]<<" "
207 // << probctl[1]<<" "
208 // << probctl[2]<<endl;
209
210 if (i==0) {
211 maxpole=prob;
212 continue;
213 }
214
215 if ( prob > maxfoundprob ) {
216 maxfoundprob = prob;
217 }
218
219 //cout << "q2,maxfoundprob:"<<q2<<" "<<maxfoundprob<<endl;
220
221 }
222 if ( EvtPDL::getWidth(meson) <= 0.0 ) {
223 //if the particle is narrow dont bother with changing the mass.
224 massiter = 4;
225 }
226
227 }
228
229 root_part->deleteTree();
230
231 poleSize=0.04*(maxpole/maxfoundprob)*4*(mass[1]*mass[1]);
232
233 //poleSize=0.002;
234
235 //cout <<"maxfoundprob,maxpole,poleSize:"<<maxfoundprob<<" "
236 // <<maxpole<<" "<<poleSize<<endl;
237
238 maxfoundprob *=1.15;
239
240 return maxfoundprob;
241
242}
double mass
EvtDiracSpinor boostTo(const EvtDiracSpinor &sp, const EvtVector4R p4)
Definition: EvtAmp.hh:30
void init(EvtId p, int ndaug, EvtId *daug)
Definition: EvtAmp.cc:70
EvtSpinDensity getSpinDensity()
Definition: EvtAmp.cc:149
Definition: EvtId.hh:27
static double getWidth(EvtId i)
Definition: EvtPDL.hh:54
static double getMeanMass(EvtId i)
Definition: EvtPDL.hh:45
static double getMinMass(EvtId i)
Definition: EvtPDL.hh:51
static double getMaxMass(EvtId i)
Definition: EvtPDL.hh:50
static double getMass(EvtId i)
Definition: EvtPDL.hh:46
void noLifeTime()
Definition: EvtParticle.hh:348
void makeDaughters(int ndaug, EvtId *id)
virtual void init(EvtId part_n, const EvtVector4R &p4)=0
int getSpinStates() const
Definition: EvtParticle.cc:118
void setDiagonalSpinDensity()
Definition: EvtParticle.cc:133
EvtParticle * getDaug(int i)
Definition: EvtParticle.cc:85
double mass() const
Definition: EvtParticle.cc:127
void deleteTree()
Definition: EvtParticle.cc:555
void init(EvtId part_n, double e, double px, double py, double pz)
double NormalizedProb(const EvtSpinDensity &d)
void SetDiag(int n)
void set(int i, double d)
Definition: EvtVector4R.hh:183
virtual void CalcAmp(EvtParticle *parent, EvtAmp &amp, EvtbTosllFF *formFactors)=0

Referenced by EvtbTosllAli::initProbMax(), and EvtbTosllBall::initProbMax().

◆ dGdsdupProb()

double EvtbTosllAmp::dGdsdupProb ( double  mb,
double  ms,
double  ml,
double  s,
double  u 
)

Definition at line 608 of file EvtbTosllAmp.cc.

610{
611 // Compute the decay probability density function given a value of s and u
612 // according to Ali's paper
613
614 double prob;
615 double f1sp, f2sp, f3sp;
616
617 double sh = s / (mb*mb);
618
619 EvtComplex c9eff = EvtbTosllAmp::GetC9Eff(sh*mb);
620 EvtComplex c7eff = EvtbTosllAmp::GetC7Eff(sh*mb);
621 EvtComplex c10eff = EvtbTosllAmp::GetC10Eff(sh*mb);
622
623 double alphas = 0.119/
624 (1 + 0.119*log(pow(4.8,2)/pow(91.1867,2))*23.0/12.0/EvtConst::pi);
625 double omega9 = - 2.0/9.0*EvtConst::pi*EvtConst::pi - 4.0/3.0*ddilog_(sh)
626 - 2.0/3.0*log(sh)*log(1.0-sh)
627 - (5.0+4.0*sh)/(3.0*(1.0+2.0*sh)) * log(1.0-sh)
628 - 2.0*sh*(1.0+sh)*(1.0-2.0*sh)
629 /(3.0*pow(1.0-sh,2)*(1.0+2.0*sh)) * log(sh)
630 + (5.0+9.0*sh-6.0*sh*sh)/(6.0*(1.0-sh)*(1.0+2.0*sh));
631 double eta9 = 1.0 + alphas*omega9/EvtConst::pi;
632 double omega7 = -8.0/3.0*log(4.8/mb)
633 -4.0/3.0*ddilog_(sh)
635 -2.0/3.0*log(sh)*log(1.0-sh)
636 -log(1-sh)*(8.0+sh)/(2.0+sh)/3.0
637 -2.0/3.0*sh*(2.0 - 2.0*sh - sh*sh)*log(sh)/pow((1.0 - sh),2)/(2.0 + sh)
638 -(16.0 - 11.0*sh - 17.0*sh*sh)/18.0/(2.0 + sh)/(1.0 - sh);
639 double eta7 = 1.0 + alphas*omega7/EvtConst::pi;
640
641 double omega79 = -4.0/3.0*log(4.8/mb)
642 -4.0/3.0*ddilog_(sh)
644 -2.0/3.0*log(sh)*log(1.0-sh)
645 -1.0/9.0*(2.0+7.0*sh)*log(1.0 - sh)/sh
646 -2.0/9.0*sh*(3.0 - 2.0*sh)*log(sh)/pow((1.0 - sh),2)
647 +1.0/18.0*(5.0 - 9.0*sh)/(1.0 - sh);
648 double eta79 = 1.0 + alphas*omega79/EvtConst::pi;
649
650 double c7c9 = abs(c7eff)*real(c9eff);
651 c7c9 *= pow(eta79,2);
652 double c7c7 = pow(abs(c7eff),2);
653 c7c7 *= pow(eta7,2);
654
655 double c9c9plusc10c10 = pow(abs(c9eff),2) + pow(abs(c10eff),2);
656 c9c9plusc10c10 *= pow(eta9,2);
657 double c9c9minusc10c10 = pow(abs(c9eff),2) - pow(abs(c10eff),2);
658 c9c9minusc10c10 *= pow(eta9,2);
659 double c7c10 = abs(c7eff)*real(c10eff);
660 c7c10 *= eta7; c7c10 *= eta9;
661 double c9c10 = real(c9eff)*real(c10eff);
662 c9c10 *= pow(eta9,2);
663
664 f1sp = ( pow(mb*mb-ms*ms,2) - s*s) * c9c9plusc10c10
665 + 4.0*( pow(mb,4) - ms*ms*mb*mb - pow(ms,4)*(1.0 - ms*ms/(mb*mb))
666 - 8.0*s*ms*ms - s*s*(1.0 + ms*ms/(mb*mb) ))*mb*mb*c7c7/s
667 // kludged mass term
668 *(1.0 + 2.0*ml*ml/s)
669 - 8.0*(s*(mb*mb + ms*ms) - pow(mb*mb-ms*ms,2)) * c7c9
670 // kludged mass term
671 *(1.0 + 2.0*ml*ml/s);
672
673 f2sp = 4.0*s*c9c10 + 8.0*(mb*mb + ms*ms)*c7c10;
674 f3sp = - (c9c9plusc10c10)
675 + 4.0*(1.0 + pow(ms/mb,4)) * mb*mb*c7c7/s
676 // kludged mass term
677 *(1.0 + 2.0*ml*ml/s);
678
679 prob = (f1sp + f2sp*u + f3sp*u*u)/ pow(mb,3);
680
681 return prob;
682}
double real(const EvtComplex &c)
Definition: EvtComplex.hh:240
double abs(const EvtComplex &c)
Definition: EvtComplex.hh:212
double ddilog_(const double &sh)
XmlRpcServer s
Definition: HelloServer.cpp:11
static const double pi
Definition: EvtConst.hh:28
EvtComplex GetC7Eff(double q2, bool nnlo=true)
EvtComplex GetC10Eff(double q2, bool nnlo=true)
EvtComplex GetC9Eff(double q2, bool nnlo=true, bool btod=false)

◆ dGdsProb()

double EvtbTosllAmp::dGdsProb ( double  mb,
double  ms,
double  ml,
double  s 
)

Definition at line 534 of file EvtbTosllAmp.cc.

536{
537 // Compute the decay probability density function given a value of s
538 // according to Ali's paper
539
540
541 double delta, lambda, prob;
542 double f1, f2, f3, f4;
543 double msh, mlh, sh;
544
545 mlh = ml / mb;
546 msh = ms / mb;
547 sh = s / (mb*mb);
548
549 EvtComplex c9eff = EvtbTosllAmp::GetC9Eff(sh*mb);
550 EvtComplex c7eff = EvtbTosllAmp::GetC7Eff(sh*mb);
551 EvtComplex c10eff = EvtbTosllAmp::GetC10Eff(sh*mb);
552
553 double alphas = 0.119/
554 (1 + 0.119*log(pow(4.8,2)/pow(91.1867,2))*23.0/12.0/EvtConst::pi);
555 double omega9 = -2.0/9.0*EvtConst::pi*EvtConst::pi - 4.0/3.0*ddilog_(sh)
556 - 2.0/3.0*log(sh)*log(1.0-sh)
557 - (5.0+4.0*sh)/(3.0*(1.0+2.0*sh)) * log(1.0-sh)
558 - 2.0*sh*(1.0+sh)*(1.0-2.0*sh)
559 /(3.0*pow(1.0-sh,2)*(1.0+2.0*sh)) * log(sh)
560 + (5.0+9.0*sh-6.0*sh*sh)/(6.0*(1.0-sh)*(1.0+2.0*sh));
561 double eta9 = 1.0 + alphas*omega9/EvtConst::pi;
562 double omega7 = -8.0/3.0*log(4.8/mb)
563 -4.0/3.0*ddilog_(sh)
565 -2.0/3.0*log(sh)*log(1.0-sh)
566 -log(1-sh)*(8.0+sh)/(2.0+sh)/3.0
567 -2.0/3.0*sh*(2.0 - 2.0*sh - sh*sh)*log(sh)/pow((1.0 - sh),2)/(2.0 + sh)
568 -(16.0 - 11.0*sh - 17.0*sh*sh)/18.0/(2.0 + sh)/(1.0 - sh);
569 double eta7 = 1.0 + alphas*omega7/EvtConst::pi;
570
571 double omega79 = -4.0/3.0*log(4.8/mb)
572 -4.0/3.0*ddilog_(sh)
574 -2.0/3.0*log(sh)*log(1.0-sh)
575 -1.0/9.0*(2.0+7.0*sh)*log(1.0 - sh)/sh
576 -2.0/9.0*sh*(3.0 - 2.0*sh)*log(sh)/pow((1.0 - sh),2)
577 +1.0/18.0*(5.0 - 9.0*sh)/(1.0 - sh);
578 double eta79 = 1.0 + alphas*omega79/EvtConst::pi;
579
580 double c7c9 = abs(c7eff)*real(c9eff);
581 c7c9 *= pow(eta79,2);
582 double c7c7 = pow(abs(c7eff),2);
583 c7c7 *= pow(eta7,2);
584
585 double c9c9plusc10c10 = pow(abs(c9eff),2) + pow(abs(c10eff),2);
586 c9c9plusc10c10 *= pow(eta9,2);
587 double c9c9minusc10c10 = pow(abs(c9eff),2) - pow(abs(c10eff),2);
588 c9c9minusc10c10 *= pow(eta9,2);
589
590 lambda = 1.0 + sh*sh + pow(msh,4) - 2.0*(sh + sh*msh*msh + msh*msh);
591
592 f1 = pow(1.0-msh*msh,2) - sh*(1.0 + msh*msh);
593 f2 = 2.0*(1.0 + msh*msh) * pow(1.0-msh*msh,2)
594 - sh*(1.0 + 14.0*msh*msh + pow(msh,4)) - sh*sh*(1.0 + msh*msh);
595 f3 = pow(1.0-msh*msh,2) + sh*(1.0 + msh*msh) - 2.0*sh*sh
596 + lambda*2.0*mlh*mlh/sh;
597 f4 = 1.0 - sh + msh*msh;
598
599 delta = ( 12.0*c7c9*f1 + 4.0*c7c7*f2/sh ) * (1.0 + 2.0*mlh*mlh/sh)
600 + c9c9plusc10c10*f3
601 + 6.0*mlh*mlh*c9c9minusc10c10*f4;
602
603 prob = sqrt(lambda*(1.0 - 4.0*mlh*mlh/sh)) * delta;
604
605 return prob;
606}
TFile * f1

◆ GetC10Eff()

EvtComplex EvtbTosllAmp::GetC10Eff ( double  q2,
bool  nnlo = true 
)

Definition at line 521 of file EvtbTosllAmp.cc.

522{
523
524 if (!nnlo) return -4.669;
525 double A10;
526 A10 = -4.592 + 0.379;
527
528 EvtComplex c10eff;
529 c10eff = A10;
530
531 return c10eff;
532}

Referenced by EvtbTosllScalarAmp::CalcAmp(), EvtbTosllVectorAmp::CalcAmp(), dGdsdupProb(), and dGdsProb().

◆ GetC7Eff()

EvtComplex EvtbTosllAmp::GetC7Eff ( double  q2,
bool  nnlo = true 
)

Definition at line 245 of file EvtbTosllAmp.cc.

246{
247
248 if (!nnlo) return -0.313;
249 double mbeff = 4.8;
250 double shat = q2/mbeff/mbeff;
251 double logshat;
252 logshat = log(shat);
253
254 double muscale;
255 muscale = 2.5;
256 double alphas;
257 alphas = 0.267;
258 double A7;
259 A7 = -0.353 + 0.023;
260 double A8;
261 A8 = -0.164;
262 double A9;
263 A9 = 4.287 + (-0.218);
264 double A10;
265 A10 = -4.592 + 0.379;
266 double C1;
267 C1 = -0.697;
268 double C2;
269 C2 = 1.046;
270 double T9;
271 T9 = 0.114 + 0.280;
272 double U9;
273 U9 = 0.045 + 0.023;
274 double W9;
275 W9 = 0.044 + 0.016;
276
277 double Lmu;
278 Lmu = log(muscale/mbeff);
279
280 EvtComplex uniti(0.0,1.0);
281
282 EvtComplex c7eff;
283 if (shat > 0.25)
284 {
285 c7eff = A7;
286 return c7eff;
287 }
288
289
290
291
292 // change energy scale to 5.0 for full NNLO calculation below shat = 0.25
293 muscale = 5.0;
294 alphas = 0.215;
295 A7 = -0.312 + 0.008;
296 A8 = -0.148;
297 A9 = 4.174 + (-0.035);
298 A10 = -4.592 + 0.379;
299 C1 = -0.487;
300 C2 = 1.024;
301 T9 = 0.374 + 0.252;
302 U9 = 0.033 + 0.015;
303 W9 = 0.032 + 0.012;
304 Lmu = log(muscale/mbeff);
305
306 EvtComplex F71;
307 EvtComplex f71;
308 EvtComplex k7100(-0.68192,-0.074998);
309 EvtComplex k7101(0.0,0.0);
310 EvtComplex k7110(-0.23935,-0.12289);
311 EvtComplex k7111(0.0027424,0.019676);
312 EvtComplex k7120(-0.0018555,-0.175);
313 EvtComplex k7121(0.022864,0.011456);
314 EvtComplex k7130(0.28248,-0.12783);
315 EvtComplex k7131(0.029027,-0.0082265);
316 f71 = k7100 + k7101*logshat + shat*(k7110 + k7111*logshat) +
317 shat*shat*(k7120 + k7121*logshat) +
318 shat*shat*shat*(k7130 + k7131*logshat);
319 F71 = (-208.0/243.0)*Lmu + f71;
320
321 EvtComplex F72;
322 EvtComplex f72;
323 EvtComplex k7200(4.0915,0.44999);
324 EvtComplex k7201(0.0,0.0);
325 EvtComplex k7210(1.4361,0.73732);
326 EvtComplex k7211(-0.016454,-0.11806);
327 EvtComplex k7220(0.011133,1.05);
328 EvtComplex k7221(-0.13718,-0.068733);
329 EvtComplex k7230(-1.6949,0.76698);
330 EvtComplex k7231(-0.17416,0.049359);
331 f72 = k7200 + k7201*logshat + shat*(k7210 + k7211*logshat) +
332 shat*shat*(k7220 + k7221*logshat) +
333 shat*shat*shat*(k7230 + k7231*logshat);
334 F72 = (416.0/81.0)*Lmu + f72;
335
336 EvtComplex F78;
337 F78 = (-32.0/9.0)*Lmu + 8.0*EvtConst::pi*EvtConst::pi/27.0 + (-44.0/9.0)
338 + (-8.0*EvtConst::pi/9.0)*uniti +
339 (4.0/3.0*EvtConst::pi*EvtConst::pi - 40.0/3.0)*shat +
340 (32.0*EvtConst::pi*EvtConst::pi/9.0 - 316.0/9.0)*shat*shat +
341 (200.0*EvtConst::pi*EvtConst::pi/27.0 - 658.0/9.0)*shat*shat*shat +
342 (-8.0*logshat/9.0)*(shat + shat*shat + shat*shat*shat);
343
344 c7eff = A7 - alphas/(4.0*EvtConst::pi)*(C1*F71 + C2*F72 + A8*F78);
345
346 return c7eff;
347}

Referenced by EvtbTosllScalarAmp::CalcAmp(), EvtbTosllVectorAmp::CalcAmp(), dGdsdupProb(), and dGdsProb().

◆ GetC9Eff()

EvtComplex EvtbTosllAmp::GetC9Eff ( double  q2,
bool  nnlo = true,
bool  btod = false 
)

Definition at line 350 of file EvtbTosllAmp.cc.

351{
352
353 if (!nnlo) return 4.344;
354 double mbeff = 4.8;
355 double shat = q2/mbeff/mbeff;
356 double logshat;
357 logshat = log(shat);
358 double mchat = 0.29;
359
360
361 double muscale;
362 muscale = 2.5;
363 double alphas;
364 alphas = 0.267;
365 double A7;
366 A7 = -0.353 + 0.023;
367 double A8;
368 A8 = -0.164;
369 double A9;
370 A9 = 4.287 + (-0.218);
371 double A10;
372 A10 = -4.592 + 0.379;
373 double C1;
374 C1 = -0.697;
375 double C2;
376 C2 = 1.046;
377 double T9;
378 T9 = 0.114 + 0.280;
379 double U9;
380 U9 = 0.045 + 0.023;
381 double W9;
382 W9 = 0.044 + 0.016;
383
384 double Lmu;
385 Lmu = log(muscale/mbeff);
386
387
388 EvtComplex uniti(0.0,1.0);
389
391 double xarg;
392 xarg = 4.0*mchat/shat;
393 hc = -4.0/9.0*log(mchat*mchat) + 8.0/27.0 + 4.0*xarg/9.0;
394
395if (xarg < 1.0)
396 {
397 hc = hc - 2.0/9.0*(2.0 + xarg)*sqrt(fabs(1.0 - xarg))*
398 (log(fabs((sqrt(1.0 - xarg)+1.0)/(sqrt(1.0 - xarg) - 1.0))) -
399 uniti*EvtConst::pi);
400 }
401 else
402 {
403 hc = hc - 2.0/9.0*(2.0 + xarg)*sqrt(fabs(1.0 - xarg))*
404 2.0*atan(1.0/sqrt(xarg - 1.0));
405 }
406
407 EvtComplex h1;
408 xarg = 4.0/shat;
409 h1 = 8.0/27.0 + 4.0*xarg/9.0;
410 if (xarg < 1.0)
411 {
412 h1 = h1 - 2.0/9.0*(2.0 + xarg)*sqrt(fabs(1.0 - xarg))*
413 (log(fabs((sqrt(1.0 - xarg)+1.0)/(sqrt(1.0 - xarg) - 1.0))) -
414 uniti*EvtConst::pi);
415 }
416 else
417 {
418 h1 = h1 - 2.0/9.0*(2.0 + xarg)*sqrt(fabs(1.0 - xarg))*
419 2.0*atan(1.0/sqrt(xarg - 1.0));
420 }
421
422
423 EvtComplex h0;
424 h0 = 8.0/27.0 - 4.0*log(2.0)/9.0 + 4.0*uniti*EvtConst::pi/9.0;
425
426
427 // X=V_{ud}^* V_ub / V_{td}^* V_tb * (4/3 C_1 +C_2) * (h(\hat m_c^2, hat s)-
428 // h(\hat m_u^2, hat s))
429 EvtComplex Vudstar(1.0 - 0.2279*0.2279/2.0, 0.0);
430 EvtComplex Vub((0.118+0.273)/2.0, -1.0*(0.305+0.393)/2.0);
431 EvtComplex Vtdstar(1.0 - (0.118+0.273)/2.0,(0.305+0.393)/2.0);
432 EvtComplex Vtb(1.0,0.0);
433
434 EvtComplex Xd;
435 Xd = (Vudstar * Vub / Vtdstar * Vtb) * (4.0/3.0*C1 + C2) * (hc - h0);
436
437
438 EvtComplex c9eff=4.344;
439 if (shat > 0.25)
440 {
441 c9eff = A9 + T9*hc + U9*h1 + W9*h0;
442 if (btod)
443 {
444 c9eff += Xd;
445 }
446
447 return c9eff;
448 }
449
450 // change energy scale to 5.0 for full NNLO calculation below shat = 0.25
451 muscale = 5.0;
452 alphas = 0.215;
453 A9 = 4.174 + (-0.035);
454 C1 = -0.487;
455 C2 = 1.024;
456 A8 = -0.148;
457 T9 = 0.374 + 0.252;
458 U9 = 0.033 + 0.015;
459 W9 = 0.032 + 0.012;
460 Lmu = log(muscale/mbeff);
461
462 EvtComplex F91;
463 EvtComplex f91;
464 EvtComplex k9100(-11.973,0.16371);
465 EvtComplex k9101(-0.081271,-0.059691);
466 EvtComplex k9110(-28.432,-0.25044);
467 EvtComplex k9111(-0.040243,0.016442);
468 EvtComplex k9120(-57.114,-0.86486);
469 EvtComplex k9121(-0.035191,0.027909);
470 EvtComplex k9130(-128.8,-2.5243);
471 EvtComplex k9131(-0.017587,0.050639);
472 f91 = k9100 + k9101*logshat + shat*(k9110 + k9111*logshat) +
473 shat*shat*(k9120 + k9121*logshat) +
474 shat*shat*shat*(k9130 + k9131*logshat);
475 F91 = (-1424.0/729.0 + 16.0*uniti*EvtConst::pi/243.0
476 + 64.0/27.0*log(mchat))*Lmu - 16.0*Lmu*logshat/243.0 +
477 (16.0/1215.0 - 32.0/135.0/mchat/mchat)*Lmu*shat +
478 (4.0/2835.0 - 8.0/315.0/mchat/mchat/mchat/mchat)*Lmu*shat*shat +
479 (16.0/76545.0 - 32.0/8505.0/mchat/mchat/mchat/mchat/mchat/mchat)*
480 Lmu*shat*shat*shat -256.0*Lmu*Lmu/243.0 + f91;
481
482 EvtComplex F92;
483 EvtComplex f92;
484 EvtComplex k9200(6.6338,-0.98225);
485 EvtComplex k9201(0.48763,0.35815);
486 EvtComplex k9210(3.3585,1.5026);
487 EvtComplex k9211(0.24146,-0.098649);
488 EvtComplex k9220(-1.1906,5.1892);
489 EvtComplex k9221(0.21115,-0.16745);
490 EvtComplex k9230(-17.12,15.146);
491 EvtComplex k9231(0.10552,-0.30383);
492 f92 = k9200 + k9201*logshat + shat*(k9210 + k9211*logshat) +
493 shat*shat*(k9220 + k9221*logshat) +
494 shat*shat*shat*(k9230 + k9231*logshat);
495 F92 = (256.0/243.0 - 32.0*uniti*EvtConst::pi/81.0
496 - 128.0/9.0*log(mchat))*Lmu + 32.0*Lmu*logshat/81.0 +
497 (-32.0/405.0 + 64.0/45.0/mchat/mchat)*Lmu*shat +
498 (-8.0/945.0 + 16.0/105.0/mchat/mchat/mchat/mchat)*Lmu*shat*shat +
499 (-32.0/25515.0 + 64.0/2835.0/mchat/mchat/mchat/mchat/mchat/mchat)*
500 Lmu*shat*shat*shat + 512.0*Lmu*Lmu/81.0 + f92;
501
502 EvtComplex F98;
503 F98 = 104.0/9.0 - 32.0*EvtConst::pi*EvtConst::pi/27.0 +
504 (1184.0/27.0 - 40.0*EvtConst::pi*EvtConst::pi/9.0)*shat +
505 (14212.0/135.0 - 32.0*EvtConst::pi*EvtConst::pi/3.0)*shat*shat +
506 (193444.0/945.0 - 560.0*EvtConst::pi*EvtConst::pi/27.0)*shat*shat*shat +
507 16.0*logshat/9.0*(1.0 + shat + shat*shat + shat*shat*shat);
508
509 Xd = (Vudstar * Vub / Vtdstar * Vtb) * (4.0/3.0*C1 + C2) * (hc - h0);
510
511 c9eff = A9 + T9*hc + U9*h1 + W9*h0 -
512 alphas/(4.0*EvtConst::pi)*(C1*F91 + C2*F92 + A8*F98);
513 if (btod)
514 {
515 c9eff += Xd;
516 }
517
518 return c9eff;
519}
const double hc
Definition: TConstant.h:41

Referenced by EvtbTosllScalarAmp::CalcAmp(), EvtbTosllVectorAmp::CalcAmp(), dGdsdupProb(), and dGdsProb().


The documentation for this class was generated from the following files: