The Richardson-Lucy Deconvolution method to Extract LAMOST 1D Spectra

TL;DR

This paper applies the Richardson-Lucy deconvolution algorithm to extract high-quality 1D spectra from LAMOST images, achieving better resolution and noise reduction compared to existing methods.

Contribution

It introduces a faster deconvolution approach that improves spectral resolution and signal-to-noise ratio while reducing ringing artifacts in LAMOST spectra.

Findings

Higher SNR and resolution in extracted spectra

Faster processing compared to other algorithms

Effective suppression of ringing artifacts

Abstract

We use the Richardson-Lucy deconvolution algorithm to extract one dimensional (1D) spectra from LAMOST spectrum images. Compared with other deconvolution algorithms, this algorithm is much more fast. The practice on a real LAMOST image illustrates that the 1D resulting spectrum of this method has a higher SNR and resolution than those extracted by the LAMOST pipeline. Furthermore, our algorithm can effectively depress the ringings that are often shown in the 1D resulting spectra of other deconvolution methods.

Tables

Table 1. Table 1: SNRs and Computation Times of Different Methods

Methods	SNR	Time (Seconds)
TDM	46.63	11.94
LIM	47.08	6.05
Equation  7 for 10 iterations	50.03	1.64

Equations

$$f(x,y)=\sum _{j=1}^N{g}_j{h}_j(x-j,y)+\eta (x,y)$$

$$a(x,y)=\sum _{j=1}^N{g}_j{h}_j(x-j,y)$$

$$P(f(x,y)|G)=\prod _{x,y}\frac{a(x,y{)}^{f(x,y)}{e}^{-a(x,y)}}{f(x,y)!}$$

$$\ln P(f(x,y)|G)=\sum _{x,y}[f(x,y)\ln a(x,y)-a(x,y)-\ln (f(x,y)!)]$$

$$\sum _{x,y}[\frac{f(x,y)h_j(x-j,y)}{a(x,y)}-h_j(x-j,y)]=0$$

$$\sum _{x,y}[\frac{f(x,y)h_j(x-j,y)}{a(x,y)}]-1=0$$

$$g_{j}^{(k+1)}=g_{j}^{(k)}\Big{\{}\sum_{x,y}\frac{f(x,y)h_{j}(x-j,y)}{a(x,y)}\Big{\}}^{p},j=1,2,3,...,N$$

$$g_{j}^{(k+1)}=g_{j}^{(k)}\Big{\{}\sum_{x,y}\frac{f(x,y)h_{j}(x-j,y)}{a(x,y)}\Big{\}},j=1,2,3,...,N$$

Full text

1 Introduction

