GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
COHDNuXSec.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  Author: Iker de Icaza <i.de-icaza-astiz \at sussex.ac.uk>
7  University of Sussex
8 
9  Costas Andreopoulos <c.andreopoulos \at cern.ch>
10  University of Liverpool
11 */
12 //____________________________________________________________________________
13 
14 #include <limits>
15 
16 #include <TMath.h>
17 #include <Math/IFunction.h>
18 #include <Math/Integrator.h>
19 
26 
27 using namespace genie;
28 using namespace genie::constants;
29 using namespace genie::utils;
30 
31 //____________________________________________________________________________
33 XSecIntegratorI("genie::COHDNuXSec")
34 {
35 
36 }
37 //____________________________________________________________________________
38 COHDNuXSec::COHDNuXSec(string config) :
39 XSecIntegratorI("genie::COHDNuXSec", config)
40 {
41 
42 }
43 //____________________________________________________________________________
45 {
46 
47 }
48 //____________________________________________________________________________
50  const XSecAlgorithmI * model, const Interaction * in) const
51 {
52  if(!model->ValidProcess(in) ) return 0.;
53  if(!in->PhaseSpace().IsAboveThreshold()) return 0.;
54 
55  Interaction interaction(*in);
56  interaction.SetBit(kISkipProcessChk);
57  interaction.SetBit(kISkipKinematicChk);
58 
59  ROOT::Math::IntegrationOneDim::Type ig_type =
61 
62  utils::gsl::dXSec_dEDNu_E func( model, & interaction, fDNuMass );
63  Range1D_t DNuEnergy = func.IntegrationRange();
64 
65  ROOT::Math::Integrator ig( func, ig_type, fGSLAbsTol, fGSLRelTol, fGSLMaxSizeOfSubintervals );
66  double xsec = ig.Integral(DNuEnergy.min, DNuEnergy.max) * (1E-38 * units::cm2); // units: GeV^-2
67 
68  const InitialState & init_state = in->InitState();
69  double Ev = init_state.ProbeE(kRfLab);
70  LOG("COHDNuXSec", pINFO)
71  << "XSec[COHDNu] (E = " << Ev << " GeV) = " << xsec/(units::cm2) << " cm^2";
72 
73  return xsec;
74 }
75 //____________________________________________________________________________
76 void COHDNuXSec::Configure(const Registry & config)
77 {
78  Algorithm::Configure(config);
79  this->LoadConfig();
80 }
81 //____________________________________________________________________________
82 void COHDNuXSec::Configure(string config)
83 {
84  Algorithm::Configure(config);
85  this->LoadConfig();
86 }
87 //____________________________________________________________________________
89 {
90  fDNuMass = 0.;
91  this->GetParam("Dark-NeutrinoMass", fDNuMass);
92 
93  // Get GSL integration type & relative tolerance
94  this->GetParamDef("gsl-integration-type", fGSLIntgType, string("adaptive"));
95  this->GetParamDef("gsl-relative-tolerance", fGSLRelTol, 1E-2);
96  this->GetParamDef("gsl-absolute-tolerance", fGSLAbsTol, std::numeric_limits<double>::max());
97 
98  int max_subint = 0 ;
99  this->GetParamDef("gsl-max-subintervals", max_subint, 500);
100  fGSLMaxSizeOfSubintervals = (unsigned int) max_subint ;
101 
102  // unused parameters from XSecIntegratorI
103  fGSLMaxEval = fGSLMinEval = 0 ;
104 }
105 //____________________________________________________________________________
Cross Section Calculation Interface.
virtual ~COHDNuXSec()
Definition: COHDNuXSec.cxx:44
const KPhaseSpace & PhaseSpace(void) const
Definition: Interaction.h:73
ROOT::Math::IntegrationOneDim::Type Integration1DimTypeFromString(string type)
Definition: GSLUtils.cxx:23
Range1D_t IntegrationRange(void) const
string fGSLIntgType
name of GSL numerical integrator
Cross Section Integrator Interface.
void LoadConfig(void)
Definition: COHDNuXSec.cxx:88
A simple [min,max] interval for doubles.
Definition: Range1.h:42
void Configure(const Registry &config)
Definition: COHDNuXSec.cxx:76
Summary information for an interaction.
Definition: Interaction.h:56
#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 cm2
Definition: Units.h:69
virtual void Configure(const Registry &config)
Definition: Algorithm.cxx:62
double func(double x, double y)
#define pINFO
Definition: Messenger.h:62
const UInt_t kISkipKinematicChk
if set, skip kinematic validity checks
Definition: Interaction.h:48
virtual bool ValidProcess(const Interaction *i) const =0
Can this cross section algorithm handle the input process?
double Integrate(const XSecAlgorithmI *model, const Interaction *i) const
Definition: COHDNuXSec.cxx:49
double max
Definition: Range1.h:53
int fGSLMaxEval
GSL max evaluations.
unsigned int fGSLMaxSizeOfSubintervals
GSL maximum number of sub-intervals for 1D integrator.
A registry. Provides the container for algorithm configuration parameters.
Definition: Registry.h:65
bool IsAboveThreshold(void) const
Checks whether the interaction is above the energy threshold.
const InitialState & InitState(void) const
Definition: Interaction.h:69
double min
Definition: Range1.h:52
bool GetParamDef(const RgKey &name, T &p, const T &def) const
bool GetParam(const RgKey &name, T &p, bool is_top_call=true) const
double ProbeE(RefFrame_t rf) const
int fGSLMinEval
GSL min evaluations. Ignored by some integrators.
const UInt_t kISkipProcessChk
if set, skip process validity checks
Definition: Interaction.h:47
Initial State information.
Definition: InitialState.h:48
double fGSLRelTol
required relative tolerance (error)