In-flight calibration of the Insight-Hard X-ray Modulation Telescope

HE adopts an array of 18 NaI(Tl)/CsI(Na) phoswich as the main detectors[4]. The diameter of each phoswich is 190 mm. The thicknesses of NaI(Tl) and CsI(Na) are about 3.5 mm and 40 mm, respectively. The working temperature is actively controlled at $18\pm 2{}^{\circ}\mathrm{C}$. The incident X-ray with most of its energy deposited in NaI(Tl) is regarded as a NaI(Tl) event. CsI(Na) is used as an active shielding detector to reject the background events from backside and events with partial energy loss in the NaI(Tl). The scintillation photons generated within the two crystals can be collected by the same photomultiplier tube (PMT). Signals from the PMT are pulse shaped to distinguish NaI(Tl) events and CsI(Na) events because of the different decay time in the two crystals. The energy loss, time of arrival and pulse width of each detected event are measured, digitized and telemetered to the ground. The CsI(Na) can also be used as a gamma-ray burst (GRB) monitor. The detected energy range in high gain mode is about 50–800 keV and is changed to 250 keV–3 MeV in low gain mode for CsI(Na) if the high voltage of PMT is decreased. In this paper, we only describe the calibration of NaI(Tl) events in high gain mode. The dead time of each detector is recorded online every second and also telemetered to the ground.

For each phoswich detector, a radioactive source ²⁴¹Am with an activity of 200 Bq is embedded into a plastic scintillator and viewed by a separate Multi-Pixel Photon Counter[4]. They are all mounted in the collimator and used as an automatic gain control (AGC) detector. A coincident measurement between the AGC detector ($\alpha$ particle of 5.5 MeV) and phoswich detector (X-ray of 59.5 keV) is labeled as a calibration event. The calibration events are saved like norm events, but with a different flag. The response of NaI(Tl) is not uniform in its large surface, so the calibration events are just used as the gain control and are not suitable to calibrate the energy scale of NaI(Tl) detectors in-orbit. Table 1 summarizes the characteristics of HE.

ME consists of three detector boxes and each box has 576 Si-PIN detector pixels read out by 18 ASIC (Application Specified Integrated Circuit)[5]. Each ASIC is responsible for the readout of 32 pixels. The thickness of the Si-PIN is 1 mm. The energy loss and arrival time of each detected event are measured, digitized and telemetered to the ground. In each detector box, two ²⁴¹Am sources are fixed in the corner of two ASIC and each one illuminates four pixels[5]. The working temperature of Si-PIN detectors in orbit is between $-40^{\circ}\mathrm{C}$ and $-10^{\circ}\mathrm{C}$ . Table 1 also summarizes the characteristics of ME.

Si-PIN detectors are fixed on the ceramic chip by silver glue. When the energy of incident X-rays is greater than 25.5 keV (K-edge of Ag), they have some probability of penetrating the Si-PIN and interact directly with silver. Fluorescence lines of silver will be generated due to the photoelectric effect with electrons in K-shell of silver and escape from the silver glue and then are detected by the Si-PIN detectors.

Although the Si-PIN signals are already screened by the onboard threshold, the events transmitted to the ground still contain low-energy thermal and electrical noise component, which varies significantly in orbit. To remove these noise events, a higher threshold must be applied to ME detectors. The new thresholds range from 6 to 10 keV for different Si-PIN detectors and we choose the new threshold at 10 keV for all of them in spectral analysis.

LE also consists of three detector boxes and each box contains 32 CCD236 which is a kind of Swept Charge Devices (SCD)[6]. CCD236 is the second-generation SCD, which has been developed by e2v company. It is built with a sensitive area of about 4 ${\rm cm}^{2}$ of each CCD236 and has four quadrants. Each quadrant has 100*100 pixels. The L-shaped readout mechanism has fast signal response but the position information is lost[6]. The readout cycle is 100 us. For events with energy above the onboard threshold, the energy loss and the readout time of each detected event are measured, digitized and telemetered to the ground. Besides this, LE also has the forced trigger events, which record the amplitude of the noise or the pedestal offset for each CCD detector every 32 ms[6]. The forced trigger events are also saved like physical events, but with a different type.

The spreading of the charge cloud over several pixels may cause split events. These split events may be read out in adjacent readout periods. In this paper, we focus on the calibration of the single events in order to reduce the effect of the charged particles. The working temperature in-orbit for CCD236 is between $-75^{\circ}\mathrm{C}$ and $-40^{\circ}\mathrm{C}$. The characteristics of LE is also summarized in Table 1.

