Simulation and photoelectron track reconstruction of soft X-ray polarimeter

Xue-Feng Huang ID Hong-Bang Liu Jin Zhang Bo Huang Wen-Jin Xie Huan-Bo Feng Xi-Chen Cai Xi-Wen Liu Zi-Li Li Jian-Yu Gu Qian Liu ID Jin Li ID En-Wei Liang

Abstract The soft X-ray polarimeter (SXP) is a detector with a wide energy range, large area, and large field of view.A SXP will be mounted on the Chinese Space Station and will mainly focus on detecting the polarization of transient soft X-ray (2–10 keV) sources, especially gamma-ray bursts (GRBs). In this work, a polarimeter detector unit is taken as an example, and Geant4 and Garfield++ software are used to simulate the detection efficiency and track production. An improved track reconstruction algorithm is proposed and used to reconstruct two-dimensional images of the tracks. In this method, the initial emission angle of photoelectrons is reconstructed from the initial part of the track by shortening or extending the initial part of the track until the remaining track is straight, and the number of pixels is within an adjustable threshold. The modulation factor of the photoelectronic tracks after reconstruction reaches approximately 57% in the photon energy range of 7–10 keV.

Keywords X-ray · Polarimetry · Track reconstruction algorithm

1 Introduction

X-ray polarization offers a new window in high-energy astrophysics. Sensitive X-ray polarimetry is a promising method for revealing crucial and unique information about the structures and physical processes of black holes, neutron stars, and all classes of X-ray sources in the universe.The measurement of X-ray polarization was first successfully achieved in an X-ray polarimeter aboard the OSO-8 satellite for Bragg peak measurement[1]although only two separate energies (2.6 keV and 5.2 keV)were measured. In 2001,Costa used a high-resolution gas detector to measure X-ray polarization in the laboratory [2]. Since then, X-ray polarization measurements have become more widespread.Various particle physics and astrophysics missions for detecting the X-ray polarization, including IXPE [3] and eXTP [4], are being prepared. Other missions, such as XIPE [5]and PRAXyS [6], have been planned but subsequently canceled.Recently,CubeSat performed new X-ray polarization measurements of the Crab Nebula [7] using gas detector technology.In these missions,the polarization of energetic astrophysical sources was mainly measured through a narrow field of view using a high-sensitivity X-ray polarization detector, which cannot respond rapidly to transient sources such as gamma-ray bursts (GRBs).Based on the radiation mechanism, it is known that the magnetic field configuration of the GRB is reflected in the radiation polarization. With this in mind, it is hoped that the magnetic field configuration of the jet can be studied by measuring the polarization of the GRB radiation to understand this extreme explosion phenomenon [8]. The SXP is currently in the conceptual design stage in which a detector with a large area and a field of view is being developed.The detector design directly affects the detector sensitivity and needs to be optimized through simulations.The obtained data must be processed by a track reconstruction algorithm to obtain the polarization of the X-ray source.A good track reconstruction algorithm can improve the accuracy of the results.

Algorithms for track reconstruction are also constantly being developed. Bellazzini et al. first proposed an algorithm for track reconstruction[9,10].Analysis of the track moments allows some properties of their distribution,such as the barycenter,the major axis,and most importantly,the projection on the pixel plane of the point of absorption of the X-ray photons, to be analyzed. Differing from the traditional moment analysis method based on the track image as well as an algorithm based on the shortest path problem in graph theory [11] and adaptive cut methods [12], the algorithm proposed here is based on the detector condition and is called the iteration cut method.

A brief introduction to SXP and the structure of the SXP detector unit is given in Sect. 2. The SXP Monte Carlo simulations are presented in Sect. 3. The algorithm for track reconstruction is proposed in Sect. 4, and the conclusions are presented in Sect. 5.

2 Detector structure

The SXP is composed of multiple detector arrays. They give an effective field of view of approximately 90° and a detection area of 300 cm2for detecting the polarization of transient X-ray transient sources. The SXP contains 96 detection modules. We use a single detector module to illustrate the algorithm in this study. A single detector module is shown in Fig. 1. The SXP will be mounted on the Chinese space station to detect X-ray polarization within the energy range of 2–10 keV. The photoelectric effect has the largest reaction cross-section in the soft Xray energy range and results in the generation of photoelectrons with lower energy.Sufficiently distinct tracks can be obtained using thick-gas electron multiplier(THGEM)[13, 14] detection. THGEM has the advantages of high detection efficiency,simplicity,and stable operation.Track detection is the key to extracting the polarization of the Xray source. The Topmetal-II-chip is chosen for charge collection and reading out the electronic signals [15]. A field-programmable gate array (FPGA) is used to store information.

