1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
|
#include <TVector3.h>
#include <vector>
#include <TTree.h>
#include <string>
struct Photon {
double t;
TVector3 pos;
TVector3 dir;
TVector3 pol;
double wavelength; // nm
unsigned int history;
int last_hit_triangle;
};
struct Track {
std::string particle;
double t;
TVector3 pos;
TVector3 dir;
double start_time;
double total_energy;
};
struct MC {
std::string particle;
TVector3 gen_pos;
TVector3 gen_dir;
double gen_total_e;
std::vector<Track> subtrack;
std::vector<Photon> photon_start;
std::vector<Photon> photon_stop;
};
struct Channel {
Channel() : channel_id(-1), t(-9999.0), q(-9999.0) { };
int channel_id;
double t;
double q;
unsigned int mc_history;
};
struct Event {
int event_id;
MC mc;
int nhit;
std::vector<Channel> channel;
// Populate arrays of length nentries with hit, time, and charge
// information, indexed by channel ID
void get_channels(unsigned int nentries, int *hit, float *time,
float *charge)
{
for (unsigned int i=0; i < nentries; i++) {
hit[i] = 0;
time[i] = -1e9f;
charge[i] = -1e9f;
}
for (unsigned int i=0; i < channel.size(); i++) {
unsigned int channel_id = channel[i].channel_id;
if (channel_id < nentries) {
hit[channel_id] = 1;
time[channel_id] = channel[i].t;
charge[channel_id] = channel[i].q;
}
}
}
};
void fill_photons(Event *ev, bool start,
unsigned int nphotons, float *pos, float *dir,
float *pol, float *wavelength, float *t0,
int *histories=0, int *last_hit_triangle=0)
{
std::vector<Photon> &photons = start ? ev->mc.photon_start : ev->mc.photon_stop;
photons.resize(nphotons);
for (unsigned int i=0; i < nphotons; i++) {
Photon &photon = photons[i];
photon.t = t0[i];
photon.pos.SetXYZ(pos[3*i], pos[3*i + 1], pos[3*i + 2]);
photon.dir.SetXYZ(dir[3*i], dir[3*i + 1], dir[3*i + 2]);
photon.pol.SetXYZ(pol[3*i], pol[3*i + 1], pol[3*i + 2]);
photon.wavelength = wavelength[i];
if (histories)
photon.history = histories[i];
else
photon.history = 0;
if (last_hit_triangle)
photon.last_hit_triangle = last_hit_triangle[i];
else
photon.last_hit_triangle = -1;
}
}
void fill_hits(Event *ev, unsigned int nchannels, float *t,
float *q, unsigned int *history)
{
ev->channel.resize(0);
ev->nhit = 0;
Channel ch;
for (unsigned int i=0; i < nchannels; i++) {
if (t[i] < 1e8) {
ev->nhit++;
ch.channel_id = i;
ch.t = t[i];
ch.q = q[i];
ch.mc_history = history[i];
ev->channel.push_back(ch);
}
}
}
#ifdef __MAKECINT__
#pragma link C++ class vector<Track>;
#pragma link C++ class vector<Photon>;
#pragma link C++ class vector<Channel>;
#endif
|
