GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions
gtestXSecCEvNS.cxx File Reference
#include <TFile.h>
#include <TGraph.h>
#include "Framework/Algorithm/Algorithm.h"
#include "Framework/Algorithm/AlgFactory.h"
#include "Framework/Conventions/Units.h"
#include "Framework/Interaction/Interaction.h"
#include "Framework/EventGen/XSecAlgorithmI.h"
#include "Framework/Messenger/Messenger.h"
#include "Framework/Utils/RunOpt.h"
Include dependency graph for gtestXSecCEvNS.cxx:

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 31 of file gtestXSecCEvNS.cxx.

References genie::RunOpt::BuildTune(), genie::Interaction::CEvNS(), genie::units::cm2, genie::AlgFactory::GetAlgorithm(), genie::Interaction::InitStatePtr(), genie::RunOpt::Instance(), genie::AlgFactory::Instance(), genie::XSecAlgorithmI::Integral(), genie::Interaction::KinePtr(), LOG, pFATAL, pINFO, pNOTICE, genie::utils::kinematics::Q2(), genie::RunOpt::ReadFromCommandLine(), genie::InitialState::SetProbeE(), and genie::Kinematics::SetQ2().

32 {
33  // <temp/>
34  RunOpt::Instance()->ReadFromCommandLine(argc,argv);
35  if ( ! RunOpt::Instance()->Tune() ) {
36  LOG("test", pFATAL) << " No TuneId in RunOption";
37  exit(-1);
38  }
39  RunOpt::Instance()->BuildTune();
40  // </temp>
41 
42  LOG("test",pINFO) << "Testing the CEvNS cross-section calculation";
43 
44  // Instantiate the coherent elastic cross-section model
45  AlgFactory * algf = AlgFactory::Instance();
46  AlgId id("genie::PattonCEvNSPXSec","Default");
47  const Algorithm * alg = algf->GetAlgorithm(id);
48  LOG("test", pINFO) << *alg;
49  const XSecAlgorithmI * xsecalg = dynamic_cast<const XSecAlgorithmI*>(alg);
50 
51  const int target = 1000180400;
52  const int probe = 12;
53  double E = 0.010; // GeV
54  double Q2 = 0.0001; //GeV^2
55 
56  Interaction * interaction = Interaction::CEvNS(target, probe, E);
57  interaction->KinePtr()->SetQ2(Q2);
58 
59 /*
60  // test differential cross-section calculation at a single point
61  double dxsec = xsecalg->XSec(interaction,kPSQ2fE);
62 
63  LOG("test", pNOTICE)
64  << "dxsec[CEvNS; target = " << target << "]/dQ2"
65  << "(E = " << E << " GeV, Q2 = " << Q2 << " GeV^2) = "
66  << dxsec/(units::cm2) << " cm^2/GeV^2";
67 */
68 
69 /*
70  // test integrated cross-section calculation at a single point
71  double xsec = xsecalg->Integral(interaction);
72 
73  LOG("test", pNOTICE)
74  << "xsec[CEvNS; target = " << target << "]"
75  << "(E = " << E << " GeV) = " << xsec/(units::cm2) << " cm^2";
76 */
77 
78  // test integrated cross-section calculation at a range of energies
79  // and compare with published predictions
80  const int n = 100;
81  const double Emin = 0.005;
82  const double Emax = 0.055;
83  const double dE = (Emax-Emin)/(n-1);
84 
85  double e_array[n] = {0};
86  double genie_xsec_array[n] = {0};
87 
88  for(int i=0; i<n; i++) {
89  double E_current = Emin + i*dE;
90  interaction->InitStatePtr()->SetProbeE(E_current);
91  double xsec_current = xsecalg->Integral(interaction)/(units::cm2);
92  e_array[i] = E_current;
93  genie_xsec_array[i] = xsec_current;
94  LOG("test", pNOTICE)
95  << "xsec[CEvNS; target = " << target << "]"
96  << "(E = " << E_current << " GeV) = " << xsec_current << " cm^2";
97  }
98  TGraph * genie_xsec = new TGraph(n, e_array, genie_xsec_array);
99  TGraph * published_xsec = new TGraph("$GENIE/src/contrib/cevns/cevns_arXiv180309183_Ar40_prediction.data");
100  published_xsec->SetMarkerStyle(20);
101  published_xsec->SetMarkerColor(kRed);
102  genie_xsec->SetLineColor(kBlue);
103  TFile f("cevns.root","recreate");
104  genie_xsec->Write("genie_xsec_nuAr40");
105  published_xsec->Write("published_xsec_nuAr40");
106  f.Close();
107  return 0;
108 }
Cross Section Calculation Interface.
double Q2(const Interaction *const i)
Definition: KineUtils.cxx:1077
void SetQ2(double Q2, bool selected=false)
Definition: Kinematics.cxx:255
Kinematics * KinePtr(void) const
Definition: Interaction.h:76
#define pFATAL
Definition: Messenger.h:56
Algorithm abstract base class.
Definition: Algorithm.h:54
void SetProbeE(double E)
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
const Algorithm * GetAlgorithm(const AlgId &algid)
Definition: AlgFactory.cxx:75
#define pINFO
Definition: Messenger.h:62
Algorithm ID (algorithm name + configuration set name)
Definition: AlgId.h:34
InitialState * InitStatePtr(void) const
Definition: Interaction.h:74
#define pNOTICE
Definition: Messenger.h:61
The GENIE Algorithm Factory.
Definition: AlgFactory.h:39
virtual double Integral(const Interaction *i) const =0