Have box_smooth accept spectral axis units instead of just pixels

Right now, `specutils.manipulation.box_smooth` only accepts pixels as the width. It would be more convenient for some cases if it instead accepted spectral axis units, because that's a lot easier to compare across spectra - e.g. if I have data and a template and want them smoothed to the same resolution, I would rather give it in angstroms than have to figure out the number of pixels by hand.

While the API per se is simple, there's a complexity here in that most spectra do *not* have a fixed pixel-to-spectral conversions - most are relatively mildly (<~10%) non-linear, although there are exceptions.

Since this is at some level a convenience given that the smoothing actually happens in pixel space, my proposal would be to do something simple and straightforward like:
1 .compute the delta-wavelength between pixels (e.g. `numpy.diff`)
2. take the median of that to get a pixel/{spectral_axis_unit}.
3. use the astropy spectral equivalency to convert the requested box size to the spectral_axis unit
4. multiple by the result of 2 to get the pixel box size

That's good enough for most "close-to-linear" wavelength axes and is at least not totally meaningless for less-linear axes. We should be sure to note in the docs something like "if your spectral axis is highly non-linear, you may be better off not using this feature" though.

Start at specutils.manipulation.box_smooth and review how its width is currently interpreted in pixels. Trace the spectral-axis unit and spectral-equivalency APIs, then define completion as accepting a spectral-axis width, converting it using the proposed median pixel spacing, and documenting the limitation for highly non-linear axes.

Written by the indexing model from the issue text.