© The Authors 2024

1 Introduction

The star β Pic is a young nearby star of about 20 million years old (Miret-Roig et al. 2020), located 19.3parsecs away from us and harbouring a gaseous and dusty debris disc seen edge-on (Smith & Terrile 1984; Vidal-Madjar et al. 1986; Kalas et al. 2000; Roberge et al. 2006; Apai et al. 2015; Brandeker et al. 2016). Two massive planets were discovered through imagery, radial velocity (RV), and mutual dynamical perturbations (Lagrange et al. 2010, 2019; Snellen & Brown 2018; Nowak et al. 2020; Lacour et al. 2021). Therefore, β Pic offers a unique opportunity to investigate the early stages of a planetary system.

For more than 30 yr, exocomets transiting the star have been detected using spectroscopy, probing the gaseous part of the cometary comas and tails (Ferlet et al. 1987; Beust et al. 1990; Vidal-Madjar et al. 1994; Kiefer et al. 2014b; Kennedy 2018). First named by the circumlocution falling evaporating bodies (FEB, see Beust et al. 1990; Vidal-Madjar et al. 1994), exocomets in β Pic are the analogous of comets in our own Solar System, that is, icy bodies whose volatiles evaporate when approaching the star (Strøm et al. 2020), producing huge gaseous and dusty tails that can be detected using transit spectroscopy and transit photometry (Zieba et al. 2019; Pavlenko et al. 2022; Lecavelier des Etangs et al. 2022). One of the most remarkable aspects of the exocomet population in the β Pic system is its unique accessibility to detailed analysis. The close proximity of this young star, coupled with orientation favouring transits of bodies with low inclination relative to the disc, allows us to scrutinise this exocomet population with a level of detail currently unattainable for most other systems. Although some systems seem to offer interesting possibilities for similar studies, such as 49 Cet (Roberge et al. 2014; Miles et al. 2016), HD 172555 (Kiefer et al. 2014a, 2023), and a few others identified by surveys (Montgomery & Welsh 2012; Welsh & Montgomery 2013; Rebollido et al. 2020), the number of comets in the β Pic system is much larger than in any other known exocometary system. β Pic also remains the only system for which the presence of exocomets has been linked to the presence of cold CO that, in the absence of H₂, must be continuously resupplied (Vidal-Madjar et al. 1994; Jolly et al. 1998; Roberge et al. 2000; Lecavelier des Etangs et al. 2001), most likely by evaporating comets.

However, despite the numerous observations of exocometary transits in the β Pic system, most previous studies focussed on modelling the observed spectroscopic events to conclude on the exocometary nature of the phenomenon (Lagrange-Henri et al. 1989; Vidal-Madjar et al. 1994), numerical simulations (e.g. Beust et al. 1990; Beust & Tagger 1993), the detection of the widest variety of species possible (e.g. Vidal-Madjar et al. 1994), or the derivation of their orbital properties (e.g. Kiefer et al. 2014b; Kennedy 2018). As a result, only very few studies concluded on the measurement of physical properties of the observed comets themselves, with notably the work by Lagrange et al. (1995) for first estimates of column densities using equivalent widths, Mouillet & Lagrange (1995) for temperature and density of the gaseous coma, Kiefer et al. (2014b) for the estimate of the nuclei evaporation efficiency, and the analysis of the photometric transits by Lecavelier des Etangs et al. (2022) for the size distribution of the cometary nucleus. This situation results from the extreme difficulty in converting the observed absorptions into intrinsic physical quantities of the transiting comets. Here, using the advantage of a large number of electronic transition lines of the same species in a few exocomets observed with the Hubble Space Telescope, we aim to provide a detailed analysis of the transit spectra in order to better characterise the physical conditions in the gaseous clouds surrounding the comets’ nuclei.

The selection of the observational data and their analysis are described in Sect. 2. The first model and the derivation of the exocomet curves of growth are discussed in Sect. 3. The estimates of the column densities, temperature, and electron density are presented in Sects. 4 and 5. We discuss, summarise, and conclude in Sects. 6 and 7.

2 Observations

2.1 Archival data

