# Analysis methods to localize and characterize X-ray sources with the Micro-
channel X-ray Telescope on board the SVOM satellite

Shaymaa Hussein  Université Paris-Saclay, CNRS/IN2P3, IJCLab, 91405 Orsay,
France  An-Najah National University, Nablus, Palestine  Florent Robinet
Université Paris-Saclay, CNRS/IN2P3, IJCLab, 91405 Orsay, France  Martin
Boutelier  Centre National d’Etudes Spatiales, Centre spatial de Toulouse, 18
avenue Edouard Belin, 31401 Toulouse Cedex 9, France  Diego Götz  Université
Paris-Saclay, Université Paris Cité, CEA, CNRS, AIM, 91191, Gif-sur-Yvette,
France  Aleksandra Gros  Université Paris-Saclay, Université Paris Cité, CEA,
CNRS, AIM, 91191, Gif-sur-Yvette, France  Benjamin Schneider  Université
Paris-Saclay, Université Paris Cité, CEA, CNRS, AIM, 91191, Gif-sur-Yvette,
France

###### Abstract

SVOM is a Sino-French space mission targeting high-energy transient
astrophysical objects such as gamma-ray bursts. The soft X-ray part of the
spectrum is covered by the Micro-channel X-ray Telescope (MXT) which is a
narrow field telescope designed to precisely localize X-ray sources. This
paper presents the method implemented on-board to characterize and localize
X-ray sources with the MXT. A specific localization method was developed to
accommodate the optical system of the MXT, which is based on “Lobster-Eye”
grazing incidence micro-pore optics. For the first time, the algorithm takes
advantage of cross-correlation techniques to achieve a localization accuracy
down to 2 arcmin with less than $200$ photons, which guarantees a rapid
follow-up for most of the gamma-ray bursts that SVOM will observe. In this
paper, we also study the limitations of the algorithm and characterize its
performance.

##  1 Introduction

The SVOM (Space-based multi-band astronomical Variable Object Monitor) is a
space mission (Wei, 2016) coordinated by the Chinese (CNSA) and French (CNES)
space agencies. The core program of SVOM is to detect and characterize gamma
ray bursts (Zhang, 2018), which are intense flashes of gamma rays lasting from
a fraction to a few hundreds of seconds. Additionnally, SVOM will also target
other astrophysical objects: active galactic nuclei (Beckmann & Shrader,
2012), black-hole binaries (McClintock & Remillard, 2006), galactic X-ray
binaries and magnetars (Revnivtsev et al., 2008), flaring stars (Kenyon &
Hartmann, 1987), or the cataclysmic variables (Robinson, 1976). Finally, SVOM
will participate to multi-messenger programs where transient sources observed
by other experiments (gravitational-wave detectors (Abbott et al., 2018) or
neutrino detectors) will be followed-up (Zhang, 2018).

The SVOM satellite will be launched in late 2023. It is equipped with four
instruments. At high energies, the ECLAIRs (Godet et al., 2014) telescope
($>4$ keV) and the Gamma-Ray Monitor (Dong et al., 2010) ($>30$ keV) will
trigger on gamma-ray sources. These sources will then be followed-up by the
Microchannel X-ray Telescope (MXT) (Götz et al., 2015) and the Visible
Telescope (VT) (Wang et al., 2010). The MXT is designed to observe X-ray
sources from 0.2 to 10 keV. For gamma-ray bursts, this energy range covers the
afterglow emission (Costa et al., 1997) following the gamma-ray burst prompt
emission. The original trigger is either detected on-board by the
ECLAIRs/Gamma-Ray Monitor or by another experiment. The SVOM platform then
slews to place the astrophysical object in the MXT $58\times 58$ arcmin2 field
of view. When stable, the telescope cumulates X-ray photons. The MXT images
are analyzed on-board to localize the X-ray source. The reconstructed position
is updated every 30 seconds and transmitted to the ground and broadcasted
world-wide to enable a follow-up of the event with ground telescopes. It will
also be used to trigger a second slew and better position the source for an
optical follow-up with the VT. One of the MXT scientific requirements is to
detect and localize 85% of gamma-ray bursts with an accuracy better than 2
arcmin after 10 minutes of stable observation.

The MXT is a light-weight ($<42$ kg) and compact (focal length $\sim 1.15$ m)
X-ray telescope. Its unique optical system is based on a “Lobster-Eye” grazing
incidence technique (Feldman et al., 2020). High-energy photons are reflected
and focused by a collection of square micropores creating a unique pattern on
the focal plane, composed of a central peak and two cross arms. The X-ray
source sky location is derived from the position of the central peak in the
focal plane.

This paper presents the algorithm developed to perform the on-board
localization analysis. A cross-correlation technique is used to best resolve
sources of X-ray photons, even with a weak intensity. This method is applied
to search for multiple sources in the MXT field of view using an iterative
subtraction process. For each identified source, the number of photons is
evaluated using the cross-correlation data, as well as the overall background
counts. This is a new approach which was specifically designed to manage the
peculiar shape of the point spread function of the telescope. Moreover, due to
space constraints, the data-processing methods had to be optimized to cope
with low on-board computing resources. This method is characterized and we
demonstrate that the scientific requirements of the MXT are achieved by the
on-board software.

This paper is organized as follows. Section 2 introduces the two data products
needed by the localization software: the camera images and the telescope point
spread function. Section 3 presents the on-board analysis steps: the source
localization, the estimation of signal and background counts, and the method
to identify multiple X-ray sources. This method is characterized in Sec. 4.
Finally, a summary and concluding remarks are given in Sec. 5.

##  2 Analysis input

###  2.1 Camera images

The MXT camera is a pn-Charge Coupled Device (pnCCD) (Schneider et al.,
submitted to Experimental Astronomy; Ceraudo et al., 2020) installed in the
telescope $yz$ focal plane. The $x$ axis is the optical axis of the telescope.
This is represented schematically in the left-hand side of Fig. 1. The camera
pixels are indexed by $i$ and $j$, running from $0$ to $N-1=255$ along the $y$
(resp. $z$) axis. We also use intrinsic coordinates, $0\leq y<1$ and $0\leq
z<1$, to measure continuous positions in the camera plane.

