Skip to content

Commit 3a2fd4f

Browse files
maureeungaromaureeungaro
committed
Added vertex information to true and digitized infos
1 parent b4c01db commit 3a2fd4f

12 files changed

Lines changed: 264 additions & 67 deletions

File tree

‎actions/generator/gPrimaryGeneratorAction.cc‎

Lines changed: 24 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -28,6 +28,28 @@ GParticleRecord make_particle_record(const GparticlePtr& particle) {
2828
vertex.z()
2929
};
3030
}
31+
32+
GParticleRecord make_particle_record(const GparticleRuntimeRecord& particle) {
33+
return {
34+
particle.name,
35+
particle.pid,
36+
particle.type,
37+
1,
38+
particle.p,
39+
particle.theta,
40+
particle.phi,
41+
particle.vertex.x(),
42+
particle.vertex.y(),
43+
particle.vertex.z()
44+
};
45+
}
46+
47+
void append_runtime_records(GParticleRecordEvent& records, const GparticlePtr& particle) {
48+
if (particle == nullptr) { return; }
49+
for (const auto& runtime_record : particle->getRuntimeRecords()) {
50+
records.emplace_back(make_particle_record(runtime_record));
51+
}
52+
}
3153
}
3254

3355
// Build the primary-generator action, load the configured particle definitions,
@@ -62,10 +84,8 @@ void GPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) {
6284
current_generated_particles.insert(current_generated_particles.end(), gparticles.begin(), gparticles.end());
6385
current_generated_tracked_particles.insert(current_generated_tracked_particles.end(), gparticles.begin(), gparticles.end());
6486
current_generated_particle_records.reserve(gparticles.size());
65-
current_generated_tracked_particle_records.reserve(gparticles.size());
6687
for (const auto& gparticle : gparticles) {
6788
current_generated_particle_records.emplace_back(make_particle_record(gparticle));
68-
current_generated_tracked_particle_records.emplace_back(make_particle_record(gparticle));
6989
}
7090

7191
const auto event_id = anEvent->GetEventID();
@@ -80,15 +100,13 @@ void GPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) {
80100
current_generated_tracked_particles.insert(current_generated_tracked_particles.end(),
81101
event_particles.begin(),
82102
event_particles.end());
83-
for (const auto& gparticle : event_particles) {
84-
current_generated_tracked_particle_records.emplace_back(make_particle_record(gparticle));
85-
}
86103
}
87104

88105
for (const auto& gparticle : gparticles) {
89106

90107
if (gparticle != nullptr) {
91108
gparticle->shootParticle(gparticleGun.get(), anEvent);
109+
append_runtime_records(current_generated_tracked_particle_records, gparticle);
92110
}
93111
}
94112

@@ -100,6 +118,7 @@ void GPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) {
100118
for (const auto& gparticle : gparticleFileEvents[static_cast<size_t>(event_id)]) {
101119
if (gparticle != nullptr) {
102120
gparticle->shootParticle(gparticleGun.get(), anEvent);
121+
append_runtime_records(current_generated_tracked_particle_records, gparticle);
103122
}
104123
}
105124
}

‎analyzer/README.md‎

Lines changed: 42 additions & 48 deletions
Original file line numberDiff line numberDiff line change
@@ -66,14 +66,20 @@ hitn, pid, tid, E, time, totEdep
6666
The true-info output includes tracking columns like:
6767

6868
```text
69-
processName, avgTime, avgx, avgy, avgz, hitn, pid, tid, totalEDeposited
69+
processName, avgTime, avgx, avgy, avgz, hitn, pid, tid, mtid, vx, vy, vz, mvx, mvy, mvz, totalEDeposited
7070
```
7171

7272
When the matching digitized CSV is available, the analyzer also adds `E` to
7373
true-info tables by matching rows on event, detector, hit, PID, and track ID.
7474
In that case `E` is the track total energy, while `totalEDeposited` remains
7575
the deposited energy.
7676

77+
The `vx`, `vy`, and `vz` columns are the current track vertex coordinates.
78+
The `mvx`, `mvy`, and `mvz` columns are the mother-track vertex coordinates
79+
when the mother track was available to GEMC hit processing; otherwise they use
80+
the GEMC uninitialized numeric sentinel. The `mtid` column stores the mother
81+
track id.
82+
7783
The ROOT streamer writes one ROOT file per worker thread. For one thread and
7884
`filename: b2`, the file is typically:
7985

