Sub-pixel precision for calculate_fwhm.

docs/analysis.rst

@@ -152,7 +152,7 @@ Each of the width analysis functions are applied to this spectrum below:
    >>> gaussian_fwhm(noisy_gaussian) # doctest: +FLOAT_CMP
    <Quantity 1.81144683 GHz>
    >>> fwhm(noisy_gaussian) # doctest: +FLOAT_CMP
-   <Quantity 1.90954774 GHz>
+   <Quantity 1.91686888 GHz>
 
 
 Reference/API

docs/fitting.rst

@@ -158,7 +158,7 @@ Then estimate the line  parameters it it for a Gaussian line profile::
               amplitude            mean             stddev
                   Jy                um                um
           ------------------ ---------------- ------------------
-          1.1845669151078486 4.57517271067525 0.3015075376884422
+          1.1845669151078486 4.57517271067525 0.3199324375933905
 
 
 If an `~astropy.modeling.Model` is used that does not have the predefined

specutils/analysis/width.py

@@ -154,8 +154,8 @@ def _compute_single_fwhm(flux, spectral_axis):
 
     argmax = flux_value.argmax()
     halfval = flux_value[argmax] / 2
-    left = flux_value[:argmax] <= halfval
-    right = flux_value[argmax + 1:] <= halfval
+    left = flux_value[:argmax] < halfval
+    right = flux_value[argmax + 1:] < halfval
 
     # Highest signal at the first point
     i0 = np.nonzero(left)[0]
@@ -176,7 +176,7 @@ def _compute_single_fwhm(flux, spectral_axis):
     if i1.size == 0:
         right_value = spectral_value[-1]
     else:
-        i1 = i1[0] + argmax
+        i1 = i1[0] + argmax + 1
         i0 = i1 - 1
         left_flux = flux_value[i0]
         left_spectral = spectral_value[i0]

specutils/tests/test_analysis.py

@@ -483,7 +483,11 @@ def test_fwhm():
 
     spectrum = Spectrum1D(spectral_axis=wavelengths, flux=flux)
     result = fwhm(spectrum)
-    assert result == 0.9090909090909092*u.um
+    # Note that this makes a little more sense than the previous version;
+    # since the maximum value occurs at wavelength=1, and the half-value of
+    # flux (0.5) occurs at exactly wavelength=2, the result should be
+    # exactly 1 (2 - 1).
+    assert result == 1.0*u.um
 
     # Highest point at the last point
     wavelengths = np.linspace(1, 10, 100) * u.um

specutils/tests/test_fitting.py

@@ -106,7 +106,7 @@ def test_single_peak_estimate():
 
     assert np.isclose(g_init.amplitude.value, 3.354169257846847)
     assert np.isclose(g_init.mean.value, 6.218588636687762)
-    assert np.isclose(g_init.stddev.value, 1.608040201005025)
+    assert np.isclose(g_init.stddev.value, 1.6339001193853715)
 
     assert g_init.amplitude.unit == u.Jy
     assert g_init.mean.unit == u.um
@@ -120,7 +120,7 @@ def test_single_peak_estimate():
 
     assert np.isclose(g_init.amplitude.value, 3.354169257846847)
     assert np.isclose(g_init.x_0.value, 6.218588636687762)
-    assert np.isclose(g_init.fwhm.value, 1.608040201005025)
+    assert np.isclose(g_init.fwhm.value, 1.6339001193853715)
 
     assert g_init.amplitude.unit == u.Jy
     assert g_init.x_0.unit == u.um
@@ -134,8 +134,8 @@ def test_single_peak_estimate():
 
     assert np.isclose(g_init.amplitude_L.value, 3.354169257846847)
     assert np.isclose(g_init.x_0.value, 6.218588636687762)
-    assert np.isclose(g_init.fwhm_L.value, 1.1370561305512321)
-    assert np.isclose(g_init.fwhm_G.value, 1.1370561305512321)
+    assert np.isclose(g_init.fwhm_L.value, 1.1553418541989058)
+    assert np.isclose(g_init.fwhm_G.value, 1.1553418541989058)
 
     assert g_init.amplitude_L.unit == u.Jy
     assert g_init.x_0.unit == u.um
@@ -158,7 +158,7 @@ def test_single_peak_estimate():
 
     assert np.isclose(g_init.amplitude.value, 3.354169257846847)
     assert np.isclose(g_init.x_0.value, 6.218588636687762)
-    assert np.isclose(g_init.stddev.value, 1.608040201005025)
+    assert np.isclose(g_init.stddev.value, 1.6339001193853715)
 
     assert g_init.amplitude.unit == u.Jy
     assert g_init.x_0.unit == u.um