M14¶
- class dust_extinction.parameter_averages.M14(*args, meta=None, name=None, **kwargs)[source]¶
Bases:
BaseExtRvModel
Maiz Apellaniz et al (2014) Milky Way & LMC R(V) dependent model
- Parameters:
- R5495: float
R5495 = A(5485)/E(4405-5495) Spectral equivalent to photometric R(V), standard value is 3.1
- Raises:
- InputParameterError
Input Rv values outside of defined range
Notes
M14 R5485-dependent model
From Maiz Apellaniz et al. (2014, A&A, 564, 63), following structure of IDL code provided in paper appendix
The published UV extinction curve is identical to Clayton, Cardelli, and Mathis (1989, CCM). Forcing the optical section to match smoothly with CCM introduces a non-physical feature at high values of R5495 around 3.9 inverse microns; see section 5 in Maiz Apellaniz et al. (2014) for more discussion. For that reason, we provide the M14 model only through 3.3 inverse microns, the limit of the optical in CCM.
R5495 = A(5485)/E(4405-5495) Spectral equivalent to photometric R(V), standard value is 3.1
Example showing M14 curves for a range of R5495 values.
import numpy as np import matplotlib.pyplot as plt import astropy.units as u from dust_extinction.parameter_averages import M14 fig, ax = plt.subplots() # temp model to get the correct x range text_model = M14() # generate the curves and plot them x = np.arange(text_model.x_range[0], text_model.x_range[1],0.1)/u.micron Rvs = ['2.0','3.1','4.0','5.0','6.0'] for cur_Rv in Rvs: ext_model = M14(Rv=cur_Rv) ax.plot(x,ext_model(x),label='R(V) = ' + str(cur_Rv)) # for 2nd x-axis with lambda values axis_xs = np.array([0.3, 0.5, 1.0, 2.0]) new_ticks = 1 / axis_xs new_ticks_labels = ["%.2f" % z for z in axis_xs] tax = ax.twiny() tax.set_xlim(ax.get_xlim()) tax.set_xticks(new_ticks) tax.set_xticklabels(new_ticks_labels) tax.set_xlabel(r"$\lambda$ [$\mu$m]") ax.set_xlabel(r'$x$ [$\mu m^{-1}$]') ax.set_ylabel(r'$A(x)/A(V)$') ax.legend(loc='best') plt.show()
(
Source code
,png
,hires.png
,pdf
)Attributes Summary
Methods Summary
evaluate
(in_x, Rv)M14 function
Attributes Documentation
- Rv_range = [2.0, 6.0]¶
- x_range = [0.3, 3.3]¶
Methods Documentation
- evaluate(in_x, Rv)[source]¶
M14 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