Hello I extracted a spectrum from the central region of the 2D spectrum of the major axis from NGC 4697 (I just took a central row in the 2D data set, see this question for more information: How to extract galaxy spectra for different radii in Python for spectra taken by long slit spectrograph?). I also have got a template spectrum from the K2-III star hd132345. Here both spectras are plotted: One can clearly see that the absorption lines in the galaxy spectrum have been broadened, due to the superposition of many stellar spectra, each one having a individual velocity along the line of sight. I now have done some preprocessing steps to keep on working with the spectra in fourier space later. These steps included: cutting off pixels at both ends of the spectrum, removing the continuum, tapering both spectra with a window function using np.blackman(). The result can be seen here: To further work with the spectra I want to remove/measure the redshift from the galaxy, to use them later in the FCQ algorithm to determine the velocity distribution along the line of sight and therefore the kinematics of NGC 4697. One can clearly see the 3 dominant lines in the template spectrum arround the x-value of 8.55, these lines should be the Magnesium triplet, which has its wavelengths at round about 5200 Angstroems ( note that the spectra are logarithmically rebinned). I wanted to use these lines to determine the redshift of the galaxy. The triplet in the galaxy has involved into a doublet, due to the broadening. I did set up vertical lines at the peaks of the absorption lines and calculated the shift graphically. One can clearly see that the galaxy absorption lines have been shifted towards the red end of the spectrum. I wanted to know if there is a more elegant/common way to tackle this problem, since it is probably one of the most basic things to do in Astronomy. I just could not find a sophisticated answer for a Python solution in the internet and did not want to use a prewritten funciton, for educational purposes. My goal is it to not have to look for the x values of the peaks manually, to make the problem easier when you are working with more spectra. The spectra in the MG-Triplet region with redshift look as follows. And after I subtracted the redshift. . Here I just calculated: loglam_gal = loglam_gal - (8.55702-8.55323). Is there a better way to do this and did I even do it correct ? Both spectra after removing the redshift from the galaxy look as follows: I have calculated the redshift from my graphically measured wavelengths to be: z = 0.003797. In Wikipedia it says z= 0.004140 for NGC 4697.
There are two ways (at least).
Do a pixel by pixel cross-correlation of the galaxy and template spectrum. The cross-correlation function will have a broad peak, and the position of the peak is the "lag" that gives you the radial velocity difference between the template and the galaxy. The lag will be in pixels, but since you have (correctly) got your data in equal bins of log wavelength, then each pixel will correspond to a fixed radial velocity increment (see https://astronomy.stackexchange.com/a/14524/2531 ). You appear in any case to be heading towards using the Tonry & Davis FFT cross-correlation route (multiplying the Fourier Transforms), which would give a similar cross-correlation function.
Direct fitting of the spectrum. Broaden and shift the template spectrum until it gives you the best match with the galaxy spectrum. The free parameters of your model will be the amount of broadening, possibly a relative normalisation and a (velocity) shift. The latter is the parameter you are looking for, so use its best-fit value. There are lots of ways to achieve the fitting in Python, from classic least squares fitting to more sophisticated Bayesian and Monte Carlo Markov Chain methods.