For our study, we used mostly unpublished archival UV spectra of β Pic, obtained with the Space Telescope Imaging Spectrograph (STIS) onboard the Hubble Space Telescope (HST) within three different guest-observer programmes (programme #7512, PI A.-M. Lagrange, #14735 and #15479, PI F. Kiefer; see Table 1). For programmes #14735 and #15479 (executed in 2017/2018), no analysis of the data has yet been published. The data of programme #7512 (December 1997) was analysed in Roberge et al. (2000) to constrain the CI and CO content of the circumstellar disc, but no study of the exocometary absorption was performed.

These STIS observations were obtained using the high-resolution (R ~ 114 000) echelle grating E230H at seven different epochs. The raw data were processed using version 3.4.2 of the calstis pipeline released on January 19, 2018, providing 1-D flux calibrated spectra that are publicly available and have been downloaded from the MAST archive. The data are wavelength calibrated on a heliocentric reference frame.

All these spectra cover the ionised iron (Fe II) line series between 2715 and 2769 Å, rising from the 3d⁶(⁵D)4s–a⁴D metastable levels, in which cometary transits are easily detected (see Sect. 2.3). In addition, the two 1997 observations (programme # 7512) also cover the Fe II line series at 2400Å and 2600 Å, which rise from a wide range of energy levels.

2.2 Data analysis

Each of the STIS/E230H spectra covers a bandwidth of about 270 Å wide, with a set of 24–40 orders of about 12Å wide, partially overlapping. For each exposure, the fluxes of all the orders have been combined into a single spectrum and resampled on a common wavelength table with a resolution of 18 m Å⁻¹. The wavelengths were also corrected for the heliocentric radial velocity of β Pic, +20 km s⁻¹ (Gontcharov 2007), so that they correspond to the β Pic reference frame.

In a second step, to correct for variation of the flux calibration from one epoch to another, we renormalised the measured flux to a common flux level, independent of exocomet activity. For that, we selected a set of 65 ~1Å-wide regions free of any variations (i.e. mostly outside spectral lines) separated by 10 to 20 Å. By measuring their flux on these different regions, we obtained, for each spectrum, a set of reference points of the form (λ_i, ƒ_i)_{i = 1,2…} The reference points of each spectrum were then fitted by cubic splines, and each spectrum was divided by its corresponding spline. This flux normalisation allowed us to obtain data that can be easily compared, unveiling variable absorption features around many lines (Fig. 1).

Finally, as no significant spectral variation was observed within each set of three exposures taken on December 6 and December 19, 1997 (both times within a single HST orbit of ~90 min), these sets of exposures were merged together, yielding a single spectrum from 2128 Å to 2888 Å for each of the two epochs. For the 2017–2018 observations, no such merging was necessary, as only one exposure was used for each visit. As a result, we ended up with seven normalised spectra: two spectra covering the 2128 Å to 2888 Å wavelength range taken 13 days apart in 1997, and five spectra covering the 2665 Å to 2931 Å wavelength range taken at five different epochs in 2017–2018.

Fig. 1 Illustration of the recalibration method. Top panel: raw spectra show significant shape disparities, probably resulting from changes in instrument calibration throughout the years. These calibration disparities are hardly fitted by polynomials. Centre panel: zoomed-in view of 2680–2820 Å region, showing the splines fitted to each spectrum. The grey areas indicate the different spectral domains where the reference fluxes were measured (black crosses), which are then used to calculate a cubic spline for each spectrum. Bottom panel: zoomed-in view of 2735–2753 Å region after renormalisation by the cubic splines. The spectra taken at different epochs are shown with different colours; they are well superimposed in the region of the stellar continuum, while variable absorptions are clearly visible around the different Fe II lines between 2737 and 2750 Å.

2.3 Detected Fe II lines

Dozens of Fe II lines are detected in the seven analysed spectra, with deep absorption originating both from the circumstellar disc (stable) and from transiting evaporating comets (variable). These exocomet signatures are detected in three Fe II main series: a first one near 2750 Å, where the high number of different visits (7) allows for distinguishing the variable absorptions from the stable stellar spectrum fairly easily (Fig. C.1), as well as two other groups near 2400 and 2600 Å, only covered by the two 1997 observations (Fig. C.2). In addition, a few lines arising from highly excited levels also showed fainter variable absorptions, as seen in Fig. C.3. The complete list of lines used in the following sections is presented in Tables A.1 (for the strongest lines) and A.2 (for the lines arising from high energy levels).

2.4 Reconstruction of the comet-free spectrum

In order to quantify the variable features observed in these many Fe II lines, particularly in the series at 2750 Å, one must first compute a stellar reference spectrum without cometary absorption. To do so, we proceeded in two main steps:

• First, we used the ~2750 Å Fe II lines to identify, for each spectrum, the radial velocity domains where cometary absorption occurs, assuming it always occurs at the same radial velocities from one line to another (this assumption is in fact easily verified; see Sect. 2.5). This identification is made possible by the large number of spectra at these wavelengths, which allows us to easily distinguish, for each observation, the spectral domains where it is affected by cometary absorption and the spectral domains where it is not (see, e.g. Fig. C.1). More quantitatively, for each radial velocity v, we considered that a given spectrum is affected by cometary absorption when its flux at the velocity υ in the 2740 and 2756 Å lines (those are the strongest lines of the 2750 Å series) is at least 5% below the brightest observation. The list of the radial velocity domains where cometary absorptions are detected is given for each observation in Table 2: we note that, for each radial velocity, there is always at least one comet-free observation, thanks to the rather high number of independent visits.

• Second, for each wavelength pixel, we used the list of the radial velocity domains where comets are detected in each spectrum (Table 2), and the list of all the Fe II lines (Tables A.1 and A.2), to determine the set of observations that can be considered comet-free at this particular wavelength. We then averaged all these un-absorbed spectra (free of cometary absorption) to obtain the reference spectrum of β Pic. For instance, for the line at 2740 Å (Fig. 2), the comet-free spectrum between +50 and +120 km s⁻¹ is obtained by averaging the 2018-03-07 and 2018-12-15 observations.

We note in Fig. 2 that it is not completely clear weather the spectrum obtained through this calculation is indeed the β Pic comet-free spectrum. For instance, shallow comets may have been missed in our identification, leading us to underestimate the reference spectrum by a small fraction (typically a few percent). On the other hand, the reference spectrum might be slightly overestimated near + 30 km s⁻¹ because of higher flux at the corresponding wavelengths in the March 7, 2018 observation (see Fig. 2, pink line). However, the impact of these biases on the calculation of the comet-free spectrum is much smaller than the typical signature of a β Pic exocomet on its host star spectrum (which can cause absorption depths as high as 50%). The conclusions presented in this paper are thus not significantly affected by such systematic errors.

It should also be noted that this comet-free spectrum can still contain interstellar and circumstellar absorption. However, because these absorptions are stable, they do not impact the measurement of cometary absorption depth. They can thus be included in the stellar reference spectrum, as if they had a stellar origin.

Fig. 2 Variable absorptions in the 2740 Å Fe II line. The black line shows the reference comet-free spectrum of β Pic, obtained by averaging the unabsorbed spectra at each individual wavelength. The deep and narrow absorption near 0 km s−1 is due to the stable component of the circumstellar disc of β Pic, which could be slightly blueshifted compared to the stellar radial velocity (see Brandeker 2011).

2.5 Cometary absorption

Cometary absorption spectra are obtained by dividing each individual observation by the comet-free spectrum previously calculated. This provides, for each observation, the absorption spectrum of the transiting gas: this spectrum is equal to 1 in spectral ranges where the transiting comets are transparent, and tends to zero in the presence of an optically thick tail. These spectra show the striking correlation between the absorption features observed in various Fe II lines, within the same observations. For example, Fig. 3 shows, for two separate visits, the superposition of their cometary absorption spectra around ten different Fe II lines, clearly revealing the transit of numerous comets within each of the observations.

The analysis developed in the next sections is focussed on the ‘absorption depth’ and ‘average absorption depth’ of exocomets. These two quantities are defined as follows:

• For a given line (l, u) (where l and u denote the lower and upper levels of the transition, respectively) and some given radial velocity υ, the absorption depth of exocomets is defined as the fraction of the stellar flux absorbed by exocomets at the velocity v in the line (l, u) in a given observation. If we denote the cometary absorption spectrum computed previously by ${\tilde{F}}_{lu}(\upsilon )$, this absorption depth writes as $${\text{abs}}_{lu}(\upsilon )=1-{\tilde{F}}_{lu}(\upsilon ).$$

  The absorption depths thus range from 0 in the absence of transiting comet to one when an optically thick cloud eclipses 100% of the stellar surface. The error bars on the measured absorption depths are estimated by accumulating all the errors on each of the used HST spectra as tabulated by the calstis pipeline.

• Then, to improve the signal-to-noise ratio (S/N) of our absorption measurements, the exocomet absorption depth can be averaged over an extended radial velocity range. For some radial velocities υ₁ < υ₂ and some line (l, u), we thus introduce the mean absorption depth of exocomets, which is defined as the average value of abs_1u(υ) over the RV domain [υ₁, υ₂]: $${\overline{\text{abs}}}_{lu}\left({\upsilon }_1,{\upsilon }_2\right)=\overline{\left\{{\text{abs}}_{lu}(\upsilon ),\upsilon \in \left[{\upsilon }_1,{\upsilon }_2\right]\right\}}.$$

3 Model description

3.1 Principle

In order to interpret the absorption features due to the transit of exocomets, we used a modified version of the cometary model that Kiefer et al. (2014b) introduced to analyse the absorptions in the 3900 Å Ca II doublet of β Pic. The gist of the model is to assume that the absorption depth (relative to the comet-free spectrum; see Sect. 2.5) observed at a radial velocity υ in a set of lines from a given chemical species is produced by a gaseous cloud of homogeneous column density, covering a fraction α(υ) of the stellar disc surface. Under these assumptions, the absorption depth of a transiting comet is given, for any of the considered lines, by $${\text{abs}}_{lu}(\upsilon )=\alpha (\upsilon )\cdot \left(1-{\text{e}}^{-{\tau }_{lu}(\upsilon )}\right),$$

with l and u being the lower and upper levels of the line (noted below (l, u)) and τ_1u(υ) the optical thickness of the cometary cloud, which can be expressed as $${\tau }_{lu}(\upsilon )=\frac{hc}{4\pi }{B}_{lu}{N}_l\varphi (\upsilon ),$$

where Bl_u is the Einstein coefficient for photon absorption of the considered line, N_l the column density of the lower energy level of the transition and ϕ the radial velocity distribution within the cloud (∫ ϕ(υ)dυ = 1). The Einstein coefficient B_lu is linked to the oscillator strength of the transition f_lu through $$\frac{hc}{4\pi }{B}_{lu}=\frac{{\text{e}}^2}{4{ϵ}_0{m}_{\text{e}}c}{f}_{lu}{\lambda }_{lu}.$$(1)

Then, if we average the cometary absorption depth over a radial velocity range [υ₁; υ₂] where τ_lu(υ) is assumed to be constant, we can re-write abs_lu and τ_lu as $$\begin{array}{c}{\overline{\text{abs}}}_{lu}\left({\upsilon }_1,{\upsilon }_2\right)=\overline{\alpha }\left({\upsilon }_1,{\upsilon }_2\right)\cdot \left(1-{\text{e}}^{-{\tau }_{lu}\left({\upsilon }_1,{\upsilon }_2\right)}\right),\\ {\tau }_{lu}\left({\upsilon }_1,{\upsilon }_2\right)=\frac{hc}{4\pi }{B}_{lu}\frac{N_l\left({\upsilon }_1,{\upsilon }_2\right)}{{\upsilon }_2-{\upsilon }_1},\end{array}$$(2)

where ${\overline{\text{abs}}}_{lu}\left({\upsilon }_1,{\upsilon }_2\right)$ is the mean absorption depth in the radial velocity range [υ₁; υ₂] in the considered line, relative to the comet-free spectrum (see Sect. 2.5), $\overline{\alpha }$[υ₁; υ₂] the mean size of the comet in this range, and N_l(υ₁, υ₂) the column density of the considered species at the lower energy level of the transition with a radial velocity between υ₁ and υ₂. This set of equations is valid for any line of the considered chemical species.

Fig. 3 Cometary absorption spectra obtained as described in Sect. 2.4 for the December 6, 1997 (top) and October 29, 2018 (bottom) observations, around 10 Fe II line around 2750 Å. Low velocities absorptions (between −10 and +10 km s−1) may be attributed to comets or circum-stellar variability.

3.2 Two-parameter model

Equation (2) can be used to fit the mean absorption depth of a comet in a set of lines rising from the same excitation level of a given species, probing its column density N_l. However, in the case where absorption features are detected in various lines rising from several excitation levels, one must find a way to link the various N_l of these levels in order to combine the information they provide and infer the global properties of the transiting comet.

The easiest way of doing this is to consider the case where the excitation levels are sufficiently close in energy to be populated according to their multiplicity, that is, N_l ∝ 𝑔_l. In particular, this hypothesis is valid at local thermodynamic equilibrium (LTE), when, for any pair of transition (l₁, u₁) and (l₂, u₂), we have $$\left|E_{l_1}-E_{l_2}\right|\ll k_{\text{B}}T,$$

where T is the temperature of the absorbing medium, $E_{l_1},E_{l_2}$ the lower level energies of the lines, and k_B the Boltzmann constant.

In these conditions, we thus have τ_lu(υ₁, υ₂) ∝ λ_lu𝑔_lf_lu (Eqs. (1) and (2)). We can therefore introduce a parameter β defined such as $${\tau }_{lu}\left({\upsilon }_1,{\upsilon }_2\right)=\beta \left({\upsilon }_1,{\upsilon }_2\right)\cdot \frac{{\lambda }_{lu}}{{\lambda }_0}{g}_l{f}_{lu},$$

where λ₀ is an arbitrary wavelength. For practicality, in the following λ₀ is chosen to be equal to the wavelength of the 2756 Å Fe II line; however, this choice has no impact on the results as it can be seen as a normalisation unit for β. Indeed, the parameter β is a unitless quantity common to all the studied lines and close to the optical thickness when g_lf_lu ~ 1.β is thus related to the column densities of the excitation levels giving rise to the observed lines through $$\beta \left({\upsilon }_1,{\upsilon }_2\right)g_l=\frac{{\text{e}}^2}{4{ϵ}_0{m}_{\text{e}}c}\frac{N_l\left({\upsilon }_1,{\upsilon }_2\right)}{\left({\upsilon }_2-{\upsilon }_1\right)}{\lambda }_0.$$

Using Eq. (3), we thus obtain the resulting curve of growth equation relating to the mean absorption depth within a radial velocity range to the line strength: $${\overline{\text{abs}}}_{lu}=\overline{\alpha }\cdot \left(1-\exp \left(-\beta \frac{{\lambda }_{lu}}{{\lambda }_0}{g}_l{f}_{lu}\right)\right).$$(3)

We note here that ${\overline{\text{abs}}}_{lu},\overline{\alpha }$ and β depend on the radial velocity range [υ₁; υ₂] where the absorption depths are measured. Thus, to use this model, one must first measure the mean absorption depths of a transiting comet in a set of lines with various absorption strengths (𝑔_lf_lu), in a fixed RV range, and then fit the curve of growth model to these measured absorptions using Eq. (3) in order to estimate the values of α and β that best fit the observations. The column densities N_l of the studied excitation levels can be then derived from the value of β using $$N_l\left({\upsilon }_1,{\upsilon }_2\right)=\frac{4{ϵ}_0{m}_{\text{e}}c}{{\text{e}}^2{\lambda }_0}\left({\upsilon }_2-{\upsilon }_1\right)\beta g_l.$$

This two-parameter model was used by Kiefer et al. (2014b) to determine $\overline{\alpha }$, the typical size of the gaseous clouds, of almost 500 comets transiting β Pic detected in the Ca II doublet. However, in this case, the model was fit to only two absorption measurements (one for each line of the doublet); it was thus not possible to check that it was realistic. To show that the model is correct and that the measured absorptions really follow a curve of growth as given by Eq. (3), a larger number of lines with various oscillator strengths is needed. The many (> 10) Fe II lines observed with STIS, which show clear cometary absorptions, are thus a very useful tool to test the curve of growth model (Sect. 4).

3.3 Thermodynamical equilibrium: A three-parameter model

In the case where the cometary absorptions are detected in lines arising from a large variety of excitation levels, the hypothesis N_l ∝ 𝑔_l may no longer be a good approximation. For instance, very excited levels may be under-populated compared to ground states. To model the relative abundances of the considered energy levels, another possibility is to assume that the absorbing cloud is at LTE, with a typical temperature of T(υ₁, υ₂) (where, again, [υ₁ ; υ₂] is the considered radial velocity range). This implies $$N_l\left({\upsilon }_1,{\upsilon }_2\right)\propto g_l{\text{e}}^{-E_l/k_{\text{B}}T\left({\upsilon }_1,{\upsilon }_2\right)}.$$

We note that this new assumption requires a sufficiently high electronic density to impose a collisional regime. At LTE, Eqs. (1) and (2) now yield ${\tau }_{lu}\left({\upsilon }_1,{\upsilon }_2\right)\propto {\lambda }_{lu}{f}_{lu}{g}_l{\text{e}}^{-E_l/k_{\text{B}}T\left({\upsilon }_1,{\upsilon }_2\right)}$. This allows the introduction of a new parameter γ(υ₁, υ₂), linking the optical thickness to the line strength through $${\tau }_{lu}\left({\upsilon }_1,{\upsilon }_2\right)=\gamma \left({\upsilon }_1,{\upsilon }_2\right)\cdot \frac{{\lambda }_{lu}}{{\lambda }_0}{g}_l{f}_{lu}{\text{e}}^{-E_l/k_{\text{B}}T\left({\upsilon }_1,{\upsilon }_2\right)}.$$

Once again, γ is a unitless parameter common to all studied lines, close to the optical thickness when $g_l{f}_{lu}{\text{e}}^{-E_l/k_{\text{B}}T}\sim 1$, and which depends on the radial velocity range where the mean absorption depths are measured. Thanks to Eqs. (1) and (2), it can be linked to the abundance of the considered excitation levels within [υ₁, υ₂], through: $$\gamma \left({\upsilon }_1,{\upsilon }_2\right)g_l{\text{e}}^{-E_l/k_{\text{B}}T\left({\upsilon }_1,{\upsilon }_2\right)}=\frac{{\text{e}}^2}{4{ϵ}_0{m}_{\text{e}}c}\frac{N_l\left({\upsilon }_1,{\upsilon }_2\right)}{\left({\upsilon }_2-{\upsilon }_1\right)}{\lambda }_0.$$

Summing the previous equation over all the excitation levels of the considered species, one then deduces $$\gamma \left({\upsilon }_1,{\upsilon }_2\right)Z\left(T\left({\upsilon }_1,{\upsilon }_2\right)\right)=\frac{{\text{e}}^2}{4{ϵ}_0{m}_{\text{e}}c}\frac{N_{\text{tot}}\left({\upsilon }_1,{\upsilon }_2\right)}{\left({\upsilon }_2-{\upsilon }_1\right)}{\lambda }_0,$$

where N_tot(υ₁, υ₂) reflects the total column density of the species with radial velocity between υ₁ and υ₂, and Z is defined by $$Z(T)={\displaystyle \underset{i=0}{\overset{+\infty }{{\sum }^{\text{}}}}}{g}_i{\text{e}}^{-E_i/k_{\text{B}}T}.$$

We therefore end up with a three-parameter curve of growth model that links the mean absorption depths of a transiting cloud in a set of lines of a given species, measured on a fixed radial velocity range, to the mean size $\left(\overline{\alpha }\right)$, temperature (T), and total column density (via the γ parameter) of the absorbing cloud, through $${\overline{\text{abs}}}_{lu}=\overline{\alpha }\cdot \left(1-\exp \left(-\gamma \frac{{\lambda }_{lu}}{{\lambda }_0}{g}_l{f}_{lu}{\text{e}}^{-E_l/k_{\text{B}}T}\right)\right),$$(4)

where, again, ${\overline{\text{abs}}}_{lu},\overline{\alpha }$, γ, and T depend on the radial velocity range [υ₁; υ₂] where the mean absorption depths are measured. Once these three quantities have been constrained by a curve of growth model fit to the measured absorption depths, one can estimate the total column density of the studied species through $$N_{\text{tot }}\left({\upsilon }_1,{\upsilon }_2\right)=\frac{4{ϵ}_0{m}_{\text{e}}c}{{\text{e}}^2{\lambda }_0}\left({\upsilon }_2-{\upsilon }_1\right)\gamma Z(T).$$

We note that this model is valid under the assumption of column density and temperature homogeneity of the absorbing gas over the considered radial velocity range. This hypothesis can be checked by fitting the curve of growth model (Eq. (4)) to the cometary absorptions, as done in Sect. 5. If the measured absorptions do indeed follow the model, then the fitted parameters ($\overline{\alpha }$, γ, T) may represent valuable estimates of the typical properties of the transiting cometary cloud.

Finally, the two-parameter model described in Sect. 3.2 can be seen as a simplified case of the three-parameter model, where the temperature T is taken to be much higher than the energy difference between the excitation level of the considered lines. Indeed, in such a case, ${\text{e}}^{-E_l/k_{\text{B}}T}$ is constant over the set of lines, thus yielding N_l(υ₁, υ₂) ∝ 𝑔l (Sect. 3.2). The β parameter introduced in the previous section is then linked to γ through $$\beta \simeq \gamma {\text{e}}^{-E_l/k_{\text{B}}T},$$

where this equation holds for any of the considered energy level l. In particular, this is the case when there is only a single excitation level for the observed set of lines, as for the Ca II doublet at 3933 and 3369 Å.

4 Validation of the two-parameter model

In order to validate our different models, we focus on two separate events in the following. The first one is the transit that occurred on December 6, 1997, between +12 and +40 km s⁻¹ (see Fig. 3). Its main advantage is that it was observed in a great variety of Fe II lines (in particular, in the three main series at 2400, 2600, and 2750 Å), which will be useful for validating the three-parameter model in the next section. The second transit is the one from the October 29, 2018 observation, between +50 and + 150 km s⁻¹ (see Fig. 3), which benefits from a very good S/N (about 100) due to the long duration of the exposure. Another advantage of these transits is their high radial velocities, making them unaffected by the circumstellar absorption near 0 km s⁻¹ (which reduces the S/N and may show long-term variability, Kiefer et al. 2019).

The two-parameter curve of growth model was validated using the Fe II line series around between 2715 and 2769 Å arising from the a⁴D Fe II term. Indeed, the four levels that make up this state share close excitation energies (corresponding to temperature between 11448 and 12771 K; see Tables A.1 and A.3), allowing us to make the hypothesis – which we verify below – that their relative abundances are proportional to their multiplicity.

Figure 4 shows a first example of a curve of growth fit (lower panel) for the redshifted comet observed in the December 6, 1997 data (upper panel). The mean absorption depths were measured in the [+12; +40] km s⁻¹ radial velocity range, corresponding to the main domain where this comet is visible. These absorption measurements were then fit with the two-parameter curve of growth model, using a Markov chain Monte Carlo (MCMC) algorithm. Remarkably, the ten absorption measurements follow the two-parameter curve of growth model very well (reduced χ² of 2.7), constraining the average size of the transiting comet to 44 ± 1% of the stellar surface.

Another example is given in Fig. 5 for the October 29, 2018 event, where the comet mean absorption depths are measured from +60 to +90 km s⁻¹ (the radial velocity range where it is the strongest). For this comet, we also took into account the reddest line of the 2747 Å doublet (see Fig. C.1). Again, the two-parameter model fits the 11 measured absorptions strikingly well (reduced χ² of 1.9), and constrains the mean size of the comet to 39 ± 1% of the stellar disc (within this velocity range).

These two examples clearly show that the variable features observed in the Fe II lines of β Pic are explained well by the transit of homogeneous cometary clouds, covering large but limited fractions of the stellar disc, and in which the sub-levels of the a⁴D Fe II term have relative abundances proportional to their multiplicity. This remarkable agreement between the measured absorptions and the two-parameter curve of growth model also validates the results of Kiefer et al. (2014b), who applied this model to the Ca II doublet at ~3900 Å to identify two families of exocomets in the β Pic system.

In addition, if we assume that the transiting cloud is at LTE (which is probably a reasonable assumption, as shown in Sect. 5) and use the three-parameter model to fit the measured absorption depths (Eq. (4)), we can put lower limits on the gas temperature around 5000 K for the December 6, 1997 comet, and 10 000 K for the October 29, 2018 comet, at a 95% confidence level. This is due to the fact that the lines rising from the most excited levels of the a⁴D term (E_l ~ 12 700 K; see Table A.3) do not show significantly lower absorption than the lines rising from the less excited levels (E_l ~ 11 400 K): the gas temperature must thus be well above the energy difference between all these states, that is, ≫ 1000 K.

Lower limits (with a 2σ confidence) can also be set on the total Fe II column density within the transiting clouds in the considered radial velocity ranges, at about 9.2 × 10¹⁴ cm⁻² (December 6, 1997) and 8.2 × 10¹⁴ cm⁻² (October 29, 2018), which would be reached for temperatures between 11 000 and 15 000 K.

The very high temperatures hinted in the tails of these two comets (similar to or above the effective temperature of β Pic) are consistent with the results that were obtained by Mouillet & Lagrange (1995) using Ca II triplet lines arising from the 3D metastable level. This shows that such objects are experiencing very strong compression, which releases heat (see Beust & Tagger 1993). This compression is also likely responsible for the production of the observed highly ionised species such as Al III or Si IV.

Fig. 4 Analysis of the December 6, 1997 comet. Top: cometary absorption spectrum of December 6, 1997 observation in the a4D Fe II lines around 2750 Å. The grey zone represents the radial velocity range [+12; +40] km s−1 over which the cometary absorptions were measured. Bottom: curve of growth of the same Fe II lines, showing the mean absorption depths measured in the [+12; +40] km s−1 range (blue dots), with one-sigma error bars computed from the STIS pipeline tabulated error bars and the two-parameter fitted model (red line) as a function of the theoretical line strength $\left(g_l{f}_{lu}\frac{{\lambda }_{lu}}{{\lambda }_0}\right)$.

Fig. 5 Analysis of the October 29, 2018 comet. Top: cometary absorption spectrum of October 29, 2018 observation in a4D Fe II lines around 2750 Å. The grey zone represents the radial velocity range [+60; +90] km s−1, which was used to measure the cometary absorptions. Bottom: curve of growth of the same Fe II lines, showing the mean absorption depths measured in the [+60; +90] km s−1 range (blue dots) and the two-parameter fitted model (red line) as a function of the line absorbing strength. The uncertainties on the absorption measurements are of the order of the marker’s size, thanks to the very good S/N (~100) of the October 29, 2018 observation.

5 Validation of the three-parameter model

5.1 Temperature measurement

The analysis in Sect. 4 yielded striking results about the shape and thermodynamical conditions in the transiting cometary clouds of β Pic. However, only lower limits could be set on their temperature and column density, as the energy difference between the different sub-levels of the FeII a⁴D term (which gives rise the line series at 2750 Å; see Tables A.1 and A.3) was too low to probe the energy distribution of Fe⁺. Luckily, the two 1997 observations also cover the 2300–2600 Å wavelength range, where a large number of Fe II lines rising from a wider range of excitation levels are present (see Table A.1). In addition to the ten previously studied lines, we thus included in our analysis 31 new lines with wavelengths between 2328 Å and 2629 Å, and with lower level energies from 0 to 12 700 K (0–8000 cm⁻¹). For the lines arising from low energy levels (≤4500 K), we limited ourselves to the ones with average values of the absorbing strength (0.08 < 𝑔_lf_lu < 0.7), avoiding both the weakest transitions (𝑔_lf_lu < 0.08), for which the S/N on the absorption depth is weak, and the strongest ones (𝑔_lf_lu > 0.7), which showed abnormally strong absorption depths (see the discussion in Sect. 6). However, wavelengths below 2600 Å are only covered by the two 1997 observations. This prevents us from building a proper reference spectrum at these wavelength (as we did in Sect. 2.4 for the 2750 Å lines), as there are wide ranges of radial velocities where both 1997 spectra seem affected by comets (see Table 2) and where, therefore, it is impossible to recover the unnoculted stellar flux. To exploit these lines and test our three-parameter curve of growth model, we thus chose to focus on the redshifted comet observed on December 6, 1997 (Fig. 4) and limited our study to the [+25; +40] km s⁻¹ radial velocity range. Indeed, in this narrower velocity range, the other 1997 observation taken 13 days later seems to be little affected by comets (Fig. 6), providing a sufficiently reliable estimate of the unnoculted β Pic spectrum and enabling us to retrieve the absorption depth of the December 6, 1997 comet in the Fe II lines between 2300 Å and 2600 Å. Here, for the sake of consistency, we also used the same December 19, 1997 observation as a reference spectrum for the 2750 Å lines, despite the fact that other observations are available. This choice, however, has very little incidence on our measurements; the results would be very similar if we had used the ’comet-free’ spectrum as a reference for these lines, as we did in the previous section.

First, we started by separately fitting the mean absorption depths of the December 6, 1997 comet in the lines rising from low-energy levels (E_l ≲ 4500 K; i.e. E_l ~ 3000 cm⁻¹) and higher energy levels (E_l ~ 12 000 K ~ 8000 cm⁻¹), using the two-parameters curve of growth model tested above. The results are provided in Fig. 7; we note that the two sets of lines constrain the size of the absorbing cloud to very close values (38.6 ± 1.8 % for the ground levels lines, 36.5 ± 1.2% for the metastable lines), hinting that the absorptions occurring in both sets are produced by the same gaseous component. We note, however, that the lines rising from ground Fe II levels are more easily saturated than the lines rising from excited levels (β ~ 3.5 for the first, β = 1.4 for the last): this reveals the significant under-population of the excited levels due to a finite temperature in the absorbing medium.

We then fitted the 41 Fe II lines all together, with the addition of a new free parameter: the temperature. Despite the very high number of measurements, three free parameters seem to be enough to fit the data all together. Indeed, our curve of growth model matches the measured absorption depths remarkably well $\left({\chi }_r^2=1.7\right)$, allowing us to constrain the temperature of the absorbing medium at 10 500 ± 500 K and its average size at 36 ± 1% of the stellar surface (Fig. 8). The total Fe II column density in the studied radial velocity range (i.e. between +25 and +40 km s⁻¹) was thus estimated to be 6.7 ± 0.4 × 10¹⁴cm⁻².

In order to better visualise the temperature of the absorbing gas, it is also possible to build an excitation diagram of Fe⁺ within the comet. Indeed, the model Eq. (4) can be re-written as $${\text{e}}^{-E_l/k_{\text{B}}T}=\frac{-\ln \left(1-\left({\overline{\text{abs}}}_{lu}/\overline{\alpha }\right)\right)}{\gamma \left({\lambda }_{lu}/{\lambda }_0\right)g_l{f}_{lu}},$$

yielding $$\frac{-E_l}{k_{\text{B}}T\ln (10)}=\log \left(\frac{N_l}{g_l{N}_{\text{tot }}}Z(T)\right)=\log \left(\frac{-\ln \left(1-\left({\overline{abs}}_{lu}/\overline{\alpha }\right)\right)}{\gamma \left({\lambda }_{lu}/{\lambda }_0\right)g_l{f}_{lu}}\right).$$(5)

The slope of log(N_l/𝑔_l) expressed as a function of E_l is therefore inversely proportional to the temperature T. On the other hand, knowing $\overline{\alpha }$, which characterises the typical size of the absorbing cloud independently of the line considered, we can derive an estimate of log(N_l/𝑔_l) from each absorption measurement abs_lu in any spectral line arising from the energy level to within one additive constant depending on γ, Z(T), and N_tot. Therefore, the temperature of the transiting cloud can be directly estimated through the variation of N_l/𝑔_l (the relative abundance of the various Fe⁺ excitation levels, normalised by their multiplicity) as a function of E_l. This is illustrated in Fig. 9, where we can clearly see that the a⁴D energy levels (E_l ~ 12 000 K) are significantly under-populated when compared to the ground levels due to the finite temperature (~ 10 500 K) within the transiting gas.

We can finally extrapolate this temperature to the whole velocity range of the December 6, 1997 comet (i.e. between +12 and +40 km s⁻¹) in order to estimate its total column density (which was found in Sect. 4 to be greater than 9.2 × 10¹⁴ cm⁻²). Indeed, fitting the 2750 Å lines’ absorption depths measured in this velocity range again with a Gaussian prior on the temperature of 10 500 ± 500 K, one can deduce a total Fe II column density in the transiting gas of 1.11 ± 0.09 × 10¹⁵ cm⁻². Given the relative size of the comet (~44% of β Pic surface, Fig. 4), the radius of β Pic (~1.8 R_⊙) and the mass of an iron atom (~9.4 × 10⁻²³ g), we can convert this column density into a mass of ~ 2.2 × 10¹⁵ g, which is equivalent to 0.28 km³ of pure iron (density ≃8). This tremendous mass shows the intense evaporation occurring at the surface of these cometary nuclei, hinting that they have very short lifetimes and may not survive to even a single periastron passage. This remark therefore indicates that a large fraction of the comets observed to transit β Pic may originate from a recent break-up of a former massive body, as was suggested in Kiefer et al. (2014b) for the D family of β Pic comets.

Fig. 6 β Pic spectra around two Fe II lines. Top: zoomed-in view of the 2740 Å Fe II line (El = 7955 cm−1), emphasising the two 1997 observations. The December 19 spectrum is close to the stellar spectrum in the [+25; +40] km s−1 range (grey zone), while the December 6 observation shows a clear cometary transit at these radial velocities. Bottom: zoomed-in view of the 2607 Å Fe II line (El = 668 cm−1), showing the only two available spectra for this line.

5.2 The electronic density

Until now, our study was focussed on Fe II lines arising from energy levels below 12 000 K, assuming that the transiting cloud is at LTE. Using the ChiantiPy Python module (Landi et al. 2012), which allows us to compute the energy distribution of many ions at given temperature and electronic density (n_e), it can be shown that this hypothesis is easily verified, roughly requiring n_e ≥ 10⁵ cm⁻³. This allowed for an estimate of the absorbing gas temperature, assuming that the excitation temperature of the low-energy Fe II levels (deduced from the excitation diagram, Fig. 9) is equal to the kinetic temperature.

In addition to the lines studied in Sects. 4 and 5, cometary absorption is also detected in a few lines arising from highly excited Fe II levels (from 30 000 to 40 000 K; see Table A.2 and Fig. C.3) at wavelengths between 2400 and 2900 Å. A key feature of these lines is that their lower levels require much higher electronic densities to be populated at LTE, around n_e ≳ 10⁸ cm⁻³ (ChiantiPy; Landi et al. 2012); at low densities, they tend to be less populated than at LTE. Thus, this set of lines can be used as a probe of the electronic density within the transiting cometary tails. To do so, we started by measuring the mean absorption depths of the December 6, 1997 comet in the [+25; +40] km s⁻¹ range for 26 high-energy Fe II lines. Then, given the cometary cloud properties measured in Sect. 5.1 (T ~ 10 500 K, $\overline{\alpha }\sim 0.356$, γ ~ 4.72), we used these absorption measurements to extract the relative abundances of the energy levels corresponding to those lines (Eq. (5)). These abundance estimates rely on the fact that the second equality in Eq. (5) remains valid for energy levels that are not populated at LTE, as long as most of the absorbing gas follows a Boltzmann distribution: $$\frac{N_l}{g_l{N}_{\text{Fe\hspace{0.17em}\hspace{0.17em}II}}}=\frac{-\ln \left(1-{\overline{\text{abs}}}_{lu}/\overline{\alpha }\right)}{\gamma \frac{{\lambda }_{lu}}{{\lambda }_0}{g}_l{f}_{lu}Z(T)}.$$

Then, we compared these relative abundances to the values predicted by LTE at the measured temperature (10 500 K) by building the complete excitation diagram of Fe II within the cometary tail (Fig. 10). This diagram clearly shows that high-energy (E_l ≳ 30 000 K) Fe II levels are under-populated compared to LTE, hinting that the electronic density is not high enough to impose a Boltzmann distribution.

Finally, in order to quantify the typical electronic density within the cometary cloud, we compared the energy distribution of Fe II measured in our spectra to the theoretical energy distribution for various electronic densities, obtained with ChiantiPy (to compute the distribution, a temperature of 10 500 K was used, as measured in Sect. 5.1). Although it is difficult to give a precise estimate, the excitation diagram of Fe II given in Fig. 10 seems to indicate an electronic density in the [2–4] × 10⁷ cm⁻³ range. These results are consistent with the work of Mouillet & Lagrange (1995), where a minimum electronic density of 10⁶ cm⁻³ in a transiting comet was inferred from its absorption in the Ca II triplet.

Fig. 7 Curve of growth of 26 Fe II lines rising from low energy levels (0–4500 K, green) and 15 lines rising from high energy levels (~12 000 K. blue) for the December 6, 1997 comet in the [+25, +40] km s−1 RV range. While both sets of lines provide very similar values for the average size of the absorbing cloud (~38%), the lines rising from the ground terms saturate a lot quicker.

Fig. 8 Curve of growth of 41 Fe II lines rising from various energy levels (from 0 to 12 700 K). The cometary absorption depths (blue dots) were measured in the [+25; +40] km s−1 velocity range in the December 6, 1997 observations.

Fig. 9 Excitation diagram of Fe II within the December 6, 1997 comet, obtained using absorption measurements in the [+25; +40] km s−1 radial velocity range. The y-axis represents the log-relative abundance of the studied Fe II energy levels, within an additive constant (see Eq. (5)), while the x-axis indicates the energy of these levels. Each point thus represents an abundance estimate of an energy level l, based on the measurement of the absorption depth in a given line (l, u). As some lines rise from the same energy level, some of these levels benefit from several independent abundance measurements. The marker sizes are inversely proportional to the vertical uncertainties.

Fig. 10 Full excitation diagram of Fe II within December 6, 1997 comet showing measured column-density of various energy levels inferred from the measured absorptions and the parameters obtained in Sect. 5.1 (blue dots), the theoretical energy distribution at LTE (red line), and the theoretical energy distribution with an electron density of 3 × 107 cm−3 for a few representative energy levels as computed using ChiantiPy (purple dots).

6 Discussion

6.1 Deviation from the curve of growth model

The analysis of the exocomet transit signature in the December 6, 1997 observation shows that the transiting gas in the [+12; +40] km s⁻¹ range is described well by an homogeneous cloud covering, on average, 44% of the stellar surface, with a total Fe II column density of ~1.1 × 10¹⁵ cm⁻², a temperature of ~10 500 K (Sect. 5.1), and an electronic density of approximately 10⁷ cm⁻³ (Sect. 5.2). These results were obtained by using most of the well-detected Fe II lines in the 2300–2800 Å wavelength domain. However, we note that for a few very strong lines (𝑔l f_lu > 0.7) rising from low energy levels (E_l ≤ 4500 K), the measured absorption depth seems to deviate from the curve of growth model. For instance, the mean absorption depth of the December 6, 1997 comet in the [+25; +40] km s⁻¹ range exceeds 50% for the lines at 2383 Å (𝑔l f_lu = 3.20), 2396 Å (𝑔l f_lu = 2.23) and 2600 Å (𝑔l f_lu = 2.63, see Fig. 11), which is well above the measured size of the cometary cloud (~36%) in this radial velocity domain.

This discrepancy between the curve of growth model fitted in Sect. 5.1, which matches the measured absorption depths on the “weak” Fe II lines (0.08 < 𝑔l f_ul < 0.7) strikingly well, and the absorption measured in the strongest ones (𝑔l f_ul > 0.7) can be revealed by plotting the curve of growth model with all the absorption measurements (Fig. 12). We see that while the cometary absorption depths in the weakest Fe II lines (blue dots) saturate towards 36%, the absorption in strong, ground-level Fe II lines (green dots) keeps increasing slowly as the line strength (≃𝑔lf_lu) increases. This shows that the complexity of the transiting object is probably underestimated by our model. For instance, the absorption measurements might be better explained by the transit of at least two gaseous components: a dense one, which quickly saturates towards a 36% absorption depth as the line strength increases, and a thinner one, saturating more slowly, explaining why the absorption increases again by about 20% when we reach the strongest lines.

Although our simple model allows for a fairly accurate characterisation of transiting comets, it thus fails to probe the complex structure that such objects probably have. Nevertheless, our study shows the tremendous potential of the many near ultraviolet Fe II lines for investigating the properties of the transiting gaseous tails. This paves the way for the development of more sophisticated models, in order to interpret the variable features observed in the spectrum of several stars in even greater detail. Such models may allow us to probe the temperature and density distribution of the transiting clouds – which are probably not perfectly homogeneous, as revealed by the study of the December 6, 1997 comet – as well as the orbital dynamics of the evaporating nucleus. In particular, one key objective of forthcoming works will be to infer the radial velocity of the nucleus from the absorption profile of its tail.

These models will also have to take into account the inhomogeneous distribution of the stellar flux at a given wavelength on the stellar disc. Indeed, this distribution is affected by limb darkening, which reduced the luminosity of the edges of the stellar disc, and by the rapid rotation velocity of stars such as β Pic (~130 km s⁻¹), which shifts the local stellar lines. Strictly speaking, occulting a fraction α of the stellar disc is thus not equivalent to reducing the stellar flux by a fraction α at each wavelength. This subtlety, which was neglected in our study, might actually be useful to constrain the 2D distribution of transiting cometary tails, by probing the location of the stellar regions they occult.

Fig. 11 Fe II line at 2600.17 Å, showing a very deep cometary absorption between +25 and +40 km s−1 in the December 6, 1997 observations. The average absorption depth is about 55 %, which is well above the size of the comet estimated using fainter lines from higher energy levels (~36 %). The deep absorption seen at −10 km s−1 relative to the β Pic radial velocity is due to the interstellar medium, which is not seen in the other Fe II lines arising from excited levels.

Fig. 12 Comparison between the curve of growth model fitted in Sect. 5 to the absorption of December 6, 1997 comet in weak Fe II lines, and the measured absorption depths of this comet in the strongest ground-level Fe II lines (green).

6.2 Optical thickness profile

The curve of growth model developed in Sect. 3 to derive the physical properties of β Pic exocomets uses the hypothesis that the optical thickness – characterised by the β parameter – is constant over the radial velocity range where the cometary absorption depths are measured. In fact, it can be shown that, even if the optical thickness profile is not flat, the model still provides very good estimates of the size, temperature, and column density of the absorbing cloud. To check this, we generated artificial absorption spectra of transiting comets with known column densities, and with quickly varying optical thickness over radial velocity. We found that, as long as the optical thickness stays below unity, or does not vary by a large factor (≤ 10) over the studied radial velocity range, the curve of growth model provides very good estimates of the true physical parameters of the absorbing gas (N_tot, α, T).

To go further, we could also subdivide the radial velocity domain [υ₁, υ₂] into elementary intervals and apply the curve of growth model to the absorption depths measured within every range [υ, υ + dυ]. This would allow us to probe the radial velocity dependency of the size (α(υ)) and optical thickness (γ(υ)) of the gaseous cometary cloud. Such an analysis could thus give a first insight into the geometry of a transiting gaseous tail, probing, in particular, how its size and optical thickness evolve as the gas is pushed away from the star by the radiation pressure. However, this analysis is far beyond the scope of the present work.

7 Conclusion

Using HST archival UV spectra of β Pic, we detected variable absorptions in dozens of Fe II lines, which are due to the transit of gaseous cometary tails. To better characterise the properties of these objects, we developed a simple model that describes the absorbing gas with a homogeneous cloud, covering a fraction α of the stellar disc. This model matches the measured absorptions in most of the Fe II lines remarkably well – with the exception of the strongest ones – showing that the variable features seen in the β Pic spectrum are well interpreted by the transit of large $(\overline{\alpha }\sim 0.4)$, dense (N_{Fe II} ~ 10¹⁵cm⁻²) and hot (T ≥ 5000 K) cometary gas tails. In the comet that transited in December 6, 1997, the population of the energy levels is consistent with LTE, which allows an estimate of the temperature within the gas tail of 10 500 ± 500 K. Using Fe II lines arising from highly excited electronic levels, we estimated the typical electronic density to be around a few 10⁷ cm⁻³ in this cometary gas cloud. All the results obtained here are consistent with previous observational studies, notably with the quantitative estimates given by Lagrange et al. (1995), Mouillet & Lagrange (1995), and Kiefer et al. (2014b). They are, however, slightly different to the theoretical predictions of Beust & Tagger (1993), for which the formation of a shock front with very high temperature (T ~ 10⁵ K) and rather low density (n_e ~ 10³−10⁶) is necessary to explain the presence of highly ionised species (such as Al III) in the transiting cometary tails. Further work will be needed to better understand this discrepancy.

Our results also confirm the huge difference between β Pic and solar comets; while the iron in β Pic comets is mainly present as Fe⁺ (only very shallow cometary absorption is detected in Fe I lines; see Welsh & Montgomery 2016), with column densities of the order of 10¹⁵cm⁻², solar comets contain mainly neutral iron, with much lower column densities (~10⁸−10¹⁰ cm⁻², Manfroid et al. 2021). This is likely to be a consequence of the very short distance between β Pic transiting comets and their star, which induces an intense sublimation of metal-rich material at the surface of the nuclei, while the tails of solar comets mainly result from the sublimation of volatile ices.

Finally, the curve of growth model appears to be an efficient tool to characterise the gaseous clouds surrounding transiting exocomets, detected through variable spectroscopic absorptions. The exocomet curve of growth is the analogue of the classical curve of growth used for interstellar clouds and stellar atmospheres. A fit to this curve provides quantitative estimates of the typical size of the absorbing cloud, and the column densities of the detected species. When several lines from various excitation energy levels are available, the temperature and possibly the electronic density can also be constrained. Further works using this new tool are already in progress.

Acknowledgements

T.V., A.L., G.H. and A.V.-M. acknowledge funding from the Centre National d’Études Spatiales (CNES).

Appendix A Line and energy level parameters

Appendix B Full studied spectra

Fig. B.1 Normalised spectra of β Pic at seven different epochs, focussing on three Fe II line series at 2400, 2600, and 1750 Å. Vertical blue lines emphasise the Fe II transitions that were used in our analysis and are listed in Table A.1. Other lines not used in this paper may also show cometary features, such as Fe II multiplets (e.g. 2344, 2750 Å) or lines from other species (e.g. Mn II at 2576, 2594 and 2606 Å).

Appendix C Line examples

Fig. C.1 Examples of β Pic spectra in Fe II lines near 2750 Å, rising from a4D metastable term. The lines are at 2737 Å (top left, El = 8680 cm−1), 2740 Å (top right, El = 7955 cm−1), 2747.8 Å (bottom left, El = 8392 cm−1) and 2756 Å (bottom right, El = 7955 cm−1). The blue lines emphasise the location of the lines in the β Pic rest frame (thus at 0 km/s).

Fig. C.2 Examples of Fell lines from low energy levels in 2300-2600 Å region: 2365Å (top-left, El = 385cm−1), 2369 Å (top right, El = 4080 cm−1), 2607 Å (bottom left, El = 668cm−1), and 2626 Å (bottom right, El = 385cm−1). Only two observations obtained in December 1997 are available on this wavelength range. Again, the blue lines emphasise the location of the lines in the stellar rest frame.

Fig. C.3 Examples of faint Fe II lines rising from highly excited levels: 2424Å (top left, El = 22637cm−1), 2512Å (top right, El = 21712 cm−1), 2754 Å (bottom left, El = 26352cm−1), and 2768 Å (bottom right, El = 26170 cm−1), with blue lines indicating the location of the lines in β Pic rest frame.

References

All Tables

All Figures

Fig. 1 Illustration of the recalibration method. Top panel: raw spectra show significant shape disparities, probably resulting from changes in instrument calibration throughout the years. These calibration disparities are hardly fitted by polynomials. Centre panel: zoomed-in view of 2680–2820 Å region, showing the splines fitted to each spectrum. The grey areas indicate the different spectral domains where the reference fluxes were measured (black crosses), which are then used to calculate a cubic spline for each spectrum. Bottom panel: zoomed-in view of 2735–2753 Å region after renormalisation by the cubic splines. The spectra taken at different epochs are shown with different colours; they are well superimposed in the region of the stellar continuum, while variable absorptions are clearly visible around the different Fe II lines between 2737 and 2750 Å.

In the text

Fig. 2 Variable absorptions in the 2740 Å Fe II line. The black line shows the reference comet-free spectrum of β Pic, obtained by averaging the unabsorbed spectra at each individual wavelength. The deep and narrow absorption near 0 km s−1 is due to the stable component of the circumstellar disc of β Pic, which could be slightly blueshifted compared to the stellar radial velocity (see Brandeker 2011).

In the text

	Fig. 3 Cometary absorption spectra obtained as described in Sect. 2.4 for the December 6, 1997 (top) and October 29, 2018 (bottom) observations, around 10 Fe II line around 2750 Å. Low velocities absorptions (between −10 and +10 km s−1) may be attributed to comets or circum-stellar variability.
In the text

Fig. 4 Analysis of the December 6, 1997 comet. Top: cometary absorption spectrum of December 6, 1997 observation in the a4D Fe II lines around 2750 Å. The grey zone represents the radial velocity range [+12; +40] km s−1 over which the cometary absorptions were measured. Bottom: curve of growth of the same Fe II lines, showing the mean absorption depths measured in the [+12; +40] km s−1 range (blue dots), with one-sigma error bars computed from the STIS pipeline tabulated error bars and the two-parameter fitted model (red line) as a function of the theoretical line strength $\left(g_l{f}_{lu}\frac{{\lambda }_{lu}}{{\lambda }_0}\right)$.

In the text

Fig. 5 Analysis of the October 29, 2018 comet. Top: cometary absorption spectrum of October 29, 2018 observation in a4D Fe II lines around 2750 Å. The grey zone represents the radial velocity range [+60; +90] km s−1, which was used to measure the cometary absorptions. Bottom: curve of growth of the same Fe II lines, showing the mean absorption depths measured in the [+60; +90] km s−1 range (blue dots) and the two-parameter fitted model (red line) as a function of the line absorbing strength. The uncertainties on the absorption measurements are of the order of the marker’s size, thanks to the very good S/N (~100) of the October 29, 2018 observation.

In the text

Fig. 6 β Pic spectra around two Fe II lines. Top: zoomed-in view of the 2740 Å Fe II line (El = 7955 cm−1), emphasising the two 1997 observations. The December 19 spectrum is close to the stellar spectrum in the [+25; +40] km s−1 range (grey zone), while the December 6 observation shows a clear cometary transit at these radial velocities. Bottom: zoomed-in view of the 2607 Å Fe II line (El = 668 cm−1), showing the only two available spectra for this line.

In the text

	Fig. 7 Curve of growth of 26 Fe II lines rising from low energy levels (0–4500 K, green) and 15 lines rising from high energy levels (~12 000 K. blue) for the December 6, 1997 comet in the [+25, +40] km s−1 RV range. While both sets of lines provide very similar values for the average size of the absorbing cloud (~38%), the lines rising from the ground terms saturate a lot quicker.
In the text

	Fig. 8 Curve of growth of 41 Fe II lines rising from various energy levels (from 0 to 12 700 K). The cometary absorption depths (blue dots) were measured in the [+25; +40] km s−1 velocity range in the December 6, 1997 observations.
In the text

Fig. 9 Excitation diagram of Fe II within the December 6, 1997 comet, obtained using absorption measurements in the [+25; +40] km s−1 radial velocity range. The y-axis represents the log-relative abundance of the studied Fe II energy levels, within an additive constant (see Eq. (5)), while the x-axis indicates the energy of these levels. Each point thus represents an abundance estimate of an energy level l, based on the measurement of the absorption depth in a given line (l, u). As some lines rise from the same energy level, some of these levels benefit from several independent abundance measurements. The marker sizes are inversely proportional to the vertical uncertainties.

In the text

Fig. 10 Full excitation diagram of Fe II within December 6, 1997 comet showing measured column-density of various energy levels inferred from the measured absorptions and the parameters obtained in Sect. 5.1 (blue dots), the theoretical energy distribution at LTE (red line), and the theoretical energy distribution with an electron density of 3 × 107 cm−3 for a few representative energy levels as computed using ChiantiPy (purple dots).

In the text

Fig. 11 Fe II line at 2600.17 Å, showing a very deep cometary absorption between +25 and +40 km s−1 in the December 6, 1997 observations. The average absorption depth is about 55 %, which is well above the size of the comet estimated using fainter lines from higher energy levels (~36 %). The deep absorption seen at −10 km s−1 relative to the β Pic radial velocity is due to the interstellar medium, which is not seen in the other Fe II lines arising from excited levels.

In the text

	Fig. 12 Comparison between the curve of growth model fitted in Sect. 5 to the absorption of December 6, 1997 comet in weak Fe II lines, and the measured absorption depths of this comet in the strongest ground-level Fe II lines (green).
In the text

Fig. B.1 Normalised spectra of β Pic at seven different epochs, focussing on three Fe II line series at 2400, 2600, and 1750 Å. Vertical blue lines emphasise the Fe II transitions that were used in our analysis and are listed in Table A.1. Other lines not used in this paper may also show cometary features, such as Fe II multiplets (e.g. 2344, 2750 Å) or lines from other species (e.g. Mn II at 2576, 2594 and 2606 Å).

In the text

	Fig. C.1 Examples of β Pic spectra in Fe II lines near 2750 Å, rising from a4D metastable term. The lines are at 2737 Å (top left, El = 8680 cm−1), 2740 Å (top right, El = 7955 cm−1), 2747.8 Å (bottom left, El = 8392 cm−1) and 2756 Å (bottom right, El = 7955 cm−1). The blue lines emphasise the location of the lines in the β Pic rest frame (thus at 0 km/s).
In the text

	Fig. C.2 Examples of Fell lines from low energy levels in 2300-2600 Å region: 2365Å (top-left, El = 385cm−1), 2369 Å (top right, El = 4080 cm−1), 2607 Å (bottom left, El = 668cm−1), and 2626 Å (bottom right, El = 385cm−1). Only two observations obtained in December 1997 are available on this wavelength range. Again, the blue lines emphasise the location of the lines in the stellar rest frame.
In the text

	Fig. C.3 Examples of faint Fe II lines rising from highly excited levels: 2424Å (top left, El = 22637cm−1), 2512Å (top right, El = 21712 cm−1), 2754 Å (bottom left, El = 26352cm−1), and 2768 Å (bottom right, El = 26170 cm−1), with blue lines indicating the location of the lines in β Pic rest frame.
In the text