sfh2exp.py 3.85 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
from collections import OrderedDict
16

17
import numpy as np
18

19
from . import SedModule
20

21

22
class Sfh2Exp(SedModule):
Yannick Roehlly's avatar
Yannick Roehlly committed
23
    """Double decreasing exponential Star Formation History
24

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

    """

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

69
70
71
72
73
    def _init_code(self):
        self.tau_main = float(self.parameters["tau_main"])
        self.tau_burst = float(self.parameters["tau_burst"])
        self.f_burst = float(self.parameters["f_burst"])
        self.burst_age = int(self.parameters["burst_age"])
74
        age = int(self.parameters["age"])
75
76
        sfr_0 = float(self.parameters["sfr_0"])
        normalise = bool(self.parameters["normalise"])
77
78

        # Time grid and age. If needed, the age is rounded to the inferior Myr
79
        time_grid = np.arange(age)
80
        time_grid_burst = np.arange(self.burst_age)
81

82
        # SFR for each component
83
        self.sfr = np.exp(-time_grid / self.tau_main)
84
        sfr_burst = np.exp(-time_grid_burst / self.tau_burst)
85
86

        # Height of the late burst to have the desired produced mass fraction
87
        sfr_burst *= self.f_burst / (1.-self.f_burst) * np.sum(self.sfr) / np.sum(sfr_burst)
88
89
90

        # We add the age burst exponential for ages superior to age -
        # burst_age
91
        self.sfr[-(time_grid_burst[-1]+1):] += sfr_burst
92

93
94
        # Compute the galaxy mass and normalise the SFH to 1 solar mass
        # produced if asked to.
95
        self.sfr_integrated = np.sum(self.sfr) * 1e6
96
        if normalise:
97
98
            self.sfr /= self.sfr_integrated
            self.sfr_integrated = 1.
99
        else:
100
            self.sfr *= sfr_0
101
            self.sfr_integrated *= sfr_0
102
103
104
105
106
107
108
109
110

    def process(self, sed):
        """Add a double decreasing exponential Star Formation History.

        Parameters
        ----------
        sed: pcigale.sed.SED object

        """
111

112
        sed.add_module(self.name, self.parameters)
113

114
        # Add the sfh and the output parameters to the SED.
115
        sed.sfh = self.sfr
116
        sed.add_info("sfh.integrated", self.sfr_integrated, True)
117
118
119
120
        sed.add_info("sfh.tau_main", self.tau_main)
        sed.add_info("sfh.tau_burst", self.tau_burst)
        sed.add_info("sfh.f_burst", self.f_burst)
        sed.add_info("sfh.burst_age", self.burst_age)
Yannick Roehlly's avatar
Yannick Roehlly committed
121

122
# SedModule to be returned by get_module
Yannick Roehlly's avatar
Yannick Roehlly committed
123
Module = Sfh2Exp