Analytic functions with blackbody

See #3077 
Ah, I think I understand why one could get `NaN`: if the wavelength is `0`, so `flux==0` and then one converts to another unit which requires dviding by the wavelength, so one has `0/0`. So, do keep! But perhaps be explicit in what you are resetting:
```
# Catch NaN (from division of 0 by 0)
flux[np.isnan(flux)] = 0.
```