GENIEGenerator
Loading...
Searching...
No Matches
Born.cxx
Go to the documentation of this file.
1//____________________________________________________________________________
2/*
3 Copyright (c) 2003-2025, The GENIE Collaboration
4 For the full text of the license visit http://copyright.genie-mc.org
5
6 Alfonso Garcia <aagarciasoto \at km3net.de>
7 IFIC & Harvard University
8*/
9//____________________________________________________________________________
10
11#include "Physics/HELepton/XSection/Born.h"
12#include "Framework/ParticleData/PDGLibrary.h"
13#include "Framework/ParticleData/PDGCodes.h"
14#include "Framework/Conventions/Constants.h"
15
16#include <TMath.h>
17
18#include <iostream>
19
20using namespace genie;
21using namespace genie::constants;
22
23//____________________________________________________________________________
24Born::Born()
25{
26
27 fGw = PDGLibrary::Instance()->Find(kPdgWM)->Width();
28 fGz = PDGLibrary::Instance()->Find(kPdgZ0)->Width();
29 fmw2c = TComplex(kMw2,fGw*kMw);
30 fmz2c = TComplex(kMz2,fGz*kMz);
31 TComplex rat = fmw2c/fmz2c;
32 fsw2 = TComplex(1.-rat.Re(),-rat.Im());
33 fcw2 = 1.-fsw2;
34 falpha = TMath::Sqrt(2.)*kGF/kPi * fmw2c * fsw2;
35
36 fgae = -1./2. + 2.*fsw2;
37 fgbe = -1./2.;
38 fgav = 1./2.;
39
40}
41//____________________________________________________________________________
42Born::~Born()
43{
44
45/*
46Make sure the p3 is always the charged lepton:
47
48nu (p1) + lp (p2) -> lp (p3) + nu (p4)
49
50 s-channel t-channel u-channel
511 \ / 3 1 --------- 3 1 ------\ / 3
52 ------ | | X
532 / \ 4 2 --------- 4 2 ------/ \ 4
54*/
55
56}
57//____________________________________________________________________________
58double Born::PXSecCCR(double s, double t, double mlin, double mlout)
59{
60
61 TComplex prop = falpha/fsw2/(s-fmw2c);
62
63 return (t-mlout*mlout)*(t-mlin*mlin) * prop.Rho2();
64
65}
66//____________________________________________________________________________
67double Born::PXSecCCV(double s, double t, double mlin, double mlout)
68{
69
70 TComplex prop = falpha/fsw2/(t-fmw2c);
71
72 return (s-mlout*mlout)*(s-mlin*mlin) * prop.Rho2();
73
74}
75//____________________________________________________________________________
76double Born::PXSecCCRNC(double s, double t, double mlin, double mlout)
77{
78
79 double u = GetU(mlin,mlout,s,t);
80
81 TComplex a = fgav*(fgae-fgbe)/(u-fmz2c)/fcw2/fsw2;
82 TComplex b = fgav*(fgae+fgbe)/(u-fmz2c)/fcw2/fsw2 + 1./(s-fmw2c)/fsw2;
83 return falpha.Rho2() * ( (t-mlout*mlout)*(t-mlin*mlin)*b.Rho2() + (s-mlout*mlout)*(s-mlin*mlin)*a.Rho2() );
84
85}
86//____________________________________________________________________________
87double Born::PXSecCCVNC(double s, double t, double mlin, double mlout)
88{
89
90 double u = GetU(mlin,mlout,s,t);
91
92 TComplex a = fgav*(fgae+fgbe)/(u-fmz2c)/fcw2/fsw2 + 1./(t-fmw2c)/fsw2;
93 TComplex b = fgav*(fgae-fgbe)/(u-fmz2c)/fcw2/fsw2;
94 return falpha.Rho2() * ( (t-mlout*mlout)*(t-mlin*mlin)*b.Rho2() + (s-mlout*mlout)*(s-mlin*mlin)*a.Rho2() );
95
96}
97//____________________________________________________________________________
98double Born::PXSecNCVnu(double s, double t, double mlin, double mlout)
99{
100
101 double u = GetU(mlin,mlout,s,t);
102
103 TComplex a = fgav*(fgae+fgbe)/(u-fmz2c)/fcw2/fsw2;
104 TComplex b = fgav*(fgae-fgbe)/(u-fmz2c)/fcw2/fsw2;
105 return falpha.Rho2() * ( (t-mlout*mlout)*(t-mlin*mlin)*b.Rho2() + (s-mlout*mlout)*(s-mlin*mlin)*a.Rho2() );
106
107}
108//____________________________________________________________________________
109double Born::PXSecNCVnubar(double s, double t, double mlin, double mlout)
110{
111
112 double u = GetU(mlin,mlout,s,t);
113
114 TComplex a = fgav*(fgae-fgbe)/(u-fmz2c)/fcw2/fsw2;
115 TComplex b = fgav*(fgae+fgbe)/(u-fmz2c)/fcw2;
116 return falpha.Rho2() * ( (t-mlout*mlout)*(t-mlin*mlin)*b.Rho2()/fsw2.Rho2() + (s-mlout*mlout)*(s-mlin*mlin)*a.Rho2() );
117
118}
119//____________________________________________________________________________
120double Born::PXSecPhoton(double s, double t, double mlout2)
121{
122
123 double u = kMw2 + mlout2 - s - t;
124
125 double ME = 8*kPi2*falpha.Rho2()*s/kMw2/fsw2.Re()/TMath::Power(mlout2 - t,2)/TMath::Power(kMw2 - u,2)*
126 ( -2*(kMw2 - u)*TMath::Power(mlout2,3) - 2*TMath::Power(mlout2,2)*(-2*kMw2*u + u*(s + u) + TMath::Power(kMw2,2))
127 + mlout2*(-(kMw2*u*(4*s + 5*u)) + (s + u)*TMath::Power(kMw2,2) + 3*TMath::Power(kMw2,3) + (s + u)*TMath::Power(u,2))
128 + 2*kMw2*((3*s + u)*TMath::Power(kMw2,2) - TMath::Power(kMw2,3) + 4*u*TMath::Power(s,2) + 2*TMath::Power(s,3) + 3*s*TMath::Power(u,2)
129 + TMath::Power(u,3) - kMw2*TMath::Power(2*s + u,2)) );
130
131 return TMath::Max(0.,ME);
132
133}
134//____________________________________________________________________________
135double Born::PXSecPhoton_T(double s12, double s13, double Q2, double ml2)
136{
137 double ME2 = 0.0;
138 ME2 = (4*falpha.Rho2()*kPi2*(TMath::Power(ml2,4)*s12*(2*TMath::Power(kMw2,2)*TMath::Power(s12,2) - 2*kMw2*Q2*s12*s13 + TMath::Power(Q2,2)*s13*(-s12 + s13)) +
139 TMath::Power(ml2,3)*(-2*TMath::Power(kMw2,3)*TMath::Power(s12,3) + TMath::Power(Q2,2)*(s12 - s13)*s13*(-(Q2*s12) + TMath::Power(s12,2) + Q2*s13 - 3*s12*s13) +
140 2*TMath::Power(kMw2,2)*TMath::Power(s12,2)*(3*Q2*s12 - 2*TMath::Power(s12,2) + 2*Q2*s13 + s12*s13) +
141 kMw2*Q2*s12*(Q2*TMath::Power(s12,2) - Q2*TMath::Power(s13,2) - 2*s12*TMath::Power(s13,2))) +
142 TMath::Power(ml2,2)*(-6*TMath::Power(kMw2,4)*TMath::Power(s12,3) + 2*kMw2*Q2*
143 (-(Q2*s12*(s12 - 3*s13)) + TMath::Power(Q2,2)*(s12 - s13) + TMath::Power(s12,2)*(s12 - s13))*(s12 - s13)*s13 +
144 TMath::Power(Q2,2)*(s12 - s13)*TMath::Power(s13,2)*(-2*Q2*s12 + 2*TMath::Power(s12,2) + 2*Q2*s13 - 3*s12*s13) +
145 2*TMath::Power(kMw2,3)*TMath::Power(s12,2)*(2*TMath::Power(s12,2) - 3*Q2*(2*s12 + s13)) +
146 TMath::Power(kMw2,2)*s12*(2*TMath::Power(s12,2)*TMath::Power(s12 - s13,2) + TMath::Power(Q2,2)*(s12 - s13)*s13 -
147 2*Q2*s12*(TMath::Power(s12,2) - 8*s12*s13 - 2*TMath::Power(s13,2)))) -
148 2*TMath::Power(kMw2,2)*TMath::Power(s12,2)*(2*TMath::Power(kMw2,4)*s12 + TMath::Power(Q2,2)*TMath::Power(s13,2)*(-s12 + s13) -
149 2*TMath::Power(kMw2,3)*(2*TMath::Power(s12,2) - Q2*s13 + s12*s13) +
150 TMath::Power(kMw2,2)*(4*Q2*s12*(s12 - 2*s13) + TMath::Power(Q2,2)*(-s12 + s13) + 2*s12*TMath::Power(s12 + s13,2)) +
151 2*kMw2*s13*(TMath::Power(Q2,2)*(s12 - s13) - s12*(TMath::Power(s12,2) + TMath::Power(s13,2)) + Q2*(-TMath::Power(s12,2) + 2*s12*s13 + TMath::Power(s13,2)))) +
152 ml2*(10*TMath::Power(kMw2,5)*TMath::Power(s12,3) - TMath::Power(Q2,2)*(Q2 - s12)*TMath::Power(s12 - s13,2)*TMath::Power(s13,3) +
153 2*TMath::Power(kMw2,4)*TMath::Power(s12,2)*(3*Q2*s12 - 4*TMath::Power(s12,2) + 4*Q2*s13 - 3*s12*s13) +
154 kMw2*Q2*(s12 - s13)*TMath::Power(s13,2)*(2*TMath::Power(Q2,2)*(s12 - s13) + 2*TMath::Power(s12,2)*(s12 - s13) + Q2*s12*(-3*s12 + 5*s13)) +
155 TMath::Power(kMw2,3)*s12*(2*Q2*s12*(5*TMath::Power(s12,2) - 16*s12*s13 - TMath::Power(s13,2)) +
156 TMath::Power(Q2,2)*(-3*TMath::Power(s12,2) + 2*s12*s13 + TMath::Power(s13,2)) + 2*TMath::Power(s12,2)*(TMath::Power(s12,2) + 6*s12*s13 + TMath::Power(s13,2))) -
157 TMath::Power(kMw2,2)*s13*(TMath::Power(Q2,3)*TMath::Power(s12 - s13,2) - 2*TMath::Power(s12,3)*TMath::Power(s12 - s13,2) +
158 TMath::Power(Q2,2)*s12*(-5*TMath::Power(s12,2) + 8*s12*s13 - 3*TMath::Power(s13,2)) +
159 2*Q2*TMath::Power(s12,2)*(3*TMath::Power(s12,2) - 9*s12*s13 + 2*TMath::Power(s13,2))))))/(TMath::Power(kMw2,3)*TMath::Power(s12,2)*TMath::Power(ml2 - kMw2 + s13,2)*TMath::Power(ml2 - kMw2 - s12 + s13,2)*fsw2.Re());
160 return TMath::Max(0.,ME2);
161}
162//____________________________________________________________________________
163double Born::PXSecPhoton_L(double s12, double s13, double Q2, double ml2)
164{
165 double ME2 = 0.0;
166 ME2 = 2*falpha.Rho2()*Q2*TMath::Power(kMw2,-3)*kPi2*TMath::Power(s12,-2)*TMath::Power(ml2 - kMw2 + s13,-2)*TMath::Power(ml2 - kMw2 - s12 + s13,-2)*
167 ((s12 - s13)*TMath::Power(ml2,5)*TMath::Power(s12,2) - 2*s12*TMath::Power(ml2,4)*(2*kMw2*TMath::Power(s12,2) - (Q2 - s12)*(-3*s12*s13 + TMath::Power(s12,2) + 2*TMath::Power(s13,2))) +
168 TMath::Power(ml2,3)*(-2*(s12 - s13)*TMath::Power(kMw2,2)*TMath::Power(s12,2) +
169 2*kMw2*s12*(Q2*(9*s12*s13 - 5*TMath::Power(s12,2) - 2*TMath::Power(s13,2)) + s12*(-11*s12*s13 + 3*TMath::Power(s12,2) + 4*TMath::Power(s13,2))) +
170 (s12 - s13)*(TMath::Power(Q2,2)*TMath::Power(s12 - 2*s13,2) + TMath::Power(s12,2)*(-6*s12*s13 + TMath::Power(s12,2) + 6*TMath::Power(s13,2)) -
171 2*Q2*s12*(-5*s12*s13 + TMath::Power(s12,2) + 6*TMath::Power(s13,2)))) +
172 2*TMath::Power(ml2,2)*(4*(2*s12 - s13)*TMath::Power(kMw2,3)*TMath::Power(s12,2) +
173 (Q2 - s12)*s13*(Q2*(s12 - 2*s13) + s12*(-s12 + s13))*(-3*s12*s13 + TMath::Power(s12,2) + 2*TMath::Power(s13,2)) +
174 s12*TMath::Power(kMw2,2)*(Q2*(-11*s12*s13 + 7*TMath::Power(s12,2) - 2*TMath::Power(s13,2)) + s12*(15*s12*s13 - 11*TMath::Power(s12,2) + 2*TMath::Power(s13,2))) -
175 kMw2*((s12 - s13)*TMath::Power(Q2,2)*TMath::Power(s12 - 2*s13,2) + TMath::Power(s12,2)*(-11*s13*TMath::Power(s12,2) + TMath::Power(s12,3) + 20*s12*TMath::Power(s13,2) - 8*TMath::Power(s13,3)) -
176 2*Q2*s12*(-10*s13*TMath::Power(s12,2) + TMath::Power(s12,3) + 17*s12*TMath::Power(s13,2) - 6*TMath::Power(s13,3)))) +
177 4*TMath::Power(kMw2,2)*TMath::Power(s12,2)*(-(s13*(Q2 - 7*s12 + 4*s13)*TMath::Power(kMw2,2)) + (s12 - 2*s13)*TMath::Power(kMw2,3) + kMw2*(2*Q2 - s12 - 2*s13)*TMath::Power(s13,2) +
178 (-Q2 + s12)*TMath::Power(s13,3)) + ml2*(-15*(s12 - s13)*TMath::Power(kMw2,4)*TMath::Power(s12,2) +
179 2*s12*TMath::Power(kMw2,3)*(Q2*(7*s12*s13 - 3*TMath::Power(s12,2) + 2*TMath::Power(s13,2)) + s12*(-21*s12*s13 + 9*TMath::Power(s12,2) + 4*TMath::Power(s13,2))) +
180 TMath::Power(kMw2,2)*((s12 - s13)*TMath::Power(Q2,2)*TMath::Power(s12 - 2*s13,2) +
181 TMath::Power(s12,2)*(-31*s13*TMath::Power(s12,2) + TMath::Power(s12,3) + 60*s12*TMath::Power(s13,2) - 14*TMath::Power(s13,3)) -
182 2*Q2*s12*(-14*s13*TMath::Power(s12,2) + TMath::Power(s12,3) + 27*s12*TMath::Power(s13,2) - 6*TMath::Power(s13,3))) -
183 2*kMw2*s13*((s12 - s13)*TMath::Power(Q2,2)*TMath::Power(s12 - 2*s13,2) +
184 TMath::Power(s12,2)*(-8*s13*TMath::Power(s12,2) + TMath::Power(s12,3) + 11*s12*TMath::Power(s13,2) - 4*TMath::Power(s13,3)) +
185 Q2*s12*(15*s13*TMath::Power(s12,2) - 2*TMath::Power(s12,3) - 25*s12*TMath::Power(s13,2) + 10*TMath::Power(s13,3))) +
186 (s12 - s13)*TMath::Power(s13,2)*TMath::Power(Q2*(s12 - 2*s13) + s12*(-s12 + s13),2)))/fsw2.Re();
187 return TMath::Max(0.,ME2);
188}
189//____________________________________________________________________________
190double Born::Lambda(double a, double b, double c)
191{
192 return a*a + b*b + c*c - 2*a*b - 2*a*c - 2*b*c;
193}
194//____________________________________________________________________________
195double Born::GetS(double mlin, double Enuin)
196{
197 return 2. * mlin * Enuin + mlin*mlin;
198}
199//____________________________________________________________________________
200double Born::GetT(double mlin, double mlout, double s, double costhCM)
201{
202 //http://edu.itp.phys.ethz.ch/hs10/ppp1/PPP1_2.pdf [Sec 2.2.1]
203 double sum = mlin*mlin+mlout*mlout;
204 return ( (TMath::Sqrt(Lambda(s,0.,mlin*mlin)*Lambda(s,mlout*mlout,0.))*costhCM+mlin*mlin*mlout*mlout)/s + sum - s ) /2.;
205}
206//____________________________________________________________________________
207double Born::GetU(double mlin, double mlout, double s, double t)
208{
209 return mlin*mlin+mlout*mlout-s-t;
210}
211//____________________________________________________________________________
212bool Born::IsInPhaseSpace(double mlin, double mlout, double Enuin, double Enuout)
213{
214
215 //https://arxiv.org/pdf/2007.14426.pdf [section 2.2]
216 double frac = Enuout/Enuin;
217 if ( frac < mlin/(mlin+2.*Enuin)+(mlout*mlout-mlin*mlin)/2./Enuin/(mlin+2.*Enuin) ) return false;
218 else if ( frac > 1.-(mlout*mlout-mlin*mlin)/2./Enuin/mlin ) return false;
219
220 return true;
221
222}
223
224
225
226
PDGCodes.h
Most commonly used PDG codes. A set of utility functions to handle PDG codes is provided in PDGUtils.
genie::Born::GetT
double GetT(double mlin, double mlout, double s, double costhCM)
Definition Born.cxx:200
genie::Born::fgbe
TComplex fgbe
Definition Born.h:59
genie::Born::PXSecPhoton_T
double PXSecPhoton_T(double s12, double s13, double Q2, double ml2)
Definition Born.cxx:135
genie::Born::fGz
double fGz
Definition Born.h:51
genie::Born::fsw2
TComplex fsw2
Definition Born.h:54
genie::Born::fgav
TComplex fgav
Definition Born.h:60
genie::Born::PXSecNCVnu
double PXSecNCVnu(double s, double t, double mlin, double mlout)
Definition Born.cxx:98
genie::Born::PXSecCCV
double PXSecCCV(double s, double t, double mlin, double mlout)
Definition Born.cxx:67
genie::Born::fcw2
TComplex fcw2
Definition Born.h:55
genie::Born::PXSecPhoton
double PXSecPhoton(double s, double t, double mlout2)
Definition Born.cxx:120
genie::Born::~Born
virtual ~Born()
Definition Born.cxx:42
genie::Born::IsInPhaseSpace
bool IsInPhaseSpace(double mlin, double mlout, double Enuin, double Enuout)
Definition Born.cxx:212
genie::Born::fmz2c
TComplex fmz2c
Definition Born.h:57
genie::Born::fmw2c
TComplex fmw2c
Definition Born.h:56
genie::Born::fGw
double fGw
Definition Born.h:50
genie::Born::falpha
TComplex falpha
Definition Born.h:53
genie::Born::fgae
TComplex fgae
Definition Born.h:58
genie::Born::PXSecNCVnubar
double PXSecNCVnubar(double s, double t, double mlin, double mlout)
Definition Born.cxx:109
genie::Born::Lambda
double Lambda(double a, double b, double c)
Definition Born.cxx:190
genie::Born::GetU
double GetU(double mlin, double mlout, double s, double t)
Definition Born.cxx:207
genie::Born::PXSecCCVNC
double PXSecCCVNC(double s, double t, double mlin, double mlout)
Definition Born.cxx:87
genie::Born::PXSecCCRNC
double PXSecCCRNC(double s, double t, double mlin, double mlout)
Definition Born.cxx:76
genie::Born::PXSecCCR
double PXSecCCR(double s, double t, double mlin, double mlout)
Definition Born.cxx:58
genie::Born::GetS
double GetS(double mlin, double Enuin)
Definition Born.cxx:195
genie::Born::PXSecPhoton_L
double PXSecPhoton_L(double s12, double s13, double Q2, double ml2)
Definition Born.cxx:163
genie::PDGLibrary::Instance
static PDGLibrary * Instance(void)
Definition PDGLibrary.cxx:68
genie::PDGLibrary::Find
TParticlePDG * Find(int pdgc, bool must_exist=true)
Definition PDGLibrary.cxx:86
a
const double a
Definition gUpMuFluxGen.cxx:164
genie::constants
Basic constants.
genie
THE MAIN GENIE PROJECT NAMESPACE
Definition AlgCmp.h:25
genie::kPdgWM
const int kPdgWM
Definition PDGCodes.h:192
genie::kPdgZ0
const int kPdgZ0
Definition PDGCodes.h:190
Generated on for GENIEGenerator by doxygen 1.14.0