An incoming X-ray photon interacts with the MXT camera pixels and deposits its
energy via photo-electric effect. This energy is distributed over one or a few
contiguous pixels, forming a pattern. The most probable pixel patterns induced
by a photon and after thresholding on the pixel amplitude are represented in
the right-hand side diagram in Fig. 1 (Godet et al., 2009).

Figure 1: An incoming photon hits the camera plane (left). From the hit
position $(y,z)$, the spherical angle $(\theta,\varphi)$ of the incoming
direction can be derived. The left-hand diagram shows the most probable pixel
patterns induced by a photon hitting the camera plane.

The MXT camera is operated in two modes. In full-frame mode, the camera front-
end electronics capture an image every 200 ms and transfer it to the MXT data
processing unit. The full-frame images are processed to estimate the pixel
noise and derive pixel-by-pixel detection thresholds. For astrophysical
observations, the camera is set to event mode. The image is integrated over
100 ms and the pixel noise component is removed by the front-end electronics
(Ceraudo et al., 2020). When the amplitude is above the detection threshold,
the pixel data is transfered to the MXT data processing unit for analysis.

In event mode, the first analysis step is to cluster adjacent pixels together.
A photon is identified if the pixel pattern matches one of those represented
in the right-hand side of Fig. 1. The other pixel patterns are attributed to
noise (e.g. cosmic rays) and are discarded. Photons are parametrized by the
energy $E_{p}=\sum_{ij}{E[i][j]}$ integrated over the $(i,j)$ pixels in the
cluster. The photon hit position $(y_{p},z_{p})$ in intrinsic coordinates is
estimated by computing the energy-weighted pixel position:

$\displaystyle y_{p}$
$\displaystyle=\frac{\sum_{ij}E[i][j]\times(i+0.5)}{N\sum_{ij}E[i][j]}\quad\mathrm{and}$
(1) $\displaystyle z_{p}$
$\displaystyle=\frac{\sum_{ij}E[i][j]\times(j+0.5)}{N\sum_{ij}E[i][j]},$ (2)

where the sums run over the pixels $(i,j)$ in the photon cluster.

As images are recorded, photons are cumulated onto a $N_{D}\times N_{D}$
counting map $D$ called the “photon map”. To limit the computational cost, the
resolution is set to $N_{D}=N/2=128$. It is shown in Sec. 4.2 that the
localization accuracy is hardly impacted by this choice of lower resolution.
The right-hand plot of Fig. 2 shows an example of a photon map.

###  2.2 Point spread function

The localization analysis presented in Sec. 3 takes the photon map $D$ as an
input and cross-correlates it with the image of a point source produced by the
MXT optical system. This expected image, called the point spread function, is
measured and modeled on the ground and uploaded to the satellite. The point
spread function of the MXT is the result of the microporous structure (Feldman
et al., 2020) of its optical system. It measures the hit position probability
of an incoming photon over the focal plane. It is composed of a central spot
(double reflection), four cross arms (single reflection) and a diffused patch
(no reflection). The MXT point spread function has been measured and modeled
by the University of Leicester team, in charge of the MXT optical system. It
is modeled by a function $P(y,z,E)$ fitted to X-ray data collected at
different energies $E$ during the MXT performance tests conducted at the
Panter facility, Neuried Germany, in November 2021. The flight-model point
spread function used in this work, $P(y,z)$, is plotted in Fig. 2. It is
weighted over the expected number of photons for an average gamma-ray burst,
taking into account a typical spectrum and the MXT photon collecting
efficiency as a function of energy. The full width at half maximum of the
central peak is about 11 arcmin ($\simeq 50$ pixels).

Figure 2: The MXT point spread function is represented on the left. An example
of a resulting photon map $D$ with 1000 photons is represented on the right.

The point spread function $P$ is discretized over a $N_{P}\times N_{P}$ grid
where $N_{P}=N=2N_{D}$. It is twice the size of the photon map in both
directions to cover all possible incoming photon directions. Moreover, the
point spread function is normalized to 1:

$\sum_{i=0}^{N_{P}-1}\sum_{j=0}^{N_{P}-1}{P[i][j]}=1.$ (3)

The right-hand side plot in Fig. 2 shows an example of a photon map generated
with 1000 photons, placing the point spread function peak at $y=0.3$ and
$z=0.3$.

##  3 Localization algorithm

The localization analysis loops continuously in the MXT data processing unit.
One cycle takes less than 2 s to complete. The analysis delivers the X-ray
source sky location in spherical coordinates $(\theta,\varphi)$, as well as a
photon count, a signal-to-noise ratio and a localization uncertainty. This
information is updated every $~{}30$ s and packaged in telemetry packets sent
to the ground. The photon map $D$, introduced in Sec. 2.1, is cross-correlated
with the point spread function presented in Sec. 2.2. As opposed to a
simplistic barycenter/centroid approach, this method, described in Sec. 3.1,
was selected as it offers a better sensitivity to faint sources. Indeed, one
can take advantage of the specific MXT’s point spread function for which the
cross-arm pattern can guide the search for source peak positions.

The cross-correlation data is then processed to search for local maxima, as
detailed in Sec. 3.2. In Sec. 3.3, we explain how to estimate the source
photon count and the background contribution. Finally, the algorithm is
designed to locate the three brightest sources in the MXT field of view. To
achieve this, an iterative subtraction process is implemented and is described
in Sec. 3.4.

###  3.1 Cross-correlation

To determine the localization of the source, the photon map $D$ is cross-
correlated with the point spread function $P$ such as:

$C[i][j]=\sum_{i^{\prime}=0}^{N_{P}-1}\sum_{j^{\prime}=0}^{N_{P}-1}D[i+i^{\prime}][j+j^{\prime}]\times
P[i^{\prime}][j^{\prime}],$ (4)

where $i$ and $j$ take values between 0 and $N_{D}-1$. In Eq. 4 the photon map
is shifted by $(i,j)$ bins over the point spread function matrix and the
cross-correlation coefficients are obtained by summing the products of all
overlapping bin values in a $N_{P}^{2}$ grid. It is worth noting that the
photon map indices may be out-of-range. As we are dealing with finite discrete
functions, we assume circular periodicity where $i+i^{\prime}$ (resp.
$j+j^{\prime}$) actually means $(i+i^{\prime})\pmod{N_{D}}$ (resp.
$(j+j^{\prime})\pmod{N_{D}}$). In the following the out-of-range issues will
no longer be addressed, and we will simply write $i+i^{\prime}$ (resp.
$j+j^{\prime}$). As a result, in the sum of Eq. 4 the photon map $D$ is used
four times.

The implementation of the cross-correlation given in Eq. 4 is computationally
expensive. Alternatively, the cross-correlation can be computed in the Fourier
space:

$\tilde{C}[k][l]=\frac{1}{N_{D}^{2}}\sum_{i=0}^{(N_{D}-1)}\sum_{j=0}^{(N_{D}-1)}C[i][j]\times
e^{-2\sqrt{-1}\pi(ik+jl)/N_{D}}.$ (5)

According to the correlation theorem (Papoulis, 1962), the cross-correlation
of two signals is equivalent to a complex conjugate multiplication of their
Fourier transforms:

$\tilde{C}[k][l]=N_{P}^{2}\times\tilde{D}[k][l]\times\tilde{P}^{*}[2k][2l],$
(6)

where we use ∗ for the complex conjugate. In Eq. 6, both the photon map
$\tilde{D}$ and the cross-correlation matrix $\tilde{C}$ are sampled at a
frequency $1/N_{D}$, while the point spread function matrix is sampled at a
frequency $1/N_{P}=1/(2N_{D})$:

$\tilde{P}[k][l]=\frac{1}{N_{P}^{2}}\sum_{i=0}^{(N_{P}-1)}\sum_{j=0}^{(N_{P}-1)}P[i][j]\times
e^{-2\sqrt{-1}\pi(ik+jl)/N_{P}}.$ (7)

As a consequence, every other coefficients of $\tilde{P}$ are considered in
Eq. 6.

Because of the limited computing power on-board and to optimize the uplink
bandwidth, the point spread function is formatted on the ground and the final
$N_{P}^{2}\times\tilde{P}^{*}[2k][2l]$ table is uploaded onboard with a
telecommand. Having one Fourier coefficient over two, the on-board point
spread function is therefore partial. This limitation has important
consequences to compute the signal and background counts: this is discussed in
Sec. 3.3. The cross-correlation map is obtained by applying an inverse Fourier
Transform to Eq. 6:

$C[i][j]=N^{2}_{P}\times\sum_{k=0}^{(N_{D}-1)}\sum_{l=0}^{(N_{D}-1)}\tilde{D}[k][l]\times\tilde{P}^{*}[2k][2l]\times
e^{+2\sqrt{-1}\pi(ik+jl)/N_{D}}.$ (8)

This cross-correlation method has two systematic biases which must be
corrected. The first bias is introduced when working with a point spread
function twice the size of the photon map. The photon map must be shifted by
half of a bin in both directions to be aligned with the point spread function.
Therefore the cross-correlation map is also shifted by half of a bin. When
localizing a peak in the cross-correlation map, an offset of $0.5/N_{D}$ must
be removed in $y$ and $z$. The second bias is caused by the Fourier transform
spectral leakage: it is more pronounced when the source peak stands near the
edges of the camera plane. This bias is discussed and corrected in Sec. 4.3

###  3.2 Source angular position

The X-ray source position is associated to a peak in the cross-correlation map
$C$. First, the global maximum $C[i_{1}][j_{1}]$ in the cross-correlation map
is identified. Then we define a window centered on this maximum:
$i_{1}-N_{w}\leq i\leq i_{1}+N_{w}$ and $j_{1}-N_{w}\leq j\leq j_{1}+N_{w}$.
As explained in Sec. 3.1, we use circular indices if the window overlaps the
edges of the cross-correlation map. The size of the window must be chosen to
fully include the central spot of the point spread function. For the MXT point
spread function, we use $N_{w}=25$. The peak position $(y_{1},z_{1})$ is
finally computed as a two-dimensional barycenter inside the window:

$y_{1}=\frac{\sum_{i,j}C[i][j]\times(i+0.5)}{N_{D}\times\sum_{i,j}C[i][j]},\quad\mathrm{and}\quad
z_{1}=\frac{\sum_{i,j}C[i][j]\times(j+0.5)}{N_{D}\times\sum_{i,j}C[i,j]}.$ (9)

The sums run over indices $(i,j)$ inside the window and where the cross-
correlation coefficients take significant values: $C[i][j]>\alpha\times
C[i_{1}][j_{1}]$. Using simulated data, we find that $\alpha=0.9$ gives the
best localization accuracy. It is worth noting that both the $\alpha$ and
$N_{w}$ parameters are configurable from the ground.

Figure 3 shows an example of a photon map simulated with a faint source (50
photons) on top of a uniform background (600 counts). The source position
cannot be identified in the photon map while it appears clearly in the cross-
correlation map. The peak position is localized reasonably well with the
method described above (circular marker). The true position is indicated with
a triangular marker. For comparison, we also add the position derived from a
simple barycenter evaluated over the entire photon map (square marker) which
fails at localizing the source.

Figure 3: Example of a photon map $D$ simulated with a flat background of 600
counts and a faint source positioned at $y=0.7$ and $z=0.7$ (triangle marker)
with 50 counts (right plot). The source position is well visible in the cross-
correlation map (right). It is localized at $y=0.73$ and $z=0.72$ by the MXT
on-board analysis (circular marker). A simple barycenter fails at localizing
the source (square marker).

As introduced in Sec. 3.1, the cross-correlation map suffers from two
systematical biases: the Fourier transform spectral leakage and the alignment
asymmetry between the photon map and the point spread function. For the first
effect we apply ad hoc corrections $\eta_{1}(y)$ and $\eta_{2}(z)$ which are
evaluated and discussed in Sec. 4.3. The second bias is corrected by applying
a shift of half of a bin:

$y^{\prime}_{1}=y_{1}+\eta_{1}(y_{1})-\frac{0.5}{N_{D}},\quad\mathrm{and}\quad
z^{\prime}_{1}=z_{1}+\eta_{2}(z_{1})-\frac{0.5}{N_{D}}.$ (10)

Finally, the source position ($y^{\prime}_{1},z^{\prime}_{1}$) is converted to
spherical angles $(\theta,\varphi)$ using the telescope focal length $F$ and
the physical size of the camera $L\times L$:

$\displaystyle\theta$
$\displaystyle=\tan^{-1}{\left(\frac{L}{F}\times\sqrt{(y^{\prime}_{1}-0.5)^{2}+(z^{\prime}_{1}-0.5)^{2}}\right)},$
$\displaystyle\varphi$
$\displaystyle=\tan^{-1}\left(\frac{z^{\prime}_{1}-0.5}{y^{\prime}_{1}-0.5}\right).$
(11)

###  3.3 Signal and noise counts

The source signal $S_{1}$ and noise $B$ counts are estimated on board the MXT
to derive the source signal-to-noise ratio $\rho_{1}=S_{1}/\sqrt{B}$. A
standard method to estimate these quantities consists of analyzing the photon
map. The signal+noise component is integrated around the main peak and the
noise-only component is estimated using a region of the image where the signal
contribution can be neglected. For the MXT, this method is not optimal due to
the cross-shaped structure of the point spread function and the possible
presence of multiple X-ray sources in the field of view. Moreover there are
large uncertainties associated to faint signals and low-statistic backgrounds.

Instead, we use the cross-correlation map $C$ which naturally integrates the
background and signal counts over the entire camera plane. Decomposing the
photon map $D$ in a linear combination of a signal component $s_{1}$ and a
background component, we can write Eq. 4 as:

$\displaystyle C[i][j]$
$\displaystyle=\sum_{i^{\prime}=0}^{N_{P}-1}\sum_{j^{\prime}=0}^{N_{P}-1}\left(s_{1}[i^{\prime}][j^{\prime}]\times
P[i^{\prime}+i][j^{\prime}+j]\right)+\frac{B_{1}}{N_{D}^{2}},$ (12)

where we assume the background $B_{1}$ to be uniformly distribued over the
camera plane and where we used the normalization of Eq. 3. The source
distribution is driven by the point spread function:

$\displaystyle s_{1}[i][j]$
$\displaystyle=\frac{S_{1}}{\sum_{i^{\prime}=0}^{N_{D}-1}\sum_{j^{\prime}=0}^{N_{D}-1}P[N_{D}+i^{\prime}-i_{1}][N_{D}+j^{\prime}-j_{1}]}\times
P[N_{D}+i-i_{1}][N_{D}+j-j_{1}],$
$\displaystyle=\frac{S_{1}}{\beta[i_{1}][j_{1}]}\times
P[N_{D}+i-i_{1}][N_{D}+j-j_{1}],$ (13)

where $i_{1}=\lfloor y_{1}^{\prime}N_{D}\rfloor$ and $j_{1}=\lfloor
z_{1}^{\prime}N_{D}\rfloor$ are the indices matching the cross-correlation
peak position $(y_{1}^{\prime},z_{1}^{\prime})$ obtained in Eq. 10 and where
we introduce:

$\beta[i_{1}][j_{1}]=\sum_{i^{\prime}=0}^{N_{D}-1}\sum_{j^{\prime}=0}^{N_{D}-1}P[N_{D}+i^{\prime}-i_{1}][N_{D}+j^{\prime}-j_{1}].$
(14)

The $\beta$ function cannot be computed exactly on board as the point spread
function is incomplete; every other Fourier coefficients
$\tilde{P}^{*}[2k][2l]$ are uploaded. Instead, the $\beta$ function is
approximated by a fit function, the parameters of which are uploaded: see the
left plot in Fig. 4.

Figure 4: For a given point spread function, the $\beta$ function (left) is
measured using Eq. 14. It is then fitted: the fit is represented by the red
contour. The function $C_{\varepsilon}$ (right) is measured with $B=0$ (Eq.
18) and a source positioned at the center.

Similarly, the signal component cannot be calculated completely so we re-write
Eq 13 as:

$s_{1}[i][j]=\frac{S_{1}}{\beta[i_{1}][j_{1}]}\times\left(P_{\mathrm{comp}}[i-i_{1}][j-j_{1}]-\varepsilon_{1}[i][j]\right),$
(15)

where we introduce a “compact” point spread function of size $N_{D}\times
N_{D}$:

$P_{\mathrm{comp}}[i][j]=P[i][j]+P[i][N_{D}+j]+P[N_{D}+i][j]+P[N_{D}+i][N_{D}+j],$
(16)

the Fourier transform of which is:

$\tilde{P}_{\mathrm{comp}}[k][l]=4\times\tilde{P}[2k][2l].$ (17)

The $\varepsilon_{1}$ function describes the missing terms of the point spread
function in Eq. 13. The cross-correlation of the signal component with the
point spread function is:

$\displaystyle C_{1}[i][j]$
$\displaystyle=\frac{S_{1}}{\beta[i_{1}][j_{1}]}\times\sum_{i^{\prime}=0}^{N_{P}-1}\sum_{j^{\prime}=0}^{N_{P}-1}P[i^{\prime}+i][j^{\prime}+j]\times\left(P_{\mathrm{comp}}[i^{\prime}-i_{1}][j^{\prime}-j_{1}]-\varepsilon_{1}[i^{\prime}][j^{\prime}]\right)$
$\displaystyle=\frac{S_{1}}{\beta[i_{1}][j_{1}]}\times(W[i-i_{1}][j-j_{1}]-C_{\varepsilon_{1}}[i][j]).$
(18)

This cross-correlation coefficient has two components. The first one, $W$, is
the cross-correlation of $P_{\mathrm{comp}}$ with the point spread function.
The $W$ function is computed on-board substituting $D$ by $P_{\mathrm{comp}}$
in Eq. 8:

$W[i][j]=4N^{2}_{P}\times\sum_{k=0}^{(N_{D}-1)}\sum_{l=0}^{(N_{D}-1)}\tilde{P}[2k][2l]\times\tilde{P}^{*}[2k][2l]\times
e^{+2\sqrt{-1}\pi(ik+jl)/N_{D}}.$ (19)

The second component, $C_{\varepsilon_{1}}$, is the cross correlation between
$\varepsilon_{1}$ and the point spread function. It is unknown.

Let us consider two points in the cross-correlation map $C$: the cross-
correlation peak at $(i_{1},j_{1})$ and the cross-correlation minimum point at
$(i_{m},j_{m})$. We get:

$\begin{cases}C[i_{1}][i_{1}]=\frac{S_{1}}{\beta[i_{1}][j_{1}]}\times(W[0][0]-C_{\varepsilon_{1}}[i_{1}][j_{1}])+\frac{B_{1}}{N_{D}^{2}}\\\
C[i_{m}][i_{m}]=\frac{S_{1}}{\beta[i_{1}][j_{1}]}\times(W[i_{m}-i_{1}][j_{m}-j_{1}]-C_{\varepsilon_{1}}[i_{m}][j_{m}])+\frac{B_{1}}{N_{D}^{2}}\end{cases}$
(20)

Solving this system of equations we get:

$\displaystyle S_{1}$
$\displaystyle=\beta[i_{1}][j_{1}]\times\frac{C[i_{1}][j_{1}]-C[i_{m}][j_{m}]}{W[0][0]-C_{\varepsilon_{1}}[i_{1}][j_{1}]-W[i_{m}-i_{1}][j_{m}-j_{1}]+C_{\varepsilon_{1}}[i_{m}][j_{m}]},$
(21) $\displaystyle B_{1}$
$\displaystyle=N_{D}^{2}\times\left(C[i_{1}][j_{1}]-\frac{S_{1}}{\beta[i_{1}][j_{1}]}\times(W[0][0]-C_{\varepsilon_{1}}[i_{1}][j_{1}])\right).$
(22)

To estimate $S_{1}$ and $B_{1}$ with Eqs. 21 and 22, we approximate the
$C_{\varepsilon_{1}}$ values by:

$\displaystyle
C_{\varepsilon_{1}}[i_{m}][j_{m}]-C_{\varepsilon_{1}}[i_{1}][j_{1}]$
$\displaystyle=\mathrm{Min}_{i,j}(C_{\varepsilon}[i][j])-\mathrm{Max}_{i,j}(C_{\varepsilon}[i][j])\quad\mathrm{and}$
(23) $\displaystyle C_{\varepsilon_{1}}[i_{1}][j_{1}]$
$\displaystyle=\mathrm{Max}_{i,j}(C_{\varepsilon}[i][j]),$ (24)

where $C_{\varepsilon}$ is a cross-correlation function computed on the ground
with a source positioned at the center of the camera plane: see the right plot
in Fig. 4.

###  3.4 Localization of multiple sources

Given the detector sensitivity (Götz et al., 2015), it is possible to have
more than one X-ray source within the field of view of the MXT. Using data
collected by the ROSAT telescope (Boller et al., 2016) which was operated in
an energy band comparable to the MXT and had a similar sensitivity, we expect
to have at most three X-ray sources within the field of view of the MXT. As a
result, the localization algorithm shall be able to manage multiple sources.
The localization analysis described in Secs. 3.1, 3.2, and 3.3 is repeated
three times, subtracting the contribution of the detected source at each
iteration.

The cross-correlation map at iteration $g=1$ is $C^{g=1}=C$, where $C$ is
computed with Eq. 8. For the next two iterations, $g=2$ and $g=3$, we use
$C^{g}=C^{g-1}-C_{g-1}$, where the contribution of source $g$, $C_{g}$, is
derived from Eq. 18. Note that the $C_{\varepsilon_{g}}$ function is
approximated by a single value:
$C_{\varepsilon_{g}}[i][j]\simeq\mathrm{Max}_{i,j}(C_{\varepsilon}[i][j])$.

After three iterations, we have three estimates for the background: $B_{1}$,
$B_{2}$, and $B_{3}$. None of them is a perfect estimate of $B$. For a single
source, the best estimator is $B_{1}$. For multiple sources, however, the
choice is not trivial: the first estimate $B_{1}$ is over-estimated by
secondary sources, while the last estimate $B_{3}$ is under-estimated by
imperfect subtractions. As a compromise, the background is updated at each
iteration:

$B=B_{g}\quad\mathrm{if}\quad S_{g}>\rho_{q}\times\sqrt{B_{g}},$ (25)

where $\rho_{q}$ is a configurable parameter.

Figure 5 illustrates the subtraction method. The photon map $D$ includes three
photon sources simulated with $S_{1}=1000$, $S_{2}=300$, and $S_{3}=300$. The
background is uniform and set to $B=600$. Running our analysis, these counts
are estimated by $S_{1}=1045$ (+4.5%), $S_{2}=386$ (+29%), $S_{3}=299$
(+0.0%), $B_{1}=824$, $B_{2}=810$, and $B_{3}=627$. The final background is
estimated at $B=B_{3}=627$ (+4.5%).

