From 6558105cd68ef53102bbf1e04a8b7043d74a479a Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 12:37:26 +0100 Subject: [PATCH 01/13] Fix wrong units being set on QL timeseries * Update plot code to use more of sunpy TimeSeries functions * Add columns keyword with useful defaults --- stixpy/timeseries/quicklook.py | 197 +++++++++------------- stixpy/timeseries/tests/test_quicklook.py | 12 +- 2 files changed, 88 insertions(+), 121 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index e92da11..f39d370 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -1,12 +1,11 @@ +import re from collections import OrderedDict import astropy.units as u -import matplotlib.dates import numpy as np from astropy.io import fits from astropy.table import QTable, Table from astropy.time.core import Time, TimeDelta -from astropy.visualization import quantity_support from sunpy.timeseries.timeseriesbase import GenericTimeSeries from sunpy.visualization import peek_show @@ -23,10 +22,12 @@ def _lfrac(trigger_rate): class QLLightCurve(GenericTimeSeries): - """ + r""" Quicklook X-ray time series. - Nominally in 5 energy bands from 4 - 150 kev + Nominally in five energy bands from 4-150 keV which are obtained by summing counts for all masked detectors and + pixels into five predefined energy ranges. They are double buffered with default integration time of 4s and depth + of 32 bits. Maximum rate is approximately 1MHz with one summed live time counter for the appropriate detectors. Examples -------- @@ -51,7 +52,7 @@ class QLLightCurve(GenericTimeSeries): @peek_show def peek(self, **kwargs): - """ + r""" Displays a graphical overview of the data in this object for user evaluation. For the creation of plots, users should instead use the `.plot` method and Matplotlib's pyplot framework. @@ -66,8 +67,8 @@ def peek(self, **kwargs): self._validate_data_for_plotting() self.plot(**kwargs) - def plot(self, axes=None, **plot_args): - """ + def plot(self, axes=None, columns='counts', **plot_args): + r""" Show a plot of the data. Parameters @@ -84,50 +85,31 @@ def plot(self, axes=None, **plot_args): axes : `~matplotlib.axes.Axes` The plot axes. """ - import matplotlib.pyplot as plt + if columns == 'counts': + count_re = re.compile(r"\d+-\d+ keV$") + columns = [column for column in self.columns if count_re.match(column)] - # Get current axes - if axes is None: - fig, axes = plt.subplots() - else: - fig = plt.gcf() + axes, columns = self._setup_axes_columns(axes, columns) - self._validate_data_for_plotting() - quantity_support() - - dates = matplotlib.dates.date2num(self.to_dataframe().index) + axes = self._data[columns].plot(ax=axes, **plot_args) - labels = [f"{col}" for col in self.columns[5:]] - - lines = [ - axes.plot_date(dates, self.to_dataframe().iloc[:, 5 + i], "-", label=labels[i], **plot_args) - for i in range(5) - ] - - axes.legend(loc="upper right") + units = set([self.units[col] for col in columns]) + if len(units) == 1: + # If units of all columns being plotted are the same, add a unit + # label to the y-axis. + unit = u.Unit(list(units)[0]) + axes.set_ylabel(unit.to_string()) + axes.set_title("STIX QL Light Curve") axes.set_yscale("log") - - axes.set_title("STIX Quick Look") - axes.set_ylabel("count s$^{-1}$ keV$^{-1}$") - - axes.yaxis.grid(True, "major") - axes.xaxis.grid(False, "major") axes.legend() - - # TODO: display better tick labels for date range (e.g. 06/01 - 06/05) - formatter = matplotlib.dates.DateFormatter("%d %H:%M") - axes.xaxis.set_major_formatter(formatter) - - axes.fmt_xdata = matplotlib.dates.DateFormatter("%d %H:%M") - fig.autofmt_xdate() - - return lines + self._setup_x_axis(axes) + return axes @classmethod def _parse_file(cls, filepath): - """ - Parses a GOES/XRS FITS file. + r""" + Parses a STIX file. Parameters ---------- @@ -139,7 +121,7 @@ def _parse_file(cls, filepath): @classmethod def _parse_hdus(cls, hdulist): - """ + r""" Parses STIX FITS data files to create TimeSeries. Parameters @@ -175,8 +157,8 @@ def _parse_hdus(cls, hdulist): units = OrderedDict( [("control_index", None), ("timedel", u.s), ("triggers", None), ("triggers_comp_err", None), ("rcr", None)] ) - units.update([(name, u.ct) for name in names]) - units.update([(f"{name}_comp_err", u.ct) for name in names]) + units.update([(name, u.ct / (u.s * u.keV)) for name in names]) + units.update([(f"{name}_comp_err", u.ct/(u.s * u.keV)) for name in names]) data["triggers"] = data["triggers"].reshape(-1) data["triggers_comp_err"] = data["triggers_comp_err"].reshape(-1) @@ -206,8 +188,12 @@ def __repr__(self): class QLBackground(GenericTimeSeries): - """ - Quicklook X-ray background detector time series. + r""" + Quicklook X-ray time series from the background detector. + + Background monitoring is done in such way that counts from background detector are summed in five specified energy + ranges. These QL data are double buffered into accumulators of 32bit depth. Maximum rate is approximately 30kHz and + one live time counter is available. Integration time is parameter with default value of 32s Examples -------- @@ -250,7 +236,7 @@ class QLBackground(GenericTimeSeries): @peek_show def peek(self, **kwargs): - """ + r""" Displays a graphical overview of the data in this object for user evaluation. For the creation of plots, users should instead use the `.plot` method and Matplotlib's pyplot framework. @@ -265,8 +251,8 @@ def peek(self, **kwargs): self._validate_data_for_plotting() self.plot(**kwargs) - def plot(self, axes=None, **plot_args): - """ + def plot(self, axes=None, columns='counts', **plot_args): + r""" Show a plot of the data. Parameters @@ -274,6 +260,8 @@ def plot(self, axes=None, **plot_args): axes : `~matplotlib.axes.Axes`, optional If provided the image will be plotted on the given axes. Defaults to `None`, so the current axes will be used. + columns : `str`, optional + Columns to plot. Defaults to 'counts'. **plot_args : `dict`, optional Additional plot keyword arguments that are handed to :meth:`pandas.DataFrame.plot`. @@ -283,45 +271,31 @@ def plot(self, axes=None, **plot_args): axes : `~matplotlib.axes.Axes` The plot axes. """ - import matplotlib.pyplot as plt + if columns == 'counts': + count_re = re.compile(r"\d+-\d+ keV$") + columns = [column for column in self.columns if count_re.match(column)] - # Get current axes - if axes is None: - fig, axes = plt.subplots() + axes, columns = self._setup_axes_columns(axes, columns) - self._validate_data_for_plotting() - quantity_support() - - dates = matplotlib.dates.date2num(self.to_dataframe().index) - - labels = [f"{col} keV" for col in self.columns[4:]] + axes = self._data[columns].plot(ax=axes, **plot_args) - [axes.plot_date(dates, self.to_dataframe().iloc[:, 4 + i], "-", label=labels[i], **plot_args) for i in range(5)] - - axes.legend(loc="upper right") + units = set([self.units[col] for col in columns]) + if len(units) == 1: + # If units of all columns being plotted are the same, add a unit + # label to the y-axis. + unit = u.Unit(list(units)[0]) + axes.set_ylabel(unit.to_string()) + axes.set_title("STIX QL Background") axes.set_yscale("log") - - axes.set_title("STIX Quick Look") - axes.set_ylabel("count s$^{-1}$ keV$^{-1}$") - - axes.yaxis.grid(True, "major") - axes.xaxis.grid(False, "major") axes.legend() - - # TODO: display better tick labels for date range (e.g. 06/01 - 06/05) - formatter = matplotlib.dates.DateFormatter("%d %H:%M") - axes.xaxis.set_major_formatter(formatter) - - axes.fmt_xdata = matplotlib.dates.DateFormatter("%d %H:%M") - fig.autofmt_xdate() - - return fig + self._setup_x_axis(axes) + return axes @classmethod def _parse_file(cls, filepath): """ - Parses a GOES/XRS FITS file. + Parses a STIX FITS file. Parameters ---------- @@ -353,7 +327,7 @@ def _parse_hdus(cls, hdulist): data["counts"] = data["counts"] / (live_time.reshape(-1, 1) * energy_delta) - names = [f'{energies["e_low"][i]}-{energies["e_high"][i]}' for i in range(5)] + names = [f'{energies["e_low"][i].astype(int)}-{energies["e_high"][i].astype(int)} keV' for i in range(5)] [data.add_column(data["counts"][:, i], name=names[i]) for i in range(5)] data.remove_column("counts") @@ -366,13 +340,11 @@ def _parse_hdus(cls, hdulist): except KeyError: data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.s) - # [f'{energies[i]["e_low"]} - {energies[i]["e_high"]} keV' for i in range(5)] - units = OrderedDict( [("control_index", None), ("timedel", u.s), ("triggers", None), ("triggers_comp_err", None), ("rcr", None)] ) - units.update([(name, u.ct) for name in names]) - units.update([(f"{name}_comp_err", u.ct) for name in names]) + units.update([(name, u.ct/(u.s * u.keV)) for name in names]) + units.update([(f"{name}_comp_err", u.ct/(u.s * u.keV)) for name in names]) data["triggers"] = data["triggers"].reshape(-1) data["triggers_comp_err"] = data["triggers_comp_err"].reshape(-1) @@ -400,7 +372,11 @@ def is_datasource_for(cls, **kwargs): class QLVariance(GenericTimeSeries): """ - Quicklook X-ray background detector time series. + Quicklook variance time series + + Variance of counts is calculated for one energy range over counts summed from selected detectors and pixels in 40 + accumulators, each accumulating for 0.1s. These accumulators are double buffered with depth of 32bits and + approximate data rate of 1 Mhz. Examples -------- @@ -440,7 +416,7 @@ class QLVariance(GenericTimeSeries): [21516 rows x 4 columns] """ - def plot(self, axes=None, **plot_args): + def plot(self, axes=None, columns='variance', **plot_args): """ Show a plot of the data. @@ -458,40 +434,26 @@ def plot(self, axes=None, **plot_args): axes : `~matplotlib.axes.Axes` The plot axes. """ - import matplotlib.pyplot as plt - - # Get current axes - if axes is None: - fig, axes = plt.subplots() - - self._validate_data_for_plotting() - quantity_support() - - dates = matplotlib.dates.date2num(self.to_dataframe().index) + if columns == 'variance': + count_re = re.compile(r"\d+-\d+ keV$") + columns = [column for column in self.columns if count_re.match(column)] - label = f"{self.columns[2]} keV" + axes, columns = self._setup_axes_columns(axes, columns) - axes.plot_date(dates, self.to_dataframe().iloc[:, 2], "-", label=label) # , **plot_args) + axes = self._data[columns].plot(ax=axes, **plot_args) - axes.legend(loc="upper right") + units = set([self.units[col] for col in columns]) + if len(units) == 1: + # If units of all columns being plotted are the same, add a unit + # label to the y-axis. + unit = u.Unit(list(units)[0]) + axes.set_ylabel(unit.to_string()) + axes.set_title("STIX QL Variance") axes.set_yscale("log") - - axes.set_title("STIX Quick Look Variance") - axes.set_ylabel("count s$^{-1}$ keV$^{-1}$") - - axes.yaxis.grid(True, "major") - axes.xaxis.grid(False, "major") axes.legend() - - # TODO: display better tick labels for date range (e.g. 06/01 - 06/05) - formatter = matplotlib.dates.DateFormatter("%d %H:%M") - axes.xaxis.set_major_formatter(formatter) - - axes.fmt_xdata = matplotlib.dates.DateFormatter("%d %H:%M") - fig.autofmt_xdate() - - return fig + self._setup_x_axis(axes) + return axes @classmethod def _parse_file(cls, filepath): @@ -526,8 +488,8 @@ def _parse_hdus(cls, hdulist): << u.keV ) name = ( - f'{energies[control["energy_bin_mask"][0]]["e_low"][0]}' - f'-{energies[control["energy_bin_mask"][0]]["e_high"][-1]}' + f'{energies[control["energy_bin_mask"][0]]["e_low"][0].astype(int)}' + f'-{energies[control["energy_bin_mask"][0]]["e_high"][-1].astype(int)} keV' ) try: @@ -542,7 +504,8 @@ def _parse_hdus(cls, hdulist): data.add_column(data["variance_comp_err"], name=f"{name}_comp_err") data.remove_column("variance_comp_err") - units = OrderedDict([("control_index", None), ("timedel", u.s), (name, u.count), (f"{name}_comp_err", u.count)]) + units = OrderedDict([("control_index", None), ("timedel", u.s), + (name, u.ct/(u.s * u.keV)), (f"{name}_comp_err", u.ct/(u.s * u.keV))]) data[name] = data[name].reshape(-1) data[f"{name}_comp_err"] = data[f"{name}_comp_err"].reshape(-1) diff --git a/stixpy/timeseries/tests/test_quicklook.py b/stixpy/timeseries/tests/test_quicklook.py index 6221877..e64e7de 100644 --- a/stixpy/timeseries/tests/test_quicklook.py +++ b/stixpy/timeseries/tests/test_quicklook.py @@ -1,3 +1,4 @@ +import astropy.units as u import pytest from sunpy.timeseries import TimeSeries @@ -10,22 +11,25 @@ def test_ql_lightcurve(): r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/QL/solo_L1_stix-ql-lightcurve_20200506_V02.fits" ) assert isinstance(ql_lc, QLLightCurve) + assert ql_lc.quantity('4-10 keV').unit == u.ct/(u.s * u.keV) @pytest.mark.remote_data def test_qlbackground(): - ql_lc = TimeSeries( + ql_bkg = TimeSeries( r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/QL/solo_L1_stix-ql-background_20200506_V02.fits" ) - assert isinstance(ql_lc, QLBackground) + assert isinstance(ql_bkg, QLBackground) + assert ql_bkg.quantity('4-10 keV').unit == u.ct / (u.s * u.keV) @pytest.mark.remote_data def test_qlvariance(): - ql_lc = TimeSeries( + ql_var = TimeSeries( r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/QL/solo_L1_stix-ql-variance_20200506_V02.fits" ) - assert isinstance(ql_lc, QLVariance) + assert isinstance(ql_var, QLVariance) + assert ql_var.quantity('4-20 keV').unit == u.ct / (u.s * u.keV) @pytest.mark.remote_data From cccf11f87dd2f8a669b02185bfd95bb873d78005 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 13:40:55 +0100 Subject: [PATCH 02/13] Add change log --- changelog/118.bugfix.rst | 1 + 1 file changed, 1 insertion(+) create mode 100644 changelog/118.bugfix.rst diff --git a/changelog/118.bugfix.rst b/changelog/118.bugfix.rst new file mode 100644 index 0000000..fd2f311 --- /dev/null +++ b/changelog/118.bugfix.rst @@ -0,0 +1 @@ +Fix a bug where the incorrect unit was set on :class:`~stixpy.timeseries.quicklook.QLLightCurve`, :class:`~stixpy.timeseries.quicklook.QLBackground` and :class:`~stixpy.timeseries.quicklook.QLVariance`. From 47650c9732e97cb5d89d8a1237487608dcc289c2 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 13:53:06 +0100 Subject: [PATCH 03/13] Remove redundant code fix columns for non default value. --- stixpy/timeseries/quicklook.py | 42 +++------------------------------- 1 file changed, 3 insertions(+), 39 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index f39d370..c54674a 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -7,7 +7,6 @@ from astropy.table import QTable, Table from astropy.time.core import Time, TimeDelta from sunpy.timeseries.timeseriesbase import GenericTimeSeries -from sunpy.visualization import peek_show __all__ = ["QLLightCurve", "QLBackground", "QLVariance", "HKMaxi"] @@ -50,23 +49,6 @@ class QLLightCurve(GenericTimeSeries): _source = "stix" - @peek_show - def peek(self, **kwargs): - r""" - Displays a graphical overview of the data in this object for user evaluation. - For the creation of plots, users should instead use the - `.plot` method and Matplotlib's pyplot framework. - - Parameters - ---------- - **kwargs : `dict` - Any additional plot arguments that should be used when plotting. - """ - - # Check we have a timeseries valid for plotting - self._validate_data_for_plotting() - self.plot(**kwargs) - def plot(self, axes=None, columns='counts', **plot_args): r""" Show a plot of the data. @@ -94,7 +76,7 @@ def plot(self, axes=None, columns='counts', **plot_args): axes = self._data[columns].plot(ax=axes, **plot_args) units = set([self.units[col] for col in columns]) - if len(units) == 1: + if len(units) == 1 and list(units)[0] is not None: # If units of all columns being plotted are the same, add a unit # label to the y-axis. unit = u.Unit(list(units)[0]) @@ -234,23 +216,6 @@ class QLBackground(GenericTimeSeries): """ - @peek_show - def peek(self, **kwargs): - r""" - Displays a graphical overview of the data in this object for user evaluation. - For the creation of plots, users should instead use the - `.plot` method and Matplotlib's pyplot framework. - - Parameters - ---------- - **kwargs : `dict` - Any additional plot arguments that should be used when plotting. - """ - - # Check we have a timeseries valid for plotting - self._validate_data_for_plotting() - self.plot(**kwargs) - def plot(self, axes=None, columns='counts', **plot_args): r""" Show a plot of the data. @@ -280,7 +245,7 @@ def plot(self, axes=None, columns='counts', **plot_args): axes = self._data[columns].plot(ax=axes, **plot_args) units = set([self.units[col] for col in columns]) - if len(units) == 1: + if len(units) == 1 and list(units)[0] is not None: # If units of all columns being plotted are the same, add a unit # label to the y-axis. unit = u.Unit(list(units)[0]) @@ -439,11 +404,10 @@ def plot(self, axes=None, columns='variance', **plot_args): columns = [column for column in self.columns if count_re.match(column)] axes, columns = self._setup_axes_columns(axes, columns) - axes = self._data[columns].plot(ax=axes, **plot_args) units = set([self.units[col] for col in columns]) - if len(units) == 1: + if len(units) == 1 and list(units)[0] is not None: # If units of all columns being plotted are the same, add a unit # label to the y-axis. unit = u.Unit(list(units)[0]) From c87e2a0a0fa994a0303a3b650eab760bffdf1db8 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 14:00:36 +0100 Subject: [PATCH 04/13] Make sure peek works as expected --- stixpy/timeseries/quicklook.py | 17 +++++++++++------ 1 file changed, 11 insertions(+), 6 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index c54674a..8be2a94 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -49,7 +49,7 @@ class QLLightCurve(GenericTimeSeries): _source = "stix" - def plot(self, axes=None, columns='counts', **plot_args): + def plot(self, axes=None, columns=None, **plot_args): r""" Show a plot of the data. @@ -58,6 +58,8 @@ def plot(self, axes=None, columns='counts', **plot_args): axes : `~matplotlib.axes.Axes`, optional If provided the image will be plotted on the given axes. Defaults to `None`, so the current axes will be used. + columns : `list`, optional + Columns to plot, defaults to 'counts'. **plot_args : `dict`, optional Additional plot keyword arguments that are handed to :meth:`pandas.DataFrame.plot`. @@ -67,7 +69,7 @@ def plot(self, axes=None, columns='counts', **plot_args): axes : `~matplotlib.axes.Axes` The plot axes. """ - if columns == 'counts': + if columns is None: count_re = re.compile(r"\d+-\d+ keV$") columns = [column for column in self.columns if count_re.match(column)] @@ -216,7 +218,7 @@ class QLBackground(GenericTimeSeries): """ - def plot(self, axes=None, columns='counts', **plot_args): + def plot(self, axes=None, columns=None, **plot_args): r""" Show a plot of the data. @@ -226,7 +228,7 @@ def plot(self, axes=None, columns='counts', **plot_args): If provided the image will be plotted on the given axes. Defaults to `None`, so the current axes will be used. columns : `str`, optional - Columns to plot. Defaults to 'counts'. + Columns to plot, defaults to 'counts'. **plot_args : `dict`, optional Additional plot keyword arguments that are handed to :meth:`pandas.DataFrame.plot`. @@ -381,15 +383,18 @@ class QLVariance(GenericTimeSeries): [21516 rows x 4 columns] """ - def plot(self, axes=None, columns='variance', **plot_args): + def plot(self, axes=None, columns=None, **plot_args): """ Show a plot of the data. Parameters ---------- + axes : `~matplotlib.axes.Axes`, optional If provided the image will be plotted on the given axes. Defaults to `None`, so the current axes will be used. + columns : + Columns to plot, defaults to 'variance'. **plot_args : `dict`, optional Additional plot keyword arguments that are handed to :meth:`pandas.DataFrame.plot`. @@ -399,7 +404,7 @@ def plot(self, axes=None, columns='variance', **plot_args): axes : `~matplotlib.axes.Axes` The plot axes. """ - if columns == 'variance': + if columns is None: count_re = re.compile(r"\d+-\d+ keV$") columns = [column for column in self.columns if count_re.match(column)] From 384323052a2f0219dc8a8a96a23b2d92c1e44658 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 16:06:04 +0100 Subject: [PATCH 05/13] Fix peek --- stixpy/timeseries/quicklook.py | 15 ++++++++------- 1 file changed, 8 insertions(+), 7 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index 8be2a94..6f5cd66 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -83,9 +83,10 @@ def plot(self, axes=None, columns=None, **plot_args): # label to the y-axis. unit = u.Unit(list(units)[0]) axes.set_ylabel(unit.to_string()) + if unit == u.ct/(u.s * u.keV): + axes.set_yscale("log") axes.set_title("STIX QL Light Curve") - axes.set_yscale("log") axes.legend() self._setup_x_axis(axes) return axes @@ -119,7 +120,7 @@ def _parse_hdus(cls, hdulist): energies = QTable(hdulist[4].data) energy_delta = energies["e_high"] - energies["e_low"] << u.keV - live_time = _lfrac(data["triggers"].reshape(-1) / (16 * data["timedel"] * u.s)) + live_time = _lfrac(data["triggers"].reshape(-1) / (16 * data["timedel"] * u.cs)) data["counts"] = data["counts"] / (live_time.reshape(-1, 1) * energy_delta) @@ -139,7 +140,7 @@ def _parse_hdus(cls, hdulist): data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.cs) units = OrderedDict( - [("control_index", None), ("timedel", u.s), ("triggers", None), ("triggers_comp_err", None), ("rcr", None)] + [("control_index", None), ("timedel", u.cs), ("triggers", None), ("triggers_comp_err", None), ("rcr", None)] ) units.update([(name, u.ct / (u.s * u.keV)) for name in names]) units.update([(f"{name}_comp_err", u.ct/(u.s * u.keV)) for name in names]) @@ -308,7 +309,7 @@ def _parse_hdus(cls, hdulist): data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.s) units = OrderedDict( - [("control_index", None), ("timedel", u.s), ("triggers", None), ("triggers_comp_err", None), ("rcr", None)] + [("control_index", None), ("timedel", u.cs), ("triggers", None), ("triggers_comp_err", None), ("rcr", None)] ) units.update([(name, u.ct/(u.s * u.keV)) for name in names]) units.update([(f"{name}_comp_err", u.ct/(u.s * u.keV)) for name in names]) @@ -473,7 +474,7 @@ def _parse_hdus(cls, hdulist): data.add_column(data["variance_comp_err"], name=f"{name}_comp_err") data.remove_column("variance_comp_err") - units = OrderedDict([("control_index", None), ("timedel", u.s), + units = OrderedDict([("control_index", None), ("timedel", u.cs), (name, u.ct/(u.s * u.keV)), (f"{name}_comp_err", u.ct/(u.s * u.keV))]) data[name] = data[name].reshape(-1) @@ -596,9 +597,9 @@ def _parse_hdus(cls, hdulist): data = Table(hdulist[2].data) try: - data["time"] = Time(header["date_obs"]) + TimeDelta(data["time"] * u.s) + data["time"] = Time(header["date_obs"]) + TimeDelta(data["time"] * u.cs) except KeyError: - data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.s) + data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.cs) data_df = data.to_pandas() data_df.index = data_df["time"] From d256b34521cfbe4dc0577ce53243d2e21f6d11e4 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 16:12:48 +0100 Subject: [PATCH 06/13] Remove peek fix doc string --- stixpy/timeseries/quicklook.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index 6f5cd66..370a495 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -178,7 +178,7 @@ class QLBackground(GenericTimeSeries): Background monitoring is done in such way that counts from background detector are summed in five specified energy ranges. These QL data are double buffered into accumulators of 32bit depth. Maximum rate is approximately 30kHz and - one live time counter is available. Integration time is parameter with default value of 32s + one live time counter is available. Integration time is parameter with default value of 32s Examples -------- From ac5c1e1f732b78077af7abf227e40a11d77d1608 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 20:54:55 +0100 Subject: [PATCH 07/13] Better tests, create QTables from hdu and get units --- stixpy/timeseries/quicklook.py | 101 ++++++++++++---------- stixpy/timeseries/tests/test_quicklook.py | 56 +++++++++--- 2 files changed, 102 insertions(+), 55 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index 370a495..4a3d17d 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -4,6 +4,8 @@ import astropy.units as u import numpy as np from astropy.io import fits +from astropy.io.fits import BinTableHDU, Header +from astropy.io.fits.connect import read_table_fits from astropy.table import QTable, Table from astropy.time.core import Time, TimeDelta from sunpy.timeseries.timeseriesbase import GenericTimeSeries @@ -15,6 +17,26 @@ tau = 12.5 * u.us +def _hdu_to_qtable(hdupair): + r""" + + + Parameters + ---------- + hdupair + HDU, header pair + + Returns + ------- + + """ + header = hdupair.header + header.pop('KEYCOMMENTS') + bintable = BinTableHDU(hdupair.data, header=Header(cards=header)) + table = read_table_fits(bintable) + qtable = QTable(table) + return qtable + def _lfrac(trigger_rate): nin = trigger_rate / (1 - (trigger_rate * (tau+eta))) return np.exp(-eta * nin) / (1 + tau * nin) @@ -114,20 +136,17 @@ def _parse_hdus(cls, hdulist): hdulist : """ header = hdulist[0].header - control = QTable(hdulist[1].data) - header["control"] = control - data = QTable(hdulist[2].data) - energies = QTable(hdulist[4].data) - energy_delta = energies["e_high"] - energies["e_low"] << u.keV + # control = _hdu_to_qtable(hdulist[1]) + data = _hdu_to_qtable(hdulist[2]) + energies = _hdu_to_qtable(hdulist[4]) + energy_delta = energies["e_high"] - energies["e_low"] - live_time = _lfrac(data["triggers"].reshape(-1) / (16 * data["timedel"] * u.cs)) + live_time = _lfrac(data["triggers"].reshape(-1) / (16 * data["timedel"])) - data["counts"] = data["counts"] / (live_time.reshape(-1, 1) * energy_delta) + data["counts"] = data["counts"] / ((data['timedel'].to(u.s)*live_time).reshape(-1, 1)*energy_delta * energy_delta) names = [ - "{:d}-{:d} keV".format(energies["e_low"][i].astype(int), energies["e_high"][i].astype(int)) - for i in range(5) - ] + f"{energies['e_low'][i].value.astype(int)}-{energies['e_high'][i].value.astype(int)} {energies['e_high'].unit}" for i in range(5)] [data.add_column(data["counts"][:, i], name=names[i]) for i in range(5)] data.remove_column("counts") @@ -135,13 +154,11 @@ def _parse_hdus(cls, hdulist): data.remove_column("counts_comp_err") try: - data["time"] = Time(header["date_obs"]) + data["time"] * u.cs + data["time"] = Time(header["date_obs"]) + data["time"] except KeyError: - data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.cs) + data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"]) - units = OrderedDict( - [("control_index", None), ("timedel", u.cs), ("triggers", None), ("triggers_comp_err", None), ("rcr", None)] - ) + units = OrderedDict((c.info.name, c.unit) for c in data.itercols() if c.info.name != 'time') units.update([(name, u.ct / (u.s * u.keV)) for name in names]) units.update([(f"{name}_comp_err", u.ct/(u.s * u.keV)) for name in names]) @@ -154,6 +171,8 @@ def _parse_hdus(cls, hdulist): return data_df, header, units + + @classmethod def is_datasource_for(cls, **kwargs): """ @@ -239,7 +258,7 @@ def plot(self, axes=None, columns=None, **plot_args): axes : `~matplotlib.axes.Axes` The plot axes. """ - if columns == 'counts': + if columns is None: count_re = re.compile(r"\d+-\d+ keV$") columns = [column for column in self.columns if count_re.match(column)] @@ -284,18 +303,16 @@ def _parse_hdus(cls, hdulist): The path to the file you want to parse. """ header = hdulist[0].header - control = Table(hdulist[1].data) - header["control"] = control - data = Table(hdulist[2].data) - energies = Table(hdulist[4].data) - - energy_delta = energies["e_high"] - energies["e_low"] << u.keV + # control = _hdu_to_qtable(hdulist[1]) + data = _hdu_to_qtable(hdulist[2]) + energies = _hdu_to_qtable(hdulist[4]) + energy_delta = energies["e_high"] - energies["e_low"] - live_time = _lfrac(data["triggers"].reshape(-1) / (16 * data["timedel"] * u.s)) + live_time = _lfrac(data["triggers"].reshape(-1) / data["timedel"]) - data["counts"] = data["counts"] / (live_time.reshape(-1, 1) * energy_delta) + data["counts"] = data["counts"] / ((data['timedel'].to(u.s)*live_time).reshape(-1, 1) * energy_delta) - names = [f'{energies["e_low"][i].astype(int)}-{energies["e_high"][i].astype(int)} keV' for i in range(5)] + names = [f'{energies["e_low"][i].value.astype(int)}-{energies["e_high"][i].value.astype(int)} {energies["e_high"].unit}' for i in range(5)] [data.add_column(data["counts"][:, i], name=names[i]) for i in range(5)] data.remove_column("counts") @@ -304,13 +321,11 @@ def _parse_hdus(cls, hdulist): data.remove_column("counts_comp_err") try: - data["time"] = Time(header["date_obs"]) + TimeDelta(data["time"] * u.s) + data["time"] = Time(header["date_obs"]) + TimeDelta(data["time"]) except KeyError: - data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.s) + data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"]) - units = OrderedDict( - [("control_index", None), ("timedel", u.cs), ("triggers", None), ("triggers_comp_err", None), ("rcr", None)] - ) + units = OrderedDict((c.info.name, c.unit) for c in data.itercols() if c.info.name != 'time') units.update([(name, u.ct/(u.s * u.keV)) for name in names]) units.update([(f"{name}_comp_err", u.ct/(u.s * u.keV)) for name in names]) @@ -418,9 +433,9 @@ def plot(self, axes=None, columns=None, **plot_args): # label to the y-axis. unit = u.Unit(list(units)[0]) axes.set_ylabel(unit.to_string()) + axes.set_yscale("log") axes.set_title("STIX QL Variance") - axes.set_yscale("log") axes.legend() self._setup_x_axis(axes) return axes @@ -449,33 +464,31 @@ def _parse_hdus(cls, hdulist): The path to the file you want to parse. """ header = hdulist[0].header - control = Table(hdulist[1].data) - header["control"] = control - data = Table(hdulist[2].data) - energies = Table(hdulist[4].data) - dE = ( + control = _hdu_to_qtable(hdulist[1]) + data = _hdu_to_qtable(hdulist[2]) + energies = _hdu_to_qtable(hdulist[4]) + delta_energy = ( energies[control["energy_bin_mask"][0]]["e_high"][-1] - energies[control["energy_bin_mask"][0]]["e_low"][0] - << u.keV ) name = ( - f'{energies[control["energy_bin_mask"][0]]["e_low"][0].astype(int)}' - f'-{energies[control["energy_bin_mask"][0]]["e_high"][-1].astype(int)} keV' + f'{energies[control["energy_bin_mask"][0]]["e_low"][0].value.astype(int)}' + f'-{energies[control["energy_bin_mask"][0]]["e_high"][-1].value.astype(int)} {energies["e_high"].unit}' ) try: - data["time"] = Time(header["date_obs"]) + TimeDelta(data["time"] * u.s) + data["time"] = Time(header["date_obs"]) + TimeDelta(data["time"]) except KeyError: - data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.s) + data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"]) - data["variance"] = data["variance"].reshape(-1) / (dE * data["timedel"]) + data["variance"] = data["variance"].reshape(-1) * u.ct / (delta_energy * data["timedel"].to(u.s)) data.add_column(data["variance"], name=name) data.remove_column("variance") data.add_column(data["variance_comp_err"], name=f"{name}_comp_err") data.remove_column("variance_comp_err") - units = OrderedDict([("control_index", None), ("timedel", u.cs), - (name, u.ct/(u.s * u.keV)), (f"{name}_comp_err", u.ct/(u.s * u.keV))]) + units = OrderedDict((c.info.name, c.unit) for c in data.itercols() if c.info.name != 'time') + units[f"{name}_comp_err"] = u.ct/(u.s * u.keV) data[name] = data[name].reshape(-1) data[f"{name}_comp_err"] = data[f"{name}_comp_err"].reshape(-1) diff --git a/stixpy/timeseries/tests/test_quicklook.py b/stixpy/timeseries/tests/test_quicklook.py index e64e7de..ea93b3c 100644 --- a/stixpy/timeseries/tests/test_quicklook.py +++ b/stixpy/timeseries/tests/test_quicklook.py @@ -6,30 +6,64 @@ @pytest.mark.remote_data -def test_ql_lightcurve(): +@pytest.fixture() +def ql_lightcurve(): ql_lc = TimeSeries( r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/QL/solo_L1_stix-ql-lightcurve_20200506_V02.fits" ) - assert isinstance(ql_lc, QLLightCurve) - assert ql_lc.quantity('4-10 keV').unit == u.ct/(u.s * u.keV) - + return ql_lc @pytest.mark.remote_data -def test_qlbackground(): +@pytest.fixture() +def ql_background(): ql_bkg = TimeSeries( r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/QL/solo_L1_stix-ql-background_20200506_V02.fits" ) - assert isinstance(ql_bkg, QLBackground) - assert ql_bkg.quantity('4-10 keV').unit == u.ct / (u.s * u.keV) - + return ql_bkg @pytest.mark.remote_data -def test_qlvariance(): +@pytest.fixture() +def ql_variance(): ql_var = TimeSeries( r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/QL/solo_L1_stix-ql-variance_20200506_V02.fits" ) - assert isinstance(ql_var, QLVariance) - assert ql_var.quantity('4-20 keV').unit == u.ct / (u.s * u.keV) + return ql_var + +@pytest.mark.remote_data +def test_ql_lightcurve(ql_lightcurve): + assert isinstance(ql_lightcurve, QLLightCurve) + assert ql_lightcurve.quantity('4-10 keV').unit == u.ct/(u.s * u.keV) + +@pytest.mark.remote_data +def test_ql_lightcurve_plot(ql_lightcurve): + ql_lightcurve.plot() + ql_lightcurve.plot(columns=['4-10 keV']) + + +@pytest.mark.remote_data +def test_ql_background(ql_background): + assert isinstance(ql_background, QLBackground) + assert ql_background.quantity('4-10 keV').unit == u.ct / (u.s * u.keV) + + +@pytest.mark.remote_data +def test_ql_background_plot(ql_background): + ql_background.plot() + ql_background.plot(columns=['4-10 keV']) + + + +@pytest.mark.remote_data +def test_qlvariance(ql_variance): + assert isinstance(ql_variance, QLVariance) + assert ql_variance.quantity('4-20 keV').unit == u.ct / (u.s * u.keV) + + +@pytest.mark.remote_data +def test_qlvariance_plot(ql_variance): + ql_variance.plot() + ql_variance.plot(columns=['4-20 keV']) + @pytest.mark.remote_data From e7cd32bbb6fbfc1b6f30f2d966671682b1f0beac Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 21:51:38 +0100 Subject: [PATCH 08/13] Fix strange but on python3.9 --- stixpy/timeseries/quicklook.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index 4a3d17d..85f60c7 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -31,7 +31,9 @@ def _hdu_to_qtable(hdupair): """ header = hdupair.header - header.pop('KEYCOMMENTS') + # TODD remove when python 3.9 and astropy 5.3.4 dropped + header.pop('keycomments') + header.pop('comment') bintable = BinTableHDU(hdupair.data, header=Header(cards=header)) table = read_table_fits(bintable) qtable = QTable(table) From 206268e36dede73aca4033bc0ccad4a347c852e5 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Thu, 16 May 2024 23:22:33 +0100 Subject: [PATCH 09/13] Replace live time correction with more accurate version. --- stixpy/timeseries/quicklook.py | 35 ++++++++++++---------------------- 1 file changed, 12 insertions(+), 23 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index 85f60c7..a1e1661 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -2,7 +2,6 @@ from collections import OrderedDict import astropy.units as u -import numpy as np from astropy.io import fits from astropy.io.fits import BinTableHDU, Header from astropy.io.fits.connect import read_table_fits @@ -10,16 +9,14 @@ from astropy.time.core import Time, TimeDelta from sunpy.timeseries.timeseriesbase import GenericTimeSeries -__all__ = ["QLLightCurve", "QLBackground", "QLVariance", "HKMaxi"] - +from stixpy.calibration.livetime import get_livetime_fraction -eta = 2.5 * u.us -tau = 12.5 * u.us +__all__ = ["QLLightCurve", "QLBackground", "QLVariance", "HKMaxi"] def _hdu_to_qtable(hdupair): r""" - + Given a HDU pair, convert it to a QTable Parameters ---------- @@ -28,10 +25,10 @@ def _hdu_to_qtable(hdupair): Returns ------- - + QTable """ header = hdupair.header - # TODD remove when python 3.9 and astropy 5.3.4 dropped + # TODD remove when python 3.9 and astropy 5.3.4 are dropped weird non-ascii error header.pop('keycomments') header.pop('comment') bintable = BinTableHDU(hdupair.data, header=Header(cards=header)) @@ -39,10 +36,6 @@ def _hdu_to_qtable(hdupair): qtable = QTable(table) return qtable -def _lfrac(trigger_rate): - nin = trigger_rate / (1 - (trigger_rate * (tau+eta))) - return np.exp(-eta * nin) / (1 + tau * nin) - class QLLightCurve(GenericTimeSeries): r""" @@ -143,9 +136,9 @@ def _parse_hdus(cls, hdulist): energies = _hdu_to_qtable(hdulist[4]) energy_delta = energies["e_high"] - energies["e_low"] - live_time = _lfrac(data["triggers"].reshape(-1) / (16 * data["timedel"])) + live_frac, *_ = get_livetime_fraction(data["triggers"].reshape(-1) / (16 * data["timedel"])) - data["counts"] = data["counts"] / ((data['timedel'].to(u.s)*live_time).reshape(-1, 1)*energy_delta * energy_delta) + data["counts"] = data["counts"] / ((data['timedel'].to(u.s) * live_frac).reshape(-1, 1) * energy_delta) names = [ f"{energies['e_low'][i].value.astype(int)}-{energies['e_high'][i].value.astype(int)} {energies['e_high'].unit}" for i in range(5)] @@ -161,8 +154,7 @@ def _parse_hdus(cls, hdulist): data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"]) units = OrderedDict((c.info.name, c.unit) for c in data.itercols() if c.info.name != 'time') - units.update([(name, u.ct / (u.s * u.keV)) for name in names]) - units.update([(f"{name}_comp_err", u.ct/(u.s * u.keV)) for name in names]) + units.update([(f"{name}_comp_err", units[name]) for name in names]) data["triggers"] = data["triggers"].reshape(-1) data["triggers_comp_err"] = data["triggers_comp_err"].reshape(-1) @@ -173,8 +165,6 @@ def _parse_hdus(cls, hdulist): return data_df, header, units - - @classmethod def is_datasource_for(cls, **kwargs): """ @@ -274,9 +264,9 @@ def plot(self, axes=None, columns=None, **plot_args): # label to the y-axis. unit = u.Unit(list(units)[0]) axes.set_ylabel(unit.to_string()) + axes.set_yscale("log") axes.set_title("STIX QL Background") - axes.set_yscale("log") axes.legend() self._setup_x_axis(axes) return axes @@ -310,9 +300,9 @@ def _parse_hdus(cls, hdulist): energies = _hdu_to_qtable(hdulist[4]) energy_delta = energies["e_high"] - energies["e_low"] - live_time = _lfrac(data["triggers"].reshape(-1) / data["timedel"]) + live_frac, *_ = get_livetime_fraction(data["triggers"].reshape(-1) / data["timedel"]) - data["counts"] = data["counts"] / ((data['timedel'].to(u.s)*live_time).reshape(-1, 1) * energy_delta) + data["counts"] = data["counts"] / ((data['timedel'].to(u.s) * live_frac).reshape(-1, 1) * energy_delta) names = [f'{energies["e_low"][i].value.astype(int)}-{energies["e_high"][i].value.astype(int)} {energies["e_high"].unit}' for i in range(5)] @@ -328,8 +318,7 @@ def _parse_hdus(cls, hdulist): data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"]) units = OrderedDict((c.info.name, c.unit) for c in data.itercols() if c.info.name != 'time') - units.update([(name, u.ct/(u.s * u.keV)) for name in names]) - units.update([(f"{name}_comp_err", u.ct/(u.s * u.keV)) for name in names]) + units.update([(f"{name}_comp_err", units[name]) for name in names]) data["triggers"] = data["triggers"].reshape(-1) data["triggers_comp_err"] = data["triggers_comp_err"].reshape(-1) From a4da0e3d99fea598b36e6428770d9bb94c4f14ca Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Fri, 17 May 2024 09:51:49 +0100 Subject: [PATCH 10/13] Update stixpy/timeseries/quicklook.py Co-authored-by: DanRyanIrish --- stixpy/timeseries/quicklook.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index a1e1661..91ded0b 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -28,7 +28,7 @@ def _hdu_to_qtable(hdupair): QTable """ header = hdupair.header - # TODD remove when python 3.9 and astropy 5.3.4 are dropped weird non-ascii error + # TODO remove when python 3.9 and astropy 5.3.4 are dropped weird non-ascii error header.pop('keycomments') header.pop('comment') bintable = BinTableHDU(hdupair.data, header=Header(cards=header)) From bb7d703790ead07cffdd1a7323a3a09ffe5bd1ea Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Fri, 17 May 2024 11:23:18 +0100 Subject: [PATCH 11/13] Apply suggestions from code review Co-authored-by: Laura Hayes --- stixpy/timeseries/quicklook.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index 91ded0b..86b7456 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -41,7 +41,7 @@ class QLLightCurve(GenericTimeSeries): r""" Quicklook X-ray time series. - Nominally in five energy bands from 4-150 keV which are obtained by summing counts for all masked detectors and + Nominally the quicklook data files contain the STIX timeseries data in five energy bands from 4-150 keV which are obtained by summing counts for all masked detectors and pixels into five predefined energy ranges. They are double buffered with default integration time of 4s and depth of 32 bits. Maximum rate is approximately 1MHz with one summed live time counter for the appropriate detectors. @@ -111,7 +111,7 @@ def plot(self, axes=None, columns=None, **plot_args): @classmethod def _parse_file(cls, filepath): r""" - Parses a STIX file. + Parses a STIX QL file. Parameters ---------- From 60298dd16fd3c7e641cf96d4ffd911dba9d449b2 Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Fri, 17 May 2024 11:43:44 +0100 Subject: [PATCH 12/13] Codre reive and update HK timeseries in the same ways as QL --- stixpy/timeseries/quicklook.py | 22 ++++++++-------- stixpy/timeseries/tests/test_quicklook.py | 31 +++++++++++++++++------ 2 files changed, 35 insertions(+), 18 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index 86b7456..5ce0152 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -5,7 +5,7 @@ from astropy.io import fits from astropy.io.fits import BinTableHDU, Header from astropy.io.fits.connect import read_table_fits -from astropy.table import QTable, Table +from astropy.table import QTable from astropy.time.core import Time, TimeDelta from sunpy.timeseries.timeseriesbase import GenericTimeSeries @@ -187,7 +187,7 @@ class QLBackground(GenericTimeSeries): r""" Quicklook X-ray time series from the background detector. - Background monitoring is done in such way that counts from background detector are summed in five specified energy + Nominally QL background files contain counts from the background detector summed over five specified energy ranges. These QL data are double buffered into accumulators of 32bit depth. Maximum rate is approximately 30kHz and one live time counter is available. Integration time is parameter with default value of 32s @@ -335,7 +335,7 @@ def is_datasource_for(cls, **kwargs): Determines if the file corresponds to a STIX QL LightCurve `~sunpy.timeseries.TimeSeries`. """ - # Check if source is explicitly assigned + # Check if a source is explicitly assigned if "source" in kwargs.keys(): if kwargs.get("source", ""): return kwargs.get("source", "").lower().startswith(cls._source) @@ -496,7 +496,7 @@ def is_datasource_for(cls, **kwargs): Determines if the file corresponds to a STIX QL LightCurve `~sunpy.timeseries.TimeSeries`. """ - # Check if source is explicitly assigned + # Check if a source is explicitly assigned if "source" in kwargs.keys(): if kwargs.get("source", ""): return kwargs.get("source", "").lower().startswith(cls._source) @@ -596,20 +596,22 @@ def _parse_hdus(cls, hdulist): The path to the file you want to parse. """ header = hdulist[0].header - control = Table(hdulist[1].data) + control = _hdu_to_qtable(hdulist[1]) header["control"] = control - data = Table(hdulist[2].data) + data = _hdu_to_qtable(hdulist[2]) try: - data["time"] = Time(header["date_obs"]) + TimeDelta(data["time"] * u.cs) + data["time"] = Time(header["date_obs"]) + TimeDelta(data["time"]) except KeyError: - data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"] * u.cs) + data["time"] = Time(header["date-obs"]) + TimeDelta(data["time"]) + + units = OrderedDict((c.info.name, c.unit) for c in data.itercols() if c.info.name != 'time') data_df = data.to_pandas() data_df.index = data_df["time"] data_df.drop(columns="time", inplace=True) - return data_df, header, None + return data_df, header, units @classmethod def is_datasource_for(cls, **kwargs): @@ -617,7 +619,7 @@ def is_datasource_for(cls, **kwargs): Determines if the file corresponds to a STIX QL LightCurve `~sunpy.timeseries.TimeSeries`. """ - # Check if source is explicitly assigned + # Check if a source is explicitly assigned if "source" in kwargs.keys(): if kwargs.get("source", ""): return kwargs.get("source", "").lower().startswith(cls._source) diff --git a/stixpy/timeseries/tests/test_quicklook.py b/stixpy/timeseries/tests/test_quicklook.py index ea93b3c..3bb7905 100644 --- a/stixpy/timeseries/tests/test_quicklook.py +++ b/stixpy/timeseries/tests/test_quicklook.py @@ -13,6 +13,7 @@ def ql_lightcurve(): ) return ql_lc + @pytest.mark.remote_data @pytest.fixture() def ql_background(): @@ -21,6 +22,7 @@ def ql_background(): ) return ql_bkg + @pytest.mark.remote_data @pytest.fixture() def ql_variance(): @@ -29,11 +31,22 @@ def ql_variance(): ) return ql_var + +@pytest.mark.remote_data +@pytest.fixture() +def hk_maxi(): + hk_maxi = TimeSeries( + r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/HK/solo_L1_stix-hk-maxi_20200506_V02U.fits" + ) + return hk_maxi + + @pytest.mark.remote_data def test_ql_lightcurve(ql_lightcurve): assert isinstance(ql_lightcurve, QLLightCurve) assert ql_lightcurve.quantity('4-10 keV').unit == u.ct/(u.s * u.keV) + @pytest.mark.remote_data def test_ql_lightcurve_plot(ql_lightcurve): ql_lightcurve.plot() @@ -52,23 +65,25 @@ def test_ql_background_plot(ql_background): ql_background.plot(columns=['4-10 keV']) - @pytest.mark.remote_data -def test_qlvariance(ql_variance): +def test_ql_variance(ql_variance): assert isinstance(ql_variance, QLVariance) assert ql_variance.quantity('4-20 keV').unit == u.ct / (u.s * u.keV) @pytest.mark.remote_data -def test_qlvariance_plot(ql_variance): +def test_ql_variance_plot(ql_variance): ql_variance.plot() ql_variance.plot(columns=['4-20 keV']) - @pytest.mark.remote_data -def test_hk_maxi(): - hk_maxi = TimeSeries( - "https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/HK/solo_L1_stix-hk-maxi_20200506_V02U.fits" - ) +def test_hk_maxi(hk_maxi): assert isinstance(hk_maxi, HKMaxi) + assert hk_maxi.quantity('hk_att_c').unit == u.mA + + +@pytest.mark.remote_data +def test_hk_maxi_plot(hk_maxi): + hk_maxi.plot() + hk_maxi.plot(columns=['hk_asp_photoa0_v']) From a3aba9db600eca0cf8fd2cbf0a187ea9cbe36fdf Mon Sep 17 00:00:00 2001 From: Shane Maloney Date: Fri, 17 May 2024 12:31:17 +0100 Subject: [PATCH 13/13] windows test fail --- stixpy/timeseries/quicklook.py | 4 ++-- stixpy/timeseries/tests/test_quicklook.py | 27 +++++++++++++---------- 2 files changed, 17 insertions(+), 14 deletions(-) diff --git a/stixpy/timeseries/quicklook.py b/stixpy/timeseries/quicklook.py index 5ce0152..e9a86e4 100644 --- a/stixpy/timeseries/quicklook.py +++ b/stixpy/timeseries/quicklook.py @@ -596,8 +596,8 @@ def _parse_hdus(cls, hdulist): The path to the file you want to parse. """ header = hdulist[0].header - control = _hdu_to_qtable(hdulist[1]) - header["control"] = control + # control = _hdu_to_qtable(hdulist[1]) + # header["control"] = control data = _hdu_to_qtable(hdulist[2]) try: diff --git a/stixpy/timeseries/tests/test_quicklook.py b/stixpy/timeseries/tests/test_quicklook.py index 3bb7905..9aa6395 100644 --- a/stixpy/timeseries/tests/test_quicklook.py +++ b/stixpy/timeseries/tests/test_quicklook.py @@ -32,13 +32,13 @@ def ql_variance(): return ql_var -@pytest.mark.remote_data -@pytest.fixture() -def hk_maxi(): - hk_maxi = TimeSeries( - r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/HK/solo_L1_stix-hk-maxi_20200506_V02U.fits" - ) - return hk_maxi +# @pytest.mark.remote_data +# @pytest.fixture() +# def hk_maxi(): +# hk_maxi = TimeSeries( +# r"https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/HK/solo_L1_stix-hk-maxi_20200506_V02U.fits" +# ) +# return hk_maxi @pytest.mark.remote_data @@ -78,12 +78,15 @@ def test_ql_variance_plot(ql_variance): @pytest.mark.remote_data -def test_hk_maxi(hk_maxi): +def test_hk_maxi(): + hk_maxi = TimeSeries( + "https://pub099.cs.technik.fhnw.ch/fits/L1/2020/05/06/HK/solo_L1_stix-hk-maxi_20200506_V02U.fits" + ) assert isinstance(hk_maxi, HKMaxi) assert hk_maxi.quantity('hk_att_c').unit == u.mA -@pytest.mark.remote_data -def test_hk_maxi_plot(hk_maxi): - hk_maxi.plot() - hk_maxi.plot(columns=['hk_asp_photoa0_v']) +# @pytest.mark.remote_data +# def test_hk_maxi_plot(hk_maxi): +# hk_maxi.plot() +# hk_maxi.plot(columns=['hk_asp_photoa0_v'])