Commit c2d19946 authored by Médéric Boquien's avatar Médéric Boquien

Allow the AGN fraction to be defined based on the luminosity between two wavelengths.

parent 1f0a652f
Pipeline #3067 skipped with stage
......@@ -133,10 +133,11 @@ class SKIRTOR2016(SedModule):
0.1
)),
('lambda_fracAGN', (
'float(min=0, max=500)',
'Wavelength in micron where to compute the AGN fraction. If zero '
'then the total dust luminosity is adopted.',
0
'string()',
'Wavelength range in microns where to compute the AGN fraction. '
'Note that it includes all the components, not just dust emission. '
'To use the the total dust luminosity set to 0/0.',
"0/0"
)),
('law', (
'cigale_list(dtype=int, options=0 & 1 & 2)',
......@@ -171,7 +172,16 @@ class SKIRTOR2016(SedModule):
self.Mcl = float(self.parameters["Mcl"])
self.i = int(self.parameters["i"])
self.fracAGN = float(self.parameters["fracAGN"])
self.lambda_fracAGN = float(self.parameters["lambda_fracAGN"]) * 1e3
if self.fracAGN == 1.:
raise ValueError("AGN fraction is exactly 1. Behaviour undefined.")
lambda_fracAGN = str(self.parameters["lambda_fracAGN"]).split('/')
self.lambdamin_fracAGN = float(lambda_fracAGN[0]) * 1e3
self.lamdamax_fracAGN = float(lambda_fracAGN[1]) * 1e3
if (self.lambdamin_fracAGN < 0 or
self.lambdamin_fracAGN > self.lamdamax_fracAGN ):
raise ValueError("lambda_fracAGN incorrect. Constrain "
f"0 < {self.lambdamin_fracAGN} < "
f"{self.lamdamax_fracAGN} not respected.")
self.law = int(self.parameters["law"])
self.EBV = float(self.parameters["EBV"])
self.temperature = float(self.parameters["temperature"])
......@@ -236,11 +246,23 @@ class SKIRTOR2016(SedModule):
# Integrate AGN luminosity for different components
self.lumin_disk = np.trapz(self.SKIRTOR2016.disk, x=self.SKIRTOR2016.wave)
if self.lambda_fracAGN > 0.:
self.AGNfluxatlambda = np.interp(self.lambda_fracAGN,
self.SKIRTOR2016.wave,
self.SKIRTOR2016.dust +
self.SKIRTOR2016.disk)
if self.lambdamin_fracAGN < self.lamdamax_fracAGN:
w = np.where((self.SKIRTOR2016.wave >= self.lambdamin_fracAGN) &
(self.SKIRTOR2016.wave <= self.lamdamax_fracAGN))
wl = np.hstack([self.lambdamin_fracAGN, self.SKIRTOR2016.wave[w],
self.lamdamax_fracAGN])
spec = np.interp(wl, self.SKIRTOR2016.wave,
self.SKIRTOR2016.dust + self.SKIRTOR2016.disk)
self.AGNlumin = np.trapz(spec, x=wl)
elif (self.lambdamin_fracAGN == 0.) & (self.lamdamax_fracAGN == 0.):
self.AGNlumin = 1.
elif self.lambdamin_fracAGN == self.lamdamax_fracAGN:
self.AGNlumin = np.interp(self.lambdamin_fracAGN,
self.SKIRTOR2016.wave,
self.SKIRTOR2016.dust +
self.SKIRTOR2016.disk)
# Store the SED wavelengths
self.wl = None
def process(self, sed):
"""Add the IR re-emission contributions
......@@ -271,25 +293,29 @@ class SKIRTOR2016(SedModule):
sed.add_info('agn.emissivity', self.emissivity)
# Compute the AGN luminosity
if self.fracAGN < 1.:
if self.lambda_fracAGN > 0.:
scale = np.interp(self.lambda_fracAGN, sed.wavelength_grid,
sed.luminosity) / self.AGNfluxatlambda
else:
scale = luminosity
agn_power = scale * (1. / (1. - self.fracAGN) - 1.)
lumin_dust = agn_power
lumin_disk = agn_power * np.trapz(self.SKIRTOR2016.disk,
x=self.SKIRTOR2016.wave)
lumin = lumin_dust + lumin_disk
else:
raise Exception("AGN fraction is exactly 1. Behaviour "
"undefined.")
if self.lambdamin_fracAGN < self.lamdamax_fracAGN:
if self.wl is None:
w = np.where((sed.wavelength_grid >= self.lambdamin_fracAGN) &
(sed.wavelength_grid <= self.lamdamax_fracAGN))
self.wl = np.hstack([self.lambdamin_fracAGN,
sed.wavelength_grid[w],
self.lamdamax_fracAGN])
spec = np.interp(self.wl, sed.wavelength_grid, sed.luminosity)
scale = np.trapz(spec, x=self.wl) / self.AGNlumin
elif (self.lambdamin_fracAGN == 0.) and (self.lamdamax_fracAGN == 0.):
scale = luminosity
elif self.lambdamin_fracAGN == self.lamdamax_fracAGN:
scale = np.interp(self.lambdamin_fracAGN, sed.wavelength_grid,
sed.luminosity) / self.AGNlumin
agn_power = scale * (1. / (1. - self.fracAGN) - 1.)
lumin_dust = agn_power
lumin_disk = agn_power * np.trapz(self.SKIRTOR2016.disk,
x=self.SKIRTOR2016.wave)
sed.add_info('agn.dust_luminosity', lumin_dust, True, unit='W')
sed.add_info('agn.disk_luminosity', lumin_disk, True, unit='W')
sed.add_info('agn.luminosity', lumin, True, unit='W')
sed.add_info('agn.luminosity', lumin_dust + lumin_disk, True, unit='W')
sed.add_contribution('agn.SKIRTOR2016_dust', self.SKIRTOR2016.wave,
agn_power * self.SKIRTOR2016.dust)
......
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