1 Introduction↩︎

Superluminous supernovae (SLSNe) are rare types of stellar explosions with energies \(\sim100\) times that of ordinary supernovae (SNe). The mechanism producing the excess energy compared to regular core-collapse SNe remains debated. Commonly proposed explanations include strong interaction between the circumstellar medium (CSM) and ejecta [1], [2], highly magnetized and rapidly rotating neutron stars (magnetars; [3]–[5]), accreting black holes [6], or a combination of different mechanisms. The former mechanism is typically invoked for hydrogen rich SLSNe (SLSNe-II), while the latter two are often proposed for hydrogen-poor SLSNe (SLSNe-I).

Using a magnetar as the driving force of SLSN-I, the observed spectra and light curves have been successfully modeled [7]–[11]. However, some SLSN-I show emerging hydrogen emission lines at later times [12]–[15], suggesting interaction between ejecta and previously expelled hydrogen layers. This transition has also been observed in ordinary stripped SNe [16]–[21]. The strong shift from stripped to interacting SNe suggests a shell-like CSM structure that could come from eruptive mass loss of the progenitor star, rather than from a steady wind. These bursts of mass loss is common in luminous blue variable (LBV) stars, which are known to be the progenitors of initially strongly interacting SNe (IIn) [22], and may also explain the superluminous population [23], [24]. Furthermore, [25] recently showed that such shell-like CSM structures are seen for a sample of SLSN-I and is consistent with mass-loss episodes from LBV-like eruptions or pulsational pair-instability events.

Despite the successful modeling of optical data, distinguishing between the power sources remains challenging as different mechanisms can produce similar optical signatures (e.g. [26], [27]). Observations at other wavelengths are therefore essential to break this degeneracy and probe the physical conditions of interaction or a central engine more directly. Some studies have addressed this by using radio observations, and detected SLSNe with support of both types of power origin [28], [29].

X-ray observations provide another important probe of both CSM interaction and central engine activity. Detecting X-rays associated with magnetar activity requires careful timing. At early epochs the dense ejecta efficiently absorb soft X-rays, but as the ejecta expand over time, the absorption and the optical depth (\(\tau\)) decrease [30], [31]. Observations must therefore be sufficiently late for the ejecta to have become partially optically thin, yet early enough for the decaying magnetar luminosity to remain detectable. In addition, magnetars, with their associated pulsar wind nebulae (PWNe), may be able to ionize the ejecta, thereby making the medium transparent to X-rays [32]. This could make the magnetar luminosity detectable at earlier epochs. Late-time detections or limits of X-ray emission, some months to years after explosion, would offer valuable constraints of magnetar properties.

While some work on X-ray emission from SLSNe have been presented [33]–[36], the detections often cover early epochs reaching only \(\sim160\) days post explosion. The majority of the X-ray observations of SLSNe result in nondetections. In most cases, the observations were obtained close to the optical peak and/or are too shallow to probe the X-ray luminosities expected from CSM-shock interaction or emerging magnetar emission (e.g., [35]). While upper X-ray luminosity limits may constrain the mass loss rates [34], late detections are required to identify the powering mechanism. To date, only one X-ray detection of an SLSN-I has been claimed, the SLSN-I SCP 06F6 at \(\sim150~\)days after explosion with an X-ray luminosity of \(\sim1\times10^{45}~\rm erg~s^{-1}\) in the 0.2–10 keV energy band [34]. [36] also report an X-ray detection which could be associated with the SLSN PTF 12dam, but state that the star formation within the host galaxy partly contributes to the X-ray luminosity and treat the result as an upper limit. Similarly, only one SLSN-II has been detected in X-rays; SN 2006gy measured to a luminosity of (0.5–2 keV) \(1.65\times10^{39}\rm~erg~s^{-1}\) at \(\sim 60~\)days post explosion [33]. This SLSN persistently show narrow lines in its spectrum, and was originally classified as a type IIn [33], [37]. The narrow lines and exceptional luminosity suggest that SN 2006gy could be considered a Type SLSN-IIn (e.g. as mentioned in [38], [39]).

In this study, we present an in depth X-ray analysis of SN 2018bsz, the to date closest observed SLSN-I at 111 Mpc. SN 2018bsz was discovered on May 17th 2018 [40], [41] by the All Sky Automated Survey for Supernovae¹ (ASAS-SN; [42]). The explosion date is estimated to be 2018 March 25 by [43], who also reclassified it from a type II SN to an SLSN-I [44]. Although SN 2018bsz shows common features with other SLSN-I, such as strong C II and O II absorption features in the early spectra [43] and a late-time spectral evolution resembling that of ordinary stripped type Ic, it also shows a number of notable deviations. Firstly, the optical light curve exhibits a long and slow rise, with a sudden steepening just before peak brightness, followed by hydrogen emission lines emerging \(\sim25~\)days after peak [45]. SN 2018bsz is also the first SLSN to show dust formation within the ejecta [46]. Although hydrogen lines have previously been detected for a small sample of SLSN-I [12]–[15], SN 2018bsz is unique in exhibiting a multi-component H\(\alpha\) profile, a characteristic more commonly associated with SN  IIn [45]. Based on this line profile and its spectral evolution, [45] proposed that this event is interacting with an asymmetric, disk-like CSM.