@@ -100,7 +106,7 @@ Read one digitized CSV file:
100106
```python
101107
from analyzer import read_output
102108

103-
output = read_output("tmp/b2_t0_digitized.csv")
109+
output = read_output("b2_t0_digitized.csv")
104110
df = output.get_frame("digitized")
105111
print(df.columns)
106112
```
@@ -110,7 +116,7 @@ Read a CSV root name when both files exist:
110116
```python
111117
from analyzer import read_output
112118

113-
output = read_output("tmp/b2_t0", kind="csv")
119+
output = read_output("b2_t0", kind="csv")
114120
print(output.summary())
115121
```
116122

@@ -119,7 +125,7 @@ Plot `totEdep` grouped by `pid`:
119125
```python
120126
from analyzer import plot_variable, read_output
121127

122-
output = read_output("tmp/b2_t0_digitized.csv")
128+
output = read_output("b2_t0_digitized.csv")
123129
plot_variable(
124130
output,
125131
"totEdep",
@@ -135,7 +141,7 @@ Read ROOT output:
135141
```python
136142
from analyzer import read_output
137143

138-
output = read_output("tmp/b2_t0.root", kind="root")
144+
output = read_output("b2_t0.root", kind="root")
139145
df = output.get_frame("digitized", detector="flux")
140146
```
141147

@@ -145,31 +151,37 @@ Run `python3 -m analyzer` from the GEMC source directory, where the `analyzer`
145151
package directory is visible to Python.
146152

147153
The `-m` flag takes a module name, not a filesystem path. Do not run
148-
`python3 -m ../analyzer`. If your shell is inside `tmp/`, either move back to
149-
the source directory or set `PYTHONPATH=..`.
154+
`python3 -m ../analyzer`. If your shell is in another directory, move back to
155+
the source directory or set `PYTHONPATH`.
150156

151157
Print a summary:
152158

153159
```sh
154-
python3 -m analyzer tmp/b2_t0_digitized.csv
160+
python3 -m analyzer b2_t0_digitized.csv
155161
```
156162

157163
Plot a digitized variable with matplotlib:
158164

159165
```sh
160-
python3 -m analyzer tmp/b2_t0_digitized.csv totEdep --kind csv --xlim 0.0 0.1
166+
python3 -m analyzer b2_t0_digitized.csv totEdep --kind csv --xlim 0.0 0.1
161167
```
162168

163169
Save a plot instead of showing it:
164170

165171
```sh
166-
python3 -m analyzer tmp/b2_t0_digitized.csv totEdep --kind csv --save tmp/b2_totEdep.png
172+
python3 -m analyzer b2_t0_digitized.csv totEdep --kind csv --save b2_totEdep.png
167173
```
168174

169175
Plot ROOT output with matplotlib:
170176

171177
```sh
172-
python3 -m analyzer tmp/b2_t0.root totEdep --kind root --detector flux --save tmp/b2_totEdep.png
178+
python3 -m analyzer b2_t0.root totEdep --kind root --detector flux --save b2_totEdep.png
179+
```
180+
181+
Plot a true-info track vertex coordinate:
182+
183+
```sh
184+
python3 -m analyzer b2_t0_true_info.csv vx --kind csv --data true_info --save b2_vertex_x.png
173185
```
174186

175187
## Dependency-Free SVG Plot
@@ -178,51 +190,47 @@ If `pandas`, `numpy`, or `matplotlib` are unavailable, create an SVG histogram
178190
directly from the CSV file:
179191

180192
```sh
181-
python3 -B analyzer/svg_plot.py tmp/b2_t0_digitized.csv totEdep --out tmp/b2_totEdep.svg --bins 30
193+
python3 -B analyzer/svg_plot.py b2_t0_digitized.csv totEdep --out b2_totEdep.svg --bins 30
182194
```
183195

184196
Add an x-axis range with:
185197

186198
```sh
187-
python3 -B analyzer/svg_plot.py tmp/b2_t0_digitized.csv totEdep --out tmp/b2_totEdep.svg --bins 30 --xlim 0.0 0.1
199+
python3 -B analyzer/svg_plot.py b2_t0_digitized.csv totEdep --out b2_totEdep.svg --bins 30 --xlim 0.0 0.1
188200
```
189201

190202
## Run the B2 Example
191203

192-
Run these commands from the GEMC source directory:
193-
194-
```sh
195-
mkdir -p tmp
196-
```
204+
Run these commands from the GEMC source directory.
197205

198206
Build the B2 geometry into a local SQLite database:
199207

200208
```sh
201209
PYTHONDONTWRITEBYTECODE=1 PYTHONPATH=/opt/projects/gemc/src/api \
202-
python3 examples/basic/b2/b2.py -f sqlite -sql tmp/gemc.db
210+
python3 examples/basic/b2/b2.py -f sqlite -sql gemc.db
203211
```
204212

205-
Run GEMC with CSV output rooted at `tmp/b2`:
213+
Run GEMC with CSV output rooted at `b2`:
206214

207215
```sh
208216
build/gemc examples/basic/b2/b2.yaml \
209217
'-gsystem=[{name: b2, factory: sqlite}]' \
210-
'-gstreamer=[{format: csv, filename: tmp/b2}]' \
211-
-sql=tmp/gemc.db \
218+
'-gstreamer=[{format: csv, filename: b2}]' \
219+
-sql=gemc.db \
212220
-n=20
213221
```
214222

215223
With one worker thread, this produces:
216224

217225
```text
218-
tmp/b2_t0_digitized.csv
219-
tmp/b2_t0_true_info.csv
226+
b2_t0_digitized.csv
227+
b2_t0_true_info.csv
220228
```
221229

222230
Inspect the digitized CSV header:
223231

224232
```sh
225-
head -1 tmp/b2_t0_digitized.csv
233+
head -1 b2_t0_digitized.csv
226234
```
227235

228236
Expected columns include:
@@ -234,31 +242,31 @@ evn, timestamp, thread_id, detector, hitn, pid, tid, E, time, totEdep
234242
Create the `totEdep` plot with the main analyzer API:
235243

236244
```sh
237-
python3 -m analyzer tmp/b2_t0_digitized.csv totEdep --kind csv --save tmp/b2_totEdep.png
245+
python3 -m analyzer b2_t0_digitized.csv totEdep --kind csv --save b2_totEdep.png
238246
```
239247

240248
Or create the same style of histogram without third-party Python packages:
241249

242250
```sh
243-
python3 -B analyzer/svg_plot.py tmp/b2_t0_digitized.csv totEdep --out tmp/b2_totEdep.svg --bins 30
251+
python3 -B analyzer/svg_plot.py b2_t0_digitized.csv totEdep --out b2_totEdep.svg --bins 30
244252
```
245253

246254
### Run B2 With ROOT Output
247255

248-
To produce ROOT output instead of CSV, keep the same `tmp/gemc.db` and run:
256+
To produce ROOT output instead of CSV, keep the same `gemc.db` and run:
249257

250258
```sh
251259
build/gemc examples/basic/b2/b2.yaml \
252260
'-gsystem=[{name: b2, factory: sqlite}]' \
253-
'-gstreamer=[{format: root, filename: tmp/b2}]' \
254-
-sql=tmp/gemc.db \
261+
'-gstreamer=[{format: root, filename: b2}]' \
262+
-sql=gemc.db \
255263
-n=20
256264
```
257265

258266
With one worker thread, this produces:
259267

260268
```text
261-
tmp/b2_t0.root
269+
b2_t0.root
262270
```
263271

264272
Read the ROOT file from Python if you want to inspect or manipulate the data
@@ -267,7 +275,7 @@ before plotting:
267275
```python
268276
from analyzer import plot_variable, read_output
269277

