Implementing AGB models

[s-goldman]

Now that we have the Padova/Colibri models incorporated into the downloadable data #357 , we need some criteria for when to use the models for a given star. At the moment we are not assuming any dust for these models. As has been discussed in #258 we may want to add dust as well as better resolution for thermal pulses.

[s-goldman]

@karllark I'm assuming we need some color or magnitude cut that may depend on the users photometric filter?

[karllark]

Thinking from the physics model standpoint, it would be based on the evolutionary phase or log(g)/Teff or something like that. Mapping the stellar atmosphere models to the tracks/isochrones is done before we have any photometry computed.

[karllark]

Not sure the isochrones we download from the website have the evolutionary phase. This needs to be checked via the file (there are lots of columns!).

[s-goldman]

The current models have ['Z', 'logz', 'logT', 'Teff', 'Reff', 'logg']

[karllark]

But there are is more info (I think) in the downloaded isochrone file.

[s-goldman]

I remember that some models that we messed with last year when we were trying to interpolate the stellar models had "evolutionary points" from Phil Rosenfield that indicated evolutionary phase.

[s-goldman]

The columns for /astro/mboyer/Science/BEAST/Aringer.AGB.grid.fits are the ones I mentioned.

[s-goldman]

The Aringer/Colibri models are also non-uniform with respect to metallicity:

[karllark]

How do the points (say in log(T) and log(g)) compare to the Kurucz stellar models? Do they overlap? Do they extend the log(g)-T(eff) coverage?

[s-goldman]

The Kurucz models are the same across each of the 8 metallicities and look like this:

[s-goldman]

So there is some overlap

[karllark]

Plot?

[s-goldman]

For example: Aringer in red, Kurucz in blue

[karllark]

Great. I think then there is an "easy" way forward. I think we can just add the models to the stellar atmosphere grid. The ordering does matter. The standard is osl = stellib.Tlusty() + stellib.Kurucz() and that means Tlusty models are used in the overlap region. We checked that Tlusty and Kurucz spectra were the same in the overlap region. This should be done for the Kurucz and Aringer overlap regions. Has this been done? Have you done this? We should at least spot check to make sure we've got the same units for the actual spectra used by the BEAST. If that goes well, then I think it would just be osl = stellib.Tlusty() + stellib.Aringer() + stellib.Kurucz(). And then checking everything of course.

[s-goldman]

When you say

We checked that Tlusty and Kurucz spectra were the same in the overlap region.

What do you mean the same?

[karllark]

If you plot the spectra for the same log(g)/Teff/metallicity from Kurucz and Aringer, do they overlap/look the same?

[s-goldman]

Got it, thanks. Is there an easy way to generate the spectra other than through stellib.Kurucz() ?

[karllark]

Not sure. It's been too long. ;-)

[s-goldman]

No worries, I think i've found a way

[s-goldman]

The flux units for Kurucz and Aringer seem different.

[s-goldman]

Two example spectra for T=3800

Kurucz:

[s-goldman]

Aringer:

[s-goldman]

Kurucz fluxes in ergs/s/cm2/A Kurucz wavelengths in nm

Aringer fluxes in (ν · Lν [erg/s]). Aringer wavelengths in Angstrom

[s-goldman]

To get spectra use:

import beast
from beast.physicsmodel.stars import stellib

# get stellib classes
kurucz = stellib.Kurucz()
aringer = stellib.Aringer()

# get spectra
kurucz_spec = kurucz.genSpectrum(np.array([[3800, 3801]]))
aringer_spec = aringer.genSpectrum(np.array([[3800, 3801]]))

# create figure
fig = plt.figure(figsize=(4, 4))
ax1 = fig.add_subplot(1, 1, 1)

# plot spectra
ax1.plot(kurucz.wavelength, kurucz_spec)
ax1.plot(aringer.wavelength / 10, aringer_spec)

# limits
ax1.set_xlim(700, 4000)

# axis labels
ax1.set_ylabel("Flux")
ax1.set_xlabel("Wavelength")

