Blind deconvolution for astronomical spectrum extraction from two-dimensional multifiber spectrum images

Although advances in astronomical observation techniques and fiber optic technology have enabled two-dimensional (2D) multifiber spectrum images, astronomers often need one-dimensional (1D) spectroscopy to study physical and chemical properties of astronomical objects. Because of faint optical flux and light pollution, determining point spread functions (PSF) for large-scale multifiber spectroscopic telescopes is difficult. We propose a new optimal extraction method that uses blind deconvolution to extract 1D astronomical spectroscopy from 2D multifiber spectrum images without knowing the exact PSF. A comparison of the performance of our blind deconvolution extraction method and those of other extraction methods showed consistent results, indicating that the blind deconvolution extraction methodology is useful in analyzing 2D multifiber spectrum images and reducing fiber-to-fiber crosstalk without degrading the spectrum’s resolution. ©2017 Optical Society of America OCIS codes: (100.1455) Blind deconvolution; (100.0100) Image processing. References and links 1. S. A. Shectman, S. D. Landy, A. Oemler, D. L. Tucker, H. Lin, R. Kirshner, and P. Schechter, “The las campanasredshift survey,” Astrophys. J. 470(1), 172–188 (1996). 2. D. G. York, J. Adelman, J. E. Anderson, Jr., S. F. Anderson, J. Annis, N. A. Bahcall, J. A. Bakken, R. Barkhouser, S. Bastian, E. Berman, W. N. Boroski, S. Bracker, C. Briegel, J. W. Briggs, J. Brinkmann, R. Brunner, S. Burles, L. Carey, M. A. Carr, F. J. Castander, B. Chen, P. L. Colestock, A. J. Connolly, J. H. Crocker, I. Csabai, P. C. Czarapata, J. E. Davis, M. Doi, T. Dombeck, D. Eisenstein, N. Ellman, B. R. Elms, M. L. Evans, X. Fan, G. R. Federwitz, L. Fiscelli, S. Friedman, J. A. Frieman, M. Fukugita, B. Gillespie, J. E. Gunn, V. K. Gurbani, E. de Haas, M. Haldeman, F. H. Harris, J. Hayes, T. M. Heckman, G. S. Hennessy, R. B. Hindsley, S. Holm, D. J. Holmgren, C. Huang, C. Hull, D. Husby, S.-I. Ichikawa, T. Ichikawa, Ž. Ivezić, S. Kent, R. S. J. Kim, E. Kinney, M. Klaene, A. N. Kleinman, S. Kleinman, G. R. Knapp, J. Korienek, R. G. Kron, P. Z. Kunszt, D. Q. Lamb, B. Lee, R. F. Leger, S. Limmongkol, C. Lindenmeyer, D. C. Long, C. Loomis, J. Loveday, R. Lucinio, R. H. Lupton, B. MacKinnon, E. J. Mannery, P. M. Mantsch, B. Margon, P. McGehee, T. A. McKay, A. Meiksin, A. Merelli, D. G. Monet, J. A. Munn, V. K. Narayanan, T. Nash, E. Neilsen, R. Neswold, H. J. Newberg, R. C. Nichol, T. Nicinski, M. Nonino, N. Okada, S. Okamura, J. P. Ostriker, R. Owen, A. G. Pauls, J. Peoples, R. L. Peterson, D. Petravick, J. R. Pier, A. Pope, R. Pordes, A. Prosapio, R. Rechenmacher, T. R. Quinn, G. T. Richards, M. W. Richmond, C. H. Rivetta, C. M. Rockosi, K. Ruthmansdorfer, D. Sandford, D. J. Schlegel, D. P. Schneider, M. Sekiguchi, G. Sergey, K. Shimasaku, W. A. Siegmund, S. Smee, J. A. Smith, S. Snedden, R. Stone, C. Stoughton, M. A. Strauss, C. Stubbs, M. SubbaRao, A. S. Szalay, I. Szapudi, G. P. Szokoly, A. R. Thakar, C. Tremonti, D. L. Tucker, A. Uomoto, D. Vanden Berk, M. S. Vogeley, P. Waddell, S. Wang, M. Watanabe, D. H. Weinberg, B. Yanny, and N. Yasuda, “The Sloan Digital Sky Survey: Technical Summary,” Astron. J. 120(3), 1579–1587 (2000). 3. R. Stone, “Astronomy. China’s LAMOST Observatory prepares for the ultimate test,” Science 320(5872), 34–35 (2008). 4. R. Stone, “China. Astronomers hope their prize telescope isn’t blinded by the light,” Science 329(5995), 1002 (2010). 5. R. Z. Zhao, F. Wang, A. L. Luo, and Y. X. Zhang, “A method for auto-extraction of spectral lines based on sparse representation,” Guangpuxue Yu Guangpu Fenxi 29(7), 2010–2013 (2009). 6. X. J. Zhao, J. H. Cai, J. F. Zhang, H. F. Yang, and Y. Ma, “Automatic classification method of star spectrum data based on classification pattern tree,” Guangpuxue Yu Guangpu Fenxi 33(10), 2875–2878 (2013). 7. K. Horne, “An optimal extraction algorithm for CCD spectroscopy,” Publ. Astron. Soc. Pac. 98, 609–617 (1986). Vol. 25, No. 5 | 6 Mar 2017 | OPTICS EXPRESS 5133 #283467 https://doi.org/10.1364/OE.25.005133 Journal © 2017 Received 12 Jan 2017; revised 23 Feb 2017; accepted 23 Feb 2017; published 27 Feb 2017 8. J. G. Robertson, “Optimal extraction of single-object spectra from observations with two-dimensional detectors,” Publ. Astron. Soc. Pac. 98, 1220–1231 (1986). 9. T. R. Marsh, “The extraction of highly distorted spectra,” Publ. Astron. Soc. Pac. 101, 1032–1037 (1989). 10. R. I. Hynes, “An optimal extraction of spatially blended spectra,” Astron. Astrophys. 382(2), 752–757 (2002). 11. S. F. Sánchez, “Techniques for reducing fiber-fed and integral-field spectroscopy data: an implementation on R3D,” Astron. Nachr. 327(9), 850–861 (2006). 12. R. Sharp and M. N. Birchall, “Optimal extraction of fibreoptic spectroscopy,” Publ. Astron. Soc. Aust. 27(01), 91–103 (2010). 13. A. S. Bolton and D. J. Schlegel, “Spectro-perfectionism: an algorithmic framework for photon noise-limited extraction of optical fiber spectroscopy,” Publ. Astron. Soc. Pac. 122, 248–257 (2010). 14. J. Yu, Q. Yin, P. Guo, and A. L. Luo, “A deconvolution extraction method for 2D multi-object fiberspectroscopy based on the regularized least squares QR-factorization algorithm,” Mon. Not. R. Astron. Soc. 443(2), 1381– 1389 (2014). 15. W. Stefan, “Image restoration by blind deconvolution,” Diploma Thesis, TechnischeUniversitätMünchen and Arizona State University (2003). 16. G. R. Ayers and J. C. Dainty, “Interative blind deconvolution method and its applications,” Opt. Lett. 13(7), 547–549 (1988). 17. L. Blanco and L. M. Mugnier, “Marginal blind deconvolution of adaptive optics retinal images,” Opt. Express 19(23), 23227–23239 (2011). 18. W. He, Z. Zhao, J. Wang, B. Zhang, F. Qian, Z. Yang, M. Shui, F. Lu, J. Teng, L. Cao, and Y. Gu, “Blind deconvolution for spatial distribution of Kα emission from ultraintense laser-plasma interaction,” Opt. Express 22(5), 5875–5882 (2014). 19. K. R. Castleman, Digital Image Processing (Prentice-Hall,1996), Chap. 9. 20. D. Krishnan, T. Tay, and R. Fergus, “Blind deconvolution using a normalized sparsity measure,” IEEE Computer Vision and Pattern Recognition (2011). 21. K. Fuerschbach, J. P. Rolland, and K. P. Thompson, “Theory of aberration fields for general optical systems with freeform surfaces,” Opt. Express 22(22), 26585–26606 (2014). 22. B. Cui, Z. F. Ye, and Z. R. Bai, “Spectral extraction algorithm for LAMOST two-dimensional fiberspectroscopic images,” Chin. Astron. Astrophys. 33(1), 99–112 (2009). 23. A.-L. Luo, H.-T. Zhang, Y.-H. Zhao, G. Zhao, X.-Q. Cui, G.-P. Li, Y.-Q. Chu, J.-R. Shi, G. Wang, J.-N. Zhang, Z.-R. Bai, X.-Y. Chen, F.-F. Wang, Y.-X. Guo, J.-J. Chen, B. Du, X. Kong, Y.-J. Lei, Y.-B. Li, Y.-H. Song, Y. Wu, Y.-X. Zhang, X.-L. Zhou, F. Zuo, P. Du, L. He, W. Hou, Y.-Q. Dong, J. Li, G.-W. Li, S. Li, J. Song, Y. Tian, M.-X. Wang, K.-F. Wu, H.-Q. Yang, H.-L. Yuan, S.-Y. Cao, H.-Y. Chen, K.-X. Chen, Y. Chen, J.-R. Chu, L. Feng, X.-F. Gong, B.-Z. Gu, Y.-H. Hou, Z.-Y. Huo, H.-Z. Hu, N.-S. Hu, Z.-W. Hu, L. Jia, F.-H. Jiang, X. Jiang, Z.-B. Jiang, G. Jin, A.-H. Li, Q. Li, X.-N. Li, Y. Li, Y.-P. Li, G.-R. Liu, G.-Q. Liu, Z.-G. Liu, Q.-S. Lu, W.-Z. Lu, Y. Luo, Y.-D. Mao, L. Men, J.-J. Ni, Y.-J. Qi, Z.-X. Qi, H.-M. Shi, D.-Q. Su, S.-W. Sun, H.-J. Su, Z.H. Tang, Q.-S. Tao, L.-P. Tu, D.-Q. Wang, D. Wang, G.-M. Wang, H. Wang, J.-N. Wang, J. Wang, J.-L. Wang, J.-P. Wang, L. Wang, S.-G. Wang, S.-Q. Wang, Y.-N. Wang, Y. Wang, Y.-F. Wang, M.-Z. Wei, X.-X. Xue, X.Z. Xing, L.-Z. Xu, X.-Q. Xu, Y. Xu, D.-H. Yang, S.-H. Yang, Z.-Q. Yao, Y. Yu, H. Yuan, C. Zhai, E.-P. Zhang, J. Zhang, L.-P. Zhang, W. Zhang, Y. Zhang, Z.-C. Zhang, M. Zhao, F. Zhou, Y.-T. Zhu, J. Zhu, and S.-C. Zou, “Data release of the LAMOST pilot survey,” Res. Astron. Astrophys. 12(9), 1243–1246 (2012). 24. Zh. Q. Zhu, J. Zhu, Sh. Wang, S. Cao, G. Bao, and Z. Ye, “Spectrum extraction for 2-D fiber spectrum images based on 2-D Gaussian model,” Publ. Astron. Soc. Aust. 28(04), 357–364 (2011). 25. R. Sharp and H. Parkinson, “Sky subtraction at the Poisson limit with fibre-optic multiobject spectroscopy,” Mon. Not. R. Astron. Soc. 408(4), 2495–2510 (2010). 26. Zh. Q. Zhu, J. Zhu, Sh. Wang, Sh. H. Cao, G. Zh. Bao, and Z. Ye, “Spectrum Extraction for 2-D fiber spectrum images based on 2-D gaussian model,” Publ. Astron. Soc. Aust. 28(04), 357–364 (2011). 27. D. Krishnan and R. Fergus, Fast image deconvolution using hyper-laplacian priors. In NIPS, 2009.


