bc03.py 5.24 KB
Newer Older
1
# -*- coding: utf-8 -*-
2 3
# Copyright (C) 2013 Centre de données Astrophysiques de Marseille
# Licensed under the CeCILL-v2 licence - see Licence_CeCILL_V2-en.txt
Yannick Roehlly's avatar
Yannick Roehlly committed
4
# Author: Yannick Roehlly
5 6 7 8

import numpy as np


9
class BC03(object):
10 11 12
    """Single Stellar Population as defined in Bruzual and Charlot (2003)

    This class holds the data associated with a single stellar population
13 14 15 16 17 18 19
    (SSP) as defined in Bruzual and Charlot (2003). Compared to the pristine
    Bruzual and Charlot SSP:
        * The time grid ranges from 1 Myr to 14 Gyr with 1 Myr steps.
        * The SSP are all interpolated on this new grid.
        * The wavelength grid is refined beyond 10 μm to avoid artefacts.
        * The wavelength is given in nm rather than Å.
        * The spectra are given in W/nm rather than Lsun.
20 21 22 23

    """

    def __init__(self, imf, metallicity, time_grid, wavelength_grid,
24
                 info_table, spec_table):
25 26 27 28 29
        """Create a new single stellar population as defined in Bruzual and
        Charlot (2003).

        Parameters
        ----------
30
        imf: string
31 32
            Initial mass function (IMF): either 'salp' for Salpeter (1955) or
            'chab' for Chabrier (2003).
33
        metallicity: float
34
            The metallicity. Possible values are 0.0001, 0.0004, 0.004, 0.008,
35
            0.02, and 0.05.
36
        time_grid: array of floats
37
            The time grid in Myr used in the info_table and the spec_table.
38
        wavelength_grid: array of floats
39 40
            The wavelength grid in nm used in spec_table.
        info_table: 2 axis array of floats
41 42
            Array containing information from some of the *.?color tables from
            Bruzual and Charlot (2003) at each time of the time_grid.
43 44 45 46 47 48 49
                * info_table[0]: Total mass in stars in solar mass
                * info_table[1]: Mass returned to the ISM by evolved stars in
                    solar mass
                * info_table[2]: rate of H-ionizing photons (s-1)
        spec_table: 2D array of floats
            Spectrum of the SSP in W/nm (first axis) every 1 Myr (second axis).

50 51 52 53 54 55 56 57 58 59
        """

        if imf in ['salp', 'chab']:
            self.imf = imf
        else:
            raise ValueError('IMF must be either sal for Salpeter or '
                             'cha for Chabrier.')
        self.metallicity = metallicity
        self.time_grid = time_grid
        self.wavelength_grid = wavelength_grid
60 61
        self.info_table = info_table
        self.spec_table = spec_table
62

63
    def convolve(self, sfh, separation_age):
64 65
        """Convolve the SSP with a Star Formation History

66
        Given an SFH, this method convolves the info table and the SSP
67
        luminosity spectrum.
68 69 70

        Parameters
        ----------
71 72
        sfh: array of floats
            Star Formation History in Msun/yr.
73 74
        separation_age: float
            Age separating the young from the old stellar populations in Myr.
75 76 77

        Returns
        -------
78 79 80 81 82 83 84
        spec_young: array of floats
            Spectrum in W/nm of the young stellar populations.
        spec_old: array of floats
            Same as spec_young but for the old stellar populations.
        info_young: dictionary
            Dictionary containing various information from the *.?color tables
            for the young stellar populations:
85 86 87
            * "m_star": Total mass in stars in Msun
            * "m_gas": Mass returned to the ISM by evolved stars in Msun
            * "n_ly": rate of H-ionizing photons (s-1)
88 89
        info_old : dictionary
            Same as info_young but for the old stellar populations.
90 91 92
        info_all: dictionary
            Same as info_young but for the entire stellar population. Also
            contains "age_mass", the stellar mass-weighted age
Yannick Roehlly's avatar
Yannick Roehlly committed
93

94
        """
95
        # We cut the SSP to the maximum age considered to simplify the
96 97 98
        # computation. We take only the first three elements from the
        # info_table as the others do not make sense when convolved with the
        # SFH (break strength).
99
        info_table = self.info_table[:, :sfh.size]
100
        spec_table = self.spec_table[:, :sfh.size]
101

102 103 104 105 106 107 108 109 110 111 112 113 114 115
        # The convolution is just a matter of reverting the SFH and computing
        # the sum of the data from the SSP one to one product. This is done
        # using the dot product. The 1e6 factor is because the SFH is in solar
        # mass per year.
        info_young = 1e6 * np.dot(info_table[:, :separation_age],
                                  sfh[-separation_age:][::-1])
        spec_young = 1e6 * np.dot(spec_table[:, :separation_age],
                                  sfh[-separation_age:][::-1])

        info_old = 1e6 * np.dot(info_table[:, separation_age:],
                                sfh[:-separation_age][::-1])
        spec_old = 1e6 * np.dot(spec_table[:, separation_age:],
                                sfh[:-separation_age][::-1])

116 117
        info_all = info_young + info_old

118 119
        info_young = dict(zip(["m_star", "m_gas", "n_ly"], info_young))
        info_old = dict(zip(["m_star", "m_gas", "n_ly"], info_old))
120 121 122 123
        info_all = dict(zip(["m_star", "m_gas", "n_ly"], info_all))

        info_all['age_mass'] = np.average(self.time_grid[:sfh.size],
                                          weights=info_table[0, :] * sfh[::-1])
124

125
        return spec_young, spec_old, info_young, info_old, info_all