astropy/specutils — pull request #55 by kbarbary (2014-02-01 23:14:57+00:00)
Title: Extinction cleanup and performance

This PR does three things:

**1. Separate functions for the various extinction models**

For example, rather than `extinction(wave, a_v=1., model='ccm89')`, you do `extinction_ccm89(wave, a_v=1)`. The old interface is still available but simply wraps the individual functions.
See the built docs here: http://www.hep.anl.gov/kbarbary/specutils-docs/specutils/extinction.html

**2. Optional class-based interface for some models with initial overheads**

This PR does away with "cacheing" used in the 'wd01' and 'd03' models, in favor of an optional class-based interface. For these two models, and to a lesser extent, 'f99', it takes a long time to initialize the function for a given `R_V`, but not very long to actually evalute it. So, if a user plans to use an extinction function with a fixed `R_V` multiple times, they can use the class-based interface:
```python
f = ExtinctionWD01(r_v=3.1) # slow as hell
f(3000., a_v=1.) # quite fast, can do multiple times.
```
The corresponding functional interface (e.g., `extinction_wd01`) simply wraps this class (create an instance and immediately use it) and does not cache anything.

In addition to making the behavior of the `wd01` and `d03` functions clearer, it also provides up to a 6x speed-up for F99:

*'f99' master*
```python
In [6]: for n in [10, 100, 1000, 10000]:
print ' n=', n,':',
wave = np.logspace(3., np.log10(30000.), n)
%timeit extinction(wave, a_v=1., model='f99')
...:
n= 10 : 1000 loops, best of 3: 664 µs per loop
n= 100 : 1000 loops, best of 3: 682 µs per loop
n= 1000 : 1000 loops, best of 3: 792 µs per loop
n= 10000 : 1000 loops, best of 3: 1.87 ms per loop
```

*'f99' this branch*
```python
In [7]: f99 = ExtinctionF99(3.1)

In [8]: for n in [10, 100, 1000, 10000]:
print ' n=', n,':',
wave = np.logspace(3., np.log10(30000.), n)
%timeit f99(wave, a_v=1.)
...:
n= 10 : 10000 loops, best of 3: 115 µs per loop
n= 100 : 10000 loops, best of 3: 125 µs per loop
n= 1000 : 1000 loops, best of 3: 212 µs per loop
n= 10000 : 1000 loops, best of 3: 1.1 ms per loop
```

**3. Cythonizes the 'ccm89', 'od94' and 'gcc09' models**

Provides a factor of ~10 speed up for small arrays. It also makes the code easier to read, in my opinion. Here are some timing comparisons for different size arrays (n=10, 100, 1000, 10000):

*'ccm89' master* (od94 and gcc09 are similar)
```python
In [3]: for n in [10, 100, 1000, 10000]:
print ' n=', n,':',
wave = np.logspace(3., np.log10(30000.), n)
%timeit extinction(wave, a_v=1., model='ccm89')
...:
n= 10 : 1000 loops, best of 3: 357 µs per loop
n= 100 : 1000 loops, best of 3: 381 µs per loop
n= 1000 : 1000 loops, best of 3: 548 µs per loop
n= 10000 : 1000 loops, best of 3: 2 ms per loop
```

*This branch*
```python
In [4]: for n in [10, 100, 1000, 10000]:
print ' n=', n,':',
wave = np.logspace(3., np.log10(30000.), n)
%timeit extinction_ccm89(wave, a_v=1.)
...:
n= 10 : 10000 loops, best of 3: 32.8 µs per loop
n= 100 : 10000 loops, best of 3: 41.8 µs per loop
n= 1000 : 10000 loops, best of 3: 126 µs per loop
n= 10000 : 1000 loops, best of 3: 943 µs per loop
```

This isn't quite ready to merge yet, still need to finish updating docstrings and check tests and accuracy.
Just wanted to share progress with you all!