EvtGen 2.2.0
Monte Carlo generator of particle decays, in particular the weak decays of heavy flavour particles such as B mesons.
Loading...
Searching...
No Matches
EvtLb2BaryonlnuFF.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
24#include "EvtGenBase/EvtId.hh"
26#include "EvtGenBase/EvtPDL.hh"
28
29#include <math.h>
30#include <stdlib.h>
31#include <string>
32using std::endl;
33
34void EvtLb2BaryonlnuFF::getdiracff( EvtId parent, EvtId daught, double q2,
35 double /* mass */, double* f1, double* f2,
36 double* f3, double* g1, double* g2,
37 double* g3 )
38{
39 // Define Event IDs for Lb and p, N+ and Lc+ states
40 static const EvtId LAMB = EvtPDL::getId( "Lambda_b0" );
41 static const EvtId LAMBB = EvtPDL::getId( "anti-Lambda_b0" );
42 static const EvtId PRO = EvtPDL::getId( "p+" );
43 static const EvtId PROB = EvtPDL::getId( "anti-p-" );
44 static const EvtId N1440 = EvtPDL::getId( "N(1440)+" );
45 static const EvtId N1440B = EvtPDL::getId( "anti-N(1440)-" );
46 static const EvtId N1535 = EvtPDL::getId( "N(1535)+" );
47 static const EvtId N1535B = EvtPDL::getId( "anti-N(1535)-" );
48 static const EvtId N1650 = EvtPDL::getId( "N(1650)+" );
49 static const EvtId N1650B = EvtPDL::getId( "anti-N(1650)-" );
50 static const EvtId N1710 = EvtPDL::getId( "N(1710)+" );
51 static const EvtId N1710B = EvtPDL::getId( "anti-N(1710)-" );
52 static const EvtId LAMCP = EvtPDL::getId( "Lambda_c+" );
53 static const EvtId LAMCM = EvtPDL::getId( "anti-Lambda_c-" );
54 static const EvtId LAMC1P = EvtPDL::getId( "Lambda_c(2593)+" );
55 static const EvtId LAMC1M = EvtPDL::getId( "anti-Lambda_c(2593)-" );
56
57 double F1, F2, F3, G1, G2, G3;
58
59 if ( ( parent == LAMB && daught == PRO ) ||
60 ( parent == LAMBB && daught == PROB ) ||
61 ( parent == LAMB && daught == LAMCP ) ||
62 ( parent == LAMBB && daught == LAMCM ) ) {
63 // Parameters needed in the calculation;
64 double mQ = 5.28;
65 double aL = 0.59;
66 double md = 0.40;
67 double MLamB = EvtPDL::getMass( parent );
68 double MLamq = EvtPDL::getMass( daught );
69 double mq = md;
70 double aLp = 0.48;
71
72 //set mq and aLp based on whether Lb->Lc* or Lb->N*
73 if ( ( parent == LAMB && daught == LAMCP ) ||
74 ( parent == LAMBB && daught == LAMCM ) ) {
75 mq = 1.89;
76 aLp = 0.55;
77 }
78
79 double aL2 = aL * aL;
80 double aLp2 = aLp * aLp;
81 double aLLp2 = 0.5 * ( aL2 + aLp2 );
82
83 // relativistic correction factor
84 double k2 = 1.0;
85 double rho2 = 3. * md * md / ( 2. * k2 * aLLp2 );
86
87 // w = scalar product of the 4 velocities of the Lb and Lc.
88 double w = 0.5 * ( MLamB * MLamB + MLamq * MLamq - q2 ) / MLamB / MLamq;
89
90 double I = pow( aL * aLp / aLLp2, 1.5 ) * exp( -rho2 * ( w * w - 1. ) );
91
92 // Calculate the form factors
93 F1 = I * ( 1.0 + ( md / aLLp2 ) * ( ( aLp2 / mq ) + ( aL2 / mQ ) ) );
94 F2 = -I * ( ( md / mq ) * ( aLp2 / aLLp2 ) -
95 aL2 * aLp2 / ( 4. * aLLp2 * mq * mQ ) );
96 F3 = -I * md * aL2 / ( mQ * aLLp2 );
97
98 G1 = I * ( 1.0 - ( aL2 * aLp2 ) / ( 12. * aLLp2 * mq * mQ ) );
99 G2 = -I * ( md * aLp2 / ( mq * aLLp2 ) +
100 ( aL2 * aLp2 ) / ( 12. * aLLp2 * mq * mQ ) *
101 ( 1. + 12. * md * md / aLLp2 ) );
102 G3 = I * ( md * aL2 / ( mQ * aLLp2 ) +
103 md * md * aL2 * aLp2 / ( mq * mQ * aLLp2 * aLLp2 ) );
104
105 // Set form factors to be passed to the amplitude calc.
106 *f1 = F1;
107 *f2 = F2;
108 *f3 = F3;
109 *g1 = G1;
110 *g2 = G2;
111 *g3 = G3;
112 } else if ( ( parent == LAMB && daught == N1440 ) ||
113 ( parent == LAMBB && daught == N1440B ) ||
114 ( parent == LAMB && daught == N1710 ) ||
115 ( parent == LAMBB && daught == N1710B ) ) {
116 // Parameters needed in the calculation;
117 double mQ = 5.28;
118 double md = 0.40;
119 double mq = md;
120 double MLamB = EvtPDL::getMass( parent );
121 double MLamq = EvtPDL::getMass( daught );
122
123 double aL = 0.59;
124 double aLp = 0.48;
125
126 double aL2 = aL * aL;
127 double aLp2 = aLp * aLp;
128 double aLLp2 = 0.5 * ( aL2 + aLp2 );
129
130 // relativistic correction factor
131 double k2 = 1.0;
132 double rho2 = 3. * md * md / ( 2. * k2 * aLLp2 );
133
134 // w = scalar product of the 4 velocities of the Lb and Lc.
135 double w = 0.5 * ( MLamB * MLamB + MLamq * MLamq - q2 ) / MLamB / MLamq;
136
137 double I = sqrt( 1.5 ) * pow( aL * aLp / aLLp2, 1.5 ) *
138 exp( -rho2 * ( w * w - 1. ) );
139
140 // Calculate the form factors
141 F1 = ( I / ( 2. * aLLp2 ) ) *
142 ( ( aL2 - aLp2 ) - ( md / ( 3. * aLLp2 ) ) *
143 ( ( aLp2 / mq ) * ( 7. * aL2 - 3. * aLp2 ) +
144 ( aL2 / mQ ) * ( 7. * aLp2 - 3. * aL2 ) ) );
145
146 F2 = -I * ( aLp2 / ( 6. * mq * aLLp2 * aLLp2 ) ) *
147 ( 7. * aL2 - 3. * aLp2 ) * ( md - ( aL2 / ( 4. * mQ ) ) );
148
149 F3 = I * ( aL2 * md / ( 6. * mQ * aLLp2 * aLLp2 ) ) *
150 ( 7. * aLp2 - 3. * aL2 );
151
152 G1 = I * ( ( aL2 - aLp2 ) / ( 2 * aLLp2 ) -
153 ( aL2 * aLp2 * ( 7. * aL2 - 3. * aLp2 ) ) /
154 ( 72. * aLLp2 * aLLp2 * mq * mQ ) );
155
156 G2 = -I * ( aLp2 / ( 6. * mq * aLLp2 * aLLp2 ) ) *
157 ( ( 7. * aL2 - 3. * aLp2 ) * ( md + ( aL2 / ( 6. * mQ ) ) ) +
158 ( 7. * md * md * aL2 * ( aL2 - aLp2 ) / ( mQ * aLLp2 ) ) );
159
160 G3 = -I * ( aL2 * md / ( 6. * mQ * aLLp2 * aLLp2 ) ) *
161 ( ( 7. * aLp2 - 3. * aL2 ) -
162 ( 7 * md * aLp2 * ( aL2 - aLp2 ) / ( mq * aLLp2 ) ) );
163
164 // Set form factors to be passed to the amplitude calc.
165 *f1 = F1;
166 *f2 = F2;
167 *f3 = F3;
168 *g1 = G1;
169 *g2 = G2;
170 *g3 = G3;
171
172 } else if ( ( parent == LAMB && daught == N1535 ) ||
173 ( parent == LAMBB && daught == N1535B ) ||
174 ( parent == LAMB && daught == N1650 ) ||
175 ( parent == LAMBB && daught == N1650B ) ||
176 ( parent == LAMB && daught == LAMC1P ) ||
177 ( parent == LAMBB && daught == LAMC1M ) ) {
178 double mQ = 5.28;
179 double md = 0.40;
180 double aL = 0.59;
181 double MLamB = EvtPDL::getMass( parent );
182 double MLamq = EvtPDL::getMass( daught );
183 double mq = md;
184 double aLp = 0.37;
185
186 //set mq and aLp based on whether Lb->Lc* or Lb->N*
187 if ( ( parent == LAMB && daught == LAMC1P ) ||
188 ( parent == LAMBB && daught == LAMC1M ) ) {
189 mq = 1.89;
190 aLp = 0.47;
191 }
192
193 double aL2 = aL * aL;
194 double aLp2 = aLp * aLp;
195 double aLLp2 = 0.5 * ( aL2 + aLp2 );
196
197 // relativistic correction factor
198 double k2 = 1.0;
199 double rho2 = 3. * md * md / ( 2. * k2 * aLLp2 );
200
201 // w = scalar product of the 4 velocities of the Lb and Lc.
202 double w = 0.5 * ( MLamB * MLamB + MLamq * MLamq - q2 ) / MLamB / MLamq;
203
204 double I = pow( aL * aLp / aLLp2, 2.5 ) * exp( -rho2 * ( w * w - 1. ) );
205
206 // Calculate the form factors
207 F1 = I * aL / 6.0 * ( 3.0 / mq - 1.0 / mQ );
208 F2 = -I * ( 2.0 * md / aL - aL / ( 2.0 * mq ) +
209 2. * md * md * aL / ( mQ * aLLp2 ) -
210 ( md * aL / ( 6. * mq * mQ * aLLp2 ) ) *
211 ( 3. * aL2 - 2. * aLp2 ) );
212 F3 = I * 2. * md * md * aL / ( mQ * aLLp2 );
213
214 G1 = I * ( 2.0 * md / aL - aL / ( 6. * mQ ) +
215 ( md * aL / ( 6. * mq * mQ * aLLp2 ) ) *
216 ( 3. * aL2 - 2. * aLp2 ) );
217 G2 = I * ( -2. * md / aL + aL / ( 2. * mq ) + aL / ( 3. * mQ ) );
218 G3 = I * aL / ( 3. * mQ ) *
219 ( 1.0 - ( md / ( 2. * mq * aLLp2 ) ) * ( 3. * aL2 - 2. * aLp2 ) );
220
221 // Set form factors to be passed to the amplitude calc.
222 *f1 = F1;
223 *f2 = F2;
224 *f3 = F3;
225 *g1 = G1;
226 *g2 = G2;
227 *g3 = G3;
228 } else {
229 EvtGenReport( EVTGEN_ERROR, "EvtGen" )
230 << "Only Lb -> N*+ transitions allowed in EvtLb2BaryonlnuFF.\n";
231 ::abort();
232 }
233
234 return;
235}
236
237void EvtLb2BaryonlnuFF::getraritaff( EvtId parent, EvtId daught, double q2,
238 double, double* f1, double* f2, double* f3,
239 double* f4, double* g1, double* g2,
240 double* g3, double* g4 )
241{
242 static const EvtId LAMB = EvtPDL::getId( "Lambda_b0" );
243 static const EvtId LAMBB = EvtPDL::getId( "anti-Lambda_b0" );
244 static const EvtId N1520 = EvtPDL::getId( "N(1520)+" );
245 static const EvtId N1520B = EvtPDL::getId( "anti-N(1520)-" );
246 static const EvtId N1720 = EvtPDL::getId( "N(1720)+" );
247 static const EvtId N1720B = EvtPDL::getId( "anti-N(1720)-" );
248 static const EvtId N1700 = EvtPDL::getId( "N(1700)+" );
249 static const EvtId N1700B = EvtPDL::getId( "anti-N(1700)-" );
250 static const EvtId N1900 = EvtPDL::getId( "N(1900)+" );
251 static const EvtId N1900B = EvtPDL::getId( "anti-N(1900)-" );
252 static const EvtId N1875 = EvtPDL::getId( "N(1875)+" );
253 static const EvtId N1875B = EvtPDL::getId( "anti-N(1875)-" );
254 static const EvtId LAMC2P = EvtPDL::getId( "Lambda_c(2625)+" );
255 static const EvtId LAMC2M = EvtPDL::getId( "anti-Lambda_c(2625)-" );
256
257 double F1, F2, F3, F4, G1, G2, G3, G4;
258
259 // 3/2 - case
260 if ( ( parent == LAMB && daught == N1520 ) ||
261 ( parent == LAMBB && daught == N1520B ) ||
262 ( parent == LAMB && daught == N1700 ) ||
263 ( parent == LAMBB && daught == N1700B ) ||
264 ( parent == LAMB && daught == N1875 ) ||
265 ( parent == LAMBB && daught == N1875B ) ||
266 ( parent == LAMB && daught == LAMC2P ) ||
267 ( parent == LAMBB && daught == LAMC2M ) ) {
268 double mQ = 5.28;
269 double md = 0.40;
270 double aL = 0.59;
271 double MLamB = EvtPDL::getMass( parent );
272 double MLamq = EvtPDL::getMass( daught );
273 double mq = md;
274 double aLp = 0.37;
275
276 //set mq and aLp based on whether Lb->Lc* or Lb->N*
277 if ( ( parent == LAMB && daught == LAMC2P ) ||
278 ( parent == LAMBB && daught == LAMC2M ) ) {
279 mq = 1.89;
280 aLp = 0.47;
281 }
282
283 double aL2 = aL * aL;
284 double aLp2 = aLp * aLp;
285 double aLLp2 = 0.5 * ( aL2 + aLp2 );
286
287 // relativistic correction factor
288 double k2 = 1.0;
289 double rho2 = 3. * md * md / ( 2. * k2 * aLLp2 );
290
291 // w = scalar product of the 4 velocities of the Lb and Lc.
292 double w = 0.5 * ( MLamB * MLamB + MLamq * MLamq - q2 ) / MLamB / MLamq;
293
294 double I = -( 1. / sqrt( 3. ) ) * pow( aL * aLp / aLLp2, 2.5 ) *
295 exp( -rho2 * ( w * w - 1. ) );
296
297 // Calculate the form factors
298 F1 = I * 3.0 * md / aL *
299 ( 1.0 + ( md / aLLp2 ) * ( ( aLp2 / mq ) + ( aL2 / mQ ) ) );
300 F2 = -I * ( ( 3. * md * md / mq ) * ( aLp2 / ( aLLp2 * aL2 ) ) -
301 5. * aL * aLp2 * md / ( 4. * aLLp2 * mq * mQ ) );
302 F3 = -I * ( 3. * md * md * aL / ( mQ * aLLp2 ) + aL / ( 2. * mQ ) );
303 F4 = I * aL / mQ;
304
305 G1 = I * ( 3.0 * md / aL -
306 ( aL / ( 2. * mQ ) ) *
307 ( 1. + 3. * md * aLp2 / ( 2. * aLLp2 * mq ) ) );
308 G2 = -I * ( ( 3. * md * md / mq ) * ( aLp2 / ( aLLp2 * aL ) ) +
309 aL * aLp2 * md / ( 4. * aLLp2 * aLLp2 * mq * mQ ) *
310 ( aLLp2 + 12. * md * md ) );
311 G3 = I * aL / ( mQ * aLLp2 ) *
312 ( aLLp2 / 2. + 3. * md * md +
313 aLp2 * md / ( mq * aLLp2 ) * ( aLLp2 + 6. * md * md ) );
314 G4 = -I * ( aL / mQ + md / ( mq * mQ ) * aLp2 * aL / aLLp2 );
315
316 // Set form factors to be passed to the amplitude calc.
317 *f1 = F1;
318 *f2 = F2;
319 *f3 = F3;
320 *f4 = F4;
321 *g1 = G1;
322 *g2 = G2;
323 *g3 = G3;
324 *g4 = G4;
325 }
326 // 3/2 + case
327 else if ( ( parent == LAMB && daught == N1720 ) ||
328 ( parent == LAMBB && daught == N1720B ) ||
329 ( parent == LAMB && daught == N1900 ) ||
330 ( parent == LAMBB && daught == N1900B )
331
332 ) {
333 double mQ = 5.28;
334 double md = 0.40;
335 double mq = md;
336 double MLamB = EvtPDL::getMass( parent );
337 double MLamq = EvtPDL::getMass( daught );
338
339 double aL = 0.59;
340 double aLp = 0.35;
341
342 double aL2 = aL * aL;
343 double aLp2 = aLp * aLp;
344 double aLLp2 = 0.5 * ( aL2 + aLp2 );
345
346 // relativistic correction factor
347 double k2 = 1.0;
348 double rho2 = 3. * md * md / ( 2. * k2 * aLLp2 );
349
350 // w = scalar product of the 4 velocities of the Lb and Lc.
351 double w = 0.5 * ( MLamB * MLamB + MLamq * MLamq - q2 ) / MLamB / MLamq;
352
353 double I = ( 1. / sqrt( 5. ) ) * pow( aL * aLp / aLLp2, 3.5 ) *
354 exp( -rho2 * ( w * w - 1. ) );
355
356 // Calculate the form factors
357 F1 = -I * ( md / 2. ) * ( ( 5. / mq ) - ( 3. / mQ ) );
358
359 F2 = I * ( md / aL ) *
360 ( ( 6. * md / aL ) - ( 5 * aL / ( 2. * mq ) ) +
361 ( 6. * md * md * aL ) / ( aLLp2 * mQ ) -
362 ( md * aL * ( aL2 - 2. * aLp2 ) ) / ( 2 * aLLp2 * mq * mQ ) );
363
364 F3 = -I * ( md / mQ ) * ( 1 + ( 6. * md * md ) / aLLp2 );
365
366 F4 = I * ( 2. * md / mQ );
367
368 G1 = -I *
369 ( ( 6. * md * md / aL2 ) - md / ( 2. * mQ ) +
370 md * md * ( 11. * aL2 - 6. * aLp2 ) / ( 6. * aLLp2 * mq * mQ ) );
371
372 G2 = I * ( 6 * md * md / aL2 - 5 * md / ( 2.0 * mq ) - ( 2 * md ) / mQ +
373 5 * aL2 / ( 12.0 * mq * mQ ) -
374 ( 2 * md * md * aL2 ) / ( 3.0 * aLLp2 * mq * mQ ) );
375
376 G3 = -I *
377 ( ( md / ( 2. * mQ ) ) - 5 * aL2 / ( 24.0 * mq * mQ ) -
378 md * md * ( 5 * aL2 - 2 * aLp2 ) / ( 4.0 * mq * mQ * aLLp2 ) );
379
380 G4 = -I * 5. * aL2 / ( 6. * mq * mQ );
381
382 // Set form factors to be passed to the amplitude calc.
383 *f1 = F1;
384 *f2 = F2;
385 *f3 = F3;
386 *f4 = F4;
387 *g1 = G1;
388 *g2 = G2;
389 *g3 = G3;
390 *g4 = G4;
391 }
392
393 else {
394 EvtGenReport( EVTGEN_ERROR, "EvtGen" )
395 << "Only Lb -> N*+ transitions allowed in EvtLb2BaryonlnuFF.\n";
396 ::abort();
397 }
398
399 return;
400}
401
402void EvtLb2BaryonlnuFF::getscalarff( EvtId, EvtId, double, double, double*,
403 double* )
404{
405 EvtGenReport( EVTGEN_ERROR, "EvtGen" )
406 << "Not implemented :getscalarff in EvtLb2BaryonlnuFF.\n";
407 ::abort();
408}
409
410void EvtLb2BaryonlnuFF::getvectorff( EvtId, EvtId, double, double, double*,
411 double*, double*, double* )
412{
413 EvtGenReport( EVTGEN_ERROR, "EvtGen" )
414 << "Not implemented :getvectorff in EvtLb2BaryonlnuFF.\n";
415 ::abort();
416}
417
418void EvtLb2BaryonlnuFF::gettensorff( EvtId, EvtId, double, double, double*,
419 double*, double*, double* )
420{
421 EvtGenReport( EVTGEN_ERROR, "EvtGen" )
422 << "Not implemented :gettensorff in EvtLb2BaryonlnuFF.\n";
423 ::abort();
424}
425
426void EvtLb2BaryonlnuFF::getbaryonff( EvtId, EvtId, double, double, double*,
427 double*, double*, double* )
428{
429 EvtGenReport( EVTGEN_ERROR, "EvtGen" )
430 << "Not implemented :getbaryonff in EvtLb2BaryonlnuFF.\n";
431 ::abort();
432}
const EvtComplex I
EvtComplex exp(const EvtComplex &c)
std::ostream & EvtGenReport(EvtGenSeverity severity, const char *facility=nullptr)
Definition EvtReport.cpp:32
@ EVTGEN_ERROR
Definition EvtReport.hh:49
Definition EvtId.hh:27
void gettensorff(EvtId parent, EvtId daught, double t, double mass, double *hf, double *kf, double *bpf, double *bmf) override
void getraritaff(EvtId parent, EvtId daught, double q2, double mass, double *f1, double *f2, double *f3, double *f4, double *g1, double *g2, double *g3, double *g4) override
void getvectorff(EvtId parent, EvtId daught, double t, double mass, double *a1f, double *a2f, double *vf, double *a0f) override
void getdiracff(EvtId parent, EvtId daught, double q2, double mass, double *f1, double *f2, double *f3, double *g1, double *g2, double *g3) override
void getscalarff(EvtId parent, EvtId daught, double t, double mass, double *fpf, double *f0f) override
void getbaryonff(EvtId, EvtId, double, double, double *, double *, double *, double *) override
static double getMass(EvtId i)
Definition EvtPDL.cpp:311
static EvtId getId(const std::string &name)
Definition EvtPDL.cpp:283