© The Authors 2024

1. Introduction

The spectral and timing behavior of hot optically thick plasma in the close neighborhood of compact objects is commonly observed by advanced X-ray satellites. This emission is interpreted as thermal radiation modified by strong scattering on free electrons, where photons interact with electrons via the Compton process. The signature of thermal Comptonization is widely observed in neutron star (NS) atmospheres and in accretion disks around NSs or black holes. The crucial parameter applied for these regions is the color-correction factor, which is the ratio of the frequency at which the best-fit spectrum displays the maximum to the frequency of the peak blackbody flux at the effective temperature of the atmosphere. Also called hardening factor, it shows how the spectrum hardens due to scattering effects and how it deviates from a blackbody shape at the source. This factor can be reduced by the effect of gravitational reddening, as seen at the Earth.

The influence of the gravitational field on the apparent radius of the hot NS in an X-ray burster was studied by van Paradijs (1979), who attempted to constrain the possible mass-radius relation for NS. This author used the observed spectra of Type I X-ray bursters during the photospheric radius expansion, assuming that the spectrum is thermalized and parameterized by the effective temperature and gravity. With this method, the NS masses and radii were estimated for particular bursters (see: Fujimoto & Taam 1986; van Paradijs & Lewin 1987), which provided insights into the details of the NS equation of state (EoS) of dense matter.

The EoS encapsulates the properties of the NS matter and allows us to close the system of equations that describes the stellar equilibrium model (see Özel & Freire 2016) by evaluating the mass and radius from the EoS because a specific EoS allows for a certain pair of values for these important NS parameters (see Haensel et al. 2007). Therefore, it is important to either have a large grid of models of continuum spectra from hot and dense NS atmospheres or to have at least a grid of color-correction factors that are correctly estimated, including all physical processes.

The color-correction factor can only be properly used when we can theoretically estimate the deviation of the actual spectrum from the blackbody spectrum, which results from the radiative transfer calculations in scattering atmospheres (Madej 1974). The exact shapes of the Comptonized continuum spectrum and the structure of hot and dense atmospheres of X-ray NS in hydrostatic and radiative equilibrium were studied by Ebisuzaki et al. (1984), London et al. (1984, 1986), Lapidus et al. (1986), Foster et al. (1986). These model atmosphere computations are especially difficult for hot and dense atmospheres, where the opacity is dominated by Compton scattering on free electrons of relativistic thermal motion (Madej 1991a,b). The angular and frequency redistribution of scattered photons was taken from Guilbert (1981). When Compton scattering was taken into account together with true absorption, the numerical computations of radiative transfer through a stratified atmosphere became very time-consuming using low-performance computers. Not many models were therefore computed by the end of the twentieth century.

The method of fitting the observed spectra to theoretical models was proposed by Madej (1997) to obtain the effective temperature T_eff, the surface gravity log g of the NS, the mass, and the radius (see also fits to Chandra spectra of thermal NS by Joss & Madej 2001; Stage et al. 2002a,b, 2004; Joss et al. 2004). The first determination of the mass and radius through fitting X-ray spectra was obtained by (Majczyna & Madej 2005) for the X-ray burster MXB1728-34. The continuum-fitting method requires using the exactest description of Compton scattering and the radiative transfer to determine the X-ray burster parameters.

There exist a number of approximate NS atmosphere models because they oversimplify the electron scattering (Heinke et al. 2006; Webb & Barret 2007; Guillot et al. 2011). These models are available in the software-fitting package XSPEC (see: Arnaud 1996; Zavlin et al. 1996; Heinke et al. 2006; Mori & Ho 2007; Ho et al. 2008, 2014; Ho & Heinke 2009). In many cases, however, observers still use color-correction factors obtained for nonray-traced spectra, which have values higher than unity (Özel & Freire 2016). The color-correction factors estimated by earlier models are mostly greater than 1.3. The influence of increased metallicities, including pure iron models, was also considered, but only for a relatively low effective temperature of NS atmospheres (Rauch et al. 2008; Nättilä et al. 2015).

Advanced models with a full treatment of the Compton scattering were presented by Suleimanov et al. (2012), Medin et al. (2016). Sulejmanov and collaborators used the exact Compton redistribution functions and Medin and collaborators used Kompaneets operators. Next, the models computed by Suleimanov et al. (2014, 2016) were introduced into XSPEC in order to measure the basic NS parameters for atmospheres with effective temperatures up to ∼1.8 × 10⁷ K. The reported color-correction values were in the range 1.4–1.5 and increased above 2.0 for more luminous atmospheres. This high increase in the hardening factor with luminosity was neither obtained in a scarce grid of hydrogen-helium models by Madej (1991a) nor in iron-rich models by Majczyna et al. (2005). The latter models also included a full treatment of the Compton scattering on free electrons in relativistic thermal velocities.

