Skip to content
Open
31 changes: 30 additions & 1 deletion stixpy/product/sources/science.py
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,10 @@
from matplotlib.widgets import Slider
from sunpy.time.timerange import TimeRange


from stixpy.calibration.livetime import get_livetime_fraction
from stixpy.io.readers import read_subc_params

from stixpy.product.product import L1Product

__all__ = [
Expand Down Expand Up @@ -854,7 +857,13 @@ def duration(self):
return self.data["timedel"]

def get_data(
self, time_indices=None, energy_indices=None, detector_indices=None, pixel_indices=None, sum_all_times=False
self,
time_indices=None,
energy_indices=None,
detector_indices=None,
pixel_indices=None,
sum_all_times=False,
livetime_correction=False,
):
"""
Return the counts, errors, times, durations and energies for selected data.
Expand Down Expand Up @@ -888,12 +897,24 @@ def get_data(
Counts, errors, times, deltatimes, energies

"""

# fmt: off
# For the moment copied from idl
trigger_to_detector = [0, 0, 7, 7, 2, 1, 1, 6, 6, 5, 2, 3, 3, 4, 4, 5, 13, 12,
12, 11, 11, 10, 13, 14, 14, 9, 9, 10, 15, 15, 8, 8,]
# fmt: on

triggers = self.data["triggers"][:, trigger_to_detector].astype(float)[...]

_, livefrac, _ = get_livetime_fraction(triggers / self.data["timedel"].to("s").reshape(-1, 1))
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So this is now in the instrument configuration namespace

STIX_INSTRUMENT = SimpleNamespace(
__doc__="""This namespace contains the instrument name, visibility information, sub-collimator ADC mapping, sub-collimator parameters, and pixel configuration.
The visibility information is obtained from the _get_uv_points_data function,
the sub-collimator ADC mapping is read from the detector ADC mapping file,
the sub-collimator parameters are read from the sub-collimator parameters file,
and the pixel configuration is read from the pixel parameters file.
The namespace is used to encapsulate all the configuration data related to the STIX instrument.
This allows for easy access and management of the instrument's configuration data.""",
name="STIX",
vis_info=_get_uv_points_data(),
subcol_adc_mapping=np.array(read_det_adc_mapping()["Adc #"]),
subcol_params=_SUBCOL_PARAMS,
pixel_config=read_pixel_params(),
)


counts = self.data["counts"]
try:
counts_var = self.data["counts_comp_err"] ** 2
except KeyError:
counts_var = self.data["counts_comp_comp_err"] ** 2
shape = counts.shape

if len(shape) < 4:
counts = counts.reshape(shape[0], 1, 1, shape[-1])
counts_var = counts_var.reshape(shape[0], 1, 1, shape[-1])
Expand All @@ -908,6 +929,7 @@ def get_data(
detector_mask[detecor_indices] = True
counts = counts[:, detector_mask, ...]
counts_var = counts_var[:, detector_mask, ...]
livefrac = livefrac[..., detector_mask]
elif detecor_indices.ndim == 2:
counts = np.hstack(
[np.sum(counts[:, dl : dh + 1, ...], axis=1, keepdims=True) for dl, dh in detecor_indices]
Expand Down Expand Up @@ -964,6 +986,7 @@ def get_data(
energies = QTable(energies * u.keV, names=["e_low", "e_high"])

t_norm = self.data["timedel"]

if time_indices is not None:
time_indices = np.asarray(time_indices)
if time_indices.ndim == 1:
Expand All @@ -972,6 +995,7 @@ def get_data(
counts = counts[time_mask, ...]
counts_var = counts_var[time_mask, ...]
t_norm = self.data["timedel"][time_mask]
livefrac = livefrac[time_mask, ...]
times = times[time_mask]
# dT = self.data['timedel'][time_mask]
elif time_indices.ndim == 2:
Expand Down Expand Up @@ -1004,6 +1028,11 @@ def get_data(
if t_norm.size != 1:
t_norm = t_norm.reshape(-1, 1, 1, 1)

livefrac = livefrac.reshape(livefrac.shape + (1, 1))

if livetime_correction:
t_norm = t_norm * livefrac

counts_err = np.sqrt(counts * u.ct + counts_var) / (e_norm * t_norm)
counts = counts / (e_norm * t_norm)

Expand Down