Python Tip: Re-sampling spectra with pysynphot

This post was inspired by a question from John Johnson.

Have you ever wanted to plot a model spectrum at lower resolution? Or compare a model spectrum with an observed spectrum? Have you ever wanted to shift several observed spectra to a common redshift to stack them up? In all these cases, you would need to resample the spectrum. Just such a scenario arose and I was asked:

What is the recommended way to resample a 1D spectrum?

I recalled using pysynphot, which is very useful for working with 1D spectra. However, the documentation is still a bit thin. After a bit of of exploring, here is an example that John Johnson and I wrote to illustrate our answer to the question. This package takes care to resample while preserving flux. We hope you find this useful! Or leave your alternative answer in the comments below.

from pysynphot import observation
from pysynphot import spectrum
 
def rebin_spec(wave, specin, wavnew):
    spec = spectrum.ArraySourceSpectrum(wave=wave, flux=specin)
    f = np.ones(len(wave))
    filt = spectrum.ArraySpectralElement(wave, f, waveunits='angstrom')
    obs = observation.Observation(spec, filt, binset=wavnew, force='taper')
 
    return obs.binflux

For those bleeding-edge python users out there, pysynphot is becoming an astropy affiliated package. You can download this new version from the pysynphot GitHub repository and test it out.

12 comments… add one
  • Christian Herenz Aug 13, 2013 @ 13:40

    I don’t get why I would want to do it this complicated, not using some simple interpolation + scaling just using numpy packages? In the above example I would have no idea in the end what comes out….

    If I know what units my input spectrum is both resampling and redshifting are straightforward mathematical expressions which can be entered as “pure” Numpy Python code.

  • Cristobal Sifon Aug 15, 2013 @ 15:14

    Sorry for the self publicity, but for an example of what Christian Herenz said, check out
    http://home.strw.leidenuniv.nl/~sifon/plotspec/
    the program is given as stand alone but it’s python code which you can dig into (comments included).

  • Jessica Lu Aug 15, 2013 @ 19:08

    Points in a spectrum are integrated fluxes over small wavelength bins, rather than flux densities. Pysynphot nicely preserves the integrated flux, while not all standard interpolation schemes do this.

  • Pey Lian Lim Apr 22, 2014 @ 12:09

    Hi. I am the current developer of the “pysynphot GitHub repository” (https://github.com/spacetelescope/pysynphot) mentioned in your blog post. It is currently still in development and not ready for general usage. However, once it is ready, it will utilize the powerful modeling package in astropy (http://astropy.org), which will greatly expand the functionality of the official release version. Until that happens, please use the pysynphot officially distributed by STScI at http://www.stsci.edu/institute/software_hardware/pyraf/stsci_python . Thank you.

    • Derek Homeier Sep 23, 2014 @ 12:29

      Hi Pey,
      is there any way to use or even just test pysynphot at this point? Building the git source fails immediately on attempts to
      from jwst_lib import modeling
      which I cannot find anywhere in stsci_python either (or any “modeling” module whatsoever). It looks like it is intended to be replaced by astropy.modeling eventually, but until that works, what jwst_lib is supposed to be used here? Thanks, Derek

  • Joe May 1, 2014 @ 16:34

    This works well, but how would you resample an associated uncertainty array? You would want to add the uncertainties in each bin in quadrature but how could you include that in this snippet using pysynphot? Thanks!

  • Pey Lian Lim May 2, 2014 @ 12:51

    Hi, Joe. Currently, pysynphot does not handle uncertainties. It might in the future, resource permitting.

  • Mike May 13, 2014 @ 7:42

    Hello

    Thank you for your nice work !

    I would like to decrease the resolution of a template-model spectrum down to the resolution of my observed spectra. Could you explain me the procedure through pysynphot or any alternative ?

    Thank you in advance !

    • Cristóbal Sifón May 13, 2014 @ 14:26

      Hi Mike,

      the code I mentioned above (http://home.strw.leidenuniv.nl/~sifon/plotspec/) can do this. It simply creates a new spectrum by averaging the signal every N pixels (where you supply N via the “smooth” parameter). This should be enough for most purposes, but I suspect you can get a more optimal smoothing with scipy.ndimage.filters.gaussian_filter and the likes.

    • Pey Lian Lim May 27, 2014 @ 10:49

      Hi, Mike. If you still wish to use pysynphot to resample your spectrum, please email all the details to help[at]stsci.edu and someone will get back to you shortly. Thank you for your interest.

  • Tugca Aug 10, 2016 @ 23:34

    Hi,
    I am trying to do the same as Mike but I couldn’t find my way…
    I wonder if there is any progress/update on the package in the last 3 years.

  • A. C. Carnall May 16, 2017 @ 5:20

    For spectral resampling, a new piece of code I’ve just released, SpectRes, may be helpful (https://github.com/ACCarnall/SpectRes). It’s a simple function which deals with uncertainties and preserves integrated flux.

Leave a Reply

Your email address will not be published. Required fields are marked *