P92¶
- class dust_extinction.shapes.P92(BKG_amp=218.57142857142858, BKG_lambda=0.047, BKG_b=90.0, BKG_n=2.0, FUV_amp=18.545454545454547, FUV_lambda=0.07, FUV_b=4.0, FUV_n=6.5, NUV_amp=0.05961038961038961, NUV_lambda=0.22, NUV_b=-1.95, NUV_n=2.0, SIL1_amp=0.0026493506493506496, SIL1_lambda=9.7, SIL1_b=-1.95, SIL1_n=2.0, SIL2_amp=0.0026493506493506496, SIL2_lambda=18.0, SIL2_b=-1.8, SIL2_n=2.0, FIR_amp=0.015896103896103898, FIR_lambda=25.0, FIR_b=0.0, FIR_n=2.0, **kwargs)[source]¶
Bases:
Fittable1DModel
Pei (1992) 24 parameter shape model
- Parameters:
- BKG_ampfloat
background term amplitude
- BKG_lambdafloat
background term central wavelength
- BKG_bfloat
background term b coefficient
- BKG_nfloat
background term n coefficient [FIXED at n = 2]
- FUV_ampfloat
far-ultraviolet term amplitude
- FUV_lambdafloat
far-ultraviolet term central wavelength
- FUV_bfloat
far-ultraviolet term b coefficent
- FUV_nfloat
far-ultraviolet term n coefficient
- NUV_ampfloat
near-ultraviolet (2175 A) term amplitude
- NUV_lambdafloat
near-ultraviolet (2175 A) term central wavelength
- NUV_bfloat
near-ultraviolet (2175 A) term b coefficent
- NUV_nfloat
near-ultraviolet (2175 A) term n coefficient [FIXED at n = 2]
- SIL1_ampfloat
1st silicate feature (~10 micron) term amplitude
- SIL1_lambdafloat
1st silicate feature (~10 micron) term central wavelength
- SIL1_bfloat
1st silicate feature (~10 micron) term b coefficent
- SIL1_nfloat
1st silicate feature (~10 micron) term n coefficient [FIXED at n = 2]
- SIL2_ampfloat
2nd silicate feature (~18 micron) term amplitude
- SIL2_lambdafloat
2nd silicate feature (~18 micron) term central wavelength
- SIL2_bfloat
2nd silicate feature (~18 micron) term b coefficient
- SIL2_nfloat
2nd silicate feature (~18 micron) term n coefficient [FIXED at n = 2]
- FIR_ampfloat
far-infrared term amplitude
- FIR_lambdafloat
far-infrared term central wavelength
- FIR_bfloat
far-infrared term b coefficent
- FIR_nfloat
far-infrared term n coefficient [FIXED at n = 2]
Notes
From Pei (1992, ApJ, 395, 130)
Applicable from the extreme UV to far-IR
Example showing a P92 curve with components identified.
import numpy as np import matplotlib.pyplot as plt import astropy.units as u from dust_extinction.shapes import P92 fig, ax = plt.subplots() # generate the curves and plot them lam = np.logspace(-3.0, 3.0, num=1000) x = (1.0/lam)/u.micron ext_model = P92() ax.plot(1/x,ext_model(x),label='total') ext_model = P92(FUV_amp=0., NUV_amp=0.0, SIL1_amp=0.0, SIL2_amp=0.0, FIR_amp=0.0) ax.plot(1./x,ext_model(x),label='BKG only') ext_model = P92(NUV_amp=0.0, SIL1_amp=0.0, SIL2_amp=0.0, FIR_amp=0.0) ax.plot(1./x,ext_model(x),label='BKG+FUV only') ext_model = P92(FUV_amp=0., SIL1_amp=0.0, SIL2_amp=0.0, FIR_amp=0.0) ax.plot(1./x,ext_model(x),label='BKG+NUV only') ext_model = P92(FUV_amp=0., NUV_amp=0.0, SIL2_amp=0.0) ax.plot(1./x,ext_model(x),label='BKG+FIR+SIL1 only') ext_model = P92(FUV_amp=0., NUV_amp=0.0, SIL1_amp=0.0) ax.plot(1./x,ext_model(x),label='BKG+FIR+SIL2 only') ext_model = P92(FUV_amp=0., NUV_amp=0.0, SIL1_amp=0.0, SIL2_amp=0.0) ax.plot(1./x,ext_model(x),label='BKG+FIR only') # Milky Way observed extinction as tabulated by Pei (1992) MW_x = [0.21, 0.29, 0.45, 0.61, 0.80, 1.11, 1.43, 1.82, 2.27, 2.50, 2.91, 3.65, 4.00, 4.17, 4.35, 4.57, 4.76, 5.00, 5.26, 5.56, 5.88, 6.25, 6.71, 7.18, 7.60, 8.00, 8.50, 9.00, 9.50, 10.00] MW_x = np.array(MW_x) MW_exvebv = [-3.02, -2.91, -2.76, -2.58, -2.23, -1.60, -0.78, 0.00, 1.00, 1.30, 1.80, 3.10, 4.19, 4.90, 5.77, 6.57, 6.23, 5.52, 4.90, 4.65, 4.60, 4.73, 4.99, 5.36, 5.91, 6.55, 7.45, 8.45, 9.80, 11.30] MW_exvebv = np.array(MW_exvebv) Rv = 3.08 MW_axav = MW_exvebv/Rv + 1.0 ax.plot(1./MW_x, MW_axav, 'o', label='MW Observed') ax.set_xscale('log') ax.set_yscale('log') ax.set_ylim(1e-3,10.) ax.set_xlabel(r'$\lambda$ [$\mu$m]') ax.set_ylabel(r'$A(x)/A(V)$') ax.legend(loc='best') plt.show()
(
Source code
,png
,hires.png
,pdf
)Attributes Summary
Function (similar to the model's
evaluate
) to compute the derivatives of the model with respect to its parameters, for use by fitting algorithms.The number of inputs.
The number of outputs.
Names of the parameters that describe models of this type.
Methods Summary
__call__
(*inputs[, model_set_axis, ...])Evaluate this model using the given input(s) and the parameter values that were specified when the model was instantiated.
evaluate
(in_x, BKG_amp, BKG_lambda, BKG_b, ...)P92 function
Attributes Documentation
- AbAv = 1.3246753246753247¶
- BKG_amp = Parameter('BKG_amp', value=218.57142857142858, bounds=(0.0, None))¶
- BKG_b = Parameter('BKG_b', value=90.0)¶
- BKG_lambda = Parameter('BKG_lambda', value=0.047)¶
- BKG_n = Parameter('BKG_n', value=2.0, fixed=True)¶
- FIR_amp = Parameter('FIR_amp', value=0.015896103896103898, bounds=(0.0, None))¶
- FIR_b = Parameter('FIR_b', value=0.0)¶
- FIR_lambda = Parameter('FIR_lambda', value=25.0, bounds=(20.0, 30.0))¶
- FIR_n = Parameter('FIR_n', value=2.0, fixed=True)¶
- FUV_amp = Parameter('FUV_amp', value=18.545454545454547, bounds=(0.0, None))¶
- FUV_b = Parameter('FUV_b', value=4.0)¶
- FUV_lambda = Parameter('FUV_lambda', value=0.07, bounds=(0.06, 0.08))¶
- FUV_n = Parameter('FUV_n', value=6.5)¶
- NUV_amp = Parameter('NUV_amp', value=0.05961038961038961, bounds=(0.0, None))¶
- NUV_b = Parameter('NUV_b', value=-1.95)¶
- NUV_lambda = Parameter('NUV_lambda', value=0.22, bounds=(0.2, 0.24))¶
- NUV_n = Parameter('NUV_n', value=2.0, fixed=True)¶
- SIL1_amp = Parameter('SIL1_amp', value=0.0026493506493506496, bounds=(0.0, None))¶
- SIL1_b = Parameter('SIL1_b', value=-1.95)¶
- SIL1_lambda = Parameter('SIL1_lambda', value=9.7, bounds=(7.0, 13.0))¶
- SIL1_n = Parameter('SIL1_n', value=2.0, fixed=True)¶
- SIL2_amp = Parameter('SIL2_amp', value=0.0026493506493506496, bounds=(0.0, None))¶
- SIL2_b = Parameter('SIL2_b', value=-1.8)¶
- SIL2_lambda = Parameter('SIL2_lambda', value=18.0, bounds=(15.0, 21.0))¶
- SIL2_n = Parameter('SIL2_n', value=2.0, fixed=True)¶
- fit_deriv = None¶
Function (similar to the model’s
evaluate
) to compute the derivatives of the model with respect to its parameters, for use by fitting algorithms. In other words, this computes the Jacobian matrix with respect to the model’s parameters.
- n_inputs = 1¶
The number of inputs.
- n_outputs = 1¶
The number of outputs.
- param_names = ('BKG_amp', 'BKG_lambda', 'BKG_b', 'BKG_n', 'FUV_amp', 'FUV_lambda', 'FUV_b', 'FUV_n', 'NUV_amp', 'NUV_lambda', 'NUV_b', 'NUV_n', 'SIL1_amp', 'SIL1_lambda', 'SIL1_b', 'SIL1_n', 'SIL2_amp', 'SIL2_lambda', 'SIL2_b', 'SIL2_n', 'FIR_amp', 'FIR_lambda', 'FIR_b', 'FIR_n')¶
Names of the parameters that describe models of this type.
The parameters in this tuple are in the same order they should be passed in when initializing a model of a specific type. Some types of models, such as polynomial models, have a different number of parameters depending on some other property of the model, such as the degree.
When defining a custom model class the value of this attribute is automatically set by the
Parameter
attributes defined in the class body.
- x_range = [0.001, 1000.0]¶
Methods Documentation
- __call__(*inputs, model_set_axis=None, with_bounding_box=False, fill_value=nan, equivalencies=None, inputs_map=None, **new_inputs)¶
Evaluate this model using the given input(s) and the parameter values that were specified when the model was instantiated.
- evaluate(in_x, BKG_amp, BKG_lambda, BKG_b, BKG_n, FUV_amp, FUV_lambda, FUV_b, FUV_n, NUV_amp, NUV_lambda, NUV_b, NUV_n, SIL1_amp, SIL1_lambda, SIL1_b, SIL1_n, SIL2_amp, SIL2_lambda, SIL2_b, SIL2_n, FIR_amp, FIR_lambda, FIR_b, FIR_n)[source]¶
P92 function
- Parameters:
- in_x: float
expects either x in units of wavelengths or frequency or assumes wavelengths in wavenumbers [1/micron]
internally wavenumbers are used
- Returns:
- axav: np array (float)
A(x)/A(V) extinction curve [mag]
- Raises:
- ValueError
Input x values outside of defined range