Skip to content

Commit 44a7611

Browse files
authored
make RhessiLoader compatible with legacy fitter args (sunpy#181)
* make RhessiLoader compatible with legacy fitter args * RHESSI lightcurves: display proper energy ranges * add time legacy properties to RHESSI loader * add back old bg subtract property * changelog * spectrum: ndcube print change (for docs to pass) * remove select_times method
1 parent 2125410 commit 44a7611

File tree

4 files changed

+85
-56
lines changed

4 files changed

+85
-56
lines changed

changelog/181.bugfix.rst

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1 @@
1+
Add back legacy parameters/functionality relating to `~sunkit_spex.extern.stix.RhessiLoader`.

examples/fitting_RHESSI_spectra.py

Lines changed: 0 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -116,20 +116,6 @@
116116
plt.show()
117117
plt.rcParams["font.size"] = default_text
118118

119-
#####################################################
120-
#
121-
# An alternative to set the event and background times like above is by using the `select_time` method. E.g.,
122-
#
123-
# .. code-block::python
124-
#
125-
# ress_spec.loaded_spec_data["spectrum1"].select_time(start="2002-10-05T10:38:32", end="2002-10-05T10:40:32", background=True)
126-
#
127-
# ress_spec.loaded_spec_data["spectrum1"].select_time(start="2002-10-05T10:41:20", end="2002-10-05T10:42:24")
128-
#
129-
# Both and end and a start time needs to be defined for the background, whereas the event time is assumed to
130-
# commence/finish at the first/final data time if the start/end time is not given.
131-
#
132-
133119
#####################################################
134120
#
135121
# Now let's get going with a model and explicitly stating a fit statistic

sunkit_spex/extern/rhessi.py

Lines changed: 83 additions & 41 deletions
Original file line numberDiff line numberDiff line change
@@ -28,16 +28,15 @@ class RhessiLoader(instruments.InstrumentBlueprint):
2828
[1] https://heasarc.gsfc.nasa.gov/xanadu/xspec/manual/XSappendixStatistics.html
2929
"""
3030

31-
def __init__(self, spectrum_fn, srm_fn, **kwargs):
31+
def __init__(self, spectrum_fn, srm_fn):
3232
"""
3333
Spectrum and SRM files are both required: attenuator state change times
3434
are in the spectrum file,
3535
and the state determines which SRM will be used.
3636
"""
37-
self._construction_string = f"RhessiLoader(spectrum_fn={spectrum_fn}, " f"srm_fn={srm_fn}," f"**{kwargs})"
37+
self._construction_string = f"RhessiLoader(spectrum_fn={spectrum_fn}, " f"srm_fn={srm_fn})"
3838
self._systematic_error = 0
3939
self.load_prepare_spectrum_srm(spectrum_fn, srm_fn)
40-
self._start_background_time, self._end_background_time = None, None
4140

4241
@property
4342
def systematic_error(self):
@@ -96,9 +95,11 @@ def load_prepare_spectrum_srm(self, spectrum_fn, srm_fn):
9695
}
9796
self._update_event_data_with_times()
9897
self._original_data = copy.deepcopy(self._loaded_spec_data)
98+
# Set these times to sensible defaults
99+
self.update_background_times(self._start_event_time, self._end_event_time)
99100

100101
@property
101-
def subtract_background(self):
102+
def subtract_background(self) -> bool:
102103
"""
103104
States whether the the data is event minus background or not.
104105
"""
@@ -411,10 +412,7 @@ def _rebin_ts(self, times, clump_bins):
411412
def lightcurve(self, energy_ranges=None, axes=None, rebin_time=1):
412413
"""Creates a RHESSI lightcurve.
413414
414-
Helps the user see the RHESSI time profile. The defined event time (defined either through `start_event_time`,
415-
`end_event_time` setters, or `select_time(...)` method) is shown with a purple shaded region and if a background
416-
time (defined either through `start_background_time`, `end_background_time` setters, or
417-
`select_time(...,background=True)` method) is defined then it is shown with an orange shaded region.
415+
The background and event times are highlighted based on what the user has set.
418416
419417
Parameters
420418
----------
@@ -445,10 +443,9 @@ def lightcurve(self, energy_ranges=None, axes=None, rebin_time=1):
445443
elif len(np.shape(energy_ranges)) == 1:
446444
energy_ranges = [energy_ranges]
447445
elif len(np.shape(energy_ranges)) == 0:
448-
print(
446+
raise ValueError(
449447
"The `energy_ranges` input should be a range of two energy values (e.g., [4,8] meaning 4-8 keV inclusive) or a list of these ranges."
450448
)
451-
return
452449

453450
ax = axes or plt.gca()
454451
font_size = plt.rcParams["font.size"]
@@ -462,11 +459,11 @@ def lightcurve(self, energy_ranges=None, axes=None, rebin_time=1):
462459

463460
# plot each energy range
464461
time_binning = (flat_times[1:] - flat_times[:-1]).to_value("s")
462+
left_edges, right_edges = self._loaded_spec_data["count_channel_bins"].T
465463
for er in energy_ranges:
466-
i = np.where(
467-
(self._loaded_spec_data["count_channel_bins"][:, 0] >= er[0])
468-
& (self._loaded_spec_data["count_channel_bins"][:, -1] <= er[-1])
469-
)
464+
# Restrict the energy range to the data energy bins
465+
ea, eb = max(er[0], left_edges[0]), min(er[1], right_edges[-1])
466+
i = np.where((left_edges >= ea) & (right_edges <= eb))
470467
e_range_cts = np.sum(self._spectrum["counts"][:, i].reshape((len(time_binning), -1)), axis=1)
471468

472469
if rebin_time > 1:
@@ -475,7 +472,7 @@ def lightcurve(self, energy_ranges=None, axes=None, rebin_time=1):
475472

476473
e_range_ctr, e_range_ctr_err = e_range_cts / time_binning, np.sqrt(e_range_cts) / time_binning
477474
lc = e_range_ctr
478-
p = ax.stairs(lc, flat_times.datetime, label=f"{er[0]}$-${er[-1]} keV")
475+
p = ax.stairs(lc, flat_times.datetime, label=f"{ea}$-${eb} keV")
479476
ax.errorbar(
480477
x=time_mids.datetime, y=e_range_ctr, yerr=e_range_ctr_err, c=p.get_edgecolor(), ls=""
481478
) # error bar in middle of the bin
@@ -515,10 +512,7 @@ def lightcurve(self, energy_ranges=None, axes=None, rebin_time=1):
515512
def spectrogram(self, axes=None, rebin_time=1, rebin_energy=1, **kwargs):
516513
"""Creates a RHESSI spectrogram.
517514
518-
Helps the user see the RHESSI time and energy evolution. The defined event time (defined either through
519-
`start_event_time`, `end_event_time` setters, or `select_time(...)` method) is shown with a violet line
520-
and if a background time (defined either through `start_background_time`, `end_background_time` setters,
521-
or `select_time(...,background=True)` method) is defined then it is shown with an orange line.
515+
The background and event times are highlighted based on what the user has set.
522516
523517
Parameters
524518
----------
@@ -604,32 +598,65 @@ def spectrogram(self, axes=None, rebin_time=1, rebin_energy=1, **kwargs):
604598

605599
return ax
606600

607-
def select_time(self, start=None, end=None, background=False):
608-
"""Provides method to set start and end time of the event or background in one line.
601+
# "retroactive" properties for setting the times like in older versions,
602+
# for legacy compatibility.
603+
@property
604+
def start_event_time(self):
605+
time_depr_warn('event')
606+
return self._start_event_time
609607

610-
Parameters
611-
----------
612-
start, end : str, `astropy.Time`, None
613-
String to be given to astropy's Time, `astropy.Time` is used directly, None sets the
614-
start/end event time to be the first time of the data. None doesn't add, or will remove,
615-
any background data in `_loaded_spec_data["extras"]` if background is set to True.
616-
Default: None
608+
@start_event_time.setter
609+
def start_event_time(self, new_time: atime.Time):
610+
time_depr_warn('event')
611+
self.update_event_times(new_time, self._end_event_time)
617612

618-
background : bool
619-
Determines whether the start and end times are for the event (False) or background
620-
time (True).
621-
Default: False
613+
@property
614+
def end_event_time(self):
615+
time_depr_warn('event')
616+
return self._end_event_time
622617

623-
Returns
624-
-------
625-
None.
626-
"""
627-
self.__warn = False if (start is not None) and (end is not None) else True
618+
@end_event_time.setter
619+
def end_event_time(self, new_time: atime.Time):
620+
time_depr_warn('event')
621+
self.update_event_times(self._start_event_time, new_time)
628622

629-
if background:
630-
self.start_background_time, self.end_background_time = start, end
631-
else:
632-
self.start_event_time, self.end_event_time = start, end
623+
@property
624+
def start_background_time(self):
625+
time_depr_warn('background')
626+
return self._start_background_time
627+
628+
@start_background_time.setter
629+
def start_background_time(self, new_time: atime.Time):
630+
time_depr_warn('background')
631+
self.update_background_times(new_time, self._end_background_time)
632+
633+
@property
634+
def end_background_time(self):
635+
time_depr_warn('background')
636+
return self._end_background_time
637+
638+
@end_background_time.setter
639+
def end_background_time(self, new_time: atime.Time):
640+
time_depr_warn('background')
641+
self.update_background_times(self._start_background_time, new_time)
642+
643+
# "retroactive" properties for specifying background subtraction,
644+
# like in older versions.
645+
@property
646+
def data2data_minus_background(self) -> bool:
647+
warnings.warn(
648+
"Please use the `subtract_background` property instead.",
649+
stacklevel=2
650+
)
651+
return self.subtract_background
652+
653+
@data2data_minus_background.setter
654+
def data2data_minus_background(self, do_subtract: bool):
655+
warnings.warn(
656+
"Please use the `subtract_background` property instead.",
657+
stacklevel=2
658+
)
659+
self.subtract_background = do_subtract
633660

634661

635662
def load_spectrum(spec_fn: str):
@@ -801,3 +828,18 @@ def load_srm(srm_file: str):
801828
}
802829

803830
return dict(channel_bins=channel_bins, photon_bins=photon_bins, srm_options=ret_srms)
831+
832+
833+
class LegacyRhessiLoader(RhessiLoader):
834+
'''Allow legacy `Fitter` to load in RHESSI data with the expected args format'''
835+
def __init__(self, spec_fn, **kw):
836+
srm_fn = kw.pop('srm_file')
837+
super().__init__(spec_fn, srm_fn)
838+
839+
840+
def time_depr_warn(partial_name: str):
841+
warnings.warn(
842+
"This time-related property is deprecated. "
843+
f"In the future, please use `RhessiLoader.update_{partial_name}_times`",
844+
stacklevel=2
845+
)

sunkit_spex/legacy/fitting/data_loader.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -136,7 +136,7 @@ def __init__(
136136
f"srm_file={srm_file},srm_custom={srm_custom},custom_channel_bins={custom_channel_bins}, **{kwargs})"
137137
)
138138

139-
self.instrument_loaders = {"NuSTAR": inst.NustarLoader, "RHESSI": rhessi.RhessiLoader, "Solar Orbiter": stix.STIXLoader}
139+
self.instrument_loaders = {"NuSTAR": inst.NustarLoader, "RHESSI": rhessi.LegacyRhessiLoader, "Solar Orbiter": stix.STIXLoader}
140140

141141
pha_file, arf_file, rmf_file, srm_file, srm_custom, custom_channel_bins, instruments = self._sort_files(
142142
pha_file=pha_file,

0 commit comments

Comments
 (0)