]>
Commit | Line | Data |
---|---|---|
da0e9ce3 | 1 | //----------------------------------------------------------------------------------------------- |
2 | // | |
3 | // Module: EvtHypNonLepton.cpp | |
4 | // | |
5 | // Desription: Routine to implement Hyperon(s=1/2) -> Baryon(s=1/2) + Scalar decays accroding to | |
6 | // Review Of Particle Physics 2004, Phys.Lett.B, Vol.592, p.864 | |
7 | // | |
8 | // Modification history: | |
9 | // | |
10 | // 09/02/2009 PR Corrected Delta sign | |
11 | // 20/02/2005 PR Module created according to PHSP and Lb2Lll model | |
12 | // | |
13 | //----------------------------------------------------------------------------------------------- | |
14 | ||
15 | #ifdef WIN32 | |
16 | #pragma warning( disable : 4786 ) | |
17 | // Disable anoying warning about symbol size | |
18 | #endif | |
19 | ||
20 | #include "EvtGenModels/EvtHypNonLepton.hh" | |
21 | #include "EvtGenBase/EvtParticle.hh" | |
22 | #include "EvtGenBase/EvtPDL.hh" | |
23 | #include "EvtGenBase/EvtDiracSpinor.hh" | |
24 | #include "EvtGenBase/EvtTensor4C.hh" | |
25 | #include "EvtGenBase/EvtVector4C.hh" | |
26 | #include "EvtGenBase/EvtVector4R.hh" | |
27 | #include "EvtGenBase/EvtComplex.hh" | |
28 | #include "EvtGenBase/EvtGammaMatrix.hh" | |
29 | ||
30 | ||
31 | EvtHypNonLepton::~EvtHypNonLepton() {} | |
32 | ||
33 | EvtDecayBase* EvtHypNonLepton::clone(){ | |
34 | return new EvtHypNonLepton; | |
35 | } | |
36 | ||
37 | std::string EvtHypNonLepton::getName(){ | |
38 | return "HypNonLepton"; | |
39 | } | |
40 | ||
41 | void EvtHypNonLepton::init(){ | |
42 | ||
43 | if(getNArg()<2 || getNArg()>3){ // alpha phi gamma delta | |
44 | report(ERROR,"EvtGen") << " ERROR: EvtHypNonLepton generator expected 2 or 3 arguments but found: " << getNArg() << std::endl; | |
45 | report(INFO ,"EvtGen") << " 1. Decay asymmetry parameter - alpha" << std::endl; | |
46 | report(INFO ,"EvtGen") << " 2. Parameter phi - in degrees (not radians)" << std::endl; | |
47 | report(INFO ,"EvtGen") << " 3. Note on every x-th decay" << std::endl; | |
48 | ::abort(); | |
49 | } | |
50 | ||
51 | if(getNDaug()!=2){ // Check that there are 2 daughters only | |
52 | report(ERROR,"EvtGen") << " ERROR: EvtHypNonLepton generator expected 2 daughters but found: " << getNDaug() << std::endl; | |
53 | ::abort(); | |
54 | } | |
55 | ||
56 | // Check particles spins | |
57 | if(EvtSpinType::getSpin2(EvtPDL::getSpinType(getParentId()))!=1){ | |
58 | report(ERROR,"EvtGen") << " ERROR: EvtHypNonLepton generator expected dirac parent particle, but found " << EvtSpinType::getSpin2(EvtPDL::getSpinType(getParentId())) << " spin degrees of freedom" << std::endl; | |
59 | ::abort(); | |
60 | } | |
61 | if(EvtSpinType::getSpin2(EvtPDL::getSpinType(getDaug(0)))!=1){ | |
62 | report(ERROR,"EvtGen") << " ERROR: EvtHypNonLepton generator expected the first child to be dirac particle, but found " << EvtSpinType::getSpin2(EvtPDL::getSpinType(getDaug(0))) << " spin degrees of freedom" << std::endl; | |
63 | ::abort(); | |
64 | } | |
65 | if(EvtSpinType::getSpin2(EvtPDL::getSpinType(getDaug(1)))!=0){ | |
66 | report(ERROR,"EvtGen") << " ERROR: EvtHypNonLepton generator expected the second child to be scalar particle, but found " << EvtSpinType::getSpin2(EvtPDL::getSpinType(getDaug(1))) << " spin degrees of freedom" << std::endl; | |
67 | ::abort(); | |
68 | } | |
69 | ||
70 | // Read all parameters | |
71 | m_alpha = getArg(0); | |
72 | m_phi = getArg(1)*EvtConst::pi/180; | |
73 | if(getNArg()==3) m_noTries = static_cast<long>(getArg(2)); | |
74 | else m_noTries = 0; | |
75 | ||
76 | // calculate additional parameters | |
77 | double p,M,m1,m2; | |
78 | double p_to_s,beta,delta,gamma; | |
79 | ||
80 | M = EvtPDL::getMass(getParentId()); | |
81 | m1 = EvtPDL::getMass(getDaug(0)); | |
82 | m2 = EvtPDL::getMass(getDaug(1)); | |
83 | ||
84 | if(m1+m2>=M){ | |
85 | report(ERROR,"EvtGen") << " ERROR: EvtHypNonLepton found impossible decay: " << M << " --> " << m1 << " + " << m2 << " GeV\n" << std::endl; | |
86 | ::abort(); | |
87 | } | |
88 | ||
89 | p = sqrt(M*M-(m1+m2)*(m1+m2))*sqrt(M*M-(m1-m2)*(m1-m2))/2./M; | |
90 | ||
91 | beta = sqrt(1.-m_alpha*m_alpha)*sin(m_phi); | |
92 | delta = -atan2(beta,m_alpha); | |
93 | gamma = sqrt(1.-m_alpha*m_alpha-beta*beta); | |
94 | p_to_s = sqrt((1.-gamma)/(1.+gamma)); | |
95 | ||
96 | m_B_to_A = p_to_s*(m1+sqrt(p*p+m1*m1))/p*EvtComplex(cos(delta),sin(delta)); | |
97 | ||
98 | } | |
99 | ||
100 | void EvtHypNonLepton::initProbMax(){ | |
101 | ||
102 | double maxProb,m1,m2,M,p; | |
103 | ||
104 | M=EvtPDL::getMass(getParentId()); | |
105 | m1=EvtPDL::getMass(getDaug(0)); | |
106 | m2=EvtPDL::getMass(getDaug(1)); | |
107 | ||
108 | if(m1+m2>=M){ | |
109 | report(ERROR,"EvtGen") << " ERROR: EvtHypNonLepton found impossible decay: " << M << " --> " << m1 << " + " << m2 << " GeV\n" << std::endl; | |
110 | ::abort(); | |
111 | } | |
112 | ||
113 | p=sqrt(M*M-(m1+m2)*(m1+m2))*sqrt(M*M-(m1-m2)*(m1-m2))/2/M; | |
114 | maxProb=16*M*(sqrt(p*p+m1*m1)+m1+abs(m_B_to_A)*abs(m_B_to_A)*(sqrt(p*p+m1*m1)-m1)); | |
115 | //maxProb *= G_F*M_pi*M_pi; | |
116 | ||
117 | setProbMax(maxProb); | |
118 | report(INFO,"EvtGen") << " EvtHypNonLepton set up maximum probability to " << maxProb << std::endl; | |
119 | ||
120 | } | |
121 | ||
122 | void EvtHypNonLepton::decay(EvtParticle* parent){ | |
123 | ||
124 | parent->initializePhaseSpace(getNDaug(),getDaugs()); | |
125 | calcAmp(&_amp2,parent); | |
126 | ||
127 | } | |
128 | ||
129 | void EvtHypNonLepton::calcAmp(EvtAmp *amp,EvtParticle *parent){ | |
130 | ||
131 | static long noTries=0; | |
132 | int i; | |
0ca57c2f | 133 | EvtComplex Matrix[2][2]; |
da0e9ce3 | 134 | |
135 | //G_F = 1.16637e-5; | |
136 | //M_pi = 0.13957; | |
137 | ||
138 | for(i=0;i<4;i++){ | |
139 | //std::cout << "--------------------------------------------------" << std::endl; | |
140 | Matrix[i/2][i%2] = EvtLeptonSCurrent(parent->sp(i/2),parent->getDaug(0)->spParent(i%2)); | |
141 | //std::cout << "Matrix = " << Matrix[i/2][i%2] << std::endl; | |
142 | Matrix[i/2][i%2] -= m_B_to_A*EvtLeptonPCurrent(parent->sp(i/2),parent->getDaug(0)->spParent(i%2)); | |
143 | //std::cout << "Matrix = " << Matrix[i/2][i%2] << std::endl; | |
144 | //Matrix[i/2][i%2] *= G_F*M_pi*M_pi; | |
145 | //std::cout << "Matrix = " << Matrix[i/2][i%2] << std::endl; | |
146 | //std::cout << "--------------------------------------------------" << std::endl; | |
147 | amp->vertex(i/2,i%2,Matrix[i/2][i%2]); | |
148 | } | |
149 | ||
150 | if(m_noTries>0) if(!((++noTries)%m_noTries)) report(DEBUG,"EvtGen") << " EvtHypNonLepton already finished " << noTries << " matrix element calculations" << std::endl; | |
151 | } |