evtgen is hosted by Hepforge, IPPP Durham
EvtGen  2.0.0
Monte Carlo generator of particle decays, in particular the weak decays of heavy flavour particles such as B mesons.
EvtSVVNONCPEIGEN.cpp
Go to the documentation of this file.
1 
2 /***********************************************************************
3 * Copyright 1998-2020 CERN for the benefit of the EvtGen authors *
4 * *
5 * This file is part of EvtGen. *
6 * *
7 * EvtGen is free software: you can redistribute it and/or modify *
8 * it under the terms of the GNU General Public License as published by *
9 * the Free Software Foundation, either version 3 of the License, or *
10 * (at your option) any later version. *
11 * *
12 * EvtGen is distributed in the hope that it will be useful, *
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15 * GNU General Public License for more details. *
16 * *
17 * You should have received a copy of the GNU General Public License *
18 * along with EvtGen. If not, see <https://www.gnu.org/licenses/>. *
19 ***********************************************************************/
20 
22 
23 #include "EvtGenBase/EvtCPUtil.hh"
24 #include "EvtGenBase/EvtConst.hh"
25 #include "EvtGenBase/EvtGenKine.hh"
26 #include "EvtGenBase/EvtPDL.hh"
28 #include "EvtGenBase/EvtPatches.hh"
29 #include "EvtGenBase/EvtRandom.hh"
30 #include "EvtGenBase/EvtReport.hh"
32 
34 
35 #include <stdlib.h>
36 #include <string>
37 
39 {
40  return "SVV_NONCPEIGEN";
41 }
42 
44 {
45  return new EvtSVVNONCPEIGEN;
46 }
47 
49 {
50  // check that there are 27 arguments
51  checkNArg( 27, 15 );
52  checkNDaug( 2 );
53 
56 
57  // The ordering of A_f is :
58  // A_f[0-2] = A_f
59  // A_f[3-5] = Abar_f
60  // A_f[6-8] = A_fbar
61  // A_f[9-11] = Abar_fbar
62  //
63  // Each of the 4 amplitudes include the 3 different helicity states in
64  // the order +, 0, -. See more about helicity amplitude ordering in ::decay
65 
66  int i = 0;
67  int j = ( getNArg() - 3 ) / 2;
68 
69  for ( i = 0; i < j; ++i ) {
70  _A_f[i] = getArg( ( 2 * i ) + 3 ) *
71  EvtComplex( cos( getArg( ( 2 * i ) + 4 ) ),
72  sin( getArg( ( 2 * i ) + 4 ) ) );
73  }
74 
75  // If only 6 amplitudes are specified, calculate the last 6 from the first 6:
76  if ( 6 == j ) {
77  for ( i = 0; i < 3; ++i ) {
78  _A_f[6 + i] = _A_f[3 + i];
79  _A_f[9 + i] = _A_f[i];
80  }
81  }
82 }
83 
85 {
86  double probMax = 0;
87  for ( int i = 0; i < 12; ++i ) {
88  double amp = abs( _A_f[i] );
89  probMax += amp * amp;
90  }
91 
92  setProbMax( probMax );
93 }
94 
96 {
97  //added by Lange Jan4,2000
98  static EvtId B0 = EvtPDL::getId( "B0" );
99  static EvtId B0B = EvtPDL::getId( "anti-B0" );
100 
101  double t;
102  EvtId other_b;
103  EvtId daugs[2];
104 
105  // MB: flip selects the final of the decay
106  int flip = ( ( p->getId() == B0 ) ? 0 : 1 );
107  daugs[0] = getDaug( 0 );
108  daugs[1] = getDaug( 1 );
109  p->initializePhaseSpace( 2, daugs );
110 
111  EvtCPUtil::getInstance()->OtherB( p, t, other_b, 0.5 );
112 
113  EvtComplex amp[3];
114 
115  double dmt2 = getArg( 0 ) * t / ( 2 * EvtConst::c );
116  double phiCKM = ( 2.0 * getArg( 1 ) + getArg( 2 ) ); // 2b+g
117  EvtComplex ePlusIPhi( cos( phiCKM ), sin( phiCKM ) );
118  EvtComplex eMinusIPhi( cos( -phiCKM ), sin( -phiCKM ) );
119 
120  // flip == 0 : D*-rho+
121  // flip == 1 : D*+rho-
122 
123  if ( !flip ) {
124  if ( other_b == B0B ) {
125  // At t=0 we have a B0
126  for ( int i = 0; i < 3; ++i ) {
127  amp[i] = _A_f[i] * cos( dmt2 ) +
128  eMinusIPhi * EvtComplex( 0.0, sin( dmt2 ) ) *
129  _A_f[i + 3];
130  }
131  }
132  if ( other_b == B0 ) {
133  // At t=0 we have a B0bar
134  for ( int i = 0; i < 3; ++i ) {
135  amp[i] = _A_f[i] * ePlusIPhi * EvtComplex( 0.0, sin( dmt2 ) ) +
136  _A_f[i + 3] * cos( dmt2 );
137  }
138  }
139  } else {
140  if ( other_b == B0B ) {
141  // At t=0 we have a B0
142 
143  // M.Baak 01/16/2004
144  // Note: \bar{H}+- = H-+
145  // If one wants to use the correct helicities for B0 and B0bar decays but the same formula-notation (as done in EvtSVV_HelAmp),
146  // count the B0bar helicities backwards. (Equivalently, one could flip the chi angle.)
147 
148  for ( int i = 0; i < 3; ++i ) {
149  amp[i] = _A_f[8 - i] * cos( dmt2 ) +
150  eMinusIPhi * EvtComplex( 0.0, sin( dmt2 ) ) *
151  _A_f[11 - i];
152  }
153  }
154  if ( other_b == B0 ) {
155  // At t=0 we have a B0bar
156  for ( int i = 0; i < 3; ++i ) {
157  amp[i] = _A_f[8 - i] * ePlusIPhi * EvtComplex( 0.0, sin( dmt2 ) ) +
158  _A_f[11 - i] * cos( dmt2 );
159  }
160  }
161  }
162 
163  EvtSVVHelAmp::SVVHel( p, _amp2, daugs[0], daugs[1], amp[0], amp[1], amp[2] );
164 
165  return;
166 }
167 
168 std::string EvtSVVNONCPEIGEN::getParamName( int i )
169 {
170  switch ( i ) {
171  case 0:
172  return "deltaM";
173  case 1:
174  return "weakPhase1";
175  case 2:
176  return "weakPhase2";
177  case 3:
178  return "AfPlusHelAmp";
179  case 4:
180  return "AfPlusHelAmpPhase";
181  case 5:
182  return "AfZeroHelAmp";
183  case 6:
184  return "AfZeroHelAmpPhase";
185  case 7:
186  return "AfMinusHelAmp";
187  case 8:
188  return "AfMinusHelAmpPhase";
189  case 9:
190  return "AbarfPlusHelAmp";
191  case 10:
192  return "AbarfPlusHelAmpPhase";
193  case 11:
194  return "AbarfZeroHelAmp";
195  case 12:
196  return "AbarfZeroHelAmpPhase";
197  case 13:
198  return "AbarfMinusHelAmp";
199  case 14:
200  return "AbarfMinusHelAmpPhase";
201  case 15:
202  return "AfbarPlusHelAmp";
203  case 16:
204  return "AfbarPlusHelAmpPhase";
205  case 17:
206  return "AfbarZeroHelAmp";
207  case 18:
208  return "AfbarZeroHelAmpPhase";
209  case 19:
210  return "AfbarMinusHelAmp";
211  case 20:
212  return "AfbarMinusHelAmpPhase";
213  case 21:
214  return "AbarfbarPlusHelAmp";
215  case 22:
216  return "AbarfbarPlusHelAmpPhase";
217  case 23:
218  return "AbarfbarZeroHelAmp";
219  case 24:
220  return "AbarfbarZeroHelAmpPhase";
221  case 25:
222  return "AbarfbarMinusHelAmp";
223  case 26:
224  return "AbarfbarMinusHelAmpPhase";
225  default:
226  return "";
227  }
228 }
229 
231 {
232  switch ( i ) {
233  case 3:
234  return "1.0";
235  case 4:
236  return "0.0";
237  case 5:
238  return "1.0";
239  case 6:
240  return "0.0";
241  case 7:
242  return "1.0";
243  case 8:
244  return "0.0";
245  case 9:
246  return "1.0";
247  case 10:
248  return "0.0";
249  case 11:
250  return "1.0";
251  case 12:
252  return "0.0";
253  case 13:
254  return "1.0";
255  case 14:
256  return "0.0";
257  default:
258  return "";
259  }
260 }
std::string getName() override
double getArg(unsigned int j)
static void SVVHel(EvtParticle *parent, EvtAmp &amp, EvtId n_v1, EvtId n_v2, const EvtComplex &hp, const EvtComplex &h0, const EvtComplex &hm)
std::string getParamDefault(int i) override
void initProbMax() override
EvtId getId() const
void decay(EvtParticle *p) override
void OtherB(EvtParticle *p, double &t, EvtId &otherb)
Definition: EvtCPUtil.cpp:372
void setProbMax(double prbmx)
std::string getParamName(int i) override
Definition: EvtId.hh:27
double initializePhaseSpace(unsigned int numdaughter, EvtId *daughters, bool forceResetMasses=false, double poleSize=-1., int whichTwo1=0, int whichTwo2=1)
void checkNDaug(int d1, int d2=-1)
void checkNArg(int a1, int a2=-1, int a3=-1, int a4=-1)
static const double c
Definition: EvtConst.hh:30
double abs(const EvtComplex &c)
Definition: EvtComplex.hh:201
static EvtId getId(const std::string &name)
Definition: EvtPDL.cpp:287
void checkSpinDaughter(int d1, EvtSpinType::spintype sp)
EvtDecayBase * clone() override
EvtComplex _A_f[12]
static EvtCPUtil * getInstance()
Definition: EvtCPUtil.cpp:43
EvtAmp _amp2
Definition: EvtDecayAmp.hh:73
void init() override
int getNArg() const
Definition: EvtDecayBase.hh:68
EvtId getDaug(int i) const
Definition: EvtDecayBase.hh:67