© The Authors 2024

1. Introduction

Cen X-3, a persistently accreting X-ray pulsar, was first identified by the Uhuru X-ray Explorer Satellite (Giacconi et al. 1971; Schreier et al. 1972). Cen X-3 is an eclipsing high-mass X-ray binary (HMXB) with a spin period of approximately 4.8 s and an orbital period of around 2.1 days (also see Klawin et al. 2023; Liu & Wang 2023). The system exhibits a decay rate of ∼1.799 × 10⁻⁶ yr⁻¹ (Falanga et al. 2015). The binary system consists of a neutron star with a mass of approximately 1.21 ± 0.21 M_⊙ and an optical companion O6-8 III supergiant star (Krzeminski 1974; Hutchings et al. 1979) with a mass of about 20.5 ± 0.7 M_⊙ and a radius of 12 R_⊙ (Ash et al. 1999). The estimated distance to the binary system is approximately 6.9 kpc based on Gaia measurements (van den Eijnden et al. 2021; Gaia Collaboration 2021). The high X-ray luminosity (∼5 × 10³⁷ erg s⁻¹, Suchy et al. 2008), a secular spin-up trend in this system, and the quasi-periodic oscillations at 40 mHz (Raichur & Paul 2008; Takeshima et al. 1991; Liu et al. 2022b) may suggest the existence of an accretion disk.

Cyclotron resonant scattering features (CRSFs, Staubert et al. 2019), also called cyclotron lines, are identified as absorption features in the continuum spectrum of X-ray pulsars. In highly magnetized accreting neutron stars, electrons perpendicular to the strong magnetic fields are quantized into discrete Landau levels. These absorption-like features form through the resonant scattering of photons on these electrons in the line-forming region. Measurements of the energy of the cyclotron line provide a direct method for determining the magnetic strength of neutron stars, because the cyclotron line energy is linked to the magnetic field, which can be expressed as follows:

$$\begin{array}{c}\hfill E_{\mathrm{cyc}}\approx \frac{11.57\mathrm{keV}}{(1+z)}{B}_{12},\end{array}$$(1)

where B₁₂ represents the magnetic field strength of the neutron star in the scattering region in units of 10¹² G, and z denotes the gravitational redshift. Santangelo et al. (1998) estimated the magnetic field of the neutron star surface for Cen X-3, finding it to be in the range of (2.4–3.0) × 10¹² G using E_cyc ∼ 28 keV.

The X-ray spectra in Cen X-3 are typically modeled by a power law with a high-energy cut-off, an iron emission line, and a soft excess below 1 keV (White et al. 1983; Burderi et al. 2000; Suchy et al. 2008). An absorption line feature around 30 keV in the continuum spectra was first observed by Nagase et al. (1992) using Ginga observations. The presence of the CRSF was later confirmed by Santangelo et al. (1998) through BeppoSAX satellite observations. Subsequently, the dependence of CRSF parameters on pulse phase in Cen X-3 was explored by Burderi et al. (2000) with BeppoSAX and by Suchy et al. (2008) with Rossi X-ray Timing Explorer (RXTE). New measurements of the cyclotron line energy were also provided by Tomar et al. (2021) using Suzaku and NuSTAR. Tamba et al. (2023) also presented a detailed time-resolved analysis with NuSTAR observation data from 2015. In the present paper, we further study the timing and spectral properties of Cen X-3 using new NuSTAR observations from 2022, which cover about two orbital periods.

The paper is organized as follows: In Sect. 2, we briefly describe the NuSTAR observations and data reduction. We show the temporal analysis in Sect. 3. We present our analysis of the spectral fitting in Sect. 4, showing the evident cyclotron lines. In Sect. 5, we focus on the dependence of CRSF energies on luminosity. In addition, the orbital dependence of CRSF and continuum parameters are discussed. We summarize our conclusions in Sect. 6.

2. Observations and data reduction

The NuSTAR observations were carried out by Harrison et al. (2013) from January 12 to January 16, 2022, with a total exposure time of 189 ks (ObsID: 30701019002). We processed the data with the standard NuSTAR analysis software in HEASoft 6.32.1. The source regions were extracted as circular regions with a radius of 180″ from the source center, and the background regions were extracted with an identical area in the off-source regions. We then applied the Solar System barycentric correction to the cleaned event files.

Light curves were extracted with a bin size of 10 ms, and spectra were generated from the barycentric corrected event files using the task nuproducts. Figure 1 shows the counts rate in 3–78 keV, which cover about two binary orbits in different intensity states. Here, we excluded the eclipse data and divided the observation into 23 different intervals based on the duration of ∼10 ks in preparation for subsequent analysis; these intervals are shown in Fig. 1. Before the temporal analysis, the shift in arrival time due to binary orbital motion was corrected for each light curve using the orbital parameters from Fermi/GBM¹. All of the spectral analyses in this paper were performed with XSPEC 12.13.1(Arnaud 1996). The uncertainties are calculated as 90% confidence levels using the Markov chain Monte Carlo (MCMC) method with ten walkers and a chain length of 10 000.

Fig. 1. FPMA and FPMB count rates at 3–78 keV of NuSTAR in 2022 with a time resolution of 100 s.

3. Temporal analysis

To find the spin period of the Cen X-3 at each time interval, we searched for the pulsations using the epoch folding technique with the efsearch tools. The best pulse period can be determined using the maximum χ² value. The pulse profile was generated by folding each light curve by the best spin period. We plotted the pulse profiles in 3–78 keV energy bands for 23 time intervals. As seen in Fig. 2, the first panel is the pulse profile from the first time interval (here the middle time at MJD 59591.2); the second panel comes from intervals 2 to 13 in the first orbital phases (∼MJD 59591.7-59593.3) with a moderate luminosity; and the third panel comes from intervals 14 to 23 over the second orbit (∼MJD 59593.8-59595.1) with high luminosity. We find that the shape of pulse profile did not evolve over the orbital phase and remained stable. These pulse profiles show obvious double-peaked shapes, sometimes with a narrow and small peak at the pulse phase of ∼0.1 in the rising phases of the main pulse. There is no significant evolution of the pulse profile within one orbit, but the second or faint peak of the pulse in the second orbit appears to be stronger than the first orbit, maybe due to the higher luminosity. We also present the pulse fraction versus the corresponding flux and find that the pulse fraction increases as the flux rises, as seen in Fig. 3. Figure 4 shows the hardness ratio in different energy bands over the 4–78 keV luminosity. The hardness ratios of different energy ranges show a negative correlation with the total luminosity, and at a higher rate (luminosities > 3 × 10³⁷ erg s⁻¹) they change, decreasing slowly before reaching a plateau where they remain stable.

Fig. 2. Pulse profiles of Cen X-3 in 3–78 keV from the summations of the FPMA and FPMB counts. These profiles are aligned by defining the minimum as phase zero. The first panel is from the first time interval. The second panel comes from intervals 2 to 13 in the first orbit with a moderate luminosity. The third panel comes from intervals 14 to 23 in the second orbit with high luminosity.

Fig. 3. Pulse fraction (PF = (Imax − Imin)/(Imax + Imin)) in energy bands of 3–78 keV of Cen X-3 as a function of luminosity.

Fig. 4. Hardness ratios of Cen X-3 as a function of the 4–78 keV luminosity in different energy bands.

To further investigate the energy-resolved pulse profiles, we focused on the time interval k of the first orbit (MJD 59593.12) at low states and another interval g of the second orbit (MJD 59594.67) at high states, and then folded the light curves into 14 different energy bands of 3–4 keV, 4–5 keV, 5–6 keV, 6–7 keV, 7–9 keV, 9–11 keV, 11–13 keV, 13–15 keV, and 15–20 keV, 20–25 keV, 25–30 keV, 30–35 keV, 35–45 keV, and 45–78 keV. Figure 5 shows the energy-dependent pulse profiles over different energy bands for the two time intervals: that is, k in the low state (left two panels) and g in the high state (right two panels). Only the energy bands below 45 keV show clear pulsations. The double-peaked profile was detected in a low energy range of 3–15 keV, and the secondary peak decreases as energy increases. At a higher energy range above 15 keV, the pulse has a single-peaked shape. An additional narrow peak in the pulse phase of ∼0.1 at the rising stages of the main pulse is seen below 20 keV. Similar energy-dependent pulse profiles have been reported (Nagase et al. 1992; Burderi et al. 2000; Suchy et al. 2008; Shirke et al. 2021; Tomar et al. 2021; Liu & Wang 2023).

Fig. 5. Pulse profiles of Cen X-3 in different energy bands. The left two panels show the energy-resolved pulse profiles for Cen X-3 from the summations of the FPMA and FPMB counts in MJD 59593.12 (the time interval k of the first orbit) at low intensity state and the right two panels show the energy-resolved pulse profiles in MJD 59594.67 (the time interval g of the second orbit) at high intensity state.

In addition, we find multiple subtle peaks in the pulse profile in the energy band of 3–5 keV at a high intensity state. Also, the main peak is divided into two subpeaks, which is particularly different from the behavior of the pulse profile seen in previous observations. The evolution from a double-peaked profile at low states to a triple-peaked profile at high states may indicate the transition from a subcritical to a critical regime. It is also worth mentioning that the maximum of the main pulse at the high states in the energy band of 3–5 keV shows an obvious shift of about 0.1 compared to other energy ranges, which may indicate the emissions come from both “pencil” and “fan” beams. Regarding the asymmetric shape of the pulse profile, this may be explained as being due to contributions from two emission regions displaced from the dipole geometry by approximately 10°, as proposed by Kraus et al. (1996). Recently, X-ray polarization research on Cen X-3 performed by the Imaging X-ray Polarimetry Explorer (IXPE) (Tsygankov et al. 2022) showed that the pulsar geometry parameters and the phase-resolved polarimetric properties are in good agreement with the pulse-profile decomposition by Kraus et al. (1996). Saathoff et al. (2024) introduced the phase-correlated variability analysis (PCVA) method and successfully decomposed the X-ray pulsar profiles of Cen X-3, but their results can produce asymmetric pulses, which are not aligned with those of Kraus et al. (1996), who assume that the emissions at the two poles are symmetrical. We expect the results of Saathoff et al. (2024) could explain the observed splitting of the main peak in our study. The complex emission profile revealed by our analysis is probably a mixture of the pencil and fan beam pattern. Tamba et al. (2023) suggested that a double-peaked shape at low energy and a single-peaked shape at high energy may indicate the presence of both pencil and fan beam emissions. Thus at higher energy states, the beam patterns probably become dominated by fan beam and a triple peak is expected, which is consistent with our results at low energy. Figure 6 shows the energy-dependent pulse fraction –defined here by the expression (I_max − I_min)/(I_max + I_min)– in energy bands from 3 to 78 keV for the two intervals above. In the range of 3–78 keV, the pulse fraction in both states increases as the energy increases, except for the band of 6–7 keV, which corresponds to the 6.4 keV iron emission line.

Fig. 6. Pulse fraction as a function of energy in 59593.12 MJD and 59594.67 MJD, respectively.

4. Spectral analysis

The X-ray spectra in Cen X-3 are usually characterized by a power law modified by an exponential cut-off, the iron emission lines, a soft excess, and the CRSF around 30 keV. For the X-ray spectral fitting, there are some analytic functions that are commonly used to describe the broadband continuum. In order to obtain a suitable or good description of the data, we used these empirical models and compared the fits among them. First, we applied a common continuum spectrum of a power law with a Fermi-Dirac cutoff (FDCUT) – also adopted by Suchy et al. (2008), Tamba et al. (2023), and Liu et al. (2023) – to fit the spectra in Cen X-3, which is characterized by

$$\begin{array}{c}\hfill C(E)=N·E^{-\mathrm{\Gamma }}{[1+exp\left(\frac{E-E_{\mathrm{cutoff}}}{E_{\mathrm{fold}}}\right)]}^{-1},\end{array}$$(2)

where Γ is the power-law index, and E_cutoff and E_fold are the cutoff and folding energies, respectively. Another choice is the power law with a high-energy cutoff (HIGHECUT) previously employed for Cen X-3 (White et al. 1983; Burderi et al. 2000), which has the form

$$\begin{array}{c}\hfill C(E)=\{\begin{array}{cc}N{E}^{-\mathrm{\Gamma }}\hfill & \mathrm{for}E\le E_{\mathrm{cutoff}},\hfill \\ N{E}^{-\mathrm{\Gamma }}\times exp[-(E-E_{\mathrm{cutoff}})/E_{\mathrm{fold}}]\hfill & \mathrm{for}E>E_{\mathrm{cut}}.\hfill \end{array}\end{array}$$(3)

Although C(E) has a discontinuous derivative around E_cutoff, which can result in wedge-shaped residuals around 15 keV, it is extensively used and adequate for modeling spectra. Finally, the sum of positive and negative power laws multiplied by an exponential cutoff (NPEX) (Mihara 1995) is expressed as

$$\begin{array}{c}\hfill C(E)=(N_1·E^{-{\mathrm{\Gamma }}_1}+N_2·E^{+{\mathrm{\Gamma }}_2})exp(-E/E_{\mathrm{fold}}),\end{array}$$(4)

where Γ₂ is usually set to 2. The photoelectric absorption $f(N_{\mathrm{H}}^{\mathrm{IM}})$ (phabs in XSPEC) component and iron line emission feature around 6.4 keV f_gauss (Gaussian in XSPEC) are also included in our models, but the soft excess is not included in the model because this contribution can be neglected. As the hydrogen column density $N_{\mathrm{H}}^{\mathrm{IM}}$ for the interstellar medium cannot be constrained well in the energy above 4 keV, we fixed it to 1.1 × 10²² atoms cm⁻² according to the galactic hydrogen column density in the same direction as the source. For the additional absorption f(N_H) for the stellar wind or the local materials, we tried the photoelectric absorption (phabs) and partial covering absorption (pcfabs in XSPEC) models using the abundance of angr (Anders & Grevesse 1989) and vern cross-sections (Verner et al. 1996). The χ² fitting statistics are almost the same for the two absorption models applied to the NPEX continuum. However, for the HIGHECUT and the FDCUT continua, the pcfabs model performs better, especially at low luminosities (<  ∼ 3 × 10³⁷ erg s⁻¹). At higher luminosities, the pcfabs parameters are not well constrained for the FDCUT model in some cases, and so we used the phabs absorption models instead.

To model the CRSF, we used a multiplicative absorption line model with a Gaussian optical depth profile (gabs in XSPEC):

$$\begin{array}{c}\hfill \mathrm{gabs}(E)=exp(-\frac{d_{\mathrm{cyc}}}{\sqrt{2\pi }{\sigma }_{\mathrm{cyc}}}{e}^{-0.5{[(E-E_{\mathrm{cyc}})/{\sigma }_{\mathrm{cyc}}]}^2}),\end{array}$$(5)

where E_cyc is the cyclotron line energy, d_cyc is the line depth, and σ_cyc is the width. We also leave the normalization constant free to vary between the FPMA and FPMB spectra. Therefore, the total model adopted here is

$$\begin{array}{c}\hfill S(E)=\mathrm{constant}\times f(N_{\mathrm{H}}^{\mathrm{IM}})\times f(N_{\mathrm{H}})\times (C(E)\ast \mathrm{gabs}(E)+f_{\mathrm{gauss}}).\end{array}$$(6)

The spectra were fitted in the energy range of 4–78 keV. We grouped the spectra from both FPMA and FPMB with a minimum of 25 counts per bin.

As an example of a specific time interval g of the second orbit (59594.67 MJD), Fig. 7 shows the best-fitting spectra and the residuals between the spectra and the model using the FDCUT continuum. The residuals around ∼28 keV after removing the CRSF from the best fit are clearly seen in panel b of Fig. 7. Panel c also shows the residuals of the spectral fitting without the CRSF, which cannot give good fits for the spectra, suggesting that the CRSF is also statistically required. The reduced χ²/d.o.f. of the fitting without and with the gabs component are 1.88/1327 and 1.09/1324, respectively. Table 1 lists all the best-fit parameters of this spectrum. The centroid energies of CRSFs for the best fitting are found at $27.{91}_{-0.44}^{+0.37}$ keV, and the iron emission line is around 6.4 keV. Table 1 also lists the fitting results for the HIGHECUT continuum and NPEX continuum models. We can see that the NPEX model is of slightly better quality (with a smaller reduced χ²) than the FDCUT and HIGHECUT models. The cyclotron line energy does not significantly change for these continuum models, except that the width and depth of the HIGHECUT model are smaller than those of the other two models.

Fig. 7. X-ray spectrum of Cen X-3 ranging from 4–78 keV for the interval g of the second orbit using the Fermi-Dirac cutoff continuum. Panel a: residuals of the best-fitting model. Panel b: residuals of the best model, setting the strength of the 30 keV cyclotron line to zero. Panel c: residuals of the spectral fitting without the CRSF.

5. Results and discussion

In this paper, we present a detailed timing and spectral analysis of Cen X-3 with observations made by NuSATR of two binary orbits covering two intensity states. Here, we focus on the dependence of the continuum and CRSF parameters on the luminosity and spectral variability during the orbital phases.

5.1. Dependence of spectral parameters on luminosity

We performed fitting on all 23 of the time intervals with the FDCUT, HIGHECUT, and NPEX continuum models. From all the fitting results, we first investigated the dependence of the continuum parameters on luminosity in Cen X-3. The luminosity is calculated from the fitted unabsorbed flux of the spectra based on the source distance of ∼6.9 kpc (van den Eijnden et al. 2021). As depicted in the left plots of Fig. 8, for the results derived from the FDCUT model (the blue points), the photon index Γ, varying from 0.8 to 1.3, exhibits a positive correlation with luminosity; the Pearson correlation coefficient is about 0.86 with a probability (the probability of the observed data given that the null hypothesis of no correlation is true) of 1.54e−07, indicating that the continuum spectrum tends to soften at a high accretion rate. Beyond the luminosity of ∼3 × 10³⁷ erg s⁻¹, the photon index increases slowly and remains flat; such spectral variability can be naturally ascribed to the changes in radiation pressure. The clustered Γ in the range >  ∼ 3 × 10³⁷ erg s⁻¹ may indicate that the accretion regime was close to the critical luminosity. Becker & Wolff (2022) suggested that their new model can explain the flatter, harder spectrum with increasing luminosity in Her X-1. Based on the amount of PdV work on the radiation field, these authors propose that in high-luminosity (supercritical) sources, the enhanced compression can lead to a flatter continuum. Conversely, a softer spectrum may indicate a subcritical source. The cut-off energy ranges from 20 to 26 keV and is positively related to luminosity; it begins to decrease around 3 × 10³⁷ erg s⁻¹. Such a high-energy quasi-exponential cutoff can be explained by the theoretical work of Becker & Wolff (2007), who successfully applied their model to reproduce the phase-averaged X-ray spectra of Cen X-3. The fold energy varies from 4 to 6 keV, which is a signature for the electron temperature in the case of nonsaturated Comptonization. The electron temperature of 4.77 keV estimated by Thalhammer et al. (2021), who successfully applied the accretion column model of Becker & Wolff (2007) to Cen X-3, is very close to the E_fold. The fold energy shows a decreasing trend at the beginning, followed by an increase, with a sudden decrease around 3 × 10³⁷ erg s⁻¹. The hydrogen column density N_H varies from 1 to 85 × 10²² atoms cm⁻², exhibits a negative correlation with the luminosity, and then remains stable, which is consistent with the results of Tamba et al. (2023), suggesting that the variations in absorption may affect the luminosity. The covering fraction decreases from ∼0.8 to ∼0.5 with large fluctuations. These spectral continuum parameters of the FDCUT model were found to show differences near 3 × 10³⁷ erg s⁻¹. The best fit for the FDCUT model yields a reduced χ² of around 1.1. Regarding the HIGHECUT model (the black points in Fig. 8), the behaviors of the photon index Γ and N_H are comparable with those seen when applying the FDCUT model. However, the fitting quality of HIGHECUT is worse at high luminosity and is not adequate to describe the data well.

Fig. 8. Spectral continuum parameters (photon index Γ, cutoff energy Ecutoff, fold energy Efold, hydrogen column density NH, and reduced χ2) as a function of X-ray unabsorbed luminosity. The left plots are the fitting results from the FDCUT model (the blue points) and the HIGHECUT model (the black points), and the right plots show the parameters from the NPEX model.

The right plots of Fig. 8 present the continuum parameters (the red points) for the NPEX model, which can also provide a good description of the data. Instead, the photon index varies between 1 and 3 and shows a negative correlation with luminosity, which is likely due to the evolution of E_fold and its effect on the continuum. The hydrogen column density decreases with increasing luminosity, consistent with the results obtained from FDCUT and HIGHECUT continuum models. At higher luminosities (> 3 × 10³⁷ erg s⁻¹), both Γ and N_H remain stable. Moreover, the folding energy also displays differences below and above this luminosity. These changes in spectral continuum parameters may indicate transitions between accretion regimes, such as from subcritical to critical. To explore this potential transition around 3 × 10³⁷ erg s⁻¹, we further investigated the dependency of CRSFs on luminosity.

The centroid energy of CRSFs in some X-ray pulsars has been observed to vary with X-ray luminosity, as discussed in a review by Staubert et al. 2019. The discovery of a positive correlation between CRSF energy and luminosity at low and moderate luminosity was first made by Staubert et al. (2007) when studying Her X-1. Similar trends have been identified in other sources, such as Vela X-1 (Fürst et al. 2014; La Parola et al. 2016; Liu et al. 2022a), Cep X-4 (Vybornov et al. 2017), and V0332+53 (Doroshenko et al. 2017; Vybornov et al. 2018). Conversely, at higher luminosities, a distinct negative correlation was observed in V 0332+53 (Makishima et al. 1990; Tsygankov et al. 2006). Furthermore, 1A 0535+262 (Kong et al. 2021) has been confirmed to exhibit an anti-correlation above the critical luminosity and a positive correlation below the critical luminosity.

As shown in Fig. 9, we present the CRSF parameters versus luminosity based on the FDCUT (the blue points) and NPEX (the red points) models. The width of the CRSF varies between 6 and 8 keV, and is correlated with the cyclotron line energy with a Pearson correlation coefficient of 0.69, and a probability of ∼0.0002 for the FDCUT model, and a Pearson correlation coefficient of 0.57 and a probability of ∼0.0043 for the NPEX model. The depth of the CRSF first shows an increasing trend and then decreases as the luminosity increases. There is probably a complex variation between X-ray luminosity and the CRSF energies for Cen X-3. For both continuum models, below ∼2 × 10³⁷ erg s⁻¹, there appears to be a decreasing trend between CRSF energy and luminosity, while above it, the CRSF energy begins to rise. However, at higher luminosities (>  ∼ 3 × 10³⁷ erg s⁻¹), a decreasing trend is seen. To further investigate the correlation between CRSF energy and luminosity around 3 × 10³⁷ erg s⁻¹, we fitted the relationship of E_cyc − L_x of the FDCUT model in the range of ∼2 − 3 × 10³⁷ erg s⁻¹ and above the luminosity of 3 × 10³⁷ erg s⁻¹ with a linear function separately. For the luminosity of 2–3 × 10³⁷ erg s⁻¹, we fitted the relation with E_cyc = (1.33 ± 0.34)×(L_x/10³⁷)+(24.50 ± 0.87) keV with a Pearson’s correlation coefficient r of 0.84 and a probability of around 0.0083, which shows a positive correlation. Above 3 × 10³⁷ erg s⁻¹, we obtained the relation of E_cyc = ( − 0.69 ± 0.24)×(L_x/10³⁷)+(30.40 ± 0.83) keV with a Pearson correlation coefficient of about −0.76 and a probability of about 0.028. For the NPEX model, the evolution of E_cyc − L_x is similar to the FDCUT model and gives comparable results, illustrating that the choice of continuum model does not significantly affect the relation between CRSF energy and luminosity. Therefore, the variation of E_cyc − L_x shows a possible positive and negative correlation below and above the transitional luminosity of ∼3 × 10³⁷ erg s⁻¹, respectively. However, we are unable to conclude that there is a distinct correlation or anticorrelation between the energy and luminosity of the CRSFs due to the large uncertainties of the centriod energy. Müller et al. (2011) reported that there is no transition in the accretion regime in Cen X-3. Recently, Liu et al. (2023) demonstrated that the CRSF energies have no correlation with luminosity in a larger range of ∼(5 − 9)×10³⁷ erg s⁻¹ based on Insight-HXMT observations obtained in 2022. Dangal et al. (2024) claimed an anticorrelation between CRSF energy and flux in their work with the same NuSTAR data as that analyzed in the present study, while using a continuum with a “newhcut” cutoff. The differences between our findings and those of Dangal et al. (2024) might stem from their choice of continuum and specifically their inclusion of a second absorption feature, for which we see no need in our analysis. Therefore, the relationship between cyclotron line energy and luminosity is still not clear. In the following, we discuss the possible critical luminosity of Cen X-3 based on different theories.

Fig. 9. CRSF parameters (energy, width, and depth) as a function of X-ray luminosity and the centroid energy versus the width of CRSF. These parameters are derived from the FDCUT model (the blue points) and the NPEX model (the red points), respectively.

A transition between subcritical and supercritical accretion regimes was clearly observed in the transient X-ray binaries V 0332+53 and 1A 0535+262 (see Mushtukov & Tsygankov 2022). As predicted by Becker et al. (2012) and Mushtukov et al. (2015a), in the supercritical case –that is, above a certain critical luminosity L_crit–, higher accretion rates lead to higher radiation pressure, and the radiation-dominated shock decelerates the falling material, forming the accretion column (Basko & Sunyaev 1976), and the emission height in the column increases with luminosity. Below the critical luminosity (subcritical), the Coulomb interaction dominates the deceleration, and the height of the line-forming region reduces with luminosity. A change of height is associated with a variation in the local magnetic field strength. Therefore, at high luminosity, the CRSF energy is expected to decrease with the luminosity. Nevertheless, at low accretion rates, a positive correlation is expected. This positive correlation below the critical luminosity may also result from the Doppler effect proposed by Mushtukov et al. (2015b). The transition from the subcritical to the critical regime may be related to the change in the emission geometry from a mixture of pencil and fan beams to only fan beams, consistent with change in the shape of the pulse profiles at low states to high states in Fig. 5. Becker et al. (2012), based on their theoretical work, show that critical luminosity can be expressed as

$$\begin{array}{c}\hfill \begin{array}{c}\hfill L_{\mathrm{crit}}\simeq 1.49\times {10}^{37}\mathrm{erg}{\mathrm{s}}^{-1}{\left(\frac{\mathrm{\Lambda }}{0.1}\right)}^{-7/5}\times B_{12}^{16/15},\end{array}\end{array}$$(7)

for M_NS = 1.4 M_⊙ and R_NS = 10 km. We assume Λ = 0.1 for the case of disk accretion, and compute the surface magnetic field strength B₁₂ using the fundamental energy as the surface value for CRSF energy. We estimate the expected L_crit to be equal to ∼4 × 10³⁷ ergs⁻¹, which is approximately the observed transitional luminosity of 3 × 10³⁷ erg s⁻¹. In our calculations of L_crit, we neglect the effect of gravitational redshift (e.g., a factor of (1 + z)⁻²), which introduces additional uncertainties in the theoretical estimates when comparing to observations. If we adopt the distance of 8.0 kpc (Krzeminski 1974), the observed transitional luminosity can be up to 4 × 10³⁷ erg s⁻¹, which is also consistent with the expected L_crit. According to the theoretical value estimated by Mushtukov et al. (2015a), the critical luminosity is around 10³⁷ erg s⁻¹ for l₀/l= 1.0 and the pure X-mode polarization, which is a little smaller than the value derived by Becker et al. (2012). Although the expected negative correlation between the CRSF energy and the luminosity above the 10³⁷ erg s⁻¹ is not seen, these theoretical results provide a good description of the accretion process, and more high-quality data are needed to verify these models in the future. Based on the theoretical critical luminosity and also considering the different behaviors of the shapes of the pulse profiles and hardness ratios from the low rate to the high rate, it is possible that the critical luminosity is contained within the range of ∼(0.5 − 4)×10³⁷ erg s⁻¹.

5.2. Spectral variability over orbital phases

Cen X-3 has a relatively short orbital period, and the spectral parameters probably depend on the orbital phases. Here we present spectral parameters over the two binary orbits observed in 2022 based on the NPEX model fitting results in Sect. 4, which cover two intensity states. Figure 10 shows the results of the evolution of CRSF energies and continuum parameters over the two orbits.

Fig. 10. NPEX fits showing the dependence of the CRSFs, the fold energy Efold, the photon index Γ, and the absorption column density NH on the orbital phases for the first binary orbit and second orbit, respectively. The orbital phases are calculated from the orbital parameters provided by Fermi/GBM.

Over the first binary orbit (∼(1.0 − 2.5)×10³⁷ erg s⁻¹, the left plots of Fig. 10), the CRSF energy varies between 27 keV and 29 keV. The evolution of the width and depth of the CSRFs is similar to the variation in their energy. The fold energy E_fold also roughly follows the CRSF energy. In addition, the behavior of the photon index Γ from 1.0 to 2.0 is opposite to luminosity. The photoelectric absorption N_H shows a large variation between 10 × 10²² atoms cm⁻² and 20 × 10²² atoms cm⁻² over the orbit. We also note that N_H slightly decreases near the beginning of the orbit and increases after that. The increase in N_H near an eclipse may be due to a larger amount of stellar material within the line of sight (Clark et al. 1988). An increase in N_H around the phase of 0.3–0.4 is observed, which is probably due to the bow shock, as mentioned by Suchy et al. (2008).

Over the second orbit (∼(2.5 − 4)×10³⁷ erg s⁻¹; see the right plots of Fig. 10), the variations in CRSF energies of 28–29 keV are almost opposite to the evolution of luminosity, except for the two points in luminosity of < 3 × 10³⁷ erg s⁻¹. Furthermore, the widths and depths of the CRSFs, and the fold energy E_fold, roughly follow the shape of the energy of the CRSFs. Here, the Γ from 0.1 to 1.0 has a smaller value compared with the first orbit and also shows an anticorrelation with luminosity. The variability of the absorption column N_H around the value of 1–10 × 10²² atoms cm⁻² is similar to the first orbit and also shows a small increase around phase 0.3.

Since the highly variable N_H probably originates from absorption by the stellar wind, we considered the possible stellar wind structure other than the accretion stream. When the majority of the material in the wind passes through our line of sight, variations in absorption are expected. In general, N_H increases after the orbital midpoint, as predicted by the simple wind model, which is consistent with our results. For example, enhanced absorption is seen in Vela X-1 (Haberl & White 1990; Diez et al. 2022), which shows a strong stellar wind and tidal stream (Blondin et al. 1991).

6. Conclusion

We carried out a detailed temporal and spectral analysis of the accreting X-ray pulsar Cen X-3 using NuSTAR observations obtained in 2022. The pulse profiles exhibit remarkable stability throughout the orbital phase. These profiles demonstrate energy dependence, displaying double peaks in the energy band below 15 keV and a single peak in the higher energy band, with the pulse fraction increasing with energy. Notably, the profile in the energy band of 3–5 keV during a high-intensity state exhibited marginal triple peaks, with the main peak further divided into two subpeaks. Additionally, we highlight a clear positive correlation between pulse fraction and flux.

From our spectral analysis, we identify a prominent cyclotron line with centroid energies ranging from 26 keV to 29 keV using different continuum models, and the spectra in Cen X-3 are well described by the FDCUT and NPEX models. Additionally, we studied the dependence of the cyclotron line energies at different luminosities and discuss the theoretical critical luminosity based on the different models. Regarding the orbital-phase evolution of the spectra, the spectral parameters display variability throughout the orbital phase, and the highly variable photoelectric absorption may suggest the presence of clumpy stellar winds.

Acknowledgments

We are grateful to the referee for the useful suggestions to improve the manuscript. This work is supported by the NSFC (12133007), the National Key Research and Development Program of China (No. 2021YFA0718503), and National Natural Science Foundation of China under grants No. 12173103. The authors thank the support from the Sino-German (CSC-DAAD) Postdoc Scholarship Program, 2022(57607866) and 2023 (57678375). This work has made use of archival data provided by the High Energy Astrophysics Science Archive Research Center (HEASARC).

References

All Tables

All Figures

	Fig. 1. FPMA and FPMB count rates at 3–78 keV of NuSTAR in 2022 with a time resolution of 100 s.
In the text

	Fig. 2. Pulse profiles of Cen X-3 in 3–78 keV from the summations of the FPMA and FPMB counts. These profiles are aligned by defining the minimum as phase zero. The first panel is from the first time interval. The second panel comes from intervals 2 to 13 in the first orbit with a moderate luminosity. The third panel comes from intervals 14 to 23 in the second orbit with high luminosity.
In the text

	Fig. 3. Pulse fraction (PF = (Imax − Imin)/(Imax + Imin)) in energy bands of 3–78 keV of Cen X-3 as a function of luminosity.
In the text

	Fig. 4. Hardness ratios of Cen X-3 as a function of the 4–78 keV luminosity in different energy bands.
In the text

Fig. 5. Pulse profiles of Cen X-3 in different energy bands. The left two panels show the energy-resolved pulse profiles for Cen X-3 from the summations of the FPMA and FPMB counts in MJD 59593.12 (the time interval k of the first orbit) at low intensity state and the right two panels show the energy-resolved pulse profiles in MJD 59594.67 (the time interval g of the second orbit) at high intensity state.

In the text

	Fig. 6. Pulse fraction as a function of energy in 59593.12 MJD and 59594.67 MJD, respectively.
In the text

	Fig. 7. X-ray spectrum of Cen X-3 ranging from 4–78 keV for the interval g of the second orbit using the Fermi-Dirac cutoff continuum. Panel a: residuals of the best-fitting model. Panel b: residuals of the best model, setting the strength of the 30 keV cyclotron line to zero. Panel c: residuals of the spectral fitting without the CRSF.
In the text

	Fig. 8. Spectral continuum parameters (photon index Γ, cutoff energy Ecutoff, fold energy Efold, hydrogen column density NH, and reduced χ2) as a function of X-ray unabsorbed luminosity. The left plots are the fitting results from the FDCUT model (the blue points) and the HIGHECUT model (the black points), and the right plots show the parameters from the NPEX model.
In the text

	Fig. 9. CRSF parameters (energy, width, and depth) as a function of X-ray luminosity and the centroid energy versus the width of CRSF. These parameters are derived from the FDCUT model (the blue points) and the NPEX model (the red points), respectively.
In the text

	Fig. 10. NPEX fits showing the dependence of the CRSFs, the fold energy Efold, the photon index Γ, and the absorption column density NH on the orbital phases for the first binary orbit and second orbit, respectively. The orbital phases are calculated from the orbital parameters provided by Fermi/GBM.
In the text