Therefore, we present spectra of hotter NSs up to 3 × 10⁷ K with an exact treatment of the Compton scattering, which is a dominant opacity source in these hot atmospheres (Madej et al. 2004). Our approach adopted the equation of transfer from Pomraning (1973) with the corresponding equations of hydrostatic and radiative equilibrium Madej (1991a, 2004). Our equations allow for a large energy exchange between free electrons and photons with an energy exceeding m_ec². We used the exact photon redistribution function by Nagirner & Poutanen (1993). Recently, we proved that the Compton redistribution functions that we used in the radiative transfer equations were computed with very high precision (see Madej et al. 2017). All our numerical models are stored at the gitlab repository¹.

We show that ray-traced intensity spectra change the value of the color-correction factor in the case of hot NS metal-rich atmospheres. We calculated this for a single value of the surface gravity, which corresponds to EoS of the NS interior that is a relativistic nonrotating polytrope computed with the code LORENE.

The code LORENE² provides spectrum-method solvers to the Einstein equations for an adiabatic index Γ = 2.7 and a pressure coefficient κ = 0.013 (for description see Gourgoulhon et al. 2016). Angle-dependent specific intensity spectra were obtained with radiative transfer calculations with the code ATM24 using the exact Compton scattering redistribution function by Nagirner & Poutanen (1993; for more details, see also Madej et al. 2004; Vincent et al. 2018). For the ray-tracing, we used the code GYOTO³ (see Vincent et al. 2011, 2012).

We note that Beloborodov (2002) and more recently Poutanen (2020) developed approximate formulae describing the light-bending effect around an NS using the Schwarzschild metric. These approximations were not used here because we prepare for the ray-tracing of emission from very rapidly rotating stars where the metric could differ from Schwarzschild.

The structure of the paper is as follows: Sect. 2 presents the most important equations and an overview of the method, including the parameter space of our calculations in Sect. 2.2. The results are described and discussed in Sect. 3. All conclusions are presented in Sect. 4.

2. Model equations

We computed model atmospheres of NSs with the following main assumptions: i) the atmosphere is static in plane-parallel geometry, ii) the gas is ideal, and the ionization populations were computed assuming local thermodynamical equilibrium (LTE), iii) hydrostatic and radiative equilibrium is assumed, iv) the NS does not rotate and therefore has a spherical shape, and v) no magnetic field was assumed.

We applied the equation of transfer with Compton-scattering terms and free-free and bound-free absorption in the following form:

$$\begin{array}{cc}& \mu \frac{\partial J_{\nu }}{\partial {\tau }_{\nu }}=I_{\nu }-{ϵ}_{\nu }{B}_{\nu }-(1-{ϵ}_{\nu })J_{\nu }\hfill \\ \hfill & +(1-{ϵ}_{\nu })J_{\nu }{\displaystyle {\int }_0^{\infty }}{\mathrm{\Phi }}_1(\nu ,{\nu }^{\prime })\mathrm{d}{\nu }^{\prime }-(1-{ϵ}_{\nu }){\displaystyle {\int }_0^{\infty }}{J}_{{\nu }^{\prime }}{\mathrm{\Phi }}_2(\nu ,{\nu }^{\prime })\mathrm{d}{\nu }^{\prime },\hfill \end{array}$$(1)

where B_ν is the Planck function, J_ν is the mean intensity, μ denotes the cosine of the zenith angle, τ_ν is the monochromatic optical depth, and ϵ_ν is the ratio of the true absorption to the total opacity. A more detailed description of the above equation of transfer and of our computing of the angle-averaged Compton redistribution functions Φ₁(ν, ν′) and Φ₂(ν, ν′) and the linearization procedure can be found in our earlier papers, (Madej 1989, 1991a; Vincent et al. 2018).

We emphasize here that the algorithms with which we computed the redistribution functions Φ₁ and Φ₂ in Eq. (1) automatically determined partial angle-averaged scattering coefficients σ(ν, ν′,T) for each pair of frequencies and temperatures in the model atmosphere. Therefore, we carefully computed the integrated scattering coefficients σ(ν, T), which always decrease below the classical Thomson-scattering coefficient for a rising photon frequency ν, similarly as in the Klein-Nishina formula.

For NS atmospheres with effective temperatures T_eff ≥ 10⁷ K, the atmosphere can be very thin (1–10 m; see Madej 1991b) for a surface gravity significantly higher than the critical gravity. In this case, the ratio of the height of the atmosphere to the NS radius is much lower than unity. This condition justifies our assumption that the atmosphere is in plane-parallel geometry.

Gravitational acceleration seen by an observer on the NS surface is

$$\begin{array}{c}\hfill g=\frac{\mathit{GM}}{R_{\star }^2}{(1-\frac{R_g}{R_{\star }})}^{-1/2}=\frac{\mathit{GM}}{R_{\star }^2}(1+z),\end{array}$$(2)

where M denotes the mass of the NS, and, R_g = 2GM/c² is the gravitational radius, G is the gravitational constant, c denotes the velocity of light, R_⋆ denotes the radius of a star measured on its surface, and 1 + z = (1 − R_g/R_⋆)^−1/2 is the gravitational redshift factor.

For each model of the local atmosphere, we assumed values of the effective temperature T_eff and the surface gravity g, where the former was used in the flux-luminosity relation, and latter to determine the NS radius, both in order to obtain the total observed luminosity,

$$\begin{array}{c}\hfill L=\frac{4\pi R_{\star }^2\sigma T_{\mathrm{eff}}^4}{{(1+z)}^2}=\frac{4\pi GM\sigma T_{\mathrm{eff}}^4}{g(1+z)}.\end{array}$$(3)

All these calculations were made by us to determine the relative luminosity l = L/L_Edd. The observed Eddington luminosity L_Edd had the usual form,

$$\begin{array}{c}\hfill L_{\mathrm{Edd}}=\frac{4\pi GMc}{{\kappa }_{\mathrm{es}}(1+z)},\end{array}$$(4)

where κ_es is the Thomson-scattering coefficient. Thus, the relative luminosity can be expressed as the ratio of the critical gravity g_crit to the actual gravity at the NS surface,

$$\begin{array}{c}\hfill l=\frac{\sigma T_{\mathrm{eff}}^4{\kappa }_{\mathrm{es}}}{\mathit{gc}}=\frac{g_{\mathrm{crit}}}{g},\end{array}$$(5)

where the critical gravity was obtained from the condition that the radiation pressure balances the downward gravitational acceleration. In the sections below, we show the dependence of our results on the relative luminosity l, fully determined numerically with a proper value of the Compton-scattering coefficient mostly lower than the classical Thomson value.

2.1. Chemical compositions

In previous years, Nättilä et al. (2017) presented the first direct fit to a bursting NS model atmosphere with the use of the models of Suleimanov et al. (2012). The authors considered several chemical abundances of atmospheres, and we followed the same approach in our computations in order to make our models more comparable. Following these authors, we consider four chemical compositions: pure hydrogen, pure helium, and a mixture of hydrogen and helium, together with different heavy element abundances scaled to the solar chemical composition (Z_⊙ and 0.01 Z_⊙). We used values of the solar elements abundances from Asplund et al. (2009). The values for the composition 0.01Z_⊙ were obtained by lowering the mass of heavy elements (by a factor of 100) and equally redistributing this mass between hydrogen and helium. The abundances of the four chemical compositions considered in this paper are presented in Table 1.

2.2. Grid of models

We calculated a grid of models for different chemical compositions (pure hydrogen, pure helium, solar composition, and 1% of solar composition). For each chemical composition, we created a grid of temperatures from T_eff = 10⁷ K to T_eff = 3.0 × 10⁷ K with ΔT = 0.02 × 10⁷ K and the surface gravity grid, which starts from the critical gravity g_crit and reaches log(g) = 15.00 (cgs), with a grid step equal to Δlog g = 0.01. The grids of the computed models (with Δg = 0.05 only for better visibility) are marked as dots on the temperature versus surface gravity plane in Fig. 1 for all considered chemical compositions.

Fig. 1. Temperature – gravity grid for atmospheres composed purely of hydrogen (upper panel), purely of helium (middle panel), and of heavy elements (bottom panel). The gray dots show the gravity grid with Δg = 0.05 in logarithmic scale, because for Δg = 0.01, the figure would not be legible. The different dashed curves represent models with constant l at the values described in the legend. The horizontal solid red line displays the value of the gravity for the best-fit mass and radius measured by Nättilä et al. (2017), and the two horizontal dashed red lines correspond to the measured errors (see Sect. 2.3).

Each individual point in Fig. 1 represents the model atmosphere computed using a grid of 1025 wavelength points in the range 0.01 − 25000 Å and a grid of 96 Thomson depths from 1e−7 to 7.5e+04. Both grids are equidistant in logarithmic scale.

The whole data set computed by us contains more than 45 000 models, which is the densest grid computed so far. For comparison, Suleimanov et al. (2011) presented a grid with 360 models. This grid was computed for only three values of the surface gravity (log(g) = 14.0, 14.3, and 14.6), for six chemical compositions, and for 20 values of the relative luminosities (starting from l = 0.001 up to l = 0.98). The grids computed by us follow the same range of relative luminosities.

The critical gravity given by Eq. (5) depends on the assumed effective temperature and on κ_es, the latter being a function of mass abundance of hydrogen, helium, and heavy elements given by X, Y, and Z, respectively (see Table 1 for the exact numbers). Thus, the value of the critical gravity is connected to the chemical abundance of the atmosphere. For atmospheres composed only of pure hydrogen (κ_es = 0.4), the value of the critical gravity is higher than the atmosphere composed only of helium (κ_es = 0.2) or an atmosphere with a solar composition (κ_es = 0.35) for the same effective temperature T_eff. The relative luminosity for which the atmosphere reaches the critical gravity equals unity, and it is marked by green solid line in each panel of Fig. 1.

2.3. Neutron star models and ray-tracing

Models computed with ATM24 are local, meaning that they describe emission from the surface of the atmosphere. When radiation is to be transported correctly from the atmosphere to a distant observer and is to include all processes that occur in close proximity to the NS, the ray-tracing procedure must be used. Following our approach (Vincent et al. 2018), we used the open-source ray-tracing code GYOTO to transport spectra to the distant observer. We obtained the spacetime solution using the LORENE/NROTSTAR library, in particular, we estimated the value of the surface gravity for a chosen EoS model, corresponding to the particular values of the NS mass and radius. In the next step of our procedure, this surface gravity, together with effective temperature, was an input parameter in the code ATM24, and after the local spectrum was computed, it was ray-traced to the observer.

The selected benchmark EoS is a relativistic polytrope model, composed of identical relativistic particles of rest-mass m₀, as defined⁴ in the LORENE library:

$$\begin{array}{c}\hfill p(n)=\kappa n^{\gamma },e(n)=p(n)/(\gamma -1)+{\mu }_{0}n,\end{array}$$(6)

where n is the numerical (baryon) density, p(n) is the pressure, e(n) is the mass-energy density, and μ₀ = m₀c² is the chemical potential at zero pressure. The values of the parameters are the particle mass m₀ = 1.66 × 10⁻³⁰ g, the polytropic index γ = 2.7, and the pressure coefficient κ = 0.013 in units of $({\rho }_{\text{nuc}}{c}^2)/n_{\text{nuc}}^{\gamma }$, with ρ_nuc = 1.66 × 10²⁰ g m⁻³, n_nuc = 0.1 fm⁻³.

With this choice, we produced a configuration consistent with the best-fit mass and radius obtained in Nättilä et al. (2017), as shown in Fig. 2 by the solid red line, intersecting the measured value marked by the central red dot, and the gray error bars in M and R correspond to the measurement errors. For comparison, we also show M(R) sequences obtained using realistic dense-matter EoSs solutions: the SLY4 EoS (Chabanat et al. 1998; Douchin & Haensel 2001), the BSK20 EoS developed by the Brussels-Montreal group (Potekhin et al. 2013), and the GM1Z0 EoS based on a relativistic mean field model (Weissenborn et al. 2012). These dense-matter EoSs were supplemented with the low-density part description of the NS crust of Haensel & Pichon (1994; for more detailed comparisons, see Bejger 2013).

Fig. 2. NS mass-radius relations for selected EoSs. The dash-dotted blue curve presentes the SLY4 EOS, the dashed yellow curve shows the BSK20 EoS, and the dotted green curve shows GM1Z0 EoS. The red curve denotes the relativistic polytropic EoS consistent with the measurements of Nättilä et al. (2017), with dots marking the center of the error box and its boundaries.

Subsequently, we ray-traced spectra for three values of the surface gravity that corresponded to the mass-radius values presented by red dots in Fig. 2, where central dot is a measurement, and the two edge dots are errors. The values of these surface gravities are marked by horizontal solid red lines (measured value), and dashed red lines (errors) in each panel of Fig. 1, meaning that ray-traced spectra span the whole effective temperature range for each chemical composition. In the case of the nonrotating spherical NS considered in this paper, the ray-traced spectra differ from the local spectra by only a factor 1 + z, where z is the gravitational redshift. Therefore, we prepared our grid of models composed of local spectra in the format required by the XSPEC fitting package. Nevertheless, a pipeline for computations of rotating NS atmospheres is ready and will be prepared in the near future.

3. Results

3.1. Temperature structure of the model atmospheres

The first direct result of our computations is an angle-dependent intensity spectrum from the stratified NS atmosphere, ray-traced in the gravitational field. The obtained spectral energy distribution may be directly fit to the observations, indicating the temperature and gravity on the NS surface. Nevertheless, in stellar atmospheres dominated by electron scattering, photons are created in the so-called thermalization layer (by thermal processes below the photosphere). In the thermalization layer, the local temperature T is higher than the effective temperature T_eff, from which we conclude that the maximum of the Planck function of B_ν(T) appears at a higher energy than the maximum of the Planck function of B_ν(T_eff). This effect has far-reaching consequences because the blackbody temperature obtained from fitting the continuum spectra is usually higher than the effective temperature T_eff, and it is indicated as the color temperature T_c.

3.2. Color-correction factor

The color-correction factor is defined as the ratio of the frequency at which the best-fit spectrum displays the maximum to the peak of the blackbody spectrum at the effective temperature of an atmosphere. It can be directly described by the ratio of the color temperature to the effective temperature of the star as f_c = T_c/T_eff. The color-correction factor can increase to very high values in atmospheres that only take Thomson scattering into account because the true absorption in these atmospheres decreases to zero (see Madej 1974). In atmospheres that include Compton scattering, this factor is larger than unity, but lower than in the case of elastic electron scattering. For a fast estimation of the value of the color temperature either from a model or from an observed spectrum, we used the following formula:

$$\begin{array}{c}\hfill T_{\mathrm{c}}=\frac{E_{\max }}{2.821k_{\mathrm{B}}},\end{array}$$(7)

where k_B is the Boltzmann constant, and E_max is the energy for which the spectrum of the NS at a given temperature T_eff reaches maximum.

Another way to obtain the color temperature is the procedure of fitting the blackbody spectrum to our simulated or observed data by the minimalization of the sum

$$\begin{array}{c}\hfill {\displaystyle \sum _{i=1}^N}(F_{{\nu }_i}-\upsilon B_{{\nu }_i}(T_j)),j=1,\cdots ,M,\end{array}$$(8)

where F_ν is the simulated or the observed flux, and υ is the normalization factor. Our procedure started with finding the maximum emission and corresponding frequency, and then we took a range of frequencies around the maximum to make our fit appropriate. After we created a table of frequencies, we set the grid of temperatures starting from T_eff of our model up to three times that temperature and divided it into M numbers of small intervals. We do not expect the color-correction factor to be higher than 3. Then we computed the renormalization factor in the following way:

$$\begin{array}{c}\hfill \upsilon =\frac{{\displaystyle {\int }_{{\nu }_{min}}^{{\nu }_{max}}}{F}_{\nu }\mathrm{d}\nu }{{\displaystyle {\int }_{{\nu }_{min}}^{{\nu }_{max}}}{B}_{\nu }(T_j)\mathrm{d}\nu },\end{array}$$(9)

where v_max and v_min correspond to the frequencies around the maximum emission of the spectral curve. We used a χ² test to determine the temperature of the blackbody spectrum that fits our simulated data best, which we attributed to the value of T_c. The frequencies around maximum were always adapted to each model atmosphere, and they might depend on different chemical composition. For instance, for a pure hydrogen atmosphere at a temperature of 1.5 × 10⁷ K, the energy range was from 8.1 × 10¹⁷ Hz up to 2.1 × 10¹⁸ Hz.

This fitting method allowed us to derive the color temperature in local spectra. In the case of a nonrotating NS, the hardening factor from ray-traced spectra was divided by 1 + z with respect to the same hardening factor for local models. In the first step, f_c, 1 is the above factor computed from the ratio of the color temperature obtained simply from Eq. (7) and the effective temperature of the model T_eff. In the second step, f_c, 2 is the ratio of the color temperature, obtained from the blackbody fitting procedure to the local spectrum to the T_eff of the model.

In Table 2 we show an example of the values of the color-correction factor for an effective temperature T_eff = 1 × 10⁷ K and surface gravity log(g) = 14.34, the latter being the result of the best-fit mass and radius obtained by Nättilä et al. (2017).

In our models, the temperature structure is determined by two sources of monochromatic opacities. In deeper layers of the atmosphere, newly generated photons undergo free-free and bound-free absorption by ions, which heats the surrounding atmosphere. As the height increases (or the opacity decreases), the absorption cross section falls. Free-free absorption dominates until the effective Thompson depth is reached (τ_eff). Now, there are very few ions to provide enough capture cross section for absorption, so that scattering starts to dominate. At T_eff ≥ 10⁷ K, Compton scattering becomes important (Madej 1991a; Madej et al. 2004; Suleimanov et al. 2016). Cold electrons in the upper atmosphere efficiently exchange energy with high-energy photons, leading to an increase in temperature. An isothermal layer develops in the photosphere with T ≥ T_eff, as presented in Figs. 3 and 4.

Fig. 3. Temperature structure for models with a fixed effective temperature Teff = 1, 2, and 3 × 107 K (upper, middle, and bottom panel, respectively). Different chemical compositions are presented by different line styles given in the legend. Two relative luminosities are presented by red lines for l = 0.95 in each figure panel and by black lines for l = 0.05, l = 0.13, l = 0.65 in the upper, middle, and bottom panel, respectively.

Fig. 4. Temperature structure for models with a fixed surface gravity equal to the measured value log(g) = 14.34, which corresponds to an NS mass 1.9 M⊙ and a radius of 12.4 km (middle panel; Nättilä et al. 2017), and the errors are log(g) = 13.70 (upper panel), and log(g) = 14.90 (bottom panel). The different chemical compositions are presented by different line styles given in the legend. Two relative luminosities are presented by red lines for l = 0.95 in the upper and middle panels, for l = 0.3 in the bottom panel, and by the black lines for l = 0.05, l = 0.13, l = 0.65 in the upper, middle, and bottom panel, respectively.

Both effects combined cause a minimum in the temperature profile below the photosphere and a significant temperature rise in the uppermost layers of the atmosphere. The depth and location of the minimum depends on the value of the surface gravity and on the chemical composition (Fig. 3). For a luminous atmosphere, the temperature minimum is less pronounced and is also shifted toward higher optical depths. The local spectrum, presented by dashed lines in all panels of Fig. 5, deviates from the blackbody curve B_v(T_eff) presented by dotted lines in the figure. All these features agree with previous papers on NS atmospheres (Madej 1991a; Madej et al. 2004; Suleimanov et al. 2016).

Fig. 5. Emergent spectra of the model atmospheres before (dashed lines) and after ray-tracing (solid lines), and emergent spectra of local blackbody emission (dotted line). The models were computed for a fixed effective temperature Teff = 1, 2, and 2.3 × 107 K (upper, middle, and bottom panel, respectively) for all chemical compositions given by colors, and for a fixed value of surface gravity log(g) = 14.34, which corresponds to an NS mass of 1.9 M⊙ and a radius of 12.4 km.

In Fig.4, the effect of the variation in T_eff on the atmosphere structure is shown at constant log(g). The temperature structure depends on the chemical composition of the NS atmosphere. This is more evident at lower T_eff. When we compare pure hydrogen NS atmospheres, their structures do not differ much from those of 0.01 Z_⊙ with the same other input parameters. For He-rich atmospheres, the characteristic temperature jump occurs at lower optical depth than in case of other chemical abundances. Since models of solar abundance are composed of a mixture of many elements, their net effect is an erratic variation. At sufficiently high T_eff, all elements are ionized and show similar variation.

3.3. Ray-traced spectra

The continuum spectra for the value of the surface gravity log(g) = 14.34 and for all chemical compositions are presented in Fig. 5 for three effective temperatures starting from T_eff = 1 × 10⁷ K. The highest possible temperature for which the above gravity does not exceed the critical gravity for all considered abundances is T_eff = 2.3 × 10⁷ K. In all figures, the spectra generated from the local ATM24 are represented by dashed lines and the ray-traced models by solid lines. The deviation from blackbody spectra (dotted green curve) is immediately apparent. This was extensively discussed previously and is summarized with the use of the color-correction factor in the section below.

In this section, we verify the correctness of our models according to the ray-tracing procedure used. Compared to the ATM24 model, all the ray-traced (RT) model spectra terminate at lower energies. This is a direct consequence of the gravitational redshift (z). The local intensity is reduced by a factor of 1+z (Vincent et al. 2018) as it reaches the observer. This has led to lower color-correction factor values.

We used a crude estimation of the redshift value from the spectrum and tried to back-calculate the mass-radius ratio as considered by the polytropic EoS in Sect. 2.3. We considered the shift in the energy that corresponds to the peak flux for both spectra computed by us, as shown in Fig. 6. Using the standard relation between an energy in which the local spectrum peaks at E′, and an energy for which the observed spectrum peaks at E, that is, E′/E = 1 + z, the z value is 0.37. Finally, for this redshift value, by using Eq. (2), we obtained an NS mass-radius ratio of 0.157. This roughly matches the mass of 1.9 M_⊙ and radius of 12.4 km of NS used in this work and presented in Fig. 2. The above exercise proves that our numerical models fully take gravitational correction factors due to the compact nature of the NS into account.

Fig. 6. Quantitative shift in peak energy due to redshift (z) from the local ATM24 model spectrum to ray-traced (RT) spectrum for Teff = 2 × 107 K and log(g) = 14.34. The red box shows the redshift calculation from the standard relation between the energy of the local intensity peak E′, marked by the vertical dash-dotted black line, and energy of the observed intensity peak E, marked by the vertical red line of the same style.

The overall spectral shape for the same effective temperature and surface gravity depends on the chemical abundance. It is interesting to note that H and He form two extremes in any type of variation in the temperature structure (Fig. 3) and continuum spectrum (Fig. 5). For the lowest effective temperature, only the He-like atmosphere spectrum stands out from the rest, especially for the high-energy end. The difference becomes larger with the rise in temperature, and it is clearly visible in the bottom panel of Fig. 5. The ray-tracing procedure does not change this trend.

In the case of a low effective temperature, Fig. 5, the Z_⊙ curve for both model spectra and RT spectra is broken at higher photon energies. These are iron absorption edges at around 8–9 kev (Suleimanov et al. 2011). Similar features were reported by Nättilä et al. (2015), which are a result of the incomplete ionization of Fe ions and other heavy elements. The spectral curve is sensitive to metal abundances at higher effective temperatures, as shown in the bottom panel of Fig. 5.

3.4. Color-correction factor versus luminosity ratio

We used the convention to plot the color-correction factor versus luminosity ratio in Figs. 7 and 8, which we did not use before. Since the temperature and surface gravity both count for the luminosity ratio value, the first plot was made for one single effective temperature, and the second plot was given for a fixed value of the surface gravity, log g = 14.43. We calculated the factors for all local spectra, f_c, 2, with the fitting procedure described in Sect. 3.2.

Fig. 7. Color-correction factor fc, 2, obtained with the fitting procedure described in Sect. 3.2, vs. luminosity for the models with a fixed value of the effective temperature Teff = 2 × 107 K for all chemical compositions.

Fig. 8. Color-correction factors for local models fc, 2, obtained with the fitting procedure described in Sect. 3.2. The plots show the luminosity dependence for a fixed value of the surface gravity log(g) = 14.22 (top panel), log(g) = 14.34 (middle panel), and log(g) = 14.43 (bottom panel), i.e., the values inferred from the NS model (see Fig. 2), and for all chemical compositions.

The color-correction factors from local spectra f_c, 2 are presented for one case of a moderate value of the temperature T_eff = 2 × 10⁷ K in Fig. 7, spanning all possible surface gravities and for all considered abundances. Fig. 8 shows color-correction factors as red lines for a reference measured value of the surface gravity log(g) = 14.34 (middle panel) and for the errors (top and bottom panels). In each case, all chemical compositions are presented. The obtained values of f_c, 2 are consistent with previous estimates made by our group (Majczyna et al. 2005). The dependence on metallicity is clearly visible here. Pure hydrogen atmospheres present the largest color shift for a wide range of luminosities, except for the case when l ≥ 0.9, where metal-rich atmospheres display a higher color correction factor. Pure helium atmospheres always have color factors lower by 0.1 than pure hydrogen models, which is consistent with Suleimanov et al. (2011). This difference is due to the variation in the number density of free electrons, which are needed for the Compton process to be important. The main difference of our models to previously published work (Suleimanov et al. 2011) is the low-luminosity end of the dependence of the color-correction factor versus luminosity ratio. In the case of our models and for all considered abundances, the value of f_c, 2 saturates at about 1.4 when the luminosity decreases. This discrepancy probably results from the differences between the Fortran procedures that compute Gaunt factors for the free-free absorption that are used by our groups.

Since our color-correction factor estimates are well derived for the wide range of l values, they may be applicable to a wide variety of systems such as the quiescent state of NS and X-ray bursters (l ≥ 0.8).

4. Conclusions

We presented a new grid of nonrotating NS models for hot stratified atmospheres resulting from radiative transfer computations. We focused on estimating the hardening factor for spectra of NS with a wide range of abundances and luminosities. The computational code ATM24 takes into account free-free and bound-free energy-dependent opacities of hydrogen, helium, and ions of heavy elements, including iron. By fitting a blackbody shape to local and ray-traced spectra near the highest emission region, the color-correction factors were calculated. They can be used to determine the spectral shape and NS parameters from observations.

The RT procedure accounts for the gravitational redshift suffered by photons emitted from an atmosphere. Photons leaving the photosphere of an NS are subjected to an immense gravitational potential that redshifts their energy for the distant observer. Additional effects can also take place, such as Doppler boosting and anisotropic emission. These effects need to be accounted for when analyzing the observed spectrum and estimating other quantities, such as the NS mass and radius (for example see: Vincent et al. 2018). The most prominent effect that contributes to a substantial difference in the color-correction factor is the gravitational redshift. Since the redshift factor z depends on the mass of the compact object, our calculations of the RT models and their use for a data analysis can place more reliable constraints on the NS mass and radius determination with a continuum-fitting method (Majczyna et al. 2005).

The composition of an NS atmosphere is a long-standing research question and is connected with the formation and evolution of an NS. We considered important elements that might comprise the NS atmosphere and showed how it affects the spectrum and temperature structure. The effect of gravity is to vary the density of lighter and heavier elements in the NS atmosphere. Depending on the physical processes involved after the formation of the NS, the atmosphere maybe dominated by H or He or by heavier elements such as Fe. Generally, accreting NSs leave behind a hydrogen-rich environment, while nonaccreting NSs are rich in iron (Ho et al. 2014). Nuclear burning processes leave behind helium or carbon (Ho et al. 2014). As shown by Zavlin et al. (1996), an iron atmosphere has a higher opacity than H or He atmospheres at a given photon energy. They also showed that at a given depth, an iron atmosphere has a higher temperature than an atmosphere that is composed of H or He.

We demonstrated that for a lower l value, the temperature minimum in the structure of an atmosphere for solar abundance (Z_⊙ = 1) is almost never as low as in case of a low metallicity or pure hydrogen or helium atmospheres. The rise in temperature due to Compton scattering is also not sharp, unlike the other abundance curves. Due to the mixture of different opacities of lighter and heavier elements in Z_⊙, the net effect is a washed-out temperature minimum.

We are not able to predict the redshift value for observed spectra of a NS since it depends on the assumed EoS. Therefore, only the direct fitting of the ray-traced spectra may give us realistic values of the mass and radius.

In a forthcoming project, we additionally plan to investigate the effect of rotation on the value of the color-correction factor. After this step is made, we can fully use our models to work with observations of hot NS emission, determine the main parameters such as the NS mass and radius, and finally, place constraints on the EoS of the dense matter under strong gravitational conditions.

Acknowledgments

This research was partially supported by the Polish National Science Center grants No. 2021/41/B/ST9/04110 and 2021/43/B/ST9/01714.

References

All Tables

All Figures

Fig. 1. Temperature – gravity grid for atmospheres composed purely of hydrogen (upper panel), purely of helium (middle panel), and of heavy elements (bottom panel). The gray dots show the gravity grid with Δg = 0.05 in logarithmic scale, because for Δg = 0.01, the figure would not be legible. The different dashed curves represent models with constant l at the values described in the legend. The horizontal solid red line displays the value of the gravity for the best-fit mass and radius measured by Nättilä et al. (2017), and the two horizontal dashed red lines correspond to the measured errors (see Sect. 2.3).

In the text

	Fig. 2. NS mass-radius relations for selected EoSs. The dash-dotted blue curve presentes the SLY4 EOS, the dashed yellow curve shows the BSK20 EoS, and the dotted green curve shows GM1Z0 EoS. The red curve denotes the relativistic polytropic EoS consistent with the measurements of Nättilä et al. (2017), with dots marking the center of the error box and its boundaries.
In the text

Fig. 3. Temperature structure for models with a fixed effective temperature Teff = 1, 2, and 3 × 107 K (upper, middle, and bottom panel, respectively). Different chemical compositions are presented by different line styles given in the legend. Two relative luminosities are presented by red lines for l = 0.95 in each figure panel and by black lines for l = 0.05, l = 0.13, l = 0.65 in the upper, middle, and bottom panel, respectively.

In the text

Fig. 4. Temperature structure for models with a fixed surface gravity equal to the measured value log(g) = 14.34, which corresponds to an NS mass 1.9 M⊙ and a radius of 12.4 km (middle panel; Nättilä et al. 2017), and the errors are log(g) = 13.70 (upper panel), and log(g) = 14.90 (bottom panel). The different chemical compositions are presented by different line styles given in the legend. Two relative luminosities are presented by red lines for l = 0.95 in the upper and middle panels, for l = 0.3 in the bottom panel, and by the black lines for l = 0.05, l = 0.13, l = 0.65 in the upper, middle, and bottom panel, respectively.

In the text

Fig. 5. Emergent spectra of the model atmospheres before (dashed lines) and after ray-tracing (solid lines), and emergent spectra of local blackbody emission (dotted line). The models were computed for a fixed effective temperature Teff = 1, 2, and 2.3 × 107 K (upper, middle, and bottom panel, respectively) for all chemical compositions given by colors, and for a fixed value of surface gravity log(g) = 14.34, which corresponds to an NS mass of 1.9 M⊙ and a radius of 12.4 km.

In the text

Fig. 6. Quantitative shift in peak energy due to redshift (z) from the local ATM24 model spectrum to ray-traced (RT) spectrum for Teff = 2 × 107 K and log(g) = 14.34. The red box shows the redshift calculation from the standard relation between the energy of the local intensity peak E′, marked by the vertical dash-dotted black line, and energy of the observed intensity peak E, marked by the vertical red line of the same style.

In the text

	Fig. 7. Color-correction factor fc, 2, obtained with the fitting procedure described in Sect. 3.2, vs. luminosity for the models with a fixed value of the effective temperature Teff = 2 × 107 K for all chemical compositions.
In the text

	Fig. 8. Color-correction factors for local models fc, 2, obtained with the fitting procedure described in Sect. 3.2. The plots show the luminosity dependence for a fixed value of the surface gravity log(g) = 14.22 (top panel), log(g) = 14.34 (middle panel), and log(g) = 14.43 (bottom panel), i.e., the values inferred from the NS model (see Fig. 2), and for all chemical compositions.
In the text