Measurement of the vertical atmospheric density profile from the X-ray Earth occultation of the Crab Nebula with Insight-HXMT

In this paper, the X-ray Earth occultation (XEO) of the Crab Nebula is investigated by using the Hard X-ray Modulation Telescope (Insight-HXMT). The pointing observation data on the 30th September, 2018 recorded by the Low Energy X-ray telescope (LE) of Insight-HXMT are selected and analyzed. The extinction lightcurves and spectra during the X-ray Earth occultation process are extracted. A forward model for the XEO lightcurve is established and the theoretical observational signal for lightcurve is predicted. The atmospheric density model is built with a scale factor to the commonly used MSIS density profile within a certain altitude range. A Bayesian data analysis method is developed for the XEO lightcurve modeling and the atmospheric density retrieval. The posterior probability distribution of the model parameters is derived through the Markov Chain Monte Carlo (MCMC) algorithm with the NRLMSISE-00 model and the NRLMSIS 2.0 model as basis functions and the best-fit density profiles are retrieved respectively. It is found that in the altitude range of 105--200 km, the retrieved density profile is 88.8% of the density of NRLMSISE-00 and 109.7% of the density of NRLMSIS 2.0 by fitting the lightcurve in the energy range of 1.0--2.5 keV based on XEOS method. In the altitude range of 95--125 km, the retrieved density profile is 81.0% of the density of NRLMSISE-00 and 92.3% of the density of NRLMSIS 2.0 by fitting the lightcurve in the energy range of 2.5--6.0 keV based on XEOS method. In the altitude range of 85--110 km, the retrieved density profile is 87.7% of the density of NRLMSISE-00 and 101.4% of the density of NRLMSIS 2.0 by fitting the lightcurve in the energy range of 6.0--10.0 keV based on XEOS method. This study demonstrates that the XEOS from the X-ray astronomical satellite Insight-HXMT can provide an approach for the study of the upper atmosphere.

obtain the spatiotemporal variation characteristics of atmospheric density on a global scale, the remote sensing by satellites is gradually developed. As a scientific instrument of Thermosphere Ionosphere Mesosphere Energetics and Dynamics (TIMED) that is the initial mission under NASA's Solar Terrestrial Probes Program, Sounding of the Atmosphere Using Broadband Emission Radiometry (SABER) can obtain vertical profiles of several atmospheric constituents, such as O 3 , H 2 O and CO 2 , as well as neutral atmospheric density in the altitude range of ∼10-140 km (Russell et al., 1999;Meier et al., 2015;Rezac et al., 2015).
In addition to the direct measurements of atmospheric density by sounding rockets, satellites and falling sphere measurements, retrieval of atmospheric density can also be carried out by an indirect method that usually refers to the method of occultation sounding. Stellar occultation has a long history as an atmospheric diagnostic method. The technique of retrieval of atmospheric density by occultation is gradually developed. There are some previous studies on the retrieval of atmospheric density of specific species by stellar occultation in the ultraviolet band. Hays and Roble (1973) obtained the nighttime vertical distribution of ozone number density in an altitude range of 60-100 km at low latitudes by analyzing the results of approximately 12 stellar occultations in the ultraviolet band near 2500 Å. Aikin et al. (1993) measured the molecular oxygen densities in the altitude range of 140-220 km based on the solar occultation data obtained from the ultraviolet spectrometer/polarimeter (UVSP) on the Solar Maximum Mission (SMM) spacecraft. The density profiles of ozone and nitrogen dioxide were inverted and evaluated by an optimal estimation algorithm using solar occultation data from SCanning Imaging Absorption spectroMeter for Atmospheric CHartographY (SCIAMACHY) in the UV-Vis wavelength range (Meyer et al., 2005). Lumpe et al. (2007) used an optimal estimation algorithm to obtain the O 2 density profiles between 110 and 240 km by analyzing solar occultation data at three nominal wavelengths (144, 161 and 171 nm). In addition to relevant studies in ultraviolet band, occultation in infrared and radio band has also been extensively studied. The water vapour number density profiles in the altitude range of 15-45 km was retrieved by the solar occultation data with SCIAMACHY in the wavelength region around 940 nm (Noël et al., 2010). The radio occultation technique can be often used to retrieve the electron density in the ionospheric by the Abel integral equations (Hajj and Romans, 1998;Lei et al., 2007;Chou et al., 2017). There are also occultation measurements that retrieve atmospheric density for specific species, such as SOFIE/AIM (McHugh et al., 2008;Rong et al., 2016), GOMOS/Envisat (Renard et al., 2008;Kyrölä et al., 2010), SAGE series (Degenstein et al., 2018;McCormick et al., 2020) and POAM series (Rusch et al., 2001;Lumpe et al., 2002).
The study of atmosphere by X-ray occultation is a new interdisciplinary study. The atmospheric extent of Titan was measured by the transit of the Crab Nebula in the X-ray band on 2003 January 5 observed by the Chandra X-Ray Observatory (Mori 2021). However, the retrieved atmospheric densities are significantly lower than the model density in some altitude ranges, and the difference between the measured values and the model values may result from the long-term accumulation of greenhouse gases, imperfect climatological estimates of solar and geomagnetic effects, temperature profile differences or gravity waves (Determan et al., 2007;Katsuda et al., 2021). Therefore, it is very important to cross-check the density structure of Earth's atmosphere by observations from other X-ray satellites like Insight-HXMT, to further verify the difference between retrieved results and model density.
The X-ray Earth occultation sounding (XEOS) has many advantages as an atmospheric diagnostic method. The X-ray photons are absorbed directly by the K-shell and L-shell electrons of atoms, including atoms within molecules, in the extinction process. Therefore, the ionized states, electronic states and chemical bonds within the molecules of atmospheric components have no effect on the absorption of X-rays. XEOS can retrieve the neutral atmospheric density in the upper mesosphere and lower thermosphere. In addition, the global distribution of neutral atmosphere in the upper mesosphere and lower thermosphere can be obtained by analyzing a large number of X-ray Earth occultation data from X-ray satellites, and the temporal evolution characteristics of neutral atmosphere density in this altitude range can also be studied. In this work, we use the Insight-HXMT observational data to study the X-ray Earth occultation of the Crab Nebula in order to demonstrate the capability of Insight-HXMT as an atmospheric diagnostic instrument. In addition, the retrieved results of Insight-HXMT and RXTE are used to cross-check the density structure of Earth's atmosphere, so as to confirm the existence of differences between retrieved results and model values. Here, the theoretical model of lightcurve is established by simulating the observations of the Low Energy X-ray telescope (LE)  to the Earth occultation of the Crab Nebula. A Bayesian data analysis framework is developed for the XEO lightcurve modeling and the atmospheric density retrieval. We use Markov Chain Monte Carlo (MCMC) algorithm to calculate the posterior probability distribution of model parameters, which is a method of inverting model parameters using Bayesian inference (Sharma, 2017). Finally, the Earth's atmospheric density in the altitude range of 85-200 km is retrieved.
The paper is structured as follows. Sect. 2 describes the observations and data reduction. Sect. 3 shows lightcurve modeling and density profile retrieval. The conclusions and discussions are given in Sect. 4.

Observations and data reduction
Insight-HXMT is the first X-ray astronomy satellite in China (Zhang et al., 2018;Li et al., 2018;. It is designed for the three main scientific objectives including Galactic Plane scanning, X-ray binaries observation, Gamma-Ray Bursts and Gravitational Wave Electromagnetic counterparts monitoring and studying . The Insight-HXMT mainly carries four scientific payloads, including the High Energy X-ray telescope (HE), the Medium Energy X-ray telescope (ME), the Low Energy X-ray telescope (LE) and the Space Environment Monitor (SEM) (Zhang et al., 2018. The quality of calibration directly determines the achievement of the three main scientific objectives. The Crab Nebula is one of the brightest X-ray sources in the sky, with its stable evolution and brightness. Therefore, the Crab Nebula is an excellent calibration source for many X-ray satellites. The Crab Nebula as a standard candle has been widely used for in-flight calibration to the X-ray source position l∞. rsat is the distance from the satellite to the center of the Earth. r is the distance from the tangent point to the center of the Earth. R is the radius of the Earth. The dotted line shows the satellite orbit. The extent of the atmosphere is marked by a double sided arrow. The location of the tangent point is marked. The Insight-HXMT and the Crab Nebula are also marked. For clarity, only three layers of the atmosphere are marked, as shown by the black solid lines. The thick, orange solid lines with arrows show the X-rays before absorption and scattering. The orange dashed lines with arrows show the absorption and scattering of X-rays by the atmosphere. The thin, orange solid lines with arrows show X-rays passing through the atmosphere. of space missions in X-ray astronomy (Kirsch et al., 2005;Meyer et al., 2010). In this work, the Crab Nebula is chosen as our observation source. The pointed observation mode is selected for the study of the X-ray Earth occultation.

Observation geometry
As Crab Nebula sets behind or rises from the limb of Earth as seen by Insight-HXMT, the X-ray flux from the Crab Nebula detected by Insight-HXMT varies due to the absorption of X-ray photons by Earth's atmosphere. The observation geometry of the X-ray Earth occultation of the Crab Nebula with Insight-HXMT is shown in Figure 1. In the process of the X-ray Earth occultation, the atmosphere reduces the flux of the X-rays detected by Insight-HXMT. As the line of sight moves closer to the Earth's surface (setting), more X-ray photons are absorbed by the atmosphere, and vice versa.

Data reduction
The observation data of the Low Energy X-ray telescope (LE) are used in this study. We select the photons observed by detectors for LE because the extinction effect at lower energy band during the occultation is obvious due to the larger cross section of the Earth's upper atmospheric compositions ( Figure 2). LE consists of three detector boxes each containing 32 pieces of CCD236  that is the second-generation Swept Charge Device (SCD) designed for X-ray spectroscopy (Holland and Pool, 2008;Zhao et al., 2019). The collimators divide each detector box of LE into four kinds of fields of view (FOV). For each detector box, 20 CCD236 have small FOVs of 1.6 • ×6 • , 6 CCD236 have wide FOVs of 4 • ×6 • , 2 CCD236   (Zhao et al., 2019;Chen et al., 2020;. The detector response matrix is generated through calibration database hxmt CALDB (v2.05) 1 . Only observations from the small FOV detectors excluding the detector ID of 29 and 87 that are damaged are used for analysis in order to get accurate background in this paper.
The information of the observation data selected in the data reduction is listed in Table 1, including the observation ID, the start and end time of the observation, the target source, and the right ascension and declination of the source in the coordinate system J2000. In the data reduction process, by using HXMTsoft(v2.04) 2 , we extract the lightcurves and spectra of the Crab Nebula during Earth occultation recorded with LE. For the LE instrument, the photon counts are recorded by CCD236 detectors.
Since we mainly study the occultation process of the Crab Nebula by the Earth atmosphere, the good time interval is screened by the following criteria: ELV (the elevation of the pointing direction above the horizon) less than 10 • .

Description of spectra and lightcurves
The comparison of the X-ray energy spectra during the occultation process is shown in Figure 3. Five X-ray spectra in different altitude ranges are shown for clarity. These five energy spectra in blue, red, orange, magenta and green in Figure 3 are derived The comparison of the X-ray energy spectra during the occultation process. These spectra cover an energy range of 1-10 keV (0.1240 -1.2398 nm). The unattenuated X-ray energy spectrum is shown in blue. The red, orange and purple data points are the attenuated X-ray energy spectra with the decreasing tangent point, respectively. The energy spectrum in green is the fully attenuated X-ray spectrum.
from the results of subsamples at altitudes of 160-170 km, 120-130 km, 100-110 km, 90-100 km and below 70 km, starting from an unattenuated energy spectrum to the partially attenuated energy spectrum, and ending with a fully attenuated energy spectrum. It is shown that the flux of the X-ray energy spectrum attenuates with reducing tangent point altitude during the occultation. Moreover, the absorption of X-ray photons by the atmosphere decreases with increasing energy.
In the process of X-ray Earth occultation of the Crab Nebula, the Crab Nebula and Earth's atmospheric disks are tangent at four contact times t I -t IV , illustrated in Figure 4. The total duration is t T = t IV − t I , the full duration is t F = t III − t II , the ingress duration is t o = t II − t I , and the egress duration is t o = t IV − t III . The occultation depth δ is the X-ray flux attenuation due to the extinction.
The egress duration is about 14 seconds during the occultation process analyzed in this study. We divide this duration into 35 bins in order to have good time resolutions and high signal-to-noise ratio (SNR). The lightcurves of three different energy bands during the occultation process are shown in Figure 5. The abscissa time is converted to tangent point altitude. In this study, the maximum height difference between two adjacent tangent points is 836 meters, and the mean height difference between two adjacent tangent points is about 673 meters. For clarity, the data points in Figure 5 are displayed by taking one point every ten points from the initial data points. The dashed lines of blue, red and green represent the modelled lightcurves.
The green and blue shadow colored regions correspond to the extinction process for the occultation in the energy bands of 1.319 keV-1.725 keV and 7.006 keV-7.412 keV, respectively. For clarity, the height range for occultation between 3.350 keV and 3.756 keV are not marked. The lightcurve in the energy range of 1.319 keV-1.725 keV starts to attenuate at 150 km and it is completely attenuated at 102 km. The lightcurve in the energy range 7.006 keV-7.412 keV starts to attenuate at 100 km and it is completely attenuated at 85 km.

Lightcurve modeling and density profile retrieval
In this section, we will describe the details of the lightcurve modeling for the X-ray Earth occultation of the Crab Nebula. In this study, we model the Earth occultation as a measurement method for atmospheric density.
X-rays can be absorbed by the photoelectric effect. The ionized states, electronic states and chemical bonds within the molecules of atmospheric components have no effect on the absorption of X-rays in the extinction process. X-ray photons are absorbed directly by the K-shell and L-shell electrons of atoms, including atoms within molecules. Therefore, the X-ray Earth occultation can work as an atmospheric diagnostics method. In this case, the source celestial coordinates and the satellite positions are known, whereupon the atmospheric density profile can be treated as the unknown. It is impossible to distinguish atoms from molecules (in the calculation process, although X-ray photons interact directly with the K-and L-shells electrons of atoms (including atoms within molecules), the O 2 (or N 2 ) counts as one absorbing "particles" in the calculation, so the total neutral atmospheric density profile can be retrieved by X-ray occultations (Katsuda et al., 2021).
The schematic of lightcurve modeling of the X-ray Earth occultation is shown in Figure 6. The attenuation process of the X-ray energy spectrum can be described by Beer-Lambert law during the occultation. The attenuation energy spectrum is convolved with the detector response matrix to obtain the forward model. Given the data and the forward model, Bayesian inference is used to estimate the model parameters. Given the prior distribution and the likelihood function, the posterior probability distribution of the model parameters is calculated by MCMC. The best fit model is obtained from the posterior probability distribution of the model parameters. The results will be compared with other measurements and models. The details of the data analysis will be described in the following subsection.

Forward model
The attenuation process of X-rays in the atmosphere can be described by Beer-Lambert law where I 0 is the unattenuated source spectrum, which is a function of energy, e −τ is the transmittance, τ is the optical depth, which has the form where s labels the gas components in the Earth's atmosphere, γ is the total correction factor, n s is the number density of each component of the atmosphere along the line of sight. Based on the spherical symmetry assumption of the Earth atmosphere, the number density is converted to column density by Abel integral. σ s is the X-ray cross section (photoelectric absorption and scattering cross section) of each component in the atmosphere. X-ray photons are absorbed or scattered by atoms, and in the energy range of interest in this paper, the scattering effect can be ignored because it is too small relative to the photoelectric absorption effect. However, the scattering cross section is still included in the calculation.
The modelled lightcurves with this forward model are shown in Figure 7 from Insight-HXMT for the X-ray Earth occultation of the Crab Nebula. The normalized flux with the orbital phase is shown. From Figure 7, we can see that the occultation depths (the difference between the highest and lowest point of the same light curve) for different energy bands are very different. Here, the atmospheric model NRLMSISE-00 (Picone et al., 2002) is chosen as our input data in the lightcurve calculations in Figure   7.
I 0 was fitted by using Xspec, a standard software package for spectrum fitting in X-ray astronomy. To fit the unattenuated spectrum of the Crab Nebula, we use the model wabs × powerlaw (Godet et al., 2009;Yan et al., 2018), where wabs is the interstellar absorption model in Xspec (Morrison and McCammon, 1983;Arnaud et al., 1999). The model powerlaw represents a simple power law shape of spectrum to fit. The data description and fitting result of unattenuated energy spectrum are listed in Table 2. The reduced chi-squared (Mighell, 1999) is 1.06 in this fitting, which indicates that the fit is good. The best-fit model and the unattenuated energy spectrum data are shown in the upper panel of Figure 8. The blue dots with the error bars are the unattenuated spectrum of the Crab Nebula observed by LE. The red solid line is the best-fit model. The lower panel of Figure 8 shows the residuals of the fit.
The main atmospheric components causing extinction are Oxygen (O, O 2 ), Nitrogen (N, N 2 ) and Argon during the occultation. The absorption of N and O to X-ray has similar characteristics because the energy dependence of the X-ray cross sections  O through X-ray occultation in our interest energy band, but their total atmospheric density (N+O+O 2 +N 2 ) distribution can be calculated. In addition, although X-ray photons interact directly with the K-and L-shells electrons of atoms (including atoms within molecules), the O 2 (or N 2 ) counts as one absorbing "particles" in the calculation. Ar is an atmospheric constituent of less content relative to N and O in the Earth's atmosphere. The atmospheric density of Ar is 0.029%-0.943% of the total density of N and O according to the NRLMSISE-00 model in the altitude range of 20 km-200 km. But the X-ray cross section of Ar is larger by almost one order of magnitude than that of N and O over the energy range of interest. Therefore, Ar is included in our model.
The number density of each atmospheric component needs to be given as input data in the process of density profile retrieval with a forward model (Determan et al., 2007). In the following, the atmospheric model NRLMSISE-00 and NRLMSIS 2.0 (Emmert et al., 2021) are chosen as our input data in the modeling, respectively. The NRLMSISE-00 model is one of the  Table 3. In Table 3, the time corresponds to the middle tangent point altitude during occultations. The geographical latitude and longitude are calculated with coordination tranformation from the coordination of the tangent point in J2000. F 10.7 and Ap are the solar activity index and the geomagnetic activity at this time.
The F 10.7 index is one of the most widely used index to characterize the level of solar activity (Tapping, 2013), and the Ap index is used to characterize the geomagnetic activity (Clúa de Gonzalez et al., 1993). These data are obtained from the Space Environment Prediction Center 3 .
The X-ray cross section σ s of each component of the gas can be obtained through the photon cross section database XCOM (Berger and Hubbell, 1987;Berger et al., 2010). The X-ray cross sections of Ar, O and N are shown in Figure 2. The specific form of the forward model of X-ray occultations is as follows (Determan et al., 2007), where R is the detector response matrix . In addition, the forward model also contains background noise B.
The background noise should be included in the forward model, as shown in eq. (3).

Density profile retrieval
Instead of subtracting the background, the background and the source counts can be modeled synchronously using Poisson statistics in the Bayesian framework (Olamaie et al., 2014). Bayes' theorem combines observation data with the prior distribution of the parameter of interest θ from a specific model to obtain the posterior probability distribution of the parameter. In this work, the posterior probability distribution of θ = {γ, B} for the forward model shown in eq. (3) applied to the data D can be given by Bayes' theorem (Bayes and Price, 1763), where D is the observation data, M is the forward model, p(θ|M ) is the prior distribution, p(D|θ, M ) is the likelihood and p(D|M ) is the Bayesian evidence.
Both the X-ray observed counts and the background data follow Poisson statistics so that the X-ray likelihood function, L, is given by where D i is the observation data in the ith bin, M i is the ith forward model value. The natural logarithm of the likelihood function can also be used for parameter estimation as the C statistic (Cash, 1979).
In this paper, we analyze the lightcurve in the energy range of 1.0-2.5 keV, 2.5-6.0 keV and 6.0-10.0 keV, it is found that these energy ranges are indeed sensitive to the altitude range of 105-200 km, 95-125 km and 85-110 km, respectively, as shown in Figure 9. The red shadow in Figure 9 indicates the occultation range and the blue shadow in Figure 9 indicates the energy range.
The Markov Chain Monte Carlo (MCMC) method is one of the parameter estimation methods used for Bayesian inference (Sharma, 2017). The density profile retrieval implements the MCMC method, which samples from a probability distribution using Markov chains (Chib and Greenberg, 1995;Dunkley et al., 2005;Hogg and Foreman-Mackey, 2018  In order to show the difference between the lightcurves, we amplified the observed lightcurve and the four model lightcurves in the altitude range of 120-124 km, 106-109 km and 95-100 km in panel (a), panel (b) and panel (c) in Figure 12, respectively.

Results testing
In this section, the Pearson's χ 2 test (Pearson, 1900;Cochran, 1952) is used to test the XEO measurements, NRLMSISE-00/NRLMSIS 2.0 model prediction for the description of the XEO lightcurve.
In this study, the following null hypotheses are proposed. The lightcurves predicted by XEO measured density profile and NRLMSISE-00/NRLMSIS 2.0 model simulated density profile fit the observed lightcurve well, respectively. It is found that the null hypothesis can not be rejected even at 84%, 90% confidence level for the two XEO measurements in the energy range of 1.0-2.5 keV, the null hypothesis can not be rejected even at 55%, 64% confidence level for the two XEO measurements in the energy range of 2.5-6.0 keV and the null hypothesis can not be rejected even at 68%, 69% confidence level for the two XEO measurements in the energy range of 6.0-10.0 keV, as shown in Table 4. It is found that the null hypothesis can not be rejected at 95% confidence level for the NRLMSISE-00/NRLMSIS 2.0 predictions, respectively, except for the NRLMSISE-00 predictions by the model lightcurve in the energy range of 1.0-2.5 keV. Goodness-of-fit testing results between the observed lightcurve and the extinction curve predictions with XEO measured density profiles, the NRLMSISE-00/NRLMSIS 2.0 model simulated density profiles are also listed in Table 4.
Goodness-of-fit testing is carried out for the observed lightcurve and four model lightcurves in the energy range of 1.0-2.5 keV. As shown in Table 4, the χ 2 /dof and p-value between the observed lightcurve and the extinction curve predicted with XEO retrieved density profile based on NRLMSISE-00 are 1.0599 and 0.1604, where, dof represents degree of freedom, i.e., the number of sample points minus the number of variables. In this paper, 551 sample points are used for fitting, with two variables of correction factor and background noise, so dof=549. The χ 2 /dof and p-value between the observed lightcurve and the extinction curve predicted with XEO retrieved density profile based on NRLMSIS 2.0 are 1.0756 and 0.1074. The χ 2 /dof and p-value between the observed lightcurve and the extinction curve predicted with NRLMSISE-00 predicted density profile are 1.1220 and 0.0249. The χ 2 /dof and p-value between the observed lightcurve and the extinction curve predicted with NRLMSIS 2.0 predicted density profile are 1.0783 and 0.0997. The results show that the lightcurves based on XEO retrieved density can better describe the observed lightcurve. Compared with the retrieved results, the atmospheric density predicted by NRLMSISE-00 model overestimates by 11.2%, and the atmospheric density predicted by NRLMSIS 2.0 model underestimates by 9.7%.
Goodness-of-fit testing is carried out for the observed lightcurve and four model lightcurves in the energy range of 2.5-6.0 keV. As shown in Table 4, the χ 2 /dof and p-value between the observed lightcurve and the extinction curve predicted with XEO retrieved density profile based on NRLMSISE-00 are 1.0091 and 0.4540. The χ 2 /dof and p-value between the observed lightcurve and the extinction curve predicted with XEO retrieved density profile based on NRLMSIS 2.0 are 1.0612 and 0.3669. The χ 2 /dof and p-value between the observed lightcurve and the extinction curve predicted with NRLMSISE-00 predicted density profile are 1.3321 and 0.0802. The χ 2 /dof and p-value between the observed lightcurve and the extinction curve predicted with NRLMSIS 2.0 predicted density profile are 1.3331 and 0.0797. It is found that the lightcurve predictions based on the XEO retrieved density profiles can describe the observed lightcurve better, and gaps between the retrieved density profiles and the model simulated ones exist. Compared with our retrieved results, the density profile of the NRLMSISE-00 model is overestimated by 19%, and the density profile of the NRLMSIS 2.0 model is overestimated by 7.7%.
Goodness-of-fit testing is carried out for the observed lightcurve and four model lightcurves in the energy range of 6.0-10.0 keV. As shown in Table 4 Table 5. AIC and BIC are used for model selection, and can also be used to compare models. Usually, we choose the model with the smallest AIC and BIC. However, the values of AIC and BIC in this paper show that solar and geomagnetic activity has a great influence on model shape. Because AIC and BIC values vary greatly under different solar and geomagnetic activity in a relatively low energy range. Goodness of fit between the observed lightcurve and the model lightcurves under the different solar activities and geomagnetic activities is also evaluated by χ 2 /dof and p-value in Table 5. It can be shown that the solar activity and the geomagnetic activity have great influence on the shape of model lightcurves. In addition, with the increase of altitude, solar and geomagnetic activities have a greater impact on the model lightcurves. The effects of solar activities and geomagnetic storms for the XEO lightcurve modeling and density retrieval will be further investigated in the future.

Comparison to the results from altitude independent method by spectrum fitting
Based on the energy spectrum fitting method during X-ray occultation (Katsuda et al., 2021;Yu et al., 2022), the altitude independent atmospheric density retrieved results can be obtained, and the overlap of the tangent point altitude can be effectively avoided. In order to prove the reliability of our retrieved results in the paper, we compare our results with the results of energy spectrum fitting (Yu et al., 2022). The optical depth of the forward model based on energy spectrum fitting is given by the following equation, where γ h represents correction factors in different altitudes ranges. By combining Eq.
(3) and Eq. (6), the forward model of the energy spectrum fitting for different altitude ranges is given. By fitting the energy spectrum data in different altitude ranges, the correction factors in corresponding altitude ranges can be obtained, namely γ h . So multiplying γ h with the input data from the NRLMSIS 2.0 model, the atmospheric density in different altitude ranges can be retrieved independently.
In the paper, by fitting the energy spectrum data in the energy range of 1-10 keV, we obtain the atmospheric density values in the altitude range of 100-200 km, and extract the energy spectrum data every 10 km. The comparison between the best-fit model and energy spectrum observational data is shown in Figure 14. The retrieved results based on energy spectrum fitting and the results with lightcurve fitting are shown in Figure 15, where the solid blue line represents the retrieved results of spectrum fitting, the solid red line represents the model density profile of NRLMSIS 2.0, and the solid green line represents the retrieved results with lightcurve fitting in the energy range of 1.0-2.5 keV. It is found that the retrieved results based on the lightcurve fitting are qualitatively consistent with the retrieved results of the energy spectrum fitting method. In the altitude range of 180-200 km, because the number of X-ray photons absorbed by the Earth's atmosphere is less than the X-ray photon counting error, the retrieved results based on energy spectrum fitting have large uncertainty. However, the reliability of the results based on lightcurve fitting is proved to be consistent with that by altitude independent method by spectrum fitting.

Conclusions and discussions
In this paper, we have studied the X-ray Earth occultation of the Crab Nebula with the pointing observation data from Insight-HXMT. We have presented a detailed Bayesian data analysis method for the extinction lightcurve modeling from the X-ray Earth occultation process. The theoretical predicted XEO observational lightcurve is calculated with the lightcurve forward model. The data recorded by the Low Energy X-ray telescope (LE) of Insight-HXMT are analyzed and the density profile is retrieved. The results are tested and validated with the measurements from the RXTE satellite and the retrieval results with the altitude independent method by spectrum fitting.
We have shown from the XEO extinction lightcurve modeling that the X-ray astronomical satellite Insight-HXMT can be used to retrieve atmospheric density by the X-ray Earth occultation of celestial sources. altitude range. In other words, an averaged scale factor for density profile in an altitude range is obtained by the lightcurve fitting method. We confirmed the measured density profile from lightcurve fitting by comparing to the ones by a standard spectrum retrieval method with an iterative inversion technique. The occultation data from larger effective areas can be used to retrieve atmospheric density to reduce the influence of energy integration. A more detailed description of this problem will be discussed in the future.
The difference between the measured and model values may result from the long-term accumulation of greenhouse gases, imperfect climatological estimates of solar and geomagnetic effects, differences in temperature profiles, and the influence of gravity waves (Determan et al., 2007;Katsuda et al., 2021). And the differences can also be due to model errors and/or retrieval errors. In order to further clarify and explain the difference, a large amount of X-ray occultation data from past, present, and future X-ray satellites will be required for further analysis. However, the atmospheric density retrieval method in this study depends on atmospheric models. The retrieved results vary with the shape of density profiles for different atmospheric models. Therefore, model independent retrieval methods needs to be developed, and we will consider this kind of XEOS method in the future.
X-ray photons with higher energy can penetrate deeper into the Earth's atmosphere. The observation data of HE and ME Cao et al., 2020) will be analyzed in our future work. In addition, we will investigate the factors for affecting the XEO lightcurve modeling and density retrieval, such as the extended X-ray source effects, the energy spectra variations, etc.. These effects for the the XEO lightcurve modeling and density retrieval will be analyzed in the next work.
Data availability. Datasets related to the paper are available from http://archive.hxmt.cn/proposal Author contributions. BL and HL were responsible for conceptualisation. HL was responsible for funding acquisition and supervision. HL was responsible for the forward model building, the retrieval algorithm and the original software development. DY and YT were responsible for the data reduction. DY, HL and YL were involved in the software update and data analysis. HL and MG were responsible for the design of the observations. XL and WX were responsible for the validation of the results. DY prepared the original draft. All co-authors reviewed and edited the paper.
Competing interests. The authors declare that no competing interests are present.
Disclaimer. Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.