In this paper, we study XMM and Chandra observations of SN 2018bsz from \(\sim3\) months to \(\sim3\) yr post explosion. The proximity and multi-epoch X-ray observations of SN 2018bsz makes it a well-suited target for investigating the high-energy emission from SLSN-I and for placing constraints on proposed powering mechanisms, including magnetar central engines and CSM interaction. The paper is organized as follows. Section 2 presents the observation and data reduction. The analysis and main results are described in Section 3. In Section 4, we discuss the results and investigate the origin of the X-ray emission. We finally present the conclusions and summary in Section 5.

2 Observations and Data Reduction↩︎

This study makes use of all available X-ray observations of SN 2018bsz. This includes four Chandra observations near in time of explosion (87–304 days post explosion, PI Margutti) and one later XMM-Newton observation (1253 days post explosion, PI: Alp) (see Table [tab:observaions]).

2.1 Chandra↩︎

SN 2018bsz was observed in four Chandra observations, each carried out with the ACIS-S instrument, specifically using chip S3. The data were downloaded from the Chandra archive and processed with CIAO ver. 4.13 [47]. Event files and images were produced through the CIAO task chandra_repro and dmcopy.

Images (bin=1) in the 0.5–8 keV energy range were used for running a source detection algorithm. Point-source identification employed the wavlet-based algorithm wavdetect, which correlates the image with wavlets of different scales (radii 1.0, 2.0 and 3.0 pixels). A significance threshold of sigthresh=\(10^{-5}\) was adopted, and standard parameters were retained otherwise. A detection of \(>3\sigma\) at the positions consistent with the optical coordinates of the SNe is seen in all Chandra images (see Table [tab:observaions] for more details) . We also detect additional sources within 10\(\arcsec\) from the SN position (see Figure 1). These additional sources could potentially contaminate the SN spectrum of the XMM observation, and are therefore analyzed as described later in Section 3.

Source and background regions were defined manually after visual examination of the 0.5–8 keV images. The circular source regions typically have radii of \(\sim3\arcsec\) (see Table [tab:observaions] for details), while background regions of \(\sim~23\arcsec\) were chosen close to the source but clear of unrelated emission. Manual definition of the source regions and sizes was preferred, as the automatic CIAO task psfsize_srcs often returned unrealistically small regions that excluded a great portion of the SN emission.

Light curves were created using dmextract to inspect and filter out intervals affected by background flares using deflare. As the light curves showed no significant flaring events, this had little to no impact on the observation and effective exposure time.

Spectral products of each observation were extracted using specextract. The auxiliary response files (arf) were generated without weighting, which is appropriate for point-source analysis, but with a point-source aperture correction applied to the unweighted arf. Owing to the limited number of counts, Poisson-based C-statistics was used for spectral analysis in XSPEC. To accommodate this, Chandra spectra were grouped with one count per spectral bin using grpmincounts=1 within spacextract.

2.2 XMM↩︎

The XMM EPIC datasets were retrieved from the XMM Science Archive and processed using the calibration files (CCF) from 2025. Data reduction followed the standard guidelines provided for the Science Analysis Software (SAS, v22.1.0; [48]). Images with a bin size of one were produced over the 0.5—10 keV energy range and used to perform source detection.

Source detection was carried out with the SAS task ewavelet, which identifies significant wavelet correlations across multiple spatial scales. This approach closely parallels the wavdetect algorithm employed for Chandra. We find a significant detection centered \(3.6\arcsec\) from the optical SN position. This offset is discussed in more detail in Section 3.2, but likely comes from the additional aforementioned X-ray sources seen in the Chandra images. Nevertheless, the optimal source region, described below, is large enough to enclose the whole source (see Figure 1).

Source regions for the XMM observations were determined using the SAS task eregionanalyse, which calculates the circular radius (between \(\sim10-13\arcsec\)) and maximizes the signal-to-noise ratio, accounting for the PSF fraction, source counts, and local background. Background regions were selected by visual inspection (see solid yellow circle in panel (f) of Figure 1), placed as close to the target as possible and always on the same CCD, while avoiding contaminating sources. They typically have radii between five and eight times larger than the associated source region. We confirm that there is negligible contamination of the source region from diffuse gas in the host galaxy by trying out different background regions (dashed yellow circles in panel (f) of Figure 1) that are smaller, but closer to the source region. The number of counts per arcsec\(^2\) in all four background regions are similar and the net count rate for the source region ranges from 24–27 in EPIC-PN, with the different backgrounds.

Source and background spectra were extracted from filtered event files in which high-background intervals had been excluded. The filtering was based on the 10—12 keV light curves, removing time periods where the count rate exceeded 4 counts/s for EPIC-PN and 3.5 for EPIC-MOS, where the EPIC-PN observation is the most affected. The total exposure time for the XMM observation is 26.6 ks, while the remaining good exposure times are \(\sim 15\) ks for PN and \(\sim 25\) ks for MOS (see Table [tab:observaions] for details). Source spectra were extracted with aperture corrections and redistribution matrix (rmf) weighting applied using the default SAS settings optimized for point sources.

