MaCh3 DUNE 1.0.0
Reference Guide
Loading...
Searching...
No Matches
samplePDFDUNEAtm Class Reference

Base class for handling atmospheric samples. More...

#include <samplePDFDUNE/samplePDFDUNEAtm.h>

Inheritance diagram for samplePDFDUNEAtm:
Collaboration diagram for samplePDFDUNEAtm:

Public Types

enum  KinematicTypes {
  kTrueNeutrinoEnergy , kRecoNeutrinoEnergy , kTrueCosZ , kRecoCosZ ,
  kOscChannel , kMode
}
 Enum to identify kinematics. More...
 

Public Member Functions

 samplePDFDUNEAtm (std::string mc_version, covarianceXsec *xsec_cov, covarianceOsc *osc_cov)
 Constructor.
 
 ~samplePDFDUNEAtm ()
 destructor
 

Protected Member Functions

void Init ()
 Initialises object.
 
int setupExperimentMC (int iSample)
 Function to setup MC from file.
 
void setupFDMC (int iSample)
 Tells FD base which variables to point to/be set to.
 
void SetupWeightPointers ()
 Sets up pointers weights for each event (oscillation/xsec/etc.)
 
void SetupSplines ()
 Sets up splines.
 
double CalcXsecWeightFunc (int iSample, int iEvent)
 NOT IMPLEMENTED: Dunder method to calculate xsec weights.
 
void applyShifts (int iSample, int iEvent)
 NOT IMPLEMENTED: Apply kinematic shifts.
 
const double * GetPointerToKinematicParameter (KinematicTypes KinPar, int iSample, int iEvent)
 Returns pointer to kinemtatic parameter for event in Structs DUNE.
 
const double * GetPointerToKinematicParameter (double KinematicVariable, int iSample, int iEvent)
 Returns pointer to kinemtatic parameter for event in Structs DUNE.
 
const double * GetPointerToKinematicParameter (std::string KinematicParameter, int iSample, int iEvent)
 Returns pointer to kinemtatic parameter for event in Structs DUNE.
 
double ReturnKinematicParameter (double KinematicVariable, int iSample, int iEvent)
 Returns pointer to kinemtatic parameter for event in Structs DUNE.
 
double ReturnKinematicParameter (std::string KinematicParameter, int iSample, int iEvent)
 Returns pointer to kinemtatic parameter for event in Structs DUNE.
 
std::vector< double > ReturnKinematicParameterBinning (std::string KinematicParameterStr)
 Gets binning for a given parameter.
 
std::vector< double > ReturnKinematicParameterBinning (KinematicTypes KinPar)
 Gets binning for a given parameter.
 
int ReturnKinematicParameterFromString (std::string KinematicStr)
 Get kinematic parameter ID from string name.
 
std::string ReturnStringFromKinematicParameter (int KinematicVariable)
 Get kinematic parameter name from ID.
 

Protected Attributes

std::vector< struct dunemc_basedunemcSamples
 Array filled with MC samples for each oscillation channel.
 
bool IsELike
 Is the sample e-like.
 

Detailed Description

Base class for handling atmospheric samples.

Definition at line 11 of file samplePDFDUNEAtm.h.

Member Enumeration Documentation

◆ KinematicTypes

Enum to identify kinematics.

Enumerator
kTrueNeutrinoEnergy 
kRecoNeutrinoEnergy 
kTrueCosZ 
kRecoCosZ 
kOscChannel 
kMode 

Definition at line 23 of file samplePDFDUNEAtm.h.

Constructor & Destructor Documentation

◆ samplePDFDUNEAtm()

samplePDFDUNEAtm::samplePDFDUNEAtm ( std::string mc_version,
covarianceXsec * xsec_cov,
covarianceOsc * osc_cov )

Constructor.

Parameters
mc_versionConfiguration file
xsec_covcross-section covariance matrix
osc_covoscillation covariance matrix

Definition at line 15 of file samplePDFDUNEAtm.cpp.

15 : samplePDFFDBase(mc_version_, xsec_cov_, osc_cov_) {
16 Initialise();
17}

◆ ~samplePDFDUNEAtm()

samplePDFDUNEAtm::~samplePDFDUNEAtm ( )

destructor

Definition at line 19 of file samplePDFDUNEAtm.cpp.

19 {
20}

Member Function Documentation

