Automatic Normalization of Spectra and Time Series of Spectra in fits-format

Normally the 1d spectra are available in fits format. These can be single spectra as well as whole time series of spectra of the same object, which shall now be subjected to a more precise evaluation (e.g. the determination of radial velocities or the equivalent width of lines).

For this purpose it is usually useful to normalize the spectra to the continuum. The Python script described in the following is a convenient way to do this. Download the script here .

For a large number of spectra to be normalized with many pixels (> 5000), sufficient computing power is required, which is why the script uses almost all of the processors (CPU’s) in your computer simultaneously (multiprocessing technique).

The program requires some Python libraries, as can be seen from the following excerpt of the script:

import numpy as np
import matplotlib.pyplot as plt
from astropy.io import fits
from scipy.stats import sigmaclip
import glob
from concurrent import futures
import os
from numba import jit
from csaps import csaps

You probably already have the first 8 installed (they are included by default in the anaconda distribution of Python 3), but You have to install the library „csasp“ separately with the command

pip install –user csaps

To start the program first create a console, change to the working directory where the spectra are located and then start the script with the command

python3 /HOME/automaticNormalization_timeseries_en.py

where /HOME/ is the path to the directory where you saved the script.

Or copy the script into your working directory where the spectra are stored and start with :

python3 automaticNormalization_timeseries_en.py

Alternatively, you can import the script into spyder and run it.

The calculations within the script run as follows:

First you will be asked for the file names of the 1d spectra, where you should use wildcards, e.g. *.fit.
Then these files are read in and a list of the spectra is printed out in the console.
The script divides an interval of 10 pixels from pixel to pixel, calculates a percentile of it and compares it with the percentiles of the neighboring intervals to find local maxima, which are counted as primary interpolation points.
The first spectrum is shown in a diagram (blue line) and the primary interpolation points are shown as red points.
Then you are asked if you want to change the interval width. With this you can change the number of primary interpolation points. Increasing the interval width also causes less primary interpolation points to be found far below the continuum in blends.
If you type „y“ (yes), you can then enter the new value for the interval width.
You can do this as often as you like. Each time the diagram is renewed and you can observe how the calculated primary interpolation points change. If you are satisfied with the result, enter anything but don’t type „y“ on the question (simply just press return).

Then the question appears in the console: Do you want to remove interpolation points in certain wavelength ranges? Then enter y:. In this example we say y, because we want to prevent interpolation points in the wavelength ranges 6510 to 6520 and 6546 to 6565 Angstrom from being considered. We enter these wavelength ranges and get the following graph. Now the interpolation points in the two wavelength ranges 6510 to 6520 and 6546 to 6565 Angstrom are removed and a spline was calculated to represent the continuum (black line).

We can influence the curvature of the spline by entering a smoothing factor. The smoothing factor can be between 0 and 1.0 simply means a straight line (regression straight line, lowest detail) through the points, 1 means that the spline goes through every point (highest detail). In the above graph you can see that the spline in the range 6540 to 6570 Angstrom is yielding downwards, it is set too soft. We want to make it stiffer, so we enter a lower smoothing factor as used to calculate the spline above. We answer the question ‚Do you want to change the smoothing factor?‘ with y and correct in this case the smoothing factor from 0.1 to 0.01. We get the new diagram with stiffer spline.

Now the course of the spline already shows the continuum quite well. But there are still points, which are clearly below the spline and influence it. We can correct this by answering the next question, If you want to delete interpolation points s1 % below and s2 % above the spline, enter y:, with ‚y‘ and then enter the two parameters s1 and s2. If we take for s1 0.5 and s2 1, then the points below and above these percent limits are deleted. The following diagram is shown in which the points and the spline has been corrected again. You can repeat this process as often as you like until the spline envelops the spectrum properly.

Now we are satisfied with the course of the spline, which now reflects the continuum shape quite well. We don’t change anything now, answer the question if the graph and the normalized spectrum and the normfunction should be saved and then answer the question if all spectra of the time series should be normalized with ‚y‘. Then all spectra are normalized with the set parameters, saved in the working directory as .fit and, if the question was answered with ‚y‘, the graphs are saved as PDF, too. The program is finished.

In this example a time series of 6 spectra of the H alpha line of the star eps Aur was normalized. The results of the normalization’s are shown below.