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