Load mass and radius from the main data table for the famous models of
    Baraffe+ (1998A&A...337..403B).

    tablelines
      An iterable yielding lines from the table data file.
      I've named the file '1998A&A...337..403B_tbl1-3.dat'
      in some repositories (it's about 150K, not too bad).
    metallicity
      The metallicity of the model to select.
    heliumfrac
      The helium fraction of the model to select.
    age_gyr
      The age of the model to select, in Gyr.
    age_tol
      The tolerance on the matched age, in Gyr.

    Returns: (mass, radius), where both are Numpy arrays.

    The ages in the data table vary slightly at fixed metallicity and helium
    fraction. Therefore, there needs to be a tolerance parameter for matching
    the age.

def load_bcah98_mass_radius (tablelines, metallicity=0, heliumfrac=0.275,
                             age_gyr=5., age_tol=0.05):
    """Load mass and radius from the main data table for the famous models of
    Baraffe+ (1998A&A...337..403B).

    tablelines
      An iterable yielding lines from the table data file.
      I've named the file '1998A&A...337..403B_tbl1-3.dat'
      in some repositories (it's about 150K, not too bad).
    metallicity
      The metallicity of the model to select.
    heliumfrac
      The helium fraction of the model to select.
    age_gyr
      The age of the model to select, in Gyr.
    age_tol
      The tolerance on the matched age, in Gyr.

    Returns: (mass, radius), where both are Numpy arrays.

    The ages in the data table vary slightly at fixed metallicity and helium
    fraction. Therefore, there needs to be a tolerance parameter for matching
    the age.

    """
    mdata, rdata = [], []

    for line in tablelines:
        a = line.strip ().split ()

        thismetallicity = float (a[0])
        if thismetallicity != metallicity:
            continue

        thisheliumfrac = float (a[1])
        if thisheliumfrac != heliumfrac:
            continue

        thisage = float (a[4])
        if abs (thisage - age_gyr) > age_tol:
            continue

        mass = float (a[3]) * cgs.msun
        teff = float (a[5])
        mbol = float (a[7])

        # XXX to check: do they specify m_bol_sun = 4.64? IIRC, yes.
        lbol = 10**(0.4 * (4.64 - mbol)) * cgs.lsun
        area = lbol / (cgs.sigma * teff**4)
        r = np.sqrt (area / (4 * np.pi))

        mdata.append (mass)
        rdata.append (r)

    return np.asarray (mdata), np.asarray (rdata)