Guoshoujing Telescope (also called Large Sky Area Multi-Object Fiber Spectroscopic Telescope, LAMOST)(Wang et al., 1996; Su & Cui, 2004; Zhao et al., 2012; Cui et al., 2012; Luo et al., 2012) can obtains 4000 spectra in one exposure, which has collected more than 10 million spectra(see http://dr7.lamost.org). So far, about 400 peer reviewed papers based on LAMOST data are published, which help us more deeply understand our Galaxy.

The spectrum extraction method is a technique to convert a two-dimensional (2D) CCD image into one-dimensional (1D) spectra, which can help astronomers explore the natures of celestial objects. The traditional extraction methods, including Aperture Extraction Method (AEM, hereafter), the Optimal Extraction Method and the Profile Fitting Method, are the most frequently used methods at present, which are discussed in detail in (Li et al., 2019). Compared to the traditional methods, the deconvolution method is a completely different method, which tries to recover the 1D spectra by eliminating the instrumental profiles (Point Spread Function, PSF) on a 2D image. The deconvolution method is first presented by Bolton & Schlegel (2010). In their work, a calibration matrix is constructed from known instrumental profiles, and then inverted to calculate the resulting 1D spectra. As a result, the method can only extract a small piece of spectrum image, and cannot be used to extract a large spectrum image (for example, a LAMOST $4K\times 4K$ spectrum image) because of huge storage for the calibration matrix, let along calculation. Furthermore, the noise in a real image would destroy the resulting 1D spectra.

To overcome the above shortcomings, Guangwei et al. (2015) presented a deconvolution algorithm based on the Tikhonov regularization method (TDM, hearafter) for practical spectrum extraction. Firstly, they gave a method to obtain all PSFs which vary with positions on a CCD image. Secondly, the big spectrum image is divided into many small subimages, whose calibration matrixes can be easily stored and inverted during calculation. Thirdly, the Tikhonov regularization item can effectively suppress the noise. This algorithm is the first practical deconvolution method to extract 1D spectra from a multi-fiber spectrum image. The signal-to-noise ratio (SNR) and resolution of the resulting 1D spectra are both higher than traditional methods. However, they did not discuss how to choose the best Tikhonov parameter. The choice of the parameter is a crucial problem, because a fixed Tikhonov parameter is not always appropriate for all fibers in a multi-fiber spectrum image.

In the last paper (Li et al., 2019), we developed a deconvolution method based on adaptive Landweber iteration (ALI, hearafter), which can extract 1D spectra with the regularization parameter adaptively selected for every fiber. The SNR and resolution of 1D resulting spectra are both as high as those of the 1D resulting spectra extracted by TDM with a deliberately selected Tikhonov regularization parameter.

In this paper, the Richardson-Lucy Iteration deconvolution method is presented. The method can not only suppress the noise and improve both SNR and resolution of the resulting 1D spectrum, but also reduce the ringings in the 1D resulting spectrum. Besides, the algorithm can be easily programmed and runs fast. The Richardson-Lucy Iteration formula is given in Section 2; Section 3 shows the extraction experiments on simulation 2D spectrum images; In Section 4, we practice our method on a real LAMOST spectrum image; Finally, the conclusion is given in Section 5.

2 Spectrum extraction method based on Richardson-Lucy iteration

Because the cross-talk on a LAMOST CCD is marginal, we only discuss how to extract 1D spectrum from an image with only one fiber. The image model can be found in (Li et al., 2019), which is:

[TABLE]

where, $N$ is the number of rows, $f(x,y)$ and ${\eta(x,y)}$ are the count and the noise at the $x$th row and $y$th column on the CCD, respectively, ${g_{j}}$ is the flux at the $j$th row, ${h_{j}(x,y)}$ is the value of the PSF at the $j$th row on $(x,y)$, and $G=(g_{1},g_{2},...,g_{N})$ is the 1D resulting spectrum. PSFs can be obtained by the method outlined in Guangwei et al. (2015) and Li et al. (2019).

The Richardson-Lucy Iteration algorithm (Richardson, 1972; Lucy, 1974) is a typical nonlinear iterative algorithm, which is widely applied in astronomical and medical image processing.

We can assume that the count of each pixel is independent and obeys Poisson distribution. We denote

[TABLE]

, which is the real spectrum image without noise. Then the likelihood function of the image is (Zou, 2001):

[TABLE]

Then,

[TABLE]

Let $\frac{\partial}{\partial g_{j}}[\ln P(f(x,y)|G)]=0,j=1,2,3,...,N$. Then

[TABLE]

or

[TABLE]

where, $j=1,2,3,...,N$.

To obtain the spectrum $G$, Meinel (1986) suggested the iteration formula:

[TABLE]

where, $k$ is the number of iteration.

If we set $p=1$, the above formula becomes Richardson-Lucy Iteration:

[TABLE]

3 experiments on simulation images

We use the similar construction method as (Li et al., 2019) to generate the simulation images. The simulation image is generated by an input 1D spectrum with 4,000 flux points given by the LAMOST pipeline convolved with PSFs with size of $13\times 15$ pixels. These PSFs are the linear interpolations on the basic PSFs from emission lines on an arc image. At last, the poisson noise is added. The resulting image is $4,000\times 15$ pixels.

We extracted the 1D spectrum from the simulation image by Equation 7 for 10 iterations. Besides, we also performed TDM and ALI for comparison. The Tikhonov regulation parameter for TDM was set to $0.02$, which is the best value in this extraction. The block sizes in TDM and ALI were both set to 20, while the computation precisions were both 100.

We ran these three algorithms on a DELL computer with CPU 3.30GHz and Operating system Windows 7 by the software MATLAB R2014a. The original Poisson noise and the 2D residuals of different deconvolution methods are shown in Fig. 1, while 1D residuals is shown in Fig. 2. These two figures illustrate that their 2D residuals are Poisson noise, and their 1D residuals are all at the similar level.

The the SNRs and computation times of different extraction methods are shown in Table 1, where the SNR is defined by the Equation 11 in (Li et al., 2019). From Table 1, we can see that Richardson-Lucy Iteration algorithm is much faster than ALI and TDM, and the SNR of the resulting spectrum extracted by Richardson-Lucy Iteration algorithm is also higher than those extracted by ALI and TDM. TDM is the slowest, because it took too much time to inverse calibration matrixes.

4 performance on a real LAMOST spectrum image

We apply Equation 7 on a real LAMOST multi-fiber spectra image. A resulting 1D spectrum after five and 10 iterations are shown in blue and red, respectively in Figure  3, while the spectra extracted by AEM, TMD and LIM are shown in black, green and grey, respectively. All spectra are overplotted together for direct comparison in the bottom of the figure. As we all known, the Gibbs phenomenon, which also called ringing artifacts in signal processing, always occurs in the results of deconvolution methods. The dashed lines in Figure 3 indicate the positions of overshoots of ringings in the spectra of TMD and LIM.

From the figure, we can see:

1. The spectra extracted by TMD, LIM and Equation  7 all have higher SNRs and resolutions than that extracted by AEM. Furthermore, their emission lines are all more symmetric, which means that these three methods can correct the distortions of PSFs on CCD.

2. The blue spectrum, which is the result of Equation  7 after five iterations, even has no overshoots. Its resolution is lower than those of the spectra extracted by TMD and LIM, but higher than the spectrum extracted by AEM.

3. After 10 iterations of Equation  7, the resolution of the resulting spectrum become similar to those of the spectra extracted by TMD and LIM, but the amplitudes of ringings are much lower.

5 Conclusions

This paper describes a deconvolution extraction algorithm based on the Richardson-Lucy Iteration to extract 1D spectra from LAMOST spectrum images. Compared with the spectrum extracted by AEM, the resolution and SNR of the spectrum extracted by the Richardson-Lucy Iteration are both higher. Compared with spectra extracted by TDM and LIM, the ringings of the spectrum extracted by the Richardson-Lucy Iteration are much weaker. Furthermore, Richardson-Lucy Iteration is the fastest deconvolution method to extract 1D spectra from a LAMOST image.

Acknowledgements.

This work is supported by the Joint Research Fund in Astronomy (U1531242) under cooperative agreement between the National Natural Science Foundation of China (NSFC) and Chinese Academy of Sciences (CAS), the National Natural Science Foundation of China (No. 11673036) and the Interdiscipline Research Funds of Beijing Normal University. Guoshoujing Telescope (the Large Sky Area Multi-Object Fiber Spectroscopic Telescope LAMOST) is a National Major Scientific Project built by the Chinese Academy of Sciences. Funding for the project has been provided by the National Development and Reform Commission. LAMOST is operated and managed by the National Astronomical Observatories, Chinese Academy of Sciences.