Spectral binning was performed with grppha, adopting one count per bin for all EPIC spectra (PN, MOS1, and MOS2), utilizing C-statistics for low count spectra.

3 Analysis and Results↩︎

In this study, we investigate the nature of the powering source of SN 2018bsz. The spectra are therefore fitted with both a power law, suggesting a millisecond (ms) pulsar origin, and a plasma model, representing an interaction driven power source. All fits include Galactic absorption that is modeled by the XSPEC component tbabs [49], where the abundance is set according to [49] and the column density of hydrogen is fixed at \(\rm NH=2.02\times 10^{21}~cm^{-2}\), obtained from the NHtot tool² [50]. Additional fits including intrinsic absorption from the host galaxy resulted in column densities consistent with zero in all model configurations described in the following sections, and where therefore excluded from the final fits.

All spectral analysis is performed using XSPEC ver. 12.15.0c [51] with the C-stat fit statistic. Chandra spectra are fitted within 0.5–8 keV, and similarly the fluxes and luminosities are derived within this interval. For the luminosities used in the pulsar analysis (Section 4.1) we extrapolate the model fits up to 10 keV and derive luminosities in the 2–10 keV energy range. XMM spectra are fitted within 0.5–10 keV and luminosities are derived in the same interval. All fluxes and luminosities presented in this work have been corrected for Galactic absorption and are presented with 90% error bars, or as \(3\sigma\) upper limits in Tables [tab:results95chandra] and [tab:results95xmm]. Below we describe the analysis methodology and presents the results for Chandra and XMM in subsections 3.1 and 3.2, respectively.

3.1 Chandra↩︎

A clear X-ray source is detected (\(> 3\sigma\)) in all four Chandra observations, located typically \(\sim0.3\arcsec\) from the updated SN coordinates [43] (up to \(\lesssim1\arcsec\) in ObsID 20289). The 1\(\sigma\) positional uncertainties from wavdetect (the statistical centroid variances) ranges between 0.11–0.23\(\arcsec\) for RA and 0.15–0.33\(\arcsec\) for DEC for the four Chandra observations. Furthermore, the systematic astrometric error for the Chandra CSC 2.0 catalog is 0.71\(\arcsec\) at 95% confidence [52]. The source is significantly detected (\(>3\sigma\), see Table [tab:observaions] for details) in the 0.5–8 keV images across 87–304 days after explosion, and we associate it with SN 2018bsz. The spectra are well described by power-law models modified by Galactic absorption, with hard spectral indices (\(\Gamma\)) in the range 0.9–1.9 (see Table [tab:results95chandra]). Attempts to fit thermal plasma models (mekal) resulted in temperatures driven to the upper bound of the model grid (78 keV), effectively reproducing a power-law shape. This behavior is expected given the scarcity of counts below 1 keV (see the binned spectra in Figure 2 (a)).

Each individual observation contains a limited number of counts, and the final fits yield in hard spectral indices with large associated uncertainties, which further show no significant temporal evolution (Table [tab:results95chandra], upper panel). This motivated a simultaneous fit of the four Chandra observations, tying \(\Gamma\) across epochs while allowing the normalization to vary. The resulting model is shown in Figure 2 (a). The corresponding luminosities, presented in Table [tab:results95chandra] (lower panel) and Figure 3, are consistent with being constant over time.

Inspection of the Chandra images and wavdetect results revealed three fainter X-ray sources near the SN position and within the XMM SN extraction region (see green boxes S1-S3 within the cyan dashed circular regions in Figure 1 (a)–(d)). These sources are detected with lower- or equivalent significance as the SN in all Chandra epochs, and with fewer net counts: S1:\(\sim\)​3–6, S2:\(\sim\)​4–6 and S3:\(\sim\)​2–9 compare to the SN in Table [tab:observaions]. Due to the lower angular resolution in XMM the emission from these sources will blend with the SN emission in the XMM spectra. We therefore analyze the combined S1–S3 Chandra spectra so that their contribution can be accounted for in the spectral model for the XMM fit.

Spectral extraction for these contaminating sources is performed using an annulus defined as the Chandra source region subtracted from the XMM EPIC-PN circular SN source region (the contamination annulus; Cont.A, illustrated by the region between the cyan dashed and white solid circles in Figure 1 (a)–(d)). These spectra are produced using the same procedures and background regions as for the SN spectra described in Section 2.1. The Cont.A spectra from all four epochs are fitted simultaneously with model parameters tied between observation, after confirming that they show no significant temporal evolution.

