A Study on the X-Ray Pulse Profile and Spectrum of the Crab Pulsar Using NICER and Insight-HXMT's Observations

We analyze the energy dependence of the X-ray pulse profile and the phase-resolved spectra (PRS) of the Crab pulsar using observations from the Neutron star Interior Composition Explorer (NICER) and the Hard X-ray Modulation Telescope (Insight-HXMT). We parameterize the pulse profiles and quantify the evolution of these parameters in the broad energy band of 0.4-250 keV. A log-parabola function is used to fit the PRS in 2-250 keV, and the curvature of the spectrum, i.e., the evolution of the photon index with energy, as represented by the parameter \{beta} of the log-parabola model, also changes with phase. The relation of \{beta} and phase has two turning points slightly later than those of the pulse intensity profile, where the values of \{beta} are the lowest, suggesting that the energy-loss rate of the particles is the lowest in the corresponding regions. A three-segment broken-power-law model is also used to fit those PRS. The differences between the hard spectral index and the soft ones have a distribution similar to that of \{beta}, confirming the fitting results of the log-parabola model, while the broken energies are generally higher in the region bridging the two pulses. We find anticorrelations between the spectral indices and the curvature of the log-parabola model fitting and a similar anticorrelation between the spectral indices and broken energies of the broken-power-law model fitting, suggesting a scenario where the highest-energy particles are produced in regions where radiation energy loss is strongest.


INTRODUCTION
The Crab pulsar (PSR B0531+21) is one of the most frequently studied objects in X-ray astrophysics. It has a period of about 33 ms, a spin down rate ofṖ = 4.2 × 10 −13 s s −1 and an energy loss power ofĖ = 5 × 10 38 erg s −1 . Although observational data on this source have been gathered for about 50 years, the physics involved in this pulsar has not yet been completely understood. The position and size of the radiation regions, the distribution of secondary pairs, the magnetic field, the inclination and viewing angles, etc., together determine the pulse profiles and spectra in different energy bands (Cheng et al. 1986a,b). Studying the pulse profile and spectrum variations will thus help to explore the links between the multi-wavelength radiation properties and the physical conditions on this pulsar.
The pulse profiles of the Crab pulsar at all wavelengths are dominated by two pulses, separated by about 0.4 pulse phase. Previous studies show that the exact pulse morphology varies as a function of photon energy (Eikenberry & Fazio 1997;Kuiper et al. 2001;Ge et al. 2012;Tuo et al. 2019): the intensities of both the second pulse and the bridge emission increase with energy relative to that of the first pulse in X-ray and followed by a downward trend above 1 MeV (Kuiper et al. 2001). However, some of these studies were either done in narrow energy bands that can not give quantitative results (Eikenberry & Fazio 1997;Tuo et al. 2019), or incomplete for only adopting the flux ratio of the two pulses to describe the evolution of the pulse profile (Kuiper et al. 2001).
Actually the information contained in the energy dependence of the pulse profile can be revealed in more detail by phase resolved X-ray spectroscopy. The overall X-ray spectrum of the Crab pulsar could be well described by a power-law function and its spectral index varies with energy (Massaro et al. 2000;Weisskopf et al. 2004;Mineo et al. 2006). The phase-averaged total spectrum of the Crab nebula and pulsar in 1-100 keV has a photon index Γ ∼ 2.1 (RXTE, BeppoSAX, EXOSAT, INTEGRAL/JEM X, Kirsch et al. (2005)), and a softer index of ∼ 2.20 − 2.25 above 100 keV (INTEGRAL/SPI/ISGRI, CGRO Mineo et al. (2006); Kuiper et al. (2001)). The pulsed spectrum has a more complex evolution (Kuiper et al. 2001) and the photon index is about 1.6 in 0.3-3.8 keV (Weisskopf et al. 2011), 1.81 in 5-60 keV, 1.91 in 15-250 keV, and 1.96 in 0.1-300 GeV (Yan et al. 2018). For the phase-resolved X-ray spectra, Pravdo et al. (1997) found a phase evolution of the spectral index across the X-ray pulse in a reverse S shape in 5-200 keV, i.e., the spectral index increases until the peak of the first pulse and then falls off, while the spectral index keeps increasing from the leading wing to the trail wing of the second pulse. Such a phase evolution trend holds albeit with different values in 0.3-3.8 keV (Weisskopf et al. 2011), 1-10 keV (Vivekanand 2021), 3-60 keV (Ge et al. 2012), 3-78 keV (Madsen et al. 2015), 11-250 keV (Tuo et al. 2019), 15-500keV (Li et al. 2019), 20-500 keV  and > 100MeV (Abdo et al. 2010).
A three-dimensional outer magnetospheric gap model was proposed to explain the pulse profile, phase-resolved spectra, and polarization of the Crab pulsar (Cheng et al. 2000;Takata et al. 2007;Tang et al. 2008). Specifically, in the model of Cheng et al. (2000), the X-ray photons at different phases are produced by the radiations from different height range of regions in the magnetosphere. Recently, particle-in-cell (PIC) simulations are developed to study the pulsar magnetosphere that is more complex than a simple dipole (Philippov & Spitkovsky 2014;Philippov et al. 2015a,b). Such PIC methods can therefore produce more complicated pulse profiles and have been successfully applied to the interpretation of Fermi observations of pulsars (Cerutti et al. 2016;Philippov & Spitkovsky 2018). As the X-ray spectrum reflects the energy distribution and radiation cooling of the high energy electrons, phase-resolved spectroscopy in a broad X-ray band can provide information for theoretical studies on the high energy particle production and cooling, as well as the magnetic field geometry of the pulsar magnetosphere.
In this work, we use observations from the Neutron star Interior Composition Explorer (NICER) and the Hard X-ray Modulation Telescope (Insight-HXMT) to analyse the pulse morphology and phase-resolved spectra (PRS) for the Crab pulsar, providing further constraints on models of pulsar X-ray emission. The organization of this paper is as follows: data processing and reduction are presented in Section 2, analysis methods in Section 3, results in Section 4, discussions on the physical implications of our results and a summary in Section 5 and Section 6, respectively. Throughout the paper, errors of the parameters are at the one standard deviation level.

Timing Ephemeris from Jodrell Bank Observatory
A 13 m radio telescope at Jodrell Bank Observatory monitors the Crab pulsar daily (Lyne et al. 1993), offering a radio ephemeris (denoted as JBE) that is used for the analyses of NICER and Insight-HXMT data. For monthly validity intervals, this ephemeris contains up-to-date parameters of the rotation frequency and its first two time derivatives. The database is available through the World Wide Web (http://www.jb.man.ac.uk/pulsar/crab.html). These JBE parameters, the pulsar coordinate R.A.=83 • .63322, decl.=22 • .01446(J2000) (Abdo et al. 2010), and solar system ephemeris JPL DE200 are used in this work.

NICER Observations and Data Reductions
The NICER is an International Space Station payload devoted to the study of neutron stars through soft X-ray timing in 0.2-12 keV. Its X-ray Timing Instrument (XTI) is an aligned collection of 56 X-ray concentrator optics (XRC) and silicon drift detector (SDD) pairs. Each XRC collects X-rays over a large geometric area from a roughly 30 arcmin 2 region of the sky and focuses them onto a small SDD. The SDD detects individual photons, recording their energies with good spectral resolution and their detection times to a ∼100 nanoseconds RMS relative to Universal Time (Prigozhin et al. 2016).
The NICER's observations used in this paper are from 2017 August 25 to 2019 November 26, but only those with available JBE are used, and the total exposure time is about 339 ks. The data reduction is carried out by using NICERDAS software version 7a included in HEAsoft version 6.27.2, with calibration database version 20200722.
Standard calibration processes as well as filtering to an entire NICER observation are performed by the command nicerl2. The cleaned photon events produced by the previous step are corrected to the solar system barycenter (SSB) and then used to fold pulse profiles in different soft X-ray energy bands. The cleaned events after corrected to SSB are also used to get the phase resolved spectrum. We extract the source spectrum with HEASoft XSELECT package (V.2.4e).

Insight-HXMT Observations and Data Reductions
The Insight-HXMT is China's first X-ray astronomy satellite. It has three main payloads: the High Energy X-ray telescope (HE, 20 − 350 keV, 5100 cm 2 ), the Medium Energy X-ray telescope (ME, 8 − 35 keV, 952 cm 2 ) and the Low Energy X-ray telescope (LE, 1 − 12 keV, 384 cm 2 ) (Zhang et al. 2020). Only the HE and ME observations of the Crab pulsar with available JBE are used in this work, because of their high time resolution (∼ 2 µs for HE and ∼ 20 µs for ME) and their complementarity to the NICER data. All data from these two instruments are reduced via data processing pipeline hpipeline included in HXMTsoft (v2.04). Events are corrected to SSB by using the command hxbary. All the data is collected in MJD 57992-59228 (2017 August 27 to 2021 January 14) and in the good time intervals (GTIs) recommended by pipeline. The selected observations have a total exposure of about 596 ks for HE and 758 ks for ME. Figure 1 shows the pulse profiles in different energy ranges, in which we denote the first pulse that is more prominent in soft X-rays as P1, and the second pulse as P2. To show the difference between the profiles more clearly, all the original profiles are first subtracted by their respective mean count rate in the off-pulse region (phase 0.6-0.8; Yan et al. (2018)) and then normalized by the peak flux of P1. After this process, the peak values of P1 for all the profiles equal to 1, as shown in Figure 1(a). The other normalization method is also used in this work: every profile in an energy band is subtracted by its off-pulse count rate and then divided by the mean photon count rate in the whole period. The profiles normalized by the second method make it easily to observe the overall energy evolution trend of the pulse profiles. The location of the maximum of P1 is aligned as phase about -0.01 by using JBE.

Normalization of Pulse Profiles
For observations from NICER, because of the high detector noise in lower energy, only profiles in 0.4-12 keV are considered, according to the suggestion at https://heasarc.gsfc.nasa.gov/docs/nicer/data analysis/nicer analysis tips.html. For the observations from Insight-HXMT, the HE profiles are in 30-250 keV and the ME profiles are in 8-30 keV.

Parameterization of the Pulse Profiles
In order to study the energy dependence of the pulses more comprehensively, we describe the shape of a pulse profile with five parameters, which are the ratio of the peak intensity of P2 to that of P1 (R I ), the separation of the two pulses (∆Φ), the full width at half maximum (FWHM) of P1 (W 1 ), the FWHM of P2 (W 2 ), and the ratio of the fluxes within the FWHM phase ranges (R F ) of the two pulses. Because the shapes including the peak phases of the two pulses change with energy, the evolution of R F with energy differs slightly from that of R I . The peak intensities and phase positions of the two pulses are obtained by using the empirical formula proposed by Nelson et al. (1970), as introduced in Ge et al. (2012).
To explore the exact profile evolution of the pulses from another aspect, we calculate the flux ratios between the left side and the right side of P1 and P2 in different energy bins. The four phase intervals from which the fluxes are calculated are: LP1 (from the left FWHM phase to the peak of P1), RP1 (from the peak of P1 to its right FWHM phase to the peak of P1), LP2 (from the left FWHM phase to the peak of P2), RP2 (from the peak of P2 to its right FWHM phase to the peak of P2).

Spectral Fitting
The XSPEC (v.12.10.0c) in the software package HEASOFT is used to fit the spectra created in Section 2. The PRS of the Crab pulsar extracted from the NICER observations were fitted by a power law model by Vivekanand (2021), and the similar spectral analyses of the Insight-HXMT data were also carried out by Tuo et al. (2019). In this work, we fit the PRS from NICER and Insight-HXMT data jointly, with the log-parabola model (Massaro et al. 2000) and the broken power law models, respectively.
The log-parabola model involves three components, wabs*constant*logpar, among which wabs accounts for the interstellar medium absorption whose column density N H is fixed as 0.36×10 22 cm −2 (Ge et al. 2012), constant is used to compensate the mis-match of the effective areas of the three instruments (i.e. XTI, ME, and HE), and logpar is a power-law with an index which varies with energy E (eq. 1) to show the detailed energy evolution of the spectrum F (E) (Massaro et al. 2000) in a broad band, where E 0 =1 keV, and α, β and the normalization value K are set as free parameters to fit. α corresponds to the photon index at 1 keV, and β is used to represent the curvature of the spectrum. The constant for the NICER spectra is fixed at 1, while for the ME and HE spectra, their constant are determined by jointly fitting the phase averaged spectra of the three instruments. As shown in Figure 2, the phase averaged spectra from these three instruments could be well fitted by the logpar model, and the obtained parameters are: constant equal to 0.623 ± 0.002 and 0.831 ± 0.004 for ME and HE, respectively, α = 1.581 ± 0.004, β = 0.109 ± 0.001, K = 0.471 ± 0.002, and the reduced χ 2 = 0.742 for degree of freedom of 1607. Specifically, the two constant values are used when we fit the PRS. We also used a two-segment broken power law function (bknpow ) and a three-segment broken function (bkn2pow ) to search for the turnover energy of a spectrum of the Crab pulsar. Again, a wabs model is used to account for the absorption of the interstellar medium and constant to adjust the mis-match between the instruments. Figure 3 shows the fitting results of a spectrum in phase 0.00-0.01 by using the wabs*constant*bknpow and wabs*constant*bkn2pow models. Compared with the bknpow model, bkn2pow gives a much better fit to the data and so the fitting results of bkn2pow model are discussed later in this work.
One spectrum is created for each phase bin with a size of 0.01 in phase -0.12 to 0.5. These spectra are fitted with the logpar model each. However, some of the spectra are merged into one when fitted with the bkn2pow model in order to improve the statistics. The background spectrum is created using all the events in phase range 0.6-0.8.  Figure 2. The joint fitting of phase average spectra for NICER and Insight-HXMT's observations by using the logpar model. The green points are the NICER data in 2-12 keV, the black ones are the HXMT-ME data in 8-35 keV, the red points are the HXMT-HE data in 27-250 keV. The upper and lower panels show the spectra and residuals in terms of sigmas with error bars, respectively.

Evolution of the Pulse Profiles with Energy
Our analyses present the evolution of the pulse profiles with energy much more quantitatively than in the previous studies. As shown in Figure 1, the peak intensity ratio and the FWHM of the two pulses all manifest an increasing trend with energy. The intensity of P1 decreases with energy and this trend reverses at phase about 0.035. The intensity of the "bridge" region between this phase and 0.4 increases with energy, and then the intensity starts to decrease with energy again, as displayed in Figure 1 Figure 4(a) shows that the intensity ratio of two pulses. R I increases monotonously and becomes higher than 1 above 120 keV. The separation between the two pulses decreases with energy in general, but might increase with energy below 0.8 keV (Figure 4(b)). The pulse widths shown in panels (c) and (d) increase with energy, albeit with different slopes. The flux ratio of two pulses (R F ) also increases with energy, from 2.08 ± 0.03 at 0.45 keV to 4.60 ± 0.14 at about 200 keV, with a slightly larger gradient than that of R I (panel(e)).
In order to quantify the evolutionary trends of the above parameters, we try to fit the data points with different models and find that the power law function Y = c 1 * E(keV ) c2 can give an acceptable fit, where Y represents different shape parameters, and c 1 and c 2 are coefficients of the evolution. All the data points from NICER and Insight-HXMT are fitted, and the parameters are listed in Table 1.
Flux ratios of the left side to the right side are also obtained for the two pulses respectively to represent the evolution of their shapes with energy, which are given in Figure 5. Combining with the results in Figure 4, we find that as the photon energy increases, (1) the flux of P1 relative to the entire flux drops, leading to a broader pulse and a lower height peak; the flux ratio between the left side and the right side of P1 decreases slightly; (2) the flux ratio of P2 relative to the entire flux rises, leading to a broader pulse and a higher peak height; the flux ratio between the left side and the right side of P2 increases obviously; (3) the flux of the bridge region relative to the entire flux rises; (4) the change trend of pulse intensity turns over at the phases about 0.035 and 0.4, which are not aligned with the peak positions of P1 and P2.
In order to compare our observational results with the theoretic prediction of the outer gap model in Zhang & Cheng (2002), flux ratios (R F,P 2/P 1 , R F,Bridge/P 1 ) are calculated according to the phase intervals used in their work, which are  (2002)'s model can produce qualitatively the observed trend, but their exact functional relationship between flux ratio and energy is not accord with the results in this work.

Phase-resolved Spectra
The fitted parameters of the PRS by the log-parabola model are displayed in Figure 7, and their values and fitting goodness are listed in Table 2. As can be seen in Figure 7(a), the phase resolved spectral properties of the two pulses are quite different: the spectral index α around P1 increases first and then falls off, with the maximum at the peak; and α keeps growing from the leading wing to the trailing wing of P2. In the bridge region, α has a down and up trend, and changes smoothly with phase. The values of the curvature parameter β for all the spectra are positive (  Figure 7(b) ), meaning that the spectral index increases with energy at all phases. β presents a "W" shape in phase and reaches the local minimum around two turning points of the pulse profiles. Both α and β in 2-250 keV show phase evolution trends similar to those in Ge et al. (2012), but with higher precision. The value of β in 2-250 keV is smaller than in 3-60 keV, indicating that the photon index changes more slowly at higher energies. The phase turning points  Figure 4. The pulse shape parameters of the Crab pulsar in X-ray. RI is the intensity ratio of two pulses. ∆Φ is the phase separation of two pulses. W1 and W2 represent the FWHM of P1 and P2, respectively. RF is used to denote the flux ratio of the two pulses within FWHM. All the shape parameters evolves with energy gradually, and the black dashed lines represent the joint fitting results of data points from the two satellites with a power law function respectively. of β in this work (gray phase regions in Figure 7(b)) are not aligned with those in Ge et al. (2012), probably due to the spectral complexity of this pulsar, because the energy ranges of the two studies are different. As shown in Figure 3, in phase interval 0-0.01, the fitting reduced χ 2 are 1.05, 1.23, 1.03 when using logpar, bknpow, and bkn2pow models, respectively. The two-segment broken power law model does not result in statistically acceptable fit to the broad band PRS, whereas a three-segment broken power law model can fit the spectra very well. The joint fitting parameters for PRS by using a three-segment broken power law are shown in Figure 8 and listed in Table  3. It can be seen that the three spectral indices (Γ 1 , Γ 2 , Γ 3 ) all change with phase and their evolution trends are similar, and distributions of the differences between them are consistent with the evolution of β with phase. E break,1 changes rather randomly, but E break,2 is positively correlated with β. Madsen et al. (2015) found that the 3-78 keV PRS extracted from the NuSTAR observations can be well fitted with a two-segment power-law function (bknpow ), with the break energy fixed at 11.7 keV or 13.1 keV. Comparing with our results, we found that these two break energies are comparable with E break,1 and the evolutions of the photon indices share the similar trend. According to the fitting results in Table 3, the reduced χ 2 is higher and closer to 1 in the phase intervals near P1 and P2 than that in the bridge intervals, it means the bkn2pow model is more suitable for the spectra fitting of P1 and P2, but the spectra in the other phase intervals with lower statistics could be over fitted with such a model. It should be noted that the NuSTAR results (Madsen et al. 2015) also support our results that the turnover or bending at high energies are around a few tens of keV. Some fitting parameters of logpar and bkn2pow models have clear correlations in certain phase intervals. As shown in Figure 9(a), the correlation between α and β are opposite on both sides of P1 before reaching the first phase turning point. Between the two phase turning points, there is a significant anti-correlation between α and β, but their functional relationship seem to be not the same for two groups of data, in phase 0.04-0.25 and 0.25-0.41 respectively (Figure 9(b)). At the right side of P2, the correlation between α and β is unclear, perhaps due to the limitation of statistics (Figure 9(c)). For the bkn2pow model, we only find a possible anti-correlation between the averaged indices ((Γ 1 + Γ 2 + Γ 3 )/3) and E break,2 in phase 0.03-0.40, and their exact functional relationship also seems to be not the same in the front and rear phase intervals of the bridge region, as shown in Figure 10.
Phase 0.03-0.18 Phase 0.18-0.40 Figure 10. The correlations between the average photon indices and E break,2 of the spectral model bkn2pow in different phase intervals. The evolution of the flux ratio of the two pulses, one of the most frequently used parameters to describe the pulse profile evolution in X-rays, is consistent with the previous work in Kuiper et al. (2001) and Tuo et al. (2019), but the current results are much more precise and cover a wider energy range of 0.4-250 keV. A multicomponent model in Massaro et al. (2006) is able to reproduce the pulse profiles in different energy bands and the evolution of the flux ratio between P2 and P1 (Kuiper et al. 2001). The flux ratio of two pulses in 0.4-250 keV in this work coincides with that in Kuiper et al. (2001) when using the same phase definitions for P1 and P2, so our results can be also explained by the multicomponent model in Massaro et al. (2006). Zhang & Cheng (2002) once fit the flux ratios of P2 to P1 and Bridge to P1 by using an outer gap model. However, as shown in Figure 6, there are obvious discrepancies between the predictions of their model and our results. The three-dimensional model in Zhang & Cheng (2002) is based on a rotating vacuum dipole, in which some processes are simplified or unreliable. For example, only photons emitting tangent to the local magnetic field lines are considered, but particles with large initial pitch angles may also emit in other directions. The magnetic inclination and viewing angles used in their work are from old theoretical works, which have not been well determined by observations and are therefore probably inaccurate. More over, the hard X-ray spectral indices they adopted are from old observations and with large uncertainties. These factors may lead to the discrepancies between model predictions and observational results. Therefore, more details need to be added and/or refined into the outer gap model of Zhang & Cheng (2002) to explain the new results in this work.
There are several studies on the spectrum of the Crab pulsar as a function of phase (Pravdo et al. 1997;Massaro et al. 2000Massaro et al. , 2006Mineo et al. 2006;Ge et al. 2012;Tuo et al. 2019;Li et al. 2019;Yeung 2020). When fitting the spectrum by using a power law model, the overall behaviors of photon indices are very similar, and indices increase clearly with energy over all the phase intervals, i.e., the higher the energy, the softer the spectra are. For the energy curvature trend of spectra which is described by β of the log-parabola model fitting, this work gives a more detailed phase evolution trend and obtains links between the spectrum curvature and the evolution of pulse profiles with energy, as shown in Figure 1

Constraints on Pulsar Radiation Models
We can use the energy dependence of the pulse profiles and the PRS to test the pulsar high energy emission models. One of the most popular models in the past is the two-gap outer gap model, proposed by Cheng, Ho, and Ruderman (Cheng et al. 1986a,b). In this model, the gap is "stratified" in energy, and lower energy radiation emanates from a region farther from the neutron star surface than higher energy radiation. The delay between the radio and X-ray pulses decreases significantly with increasing X-ray energy, which indicates that the softer X-ray photons originate from higher altitudes (Molkov et al. 2010). Qualitatively, the energy evolution trends of the pulse shape and spectral index of the Crab pulsar obtained in this paper all conform to the "stratified" structure of this outer gap model. The gap closer to the neutron star surface may cover a larger angular size (relative to the neutron star center) as indicated by the increasing widths of the two pulses with energy. Because the magnetic field in deeper magnetosphere is stronger and P2 becomes more prominent at higher energy, the particle generation (and thus synchrotron emission) in the region close to the neutron star surface and corresponds to P2 maybe more intense (Cheng et al. 2000). According to the energy evolution trend of the whole pulse profiles, the pulse intensity of P1 changes more distinctly than P2, and corresponding radiation area of P1 may be longer and narrower than P2.
In recent years, global PIC simulations were used to model the pulsar magnetosphere (such as Philippov & Spitkovsky (2014); Cerutti et al. (2016); Philippov & Spitkovsky (2018)). In these simulations, the high energy emissions are mainly produced in current sheets just beyond the light cylinder. Cerutti et al. (2016) modeled the pulse profiles and phase-averaged spectra produced by synchrotron radiation process at different viewing and magnetic inclination angles. When the viewing angle is ∼ 90-100 • and the magnetic inclination is ∼ 30-45 • , the model predicted spectra and pulse profiles of the Crab pulsar are well consistent with the observational results in energy band from soft X-ray to γ-ray as presented in Kuiper et al. (2001). The PIC simulations and modeling show that electrons and positrons are generated mainly within the magnetosphere. Near the pulsar surface, the electric potential is highest in the equatorial plane and decreases toward the poles for a quasi-aligned rotator (Goldreich & Julian 1969). The freshly generated positrons in the magnetosphere are accelerated outward to high energy in the reconnection regions outside the light cylinder, where electrons coming from farther in the equatorial plane are accelerated inward. Because the supply of outgoing positrons is much higher than the ingoing electrons, high energy radiation is mainly generated by positrons for a quasi-aligned rotator (Cerutti et al. 2016), and on the contrary, for an anti-aligned rotator, it is expected that the high energy radiation is generated mainly by electrons. The highest energy that a positron can acquire is determined by the accelerating electric field in the reconnection region, multiplied by the acceleration distance that is available in the undulating current sheet. As a consequence, the spectrum will be harder for a lower magnetic inclination angle. On the other hand, charged particles in the current sheet can screen the accelerating electric field (Philippov & Spitkovsky 2018). These combined factors imply that the high intensity of P1 in low energy band is generated by regions where positrons are more numerous and therefore electric yielding is more significant. In higher energy band, the intensity of P2 becomes higher because P2 is likely produced by regions where electric yielding is weak and positrons can be accelerated to higher energy. The spectra of the bridge component are harder than those of both P1 and P2, while its intensity is lower than P1 and P2. This behavior can also be explained if the bridge component is generated by regions where positron density is lower and the fraction of high energy positron is higher. Future comparison between simulations of the current sheet in different energy bands with the results in our paper could provide more detailed information about the particle acceleration and radiation in the reconnection regions. Weisskopf et al. (2011) reported that the soft X-ray spectral index has an abrupt change in phase 0.83-0.95, where transition structures of the optical and X-ray polarization properties were also detected (S lowikowska et al. 2009;Vadawale et al. 2018). We find that around that phase the break energies of the bkn2pow are the lowest. It is yet unknown whether these properties are physically associated. Besides this, just behind the two pulses the polarization degree has two local minimum points, where the values of β are also the lowest. We speculate that in the two corresponding regions the particles evolve the least due to the lowest magnetic field, which could be less regulated and thus the radiations are less polarized.

SUMMARY
In this work, we adopt the observations from NICER and Insight-HXMT, and study the pulse profile and spectrum of the Crab pulsar systematically. The main results and conclusions are: (1) in 0.4-250 keV, the intensity/flux ratio of two pulses show an increasing trend with energy, so do the widths of the two pulses; the separation of the two pulses decreases with energy; the energy evolution trends of the two sides of P1 and P2 are asymmetric; (2) in 2-250 keV, where the spectral shape is not significantly affected by the chemical composition of the absorption interstellar medium, the spectral index increases with energy over all the phase intervals; the spectral index changes more slowly at higher energies; the change of spectral index reaches its local minimum around phase 0.035 and 0.4 that are later than the peak phase of P1 and P2; there are obvious energy breaks on the spectra near the two peaks of the Crab pulsar, and the break energies are generally higher in the region bridging the two pulses; (3) the energy evolution trends of pulse profile and spectrum conform to results in previous studies qualitatively, but more accurate to allow for a detailed analysis; (4) the anti-correlations among spectral fitting parameters (α and β) indicate that the highest energy particles are produced in regions where the radiation energy loss is also the strongest. We simply compared our results with the predictions of the outer gap model and the PIC simulations of pulsar magnetosphere. It turns out that the theoretical models need to be further refined to explain the details of the observational results.