All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Public Attributes | Private Member Functions | List of all members
ana::SBNOsc::InteractionHistos Struct Reference

#include <InteractionHisto.h>

Inheritance diagram for ana::SBNOsc::InteractionHistos:
ana::SBNOsc::HistoList

Public Member Functions

void Initialize (const std::string &prefix, const geo::BoxBoundedGeo &detector_volume, const std::vector< double > &tagger_volume)
 
void Fill (const numu::RecoInteraction &vertex, const numu::RecoEvent &event, const std::vector< event::Interaction > &core_truth)
 
void Fill (const event::Interaction &interaction, unsigned mctruth_id, const numu::RecoEvent &event)
 
- Public Member Functions inherited from ana::SBNOsc::HistoList
void Scale (double scale)
 
void Add (const HistoList &other)
 
void Write ()
 
void StoreHisto (TH1 *histo)
 
void Merge (const HistoList &merge)
 

Public Attributes

TH1D * track_length
 Length of the reconstructed primary track. More...
 
TH1D * nuE
 Neutrino energy. More...
 
TH1D * track_p
 Primary track momentum. More...
 
TH1D * beam_center_distance
 Distance of the neutrino interaction to the beam center. More...
 
TH1D * Q2
 Q2 of the interaction. More...
 
TH1D * true_contained_length
 True contained length of primary track. More...
 
TH1D * true_track_multiplicity
 True particle multiplicity of the interaction. More...
 
TH1D * crosses_tpc
 Whether the primary track crosses a TPC boundary. More...
 
TH1D * dist_to_match
 Distance from this vertex to the closest matching vertex reco->truth and truth->reco. More...
 
TH1D * primary_track_completion
 Completion of the primary track. More...
 
TH1D * n_reco_vertices
 Number of reconstructed vertices in the event with this vertex. More...
 
TH1D * maxpe_crt_intime_hit
 Maximum number of PE's in a single CRT hit in time with the beam. More...
 
TH1D * crt_hit_times
 
TH1D * closest_crt_hit_time
 
TH1D * crt_pes
 
TH1D * fmatch_score
 
TH2D * fmatch_score_true_time
 
TH2D * fmatch_score_true_time_zoom
 
TH1D * fmatch_score_outtime
 
TH1D * fmatch_score_intime
 
TH2D * fmatch_time_true_time_zoom
 
TH1D * fmatch_time
 
TH1D * fmatch_time_real_time
 
TH2D * intime_crt_hits_xy
 
TH2D * intime_crt_hits_xz
 
TH2D * intime_crt_hits_yz
 
TH2D * vertex_xy
 
TH2D * vertex_yz
 
TH2D * vertex_xz
 
TH1D * vertex_x
 
TH1D * vertex_y
 
TH1D * vertex_z
 
TH2D * light_trigger
 
- Public Attributes inherited from ana::SBNOsc::HistoList
std::vector< TH1 * > fAllHistos
 
std::vector< TDirectory * > fLocations
 

Private Member Functions

void FillEvent (const numu::RecoEvent &event)
 

Detailed Description

Histograms associated with neutrino interactions. Filled for the list of all true and reco vertices. These histograms are constructed per interaction mode per cut.

Definition at line 25 of file InteractionHisto.h.

Member Function Documentation

void ana::SBNOsc::InteractionHistos::Fill ( const numu::RecoInteraction vertex,
const numu::RecoEvent event,
const std::vector< event::Interaction > &  core_truth 
)

Fill the histograms with a single reconstructed neutrino candidate

Parameters
vertexThe reconstructed vertex being filled
eventThe reco event object
core_truthThe list of true interactions from sbncode core

Definition at line 90 of file InteractionHisto.cc.

93  {
94 
95  const numu::RecoTrack &primary_track = event.tracks.at(vertex.slice.primary_track_index);
96 
98  primary_track_completion->Fill(primary_track.match.completion);
99 
100  vertex_xy->Fill(vertex.position.X(), vertex.position.Y());
101  vertex_yz->Fill(vertex.position.Y(), vertex.position.Z());
102  vertex_xz->Fill(vertex.position.X(), vertex.position.Z());
103 
104  vertex_x->Fill(vertex.position.X());
105  vertex_y->Fill(vertex.position.Y());
106  vertex_z->Fill(vertex.position.Z());
107 
108 
109  if (vertex.slice.flash_match.present) {
110  fmatch_score->Fill(vertex.slice.flash_match.score);
111  fmatch_time->Fill(vertex.slice.flash_match.time);
112  if (vertex.slice.primary_track_index >= 0 && event.tracks.at(vertex.slice.primary_track_index).match.has_match) {
113  int mcparticle_id = event.tracks.at(vertex.slice.primary_track_index).match.mcparticle_id;
114  double true_time = event.particles.at(mcparticle_id).start_time;
115  fmatch_score_true_time->Fill(vertex.slice.flash_match.score, true_time);
116  fmatch_score_true_time_zoom->Fill(vertex.slice.flash_match.score, true_time);
117  fmatch_time_true_time_zoom->Fill(vertex.slice.flash_match.time, true_time);
118  if (true_time < 0. || true_time > 1.6) {
120  }
121  else {
123  }
124 
125  fmatch_time_real_time->Fill(vertex.slice.flash_match.time - true_time);
126  }
127  }
128 
129  track_length->Fill(primary_track.length);
130  if (primary_track.match.has_match) {
131  int mcparticle_id = primary_track.match.mcparticle_id;
132 
133  double true_track_momentum = event.particles.at(mcparticle_id).start_momentum.Mag();
134  track_p->Fill(true_track_momentum);
135 
136  int crosses_tpc_val = event.particles.at(mcparticle_id).crosses_tpc;
137  crosses_tpc->Fill(crosses_tpc_val);
138 
139  double length = event.particles.at(mcparticle_id).length;
140  true_contained_length->Fill(length);
141  }
142 
143  if (vertex.match.mctruth_track_id >= 0) {
144  int mctruth_id = vertex.match.mctruth_track_id;
145 
146  true_track_multiplicity->Fill(truth[mctruth_id].nfinalstate);
147  nuE->Fill(truth[mctruth_id].neutrino.energy);
148  Q2->Fill(truth[mctruth_id].neutrino.Q2);
149  // get the distance from the beam center
150  /*
151  float beam_center_distance = sqrt( (truth[mctruth_id].neutrino.position.X() - _config.beamCenterX) *
152  (truth[mctruth_id].neutrino.position.X() - _config.beamCenterX) +
153  (truth[mctruth_id].neutrino.position.Y() - _config.beamCenterY) *
154  (truth[mctruth_id].neutrino.position.Y() - _config.beamCenterY));
155 
156  beam_center_distance->Fill(beam_center_distance);
157  */
158  }
159 
160  // fill the event-info based histograms
161  FillEvent(event);
162 }
TH1D * true_contained_length
True contained length of primary track.
void FillEvent(const numu::RecoEvent &event)
TH1D * track_p
Primary track momentum.
TH1D * dist_to_match
Distance from this vertex to the closest matching vertex reco-&gt;truth and truth-&gt;reco.
TruthMatch match
Info for mathing to truth.
Definition: RecoEvent.h:42
TH1D * nuE
Neutrino energy.
int mcparticle_id
MCParticle ID of the particle this track matches to (same as the ID of the RecoTrack of that particle...
int primary_track_index
Index of the primary track.
Definition: RecoEvent.h:27
TVector3 position
location of the vertex
Definition: RecoEvent.h:40
std::map< size_t, RecoTrack > tracks
Map of track indices to Track information.
Definition: RecoEvent.h:48
FlashMatch flash_match
Result of flash matching algorithm on this slice.
Definition: RecoEvent.h:30
TH1D * Q2
Q2 of the interaction.
float length
Length of track.
Definition: RecoTrack.h:48
float completion
Fraction of energy deposits by true particle matched by this track.
bool has_match
Whether a track match exists.
TrackTruthMatch match
Truth matching information.
Definition: RecoTrack.h:57
RecoSlice slice
Particle content of the interaction.
Definition: RecoEvent.h:39
float truth_vertex_distance
Distance from truth interaction vertex to this vertex.
TH1D * track_length
Length of the reconstructed primary track.
int mctruth_track_id
index of the primary track into the list of MCTruths. -1 if no match.
TH1D * crosses_tpc
Whether the primary track crosses a TPC boundary.
TH1D * primary_track_completion
Completion of the primary track.
TH1D * true_track_multiplicity
True particle multiplicity of the interaction.
void ana::SBNOsc::InteractionHistos::Fill ( const event::Interaction interaction,
unsigned  mctruth_id,
const numu::RecoEvent event 
)

Fill the histograms with a single true interaction

Parameters
interactionThe true interaction
mctruth_idID of the interaction (index into the list of Interactions)
eventThe reco event object

Definition at line 60 of file InteractionHisto.cc.

60  {
61  double dist = numu::dist2Match(interaction, event.reco);
62  dist_to_match->Fill(dist);
63  primary_track_completion->Fill(numu::trackMatchCompletion(mctruth_id, event));
64 
65  vertex_xy->Fill(interaction.neutrino.position.X(), interaction.neutrino.position.Y());
66  vertex_yz->Fill(interaction.neutrino.position.Y(), interaction.neutrino.position.Z());
67  vertex_xz->Fill(interaction.neutrino.position.X(), interaction.neutrino.position.Z());
68 
69  vertex_x->Fill(interaction.neutrino.position.X());
70  vertex_y->Fill(interaction.neutrino.position.Y());
71  vertex_z->Fill(interaction.neutrino.position.Z());
72 
73 
74  // TODO: how to fix these histograms?
75  // crosses_tpc
76 
77  // only fill for muon tracks (CC interactions)
78  if (abs(interaction.lepton.pdg) == 13) {
79  track_p->Fill(interaction.lepton.momentum.Mag());
81  }
82 
83  true_track_multiplicity->Fill(interaction.nfinalstate);
84  nuE->Fill(interaction.neutrino.energy);
85  Q2->Fill(interaction.neutrino.Q2);
86 
87  FillEvent(event);
88 }
TVector3 momentum
Three-momentum.
Definition: Event.hh:129
TH1D * true_contained_length
True contained length of primary track.
void FillEvent(const numu::RecoEvent &event)
int pdg
PDG Code.
Definition: Event.hh:127
FinalStateParticle lepton
The primary final state lepton.
Definition: Event.hh:153
TH1D * track_p
Primary track momentum.
TH1D * dist_to_match
Distance from this vertex to the closest matching vertex reco-&gt;truth and truth-&gt;reco.
TH1D * nuE
Neutrino energy.
Neutrino neutrino
The neutrino.
Definition: Event.hh:152
float Q2
Q squared.
Definition: Event.hh:82
std::vector< RecoInteraction > reco
List of reconstructed vertices.
Definition: RecoEvent.h:50
T abs(T value)
TH1D * Q2
Q2 of the interaction.
float dist2Match(const event::Interaction &truth, const std::vector< numu::RecoInteraction > &candidates)
Definition: Derived.cc:3
size_t nfinalstate
Size of finalstate.
Definition: Event.hh:155
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
float energy
Neutrino energy (GeV)
Definition: Event.hh:90
TVector3 position
Neutrino interaction position.
Definition: Event.hh:92
float trackMatchCompletion(unsigned truth_index, const numu::RecoEvent &event)
Definition: Derived.cc:13
TH1D * primary_track_completion
Completion of the primary track.
TH1D * true_track_multiplicity
True particle multiplicity of the interaction.
void ana::SBNOsc::InteractionHistos::FillEvent ( const numu::RecoEvent event)
private

Definition at line 164 of file InteractionHisto.cc.

164  {
165  n_reco_vertices->Fill(event.reco.size());
166 
167  double maxpe = 0.;
168  double closest_time_dist = -1;
169  double closest_time = 0.;
170  for (const numu::CRTHit &hit: event.in_time_crt_hits) {
171  intime_crt_hits_xy->Fill(hit.location.X(), hit.location.Y());
172  intime_crt_hits_xz->Fill(hit.location.X(), hit.location.Z());
173  intime_crt_hits_yz->Fill(hit.location.Y(), hit.location.Z());
174 
175  crt_pes->Fill(hit.pes);
176  if (hit.pes > maxpe) maxpe = hit.pes;
177 
178  //if (hit.pes < 100.) continue;
179 
180  crt_hit_times->Fill(hit.time);
181 
182  if (closest_time_dist < 0. || closest_time_dist > 1e-3) {
183  double this_time_dist = -1;
184  if (hit.time > -0.2 && hit.time < 1.8) this_time_dist = 0.;
185  else if (hit.time < 0.) this_time_dist = -hit.time -0.2;
186  if (this_time_dist >= 0. && (this_time_dist < closest_time_dist || closest_time_dist < 0.)) {
187  closest_time = hit.time;
188  closest_time_dist = this_time_dist;
189  }
190  }
191  }
192  if (closest_time_dist >= 0.) {
193  closest_crt_hit_time->Fill(closest_time);
194  }
195  else {
196  closest_crt_hit_time->Fill(30.);
197  }
198 
199  std::cout << "Max PE: " << maxpe << std::endl;
200  if (event.in_time_crt_hits.size() == 0) {
201  maxpe_crt_intime_hit->Fill(-1);
202  }
203  else {
204  maxpe_crt_intime_hit->Fill(maxpe);
205  }
206 }
TH1D * n_reco_vertices
Number of reconstructed vertices in the event with this vertex.
process_name hit
Definition: cheaterreco.fcl:51
std::vector< CRTHit > in_time_crt_hits
List of crt hits in time with the beam spill.
Definition: RecoEvent.h:51
TVector3 location
Location of the hit.
Definition: DetInfo.h:16
std::vector< RecoInteraction > reco
List of reconstructed vertices.
Definition: RecoEvent.h:50
TH1D * maxpe_crt_intime_hit
Maximum number of PE&#39;s in a single CRT hit in time with the beam.
float time
CRT Hit time.
Definition: DetInfo.h:13
float pes
Number of PE&#39;s in hit.
Definition: DetInfo.h:15
do i e
BEGIN_PROLOG could also be cout
void ana::SBNOsc::InteractionHistos::Initialize ( const std::string &  prefix,
const geo::BoxBoundedGeo detector_volume,
const std::vector< double > &  tagger_volume 
)

Intialize the histograms

Parameters
prefixA prefix to be added to each histogram name
modeThe mode of interaction for these histograms
indexThe cut index for these histograms.

Definition at line 11 of file InteractionHisto.cc.

11  {
12 #define INT_HISTO(name, n_bins, lo, hi) name = new TH1D((#name + postfix).c_str(), #name, n_bins, lo, hi); StoreHisto(name)
13 #define INT_HISTO2D(name, n_binsx, xlo, xhi, n_binsy, ylo, yhi) name = new TH2D((#name + postfix).c_str(), #name, n_binsx, xlo, xhi, n_binsy, ylo, yhi); StoreHisto(name)
14 #define INT_HISTO2D_BINSY(name, n_binsx, xlo, xhi, n_binsy, binsy) name = new TH2D((#name + postfix).c_str(), #name, n_binsx, xlo, xhi, n_binsy, binsy); StoreHisto(name)
15 
16  INT_HISTO(track_length, 100, 0., 600.);
17  INT_HISTO(track_p, 50, 0., 5.);
18  INT_HISTO(nuE, 50, 0., 5.);
19  INT_HISTO(beam_center_distance, 60, 0., 300.);
20  INT_HISTO(Q2, 50, 0., 10.);
21  INT_HISTO(true_contained_length, 100, 0., 600.);
23  INT_HISTO(crosses_tpc, 2, -0.5, 1.5);
24  INT_HISTO(dist_to_match, 101, -1., 100.);
25  INT_HISTO(primary_track_completion, 100, 0., 1.);
26  INT_HISTO(n_reco_vertices, 10, -0.5, 9.5);
27  INT_HISTO(maxpe_crt_intime_hit, 1000, 0., 1000.);
28  INT_HISTO(crt_pes, 1000, 0., 1000.);
29 
30  INT_HISTO(crt_hit_times, 300, -20., 10.);
31  INT_HISTO(closest_crt_hit_time, 300, -20., 10.);
32 
33  INT_HISTO(fmatch_score, 200, 0., 200.);
34  INT_HISTO2D(fmatch_score_true_time, 100, 0., 100., 1400, -4000., 3000.);
35  INT_HISTO2D(fmatch_score_true_time_zoom, 100, 0., 100., 100, -5., 5.);
36  double intout_times[4] = {-3000., 0., 1.6, 3000.};
37  INT_HISTO(fmatch_score_intime, 100, 0., 100.);
38  INT_HISTO(fmatch_score_outtime, 100, 0., 100.);
39  INT_HISTO2D(fmatch_time_true_time_zoom, 100, -5., 5., 100, -5., 5.);
40  INT_HISTO(fmatch_time_real_time, 4000, -2, 2);
41  INT_HISTO(fmatch_time, 700, -5., 2.);
42 
43  // INT_HISTO2D(light_trigger, 20, 0.5, 20.5, 200, 6000, 8000);
44 
45  INT_HISTO2D(intime_crt_hits_xy, 100, tagger_volume[0], tagger_volume[3], 100, tagger_volume[1], tagger_volume[4]);
46  INT_HISTO2D(intime_crt_hits_xz, 100, tagger_volume[0], tagger_volume[3], 100, tagger_volume[2], tagger_volume[5]);
47  INT_HISTO2D(intime_crt_hits_yz, 100, tagger_volume[1], tagger_volume[4], 100, tagger_volume[2], tagger_volume[5]);
48 
49  INT_HISTO2D(vertex_xy, 100, detector_volume.MinX(), detector_volume.MaxX(), 100, detector_volume.MinY(), detector_volume.MaxY());
50  INT_HISTO2D(vertex_yz, 100, detector_volume.MinY(), detector_volume.MaxY(), 100, detector_volume.MinZ(), detector_volume.MaxZ());
51  INT_HISTO2D(vertex_xz, 100, detector_volume.MinX(), detector_volume.MaxX(), 100, detector_volume.MinZ(), detector_volume.MaxZ());
52 
53  INT_HISTO(vertex_x, 100, detector_volume.MinX(), detector_volume.MaxX());
54  INT_HISTO(vertex_y, 100, detector_volume.MinY(), detector_volume.MaxY());
55  INT_HISTO(vertex_z, 100, detector_volume.MinZ(), detector_volume.MaxZ());
56 
57 #undef INT_HISTO
58 }
TH1D * true_contained_length
True contained length of primary track.
TH1D * track_p
Primary track momentum.
TH1D * dist_to_match
Distance from this vertex to the closest matching vertex reco-&gt;truth and truth-&gt;reco.
TH1D * n_reco_vertices
Number of reconstructed vertices in the event with this vertex.
TH1D * nuE
Neutrino energy.
double MinX() const
Returns the world x coordinate of the start of the box.
Definition: BoxBoundedGeo.h:88
double MaxX() const
Returns the world x coordinate of the end of the box.
Definition: BoxBoundedGeo.h:91
TH1D * beam_center_distance
Distance of the neutrino interaction to the beam center.
TH1D * maxpe_crt_intime_hit
Maximum number of PE&#39;s in a single CRT hit in time with the beam.
TH1D * Q2
Q2 of the interaction.
double MinZ() const
Returns the world z coordinate of the start of the box.
double MaxY() const
Returns the world y coordinate of the end of the box.
double MaxZ() const
Returns the world z coordinate of the end of the box.
#define INT_HISTO(name, n_bins, lo, hi)
TH1D * track_length
Length of the reconstructed primary track.
TH1D * crosses_tpc
Whether the primary track crosses a TPC boundary.
TH1D * primary_track_completion
Completion of the primary track.
double MinY() const
Returns the world y coordinate of the start of the box.
TH1D * true_track_multiplicity
True particle multiplicity of the interaction.
#define INT_HISTO2D(name, n_binsx, xlo, xhi, n_binsy, ylo, yhi)

Member Data Documentation

TH1D* ana::SBNOsc::InteractionHistos::beam_center_distance

Distance of the neutrino interaction to the beam center.

Definition at line 30 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::closest_crt_hit_time

Definition at line 40 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::crosses_tpc

Whether the primary track crosses a TPC boundary.

Definition at line 34 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::crt_hit_times

Definition at line 39 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::crt_pes

Definition at line 41 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::dist_to_match

Distance from this vertex to the closest matching vertex reco->truth and truth->reco.

Definition at line 35 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::fmatch_score

Definition at line 42 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::fmatch_score_intime

Definition at line 46 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::fmatch_score_outtime

Definition at line 45 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::fmatch_score_true_time

Definition at line 43 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::fmatch_score_true_time_zoom

Definition at line 44 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::fmatch_time

Definition at line 48 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::fmatch_time_real_time

Definition at line 49 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::fmatch_time_true_time_zoom

Definition at line 47 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::intime_crt_hits_xy

Definition at line 50 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::intime_crt_hits_xz

Definition at line 51 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::intime_crt_hits_yz

Definition at line 52 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::light_trigger

Definition at line 59 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::maxpe_crt_intime_hit

Maximum number of PE's in a single CRT hit in time with the beam.

Definition at line 38 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::n_reco_vertices

Number of reconstructed vertices in the event with this vertex.

Definition at line 37 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::nuE

Neutrino energy.

Definition at line 28 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::primary_track_completion

Completion of the primary track.

Definition at line 36 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::Q2

Q2 of the interaction.

Definition at line 31 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::track_length

Length of the reconstructed primary track.

Definition at line 27 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::track_p

Primary track momentum.

Definition at line 29 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::true_contained_length

True contained length of primary track.

Definition at line 32 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::true_track_multiplicity

True particle multiplicity of the interaction.

Definition at line 33 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::vertex_x

Definition at line 56 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::vertex_xy

Definition at line 53 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::vertex_xz

Definition at line 55 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::vertex_y

Definition at line 57 of file InteractionHisto.h.

TH2D* ana::SBNOsc::InteractionHistos::vertex_yz

Definition at line 54 of file InteractionHisto.h.

TH1D* ana::SBNOsc::InteractionHistos::vertex_z

Definition at line 58 of file InteractionHisto.h.


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