We find that the spectrum of the contamination region is well described by a power-law model with best fit parameters \(\Gamma_{\rm cont. A}=1.07_{-0.53}^{+0.58}\) and normalization \(N_{\rm Cont. A}=8.77_{-3.97}^{+5.73}\times10^{-7}\). The true nature of S1–S3 is difficult to asses as there are limited number of multi-wavelength observations in this field. It is possible that the emission comes from the hot gas in the host galaxy, as Pan-STARRS pre-explosion images show that the system is an ongoing merger of two dwarf galaxies [43]. A caveat is that such emission is typically well described by a low-temperature thermal model [53]–[55], while S1–S3 is best fit with a hard power law. Alternatively, the emission could come from high mass X-ray binaries or similar.

3.2 XMM↩︎

Running ewavelet on the XMM images, we detect (at \(\ge3\sigma\)) a source at \(\sim3.6\arcsec\) from the updated SN position [43]. The absolute astrometric uncertainty of XMM EPIC is \(1.2\arcsec(1\sigma)\)³, making the optical SN position and the detected X-ray source consistent within \(3\sigma\). As noted above, additional X-ray sources were identified near the SN position in the Chandra images. Due to the lower angular resolution of XMM, these different sources are unresolved and identified as a single source by ewavelet, explaining the slightly larger offset from the SN coordinate. In particular, the algorithm places the center of the source centroid between the SN coordinate and the closest contaminating source, S2 (Figure 1e). S2 is one of the brightest nearby X-ray sources and is likely to significantly contaminate the XMM spectrum, especially if the SN emission has faded since the Chandra observations. The best fit Cont.A model described at the end of Section 3.1, predicts 15 net counts in the XMM EPIC-PN spectrum, compared to the total number of 26 net counts in the spectrum for the energy range 0.5–10 keV.

To isolate the SN contribution and account for blended emission from the unresolved nearby sources, we introduce an additional power-law component as described at the end of Section 3.1. This component is included in the XMM SN fit with \(\Gamma_{\rm Cont. A}\) and \(N_{\rm Cont.A}\) frozen at previously determined values. Three scenarios are considered using the best fit \(N_{\rm Cont. A}\) and its lower and upper 90% error bounds, corresponding to the average, minimal, and maximal levels of contamination from S1–S3. The resulting model for the XMM fits is thus tbabs(pow\(_{\rm Cont.A}\)+ pow), with the unabsorbed flux extracted from the second power-law component representing the SN, having free \(\Gamma\) and normalization.

The PN, MOS1 and MOS2 datasets are fitted simultaneously, with both \(\Gamma\) and normalizations tied across instruments but allowing for 8% relative variation [56] through a multiplicative constant (fixed to 1 for PN, free between 0.92–1.08 for MOS1 and MOS2). Table [tab:results95xmm] includes the resulting luminosities and limits for the power-law model in the upper panel, while Figure 2 (b) shows the spectral fit corresponding to the top row of Table [tab:results95xmm]. For the best-fit level of contamination, the unabsorbed SN luminosity is constrained at \(3\sigma\) in the 0.5–10 keV band (but not in 2–10 keV), implying a detection in the full band. The SN is also detected in the case of minimal contamination as expected, while we only obtain an upper limit on the SN luminosity in the maximal contamination scenario. The soft \(\Gamma\) for all contamination scenarios suggest a significant softening of the SN spectrum since the early Chandra observations, with luminosities comparable or slightly lower at the time of the XMM observation.

Although CSM-ejecta interaction is often invoked for SLSN-II and SLSN-IIn in particular, it may also be relevant for SLSN-I [15], [45], [57]. Motivated by the soft \(\Gamma\) of the late XMM spectra, we additionally fit a mekal model alongside the contamination component, giving a total model of tbabs(pow\(_{\rm Cont.A}\)+ mekal). Unlike the Chandra spectra, the late-time XMM spectra allow us to constrain the plasma temperature when adopting the best-fit \(N_{\rm Cont.A}\) (Table [tab:results95xmm], first row of the lower panel), yielding in a soft spectrum with \(kT\sim0.2~\)keV. We note, however, that only in the minimal contamination scenario can we constrain the 0.5–10 keV luminosity, but with a large upper bound for \(kT\), which is outside the spectral range and therefore equivalent to unconstrained.

4 Discussion↩︎

We detect SN 2018bsz in five observations spanning 87 to 1253 days after explosion, with the final XMM observation requiring additional caution in its interpretation due to the nearby contaminating sources. This makes SN 2018bsz the only other X-ray detected SLSN-I besides the unusually bright SLSN-I SCP06F6 [34]. The SN luminosity remains nearly constant as seen in Figure 3, while the spectrum evolves from the very hard power law seen in the early Chandra data to a softer spectrum at later times. Although the \(\Gamma\) of the late XMM spectrum is weakly constrained, the low plasma temperature obtained from alternative fits likewise supports a softer spectral shape (Table [tab:results95xmm]).

In the following sections, we explore the two commonly invoked powering mechanisms for SLSNe: spin-down of a rapidly rotating magnetar and interaction between the ejecta and CSM, and assess their ability to account for the X-ray emission in SN 2018bsz.

4.1 Magnetar origin↩︎

Millisecond magnetars are often invoked as the power source of SLSN-I (e.g. [4], [5]). We therefore test whether the X-ray results of SN 2018bsz are compatible with such an engine.