Many experiments have indicated that the light output response of the NaI(Tl) crystal is not proportional to the deposited energies of X-rays and electrons [9],[10]. The light output is also non-monotonic at iodine K-edge at 33.17 keV and 50.2 keV for HE detectors [11]. On ground, three quadratic functions at three energy ranges (below iodine K-edge, from iodine K-edge to 50.2 keV and above the 50.2 keV) were used to parameterize the energy-channel (E-C) relationship of HE detectors[11].

The background of HE in-orbit is dominated by internal activation effects. Prominent background lines due to activation of iodine by cosmic and SAA protons are evident at 31, 56, 67 and 191 keV [8]. These four lines could be used to calibrate the E-C relation or energy scale of HE. Statistics sensitive Nonlinear Iterative Peak clipping (SNIP,[12]) method, which is also widely used in the gamma-ray spectrometry, is utilized to determine the continuum of the observed background spectrum. After subtracting the continuum component, the peak centroids of the four lines can be fitted with Gaussian functions. As shown in Figure 2, it illustrates the process of how to obtain the peak centroids of the four lines in the spectrum of blank sky observation.

According to the E-C relationship on ground, the peak centroids of the four lines are different with the expected values. As shown in Figure 3, the peak ratios of the net channels (where the offset of the electronics has been subtracted) between pre-launch and post-launch of the four background lines are not equal to 1 and here we just show the ratios for detID 0-5. Therefore, the E-C of HE detectors in orbit are different with the pre-launch calibration results, and thus must be re-calibrated.

In pre-launch calibration experiments, to model the E-C relation [11], we defined

where x is the net channel of the full energy peak, $E$ is the energy of incident photons from the monochromatic radiation, $A_{1},A_{2}$ and $A_{3}$ are the fitting parameters. In-flight, the parameter $\kappa$ that lessens the difference among the four background lines could be adopted to correct the pre-launch E-C and $\kappa$ is affected by the efficiency of the light collection and high voltage of the PMT. We construct the following model to represent the in-flight EC,

where $\kappa$ is the ratio of net channels between pre-launch and orbit. The half of the maximum and minimum of $\kappa$ of the four background lines as plotted in Figure 3 is regarded as the ratio of different detectors of HE. The accuracy of this gain calibration is estimated to be about 1% based on the deviations of the four calibration points. In order to validate the energy scale results, we jointly fitted the energy of cyclotron resonance scattering features (CRSF) of Her X-1 using the simultaneous observation with NuSTAR telescope and obtained consistent results at around 37.5 keV [13].

In order to monitor the energy scale of HE, we make full use of regular blank sky observations which are used to estimate the background model. Figure 4 shows the peak centriods of 191 keV line versus time. The standard deviations of the peak centroids range from 0.23 to 0.49 channel for all the HE detectors and most of the detectors are at 0.3 channel. We also utilize the data of ²⁴¹Am, which is located in the edge of each detector and permanently illuminates a small area of one HE detector, to monitor the energy scale. The top panel of Figure 5 shows the variations of the peak centroid of the 59.5 keV versus time. The variations are less than 0.01 channel for all the detectors after about 90 days. From Figure 4 and the top panel of Figure 5, the energy scale of HE remains stable after about three months in-orbit.

In the pre-launch calibration experiments of HE, we modeled the resolution of HE in channel space as

where $R(x)$ is the Full Width at Half Maximum (FWHM) in channel space and $x$ is the channel value, $B_{1},B_{2}$ and $B_{3}$ are the fitting parameters [11]. In-orbit, we have introduced another parameter $\lambda$ to estimate the resolution in channel space as

Before launch, we utilized an object-oriented toolkit, Geant4 ([14], [15], [16], version 4.9.4), to estimate the background of HE in orbit and found some lines can contribute to the profile of 31 keV line as shown in Table 2. The energies and the corresponding intensities which come from the Geant4 simulation, are also displayed in Table 2. We use the in-flight E-C to get the peak channel of the blended lines and adjust the parameter $\lambda$ from 0.9 to 3 with bisection method to broaden the lines. Till the difference between the width of simulated profile of 31 keV and that of the observed is less than 0.001 channel, the parameter $\lambda$ of each detector can be obtained. This process is shown in Figure 6. A Gaussian function with an exponential function and a small constant provides an excellent fit to the profile of 31 keV by use of the blank sky data. The width of 191 keV is also used to verify the parameter $\lambda$.

Except background line of 31 keV, we also monitor the width of 59.5 keV emitted by the ²⁴¹Am. A Gaussian function is used to fit the full energy peak profile of 59.5 keV, the sigma of the Gaussian versus time is shown in the bottom panel of Figure 5. The variations are less than 0.02 channel after about 90 days in orbit. The resolution of the 18 HE detectors also remains stable after about 90 days launch.

Each ME detector box contains two ²⁴¹Am sources which can continuously illuminate 8 pixels. The spectrum of pixels carried with ²⁴¹Am can be accumulated during the blank sky observations as shown in Figure 7. The energies and intensities about each line are from the database of National Nuclear Data Center (NNDC, [17].) After taking the efficiencies and resolution of ME detectors into account, the model energies of the four lines can be determined as 13.94 keV, 17.58 keV, 21.30 keV, and 26.34 keV.

In the calibration experiments on ground, the E-C relation of each pixel was linear from 11 keV to 30 keV using the spectrum of ²⁴¹Am source. The slopes and intercepts of E-C for all pixels were not a constant at different temperatures [8]. The slopes of E-C increased with the temperature whereas the intercepts remained almost the same at different temperatures. The slopes and intercepts at seven temperatures ranging from $-30^{\circ}\mathrm{C}$ to $-2^{\circ}\mathrm{C}$ were measured on ground and stored in CALDB as the primary E-C calculation. The slopes and intercepts at a given temperature can be obtained from the linear interpolation at two adjacent temperatures for each pixel.

In orbit, we use the pre-launch E-C model of each pixel to convert primary ADC channel to energy, then utilize the same linear function of energy versus pulse invariant (PI) channel to get the PI channel spectrum of ²⁴¹Am as x-axis plotted in Figure 7. A model with four Gaussians provides an excellent fit to the ²⁴¹Am spectrum between 11 keV and 27 keV. The spectrum of other pixels not illuminated with ²⁴¹Am are also generated from the blank sky data to fit the Ag peak (22.5 keV) using the pre-launch E-C. A model with a power-law and a Gaussian line is used to fit the data near Ag peak between 20 keV and 24 keV. Over this bandpass, this simple model provides an acceptable fit. The peak values of the five lines increase slowly with time as shown in Figure 8. A linear function is used to fit the evolution of the five lines. Till now, they increase by less than 1%. In the current status of software (HXMTDAS V2.02), this effect is not corrected in the gain for this small and slow change. But in the next software update, this effect will be added.

The FWHM of Ag line and the four nuclear lines from ²⁴¹Am are utilized to monitor the resolution of ME in orbit. As shown in Figure 9, the FWHM of ²⁴¹Am, which is used as the same parameter for the four lines in the spectral fitting, is plotted as a function of time. Although, there are some structures due to the effect of temperature, the change is very small. The resolution at 13.94 keV has increased by 0.5% till 700 days in orbit. Therefore, the resolution measured on ground can also be used for the in orbit data for all the pixels.

For LE, the data suitable for parameterizing the energy scale and monitoring its variations come from three regular observed sources: the internal background lines (Ni, Cu, and Zn), observations of Cas A with rich lines (Si, S, Fe, et al) and regular monitoring observations of the Crab Nebula providing an opportunity to measure the location of the K-edge of silicon at 1.839 keV.

The supernova remnant Cassiopeia A produces several strong lines easily visible in the LE spectrum [1]. Cas A is observed almost every month as a calibration source if the solar avoidance angle is allowed. We have fitted a model with a power-law continuum and several Gaussian lines to the data between different ranges as listed in Table 3 to obtain the energy peaks of these lines using the pre-launch E-C calibration. The emission lines of Ni, Cu and Zn generated by materials near the CCD detectors during the blank sky observations are also used to calibrate the energy gain of LE. As shown in Figure 10. the three lines are prominent in the measured spectrum of blank sky.

In the calibration experiments of LE on ground, the E-C relation of each CCD was also linear. The slopes and intercepts of E-C for CCD were not constant at different temperatures [8]. The same method for ME was also used to obtain the slopes and intercepts at different temperatures. In orbit, the forced trigger events of LE are accumulated every minute to fit the peak of the pedestal, then the physical events in this minute will subtract the pedestal value to obtain the pure pulse height. As for LE, the channel is referred to as the pedestal subtracted one. We use the pre-launch E-C results at different temperatures stored in CALDB to fit the peak energies of the lines listed in Table 3 with all the data of in-flight observations of Cas A and blank sky.

First, we check the effect of temperature and the differences between different CCDs in each observed spectrum of Cas A. The differences of Si and S energies between different small FOV CCDs is less than 8 eV at different observation time. So the data of small FOV CCDs are used to accumulate the spectrum of Cas A lines.

The observed spectrum of XMM/MOS is also used to fit the energies and intrinsic width of Cas A lines together with LE spectrum at different time. The fit results of energy peaks by XMM/MOS are regarded as the model energies of the lines as shown in the middle column of Table 3. If the pre-launch E-C of LE is suitable for the in-flight data, the observed energy minus the model energy of the lines would be expected to be around zero. But the gain has shifted in a non-linear way as shown in Figure 11 on different observation date. From Figure 12, all the peaks of Cas A decrease gradually and this phenomenon may be caused by the decrease of the charge transfer efficiency of the CCDs. In order to describe the evolution, a quadratic polynomial function is used to fit the change. From the fit result as shown in Figure 12, the peak values could be obtained on any day even the Cas A is invisible due to the observation constraint.

According to the residual shape shown in Figure 11, a cubic function is used to describe the gain model of LE as shown in Figure 13. The x-axis is the energy from the polynomial fit result of different peaks as show in Figure 12 on the same day and y-axis is the model energy fit with XMM/MOS. We compare the cubic function on different days and decide to supply the gain of LE every month to get the accurate gain results. The parameters of the cubic function in each month are stored in CALDB.

We also use the the following model to fit the K-edge of silicon to verify the gain of LE when Crab is observed,

where $E_{0}$ is the fit energy of Si K-edge, $A_{0}$ is the normalization coefficient, $C_{0}$ is the constant term and $\sigma_{0}$ is the width of $E_{0}$. $\rm{erfc}$ is the complementary error function. $Y$ is the rate of Cas A spectrum. The differences between the fit energy and K-edge energy at 1.839 keV are shown in Figure 14 and most of them are within $\pm 5$ eV. Therefore, the gain calibration of LE is effective.

After the calibration of LE gain, we generate the Cas A spectrum using the new gain in-orbit again. The widths of Si, S, and Fe using the pre-launch response file of LE are jointly fitted with XMM/MOS using Xspec. If the energy resolution of LE keeps the same as the measurement on ground, the intrinsic width of Si, S, and Fe will be same as the fit result of XMM/MOS. Actually, the fitted intrinsic widths of Si, S, and Fe for LE are larger than that of XMM/MOS, so the energy resolution is also changed, compared with the pre-launch calibration results. After subtracting the intrinsic width of Si, S, and Fe from the fit result of XMM/MOS, the extra broadening of the resolution of LE is obtained and it evolves with time and temperature. In order to parameterise the extra broadening of LE, we define the following two dimensional function,

where $W(t,T)$ is the extra broadening of LE except the pre-launch spreading, $t$ is the observation time and $T$ is the temperature, $c_{0}-c_{5}$ are the fitting parameters.

The results are plotted for each line in Figure 15. The residuals for the Si, S, and Fe are plotted in the bottom panel of Figure 15. For an observation, we calculate the mean temperature in its good time interval (GTI) and also the mean time. From the Equation (6), the extra broadening of Si, S, and Fe can be derived for this observation. Once the extra boarding widths of Si, S, and Fe have been derived, a liner function is used to fit the extra broadening of LE as shown in Figure 16.

Extensive pre-launch calibration experiments and modeling of the response were performed for HE detectors [11]. About 32 discrete energies covering the band from 20 keV to 356 keV were used to calibrate the energy scale and resolution of HE detectors. The 32 energies and their corresponding full peak channels were used to get the non-proportional response (NPR) of NaI(Tl) to electrons.

Geant4 (version 4.10.2) is utilized to perform the simulations. The model of low energy electromagnetic interaction, G4EmLivermorePhysics, is invoked in our simulation, and the fluorescence line and Auger process are activated. The primary events, including the type, energy, direction, and position of the initial particles are generated according to the calibration experiment setup on ground. Geant4 treats the particles one by one and tracks the trajectories of primary particles and secondary particles step by step.

At each step, we record the information of particles, particle type, kinetic energy, deposited energy and physical process involved and so on. Then the recorded information of every step can be used to select the particles of interest. The following model is used to represent the NPR of NaI(Tl) crystal to electrons,

where $E_{\rm{e}}$ is the kinetic energy of electrons generated by the incident X-ray interacting with NaI(Tl) crystal through the phot-electric effect and Compton scattering. Figure 17 is an example showing the simulated distribution of $E_{\rm{e}}$, which is generated in NaI(Tl) by absorption of 90 keV X-ray photons. $NPR(E_{e})$ is the pulse height generated by electrons with energy of $E_{\rm{e}}$ per keV with unit of channel/keV. $p_{1}-p_{5}$ are the fitting parameters.

We convolve the distribution of kinetic energy of electrons with NRP model to get the theoretical channel of full-energy peak. A minimization method called fmincon in MATLAB [18] is used to find the minimum of the nonlinear multi-variable function. After minimization, the parameters of NPR can be derived and the NPR function is shown in Figure 18. The NPR model on ground and in-orbit are both plotted in Figure 18 using the different input of E-C relation. Before launch, we compared the differences of full peak channel of 32 energies between measured and predicted with NPR model; most of the differences are within $\pm 0.2$ channel as plotted in Figure 19.

The pre-launch energy resolution in channel space was used to spread the simulated spectrum and we also compared it with the measured spectrum as shown in Figure 20. The centroids of full energy peak and escape peak are the same as the data. As a result, the NPR model can predict the response function very well. If we do not take the NPR of electrons in NaI(Tl) in the Geant4 simulation into account and only use the deposited energy to get the simulated spectrum, the escape peak shifts by about 3 keV in the data when the incident photons have energy of 90 keV. Once we have verified the simulated response with the measured 32 energies on ground, we generate the response function at other energies.

In orbit, the E-C model and resolution in channel space can be utilized to produce the NPR model and response matrix files as same as we have done on ground.

The pre-launch RMFs was calibrated across the ME energy range (9-30 keV) at the calibration facility in the Institute of High Energy Physics (IHEP). The double crystal monochromator was used to generate mono-energetic photons [5]. The fluorescence lines of Ag was found in the measured spectrum when the photon energies were above the K-edge of silver (25.5 keV), and we can use the measured spectrum and Monte Carlo (MC) simulation to confirm the thickness of Ag glue under the Si-PIN detectors of ME.

We also utilize Geant4 (version 4.9.4) to perform the simulations. We invoke the low-energy electromagnetic process and consider the photoelectric effect, Compton scattering, Rayleigh scattering in the simulation. Fluorescence and Auger processes are also loaded in the photoelectric effect. The X-ray events spectrum resulting from monochochromatic radiation with an energy $E$ shows a simple Gaussian for photon energies below the Ag K-shell edge (25.5 keV). But for photons with energy above the K-edge of Ag, it consists of a Gaussian photo-peak and Ag fluorescence peak at 22.5 keV as displayed in Figure 21. We have adjusted the thickness of Ag glue under the SI-PIN till the differences of MC simulation spectra and measured spectra are the minimum as shown in Figure 21. Finally, the thickness of Ag glue under the Si-PIN is estimated as 14 um.

As the resolution of ME has changed only slightly in orbit as shown in Figure 9, so the redistribution matrix of ME is still taken as the pre-launch one.

As LE is a type of CCD detector, the spreading of the charge cloud over several pixels produced by a photon or a charged particle can be read out by several adjacent periods. To eliminate events due to charged particles, and to obtain good energy resolution, we only consider single events, which are not split as two or more adjacent readout, as valid X-ray events.

The pre-launch RMF calibration of LE was done using the calibration facility in IHEP, using twenty discrete energies covering the 0.9-12 keV energy range. The X-ray single events spectrum resulting from monochochromatic radiation with an energy $E$ significantly differs from a simple Gaussian. It consists of multiple components: a Gaussian photo-peak with a shoulder on the low energy side, a plateau extending to low energies. For photon energies above the Si K-shell edge (1.839 keV) two additional features appear: an escape peak of energy $E-E_{\rm{Si}}$ and a Si $K_{\alpha}$ fluorescence peak at 1.74 keV. Displayed in Figure 22 is a measured spectrum with energy of 7.48 keV.

At first, we also used Geant4 to simulate the response function of LE and considered the charge transfer process [19] in the simulation. It is very difficult to use the same parameters, like the coefficient of the primary cloud radius, diffusion length in depletion region and so on to make all the MC spectra the same as the measured spectra. Finally it was decided to use the two dimensional (one is the energy of incident photons and the other is the ADC channel of LE) interpolation to generate the response function of LE, since we have twenty discrete energies covering from 0.93 keV to 11.9 keV, and almost for every 0.5 keV interval we have a measured spectrum below 8 keV. After interpolation, the probability in every ADC channel at energy $E$ will be calculated, and we smooth it to get the final model of response function. As a verification, we compare the measured spectrum and model spectrum for photons with energy of 7.48 keV and they show good agreement with each other as plotted in Figure 22. The model spectrum of 7.48 keV is interpolated by photons with energy of 7.01 keV and 8.02 keV. Therefore, the pre-launch RMF is generated from 0.855 keV to 11.915 keV and step is 0.01 keV.

When we use the pre-launch RMF to fit the in-flight spectrum of Cas A, extra broadening is needed to fit the line profiles well as described in section 3.5. The extra broadening of the resolution is convolved with the pre-launch RMF to generate the RMF in-flight.