-
Notifications
You must be signed in to change notification settings - Fork 1
/
corrected_areas.py
68 lines (51 loc) · 2.08 KB
/
corrected_areas.py
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
from typing import Tuple
import numpy as np
import strax
import amstrax
export, __all__ = strax.exporter()
@export
@strax.takes_config(
strax.Option(
"elife",
default=30000,
help="electron lifetime in [ns] (should be implemented in db soon)",
),
)
class CorrectedAreas(strax.Plugin):
"""Plugin which applies light collection efficiency maps and electron life time to the data.
Computes the cS1/cS2 for the main/alternative S1/S2 as well as the
corrected life time.
Note:
Please be aware that for both, the main and alternative S1, the
area is corrected according to the xy-position of the main S2.
There are now 3 components of cS2s: cs2_top, cS2_bottom and cs2.
cs2_top and cs2_bottom are corrected by the corresponding maps,
and cs2 is the sum of the two.
"""
__version__ = "0.5.1"
depends_on = ("event_basics", "event_positions")
def infer_dtype(self):
dtype = []
dtype += strax.time_fields
for peak_type, peak_name in zip(["", "alt_"], ["main", "alternate"]):
# Only apply
dtype += [
(f"{peak_type}cs1", np.float32, f"Corrected area of {peak_name} S1 [PE]"),
]
dtype += [
(f"{peak_type}cs2", np.float32, f"Corrected area of {peak_name} S2 [PE]"),
]
return dtype
def compute(self, events):
result = np.zeros(len(events), self.dtype)
result["time"] = events["time"]
result["endtime"] = events["endtime"]
# S1 corrections depend on the actual corrected event position.
# We use this also for the alternate S1; for e.g. Kr this is
# fine as the S1 correction varies slowly.
# event_positions = np.vstack([events["x"], events["y"], events["z"]]).T
elife = self.config["elife"]
for peak_type in ["", "alt_"]:
result[f"{peak_type}cs1"] = events[f"{peak_type}s1_area"]
result[f"{peak_type}cs2"] = events[f"{peak_type}s2_area"]*np.exp(events["drift_time"]/elife)
return result