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: EvtTauScalarnu.cc
12 : //
13 : // Description: The leptonic decay of the tau meson.
14 : // E.g., tau- -> e- nueb nut
15 : //
16 : // Modification history:
17 : //
18 : // RYD January 17, 1997 Module created
19 : //
20 : //------------------------------------------------------------------------
21 : //
22 : #include "EvtGenBase/EvtPatches.hh"
23 : #include "EvtGenBase/EvtPatches.hh"
24 : #include <iostream>
25 : #include <string>
26 : #include "EvtGenBase/EvtParticle.hh"
27 : #include "EvtGenBase/EvtDiracParticle.hh"
28 : #include "EvtGenBase/EvtPDL.hh"
29 : #include "EvtGenBase/EvtIdSet.hh"
30 : #include "EvtGenBase/EvtGenKine.hh"
31 : #include "EvtGenModels/EvtBsquark.hh"
32 : #include "EvtGenBase/EvtDiracSpinor.hh"
33 : #include "EvtGenBase/EvtGammaMatrix.hh"
34 : #include "EvtGenBase/EvtReport.hh"
35 :
36 0 : EvtBsquark::~EvtBsquark() {}
37 :
38 : std::string EvtBsquark::getName(){
39 :
40 0 : return "BSQUARK";
41 :
42 : }
43 :
44 :
45 : EvtDecayBase* EvtBsquark::clone(){
46 :
47 0 : return new EvtBsquark;
48 :
49 0 : }
50 :
51 : void EvtBsquark::init(){
52 :
53 : // check that there are 5 arguments
54 0 : checkNArg(5);
55 0 : }
56 :
57 : void EvtBsquark::initProbMax(){
58 :
59 : //For now do not set a maximum.
60 :
61 : //SetProbMax(0.000000000005);
62 :
63 0 : }
64 :
65 : void EvtBsquark::decay(EvtParticle *p){
66 :
67 0 : static EvtId cquark=EvtPDL::getId("c");
68 0 : static EvtId anticquark=EvtPDL::getId("anti-c");
69 :
70 0 : static EvtIdSet leptons("e-","mu-","tau-");
71 :
72 0 : p->initializePhaseSpace(getNDaug(),getDaugs());
73 :
74 :
75 : int charge=1;
76 :
77 : EvtParticle* lepton;
78 0 : lepton = p->getDaug(1);
79 0 : if (leptons.contains(lepton->getId())){
80 : charge=-1;
81 0 : }
82 :
83 :
84 0 : EvtDiracParticle charmquark;
85 :
86 : //this is a very crude approximation...
87 0 : if (charge==-1){
88 0 : charmquark.init(cquark,p->getDaug(0)->getP4());
89 : }
90 : else{
91 0 : charmquark.init(anticquark,p->getDaug(0)->getP4());
92 : }
93 :
94 0 : EvtVector4R p4c = p->getDaug(0)->getP4();
95 :
96 0 : EvtVector4R p4sn = p->getDaug(2)->getP4();
97 :
98 0 : EvtVector4R p4b(p->mass(),0.0,0.0,0.0);
99 :
100 0 : EvtComplex M[2][2];
101 :
102 : int il,ic;
103 :
104 : //project out the right handed current
105 0 : EvtGammaMatrix PR=0.5*(EvtGammaMatrix::id()+EvtGammaMatrix::g5());
106 :
107 0 : double tanbeta=getArg(1);
108 0 : double cosbeta=cos(atan(tanbeta));
109 0 : double sinbeta=sin(atan(tanbeta));
110 :
111 : double mb=4.9;
112 : double mc=1.3;
113 : double mw=80.4;
114 :
115 0 : double Mass=getArg(2);
116 0 : double mu=getArg(3);
117 0 : double mchargino=getArg(4);
118 :
119 :
120 0 : double tan2phim=2*sqrt(2.0)*mw*(mu*cosbeta+Mass*sinbeta)/
121 0 : (Mass*Mass-mu*mu+2*mw*mw*cos(2*atan(tanbeta)));
122 :
123 0 : double phim=0.5*atan(tan2phim);
124 :
125 0 : EvtComplex U11=cos(phim);
126 0 : EvtComplex U12=sin(phim);
127 0 : EvtComplex U21=-sin(phim);
128 0 : EvtComplex U22=cos(phim);
129 :
130 0 : double tan2phip=2*sqrt(2.0)*mw*(mu*cosbeta+Mass*sinbeta)/
131 0 : (Mass*Mass-mu*mu-2*mw*mw*cos(2*atan(tanbeta)));
132 :
133 0 : double phip=0.5*atan(tan2phip);
134 :
135 0 : EvtComplex V11=cos(phip);
136 0 : EvtComplex V12=sin(phip);
137 0 : EvtComplex V21=-sin(phip);
138 0 : EvtComplex V22=cos(phip);
139 :
140 :
141 0 : double theta=getArg(0);
142 0 : double ctheta=cos(theta);
143 0 : double stheta=sin(theta);
144 :
145 : double vcsb=0.08;
146 : double mchi1=mchargino;
147 : double mchi2=mchargino;
148 :
149 : //overall scale factor
150 : double g=1.0;
151 :
152 0 : EvtComplex a1=mchi1*(U11*ctheta-mb*U12*stheta/(sqrt(2.0)*mw*cosbeta));
153 0 : EvtComplex a2=mchi2*(U21*ctheta-mb*U22*stheta/(sqrt(2.0)*mw*cosbeta));
154 :
155 0 : EvtComplex b1=mc*conj(V12)*ctheta/(sqrt(2.0)*mw*sinbeta);
156 0 : EvtComplex b2=mc*conj(V22)*ctheta/(sqrt(2.0)*mw*sinbeta);
157 :
158 0 : EvtComplex f1=-(g*g*V11*vcsb)/((p4b-p4c).mass2()-mchi1*mchi1);
159 0 : EvtComplex f2=-(g*g*V21*vcsb)/((p4b-p4c).mass2()-mchi1*mchi2);
160 :
161 : //report(Severity::Info,"EvtGen") <<g<<" "<<V11<<" "<<FL<<" "<<vcsb<<" "<<mchi1<<endl;
162 : //report(Severity::Info,"EvtGen") << "f1:"<<f1<<" "<<(p4b-p4c).mass2()<<endl;
163 : //report(Severity::Info,"EvtGen") << "f2:"<<f2<<" "<<(p4b-p4c).mass2()<<endl;
164 :
165 : //report(Severity::Info,"EvtGen") << "p4sn:"<<p4sn<<endl;
166 :
167 0 : EvtGammaMatrix pslash=p4sn.get(0)*EvtGammaMatrix::g0()
168 0 : -p4sn.get(1)*EvtGammaMatrix::g1()
169 0 : -p4sn.get(2)*EvtGammaMatrix::g2()
170 0 : -p4sn.get(3)*EvtGammaMatrix::g3();
171 :
172 : //report(Severity::Info,"EvtGen") << "pslash:"<<pslash<<endl;
173 :
174 :
175 :
176 0 : for(il=0;il<2;il++){
177 0 : for(ic=0;ic<2;ic++){
178 :
179 0 : EvtComplex a=0.0;
180 0 : EvtComplex b=0.0;
181 :
182 0 : if (charge==-1){
183 0 : a=charmquark.spParent(ic)*(PR*lepton->spParent(il));
184 0 : b=charmquark.spParent(ic)*((pslash*PR)*lepton->spParent(il));
185 0 : }
186 : else{
187 0 : a=lepton->spParent(il)*(PR*charmquark.spParent(ic));
188 0 : b=lepton->spParent(il)*((pslash*PR)*charmquark.spParent(ic));
189 : }
190 :
191 : //report(Severity::Info,"EvtGen") <<"pslash*PR:"<<pslash*PR<<endl;
192 : //report(Severity::Info,"EvtGen") <<"sp charm:"<<charmquark.spParent(ic)<<endl;
193 : //report(Severity::Info,"EvtGen") <<"sp lepton:"<<lepton->spParent(il)<<endl;
194 :
195 0 : M[ic][il]=f1*(a1*a+b1*b)+f2*(a2*a+b2*b);
196 :
197 : //report(Severity::Info,"EvtGen") << "Contr1:"<<a1<<" "<<a<<" "<<b1<<" "<<b<<endl;
198 : //report(Severity::Info,"EvtGen") << "Contr2:"<<a2<<" "<<a<<" "<<b2<<" "<<b<<endl;
199 :
200 : //report(Severity::Info,"EvtGen") <<"case1:"<<f1<<" "<<a1<<" "<<b1<<" "<<a<<" "<<b<<endl;
201 : //report(Severity::Info,"EvtGen") <<"case2:"<<f2<<" "<<a2<<" "<<b2<<" "<<a<<" "<<b<<endl;
202 :
203 0 : }
204 : }
205 :
206 0 : double prob=real(M[0][0]*conj(M[0][0])+
207 0 : M[1][0]*conj(M[1][0])+
208 0 : M[0][1]*conj(M[0][1])+
209 0 : M[1][1]*conj(M[1][1]));
210 :
211 : //report(Severity::Info,"EvtGen") <<"prob:"<<prob<<endl;
212 :
213 0 : setProb(prob);
214 :
215 : return;
216 :
217 0 : }
218 :
|