We use the common pulsar spin-down model where pulsars are modeled as rotating dipoles in vacuum [58], [59]. The spin-down energy is given by \[\label{eqn:Erot95dot} \dot{E}_{\rm rot}=-4\pi^2I\frac{\dot{P}}{P^3} ,\tag{1}\] where \(I ~(=10^{45}~\rm g~cm^{2})\) is the moment of inertia, \(P\) is the spin period and \(\dot{P}\) its time derivative. \(P\) evolves as

\[\label{eqn:p} P=\left[ P_0^2+\left(\frac{16\pi^2R^6B^2}{3Ic^3}\right)t\right]^{1/2},\tag{2}\] with \(B\) being the magnetic field which is assumed constant with time, \(P_0\) is the pulsar spin down at birth, \(c\) is the speed of light, and \(R ~(=12~\rm km)\) is the radius of the neutron star. We adopt the magnetic field prescription from [60],

\[\label{eqn:b95field} B^2= \frac{3}{2}\frac{Ic^3}{(2\pi)^2R^6}P\dot{P},\tag{3}\] which allows us to express \(\dot{E}_{\rm rot}\) directly as a function of \(P_0\) and \(B\). A fraction \(\eta=L_{\rm X}/\dot{E}_{\rm rot}\) of the spin-down power is assumed to emerge as X-ray emission in the 2–10 keV range. Since the X-ray efficiency of young pulsars and magnetars are poorly constrained, we adopt the empirically derived Crab value, \(\eta=0.022\). This value was obtained from \(\dot{E}_{\rm rot}=4.6\times10^{38}~\rm erg~s^{-1}\) [61] and the X-ray luminosity, which was derived from an extrapolated power law with \(\Gamma=2.106\) [62], extrapolated from the energy range from 3–50 keV into 2–10 keV.

To constrain \(P_0\) and \(B\), we run a Markov Chain Monte Carlo (MCMC) fitter from the emcee package [63], using the pulsar model above. The fit accounts for asymmetric luminosity uncertainties by using a modified log-likelihood function, which converts the reported 90% confidence intervals and 3\(\sigma\) upper limit to equivalent 1\(\sigma\) values. We run the fit over the 2–10 keV Chandra luminosities and XMM limit (black solid markers in Figure 3), where the XMM limit corresponds to the average Cont.A results. We assume no ejecta absorption in these fits, which could be due to either ejecta asymmetries or ionization of the ejecta from the magnetar, discussed further below. Furthermore, we assume flat initial priors in log space (\(\log_{10}(P_0)\) and \(\log_{10}(B)\)), with \(\rm0.5~ms< P_0<1~s\) and \(\rm 10^{11}~G<B<10^{16}~G\). The MCMC fit is initialized with 1000 walkers, evolved over 10,000 steps, with the first 200 discarded as burn-in.

The MCMC fit results in \(\log B=13.1^{+0.6}_{-1.2}\) and \(\log P_0=-1.8^{+0.3}_{-0.6}\) (corresponding to \(P_0\sim16~\rm ms\)), as shown in the corner plot in Figure 4. The fit is also illustrated by the solid black line in Figure 3. If we instead use the XMM limit from min Cont.A, we get a similar result with \(\log B\sim13.5\) and \(\log P_0\sim-1.6\) and a good fit to the data. Using the max Cont.A value yields a fit that deviates substantially from the observed data, resulting in parameter estimates of \(\log B\sim 14.9\) and \(\log P_0\sim-1.0\).

We compare these results to those inferred from fits of the optical light curve of SN 2018bsz by [11], who found best-fit values of \(P_0=0.9_{-0.2}^{+0.3}\) ms and \(\log B=14.4_{-0.2}^{+0.3}\) (from here on we refer to these best-fit values values as \(BP_{\rm optical}\)), with corresponding \(1\sigma\) confidence interval. We use \(BP_{\rm optical}\) in our pulsar model and include the results in Figure 3 and Figure 4 as a green line for the best-fit values, and green-shaded regions for the 1\(\sigma\) confidence interval. From the corner plot in Figure 4, we see that the MCMC results are generally not in line with \(BP_{\rm optical}\) (green region) nor with the sample results from light curve fits of 262 SLSNe from [11] (pink region), including SN 2018bsz. The corner plot also reveals a strong degeneracy between \(P_0\) and \(B\) (bottom left panel of Figure 4). Although the \(B-P_0\) degeneracy allows for some overlap with \(BP_{\rm optical}\) at larger \(B\), these solutions lie outside the 2\(\sigma\) contours and ar disfavored (see bottom left panel of Figure 4). A larger \(B\), such as suggested by \(BP_{\rm optical}\), is clearly not compatible with the temporal evolution of our observational results, as seen in Figure 3. The observed light curve requires a lower \(B\) for a more flat evolution.

We also explore whether a lower \(\eta\) could reconcile the X-ray light curve with \(BP_{\rm optical}\) by repeating the MCMC with fixed \(B\) and \(P_0\) at \(BP_{\rm optical}\) and fitting for \(\eta\) with flat priors (\(10^{-6}<\eta<1\)), we obtain \(\eta=0.01\). This value is below that of the Crab and smaller than typical efficiencies for pulsars in X-rays (e.g. [64], [65]), but foremost, it does not change the temporal luminosity evolution which is the main discrepancy between our MCMC results and \(BP_{\rm optical}\).

4.1.1 Impact of ejecta absorption↩︎

The X-ray luminosities derived from our spectral fits assume negligible ejecta absorption, since the low spectral count rates prevent us from fitting for intrinsic absorption directly and the spectra do not show strong absorption features (see Figure 2). Such an assumption of no ejecta absorption, is not necessarily expected a priori. The ejecta masses inferred for SLSN-I are typically in the range of \(2-30~M_\Sun\) [10], [66], [67], where such ejecta are generally expected to remain optically thick to soft X-rays for a significant time after explosion. For SN 2018bsz specifically, [11] inferred an ejecta mass of \(M_{\rm ej}=9.5_{-4.4}^{+7.5}~\rm M_\Sun\) based on optical light-curve modeling.

The apparent lack of strong absorption signatures may indicate that the magnetar is in fact fully absorbed and that the X-ray emission is due to CSM interaction (a scenario discussed in Section 4.2). This could also explain the discrepancy between the \(B,P_0\) inferred from optical data and our X-ray data. Alternatively, the absorption of the magnetar emission is very low, either due to geometric asymmetries of the ejecta, or the ejecta have become highly ionized by radiation from the magnetar. In the latter case, the opacity to soft X-rays can be reduced as the PWN progressively ionizes the surrounding ejecta [32]. Once the ionization front reaches the outer ejecta surface, X-rays can escape freely, a process referred to ar ionization breakout. A related possibility is that the emission comes from an extended PWN, in which case the observed spectrum may represent a superposition of emission components with different \(\Gamma\) and absorption, which could further obscure clear absorption signatures.

Following [36], we estimate the ionization breakout time (\(t_{\rm ion}\)) and the corresponding X-ray luminosity at \(t_{\rm ion}\) using the magnetar parameters \(BP_{\rm optical}\) together with the ejecta mass \(M_{\rm ej}=9.5_{-4.4}^{+7.5}~\rm M_\Sun\) and velocity \(v_{\rm ej}=15,127^{+1873}_{-2345}\rm km~s^{-1}\) inferred for SN 2018bsz by [11]. For this parameter set, the ionization breakout occurs at \(t_{\rm ion}\approx 72~\rm yr\) with \(L_{\rm x}(t_{\rm ion})\approx1.6\times10^{37}~\rm erg~s^{-1}\). Adopting the most favorable combination within the quoted \(1\sigma\) confidence intervals of \(BP_{\rm optical}\), \(M_{\rm ej}\) and \(v_{\rm ej}\), corresponding to the lowest \(B\) and lowest ejecta absorption, yields in \(t_{\rm ion}\approx 5~\rm yr\) with \(L_{\rm x}(t_{\rm ion}) \approx 9.3\times10^{39} ~\rm erg~s^{-1}\). In both cases, the predicted ionization breakout occurs beyond the temporal window of our observations, though with a detectable flux of \(6.3\times10^{-15}~\rm erg~s^{-1}~cm^{-2}\) for the most favorable scenario, at a distance of 111 Mpc. Given the late \(t_{\rm ion}\), it is unlikely that the X-ray emission is due to ionization breakout from a magnetar with \(BP_{\rm optical}\) properties.

We note, however, that the ionization breakout models are one-dimensional and assume spherical symmetry. In more realistic scenarios, asymmetric ejecta structures may lead to low-density channels with reduced absorption, giving an asymmetric ionization breakout [68]. Such effects could, in principle, relax the constraints derived from spherically symmetric models and partially reconcile a magnetar-powered scenario with the observed luminosities. However, reproducing both the high X-ray luminosities and the relatively flat temporal evolution would still require a finely tuned combination of magnetar parameters and ejecta properties. Even when accounting for asymmetric ionization and reduced ejecta absorption, it remains challenging for a magnetar with \(BP_{\rm optical}\)-like properties to account for the observed X-ray evolution.

Figure 5 illustrates the excluded parameter space in the \(B-P_0\) plane based on the 2–10 keV luminosity limit from the XMM observation (best-fit Cont.A), which provides the strongest constraints due to the late epoch. The figure shows excluded regions for the case of negligible ejecta absorption and for two examples of absorption by ejecta with \(M_{\rm ej}\) and \(v_{\rm ej}\) within the constraints from [11]. The curve corresponding to no ejecta absorption (black) is derived using Equations 1 –3 , while the curves including ejecta absorption are derived using the ionization-breakout formalism (Equations 1 and 4 of [36]). We note that our adopted pulsar models differs somewhat from the setup used by [36], resulting in spin-down luminosities higher by a factor \(\sim 3\).⁴ This analysis excludes pulsars born with low \(B\) and short \(P_0\), which would have fully ionized the ejecta at earlier times and produced significantly higher luminosities than observed. Magnetar properties inferred for SN 2018bsz by [11] (green cross and shaded region in Figure 5) are consistent with the observations for expected levels of X-ray absorption.

We note, however, that the complex optical light curve implies that there are likely systematic uncertainties on \(BP_{\rm optical}\) that are not fully captured by the shaded region in Figure 5. Using the online data provided by [11], we reconstruct the light curves since the figure included in the paper is too small for a detailed evaluation. We find that the post-peak plateau is not fully captured by the fit, which may indicate that additional physical processes beyond magnetar powering contribute to the observed the light curve. This interpretation is consistent with [45], who suggested that the post-peak plateau could originate from an optically thick, H-rich CSM.

In summary, the temporal evolution of the observed X-ray emission is not consistent with a magnetar central engine characterized by \(BP_{\rm optical}\), but instead requires a lower \(B\) to reproduce the flatter luminosity evolution. While the earliest Chandra observation could be marginally compatible with \(BP_{\rm optical}\) in the presence of little absorption, the later observed luminosities exceed the predictions of the magnetar model. A more plausible explanation is that the X-ray emission is dominated by interaction between the ejecta and CSM.

4.2 CSM interaction↩︎

An alterative energy source proposed to explain the brightness of some SLSN-I (e.g. [15], [69], [70]), including SN 2018bsz [45], is the ejecta-CSM interaction. For SN 2018bsz, [45] identified hydrogen emission lines emerging \(\sim10-40~\rm days\) after peak optical luminosity, with spectral evolution more reminiscent of common SNe IIn than other SLSNe-I, interpreted as interaction with an asymmetric, disk-like CSM. There are also some SLSN-I which show similar behaviors, where the spectra develop H emission at \(\gtrsim100~\)days, indicating presence of CSM interaction [13], [15], [29]. SN 2017ens is such an example, which showed H-line emission and CSM interaction that was detected in radio at \(\sim1250~\)days [18], [29]. These signatures suggested a dense shell-like and H-rich CSM structure.

Our results for SN 2018bsz could likewise be interpreted as signatures of strong early CSM interaction. The early Chandra spectra are best described by a hard power law with \(\Gamma=1.32\) (see Table [tab:results95chandra] for details), whereas by the XMM epoch the spectrum has softened to \(\Gamma<4.0\) and is equally well fit by a thermal plasma model with \(kT=0.20~\rm keV\) (Table [tab:results95xmm]). The spectral softening is significant at \(2\sigma\) (assuming average Cont.A). This behavior is consistent with an initially strong interaction that weakens with time (e.g. [71], [72]). The little to no decline in luminosity further favors a CSM interaction scenario over a magnetar powered one, which would be expected to fade more rapidly. We note, however, that the inferred XMM luminosity depends on the assumed level of contamination; in the maximal-contamination case, a significant decrease is observed. Although our \(\Gamma\) uncertainties are large or in some cases only lower limits (Table [tab:results95xmm]), we also find constrained plasma temperatures indicative of a soft spectrum at late times. Similar spectral softening in X-rays and late-time hydrogen-line emergence have been observed in ordinary stripped type SNe transitioning into CSM-interacting systems (e.g. [16]–[18], [20], [73]).

The inferred X-ray luminosities and limits of SN 2018bsz are consistent with the highest luminosities expected from CSM interaction as suggested by [34], assuming a high mass-loss rate of \(\sim10^{-4}~M_\Sun~\rm yr^{-1}\). Higher X-ray luminosities would require larger ejecta masses and higher mass-loss rates. Using the mass-loss-luminosity relation from [72], we estimate the progenitor mass-loss rate based on each observation, using the formula:

\[Lx\approx 3\times 10^{39}g_{ff}C_n\left(\frac{\dot{M}_{-5}}{v_{w1}}\right)^2t_{10}^{-1},\]

where \(C_n=(n-3)(n-4)^2/[4(n-2)]\) for the reverse shock (adopting \(\rm n=7\)), \(g_{ff}\) is the Gaunt factor in order unity, and \(\dot{M_{5}}\) is the mass-loss rate in units of \(10^{-5}~M_\Sun~\rm yr^{-1}\). The wind velocity, \(v_{w1}\), is in units of 10 km s\(^{-1}\), which we adopt for stripped type SNe progenitors with typical value of \(\sim 1000~\rm km~s^{-1}\) [22], [72]. Finally, \(t_{10}\) is time in units of 10 days. We stress that the estimated \(\dot{M_{5}}\) are only indicative as these calculations assume a steady wind and rely on additional assumptions, such as \(v_{w1}\) and \(n\).

The resulting mass-loss rates for the Chandra epochs averages to \(\sim 5.4\times10^{-3}~M_\Sun~\rm yr^{-1}\) and for the XMM epoch \(\sim 1.4\times10^{-2}~M_\Sun~\rm yr^{-1}\) adopting the luminosity from the thermal fit with the average Cont.A. These values exceed those of typical ordinary stripped SNe [22], [72], most SLSNe-I [36] and that of SN 2017ens [18], [29], and are closer to some mass-loss rates for SLSN-II [33], though the ranges for SLSN-II often reaches \(\gtrsim 0.01~M_\Sun~\rm yr^{-1}\) [74], [75]. Our mass-loss rates are most comparable to those of the subset of stripped SNe that evolve into strong interacting systems [73]. It is noteworthy that multiwavelength observations of SN 2018bsz already indicate atypical behavior for an SLSN-I, including an early optical plateau and hydrogen emission [43], [45]. This agrees well with the slightly larger mass-loss rates found compared to other SLSN-I, and could have been ejected as previously discussed, in a violent eruption of the progenitor star.

Given the derived mass-loss rates, we estimate the time of the progenitor mass loss to assess whether the X-ray detections from Chandra could originate from the disk-like CSM structure proposed in [45]. We use \(t_{\rm CSM}=t_{\rm Xray,obs}v_{\rm shock}/v_w\) and estimate the shock velocity (\(v_{\rm shock}\)) for a spherically symmetric wind following [72]. We adopt the same assumptions for \(n\) and \(v_w\) as previously, as well as the values of \(E_{\rm ej}\) and \(M_{\rm ej}\) for SN 2018bsz from [11]. This yields \(t_{\rm CSM}\approx 5-18\) yr before explosion, implying a mass-loss episode lasting for \(\sim13~\rm yr\). The corresponding CSM radius is \(R_{\rm CSM\approx (1.7-5.8)\times10^{16}~\rm cm}\). For comparison, [45] derived an upper limit for the disk-like CSM radius of \(<6.5\times10^{15}~\rm cm\). While our estimated radius is larger, there are many simplifying assumptions, including considering spherical symmetry rather than a disk-like geometry, and we conclude that the X-ray emission is compatible with an origin in the CSM disk proposed by [45].

5 Summary & Conclusion↩︎

We detect X-ray emission from SN 2018bsz in Chandra observations spanning 87 to 304 days after explosion at \(3\sigma\) significance, with an average luminosity of \(L=1.1\times 10^{40}~\rm erg~s^{-1}\) in the 0.5–8 keV energy range. A possible detection is also seen in the 0.5–10 keV XMM observation at 1253 days post explosion; however, this result is sensitive to assumptions about nearby contaminating sources and should be interpreted with caution. Independent of the XMM detection, SN 2018bsz is the second SLSN-I detected in X-rays, and the third of all SLSN detected in X-rays. It exhibits the lowest X-ray luminosity among this small sample.

Interpreted within a magnetar-powered framework and assuming negligible ejecta absorption, the observed X-ray luminosities favor magnetar parameters with a comparatively weak magnetic field (\(\log B\sim\) 13.1) and a long initial spin period (\(P_0\sim16~\rm ms\)). These values lie outside the range typically inferred for ms magnetars (see e.g. [10], [11]), and are difficult to reconcile with the energetics required to power the superluminous optical emission.

Adopting instead the magnetar parameters inferred from the optical light curve by [11], we find that the predicted X-ray luminosities are lower than the observed values at all but the earliest epoch, and fail to reproduce the relatively flat temporal evolution. Accounting for ejecta absorption further increases this discrepancy. While asymmetric ejecta and ionization effects may reduce the effective absorption along the line of sight, our estimated timescale of ionization breakout occurs at best \(\sim2\) yr after our observational window assuming \(BP_{\rm optical}\), though at an observable flux (\(6.3\times10^{-15}~\rm erg~s^{-1}~cm^{-2}\)). Although such effects could partially alleviate the tension, reproducing both the X-ray luminosity level and temporal evolution would be difficult given the ms magnetar properties of \(BP_{\rm optical}\).

Taken together, these results disfavor a standard magnetar-powered origin for the observed X-ray emission. Instead, the data are more readily explained by strong early ejecta-CSM interaction. The X-ray spectrum is initially hard at epochs 87–304 days, with a photon index of \(\Gamma=1.32\), and softens substantially by the XMM observation at 1253 days, where it is consistent with either \(\Gamma<4.0\) or a thermal plasma temperature of \(kT=0.20~\rm keV\). This spectral evolution, combined with the relatively slow luminosity decline, is consistent with a scenario in which strong early interaction weakens over time.

From the observed X-ray luminosities, we estimate progenitor mass-loss rates of order a few \(10^{-3}~\rm M_\Sun~yr^{-1}\). These values exceed those typically inferred for ordinary stripped SNe and SLSNe-I, and are more comparable to the subset of stripped SNe which evolve into strongly interacting systems following eruptive mass-loss episodes that produce dense, shell-like CSM.

The presence of CSM interaction and elevated mass-loss rates is consistent with previous multi-wavelength studies of SN 2018bsz [43], [45], [46], which revealed hydrogen emission lines and evidence for interaction with a disk-like CSM [45]. Together, these results suggest that SN 2018bsz represents a distinct subclass within the SLSN-I population, where CSM interaction plays a significant role in shaping the high-energy emission, and where a magnetar central engine, if present, is unlikely to dominate the observed X-ray output.

References↩︎