False PHOENIX model spectra for Sun-like stars

[MNGuenther]

Either Pysynphot reads the wrong PHOENIX spectra, interpolates them incorrectly, or there is something wrong with PHOENIX itself. The spectra do not match what is expected at the given temperatures. This seems to only affect 'even-hundred' temperatures, the 'odd-hundred' ones seem to work correctly.

pysynphot_phoenix_error_even_temperatures.pdf

pysynphot_phoenix_error_odd_temperatures.pdf

I used the following code to generate these plots:

import numpy as np
import matplotlib.pyplot as plt
import pysynphot as S

#::: odd numbers
plt.figure()
for teff in np.arange(5900,3900,-200):
    sp = S.Icat('phoenix', teff, 0., 4.5)
    plt.plot(sp.wave, sp.flux, label=str(teff))
plt.xlim(3000,10000)
plt.xlabel(r'$T_{eff}$ (K)')
plt.ylabel('flam')
plt.legend(bbox_to_anchor=(1.04, 0.2),  borderaxespad=0)
plt.subplots_adjust(right=0.7)
plt.savefig('pysynphot_phoenix_error_odd_temperatures.pdf')

#::: even numbers
plt.figure()
for teff in np.arange(6000,3900,-200):
    sp = S.Icat('phoenix', teff, 0., 4.5)
    plt.plot(sp.wave, sp.flux, label=str(teff))
plt.xlim(3000,10000)
plt.xlabel(r'$T_{eff}$ (K)')
plt.ylabel('flam')
plt.legend(bbox_to_anchor=(1.04, 0.2),  borderaxespad=0)
plt.subplots_adjust(right=0.7)
plt.savefig('pysynphot_phoenix_error_even_temperatures.pdf')

[pllim]

Just to make sure, what is your version of PySynphot?

[pllim]

And to clarify, what do you expect for 'even-hundred' temperatures that you are not seeing?

[MNGuenther]

Hi pllim, thanks for your swift response!

pysynphot version 0.9.7

I expect the 'even-hundred' to smoothly increase in total flux with increasing temperature, just as the 'odd-hundred' do. But particularly 5600K and 5800K look significantly different. Best visible here:

pysynphot_phoenix_error_total_flux.pdf

Figure created with this code:

Teff = np.arange(6000,3900,-100) 
bp = S.UniformTransmission(1.)

def get_total_flux(teff):
    sp = S.Icat('phoenix', teff, 0., 4.5)
    obs = S.Observation(sp, bp, binset=sp.wave)
    return obs.countrate()

total_flux = [ get_total_flux(teff) for teff in Teff ]

plt.figure()
plt.plot(Teff, total_flux, 'k.')
plt.xlabel(r'$T_{eff}$ (K)')
plt.ylabel('Total counts')
plt.xticks(range(4000,6001,200))
plt.savefig('pysynphot_phoenix_error_total_flux.pdf')

[pllim]

Thank you for the info. 0.9.7 is not the latest release but it should not matter for this case. Just making sure you are not too far behind. I'll investigate and get back to you.

[pllim]

@MNGuenther , after some investigation, this turns out to be a problem with the data itself and not the software. In fact, the parameters you gave did not require interpolation at all, which means the flux you are seeing are the exact values from the data grid itself. When I play the flux through your range of temperatures, the decline only affects Teff at 5600 K and 5800 K.

Here is how you can visualize it as I did:

# This example is only compatible with Python 3, not Python 2.

import os
import matplotlib.pyplot as plt
import pysynphot as S
from astropy.table import Table  # You have to install astropy

def play():
    fig = plt.figure()
    ax = fig.add_subplot(111)

    # This shows flux decline at Teff 5600 and 5800
    for teff in range(3900, 6201, 100):
        filename = os.path.join(
            S.locations.rootdir, 'grid', 'phoenix', 'phoenixm00',
            'phoenixm00_{}.fits'.format(teff))
        tab = Table.read(filename)
        wave = tab['WAVELENGTH']
        flux = tab['g45']

        ax.plot(wave, flux, label='Teff = {}'.format(teff))
        ax.set_xlim(3000, 10000)
        ax.set_ylim(0, 1.8e7)
        ax.set_xlabel('Angstrom')
        ax.set_ylabel('FLAM')
        ax.legend()
        plt.show()

        input('Enter to continue: ')  # Pause

        ax.clear()

play()

@rizeladiaz -- Can you please clarify why the data is this way?

@mgennaro -- Another Phoenix models quirk you might want to know about.

[pllim]

@MNGuenther , this problem might be caused by data obtained upstream from Allard et al. @rizeladiaz is investigating. Thank you for your patience.

[MNGuenther]

@pllim @rizeladiaz many thanks.

Could you clarify from which publication the PHOENIX models are drawn? If looking at the raw data of the recent publication in Husser et al. 2013, it seems that their data for Sun-like stars is consistent with expectations/observations.

Might it hence be an issue in older PHOENIX models?

[dhomeier]

Summarising some findings from an email conversation, I take it the data used for the included ‘phoenix’ grid are pulled from ftp://ftp.stsci.edu/cdbs/grid/phoenix/phoenixm00

Their origin in turn is documented somewhat fuzzily, but they were clearly pulled from the Lyon simulator; from the date and metallicity coverage they can only be from the PHOENIX AGSS2009 grid which is actually the first ‘GAIA-Dusty’ grid I had calculated at the University of Göttingen. Since the files were renamed when uploaded to the Lyon server, and the CDBS FITS headers don’t give any more specific file identification, this is to some extent inference, but I can’t think of any other possible source. However the spectra themselves still available on the Lyon simulator don’t show any such discrepancies - see attached for a comparison of the respective CDBS log g = 4.5 columns and the Lyon/Göttingen spectra, the latter smoothed with a 100-pixel boxcar: sol_5700_pysynphot.pdf sol_5800_pysynphot.pdf

I haven't inspected any of the other "even Teff" data, but in the phoenixm00_5800.fits file, the g40 and g50 columns appear OK in comparison.

I recommend that all spectra at least be checked by comparing their integrated fluxes against constants.sigma_sb * T_eff**4.

[pllim]

Thanks, @dhomeier ! The last I heard, RedCat Team (website: http://www.stsci.edu/hst/observatory/crds) is working on updating the Phoenix models that STScI distributes.

[dhomeier]

Thanks for the link, @pllim - I couldn't find anything related to stellar models on a quick glance, but will keep this in mind. Is @rizeladiaz still involved with the archiving? On a closer look at higher resolution, there still seem to be some discrepancies even with models that overall match the SED (like 5700_g45), and I am not sure if they can be explained by downsampling artefacts.

[pllim]

@dhomeier , AFAIK yes. Unfortunately, PySynphot package is not responsible for the data themselves. Please contact RedCat Team (contact info on the given website) for further discussions (feel free to c/c me if you have my work email).