plt.subplots_adjust(wspace=0, hspace=0)
fig.savefig("Model_spectra_comparison.png", dpi=200, bbox_inches="tight")
plt.close()

[s-goldman]

@karllark How should we proceed in converting the fluxes between the spectra and interpolating them to overlapping wavelengths. Do we have functions already in place for this?

[karllark]

Not sure. You can look in the kurucz subdir to see (physicsmodel/stars). The spectra look very different. Are these for the same log(g), Teff?

[s-goldman]

fluxes are normalized and are in ergs/s/cm2/A versus ergs/s

[karllark]

Putting the script into somewhere in physicsmodel/stars would be great. Hopefully named something descriptive.

[karllark]

Even normalized, the shapes are different. And you didn't say if they were for the same log(g),Teff.

[s-goldman]

That might be why they look different. I didn't specify.

[s-goldman]

Specifying g0 doesn't seem to do anything:

kurucz_spec = kurucz.genSpectrum(np.array([[3800, 3801]]), g0=0)

vs

kurucz_spec = kurucz.genSpectrum(np.array([[3800, 3801]]), g0=5)

[s-goldman]

Need to divide aringer fluxes by the associated frequency (specific luminosity * frequency -- > specific luminosity) and then get them in the form ergs/s/cm2/A. THEN need to interpolate those values to the wavelength grid of the Kurucz models AND extrapolate them toward bluer wavelengths (maybe using a blackbody).

These two steps are important otherwise the composite spectral grid will include both sets of wavelength points. @lcjohnso did something similar for BOSZ in IDL in #180 .

A function for doing this should be created and added to a file in physicsmodel/stars.

[s-goldman]

Aringer models here

[s-goldman]

I've parsed the COMARCS models into a useable form, however the flux conversion is tricky.

The COMARCS (AGB) models are in specific luminosity νL(ν) (ergs/s).

The models need to be in (erg/s/cm2/Å) for the BEAST.

I can use the specific luminosity to get a flux (L=4πr^2) and doing some conversions, but for this we need to assume a distance.

How does the BEAST currently scale model fluxes by distance?

[karllark]

Distance is handled by the BEAST in a different part. If you can get the COMARCS models into surface flux units like the other stellar atmosphere models (e.g., ergs/s/cm2/A), then this will be handled. In this case, the cm2 is based on the radius of the star. At least if I understand things correctly.

[karllark]

The conversion from the stellar atmosphere surface flux is done when the stellar evolutionary results are combined with the stellar atmospheres. This is where the radius is used to compute the luminosity and place star at a default distance. Then the actual distance is applied to get fluxes. Just a FYI.

[s-goldman]

Okay great, and we can get radius with logg. Do you know how or in what file the surface flux is scaled by distance (just curious at this point).

[karllark]

Somewhere in the creation of the spectral grid (not SED grid). To my memory.

[s-goldman]

Okay it looks like the SEDs are scaled to 10pc and then again by the distance in apply_distance_and_spectral_props.

[karllark]

That sound like what I remember. ;-)

[marthaboyer]

Here is the plot I made when I added the models. Looks the like overlap between Tlusty, Kurucz, and Aringer is as expected, so units should be OK... I think.

check.pdf

[karllark]

@marthaboyer : your plot of the overlap of between Aringer and Kurucz models in Teff/logg space is different than Steve's earlier this this issue. Maybe "cleaning" the Aringer models has already been done?

[marthaboyer]

It looks like Steve plotted just one metallicity at a time. The weird extra evolutionary tracks are for just 2 metallicities, so they'll only show up if you plot one of those. I'm still working on cleaning them out. It might take a week or two (in the middle of a JWST rehearsal right now), but here's my todo list:

1. finish extrapolating to lower wavelengths (just about complete)
2. take out those extra evolutionary tracks in the grid.
3. fix the metallicities listed in the headers (I found an error in the metallicity currently included).

[karllark]

Totally awesome. Thanks!