Script to convert Kurucz ODFNEW model packages to HDF5

extract_kurucz.py

import io

import numpy as np
import atpy


def chunk(string, n):
    return [string[i*n:i*n+n] for i in range(len(string) / n)]

# Open file
f = open('fp00k2odfnew.pck', 'rb')

# Replace strange line endings
f = io.StringIO(f.read().replace('\r', '\n'))

# Skip header
while f.readline().strip() != "END":
    pass

# Read in wavelengths
n_wav = 1221
wav = []
while len(wav) < n_wav:
    wav += f.readline().strip().split()
wav = np.array(wav, dtype=np.float32)

models = []

for i in range(476):

    # Read in header line
    header = f.readline()

    # Extract information from header
    _, teff, _, logg, _ = header.split(None, 4)
    teff, logg = float(teff), float(logg)

    # Read model flux
    fnu = []
    while len(fnu) < n_wav:
        fnu += chunk(f.readline().strip(), 10)
    fnu = np.array(fnu, dtype=np.float32)

    # Read continuum flux
    fnu_cont = []
    while len(fnu_cont) < n_wav:
        fnu_cont += chunk(f.readline().strip(), 10)
    fnu_cont = np.array(fnu_cont, dtype=np.float32)

    # Create table to contain spectrum
    t = atpy.Table()

    # Write parameters to header
    t.keywords['teff'] = teff
    t.keywords['logg'] = logg

    # Add wavelength/frequency info
    t.add_column('wav', wav)
    t.add_column('nu', 3e8 / (wav * 1.e-9))

    # Add flux columns
    t.add_column('fnu', fnu)
    t.add_column('fnu_cont', fnu_cont)

    # Write to HDF5 file
    t.write('models/kt%05ig+%3.1f.hdf5' % (teff, logg), verbose=False)