270-
output = read_output("tmp/b2_t0.root", kind="root")
278+
output = read_output("b2_t0.root", kind="root")
271279
print(output.summary())
272280

273281
df = output.get_frame("digitized", detector="flux")
@@ -287,28 +295,14 @@ The Python inspection step is not required for plotting. To plot directly from
287295
the command line, use:
288296

289297
```sh
290-
python3 -m analyzer tmp/b2_t0.root totEdep --kind root --detector flux --save tmp/b2_root_totEdep.png
291-
```
292-
293-
If your shell is inside `tmp/`, do not use `python3 -m ../analyzer`. The `-m`
294-
flag accepts a module name, not a relative path. Use either:
295-
296-
```sh
297-
cd ..
298-
python3 -m analyzer tmp/b2_t0.root totEdep --kind root --detector flux --save tmp/b2_root_totEdep.png
299-
```
300-
301-
or:
302-
303-
```sh
304-
PYTHONPATH=.. python3 -m analyzer b2_t0.root totEdep --kind root --detector flux --save b2_root_totEdep.png
298+
python3 -m analyzer b2_t0.root totEdep --kind root --detector flux --save b2_root_totEdep.png
305299
```
306300

307301
If matplotlib reports that its default cache directory is not writable, set a
308302
writable `MPLCONFIGDIR`:
309303

