From ee87cdbeb1c19f66b1a7260010ba513c14ae0686 Mon Sep 17 00:00:00 2001 From: brechmos Date: Thu, 9 Aug 2018 14:32:16 -0400 Subject: [PATCH 01/18] intiial pass of centroid --- specutils/analysis/__init__.py | 1 + specutils/analysis/centroid.py | 84 ++++++++++++++++++++ specutils/tests/test_analysis.py | 127 ++++++++++++++++++++++++++++++- 3 files changed, 210 insertions(+), 2 deletions(-) create mode 100644 specutils/analysis/centroid.py diff --git a/specutils/analysis/__init__.py b/specutils/analysis/__init__.py index 47dd37e9d..2514ef9df 100644 --- a/specutils/analysis/__init__.py +++ b/specutils/analysis/__init__.py @@ -2,3 +2,4 @@ from .equivalent_width import * # noqa from .snr import snr # noqa +from .centroid import centroid # noqa diff --git a/specutils/analysis/centroid.py b/specutils/analysis/centroid.py new file mode 100644 index 000000000..edcc3a11a --- /dev/null +++ b/specutils/analysis/centroid.py @@ -0,0 +1,84 @@ +import numpy as np +from ..spectra import SpectralRegion + +__all__ = ['centroid'] + + +def centroid(spectrum, region=None): + """ + Calculate the centroid of the spectrum based on the flux and uncertainty + in the spectrum. This will be calculated over the regions, if they + are specified. + + Parameters + ---------- + spectrum : `~specutils.spectra.spectrum1d.Spectrum1D` + The spectrum object overwhich the equivalent width will be calculated. + + region: `~specutils.utils.SpectralRegion` or list of `~specutils.utils.SpectralRegion` + Region within the spectrum to calculate the centroid. + + Returns + ------- + centroid : float or list (based on region input) + Signal to noise ratio of the spectrum or within the regions + + Notes + ----- + The spectrum will need to be continuum subtracted before calling + this method. + + """ + + # No region, therefore whole spectrum. + if region is None: + return _centroid_single_region(spectrum) + + # Single region + elif isinstance(region, SpectralRegion): + return _ntroid(spectrum, region=region) + + # List of regions + elif isinstance(region, list): + return [_centroid_single_region(spectrum, region=reg) + for reg in region] + + +def _centroid_single_region(spectrum, region=None): + """ + Calculate the centroid of the spectrum based on the flux and uncertainty + in the spectrum. + + Parameters + ---------- + spectrum : `~specutils.spectra.spectrum1d.Spectrum1D` + The spectrum object overwhich the equivalent width will be calculated. + + region: `~specutils.utils.SpectralRegion` + Region within the spectrum to calculate the centroid. + + Returns + ------- + centroid : float or list (based on region input) + Centroid of the spectrum or within the regions + + Notes + ----- + This is a helper function for the above `centroid()` method. + + """ + + if region is not None: + calc_spectrum = region.extract(spectrum) + else: + calc_spectrum = spectrum + + flux = calc_spectrum.flux + dispersion = calc_spectrum.spectral_axis + + if len(flux.shape) > 1: + dispersion = np.tile(dispersion[np.newaxis].T, [1, flux.shape[1]]) + + # the axis=0 will enable this to run on single-dispersion, single-flux + # and single-dispersion, multiple-flux + return np.sum(flux * dispersion, axis=0) / np.sum(calc_spectrum.spectral_axis, axis=0) diff --git a/specutils/tests/test_analysis.py b/specutils/tests/test_analysis.py index 9aaa771ce..27be99b29 100644 --- a/specutils/tests/test_analysis.py +++ b/specutils/tests/test_analysis.py @@ -7,7 +7,7 @@ from ..tests.spectral_examples import simulated_spectra from ..spectra.spectrum1d import Spectrum1D import astropy.units as u -from ..analysis import equivalent_width, snr +from ..analysis import equivalent_width, snr, centroid from ..manipulation import noise_region_uncertainty from ..spectra import SpectralRegion from astropy.nddata import StdDevUncertainty @@ -178,7 +178,7 @@ def test_snr_single_region_with_noise_region(simulated_spectra): region = SpectralRegion(0.52*u.um, 0.59*u.um) noise_region = SpectralRegion(0.40*u.um, 0.45*u.um) - + # # Set up the data # @@ -205,4 +205,127 @@ def test_snr_single_region_with_noise_region(simulated_spectra): assert np.allclose(spec_snr.value, spec_snr_expected.value) +# -- + +def test_centroid(simulated_spectra): + """ + Test the simple version of the spectral centroid. + """ + + np.random.seed(42) + + # + # Set up the data and add the uncertainty and calculate the expected SNR + # + + spectrum = simulated_spectra.s1_um_mJy_e1 + uncertainty = StdDevUncertainty(0.1*np.random.random(len(spectrum.flux))*u.mJy) + spectrum.uncertainty = uncertainty + + wavelengths = spectrum.spectral_axis + flux = spectrum.flux + + spec_centroid_expected = np.sum(flux * wavelengths) / np.sum(wavelengths) + + # + # SNR of the whole spectrum + # + + spec_centroid = centroid(spectrum) + + assert isinstance(spec_centroid, u.Quantity) + assert np.allclose(spec_centroid.value, spec_centroid_expected.value) + + +def test_snr_multiple_flux(simulated_spectra): + """ + Test the simple version of the spectral SNR, with multiple flux per single dispersion. + """ + + # + # Set up the data and add the uncertainty and calculate the expected SNR + # + + np.random.seed(42) + spec = Spectrum1D(spectral_axis=np.arange(10) * u.um, + flux=np.random.sample((10, 5)) * u.Jy) + + centroid_spec = centroid(spec) + + assert np.allclose(centroid_spec.value, np.array([0.5191939 , 0.31976068, 0.30832925, 0.65393441, 0.36912737])) + + +# def test_snr_single_region(simulated_spectra): +# """ +# Test the simple version of the spectral SNR over a region of the spectrum. +# """ +# +# np.random.seed(42) +# +# region = SpectralRegion(0.52*u.um, 0.59*u.um) +# +# # +# # Set up the data +# # +# +# spectrum = simulated_spectra.s1_um_mJy_e1 +# uncertainty = StdDevUncertainty(0.1*np.random.random(len(spectrum.flux))*u.mJy) +# spectrum.uncertainty = uncertainty +# +# wavelengths = spectrum.spectral_axis +# flux = spectrum.flux +# +# l = np.nonzero(wavelengths>region.lower)[0][0] +# r = np.nonzero(wavelengthsregion.lower)[0][0] +# r = np.nonzero(wavelengths Date: Fri, 10 Aug 2018 08:49:09 -0400 Subject: [PATCH 02/18] fixed incorrect definition of centroid --- specutils/analysis/centroid.py | 2 +- specutils/tests/test_analysis.py | 5 +++-- 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/specutils/analysis/centroid.py b/specutils/analysis/centroid.py index edcc3a11a..8756630a5 100644 --- a/specutils/analysis/centroid.py +++ b/specutils/analysis/centroid.py @@ -81,4 +81,4 @@ def _centroid_single_region(spectrum, region=None): # the axis=0 will enable this to run on single-dispersion, single-flux # and single-dispersion, multiple-flux - return np.sum(flux * dispersion, axis=0) / np.sum(calc_spectrum.spectral_axis, axis=0) + return np.sum(flux * dispersion, axis=0) / np.sum(flux, axis=0) diff --git a/specutils/tests/test_analysis.py b/specutils/tests/test_analysis.py index 27be99b29..9be47b6ad 100644 --- a/specutils/tests/test_analysis.py +++ b/specutils/tests/test_analysis.py @@ -225,7 +225,7 @@ def test_centroid(simulated_spectra): wavelengths = spectrum.spectral_axis flux = spectrum.flux - spec_centroid_expected = np.sum(flux * wavelengths) / np.sum(wavelengths) + spec_centroid_expected = np.sum(flux * wavelengths) / np.sum(flux) # # SNR of the whole spectrum @@ -253,7 +253,8 @@ def test_snr_multiple_flux(simulated_spectra): centroid_spec = centroid(spec) - assert np.allclose(centroid_spec.value, np.array([0.5191939 , 0.31976068, 0.30832925, 0.65393441, 0.36912737])) + assert np.allclose(centroid_spec.value, np.array([5.39321967, 3.6856305 , 3.09779811, 4.99442161, 4.50267016])) + assert centroid_spec.unit == u.um # def test_snr_single_region(simulated_spectra): From 3648a007a6d4496e1ca0209a832bbdfc9ae8ae4f Mon Sep 17 00:00:00 2001 From: brechmos Date: Fri, 10 Aug 2018 08:59:44 -0400 Subject: [PATCH 03/18] added documentation for centroid --- docs/basic_analysis.rst | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) diff --git a/docs/basic_analysis.rst b/docs/basic_analysis.rst index 3d0a30a12..524d4cc7e 100644 --- a/docs/basic_analysis.rst +++ b/docs/basic_analysis.rst @@ -42,6 +42,23 @@ Currently, specutils supports basic signal-to-noise ratio calculations. >>> snr(spec) #doctest:+SKIP +Centroid +-------- + +Currently, specutils supports basic centroid calculations. + +.. code-block:: python + + >>> import numpy as np + >>> import astropy.units as u + >>> from specutils.spectra import Spectrum1D + >>> from astropy.nddata import StdDevUncertainty + >>> from specutils.analysis import centroid + + >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(3+np.random.randn(50))*u.Jy) + >>> centroid(spec) #doctest:+SKIP + + Reference/API ------------- .. automodapi:: specutils.analysis From 255dd5103debeb8c4d53601ca29685cdd54e8124 Mon Sep 17 00:00:00 2001 From: brechmos Date: Fri, 10 Aug 2018 12:38:51 -0400 Subject: [PATCH 04/18] intermediate --- specutils/analysis/centroid.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/specutils/analysis/centroid.py b/specutils/analysis/centroid.py index 8756630a5..282b783d5 100644 --- a/specutils/analysis/centroid.py +++ b/specutils/analysis/centroid.py @@ -36,7 +36,7 @@ def centroid(spectrum, region=None): # Single region elif isinstance(region, SpectralRegion): - return _ntroid(spectrum, region=region) + return _centroid_single_region(spectrum, region=region) # List of regions elif isinstance(region, list): @@ -69,7 +69,10 @@ def _centroid_single_region(spectrum, region=None): """ if region is not None: + print('original spectrum is {}'.format(spectrum.spectral_axis)) calc_spectrum = region.extract(spectrum) + print('extracted region is {}'.format(calc_spectrum)) + print('extracted region is {}'.format(calc_spectrum.spectral_axis)) else: calc_spectrum = spectrum From 234e555208973afed4138d8436517cd6702b5cef Mon Sep 17 00:00:00 2001 From: brechmos Date: Mon, 27 Aug 2018 07:30:39 -0400 Subject: [PATCH 05/18] small documentation fixes --- specutils/analysis/centroid.py | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/specutils/analysis/centroid.py b/specutils/analysis/centroid.py index 282b783d5..7228eb4dd 100644 --- a/specutils/analysis/centroid.py +++ b/specutils/analysis/centroid.py @@ -13,7 +13,7 @@ def centroid(spectrum, region=None): Parameters ---------- spectrum : `~specutils.spectra.spectrum1d.Spectrum1D` - The spectrum object overwhich the equivalent width will be calculated. + The spectrum object overwhich the centroid will be calculated. region: `~specutils.utils.SpectralRegion` or list of `~specutils.utils.SpectralRegion` Region within the spectrum to calculate the centroid. @@ -21,12 +21,13 @@ def centroid(spectrum, region=None): Returns ------- centroid : float or list (based on region input) - Signal to noise ratio of the spectrum or within the regions + Centroid of the spectrum or within the regions Notes ----- The spectrum will need to be continuum subtracted before calling - this method. + this method. See the + `analysis documentation `_ for more information. """ @@ -52,7 +53,7 @@ def _centroid_single_region(spectrum, region=None): Parameters ---------- spectrum : `~specutils.spectra.spectrum1d.Spectrum1D` - The spectrum object overwhich the equivalent width will be calculated. + The spectrum object overwhich the centroid will be calculated. region: `~specutils.utils.SpectralRegion` Region within the spectrum to calculate the centroid. From 569cf4ae2192e3535f136b55e9af7a35a8715eb2 Mon Sep 17 00:00:00 2001 From: brechmos Date: Mon, 27 Aug 2018 07:32:48 -0400 Subject: [PATCH 06/18] uncommented and added pytest xfail --- specutils/tests/test_analysis.py | 148 ++++++++++++++++--------------- 1 file changed, 75 insertions(+), 73 deletions(-) diff --git a/specutils/tests/test_analysis.py b/specutils/tests/test_analysis.py index 9be47b6ad..2a01cf1ee 100644 --- a/specutils/tests/test_analysis.py +++ b/specutils/tests/test_analysis.py @@ -257,76 +257,78 @@ def test_snr_multiple_flux(simulated_spectra): assert centroid_spec.unit == u.um -# def test_snr_single_region(simulated_spectra): -# """ -# Test the simple version of the spectral SNR over a region of the spectrum. -# """ -# -# np.random.seed(42) -# -# region = SpectralRegion(0.52*u.um, 0.59*u.um) -# -# # -# # Set up the data -# # -# -# spectrum = simulated_spectra.s1_um_mJy_e1 -# uncertainty = StdDevUncertainty(0.1*np.random.random(len(spectrum.flux))*u.mJy) -# spectrum.uncertainty = uncertainty -# -# wavelengths = spectrum.spectral_axis -# flux = spectrum.flux -# -# l = np.nonzero(wavelengths>region.lower)[0][0] -# r = np.nonzero(wavelengthsregion.lower)[0][0] -# r = np.nonzero(wavelengthsregion.lower)[0][0] + r = np.nonzero(wavelengthsregion.lower)[0][0] + r = np.nonzero(wavelengths Date: Mon, 27 Aug 2018 07:39:50 -0400 Subject: [PATCH 07/18] more documentation clarification --- specutils/analysis/centroid.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/specutils/analysis/centroid.py b/specutils/analysis/centroid.py index 7228eb4dd..b7fbf533e 100644 --- a/specutils/analysis/centroid.py +++ b/specutils/analysis/centroid.py @@ -4,11 +4,9 @@ __all__ = ['centroid'] -def centroid(spectrum, region=None): +def centroid(spectrum, region): """ - Calculate the centroid of the spectrum based on the flux and uncertainty - in the spectrum. This will be calculated over the regions, if they - are specified. + Calculate the centroid of a region, or regions, of the spectrum. Parameters ---------- From b06ce1ef2fff335180636503418546de686ce86b Mon Sep 17 00:00:00 2001 From: brechmos Date: Mon, 27 Aug 2018 07:45:36 -0400 Subject: [PATCH 08/18] added in continuum documentation --- docs/basic_analysis.rst | 15 +++++++++++++++ 1 file changed, 15 insertions(+) diff --git a/docs/basic_analysis.rst b/docs/basic_analysis.rst index 524d4cc7e..0672ebded 100644 --- a/docs/basic_analysis.rst +++ b/docs/basic_analysis.rst @@ -59,6 +59,21 @@ Currently, specutils supports basic centroid calculations. >>> centroid(spec) #doctest:+SKIP +And if the spectrum contains a continuum, then it should be subtracted first: +.. code-block:: python + + >>> import numpy as np + >>> import astropy.units as u + >>> from specutils.spectra import Spectrum1D + >>> from astropy.nddata import StdDevUncertainty + >>> from specutils.fitting import continuum + >>> from specutils.analysis import centroid + + >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(3+np.random.randn(50))*u.Jy) + >>> continuum_baseline = continuum(spec) + >>> centroid(spec-continuum_baseline) #doctest:+SKIP + + Reference/API ------------- .. automodapi:: specutils.analysis From 2cb7ce4100eb0f7b29d7c6713828aaddcc6ed5f5 Mon Sep 17 00:00:00 2001 From: brechmos Date: Mon, 27 Aug 2018 07:46:09 -0400 Subject: [PATCH 09/18] added in continuum subtraction example --- docs/basic_analysis.rst | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/docs/basic_analysis.rst b/docs/basic_analysis.rst index 0672ebded..f28c14395 100644 --- a/docs/basic_analysis.rst +++ b/docs/basic_analysis.rst @@ -71,8 +71,7 @@ And if the spectrum contains a continuum, then it should be subtracted first: >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(3+np.random.randn(50))*u.Jy) >>> continuum_baseline = continuum(spec) - >>> centroid(spec-continuum_baseline) #doctest:+SKIP - + >>> c = centroid(spec-continuum_baseline) #doctest:+SKIP Reference/API ------------- From 0d0d8b36a3d6363470f795fc4d2481adc642cb15 Mon Sep 17 00:00:00 2001 From: Erik Tollerud Date: Wed, 5 Sep 2018 17:05:14 -0400 Subject: [PATCH 10/18] fix failing test --- specutils/tests/test_analysis.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/specutils/tests/test_analysis.py b/specutils/tests/test_analysis.py index 2a01cf1ee..fec19f44c 100644 --- a/specutils/tests/test_analysis.py +++ b/specutils/tests/test_analysis.py @@ -43,7 +43,7 @@ def test_snr(simulated_spectra): """ np.random.seed(42) - + # # Set up the data and add the uncertainty and calculate the expected SNR # @@ -231,7 +231,7 @@ def test_centroid(simulated_spectra): # SNR of the whole spectrum # - spec_centroid = centroid(spectrum) + spec_centroid = centroid(spectrum, None) assert isinstance(spec_centroid, u.Quantity) assert np.allclose(spec_centroid.value, spec_centroid_expected.value) @@ -251,7 +251,7 @@ def test_snr_multiple_flux(simulated_spectra): spec = Spectrum1D(spectral_axis=np.arange(10) * u.um, flux=np.random.sample((10, 5)) * u.Jy) - centroid_spec = centroid(spec) + centroid_spec = centroid(spec, None) assert np.allclose(centroid_spec.value, np.array([5.39321967, 3.6856305 , 3.09779811, 4.99442161, 4.50267016])) assert centroid_spec.unit == u.um @@ -305,7 +305,7 @@ def test_snr_two_regions(simulated_spectra): # regions = [SpectralRegion(0.52*u.um, 0.59*u.um), SpectralRegion(0.8*u.um, 0.9*u.um)] - + # # Set up the data # From a61f38c8403ccd39bab14f6b78c5b45a7c6c5ed2 Mon Sep 17 00:00:00 2001 From: Erik Tollerud Date: Wed, 5 Sep 2018 17:20:59 -0400 Subject: [PATCH 11/18] add invert test --- specutils/tests/test_analysis.py | 21 +++++++++++++++++++++ 1 file changed, 21 insertions(+) diff --git a/specutils/tests/test_analysis.py b/specutils/tests/test_analysis.py index fec19f44c..a27808a2c 100644 --- a/specutils/tests/test_analysis.py +++ b/specutils/tests/test_analysis.py @@ -237,6 +237,21 @@ def test_centroid(simulated_spectra): assert np.allclose(spec_centroid.value, spec_centroid_expected.value) +def test_inverted_centroid(simulated_spectra): + """ + Ensures the centroid calculation also works for *inverted* spectra - i.e. + continuum-subtracted absorption lines. + """ + spectrum = simulated_spectra.s1_um_mJy_e1 + spec_centroid_expected = (np.sum(spectrum.flux * spectrum.spectral_axis) / + np.sum(spectrum.flux)) + + spectrum_inverted = Spectrum1D(spectral_axis=spectrum.spectral_axis, + flux=-spectrum.flux) + spec_centroid_inverted = centroid(spectrum_inverted, None) + assert np.allclose(spec_centroid_inverted.value, spec_centroid_expected.value) + + def test_snr_multiple_flux(simulated_spectra): """ Test the simple version of the spectral SNR, with multiple flux per single dispersion. @@ -256,6 +271,12 @@ def test_snr_multiple_flux(simulated_spectra): assert np.allclose(centroid_spec.value, np.array([5.39321967, 3.6856305 , 3.09779811, 4.99442161, 4.50267016])) assert centroid_spec.unit == u.um + # spec_inverted = Spectrum1D(spectral_axis=spec.spectral_axis, + # flux=-spec.flux) + # spec_centroid_inverted = centroid(spec_inverted, None, invert=True) + # assert np.allclose(spec_centroid_inverted.value, spec_centroid_expected_inverted.value) + + @pytest.mark.xfail def test_snr_single_region(simulated_spectra): From abc67441a09991bb81c7ba214f3f2699822ed63d Mon Sep 17 00:00:00 2001 From: brechmos Date: Thu, 6 Sep 2018 13:44:14 -0400 Subject: [PATCH 12/18] noop for tests --- specutils/tests/test_analysis.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/specutils/tests/test_analysis.py b/specutils/tests/test_analysis.py index a27808a2c..477fc8f76 100644 --- a/specutils/tests/test_analysis.py +++ b/specutils/tests/test_analysis.py @@ -299,8 +299,8 @@ def test_snr_single_region(simulated_spectra): wavelengths = spectrum.spectral_axis flux = spectrum.flux - l = np.nonzero(wavelengths>region.lower)[0][0] - r = np.nonzero(wavelengths region.lower)[0][0] + r = np.nonzero(wavelengths < region.upper)[0][-1] spec_snr_expected = np.mean(flux[l:r] / (uncertainty.array[l:r]*uncertainty.unit)) From fd7673db0b23775fc308e7b3f229e750ccd4fdf2 Mon Sep 17 00:00:00 2001 From: brechmos Date: Thu, 6 Sep 2018 14:21:34 -0400 Subject: [PATCH 13/18] fixed for incorrect transpose on multi-flux, single dispersion --- specutils/analysis/centroid.py | 9 +++------ specutils/tests/test_analysis.py | 4 ++-- 2 files changed, 5 insertions(+), 8 deletions(-) diff --git a/specutils/analysis/centroid.py b/specutils/analysis/centroid.py index b7fbf533e..736269394 100644 --- a/specutils/analysis/centroid.py +++ b/specutils/analysis/centroid.py @@ -68,10 +68,7 @@ def _centroid_single_region(spectrum, region=None): """ if region is not None: - print('original spectrum is {}'.format(spectrum.spectral_axis)) calc_spectrum = region.extract(spectrum) - print('extracted region is {}'.format(calc_spectrum)) - print('extracted region is {}'.format(calc_spectrum.spectral_axis)) else: calc_spectrum = spectrum @@ -79,8 +76,8 @@ def _centroid_single_region(spectrum, region=None): dispersion = calc_spectrum.spectral_axis if len(flux.shape) > 1: - dispersion = np.tile(dispersion[np.newaxis].T, [1, flux.shape[1]]) + dispersion = np.tile(dispersion, [flux.shape[0], 1]) - # the axis=0 will enable this to run on single-dispersion, single-flux + # the axis=-1 will enable this to run on single-dispersion, single-flux # and single-dispersion, multiple-flux - return np.sum(flux * dispersion, axis=0) / np.sum(flux, axis=0) + return np.sum(flux * dispersion, axis=-1) / np.sum(flux, axis=-1) diff --git a/specutils/tests/test_analysis.py b/specutils/tests/test_analysis.py index 477fc8f76..0840e1598 100644 --- a/specutils/tests/test_analysis.py +++ b/specutils/tests/test_analysis.py @@ -264,11 +264,11 @@ def test_snr_multiple_flux(simulated_spectra): np.random.seed(42) spec = Spectrum1D(spectral_axis=np.arange(10) * u.um, - flux=np.random.sample((10, 5)) * u.Jy) + flux=np.random.sample((5, 10)) * u.Jy) centroid_spec = centroid(spec, None) - assert np.allclose(centroid_spec.value, np.array([5.39321967, 3.6856305 , 3.09779811, 4.99442161, 4.50267016])) + assert np.allclose(centroid_spec.value, np.array([4.46190995, 4.17223565, 4.37778249, 4.51595259, 4.7429066])) assert centroid_spec.unit == u.um # spec_inverted = Spectrum1D(spectral_axis=spec.spectral_axis, From a8bbf128f800644c5745d62bdd70370c5de315a5 Mon Sep 17 00:00:00 2001 From: brechmos Date: Thu, 6 Sep 2018 14:50:51 -0400 Subject: [PATCH 14/18] fixed copy-paste error (sigh) --- docs/basic_analysis.rst | 1 - 1 file changed, 1 deletion(-) diff --git a/docs/basic_analysis.rst b/docs/basic_analysis.rst index f28c14395..24b9599dd 100644 --- a/docs/basic_analysis.rst +++ b/docs/basic_analysis.rst @@ -66,7 +66,6 @@ And if the spectrum contains a continuum, then it should be subtracted first: >>> import astropy.units as u >>> from specutils.spectra import Spectrum1D >>> from astropy.nddata import StdDevUncertainty - >>> from specutils.fitting import continuum >>> from specutils.analysis import centroid >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(3+np.random.randn(50))*u.Jy) From 97472d5b146f27a51533a9fc23d708537829c3f0 Mon Sep 17 00:00:00 2001 From: brechmos Date: Fri, 7 Sep 2018 09:21:45 -0400 Subject: [PATCH 15/18] fixed up tests --- specutils/tests/test_analysis.py | 86 +------------------------------- 1 file changed, 1 insertion(+), 85 deletions(-) diff --git a/specutils/tests/test_analysis.py b/specutils/tests/test_analysis.py index 0840e1598..4b84d70db 100644 --- a/specutils/tests/test_analysis.py +++ b/specutils/tests/test_analysis.py @@ -205,7 +205,6 @@ def test_snr_single_region_with_noise_region(simulated_spectra): assert np.allclose(spec_snr.value, spec_snr_expected.value) -# -- def test_centroid(simulated_spectra): """ @@ -252,7 +251,7 @@ def test_inverted_centroid(simulated_spectra): assert np.allclose(spec_centroid_inverted.value, spec_centroid_expected.value) -def test_snr_multiple_flux(simulated_spectra): +def test_centroid_multiple_flux(simulated_spectra): """ Test the simple version of the spectral SNR, with multiple flux per single dispersion. """ @@ -270,86 +269,3 @@ def test_snr_multiple_flux(simulated_spectra): assert np.allclose(centroid_spec.value, np.array([4.46190995, 4.17223565, 4.37778249, 4.51595259, 4.7429066])) assert centroid_spec.unit == u.um - - # spec_inverted = Spectrum1D(spectral_axis=spec.spectral_axis, - # flux=-spec.flux) - # spec_centroid_inverted = centroid(spec_inverted, None, invert=True) - # assert np.allclose(spec_centroid_inverted.value, spec_centroid_expected_inverted.value) - - - -@pytest.mark.xfail -def test_snr_single_region(simulated_spectra): - """ - Test the simple version of the spectral SNR over a region of the spectrum. - """ - - np.random.seed(42) - - region = SpectralRegion(0.52*u.um, 0.59*u.um) - - # - # Set up the data - # - - spectrum = simulated_spectra.s1_um_mJy_e1 - uncertainty = StdDevUncertainty(0.1*np.random.random(len(spectrum.flux))*u.mJy) - spectrum.uncertainty = uncertainty - - wavelengths = spectrum.spectral_axis - flux = spectrum.flux - - l = np.nonzero(wavelengths > region.lower)[0][0] - r = np.nonzero(wavelengths < region.upper)[0][-1] - - spec_snr_expected = np.mean(flux[l:r] / (uncertainty.array[l:r]*uncertainty.unit)) - - # - # SNR of the whole spectrum - # - - spec_snr = snr(spectrum, region) - - assert np.allclose(spec_snr.value, spec_snr_expected.value) - - -@pytest.mark.xfail -def test_snr_two_regions(simulated_spectra): - """ - Test the simple version of the spectral SNR within two regions. - """ - - np.random.seed(42) - - # - # Set the regions over which the SNR is calculated - # - - regions = [SpectralRegion(0.52*u.um, 0.59*u.um), SpectralRegion(0.8*u.um, 0.9*u.um)] - - # - # Set up the data - # - - spectrum = simulated_spectra.s1_um_mJy_e1 - uncertainty = StdDevUncertainty(0.1*np.random.random(len(spectrum.flux))*u.Jy) - spectrum.uncertainty = uncertainty - - wavelengths = spectrum.spectral_axis - flux = spectrum.flux - - spec_snr_expected = [] - for region in regions: - - l = np.nonzero(wavelengths>region.lower)[0][0] - r = np.nonzero(wavelengths Date: Fri, 7 Sep 2018 14:06:50 -0400 Subject: [PATCH 16/18] added in region for centroid --- docs/basic_analysis.rst | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/docs/basic_analysis.rst b/docs/basic_analysis.rst index 24b9599dd..632e217ea 100644 --- a/docs/basic_analysis.rst +++ b/docs/basic_analysis.rst @@ -56,7 +56,7 @@ Currently, specutils supports basic centroid calculations. >>> from specutils.analysis import centroid >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(3+np.random.randn(50))*u.Jy) - >>> centroid(spec) #doctest:+SKIP + >>> centroid(spec, region=None) #doctest:+SKIP And if the spectrum contains a continuum, then it should be subtracted first: @@ -66,11 +66,14 @@ And if the spectrum contains a continuum, then it should be subtracted first: >>> import astropy.units as u >>> from specutils.spectra import Spectrum1D >>> from astropy.nddata import StdDevUncertainty + >>> from specutils.fitting import fit_generic_continuum >>> from specutils.analysis import centroid - >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(3+np.random.randn(50))*u.Jy) - >>> continuum_baseline = continuum(spec) - >>> c = centroid(spec-continuum_baseline) #doctest:+SKIP + >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(10+np.random.randn(50))*u.Jy) + >>> continuum_baseline = fit_generic_continuum(spec) + >>> continuum_flux = continuum_baseline(spec.spectral_axis.value) + >>> continuum = Spectrum1D(spectral_axis=spec.spectral_axis, flux=continuum_flux) + >>> c = centroid(spec-continuum, region=None) #doctest:+SKIP Reference/API ------------- From 4c701fa011f42e5833d4a277b218043450f09073 Mon Sep 17 00:00:00 2001 From: brechmos Date: Mon, 10 Sep 2018 07:54:36 -0400 Subject: [PATCH 17/18] doc commented the results for the continuum fit --- docs/basic_analysis.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/basic_analysis.rst b/docs/basic_analysis.rst index 632e217ea..c21bd3cf9 100644 --- a/docs/basic_analysis.rst +++ b/docs/basic_analysis.rst @@ -70,7 +70,7 @@ And if the spectrum contains a continuum, then it should be subtracted first: >>> from specutils.analysis import centroid >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(10+np.random.randn(50))*u.Jy) - >>> continuum_baseline = fit_generic_continuum(spec) + >>> continuum_baseline = fit_generic_continuum(spec) #doctest:+SKIP >>> continuum_flux = continuum_baseline(spec.spectral_axis.value) >>> continuum = Spectrum1D(spectral_axis=spec.spectral_axis, flux=continuum_flux) >>> c = centroid(spec-continuum, region=None) #doctest:+SKIP From 5e455d2535b309c81df4ebf243cf9bf17fcf04dd Mon Sep 17 00:00:00 2001 From: brechmos Date: Mon, 10 Sep 2018 08:19:26 -0400 Subject: [PATCH 18/18] doctest ignores --- docs/basic_analysis.rst | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/docs/basic_analysis.rst b/docs/basic_analysis.rst index c21bd3cf9..e88d7a47b 100644 --- a/docs/basic_analysis.rst +++ b/docs/basic_analysis.rst @@ -71,8 +71,8 @@ And if the spectrum contains a continuum, then it should be subtracted first: >>> spec = Spectrum1D(spectral_axis=np.arange(50), flux=(10+np.random.randn(50))*u.Jy) >>> continuum_baseline = fit_generic_continuum(spec) #doctest:+SKIP - >>> continuum_flux = continuum_baseline(spec.spectral_axis.value) - >>> continuum = Spectrum1D(spectral_axis=spec.spectral_axis, flux=continuum_flux) + >>> continuum_flux = continuum_baseline(spec.spectral_axis.value) #doctest:+SKIP + >>> continuum = Spectrum1D(spectral_axis=spec.spectral_axis, flux=continuum_flux) #doctest:+SKIP >>> c = centroid(spec-continuum, region=None) #doctest:+SKIP Reference/API