sfh2exp.py 4.21 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

Yannick Roehlly's avatar
Yannick Roehlly committed
6
7
8
9
10
11
12
13
14
"""
Double decreasing exponential star formation history module
===========================================================

This module implements a star formation history (SFH) composed of two
decreasing exponentials.

"""

15
import numpy as np
16
from collections import OrderedDict
17
from . import CreationModule
18

Yannick Roehlly's avatar
Yannick Roehlly committed
19
# Time lapse used in the age grid in Myr. If should be consistent with the
20
# time lapse in the SSP modules.
Yannick Roehlly's avatar
Yannick Roehlly committed
21
AGE_LAPSE = 1
22
23


Yannick Roehlly's avatar
Yannick Roehlly committed
24
class Sfh2Exp(CreationModule):
Yannick Roehlly's avatar
Yannick Roehlly committed
25
    """Double decreasing exponential Star Formation History
26

Yannick Roehlly's avatar
Yannick Roehlly committed
27
28
    This module sets the SED star formation history (SFH) as a combination of
    two exp(-t/τ) exponentials.
29
30
31

    """

32
33
    parameter_list = OrderedDict([
        ("tau_main", (
34
            "float",
Yannick Roehlly's avatar
Yannick Roehlly committed
35
            "e-folding time of the main stellar population model in Myr.",
36
            6000.
37
38
        )),
        ("tau_burst", (
39
            "float",
Yannick Roehlly's avatar
Yannick Roehlly committed
40
            "e-folding time of the late starburst population model in Myr.",
41
            50.
42
43
        )),
        ("f_burst", (
44
45
            "float",
            "Mass fraction of the late burst population.",
46
            0.01
47
48
        )),
        ("age", (
49
            "integer",
50
51
52
            "Age of the main stellar population in the galaxy in Myr."
            "The precision is 1 Myr.",
            5000.
53
54
        )),
        ("burst_age", (
55
            "integer",
Yannick Roehlly's avatar
Yannick Roehlly committed
56
            "Age of the late burst in Myr. Precision is 1 Myr.",
57
            20.
Denis's avatar
Denis committed
58
59
        )),
         ("sfr_0", (
60
            "float",
61
            "Value of SFR at t = 0 in M_sun/yr.",
62
63
64
65
66
67
            1.
        )),
        ("normalise", (
            "boolean",
            "Normalise the SFH to produce one solar mass.",
            "False"
68
        )),
69
    ])
70

71
72
73
74
75
76
    out_parameter_list = OrderedDict([
        ("tau_main", "e-folding time of the main stellar population model "
                     "in Myr."),
        ("tau_burst", "e-folding time of the late starburst population model "
                      "in Myr."),
        ("f_burst", "Produced mass fraction of the late burst population."),
77
        ("age", "Age of the main stellar population in the galaxy in Myr."),
Denis's avatar
Denis committed
78
        ("burst_age", "Age of the late burst in Myr."),
79
        ("sfr_0", "SFR at t = 0 in M_sun/yr.")
80
    ])
81

82
    def process(self, sed):
83
84
85
86
        """Add a double decreasing exponential Star Formation History.

        Parameters
        ----------
87
        sed: pcigale.sed.SED object
88
89

        """
90
91
92
93
94
        tau_main = float(self.parameters["tau_main"])
        tau_burst = float(self.parameters["tau_burst"])
        f_burst = float(self.parameters["f_burst"])
        age = int(self.parameters["age"])
        burst_age = int(self.parameters["burst_age"])
95
96
        sfr_0 = int(self.parameters["sfr_0"])
        normalise = (self.parameters["normalise"].lower() == "true")
97
98
99
100
101
102

        # Time grid and age. If needed, the age is rounded to the inferior Myr
        time_grid = np.arange(AGE_LAPSE, age + AGE_LAPSE, AGE_LAPSE)
        age = np.max(time_grid)

        # Main exponential
103
        sfr = sfr_0 * np.exp(-time_grid / tau_main)
104
105
106

        # Height of the late burst to have the desired produced mass fraction
        # (assuming that the main burst as a height of 1).
107
108
        burst_height = (f_burst/(1-f_burst) * tau_main/tau_burst *
            (1-np.exp(-age/tau_main))/(1-np.exp(-burst_age/tau_burst)))
109
110
111
112
113
114
115

        # We add the age burst exponential for ages superior to age -
        # burst_age
        mask = (time_grid >= (age - burst_age))
        sfr[mask] = sfr[mask] + burst_height * np.exp(
            (-time_grid[mask] + age - burst_age) / tau_burst)

116
117
118
        # Normalise the SFH to 1 solar mass produced if asked to.
        if normalise:
            sfr = sfr / np.trapz(sfr * 1e6, time_grid)
119

120
        sed.add_module(self.name, self.parameters)
121

122
        # Add the sfh and the output parameters to the SED.
Yannick Roehlly's avatar
Yannick Roehlly committed
123
        sed.sfh = (time_grid, sfr)
124
125
126
127
        sed.add_info("sfh.tau_main", tau_main)
        sed.add_info("sfh.tau_burst", tau_burst)
        sed.add_info("sfh.f_burst", f_burst)
        sed.add_info("sfh.burst_age", burst_age)
Yannick Roehlly's avatar
Yannick Roehlly committed
128
129
130

# CreationModule to be returned by get_module
Module = Sfh2Exp