00001 //-------------------------------------------------------------------------- 00002 // 00003 // Environment: 00004 // This software is part of the EvtGen package developed jointly 00005 // for the BaBar and CLEO collaborations. If you use all or part 00006 // of it, please give an appropriate acknowledgement. 00007 // 00008 // Copyright Information: See EvtGen/COPYRIGHT 00009 // Copyright (C) 1998 Caltech, UCSB 00010 // 00011 // Module: EvtTauScalarnu.cc 00012 // 00013 // Description: The leptonic decay of the tau meson. 00014 // E.g., tau- -> e- nueb nut 00015 // 00016 // Modification history: 00017 // 00018 // RYD January 17, 1997 Module created 00019 // 00020 //------------------------------------------------------------------------ 00021 // 00022 #include "EvtGenBase/EvtPatches.hh" 00023 #include "EvtGenBase/EvtPatches.hh" 00024 #include <iostream> 00025 #include <string> 00026 #include "EvtGenBase/EvtParticle.hh" 00027 #include "EvtGenBase/EvtDiracParticle.hh" 00028 #include "EvtGenBase/EvtPDL.hh" 00029 #include "EvtGenBase/EvtIdSet.hh" 00030 #include "EvtGenBase/EvtGenKine.hh" 00031 #include "EvtGenModels/EvtBsquark.hh" 00032 #include "EvtGenBase/EvtDiracSpinor.hh" 00033 #include "EvtGenBase/EvtGammaMatrix.hh" 00034 #include "EvtGenBase/EvtReport.hh" 00035 00036 EvtBsquark::~EvtBsquark() {} 00037 00038 void EvtBsquark::getName(std::string& model_name){ 00039 00040 model_name="BSQUARK"; 00041 00042 } 00043 00044 00045 EvtDecayBase* EvtBsquark::clone(){ 00046 00047 return new EvtBsquark; 00048 00049 } 00050 00051 void EvtBsquark::init(){ 00052 00053 // check that there are 5 arguments 00054 checkNArg(5); 00055 } 00056 00057 void EvtBsquark::initProbMax(){ 00058 00059 //For now do not set a maximum. 00060 00061 //SetProbMax(0.000000000005); 00062 00063 } 00064 00065 void EvtBsquark::decay(EvtParticle *p){ 00066 00067 static EvtId cquark=EvtPDL::getId("c"); 00068 static EvtId anticquark=EvtPDL::getId("anti-c"); 00069 00070 static EvtIdSet leptons("e-","mu-","tau-"); 00071 00072 p->initializePhaseSpace(getNDaug(),getDaugs()); 00073 00074 00075 int charge=1; 00076 00077 EvtParticle* lepton; 00078 lepton = p->getDaug(1); 00079 if (leptons.contains(lepton->getId())){ 00080 charge=-1; 00081 } 00082 00083 00084 EvtDiracParticle charmquark; 00085 00086 //this is a very crude approximation... 00087 if (charge==-1){ 00088 charmquark.init(cquark,p->getDaug(0)->getP4()); 00089 } 00090 else{ 00091 charmquark.init(anticquark,p->getDaug(0)->getP4()); 00092 } 00093 00094 EvtVector4R p4c = p->getDaug(0)->getP4(); 00095 00096 EvtVector4R p4sn = p->getDaug(2)->getP4(); 00097 00098 EvtVector4R p4b(p->mass(),0.0,0.0,0.0); 00099 00100 EvtComplex M[2][2]; 00101 00102 int il,ic; 00103 00104 //project out the right handed current 00105 EvtGammaMatrix PR=0.5*(EvtGammaMatrix::id()+EvtGammaMatrix::g5()); 00106 00107 double tanbeta=getArg(1); 00108 double cosbeta=cos(atan(tanbeta)); 00109 double sinbeta=sin(atan(tanbeta)); 00110 00111 double mb=4.9; 00112 double mc=1.3; 00113 double mw=80.4; 00114 00115 double Mass=getArg(2); 00116 double mu=getArg(3); 00117 double mchargino=getArg(4); 00118 00119 00120 double tan2phim=2*sqrt(2.0)*mw*(mu*cosbeta+Mass*sinbeta)/ 00121 (Mass*Mass-mu*mu+2*mw*mw*cos(2*atan(tanbeta))); 00122 00123 double phim=0.5*atan(tan2phim); 00124 00125 EvtComplex U11=cos(phim); 00126 EvtComplex U12=sin(phim); 00127 EvtComplex U21=-sin(phim); 00128 EvtComplex U22=cos(phim); 00129 00130 double tan2phip=2*sqrt(2.0)*mw*(mu*cosbeta+Mass*sinbeta)/ 00131 (Mass*Mass-mu*mu-2*mw*mw*cos(2*atan(tanbeta))); 00132 00133 double phip=0.5*atan(tan2phip); 00134 00135 EvtComplex V11=cos(phip); 00136 EvtComplex V12=sin(phip); 00137 EvtComplex V21=-sin(phip); 00138 EvtComplex V22=cos(phip); 00139 00140 00141 double theta=getArg(0); 00142 double ctheta=cos(theta); 00143 double stheta=sin(theta); 00144 00145 double vcsb=0.08; 00146 double mchi1=mchargino; 00147 double mchi2=mchargino; 00148 00149 //overall scale factor 00150 double g=1.0; 00151 00152 EvtComplex a1=mchi1*(U11*ctheta-mb*U12*stheta/(sqrt(2.0)*mw*cosbeta)); 00153 EvtComplex a2=mchi2*(U21*ctheta-mb*U22*stheta/(sqrt(2.0)*mw*cosbeta)); 00154 00155 EvtComplex b1=mc*conj(V12)*ctheta/(sqrt(2.0)*mw*sinbeta); 00156 EvtComplex b2=mc*conj(V22)*ctheta/(sqrt(2.0)*mw*sinbeta); 00157 00158 EvtComplex f1=-(g*g*V11*vcsb)/((p4b-p4c).mass2()-mchi1*mchi1); 00159 EvtComplex f2=-(g*g*V21*vcsb)/((p4b-p4c).mass2()-mchi1*mchi2); 00160 00161 //report(INFO,"EvtGen") <<g<<" "<<V11<<" "<<FL<<" "<<vcsb<<" "<<mchi1<<endl; 00162 //report(INFO,"EvtGen") << "f1:"<<f1<<" "<<(p4b-p4c).mass2()<<endl; 00163 //report(INFO,"EvtGen") << "f2:"<<f2<<" "<<(p4b-p4c).mass2()<<endl; 00164 00165 //report(INFO,"EvtGen") << "p4sn:"<<p4sn<<endl; 00166 00167 EvtGammaMatrix pslash=p4sn.get(0)*EvtGammaMatrix::g0() 00168 -p4sn.get(1)*EvtGammaMatrix::g1() 00169 -p4sn.get(2)*EvtGammaMatrix::g2() 00170 -p4sn.get(3)*EvtGammaMatrix::g3(); 00171 00172 //report(INFO,"EvtGen") << "pslash:"<<pslash<<endl; 00173 00174 00175 00176 for(il=0;il<2;il++){ 00177 for(ic=0;ic<2;ic++){ 00178 00179 EvtComplex a=0.0; 00180 EvtComplex b=0.0; 00181 00182 if (charge==-1){ 00183 a=charmquark.spParent(ic)*(PR*lepton->spParent(il)); 00184 b=charmquark.spParent(ic)*((pslash*PR)*lepton->spParent(il)); 00185 } 00186 else{ 00187 a=lepton->spParent(il)*(PR*charmquark.spParent(ic)); 00188 b=lepton->spParent(il)*((pslash*PR)*charmquark.spParent(ic)); 00189 } 00190 00191 //report(INFO,"EvtGen") <<"pslash*PR:"<<pslash*PR<<endl; 00192 //report(INFO,"EvtGen") <<"sp charm:"<<charmquark.spParent(ic)<<endl; 00193 //report(INFO,"EvtGen") <<"sp lepton:"<<lepton->spParent(il)<<endl; 00194 00195 M[ic][il]=f1*(a1*a+b1*b)+f2*(a2*a+b2*b); 00196 00197 //report(INFO,"EvtGen") << "Contr1:"<<a1<<" "<<a<<" "<<b1<<" "<<b<<endl; 00198 //report(INFO,"EvtGen") << "Contr2:"<<a2<<" "<<a<<" "<<b2<<" "<<b<<endl; 00199 00200 //report(INFO,"EvtGen") <<"case1:"<<f1<<" "<<a1<<" "<<b1<<" "<<a<<" "<<b<<endl; 00201 //report(INFO,"EvtGen") <<"case2:"<<f2<<" "<<a2<<" "<<b2<<" "<<a<<" "<<b<<endl; 00202 00203 } 00204 } 00205 00206 double prob=real(M[0][0]*conj(M[0][0])+ 00207 M[1][0]*conj(M[1][0])+ 00208 M[0][1]*conj(M[0][1])+ 00209 M[1][1]*conj(M[1][1])); 00210 00211 //report(INFO,"EvtGen") <<"prob:"<<prob<<endl; 00212 00213 setProb(prob); 00214 00215 return; 00216 00217 } 00218