Point source detection performance of Hard X-ray Modulation Telescope imaging observation
The Hard X-ray Modulation Telescope (HXMT) will perform an all-sky survey in hard X-ray band as well as deep imaging of a series of small sky regions. We expect various compact objects to be detected in these imaging observations. Point source detection performance of HXMT imaging observation depends not only on the instrument but also on its data analysis since images are reconstructed from HXMT observed data with numeric methods. Denoising technique plays an import part in HXMT imaging data analysis pipeline alongside with demodulation and source detection. In this paper we have implemented several methods for denoising HXMT data and evaluated the point source detection performances in terms of sensitivities and location accuracies. The results show that direct demodulation with -fold cross correlation should be the default reconstruction and regularization methods, although both sensitivity and location accuracy could be further imporved by selecting and tuning numerical methods in data analysis of HXMT imaging observations.
Vol.0 (200x) No.0, 000–000
Hard X-ray Modulation Telescope (HXMT) is a planned Chinese space X-ray telescope. The telescope will perform an all-sky survey in hard X-ray band ( – ), a series of deep imaging observations at small sky regions and pointed observations.
We expect a large number of X-ray sources, e.g., AGNs, to be detected in its all-sky survey. We also expect through a series of deep imaging observations at the Galactic plane X-ray transients can be detected(Li, 2007; Lu, 2012). Therefore the point source detection performance is one of our concerns on HXMT data analysis.
Methods and corresponding sensitivities of pointed observation have been discussed by Jin et al. (2010). In imaging observations such as all-sky survey and deep imaging at small sky regions, a variety of data analysis aspects and methods are involved. First, images are computed instead of recorded directly by the optical instrument. Mapping as well as image reconstruction methods are useful. The direct scientific data products from imaging observations of HXMT are still scientific events, more specifically, X-ray photon arrival events. The attitude control system (ACS) of the spacecraft reports the attitudes periodically. We take these reported attitudes as nodes to perform certain interpolations to determine the instantaneous attitude for each scientific event. In this way, a set of parameters are assigned to each event, including the coordinates on the celestial sphere where the telescope is pointing at. Hence we call this process events mapping, where scientific events are mapped from time domain to the celestial sphere. The product of this process is refered to as the observed image, which implies the dimensionality of the data. Second, the exposure to a specific source is limited more strictly thus the signal-to-noise ratio (SNR) is tightly restricted. Hence the importance of regularization methods becomes prominent. Finally, a picture is worth a thousand words. Various information can be extracted from an image by numerical methods.
In this paper we investigate the point source detection performance of the imaging and detecting system synthesised from the telescope as well as diverse combinations of data analysis methods, especially the regularization methods.
1.2 Modulation in HXMT imaging observation
The PSF of HXMT HE, which is a composite telescope consisting of 18 collimated detectors, describes the response of the telescope to a point source when the telescope is pointing at the source as well as different locations around it. In other word, the PSF is a density function of the distribution of responses occur in different observation states, which is denoted by the instantaneous attitude of the telescope as a spacecraft. So the PSF takes the attitude of the telescope as its input.
We use the proper Euler’s angles to describe the attitude of the telescope, i.e., and the longitude and latitude of the pointing, as well as denoting the rotation angle of the telescope around its own pointing axis, namely, the position angle. The modulation equation that corresponds to the imaging observation over the sphere surface is
where is the object function (i.e., the image) defined in a compact subset of the sphere surface , is the modulation kernel function that relates the value of object function defined on a neighbourhood of the point of the subset to the instantaneous response of the telescope, while its status is .
The modulation kernel function is determined by the point spread function (PSF) of the telescope. Suppose a unit point source is located at the zenith of the celestial sphere, i.e., the point in the corresponding Cartesian coordinate system. Fix the position angle , while slew the telescope across the polar cap, and assign responses of the telescope to the unit source to pixels on the celestial sphere. Then we obtain a map represents the PSF on the celestial sphere with fixed rotation angle .
The map is then projected to a tangent plane of the celestial sphere by gnomonic mapping, i.e.,
where and are local Cartesian coordinates on the tangent plane. Now we have
i.e., the PSF defined on the tangent plane.
To provide for data analysis we have two discrete forms of ,
where is the neighbourhood of the pixel of the 2-D pixel grid, and the normalized form
Given the discrete image , the detection area and the duration of exposure on each pixel , the observed data is
where is the constant background count rate, denotes the convolution, and is the normalized PSF on the tangent plane. Eq. 6 approximates modulations in HXMT imaging observations.
Distortion occurs when projecting a set of points from (or to) a sphere surface to (or from) a plane. For example, distance between any two points, area of any continuous subset, angle between any two crossing lines (or tangents of curves) are altered non-uniformly. On the other hand, the rotation angle is not always fixed during HXMT imaging observations. Both distortions in image reconstruction from HXMT observed data and rotations during imaging observations have been discussed by Huo & Zhou (2013). Here for the sake of simplicity, we ignore them in this article.
2 Numerical methods
2.1 Single point source detection performance estimation
We estimate the single point source performance in terms of sensitivity and position accuracy through the following procedures.
Determine flux threshold for point source detection.
Simulate a frame of observed data contains only background counts.
Run denoise program on the simulated data to try to increase the signal-to-noise ratio.
Demodulate the denoised data.
Run SExtractor, a source detection program, by Bertin & Arnouts (1996) on the demodulated image to detect point sources and extract their intensities, coordinates and other parameters. At this point a catalog of point sources is compiled from the simulated data. Point sources detected here, i.e., from images demodulated from background data are false detections.
Choose a cut from the histogram as the flux threshold so that a certain percent of the false detections are rejected and the rejection percentage is precise enough. The rejection percentage, e.g., or etc., reflects the significance of detections above the corresponding threshold.
Estimate detection efficiency and position accuracy of a point source of specific flux intensity.
Simulate observed data contains a single point source of given flux intensity located at in the model image.
Examine each detection is the catalog. Provided the th source in the catalog is detected at in the demodulated image, the distance between the extracted source and the true point source
as well as the flux of the extracted source are investigated to determine whether the th source is a true source or not. We define the score of the current catalog in detecting the single point source as
where is the index of the current catalog, and are position accuracy and flux thresholds, therefore the score reveals whether the th catalog contains the true source or not, in other word, through the previous steps (simulated observation, denoising, demodulation, source extraction and thresholding) whether we have detected the true source effectively. If we have, the outcome of these steps is counted as an effective detection of the true source, otherwise it is ineffective.
which is defined as the detection efficiency.
Let be the position of the brightest source in the th catalog, the position accuracy is calculated as
Find an flux intensity so that the corresponding detection effeciency . The intensity marks the point source sensitivity of the detecting system synthesised from both the telescope in specific status and the data analysis program chain.
2.2 Imaging and mapping
Direct demodulation (DD) method (Li & Wu, 1994) is used to estimate the true image from observed data. Residual map calculated with CLEAN algorithm (Högbom, 1974) is used as lower limit constraint in DD method. The skewness of the residual map is calculated in each CLEAN iteration and its minimum absolute value serves as the main stopping criterion of iterations.
In contrary to detecting and extracting point sources from demodulated images, it is also feasible to do so from cross-correlated maps as long as the point sources are isolated with each other, compared to the FWHM of the PSF, since the position of a peak of the expected value of such a map coincides with the position of a source regardless of the PSF. Cross-correlating the observed data and the PSF yields the correlated map
where denotes cross correlation, and denote the observed data and the PSF respectively. The peak of that coincides with a point source of flux intensity is
according to Eq. 6. On the other hand, the background variance of the correlated map is
since follows Poisson distribution. Hence the significance of the peak in term of numbers of is
where is total duration of exposure on the 2-D pixel grid, is the square root of the arithme mean of over the 2-D pixel grid, which is determined by the PSF as well as the range of the pixel grid only, provided the pixel grid is fine enough (see Eq. 5).
Cross-correlation significance of isolated point source defined Eq. 14 can be evaluated directly, given only the flux intensity of the source, the background count rate, the duration of exposure, the detection area and the PSF. Hence it is determined by the object (i.e., the point source), the telescope and the status of observation thus effects from data analysis programs are minimized.
We use the cross-correlation significance as a reference. For example, in our simulations an isolated point source of flux has significance.
2.3.1 Linear methods
Gaussian smoothing is often used in digital image processing to suppress the noise at the cost of reduction in resolution. Trade-off between noise suppression and resolution conservation is adjusted through the standard deviation of the Gaussian distribution serving as the smoothing kernel function. The best resolution of HXMT HE observed data is about , limited by the FWHM of its narrow-field collimator. We set to so that the FWHM of the Gaussian kernel is also .
fold cross correlation transform () can be used in DD method to regularize the ill-posed problem, more specifically speaking, to ensure the convergence as well as stability of the solution(Li, 2003). Here we put this technique in the denoising category. We have tested fold and fold cross correlated DD methods respectively in this article.
2.3.2 Non-linear methods
Non-local means denoising (Buades et al., 2005) is an edge-preserving non-linear denoising method. To increase its performance we have implemented this method with fast fourier transforms (FFTs). The pixel-wisely evaluation of the general Euclidean distance between the th pixel and other pixels of an image is replaced by
where is neighbourhood of the th pixel, is the pixel in the neighbourhood, is the image, is weight coefficient of the distance function. We use pixels Gaussian kernel with standard deviation (in pixels) as the weight coefficient in our simulations. We reduce the complexity of NLMeans method by computing convolutions in Eq. 15 with FFTs instead of using searching windows. (Buades et al., 2008).
Median filter is another non-linear edge-preserving denoising method. This method is effective in removing salt-and-pepper noise in digital images. In HXMT observed data such noise is incurred typically by missing data or charged particles of cosmic rays. We fixed the size of the filter at (about pixels) in HXMT HE data denoising.
The last non-linear denoising method included in this article is the adaptive wavelet thresholding with multiresolution support (Starck & Murtagh, 2006). The multiresolution support of a noisy image is a subset that contains significant coefficients only. So wavelet coefficients that dominated by noise are discarded. In this article we implemented the non-iterative algorithm. The spline wavelet is used for multiresolution decomposition.
3 Simulation and results
3.1 In-orbit background simulation
HXMT HE in-orbit background count rate ranges from to (Li et al., 2009). We use a constant count rate to simulate the average in-orbit background of HXMT HE in this article.
3.2 Source energy spectrum and telescope detection efficiency
to model energy spectra of Crab-like sources from to , where is in and the flux is in .
The detector efficiency of HXMT HE is derived from its simulated energy response, as shown in Fig. 1.
Detection efficiency of HXMT HE to a Crab-like source is . Count rate of HXMT HE corresponds to the intensity is , given the detection area of HXMT HE is about .
3.3 PSF and modulation
We use the diagram shown in Fig. 2 to simulate the PSF on the tangent plane.
We use the concentric average of , namely,
and the cumulative sum of ,
to characterize the radial fade-out of the PSF and the concentration of the PSF respectively, as plotted in Fig. 3.
From Fig. 3 we see that the FWHM of the simulated PSF is about while responses occur in a diameter of .
Despite the fact that the direct observed data is scientific events instead of 2-dimensional images, we start our simulation from simulated observed data in the form of images defined on 2-dimensional cartesian pixel grid. We use a model image for simulations. Given the diameter of the PSF, the central region is fully modulated, that is, all contributions to observed data in this region are from the model image only. The surrounding region is partially modulated, i.e., part of the contributions to observed data is this region are from the model image.
The average exposure per unit solid angle is in HXMT half-year all-sky survey. The partially modulated region is discretized by pixel grid, so . The detection area of each HXMT HE detector is approximately so the total area of all the 17 HXMT HE detectors is .
We have simulated frames of observed data that contains the in-orbit background counts only for each methods to estimate the corresponding flux thresholds by the method specified in Procedure 1 of Section 2.1. From false detections we have obtained and flux thresholds, see Table 1 for results.
|Denoising||thres., in||thres., in|
We have simulated frames of observd data that contains a Crab-like point source of , , , , , , , , and respectively, i.e., frames of observed data in total.
With these simulated data we have estimated location accuracies as well as detection efficiencies by the method described in Procedure 2 of the following methods:
DD without denoising,
DD with -fold cross correlation,
DD with -fold cross correlation,
DD with Gaussian smoothing, where ,
DD with NLMeans filtering, the size of the filter is and (both parameters are in pixels),
DD with median filter,
DD with median filter, and
DD with adaptive wavelet thresholding.
Implementation details of the above methods are in Sect. 2.3.
Table 2 shows the location accuracies on simulated data of , , and point sources.
Table 3 shows the detection efficiencies on simulated data of , , , and point sources.
Although the RL iteration itself we employed in DD method is total-counts-conservative, i.e., the sum of counts of each pixel is conserved after the iteration (Richardson, 1972), the regularizations including the background constrains as well as various denoising techniques are not necessarily counts-conservative. As a result, the absolute flux threshold (Table 1) for rejecting false detections doesn’t reflect the sensitivity directly but the detection efficiency (Table 3) does.
Errors of flux thresholds, location accuracies and detection efficiencies in Table 1, Table 2 and Table 3 are calculated by bootstrapping. For example, following Procedure 1 a set of frames of demodulated images are obtained, from which false detections are calculated and a histogram is plotted, where both the and thresholds are determined. Now let’s generate a new set with the same volume by resampling from the original set with replacement in order to calculate both the thresholds again. We repeat this resampling process until we get enough thresholds calculated to estimate their standard deviations. In Table 1, Table 2 and Table 3 each of the errors is a standard deviation that calculated from resampled sets.
A comprehensive summary of all the tested methods on all simulated data is shown in Fig. 5.
According to the results from our tests no denoising method shows significant advantage over -fold cross correlated DD method in single point source detection efficiency. Therefore it’s suggested that -fold cross correlation should be the default regularization method for single point source detection in HXMT imaging data analysis.
In the other hand, the location accuracy can be improved with alternative denoising methods, such as median filter, wavelet thresholding, Gaussian smoothing with little kernel, or without denoising, according to the results in this work.
This article is focused on the single point source detection of HXMT imaging data analysis, where other interesting topics can not be all covered. Although the alternative denoising methods have been out-performed more or less by the default -fold cross correlation in their contributions to the detection efficiency, they have shown certain advantages in location accuracy, and these features are promising for locating bright transients, multiple sources resolving, and so on.
In this work we made use of SciPy(Jones et al., 2001–), PyFITS and AIRE. PyFITS is a product of the Space Telescope Science Institute, which is operated by AURA for NASA. AIRE is a set of computing facilities initiated by Tsinghua Centre for Astrophysics. This work was supported by National Natrual Science Foundation of China (NSFC) under grants No. 11373025, No. 11173038 and No. 11403014 as well as Tsinghua University Initiative Scientific Research Program under grant No. 20111081102.
- Bertin & Arnouts (1996) Bertin, E., & Arnouts, S. 1996, Astronomy and Astrophysics Supplement Series, 117, 393
- Buades et al. (2005) Buades, A., Coll, B., & Morel, J.-M. 2005, in Computer Vision and Pattern Recognition, 2005. CVPR 2005. IEEE Computer Society Conference on, vol. 2, 60–65 (IEEE)
- Buades et al. (2008) Buades, A., Coll, B., & Morel, J.-M. 2008, International Journal of Computer Vision, 76, 123
- Högbom (1974) Högbom, J. A. 1974, \aaps, 15, 417
- Huo & Zhou (2013) Huo, Z.-X., & Zhou, J.-F. 2013, Research in Astronomy and Astrophysics, 13, 991
- Jin et al. (2010) Jin, J., Chen, Y., Zhang, S.-N., et al. 2010, Chinese Physics C, 34, 66
- Jones et al. (2001–) Jones, E., Oliphant, T., Peterson, P., et al. 2001–, SciPy: Open source scientific tools for Python
- Jourdain & Roques (2009) Jourdain, E., & Roques, J. 2009, The Astrophysical Journal, 704, 17
- Li et al. (2009) Li, G., Wu, M., Zhang, S., & Jin, Y.-K. 2009, Chinese Astronomy and Astrophysics, 33, 333
- Li (2003) Li, T.-P. 2003, in High Energy Processes and Phenomena in Astrophysics, IAU Symposium, vol. 214, edited by X. D. Li, V. Trimble, & Z. R. Wang, 70
- Li (2007) Li, T.-P. 2007, Nuclear Physics B - Proceedings Supplements, 166, 131, proceedings of the Third International Conference on Particle and Fundamental Physics in Space
- Li & Wu (1994) Li, T.-P., & Wu, M. 1994, Astrophysics and Space Science, 215, 213
- Lu (2012) Lu, F. 2012, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, vol. 8443
- Massaro et al. (2000) Massaro, E., Cusumano, G., Litterio, M., & Mineo, T. 2000, \aap, 361, 695
- Richardson (1972) Richardson, W. H. 1972, JOSA, 62, 55
- Starck & Murtagh (2006) Starck, J.-L., & Murtagh, F. 2006, Astronomical Image and Data Analysis (Berlin: Springer-Verlag)