Fig.1 (Color online)Schematic diagram of a single detector module of the SXP,including the gas detector chamber,the THGEM,and the signal collection system

The internal structure of the detector is divided into three parts, namely the drift region, multiplier region, and induction region.As shown in Fig. 2,a photon incident in a direction perpendicular to the drift region of the detector is converted into a photoelectron through the photoelectric effect in the working gas. When the photoelectron passes through the working gas, it ionizes the gas, giving rise to electron-ion pairs.The electric field causes the electrons to drift toward the multiplier region,where they are amplified through avalanche generation. The photoelectron track is finally detected via two-dimensional imaging. The initial direction of the photoelectron is related to the polarization of the incident photons. The ejection probability is the highest when the initial direction has the same angle as the polarization. The polarization probability is given by the distribution of cos2θ [16], where θ is the polar angle with respect to the incident X-ray direction.

Fig. 2 (Color online) Photoelectric effect detection principle

The THGEMs are fabricated via standard PCB mechanical drilling and global etching processes, which allow for large-scale economic production of large-size detectors [17]. When two different negative voltages are applied to the extremities of the metal layer, electrons can drift through the etched hole and avalanche multiplication can occur, thus amplifying the charge for detection. The THGEM has several advantages, including discharge prevention,sub-millimeter position resolution,large gain,and ease of manufacture and maintenance. A thin THGEM[18–20]with 200 μm thickness,200 μm hole diameter,500 μm hole pitch,30 μm rim width,and 10 μm thickness was produced.

The Topmetal-II-is a direct charge sensor that collects drift charges directly without post-processing and reads them out into a 2D projection image. The pixel array is patterned in a matrix of 72×72 pixels with an 83 μm pitch, resulting in a charge-sensitive area of 6×6 mm2.Each pixel contains a low-noise charge-sensitive preamplifier (CSA) to establish an analog signal and a discriminator with a tunable threshold to generate hits. The Topmetal-II-is characterized by its small volume with the overall geometric size of 8 mm×9 mm×0.7 mm. The average equivalent noise charge measured is less than 15 e-. The low noise allows for improved trail detection.Furthermore, the sensor is resistant to radiation and remains functional even after 28 h of irradiation with a60Co source at the dose rate of 3.6 Gy h-1.

3 Monte Carlo simulations

The expected track was simulated using the Geant4[21],Garfield++ [22], and COMSOL Multiphysics [23] software packages. Geant4 is a toolkit for simulating particle–matter interactions. The detector model was constructed in Geant4 and used to simulate the initial track.The structure of the SXP is shown in Fig. 3. Here, photons enter the detector through a collimated aperture and are detected by the photoelectric effect in the gas. The effect of the THGEM is modeled in the following simulations.

The main criteria for selecting the gas are high detection efficiency and long track length. High detection efficiency is beneficial for the detection of weak sources, whereas long tracks are useful for reconstruction. As shown in Fig. 4,the detection efficiency and track length of different gases with the same thickness are simulated (panel(a) and panel (b)). It can be seen that the cross-section of the photoelectric effect decreases with increasing energy of the incident photons, thus reducing the detection efficiency.The efficiency at 2 keV is lower than that at 3 keV because of the Be window, which acts as a shield against lowenergy photons. Another factor worth considering is the gas diffusion coefficient. The diffusion coefficient of the gas should be as small as possible.Through experiments,it was found that the diffusion coefficients of the candidate gases were all within a certain range which did not affect the reconstruction algorithm.