◆ applyShifts()

void samplePDFDUNEAtm::applyShifts ( int iSample,
int iEvent )
inlineprotected

NOT IMPLEMENTED: Apply kinematic shifts.

Parameters
iSampleSample Number
iEventEvent number

Definition at line 54 of file samplePDFDUNEAtm.h.

54{}

◆ CalcXsecWeightFunc()

double samplePDFDUNEAtm::CalcXsecWeightFunc ( int iSample,
int iEvent )
inlineprotected

NOT IMPLEMENTED: Dunder method to calculate xsec weights.

Parameters
iSamplesample ID
iEventEvent number

Definition at line 49 of file samplePDFDUNEAtm.h.

49{return 1.;}

◆ GetPointerToKinematicParameter() [1/3]

const double * samplePDFDUNEAtm::GetPointerToKinematicParameter ( double KinematicVariable,
int iSample,
int iEvent )
protected

Returns pointer to kinemtatic parameter for event in Structs DUNE.

Parameters
KinematicVariableKinematic parameter as double (gets cast -> int)
iSampleSample ID
iEventEvent ID
Returns
Pointer to KinPar for a given event

Definition at line 174 of file samplePDFDUNEAtm.cpp.

174 {
175 KinematicTypes KinPar = (KinematicTypes) std::round(KinematicVariable);
176 return GetPointerToKinematicParameter(KinPar,iSample,iEvent);
177}
KinematicTypes
Enum to identify kinematics.
const double * GetPointerToKinematicParameter(KinematicTypes KinPar, int iSample, int iEvent)
Returns pointer to kinemtatic parameter for event in Structs DUNE.

◆ GetPointerToKinematicParameter() [2/3]

const double * samplePDFDUNEAtm::GetPointerToKinematicParameter ( KinematicTypes KinPar,
int iSample,
int iEvent )
protected

Returns pointer to kinemtatic parameter for event in Structs DUNE.

Parameters
KinParKinematic parameter enum val
iSampleSample ID
iEventEvent ID
Returns
Pointer to KinPar for a given event

Definition at line 144 of file samplePDFDUNEAtm.cpp.

144 {
145 double* KinematicValue;
146
147 switch (KinPar) {
149 KinematicValue = &(dunemcSamples[iSample].rw_etru[iEvent]);
150 break;
152 KinematicValue = &(dunemcSamples[iSample].rw_erec[iEvent]);
153 break;
154 case kTrueCosZ:
155 KinematicValue = &(dunemcSamples[iSample].rw_truecz[iEvent]);
156 break;
157 case kRecoCosZ:
158 KinematicValue = &(dunemcSamples[iSample].rw_theta[iEvent]);
159 break;
160 case kOscChannel:
161 KinematicValue = &(MCSamples[iSample].ChannelIndex);
162 break;
163 case kMode:
164 KinematicValue = &(dunemcSamples[iSample].mode[iEvent]);
165 break;
166 default:
167 MACH3LOG_ERROR("Unknown KinPar: {}",KinPar);
168 throw MaCh3Exception(__FILE__, __LINE__);
169 }
170
171 return KinematicValue;
172}
std::vector< struct dunemc_base > dunemcSamples
Array filled with MC samples for each oscillation channel.

◆ GetPointerToKinematicParameter() [3/3]

const double * samplePDFDUNEAtm::GetPointerToKinematicParameter ( std::string KinematicParameter,
int iSample,
int iEvent )
protected

Returns pointer to kinemtatic parameter for event in Structs DUNE.

Parameters
KinematicParameterKinematic parameter name as string (gets cast -> int)
iSampleSample ID
iEventEvent ID
Returns
Pointer to KinPar for a given event

Definition at line 179 of file samplePDFDUNEAtm.cpp.

179 {
180 KinematicTypes KinPar = static_cast<KinematicTypes>(ReturnKinematicParameterFromString(KinematicParameter));
181 return GetPointerToKinematicParameter(KinPar,iSample,iEvent);
182}
int ReturnKinematicParameterFromString(std::string KinematicStr)
Get kinematic parameter ID from string name.

◆ Init()

void samplePDFDUNEAtm::Init ( )
protected

Initialises object.

Definition at line 22 of file samplePDFDUNEAtm.cpp.

