GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ReinSehgalRESPXSec.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*
3  Copyright (c) 2003-2024, The GENIE Collaboration
4  For the full text of the license visit http://copyright.genie-mc.org
5 
6  Costas Andreopoulos <c.andreopoulos \at cern.ch>
7  University of Liverpool
8 */
9 //____________________________________________________________________________
10 
11 #include <TMath.h>
12 #include <TSystem.h>
13 
17 #include "Framework/Conventions/GBuild.h"
28 #include "Framework/Utils/Range1.h"
29 #include "Framework/Utils/BWFunc.h"
37 
38 using namespace genie;
39 using namespace genie::constants;
40 
41 //____________________________________________________________________________
43 XSecAlgorithmI("genie::ReinSehgalRESPXSec")
44 {
45  fNuTauRdSpl = 0;
46  fNuTauBarRdSpl = 0;
47 }
48 //____________________________________________________________________________
50 XSecAlgorithmI("genie::ReinSehgalRESPXSec", config)
51 {
52  fNuTauRdSpl = 0;
53  fNuTauBarRdSpl = 0;
54 }
55 //____________________________________________________________________________
57 {
58  if(fNuTauRdSpl) delete fNuTauRdSpl;
59  if(fNuTauBarRdSpl) delete fNuTauBarRdSpl;
60 }
61 //____________________________________________________________________________
63  const Interaction * interaction, KinePhaseSpace_t kps) const
64 {
65  if(! this -> ValidProcess (interaction) ) return 0.;
66  if(! this -> ValidKinematics (interaction) ) return 0.;
67 
68  const InitialState & init_state = interaction -> InitState();
69  const ProcessInfo & proc_info = interaction -> ProcInfo();
70  const Target & target = init_state.Tgt();
71 
72  // Get kinematical parameters
73  const Kinematics & kinematics = interaction -> Kine();
74  double W = kinematics.W();
75  double q2 = kinematics.q2();
76 
77  // Under the DIS/RES joining scheme, xsec(RES)=0 for W>=Wcut
78  if(fUsingDisResJoin) {
79  if(W>=fWcut) {
80 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
81  LOG("ReinSehgalRes", pDEBUG)
82  << "RES/DIS Join Scheme: XSec[RES, W=" << W
83  << " >= Wcut=" << fWcut << "] = 0";
84 #endif
85  return 0;
86  }
87  }
88 
89  // Get the input baryon resonance
90  Resonance_t resonance = interaction->ExclTag().Resonance();
91  string resname = utils::res::AsString(resonance);
92  bool is_delta = utils::res::IsDelta (resonance);
93 
94  // Get the neutrino, hit nucleon & weak current
95  int nucpdgc = target.HitNucPdg();
96  int probepdgc = init_state.ProbePdg();
97  bool is_nu = pdg::IsNeutrino (probepdgc);
98  bool is_nubar = pdg::IsAntiNeutrino (probepdgc);
99  bool is_lplus = pdg::IsPosChargedLepton (probepdgc);
100  bool is_lminus = pdg::IsNegChargedLepton (probepdgc);
101  bool is_p = pdg::IsProton (nucpdgc);
102  bool is_n = pdg::IsNeutron (nucpdgc);
103  bool is_CC = proc_info.IsWeakCC();
104  bool is_NC = proc_info.IsWeakNC();
105  bool is_EM = proc_info.IsEM();
106 
107  if(is_CC && !is_delta) {
108  if((is_nu && is_p) || (is_nubar && is_n)) return 0;
109  }
110 
111  // Get baryon resonance parameters
112  int IR = utils::res::ResonanceIndex (resonance);
113  int LR = utils::res::OrbitalAngularMom (resonance);
114  double MR = utils::res::Mass (resonance);
115  double WR = utils::res::Width (resonance);
117 
118  // Following NeuGEN, avoid problems with underlying unphysical
119  // model assumptions by restricting the allowed W phase space
120  // around the resonance peak
121  if (fNormBW) {
122  if (W > MR + fN0ResMaxNWidths * WR && IR==0) return 0.;
123  else if (W > MR + fN2ResMaxNWidths * WR && IR==2) return 0.;
124  else if (W > MR + fGnResMaxNWidths * WR) return 0.;
125  }
126 
127  // Compute auxiliary & kinematical factors
128  double E = init_state.ProbeE(kRfHitNucRest);
129  double Mnuc = target.HitNucMass();
130  double W2 = TMath::Power(W, 2);
131  double Mnuc2 = TMath::Power(Mnuc, 2);
132  double k = 0.5 * (W2 - Mnuc2)/Mnuc;
133  double v = k - 0.5 * q2/Mnuc;
134  double v2 = TMath::Power(v, 2);
135  double Q2 = v2 - q2;
136  double Q = TMath::Sqrt(Q2);
137  double Eprime = E - v;
138  double U = 0.5 * (E + Eprime + Q) / E;
139  double V = 0.5 * (E + Eprime - Q) / E;
140  double U2 = TMath::Power(U, 2);
141  double V2 = TMath::Power(V, 2);
142  double UV = U*V;
143 
144 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
145  LOG("ReinSehgalRes", pDEBUG)
146  << "Kinematical params V = " << V << ", U = " << U;
147 #endif
148 
149  // Calculate the Feynman-Kislinger-Ravndall parameters
150 
151  double Go = TMath::Power(1 - 0.25 * q2/Mnuc2, 0.5-IR);
152  double GV = Go * TMath::Power( 1./(1-q2/fMv2), 2);
153  double GA = Go * TMath::Power( 1./(1-q2/fMa2), 2);
154 
155  if(is_EM) {
156  GA = 0.; // zero the axial term for EM scattering
157  }
158 
159  double d = TMath::Power(W+Mnuc,2.) - q2;
160  double sq2omg = TMath::Sqrt(2./fOmega);
161  double nomg = IR * fOmega;
162  double mq_w = Mnuc*Q/W;
163 
164  fFKR.Lamda = sq2omg * mq_w;
165  fFKR.Tv = GV / (3.*W*sq2omg);
166  fFKR.Rv = kSqrt2 * mq_w*(W+Mnuc)*GV / d;
167  fFKR.S = (-q2/Q2) * (3*W*Mnuc + q2 - Mnuc2) * GV / (6*Mnuc2);
168  fFKR.Ta = (2./3.) * (fZeta/sq2omg) * mq_w * GA / d;
169  fFKR.Ra = (kSqrt2/6.) * fZeta * (GA/W) * (W+Mnuc + 2*nomg*W/d );
170  fFKR.B = fZeta/(3.*W*sq2omg) * (1 + (W2-Mnuc2+q2)/ d) * GA;
171  fFKR.C = fZeta/(6.*Q) * (W2 - Mnuc2 + nomg*(W2-Mnuc2+q2)/d) * (GA/Mnuc);
172  fFKR.R = fFKR.Rv;
173  fFKR.Rplus = - (fFKR.Rv + fFKR.Ra);
174  fFKR.Rminus = - (fFKR.Rv - fFKR.Ra);
175  fFKR.T = fFKR.Tv;
176  fFKR.Tplus = - (fFKR.Tv + fFKR.Ta);
177  fFKR.Tminus = - (fFKR.Tv - fFKR.Ta);
178 
179 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
180  LOG("FKR", pDEBUG)
181  << "FKR params for RES = " << resname << " : " << fFKR;
182 #endif
183 
184  // Calculate the Rein-Sehgal Helicity Amplitudes
185 
186  const RSHelicityAmplModelI * hamplmod = 0;
187  if(is_CC) {
188  hamplmod = fHAmplModelCC;
189  }
190  else
191  if(is_NC) {
192  if (is_p) { hamplmod = fHAmplModelNCp;}
193  else { hamplmod = fHAmplModelNCn;}
194  }
195  else
196  if(is_EM) {
197  if (is_p) { hamplmod = fHAmplModelEMp;}
198  else { hamplmod = fHAmplModelEMn;}
199  }
200  assert(hamplmod);
201 
202  const RSHelicityAmpl & hampl = hamplmod->Compute(resonance, fFKR);
203 
204 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
205  LOG("RSHAmpl", pDEBUG)
206  << "Helicity Amplitudes for RES = " << resname << " : " << hampl;
207 #endif
208 
209  double g2 = kGF2;
210  if(is_CC) g2 = kGF2*fVud2;
211  // For EM interaction replace G_{Fermi} with :
212  // a_{em} * pi / ( sqrt(2) * sin^2(theta_weinberg) * Mass_{W}^2 }
213  // See C.Quigg, Gauge Theories of the Strong, Weak and E/M Interactions,
214  // ISBN 0-8053-6021-2, p.112 (6.3.57)
215  // Also, take int account that the photon propagator is 1/p^2 but the
216  // W propagator is 1/(p^2-Mass_{W}^2), so weight the EM case with
217  // Mass_{W}^4 / q^4
218  // So, overall:
219  // G_{Fermi}^2 --> a_{em}^2 * pi^2 / (2 * sin^4(theta_weinberg) * q^{4})
220  //
221  if(is_EM) {
222  double q4 = q2*q2;
223  g2 = kAem2 * kPi2 / (2.0 * fSin48w * q4);
224  }
225 
226  // Compute the cross section
227 
228  double sig0 = 0.125*(g2/kPi)*(-q2/Q2)*(W/Mnuc);
229  double scLR = W/Mnuc;
230  double scS = (Mnuc/W)*(-Q2/q2);
231  double sigL = scLR* (hampl.Amp2Plus3 () + hampl.Amp2Plus1 ());
232  double sigR = scLR* (hampl.Amp2Minus3() + hampl.Amp2Minus1());
233  double sigS = scS * (hampl.Amp20Plus () + hampl.Amp20Minus());
234 
235 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
236  LOG("ReinSehgalRes", pDEBUG) << "sig_{0} = " << sig0;
237  LOG("ReinSehgalRes", pDEBUG) << "sig_{L} = " << sigL;
238  LOG("ReinSehgalRes", pDEBUG) << "sig_{R} = " << sigR;
239  LOG("ReinSehgalRes", pDEBUG) << "sig_{S} = " << sigS;
240 #endif
241 
242  double xsec = 0.0;
243  if (is_nu || is_lminus) {
244  xsec = sig0*(V2*sigR + U2*sigL + 2*UV*sigS);
245  }
246  else
247  if (is_nubar || is_lplus) {
248  xsec = sig0*(U2*sigR + V2*sigL + 2*UV*sigS);
249  }
250  xsec = TMath::Max(0.,xsec);
251 
252  double mult = 1.0;
253  if(is_CC && is_delta) {
254  if((is_nu && is_p) || (is_nubar && is_n)) mult=3.0;
255  }
256  xsec *= mult;
257 
258  // Check whether the cross section is to be weighted with a
259  // Breit-Wigner distribution (default: true)
260  double bw = 1.0;
261  if(fWghtBW) {
262  //different Delta photon decay branch
263  if(is_delta){
264  bw = utils::bwfunc::BreitWignerLGamma(W,LR,MR,WR,NR);
265  }
266  else{
267  bw = utils::bwfunc::BreitWignerL(W,LR,MR,WR,NR);
268  }
269  }
270 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
271  LOG("ReinSehgalRes", pDEBUG)
272  << "BreitWigner(RES=" << resname << ", W=" << W << ") = " << bw;
273 #endif
274  xsec *= bw;
275 
276  // Apply NeuGEN nutau cross section reduction factors
277  double rf = 1.0;
278  Spline * spl = 0;
279  if (is_CC && fUsingNuTauScaling) {
280  if (pdg::IsNuTau (probepdgc)) spl = fNuTauRdSpl;
281  else if (pdg::IsAntiNuTau(probepdgc)) spl = fNuTauBarRdSpl;
282 
283  if(spl) {
284  if(E <spl->XMax()) rf = spl->Evaluate(E);
285  }
286  }
287  xsec *= rf;
288 
289  // Apply given scaling factor
290  double xsec_scale = 1.;
291  if (is_CC) { xsec_scale = fXSecScaleCC; }
292  else if (is_NC) { xsec_scale = fXSecScaleNC; }
293  else if (is_EM) { xsec_scale = fXSecScaleEM; }
294  xsec *= xsec_scale;
295 
296 #ifdef __GENIE_LOW_LEVEL_MESG_ENABLED__
297  LOG("ReinSehgalRes", pINFO)
298  << "\n d2xsec/dQ2dW" << "[" << interaction->AsString()
299  << "](W=" << W << ", q2=" << q2 << ", E=" << E << ") = " << xsec;
300 #endif
301 
302  // The algorithm computes d^2xsec/dWdQ2
303  // Check whether variable tranformation is needed
304  if(kps!=kPSWQ2fE) {
305  double J = utils::kinematics::Jacobian(interaction,kPSWQ2fE,kps);
306  xsec *= J;
307  }
308 
309  // If requested return the free nucleon xsec even for input nuclear tgt
310  if( interaction->TestBit(kIAssumeFreeNucleon) ) return xsec;
311 
312 
313  int Z = target.Z();
314  int A = target.A();
315  int N = A-Z;
316 
317  // Take into account the number of scattering centers in the target
318  int NNucl = (is_p) ? Z : N;
319 
320  xsec*=NNucl; // nuclear xsec (no nuclear suppression factor)
321 
322  if (fUsePauliBlocking && A!=1)
323  {
324  // Calculation of Pauli blocking according references:
325  //
326  // [1] S.L. Adler, S. Nussinov, and E.A. Paschos, "Nuclear
327  // charge exchange corrections to leptonic pion production
328  // in the (3,3) resonance region," Phys. Rev. D 9 (1974)
329  // 2125-2143 [Erratum Phys. Rev. D 10 (1974) 1669].
330  // [2] J.Y. Yu, "Neutrino interactions and nuclear effects in
331  // oscillation experiments and the nonperturbative disper-
332  // sive sector in strong (quasi-)abelian fields," Ph. D.
333  // Thesis, Dortmund U., Dortmund, 2002 (unpublished).
334  // [3] E.A. Paschos, J.Y. Yu, and M. Sakuda, "Neutrino pro-
335  // duction of resonances," Phys. Rev. D 69 (2004) 014013
336  // [arXiv: hep-ph/0308130].
337 
338  double P_Fermi = 0.0;
339 
340  // Maximum value of Fermi momentum of target nucleon (GeV)
341  if (A<6 || !fUseRFGParametrization)
342  {
343  // Look up the Fermi momentum for this target
345  const FermiMomentumTable * kft = kftp->GetTable(fKFTable);
346  P_Fermi = kft->FindClosestKF(pdg::IonPdgCode(A, Z), nucpdgc);
347  }
348  else {
349  // Define the Fermi momentum for this target
351  // Correct the Fermi momentum for the struck nucleon
352  if(is_p) { P_Fermi *= TMath::Power( 2.*Z/A, 1./3); }
353  else { P_Fermi *= TMath::Power( 2.*N/A, 1./3); }
354  }
355 
356  double FactorPauli_RES = 1.0;
357 
358  double k0 = 0., q = 0., q0 = 0.;
359 
360  if (P_Fermi > 0.)
361  {
362  k0 = (W2-Mnuc2-Q2)/(2*W);
363  k = TMath::Sqrt(k0*k0+Q2); // previous value of k is overridden
364  q0 = (W2-Mnuc2+kPionMass2)/(2*W);
365  q = TMath::Sqrt(q0*q0-kPionMass2);
366  }
367 
368  if (2*P_Fermi < k-q)
369  FactorPauli_RES = 1.0;
370  if (2*P_Fermi >= k+q)
371  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);
372  if (2*P_Fermi >= k-q && 2*P_Fermi <= k+q)
373  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);
374 
375  xsec *= FactorPauli_RES;
376  }
377 
378  return xsec;
379 }
380 //____________________________________________________________________________
381 double ReinSehgalRESPXSec::Integral(const Interaction * interaction) const
382 {
383  double xsec = fXSecIntegrator->Integrate(this,interaction);
384  return xsec;
385 }
386 //____________________________________________________________________________
387 bool ReinSehgalRESPXSec::ValidProcess(const Interaction * interaction) const
388 {
389  if(interaction->TestBit(kISkipProcessChk)) return true;
390 
391  const InitialState & init_state = interaction->InitState();
392  const ProcessInfo & proc_info = interaction->ProcInfo();
393  const XclsTag & xcls = interaction->ExclTag();
394 
395  if(!proc_info.IsResonant()) return false;
396  if(!xcls.KnownResonance()) return false;
397 
398  int hitnuc = init_state.Tgt().HitNucPdg();
399  bool is_pn = (pdg::IsProton(hitnuc) || pdg::IsNeutron(hitnuc));
400 
401  if (!is_pn) return false;
402 
403  int probe = init_state.ProbePdg();
404  bool is_weak = proc_info.IsWeak();
405  bool is_em = proc_info.IsEM();
406  bool nu_weak = (pdg::IsNeutralLepton(probe) && is_weak);
407  bool l_em = (pdg::IsChargedLepton(probe) && is_em );
408 
409  if (!nu_weak && !l_em) return false;
410 
411  return true;
412 }
413 //____________________________________________________________________________
415 {
416  Algorithm::Configure(config);
417  this->LoadConfig();
418 }
419 //____________________________________________________________________________
420 void ReinSehgalRESPXSec::Configure(string config)
421 {
422  Algorithm::Configure(config);
423  this->LoadConfig();
424 }
425 //____________________________________________________________________________
427 {
428  // Cross section scaling factors
429  this->GetParam( "RES-CC-XSecScale", fXSecScaleCC ) ;
430  this->GetParam( "RES-NC-XSecScale", fXSecScaleNC ) ;
431  this->GetParam( "RES-EM-XSecScale", fXSecScaleEM ) ;
432 
433  this->GetParam( "RES-Zeta", fZeta ) ;
434  this->GetParam( "RES-Omega", fOmega ) ;
435 
436  double ma, mv ;
437  this->GetParam( "RES-Ma", ma ) ;
438  this->GetParam( "RES-Mv", mv ) ;
439  fMa2 = TMath::Power(ma,2);
440  fMv2 = TMath::Power(mv,2);
441 
442  this->GetParamDef( "BreitWignerWeight", fWghtBW, true ) ;
443  this->GetParamDef( "BreitWignerNorm", fNormBW, true);
444 
445  double thw ;
446  this->GetParam( "WeinbergAngle", thw ) ;
447  fSin48w = TMath::Power( TMath::Sin(thw), 4 );
448  double Vud;
449  this->GetParam("CKM-Vud", Vud );
450  fVud2 = TMath::Power( Vud, 2 );
451  this->GetParam("FermiMomentumTable", fKFTable);
452  this->GetParam("RFG-UseParametrization", fUseRFGParametrization);
453  this->GetParam("UsePauliBlockingForRES", fUsePauliBlocking);
454 
455  // Load all the sub-algorithms needed
456 
457  fHAmplModelCC = 0;
458  fHAmplModelNCp = 0;
459  fHAmplModelNCn = 0;
460  fHAmplModelEMp = 0;
461  fHAmplModelEMn = 0;
462 
463  AlgFactory * algf = AlgFactory::Instance();
464 
465  fHAmplModelCC = dynamic_cast<const RSHelicityAmplModelI *> (
466  algf->GetAlgorithm("genie::RSHelicityAmplModelCC","Default"));
467  fHAmplModelNCp = dynamic_cast<const RSHelicityAmplModelI *> (
468  algf->GetAlgorithm("genie::RSHelicityAmplModelNCp","Default"));
469  fHAmplModelNCn = dynamic_cast<const RSHelicityAmplModelI *> (
470  algf->GetAlgorithm("genie::RSHelicityAmplModelNCn","Default"));
471  fHAmplModelEMp = dynamic_cast<const RSHelicityAmplModelI *> (
472  algf->GetAlgorithm("genie::RSHelicityAmplModelEMp","Default"));
473  fHAmplModelEMn = dynamic_cast<const RSHelicityAmplModelI *> (
474  algf->GetAlgorithm("genie::RSHelicityAmplModelEMn","Default"));
475 
476  assert( fHAmplModelCC );
477  assert( fHAmplModelNCp );
478  assert( fHAmplModelNCn );
479  assert( fHAmplModelEMp );
480  assert( fHAmplModelEMn );
481 
482  // Use algorithm within a DIS/RES join scheme. If yes get Wcut
483  this->GetParam( "UseDRJoinScheme", fUsingDisResJoin ) ;
484  fWcut = 999999;
485  if(fUsingDisResJoin) {
486  this->GetParam( "Wcut", fWcut ) ;
487  }
488 
489  // NeuGEN limits in the allowed resonance phase space:
490  // W < min{ Wmin(physical), (res mass) + x * (res width) }
491  // It limits the integration area around the peak and avoids the
492  // problem with huge xsec increase at low Q2 and high W.
493  // In correspondence with Hugh, Rein said that the underlying problem
494  // are unphysical assumptions in the model.
495  this->GetParamDef( "MaxNWidthForN2Res", fN2ResMaxNWidths, 2.0 ) ;
496  this->GetParamDef( "MaxNWidthForN0Res", fN0ResMaxNWidths, 6.0 ) ;
497  this->GetParamDef( "MaxNWidthForGNRes", fGnResMaxNWidths, 4.0 ) ;
498 
499  // NeuGEN reduction factors for nu_tau: a gross estimate of the effect of
500  // neglected form factors in the R/S model
501  this->GetParamDef( "UseNuTauScalingFactors", fUsingNuTauScaling, true ) ;
502  if(fUsingNuTauScaling) {
503  if(fNuTauRdSpl) delete fNuTauRdSpl;
504  if(fNuTauBarRdSpl) delete fNuTauBarRdSpl;
505 
506  assert( std::getenv( "GENIE") );
507  string base = std::getenv( "GENIE") ;
508 
509  string filename = base + "/data/evgen/rein_sehgal/res/nutau_xsec_scaling_factors.dat";
510  LOG("ReinSehgalRes", pNOTICE)
511  << "Loading nu_tau xsec reduction spline from: " << filename;
512  fNuTauRdSpl = new Spline(filename);
513 
514  filename = base + "/data/evgen/rein_sehgal/res/nutaubar_xsec_scaling_factors.dat";
515  LOG("ReinSehgalRes", pNOTICE)
516  << "Loading bar{nu_tau} xsec reduction spline from: " << filename;
517  fNuTauBarRdSpl = new Spline(filename);
518  }
519 
520  // load the differential cross section integrator
522  dynamic_cast<const XSecIntegratorI *> (this->SubAlg("XSec-Integrator"));
523  assert(fXSecIntegrator);
524 }
525 //____________________________________________________________________________
bool IsDelta(Resonance_t res)
is it a Delta resonance?
bool IsResonant(void) const
Definition: ProcessInfo.cxx:99
Cross Section Calculation Interface.
double W(bool selected=false) const
Definition: Kinematics.cxx:157
bool ValidProcess(const Interaction *i) const
Can this cross section algorithm handle the input process?
bool IsNuTau(int pdgc)
Definition: PDGUtils.cxx:168
bool fNormBW
normalize resonance breit-wigner to 1?
bool IsWeak(void) const
bool IsWeakCC(void) const
bool IsNeutrino(int pdgc)
Definition: PDGUtils.cxx:110
double J(double q0, double q3, double Enu, double ml)
Definition: MECUtils.cxx:147
double Rminus
Definition: FKR.h:50
Spline * fNuTauBarRdSpl
xsec reduction spline for nu_tau_bar
Cross Section Integrator Interface.
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1077
bool fWghtBW
weight with resonance breit-wigner?
int HitNucPdg(void) const
Definition: Target.cxx:304
const RSHelicityAmplModelI * fHAmplModelEMn
const RSHelicityAmplModelI * fHAmplModelNCp
double Ra
Definition: FKR.h:42
double fVud2
|Vud|^2(square of magnitude ud-element of CKM-matrix)
int A(void) const
Definition: Target.h:70
virtual const RSHelicityAmpl & Compute(Resonance_t res, const FKR &fkr) const =0
bool KnownResonance(void) const
Definition: XclsTag.h:68
double HitNucMass(void) const
Definition: Target.cxx:233
A numeric analysis tool class for interpolating 1-D functions.
Definition: Spline.h:58
static FermiMomentumTablePool * Instance(void)
Generated/set kinematical variables for an event.
Definition: Kinematics.h:39
bool IsAntiNuTau(int pdgc)
Definition: PDGUtils.cxx:183
double Lamda
Definition: FKR.h:37
bool IsChargedLepton(int pdgc)
Definition: PDGUtils.cxx:101
double Mass(Resonance_t res)
resonance mass (GeV)
double R
Definition: FKR.h:45
A table of Fermi momentum constants.
double XSec(const Interaction *i, KinePhaseSpace_t k) const
Compute the cross section for the input interaction.
double Width(Resonance_t res)
resonance width (GeV)
double fMv2
(vector mass)^2
double BreitWignerLGamma(double W, int L, double mass, double width0, double norm)
Definition: BWFunc.cxx:22
double Evaluate(double x) const
Definition: Spline.cxx:363
double BreitWignerL(double W, int L, double mass, double width0, double norm)
Definition: BWFunc.cxx:99
enum genie::EKinePhaseSpace KinePhaseSpace_t
const RSHelicityAmplModelI * fHAmplModelCC
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 * fHAmplModelNCn
string AsString(void) const
Contains minimal information for tagging exclusive processes.
Definition: XclsTag.h:39
double Integral(const Interaction *i) const
const RSHelicityAmplModelI * fHAmplModelEMp
double W(const Interaction *const i)
Definition: KineUtils.cxx:1101
bool IsNeutron(int pdgc)
Definition: PDGUtils.cxx:341
double fZeta
FKR parameter Zeta.
bool IsPosChargedLepton(int pdgc)
Definition: PDGUtils.cxx:148
Summary information for an interaction.
Definition: Interaction.h:56
double Tv
Definition: FKR.h:38
virtual bool ValidKinematics(const Interaction *i) const
Is the input kinematical point a physically allowed one?
double q2(bool selected=false) const
Definition: Kinematics.cxx:141
A class holding the Rein-Sehgal&#39;s helicity amplitudes.
bool IsProton(int pdgc)
Definition: PDGUtils.cxx:336
bool IsWeakNC(void) const
Singleton class to load &amp; serve tables of Fermi momentum constants.
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
static constexpr double A
Definition: Units.h:74
const FermiMomentumTable * GetTable(string name)
A class encapsulating an enumeration of interaction types (EM, Weak-CC, Weak-NC) and scattering types...
Definition: ProcessInfo.h:46
double fWcut
apply DIS/RES joining scheme &lt; Wcut
double T
Definition: FKR.h:46
double Rv
Definition: FKR.h:39
bool IsAntiNeutrino(int pdgc)
Definition: PDGUtils.cxx:118
A Neutrino Interaction Target. Is a transparent encapsulation of quite different physical systems suc...
Definition: Target.h:40
virtual void Configure(const Registry &config)
Definition: Algorithm.cxx:62
int ProbePdg(void) const
Definition: InitialState.h:64
const Algorithm * GetAlgorithm(const AlgId &algid)
Definition: AlgFactory.cxx:75
int OrbitalAngularMom(Resonance_t res)
orbital angular momentum
int Z(void) const
Definition: Target.h:68
#define pINFO
Definition: Messenger.h:62
Pure abstract base class. Defines the RSHelicityAmplModelI interface.
Resonance_t Resonance(void) const
Definition: XclsTag.h:69
bool IsEM(void) const
bool IsNeutralLepton(int pdgc)
Definition: PDGUtils.cxx:95
double C
Definition: FKR.h:44
double FermiMomentumForIsoscalarNucleonParametrization(const Target &target)
double Tplus
Definition: FKR.h:47
const XSecIntegratorI * fXSecIntegrator
double fXSecScaleCC
external CC xsec scaling factor
static AlgFactory * Instance()
Definition: AlgFactory.cxx:64
double B
Definition: FKR.h:43
double fOmega
FKR parameter Omega.
double Rplus
Definition: FKR.h:49
bool fUseRFGParametrization
use parametrization for fermi momentum insted of table?
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:65
bool fUsingDisResJoin
use a DIS/RES joining scheme?
const UInt_t kIAssumeFreeNucleon
Definition: Interaction.h:49
const XclsTag & ExclTag(void) const
Definition: Interaction.h:72
double fXSecScaleEM
external EM xsec scaling factor
double Tminus
Definition: FKR.h:48
int IonPdgCode(int A, int Z)
Definition: PDGUtils.cxx:71
double fSin48w
sin^4(Weingberg angle)
double fN0ResMaxNWidths
limits allowed phase space for n=0 res
double Jacobian(const Interaction *const i, KinePhaseSpace_t f, KinePhaseSpace_t t)
Definition: KineUtils.cxx:130
const InitialState & InitState(void) const
Definition: Interaction.h:69
const char * AsString(Resonance_t res)
resonance id -&gt; string
double fGnResMaxNWidths
limits allowed phase space for other res
const ProcessInfo & ProcInfo(void) const
Definition: Interaction.h:70
double FindClosestKF(int target_pdgc, int nucleon_pdgc) const
#define pNOTICE
Definition: Messenger.h:61
bool GetParamDef(const RgKey &name, T &p, const T &def) const
bool GetParam(const RgKey &name, T &p, bool is_top_call=true) const
const Target & Tgt(void) const
Definition: InitialState.h:66
void Configure(const Registry &config)
The GENIE Algorithm Factory.
Definition: AlgFactory.h:39
double fXSecScaleNC
external NC xsec scaling factor
double fMa2
(axial mass)^2
bool fUsingNuTauScaling
use NeuGEN nutau xsec reduction factors?
string fKFTable
table of Fermi momentum (kF) constants for various nuclei
virtual double Integrate(const XSecAlgorithmI *model, const Interaction *interaction) const =0
double ProbeE(RefFrame_t rf) const
Most commonly used PDG codes. A set of utility functions to handle PDG codes is provided in PDGUtils...
bool IsNegChargedLepton(int pdgc)
Definition: PDGUtils.cxx:139
double S
Definition: FKR.h:40
double Ta
Definition: FKR.h:41
const UInt_t kISkipProcessChk
if set, skip process validity checks
Definition: Interaction.h:47
int ResonanceIndex(Resonance_t res)
resonance idx, quark model / SU(6)
bool fUsePauliBlocking
account for Pauli blocking?
Initial State information.
Definition: InitialState.h:48
#define pDEBUG
Definition: Messenger.h:63
double fN2ResMaxNWidths
limits allowed phase space for n=2 res
Spline * fNuTauRdSpl
xsec reduction spline for nu_tau
const Algorithm * SubAlg(const RgKey &registry_key) const
Definition: Algorithm.cxx:345