Skip to content

Commit

Permalink
Merge 0b6d17c into 638a2af
Browse files Browse the repository at this point in the history
  • Loading branch information
ibusko committed Jul 24, 2019
2 parents 638a2af + 0b6d17c commit 8d6326b
Show file tree
Hide file tree
Showing 3 changed files with 59 additions and 3 deletions.
9 changes: 8 additions & 1 deletion specutils/spectra/spectrum1d.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,9 @@ class Spectrum1D(OneDSpectrumMixin, NDDataRef):
Any quantity supported by the standard spectral equivalencies
(wavelength, energy, frequency, wave number). Describes the rest value
of the spectral axis for use with velocity conversions.
redshift : float
The redshift correction to be applied to the spectral axis for use
with velocity conversions.
uncertainty : `~astropy.nddata.NDUncertainty`
Contains uncertainty information along with propagation rules for
spectrum arithmetic. Can take a unit, but if none is given, will use
Expand All @@ -41,7 +44,7 @@ class Spectrum1D(OneDSpectrumMixin, NDDataRef):
around with the spectrum container object.
"""
def __init__(self, flux=None, spectral_axis=None, wcs=None,
velocity_convention=None, rest_value=None, *args, **kwargs):
velocity_convention=None, rest_value=None, redshift=None, *args, **kwargs):
# If the flux (data) argument is a subclass of nddataref (as it would
# be for internal arithmetic operations), avoid setup entirely.
if isinstance(flux, NDDataRef):
Expand Down Expand Up @@ -108,6 +111,10 @@ def __init__(self, flux=None, spectral_axis=None, wcs=None,
elif not self._rest_value.unit.is_equivalent(u.AA) and not self._rest_value.unit.is_equivalent(u.Hz):
raise u.UnitsError("Rest value must be energy/wavelength/frequency equivalent.")

self._redshift = redshift
if redshift is None:
self._redshift = 0.

super(Spectrum1D, self).__init__(
data=flux.value if isinstance(flux, u.Quantity) else flux,
wcs=wcs, **kwargs)
Expand Down
7 changes: 5 additions & 2 deletions specutils/spectra/spectrum_mixin.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,8 +2,8 @@
from copy import deepcopy

import numpy as np
from astropy.wcs import WCSSUB_SPECTRAL
from astropy.units import Unit
from astropy.units import Quantity
from astropy.constants import c as C
from astropy import units as u
import astropy.units.equivalencies as eq
from astropy.utils.decorators import lazyproperty
Expand Down Expand Up @@ -185,6 +185,9 @@ def velocity(self):

new_data = self.spectral_axis.to(u.km/u.s, equivalencies=equiv)

# Redshift-correct velocities using approximate low-z correction.
new_data += Quantity((self._redshift * C), new_data.unit)

return new_data

def with_spectral_unit(self, unit, velocity_convention=None,
Expand Down
46 changes: 46 additions & 0 deletions specutils/tests/test_spectrum1d.py
Original file line number Diff line number Diff line change
Expand Up @@ -106,6 +106,52 @@ def test_spectral_axis_conversions():
new_spec = spec.with_spectral_unit(u.GHz)


def test_redshift():
spec = Spectrum1D(flux=np.array([26.0, 30., 44.5]) * u.Jy,
spectral_axis=np.array([4000, 6000, 8000]) * u.AA,
velocity_convention='optical',
rest_value=6000 * u.AA)

assert spec.velocity.unit == u.Unit('km/s')
assert spec.velocity[0].value == pytest.approx(-99930.8, abs=0.5)
assert spec.velocity[1].value == pytest.approx(0.0, abs=0.5)
assert spec.velocity[2].value == pytest.approx(99930.8, abs=0.5)

spec = Spectrum1D(flux=np.array([26.0, 30., 44.5]) * u.Jy,
spectral_axis=np.array([4000, 6000, 8000]) * u.AA,
velocity_convention='optical',
rest_value=6000 * u.AA,
redshift = 0.1)

assert spec.velocity.unit == u.Unit('km/s')
assert spec.velocity[0].value == pytest.approx(-69951.3, abs=0.5)
assert spec.velocity[1].value == pytest.approx(29979.2, abs=0.5)
assert spec.velocity[2].value == pytest.approx(129910.1, abs=0.5)

#-------------------------

spec = Spectrum1D(flux=np.array([26.0, 30.0, 44.5]) * u.Jy,
spectral_axis=np.array([10.5, 11.0, 11.5]) * u.GHz,
velocity_convention='radio',
rest_value=11.0 * u.GHz)

assert spec.velocity.unit == u.Unit('km/s')
assert spec.velocity[0].value == pytest.approx(13626., abs=1.)
assert spec.velocity[1].value == pytest.approx(0., abs=1.)
assert spec.velocity[2].value == pytest.approx(-13626., abs=1.)

spec = Spectrum1D(flux=np.array([26.0, 30.0, 44.5]) * u.Jy,
spectral_axis=np.array([10.5, 11.0, 11.5]) * u.GHz,
velocity_convention='radio',
rest_value=11.0 * u.GHz,
redshift = 0.1)

assert spec.velocity.unit == u.Unit('km/s')
assert spec.velocity[0].value == pytest.approx(43606., abs=1.)
assert spec.velocity[1].value == pytest.approx(29979., abs=1.)
assert spec.velocity[2].value == pytest.approx(16352., abs=1.)


def test_flux_unit_conversion():
# By default the flux units should be set to Jy
s = Spectrum1D(flux=np.array([26.0, 44.5]) * u.Jy,
Expand Down

0 comments on commit 8d6326b

Please sign in to comment.