Introduction
With the development of astronomical observation techniques and fiber optic technology, a variety of large-scale fiber-based spectroscopic telescopes with a large field of view, deep depth and high efficiency are used worldwide.These telescopes include the Las Campanas Red shift Survey (LCRS) [1], Sloan Digital Sky Survey (SDSS) [2], and GuoShouJing Telescope (GSJT) in Xinglong, Hebei Province, China [3,4].The GSJT was formerly known as the Large Area Multi-Object Spectroscopic Telescope (LAMOST), which consists of "1meterwide hexagonal sections of its two mirrors and the 4000 optical fibers on its focal surface that will feed starlight into a battalion of spectrographs" [3].The LAMOST makes it possible to collect spectral deluge of tens of millions of galaxies and stars, and thus to offer new insights into galaxy formation, including our own Milky Way [3].Specifically, the LAMOST exceeds all other sky survey telescopes combined, and it has 4000 optical fibers that simultaneously sample 4000 astronomical objects and present 4000 respective two-dimensional (2D) fiber spectroscopies on a set of charge-coupled devices(CCDs) [3].However, astronomers need one-dimensional (1D) wavelength-dependent spectroscopy to study a variety of physical and chemical properties of astronomical objects [5,6].Extracting 1D astronomical spectra from 2D 4000-fiber spectrum spectrographs is one of the most important and challenging steps in LAMOST data analysis, since accurate extraction affects directly the quality of the sky survey observations and subsequent interpretations.
Four extraction methods for 2D multifiber spectrum images have been studied; they are the aperture extraction method [7,8], optimal aperture extraction method [9,10], profile fitting method [11,12], and deconvolution method [13].The current, widely adopted extraction algorithms are the first three methods, but they are inaccurate with reduced resolutions for faint objects.However, the deconvolution method does not suffer from this drawback.It can characterize the resolution accurately to upper limits of native instruments.In a previous work on the point spread function (PSF), we developed a deconvolution algorithm that is based on the principles of optical imaging to extract the spectra from a multifiber spectroscopic image with graphics processing unit acceleration [14].We found that the method is feasible for solving computational problems that involve a large amount of data [14].
In theory, we could utilize a direct deconvolution operation/method to extract 1D spectroscopy signals from 2D multifiber spectrum images if given a known PSF [15].For a laboratory-based spectroscopic/imaging device, a point source can be used along with the imaging device to experimentally measure the corresponding PSF, which should be estimated in advance and incorporated into the deconvolution formula.Currently, a common method used to estimate a PSF is by quantifying or calibrating an arc lamp spectrum image.In practice, however, PSF is very difficult to obtain for the LAMOST system.The first reason is that there are several types of analytical expressions for PSFs, such as the "disc,""Gaussian," and exponential families, whose parameters are often unspecified for an arc lamp spectrum.The second reason is that each light source at the slit of the spectroscope of the LAMOST system is not a point source.Each fiber has a diameter of 320 μm.If we treat each fiber as a point source, this approximation will cause large computing errors.Meanwhile, the LAMOST system was built a decade ago with 4000 fibers on the focal plane; the direction of each fiber needs to be adjusted every time to locate a sky object during the LAMOST survey.Therefore, the PSF of a fiber may vary from time to time and cannot be estimated for every observation night.To solve these problems, in this paper, we propose a novel method for building a spectrum extraction model using the blind deconvolution principle to extract 1D spectroscopic information from 2D fiber-based images.
To date, the blind deconvolution principle has been used in a variety of research studies [16][17][18], but it has not been used in spectroscopic extraction.In this paper, we present an algorithmic framework based on the blind deconvolution principle for the extraction of optical fiber spectroscopy.This algorithm performs a simultaneous estimation progress for both the object and PSF of the spectroscopic/image device.Consequently, we can obtain an optimal estimation of the PSF and the extracted spectrum simultaneously.
The remainder of this paper is organized as follows.Section 2 describes the imaging model, the blind deconvolution algorithm, and the proposed spectrum extraction framework.Section 3 reports the implementation and analysis of results, followed by a discussion in Section 4 and conclusions in Section 5.

Imaging Model
In principle, a spectroscopic measurement of an astronomical object collects light passing through a set of optical instruments, such as a collimator mirror, dichroic beam splitter, and diffraction grating.These optical instruments will disperse the detected light into spectral components such that it is recorded by a CCD.It is unavoidable that all of the optical devices give rise to broadening and distortion effects on the measured spectroscopic signals of the objects.Such an effect on a point light source is called the PSF [14].
Let us first consider a general imaging model.Monochromatic light at a wavelength λ imaged on an observation plane can be represented generally by the following integral equation [13]: where ( ) For a CCD image, the corresponding discrete form can be written as [19] ( ) ( ) ( ) Following [14], Eq. ( 2) can be expressed in matrix form as where G represents the CCD-recorded image, H is a kernel of the PSF, ⊗ denotes the convolution operator, and F is a matrix in which the rows are related to the wave lengths and columns are composed of multifiber spectrum lines.For the LAMOST, F contains 250 spectrum lines with a wavelength range from 570 to 900 nm (red channel).At each wavelength position, point source will spread over 15 pixels on the detector.Because each fiber spot has a diameter of 320μm, in our imaging model we assume that each fiber spot includes at least 15 pixels.This mean that a spatial size of about 21.3 μm is covered by each pixel (i.e., 320 μm/15 pixel), which is similar to the double physical size of the pixel.
The original optical source should be point source.However, in the imaging process, the point source is spread by the PSF and changed into a spot.Therefore, we can assume that the original image is sharper than the CCD image, and, according to the above imaging model, spectrum extraction of 2D multifiber spectrographs can be considered an image restoration problem.

Blind deconvolution algorithm
Blind image deconvolution, an ill-posed problem that requires a regularization to solve, is often adopted to resolve image restoration problems.A large number of blind deconvolution algorithms have been proposed in the literature.In this work, we adopted Krishnan's method [20] because it is fast and robust when working with small kernel errors, and also robust in the choice of parameters.
In this algorithm, Krishnan first conducted blind kernel estimation for the high frequencies of an image y, a concatenation of the gradient images, x ∇ g and y ∇ g .The cost function is Here, x is the unknown sharp image in the high-frequency space, k is the unknown blurring (PSF) kernel, || || is the alpha norm, i.e., the log-prior function with the parameter alpha, and ⊗ denotes the convolution operator.In Eq. ( 4), the first term is the likelihood that accounts for the model.The subscripts 1 and 2 in Eq. ( 4) denote l 1 and l 2 norm regularizers, respectively, i.e., the respective values of the alpha parameter are 1and 2.The l 1 norm is widely used to impose signal sparsity, but it is scale variant so the norm can be minimized by simply reducing the signal.In a blind deconvolution setting where the kernel is only loosely constrained, minimizing the l 1 norm for the high frequencies of the image will result in a blurry image.Therefore, the l 1 /l 2 function is used as a normalized version of the l 1 norm, making it scale invariant.The third term is a regularization of k for the reduction of noise in the kernel.The scalar weights μ and ψ control the relative strengths of the kernel and image regularization terms.
The standard approach to optimizing such a problem is to start with an initialization of x and k, and then alternate between x and k updates.The alternating update progress of x and k can be found in [20].Once the kernel k has been estimated, a nonblind deconvolution method is adopted to obtain f from g.This deconvolution algorithm uses a continuation method to solve the following cost function: where α stands for the α norm.With these algorithms, from the observed CCD image g, we can obtain the deconvolved image f.

Proposed spectrum extraction framework
However, the blind deconvolution algorithm cannot be used directly to extract the spectrum.The details of spectrum imaging and extraction progress need to be considered.The nature of spectrum extraction is to obtain the 1D spectrum from a 2D image.The result of blind deconvolution is narrow and sharp, similar to the spectrum extraction progress.However, it is still an image of 2D fiber spots.For a multifiber spectrum imaging system, each fiber cannot be approximated as a point source because it has a certain size (320μm) and covers about 15 pixels.Therefore, an integral matrix T is used in our method to perform the integral in the slit height direction to obtain the spectra of 250 fibers.Furthermore, to obtain the "undistorted" spectra, we need to not only conduct sky subtraction and noise reduction, but also run a deconvolution computation for Eq.(3).To do so, several significant challenges exist: (1) the 2D PSF for the multifiber spectral measurements must be known, (2) the measured optical data have no or minimal noise, (3) the spectroscope has no aberration.However, because of optical noise effects, imperfection of spectroscopic systems, and unknown spectral energy distributions of the imaged objects, it is a complex process to determine the PSF either analytically or experimentally.Moreover, for real-world application of 2D multifiber spectrum imaging, the computational complexity should also be considered.This consideration is especially true when working with data from high-tech surveys such as the LAMOST, which has recorded a total of 4096 × 4096 pixels on a CCD.
From a computational point of view, in a common spectroscopy system the PSF is spatially varied because of system aberration, and hence the matrix that represents the convolution kernel H is not a Toeplitz matrix.This circumstance makes computational complexity very high.If the matrix is m × n and non-Toeplitz, it will require an enormous amount of computational memory to store, leading to a very high spatial complexity.Furthermore, the time complexity will depend on the algorithm.If we use the least squares method to compute the generalized inversion, the time complexity is approximately of the order of O(n3), which cannot be implemented on an ordinary workstation.If we use the least squares QR factorization algorithm to perform the deconvolution operation as we did before [14], the matrix must be sparse, and the time complexity is of the order of O(m + n).If the matrix is not sparse, we must consider some other methods to address the high-complexity situation in order to solve the real-world spectrum extraction problem.
Thus, to reduce computational complexity for practical implementation, in this paper we propose to use an aberration matrix for 1D astronomical spectral extraction from a 2D multifiber spectrograph.The aberration matrix a H is instrument-dependent. Precise computation of a H would require the use of nodal aberration theory [21], but we can empirically obtain it with an ad hoc method.For example, spectrum contour peak trace is a simple method that can be applied for aberration correction.
In our proposed methodology, we mainly focus on the spectrum extraction method, not the sky subtraction or noise reduction problems.Our blind deconvolution approach is an alternating variable optimization method that was developed to estimate both the 2Dmultifiber spectrograph and PSF.The procedure of the proposed spectrum extraction framework is as follows.
Step 1.For a given multifiber CCD image g , compute gradient images x ∇ g and y ∇ g to construct the image y.
Step 2. Use y to estimate the PSF kernel k using Eq. ( 4).
Step 3.With the estimated kernel k, obtain a deconvolved image f using Eq. ( 5).
Step 4. With the precomputed aberration matrix a H , let a = f H f a and compute to obtain the aberration-corrected image f a .
Step 5.With a predefined matrix T, which is constructed using a fiber spot kernel t, compute the convolution = f Tf b a by performing an integral in the slit height direction.Step 6.Finally, obtain the multifiber spectrum matrix f b , in which each column represents one spectrum line.

Implementation
In our experiments, the 2D fiber spectrograph of the LAMOST project [22,23] was used to verify our proposed method.The full size of a CCD image from the LAMOST was 4096 × 4096 pixels, and each CCD image contained 250 fiber spectra.Initial value setting of the kernel k has no accepted initial method up to now.The method generally used is to choose the value by prior knowledge.In our experiments, the size of the initial kernel k is 3 × 3. We also analyzed the impact of scale on the other parameters in spectrum extraction, which will be discussed in Section 3.2.
We perform an integral with a predefined matrix T, which is constructed using a fiber spot kernel t, in the slit height direction according to Eq. ( 1).Then, we were able to obtain the convolution kernel for this fiber spot, which was designed with a size of 15 × 15.The kernel t takes the form .
a a a a In this kernel t, each column corresponds to the spatial information at a given wavelength, and each row provides the spectral information at a corresponding astronomical location given by the respective fiber.All of the elements are zero except for the elements in the first column.The elements in the first column are 0,1,…,1,0.
We only perform a simple optical fiber tracing to straighten the curved spectrum lines, not to address the distortion.When the observed celestial objects have high magnitudes, the observed targets appear as dark celestial bodies, and the obtained 2D fiber spectrum data have a very low signal-to-noise ratio (SNR).For this type of spectrum data, the spectrum profiles are polluted so heavily that the centers of the profiles deviate from the ideal locations, and the shapes of the contours are no longer smooth and regular, meaning that the profile fitting method cannot be performed directly in the spectrum extraction process [24].

Analysis of results
This section presents the results of spectra extracted using the proposed blind deconvolution method.First, we applied the proposed blind deconvolution method to extract real 2D arcs and object spectra from the data obtained by the LAMOST.We then compared these results with the spectra that were extracted using the aperture and non-blind deconvolution methods reported in our previous work [14].Note that we have not performed sky subtraction [25] in this study because it can be performed after we validate this proposed method.To avoid diffusing the effects of different flux base values for different extraction methods, the flux magnitude of the extraction methods has been changed to the same scale.
Figure 1 shows the estimated PSF in a 3D view after implementing and utilizing the proposed blind deconvolution algorithm.The regularization parameters need to be selected properly.Figure 2(a) is the image constructed by one simulated arc spectrum image with added noise. Figure 2(b) is the result obtained by utilizing the proposed blind deconvolution method to extract from Fig. 2(a) and reconstructing with the same estimated 2D PSF of blind deconvolution above is shown in Fig. 2(c).The residual image obtained by subtracting Fig. 2(b) from Fig. 2(a), and is shown in Fig. 2(d).The Fig. 3 is the 3Dview of the Fig. 2. Using the blind deconvolution method, crosstalk [13] between fibers can be improved significantly (see Figs. 2 and 3).  Figure 4 shows a comparison of spectrum wave forms extracted using the aperture and proposed blind deconvolution methods.The spectra extraction results and detailed results are shown in Figs. 5 and 6, which show the practical LAMOST object spectra extracted using the aperture, deconvolution, and proposed blind deconvolution methods.Fig. 6.Zoomed detailed object spectra extracted using the aperture, deconvolution, and blind deconvolution methods.
One advantage of using the blind deconvolution method is that it can resolve crosstalk problems that were not solved in the aperture method.From Figs. 2 and 3, we can see that crosstalk between the fibers can be improved significantly and that only a very small residual error remained because of the noise effect.
Another advantage is that the blind deconvolution method does not degrade the resolution compared with the other methods.From Fig. 4, we can see that the spectra waveforms extracted using the proposed method are narrower than those extracted using the aperture method.Therefore, we can draw the conclusion that spectra extracted using the proposed blind deconvolution method can characterize resolution accurately up to the upper limit of the native instrument by adjusting the appropriate regularization parameters.If the regularization parameters become small enough, then noise can be suppressed.If they become too large, then we know that too much noise was recovered.Therefore, by selecting appropriate regularization parameters, we can make a tradeoff between spectrum resolution and noise suppression.
Not only does the blind deconvolution method have the aforementioned advantages, but the extraction results also prove to be better than those extracted using the aperture method.From Figs. 5 and 6, we can see that resolution of the spectrum extracted using the proposed method is almost the same as that of the spectra extracted using the aperture and deconvolution methods [14], meaning that the proposed method can obtain satisfactory results as well.
To demonstrate the accuracy and effectiveness of the proposed method, we adopted a SNR for evaluating the results of different extraction methods.The SNR of the extracted spectrum is calculated from the ratio of the spectra extracted from the realistic simulated spectral detector data minus instrument noise of the detector and instrument noise of the detector.
We have simulated an image with 250 input real spectra, which come from the LAMOST system and which are the effective spectra identified by the astronomer.The simulated process is as follows.
a. Trace the peak of spectra; b.Use the results of tracing find the value of pixels nearby the spectra from the original image (above 7 pixels and below 7 pixels); c.Use these values and 250 spectra to build one simulated real spectrum data.We generate realistic simulated spectral detector data with known realistic input object spectra obtained from LAMOST, and kernels and noise of the instruments.Different source spectra are simulated to compare the SNR.The comparison results of the SNR are shown in Table 1.f stands for an object spectrum with a 10dB white noise convoluted with a PSF, whose full wave at half maximum is 2 pixels, and 3 o f stands for an object spectrum with a 10dB white noise convoluted with a PSF, whose full wave at half maximum is 4 pixels.Even for the input spectra with noise, the results show that the proposed method is better than that in our previous work [14].In [14], we have done a comparison study among aperture extraction, profile fitting method and deconvolution method, and obtain the results that deconvolution method is better than other two methods.Now we focus on the work to compare the proposed method with the state of the art extraction method (deconvolution method) only in this paper.From Table 1, we also find that even though the spectra 2 o f and 3 o f are not narrow compared with 1 o f , we can still obtain better results.This can confirm that performance of the proposed method does not depend on the presence of sharp, narrow emission features in the input spectrum.