Figure 5: The top-left plot shows the distribution of photons from three
separate sources ($S_{1}=1000$, $S_{2}=300$, and $S_{3}=300$) on top of a
uniform background ($B=600$). The three other plots show the initial cross-
correlation map $C^{1}$ (top-right), after subtracting the first source
$C^{2}$ (bottom-left), and after subtracting the second source $C^{3}$
(bottom-right). The true source positions are indicated by the circles.

##  4 Characterization of the localization method

The localization algorithm presented in Sec. 3 is characterized using
different test scenarios. For each scenario we generate camera images with the
SatAndLight simulation toolkit (Robinet & Hussein, 2021), which is designed to
include many effects: photons from astrophysical and background sources,
cosmic rays, corrupted camera pixels, etc. Astrophysical sources generate
photons at a given energy, or following a given energy distribution. Photons
are reflected by the MXT optical system following the point spread function
(see Sec. 2.2) and interact with the camera CCD pixels. A cosmic X-ray
background is simulated by injecting high-energy photons uniformly on the
camera plane. For the MXT, we expect to detect one background photon every
second. In the following we will often consider a canonical background of
$B=600$ counts, corresponding to an observation of 10 minutes.

###  4.1 The source intensity

In this section, the localization accuracy and its associated uncertainty is
evaluated for different source and noise intensities. We use the $r_{90}$
quantity,

$r_{90}=\left[\angle(\vec{r}_{\mathrm{meas}},\vec{r}_{\mathrm{true}})\right]_{90\%},$
(26)

which measures the angular distance between the source position recontructed
by the localization algorithm $\vec{r}_{\mathrm{meas}}$ (Eq. 10) and the true
source position defined in the simulation $\vec{r}_{\mathrm{true}}$. This
angular distance is measured with 1000 simulations of different positions
$\vec{r}_{\mathrm{true}}$ randomly drawn in the field of view of the
telescope. The $r_{90}$ quantity is obtained from the 90th percentile of the
distribution of $\angle(\vec{r}_{\mathrm{meas}},\vec{r}_{\mathrm{true}})$.

The MXT was designed to achieve a source localization accuracy of $r_{90}<2$
arcmin for 85% of the detected gamma-ray bursts. In Fig. 6, $r_{90}$ is
evaluated for an X-ray source with an intensity ranging from 10 to 10000
photons. The source photons are cumulated on top of a background noise of
different intensities: 30, 300, and 600 counts. For these three background
levels, the localization uncertainty is below 2 arcmin for 120, 164, and 197
source photons respectively.

Figure 6: Localization uncertainty $r_{90}$ evaluated for a source intensity
ranging from 10 to 10000 photons and with different background levels. Each
point in this graph results from 1000 source positions randomly drawn in a
$17\times 17$ arcmin2 region around the camera center. Here, we use
$N_{D}=N=256$.

###  4.2 The photon map resolution

After photons are reconstructed, as explained in Sec. 2.1, they are cumulated
onto a $N_{D}\times N_{D}$ two-dimensional grid called the photon map. The
resolution of the photon map must be kept at minimum to limit the computing
cost associated to the Fourier transform used to cross-correlate the data with
the point spread function (see Sec. 3.1). However, it should be chosen high
enough to fully resolve the peak of the point spread function and maximize the
localization performance.

To drive this choice, we measure the impact of the photon map resolution on
the localization performance. The results are presented in Fig. 7. There is no
significant impact between the full resolution ($N_{D}=N=256$) and using a
$128\times 128$ grid. We observe a $\sim 8\%$ loss when the resolution is
further reduced ($64\times 64$). We only consider power-of-two values for
$N_{D}$ to perform optimal Fourier transforms. For the MXT localization
algorithm, and for the rest of this paper, we fix $N_{D}=128$. This resolution
is sufficiently high to completely sample the shape of the point spread
function and to offer an optimal localization accuracy. To perform the cross-
correlation, we developed a two-dimensional $128\times 128$ complex-to-complex
Fourier transform algorithm. It takes approximately 300 ms to complete one
Fourier transform with the on-board computer.

Figure 7: Localization uncertainty $r_{90}$ vs. number of photons evaluated
for various photon map resolutions: $64\times 64$ (light blue), $128\times
128$ (bright blue), and $256\times 256$ (dark blue). The background level is
fixed to 600 counts.

###  4.3 Detector edge bias

The cross-correlation between the photon map and the point spread function is
performed in the Fourier domain (Eq. 6). The discrete Fourier transform of
size $N_{D}$ assumes periodicity: $D[0]\simeq D[N_{D}-1]$ in both the $y$ and
$z$ directions. This assumption is reasonable when the X-ray source is
positioned around the center of the field of view. It is highly violated when
the source peak gets closer to the camera edges, when a fraction of the peak
is chopped off. In this situation, spectral-leakage effects can be large and
they bias the peak position in the cross-correlation map $C$.

To evaluate this bias along the $y$ direction, we move the source position
horizontally, from $(y=0,z=0.5)$ to $(y=1,z=0.5)$. Using discrete positions
along this axis, we measure the deviation between the measured and the true
positions: $dy=y_{\mathrm{true}}-y_{\mathrm{meas}}$. This deviation is
represented as a function of $y_{\mathrm{meas}}$ in Fig. 8. Usually, this bias
can be mitigated using padding or windowing methods. Padding methods would
increase data sizes and would lead to longer Fourier transform calculations.
Windowing the data would alter the signal and noise counts used to compute the
signal-to-noise ratio. Instead, the bias $dy(y_{\mathrm{meas}})$ is fitted by
an empirical function:

$\eta_{1}(y_{\mathrm{meas}})=\begin{cases}H_{0}\times e^{H_{1}\times
y_{\mathrm{meas}}}&\mathrm{if}\quad y_{\mathrm{meas}}>0.5\\\ H_{2}\times
e^{H_{3}\times y_{\mathrm{meas}}}&\mathrm{otherwise},\end{cases}$ (27)

where $H_{0}$, $H_{1}$, $H_{2}$, and $H_{3}$ are the fit parameters. The bias
in the $z$ direction is corrected using the same method. Indeed the edge
corrections in the $y$ and $z$ directions can be be addressed separately
except for a small region near the corners of the camera plane where the
corrections should not be decoupled. Given the low probability of finding a
source in these positions, we neglect this effect.

Figure 8: The edge bias $dy=y_{\mathrm{true}}-y_{\mathrm{meas}}$ is measured
in the $y$ direction. The red line shows the best fit using the function in
Eq. 27.

The $\eta_{1}$ and $\eta_{2}$ corrections are applied in Eq. 10 to obtain the
final source position. Figure 9 shows the localization uncertainty when the
source position moves across from the camera plane with and without the edge-
bias correction. The edge bias appears when the source is less than 30 pixels
away from the camera edge. After correction, the bias is still visible but it
is below the MXT design requirement of 2 arcmin. This comparison is done with
500 and 1000 photons showing that the edge correction method applies equally
to all source intensities.

Figure 9: Localization uncertainty $r_{90}$ when the source position moves
across the camera plane along the $y$ axis ($z=0.5$). The dashed lines do not
include the edge-bias correction while the continuous lines do. The background
count is fixed to 600 and two source intensities are tested: 1000 photons (red
curves) and 500 photons (blue curves).

###  4.4 Multiple sources

As explained in Sec. 3.4, three X-ray sources are identified and localized in
the field of view of the MXT. To characterize the algorithm with additional
sources, we consider the worst-case scenario where two X-ray sources share the
same $z$ coordinate; given the shape of the point spread function, the photon
patterns of the two sources overlap. We fix the brightest source at center of
the field of view and we move the second source from the center to the camera
edge in the $y$ direction. We also vary the intensity of the second source
relatively to the first source. The two sources are localized running the
algorithm described in Sec. 3.4 and we measure the localization bias,
$dr=\angle(\vec{r}_{\mathrm{meas}},\vec{r}_{\mathrm{true}})$, for each source.

Figure 10 shows the variation of $dr$ when the second source is moved away
from the first source. The two sources cannot be separated as long as their
angular distance is below 10 arcmin, which is approximately the size of the
point spread function central peak, as expected. When the separation between
the two sources is between 10 and 20 arcmin, the localization is biased by the
presence of another source. Above 20 arcmin, the interference between the two
sources can be neglected. We also test the impact of the relative intensity
between the two sources. When the second source is less intense by a factor
two, the localization performance for the first source is improved roughly by
the same factor. The localization of the second source also improves thanks to
a better subtraction of the first source.

Figure 10: Localization bias ($dr$) when two co-aligned X-ray sources are
present in the MXT field of view as a function of the angular distance between
the two sources ($r_{2}-r_{1}$). We have masked the region where
$r_{2}-r_{1}<10$ arcmin (in green) for which the peaks of the point spread
function are overlapping and cannot be separated. The intensity of the first
source is fixed to 1000 photons. Two intensities for the second source are
tested: 500 photons (left plot) and 250 photons (right plot).

###  4.5 Signal and noise counts

In Sec. 3.3, we developed an inovative but non-trivial method to estimate the
signal and background counts. To test the signal and background estimators, we
vary the source intensity ($S_{1}$) and the background level ($B$). Both
quantities are then estimated using Eqs. 21 and 22. The result of this study
is presented in Fig. 11. The signal count is systematically over-estimated
because of the simplification used to estimate the $C_{\varepsilon_{1}}$
matrix (see Sec. 3.3). For instance, when $S_{1}\simeq 200$, the signal is
over-estinated by $\sim 10\%$ when $B=600$. The same effect leads to a
systematic bias when estimating the background count: it is under-estimated
for low signal fluxes and over-estimated for large signal fluxes.

Figure 11: The signal (smooth curves) and background (curves with markers)
counts are estimated as a function of the signal intensity. For each point,
1000 source positions are randomly drawn in $0.3<y<0.7$ and $0.3<z<0.7$. Three
background levels are considered: 60 counts (blue), 600 counts (magenta) and
3000 (dark red). They are represented by horizontal dashed lines.

Figure 12 shows the resulting signal-to-noise ratio,
$\rho_{1}=S_{1}/\sqrt{B}$, estimated when fixing $B=600$. It is slightly over-
estimated: for $10<\rho_{1}<100$, the signal-to-noise bias does not exceed
20%.

Figure 12: Estimated signal-to-noise ratio as a function of the true signal-
to-noise ratio, fixing $B=600$. The diagonal is indicated by a dashed line.

The method to estimate $S$ and $B$ must be robust when several sources are
present in the MXT field of view. The particular example of Fig. 5, and
commented at the end of Sec. 3.4, shows that our estimators give satifactory
results. This method to estimate counts is iterative and is, therefore,
intrinsically biased by the presence of less intense sources of X-ray photons.
This bias is acceptable for on-board calculations. More advanced analyses can
be conducted offline using X-ray source catalogs and multiple-source fitting
techniques.

##  5 Summary

We presented the localization algorithm developed to process the images
captured by the Microchannel X-ray Telescope mounted on the SVOM satellite.
The images are analyzed on-board to rapidly detect and localize X-ray sources
in the telescope field of view. Given the Lobster-Eye optical design of the
telescope and the resulting cross-shaped structure of the point spread
function, we developed specific analysis methods. Cross-correlation techniques
were selected to maximize the sensitivity to faint sources. This paper shows
that a localization accuracy better than 2 arcmin can be achieved when the
telescope cumulates $\sim 150$ X-ray photons. This number guarantees that more
than 85% of gamma-ray bursts will be localized with an uncertainty better than
2 arcmin after 30 minutes (Götz et al., submitted to Experimental Astronomy).
Moreover, we developed a new and reliable method to estimate the signal and
background counts from the cross-correlation map from which we derive the
signal-to-noise ratio for multiple sources in the field of view.