For the above reasons, a mix of neon (80%) and dimethyl ether(DME 20%)was chosen as the detection gas.Nebased gas mixtures have longer photoelectron tracks in the lower energy bands while retaining reasonable efficiency,while the DME absorbs UV photons produced by charge propagation in the gas with small transversal diffusion and relatively high gain [24]. After completing the track simulation in Geant4,Garfield++and COMSOL Multiphysics were used to perform the track simulations in more detail.Garfield++ is a particle detector tool that is useful for detailed simulations in which a gas is used as a sensitive medium,while COMSOL Multiphysics is a multi-physical field simulation software based on finite element calculations. The gas components and the THGEM model of the detector were input into COMSOL Multiphysics to calculate the electric field distribution.The potential distribution calculated via COMSOL Multiphysics is shown in Fig. 5.The low-energy physics model in Geant4 was used to generate the secondary electrons and then reconstruct them into delta electrons in Garfield++ for simulations [25].The Topmetal-II-chip senses the charge signal through the exposed copper on its surface;it not only generates a signal at the incident location of each electron but also in the area before impact. To simulate the actual electronic signal collection, the charge value at each pixel position in the Garfield++simulation was redistributed between the pixel and the surrounding pixels with a Gaussian distribution,and the appropriate σ value was selected for comparison with the experiment to obtain more accurate track simulations.

Using Garfield++ simulations, which are intrinsically slow,a Gaussian distribution for the diffusion can be added to the slender electron track information produced by Geant4 to make the tracks thicker to simulate the signal sensed by the electrons on the Topmetal-II-chip. The width of the Garfield++track is first fitted with a Gaussian distribution to obtain the corresponding σ value. The fitted Gaussian distribution is then added to the Geant4 track to obtain a track width similar to that simulated with Garfield++. At the same time, different apertures in the THGEM directly affect the transverse diffusion of electrons, leading to different track thicknesses. The effect of different pore sizes on the track can be easily simulated by varying the σ value of the Gaussian distribution. Figure 6 shows an example of a track obtained through simulations.The size of each histogram bin is equivalent to the size of each pixel in the Topmetal-II-. Panel (a) shows the track used in the reconstruction.The size of the square represents the charge value at the coordinate position. The label numbers at the horizontal and vertical axes represent the detector readout positions, which coincide with the Topmetal-II-pixels.

Fig.3 (Color online)Detector structure used in Geant4.The labels in the figure denote (1) the 100-μm-thick beryllium window that transmits photons and are subjected to the pressure difference between the inside and outside of the detector chamber; (2) the vacuum collimation aperture in the kovar lid, which acts as a directional guide and has an area of 1.6 mm× 1.6 mm and height of 55 mm;(3)the aluminum frame layer,which is mainly used to shield low-energy electrons and background photons;(4)the ceramic,which encases the detection gas; (5) the gas required for detection; and (6)the kovar lid and base,which represent the entire shell of the detector

Fig. 4 (Color online) a Photon detection efficiency for different photon energies. b Track length distribution for different energies. The simulations were carried out in Geant4 using different detection gases with the gas thickness of 1 cm. The other conditions were the same

Fig. 5 (Color online) Distribution of cross-sectional potential in the center of the THGEM. The vertical and horizontal coordinates represent the distances. The structure of the hole is clearly visible

4 Track reconstruction algorithm

4.1 Event selection

Fig. 6 (Color online) A detailed diagram of the reconstruction process(a)and a typical track image(b).This track was produced by a 8 keV photon in neon (80%) mixed with dimethyl ether (DME 20%) at 0.8 atmosphere. The horizontal and vertical coordinates correspond to pixels with the same size as the actual Topmetal-IIpixels. In b, the photoelectron hits in the figure correspond to the actual positions of the chip pixels.The details of Fig. 6a are described in Sect. 4

Several types of tracks can be difficult to reconstruct even under ideal simulation conditions. These include certain completely symmetrical tracks and certain elliptic or circular tracks. The experimental conditions are more complicated than the simulated ones.The presence of both electronic noise and multiple tracks in a frame may affect the reconstruction results. The isolated noise points can be removed through image erosion and dilation algorithm processing of the tracks. To estimate the track range, the coordinate values of the leftmost point (xmin), lowest point(ymin),rightmost point(xmax),and highest point(ymax)need to be considered. When the conditions xmax-xmin＜6 and ymax-ymin＜6 are satisfied,the track range is limited to six pixels, which may correspond to circular or elliptical tracks. By replacing the longer and shorter axes of x and y by L and S, respectively, a track that satisfies （Lmax-Lmin）/（Smax-Smin）＜1.5 is considered as an elliptical track and eliminated.The thresholds of 1.5 and 6 used here can be adjusted in the future. In addition, this process removes several difficult parts in the reconstruction process of symmetric charge distributions, such as finding the Bragg peak track, as explained in detail below. These tracks should be screened out in advance, while the remaining tracks are reconstructed in the following steps.

4.2 Event reconstruction

