]> git.uio.no Git - u/mrichter/AliRoot.git/blob - TEvtGen/EvtGenModels/EvtSVVNONCPEIGEN.cxx
AliDecayer realisation for the EvtGen code and EvtGen itself.
[u/mrichter/AliRoot.git] / TEvtGen / EvtGenModels / EvtSVVNONCPEIGEN.cxx
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) 2001      Royal Holloway, University of London
10 //
11 // Module: EvtSVVNONCPEIGEN.cc
12 //
13 // Description: Routine to decay scalar -> vector vector 
14 //              and has CP violation.
15 //
16 //              This model does all the ckm-suppressed decays and mixing for you. It randomly 'overwrites' 
17 //              any reco or tagging state as set in the Y(4S) decay model (VSS_(B)MIX) with its own generated states.
18 //
19 //              As such, the corresponding dec file requires only one decay-mode description, for example:
20 //              Decay MyB0
21 //              1.000    rho+ MyD*-       SVV_NONCPEIGEN dm beta gamma 0.322 0.31 0.941 0 0.107 1.42 0.02 0 0.02 0 0.02 0 ;
22 //              EndDecay
23 //              and furthermore Y(4S) only needs to decay to B0's (or B0bar's).
24 //              The decay above should be a CKM-favored mode (eg. B0->D*-rho+ or B0bar->D*+rho-).
25 //              All ckm-suppressed decays and the mixing are derived from this line in the ::Decay function.
26 //
27 //              There are 15 or 27 arguments. The first three are dm, phase1
28 //              and phase2. dm is the B0-B0bar mass difference. Phases 1
29 //              and 2 are the CKM weak phases relevant for the particular mode, 
30 //              eg for B-->DstRho phase1 is beta and phase2 is gamma.
31 //
32 //              The next arguments are the 2 amplitudes (= 12 input parameters) 
33 //              in the order: A_f, Abar_f. In the example above, the 'A_f' amplitude now 
34 //              stands for the ckm-favored decay 'B0->D*-rho+', and 'Abar_f' stands for 'B0bar->D*-rho+'
35 //
36 //              Each amplitude has its 3 helicity states in the order +, 0, -, which are each 
37 //              specified by a magnitude and a strong phase.
38 //
39 //              The last 2 arguments A_fbar and Abar_fbar (=12 input parameters) are not necessary, 
40 //              but can included if one wants to set them differently from A_f, Abar_f.
41 //
42 //              Mind you that Hbar_+- = H_-+ (ignoring the weak phase, which flips sign).
43 //              It is custumary to select one set of helicity states (eg H_+-) and to adopt these for
44 //              the CP-conjugate decays as well (ie. depict Hbar_-+ with H_+-), which is the interpretation
45 //              we use for the input-parameters above. 
46 //              However, the angular decay in EvtGen is just a formula in which helicity amplitudes are 'plugged' in,
47 //              making no difference between B0 or B0bar decays. In the model below we (thus) account for the +- 
48 //              flipping between B0 and B0bar.
49 //              
50 //
51 // Modification history:
52 //    Ajit Kurup 9 March 2001        Module created (from EvtSVSNONCPEIGEN)
53 //    Max Baak 01/16/2004            Fix of Helicity amplitude ordering.
54 //                                   Decay also works for B0bar decays.
55 //------------------------------------------------------------------------
56 // 
57 #include "EvtGenBase/EvtPatches.hh"
58 #include <stdlib.h>
59 #include "EvtGenBase/EvtParticle.hh"
60 #include "EvtGenBase/EvtRandom.hh"
61 #include "EvtGenBase/EvtGenKine.hh"
62 #include "EvtGenBase/EvtCPUtil.hh"
63 #include "EvtGenBase/EvtPDL.hh"
64 #include "EvtGenBase/EvtReport.hh"
65 #include "EvtGenBase/EvtVector4C.hh"
66 #include "EvtGenModels/EvtSVVNONCPEIGEN.hh"
67 #include <string>
68 #include "EvtGenModels/EvtSVVHelAmp.hh"
69 #include "EvtGenBase/EvtConst.hh"
70
71 EvtSVVNONCPEIGEN::~EvtSVVNONCPEIGEN() {}
72
73 std::string EvtSVVNONCPEIGEN::getName(){
74
75   return "SVV_NONCPEIGEN";     
76
77 }
78
79
80 EvtDecayBase* EvtSVVNONCPEIGEN::clone(){
81
82   return new EvtSVVNONCPEIGEN;
83
84 }
85
86 void EvtSVVNONCPEIGEN::init(){
87
88   // check that there are 27 arguments
89   checkNArg(27,15);
90   checkNDaug(2);
91
92   checkSpinDaughter(0,EvtSpinType::VECTOR);
93   checkSpinDaughter(1,EvtSpinType::VECTOR);
94
95   //  The ordering of A_f is :
96   //  A_f[0-2] = A_f
97   //  A_f[3-5] = Abar_f
98   //  A_f[6-8] = A_fbar 
99   //  A_f[9-11] = Abar_fbar
100   //  
101   //  Each of the 4 amplitudes include the 3 different helicity states in 
102   //  the order +, 0, -. See more about helicity amplitude ordering in ::decay
103
104   int i=0;
105   int j=(getNArg()-3)/2;
106
107   for(i=0; i<j; ++i){
108     _A_f[i] = getArg((2*i)+3) * EvtComplex( cos(getArg((2*i)+4)),sin(getArg((2*i)+4)) );
109   }
110
111   //  If only 6 amplitudes are specified, calculate the last 6 from the first 6:
112   if(6 == j){
113     for(i = 0; i < 3; ++i){
114       _A_f[6+i] = _A_f[3+i];
115       _A_f[9+i] = _A_f[i];
116     }
117   }
118 }
119
120 void EvtSVVNONCPEIGEN::initProbMax() {
121   double probMax = 0;
122   for (int i = 0; i < 12; ++i){
123     double amp = abs(_A_f[i]);
124     probMax += amp * amp;
125   }
126
127   setProbMax(probMax); 
128 }
129
130 void EvtSVVNONCPEIGEN::decay( EvtParticle *p){
131
132   //added by Lange Jan4,2000
133   static EvtId B0=EvtPDL::getId("B0");
134   static EvtId B0B=EvtPDL::getId("anti-B0");
135
136   double t;
137   EvtId other_b;
138   EvtId daugs[2];
139
140
141   // MB: flip selects the final of the decay
142   int flip = ((p->getId() == B0) ? 0 : 1);
143   daugs[0]=getDaug(0);
144   daugs[1]=getDaug(1);
145   p->initializePhaseSpace(2,daugs);
146
147   EvtCPUtil::OtherB(p,t,other_b,0.5);
148
149   EvtComplex amp[3];
150
151   double dmt2 = getArg(0) * t / (2 * EvtConst::c);
152   double phiCKM = (2.0 * getArg(1) + getArg(2));   // 2b+g
153   EvtComplex ePlusIPhi(cos(phiCKM), sin(phiCKM));
154   EvtComplex eMinusIPhi(cos(-phiCKM), sin(-phiCKM));
155
156   // flip == 0 : D*-rho+
157   // flip == 1 : D*+rho-
158
159   if (!flip) {
160     if (other_b==B0B){
161       // At t=0 we have a B0
162       for (int i=0; i<3; ++i) {
163         amp[i] = _A_f[i]*cos(dmt2) + eMinusIPhi*EvtComplex(0.0,sin(dmt2))*_A_f[i+3];
164       }
165     }
166     if (other_b==B0){
167       // At t=0 we have a B0bar
168       for(int i=0; i<3; ++i) {
169         amp[i] = _A_f[i]*ePlusIPhi*EvtComplex(0.0,sin(dmt2)) + _A_f[i+3]*cos(dmt2);
170       }
171     }
172   } else{
173     if (other_b==B0B){
174       // At t=0 we have a B0
175
176       // M.Baak 01/16/2004
177       // Note: \bar{H}+- = H-+ 
178       // If one wants to use the correct helicities for B0 and B0bar decays but the same formula-notation (as done in EvtSVV_HelAmp), 
179       // count the B0bar helicities backwards. (Equivalently, one could flip the chi angle.)
180
181       for(int i=0; i<3; ++i) { 
182         amp[i] = _A_f[8-i]*cos(dmt2) + eMinusIPhi*EvtComplex(0.0,sin(dmt2))*_A_f[11-i];
183       }
184     }
185     if (other_b==B0){
186       // At t=0 we have a B0bar
187       for(int i=0; i<3; ++i) {
188         amp[i] = _A_f[8-i] * ePlusIPhi * EvtComplex(0.0,sin(dmt2)) + _A_f[11-i]*cos(dmt2);
189       }
190     }
191   }
192   
193   EvtSVVHelAmp::SVVHel(p,_amp2,daugs[0],daugs[1],amp[0],amp[1],amp[2]);
194
195   return ;
196 }
197