Comments on the procedure

[nsecrest]

1. You can get the variance from the same header data unit (HDU) the flux column comes from. flux_ivar = hdu[1].data['ivar']. Then the standard deviation of the flux is ferr = np.sqrt(1/fivar). You may have certain columns with bad data. In this case, ivar might be zero. I usually do:

wav, flux, flux_ivar = wav[flux_ivar != 0], flux[flux_ivar != 0], flux_ivar[flux_ivar != 0]
flux_err = np.sqrt(1/flux_ivar)

3. Instead of averaging the endpoints, it might be worthwhile to try using the median, which is less sensitive to outliers (like emission lines). Or a weighted average using the above inverse variance.

For the Todos:

I'm not sure anything else in the FITS file, besides the redshift and the variance, will help train the classifier. Most of the other information is instrumental stuff.

Other features we can extract from the data are the strong line fluxes (H-beta, [OIII]5007, H-alpha, [NII]6584, and the [SII] doublet. These lines may be correlated with the presence or non-presence of WR activity. Maybe an ANN can pick out such correlations?