Discussion
As mentioned above, researchers must determine the PSF before using the deconvolution method to extract the spectra.In a laboratory, one can shine an arc lamp to estimate the PSF by measuring the image spread.This procedure is similar to the profile fitting method [11,12,26], and an analytical form for the PSF must be provided in the first place.However, there is no theory to guide the selection of a suitable PSF form from several types of profiles, such as the "disc,""Gaussian," and exponential families.In most cases, the PSF is selected empirically to fit the shined light spot in the spectrum image.Because of the existence of instrumental aberration, the shape of the shined light in the spectrum image varies at different positions, and it becomes difficult to determine which spot is selected for the estimation.
In the LAMOST system, there are 4000 fibers on the focal plane, and these fibers are controlled to move to the positions where the observed astronomical objects are imaged well.Figure 7(a) shows the focal plane of the LAMOST, and Fig. 7(b) is the enlarged focal plane, which shows the fiber movement control module.At astronomical observation night, the position of each fiber will be adjusted to locate an astronomical object.This movement of the fiber position will lead to a possible change of the PSF.Furthermore, in the real LAMOST environment, limited observation time restricts us from calibrating the PSF with the arc lamp line during each astronomical observation night.In the faint light environment, it is even more difficult to estimate the PSF accurately [26].In contrast to the previous PSF estimate method, the largest advantage of blind deconvolution is that it does not require prior knowledge of the PSF, nor does it need to be given an analytical form of the PSF.From Fig. 1, we can see that the estimated PSF in the blind deconvolution framework looks like a bell function.The long axis is approximately 8 pixels in length, and the short axis is approximately 4 pixels in width.
It is worthwhile to address that in practice, PSF estimation using the arc-lamp line method is not practical to integrate into the LAMOST data processing pipeline because it requires interactively selecting the PSF form and estimating the model parameters.With the blind deconvolution method, it becomes possible to integrate the data processing pipeline into the LAMOST because it automatically calculates the PSF and extracts the spectra, and hence no human interactive input is required when running the program.
In our previous work [14], the situation of a spatially varied PSF was not considered.In this paper, aberration correction is used to resolve the problem of spatial complexity, which can be reduced by the proposed method.The runtime of the deconvolution extraction method [14] for a full CCD image was approximately 2 to 3 minutes when running on a computer that was equipped with an Intel Xeon X5680 3.33 GHz processor and 12 GB RAM.However, the time duration of the PSF estimation before the spectrum extraction is at least a half hour, which includes establishing the arc lamp environment, obtaining the CCD image, and estimating the PSF.Therefore, the method of [14] cannot be used in the practical LAMOST environment.The method proposed in this paper only requires approximately 8 minutes for the PSF estimation and spectrum extraction, and hence can be used in the practical LAMOST environment.

Conclusions
We presented a new algorithm for the extraction of 1D spectra from 2D multifiber spectrographs based on blind deconvolution.This framework can extract 1D astronomical spectra and estimate the telescopic imaging system's PSF simultaneously from 2D multifiber spectroscope images.Experimental results demonstrated that the proposed algorithm is more efficient and cost effective than a previous algorithm [14].As a more precise spectrum extraction method, the proposed algorithm has practical applications in 2D multifiber spectrum imaging telescope system.

Appendix A
Here shows the core algorithm of the Krishnan's method [20].

−
represents the image intensity at the position (x, y) in the observation plane, represents the PSF value of the monochromatic light, and ε denotes an additive noise at the position (x, y).This means that the observed image is obtained by the convolution of the source image and a PSF, plus an additive noise term.

Fig. 1 .Fig. 2 .
Fig. 1.PSF of an object spectrum in a 3D view in the proposed blind deconvolution framework.

Fig. 4 .
Fig. 4. Comparison of spectrum waveforms extracted using the aperture and proposed blind deconvolution methods.