From bda741f8d348fff3d275845cd07d7ede2cf26676 Mon Sep 17 00:00:00 2001 From: ojustino Date: Fri, 30 Sep 2022 09:57:35 -0400 Subject: [PATCH 01/12] Better handled Spectrum1D images across classes --- specreduce/background.py | 14 +++++++++----- specreduce/extract.py | 41 ++++++++++++++++++++++++---------------- specreduce/tracing.py | 20 ++++++++++++++++---- 3 files changed, 50 insertions(+), 25 deletions(-) diff --git a/specreduce/background.py b/specreduce/background.py index 23d34018..328df0a9 100644 --- a/specreduce/background.py +++ b/specreduce/background.py @@ -27,7 +27,7 @@ class Background: Parameters ---------- - image : `~astropy.nddata.NDData` or array-like + image : `~astropy.nddata.NDData`-like or array-like image with 2-D spectral image data traces : List list of trace objects (or integers to define FlatTraces) to @@ -60,7 +60,7 @@ def __post_init__(self): Parameters ---------- - image : `~astropy.nddata.NDData` or array-like + image : `~astropy.nddata.NDData`-like or array-like image with 2-D spectral image data traces : List list of trace objects (or integers to define FlatTraces) to @@ -88,7 +88,6 @@ def _to_trace(trace): if self.width < 0: raise ValueError("width must be positive") - if self.width == 0: self.bkg_array = np.zeros(self.image.shape[self.disp_axis]) return @@ -96,6 +95,11 @@ def _to_trace(trace): if isinstance(self.traces, Trace): self.traces = [self.traces] + if isinstance(self.image, NDData): + # NOTE: should the NDData structure instead be preserved? + # (NDData includes Spectrum1D under its umbrella) + self.image = self.image.data + bkg_wimage = np.zeros_like(self.image, dtype=np.float64) for trace in self.traces: trace = _to_trace(trace) @@ -150,7 +154,7 @@ def two_sided(cls, image, trace_object, separation, **kwargs): Parameters ---------- - image : nddata-compatible image + image : `~astropy.nddata.NDData`-like or array-like image with 2-D spectral image data trace_object: Trace estimated trace of the spectrum to center the background traces @@ -183,7 +187,7 @@ def one_sided(cls, image, trace_object, separation, **kwargs): Parameters ---------- - image : nddata-compatible image + image : `~astropy.nddata.NDData`-like or array-like image with 2-D spectral image data trace_object: Trace estimated trace of the spectrum to center the background traces diff --git a/specreduce/extract.py b/specreduce/extract.py index 9965794f..51ce09a7 100644 --- a/specreduce/extract.py +++ b/specreduce/extract.py @@ -137,15 +137,15 @@ class BoxcarExtract(SpecreduceOperation): Parameters ---------- - image : nddata-compatible image + image : `~astropy.nddata.NDData`-like or array-like, required image with 2-D spectral image data - trace_object : Trace + trace_object : Trace, required trace object - width : float + width : float, optional width of extraction aperture in pixels - disp_axis : int + disp_axis : int, optional dispersion axis - crossdisp_axis : int + crossdisp_axis : int, optional cross-dispersion axis Returns @@ -171,15 +171,15 @@ def __call__(self, image=None, trace_object=None, width=None, Parameters ---------- - image : nddata-compatible image + image : `~astropy.nddata.NDData`-like or array-like, required image with 2-D spectral image data - trace_object : Trace + trace_object : Trace, required trace object - width : float + width : float, optional width of extraction aperture in pixels [default: 5] - disp_axis : int + disp_axis : int, optional dispersion axis [default: 1] - crossdisp_axis : int + crossdisp_axis : int, optional cross-dispersion axis [default: 0] @@ -195,6 +195,14 @@ def __call__(self, image=None, trace_object=None, width=None, disp_axis = disp_axis if disp_axis is not None else self.disp_axis crossdisp_axis = crossdisp_axis if crossdisp_axis is not None else self.crossdisp_axis + # handle image processing based on its type + if isinstance(image, Spectrum1D): + img = image.data + unit = image.unit + else: + img = image + unit = getattr(image, 'unit', u.DN) + # TODO: this check can be removed if/when implemented as a check in FlatTrace if isinstance(trace_object, FlatTrace): if trace_object.trace_pos < 1: @@ -204,16 +212,16 @@ def __call__(self, image=None, trace_object=None, width=None, raise ValueError("width must be positive") # weight image to use for extraction - wimage = _ap_weight_image( + wimg = _ap_weight_image( trace_object, width, disp_axis, crossdisp_axis, - image.shape) + img.shape) # extract - ext1d = np.sum(image * wimage, axis=crossdisp_axis) - return _to_spectrum1d_pixels(ext1d * getattr(image, 'unit', u.DN)) + ext1d = np.sum(img * wimg, axis=crossdisp_axis) * unit + return _to_spectrum1d_pixels(ext1d) @dataclass @@ -225,7 +233,7 @@ class HorneExtract(SpecreduceOperation): Parameters ---------- - image : `~astropy.nddata.NDData` or array-like, required + image : `~astropy.nddata.NDData`-like or array-like, required The input 2D spectrum from which to extract a source. An NDData object must specify uncertainty and a mask. An array requires use of the ``variance``, ``mask``, & ``unit`` arguments. @@ -288,7 +296,7 @@ def __call__(self, image=None, trace_object=None, Parameters ---------- - image : `~astropy.nddata.NDData` or array-like, required + image : `~astropy.nddata.NDData`-like or array-like, required The input 2D spectrum from which to extract a source. An NDData object must specify uncertainty and a mask. An array requires use of the ``variance``, ``mask``, & ``unit`` arguments. @@ -341,6 +349,7 @@ def __call__(self, image=None, trace_object=None, # handle image and associated data based on image's type if isinstance(image, NDData): + # (NDData includes Spectrum1D under its umbrella) img = np.ma.array(image.data, mask=image.mask) unit = image.unit if image.unit is not None else u.Unit() diff --git a/specreduce/tracing.py b/specreduce/tracing.py index d22d01e6..687beba8 100644 --- a/specreduce/tracing.py +++ b/specreduce/tracing.py @@ -5,9 +5,10 @@ import warnings from astropy.modeling import Model, fitting, models -from astropy.nddata import CCDData, NDData +from astropy.nddata import NDData from astropy.stats import gaussian_sigma_to_fwhm from astropy.utils.decorators import deprecated +from specutils import Spectrum1D import numpy as np __all__ = ['Trace', 'FlatTrace', 'ArrayTrace', 'FitTrace'] @@ -20,7 +21,7 @@ class Trace: Parameters ---------- - image : `~astropy.nddata.CCDData` + image : `~astropy.nddata.NDData`-like or array-like, required Image to be traced Properties @@ -28,7 +29,7 @@ class Trace: shape : tuple Shape of the array describing the trace """ - image: CCDData + image: NDData def __post_init__(self): self.trace_pos = self.image.shape[0] / 2 @@ -37,6 +38,11 @@ def __post_init__(self): def __getitem__(self, i): return self.trace[i] + def _parse_image(self): + if isinstance(self.image, Spectrum1D): + # NOTE: should the Spectrum1D structure instead be preserved? + self.image = self.image.data + @property def shape(self): return self.trace.shape @@ -95,6 +101,8 @@ class FlatTrace(Trace): trace_pos: float def __post_init__(self): + super()._parse_image() + self.set_position(self.trace_pos) def set_position(self, trace_pos): @@ -124,6 +132,8 @@ class ArrayTrace(Trace): trace: np.ndarray def __post_init__(self): + super()._parse_image() + nx = self.image.shape[1] nt = len(self.trace) if nt != nx: @@ -156,7 +166,7 @@ class FitTrace(Trace): Parameters ---------- - image : `~astropy.nddata.NDData` or array-like, required + image : `~astropy.nddata.NDData`-like or array-like, required The image over which to run the trace. Assumes cross-dispersion (spatial) direction is axis 0 and dispersion (wavelength) direction is axis 1. @@ -203,6 +213,8 @@ class FitTrace(Trace): _disp_axis = 1 def __post_init__(self): + super()._parse_image() + # handle multiple image types and mask uncaught invalid values if isinstance(self.image, NDData): img = np.ma.masked_invalid(np.ma.masked_array(self.image.data, From 4b73bc1ba307e539f511f634e563cdafdcc89dc4 Mon Sep 17 00:00:00 2001 From: ojustino Date: Tue, 8 Nov 2022 12:30:16 -0500 Subject: [PATCH 02/12] Added image parsers to all classes that take images Quantity arrays are now acceptable inputs. After parsing, the attribute is now of type instead of as before. --- specreduce/background.py | 52 +++++++++-- specreduce/extract.py | 190 +++++++++++++++++++++++++-------------- specreduce/tracing.py | 51 ++++++++--- 3 files changed, 208 insertions(+), 85 deletions(-) diff --git a/specreduce/background.py b/specreduce/background.py index 328df0a9..31884adc 100644 --- a/specreduce/background.py +++ b/specreduce/background.py @@ -4,8 +4,9 @@ from dataclasses import dataclass, field import numpy as np -from astropy.nddata import NDData +from astropy.nddata import NDData, VarianceUncertainty from astropy import units as u +from specutils import Spectrum1D from specreduce.extract import _ap_weight_image, _to_spectrum1d_pixels from specreduce.tracing import Trace, FlatTrace @@ -54,6 +55,41 @@ class Background: disp_axis: int = 1 crossdisp_axis: int = 0 + def _parse_image(self): + """ + Convert all accepted image types to a consistently formatted Spectrum1D. + """ + + if isinstance(self.image, np.ndarray): + img = self.image + elif isinstance(self.image, u.quantity.Quantity): + img = self.image.value + else: # NDData, including CCDData and Spectrum1D + img = self.image.data + + # mask and uncertainty are set as None when they aren't specified upon + # creating a Spectrum1D object, so we must check whether these + # attributes are absent *and* whether they are present but set as None + if getattr(self.image, 'mask', None) is not None: + mask = self.image.mask + else: + mask = np.ma.masked_invalid(img).mask + + if getattr(self.image, 'uncertainty', None) is not None: + uncertainty = self.image.uncertainty + else: + uncertainty = VarianceUncertainty(np.ones(img.shape)) + + unit = getattr(self.image, 'unit', u.Unit('DN')) # or u.Unit()? + + spectral_axis = getattr(self.image, 'spectral_axis', + (np.arange(img.shape[self.disp_axis]) + if hasattr(self, 'disp_axis') + else np.arange(img.shape[1])) * u.pix) + + self.image = Spectrum1D(img * unit, spectral_axis=spectral_axis, + uncertainty=uncertainty, mask=mask) + def __post_init__(self): """ Determine the background from an image for subtraction. @@ -86,6 +122,8 @@ def _to_trace(trace): raise ValueError('trace_object.trace_pos must be >= 1') return trace + self._parse_image() + if self.width < 0: raise ValueError("width must be positive") if self.width == 0: @@ -95,12 +133,7 @@ def _to_trace(trace): if isinstance(self.traces, Trace): self.traces = [self.traces] - if isinstance(self.image, NDData): - # NOTE: should the NDData structure instead be preserved? - # (NDData includes Spectrum1D under its umbrella) - self.image = self.image.data - - bkg_wimage = np.zeros_like(self.image, dtype=np.float64) + bkg_wimage = np.zeros_like(self.image.data, dtype=np.float64) for trace in self.traces: trace = _to_trace(trace) windows_max = trace.trace.data.max() + self.width/2 @@ -131,10 +164,11 @@ def _to_trace(trace): self.bkg_wimage = bkg_wimage if self.statistic == 'average': - self.bkg_array = np.average(self.image, weights=self.bkg_wimage, + self.bkg_array = np.average(self.image.data, + weights=self.bkg_wimage, axis=self.crossdisp_axis) elif self.statistic == 'median': - med_image = self.image.copy() + med_image = self.image.data.copy() med_image[np.where(self.bkg_wimage) == 0] = np.nan self.bkg_array = np.nanmedian(med_image, axis=self.crossdisp_axis) else: diff --git a/specreduce/extract.py b/specreduce/extract.py index 51ce09a7..69ed41b4 100644 --- a/specreduce/extract.py +++ b/specreduce/extract.py @@ -7,7 +7,7 @@ from astropy import units as u from astropy.modeling import Model, models, fitting -from astropy.nddata import NDData +from astropy.nddata import NDData, VarianceUncertainty from specreduce.core import SpecreduceOperation from specreduce.tracing import Trace, FlatTrace @@ -164,6 +164,41 @@ class BoxcarExtract(SpecreduceOperation): def spectrum(self): return self.__call__() + def _parse_image(self): + """ + Convert all accepted image types to a consistently formatted Spectrum1D. + """ + + if isinstance(self.image, np.ndarray): + img = self.image + elif isinstance(self.image, u.quantity.Quantity): + img = self.image.value + else: # NDData, including CCDData and Spectrum1D + img = self.image.data + + # mask and uncertainty are set as None when they aren't specified upon + # creating a Spectrum1D object, so we must check whether these + # attributes are absent *and* whether they are present but set as None + if getattr(self.image, 'mask', None) is not None: + mask = self.image.mask + else: + mask = np.ma.masked_invalid(img).mask + + if getattr(self.image, 'uncertainty', None) is not None: + uncertainty = self.image.uncertainty + else: + uncertainty = VarianceUncertainty(np.ones(img.shape)) + + unit = getattr(self.image, 'unit', u.Unit('DN')) # or u.Unit()? + + spectral_axis = getattr(self.image, 'spectral_axis', + (np.arange(img.shape[self.disp_axis]) + if hasattr(self, 'disp_axis') + else np.arange(img.shape[1])) * u.pix) + + self.image = Spectrum1D(img * unit, spectral_axis=spectral_axis, + uncertainty=uncertainty, mask=mask) + def __call__(self, image=None, trace_object=None, width=None, disp_axis=None, crossdisp_axis=None): """ @@ -285,6 +320,86 @@ class HorneExtract(SpecreduceOperation): def spectrum(self): return self.__call__() + def _parse_image(self, variance=None, mask=None, unit=None): + """ + Convert all accepted image types to a consistently formatted Spectrum1D. + Takes some extra arguments exactly as they come from self.__call__() to + handle cases where users specify them as arguments instead of as + attributes of their image object. + """ + + if isinstance(self.image, np.ndarray): + img = self.image + elif isinstance(self.image, u.quantity.Quantity): + img = self.image.value + else: # NDData, including CCDData and Spectrum1D + img = self.image.data + + # mask is set as None when not specified upon creating a Spectrum1D + # object, so we must check whether it is absent *and* whether it's + # present but set as None + if getattr(self.image, 'mask', None) is not None: + mask = self.image.mask + else: + mask = np.ma.masked_invalid(img).mask + + # Process uncertainties, converting to variances when able and throwing + # an error when uncertainties are missing or less easily converted + if (hasattr(self.image, 'uncertainty') + and self.image.uncertainty is not None): + if self.image.uncertainty.uncertainty_type == 'var': + variance = self.image.uncertainty.array + elif self.image.uncertainty.uncertainty_type == 'std': + warnings.warn("image NDData object's uncertainty " + "interpreted as standard deviation. if " + "incorrect, use VarianceUncertainty when " + "assigning image object's uncertainty.") + variance = self.image.uncertainty.array**2 + elif self.image.uncertainty.uncertainty_type == 'ivar': + variance = 1 / self.image.uncertainty.array + else: + # other options are InverseVariance and UnknownVariance + raise ValueError("image NDData object has unexpected " + "uncertainty type. instead, try " + "VarianceUncertainty or StdDevUncertainty.") + elif (hasattr(self.image, 'uncertainty') + and self.image.uncertainty is None): + # ignore variance arg to focus on updating NDData object + raise ValueError('image NDData object lacks uncertainty') + else: + if variance is None: + raise ValueError("if image is a numpy or Quantity array, a " + "variance must be specified. consider " + "wrapping it into one object by instead " + "passing an NDData image.") + elif self.image.shape != variance.shape: + raise ValueError("image and variance shapes must match") + + if np.any(variance < 0): + raise ValueError("variance must be fully positive") + if np.all(variance == 0): + # technically would result in infinities, but since they're all + # zeros, we can override ones to simulate an unweighted case + variance = np.ones_like(variance) + if np.any(variance == 0): + # exclude such elements by editing the input mask + mask[variance == 0] = True + # replace the variances to avoid a divide by zero warning + variance[variance == 0] = np.nan + + variance = VarianceUncertainty(variance) + + unit = getattr(self.image, 'unit', + u.Unit(self.unit) if self.unit is not None else u.Unit()) + + spectral_axis = getattr(self.image, 'spectral_axis', + (np.arange(img.shape[self.disp_axis]) + if hasattr(self, 'disp_axis') + else np.arange(img.shape[1])) * u.pix) + + self.image = Spectrum1D(img * unit, spectral_axis=spectral_axis, + uncertainty=variance, mask=mask) + def __call__(self, image=None, trace_object=None, disp_axis=None, crossdisp_axis=None, bkgrd_prof=None, @@ -347,71 +462,16 @@ def __call__(self, image=None, trace_object=None, mask = mask if mask is not None else self.mask unit = unit if unit is not None else self.unit - # handle image and associated data based on image's type - if isinstance(image, NDData): - # (NDData includes Spectrum1D under its umbrella) - img = np.ma.array(image.data, mask=image.mask) - unit = image.unit if image.unit is not None else u.Unit() - - if image.uncertainty is not None: - # prioritize NDData's uncertainty over variance argument - if image.uncertainty.uncertainty_type == 'var': - variance = image.uncertainty.array - elif image.uncertainty.uncertainty_type == 'std': - # NOTE: CCDData defaults uncertainties given as pure arrays - # to std and logs a warning saying so upon object creation. - # should we remind users again here? - warnings.warn("image NDData object's uncertainty " - "interpreted as standard deviation. if " - "incorrect, use VarianceUncertainty when " - "assigning image object's uncertainty.") - variance = image.uncertainty.array**2 - elif image.uncertainty.uncertainty_type == 'ivar': - variance = 1 / image.uncertainty.array - else: - # other options are InverseVariance and UnknownVariance - raise ValueError("image NDData object has unexpected " - "uncertainty type. instead, try " - "VarianceUncertainty or StdDevUncertainty.") - else: - # ignore variance arg to focus on updating NDData object - raise ValueError('image NDData object lacks uncertainty') - - else: - if variance is None: - raise ValueError('if image is a numpy array, a variance must ' - 'be specified. consider wrapping it into one ' - 'object by instead passing an NDData image.') - elif image.shape != variance.shape: - raise ValueError('image and variance shapes must match') - - # check optional arguments, filling them in if absent - if mask is None: - mask = np.ma.masked_invalid(image).mask - elif image.shape != mask.shape: - raise ValueError('image and mask shapes must match.') - - if isinstance(unit, str): - unit = u.Unit(unit) - else: - unit = unit if unit is not None else u.Unit() - - # create image - img = np.ma.array(image, mask=mask) - - if np.all(variance == 0): - # technically would result in infinities, but since they're all zeros - # we can just do the unweighted case by overriding with all ones - variance = np.ones_like(variance) + # parse image and replace optional arguments with updated values + self._parse_image(variance, mask, unit) + variance = self.image.uncertainty.array + unit = self.image.unit - if np.any(variance < 0): - raise ValueError("variance must be fully positive") - - if np.any(variance == 0): - # exclude these elements by editing the input mask - img.mask[variance == 0] = True - # replace the variances to avoid a divide by zero warning - variance[variance == 0] = np.nan + # mask any previously uncaught invalid values + or_mask = np.logical_or(mask, + np.ma.masked_invalid(self.image.data).mask) + img = np.ma.masked_array(self.image.data, or_mask) + mask = img.mask # co-add signal in each image column ncols = img.shape[crossdisp_axis] diff --git a/specreduce/tracing.py b/specreduce/tracing.py index 687beba8..7865bda2 100644 --- a/specreduce/tracing.py +++ b/specreduce/tracing.py @@ -5,9 +5,10 @@ import warnings from astropy.modeling import Model, fitting, models -from astropy.nddata import NDData +from astropy.nddata import NDData, VarianceUncertainty from astropy.stats import gaussian_sigma_to_fwhm from astropy.utils.decorators import deprecated +from astropy import units as u from specutils import Spectrum1D import numpy as np @@ -39,9 +40,39 @@ def __getitem__(self, i): return self.trace[i] def _parse_image(self): - if isinstance(self.image, Spectrum1D): - # NOTE: should the Spectrum1D structure instead be preserved? - self.image = self.image.data + """ + Convert all accepted image types to a consistently formatted Spectrum1D. + """ + + if isinstance(self.image, np.ndarray): + img = self.image + elif isinstance(self.image, u.quantity.Quantity): + img = self.image.value + else: # NDData, including CCDData and Spectrum1D + img = self.image.data + + # mask and uncertainty are set as None when they aren't specified upon + # creating a Spectrum1D object, so we must check whether these + # attributes are absent *and* whether they are present but set as None + if getattr(self.image, 'mask', None) is not None: + mask = self.image.mask + else: + mask = np.ma.masked_invalid(img).mask + + if getattr(self.image, 'uncertainty', None) is not None: + uncertainty = self.image.uncertainty + else: + uncertainty = VarianceUncertainty(np.ones(img.shape)) + + unit = getattr(self.image, 'unit', u.Unit('DN')) # or u.Unit()? + + spectral_axis = getattr(self.image, 'spectral_axis', + (np.arange(img.shape[self._disp_axis]) + if hasattr(self, '_disp_axis') + else np.arange(img.shape[1])) * u.pix) + + self.image = Spectrum1D(img * unit, spectral_axis=spectral_axis, + uncertainty=uncertainty, mask=mask) @property def shape(self): @@ -115,7 +146,7 @@ def set_position(self, trace_pos): Position of the trace """ self.trace_pos = trace_pos - self.trace = np.ones_like(self.image[0]) * self.trace_pos + self.trace = np.ones_like(self.image.data[0]) * self.trace_pos self._bound_trace() @@ -215,12 +246,10 @@ class FitTrace(Trace): def __post_init__(self): super()._parse_image() - # handle multiple image types and mask uncaught invalid values - if isinstance(self.image, NDData): - img = np.ma.masked_invalid(np.ma.masked_array(self.image.data, - mask=self.image.mask)) - else: - img = np.ma.masked_invalid(self.image) + # mask any previously uncaught invalid values + or_mask = np.logical_or(self.image.mask, + np.ma.masked_invalid(self.image.data).mask) + img = np.ma.masked_array(self.image.data, or_mask) # validate arguments valid_peak_methods = ('gaussian', 'centroid', 'max') From 3eb1527968c3e68923ac57fe2b7ea2e9559deb26 Mon Sep 17 00:00:00 2001 From: ojustino Date: Mon, 14 Nov 2022 23:34:47 -0500 Subject: [PATCH 03/12] Consolidate image parsing in _ImageParser mixin --- CHANGES.rst | 4 ++- specreduce/background.py | 42 +++------------------- specreduce/core.py | 75 +++++++++++++++++++++++++++++++++++++++- specreduce/extract.py | 50 +++------------------------ specreduce/tracing.py | 49 +++++--------------------- 5 files changed, 94 insertions(+), 126 deletions(-) diff --git a/CHANGES.rst b/CHANGES.rst index 31515802..4b14c559 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -13,7 +13,9 @@ API Changes - Renamed KosmosTrace as FitTrace, a conglomerate class for traces that are fit to images instead of predetermined [#128] - The default number of bins for FitTrace is now its associated image's number -of dispersion pixels instead of 20. Its default peak_method is now 'max'. [#128] +of dispersion pixels instead of 20. Its default peak_method is now 'max' [#128] +- All operations now accept Spectrum1D and Quantity-type images. All accepted +image types are now processed internally as Spectrum1D objects [#146] Bug Fixes ^^^^^^^^^ diff --git a/specreduce/background.py b/specreduce/background.py index 31884adc..aebecdb4 100644 --- a/specreduce/background.py +++ b/specreduce/background.py @@ -4,10 +4,11 @@ from dataclasses import dataclass, field import numpy as np -from astropy.nddata import NDData, VarianceUncertainty +from astropy.nddata import NDData from astropy import units as u from specutils import Spectrum1D +from specreduce.core import _ImageParser from specreduce.extract import _ap_weight_image, _to_spectrum1d_pixels from specreduce.tracing import Trace, FlatTrace @@ -15,7 +16,7 @@ @dataclass -class Background: +class Background(_ImageParser): """ Determine the background from an image for subtraction. @@ -55,41 +56,6 @@ class Background: disp_axis: int = 1 crossdisp_axis: int = 0 - def _parse_image(self): - """ - Convert all accepted image types to a consistently formatted Spectrum1D. - """ - - if isinstance(self.image, np.ndarray): - img = self.image - elif isinstance(self.image, u.quantity.Quantity): - img = self.image.value - else: # NDData, including CCDData and Spectrum1D - img = self.image.data - - # mask and uncertainty are set as None when they aren't specified upon - # creating a Spectrum1D object, so we must check whether these - # attributes are absent *and* whether they are present but set as None - if getattr(self.image, 'mask', None) is not None: - mask = self.image.mask - else: - mask = np.ma.masked_invalid(img).mask - - if getattr(self.image, 'uncertainty', None) is not None: - uncertainty = self.image.uncertainty - else: - uncertainty = VarianceUncertainty(np.ones(img.shape)) - - unit = getattr(self.image, 'unit', u.Unit('DN')) # or u.Unit()? - - spectral_axis = getattr(self.image, 'spectral_axis', - (np.arange(img.shape[self.disp_axis]) - if hasattr(self, 'disp_axis') - else np.arange(img.shape[1])) * u.pix) - - self.image = Spectrum1D(img * unit, spectral_axis=spectral_axis, - uncertainty=uncertainty, mask=mask) - def __post_init__(self): """ Determine the background from an image for subtraction. @@ -122,7 +88,7 @@ def _to_trace(trace): raise ValueError('trace_object.trace_pos must be >= 1') return trace - self._parse_image() + self.image = self._parse_image(self.image) if self.width < 0: raise ValueError("width must be positive") diff --git a/specreduce/core.py b/specreduce/core.py index de7baa34..4dad59aa 100644 --- a/specreduce/core.py +++ b/specreduce/core.py @@ -1,13 +1,86 @@ # Licensed under a 3-clause BSD style license - see LICENSE.rst import inspect +import numpy as np + +from astropy import units as u +from astropy.nddata import VarianceUncertainty from dataclasses import dataclass +from specutils import Spectrum1D __all__ = ['SpecreduceOperation'] +class _ImageParser: + """ + Coerces images from accepted formats to Spectrum1D objects for + internal use in specreduce's operation classes. + + Fills any and all of uncertainty, mask, units, and spectral axis + that are missing in the provided image with generic values. + Accepted image types are: + + - `~specutils.spectra.spectrum1d.Spectrum1D` (preferred) + - `~astropy.nddata.ccddata.CCDData` + - `~astropy.nddata.ndddata.NDDData` + - `~astropy.units.quantity.Quantity` + - `~numpy.ndarray` + """ + def _parse_image(self, image, disp_axis=1): + """ + Convert all accepted image types to a consistently formatted + Spectrum1D object. + + Parameters + ---------- + image : `~astropy.nddata.NDData`-like or array-like, required + The image to be parsed. If None, defaults to class' own + image attribute. + disp_axis : int, optional + The index of the image's dispersion axis. Should not be + changed until operations can handle variable image + orientations. [default: 1] + """ + + # would be nice to handle (cross)disp_axis consistently across + # operations (public attribute? private attribute? argument only?) so + # it can be called from self instead of via kwargs... + + if image is None: + # useful for Background's instance methods + return self.image + + if isinstance(image, np.ndarray): + img = image + elif isinstance(image, u.quantity.Quantity): + img = image.value + else: # NDData, including CCDData and Spectrum1D + img = image.data + + # mask and uncertainty are set as None when they aren't specified upon + # creating a Spectrum1D object, so we must check whether these + # attributes are absent *and* whether they are present but set as None + if getattr(image, 'mask', None) is not None: + mask = image.mask + else: + mask = np.ma.masked_invalid(img).mask + + if getattr(image, 'uncertainty', None) is not None: + uncertainty = image.uncertainty + else: + uncertainty = VarianceUncertainty(np.ones(img.shape)) + + unit = getattr(image, 'unit', u.Unit('DN')) # or u.Unit()? + + spectral_axis = getattr(image, 'spectral_axis', + np.arange(img.shape[disp_axis]) * u.pix) + + return Spectrum1D(img * unit, spectral_axis=spectral_axis, + uncertainty=uncertainty, mask=mask) + + @dataclass -class SpecreduceOperation: +class SpecreduceOperation(_ImageParser): """ An operation to perform as part of a spectroscopic reduction pipeline. diff --git a/specreduce/extract.py b/specreduce/extract.py index 69ed41b4..1aab8cb4 100644 --- a/specreduce/extract.py +++ b/specreduce/extract.py @@ -9,7 +9,7 @@ from astropy.modeling import Model, models, fitting from astropy.nddata import NDData, VarianceUncertainty -from specreduce.core import SpecreduceOperation +from specreduce.core import _ImageParser, SpecreduceOperation from specreduce.tracing import Trace, FlatTrace from specutils import Spectrum1D @@ -164,41 +164,6 @@ class BoxcarExtract(SpecreduceOperation): def spectrum(self): return self.__call__() - def _parse_image(self): - """ - Convert all accepted image types to a consistently formatted Spectrum1D. - """ - - if isinstance(self.image, np.ndarray): - img = self.image - elif isinstance(self.image, u.quantity.Quantity): - img = self.image.value - else: # NDData, including CCDData and Spectrum1D - img = self.image.data - - # mask and uncertainty are set as None when they aren't specified upon - # creating a Spectrum1D object, so we must check whether these - # attributes are absent *and* whether they are present but set as None - if getattr(self.image, 'mask', None) is not None: - mask = self.image.mask - else: - mask = np.ma.masked_invalid(img).mask - - if getattr(self.image, 'uncertainty', None) is not None: - uncertainty = self.image.uncertainty - else: - uncertainty = VarianceUncertainty(np.ones(img.shape)) - - unit = getattr(self.image, 'unit', u.Unit('DN')) # or u.Unit()? - - spectral_axis = getattr(self.image, 'spectral_axis', - (np.arange(img.shape[self.disp_axis]) - if hasattr(self, 'disp_axis') - else np.arange(img.shape[1])) * u.pix) - - self.image = Spectrum1D(img * unit, spectral_axis=spectral_axis, - uncertainty=uncertainty, mask=mask) - def __call__(self, image=None, trace_object=None, width=None, disp_axis=None, crossdisp_axis=None): """ @@ -231,12 +196,7 @@ def __call__(self, image=None, trace_object=None, width=None, crossdisp_axis = crossdisp_axis if crossdisp_axis is not None else self.crossdisp_axis # handle image processing based on its type - if isinstance(image, Spectrum1D): - img = image.data - unit = image.unit - else: - img = image - unit = getattr(image, 'unit', u.DN) + self.image = self._parse_image(image) # TODO: this check can be removed if/when implemented as a check in FlatTrace if isinstance(trace_object, FlatTrace): @@ -252,11 +212,11 @@ def __call__(self, image=None, trace_object=None, width=None, width, disp_axis, crossdisp_axis, - img.shape) + self.image.shape) # extract - ext1d = np.sum(img * wimg, axis=crossdisp_axis) * unit - return _to_spectrum1d_pixels(ext1d) + ext1d = np.sum(self.image.data * wimg, axis=crossdisp_axis) + return _to_spectrum1d_pixels(ext1d * self.image.unit) @dataclass diff --git a/specreduce/tracing.py b/specreduce/tracing.py index 7865bda2..27c6a448 100644 --- a/specreduce/tracing.py +++ b/specreduce/tracing.py @@ -12,6 +12,8 @@ from specutils import Spectrum1D import numpy as np +from specreduce.core import _ImageParser + __all__ = ['Trace', 'FlatTrace', 'ArrayTrace', 'FitTrace'] @@ -39,41 +41,6 @@ def __post_init__(self): def __getitem__(self, i): return self.trace[i] - def _parse_image(self): - """ - Convert all accepted image types to a consistently formatted Spectrum1D. - """ - - if isinstance(self.image, np.ndarray): - img = self.image - elif isinstance(self.image, u.quantity.Quantity): - img = self.image.value - else: # NDData, including CCDData and Spectrum1D - img = self.image.data - - # mask and uncertainty are set as None when they aren't specified upon - # creating a Spectrum1D object, so we must check whether these - # attributes are absent *and* whether they are present but set as None - if getattr(self.image, 'mask', None) is not None: - mask = self.image.mask - else: - mask = np.ma.masked_invalid(img).mask - - if getattr(self.image, 'uncertainty', None) is not None: - uncertainty = self.image.uncertainty - else: - uncertainty = VarianceUncertainty(np.ones(img.shape)) - - unit = getattr(self.image, 'unit', u.Unit('DN')) # or u.Unit()? - - spectral_axis = getattr(self.image, 'spectral_axis', - (np.arange(img.shape[self._disp_axis]) - if hasattr(self, '_disp_axis') - else np.arange(img.shape[1])) * u.pix) - - self.image = Spectrum1D(img * unit, spectral_axis=spectral_axis, - uncertainty=uncertainty, mask=mask) - @property def shape(self): return self.trace.shape @@ -116,7 +83,7 @@ def __sub__(self, delta): @dataclass -class FlatTrace(Trace): +class FlatTrace(Trace, _ImageParser): """ Trace that is constant along the axis being traced. @@ -132,7 +99,7 @@ class FlatTrace(Trace): trace_pos: float def __post_init__(self): - super()._parse_image() + self.image = self._parse_image(self.image) self.set_position(self.trace_pos) @@ -151,7 +118,7 @@ def set_position(self, trace_pos): @dataclass -class ArrayTrace(Trace): +class ArrayTrace(Trace, _ImageParser): """ Define a trace given an array of trace positions. @@ -163,7 +130,7 @@ class ArrayTrace(Trace): trace: np.ndarray def __post_init__(self): - super()._parse_image() + self.image = self._parse_image(self.image) nx = self.image.shape[1] nt = len(self.trace) @@ -180,7 +147,7 @@ def __post_init__(self): @dataclass -class FitTrace(Trace): +class FitTrace(Trace, _ImageParser): """ Trace the spectrum aperture in an image. @@ -244,7 +211,7 @@ class FitTrace(Trace): _disp_axis = 1 def __post_init__(self): - super()._parse_image() + self.image = self._parse_image(self.image) # mask any previously uncaught invalid values or_mask = np.logical_or(self.image.mask, From aca0ec70343937355ffef6e76918f9fa476b4dc0 Mon Sep 17 00:00:00 2001 From: ojustino Date: Mon, 14 Nov 2022 23:47:06 -0500 Subject: [PATCH 04/12] Fixed Spectrum1D parsing in Horne and its tests --- specreduce/extract.py | 103 +++++++++++++++++++++---------- specreduce/tests/test_extract.py | 56 +++++++++++------ 2 files changed, 106 insertions(+), 53 deletions(-) diff --git a/specreduce/extract.py b/specreduce/extract.py index 1aab8cb4..a3b1cf88 100644 --- a/specreduce/extract.py +++ b/specreduce/extract.py @@ -216,7 +216,8 @@ def __call__(self, image=None, trace_object=None, width=None, # extract ext1d = np.sum(self.image.data * wimg, axis=crossdisp_axis) - return _to_spectrum1d_pixels(ext1d * self.image.unit) + return Spectrum1D(ext1d * self.image.unit, + spectral_axis=self.image.spectral_axis) @dataclass @@ -280,50 +281,86 @@ class HorneExtract(SpecreduceOperation): def spectrum(self): return self.__call__() - def _parse_image(self, variance=None, mask=None, unit=None): + def _parse_image(self, image, + variance=None, mask=None, unit=None, disp_axis=1): """ - Convert all accepted image types to a consistently formatted Spectrum1D. - Takes some extra arguments exactly as they come from self.__call__() to - handle cases where users specify them as arguments instead of as - attributes of their image object. + Convert all accepted image types to a consistently formatted + Spectrum1D object. + + HorneExtract needs its own version of this method because it is + more stringent in its requirements for input images. The extra + arguments are needed to handle cases where these parameters were + specified as arguments and those where they came as attributes + of the image object. + + Parameters + ---------- + image : `~astropy.nddata.NDData`-like or array-like, required + The image to be parsed. If None, defaults to class' own + image attribute. + variance : `~numpy.ndarray`, optional + (Only used if ``image`` is not an NDData object.) + The associated variances for each pixel in the image. Must + have the same dimensions as ``image``. If all zeros, the variance + will be ignored and treated as all ones. If any zeros, those + elements will be excluded via masking. If any negative values, + an error will be raised. + mask : `~numpy.ndarray`, optional + (Only used if ``image`` is not an NDData object.) + Whether to mask each pixel in the image. Must have the same + dimensions as ``image``. If blank, all non-NaN pixels are + unmasked. + unit : `~astropy.units.Unit` or str, optional + (Only used if ``image`` is not an NDData object.) + The associated unit for the data in ``image``. If blank, + fluxes are interpreted as unitless. + disp_axis : int, optional + The index of the image's dispersion axis. Should not be + changed until operations can handle variable image + orientations. [default: 1] """ - if isinstance(self.image, np.ndarray): - img = self.image - elif isinstance(self.image, u.quantity.Quantity): - img = self.image.value + if isinstance(image, np.ndarray): + img = image + elif isinstance(image, u.quantity.Quantity): + img = image.value else: # NDData, including CCDData and Spectrum1D - img = self.image.data + img = image.data # mask is set as None when not specified upon creating a Spectrum1D # object, so we must check whether it is absent *and* whether it's # present but set as None - if getattr(self.image, 'mask', None) is not None: - mask = self.image.mask + if getattr(image, 'mask', None) is not None: + mask = image.mask + elif mask is not None: + pass else: mask = np.ma.masked_invalid(img).mask + if img.shape != mask.shape: + raise ValueError('image and mask shapes must match.') + # Process uncertainties, converting to variances when able and throwing # an error when uncertainties are missing or less easily converted - if (hasattr(self.image, 'uncertainty') - and self.image.uncertainty is not None): - if self.image.uncertainty.uncertainty_type == 'var': - variance = self.image.uncertainty.array - elif self.image.uncertainty.uncertainty_type == 'std': + if (hasattr(image, 'uncertainty') + and image.uncertainty is not None): + if image.uncertainty.uncertainty_type == 'var': + variance = image.uncertainty.array + elif image.uncertainty.uncertainty_type == 'std': warnings.warn("image NDData object's uncertainty " "interpreted as standard deviation. if " "incorrect, use VarianceUncertainty when " "assigning image object's uncertainty.") - variance = self.image.uncertainty.array**2 - elif self.image.uncertainty.uncertainty_type == 'ivar': - variance = 1 / self.image.uncertainty.array + variance = image.uncertainty.array**2 + elif image.uncertainty.uncertainty_type == 'ivar': + variance = 1 / image.uncertainty.array else: - # other options are InverseVariance and UnknownVariance + # other options are InverseUncertainty and UnknownUncertainty raise ValueError("image NDData object has unexpected " "uncertainty type. instead, try " "VarianceUncertainty or StdDevUncertainty.") - elif (hasattr(self.image, 'uncertainty') - and self.image.uncertainty is None): + elif (hasattr(image, 'uncertainty') + and image.uncertainty is None): # ignore variance arg to focus on updating NDData object raise ValueError('image NDData object lacks uncertainty') else: @@ -332,7 +369,7 @@ def _parse_image(self, variance=None, mask=None, unit=None): "variance must be specified. consider " "wrapping it into one object by instead " "passing an NDData image.") - elif self.image.shape != variance.shape: + elif image.shape != variance.shape: raise ValueError("image and variance shapes must match") if np.any(variance < 0): @@ -349,16 +386,14 @@ def _parse_image(self, variance=None, mask=None, unit=None): variance = VarianceUncertainty(variance) - unit = getattr(self.image, 'unit', - u.Unit(self.unit) if self.unit is not None else u.Unit()) + unit = getattr(image, 'unit', + u.Unit(unit) if unit is not None else u.Unit()) - spectral_axis = getattr(self.image, 'spectral_axis', - (np.arange(img.shape[self.disp_axis]) - if hasattr(self, 'disp_axis') - else np.arange(img.shape[1])) * u.pix) + spectral_axis = getattr(image, 'spectral_axis', + np.arange(img.shape[disp_axis]) * u.pix) - self.image = Spectrum1D(img * unit, spectral_axis=spectral_axis, - uncertainty=variance, mask=mask) + return Spectrum1D(img * unit, spectral_axis=spectral_axis, + uncertainty=variance, mask=mask) def __call__(self, image=None, trace_object=None, disp_axis=None, crossdisp_axis=None, @@ -423,7 +458,7 @@ def __call__(self, image=None, trace_object=None, unit = unit if unit is not None else self.unit # parse image and replace optional arguments with updated values - self._parse_image(variance, mask, unit) + self.image = self._parse_image(image, variance, mask, unit, disp_axis) variance = self.image.uncertainty.array unit = self.image.unit diff --git a/specreduce/tests/test_extract.py b/specreduce/tests/test_extract.py index d378bb35..5238045d 100644 --- a/specreduce/tests/test_extract.py +++ b/specreduce/tests/test_extract.py @@ -2,7 +2,7 @@ import pytest import astropy.units as u -from astropy.nddata import CCDData +from astropy.nddata import CCDData, VarianceUncertainty, UnknownUncertainty from specreduce.extract import BoxcarExtract, HorneExtract, OptimalExtract from specreduce.tracing import FlatTrace, ArrayTrace @@ -79,7 +79,7 @@ def test_boxcar_array_trace(): assert np.allclose(spectrum.flux.value, np.full_like(spectrum.flux.value, 75.)) -def test_horne_array_validation(): +def test_horne_image_validation(): # # Test HorneExtract scenarios specific to its use with an image of # type `~numpy.ndarray` (instead of the default `~astropy.nddata.NDData`). @@ -88,47 +88,65 @@ def test_horne_array_validation(): extract = OptimalExtract(image.data, trace) # equivalent to HorneExtract # an array-type image must come with a variance argument - with pytest.raises(ValueError, match=r'.*array.*variance.*specified.*'): + with pytest.raises(ValueError, match=r'.*variance must be specified.*'): ext = extract() + # an NDData-type image can't have an empty uncertainty attribute + with pytest.raises(ValueError, match=r'.*NDData object lacks uncertainty'): + ext = extract(image=image) + + # an NDData-type image's uncertainty must be of type VarianceUncertainty + # or type StdDevUncertainty + with pytest.raises(ValueError, match=r'.*unexpected uncertainty type.*'): + err = UnknownUncertainty(np.ones_like(image)) + image.uncertainty = err + ext = extract(image=image) + # an array-type image must have the same dimensions as its variance argument with pytest.raises(ValueError, match=r'.*shapes must match.*'): + # remember variance, mask, and unit args are only checked if image + # object doesn't have those attributes (e.g., numpy and Quantity arrays) err = np.ones_like(image[0]) - ext = extract(variance=err) + ext = extract(image=image.data, variance=err) # an array-type image must have the same dimensions as its mask argument with pytest.raises(ValueError, match=r'.*shapes must match.*'): err = np.ones_like(image) mask = np.zeros_like(image[0]) - ext = extract(variance=err, mask=mask) + ext = extract(image=image.data, variance=err, mask=mask) # an array-type image given without mask and unit arguments is fine - # and produces a unitless result + # and produces an extraction with unitless flux and spectral axis in pixels err = np.ones_like(image) - ext = extract(variance=err) + ext = extract(image=image.data, variance=err, mask=None, unit=None) assert ext.unit == u.Unit() + assert np.all(ext.spectral_axis + == np.arange(image.shape[extract.disp_axis]) * u.pix) def test_horne_variance_errors(): trace = FlatTrace(image, 3.0) - # all zeros are treated as non-weighted (give non-zero fluxes) - err = np.zeros_like(image) - mask = np.zeros_like(image) - extract = HorneExtract(image.data, trace, variance=err, mask=mask, unit=u.Jy) + # all zeros are treated as non-weighted (i.e., as non-zero fluxes) + image.uncertainty = VarianceUncertainty(np.zeros_like(image)) + image.mask = np.zeros_like(image) + extract = HorneExtract(image, trace) ext = extract.spectrum assert not np.all(ext == 0) - # single zero value adjusts mask (does not raise error) + # single zero value adjusts mask and does not raise error err = np.ones_like(image) - err[0] = 0 - mask = np.zeros_like(image) - ext = extract(variance=err, mask=mask, unit=u.Jy) - assert not np.all(ext == 0) + err[0][0] = 0 + image.uncertainty = VarianceUncertainty(err) + ext = extract(image) + assert not np.all(ext == 1) # single negative value raises error - err = np.ones_like(image) - err[0] = -1 + err = image.uncertainty.array + err[0][0] = -1 mask = np.zeros_like(image) with pytest.raises(ValueError, match='variance must be fully positive'): - ext = extract(variance=err, mask=mask, unit=u.Jy) + # remember variance, mask, and unit args are only checked if image + # object doesn't have those attributes (e.g., numpy and Quantity arrays) + ext = extract(image=image.data, variance=err, + mask=image.mask, unit=u.Jy) From 98fff7fc2432f9e511613af93f5ef69eda7060ea Mon Sep 17 00:00:00 2001 From: ojustino Date: Mon, 14 Nov 2022 23:58:12 -0500 Subject: [PATCH 05/12] Fixed Spectrum1D parsing in Background and tests --- setup.cfg | 2 +- specreduce/background.py | 75 ++++++++++++++++++----------- specreduce/tests/test_background.py | 11 +++-- 3 files changed, 54 insertions(+), 34 deletions(-) diff --git a/setup.cfg b/setup.cfg index 756e81a3..b7d37a2a 100644 --- a/setup.cfg +++ b/setup.cfg @@ -18,7 +18,7 @@ python_requires = >=3.7 setup_requires = setuptools_scm install_requires = astropy - specutils + specutils>=1.7 synphot matplotlib photutils diff --git a/specreduce/background.py b/specreduce/background.py index aebecdb4..5c91f642 100644 --- a/specreduce/background.py +++ b/specreduce/background.py @@ -5,11 +5,12 @@ import numpy as np from astropy.nddata import NDData +from astropy.units import UnitTypeError from astropy import units as u from specutils import Spectrum1D from specreduce.core import _ImageParser -from specreduce.extract import _ap_weight_image, _to_spectrum1d_pixels +from specreduce.extract import _ap_weight_image from specreduce.tracing import Trace, FlatTrace __all__ = ['Background'] @@ -155,8 +156,10 @@ def two_sided(cls, image, trace_object, separation, **kwargs): Parameters ---------- image : `~astropy.nddata.NDData`-like or array-like - image with 2-D spectral image data - trace_object: Trace + Image with 2-D spectral image data. Assumes cross-dispersion + (spatial) direction is axis 0 and dispersion (wavelength) + direction is axis 1. + trace_object: `~specreduce.tracing.Trace` estimated trace of the spectrum to center the background traces separation: float separation from ``trace_object`` for the background regions @@ -171,6 +174,7 @@ def two_sided(cls, image, trace_object, separation, **kwargs): crossdisp_axis : int cross-dispersion axis """ + image = cls._parse_image(cls, image) kwargs['traces'] = [trace_object-separation, trace_object+separation] return cls(image=image, **kwargs) @@ -188,8 +192,10 @@ def one_sided(cls, image, trace_object, separation, **kwargs): Parameters ---------- image : `~astropy.nddata.NDData`-like or array-like - image with 2-D spectral image data - trace_object: Trace + Image with 2-D spectral image data. Assumes cross-dispersion + (spatial) direction is axis 0 and dispersion (wavelength) + direction is axis 1. + trace_object: `~specreduce.tracing.Trace` estimated trace of the spectrum to center the background traces separation: float separation from ``trace_object`` for the background, positive will be @@ -205,6 +211,7 @@ def one_sided(cls, image, trace_object, separation, **kwargs): crossdisp_axis : int cross-dispersion axis """ + image = cls._parse_image(cls, image) kwargs['traces'] = [trace_object+separation] return cls(image=image, **kwargs) @@ -214,18 +221,20 @@ def bkg_image(self, image=None): Parameters ---------- - image : nddata-compatible image or None - image with 2-D spectral image data. If None, will extract - the background from ``image`` used to initialize the class. + image : `~astropy.nddata.NDData`-like or array-like, optional + Image with 2-D spectral image data. Assumes cross-dispersion + (spatial) direction is axis 0 and dispersion (wavelength) + direction is axis 1. If None, will extract the background + from ``image`` used to initialize the class. [default: None] Returns ------- - array with same shape as ``image``. + Spectrum1D object with same shape as ``image``. """ - if image is None: - image = self.image - - return np.tile(self.bkg_array, (image.shape[0], 1)) + image = self._parse_image(image) + return Spectrum1D(np.tile(self.bkg_array, + (image.shape[0], 1)) * image.unit, + spectral_axis=image.spectral_axis) def bkg_spectrum(self, image=None): """ @@ -233,9 +242,11 @@ def bkg_spectrum(self, image=None): Parameters ---------- - image : nddata-compatible image or None - image with 2-D spectral image data. If None, will extract - the background from ``image`` used to initialize the class. + image : `~astropy.nddata.NDData`-like or array-like, optional + Image with 2-D spectral image data. Assumes cross-dispersion + (spatial) direction is axis 0 and dispersion (wavelength) + direction is axis 1. If None, will extract the background + from ``image`` used to initialize the class. [default: None] Returns ------- @@ -244,10 +255,15 @@ def bkg_spectrum(self, image=None): units as the input image (or u.DN if none were provided) and the spectral axis expressed in pixel units. """ - bkg_image = self.bkg_image(image=image) + bkg_image = self.bkg_image(image) - ext1d = np.sum(bkg_image, axis=self.crossdisp_axis) - return _to_spectrum1d_pixels(ext1d * getattr(image, 'unit', u.DN)) + try: + return bkg_image.collapse(np.sum, axis=self.crossdisp_axis) + except UnitTypeError: + # can't collapse with a spectral axis in pixels because + # SpectralCoord only allows frequency/wavelength equivalent units... + ext1d = np.sum(bkg_image.flux, axis=self.crossdisp_axis) + return Spectrum1D(ext1d, bkg_image.spectral_axis) def sub_image(self, image=None): """ @@ -263,14 +279,12 @@ def sub_image(self, image=None): ------- array with same shape as ``image`` """ - if image is None: - image = self.image + image = self._parse_image(image) - if isinstance(image, NDData): - # https://docs.astropy.org/en/stable/nddata/mixins/ndarithmetic.html - return image.subtract(self.bkg_image(image)*image.unit) - else: - return image - self.bkg_image(image) + # https://docs.astropy.org/en/stable/nddata/mixins/ndarithmetic.html + # (compare_wcs argument needed to avoid TypeError from SpectralCoord + # when image's spectral axis is in pixels) + return image.subtract(self.bkg_image(image), compare_wcs=None) def sub_spectrum(self, image=None): """ @@ -291,8 +305,13 @@ def sub_spectrum(self, image=None): """ sub_image = self.sub_image(image=image) - ext1d = np.sum(sub_image, axis=self.crossdisp_axis) - return _to_spectrum1d_pixels(ext1d * getattr(image, 'unit', u.DN)) + try: + return sub_image.collapse(np.sum, axis=self.crossdisp_axis) + except UnitTypeError: + # can't collapse with a spectral axis in pixels because + # SpectralCoord only allows frequency/wavelength equivalent units... + ext1d = np.sum(sub_image.flux, axis=self.crossdisp_axis) + return Spectrum1D(ext1d, spectral_axis=sub_image.spectral_axis) def __rsub__(self, image): """ diff --git a/specreduce/tests/test_background.py b/specreduce/tests/test_background.py index 53d13892..f0e98099 100644 --- a/specreduce/tests/test_background.py +++ b/specreduce/tests/test_background.py @@ -2,7 +2,7 @@ import numpy as np import astropy.units as u -from astropy.nddata import CCDData +from astropy.nddata import CCDData, VarianceUncertainty from specutils import Spectrum1D from specreduce.background import Background @@ -16,7 +16,8 @@ image = np.ones(shape=(30, 10)) for j in range(image.shape[0]): image[j, ::] *= j -image = CCDData(image, unit=u.Jy) +image = Spectrum1D(image * u.DN, + uncertainty=VarianceUncertainty(np.ones_like(image))) def test_background(): @@ -45,8 +46,8 @@ def test_background(): sub1 = image - bg1 sub2 = bg1.sub_image(image) sub3 = bg1.sub_image() - assert np.allclose(sub1, sub2) - assert np.allclose(sub1, sub3) + assert np.allclose(sub1.flux, sub2.flux) + assert np.allclose(sub1.flux, sub3.flux) bkg_spec = bg1.bkg_spectrum() assert isinstance(bkg_spec, Spectrum1D) @@ -54,7 +55,7 @@ def test_background(): assert isinstance(sub_spec, Spectrum1D) # test that width==0 results in no background bg = Background.two_sided(image, trace, bkg_sep, width=0) - assert np.all(bg.bkg_image() == 0) + assert np.all(bg.bkg_image().flux == 0) def test_warnings_errors(): From 68a5766460393f79ce8a08a2728716e785fcd572 Mon Sep 17 00:00:00 2001 From: ojustino Date: Tue, 15 Nov 2022 00:02:00 -0500 Subject: [PATCH 06/12] Drop _to_spectrum1d_pixels to allow other axes --- specreduce/extract.py | 9 ++------- 1 file changed, 2 insertions(+), 7 deletions(-) diff --git a/specreduce/extract.py b/specreduce/extract.py index a3b1cf88..a4bd9a9c 100644 --- a/specreduce/extract.py +++ b/specreduce/extract.py @@ -117,12 +117,6 @@ def _ap_weight_image(trace, width, disp_axis, crossdisp_axis, image_shape): return wimage -def _to_spectrum1d_pixels(fluxes): - # TODO: add wavelength units, uncertainty and mask to spectrum1D object - return Spectrum1D(spectral_axis=np.arange(len(fluxes)) * u.pixel, - flux=fluxes) - - @dataclass class BoxcarExtract(SpecreduceOperation): """ @@ -515,7 +509,8 @@ def __call__(self, image=None, trace_object=None, extraction = result * norms # convert the extraction to a Spectrum1D object - return _to_spectrum1d_pixels(extraction * unit) + return Spectrum1D(extraction * unit, + spectral_axis=self.image.spectral_axis) @dataclass From 0ec792d35708d6ba50614587ae9edd6beabc9b53 Mon Sep 17 00:00:00 2001 From: ojustino Date: Tue, 15 Nov 2022 00:16:00 -0500 Subject: [PATCH 07/12] Wrote tests for image parsing methods + codestyle --- specreduce/background.py | 1 - specreduce/extract.py | 2 +- specreduce/tests/test_background.py | 10 ++- specreduce/tests/test_extract.py | 1 - specreduce/tests/test_image_parsing.py | 104 +++++++++++++++++++++++++ specreduce/tracing.py | 4 +- 6 files changed, 112 insertions(+), 10 deletions(-) create mode 100644 specreduce/tests/test_image_parsing.py diff --git a/specreduce/background.py b/specreduce/background.py index 5c91f642..b1887db9 100644 --- a/specreduce/background.py +++ b/specreduce/background.py @@ -6,7 +6,6 @@ import numpy as np from astropy.nddata import NDData from astropy.units import UnitTypeError -from astropy import units as u from specutils import Spectrum1D from specreduce.core import _ImageParser diff --git a/specreduce/extract.py b/specreduce/extract.py index a4bd9a9c..64dfc0a7 100644 --- a/specreduce/extract.py +++ b/specreduce/extract.py @@ -9,7 +9,7 @@ from astropy.modeling import Model, models, fitting from astropy.nddata import NDData, VarianceUncertainty -from specreduce.core import _ImageParser, SpecreduceOperation +from specreduce.core import SpecreduceOperation from specreduce.tracing import Trace, FlatTrace from specutils import Spectrum1D diff --git a/specreduce/tests/test_background.py b/specreduce/tests/test_background.py index f0e98099..e85fee69 100644 --- a/specreduce/tests/test_background.py +++ b/specreduce/tests/test_background.py @@ -2,7 +2,7 @@ import numpy as np import astropy.units as u -from astropy.nddata import CCDData, VarianceUncertainty +from astropy.nddata import VarianceUncertainty from specutils import Spectrum1D from specreduce.background import Background @@ -43,11 +43,13 @@ def test_background(): bg = Background(image, trace, width=bkg_width) # test that image subtraction works - sub1 = image - bg1 + # NOTE: uncomment sub1 test once Spectrum1D and Background subtraction works + # (meaning specutils PR #988 is merged, released, and pinned here) + # sub1 = image - bg1 sub2 = bg1.sub_image(image) sub3 = bg1.sub_image() - assert np.allclose(sub1.flux, sub2.flux) - assert np.allclose(sub1.flux, sub3.flux) + # assert np.allclose(sub1.flux, sub2.flux) + assert np.allclose(sub2.flux, sub3.flux) bkg_spec = bg1.bkg_spectrum() assert isinstance(bkg_spec, Spectrum1D) diff --git a/specreduce/tests/test_extract.py b/specreduce/tests/test_extract.py index 5238045d..17552922 100644 --- a/specreduce/tests/test_extract.py +++ b/specreduce/tests/test_extract.py @@ -144,7 +144,6 @@ def test_horne_variance_errors(): # single negative value raises error err = image.uncertainty.array err[0][0] = -1 - mask = np.zeros_like(image) with pytest.raises(ValueError, match='variance must be fully positive'): # remember variance, mask, and unit args are only checked if image # object doesn't have those attributes (e.g., numpy and Quantity arrays) diff --git a/specreduce/tests/test_image_parsing.py b/specreduce/tests/test_image_parsing.py new file mode 100644 index 00000000..9cbb89fc --- /dev/null +++ b/specreduce/tests/test_image_parsing.py @@ -0,0 +1,104 @@ +import numpy as np + +from astropy import units as u +from astropy.io import fits +from astropy.nddata import CCDData, NDData, VarianceUncertainty +from astropy.utils.data import download_file + +from specreduce.extract import HorneExtract +from specreduce.tracing import FlatTrace +from specutils import Spectrum1D, SpectralAxis + +# fetch test image +fn = download_file('https://stsci.box.com/shared/static/exnkul627fcuhy5akf2gswytud5tazmw.fits', + cache=True) + +# duplicate image in all accepted formats +# (one Spectrum1D variant has a physical spectral axis; the other is in pixels) +img = fits.getdata(fn).T +flux = img * u.MJy / u.sr +sax = SpectralAxis(np.linspace(14.377, 3.677, flux.shape[-1]) * u.um) +unc = VarianceUncertainty(np.random.rand(*flux.shape)) + +all_images = {} +all_images['arr'] = img +all_images['s1d'] = Spectrum1D(flux, spectral_axis=sax, uncertainty=unc) +all_images['s1d_pix'] = Spectrum1D(flux, uncertainty=unc) +all_images['ccd'] = CCDData(img, uncertainty=unc, unit=flux.unit) +all_images['ndd'] = NDData(img, uncertainty=unc, unit=flux.unit) +all_images['qnt'] = img * flux.unit + +# save default values used for spectral axis and uncertainty when they are not +# available from the image object or provided by the user +sax_def = np.arange(img.shape[1]) * u.pix +unc_def = np.ones_like(img) + + +# (for use inside tests) +def compare_images(key, collection, compare='s1d'): + # was input converted to Spectrum1D? + assert isinstance(collection[key], Spectrum1D), (f"image '{key}' not " + "of type Spectrum1D") + + # do key's fluxes match its comparison's fluxes? + assert np.allclose(collection[key].data, + collection[compare].data), (f"images '{key}' and " + f"'{compare}' have unequal " + "flux values") + + # if the image came with a spectral axis, was it kept? if not, was the + # default spectral axis in pixels applied? + sax_provided = hasattr(all_images[key], 'spectral_axis') + assert np.allclose(collection[key].spectral_axis, + (all_images[key].spectral_axis if sax_provided + else sax_def)), (f"spectral axis of image '{key}' does " + f"not match {'input' if sax_provided else 'default'}") + + # if the image came with an uncertainty, was it kept? if not, was the + # default uncertainty created? + unc_provided = hasattr(all_images[key], 'uncertainty') + assert np.allclose(collection[key].uncertainty.array, + (all_images[key].uncertainty.array if unc_provided + else unc_def)), (f"uncertainty of image '{key}' does " + f"not match {'input' if unc_provided else 'default'}") + + # were masks created despite none being given? (all indices should be False) + assert (getattr(collection[key], 'mask', None) + is not None), f"no mask was created for image '{key}'" + assert np.all(collection[key].mask == 0), ("mask not all False " + f"for image '{key}'") + + +# test consistency of general image parser results +def test_parse_general(): + all_images_parsed = {k: FlatTrace._parse_image(object, im) + for k, im in all_images.items()} + + for key in all_images_parsed.keys(): + compare_images(key, all_images_parsed) + + +# use verified general image parser results to check HorneExtract's image parser +def test_parse_horne(): + # HorneExtract's parser is more stringent than the general one, hence the + # separate test. Given proper inputs, both should produce the same results. + images_collection = {k: {} for k in all_images.keys()} + + for key, col in images_collection.items(): + img = all_images[key] + col['general'] = FlatTrace._parse_image(object, img) + + if hasattr(all_images[key], 'uncertainty'): + defaults = {} + else: + # save default values of attributes used in general parser when + # they are not available from the image object. HorneExtract always + # requires a variance, so it's chosen here to be on equal footing + # with the general case + defaults = {'variance': unc_def, + 'mask': np.ma.masked_invalid(img).mask, + 'unit': getattr(img, 'unit', u.DN)} + + col[key] = HorneExtract._parse_image(object, img, **defaults) + + compare_images(key, col, compare='general') diff --git a/specreduce/tracing.py b/specreduce/tracing.py index 27c6a448..6a4b8d12 100644 --- a/specreduce/tracing.py +++ b/specreduce/tracing.py @@ -5,11 +5,9 @@ import warnings from astropy.modeling import Model, fitting, models -from astropy.nddata import NDData, VarianceUncertainty +from astropy.nddata import NDData from astropy.stats import gaussian_sigma_to_fwhm from astropy.utils.decorators import deprecated -from astropy import units as u -from specutils import Spectrum1D import numpy as np from specreduce.core import _ImageParser From 4b0bdd54cb8ccfaa8c7ca9a277b4e834d3988126 Mon Sep 17 00:00:00 2001 From: ojustino Date: Wed, 16 Nov 2022 13:33:37 -0500 Subject: [PATCH 08/12] Fix Background.sub_image for physical wavelengths --- specreduce/background.py | 18 ++++++++++++------ specreduce/tests/test_background.py | 22 +++++++++++++++++++++- 2 files changed, 33 insertions(+), 7 deletions(-) diff --git a/specreduce/background.py b/specreduce/background.py index b1887db9..b6ca2586 100644 --- a/specreduce/background.py +++ b/specreduce/background.py @@ -5,7 +5,7 @@ import numpy as np from astropy.nddata import NDData -from astropy.units import UnitTypeError +from astropy import units as u from specutils import Spectrum1D from specreduce.core import _ImageParser @@ -258,7 +258,7 @@ def bkg_spectrum(self, image=None): try: return bkg_image.collapse(np.sum, axis=self.crossdisp_axis) - except UnitTypeError: + except u.UnitTypeError: # can't collapse with a spectral axis in pixels because # SpectralCoord only allows frequency/wavelength equivalent units... ext1d = np.sum(bkg_image.flux, axis=self.crossdisp_axis) @@ -280,10 +280,14 @@ def sub_image(self, image=None): """ image = self._parse_image(image) + # a compare_wcs argument is needed for Spectrum1D.subtract() in order to + # avoid a TypeError from SpectralCoord when image's spectral axis is in + # pixels. it is not needed when image's spectral axis has physical units + kwargs = ({'compare_wcs': None} if image.spectral_axis.unit == u.pix + else {}) + # https://docs.astropy.org/en/stable/nddata/mixins/ndarithmetic.html - # (compare_wcs argument needed to avoid TypeError from SpectralCoord - # when image's spectral axis is in pixels) - return image.subtract(self.bkg_image(image), compare_wcs=None) + return image.subtract(self.bkg_image(image), **kwargs) def sub_spectrum(self, image=None): """ @@ -306,7 +310,7 @@ def sub_spectrum(self, image=None): try: return sub_image.collapse(np.sum, axis=self.crossdisp_axis) - except UnitTypeError: + except u.UnitTypeError: # can't collapse with a spectral axis in pixels because # SpectralCoord only allows frequency/wavelength equivalent units... ext1d = np.sum(sub_image.flux, axis=self.crossdisp_axis) @@ -316,4 +320,6 @@ def __rsub__(self, image): """ Subtract the background from an image. """ + # NOTE: will not be called until specutils PR #988 is merged, released, + # and pinned here return self.sub_image(image) diff --git a/specreduce/tests/test_background.py b/specreduce/tests/test_background.py index e85fee69..f36cec71 100644 --- a/specreduce/tests/test_background.py +++ b/specreduce/tests/test_background.py @@ -18,6 +18,9 @@ image[j, ::] *= j image = Spectrum1D(image * u.DN, uncertainty=VarianceUncertainty(np.ones_like(image))) +image_um = Spectrum1D(image.flux, + spectral_axis=np.arange(image.data.shape[1]) * u.um, + uncertainty=VarianceUncertainty(np.ones_like(image.data))) def test_background(): @@ -29,13 +32,22 @@ def test_background(): trace = FlatTrace(image, trace_pos) bkg_sep = 5 bkg_width = 2 - # all the following should be equivalent: + + # all the following should be equivalent, whether image's spectral axis + # is in pixels or physical units: bg1 = Background(image, [trace-bkg_sep, trace+bkg_sep], width=bkg_width) bg2 = Background.two_sided(image, trace, bkg_sep, width=bkg_width) bg3 = Background.two_sided(image, trace_pos, bkg_sep, width=bkg_width) assert np.allclose(bg1.bkg_array, bg2.bkg_array) assert np.allclose(bg1.bkg_array, bg3.bkg_array) + bg4 = Background(image_um, [trace-bkg_sep, trace+bkg_sep], width=bkg_width) + bg5 = Background.two_sided(image_um, trace, bkg_sep, width=bkg_width) + bg6 = Background.two_sided(image_um, trace_pos, bkg_sep, width=bkg_width) + assert np.allclose(bg1.bkg_array, bg4.bkg_array) + assert np.allclose(bg1.bkg_array, bg5.bkg_array) + assert np.allclose(bg1.bkg_array, bg6.bkg_array) + # test that creating a one_sided background works Background.one_sided(image, trace, bkg_sep, width=bkg_width) @@ -51,6 +63,14 @@ def test_background(): # assert np.allclose(sub1.flux, sub2.flux) assert np.allclose(sub2.flux, sub3.flux) + # NOTE: uncomment sub4 test once Spectrum1D and Background subtraction works + # sub4 = image_um - bg4 + sub5 = bg4.sub_image(image_um) + sub6 = bg4.sub_image() + assert np.allclose(sub2.flux, sub5.flux) + # assert np.allclose(sub4.flux, sub5.flux) + assert np.allclose(sub5.flux, sub6.flux) + bkg_spec = bg1.bkg_spectrum() assert isinstance(bkg_spec, Spectrum1D) sub_spec = bg1.sub_spectrum() From c9e076ef47b6149fe023f50e3b4ee9aa346daac9 Mon Sep 17 00:00:00 2001 From: ojustino Date: Mon, 21 Nov 2022 16:47:08 -0500 Subject: [PATCH 09/12] Make default flux unit DN; bump minimum Python req --- setup.cfg | 2 +- specreduce/core.py | 2 +- specreduce/extract.py | 8 ++++---- specreduce/tests/test_extract.py | 4 ++-- 4 files changed, 8 insertions(+), 8 deletions(-) diff --git a/setup.cfg b/setup.cfg index b7d37a2a..673d111b 100644 --- a/setup.cfg +++ b/setup.cfg @@ -14,7 +14,7 @@ github_project = astropy/specreduce [options] zip_safe = False packages = find: -python_requires = >=3.7 +python_requires = >=3.8 setup_requires = setuptools_scm install_requires = astropy diff --git a/specreduce/core.py b/specreduce/core.py index 4dad59aa..c8fa00b4 100644 --- a/specreduce/core.py +++ b/specreduce/core.py @@ -70,7 +70,7 @@ def _parse_image(self, image, disp_axis=1): else: uncertainty = VarianceUncertainty(np.ones(img.shape)) - unit = getattr(image, 'unit', u.Unit('DN')) # or u.Unit()? + unit = getattr(image, 'unit', u.Unit('DN')) spectral_axis = getattr(image, 'spectral_axis', np.arange(img.shape[disp_axis]) * u.pix) diff --git a/specreduce/extract.py b/specreduce/extract.py index 64dfc0a7..b792be88 100644 --- a/specreduce/extract.py +++ b/specreduce/extract.py @@ -258,7 +258,7 @@ class HorneExtract(SpecreduceOperation): unit : `~astropy.units.Unit` or str, optional (Only used if ``image`` is not an NDData object.) The associated unit for the data in ``image``. If blank, - fluxes are interpreted as unitless. [default: None] + fluxes are interpreted in DN. [default: None] """ image: NDData @@ -307,7 +307,7 @@ def _parse_image(self, image, unit : `~astropy.units.Unit` or str, optional (Only used if ``image`` is not an NDData object.) The associated unit for the data in ``image``. If blank, - fluxes are interpreted as unitless. + fluxes are interpreted in DN. disp_axis : int, optional The index of the image's dispersion axis. Should not be changed until operations can handle variable image @@ -381,7 +381,7 @@ def _parse_image(self, image, variance = VarianceUncertainty(variance) unit = getattr(image, 'unit', - u.Unit(unit) if unit is not None else u.Unit()) + u.Unit(unit) if unit is not None else u.Unit('DN')) spectral_axis = getattr(image, 'spectral_axis', np.arange(img.shape[disp_axis]) * u.pix) @@ -434,7 +434,7 @@ def __call__(self, image=None, trace_object=None, unit : `~astropy.units.Unit` or str, optional (Only used if ``image`` is not an NDData object.) The associated unit for the data in ``image``. If blank, - fluxes are interpreted as unitless. + fluxes are interpreted in DN. Returns diff --git a/specreduce/tests/test_extract.py b/specreduce/tests/test_extract.py index 17552922..c711f2e3 100644 --- a/specreduce/tests/test_extract.py +++ b/specreduce/tests/test_extract.py @@ -116,10 +116,10 @@ def test_horne_image_validation(): ext = extract(image=image.data, variance=err, mask=mask) # an array-type image given without mask and unit arguments is fine - # and produces an extraction with unitless flux and spectral axis in pixels + # and produces an extraction with flux in DN and spectral axis in pixels err = np.ones_like(image) ext = extract(image=image.data, variance=err, mask=None, unit=None) - assert ext.unit == u.Unit() + assert ext.unit == u.Unit('DN') assert np.all(ext.spectral_axis == np.arange(image.shape[extract.disp_axis]) * u.pix) From 53fd1a456d33231cf101b32d57f7e47d46e3cc33 Mon Sep 17 00:00:00 2001 From: ojustino Date: Mon, 21 Nov 2022 18:21:30 -0500 Subject: [PATCH 10/12] Begin deprecation process for bkg_array --- specreduce/background.py | 16 ++++++++++------ specreduce/tests/test_background.py | 10 +++++----- 2 files changed, 15 insertions(+), 11 deletions(-) diff --git a/specreduce/background.py b/specreduce/background.py index b6ca2586..e6aac876 100644 --- a/specreduce/background.py +++ b/specreduce/background.py @@ -5,6 +5,7 @@ import numpy as np from astropy.nddata import NDData +from astropy.utils.decorators import deprecated_attribute from astropy import units as u from specutils import Spectrum1D @@ -56,6 +57,9 @@ class Background(_ImageParser): disp_axis: int = 1 crossdisp_axis: int = 0 + # TO-DO: update bkg_array with Spectrum1D alternative (is bkg_image enough?) + bkg_array = deprecated_attribute('bkg_array', '1.3') + def __post_init__(self): """ Determine the background from an image for subtraction. @@ -93,7 +97,7 @@ def _to_trace(trace): if self.width < 0: raise ValueError("width must be positive") if self.width == 0: - self.bkg_array = np.zeros(self.image.shape[self.disp_axis]) + self._bkg_array = np.zeros(self.image.shape[self.disp_axis]) return if isinstance(self.traces, Trace): @@ -130,13 +134,13 @@ def _to_trace(trace): self.bkg_wimage = bkg_wimage if self.statistic == 'average': - self.bkg_array = np.average(self.image.data, - weights=self.bkg_wimage, - axis=self.crossdisp_axis) + self._bkg_array = np.average(self.image.data, + weights=self.bkg_wimage, + axis=self.crossdisp_axis) elif self.statistic == 'median': med_image = self.image.data.copy() med_image[np.where(self.bkg_wimage) == 0] = np.nan - self.bkg_array = np.nanmedian(med_image, axis=self.crossdisp_axis) + self._bkg_array = np.nanmedian(med_image, axis=self.crossdisp_axis) else: raise ValueError("statistic must be 'average' or 'median'") @@ -231,7 +235,7 @@ def bkg_image(self, image=None): Spectrum1D object with same shape as ``image``. """ image = self._parse_image(image) - return Spectrum1D(np.tile(self.bkg_array, + return Spectrum1D(np.tile(self._bkg_array, (image.shape[0], 1)) * image.unit, spectral_axis=image.spectral_axis) diff --git a/specreduce/tests/test_background.py b/specreduce/tests/test_background.py index f36cec71..2182ae4a 100644 --- a/specreduce/tests/test_background.py +++ b/specreduce/tests/test_background.py @@ -38,15 +38,15 @@ def test_background(): bg1 = Background(image, [trace-bkg_sep, trace+bkg_sep], width=bkg_width) bg2 = Background.two_sided(image, trace, bkg_sep, width=bkg_width) bg3 = Background.two_sided(image, trace_pos, bkg_sep, width=bkg_width) - assert np.allclose(bg1.bkg_array, bg2.bkg_array) - assert np.allclose(bg1.bkg_array, bg3.bkg_array) + assert np.allclose(bg1.bkg_image().flux, bg2.bkg_image().flux) + assert np.allclose(bg1.bkg_image().flux, bg3.bkg_image().flux) bg4 = Background(image_um, [trace-bkg_sep, trace+bkg_sep], width=bkg_width) bg5 = Background.two_sided(image_um, trace, bkg_sep, width=bkg_width) bg6 = Background.two_sided(image_um, trace_pos, bkg_sep, width=bkg_width) - assert np.allclose(bg1.bkg_array, bg4.bkg_array) - assert np.allclose(bg1.bkg_array, bg5.bkg_array) - assert np.allclose(bg1.bkg_array, bg6.bkg_array) + assert np.allclose(bg1.bkg_image().flux, bg4.bkg_image().flux) + assert np.allclose(bg1.bkg_image().flux, bg5.bkg_image().flux) + assert np.allclose(bg1.bkg_image().flux, bg6.bkg_image().flux) # test that creating a one_sided background works Background.one_sided(image, trace, bkg_sep, width=bkg_width) From 14f48c65d1b1016c6d92b58fa25b13026cc55471 Mon Sep 17 00:00:00 2001 From: ojustino Date: Tue, 22 Nov 2022 11:46:10 -0500 Subject: [PATCH 11/12] Bump specutils requirement to v1.9.1 --- setup.cfg | 2 +- specreduce/background.py | 2 -- specreduce/tests/test_background.py | 13 +++++-------- 3 files changed, 6 insertions(+), 11 deletions(-) diff --git a/setup.cfg b/setup.cfg index 673d111b..2bed71fe 100644 --- a/setup.cfg +++ b/setup.cfg @@ -18,7 +18,7 @@ python_requires = >=3.8 setup_requires = setuptools_scm install_requires = astropy - specutils>=1.7 + specutils>=1.9.1 synphot matplotlib photutils diff --git a/specreduce/background.py b/specreduce/background.py index e6aac876..59395b00 100644 --- a/specreduce/background.py +++ b/specreduce/background.py @@ -324,6 +324,4 @@ def __rsub__(self, image): """ Subtract the background from an image. """ - # NOTE: will not be called until specutils PR #988 is merged, released, - # and pinned here return self.sub_image(image) diff --git a/specreduce/tests/test_background.py b/specreduce/tests/test_background.py index 2182ae4a..558969dd 100644 --- a/specreduce/tests/test_background.py +++ b/specreduce/tests/test_background.py @@ -55,20 +55,17 @@ def test_background(): bg = Background(image, trace, width=bkg_width) # test that image subtraction works - # NOTE: uncomment sub1 test once Spectrum1D and Background subtraction works - # (meaning specutils PR #988 is merged, released, and pinned here) - # sub1 = image - bg1 + sub1 = image - bg1 sub2 = bg1.sub_image(image) sub3 = bg1.sub_image() - # assert np.allclose(sub1.flux, sub2.flux) + assert np.allclose(sub1.flux, sub2.flux) assert np.allclose(sub2.flux, sub3.flux) - # NOTE: uncomment sub4 test once Spectrum1D and Background subtraction works - # sub4 = image_um - bg4 + sub4 = image_um - bg4 sub5 = bg4.sub_image(image_um) sub6 = bg4.sub_image() - assert np.allclose(sub2.flux, sub5.flux) - # assert np.allclose(sub4.flux, sub5.flux) + assert np.allclose(sub1.flux, sub4.flux) + assert np.allclose(sub4.flux, sub5.flux) assert np.allclose(sub5.flux, sub6.flux) bkg_spec = bg1.bkg_spectrum() From 90b6460f04150755469327302a8303b7de091eb6 Mon Sep 17 00:00:00 2001 From: ojustino Date: Tue, 22 Nov 2022 11:57:38 -0500 Subject: [PATCH 12/12] Update changelog --- CHANGES.rst | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/CHANGES.rst b/CHANGES.rst index 4b14c559..25cb6693 100644 --- a/CHANGES.rst +++ b/CHANGES.rst @@ -15,7 +15,8 @@ to images instead of predetermined [#128] - The default number of bins for FitTrace is now its associated image's number of dispersion pixels instead of 20. Its default peak_method is now 'max' [#128] - All operations now accept Spectrum1D and Quantity-type images. All accepted -image types are now processed internally as Spectrum1D objects [#146] +image types are now processed internally as Spectrum1D objects [#144] +- All operations' ``image`` attributes are now coerced Spectrum1D objects [#144] Bug Fixes ^^^^^^^^^