Electrons may be scattered multiple times as they pass through the gas, thus deviating from the initial ejection direction. When electrons are close to stopping, they release most of their energy in the form of a Bragg peak.Thus, only a modest charge is accumulated at the initial interaction point. In this case, monitoring only the start of the track is sufficient to achieve full reconstruction. The repeated iteration method is divided into the following four main steps,of which the first two are similar to those in the traditional method.

(1) The barycenter and gravity line are calculated from the location and charge (as the weight) of the track in Fig. 6a. This line roughly marks the direction of the track trend.

(2) A line represented by the green line in Fig. 6a is traced at the barycenter perpendicular to the gravity line.Because Bragg peaks are accumulated at the tail of the track, the head and tail can be separated using third-order moment calculations to determine the half-zone where the photoelectron interaction point is located. The third-order moment Mdis calculated as

where d is the distance between the pixel and the line perpendicular to the gravity line (d may be negative), q is the charge of the pixel, and ∑pixelsindicates summation over all the pixels. Thus, by comparing the d and Mdof each pixel with similar or dissimilar signs, it can be determined whether the pixel is within the region where the photoelectron Bragg peak is located,or in the region where the photoelectrons are produced, which is indicated by the semicircle in the figure. Tracks with Md=0 are symmetrical in shape and difficult to reconstruct; thus, such tracks are rejected.

Fig. 7 (Color online) Superposition distribution of track reconstruction angle.The selected photons were 100%polarized with the energy of 8 keV. The photons entered the detector vertically with the angles θ and φ fixed to π/4 and π/2,respectively. θ is the azimuth angle of the X-ray electric vector, and φ is the polar angle with respect to the incident X-ray direction

where xnand ynrepresent the coordinates of each pixel,and x0and y0are the coordinates of the barycenter. The track length value is calculated at the best applicable energy point of the traditional algorithm,which is approximately a Gaussian distribution, and the mean value of the Gaussian distribution is used as the initial threshold for determining if the track length is long or short. If the track length is smaller than this value, the track is a short track. Otherwise,it is a long track.At the same time,the mean value of the Gaussian distribution is used as one of the parameters that can be adjusted later. Furthermore, for tracks that are too short, a re-division needs to be performed along a direction that makes the half-zone of the starting track longer. The spacing between the green lines is set as a modifiable parameter. The degree of curvature of the track can be determined by calculating ∑pixels Mmin

2 .Unlike the third-order moment calculation,the influence of the charge value is removed in the second-order moment calculation and only the shape of the track is considered, thus preventing the uneven charge distribution of the track from affecting the results. The variation in the second-order moment is expressed by

2 ).The variation is summed over the initial second-order moment to normalize the variation of the second-order moments for different tracks. The second-order moments of the green lines are computed in sequence, and the final split position chosen is the position with the most dramatic change. The above equation provides the maximum value at which a split will occur. To prevent the green line from dividing the outside of the track, a minimum track head length was set to preserve the first part of the track and the green line was forbidden from crossing this position. The minimum track head length is also a variable that can be adjusted. The actual parameters used were obtained after several rounds of empirical testing.

(4) The barycenter and center of gravity line are calculated in the remaining half of the partition.The center of gravity is identified in the diagram by the interaction point(marked by an asterisk in Fig. 6a)generated by the algorithm. The remaining track gravity line is illustrated in red in the diagram, while the arrow indicates the polarization direction obtained via the algorithm because the initial direction of photoelectron emission is also the polarization angle.

Figure 6 shows the tracks obtained from the simulation and the reconstruction process. Some of these parameters can be modified at the beginning of the reconstruction,including the pixel blocks used for filtering circular tracks,and the value of the minimum preserved track head length.These parameters have a direct impact on all aspects of the reconstruction, including the selection of pre-reconstruction events.Monte Carlo simulations were used to generate a certain number of photoelectron tracks. The above algorithm was then used to analyze the tracks line by line.The initial emission angle of each track was obtained through parameter optimization. The overall angle distribution is shown in Fig. 7, where the horizontal coordinate represents the angle and the vertical coordinate represents the number of instances.

Fig. 8 (Color online) Modulation factor as a function of photon energy. The selected photons were 100% polarized and entered the detector vertically at the same angle as that in Fig. 7. The signal acceptance rate is the proportion of events remaining after event selection

Fig.9 (Color online)The actual detector in the laboratory.The upper part is the detector chamber and the THGEM,and the bottom part is the Topmetal-II-

4.3 Event fitting analysis

The angular distribution of the K-shell photoelectron emission in the non-relativistic region is theoretically given by [26]:

Fig. 10 (Color online) Experimentally obtained tracks. a The THGEM structures were used to multiply the resulting track. Some background noise was present. b The tracks were processed with an image erosion and dilation algorithm for reconstruction,which almost completely removed the effect of the background

where β is the ratio of the photoelectron velocity to the speed of light, θ is the azimuth angle of the X-ray electric vector,and φ is the polar angle with respect to the incident X-ray direction.

The corresponding equation for the response amplitude after projecting the track onto a two-dimensional plane for processing is

The emission direction of polarized photoelectrons obeys such a distribution in theory.This equation was also used to fit the data shown in Fig. 7.In practice,modulation factors are often used to describe the polarization of the resulting photons because of slight deviations due to limited detection capability and the precision of the algorithm itself.The modulation factor is calculated as

where Mmaxand Mmincorrespond to the maximum and minimum values of the angle distribution in Fig. 7. The track reconstruction efficiency and the corresponding prescreening efficiency of the SXP are shown in Fig. 8. The low-energy portion is difficult to reconstruct because the track lengths are too short. The middle energy portion is the most suitable area for the algorithm.At higher energies and longer track lengths, the modulation factor decreases because the parameters of the algorithm need to be adjusted. The circles show the efficiency of the pre-selection,and the ordinate the proportion of tracks retained for reconstruction. It can be observed that only a few cases have been screened out and that the vast majority has been retained for reconstruction.

Based on the conditions set in the simulations,a 4.5 keV X-ray polarization source was used in the experiment to reconstruct the measurements. Figure 9 shows the actual structure of the detector. The upper part of the structure is the SXP detector, while the center is the THGEM. The bottom part is the Topmetal-II-.Both sides were equipped with gas pipes for ventilation and power supply connection wires for the application of high voltages. The signal collected by the Topmetal-II-is shown in Fig. 10.

Fig.11 (Color online)Angular distribution of experimentally obtained tracks after reconstruction.A 4.5 keV X-ray polarization source was used in the left figure. A 4.5 keV unpolarized X-ray source was used in the right figure

The experimental tracks were reconstructed, and their angular distributions are shown in Fig. 11. The experimental tracks yielded the modulation factor of ～28.60±0.86%,while the event rejection rate was 46.68%.The modulation factor of the unpolarized X-ray source was 2.057%.The tracks obtained from the THGEM experiment were very rough and poorly reconstructed because of the large pore size. There were very few tracks that could be easily reconstructed,as shown in Fig. 10.In addition to the event selection rules described in Sect. 4.1 (the range for the circular track judgment was nine pixels), there were tracks distributed at the edge of the image in the experiment that were not involved in the reconstruction. For a track to be reconstructed, the track must have a sufficient number of pixels in the central area. Both the size of the central area and the number of pixels used as the threshold are adjustable input parameters. The main reason for the inconsistency between the experimental and simulation results is the influence of the experimental background;for example, the X-ray sources were not 100% polarized, and readout noise or dark current was present. Meanwhile, the spatial resolution was insufficient, resulting in possible inaccuracies in the reconstructed photoelectron emission direction.The spatial resolution will be optimized in future studies.

5 Discussion and summary

The SXP adopts the most popular soft X-ray polarization detection method at present. The appropriate experimental gas and structure were selected through experiments and simulations. The iteration cut algorithm can provide reconstruction angles that are important for measuring the polarization and relevant for detection sources in the 2–10 keV energy range. The modulation factor reaches as high as ～57% above 7 keV. At low energies,the modulation factor is relatively low because of the short photoelectron tracks and the effect of electron diffusion in the THGEM avalanche. Experiments have shown that reducing the hole size of the THGEM structure can yield better results. The next step will be to optimize the hardware structure of the SXP and the algorithm, for example,by using deep learning. It is envisaged that highprecision X-ray polarization detection will have a significant impact in revealing the physical mechanisms of astrophysical sources in the future.

Author ContributionsAll authors contributed to the study conception and design. Material preparation, data collection, and analysis were performed by Xue-Feng Huang,Jin Zhang,Huan-Bo Feng,and Xi-Chen Cai. The first draft of the manuscript was written by Xue-Feng Huang and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript.