Compute the cross section for the input interaction.
86 double W = kinematics.
W();
87 double q2 = kinematics.
q2();
88 double costh = kinematics.
FSLeptonP4().CosTheta();
93 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 95 <<
"RES/DIS Join Scheme: XSec[RES, W=" << W
96 <<
" >= Wcut=" <<
fWcut <<
"] = 0";
109 int probepdgc = init_state.
ProbePdg();
118 bool is_EM = proc_info.
IsEM();
124 if(is_CC && !is_delta) {
125 if((is_nu && is_p) || (is_nubar && is_n))
return 0;
147 double W2 = TMath::Power(W, 2);
148 double Mnuc2 = TMath::Power(Mnuc, 2);
149 double k = 0.5 * (W2 - Mnuc2)/Mnuc;
150 double v = k - 0.5 * q2/Mnuc;
151 double v2 = TMath::Power(v, 2);
153 double Q = TMath::Sqrt(Q2);
154 double Eprime = E - v;
155 double U = 0.5 * (E + Eprime + Q) / E;
156 double V = 0.5 * (E + Eprime - Q) / E;
157 double U2 = TMath::Power(U, 2);
158 double V2 = TMath::Power(V, 2);
166 if(
fKLN && is_CC) is_KLN=
true;
169 if(
fBRS && is_CC) is_BRS=
true;
172 double Pl = TMath::Sqrt(Eprime*Eprime - ml*ml);
174 double vstar = (Mnuc*v + q2)/W;
175 double Qstar = TMath::Sqrt(-q2 + vstar*vstar);
176 double sqrtq2 = TMath::Sqrt(-q2);
177 double a = 1. + 0.5*(W2-q2+Mnuc2)/Mnuc/W;
179 double KNL_Alambda_plus = 0;
180 double KNL_Alambda_minus = 0;
181 double KNL_j0_plus = 0;
182 double KNL_j0_minus = 0;
183 double KNL_jx_plus = 0;
184 double KNL_jx_minus = 0;
185 double KNL_jy_plus = 0;
186 double KNL_jy_minus = 0;
187 double KNL_jz_plus = 0;
188 double KNL_jz_minus = 0;
189 double KNL_Qstar_plus =0;
190 double KNL_Qstar_minus =0;
192 double KNL_K = Q/E/TMath::Sqrt(2*(-q2));
194 double KNL_cL_plus = 0;
195 double KNL_cL_minus = 0;
197 double KNL_cR_plus = 0;
198 double KNL_cR_minus = 0;
200 double KNL_cS_plus = 0;
201 double KNL_cS_minus = 0;
203 double KNL_vstar_plus = 0;
204 double KNL_vstar_minus = 0;
206 if(is_CC && (is_KLN || is_BRS)){
208 LOG(
"BSKLNBaseRESPXSec2014",
pINFO)
"costh1="<<costh;
209 costh = (q2 - ml*ml + 2.*E*Eprime)/2./E/Pl;
211 LOG(
"BSKLNBaseRESPXSec2014",
pINFO)
"q2="<<q2<<
"m2="<<ml*ml<<
" 2.*E*Eprime="<<2.*E*Eprime<<
" nom="<< (q2 - ml*ml + 2.*E*Eprime)<<
" den="<<2.*E*Pl;
212 LOG(
"BSKLNBaseRESPXSec2014",
pINFO)
"costh2="<<costh;
214 KNL_Alambda_plus = TMath::Sqrt(E*(Eprime - Pl));
215 KNL_Alambda_minus = TMath::Sqrt(E*(Eprime + Pl));
217 <<
"\n+++++++++++++++++++++++ \n" 218 <<
"E="<<E <<
" K= "<<KNL_K <<
"\n" 219 <<
"El="<<Eprime<<
" Pl="<<Pl<<
" ml="<<ml <<
"\n" 220 <<
"W="<<W<<
" Q="<<Q<<
" q2="<<q2 <<
"\n" 221 <<
"A-="<<KNL_Alambda_minus<<
" A+="<<KNL_Alambda_plus <<
"\n" 222 <<
"xxxxxxxxxxxxxxxxxxxxxxx";
224 KNL_j0_plus = KNL_Alambda_plus /W * TMath::Sqrt(1 - costh) * (Mnuc - Eprime - Pl);
225 KNL_j0_minus = KNL_Alambda_minus/W * TMath::Sqrt(1 + costh) * (Mnuc - Eprime + Pl);
227 KNL_jx_plus = KNL_Alambda_plus/ Q * TMath::Sqrt(1 + costh) * (Pl -
E);
228 KNL_jx_minus = KNL_Alambda_minus/Q * TMath::Sqrt(1 - costh) * (Pl +
E);
230 KNL_jy_plus = KNL_Alambda_plus * TMath::Sqrt(1 + costh);
231 KNL_jy_minus = -KNL_Alambda_minus * TMath::Sqrt(1 - costh);
233 KNL_jz_plus = KNL_Alambda_plus /W/Q * TMath::Sqrt(1 - costh) * ( (E + Pl)*(Mnuc -Eprime) + Pl*( E + 2*E*costh -Pl) );
234 KNL_jz_minus = KNL_Alambda_minus/W/Q * TMath::Sqrt(1 + costh) * ( (E - Pl)*(Mnuc -Eprime) + Pl*( -E + 2*E*costh -Pl) );
236 if (is_nu || is_lminus) {
237 KNL_Qstar_plus = sqrtq2 * KNL_j0_plus / TMath::Sqrt(TMath::Abs(KNL_j0_plus*KNL_j0_plus - KNL_jz_plus*KNL_jz_plus) );
238 KNL_Qstar_minus = sqrtq2 * KNL_j0_minus / TMath::Sqrt(TMath::Abs(KNL_j0_minus*KNL_j0_minus - KNL_jz_minus*KNL_jz_minus) );
241 else if (is_nubar || is_lplus){
242 KNL_Qstar_plus = sqrtq2 * KNL_j0_minus / TMath::Sqrt(TMath::Abs(KNL_j0_minus*KNL_j0_minus - KNL_jz_minus*KNL_jz_minus) );
243 KNL_Qstar_minus = sqrtq2 * KNL_j0_plus / TMath::Sqrt(TMath::Abs(KNL_j0_plus*KNL_j0_plus - KNL_jz_plus*KNL_jz_plus) );
246 if (is_nu || is_lminus) {
247 KNL_vstar_plus = sqrtq2 * KNL_jz_plus / TMath::Sqrt(TMath::Abs(KNL_j0_plus*KNL_j0_plus - KNL_jz_plus*KNL_jz_plus) );
248 KNL_vstar_minus = sqrtq2 * KNL_jz_minus / TMath::Sqrt(TMath::Abs(KNL_j0_minus*KNL_j0_minus - KNL_jz_minus*KNL_jz_minus) );
250 else if (is_nubar || is_lplus) {
251 KNL_vstar_minus = sqrtq2 * KNL_jz_plus / TMath::Sqrt(TMath::Abs(KNL_j0_plus*KNL_j0_plus - KNL_jz_plus*KNL_jz_plus) );
252 KNL_vstar_plus = sqrtq2 * KNL_jz_minus / TMath::Sqrt(TMath::Abs(KNL_j0_minus*KNL_j0_minus - KNL_jz_minus*KNL_jz_minus) );
255 if(is_nu || is_lminus){
256 KNL_cL_plus = TMath::Sqrt(0.5)* KNL_K * (KNL_jx_plus - KNL_jy_plus);
257 KNL_cL_minus = TMath::Sqrt(0.5)* KNL_K * (KNL_jx_minus - KNL_jy_minus);
259 KNL_cR_plus = -TMath::Sqrt(0.5)* KNL_K * (KNL_jx_plus + KNL_jy_plus);
260 KNL_cR_minus = -TMath::Sqrt(0.5)* KNL_K * (KNL_jx_minus + KNL_jy_minus);
262 KNL_cS_plus = KNL_K * TMath::Sqrt(TMath::Abs(KNL_j0_plus *KNL_j0_plus - KNL_jz_plus *KNL_jz_plus ) );
263 KNL_cS_minus = KNL_K * TMath::Sqrt(TMath::Abs(KNL_j0_minus*KNL_j0_minus - KNL_jz_minus*KNL_jz_minus) );
266 if (is_nubar || is_lplus) {
267 KNL_cL_plus = -1 * TMath::Sqrt(0.5)* KNL_K * (KNL_jx_minus + KNL_jy_minus);
268 KNL_cL_minus = 1 * TMath::Sqrt(0.5)* KNL_K * (KNL_jx_plus + KNL_jy_plus);
270 KNL_cR_plus = 1 * TMath::Sqrt(0.5)* KNL_K * (KNL_jx_minus - KNL_jy_minus);
271 KNL_cR_minus = -1 * TMath::Sqrt(0.5)* KNL_K * (KNL_jx_plus - KNL_jy_plus);
273 KNL_cS_plus = -1 * KNL_K * TMath::Sqrt(TMath::Abs(KNL_j0_minus*KNL_j0_minus - KNL_jz_minus*KNL_jz_minus) );
274 KNL_cS_minus = 1 * KNL_K * TMath::Sqrt(TMath::Abs(KNL_j0_plus*KNL_j0_plus - KNL_jz_plus*KNL_jz_plus) );
278 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"j0-="<<KNL_j0_minus<<
" j0+="<<KNL_j0_plus;
279 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"jx-="<<KNL_jx_minus<<
" jx+="<<KNL_jx_plus;
280 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"jy-="<<KNL_jy_minus<<
" jy+="<<KNL_jy_plus;
281 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"jz-="<<KNL_jz_minus<<
" jz+="<<KNL_jz_plus;
283 LOG(
"BSKLNBaseRESPXSec2014",
pINFO)
"sqrt2="<<sqrtq2<<
" jz+=:"<<KNL_jz_plus<<
" j0+="<<KNL_j0_plus<<
" denom="<<TMath::Sqrt(TMath::Abs(KNL_j0_plus*KNL_j0_plus - KNL_jz_plus*KNL_jz_plus) );
285 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"vstar-="<<KNL_vstar_minus<<
" vstar+="<<KNL_vstar_plus;
286 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"Qstar-="<<KNL_Qstar_minus<<
" Qstar+="<<KNL_Qstar_plus;
288 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 290 <<
"Kinematical params V = " << V <<
", U = " << U;
295 double Go = TMath::Power(1 - 0.25 * q2/Mnuc2, 0.5-IR);
296 double GV = Go * TMath::Power( 1./(1-q2/
fMv2), 2);
297 double GA = Go * TMath::Power( 1./(1-q2/
fMa2), 2);
301 LOG(
"BSKLNBaseRESPXSec2014",
pDEBUG) <<
"Using new GV";
302 double CV0 = 1./(1-q2/
fMv2/4.);
303 double CV3 = 2.13 * CV0 * TMath::Power( 1-q2/
fMv2,-2);
304 double CV4 = -1.51 * CV0 * TMath::Power( 1-q2/
fMv2,-2);
305 double CV5 = 0.48 * CV0 * TMath::Power( 1-q2/
fMv2/0.766, -2);
307 double GV3 = 0.5 / TMath::Sqrt(3) * ( CV3 * (W + Mnuc)/Mnuc
308 + CV4 * (W2 + q2 -Mnuc2)/2./Mnuc2
309 + CV5 * (W2 - q2 -Mnuc2)/2./Mnuc2 );
311 double GV1 = - 0.5 / TMath::Sqrt(3) * ( CV3 * (Mnuc2 -q2 +Mnuc*
W)/W/Mnuc
312 + CV4 * (W2 +q2 - Mnuc2)/2./Mnuc2
313 + CV5 * (W2 -q2 - Mnuc2)/2./Mnuc2 );
315 GV = 0.5 * TMath::Power( 1 - q2/(Mnuc + W)/(Mnuc + W), 0.5-IR)
316 * TMath::Sqrt( 3 * GV3*GV3 + GV1*GV1);
320 LOG(
"BSKLNBaseRESPXSec2014",
pDEBUG) <<
"Using new GA";
323 double CA5 = CA5_0 * TMath::Power( 1./(1-q2/
fMa2), 2);
325 GA = 0.5 * TMath::Sqrt(3.) * TMath::Power( 1 - q2/(Mnuc + W)/(Mnuc + W), 0.5-IR) * (1- (W2 +q2 -Mnuc2)/8./Mnuc2) * CA5;
327 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"GA= " <<GA <<
" C5A= " <<CA5;
335 double d = TMath::Power(W+Mnuc,2.) - q2;
336 double sq2omg = TMath::Sqrt(2./
fOmega);
337 double nomg = IR *
fOmega;
338 double mq_w = Mnuc*Q/
W;
341 fFKR.
Tv = GV / (3.*W*sq2omg);
343 fFKR.
S = (-q2/
Q2) * (3*W*Mnuc + q2 - Mnuc2) * GV / (6*Mnuc2);
344 fFKR.
Ta = (2./3.) * (
fZeta/sq2omg) * mq_w * GA / d;
346 fFKR.
B =
fZeta/(3.*W*sq2omg) * (1 + (W2-Mnuc2+q2)/
d) * GA;
347 fFKR.
C =
fZeta/(6.*Q) * (W2 - Mnuc2 + nomg*(W2-Mnuc2+q2)/
d) * (GA/Mnuc);
356 double KNL_S_plus = 0;
357 double KNL_S_minus = 0;
358 double KNL_B_plus = 0;
359 double KNL_B_minus = 0;
360 double KNL_C_plus = 0;
361 double KNL_C_minus = 0;
364 KNL_S_plus = (KNL_vstar_plus*vstar - KNL_Qstar_plus *Qstar )* (Mnuc2 -q2 - 3*W*Mnuc ) * GV / (6*Mnuc2)/Q2;
365 KNL_S_minus = (KNL_vstar_minus*vstar - KNL_Qstar_minus*Qstar )* (Mnuc2 -q2 - 3*W*Mnuc ) * GV / (6*Mnuc2)/Q2;
367 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"KNL S= " <<KNL_S_plus<<
"\t"<<KNL_S_minus<<
"\t"<<
fFKR.
S;
369 KNL_B_plus =
fZeta/(3.*W*sq2omg)/Qstar * (KNL_Qstar_plus + KNL_vstar_plus *Qstar/a/Mnuc ) * GA;
370 KNL_B_minus =
fZeta/(3.*W*sq2omg)/Qstar * (KNL_Qstar_minus + KNL_vstar_minus*Qstar/a/Mnuc ) * GA;
371 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"KNL B= " <<KNL_B_plus<<
"\t"<<KNL_B_minus<<
"\t"<<
fFKR.
B;
373 KNL_C_plus = ( (KNL_Qstar_plus*Qstar - KNL_vstar_plus*vstar ) * ( 1./3. + vstar/a/Mnuc)
374 + KNL_vstar_plus*(2./3.*W +q2/a/Mnuc + nomg/3./a/Mnuc) )*
fZeta * (GA/2./W/Qstar);
376 KNL_C_minus = ( (KNL_Qstar_minus*Qstar - KNL_vstar_minus*vstar ) * ( 1./3. + vstar/a/Mnuc)
377 + KNL_vstar_minus*(2./3.*W +q2/a/Mnuc + nomg/3./a/Mnuc) )*
fZeta * (GA/2./W/Qstar);
379 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"KNL C= "<<KNL_C_plus<<
"\t"<<KNL_C_minus<<
"\t"<<
fFKR.
C;
381 double BRS_S_plus = 0;
382 double BRS_S_minus = 0;
383 double BRS_B_plus = 0;
384 double BRS_B_minus = 0;
385 double BRS_C_plus = 0;
386 double BRS_C_minus = 0;
391 KNL_S_plus = (KNL_vstar_plus*vstar - KNL_Qstar_plus *Qstar )* (Mnuc2 -q2 - 3*W*Mnuc ) * GV / (6*Mnuc2)/Q2;
392 KNL_S_minus = (KNL_vstar_minus*vstar - KNL_Qstar_minus*Qstar )* (Mnuc2 -q2 - 3*W*Mnuc ) * GV / (6*Mnuc2)/Q2;
395 KNL_B_plus =
fZeta/(3.*W*sq2omg)/Qstar * (KNL_Qstar_plus + KNL_vstar_plus *Qstar/a/Mnuc ) * GA;
396 KNL_B_minus =
fZeta/(3.*W*sq2omg)/Qstar * (KNL_Qstar_minus + KNL_vstar_minus*Qstar/a/Mnuc ) * GA;
399 KNL_C_plus = ( (KNL_Qstar_plus*Qstar - KNL_vstar_plus*vstar ) * ( 1./3. + vstar/a/Mnuc)
400 + KNL_vstar_plus*(2./3.*W +q2/a/Mnuc + nomg/3./a/Mnuc) )*
fZeta * (GA/2./W/Qstar);
402 KNL_C_minus = ( (KNL_Qstar_minus*Qstar - KNL_vstar_minus*vstar ) * ( 1./3. + vstar/a/Mnuc)
403 + KNL_vstar_minus*(2./3.*W +q2/a/Mnuc + nomg/3./a/Mnuc) )*
fZeta * (GA/2./W/Qstar);
405 BRS_S_plus = KNL_S_plus;
406 BRS_S_minus = KNL_S_minus;
407 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"BRS S= " <<KNL_S_plus<<
"\t"<<KNL_S_minus<<
"\t"<<
fFKR.
S;
409 BRS_B_plus = KNL_B_plus +
fZeta*GA/2./W/Qstar*( KNL_Qstar_plus*vstar - KNL_vstar_plus*Qstar)
410 *( 2./3 /sq2omg *(vstar + Qstar*Qstar/Mnuc/a))/(
kPionMass2 -q2);
412 BRS_B_minus = KNL_B_minus +
fZeta*GA/2./W/Qstar*( KNL_Qstar_minus*vstar - KNL_vstar_minus*Qstar)
413 *( 2./3 /sq2omg *(vstar + Qstar*Qstar/Mnuc/a))/(
kPionMass2 -q2);
414 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"BRS B= " <<KNL_B_plus<<
"\t"<<KNL_B_minus<<
"\t"<<
fFKR.
B;
416 BRS_C_plus = KNL_C_plus +
fZeta*GA/2./W/Qstar*( KNL_Qstar_plus*vstar - KNL_vstar_plus*Qstar)
417 * Qstar*(2./3.*W +q2/Mnuc/a +nomg/3./a/Mnuc)/(
kPionMass2 -q2);
419 BRS_C_minus = KNL_C_minus +
fZeta*GA/2./W/Qstar*( KNL_Qstar_minus*vstar - KNL_vstar_minus*Qstar)
420 * Qstar*(2./3.*W +q2/Mnuc/a +nomg/3./a/Mnuc)/(
kPionMass2 -q2);
421 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"BRS C= " <<KNL_C_plus<<
"\t"<<KNL_C_minus<<
"\t"<<
fFKR.
C;
424 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 426 <<
"FKR params for RES = " << resname <<
" : " <<
fFKR;
430 double sigL_minus = 0;
431 double sigR_minus = 0;
432 double sigS_minus = 0;
434 double sigL_plus = 0;
435 double sigR_plus = 0;
436 double sigS_plus = 0;
465 double sig0 = 0.125*(g2/
kPi)*(-q2/Q2)*(W/Mnuc);
466 double scLR = W/Mnuc;
467 double scS = (Mnuc/
W)*(-Q2/q2);
477 if(is_CC && !(is_KLN || is_BRS) ) {
492 if(is_CC && is_KLN ){
493 fFKR.S = KNL_S_minus;
494 fFKR.B = KNL_B_minus;
495 fFKR.C = KNL_C_minus;
499 assert(hamplmod_KNL_minus);
512 assert(hamplmod_KNL_plus);
522 if(is_CC && is_BRS ){
523 fFKR.S = BRS_S_minus;
524 fFKR.B = BRS_B_minus;
525 fFKR.C = BRS_C_minus;
528 assert(hamplmod_BRS_minus);
540 assert(hamplmod_BRS_plus);
550 if(is_KLN || is_BRS) {
560 <<
"sL,R,S minus = " << sigL_minus <<
"," << sigR_minus <<
"," << sigS_minus;
562 <<
"sL,R,S plus = " << sigL_plus <<
"," << sigR_plus <<
"," << sigS_plus;
574 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 575 LOG(
"BSKLNBaseRESPXSec2014",
pDEBUG) <<
"sig_{0} = " << sig0;
576 LOG(
"BSKLNBaseRESPXSec2014",
pDEBUG) <<
"sig_{L} = " << sigL;
577 LOG(
"BSKLNBaseRESPXSec2014",
pDEBUG) <<
"sig_{R} = " << sigR;
578 LOG(
"BSKLNBaseRESPXSec2014",
pDEBUG) <<
"sig_{S} = " << sigS;
583 if(is_KLN || is_BRS) {
584 xsec = TMath::Power(KNL_cL_minus,2)*sigL_minus + TMath::Power(KNL_cL_plus,2)*sigL_plus
585 + TMath::Power(KNL_cR_minus,2)*sigR_minus + TMath::Power(KNL_cR_plus,2)*sigR_plus
586 + TMath::Power(KNL_cS_minus,2)*sigS_minus + TMath::Power(KNL_cS_plus,2)*sigS_plus;
589 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"A-="<<KNL_Alambda_minus<<
" A+="<<KNL_Alambda_plus;
591 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<q2<<
"\t"<<xsec<<
"\t"<<sig0*(V2*sigR + U2*sigL + 2*UV*sigS)<<
"\t"<<xsec/TMath::Max(sig0*(V2*sigRSR + U2*sigRSL + 2*UV*sigRSS),1.0e-100);
592 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"fFKR.B="<<fFKR.B<<
" fFKR.C="<<fFKR.C<<
" fFKR.S="<<fFKR.S;
593 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"CL-="<<TMath::Power(KNL_cL_minus,2)<<
" CL+="<<TMath::Power(KNL_cL_plus,2)<<
" U2="<<U2;
594 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"SL-="<<sigL_minus<<
" SL+="<<sigL_plus<<
" SL="<<sigRSL;
596 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"CR-="<<TMath::Power(KNL_cR_minus,2)<<
" CR+="<<TMath::Power(KNL_cR_plus,2)<<
" V2="<<V2;
597 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"SR-="<<sigR_minus<<
" SR+="<<sigR_plus<<
" sR="<<sigRSR;
599 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"CS-="<<TMath::Power(KNL_cS_minus,2)<<
" CS+="<<TMath::Power(KNL_cS_plus,2)<<
" UV="<<UV;
600 LOG(
"BSKLNBaseRESPXSec2014",
pINFO) <<
"SS-="<<sigL_minus<<
" SS+="<<sigS_plus<<
" sS="<<sigRSS;
603 if (is_nu || is_lminus) {
604 xsec = sig0*(V2*sigR + U2*sigL + 2*UV*sigS);
607 if (is_nubar || is_lplus) {
608 xsec = sig0*(U2*sigR + V2*sigL + 2*UV*sigS);
610 xsec = TMath::Max(0.,xsec);
613 if ( is_CC && is_delta ) {
614 if ( (is_nu && is_p) || (is_nubar && is_n) ) mult=3.0;
624 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 626 <<
"BreitWigner(RES=" << resname <<
", W=" << W <<
") = " << bw;
630 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__ 632 <<
"\n d2xsec/dQ2dW" <<
"[" <<
interaction->AsString()
633 <<
"](W=" << W <<
", q2=" << q2 <<
", E=" << E <<
") = " << xsec;
655 int NNucl = (is_p) ? Z : N;
674 double P_Fermi = 0.0;
688 if(is_p) { P_Fermi *= TMath::Power( 2.*Z/A, 1./3); }
689 else { P_Fermi *= TMath::Power( 2.*N/A, 1./3); }
692 double FactorPauli_RES = 1.0;
694 double k0 = 0., q = 0., q0 = 0.;
698 k0 = (W2-Mnuc2-
Q2)/(2*W);
699 k = TMath::Sqrt(k0*k0+Q2);
704 if ( 2*P_Fermi < k-q )
705 FactorPauli_RES = 1.0;
706 if ( 2*P_Fermi >= k+q )
707 FactorPauli_RES = ((3*k*k+q*q)/(2*P_Fermi)-(5*TMath::Power(k,4)+TMath::Power(q,4)+10*k*k*q*q)/(40*TMath::Power(P_Fermi,3)))/(2*
k);
708 if ( 2*P_Fermi >= k-q && 2*P_Fermi <= k+q )
709 FactorPauli_RES = ((q+
k)*(q+k)-4*P_Fermi*P_Fermi/5-TMath::Power(k-q, 3)/(2*P_Fermi)+TMath::Power(k-q, 5)/(40*TMath::Power(P_Fermi, 3)))/(4*q*
k);
711 xsec *= FactorPauli_RES;
bool IsDelta(Resonance_t res)
is it a Delta resonance?
bool fNormBW
normalize resonance breit-wigner to 1?
virtual const RSHelicityAmpl & Compute(Resonance_t res, const FKR &fkr) const =0
string fKFTable
table of Fermi momentum (kF) constants for various nuclei
double fOmega
FKR parameter Omega.
double W(bool selected=false) const
bool IsWeakCC(void) const
static const double kSqrt2
bool IsNeutrino(int pdgc)
double fMv2
(vector mass)^2
bool fUsingDisResJoin
use a DIS/RES joining scheme?
double fXSecScaleNC
external NC xsec scaling factor
double J(double q0, double q3, double Enu, double ml)
double Q2(const Interaction *const i)
int HitNucPdg(void) const
double Amp2Plus3(void) const
double Amp2Minus3(void) const
double HitNucMass(void) const
double fN0ResMaxNWidths
limits allowed phase space for n=0 res
static FermiMomentumTablePool * Instance(void)
Generated/set kinematical variables for an event.
double Mass(Resonance_t res)
resonance mass (GeV)
A table of Fermi momentum constants.
double Width(Resonance_t res)
resonance width (GeV)
double Amp2Plus1(void) const
double Amp2Minus1(void) const
return |helicity amplitude|^2
double BreitWignerL(double W, int L, double mass, double width0, double norm)
double BWNorm(Resonance_t res, double N0ResMaxNWidths=6, double N2ResMaxNWidths=2, double GnResMaxNWidths=4)
breit-wigner normalization factor
enum genie::EResonance Resonance_t
const RSHelicityAmplModelI * fHAmplModelEMp
const RSHelicityAmplModelI * fHAmplModelCC
double fVud2
|Vud|^2(square of magnitude ud-element of CKM-matrix)
bool IsPosChargedLepton(int pdgc)
virtual bool ValidKinematics(const Interaction *i) const
Is the input kinematical point a physically allowed one?
double q2(bool selected=false) const
A class holding the Rein-Sehgal's helicity amplitudes.
bool IsWeakNC(void) const
const TLorentzVector & FSLeptonP4(void) const
Singleton class to load & serve tables of Fermi momentum constants.
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
bool fWghtBW
weight with resonance breit-wigner?
const FermiMomentumTable * GetTable(string name)
static const double kAem2
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
bool IsAntiNeutrino(int pdgc)
double fXSecScaleCC
external CC xsec scaling factor
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
bool fUsePauliBlocking
account for Pauli blocking?
double fWcut
apply DIS/RES joining scheme < Wcut
int OrbitalAngularMom(Resonance_t res)
orbital angular momentum
const RSHelicityAmplModelI * fHAmplModelEMn
Pure abstract base class. Defines the RSHelicityAmplModelI interface.
double fMa2
(axial mass)^2
double Amp20Minus(void) const
double fGnResMaxNWidths
limits allowed phase space for other res
double FermiMomentumForIsoscalarNucleonParametrization(const Target &target)
bool fUseRFGParametrization
use parametrization for fermi momentum insted of table?
const RSHelicityAmplModelI * fHAmplModelNCp
const UInt_t kIAssumeFreeNucleon
int IonPdgCode(int A, int Z)
double fSin48w
sin^4(Weingberg angle)
double Jacobian(const Interaction *const i, KinePhaseSpace_t f, KinePhaseSpace_t t)
double Amp20Plus(void) const
const char * AsString(Resonance_t res)
resonance id -> string
bool ValidProcess(const Interaction *i) const
Can this cross section algorithm handle the input process?
double FindClosestKF(int target_pdgc, int nucleon_pdgc) const
const RSHelicityAmplModelI * fHAmplModelNCn
const Target & Tgt(void) const
double fN2ResMaxNWidths
limits allowed phase space for n=2 res
double fZeta
FKR parameter Zeta.
double ProbeE(RefFrame_t rf) const
bool IsNegChargedLepton(int pdgc)
int ResonanceIndex(Resonance_t res)
resonance idx, quark model / SU(6)
Initial State information.
static const double kPionMass2