22 {
23 dunemcSamples.resize(nSamples,dunemc_base());
24
25 IsELike = SampleManager->raw()["SampleBools"]["IsELike"].as<bool>();
26}
bool IsELike
Is the sample e-like.

◆ ReturnKinematicParameter() [1/2]

double samplePDFDUNEAtm::ReturnKinematicParameter ( double KinematicVariable,
int iSample,
int iEvent )
protected

Returns pointer to kinemtatic parameter for event in Structs DUNE.

Parameters
KinematicVariableKinematic parameter ID as double (gets cast -> int)
iSampleSample ID
iEventEvent ID
Returns
Value of kinematic parameter corresponding for a given event

Definition at line 184 of file samplePDFDUNEAtm.cpp.

184 {
185 return *GetPointerToKinematicParameter(KinematicVariable, iSample, iEvent);
186}

◆ ReturnKinematicParameter() [2/2]

double samplePDFDUNEAtm::ReturnKinematicParameter ( std::string KinematicParameter,
int iSample,
int iEvent )
protected

Returns pointer to kinemtatic parameter for event in Structs DUNE.

Parameters
KinematicParameterKinematic parameter name as string (gets cast -> int)
iSampleSample ID
iEventEvent ID
Returns
Value of kinematic parameter corresponding for a given event

Definition at line 188 of file samplePDFDUNEAtm.cpp.

188 {
189 return *GetPointerToKinematicParameter(KinematicParameter, iSample, iEvent);
190}

◆ ReturnKinematicParameterBinning() [1/2]

std::vector< double > samplePDFDUNEAtm::ReturnKinematicParameterBinning ( KinematicTypes KinPar)
protected

Gets binning for a given parameter.

Parameters
KinParParameter ID
Returns
Vector containing parameter bins

Definition at line 197 of file samplePDFDUNEAtm.cpp.

197 {
198 std::vector<double> binningVector;
199
200 switch (KinPar) {
201
204 for (int i=0;i<20;i++) {
205 binningVector.emplace_back((double)i);
206 }
207 binningVector.emplace_back(100.);
208 binningVector.emplace_back(1000.);
209 break;
210
211 default:
212 MACH3LOG_ERROR("Unknown KinPar: {}",KinPar);
213 throw MaCh3Exception(__FILE__, __LINE__);
214 }
215
216 return binningVector;
217}

◆ ReturnKinematicParameterBinning() [2/2]

std::vector< double > samplePDFDUNEAtm::ReturnKinematicParameterBinning ( std::string KinematicParameterStr)
protected

Gets binning for a given parameter.

Parameters
KinematicParameterStrParameter name
Returns
Vector containing parameter bins

Definition at line 192 of file samplePDFDUNEAtm.cpp.

192 {
193 KinematicTypes KinPar = static_cast<KinematicTypes>(ReturnKinematicParameterFromString(KinematicParameterStr));
194 return ReturnKinematicParameterBinning(KinPar);
195}
std::vector< double > ReturnKinematicParameterBinning(std::string KinematicParameterStr)
Gets binning for a given parameter.

◆ ReturnKinematicParameterFromString()

int samplePDFDUNEAtm::ReturnKinematicParameterFromString ( std::string KinematicStr)
inlineprotected

Get kinematic parameter ID from string name.

Parameters
KinematicStr
Returns
Parameter ID

Definition at line 219 of file samplePDFDUNEAtm.cpp.

219 {
220 int ReturnVal;
221
222 if (KinematicParameterStr == "TrueNeutrinoEnergy") {ReturnVal = kTrueNeutrinoEnergy;}
223 else if (KinematicParameterStr == "RecoNeutrinoEnergy") {ReturnVal = kRecoNeutrinoEnergy;}
224 else if (KinematicParameterStr == "TrueCosineZ") {ReturnVal = kTrueCosZ;}
225 else if (KinematicParameterStr == "RecoCosineZ") {ReturnVal = kRecoCosZ;}
226 else if (KinematicParameterStr == "OscChannel") {ReturnVal = kOscChannel;}
227 else if (KinematicParameterStr == "Mode") {ReturnVal = kMode;}
228 else {
229 MACH3LOG_ERROR("KinematicParameterStr: {} not found",KinematicParameterStr);
230 throw MaCh3Exception(__FILE__, __LINE__);
231 }
232
233 return ReturnVal;
234}

◆ ReturnStringFromKinematicParameter()

std::string samplePDFDUNEAtm::ReturnStringFromKinematicParameter ( int KinematicVariable)
inlineprotected

Get kinematic parameter name from ID.

Parameters
KinematicVariableParameter ID
Returns
Parameter Name

Definition at line 236 of file samplePDFDUNEAtm.cpp.

236 {
237 return "";
238}

◆ setupExperimentMC()

int samplePDFDUNEAtm::setupExperimentMC ( int iSample)
protected

Function to setup MC from file.

Parameters
iSamplesample ID
Returns
Total number of events

Definition at line 44 of file samplePDFDUNEAtm.cpp.

44 {
45 int CurrErrorLevel = gErrorIgnoreLevel;
46 gErrorIgnoreLevel = kFatal;
47
48 caf::StandardRecord* sr = new caf::StandardRecord();
49 dunemc_base* duneobj = &dunemcSamples[iSample];
50
51 std::string FileName = mc_files[iSample];
52 MACH3LOG_INFO("Reading File: {}",FileName);
53 TFile* File = TFile::Open(FileName.c_str());
54 if (!File || File->IsZombie()) {
55 MACH3LOG_ERROR("Did not find File: {}",FileName);
56 throw MaCh3Exception(__FILE__, __LINE__);
57 }
58 TTree* Tree = File->Get<TTree>("cafTree");
59 if (!Tree){
60 MACH3LOG_ERROR("Did not find Tree::cafTree in File: {}",FileName);
61 throw MaCh3Exception(__FILE__, __LINE__);
62 }
63
64 Tree->SetBranchStatus("*", 1);
65 Tree->SetBranchAddress("rec", &sr);
66
67 duneobj->nEvents = Tree->GetEntries();
68 duneobj->norm_s = 1;
69 duneobj->pot_s = 1;
70
71 duneobj->nupdg = new int[duneobj->nEvents];
72 duneobj->nupdgUnosc = new int[duneobj->nEvents];
73
74 duneobj->mode = new double[duneobj->nEvents];
75 duneobj->rw_isCC = new int[duneobj->nEvents];
76 duneobj->Target = new int[duneobj->nEvents];
77
78 duneobj->rw_etru = new double[duneobj->nEvents];
79 duneobj->rw_truecz = new double[duneobj->nEvents];
80 duneobj->flux_w = new double[duneobj->nEvents];
81 duneobj->rw_erec = new double[duneobj->nEvents];
82 duneobj->rw_theta = new double[duneobj->nEvents];
83
84 for (int iEvent=0;iEvent<duneobj->nEvents;iEvent++) {
85 Tree->GetEntry(iEvent);
86
87 if ((iEvent % (duneobj->nEvents/10))==0) {
88 MACH3LOG_INFO("\tProcessing event: {}/{}",iEvent,duneobj->nEvents);
89 }
90
91 duneobj->nupdg[iEvent] = sample_nupdg[iSample];
92 duneobj->nupdgUnosc[iEvent] = sample_nupdgunosc[iSample];
93
94 duneobj->mode[iEvent] = (double)SIMBMode_ToMaCh3Mode(TMath::Abs(sr->mc.nu[0].mode),(int)sr->mc.nu[0].iscc);
95 duneobj->rw_isCC[iEvent] = sr->mc.nu[0].iscc;
96 duneobj->Target[iEvent] = 40;
97
98 duneobj->rw_etru[iEvent] = (double)(sr->mc.nu[0].E);
99
100 TVector3 TrueNuMomentumVector = (TVector3(sr->mc.nu[0].momentum.X(),sr->mc.nu[0].momentum.Y(),sr->mc.nu[0].momentum.Z())).Unit();
101 duneobj->rw_truecz[iEvent] = TrueNuMomentumVector.Y();
102
103 duneobj->flux_w[iEvent] = sr->mc.nu[0].genweight;
104
105 TVector3 RecoNuMomentumVector;
106 if (IsELike) {
107 duneobj->rw_erec[iEvent] = sr->common.ixn.pandora[0].Enu.e_calo;
108 RecoNuMomentumVector = (TVector3(sr->common.ixn.pandora[0].dir.heshw.X(),sr->common.ixn.pandora[0].dir.heshw.Y(),sr->common.ixn.pandora[0].dir.heshw.Z())).Unit();
109 } else {
110 duneobj->rw_erec[iEvent] = sr->common.ixn.pandora[0].Enu.lep_calo;
111 RecoNuMomentumVector = (TVector3(sr->common.ixn.pandora[0].dir.lngtrk.X(),sr->common.ixn.pandora[0].dir.lngtrk.Y(),sr->common.ixn.pandora[0].dir.lngtrk.Z())).Unit();
112 }
113 duneobj->rw_theta[iEvent] = RecoNuMomentumVector.Y();
114
115 }
116
117 delete Tree;
118 delete File;
119
120 gErrorIgnoreLevel = CurrErrorLevel;
121
122 return duneobj->nEvents;
123}
int SIMBMode_ToMaCh3Mode(int SIMB_mode, int isCC)
int * nupdgUnosc
Definition StructsDUNE.h:15
double * mode
Definition StructsDUNE.h:65
double pot_s
Definition StructsDUNE.h:59
double * rw_erec
Definition StructsDUNE.h:17
int * rw_isCC
Definition StructsDUNE.h:47
double * flux_w
Definition StructsDUNE.h:63
double * rw_etru
Definition StructsDUNE.h:35
double * rw_theta
Definition StructsDUNE.h:37
double * rw_truecz
Definition StructsDUNE.h:68
double norm_s
Definition StructsDUNE.h:60

◆ setupFDMC()

void samplePDFDUNEAtm::setupFDMC ( int iSample)
protected

Tells FD base which variables to point to/be set to.

Parameters
iSampleSample ID

Definition at line 125 of file samplePDFDUNEAtm.cpp.

125 {
126 dunemc_base *duneobj = &(dunemcSamples[iSample]);
127 FarDetectorCoreInfo *fdobj = &(MCSamples[iSample]);
128
129 //Make sure that this is only set if you're an atmoshperic object
130 fdobj->rw_truecz.resize(fdobj->nEvents);
131
132 for(int iEvent = 0 ;iEvent < fdobj->nEvents ; ++iEvent) {
133 fdobj->rw_etru[iEvent] = &(duneobj->rw_etru[iEvent]);
134 fdobj->mode[iEvent] = &(duneobj->mode[iEvent]);
135 fdobj->Target[iEvent] = &(duneobj->Target[iEvent]);
136 fdobj->isNC[iEvent] = !duneobj->rw_isCC[iEvent];
137 fdobj->nupdg[iEvent] = &(duneobj->nupdg[iEvent]);
138 fdobj->nupdgUnosc[iEvent] = &(duneobj->nupdgUnosc[iEvent]);
139
140 fdobj->rw_truecz[iEvent] = &(duneobj->rw_truecz[iEvent]);
141 }
142}

◆ SetupSplines()

void samplePDFDUNEAtm::SetupSplines ( )
protected

Sets up splines.

Definition at line 28 of file samplePDFDUNEAtm.cpp.

28 {
29 SplineHandler = nullptr;
30}

◆ SetupWeightPointers()

void samplePDFDUNEAtm::SetupWeightPointers ( )
protected

Sets up pointers weights for each event (oscillation/xsec/etc.)

Definition at line 32 of file samplePDFDUNEAtm.cpp.

32 {
33 for (int i = 0; i < (int)dunemcSamples.size(); ++i) {
34 for (int j = 0; j < dunemcSamples[i].nEvents; ++j) {
35 MCSamples[i].ntotal_weight_pointers[j] = 3;
36 MCSamples[i].total_weight_pointers[j].resize(MCSamples[i].ntotal_weight_pointers[j]);
37 MCSamples[i].total_weight_pointers[j][0] = &(dunemcSamples[i].flux_w[j]);
38 MCSamples[i].total_weight_pointers[j][1] = MCSamples[i].osc_w_pointer[j];
39 MCSamples[i].total_weight_pointers[j][2] = &(MCSamples[i].xsec_w[j]);
40 }
41 }
42}

Member Data Documentation

◆ dunemcSamples

std::vector<struct dunemc_base> samplePDFDUNEAtm::dunemcSamples
protected

Array filled with MC samples for each oscillation channel.

Definition at line 112 of file samplePDFDUNEAtm.h.

◆ IsELike

bool samplePDFDUNEAtm::IsELike
protected

Is the sample e-like.

Definition at line 115 of file samplePDFDUNEAtm.h.


The documentation for this class was generated from the following files: