# First are the imports we will need as well as loading in the example data:
#
from astropy.io import fits
from astropy import units as u
import numpy as np
from matplotlib import pyplot as plt
from astropy.visualization import quantity_support
quantity_support()  # for getting units on the axes below  # doctest: +IGNORE_OUTPUT
#
filename = 'https://data.sdss.org/sas/dr16/sdss/spectro/redux/26/spectra/1323/spec-1323-52797-0012.fits'
with fits.open(filename) as f:  # doctest: +IGNORE_OUTPUT +REMOTE_DATA
    specdata = f[1].data[1020:1250]  # doctest: +REMOTE_DATA
#
# Then we re-format this dataset into astropy quantities, and create a
# `~specutils.Spectrum1D` object:
#
from specutils import Spectrum1D
lamb = 10**specdata['loglam'] * u.AA # doctest: +REMOTE_DATA
flux = specdata['flux'] * 10**-17 * u.Unit('erg cm-2 s-1 AA-1') # doctest: +REMOTE_DATA
input_spec = Spectrum1D(spectral_axis=lamb, flux=flux) # doctest: +REMOTE_DATA
#
f, ax = plt.subplots()  # doctest: +IGNORE_OUTPUT +REMOTE_DATA
ax.step(input_spec.spectral_axis, input_spec.flux) # doctest: +IGNORE_OUTPUT +REMOTE_DATA