Flux with electron/s cannot be read back in from file (roundtrip problem)

Consider this use case:

```python
import numpy as np
from astropy import units as u
from specutils import Spectrum1D

flux = np.random.random(10) * (u.electron / u.s)
wave = (np.arange(10) + 1) * u.m
sp = Spectrum1D(flux=flux, spectral_axis=wave) # Works fine
```

But then when writing it back out into a FITS file:

```python
>>> sp.write('flux_es.fits')
WARNING: The unit 'electron / s' could not be saved in native FITS format
and cannot be recovered in reading. It can roundtrip within astropy by using
QTable both to write and read back, though one has to enable the unit before
reading. [astropy.io.fits.convenience]
```

Now, that file cannot be read back in, as the warning indicated. But why can't `specutils` write it out in a way that can roundtrip?

```
>>> sp2 = Spectrum1D.read('flux_es.fits')
.../astropy/nddata/mixins/ndio.py:59, in NDDataRead.__call__(self, *args, **kwargs)
58 def __call__(self, *args, **kwargs):
---> 59 return self.registry.read(self._cls, *args, **kwargs)

.../astropy/io/registry/core.py:219, in UnifiedInputRegistry.read(self, cls, format, cache, *args, **kwargs)
214 format = self._get_valid_format(
215 "read", cls, path, fileobj, args, kwargs
216 )
218 reader = self.get_reader(format, cls)
--> 219 data = reader(*args, **kwargs)
221 if not isinstance(data, cls):
222 # User has read with a subclass where only the parent class is
223 # registered. This returns the parent class, so try coercing
224 # to desired subclass.
225 try:

.../specutils/io/default_loaders/tabular_fits.py:86, in tabular_fits_loader(file_obj, column_mapping, hdu, **kwargs)
83 # If no column mapping is given, attempt to parse the file using
84 # unit information
85 if column_mapping is None:
---> 86 return generic_spectrum_from_table(tab, wcs=wcs, **kwargs)
88 return spectrum_from_column_mapping(tab, column_mapping, wcs=wcs)

.../specutils/io/parsing_utils.py:258, in generic_spectrum_from_table(table, wcs, **kwargs)
256 flux_column = _find_spectral_column(table, colnames, spectral_axis)
257 if flux_column is None:
--> 258 raise IOError("Could not identify column containing the flux")
259 flux = table[flux_column].to(table[flux_column].unit)
260 colnames.remove(flux_column)

OSError: Could not identify column containing the flux
```

Desired behavior: `specutils` should know how to workaround this limitation for proper I/O roundtripping of its own product.

[🐱](https://jira.stsci.edu/browse/JDAT-1651)

Reproduce the example with electron-per-second flux, then inspect specutils/io/default_loaders/tabular_fits.py and specutils/io/parsing_utils.py along the shown read path. Compare the write warning with the generated table and determine how the reader identifies flux; done means a Spectrum1D written by specutils can be read back with its flux column and unit preserved.

Written by the indexing model from the issue text.