The MXT localization algorithm relies on many parameters, some of which were
characterized and presented in this paper. All of them are configurable from
the ground. During the first months of the SVOM mission, these parameters will
be adjusted using real data collected in space. They will then be uploaded on
board to maximize the performance of the localization algorithm. These
parameters will be updated over time to follow aging effects.

The SVOM satelite will be launched in 2023 and will observe gamma-ray bursts
shortly after that. Moreover, SVOM will continue X-ray observations and step
over when other satellites, such as Swift (Gehrels, 2004) or Fermi (Michelson
et al., 2010), will terminate their missions. SVOM will play an important role
in the multi-messenger astronomy, in particular in coincidence with
gravitational-wave detectors. The localization algorithm described in this
paper will provide a rapid follow-up of gravitational-wave alerts to identify
an electromagnetic counterpart.

This work was developed in the context of the SVOM/MXT science group. We
gratefully acknowledge the support from the group and we value the fruitful
discussions that took place. In particular, we are grateful to the MXT team at
the University of Leicester for providing the flight-model point spread
function used in this paper, with special thanks to Richard “PSF wizard”
Willingale. We also thank the MXT team at the Max Planck Institute for
Extraterrestrial Physics for hosting the final performance tests of the MXT
conducted at the Panter facility, Neuried Germany, in November 2021. The PhD
thesis of Shaymaa Hussein was supported by the Centre National d’Etudes
Spatiales in France.

## References

  * Wei (2016) 2016, The Deep and Transient Universe in the SVOM Era: New Challenges and Opportunities - Scientific prospects of the SVOM mission.  https://arxiv.org/abs/1610.06892
  * Abbott et al. (2018) Abbott, B. P., et al. 2018, Living Rev. Rel., 21, 3, doi: 10.1007/s41114-020-00026-9
  * Beckmann & Shrader (2012) Beckmann, V., & Shrader, C. 2012, The Central Engine (John Wiley & Sons, Ltd), 35–88, doi: https://doi.org/10.1002/9783527666829.ch3
  * Boller et al. (2016) Boller, T., Freyberg, M. J., Truemper, J., et al. 2016, Astron. Astrophys., 588, 103, doi: 10.1051/0004-6361/201525648
  * Ceraudo et al. (2020) Ceraudo, F., Meuris, A., Doumayrou, E., et al. 2020, Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment, 973, 164164, doi: https://doi.org/10.1016/j.nima.2020.164164
  * Costa et al. (1997) Costa, E., et al. 1997, Nature, 387, 783, doi: 10.1038/42885
  * Dong et al. (2010) Dong, Y., Wu, B., Li, Y., Zhang, Y., & Zhang, S. 2010, Sci. China Phys. Mech. Astron., 53, 40, doi: 10.1007/s11433-010-0011-7
  * Feldman et al. (2020) Feldman, C. H., Willingale, R., Pearson, J., et al. 2020, in Space Telescopes and Instrumentation 2020: Ultraviolet to Gamma Ray, ed. J.-W. A. den Herder, S. Nikzad, & K. Nakazawa, Vol. 11444, International Society for Optics and Photonics (SPIE), 114441K, doi: 10.1117/12.2561739
  * Gehrels (2004) Gehrels, N. 2004, AIP Conf. Proc., 727, 637, doi: 10.1063/1.1810924
  * Godet et al. (2009) Godet, O., et al. 2009, Astron. Astrophys., 494, 775, doi: 10.1051/0004-6361:200811157
  * Godet et al. (2014) —. 2014, Proc. SPIE Int. Soc. Opt. Eng., 9144, 914424, doi: 10.1117/12.2055507
  * Götz et al. (2015) Götz, D., et al. 2015, PoS, SWIFT10, 074, doi: 10.22323/1.233.0074
  * Götz et al. (submitted to Experimental Astronomy) Götz, D., et al. submitted to Experimental Astronomy, Experimental Astronomy
  * Kenyon & Hartmann (1987) Kenyon, S. J., & Hartmann, L. 1987, Astrophys. J., 323, 714, doi: 10.1086/165866
  * McClintock & Remillard (2006) McClintock, J. E., & Remillard, R. A. 2006, Black hole binaries, ed. W. Lewin & M. van der Klis, Cambridge Astrophysics (Cambridge University Press), doi: 10.1017/CBO9780511536281
  * Michelson et al. (2010) Michelson, P. F., Atwood, W. B., & Ritz, S. 2010, Rept. Prog. Phys., 73, 074901, doi: 10.1088/0034-4885/73/7/074901
  * Papoulis (1962) Papoulis, A. 1962, The Fourier Integral and its Applications (New York: McGraw-Hill)
  * Revnivtsev et al. (2008) Revnivtsev, M., Sazonov, S., Krivonos, R., Ritter, H., & Sunyaev, R. 2008, Astron. Astrophys., 489, 1121, doi: 10.1051/0004-6361:200810213
  * Robinet & Hussein (2021) Robinet, F., & Hussein, S. 2021, SatAndLight: a simulation toolkit for space telescopes.  {https://svom.pages.in2p3.fr/SatAndLight/}
  * Robinson (1976) Robinson, E. L. 1976, Ann. Rev. Astron. Astrophys., 14, 119, doi: 10.1146/annurev.aa.14.090176.001003
  * Schneider et al. (submitted to Experimental Astronomy) Schneider, B., Renault-Tinacci, N., Götz, D., et al. submitted to Experimental Astronomy, Experimental Astronomy, 20
  * Wang et al. (2010) Wang, T., Qiu, Y., Cai, H., & Deng, J. 2010, Sci. China Phys. Mech. Astron., 53, 51, doi: 10.1007/s11433-010-0010-8
  * Zhang (2018) Zhang, B. 2018, The Physics of Gamma-Ray Bursts (Cambridge University Press), doi: 10.1017/9781139226530