Commit bfcea745 authored by Laure Ciesla's avatar Laure Ciesla Committed by Yannick Roehlly
Browse files

Implementaiton of the Dale+14 model (AGN+SF)

parent b5622495
...@@ -20,7 +20,7 @@ import numpy as np ...@@ -20,7 +20,7 @@ import numpy as np
from scipy import interpolate from scipy import interpolate
import scipy.constants as cst import scipy.constants as cst
from pcigale.data import (Database, Filter, SspM2005, SspBC03, AgnFritz2006, from pcigale.data import (Database, Filter, SspM2005, SspBC03, AgnFritz2006,
DL2007) DALE2014, DL2007)
def read_bc03_ssp(filename): def read_bc03_ssp(filename):
...@@ -319,6 +319,67 @@ def build_dh2002(base): ...@@ -319,6 +319,67 @@ def build_dh2002(base):
base.add_dh2002_infrared_templates(data) base.add_dh2002_infrared_templates(data)
def build_dale2014(base):
dh2002_dir = os.path.join(os.path.dirname(__file__), 'dh2002/')
dale2014_dir = os.path.join(os.path.dirname(__file__), 'dale2014/')
# Getting the alpha grid for the templates
d14cal = np.genfromtxt(dh2002_dir + 'dhcal.dat')
alpha_grid = d14cal[:, 1]
# Getting the lambda grid for the templates and convert from microns to nm.
first_template = np.genfromtxt(dale2014_dir + 'spectra.0.00AGN.dat')
wave = first_template[:, 0] * 1E3
# Getting the stellar emission and interpolate it at the same wavelength grid
stell_emission_file = np.genfromtxt(dale2014_dir + 'stellar_SED_age13Gyr_tau10Gyr.spec')
# A -> to nm
wave_stell = stell_emission_file[:,0] * 0.1
# W/A -> W/nm
stell_emission = stell_emission_file[:,1] * 10
stell_emission_interp = np.interp(wave,wave_stell,stell_emission)
# The models are in nuFnu and contain stellar emission.
# We convert this to W/nm and remove the stellar emission.
# Emission from dust heated by SB
fraction = 0.0
filename = dale2014_dir + "spectra.0.00AGN.dat"
print("Importing {}...".format(filename))
datafile = open(filename)
data = "".join(datafile.readlines())
datafile.close()
for al in range(1,len(alpha_grid),1):
lumin_with_stell = np.genfromtxt(io.BytesIO(data.encode()), usecols=(al))
lumin_with_stell = pow(10,lumin_with_stell) / wave
constant = lumin_with_stell[7] / stell_emission_interp[7]
lumin = lumin_with_stell - stell_emission_interp * constant
lumin[lumin<0] = 0
lumin[wave<2E3] = 0
norm = np.trapz(lumin, x = wave)
lumin = lumin/norm
base.add_dale2014(DALE2014(fraction, alpha_grid[al-1], wave, lumin))
# Emission from dust heated by AGN - Quasar template
fraction = 1.0
filename = dale2014_dir + "spectra.1.00AGN.dat"
print("Importing {}...".format(filename))
datafile = open(filename)
data = "".join(datafile.readlines())
datafile.close()
for al in range(1,len(alpha_grid),1):
lumin_quasar = np.genfromtxt(io.BytesIO(data.encode()), usecols=(al))
lumin_quasar = pow(10,lumin_quasar) / wave
lumin_quasar[lumin_quasar<0] = 0
norm = np.trapz(lumin_quasar, x = wave)
lumin_quasar = lumin_quasar/norm
base.add_dale2014(DALE2014(fraction, alpha_grid[al-1], wave, lumin_quasar))
def build_dl2007(base): def build_dl2007(base):
dl2007_dir = os.path.join(os.path.dirname(__file__), 'dl2007/') dl2007_dir = os.path.join(os.path.dirname(__file__), 'dl2007/')
...@@ -457,6 +518,11 @@ def build_base(): ...@@ -457,6 +518,11 @@ def build_base():
print("\nDONE\n") print("\nDONE\n")
print('#' * 78) print('#' * 78)
print("7- Importing Dale et al (2014) templates\n")
build_dale2014(base)
print("\nDONE\n")
print('#' * 78)
base.session.close_all() base.session.close_all()
......
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
...@@ -24,6 +24,7 @@ from .filters import Filter ...@@ -24,6 +24,7 @@ from .filters import Filter
from .ssp_m2005 import SspM2005 from .ssp_m2005 import SspM2005
from .ssp_bc03 import SspBC03 from .ssp_bc03 import SspBC03
from .ir_templates_dh2002 import IrTemplatesDH2002 from .ir_templates_dh2002 import IrTemplatesDH2002
from .ir_agn_templates_dale2014 import DALE2014
from .ir_models_dl2007 import DL2007 from .ir_models_dl2007 import DL2007
from .agn_fritz2006 import AgnFritz2006 from .agn_fritz2006 import AgnFritz2006
...@@ -118,6 +119,21 @@ class _DH2002InfraredTemplates(BASE): ...@@ -118,6 +119,21 @@ class _DH2002InfraredTemplates(BASE):
self.description = description self.description = description
self.data = data self.data = data
class _DALE2014(BASE):
"""Storage for Dale et al (2014) infra-red templates
"""
__tablename__ = 'dale2014_templates'
fracAGN = Column(Float, primary_key=True)
alpha = Column(String, primary_key=True)
wave = Column(PickleType)
lumin = Column(PickleType)
def __init__(self, iragn):
self.fracAGN = iragn.fracAGN
self.alpha = iragn.alpha
self.wave = iragn.wave
self.lumin = iragn.lumin
class _DL2007(BASE): class _DL2007(BASE):
"""Storage for Draine and Li (2007) IR models """Storage for Draine and Li (2007) IR models
...@@ -288,6 +304,28 @@ class Database(object): ...@@ -288,6 +304,28 @@ class Database(object):
raise StandardError('The template is already in the base.') raise StandardError('The template is already in the base.')
else: else:
raise StandardError('The database is not writable.') raise StandardError('The database is not writable.')
def add_dale2014(self, iragn):
"""
Add Dale et al (2014) templates the collection.
Parameters
----------
iragn : pcigale.data.DALE2014
"""
if self.is_writable:
template = _DALE2014(iragn)
self.session.add(template)
try:
self.session.commit()
except exc.IntegrityError:
self.session.rollback()
raise StandardError('The DALE2014 template is already in the base.')
else:
raise StandardError('The database is not writable.')
def add_dl2007(self, model): def add_dl2007(self, model):
""" """
...@@ -431,6 +469,29 @@ class Database(object): ...@@ -431,6 +469,29 @@ class Database(object):
else: else:
return None return None
def get_dale2014(self, fracAGN, alpha):
"""
Get the Dale et al (2014) template corresponding to the given set of
parameters.
Parameters
----------
fracAGN: float
contribution of the AGN to the IR luminosity
alpha: float
alpha corresponding to the updated Dale & Helou (2002) star forming template.
"""
result = (self.session.query(_DALE2014).
filter(_DALE2014.fracAGN == fracAGN).
filter(_DALE2014.alpha == alpha).
first())
if result:
return DALE2014(result.fracAGN, result.alpha, result.wave,
result.lumin)
else:
return None
def get_agn_fritz2006(self, model_nb): def get_agn_fritz2006(self, model_nb):
""" """
Get the Fritz et al. (2006) AGN model corresponding to the number. Get the Fritz et al. (2006) AGN model corresponding to the number.
......
## -*- coding: utf-8 -*-
# Copyright (C) 2013 Centre de données Astrophysiques de Marseille
# Licensed under the CeCILL-v2 licence - see Licence_CeCILL_V2-en.txt
# Author: Yannick Roehlly <yannick.roehlly@oamp.fr>, Laure Ciesla <ciesla@physics.uoc.gr>
class DALE2014(object):
"""Dale et al (2014) IR templates containing an AGN component.
This class holds the data associated with the Dale et al (2014)
dust templates.
"""
def __init__(self, fracAGN, alpha, wave, lumin):
"""Create a new IR model
Parameters
----------
fracAGN: float
Contribution of the AGN
alpha: float
Dale & Helou (2002) alpha slope.
lir: float
Total IR luminosity between 8 and 1000 microns (AGN + SB)
wave: array
Vector of the λ grid used in the templates [nm]
lumin: array
Model data in an array containing the luminosity density
versus the wavelength λ
"""
self.fracAGN = fracAGN
self.alpha = alpha
self.wave = wave
self.lumin = lumin
\ No newline at end of file
# -*- coding: utf-8 -*-
# Copyright (C) 2013 Centre de données Astrophysiques de Marseille
# Licensed under the CeCILL-v2 licence - see Licence_CeCILL_V2-en.txt
# Author: Yannick Roehlly <yannick.roehlly@oamp.fr>, Laure Ciesla <ciesla@physics.uoc.gr>
from . import common
from collections import OrderedDict
import numpy as np
from pcigale.data import Database
class Module(common.SEDCreationModule):
"""
Module computing the infra-red re-emission corresponding to an amount of
attenuation using the Dale et al (2014) models.
Given an amount of attenuation (e.g. resulting from the action of a dust
attenuation module) this module normalises the Dale et al (2014)
template corresponding to a given α to this amount of energy and add it
to the SED.
Information added to the SED: NAME_fracAGN, NAME_alpha.
"""
parameter_list = OrderedDict([
('fracAGN', (
'float',
"Contribution of the AGN",
None
)),
('alpha', (
'float',
"Alpha slope",
None
)),
('attenuation_value_keys', (
'string',
"Keys of the SED information dictionary where the module will "
"look for the attenuation (in W) to re-emit. You can give several "
"keys separated with a & (don't use commas), a re-emission "
"contribution will be added for each key.",
"attenuation"
))
])
out_parameter_list = OrderedDict([
('fracAGN', 'Contribution of the AGN'),
('alpha', 'Alpha slope'),
('lir', 'Total IR luminosity between 8 and 1000 microns (AGN + SB)')
])
def _init_code(self):
"""
Get the models out of the database
model_sb corresponds to an AGN fraction of 0%: only the dust heated by star formation
model_quasar corresponds to an AGN fraction of 100%: only the SED of the quasar
The energy attenuated is re-injected in model_sb only.
"""
alpha = self.parameters["alpha"]
database = Database()
self.model_sb = database.get_dale2014(0.00, alpha)
self.model_quasar = database.get_dale2014(1.00, alpha)
database.session.close_all()
def process(self, sed):
"""Add the IR re-emission contributions
Parameters
----------
sed : pcigale.sed.SED object
parameters : dictionary containing the parameters
"""
attenuation_value_keys = [
item.strip() for item in
self.parameters["attenuation_value_keys"].split("&")]
fracAGN = self.parameters["fracAGN"]
# Base name for adding information to the SED.
name = self.name or 'dale2014'
sed.add_module(name, self.parameters)
sed.add_info("fracAGN" + self.postfix, self.parameters["fracAGN"])
sed.add_info("alpha" + self.postfix, self.parameters["alpha"])
for attenuation in attenuation_value_keys:
sed.add_contribution(
name + '_sb_' + attenuation,
self.model_sb.wave,
sed.info[attenuation] * self.model_sb.lumin
)
sed.add_contribution(
name + '_quasar',
self.model_quasar.wave,
fracAGN * sed.info[attenuation] * self.model_quasar.lumin
)
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment