All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Private Member Functions | Private Attributes | List of all members
FluxReaderAna Class Reference
Inheritance diagram for FluxReaderAna:

Public Member Functions

 FluxReaderAna (fhicl::ParameterSet const &p)
 
 FluxReaderAna (FluxReaderAna const &)=delete
 
 FluxReaderAna (FluxReaderAna &&)=delete
 
FluxReaderAnaoperator= (FluxReaderAna const &)=delete
 
FluxReaderAnaoperator= (FluxReaderAna &&)=delete
 
void analyze (art::Event const &e) override
 

Private Member Functions

TVector3 GetIntersection (TVector3 nu_pos, TVector3 nu_dir, float z_location=0)
 Returns the intersection point with the front face of the TPC. More...
 

Private Attributes

std::string _flux_label
 Label for flux dataproduct (to be set via fcl) More...
 
float _x_shift
 
float _baseline
 
float _nu_intersection_z
 Z position where to evaluate the flux (to be set via fcl) More...
 
float _nu_other_intersection_z
 Additional Z position where to evaluate the flux (to be set via fcl) More...
 
bool _apply_position_cuts
 
float _x_cut
 
float _y_cut
 
const TDatabasePDG * _pdg_database = TDatabasePDG::Instance()
 
TTree * _tree
 
bool _nu_hit
 
int _nu_pdg
 True if the neutrino hit the requested volumes. More...
 
float _nu_e
 PDG of neutrino. More...
 
float _nu_t
 Energy of the neutrino. More...
 
float _nu_w
 Time of the neutrino. More...
 
float _nu_x
 Neutrino weight. More...
 
float _nu_y
 X poisition of neutrino at the front face of the TPC. More...
 
float _nu_z
 Y poisition of neutrino at the front face of the TPC. More...
 
float _nu_vx
 Z poisition of neutrino at the front face of the TPC. More...
 
float _nu_vy
 X poisition of neutrino at neutrino production point. More...
 
float _nu_vz
 Y poisition of neutrino at neutrino production point. More...
 
float _nu_px
 Z poisition of neutrino at neutrino production point. More...
 
float _nu_py
 X momentum of neutrino. More...
 
float _nu_pz
 Y momentum of neutrino. More...
 
int _nu_decay
 Z momentum of neutrino. More...
 
float _nu_dk2gen
 Neutrino parent decay code. More...
 
float _nu_p_angle
 Distance from decay to ray origin. More...
 
int _nu_p_type
 Angle between neutrino and parent direction. More...
 
float _nu_p_dpx
 Neutrino parent pdg. More...
 
float _nu_p_dpy
 Neutrino parent momentum x. More...
 
float _nu_p_dpz
 Neutrino parent momentum x. More...
 
float _nu_p_e
 Neutrino parent momentum x. More...
 
float _nu_r
 Neutrino parent energy. More...
 
float _nu_oaa
 Neutrino r. More...
 
float _nu_other_x
 Neutrino off axis angle. More...
 
float _nu_other_y
 X poisition of neutrino at the front face of the TPC (other) More...
 
float _nu_other_z
 Y poisition of neutrino at the front face of the TPC (other) More...
 
float _nu_other_t
 Z poisition of neutrino at the front face of the TPC (other) More...
 

Detailed Description

Definition at line 41 of file FluxReaderAna_module.cc.

Constructor & Destructor Documentation

FluxReaderAna::FluxReaderAna ( fhicl::ParameterSet const &  p)
explicit

Definition at line 107 of file FluxReaderAna_module.cc.

108  : EDAnalyzer{p}
109 {
110 
111  _flux_label = p.get<std::string>("FluxLabel", "flux");
112 
113  _baseline = p.get<float>("Baseline"); // cm, distance from detector to target position
114  _x_shift = p.get<float>("XShift", 0.); // cm, detector shift along X w.r.t. beamline coordinate system
115 
116  _nu_intersection_z = p.get<float>("NuIntersectionZ", 0.);
117  _nu_other_intersection_z = p.get<float>("NuOtherIntersectionZ");
118  // 49000 is ICARUS location in SBND coordinate system (600 - 110)
119  // 36000 is MicroBooNE location in SBND coordinate system (470 - 110)
120 
121  _apply_position_cuts = p.get<bool>("ApplyPositionCuts", false);
122  _x_cut = p.get<float>("XCut"); // cm
123  _y_cut = p.get<float>("YCut"); // cm
124 
125 
126  art::ServiceHandle<art::TFileService> fs;
127  _tree = fs->make<TTree>("tree", "");
128  _tree->Branch("nu_hit", &_nu_hit, "nu_hit/O");
129  _tree->Branch("nu_pdg", &_nu_pdg, "nu_pdg/I");
130  _tree->Branch("nu_e", &_nu_e, "nu_e/F");
131  _tree->Branch("nu_t", &_nu_t, "nu_t/F");
132  _tree->Branch("nu_w", &_nu_w, "nu_w/F");
133  _tree->Branch("nu_x", &_nu_x, "nu_x/F");
134  _tree->Branch("nu_y", &_nu_y, "nu_y/F");
135  _tree->Branch("nu_z", &_nu_z, "nu_z/F");
136  _tree->Branch("nu_vx", &_nu_vx, "nu_vx/F");
137  _tree->Branch("nu_vy", &_nu_vy, "nu_vy/F");
138  _tree->Branch("nu_vz", &_nu_vz, "nu_vz/F");
139  _tree->Branch("nu_px", &_nu_px, "nu_px/F");
140  _tree->Branch("nu_py", &_nu_py, "nu_py/F");
141  _tree->Branch("nu_pz", &_nu_pz, "nu_pz/F");
142  _tree->Branch("nu_decay", &_nu_decay, "nu_decay/I");
143  _tree->Branch("nu_dk2gen", &_nu_dk2gen, "nu_dk2gen/F");
144  _tree->Branch("nu_p_angle", &_nu_p_angle, "nu_p_angle/F");
145  _tree->Branch("nu_p_type", &_nu_p_type, "nu_p_type/I");
146  _tree->Branch("nu_p_dpx", &_nu_p_dpx, "nu_p_dpx/F");
147  _tree->Branch("nu_p_dpy", &_nu_p_dpy, "nu_p_dpy/F");
148  _tree->Branch("nu_p_dpz", &_nu_p_dpz, "nu_p_dpz/F");
149  _tree->Branch("nu_p_e", &_nu_p_e, "nu_p_e/F");
150  _tree->Branch("nu_r", &_nu_r, "nu_r/F");
151  _tree->Branch("nu_oaa", &_nu_oaa, "nu_oaa/F");
152 
153  _tree->Branch("nu_other_x", &_nu_other_x, "nu_other_x/F");
154  _tree->Branch("nu_other_y", &_nu_other_y, "nu_other_y/F");
155  _tree->Branch("nu_other_z", &_nu_other_z, "nu_other_z/F");
156  _tree->Branch("nu_other_t", &_nu_other_t, "nu_other_t/F");
157 }
float _nu_vz
Y poisition of neutrino at neutrino production point.
pdgs p
Definition: selectors.fcl:22
std::string _flux_label
Label for flux dataproduct (to be set via fcl)
float _nu_p_dpy
Neutrino parent momentum x.
float _nu_other_x
Neutrino off axis angle.
float _nu_p_dpx
Neutrino parent pdg.
float _nu_py
X momentum of neutrino.
float _nu_vy
X poisition of neutrino at neutrino production point.
float _nu_dk2gen
Neutrino parent decay code.
float _nu_other_intersection_z
Additional Z position where to evaluate the flux (to be set via fcl)
int _nu_p_type
Angle between neutrino and parent direction.
float _nu_other_t
Z poisition of neutrino at the front face of the TPC (other)
float _nu_other_y
X poisition of neutrino at the front face of the TPC (other)
float _nu_oaa
Neutrino r.
float _nu_r
Neutrino parent energy.
float _nu_px
Z poisition of neutrino at neutrino production point.
float _nu_y
X poisition of neutrino at the front face of the TPC.
float _nu_z
Y poisition of neutrino at the front face of the TPC.
float _nu_intersection_z
Z position where to evaluate the flux (to be set via fcl)
float _nu_x
Neutrino weight.
int _nu_pdg
True if the neutrino hit the requested volumes.
float _nu_p_e
Neutrino parent momentum x.
float _nu_w
Time of the neutrino.
float _nu_t
Energy of the neutrino.
float _nu_vx
Z poisition of neutrino at the front face of the TPC.
int _nu_decay
Z momentum of neutrino.
float _nu_p_angle
Distance from decay to ray origin.
float _nu_p_dpz
Neutrino parent momentum x.
float _nu_pz
Y momentum of neutrino.
float _nu_e
PDG of neutrino.
float _nu_other_z
Y poisition of neutrino at the front face of the TPC (other)
FluxReaderAna::FluxReaderAna ( FluxReaderAna const &  )
delete
FluxReaderAna::FluxReaderAna ( FluxReaderAna &&  )
delete