310304
```sh
311-
PYTHONPATH=.. MPLCONFIGDIR=. python3 -m analyzer b2_t0.root totEdep --kind root --detector flux --save b2_root_totEdep.png
305+
MPLCONFIGDIR=. python3 -m analyzer b2_t0.root totEdep --kind root --detector flux --save b2_root_totEdep.png
312306
```
313307

314308
## Extending Readers

‎analyzer/plotting.py‎

Lines changed: 16 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -22,12 +22,28 @@
2222
"avgTime": "Average Time",
2323
"E": "Energy (MeV)",
2424
"totalE": "Total Energy (MeV)",
25+
"p": "Momentum (MeV)",
26+
"theta": "Theta (rad)",
27+
"phi": "Phi (rad)",
28+
"vx": "Track Vertex X (mm)",
29+
"vy": "Track Vertex Y (mm)",
30+
"vz": "Track Vertex Z (mm)",
31+
"mvx": "Mother Track Vertex X (mm)",
32+
"mvy": "Mother Track Vertex Y (mm)",
33+
"mvz": "Mother Track Vertex Z (mm)",
34+
"mtid": "Mother Track ID",
2535
}
2636

2737
VARIABLE_ALIASES: Mapping[str, tuple[str, ...]] = {
2838
"E": ("totalE", "etot"),
2939
"totalEDeposited": ("totEdep",),
3040
"etot": ("totalE",),
41+
"track_vx": ("vx",),
42+
"track_vy": ("vy",),
43+
"track_vz": ("vz",),
44+
"mother_vx": ("mvx",),
45+
"mother_vy": ("mvy",),
46+
"mother_vz": ("mvz",),
3147
}
3248

3349

‎examples/optical/cherenkov/cherenkov.py‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -109,7 +109,7 @@
109109
# by crossing into a material without matching optical properties.
110110
backplate.material = material
111111

112-
backplate.color = "slategray"
112+
backplate.color = "lightpink"
113113
backplate.digitization = "flux"
114114
backplate.set_identifier("detector", panel_id)
115115
backplate.publish(cfg)

‎gdynamicDigitization/gdynamicdigitization.cc‎

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -32,9 +32,12 @@ std::unique_ptr<GTrueInfoData> GDynamicDigitization::collectTrueInformationImpl(
3232
// Average positions are computed at the hit level by GHit and returned here.
3333
G4ThreeVector avgGlobalPos = ghit->getAvgGlobaPosition();
3434
G4ThreeVector avgLocalPos = ghit->getAvgLocalPosition();
35+
G4ThreeVector trackVertex = ghit->getTrackVertexPosition();
36+
G4ThreeVector motherTrackVertex = ghit->getMotherTrackVertexPosition();
3537

3638
trueInfoData->includeVariable("pid", ghit->getPid());
3739
trueInfoData->includeVariable("tid", ghit->getTid());
40+
trueInfoData->includeVariable("mtid", ghit->getMotherTid());
3841
trueInfoData->includeVariable("totalEDeposited", ghit->getTotalEnergyDeposited());
3942
trueInfoData->includeVariable("avgTime", ghit->getAverageTime());
4043
trueInfoData->includeVariable("avgx", avgGlobalPos.getX());
@@ -43,6 +46,12 @@ std::unique_ptr<GTrueInfoData> GDynamicDigitization::collectTrueInformationImpl(
4346
trueInfoData->includeVariable("avglx", avgLocalPos.getX());
4447
trueInfoData->includeVariable("avgly", avgLocalPos.getY());
4548
trueInfoData->includeVariable("avglz", avgLocalPos.getZ());
49+
trueInfoData->includeVariable("vx", trackVertex.getX());
50+
trueInfoData->includeVariable("vy", trackVertex.getY());
51+
trueInfoData->includeVariable("vz", trackVertex.getZ());
52+
trueInfoData->includeVariable("mvx", motherTrackVertex.getX());
53+
trueInfoData->includeVariable("mvy", motherTrackVertex.getY());
54+
trueInfoData->includeVariable("mvz", motherTrackVertex.getZ());
4655
trueInfoData->includeVariable("hitn", static_cast<int>(hitn)); // assume hitn < INT_MAX
4756

4857
// Bit 1 typically includes metadata like the process name.

0 commit comments

Comments
 (0)