#include <EvtSemiLeptonicBaryonAmp.hh>
Inheritance diagram for EvtSemiLeptonicBaryonAmp:
Public Member Functions | |
| void | CalcAmp (EvtParticle *parent, EvtAmp &, EvtSemiLeptonicFF *FormFactors) |
| double | CalcMaxProb (EvtId parent, EvtId meson, EvtId lepton, EvtId nudaug, EvtSemiLeptonicFF *FormFactors) |
Definition at line 30 of file EvtSemiLeptonicBaryonAmp.hh.
| void EvtSemiLeptonicBaryonAmp::CalcAmp | ( | EvtParticle * | parent, | |
| EvtAmp & | amp, | |||
| EvtSemiLeptonicFF * | FormFactors | |||
| ) | [virtual] |
Reimplemented from EvtSemiLeptonicAmp.
Definition at line 38 of file EvtSemiLeptonicBaryonAmp.cc.
References EvtVector4C::cont(), calibUtil::ERROR, EvtLeptonVACurrent(), EvtGammaMatrix::g0(), EvtGammaMatrix::g1(), EvtGammaMatrix::g2(), EvtGammaMatrix::g3(), EvtGammaMatrix::g5(), EvtSemiLeptonicFF::getbaryonff(), EvtParticle::getDaug(), EvtParticle::getId(), EvtPDL::getId(), EvtParticle::getP4(), EvtParticle::mass(), q, report(), EvtVector4C::set(), EvtVector4R::set(), EvtParticle::sp(), EvtParticle::spParent(), EvtParticle::spParentNeutrino(), and EvtAmp::vertex().
00040 { 00041 00042 static EvtId EM=EvtPDL::getId("e-"); 00043 static EvtId MUM=EvtPDL::getId("mu-"); 00044 static EvtId TAUM=EvtPDL::getId("tau-"); 00045 static EvtId EP=EvtPDL::getId("e+"); 00046 static EvtId MUP=EvtPDL::getId("mu+"); 00047 static EvtId TAUP=EvtPDL::getId("tau+"); 00048 00049 00050 //Add the lepton and neutrino 4 momenta to find q2 00051 00052 EvtVector4R q = parent->getDaug(1)->getP4() 00053 + parent->getDaug(2)->getP4(); 00054 double q2 = (q.mass2()); 00055 00056 double f1v,f1a,f2v,f2a; 00057 double m_meson = parent->getDaug(0)->mass(); 00058 00059 FormFactors->getbaryonff(parent->getId(), 00060 parent->getDaug(0)->getId(), 00061 q2, 00062 m_meson, 00063 &f1v, 00064 &f1a, 00065 &f2v, 00066 &f2a); 00067 00068 EvtVector4R p4b; 00069 p4b.set(parent->mass(),0.0,0.0,0.0); 00070 00071 EvtVector4C temp_00_term1; 00072 EvtVector4C temp_00_term2; 00073 00074 EvtVector4C temp_01_term1; 00075 EvtVector4C temp_01_term2; 00076 00077 EvtVector4C temp_10_term1; 00078 EvtVector4C temp_10_term2; 00079 00080 EvtVector4C temp_11_term1; 00081 EvtVector4C temp_11_term2; 00082 00083 EvtDiracSpinor p0=parent->sp(0); 00084 EvtDiracSpinor p1=parent->sp(1); 00085 00086 EvtDiracSpinor d0=parent->getDaug(0)->spParent(0); 00087 EvtDiracSpinor d1=parent->getDaug(0)->spParent(1); 00088 00089 temp_00_term1.set(0,f1v*(d0*(EvtGammaMatrix::g0()*p0))); 00090 temp_00_term2.set(0,f1a*(d0*((EvtGammaMatrix::g0()*EvtGammaMatrix::g5())*p0))); 00091 temp_01_term1.set(0,f1v*(d0*(EvtGammaMatrix::g0()*p1))); 00092 temp_01_term2.set(0,f1a*(d0*((EvtGammaMatrix::g0()*EvtGammaMatrix::g5())*p1))); 00093 temp_10_term1.set(0,f1v*(d1*(EvtGammaMatrix::g0()*p0))); 00094 temp_10_term2.set(0,f1a*(d1*((EvtGammaMatrix::g0()*EvtGammaMatrix::g5())*p0))); 00095 temp_11_term1.set(0,f1v*(d1*(EvtGammaMatrix::g0()*p1))); 00096 temp_11_term2.set(0,f1a*(d1*((EvtGammaMatrix::g0()*EvtGammaMatrix::g5())*p1))); 00097 00098 temp_00_term1.set(1,f1v*(d0*(EvtGammaMatrix::g1()*p0))); 00099 temp_00_term2.set(1,f1a*(d0*((EvtGammaMatrix::g1()*EvtGammaMatrix::g5())*p0))); 00100 temp_01_term1.set(1,f1v*(d0*(EvtGammaMatrix::g1()*p1))); 00101 temp_01_term2.set(1,f1a*(d0*((EvtGammaMatrix::g1()*EvtGammaMatrix::g5())*p1))); 00102 temp_10_term1.set(1,f1v*(d1*(EvtGammaMatrix::g1()*p0))); 00103 temp_10_term2.set(1,f1a*(d1*((EvtGammaMatrix::g1()*EvtGammaMatrix::g5())*p0))); 00104 temp_11_term1.set(1,f1v*(d1*(EvtGammaMatrix::g1()*p1))); 00105 temp_11_term2.set(1,f1a*(d1*((EvtGammaMatrix::g1()*EvtGammaMatrix::g5())*p1))); 00106 00107 temp_00_term1.set(2,f1v*(d0*(EvtGammaMatrix::g2()*p0))); 00108 temp_00_term2.set(2,f1a*(d0*((EvtGammaMatrix::g2()*EvtGammaMatrix::g5())*p0))); 00109 temp_01_term1.set(2,f1v*(d0*(EvtGammaMatrix::g2()*p1))); 00110 temp_01_term2.set(2,f1a*(d0*((EvtGammaMatrix::g2()*EvtGammaMatrix::g5())*p1))); 00111 temp_10_term1.set(2,f1v*(d1*(EvtGammaMatrix::g2()*p0))); 00112 temp_10_term2.set(2,f1a*(d1*((EvtGammaMatrix::g2()*EvtGammaMatrix::g5())*p0))); 00113 temp_11_term1.set(2,f1v*(d1*(EvtGammaMatrix::g2()*p1))); 00114 temp_11_term2.set(2,f1a*(d1*((EvtGammaMatrix::g2()*EvtGammaMatrix::g5())*p1))); 00115 00116 temp_00_term1.set(3,f1v*(d0*(EvtGammaMatrix::g3()*p0))); 00117 temp_00_term2.set(3,f1a*(d0*((EvtGammaMatrix::g3()*EvtGammaMatrix::g5())*p0))); 00118 temp_01_term1.set(3,f1v*(d0*(EvtGammaMatrix::g3()*p1))); 00119 temp_01_term2.set(3,f1a*(d0*((EvtGammaMatrix::g3()*EvtGammaMatrix::g5())*p1))); 00120 temp_10_term1.set(3,f1v*(d1*(EvtGammaMatrix::g3()*p0))); 00121 temp_10_term2.set(3,f1a*(d1*((EvtGammaMatrix::g3()*EvtGammaMatrix::g5())*p0))); 00122 temp_11_term1.set(3,f1v*(d1*(EvtGammaMatrix::g3()*p1))); 00123 temp_11_term2.set(3,f1a*(d1*((EvtGammaMatrix::g3()*EvtGammaMatrix::g5())*p1))); 00124 00125 00126 00127 EvtVector4C l1,l2; 00128 00129 EvtId l_num = parent->getDaug(1)->getId(); 00130 if (l_num==EM||l_num==MUM||l_num==TAUM){ 00131 00132 l1=EvtLeptonVACurrent(parent->getDaug(1)->spParent(0), 00133 parent->getDaug(2)->spParentNeutrino()); 00134 l2=EvtLeptonVACurrent(parent->getDaug(1)->spParent(1), 00135 parent->getDaug(2)->spParentNeutrino()); 00136 } 00137 else{ 00138 if (l_num==EP||l_num==MUP||l_num==TAUP){ 00139 l1=EvtLeptonVACurrent(parent->getDaug(2)->spParentNeutrino(), 00140 parent->getDaug(1)->spParent(0)); 00141 l2=EvtLeptonVACurrent(parent->getDaug(2)->spParentNeutrino(), 00142 parent->getDaug(1)->spParent(1)); 00143 } 00144 else{ 00145 report(ERROR,"EvtGen") << "Wrong lepton number"<<endl; 00146 } 00147 } 00148 00149 amp.vertex(0,0,0,l1.cont(temp_00_term1+temp_00_term2)); 00150 amp.vertex(0,0,1,l2.cont(temp_00_term1+temp_00_term2)); 00151 00152 amp.vertex(0,1,0,l1.cont(temp_01_term1+temp_01_term2)); 00153 amp.vertex(0,1,1,l2.cont(temp_01_term1+temp_01_term2)); 00154 00155 amp.vertex(1,0,0,l1.cont(temp_10_term1+temp_10_term2)); 00156 amp.vertex(1,0,1,l2.cont(temp_10_term1+temp_10_term2)); 00157 00158 amp.vertex(1,1,0,l1.cont(temp_11_term1+temp_11_term2)); 00159 amp.vertex(1,1,1,l2.cont(temp_11_term1+temp_11_term2)); 00160 00161 return; 00162 }
| double EvtSemiLeptonicAmp::CalcMaxProb | ( | EvtId | parent, | |
| EvtId | meson, | |||
| EvtId | lepton, | |||
| EvtId | nudaug, | |||
| EvtSemiLeptonicFF * | FormFactors | |||
| ) | [inherited] |
Definition at line 38 of file EvtSemiLeptonicAmp.cc.
References boostTo(), EvtSemiLeptonicAmp::CalcAmp(), EvtParticle::deleteTree(), EvtParticle::getDaug(), EvtPDL::getMass(), EvtPDL::getMaxMass(), EvtPDL::getMeanMass(), EvtPDL::getMinMass(), EvtAmp::getSpinDensity(), EvtParticle::getSpinStates(), EvtPDL::getWidth(), genRecEmupikp::i, EvtParticle::init(), EvtAmp::init(), EvtScalarParticle::init(), ganga-rec::j, EvtParticle::makeDaughters(), EvtParticle::mass(), mass, EvtParticle::noLifeTime(), EvtSpinDensity::NormalizedProb(), EvtVector4R::set(), EvtSpinDensity::SetDiag(), and EvtParticle::setDiagonalSpinDensity().
Referenced by EvtSLPole::initProbMax(), EvtSLBKPole::initProbMax(), EvtHQET2::initProbMax(), and EvtHQET::initProbMax().
00040 { 00041 00042 //This routine takes the arguements parent, meson, and lepton 00043 //number, and a form factor model, and returns a maximum 00044 //probability for this semileptonic form factor model. A 00045 //brute force method is used. The 2D cos theta lepton and 00046 //q2 phase space is probed. 00047 00048 //Start by declaring a particle at rest. 00049 00050 //It only makes sense to have a scalar parent. For now. 00051 //This should be generalized later. 00052 00053 EvtScalarParticle *scalar_part; 00054 EvtParticle *root_part; 00055 00056 scalar_part=new EvtScalarParticle; 00057 00058 //cludge to avoid generating random numbers! 00059 scalar_part->noLifeTime(); 00060 00061 EvtVector4R p_init; 00062 00063 p_init.set(EvtPDL::getMass(parent),0.0,0.0,0.0); 00064 scalar_part->init(parent,p_init); 00065 root_part=(EvtParticle *)scalar_part; 00066 // root_part->set_type(EvtSpinType::SCALAR); 00067 root_part->setDiagonalSpinDensity(); 00068 00069 EvtParticle *daughter, *lep, *trino; 00070 00071 EvtAmp amp; 00072 00073 EvtId listdaug[3]; 00074 listdaug[0] = meson; 00075 listdaug[1] = lepton; 00076 listdaug[2] = nudaug; 00077 00078 amp.init(parent,3,listdaug); 00079 00080 root_part->makeDaughters(3,listdaug); 00081 daughter=root_part->getDaug(0); 00082 lep=root_part->getDaug(1); 00083 trino=root_part->getDaug(2); 00084 00085 //cludge to avoid generating random numbers! 00086 daughter->noLifeTime(); 00087 lep->noLifeTime(); 00088 trino->noLifeTime(); 00089 00090 00091 //Initial particle is unpolarized, well it is a scalar so it is 00092 //trivial 00093 EvtSpinDensity rho; 00094 rho.SetDiag(root_part->getSpinStates()); 00095 00096 double mass[3]; 00097 00098 double m = root_part->mass(); 00099 00100 EvtVector4R p4meson, p4lepton, p4nu, p4w; 00101 double q2max; 00102 00103 double q2, elepton, plepton; 00104 int i,j; 00105 double erho,prho,costl; 00106 00107 double maxfoundprob = 0.0; 00108 double prob = -10.0; 00109 int massiter; 00110 00111 for (massiter=0;massiter<3;massiter++){ 00112 00113 mass[0] = EvtPDL::getMeanMass(meson); 00114 mass[1] = EvtPDL::getMeanMass(lepton); 00115 mass[2] = EvtPDL::getMeanMass(nudaug); 00116 if ( massiter==1 ) { 00117 mass[0] = EvtPDL::getMinMass(meson); 00118 } 00119 if ( massiter==2 ) { 00120 mass[0] = EvtPDL::getMaxMass(meson); 00121 if ( (mass[0]+mass[1]+mass[2])>m) mass[0]=m-mass[1]-mass[2]-0.00001; 00122 } 00123 00124 q2max = (m-mass[0])*(m-mass[0]); 00125 00126 //loop over q2 00127 00128 for (i=0;i<25;i++) { 00129 q2 = ((i+0.5)*q2max)/25.0; 00130 00131 erho = ( m*m + mass[0]*mass[0] - q2 )/(2.0*m); 00132 00133 prho = sqrt(erho*erho-mass[0]*mass[0]); 00134 00135 p4meson.set(erho,0.0,0.0,-1.0*prho); 00136 p4w.set(m-erho,0.0,0.0,prho); 00137 00138 //This is in the W rest frame 00139 elepton = (q2+mass[1]*mass[1])/(2.0*sqrt(q2)); 00140 plepton = sqrt(elepton*elepton-mass[1]*mass[1]); 00141 00142 double probctl[3]; 00143 00144 for (j=0;j<3;j++) { 00145 00146 costl = 0.99*(j - 1.0); 00147 00148 //These are in the W rest frame. Need to boost out into 00149 //the B frame. 00150 p4lepton.set(elepton,0.0, 00151 plepton*sqrt(1.0-costl*costl),plepton*costl); 00152 p4nu.set(plepton,0.0, 00153 -1.0*plepton*sqrt(1.0-costl*costl),-1.0*plepton*costl); 00154 00155 EvtVector4R boost((m-erho),0.0,0.0,1.0*prho); 00156 p4lepton=boostTo(p4lepton,boost); 00157 p4nu=boostTo(p4nu,boost); 00158 00159 //Now initialize the daughters... 00160 00161 daughter->init(meson,p4meson); 00162 lep->init(lepton,p4lepton); 00163 trino->init(nudaug,p4nu); 00164 00165 CalcAmp(root_part,amp,FormFactors); 00166 00167 //Now find the probability at this q2 and cos theta lepton point 00168 //and compare to maxfoundprob. 00169 00170 //Do a little magic to get the probability!! 00171 prob = rho.NormalizedProb(amp.getSpinDensity()); 00172 00173 probctl[j]=prob; 00174 } 00175 00176 //probclt contains prob at ctl=-1,0,1. 00177 //prob=a+b*ctl+c*ctl^2 00178 00179 double a=probctl[1]; 00180 double b=0.5*(probctl[2]-probctl[0]); 00181 double c=0.5*(probctl[2]+probctl[0])-probctl[1]; 00182 00183 prob=probctl[0]; 00184 if (probctl[1]>prob) prob=probctl[1]; 00185 if (probctl[2]>prob) prob=probctl[2]; 00186 00187 if (fabs(c)>1e-20){ 00188 double ctlx=-0.5*b/c; 00189 if (fabs(ctlx)<1.0){ 00190 double probtmp=a+b*ctlx+c*ctlx*ctlx; 00191 if (probtmp>prob) prob=probtmp; 00192 } 00193 00194 } 00195 00196 //report(DEBUG,"EvtGen") << "prob,probctl:"<<prob<<" " 00197 // << probctl[0]<<" " 00198 // << probctl[1]<<" " 00199 // << probctl[2]<<endl; 00200 00201 if ( prob > maxfoundprob ) { 00202 maxfoundprob = prob; 00203 } 00204 00205 } 00206 if ( EvtPDL::getWidth(meson) <= 0.0 ) { 00207 //if the particle is narrow dont bother with changing the mass. 00208 massiter = 4; 00209 } 00210 00211 } 00212 root_part->deleteTree(); 00213 00214 maxfoundprob *=1.1; 00215 return maxfoundprob; 00216 00217 }
1.4.7