From 49c3a298e8e252103403abbd9149cc75252d2cf0 Mon Sep 17 00:00:00 2001 From: Karl Gordon Date: Fri, 12 Jan 2018 17:22:07 -0500 Subject: [PATCH 1/6] Adding in the AlvAvtoElv helper function and PEP8 changes --- docs/dust_extinction/choose_model.rst | 33 +- docs/dust_extinction/fit_extinction.rst | 20 +- docs/dust_extinction/model_flavors.rst | 47 ++- dust_extinction/dust_extinction.py | 337 ++++++++++-------- .../{test_g09_ave.py => test_gcc09_ave.py} | 61 ++-- 5 files changed, 275 insertions(+), 223 deletions(-) rename dust_extinction/tests/{test_g09_ave.py => test_gcc09_ave.py} (76%) diff --git a/docs/dust_extinction/choose_model.rst b/docs/dust_extinction/choose_model.rst index 0b708fc..0d9155c 100644 --- a/docs/dust_extinction/choose_model.rst +++ b/docs/dust_extinction/choose_model.rst @@ -2,9 +2,9 @@ How to Choose a Model ##################### -The ``dust_extinction`` package provides a suite of dust extinction models. +The ``dust_extinction`` package provides a suite of dust extinction models. Which model to use can depend on the wavelength range of interest, the expected -type of extinction, or some other property. +type of extinction, or some other property. Average Models ============== @@ -14,7 +14,7 @@ Simple Average Curves These are straightforward averages of observed extinction curves. They are the simpliest models and include models for the MW -(:class:`~dust_extinction.dust_extinction.G09_MWAvg`), the LMC +(:class:`~dust_extinction.dust_extinction.GCC09_MWAvg`), the LMC (:class:`~dust_extinction.dust_extinction.G03_LMCAvg`, :class:`~dust_extinction.dust_extinction.G03_LMC2`) and the SMC (:class:`~dust_extinction.dust_extinction.G03_SMCBar`). @@ -23,17 +23,17 @@ One often used alternative to these straight average models is to use one of the parameter dependent models with the average R(V) value. For the Milky Way, the usual average used is R(V) = 3.1. -+------------+-------------+------------------+--------------+ -| Model | x range | wavelength range | galaxy | -+============+=============+==================+==============+ -| G09_MWAvg | 0.3 - 10.96 | 0.0912 - 3.3 | MW | -+------------+-------------+------------------+--------------+ -| G03_LMCAvg | 0.3 - 10.0 | 0.1 - 3.3 | LMC | -+------------+-------------+------------------+--------------+ -| G03_LMC2 | 0.3 - 10.0 | 0.1 - 3.3 | LMC (30 Dor) | -+------------+-------------+------------------+--------------+ -| G03_SMCBar | 0.3 - 10.0 | 0.1 - 3.3 | SMC | -+------------+-------------+------------------+--------------+ ++--------------+-------------+------------------+--------------+ +| Model | x range | wavelength range | galaxy | ++==============+=============+==================+==============+ +| GCC09_MWAvg | 0.3 - 10.96 | 0.0912 - 3.3 | MW | ++--------------+-------------+------------------+--------------+ +| G03_LMCAvg | 0.3 - 10.0 | 0.1 - 3.3 | LMC | ++--------------+-------------+------------------+--------------+ +| G03_LMC2 | 0.3 - 10.0 | 0.1 - 3.3 | LMC (30 Dor) | ++--------------+-------------+------------------+--------------+ +| G03_SMCBar | 0.3 - 10.0 | 0.1 - 3.3 | SMC | ++--------------+-------------+------------------+--------------+ Parameter Dependent Average Curves @@ -55,14 +55,14 @@ is based on signifincatly more extinction curves than the :class:`~dust_extinction.dust_extinction.O94` models. +---------+-------------+-------------+------------------+--------------+ -| Model | Parameters | x range | wavelength range | galaxy | +| Model | Parameters | x range | wavelength range | galaxy | | | | [1/micron] | [micron] | | +=========+=============+=============+==================+==============+ | CCM89 | R(V) | 0.3 - 10.0 | 0.1 - 3.3 | MW | +---------+-------------+-------------+------------------+--------------+ | O94 | R(V) | 0.3 - 10.0 | 0.1 - 3.3 | MW | +---------+-------------+-------------+------------------+--------------+ -| F99 | R(V) | 0.3 - 10.0 | 0.1 - 3.3 | MW | +| F99 | R(V) | 0.3 - 10.0 | 0.1 - 3.3 | MW | +---------+-------------+-------------+------------------+--------------+ | G16 | R(V)_A, f_A | 0.3 - 10.0 | 0.1 - 3.3 | MW, LMC, SMC | +---------+-------------+-------------+------------------+--------------+ @@ -85,4 +85,3 @@ but only covers the UV wavelength range. +------------+--------------+------------------+-------------------+ | P92 | 0.001 - 1000 | 0.001 - 1000 | 19 (24 possible) | +------------+--------------+------------------+-------------------+ - diff --git a/docs/dust_extinction/fit_extinction.rst b/docs/dust_extinction/fit_extinction.rst index 11cbbea..9464b13 100644 --- a/docs/dust_extinction/fit_extinction.rst +++ b/docs/dust_extinction/fit_extinction.rst @@ -6,14 +6,14 @@ The ``dust_extinction`` package is built on the `astropy.modeling `_ package. Fitting is done in the standard way for this package where the model is initialized with a starting point (either the default or user input), the fitter -is choosen, and the fit performed. +is choosen, and the fit performed. Example: FM90 Fit ================= In this example, the FM90 model is used to fit the observed average -extinction curve for the LMC outside of the LMC2 supershell region -(G03_LMCAvg ``dust_extinction`` model). +extinction curve for the LMC outside of the LMC2 supershell region +(G03_LMCAvg ``dust_extinction`` model). .. plot:: :include-source: @@ -50,7 +50,7 @@ extinction curve for the LMC outside of the LMC2 supershell region ax.plot(x, y, 'ko', label='Observed Curve') ax.plot(x[gindxs], fm90_init(x[gindxs]), label='Initial guess') ax.plot(x[gindxs], g03_fit(x[gindxs]), label='Fitted model') - + ax.set_xlabel('$x$ [$\mu m^{-1}$]') ax.set_ylabel('$E(x-V)/E(B-V)$') @@ -64,7 +64,7 @@ Example: P92 Fit ================ In this example, the P92 model is used to fit the observed average -extinction curve for the MW (G09_MWAvg ``dust_extinction`` model). +extinction curve for the MW (GCC09_MWAvg ``dust_extinction`` model). The fit is done using the observed uncertainties that are passed as weights. The weights assume the noise is Gaussian and not correlated between data points. @@ -80,10 +80,10 @@ between data points. from astropy.modeling.fitting import LevMarLSQFitter - from dust_extinction.dust_extinction import (P92, G09_MWAvg) + from dust_extinction.dust_extinction import (P92, GCC09_MWAvg) # get an observed extinction curve to fit - g09_model = G09_MWAvg() + g09_model = GCC09_MWAvg() # get an observed extinction curve to fit x = g09_model.obsdata_x @@ -106,14 +106,14 @@ between data points. p92_init.FIR_amp.fixed = True p92_init.FIR_lambda.fixed = True p92_init.FIR_b.fixed = True - + # pick the fitter fit = LevMarLSQFitter() # set to avoid the "fit may have been unsuccessful" warning # fit is fine, but this means the build of the docs fails warnings.simplefilter('ignore', category=AstropyWarning) - + # fit the data to the P92 model using the fitter # use the initialized model as the starting point # accuracy set to avoid warning the fit may have failed @@ -129,7 +129,7 @@ between data points. ax.set_xlabel('$x$ [$\mu m^{-1}$]') ax.set_ylabel('$A(x)/A(V)$') - ax.set_title('Example P92 Fit to MW average curve') + ax.set_title('Example P92 Fit to GCC09_MWAvg average curve') ax.legend(loc='best') plt.tight_layout() diff --git a/docs/dust_extinction/model_flavors.rst b/docs/dust_extinction/model_flavors.rst index cfa1881..97bc401 100644 --- a/docs/dust_extinction/model_flavors.rst +++ b/docs/dust_extinction/model_flavors.rst @@ -2,7 +2,7 @@ Model Flavors ############# -There are three differnet types of models: average, R(V)+ dependent prediction, +There are three different types of models: average, R(V)+ dependent prediction, and shape fitting. Average models @@ -21,41 +21,41 @@ Average models import numpy as np import matplotlib.pyplot as plt import astropy.units as u - - from dust_extinction.dust_extinction import (G09_MWAvg, + + from dust_extinction.dust_extinction import (GCC09_MWAvg, G03_SMCBar, G03_LMCAvg, - G03_LMC2) - + G03_LMC2) + fig, ax = plt.subplots() - + # generate the curves and plot them x = np.arange(0.3,10.0,0.1)/u.micron ext_model = G09_MWAvg() - ax.plot(x,ext_model(x),label='G09 MWAvg') - + ax.plot(x,ext_model(x),label='GCC09 MWAvg') + ext_model = G03_SMCBar() ax.plot(x,ext_model(x),label='G03 SMCBar') ext_model = G03_LMCAvg() ax.plot(x,ext_model(x),label='G03 LMCAvg') - + ext_model = G03_LMC2() ax.plot(x,ext_model(x),label='G03 LMC2') - + ax.set_xlabel('$x$ [$\mu m^{-1}$]') ax.set_ylabel('$A(x)/A(V)$') - + ax.legend(loc='best') plt.tight_layout() plt.show() - + R(V) (+ other variables) dependent prediction models ==================================================== These models provide predictions of the shape of the dust extinction - given input variable(s). + given input variable(s). The R(V) dependent models include CCM89 the original such model (Cardelli, Clayton, and Mathis 1989), the O94 model that updates the @@ -63,9 +63,9 @@ R(V) (+ other variables) dependent prediction models (Fitzpatrick 1999). These models are based on the average behavior of extinction in the Milky Way as a function of R(V). - In addition, the (R(V), f_A) two parameter relationship from + In addition, the (R(V), f_A) two parameter relationship from Gordon et al. (2016) is included. This model is based on the average - behavior of extinction in the Milky Way, Large Magellanic Cloud, and + behavior of extinction in the Milky Way, Large Magellanic Cloud, and Small Magellanic Cloud. .. plot:: @@ -121,7 +121,7 @@ R(V) (+ other variables) dependent prediction models ax.legend(loc='best') plt.tight_layout() plt.show() - + .. plot:: import numpy as np @@ -219,7 +219,7 @@ Shape fitting models These models are used to fit the detailed shape of dust extinction curves. The FM90 (Fitzpatrick & Mass 1990) model uses 6 parameters to fit the shape of the ultraviolet extinction. - The P92 (Pei 1992) uses 19 parameters to fit the shape of the X-ray to + The P92 (Pei 1992) uses 19 parameters to fit the shape of the X-ray to far-infrared extinction. .. plot:: @@ -273,11 +273,11 @@ Shape fitting models ext_model = P92() ax.plot(1/x,ext_model(x),label='total') - ext_model = P92(FUV_amp=0., NUV_amp=0.0, + 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, + 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') @@ -285,15 +285,15 @@ Shape fitting models 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, + 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, + 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, + 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') @@ -310,6 +310,3 @@ Shape fitting models ax.legend(loc='best') plt.tight_layout() plt.show() - - - diff --git a/dust_extinction/dust_extinction.py b/dust_extinction/dust_extinction.py index 37ae20c..f435135 100644 --- a/dust_extinction/dust_extinction.py +++ b/dust_extinction/dust_extinction.py @@ -1,9 +1,6 @@ from __future__ import (absolute_import, print_function, division) -# STDLIB -import warnings - import numpy as np from scipy import interpolate @@ -13,17 +10,19 @@ __all__ = ['CCM89', 'FM90', 'P92', 'O94', 'F99', 'G03_SMCBar', 'G03_LMCAvg', 'G03_LMC2', - 'G09_MWAvg', 'G16'] + 'GCC09_MWAvg', 'G16', + 'AlAvToElv'] -x_range_CCM89 = [0.3,10.0] -x_range_FM90 = [1.0/0.32,1.0/0.0912] -x_range_P92 = [1.0/1e3,1.0/1e-3] +x_range_CCM89 = [0.3, 10.0] +x_range_FM90 = [1.0/0.32, 1.0/0.0912] +x_range_P92 = [1.0/1e3, 1.0/1e-3] x_range_O94 = x_range_CCM89 -x_range_F99 = [0.3,10.0] -x_range_G03 = [0.3,10.0] -x_range_G09 = [0.3,1.0/0.0912] +x_range_F99 = [0.3, 10.0] +x_range_G03 = [0.3, 10.0] +x_range_GCC09 = [0.3, 1.0/0.0912] x_range_G16 = x_range_G03 + def _test_valid_x_range(x, x_range, outname): """ Test if any of the x values are outside of the valid range @@ -41,13 +40,14 @@ def _test_valid_x_range(x, x_range, outname): """ if np.logical_or(np.any(x < x_range[0]), np.any(x > x_range[1])): - raise ValueError('Input x outside of range defined for ' + outname \ + raise ValueError('Input x outside of range defined for ' + outname + ' [' + str(x_range[0]) - + ' <= x <= ' + + ' <= x <= ' + str(x_range[1]) + ', x has units 1/micron]') + def _curve_F99_method(in_x, Rv, C1, C2, C3, C4, xo, gamma, optnir_axav_x, optnir_axav_y, @@ -119,7 +119,7 @@ def _curve_F99_method(in_x, Rv, x_cutval_uv = 10000.0/2700.0 # required UV points for spline interpolation - x_splineval_uv = 10000.0/np.array([2700.0,2600.0]) + x_splineval_uv = 10000.0/np.array([2700.0, 2600.0]) # UV points in input x indxs_uv, = np.where(x >= x_cutval_uv) @@ -153,7 +153,7 @@ def _curve_F99_method(in_x, Rv, x_splineval_optir = optnir_axav_x # determine optical/IR values at spline points - y_splineval_optir = optnir_axav_y + y_splineval_optir = optnir_axav_y # add in zero extinction at infinite wavelength x_splineval_optir = np.insert(x_splineval_optir, 0, 0.0) @@ -168,6 +168,7 @@ def _curve_F99_method(in_x, Rv, # return A(x)/A(V) return axav + class BaseExtModel(Model): """ Base Extinction Model. Do not use. @@ -212,7 +213,8 @@ def extinguish(self, x, Av=None, Ebv=None): Av = self.Rv*Ebv # return fractional extinction - return np.power(10.0,-0.4*axav*Av) + return np.power(10.0, -0.4*axav*Av) + class BaseExtAve(Model): """ @@ -258,14 +260,14 @@ def extinguish(self, x, Av=None, Ebv=None): Av = self.Rv*Ebv # return fractional extinction - return np.power(10.0,-0.4*axav*Av) + return np.power(10.0, -0.4*axav*Av) + class BaseExtRvModel(BaseExtModel): """ Base Extinction R(V)-dependent Model. Do not use. """ - - Rv = Parameter(description="R(V) = A(V)/E(B-V) = " \ + Rv = Parameter(description="R(V) = A(V)/E(B-V) = " + "total-to-selective extinction", default=3.1) @@ -290,14 +292,15 @@ def Rv(self, value): + " and " + str(self.Rv_range[1])) + class BaseExtRvAfAModel(BaseExtModel): """ Base Extinction R(V)_A, f_A -dependent Model. Do not use. """ - RvA = Parameter(description="R_A(V) = A(V)/E(B-V) = " \ - + "total-to-selective extinction of component A", - default=3.1) + RvA = Parameter(description="R_A(V) = A(V)/E(B-V) = " + + "total-to-selective extinction of component A", + default=3.1) fA = Parameter(description="f_A = mixture coefficent of component A", default=1.0) @@ -343,6 +346,43 @@ def fA(self, value): + " and " + str(self.fA_range[1])) + +class AlAvToElv(Fittable1DModel): + """ + Model to convert from A(lambda)/A(V) to E(lambda-V) + + Paramters + --------- + Av : float + dust column in A(V) [mag] + """ + inputs = ('alav',) + outputs = ('elv',) + + Av = Parameter(description="A(V)", + default=1.0, min=0.0) + + @staticmethod + def evaluate(alav, Av): + """ + AlAvToElv function + + Paramters + --------- + alav : np array (float) + E(lambda-V)/E(B-V) values + + Returns + ------- + elv : np array (float) + E(lamda - V) + """ + return (alav - 1.0)*Av + + # use numerical derivaties (need to add analytic) + fit_deriv = None + + class CCM89(BaseExtRvModel): """ CCM89 extinction model calculation @@ -390,8 +430,7 @@ class CCM89(BaseExtRvModel): ax.legend(loc='best') plt.show() """ - - Rv_range = [2.0,6.0] + Rv_range = [2.0, 6.0] x_range = x_range_CCM89 @staticmethod @@ -435,11 +474,11 @@ def evaluate(in_x, Rv): b = np.zeros(n_x) # define the ranges - ir_indxs = np.where(np.logical_and(0.3 <= x,x < 1.1)) - opt_indxs = np.where(np.logical_and(1.1 <= x,x < 3.3)) - nuv_indxs = np.where(np.logical_and(3.3 <= x,x <= 8.0)) - fnuv_indxs = np.where(np.logical_and(5.9 <= x,x <= 8)) - fuv_indxs = np.where(np.logical_and(8 < x,x <= 10)) + ir_indxs = np.where(np.logical_and(0.3 <= x, x < 1.1)) + opt_indxs = np.where(np.logical_and(1.1 <= x, x < 3.3)) + nuv_indxs = np.where(np.logical_and(3.3 <= x, x <= 8.0)) + fnuv_indxs = np.where(np.logical_and(5.9 <= x, x <= 8)) + fuv_indxs = np.where(np.logical_and(8 < x, x <= 10)) # Infrared y = x[ir_indxs]**1.61 @@ -454,11 +493,11 @@ def evaluate(in_x, Rv): 1.07233, 2.28305, 1.41338, 0), y) # NUV - a[nuv_indxs] = 1.752-.316*x[nuv_indxs] \ - - 0.104/((x[nuv_indxs] - 4.67)**2 + .341) - b[nuv_indxs] = -3.09 + \ - 1.825*x[nuv_indxs] \ - + 1.206/((x[nuv_indxs] - 4.62)**2 + .263) + a[nuv_indxs] = (1.752-.316*x[nuv_indxs] + - 0.104/((x[nuv_indxs] - 4.67)**2 + .341)) + b[nuv_indxs] = (-3.09 + + 1.825*x[nuv_indxs] + + 1.206/((x[nuv_indxs] - 4.62)**2 + .263)) # far-NUV y = x[fnuv_indxs] - 5.9 @@ -473,6 +512,7 @@ def evaluate(in_x, Rv): # return A(x)/A(V) return a + b/Rv + class FM90(Fittable1DModel): """ FM90 extinction model calculation @@ -623,7 +663,7 @@ def fit_deriv(in_x, C1, C2, C3, C4, xo, gamma): denom = (x2mxo2_2 - x2*g2)**2 # derivatives - d_C1 = np.full((len(x)),1.) + d_C1 = np.full((len(x)), 1.) d_C2 = x d_C3 = (x2/(x2mxo2_2 + x2*g2)) @@ -640,7 +680,6 @@ def fit_deriv(in_x, C1, C2, C3, C4, xo, gamma): return [d_C1, d_C2, d_C3, d_C4, d_xo, d_gamma] - #fit_deriv = None class P92(Fittable1DModel): """ @@ -791,32 +830,32 @@ class P92(Fittable1DModel): BKG_lambda = Parameter(description="BKG term: center wavelength", default=0.047) BKG_b = Parameter(description="BKG term: b coefficient", - default=90.) + default=90.) BKG_n = Parameter(description="BKG term: n coefficient", - default=2.0, fixed=True) + default=2.0, fixed=True) FUV_amp = Parameter(description="FUV term: amplitude", default=14.*AbAv, min=0.0) FUV_lambda = Parameter(description="FUV term: center wavelength", - default=0.08, bounds=(0.07,0.09)) + default=0.08, bounds=(0.07, 0.09)) FUV_b = Parameter(description="FUV term: b coefficient", - default=4.0) + default=4.0) FUV_n = Parameter(description="FUV term: n coefficient", - default=6.5) + default=6.5) NUV_amp = Parameter(description="NUV term: amplitude", default=0.045*AbAv, min=0.0) NUV_lambda = Parameter(description="NUV term: center wavelength", - default=0.22, bounds=(0.20,0.24)) + default=0.22, bounds=(0.20, 0.24)) NUV_b = Parameter(description="NUV term: b coefficient", - default=-1.95) + default=-1.95) NUV_n = Parameter(description="NUV term: n coefficient", - default=2.0, fixed=True) + default=2.0, fixed=True) SIL1_amp = Parameter(description="SIL1 term: amplitude", default=0.002*AbAv, min=0.0) SIL1_lambda = Parameter(description="SIL1 term: center wavelength", - default=9.7, bounds=(7.0,13.0)) + default=9.7, bounds=(7.0, 13.0)) SIL1_b = Parameter(description="SIL1 term: b coefficient", default=-1.95) SIL1_n = Parameter(description="SIL1 term: n coefficient", @@ -825,20 +864,20 @@ class P92(Fittable1DModel): SIL2_amp = Parameter(description="SIL2 term: amplitude", default=0.002*AbAv, min=0.0) SIL2_lambda = Parameter(description="SIL2 term: center wavelength", - default=18.0, bounds=(15.0,21.0)) + default=18.0, bounds=(15.0, 21.0)) SIL2_b = Parameter(description="SIL2 term: b coefficient", - default=-1.80) + default=-1.80) SIL2_n = Parameter(description="SIL2 term: n coefficient", - default=2.0, fixed=True) + default=2.0, fixed=True) FIR_amp = Parameter(description="FIR term: amplitude", default=0.012*AbAv, min=0.0) FIR_lambda = Parameter(description="FIR term: center wavelength", - default=25.0, bounds=(20.0,30.0)) + default=25.0, bounds=(20.0, 30.0)) FIR_b = Parameter(description="FIR term: b coefficient", - default=0.00) + default=0.00) FIR_n = Parameter(description="FIR term: n coefficient", - default=2.0, fixed=True) + default=2.0, fixed=True) x_range = x_range_P92 @@ -872,7 +911,7 @@ def _p92_single_term(in_lambda, amplitude, cen_wave, b, n): """ l_norm = in_lambda/cen_wave - return amplitude/(np.power(l_norm,n) + np.power(l_norm,-1*n) + b) + return amplitude/(np.power(l_norm, n) + np.power(l_norm, -1*n) + b) def evaluate(self, in_x, BKG_amp, BKG_lambda, BKG_b, BKG_n, @@ -917,11 +956,14 @@ def evaluate(self, in_x, # calculate the terms lam = 1.0/x axav = (self._p92_single_term(lam, BKG_amp, BKG_lambda, BKG_b, BKG_n) - + self._p92_single_term(lam, FUV_amp, FUV_lambda, FUV_b, FUV_n) - + self._p92_single_term(lam, NUV_amp, NUV_lambda, NUV_b, NUV_n) - + self._p92_single_term(lam, SIL1_amp, SIL1_lambda, SIL1_b, SIL1_n) - + self._p92_single_term(lam, SIL2_amp, SIL2_lambda, SIL2_b, SIL2_n) - + self._p92_single_term(lam, FIR_amp, FIR_lambda, FIR_b, FIR_n)) + + self._p92_single_term(lam, FUV_amp, FUV_lambda, FUV_b, FUV_n) + + self._p92_single_term(lam, NUV_amp, NUV_lambda, NUV_b, NUV_n) + + self._p92_single_term(lam, SIL1_amp, SIL1_lambda, + SIL1_b, SIL1_n) + + self._p92_single_term(lam, SIL2_amp, SIL2_lambda, + SIL2_b, SIL2_n) + + self._p92_single_term(lam, FIR_amp, FIR_lambda, + FIR_b, FIR_n)) # return A(x)/A(V) return axav @@ -929,6 +971,7 @@ def evaluate(self, in_x, # use numerical derivaties (need to add analytic) fit_deriv = None + class O94(BaseExtRvModel): """ O94 extinction model calculation @@ -977,8 +1020,7 @@ class O94(BaseExtRvModel): ax.legend(loc='best') plt.show() """ - - Rv_range = [2.0,6.0] + Rv_range = [2.0, 6.0] x_range = x_range_O94 @staticmethod @@ -1022,11 +1064,11 @@ def evaluate(in_x, Rv): b = np.zeros(n_x) # define the ranges - ir_indxs = np.where(np.logical_and(0.3 <= x,x < 1.1)) - opt_indxs = np.where(np.logical_and(1.1 <= x,x < 3.3)) - nuv_indxs = np.where(np.logical_and(3.3 <= x,x <= 8.0)) - fnuv_indxs = np.where(np.logical_and(5.9 <= x,x <= 8)) - fuv_indxs = np.where(np.logical_and(8 < x,x <= 10)) + ir_indxs = np.where(np.logical_and(0.3 <= x, x < 1.1)) + opt_indxs = np.where(np.logical_and(1.1 <= x, x < 3.3)) + nuv_indxs = np.where(np.logical_and(3.3 <= x, x <= 8.0)) + fnuv_indxs = np.where(np.logical_and(5.9 <= x, x <= 8)) + fuv_indxs = np.where(np.logical_and(8 < x, x <= 10)) # Infrared y = x[ir_indxs]**1.61 @@ -1041,11 +1083,11 @@ def evaluate(in_x, Rv): -7.985, -3.989, 2.908, 1.952, 0), y) # NUV - a[nuv_indxs] = 1.752-.316*x[nuv_indxs] \ - - 0.104/((x[nuv_indxs] - 4.67)**2 + .341) - b[nuv_indxs] = -3.09 + \ - 1.825*x[nuv_indxs] \ - + 1.206/((x[nuv_indxs] - 4.62)**2 + .263) + a[nuv_indxs] = (1.752-.316*x[nuv_indxs] + - 0.104/((x[nuv_indxs] - 4.67)**2 + .341)) + b[nuv_indxs] = (-3.09 + + 1.825*x[nuv_indxs] + + 1.206/((x[nuv_indxs] - 4.62)**2 + .263)) # far-NUV y = x[fnuv_indxs] - 5.9 @@ -1060,6 +1102,7 @@ def evaluate(in_x, Rv): # return A(x)/A(V) return a + b/Rv + class F99(BaseExtRvModel): """ F99 extinction model calculation @@ -1100,7 +1143,8 @@ class F99(BaseExtRvModel): text_model = F99() # generate the curves and plot them - x = np.arange(text_model.x_range[0], text_model.x_range[1],0.1)/u.micron + x = np.arange(text_model.x_range[0], + text_model.x_range[1],0.1)/u.micron Rvs = ['2.0','3.0','4.0','5.0','6.0'] for cur_Rv in Rvs: @@ -1113,8 +1157,7 @@ class F99(BaseExtRvModel): ax.legend(loc='best') plt.show() """ - - Rv_range = [2.0,6.0] + Rv_range = [2.0, 6.0] x_range = x_range_F99 def evaluate(self, in_x, Rv): @@ -1154,8 +1197,8 @@ def evaluate(self, in_x, Rv): C1 = 2.030 - 3.007*C2 # spline points - optnir_axav_x = 10000./np.array([26500.0,12200.0,6000.0, - 5470.0,4670.0,4110.0]) + optnir_axav_x = 10000./np.array([26500.0, 12200.0, 6000.0, + 5470.0, 4670.0, 4110.0]) # determine optical/IR values at spline points # Final term has a "-1.208" in Table 4 of F99, but then does @@ -1170,14 +1213,15 @@ def evaluate(self, in_x, Rv): -0.050 + 1.0016*Rv, 0.701 + 1.0016*Rv, 1.208 + 1.0032*Rv - 0.00033*(Rv**2)]) - nir_axebv_y = np.array([0.265,0.829])*Rv/3.1 - optnir_axebv_y = np.concatenate([nir_axebv_y,opt_axebv_y]) + nir_axebv_y = np.array([0.265, 0.829])*Rv/3.1 + optnir_axebv_y = np.concatenate([nir_axebv_y, opt_axebv_y]) # return A(x)/A(V) return _curve_F99_method(in_x, Rv, C1, C2, C3, C4, xo, gamma, optnir_axav_x, optnir_axebv_y/Rv, self.x_range, 'F99') + class G03_SMCBar(BaseExtAve): """ G03 SMCBar Average Extinction Curve @@ -1292,6 +1336,7 @@ def evaluate(self, in_x): optnir_axav_x, optnir_axav_y, self.x_range, 'G03') + class G03_LMCAvg(BaseExtAve): """ G03 LMCAvg Average Extinction Curve @@ -1521,7 +1566,7 @@ def evaluate(self, in_x): self.x_range, 'G03') -class G09_MWAvg(BaseExtAve): +class GCC09_MWAvg(BaseExtAve): """ G09 MW Average Extinction Curve @@ -1548,7 +1593,7 @@ class G09_MWAvg(BaseExtAve): import matplotlib.pyplot as plt import astropy.units as u - from dust_extinction.dust_extinction import G09_MWAvg + from dust_extinction.dust_extinction import GCC09_MWAvg fig, ax = plt.subplots() @@ -1556,19 +1601,19 @@ class G09_MWAvg(BaseExtAve): x = np.arange(0.3,1.0/0.0912,0.1)/u.micron # define the extinction model - ext_model = G09_MWAvg() + ext_model = GCC09_MWAvg() # generate the curves and plot them x = np.arange(ext_model.x_range[0], ext_model.x_range[1],0.1)/u.micron - ax.plot(x,ext_model(x),label='G09_MWAvg') - ax.errorbar(ext_model.obsdata_x_fuse, ext_model.obsdata_axav_fuse, + ax.plot(x,ext_model(x),label='GCC09_MWAvg') + ax.errorbar(ext_model.obsdata_x_fuse, ext_model.obsdata_axav_fuse, yerr=ext_model.obsdata_axav_unc_fuse, fmt='ko', label='obsdata (FUSE)') - ax.errorbar(ext_model.obsdata_x_iue, ext_model.obsdata_axav_iue, + ax.errorbar(ext_model.obsdata_x_iue, ext_model.obsdata_axav_iue, yerr=ext_model.obsdata_axav_unc_iue, fmt='bs', label='obsdata (IUE)') - ax.errorbar(ext_model.obsdata_x_bands, ext_model.obsdata_axav_bands, + ax.errorbar(ext_model.obsdata_x_bands, ext_model.obsdata_axav_bands, yerr=ext_model.obsdata_axav_unc_bands, fmt='g^', label='obsdata (Opt/NIR)') @@ -1579,11 +1624,11 @@ class G09_MWAvg(BaseExtAve): plt.show() """ - x_range = x_range_G09 + x_range = x_range_GCC09 Rv = 3.1 - # G09 sigma clipped average of 75 sightlines + # GCC09 sigma clipped average of 75 sightlines # FUSE range obsdata_x_fuse = [10.9546, 10.9103, 10.8662, 10.8223, 10.7785, 10.7349, 10.6915, 10.6483, 10.6052, 10.5623, @@ -1618,7 +1663,7 @@ class G09_MWAvg(BaseExtAve): 0.135772, 0.133049, 0.113491, 0.112876, 0.130662, 0.112918, 0.112409, 0.134404, 0.108662, 0.110052, 0.127927, 0.105854, 0.101012, 0.101105, 0.103142, - 0.100781, 0.0992076, 0.100922, 0.0980583, 0.095052, + 0.100781, 0.0992076, 0.100922, 0.0980583, .095052, 0.0977723, 0.0999428, 0.0968485, 0.0935959, 0.0962254, 0.0934298, 0.0901496, 0.0914562, 0.0869283, 0.0827114, 0.0839611, 0.0824688, @@ -1709,66 +1754,68 @@ class G09_MWAvg(BaseExtAve): 1.74919, 1.73011, 1.7181, 1.70523, 1.69983, 1.70104, 1.67963, 1.69661, 1.69288, 1.64255, 1.69986, 1.64438] obsdata_axav_iue = np.array(obsdata_axav_iue) - obsdata_axav_unc_iue = [0.0872464, 0.0870429, 0.0884721, 0.0856429, 0.085733, - 0.0870081, 0.0868195, 0.0824011, 0.0831442, 0.0700169, - 0.0753427, 0.0706412, 0.0708828, 0.0725433, 0.0718813, - 0.0733914, 0.0746335, 0.070549, 0.0701482, 0.0751515, - 0.0703968, 0.0704839, 0.0716527, 0.0665725, 0.0645906, - 0.0646839, 0.0667977, 0.0638441, 0.0644936, 0.0646092, - 0.0632691, 0.062183, 0.0614076, 0.0609257, 0.062895, - 0.0621621, 0.0616836, 0.0667685, 0.0637902, 0.0628244, - 0.06145, 0.0599388, 0.0582835, 0.0588209, 0.0552912, - 0.0570526, 0.0570427, 0.0551712, 0.0551769, 0.0532115, - 0.0537545, 0.05434, 0.0542661, 0.0532496, 0.052242, - 0.0526348, 0.0528077, 0.0543117, 0.0530651, 0.0543711, - 0.0535841, 0.0527281, 0.0615197, 0.0744345, 0.0604337, - 0.0540828, 0.0502652, 0.0493793, 0.0493413, 0.0494069, - 0.0490926, 0.0485171, 0.0476265, 0.0488545, 0.048111, - 0.0465552, 0.0456084, 0.0465916, 0.0438029, 0.0443429, - 0.0436129, 0.0431271, 0.0438585, 0.0437505, 0.0432799, - 0.0428038, 0.0418432, 0.0427744, 0.0426111, 0.0430678, - 0.043576, 0.0449369, 0.0427624, 0.0420389, 0.0433408, - 0.0431901, 0.0408719, 0.0410534, 0.0432204, 0.0428762, - 0.0412307, 0.0434997, 0.0390442, 0.0404766, 0.0409797, - 0.0407884, 0.0407895, 0.0418381, 0.0441335, 0.0435533, - 0.0417593, 0.0411644, 0.0394862, 0.0421479, 0.0443296, - 0.0444373, 0.0450305, 0.0445311, 0.0444786, 0.0454764, - 0.0437134, 0.0434556, 0.0441117, 0.045317, 0.0426919, - 0.0407817, 0.0415202, 0.0436795, 0.0460842, 0.0441931, - 0.0429678, 0.0416575, 0.0405409, 0.0449516, 0.0422263, - 0.04354, 0.0447088, 0.0427848, 0.0427436, 0.0454757, - 0.0459811, 0.0452352, 0.0473206, 0.0440146, 0.0460769, - 0.0476268, 0.0476545, 0.0466497, 0.0461933, 0.0468001, - 0.0456869, 0.0484227, 0.0460075, 0.0460413, 0.0450113, - 0.0462236, 0.0468239, 0.0438579, 0.0437953, 0.0440855, - 0.0427247, 0.0428557, 0.040375, 0.0404052, 0.038604, - 0.0383919, 0.0376932, 0.0359274, 0.0366497, 0.0347373, - 0.0364086, 0.0331321, 0.0339372, 0.0342494, 0.0340913, - 0.0320641, 0.0325347, 0.03207, 0.0303993, 0.0301872, - 0.0304542, 0.0295561, 0.0285186, 0.0300075, 0.0277648, - 0.0288373, 0.0272505, 0.0259275, 0.026904, 0.0273099, - 0.0270734, 0.0255038, 0.0250177, 0.0253357, 0.0243852, - 0.0235365, 0.0243051, 0.0234035, 0.0243984, 0.0238785, - 0.0233703, 0.0230637, 0.0237698, 0.0217698, 0.0231451, - 0.0234694, 0.0221487, 0.0220043, 0.0225651, 0.0218849, - 0.0214238, 0.0226075, 0.0199234, 0.0210254, 0.0208846, - 0.0208468, 0.021534, 0.0201173, 0.0219919, 0.019761, - 0.0205099, 0.0203457, 0.0206186, 0.0210671, 0.0209995, - 0.0201318, 0.0199858, 0.0191492, 0.0185285, 0.0190833, - 0.0181813, 0.0176767, 0.0188425, 0.0177012, 0.0170917, - 0.0177921, 0.0158575, 0.0162902, 0.0180801, 0.0178876, - 0.0183926, 0.0195394, 0.0186493, 0.0179556, 0.0174061, - 0.017488] - obsdata_axav_unc_iue = np.array(obsdata_axav_unc_iue) + axav_unc_iue = [0.0872464, 0.0870429, 0.0884721, 0.0856429, 0.085733, + 0.0870081, 0.0868195, 0.0824011, 0.0831442, 0.0700169, + 0.0753427, 0.0706412, 0.0708828, 0.0725433, 0.0718813, + 0.0733914, 0.0746335, 0.070549, 0.0701482, 0.0751515, + 0.0703968, 0.0704839, 0.0716527, 0.0665725, 0.0645906, + 0.0646839, 0.0667977, 0.0638441, 0.0644936, 0.0646092, + 0.0632691, 0.062183, 0.0614076, 0.0609257, 0.062895, + 0.0621621, 0.0616836, 0.0667685, 0.0637902, 0.0628244, + 0.06145, 0.0599388, 0.0582835, 0.0588209, 0.0552912, + 0.0570526, 0.0570427, 0.0551712, 0.0551769, 0.0532115, + 0.0537545, 0.05434, 0.0542661, 0.0532496, 0.052242, + 0.0526348, 0.0528077, 0.0543117, 0.0530651, 0.0543711, + 0.0535841, 0.0527281, 0.0615197, 0.0744345, 0.0604337, + 0.0540828, 0.0502652, 0.0493793, 0.0493413, 0.0494069, + 0.0490926, 0.0485171, 0.0476265, 0.0488545, 0.048111, + 0.0465552, 0.0456084, 0.0465916, 0.0438029, 0.0443429, + 0.0436129, 0.0431271, 0.0438585, 0.0437505, 0.0432799, + 0.0428038, 0.0418432, 0.0427744, 0.0426111, 0.0430678, + 0.043576, 0.0449369, 0.0427624, 0.0420389, 0.0433408, + 0.0431901, 0.0408719, 0.0410534, 0.0432204, 0.0428762, + 0.0412307, 0.0434997, 0.0390442, 0.0404766, 0.0409797, + 0.0407884, 0.0407895, 0.0418381, 0.0441335, 0.0435533, + 0.0417593, 0.0411644, 0.0394862, 0.0421479, 0.0443296, + 0.0444373, 0.0450305, 0.0445311, 0.0444786, 0.0454764, + 0.0437134, 0.0434556, 0.0441117, 0.045317, 0.0426919, + 0.0407817, 0.0415202, 0.0436795, 0.0460842, 0.0441931, + 0.0429678, 0.0416575, 0.0405409, 0.0449516, 0.0422263, + 0.04354, 0.0447088, 0.0427848, 0.0427436, 0.0454757, + 0.0459811, 0.0452352, 0.0473206, 0.0440146, 0.0460769, + 0.0476268, 0.0476545, 0.0466497, 0.0461933, 0.0468001, + 0.0456869, 0.0484227, 0.0460075, 0.0460413, 0.0450113, + 0.0462236, 0.0468239, 0.0438579, 0.0437953, 0.0440855, + 0.0427247, 0.0428557, 0.040375, 0.0404052, 0.038604, + 0.0383919, 0.0376932, 0.0359274, 0.0366497, 0.0347373, + 0.0364086, 0.0331321, 0.0339372, 0.0342494, 0.0340913, + 0.0320641, 0.0325347, 0.03207, 0.0303993, 0.0301872, + 0.0304542, 0.0295561, 0.0285186, 0.0300075, 0.0277648, + 0.0288373, 0.0272505, 0.0259275, 0.026904, 0.0273099, + 0.0270734, 0.0255038, 0.0250177, 0.0253357, 0.0243852, + 0.0235365, 0.0243051, 0.0234035, 0.0243984, 0.0238785, + 0.0233703, 0.0230637, 0.0237698, 0.0217698, 0.0231451, + 0.0234694, 0.0221487, 0.0220043, 0.0225651, 0.0218849, + 0.0214238, 0.0226075, 0.0199234, 0.0210254, 0.0208846, + 0.0208468, 0.021534, 0.0201173, 0.0219919, 0.019761, + 0.0205099, 0.0203457, 0.0206186, 0.0210671, 0.0209995, + 0.0201318, 0.0199858, 0.0191492, 0.0185285, 0.0190833, + 0.0181813, 0.0176767, 0.0188425, 0.0177012, 0.0170917, + 0.0177921, 0.0158575, 0.0162902, 0.0180801, 0.0178876, + 0.0183926, 0.0195394, 0.0186493, 0.0179556, 0.0174061, + 0.017488] + obsdata_axav_unc_iue = np.array(axav_unc_iue) # Opt/NIR range - obsdata_x_bands = np.array([2.73224, 2.28311, 0.819672, 0.613497, 0.456621]) + obsdata_x_bands = np.array([2.73224, 2.28311, 0.819672, 0.613497, + 0.456621]) obsdata_axav_bands = np.array([1.53296, 1.30791, 0.291042, 0.188455, 0.095588]) obsdata_axav_unc_bands = np.array([0.0105681, 0.00506663, 0.00407895, - 0.00307513, 0.00371036]) + 0.00307513, 0.00371036]) # put them together - obsdata_x = np.concatenate((obsdata_x_fuse, obsdata_x_iue, obsdata_x_bands)) + obsdata_x = np.concatenate((obsdata_x_fuse, obsdata_x_iue, + obsdata_x_bands)) obsdata_axav = np.concatenate((obsdata_axav_fuse, obsdata_axav_iue, obsdata_axav_bands)) obsdata_axav_unc = np.concatenate((obsdata_axav_unc_fuse, @@ -1780,7 +1827,7 @@ class G09_MWAvg(BaseExtAve): def evaluate(self, in_x): """ - G09_MWAvg function + GCC09_MWAvg function Parameters ---------- @@ -1801,16 +1848,14 @@ def evaluate(self, in_x): Input x values outside of defined range """ # convert to wavenumbers (1/micron) if x input in units - # otherwise, assume x in appropriate wavenumber units with u.add_enabled_equivalencies(u.spectral()): x_quant = u.Quantity(in_x, 1.0/u.micron, dtype=np.float64) # strip the quantity to avoid needing to add units to all the # polynomical coefficients x = x_quant.value - # check that the wavenumbers are within the defined range - _test_valid_x_range(x, x_range_G09, 'G09') + _test_valid_x_range(x, x_range_GCC09, 'GCC09') # P92 parameters fit to the data using uncs as weights p92_fit = P92(BKG_amp=203.805939127, BKG_lambda=0.0508199427208, @@ -1876,7 +1921,8 @@ class G16(BaseExtRvAfAModel): text_model = G16() # generate the curves and plot them - x = np.arange(text_model.x_range[0], text_model.x_range[1],0.1)/u.micron + x = np.arange(text_model.x_range[0], + text_model.x_range[1],0.1)/u.micron Rvs = ['2.0','3.0','4.0','5.0','6.0'] for cur_Rv in Rvs: @@ -1904,7 +1950,8 @@ class G16(BaseExtRvAfAModel): text_model = G16() # generate the curves and plot them - x = np.arange(text_model.x_range[0], text_model.x_range[1],0.1)/u.micron + x = np.arange(text_model.x_range[0], + text_model.x_range[1],0.1)/u.micron fAs = [0.0, 0.2, 0.4, 0.6, 0.8, 1.0] for cur_fA in fAs: diff --git a/dust_extinction/tests/test_g09_ave.py b/dust_extinction/tests/test_gcc09_ave.py similarity index 76% rename from dust_extinction/tests/test_g09_ave.py rename to dust_extinction/tests/test_gcc09_ave.py index ba3392e..5bd53b2 100644 --- a/dust_extinction/tests/test_g09_ave.py +++ b/dust_extinction/tests/test_gcc09_ave.py @@ -4,66 +4,72 @@ import astropy.units as u from astropy.modeling import InputParameterError -from ..dust_extinction import G09_MWAvg +from ..dust_extinction import GCC09_MWAvg x_bad = [-1.0, 0.1, 11., 100.] + @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): - tmodel = G09_MWAvg() + tmodel = GCC09_MWAvg() with pytest.raises(ValueError) as exc: tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for G09' \ + assert exc.value.args[0] == 'Input x outside of range defined for GCC09' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): - tmodel = G09_MWAvg() + tmodel = GCC09_MWAvg() with pytest.raises(ValueError) as exc: tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for G09' \ + assert exc.value.args[0] == 'Input x outside of range defined for GCC09' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): - tmodel = G09_MWAvg() + tmodel = GCC09_MWAvg() with pytest.raises(ValueError) as exc: tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for G09' \ + assert exc.value.args[0] == 'Input x outside of range defined for GCC09' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): - tmodel = G09_MWAvg() + tmodel = GCC09_MWAvg() with pytest.raises(ValueError) as exc: tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for G09' \ + assert exc.value.args[0] == 'Input x outside of range defined for GCC09' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + def test_extinguish_no_av_or_ebv(): - tmodel = G09_MWAvg() + tmodel = GCC09_MWAvg() with pytest.raises(InputParameterError) as exc: tmodel.extinguish([1.0]) assert exc.value.args[0] == 'neither Av or Ebv passed, one required' -def test_extinction_G09_values(): - tmodel = G09_MWAvg() + +def test_extinction_GCC09_values(): + tmodel = GCC09_MWAvg() # test # not to numerical precision as we are using the FM90 fits # and spline functions and the correct values are the data @@ -71,36 +77,39 @@ def test_extinction_G09_values(): tmodel.obsdata_axav, rtol=tmodel.obsdata_tolerance) -def test_extinction_G09_single_values(): - tmodel = G09_MWAvg() + +def test_extinction_GCC09_single_values(): + tmodel = GCC09_MWAvg() # test for x, cor_val in zip(tmodel.obsdata_x, tmodel.obsdata_axav): np.testing.assert_allclose(tmodel(x), cor_val, rtol=tmodel.obsdata_tolerance) -def test_extinction_G09_extinguish_values_Av(): - tmodel = G09_MWAvg() - x = np.arange(0.3,1.0/0.0912,0.1)/u.micron +def test_extinction_GCC09_extinguish_values_Av(): + tmodel = GCC09_MWAvg() + + x = np.arange(0.3, 1.0/0.0912, 0.1)/u.micron cor_vals = tmodel(x) # calculate the cor_vals in fractional units Av = 1.0 - cor_vals = np.power(10.0,-0.4*cor_vals*Av) + cor_vals = np.power(10.0, -0.4*cor_vals*Av) # test np.testing.assert_equal(tmodel.extinguish(x, Av=Av), cor_vals) -def test_extinction_G09_extinguish_values_Ebv(): - tmodel = G09_MWAvg() - x = np.arange(0.3,1.0/0.0912,0.1)/u.micron +def test_extinction_GCC09_extinguish_values_Ebv(): + tmodel = GCC09_MWAvg() + + x = np.arange(0.3, 1.0/0.0912, 0.1)/u.micron cor_vals = tmodel(x) # calculate the cor_vals in fractional units Ebv = 1.0 Av = Ebv*tmodel.Rv - cor_vals = np.power(10.0,-0.4*cor_vals*Av) + cor_vals = np.power(10.0, -0.4*cor_vals*Av) # test np.testing.assert_equal(tmodel.extinguish(x, Ebv=Ebv), cor_vals) From a0fc22a7bb7efb41cadb72189059d360434747ac Mon Sep 17 00:00:00 2001 From: Karl Gordon Date: Fri, 12 Jan 2018 17:46:25 -0500 Subject: [PATCH 2/6] Test for AlAvtoElv model --- docs/dust_extinction/model_flavors.rst | 3 +-- dust_extinction/tests/test_p92.py | 1 - 2 files changed, 1 insertion(+), 3 deletions(-) diff --git a/docs/dust_extinction/model_flavors.rst b/docs/dust_extinction/model_flavors.rst index 97bc401..e6cf93f 100644 --- a/docs/dust_extinction/model_flavors.rst +++ b/docs/dust_extinction/model_flavors.rst @@ -31,8 +31,7 @@ Average models # generate the curves and plot them x = np.arange(0.3,10.0,0.1)/u.micron - - ext_model = G09_MWAvg() + ext_model = GCC09_MWAvg() ax.plot(x,ext_model(x),label='GCC09 MWAvg') ext_model = G03_SMCBar() diff --git a/dust_extinction/tests/test_p92.py b/dust_extinction/tests/test_p92.py index a22d654..8e4478e 100644 --- a/dust_extinction/tests/test_p92.py +++ b/dust_extinction/tests/test_p92.py @@ -72,7 +72,6 @@ def get_axav_cor_vals(): Rv = 3.08 MW_axav = MW_exvebv/Rv + 1.0 - # add units x = MW_x/u.micron From c89c64ea6558facaaa78523e6bd1915ef0f30a40 Mon Sep 17 00:00:00 2001 From: Karl Gordon Date: Fri, 12 Jan 2018 17:46:50 -0500 Subject: [PATCH 3/6] Actual test added now --- dust_extinction/tests/test_alavtoelv.py | 63 +++++++++++++++++++++++++ 1 file changed, 63 insertions(+) create mode 100644 dust_extinction/tests/test_alavtoelv.py diff --git a/dust_extinction/tests/test_alavtoelv.py b/dust_extinction/tests/test_alavtoelv.py new file mode 100644 index 0000000..88a548a --- /dev/null +++ b/dust_extinction/tests/test_alavtoelv.py @@ -0,0 +1,63 @@ +import numpy as np +import astropy.units as u + +from astropy.modeling.fitting import LevMarLSQFitter + +from ..dust_extinction import (P92, AlAvToElv) + + +def get_axav_cor_vals(): + + # 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 + + return (MW_x/u.micron, MW_axav) + + +class P92_Elv(P92 | AlAvToElv): + """ + Evalute P92 on E(l-V) data including solving for A(V) + """ + + +def test_P92_fitting(): + + # get an observed extinction curve to fit + x_quant, y_axav = get_axav_cor_vals() + x = x_quant.value + + # pick an A(V) value + av = 1.3 + + # transform the extinction from AxAv to Exv + y = (y_axav - 1.0)*av + + # change from defaults to make the best fit harder to find + p92_init = P92_Elv() + + fit = LevMarLSQFitter() + # accuracy set to avoid warning the fit may have failed + p92_fit = fit(p92_init, x, y, acc=1e-3) + + fit_vals = p92_fit._parameters + + good_vals = [220.806743348, 0.0474488526572, 88.4203515994, 2.0, + 17.0089730125, 0.0729612929387, 1.38876455935, 5.63200040759, + 0.052051591862, 0.218617616831, -1.95160129592, 2.0, + 0.0, 7.0, -278.831631235, 2.0, + 0.0, 21.0, -6617.22039032, 2.0, + 0.0, 30.0, -4992.15823779, 2.0, + 1.29867067922] + + np.testing.assert_allclose(good_vals, fit_vals) From 161ce4ca1f8b8489d28ad6b1526670692d79e3b5 Mon Sep 17 00:00:00 2001 From: Karl Gordon Date: Fri, 12 Jan 2018 18:17:45 -0500 Subject: [PATCH 4/6] Updates for PEP8 --- dust_extinction/tests/test_ccm89.py | 96 ++++++++++++++++------------- dust_extinction/tests/test_f99.py | 25 +++++--- dust_extinction/tests/test_fm90.py | 15 +++-- dust_extinction/tests/test_g03.py | 24 +++++--- dust_extinction/tests/test_g16.py | 26 +++++--- dust_extinction/tests/test_o94.py | 85 ++++++++++++++----------- dust_extinction/tests/test_p92.py | 30 ++++++--- 7 files changed, 183 insertions(+), 118 deletions(-) diff --git a/dust_extinction/tests/test_ccm89.py b/dust_extinction/tests/test_ccm89.py index 5f0c589..30bdd5e 100644 --- a/dust_extinction/tests/test_ccm89.py +++ b/dust_extinction/tests/test_ccm89.py @@ -6,12 +6,14 @@ from ..dust_extinction import CCM89 -@pytest.mark.parametrize("Rv_invalid", [-1.0,0.0,1.9,6.1,10.]) + +@pytest.mark.parametrize("Rv_invalid", [-1.0, 0.0, 1.9, 6.1, 10.]) def test_invalid_Rv_input(Rv_invalid): with pytest.raises(InputParameterError) as exc: tmodel = CCM89(Rv=Rv_invalid) assert exc.value.args[0] == 'parameter Rv must be between 2.0 and 6.0' + @pytest.mark.parametrize("x_invalid", [-1.0, 0.2, 10.1, 100.]) def test_invalid_wavenumbers(x_invalid): tmodel = CCM89(Rv=3.1) @@ -20,10 +22,11 @@ def test_invalid_wavenumbers(x_invalid): assert exc.value.args[0] == 'Input x outside of range defined for CCM89' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_wavenumber", [-1.0, 0.2, 10.1, 100.]/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): @@ -33,10 +36,11 @@ def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): assert exc.value.args[0] == 'Input x outside of range defined for CCM89' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_micron", u.micron/[-1.0, 0.2, 10.1, 100.]) def test_invalid_micron(x_invalid_micron): @@ -46,10 +50,11 @@ def test_invalid_micron(x_invalid_micron): assert exc.value.args[0] == 'Input x outside of range defined for CCM89' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/[-1.0, 0.2, 10.1, 100.]) def test_invalid_micron(x_invalid_angstrom): @@ -59,10 +64,11 @@ def test_invalid_micron(x_invalid_angstrom): assert exc.value.args[0] == 'Input x outside of range defined for CCM89' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + def test_axav_ccm89_table3(): # values from Table 3 of Cardelli et al. (1989) # ignoring the last value at L band as it is outside the @@ -80,56 +86,57 @@ def test_axav_ccm89_table3(): # test (table in paper has limited precision) np.testing.assert_allclose(tmodel(x), cor_vals, atol=1e-2) + def get_axav_cor_vals(Rv): # testing wavenumbers - x = np.array([ 10. , 9. , 8. , 7. , - 6. , 5. , 4.6 , 4. , - 3. , 2.78, 2.27, 1.82, - 1.43, 1.11, 0.8 , 0.63, - 0.46]) + x = np.array([10., 9., 8., 7., + 6., 5., 4.6, 4., + 3., 2.78, 2.27, 1.82, + 1.43, 1.11, 0.8, 0.63, + 0.46]) # add units x = x/u.micron # correct values if Rv == 3.1: - cor_vals = np.array([ 5.23835484, 4.13406452, 3.33685933, 2.77962453, - 2.52195399, 2.84252644, 3.18598916, 2.31531711, - 1.64254927, 1.56880904, 1.32257836, 1. , - 0.75125994, 0.4780346 , 0.28206957, 0.19200814, - 0.11572348]) + cor_vals = np.array([5.23835484, 4.13406452, 3.33685933, 2.77962453, + 2.52195399, 2.84252644, 3.18598916, 2.31531711, + 1.64254927, 1.56880904, 1.32257836, 1., + 0.75125994, 0.47803460, 0.28206957, 0.19200814, + 0.11572348]) elif Rv == 2.0: - cor_vals = np.array([ 9.407 , 7.3065 , 5.76223881, 4.60825807, - 4.01559036, 4.43845534, 4.93952892, 3.39275574, - 2.068771 , 1.9075018 , 1.49999733, 1. , - 0.68650255, 0.36750326, 0.21678862, 0.14757062, - 0.08894094]) + cor_vals = np.array([9.407, 7.3065, 5.76223881, 4.60825807, + 4.01559036, 4.43845534, 4.93952892, 3.39275574, + 2.068771, 1.9075018, 1.49999733, 1., + 0.68650255, 0.36750326, 0.21678862, 0.14757062, + 0.08894094]) elif Rv == 3.0: - cor_vals = np.array([ 5.491 , 4.32633333, 3.48385202, 2.8904508 , - 2.6124774 , 2.9392494 , 3.2922643 , 2.38061642, - 1.66838089, 1.58933588, 1.33333103, 1. , - 0.74733525, 0.47133573, 0.27811315, 0.18931496, - 0.11410029]) + cor_vals = np.array([5.491, 4.32633333, 3.48385202, 2.8904508, + 2.6124774, 2.9392494, 3.2922643, 2.38061642, + 1.66838089, 1.58933588, 1.33333103, 1., + 0.74733525, 0.47133573, 0.27811315, 0.18931496, + 0.11410029]) elif Rv == 4.0: - cor_vals = np.array([ 3.533 , 2.83625 , 2.34465863, 2.03154717, - 1.91092092, 2.18964643, 2.46863199, 1.87454675, - 1.46818583, 1.43025292, 1.24999788, 1. , - 0.7777516 , 0.52325196, 0.30877542, 0.21018713, - 0.12667997]) + cor_vals = np.array([3.533, 2.83625, 2.34465863, 2.03154717, + 1.91092092, 2.18964643, 2.46863199, 1.87454675, + 1.46818583, 1.43025292, 1.24999788, 1., + 0.7777516, 0.52325196, 0.30877542, 0.21018713, + 0.12667997]) elif Rv == 5.0: - cor_vals = np.array([ 2.3582 , 1.9422 , 1.66114259, 1.51620499, - 1.48998704, 1.73988465, 1.97445261, 1.57090496, - 1.3480688 , 1.33480314, 1.19999799, 1. , - 0.79600141, 0.5544017 , 0.32717278, 0.22271044, - 0.13422778]) + cor_vals = np.array([2.3582, 1.9422, 1.66114259, 1.51620499, + 1.48998704, 1.73988465, 1.97445261, 1.57090496, + 1.3480688, 1.33480314, 1.19999799, 1., + 0.79600141, 0.5544017, 0.32717278, 0.22271044, + 0.13422778]) elif Rv == 6.0: - cor_vals = np.array([ 1.575 , 1.34616667, 1.20546523, 1.17264354, - 1.20936444, 1.44004346, 1.64499968, 1.36847709, - 1.26799077, 1.27116996, 1.16666472, 1. , - 0.80816794, 0.5751682 , 0.33943769, 0.23105931, - 0.13925965]) + cor_vals = np.array([1.575, 1.34616667, 1.20546523, 1.17264354, + 1.20936444, 1.44004346, 1.64499968, 1.36847709, + 1.26799077, 1.27116996, 1.16666472, 1., + 0.80816794, 0.5751682, 0.33943769, 0.23105931, + 0.13925965]) else: - cor_vals = np.array([ 0.0 ]) + cor_vals = np.array([0.0]) return (x, cor_vals) @@ -145,12 +152,14 @@ def test_extinction_CCM89_values(Rv): # test np.testing.assert_allclose(tmodel(x), cor_vals) + def test_extinguish_no_av_or_ebv(): tmodel = CCM89() with pytest.raises(InputParameterError) as exc: tmodel.extinguish([1.0]) assert exc.value.args[0] == 'neither Av or Ebv passed, one required' + @pytest.mark.parametrize("Rv", [2.0, 3.0, 3.1, 4.0, 5.0, 6.0]) def test_extinction_CCM89_extinguish_values_Av(Rv): # get the correct values @@ -158,7 +167,7 @@ def test_extinction_CCM89_extinguish_values_Av(Rv): # calculate the cor_vals in fractional units Av = 1.0 - cor_vals = np.power(10.0,-0.4*(cor_vals*Av)) + cor_vals = np.power(10.0, -0.4*(cor_vals*Av)) # initialize extinction model tmodel = CCM89(Rv=Rv) @@ -166,6 +175,7 @@ def test_extinction_CCM89_extinguish_values_Av(Rv): # test np.testing.assert_allclose(tmodel.extinguish(x, Av=Av), cor_vals) + @pytest.mark.parametrize("Rv", [2.0, 3.0, 3.1, 4.0, 5.0, 6.0]) def test_extinction_CCM89_extinguish_values_Ebv(Rv): # get the correct values @@ -174,7 +184,7 @@ def test_extinction_CCM89_extinguish_values_Ebv(Rv): # calculate the cor_vals in fractional units Ebv = 1.0 Av = Ebv*Rv - cor_vals = np.power(10.0,-0.4*(cor_vals*Av)) + cor_vals = np.power(10.0, -0.4*(cor_vals*Av)) # initialize extinction model tmodel = CCM89(Rv=Rv) diff --git a/dust_extinction/tests/test_f99.py b/dust_extinction/tests/test_f99.py index aaba2df..e545cd3 100644 --- a/dust_extinction/tests/test_f99.py +++ b/dust_extinction/tests/test_f99.py @@ -2,11 +2,13 @@ import pytest import astropy.units as u -from astropy.modeling import InputParameterError from ..dust_extinction import F99 + x_bad = [-1.0, 0.1, 12.0, 100.] + + @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): tmodel = F99() @@ -15,10 +17,11 @@ def test_invalid_wavenumbers(x_invalid): assert exc.value.args[0] == 'Input x outside of range defined for F99' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): tmodel = F99() @@ -27,10 +30,11 @@ def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): assert exc.value.args[0] == 'Input x outside of range defined for F99' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): tmodel = F99() @@ -39,10 +43,11 @@ def test_invalid_micron(x_invalid_micron): assert exc.value.args[0] == 'Input x outside of range defined for F99' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): tmodel = F99() @@ -51,17 +56,19 @@ def test_invalid_micron(x_invalid_angstrom): assert exc.value.args[0] == 'Input x outside of range defined for F99' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + def get_axav_cor_vals(): # testing wavenumbers # from previous version of code - #x = np.array([0.5, 1.0, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, 8.0]) + # x = np.array([0.5, 1.0, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, + # 8.0]) # correct values (not quite right) - #cor_vals = np.array([0.124997, 0.377073, 1.130636, 1.419644, + # cor_vals = np.array([0.124997, 0.377073, 1.130636, 1.419644, # 1.630377, 1.888546, 2.275900, 3.014577, # 2.762256, 2.475272, 2.711508, 3.197144]) @@ -94,7 +101,9 @@ def test_extinction_F99_values(): x_vals, axav_vals, tolerance = get_axav_cor_vals() -test_vals = zip(x_vals, axav_vals, np.full(len(x_vals),tolerance)) +test_vals = zip(x_vals, axav_vals, np.full(len(x_vals), tolerance)) + + @pytest.mark.parametrize("test_vals", test_vals) def test_extinction_F99_single_values(test_vals): x, cor_val, tolerance = test_vals diff --git a/dust_extinction/tests/test_fm90.py b/dust_extinction/tests/test_fm90.py index 895e76a..ba1dfa4 100644 --- a/dust_extinction/tests/test_fm90.py +++ b/dust_extinction/tests/test_fm90.py @@ -2,12 +2,13 @@ import pytest import astropy.units as u -from astropy.modeling import InputParameterError from astropy.modeling.fitting import LevMarLSQFitter from ..dust_extinction import FM90, G03_LMCAvg x_bad = [-1.0, 0.2, 3.0, 11.0, 100.] + + @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): tmodel = FM90() @@ -16,10 +17,11 @@ def test_invalid_wavenumbers(x_invalid): assert exc.value.args[0] == 'Input x outside of range defined for FM90' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): tmodel = FM90() @@ -28,10 +30,11 @@ def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): assert exc.value.args[0] == 'Input x outside of range defined for FM90' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): tmodel = FM90() @@ -40,10 +43,11 @@ def test_invalid_micron(x_invalid_micron): assert exc.value.args[0] == 'Input x outside of range defined for FM90' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): tmodel = FM90() @@ -52,10 +56,11 @@ def test_invalid_micron(x_invalid_angstrom): assert exc.value.args[0] == 'Input x outside of range defined for FM90' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + def get_elvebv_cor_vals(): # testing wavenumbers x = np.array([3.5, 4.0, 4.5, 5.0, 6.0, 7.0, 9.0, 10.0]) diff --git a/dust_extinction/tests/test_g03.py b/dust_extinction/tests/test_g03.py index 6fbcee7..42bf070 100644 --- a/dust_extinction/tests/test_g03.py +++ b/dust_extinction/tests/test_g03.py @@ -18,10 +18,11 @@ def test_invalid_wavenumbers(x_invalid, tmodel): assert exc.value.args[0] == 'Input x outside of range defined for G03' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) @pytest.mark.parametrize("tmodel", models) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber, tmodel): @@ -30,10 +31,11 @@ def test_invalid_wavenumbers_imicron(x_invalid_wavenumber, tmodel): assert exc.value.args[0] == 'Input x outside of range defined for G03' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) @pytest.mark.parametrize("tmodel", models) def test_invalid_micron(x_invalid_micron, tmodel): @@ -42,10 +44,11 @@ def test_invalid_micron(x_invalid_micron, tmodel): assert exc.value.args[0] == 'Input x outside of range defined for G03' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) @pytest.mark.parametrize("tmodel", models) def test_invalid_micron(x_invalid_angstrom, tmodel): @@ -54,16 +57,18 @@ def test_invalid_micron(x_invalid_angstrom, tmodel): assert exc.value.args[0] == 'Input x outside of range defined for G03' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("tmodel", models) def test_extinguish_no_av_or_ebv(tmodel): with pytest.raises(InputParameterError) as exc: tmodel.extinguish([1.0]) assert exc.value.args[0] == 'neither Av or Ebv passed, one required' + @pytest.mark.parametrize("tmodel", models) def test_extinction_G03_values(tmodel): # test @@ -73,6 +78,7 @@ def test_extinction_G03_values(tmodel): tmodel.obsdata_axav, rtol=tmodel.obsdata_tolerance) + @pytest.mark.parametrize("tmodel", models) def test_extinction_G03_single_values(tmodel): # test @@ -80,27 +86,29 @@ def test_extinction_G03_single_values(tmodel): np.testing.assert_allclose(tmodel(x), cor_val, rtol=tmodel.obsdata_tolerance) + @pytest.mark.parametrize("tmodel", models) def test_extinction_G03_extinguish_values_Av(tmodel): - x = np.arange(0.3,10.0,0.1)/u.micron + x = np.arange(0.3, 10.0, 0.1)/u.micron cor_vals = tmodel(x) # calculate the cor_vals in fractional units Av = 1.0 - cor_vals = np.power(10.0,-0.4*(cor_vals*Av)) + cor_vals = np.power(10.0, -0.4*(cor_vals*Av)) # test np.testing.assert_equal(tmodel.extinguish(x, Av=Av), cor_vals) + @pytest.mark.parametrize("tmodel", models) def test_extinction_G03_extinguish_values_Ebv(tmodel): - x = np.arange(0.3,10.0,0.1)/u.micron + x = np.arange(0.3, 10.0, 0.1)/u.micron cor_vals = tmodel(x) # calculate the cor_vals in fractional units Ebv = 1.0 Av = Ebv*tmodel.Rv - cor_vals = np.power(10.0,-0.4*cor_vals*Av) + cor_vals = np.power(10.0, -0.4*cor_vals*Av) # test np.testing.assert_equal(tmodel.extinguish(x, Ebv=Ebv), cor_vals) diff --git a/dust_extinction/tests/test_g16.py b/dust_extinction/tests/test_g16.py index 9365810..bdd2310 100644 --- a/dust_extinction/tests/test_g16.py +++ b/dust_extinction/tests/test_g16.py @@ -8,19 +8,24 @@ from ..dust_extinction import G03_SMCBar from .test_f99 import get_axav_cor_vals as get_axav_cor_vals_fA_1 -@pytest.mark.parametrize("RvA_invalid", [-1.0,0.0,1.9,6.1,10.]) + +@pytest.mark.parametrize("RvA_invalid", [-1.0, 0.0, 1.9, 6.1, 10.]) def test_invalid_RvA_input(RvA_invalid): with pytest.raises(InputParameterError) as exc: tmodel = G16(RvA=RvA_invalid) assert exc.value.args[0] == 'parameter RvA must be between 2.0 and 6.0' -@pytest.mark.parametrize("fA_invalid", [-1.0,-0.1,1.1,10.0]) + +@pytest.mark.parametrize("fA_invalid", [-1.0, -0.1, 1.1, 10.0]) def test_invalid_fA_input(fA_invalid): with pytest.raises(InputParameterError) as exc: tmodel = G16(fA=fA_invalid) assert exc.value.args[0] == 'parameter fA must be between 0.0 and 1.0' + x_bad = [-1.0, 0.1, 12.0, 100.] + + @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): tmodel = G16() @@ -29,10 +34,11 @@ def test_invalid_wavenumbers(x_invalid): assert exc.value.args[0] == 'Input x outside of range defined for G16' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): tmodel = G16() @@ -41,10 +47,11 @@ def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): assert exc.value.args[0] == 'Input x outside of range defined for G16' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): tmodel = G16() @@ -53,10 +60,11 @@ def test_invalid_micron(x_invalid_micron): assert exc.value.args[0] == 'Input x outside of range defined for G16' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): tmodel = G16() @@ -65,7 +73,7 @@ def test_invalid_micron(x_invalid_angstrom): assert exc.value.args[0] == 'Input x outside of range defined for G16' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' @@ -80,6 +88,7 @@ def test_extinction_G16_fA_1_values(): # test np.testing.assert_allclose(tmodel(x), cor_vals, rtol=tolerance) + def test_extinction_G16_fA_0_values(): # initialize the model tmodel = G16(fA=0.0) @@ -93,8 +102,11 @@ def test_extinction_G16_fA_0_values(): # test np.testing.assert_allclose(tmodel(x), cor_vals, rtol=tolerance) + x_vals, axav_vals, tolerance = get_axav_cor_vals_fA_1() -test_vals = zip(x_vals, axav_vals, np.full(len(x_vals),tolerance)) +test_vals = zip(x_vals, axav_vals, np.full(len(x_vals), tolerance)) + + @pytest.mark.parametrize("test_vals", test_vals) def test_extinction_G16_fA_1_single_values(test_vals): x, cor_val, tolerance = test_vals diff --git a/dust_extinction/tests/test_o94.py b/dust_extinction/tests/test_o94.py index 00c1184..b8c7108 100644 --- a/dust_extinction/tests/test_o94.py +++ b/dust_extinction/tests/test_o94.py @@ -6,12 +6,14 @@ from ..dust_extinction import O94 -@pytest.mark.parametrize("Rv_invalid", [-1.0,0.0,1.9,6.1,10.]) + +@pytest.mark.parametrize("Rv_invalid", [-1.0, 0.0, 1.9, 6.1, 10.]) def test_invalid_Rv_input(Rv_invalid): with pytest.raises(InputParameterError) as exc: tmodel = O94(Rv=Rv_invalid) assert exc.value.args[0] == 'parameter Rv must be between 2.0 and 6.0' + @pytest.mark.parametrize("x_invalid", [-1.0, 0.2, 10.1, 100.]) def test_invalid_wavenumbers(x_invalid): tmodel = O94(Rv=3.1) @@ -20,10 +22,11 @@ def test_invalid_wavenumbers(x_invalid): assert exc.value.args[0] == 'Input x outside of range defined for O94' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_wavenumber", [-1.0, 0.2, 10.1, 100.]/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): @@ -33,10 +36,11 @@ def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): assert exc.value.args[0] == 'Input x outside of range defined for O94' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_micron", u.micron/[-1.0, 0.2, 10.1, 100.]) def test_invalid_micron(x_invalid_micron): @@ -46,10 +50,11 @@ def test_invalid_micron(x_invalid_micron): assert exc.value.args[0] == 'Input x outside of range defined for O94' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/[-1.0, 0.2, 10.1, 100.]) def test_invalid_micron(x_invalid_angstrom): @@ -59,10 +64,11 @@ def test_invalid_micron(x_invalid_angstrom): assert exc.value.args[0] == 'Input x outside of range defined for O94' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + def test_axav_o94_rv31(): # values from Bastiaasen (1992) Table 6 x = np.array([2.939, 2.863, 2.778, 2.642, 2.476, 2.385, @@ -81,54 +87,56 @@ def test_axav_o94_rv31(): mod_vals = tmodel(x) norm_vals = tmodel([1.5, 2.2]) mod_vals = (mod_vals - norm_vals[0])/(norm_vals[1] - norm_vals[0]) - + # test (table in paper has limited precision) np.testing.assert_allclose(mod_vals, cor_vals, atol=6e-2) + def get_axav_cor_vals(Rv): # testing only NIR or UV wavenumbers (optical tested in previous test) # O94 is the same as CCM89 for these wavelengths - x = np.array([ 10. , 9. , 8. , 7. , - 6. , 5. , 4.6 , 4. , - 0.8 , 0.63, - 0.46]) + + x = np.array([10., 9., 8., 7., + 6., 5., 4.6, 4., + 0.8, 0.63, + 0.46]) # add units x = x/u.micron # correct values if Rv == 3.1: - cor_vals = np.array([ 5.23835484, 4.13406452, 3.33685933, 2.77962453, - 2.52195399, 2.84252644, 3.18598916, 2.31531711, - 0.28206957, 0.19200814, - 0.11572348]) + cor_vals = np.array([5.23835484, 4.13406452, 3.33685933, 2.77962453, + 2.52195399, 2.84252644, 3.18598916, 2.31531711, + 0.28206957, 0.19200814, + 0.11572348]) elif Rv == 2.0: - cor_vals = np.array([ 9.407 , 7.3065 , 5.76223881, 4.60825807, - 4.01559036, 4.43845534, 4.93952892, 3.39275574, - 0.21678862, 0.14757062, - 0.08894094]) + cor_vals = np.array([9.407, 7.3065, 5.76223881, 4.60825807, + 4.01559036, 4.43845534, 4.93952892, 3.39275574, + 0.21678862, 0.14757062, + 0.08894094]) elif Rv == 3.0: - cor_vals = np.array([ 5.491 , 4.32633333, 3.48385202, 2.8904508 , - 2.6124774 , 2.9392494 , 3.2922643 , 2.38061642, - 0.27811315, 0.18931496, - 0.11410029]) + cor_vals = np.array([5.491, 4.32633333, 3.48385202, 2.8904508, + 2.6124774, 2.9392494, 3.2922643, 2.38061642, + 0.27811315, 0.18931496, + 0.11410029]) elif Rv == 4.0: - cor_vals = np.array([ 3.533 , 2.83625 , 2.34465863, 2.03154717, - 1.91092092, 2.18964643, 2.46863199, 1.87454675, - 0.30877542, 0.21018713, - 0.12667997]) + cor_vals = np.array([3.533, 2.83625, 2.34465863, 2.03154717, + 1.91092092, 2.18964643, 2.46863199, 1.87454675, + 0.30877542, 0.21018713, + 0.12667997]) elif Rv == 5.0: - cor_vals = np.array([ 2.3582 , 1.9422 , 1.66114259, 1.51620499, - 1.48998704, 1.73988465, 1.97445261, 1.57090496, - 0.32717278, 0.22271044, - 0.13422778]) + cor_vals = np.array([2.3582, 1.9422, 1.66114259, 1.51620499, + 1.48998704, 1.73988465, 1.97445261, 1.57090496, + 0.32717278, 0.22271044, + 0.13422778]) elif Rv == 6.0: - cor_vals = np.array([ 1.575 , 1.34616667, 1.20546523, 1.17264354, - 1.20936444, 1.44004346, 1.64499968, 1.36847709, - 0.33943769, 0.23105931, - 0.13925965]) + cor_vals = np.array([1.575, 1.34616667, 1.20546523, 1.17264354, + 1.20936444, 1.44004346, 1.64499968, 1.36847709, + 0.33943769, 0.23105931, + 0.13925965]) else: - cor_vals = np.array([ 0.0 ]) + cor_vals = np.array([0.0]) return (x, cor_vals) @@ -144,12 +152,14 @@ def test_extinction_O94_values(Rv): # test np.testing.assert_allclose(tmodel(x), cor_vals) + def test_extinguish_no_av_or_ebv(): tmodel = O94() with pytest.raises(InputParameterError) as exc: tmodel.extinguish([1.0]) assert exc.value.args[0] == 'neither Av or Ebv passed, one required' + @pytest.mark.parametrize("Rv", [2.0, 3.0, 3.1, 4.0, 5.0, 6.0]) def test_extinction_O94_extinguish_values_Av(Rv): # get the correct values @@ -157,7 +167,7 @@ def test_extinction_O94_extinguish_values_Av(Rv): # calculate the cor_vals in fractional units Av = 1.0 - cor_vals = np.power(10.0,-0.4*(cor_vals*Av)) + cor_vals = np.power(10.0, -0.4*(cor_vals*Av)) # initialize extinction model tmodel = O94(Rv=Rv) @@ -165,6 +175,7 @@ def test_extinction_O94_extinguish_values_Av(Rv): # test np.testing.assert_allclose(tmodel.extinguish(x, Av=Av), cor_vals) + @pytest.mark.parametrize("Rv", [2.0, 3.0, 3.1, 4.0, 5.0, 6.0]) def test_extinction_O94_extinguish_values_Ebv(Rv): # get the correct values @@ -173,7 +184,7 @@ def test_extinction_O94_extinguish_values_Ebv(Rv): # calculate the cor_vals in fractional units Ebv = 1.0 Av = Ebv*Rv - cor_vals = np.power(10.0,-0.4*(cor_vals*Av)) + cor_vals = np.power(10.0, -0.4*(cor_vals*Av)) # initialize extinction model tmodel = O94(Rv=Rv) diff --git a/dust_extinction/tests/test_p92.py b/dust_extinction/tests/test_p92.py index 8e4478e..ec4a9a2 100644 --- a/dust_extinction/tests/test_p92.py +++ b/dust_extinction/tests/test_p92.py @@ -2,24 +2,30 @@ import pytest import astropy.units as u -from astropy.modeling import InputParameterError from astropy.modeling.fitting import LevMarLSQFitter from ..dust_extinction import P92 x_bad = [-1.0, 1001.] -@pytest.mark.parametrize("x_invalid", x_bad) -def test_invalid_wavenumbers(x_invalid): - tmodel = P92() + + +def _invalid_x_range(x, tmodel, modname): with pytest.raises(ValueError) as exc: - tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for P92' \ + tmodel(x) + assert exc.value.args[0] == 'Input x outside of range defined for ' \ + + modname \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + +@pytest.mark.parametrize("x_invalid", x_bad) +def test_invalid_wavenumbers(x_invalid): + _invalid_x_range(x_invalid, P92(), 'P92') + + @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): tmodel = P92() @@ -28,10 +34,11 @@ def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): assert exc.value.args[0] == 'Input x outside of range defined for P92' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): tmodel = P92() @@ -40,10 +47,11 @@ def test_invalid_micron(x_invalid_micron): assert exc.value.args[0] == 'Input x outside of range defined for P92' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): tmodel = P92() @@ -52,10 +60,11 @@ def test_invalid_micron(x_invalid_angstrom): assert exc.value.args[0] == 'Input x outside of range defined for P92' \ + ' [' \ + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ + + ' <= x <= ' \ + str(tmodel.x_range[1]) \ + ', x has units 1/micron]' + def get_axav_cor_vals(): # Milky Way observed extinction as tabulated by Pei (1992) @@ -80,6 +89,7 @@ def get_axav_cor_vals(): return (x, cor_vals) + def test_extinction_P92_values(): # get the correct values x, cor_vals = get_axav_cor_vals() From 4c59c3b8dacf216a19b1128cf82a4ee9b0a24464 Mon Sep 17 00:00:00 2001 From: Karl Gordon Date: Sat, 13 Jan 2018 09:10:21 -0500 Subject: [PATCH 5/6] Update tests to use helper function instead of 4 repeats in each file --- dust_extinction/tests/helpers.py | 13 +++++++ dust_extinction/tests/test_ccm89.py | 41 +++----------------- dust_extinction/tests/test_f99.py | 51 +++---------------------- dust_extinction/tests/test_fm90.py | 44 ++++----------------- dust_extinction/tests/test_g03.py | 41 ++++---------------- dust_extinction/tests/test_g16.py | 41 +++----------------- dust_extinction/tests/test_gcc09_ave.py | 42 +++----------------- dust_extinction/tests/test_o94.py | 41 +++----------------- dust_extinction/tests/test_p92.py | 43 ++------------------- 9 files changed, 57 insertions(+), 300 deletions(-) create mode 100644 dust_extinction/tests/helpers.py diff --git a/dust_extinction/tests/helpers.py b/dust_extinction/tests/helpers.py new file mode 100644 index 0000000..0483d1e --- /dev/null +++ b/dust_extinction/tests/helpers.py @@ -0,0 +1,13 @@ +import pytest + + +def _invalid_x_range(x, tmodel, modname): + with pytest.raises(ValueError) as exc: + tmodel(x) + assert exc.value.args[0] == 'Input x outside of range defined for ' \ + + modname \ + + ' [' \ + + str(tmodel.x_range[0]) \ + + ' <= x <= ' \ + + str(tmodel.x_range[1]) \ + + ', x has units 1/micron]' diff --git a/dust_extinction/tests/test_ccm89.py b/dust_extinction/tests/test_ccm89.py index 30bdd5e..621ea24 100644 --- a/dust_extinction/tests/test_ccm89.py +++ b/dust_extinction/tests/test_ccm89.py @@ -5,6 +5,7 @@ from astropy.modeling import InputParameterError from ..dust_extinction import CCM89 +from .helpers import _invalid_x_range @pytest.mark.parametrize("Rv_invalid", [-1.0, 0.0, 1.9, 6.1, 10.]) @@ -16,57 +17,25 @@ def test_invalid_Rv_input(Rv_invalid): @pytest.mark.parametrize("x_invalid", [-1.0, 0.2, 10.1, 100.]) def test_invalid_wavenumbers(x_invalid): - tmodel = CCM89(Rv=3.1) - with pytest.raises(ValueError) as exc: - tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for CCM89' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid, CCM89(Rv=3.1), 'CCM89') @pytest.mark.parametrize("x_invalid_wavenumber", [-1.0, 0.2, 10.1, 100.]/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): - tmodel = CCM89(Rv=3.1) - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for CCM89' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_wavenumber, CCM89(Rv=3.1), 'CCM89') @pytest.mark.parametrize("x_invalid_micron", u.micron/[-1.0, 0.2, 10.1, 100.]) def test_invalid_micron(x_invalid_micron): - tmodel = CCM89(Rv=3.1) - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for CCM89' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_micron, CCM89(Rv=3.1), 'CCM89') @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/[-1.0, 0.2, 10.1, 100.]) def test_invalid_micron(x_invalid_angstrom): - tmodel = CCM89(Rv=3.1) - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for CCM89' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_angstrom, CCM89(Rv=3.1), 'CCM89') def test_axav_ccm89_table3(): diff --git a/dust_extinction/tests/test_f99.py b/dust_extinction/tests/test_f99.py index e545cd3..c2298f4 100644 --- a/dust_extinction/tests/test_f99.py +++ b/dust_extinction/tests/test_f99.py @@ -4,6 +4,7 @@ import astropy.units as u from ..dust_extinction import F99 +from .helpers import _invalid_x_range x_bad = [-1.0, 0.1, 12.0, 100.] @@ -11,67 +12,25 @@ @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): - tmodel = F99() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for F99' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid, F99(), 'F99') @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): - tmodel = F99() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for F99' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_wavenumber, F99(), 'F99') @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): - tmodel = F99() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for F99' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_micron, F99(), 'F99') @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): - tmodel = F99() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for F99' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_angstrom, F99(), 'F99') def get_axav_cor_vals(): - # testing wavenumbers - # from previous version of code - # x = np.array([0.5, 1.0, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, - # 8.0]) - - # correct values (not quite right) - # cor_vals = np.array([0.124997, 0.377073, 1.130636, 1.419644, - # 1.630377, 1.888546, 2.275900, 3.014577, - # 2.762256, 2.475272, 2.711508, 3.197144]) - # from Fitzpatrick (1999) Table 3 x = np.array([0.377, 0.820, 1.667, 1.828, 2.141, 2.433, 3.704, 3.846]) diff --git a/dust_extinction/tests/test_fm90.py b/dust_extinction/tests/test_fm90.py index ba1dfa4..acb4506 100644 --- a/dust_extinction/tests/test_fm90.py +++ b/dust_extinction/tests/test_fm90.py @@ -4,61 +4,31 @@ import astropy.units as u from astropy.modeling.fitting import LevMarLSQFitter -from ..dust_extinction import FM90, G03_LMCAvg +from ..dust_extinction import (FM90, G03_LMCAvg) +from .helpers import _invalid_x_range + x_bad = [-1.0, 0.2, 3.0, 11.0, 100.] @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): - tmodel = FM90() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for FM90' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid, FM90(), 'FM90') @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): - tmodel = FM90() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for FM90' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_wavenumber, FM90(), 'FM90') @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): - tmodel = FM90() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for FM90' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_micron, FM90(), 'FM90') @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): - tmodel = FM90() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for FM90' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_angstrom, FM90(), 'FM90') def get_elvebv_cor_vals(): diff --git a/dust_extinction/tests/test_g03.py b/dust_extinction/tests/test_g03.py index 42bf070..d20f71a 100644 --- a/dust_extinction/tests/test_g03.py +++ b/dust_extinction/tests/test_g03.py @@ -4,62 +4,35 @@ import astropy.units as u from astropy.modeling import InputParameterError -from ..dust_extinction import G03_SMCBar, G03_LMCAvg, G03_LMC2 +from ..dust_extinction import (G03_SMCBar, G03_LMCAvg, G03_LMC2) +from .helpers import _invalid_x_range x_bad = [-1.0, 0.1, 10.1, 100.] models = [G03_SMCBar(), G03_LMCAvg(), G03_LMC2()] + @pytest.mark.parametrize("x_invalid", x_bad) @pytest.mark.parametrize("tmodel", models) def test_invalid_wavenumbers(x_invalid, tmodel): - - with pytest.raises(ValueError) as exc: - tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for G03' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid, tmodel, 'G03') @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) @pytest.mark.parametrize("tmodel", models) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber, tmodel): - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for G03' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_wavenumber, tmodel, 'G03') @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) @pytest.mark.parametrize("tmodel", models) def test_invalid_micron(x_invalid_micron, tmodel): - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for G03' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_micron, tmodel, 'G03') @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) @pytest.mark.parametrize("tmodel", models) def test_invalid_micron(x_invalid_angstrom, tmodel): - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for G03' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_angstrom, tmodel, 'G03') @pytest.mark.parametrize("tmodel", models) diff --git a/dust_extinction/tests/test_g16.py b/dust_extinction/tests/test_g16.py index bdd2310..53c8cc6 100644 --- a/dust_extinction/tests/test_g16.py +++ b/dust_extinction/tests/test_g16.py @@ -7,6 +7,7 @@ from ..dust_extinction import G16 from ..dust_extinction import G03_SMCBar from .test_f99 import get_axav_cor_vals as get_axav_cor_vals_fA_1 +from .helpers import _invalid_x_range @pytest.mark.parametrize("RvA_invalid", [-1.0, 0.0, 1.9, 6.1, 10.]) @@ -28,54 +29,22 @@ def test_invalid_fA_input(fA_invalid): @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): - tmodel = G16() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for G16' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid, G16(), 'G16') @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): - tmodel = G16() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for G16' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_wavenumber, G16(), 'G16') @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): - tmodel = G16() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for G16' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_micron, G16(), 'G16') @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): - tmodel = G16() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for G16' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_angstrom, G16(), 'G16') def test_extinction_G16_fA_1_values(): diff --git a/dust_extinction/tests/test_gcc09_ave.py b/dust_extinction/tests/test_gcc09_ave.py index 5bd53b2..2a2b7c5 100644 --- a/dust_extinction/tests/test_gcc09_ave.py +++ b/dust_extinction/tests/test_gcc09_ave.py @@ -5,60 +5,30 @@ from astropy.modeling import InputParameterError from ..dust_extinction import GCC09_MWAvg +from .helpers import _invalid_x_range + x_bad = [-1.0, 0.1, 11., 100.] @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): - tmodel = GCC09_MWAvg() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for GCC09' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid, GCC09_MWAvg(), 'GCC09') @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): - tmodel = GCC09_MWAvg() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for GCC09' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_wavenumber, GCC09_MWAvg(), 'GCC09') @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): - tmodel = GCC09_MWAvg() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for GCC09' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_micron, GCC09_MWAvg(), 'GCC09') @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): - tmodel = GCC09_MWAvg() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for GCC09' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_angstrom, GCC09_MWAvg(), 'GCC09') def test_extinguish_no_av_or_ebv(): diff --git a/dust_extinction/tests/test_o94.py b/dust_extinction/tests/test_o94.py index b8c7108..1a8b5fe 100644 --- a/dust_extinction/tests/test_o94.py +++ b/dust_extinction/tests/test_o94.py @@ -5,6 +5,7 @@ from astropy.modeling import InputParameterError from ..dust_extinction import O94 +from .helpers import _invalid_x_range @pytest.mark.parametrize("Rv_invalid", [-1.0, 0.0, 1.9, 6.1, 10.]) @@ -16,57 +17,25 @@ def test_invalid_Rv_input(Rv_invalid): @pytest.mark.parametrize("x_invalid", [-1.0, 0.2, 10.1, 100.]) def test_invalid_wavenumbers(x_invalid): - tmodel = O94(Rv=3.1) - with pytest.raises(ValueError) as exc: - tmodel(x_invalid) - assert exc.value.args[0] == 'Input x outside of range defined for O94' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid, O94(Rv=3.1), 'O94') @pytest.mark.parametrize("x_invalid_wavenumber", [-1.0, 0.2, 10.1, 100.]/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): - tmodel = O94(Rv=3.1) - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for O94' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_wavenumber, O94(Rv=3.1), 'O94') @pytest.mark.parametrize("x_invalid_micron", u.micron/[-1.0, 0.2, 10.1, 100.]) def test_invalid_micron(x_invalid_micron): - tmodel = O94(Rv=3.1) - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for O94' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_micron, O94(Rv=3.1), 'O94') @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/[-1.0, 0.2, 10.1, 100.]) def test_invalid_micron(x_invalid_angstrom): - tmodel = O94(Rv=3.1) - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for O94' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_angstrom, O94(Rv=3.1), 'O94') def test_axav_o94_rv31(): diff --git a/dust_extinction/tests/test_p92.py b/dust_extinction/tests/test_p92.py index ec4a9a2..a7c46b4 100644 --- a/dust_extinction/tests/test_p92.py +++ b/dust_extinction/tests/test_p92.py @@ -5,22 +5,11 @@ from astropy.modeling.fitting import LevMarLSQFitter from ..dust_extinction import P92 +from .helpers import _invalid_x_range x_bad = [-1.0, 1001.] -def _invalid_x_range(x, tmodel, modname): - with pytest.raises(ValueError) as exc: - tmodel(x) - assert exc.value.args[0] == 'Input x outside of range defined for ' \ - + modname \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' - - @pytest.mark.parametrize("x_invalid", x_bad) def test_invalid_wavenumbers(x_invalid): _invalid_x_range(x_invalid, P92(), 'P92') @@ -28,41 +17,17 @@ def test_invalid_wavenumbers(x_invalid): @pytest.mark.parametrize("x_invalid_wavenumber", x_bad/u.micron) def test_invalid_wavenumbers_imicron(x_invalid_wavenumber): - tmodel = P92() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_wavenumber) - assert exc.value.args[0] == 'Input x outside of range defined for P92' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_wavenumber, P92(), 'P92') @pytest.mark.parametrize("x_invalid_micron", u.micron/x_bad) def test_invalid_micron(x_invalid_micron): - tmodel = P92() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_micron) - assert exc.value.args[0] == 'Input x outside of range defined for P92' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_micron, P92(), 'P92') @pytest.mark.parametrize("x_invalid_angstrom", u.angstrom*1e4/x_bad) def test_invalid_micron(x_invalid_angstrom): - tmodel = P92() - with pytest.raises(ValueError) as exc: - tmodel(x_invalid_angstrom) - assert exc.value.args[0] == 'Input x outside of range defined for P92' \ - + ' [' \ - + str(tmodel.x_range[0]) \ - + ' <= x <= ' \ - + str(tmodel.x_range[1]) \ - + ', x has units 1/micron]' + _invalid_x_range(x_invalid_angstrom, P92(), 'P92') def get_axav_cor_vals(): From eec1bcbbd2aa2f4140b7d2b453dfb55a0c49feef Mon Sep 17 00:00:00 2001 From: Karl Gordon Date: Sat, 13 Jan 2018 09:14:40 -0500 Subject: [PATCH 6/6] Change alavtoelv to axavtoexv to be consist with rest of code --- dust_extinction/dust_extinction.py | 22 +++++++++---------- .../{test_alavtoelv.py => test_axavtoexv.py} | 8 +++---- 2 files changed, 15 insertions(+), 15 deletions(-) rename dust_extinction/tests/{test_alavtoelv.py => test_axavtoexv.py} (91%) diff --git a/dust_extinction/dust_extinction.py b/dust_extinction/dust_extinction.py index f435135..a064319 100644 --- a/dust_extinction/dust_extinction.py +++ b/dust_extinction/dust_extinction.py @@ -11,7 +11,7 @@ __all__ = ['CCM89', 'FM90', 'P92', 'O94', 'F99', 'G03_SMCBar', 'G03_LMCAvg', 'G03_LMC2', 'GCC09_MWAvg', 'G16', - 'AlAvToElv'] + 'AxAvToExv'] x_range_CCM89 = [0.3, 10.0] x_range_FM90 = [1.0/0.32, 1.0/0.0912] @@ -347,37 +347,37 @@ def fA(self, value): + str(self.fA_range[1])) -class AlAvToElv(Fittable1DModel): +class AxAvToExv(Fittable1DModel): """ - Model to convert from A(lambda)/A(V) to E(lambda-V) + Model to convert from A(x)/A(V) to E(x-V) Paramters --------- Av : float dust column in A(V) [mag] """ - inputs = ('alav',) - outputs = ('elv',) + inputs = ('axav',) + outputs = ('exv',) Av = Parameter(description="A(V)", default=1.0, min=0.0) @staticmethod - def evaluate(alav, Av): + def evaluate(axav, Av): """ AlAvToElv function Paramters --------- - alav : np array (float) - E(lambda-V)/E(B-V) values + axav : np array (float) + E(x-V)/E(B-V) values Returns ------- - elv : np array (float) - E(lamda - V) + exv : np array (float) + E(x - V) """ - return (alav - 1.0)*Av + return (axav - 1.0)*Av # use numerical derivaties (need to add analytic) fit_deriv = None diff --git a/dust_extinction/tests/test_alavtoelv.py b/dust_extinction/tests/test_axavtoexv.py similarity index 91% rename from dust_extinction/tests/test_alavtoelv.py rename to dust_extinction/tests/test_axavtoexv.py index 88a548a..1b52f0e 100644 --- a/dust_extinction/tests/test_alavtoelv.py +++ b/dust_extinction/tests/test_axavtoexv.py @@ -3,7 +3,7 @@ from astropy.modeling.fitting import LevMarLSQFitter -from ..dust_extinction import (P92, AlAvToElv) +from ..dust_extinction import (P92, AxAvToExv) def get_axav_cor_vals(): @@ -25,9 +25,9 @@ def get_axav_cor_vals(): return (MW_x/u.micron, MW_axav) -class P92_Elv(P92 | AlAvToElv): +class P92_Exv(P92 | AxAvToExv): """ - Evalute P92 on E(l-V) data including solving for A(V) + Evalute P92 on E(x-V) data including solving for A(V) """ @@ -44,7 +44,7 @@ def test_P92_fitting(): y = (y_axav - 1.0)*av # change from defaults to make the best fit harder to find - p92_init = P92_Elv() + p92_init = P92_Exv() fit = LevMarLSQFitter() # accuracy set to avoid warning the fit may have failed