Member Function Documentation

void FluxReaderAna::analyze ( art::Event const &  e)
override

Definition at line 159 of file FluxReaderAna_module.cc.

160 {
161 
162  art::Handle< std::vector<simb::MCFlux> > mcFluxHandle;
163  e.getByLabel(_flux_label, mcFluxHandle);
164  std::vector<simb::MCFlux> const& fluxlist = *mcFluxHandle;
165 
166  art::Handle< std::vector<simb::MCTruth> > mctruthHandle;
167  e.getByLabel(_flux_label, mctruthHandle);
168  std::vector<simb::MCTruth> const& mclist = *mctruthHandle;
169 
170  for(unsigned int inu = 0; inu < mclist.size(); inu++) {
171  simb::MCParticle nu = mclist[inu].GetNeutrino().Nu();
172  simb::MCFlux flux = fluxlist[inu];
173 
174  // std::cout << "This neutrino has vtx " << nu.Vx() << ", " << nu.Vy() << ", " << nu.Vz() << std::endl;
175  // std::cout << "This neutrino has dir " << nu.Px() << ", " << nu.Py() << ", " << nu.Pz() << std::endl;
176 
177  // Calculate the position of the neutrino at Z = _nu_intersection_z
178  TVector3 intersection = GetIntersection(TVector3(nu.Vx(),nu.Vy(),nu.Vz()),
179  TVector3(nu.Px(),nu.Py(),nu.Pz()),
181  // Calculate the distance from the original to the new position
182  float dist = (TVector3(nu.Vx(),nu.Vy(),nu.Vz()) - intersection).Mag();
183 
184  // Do the same for an additional Z
185  TVector3 intersection_other = GetIntersection(TVector3(nu.Vx(),nu.Vy(),nu.Vz()),
186  TVector3(nu.Px(),nu.Py(),nu.Pz()),
188  float other_dist = (TVector3(nu.Vx(),nu.Vy(),nu.Vz()) - intersection_other).Mag();
189 
190  float light_speed = TMath::C() / 1e7; // cm / ns
191 
192  _nu_pdg = nu.PdgCode();
193  _nu_e = nu.E();
194  _nu_t = nu.T() + dist / light_speed;
195  _nu_w = flux.fnimpwt;
196  _nu_x = intersection.X();
197  _nu_y = intersection.Y();
198  _nu_z = intersection.Z();
199  _nu_vx = nu.Vx();
200  _nu_vy = nu.Vy();
201  _nu_vz = nu.Vz();
202  _nu_px = nu.Px();
203  _nu_py = nu.Py();
204  _nu_pz = nu.Pz();
205  _nu_decay = flux.fndecay;
206  _nu_dk2gen = flux.fdk2gen;
207  _nu_p_type = flux.fptype;
208  _nu_p_dpx = flux.fpdpx;
209  _nu_p_dpy = flux.fpdpy;
210  _nu_p_dpz = flux.fpdpz;
211  float parent_mass = _pdg_database->GetParticle(_nu_p_type)->Mass();
212  _nu_p_e = std::sqrt(flux.fpdpx * flux.fpdpx +
213  flux.fpdpy * flux.fpdpy +
214  flux.fpdpz * flux.fpdpz +
215  parent_mass * parent_mass);
216  _nu_p_angle = TVector3(_nu_px, _nu_py, _nu_pz).Angle(TVector3(_nu_p_dpx, _nu_p_dpy, _nu_p_dpz));
217  _nu_r = std::sqrt((_nu_x - _x_shift) * (_nu_x - _x_shift) + _nu_y * _nu_y);
218  _nu_oaa = std::atan(_nu_r * _nu_r / _baseline);
219 
220  _nu_other_x = intersection_other.X();
221  _nu_other_y = intersection_other.Y();
222  _nu_other_z = intersection_other.Z();
223  _nu_other_t = nu.T() + other_dist / light_speed;
224 
225  if (_apply_position_cuts) {
226  if (_nu_x < -_x_cut || _nu_x > _x_cut || _nu_y < -_y_cut || _nu_y > _y_cut) {
227  continue;
228  }
229  }
230 
231  _tree->Fill();
232  }
233 }
float _nu_vz
Y poisition of neutrino at neutrino production point.
std::string _flux_label
Label for flux dataproduct (to be set via fcl)
float _nu_p_dpy
Neutrino parent momentum x.
float _nu_other_x
Neutrino off axis angle.
float _nu_p_dpx
Neutrino parent pdg.
float _nu_py
X momentum of neutrino.
float _nu_vy
X poisition of neutrino at neutrino production point.
float _nu_dk2gen
Neutrino parent decay code.
float _nu_other_intersection_z
Additional Z position where to evaluate the flux (to be set via fcl)
int _nu_p_type
Angle between neutrino and parent direction.
float _nu_other_t
Z poisition of neutrino at the front face of the TPC (other)
float _nu_other_y
X poisition of neutrino at the front face of the TPC (other)
float _nu_oaa
Neutrino r.
float _nu_r
Neutrino parent energy.
float _nu_px
Z poisition of neutrino at neutrino production point.
TVector3 GetIntersection(TVector3 nu_pos, TVector3 nu_dir, float z_location=0)
Returns the intersection point with the front face of the TPC.
float _nu_y
X poisition of neutrino at the front face of the TPC.
float _nu_z
Y poisition of neutrino at the front face of the TPC.
float _nu_intersection_z
Z position where to evaluate the flux (to be set via fcl)
float _nu_x
Neutrino weight.
int _nu_pdg
True if the neutrino hit the requested volumes.
float _nu_p_e
Neutrino parent momentum x.
float _nu_w
Time of the neutrino.
float _nu_t
Energy of the neutrino.
const TDatabasePDG * _pdg_database
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
float _nu_vx
Z poisition of neutrino at the front face of the TPC.
int _nu_decay
Z momentum of neutrino.
do i e
float _nu_p_angle
Distance from decay to ray origin.
float _nu_p_dpz
Neutrino parent momentum x.
float _nu_pz
Y momentum of neutrino.
BEGIN_PROLOG SN nu
float _nu_e
PDG of neutrino.
float _nu_other_z
Y poisition of neutrino at the front face of the TPC (other)
TVector3 FluxReaderAna::GetIntersection ( TVector3  nu_pos,
TVector3  nu_dir,
float  z_location = 0 
)
private

Returns the intersection point with the front face of the TPC.

Definition at line 235 of file FluxReaderAna_module.cc.

235  {
236 
237  TVector3 plane_point(0, 0, z_location);
238  TVector3 plane_normal(0, 0, 1);
239 
240  TVector3 diff = nu_pos - plane_point;
241  double prod1 = diff.Dot(plane_normal);
242  double prod2 = nu_dir.Dot(plane_normal);
243  double prod3 = prod1 / prod2;
244  return nu_pos - nu_dir * prod3;
245 
246 }
FluxReaderAna& FluxReaderAna::operator= ( FluxReaderAna const &  )
delete
FluxReaderAna& FluxReaderAna::operator= ( FluxReaderAna &&  )
delete

Member Data Documentation

bool FluxReaderAna::_apply_position_cuts
private

Definition at line 64 of file FluxReaderAna_module.cc.

float FluxReaderAna::_baseline
private

Definition at line 60 of file FluxReaderAna_module.cc.

std::string FluxReaderAna::_flux_label
private

Label for flux dataproduct (to be set via fcl)

Definition at line 58 of file FluxReaderAna_module.cc.

int FluxReaderAna::_nu_decay
private

Z momentum of neutrino.

Definition at line 89 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_dk2gen
private

Neutrino parent decay code.

Definition at line 90 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_e
private

PDG of neutrino.

Definition at line 77 of file FluxReaderAna_module.cc.

bool FluxReaderAna::_nu_hit
private

Definition at line 75 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_intersection_z
private

Z position where to evaluate the flux (to be set via fcl)

Definition at line 61 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_oaa
private

Neutrino r.

Definition at line 98 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_other_intersection_z
private

Additional Z position where to evaluate the flux (to be set via fcl)

Definition at line 62 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_other_t
private

Z poisition of neutrino at the front face of the TPC (other)

Definition at line 103 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_other_x
private

Neutrino off axis angle.

Definition at line 100 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_other_y
private

X poisition of neutrino at the front face of the TPC (other)

Definition at line 101 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_other_z
private

Y poisition of neutrino at the front face of the TPC (other)

Definition at line 102 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_p_angle
private

Distance from decay to ray origin.

Definition at line 91 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_p_dpx
private

Neutrino parent pdg.

Definition at line 93 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_p_dpy
private

Neutrino parent momentum x.

Definition at line 94 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_p_dpz
private

Neutrino parent momentum x.

Definition at line 95 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_p_e
private

Neutrino parent momentum x.

Definition at line 96 of file FluxReaderAna_module.cc.

int FluxReaderAna::_nu_p_type
private

Angle between neutrino and parent direction.

Definition at line 92 of file FluxReaderAna_module.cc.

int FluxReaderAna::_nu_pdg
private

True if the neutrino hit the requested volumes.

Definition at line 76 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_px
private

Z poisition of neutrino at neutrino production point.

Definition at line 86 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_py
private

X momentum of neutrino.

Definition at line 87 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_pz
private

Y momentum of neutrino.

Definition at line 88 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_r
private

Neutrino parent energy.

Definition at line 97 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_t
private

Energy of the neutrino.

Definition at line 78 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_vx
private

Z poisition of neutrino at the front face of the TPC.

Definition at line 83 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_vy
private

X poisition of neutrino at neutrino production point.

Definition at line 84 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_vz
private

Y poisition of neutrino at neutrino production point.

Definition at line 85 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_w
private

Time of the neutrino.

Definition at line 79 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_x
private

Neutrino weight.

Definition at line 80 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_y
private

X poisition of neutrino at the front face of the TPC.

Definition at line 81 of file FluxReaderAna_module.cc.

float FluxReaderAna::_nu_z
private

Y poisition of neutrino at the front face of the TPC.

Definition at line 82 of file FluxReaderAna_module.cc.

const TDatabasePDG* FluxReaderAna::_pdg_database = TDatabasePDG::Instance()
private

Definition at line 69 of file FluxReaderAna_module.cc.

TTree* FluxReaderAna::_tree
private

Definition at line 74 of file FluxReaderAna_module.cc.

float FluxReaderAna::_x_cut
private

Definition at line 65 of file FluxReaderAna_module.cc.

float FluxReaderAna::_x_shift
private

Definition at line 59 of file FluxReaderAna_module.cc.

float FluxReaderAna::_y_cut
private

Definition at line 66 of file FluxReaderAna_module.cc.


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