This is an early version, but it seems to be my first succesful use of both pyfits and looping.

% getspec_from_fits.py FeSpesvp8303.fits sp.txt   

FeSpesvp8303.fits = the local FITS file with a spectrum.

% head sp.txt
3553.875       -5.453555
3556.099       -5.984172
3558.322      -19.443817
3560.546        8.104003
3562.770      -11.315107
3564.993      -17.651787
3567.217       -3.473687
3569.440        4.101283
3571.664       15.787663
3573.888       -7.787524

% tail sp.txt
5808.637       -6.461044
5810.860      -35.739651
5813.084       -8.235651
5815.308        0.000000
5817.531        0.000000
5819.755        0.000000
5821.979        0.000000
5824.202        0.000000
5826.426        0.000000
5828.649        0.000000

This appears to be correct!

Since this represents a new era in my python coding, I include the source code for my early version.

#!/usr/bin/python
# 
# Usage: getspec_from_fits.py a_spec_file.fits out.txt

from getRAdeg import convHMS,convDMS
from sys import argv
script_name, fitsname, fout = argv

import numpy as np
from astropy.io import fits as pyfits

# compute the RA,DEC as floating point in DEGREES 
print "Input image = %s" % (fitsname) 

# get header info 
data, header = pyfits.getdata(fitsname, header=True)
naxis1 = header["NAXIS1"]
naxis2 = header["NAXIS2"]
print "NAXIS1,NAXIS2 = %s, %s" % (naxis1,naxis2) 

ref_pixel = header["CRPIX1"]
coord_ref_pixel = header["CRVAL1"]
wave_pixel = header["CDELT1"]
print "ref_pixel,coord_ref_pixel,wave_pixel = %s, %s, %s" % (ref_pixel,coord_ref_pixel,wave_pixel) 

# Open the spectrum FITS file and load data values 
flux = pyfits.open(fitsname)
scidata = flux[0].data
print scidata
print scidata[0,100]

del_wave = float(wave_pixel)
wave0 = float(coord_ref_pixel) 

wave_dat = []
flux_dat = []
# Steve does a loop 
i = 0
while True:
 if i < naxis1:
   aaa = scidata[0,i]
   flux_dat.append(float(aaa)) 
   ww = wave0 + i*del_wave 
   wave_dat.append(ww) 
#  print ww,aaa
 else: 
   break
 i=i+1

# Create the numpy arrays 
xv1 = np.array(wave_dat)
yv1 = np.array(flux_dat)

# Open the for reading
f2 = open(fout,'w')

# Steve does another loop 
i = 0
while True:
 if i < naxis1:
   ww = xv1[i]
   ff = yv1[i] 
   lino = "%.3f  %14.6f\n" % (ww,ff) 
#   print lino
   f2.write(lino)
 else: 
   break
 i=i+1

# close the output file
f2.close()