© The Authors 2025

1. Introduction

Fast radio bursts (FRBs) are extremely bright radio pulses that typically occur within milliseconds. Their dispersion measure (DM) – the spectro-temporal pattern imprinted in the signal by the interstellar medium – can be used to infer a distant extragalactic origin, which has been confirmed by the localisation of some host galaxies (e.g. Petroff et al. 2022).

Two broad categories stand out (e.g. Pleunis et al. 2021): first, the one-off events that have been seen only a single time; second, the so-called repeaters that have been seen to repeat in a random, albeit clustered, fashion (e.g. Cruces et al. 2021). Two repeaters have been shown to follow a periodic activity window (Rajwade et al. 2020; The CHIME/FRB Collaboration 2020).

Spectro-temporal burst morphologies appear to be statistically different between the two categories (e.g. Pleunis et al. 2021). One-offs appear to be shorter and broad-band, and usually made of a single component, or sub-burst. Repeaters, on the other hand, are narrow-band events that tend to last longer, with successive components drifting downwards in frequency (Hessels et al. 2019). This difference in morphology leads one to wonder if indeed there might be two different types of progenitors.

Many theories have been proposed to explain FRBs, an increasing number of which are focusing on neutron stars, and particularly variants of magnetars (e.g. Zhang 2020; Voisin 2021). Magnetars are appealing because of their ample magnetic energy reservoir, and their known randomly eruptive behaviour, which appears to be statistically consistent with FRBs (Wadiasingh & Timokhin 2019) although it occurs primarily in the X-ray and gamma-ray bands with occasional radio counterparts. These radio counterparts are, however, much weaker than known FRBs and, conversely, high-energy counterparts have so far not been detected at extragalactic distances (Zhang 2024). Thus, if FRBs are caused by magnetars, then these are likely unusual by Galactic standards. One putative FRB has been seen from the Galactic magnetar SGR1935+2154 (Bochenek et al. 2020; Andersen et al. 2020), with a burst both a thousand times brighter than other known magnetar emissions, and a hundred times weaker than extragalactic FRBs. A range of models involve magnetars, including FRBs, which are caused by very young magnetars (Metzger et al. 2017) or low-twist ones (Wadiasingh & Timokhin 2019).

Inferring source properties from observed FRBs proves to be very challenging. So far, to the best of our knowledge, events have been fitted independently using empirical functions (e.g. Hessels et al. 2019; Fonseca et al. 2024). This allows one to extract burst-specific quantities such as the bandwidth, duration, or frequency drift rate.

In Voisin (2023), we proposed a geometrical model that can constrain the spectro-temporal properties of an event as a function of its time of arrival and source parameters. The model assumes that emission is due to a localised ultra-relativistic plasma flowing along curved streamlines. Tangential emission along magnetic field lines is notably at the core of the rotating vector model (Radhakrishnan & Cooke 1969; Pétri 2017) successfully applied to explain the polarisation swing of many pulsars. Similar swings have been observed in a couple of FRBs (Luo et al. 2020; Mckinven et al. 2025). This suggests that this type of geometrical model is also relevant for FRBs, albeit with possibly different magnetic configurations or lines of sight.

In our model the emission mechanism need not be specified, and it is enough to say that emission is beamed forwards as a result of relativistic beaming. The model shows that the spectro-temporal morphology of a burst is constrained by a geometrical envelope that depends only on the local geometry of the streamlines, the spin period of the star, and the degree of beaming. Thus, the morphology of events from a given source are constrained by its spin and global magnetic geometry. Said differently, global source properties can be inferred from burstmorphologies.

In this work, we assume that FRBs are produced by plasma flowing along magnetic field lines inside the co-rotating region of a magnetar magnetosphere. Such magnetospheres are thought to be twisted; that is, to possess a toroidal magnetic component. This can result from magnetic footpoint motion during eruptions (Thompson et al. 2002), or from the formation and structure of the star (e.g. Barrère et al. 2022; Uryu et al. 2023). We used a very simplified phenomenological model of the magnetic structure, considering it as the sum of a dipole aligned with the spin axis, and of a uniform toroidal component decaying radially as a power law (see Sect. 2).

Using our model, we performed a fit of the data of the two most frequent repeaters in the first CHIME/FRB catalogue (The CHIME/FRB Collaboration 2021): FRB20180814A and FRB20180916B, which in the following we call A and B, respectively. B is famous for being one of the two repeaters exhibiting a periodic activity window (The CHIME/FRB Collaboration 2020).

All the events of a given source are fitted jointly, since the model predicts common morphological constraints depending mostly on spin and magnetic geometry. This is also why repeaters are more relevant than one-off events for the practical tests of our model. To our knowledge, it is the first time such a joint fit has been performed. Although the number of bursts is relatively low, the fit has the advantage of keeping computational needs manageable for this exploratory work.

2. Methods

2.1. Data

The data for the two sources with the largest number of repetitions in the CHIME/FRB catalogue (The CHIME/FRB Collaboration 2021), namely FRB 20180814A (A) and FRB 20180916B (B), was retrieved from the dedicated website¹. The number of events is, respectively, 11 for A and 19 for B. The catalogue itself was retrieved thanks to the python package cfod².

The data spans from August 2018 to June 2019 for A, and from September 2018 to June 2019 for B. In total, the public data in the catalogue contains 11 events for A and 19 events for B as determined by the CHIME/FRB pipeline. Events correspond to one or several components forming a continuous signal, with the notable exception of FRB20181019A from source B for which the two components appear to be fully separated by about 60 ms. In the case of A, seven out of eleven events are composed of a single component (also called sub-burst below) and the four others have multiple components. For B, there are eleven singles and eight multiples. Both sources have events across the entire bandwidth of 400−800 MHz, with durations shorter than 100 ms. In the catalogue, each event is de-dispersed individually, with DMs of ∼190 pc cm⁻³ for A and ∼349 pc cm⁻³ for B (see also Sect. 2.4).

Basic pre-processing such as interference removal and generation of dynamic spectra (waterfall plots) was adapted from the documentation kindly provided by the CHIME/FRB collaboration³. The original time resolution of 0.983 ms was kept, while the frequency resolution was reduced to 128 bins corresponding to bin widths of 3.125 MHz in order to reduce the computational load.

For each event, quantities provided in the catalogue, such as the burst arrival time, width, and fluxes, were retrieved in order to be used as first guesses for the fitting procedure. Gaussian profiles were then fitted to the intensity profile, requiring manual interventions in a select number of cases in which many components were needed. These fits finally provided initial estimates for the peak times, Gaussian width, and fluence of each individual component within events.

2.2. Model

The model presented in Voisin (2023) constrains burst spectro-temporal properties using the geometry of streamlines followed by the emitting plasma. Its key assumptions are (i) the relativistic nature of the emitting plasma and (ii) the locality of the emission region. The main consequence of assumption (i) is that the emission is relativistically beamed forwards in the direction of motion. Assumption (ii) ensures that streamlines within the emission region can be approximated to a single streamline endowed with an effective emission angle. Thus, this emission angle accounts for the local field line divergence as well as the intrinsic beaming of the emission process. Locally, streamlines are approximated by their tangent direction and curvature. The fact that a line is curved together with the narrowness of the emission cone ensures that only a small segment of the line is visible by an observer at any time. The small lateral extent of the emission region effectively means that only a small bundle of lines radiates. In addition, the model accounts for streamlines defined in a rotating frame with period P_*.

In this work, we assume that streamlines follow the magnetic field lines of a co-rotating neutron-star magnetosphere. This implies that emission regions are located sufficiently close to the star in order to avoid relativistic co-rotation speeds that are not accounted for by the model as yet. This condition was not imposed a priori, but we checked a posteriori that this is indeed the case that r ≪ cP_*/2π, where r is the distance to the stellar centre, c is the speed of light, and the right-hand side is the so-called light-cylinder radius, which is the distance at which the co-rotation speed equals the speed of light.

The emission frequency, f, was mapped to geometry through the radius-to-frequency mapping phenomenological law (e.g. Lyutikov 2020), f = f_*R_*/r, where R_* is the stellar radius and f_* is the emission frequency at the stellar surface. Thanks to this mapping, the emission region maps directly onto an area in the dynamic spectrum of the observer. Indeed it is clear that, as the star rotates, the observer moves out of the emission beam, and thus defines a finite visible segment of time, while the radial extent of the emission region maps to a frequency interval (see Fig. 1). The radial extent of the region within the line of sight of the observer changes in time, which maps to a non-trivial spectro-temporal envelope⁴. The centre of each spectro-temporal envelope (T₀, f₀) was mapped to a point in the magnetosphere. This point was numerically computed as the point at which the tangent to the field line points towards the observer at time T₀ at the altitude corresponding to f₀.

Fig. 1. Illustration of the model, from Voisin (2023). The emitter travels from left to right, emitting in a forward cone (the cone is illustrative, in practice we use a smooth angular profile in this work). Emission frequency varies by Δf as it travels radially by Δr. Dotted lines show the observer’s direction, with γ the angle with respect to the path. The two thick ticks delimit the visible segment beyond which no emission can be seen from within the cone. The characteristic emission frequency, fc, maps to the radius, rc, at the centre of the segment. The visible segment being in a rotating frame, it varies with time until it vanishes.

The magnetic field was constructed in the simplest possible way: a dipolar field, B_d, superposed with a toroidal field (relative to the rotation axis) defined by B_t = −αB_d*(R_*/r)^pe_ϕ, where B_d* is the dipolar field strength at the pole, and e_ϕ is the azimuthal unit vector. Magnetic geometry is thus parametrised by the toroidal-to-dipolar ratio, α, and by the exponent p. The ratio, α, can be either positive or negative. The strength, B_d*, is irrelevant for geometrical purposes and can be set to 1. In order to simplify further, the dipole is assumed to be aligned with the spin axis of the star, making the geometry cylindrically symmetric. It follows that the orientation of the star with respect to the observer is only determined by the inclination, i, of its axis with respect to the normal to the plane of the sky (oriented away from the observer), and the rotational phase is irrelevant in this particular case. If the toroidal field is important, or dominant, axial symmetry is partially realised, which makes it a more reasonable approximation.

We also tried a model where f = f_*R_*^p/r^p, such that the radius-to-frequency mapping follows the same scaling as the toroidal component of the magnetic field (which, in our model, is expected to dominate in the inner magnetosphere). However we were unable to draw conclusions from this model, as the MCMC could not converge and was mostly constrained by priors on magnetic parameters and frequency parameters. Thus we focused only on the inverse law mentioned above.

So far, we have reasoned with a finite emission cone emerging from a radiating moving point. In practice, bursts were modelled using a pseudo-Gaussian model (Voisin 2023) whereby a radiating blob (it could be a wave packet or an actual blob of matter) is injected along a field line following a Gaussian time profile in intensity with a characteristic width, w, fluence, F, and injection time, t_i. A blob radiates across a Gaussian angular profile of characteristic opening angle Ω. This angle is assumed to be common to all events of a given source, although this is probably not the case in general.

Under the pseudo-Gaussian model (with f ∝ r⁻¹), an event is the sum of sub-bursts of the form (see Eq. (63) of Voisin 2023)

$$\begin{array}{c}\hfill I(t_a,f)\propto \underset{\mathcal{E}(t_a,f)}{\underbrace{\frac{1}{f^2{\mathrm{\Omega }}^{2}cos\rho }exp(-\frac{{\gamma }^2(t_a,f)}{2{\mathrm{\Omega }}^2})}}Fexp(-\frac{{(t_a-t_i)}^2}{2w^2}),\end{array}$$(1)

where I is the intensity per unit time of arrival, t_a, per unit frequency, f. The function γ characterises the emission beam, and is defined as the angle between the line of sight and the tangent to the magnetic field at the magnetospheric location that produces the emission received at (t_a, f) in the dynamic spectrum (see Fig. 1). This function depends on the global geometry, orientation, and spin; that is, it depends on the global parameters. It also depends on the location of the centre of the emission region, x₀, which is given by the co-ordinates (T₀, f₀), it being understood that the model maps time and frequency to the magnetospheric location. The angle ρ is the angle of the magnetic field with the radial direction at x₀.

As a result, each sub-burst defined by a triplet (F, t_i, w) is a Gaussian modulated by what may be called the envelope function, ℰ, of the emission region at x₀. Importantly, it is independent of the sub-burst parameters (F, t_i, w) but rather delimits where in the dynamic spectrum sub-bursts can be seen whenever they occur. The envelope only depends on the global parameters and the location of the emission region. It is thus the key element that links all sub-bursts of all events to the properties of the source.

The one-Ω envelope is the contour in the dynamic spectrum defined by γ(t_a, f) = Ω. It is represented by a blue line in the residual plots of Fig. 2 and in the additional material available on the Zenodo repository. It is important to keep in mind that sub-bursts can occur outside of this contour, both in time and frequency, indicating that the observer is seeing the outer part of the beam. We could also have used an iso-contour of ℰ(t_a, f), which has an additional dependency on f, but using γ presents the advantage of having an analytical expression (Voisin 2023), while the exponential makes it a good approximation.

Fig. 2. Residuals of a particular event from B, FRB20190519B, using best-fit parameters. Left-hand side: CHIME/FRB dynamic spectrum (middle), intensity curve (top), and spectrum (left). The one-Ω envelope (or emission angle iso-contour) of the emission beam is shown in blue, with a dashed purple line showing its characteristic frequency as a function of time. The vertical dashed red lines mark the injection time for each sub-burst. In the spectrum and intensity plots the blue lines represent the data, and the orange lines the model. Right-hand side: Same but showing the residuals between the data and the model.

In this work, we assume that each event in the CHIME/FRB catalogue is caused by a single emission region. This means that a single envelope function is associated with each event such that the two words can be used interchangeably. One could, however, imagine an event with a subset of sub-bursts associated with one envelope, and another subset associated with another envelope. In other words, emission would come from several regions simultaneously.

Propagation effects, the DM, and scattering are taken into account. In the CHIME/FRB catalogue, each event is de-dispersed with a specific DM and fitted for a specific scattering time. In this work, we fit only one DM value common to all events. It follows that we need to disperse the modelled bursts by the difference between the fitted DM and the CHIME/FRB DM. We fitted for a single common DM because it is the simplest approach, but also because our model predicts intrinsic effects that can mimic local DM variations (Voisin 2023) and that might otherwise be difficult to disentangle from local DM variations. Pseudo-Gaussian bursts have a Gaussian temporal profile at a given frequency (but the parameters of the Gaussian vary in frequency). This is because the γ²(t_a, f) can be expanded to second order in t_a. As a result, the effect of scattering can be computed using the formula in McKinnon (2014). In order to reduce the complexity of the fit, we kept the values provided by the CHIME/FRB catalogue fixed.

Altogether our model has seven global parameters (f_*,i,P_*,α,p,Ω,DM); that is, parameters that are common to all events of a given source. On the other hand, we have a relatively large number of local parameters: the two parameters (T₀, f₀) for each envelope, and within it the triplet (F, t_i, w) for each sub-burst (see also Fig. 3).

Fig. 3. Diagram of the two-stage fitting process. Global parameters are iterated using nested sampling. For a fixed set of global parameters, the ne envelopes (one per event) each have a set of local parameters comprising two envelope parameters, f0, T0 (see main text), and parameters for each of the n(k) sub-bursts (also denoted as nb/e in the main text), with k the envelope index. Local parameters were fitted independently to the dynamic spectrum of each event using least-square minimisation. Initial guesses were obtained thanks to pre-processing.

Nonetheless, the predictive power of this model is superior to an empirical model that fits each sub-burst independently from one another, in the sense that the number of parameters required is smaller. In order to demonstrate that assertion, we compare the number of parameters of our model with a fully empirical model whereby each sub-burst is modelled with a two-dimensional Gaussian of the form Aexp(^T(y−m)C(y−m)), where y = (t_a, f). Each sub-burst then requires six independent parameters: one amplitude, A, two means, m, and three parameters for the correlation matrix, C. Thus, one has $6{\overline{n}}_{\mathrm{b}/\mathrm{e}}{n}_{\mathrm{e}}$ parameters in the empirical model, where n_e is the number of events and ${\overline{n}}_{b/e}$ is the average number of sub-burst per event. On the other hand, our model has $(2+3{\overline{n}}_{\mathrm{b}/\mathrm{e}})n_{\mathrm{e}}+7$ parameters: seven global parameters, two parameters (T₀, f₀) per envelope (equivalent to events here), and three parameters (F, t_i, w) per sub-burst. Thus, it is enough to have n_e > 7 in order for our model to require fewer parameters than the empirical model, and a single event is sufficient if it has at least three sub-bursts.

However, fitting such a large parameter space remains challenging. Since the total fluence of an event can be measured independently, we further reduced the parameter space by fixing the sum of the sub-burst fluences of a given event to be equal to the measured total fluence of that event. As a result, only n_b/e − 1 fluence parameters, F, were fitted per event (instead of n_b/e).

For A, $n_{\mathrm{e}}=11,{\overline{n}}_{\mathrm{b}/\mathrm{e}}=22/11\simeq 2$ such that we have 84 parameters (77 local ones), compared to 123 for an empirical Gaussian model (assuming total fluences to be fixed in both cases). For B, $n_{\mathrm{e}}=19,{\overline{n}}_{\mathrm{b}/\mathrm{e}}=34/19\simeq 1.8$ such that we have 128 parameters (121 local ones), compared to 185 for an empirical Gaussian model.

We developed an implementation of the model as a python package called frbgeom. We release it together with posterior samples alongside this paper⁵.

2.3. Bayesian inference

We sampled the posterior distribution functions of global parameters using Multinest (Feroz et al. 2009) through the PyMultinest python interface Buchner et al. (2014)⁶. The main advantage of the nested sampling method is the fact that it explores a bounded region of the parameter space without being a priori biased by a local seed initialisation, as is the case in a typical Metropolis-Hastings Markov chain Monte Carlo. Related to this is its capacity to capture efficiently local modes without remaining trapped in one, which proved particularly helpful in this work.

Due to the large number of parameters involved, 84 for A and 128 for B, it is not possible to explore the full parameter space with reasonable computing resources. However, we used the fact that local parameters are mostly uncorrelated with global ones. Indeed, local parameters only influence their respective events, while global parameters affect the whole dataset (see Sect. 2.2). Thus, the data were fitted following two nested stages Fig. 3. First, given a set of global parameters, local parameters were fitted independently on each event’s data using a least-square minimiser⁷. The global log-likelihood was then obtained as −χ²/2 = −∑_imin(χ_i²)/2, where min(χ_i²) are the locally minimised chi-square values of each event. If this can be seen as a form of marginalisation over local parameters, we caution that it is not strictly speaking the case.

Overall, each of the two runs took ∼20 000 CPU hours, running on 96 Intel Xeon Gold 6126 cores. Since the code is currently written in python, there is a potential for a significant speed-up by converting at least some key parts of the code into a compiled language.

The local optimisation, by abruptly reducing the number of dimensions from ∼100 down to k = 7, tends to create a more disconnected χ² landscape, which is therefore harder to sample. Apart from optimisation itself, the main reason for this disconnection is the systematic noise created by the local fitting procedure. Indeed, there is a significant variance of the optimal local χ_i² depending on the initial conditions of the fit. By significant, we mean much larger than the expected standard deviation of the χ² distribution for a linear model with k degrees of freedom, which is ${\sigma }_k=\sqrt{2k}\simeq 3.7$. Any variation ≫σ_k is exponentially penalised in the posterior estimations, which can create huge difficulties in the sampling process⁸. These systematics must be accounted for, both to ensure convergence and to obtain estimates of posterior uncertainties that are sufficiently conservative.

Typically we observe a characteristic dispersion due to the fitting noise, σ_fn < 0.01χ², considering that χ² ∼ n_d.o.f. ∼ 10⁵ with n_d.o.f. the number of likelihood degrees of freedom. This is equivalent to an identical proportion of bins changing their value by one standard deviation of the background noise. Thus, this is not something that can usually be easily noticeable by eye, and the fit remains qualitatively similar. A more quantitative way of estimating the impact of this systematic noise is to compare it to the standard deviation of the likelihood’s χ², ${\sigma }_n=\sqrt{2n}\sim 500$, where n ≃ n_d.o.f. is the number of data bins. Thus we see that σ_n ∼ σ_fn, meaning that the impact of the fitting noise has a similar amplitude to that of having a different realisation of the background noise.

The solution we adopted to smooth the landscape was to include a temperature, T, such that the log-likelihood is −χ²/2T. The temperature was chosen such that var(χ²/T)≃2k, which is the variance of a Gaussian mode with k degrees of freedom. We find this value to smooth the landscape enough for the sampling to perform well. We chose T = 8 for A and T = 22 for B. Figure 4 shows that with these temperatures we recover posterior χ² samples broadly similar to the χ² distribution that would be obtained if the posterior were Gaussian (with seven degrees of freedom). Thanks to this, we expect the sampling to be reasonably well behaved. Due to large uncertainties on the surface frequency parameter, f_*, and the toroidal-to-poloidal ratio, α, sampling was performed against log₁₀(f_*/MHz) and log₁₀(α), and thus effectively used log priors.

Fig. 4. Posterior distribution of χ2 for source A (left) and B (right). As a reference, the orange lines represent the theoretical χ2 distribution for a Gaussian posterior with 7 degrees of freedom.

2.4. Priors

Nested sampling requires a finite parameter space to sample. All our priors are flat, once the logarithmic transformation mentioned above is considered; they are summarised in Table 1.

The lower bound on f_* simply coincides with the highest observed frequencies in the CHIME/FRB data, which is 800 MHz. The upper bound at 10⁴ GHz is much higher than any FRB observation, and not reached in any of our runs.

The lower bound on the spin period is set to approximately three times the duration of the longest event. This is dictated by the locality hypothesis, whereby an emission region can only occupy a fraction of the magnetosphere. The upper bound of 15 s is somewhat larger than most known confirmed magnetars (see footnote 11), albeit faster than the slowest-spinning neutron stars known (Agar et al. 2021; Caleb et al. 2022). Boundaries were not reached in any run.

The prior on α is only positive such that the orientation of the toroidal field can only be counter-rotating. Since we sample against log₁₀(α), it is impossible to account for both signs simultaneously. However, test runs with linear scaling showed that negative values were thoroughly dismissed. This is readily understood considering that a negative value cannot possibly account for the frequency downward-drifting patterns, as is demonstrated in Voisin (2023). Therefore we decided to set a lower bound at α > 0.1 (0 being excluded by the logarithmic scaling). The upper bound of 10⁴ far exceeds the magnitude of any inferred none-dipolar magnetic component in magnetars, and was never reached in our runs.

The lower bound on the toroidal power-law index, p, includes the monopole limit, p = 1, of self-similar force-free solutions Thompson et al. (2002). On the other hand, the upper bound of 8 implies an extremely steep decline of the toroidal field with respect to the dipolar one that would confine it very close to the neutron star unless extreme ratios α were also given. These boundaries were not reached in any run.

The characteristic opening angle, Ω, must be small by model assumption. An upper bound of 0.5 rad fulfils this criterion without being too constraining. It was not reached in any run.

The DM measurements for individual events in the CHIME/FRB catalogue range between [188.5, 192.5] pc cm⁻³ for A, and [348.7, 350.2] pc cm⁻³ for B. We therefore chose ranges wider than these by a couple of DM units, which appears to be sufficient as these boundaries were not reached in either run.

2.5. Mode selection

We used a Gaussian mixture model in order to characterise the modes that appear in the parameter sample. In a marginalised distribution, modes are mostly visible for the surface frequency parameter, f_*. However since this one correlates with the magnetic parameters α and p, including these generates much better results. At the same time, including the other parameters of the model creates difficulties due to the increasing dimensionality of the problem without significant improvement. The Gaussian mixture was fitted using the expectation-maximisation algorithm⁹, the result of which is shown in Fig. 5. The number of components was kept limited such that the weight of an extra component would be less than 2%. This led to five and two components for A and B, respectively.

Fig. 5. Correlation plot of the posterior sample of source A (a) and B (b) marginalised over all but three parameters: the surface emission frequency, log10(f*), the surface toroidal-to-poloidal magnetic-field ratio, log10(α), and the radial power-law exponent of the toroidal field, p. A Gaussian mixture model has been fitted to this marginalised distribution. On the 2D correlation plots, grey iso-contours of the model are shown together with coloured ellipses representing the one-standard-deviation area of each mode. The same colours are used on the 1D plots to represent the parts of the histogram attributed to each mode, and the line represents the model. In particular, the low mode is blue, the medium mode is orange, and the high mode is green.

The three main modes of A together represent 97% of the distribution (see Appendix B), and are labelled as low, medium, and high modes according to the surface frequency, f_*, they correspond to (see Fig. 5a). For B, only the low and medium mode are present, representing 100% of the sample (see Fig. 5b). Each element of the sample can be assigned to a mode based on the probability derived from the Gaussian mixture model. In the case of A, the medium mode (in orange in Fig. 5a) is actually the union of two Gaussian components that share the same qualitative properties, such that discussing them separately would be redundant. These two sub-components are visible as the two orange ellipses in the 2D histograms in Fig. 5a.

2.6. Maximum transverse size of emission regions

For a given observer direction, u, and a unit vector indicating the direction of the magnetic field, n = B/B, the characteristic region that can be seen by the observer is characterised by u ⋅ n = cosΩ, where Ω is the characteristic opening of the emission beam. The visible transverse region is the cross-section of the visible region through the plane orthogonal to the magnetic field at the emission site. Assuming that Ω is fully due to the divergence of field lines, that is that emission is infinitely beamed around the magnetic field direction, then one can define the maximum transverse region as above.

Expanding n to quadratic order in the orthogonal plane, the above equation, u ⋅ n = cosΩ, becomes that of an ellipse of the semi-major axis, a, and semi-minor axis, b. As long as the ellipse is small compared to the typical scale of variation in the magnetic field, which scales as the distance to the centre of the star, the second-order approximation is reasonable. As is shown in Fig. 6, this is generally the case in our study except for a few outliers.

Fig. 6. Posterior distributions of transverse sizes: semi-major axis, a, and semi-minor axis, b, of the approximating ellipses, relative to the radius, r. The first row represents values for A, and the second row values for B. For better visibility the left-hand histograms are limited to values of a/r < 0.25, which excludes from the plot about 0.1% of all samples for A and 0.004% for B. The colours represent modes according the convention of Fig. 5.

3. Results

3.1. General parameters

Table 2 summarises the results of our MCMC sampling (see Sect. 2), and in particular the seven global parameters described in Sect. 2.2. Correlation plots of posterior probabilities are available in Appendix A and best-fit residuals are available on the Zenodo repository. Interestingly, both sources give broadly similar results.

An example of a best-fit residual for a particular event is given in Fig. 2 (see the Zenodo repository for the others). We caution that the fit was performed on the set of all events, and therefore it cannot be judged by looking at a singe event. In the particular case of Fig. 2, one can see that the four components are following a clear downward-drifting pattern that is well within the one-Ω envelope. This is not always the case, as a burst can be seen from the outer part of the beam (see Sect. 2.2, and the Zenodo repository for examples). Moreover, if the fit appears here to be visually quite satisfactory, there are two events (out of 19) from B for which this is clearly not the case. This suggests that a more complex geometry is involved. One of them is FRB20181019A, which is the only one where the two components are fully separated by a relatively long 60 ms emission gap (see also Sect. 2.1). One may think that two emission regions are at play, in which case this event could be split in two with two different envelopes. The other ill-fitting event is FRB20190605B.

In both cases, the source DM is well determined and lies within the interval of per-event DMs calculated in the CHIME/FRB catalogue. Inclinations, i, are distributed nearly symmetrically around 90° (see Appendix A). The opening angle, Ω, is well determined with values of $5_{-1}^{+1}\times {10}^{-2}$ rad (∼3°) in both cases, which is compatible with the model’s assumption that Ω ≪ 1 rad.

The spin periods are the main difference between A and B, with $P_{\ast }=2.3_{-0.5}^{+0.5}$ s and $P_{\ast }=0.8_{-0.1}^{+0.2}$ s respectively. This parameter moderately correlates with Ω and p (see Appendix A), but the shorter period of B is caused by the stronger rate of downward drifting observed from this source (Voisin 2023).

The magnetic parameters, α, p, and surface emission frequency, f_*, are correlated (Fig. 5). In particular, the toroidal-to-poloidal ratio, α, visibly correlates with p and f_*. Its values are ${log}_{10}\alpha =2_{-1}^{+1}$ and ${log}_{10}\alpha =0.7_{-0.6}^{+0.6}$ for A and B, respectively, compatible with the two sources thanks to the relatively wide uncertainties. The toroidal index, p, is well determined and similar for the two sources, with values of $p=4.{(1)}_{-2}^{+2}$ and $p=3.{(7)}_{-1}^{+1}$ for A and B, respectively. The surface emission frequency, f_*, is in the range of 10 − 40 GHz for A and 15 − 100 GHz for B.

3.2. Emission sites

Figure 7 shows the distribution of emission sites in the magnetosphere averaged over all events (each event has its own posterior distribution). The altitude is broadly distributed between a few stellar radii up to about 1000, which, given the spin periods involved, is compatible with the model assumption of a distance, r, sufficiently low to assume co-rotation. Colatitudes are relatively well localised, respectively, around 45° and 15° for A and B. The ratio between the toroidal and poloidal components of the magnetic field at the emission site is narrowly distributed between 0.3 and 0.8. Since this ratio is directly related to the narrow-band downward-drifting pattern, we indeed expect a relatively stable value across the parameter sample. It is notably in line with Voisin (2023), which proposed that a value of order one was necessary to reproduce the observed features of repeaters.

Fig. 7. Posterior distributions of radius r, colatitude θ, and toroidal-to-poloidal magnetic ratio of emission site for source A, top row, and source B, bottom row. The colours represent modes according the convention of Fig. 5.

The characteristic transverse sizes of the emission regions (see Sect. 2.6) appear to be small, with ratios of the semi-major axis to the radius, a/r, and semi-minor axis to the radius, b/r, that are typically of the order of ∼a few × 10⁻² (Table 2 and Fig. 6). As before, we considered size distributions averaged over all events. Such a limited extent is compatible with the model assumption of locality. The semi-major axis a typically is not larger than a few times the semi-minor axis b, indicating that the emission ellipse can be quite flat. Mean semi-major axis are between ${13}_{-13}^{+16}{R}_{\ast }$ for A and $4_{-4}^{+4}{R}_{\ast }$ for B (Table 2), with the 99th percentile of the distribution reaching 1.4 × 10²R_* for A and 33R_* for B.

3.3. Modes

Posterior distribution functions can be decomposed into three main modes corresponding to low, medium, and high surface emission frequencies, f_*, corresponding, respectively, to f_* ∼ 3.5 GHz, ∼10 − 100 GHz, and ∼500 GHz (see Fig. 5). Source A shows all three modes, detailed in the extended data in Table B.1, while B shows only the low and medium ones (Table B.2). These modes also differ by their magnetic surface ratios, α. Orders of magnitudes are α ≲ 10, 10 ≲ α ≲ 100, and 100 ≲ α ≲ 10⁴ for the low, medium, and high modes, respectively. Somewhat larger values for A are found relative to B. The spin period, in particular, is similar for each mode (see Fig. 8).

Fig. 8. Posterior distributions of the spin periods for FRB20180814A (A), left-hand side, and FRB20180916B (B), right-hand side. The colours represent modes according to the convention of Fig. 5.

4. Discussion

4.1. Preferred mode

The presence of multiple modes is at least partly a consequence of the high degree of symmetry of the model. This allows distinct regions of the magnetosphere, and in particular different altitudes, to match the line of sight. We thus have three main modes corresponding to low, medium, and high altitudes (see Fig. 7).

Prior knowledge can be mobilised to select a preferred mode. So far one source, FRB 20121102, has been observed almost up to 8 GHz (Gajjar et al. 2018). If this result could be generalised to A and B, then it would discard the lower mode by setting a lower limit on f_*. A second argument for discarding this mode lies in the correspondingly low altitudes of the emission regions, less than ∼10R_*, which reduces the chances of radiation escaping the magnetosphere (Beloborodov 2021; Huang & Dai 2024; Qu et al. 2022).

The medium mode is compatible with observational limits in surface frequency with an overall interval comprised between 10 GHz and 100 GHz, and altitudes broadly distributed around ∼100R_* (Fig. 7, Appendix B). Interestingly, recent numerical work shows that FRBs might be sourced by the dissipation of shocks in the same altitude range, assuming typical magnetar settings (Beloborodov 2023; Vanthieghem & Levinson 2025). The medium mode also has the largest weight in the Gaussian mixture of both sources (see Tables B.1 and B.2).

The high mode of source A corresponds to emission altitudes of ∼1000R_* (Fig. 7, Table B.1), and predicts an upper frequency much higher than currently observed unless a cut-off at a fairly high altitude is additionally assumed. In addition, this mode corresponds to a surface magnetic-field ratio of ${log}_{10}(\alpha )=3.2_{-0.4}^{+0.4}$ (6 × 10² < α < 4.0 × 10³, Appendix B), the central value of which is almost an order of magnitude larger than the largest inferred ratios between non-dipolar and dipolar components from observations of two low-field magnetars (Tiengo et al. 2013; Rodríguez Castillo et al. 2016), while the medium and low mode are both similar to or under this value. Additionally the fact that the medium mode is common to both sources is somewhat expected if one assumes that their similar behaviour arises from a similar physical configuration. Indeed all parameters are compatible within uncertainties except for P_* and DM, meaning that both magnetic configurations are similar. For all these reasons, we suggest that the medium mode is more likely than the two others given prior knowledge, and report its detailed parameters in Table 3. Nonetheless, in the remainder of this section we conservatively quote the general values of the parameters (Table 2), unless otherwise stated. This does not affect the conclusions, and the reader may compare these values with the more accurate ones in Table 3 that are associated with the common mode.

4.2. Magnetic configuration

It is quite remarkable that the toroidal magnetic index is p ∼ 4 for both sources. This value is evocative of the quadrupolar order of a multi-polar decomposition of the magnetic field in vacuum. Vacuum configurations are relevant because the magnetic field is expected to strongly dominate the plasma at low altitudes in most of the magnetosphere. In addition, the quadrupole order is more likely to be detected at intermediate altitudes than other multipoles due to its slower decay with distance. However, vacuum quadrupolar components cannot create a homogeneous toroidal field as assumed here, but could contribute to it under suitable lines of sight.

In known semi-analytical self-similar solutions of twisted axisymmetric magnetospheres (Thompson et al. 2002; Wolfson 1995), the toroidal and poloidal components share the same radial dependence, whereby the toroidal component results from magnetospheric currents. In these particular solutions, integer exponents, p, correspond to vacuum solutions, and are thus without a toroidal component. Non-integer values correspond to twisted magnetospheres (Wolfson 1995). A recent solution (Voisin 2025) does possess the satisfying combination of an aligned dipole with a toroidal component decaying as 1/r⁴. However, the surface toroidal-to-dipolar ratio is about 1 in these solutions.

Recent numerical approaches (Viganò et al. 2011; Fujisawa & Kisaka 2014; Pili et al. 2015; Glampedakis et al. 2014; Akgün et al. 2017; Uryu et al. 2023) simultaneously solve for the stellar interior and the magnetosphere of extremely magnetised stars, and show a mix of toroidal and poloidal components concentrated near the stellar surface in the equatorial region, as well as a dipole-like poloidal field on larger scales. Indeed, except in Uryu et al. (2023), the twist is confined within a given magnetic field line. In all these cases, the ratio between poloidal and toroidal components does not exceed ∼3 at the stellar surface (Uryu, private communication). Moreover, the surface ratio remains somewhat smaller than the lower end of the range found in this paper (10 < α < 1000 for A and 1.3 < α < 20 for B at a 68% confidence level, but, respectively, 25 < α < 158 and 4 < α < 25 for the common mode). More generally, the limited surface toroidal-to-dipolar ratio is intrinsically related to the existence of a maximum twist in force-free solutions (e.g. Akgün et al. 2016; Mahlmann et al. 2019).

Thus, if the results of our phenomenological model are evocative of what may be called a twisted quadrupole, the discrepancy with available self-consistent magnetospheric solutions is significant and work is needed to reconcile both approaches. Notably, fitting multi-polar magnetic configurations or going beyond force-free solutions might be necessary.

More generally the large ratio between dipolar and non-dipolar magnetic components suggests a relatively low amplitude of the dipolar component. Indeed, assuming a maximum surface field of ∼10¹¹ T (10¹⁵ G), the dipolar component should then be B_dip ≲ 10¹⁰ T (10¹⁴ G). This is evocative of low-field magnetars, which are characterised by a relatively small dipolar component but large multi-polar ones. In particular, the three known low-field magnetars are SGR 0418+5729 with a surface dipolar strength of B_dip ≃ 6 × 10⁸ T (6 × 10¹² G) (Rea et al. 2010, 2013), SWIFT J1822.3−1606 with B_dip ≃ 3 × 10⁹ T (3 × 10¹³ G) (Scholz et al. 2012; Rea et al. 2012), and 3XMM J185246.6+003317 with B_dip ≃ 4.1 × 10⁹ T (4.1 × 10¹³ G) (Rea et al. 2014; Zhou et al. 2014). The first two show spectral lines that, if interpreted as cyclotron resonances, indicate multi-polar components that are one to two orders of magnitude stronger than the dipolar one.

4.3. Emission regions, beaming, and minimum Lorentz factor

Interestingly, the emission altitude is consistent with a pulsar-like emission mechanism from low-field magnetars or some pulsars. Indeed, Gil & Kijak (1993) have shown that radio emission starts at altitudes where the magnetic field drops below 1000 T (10⁷ G) across all types of pulsars. Thus, assuming that B < 1000 T at the emission site at r ∼ 100R_* (Table 2 and 3 and Fig. 7), we obtain a surface field of B_dip < 10⁹ T (10¹³ G). The toroidal field has been neglected in this order-of-magnitude estimate as B_ϕ < B_pol at the emission site (Fig. 7). Such dipolar strength is compatible with low-field magnetars (see also above) in a way that is noticeably consistent with the above discussion on the toroidal-to-dipolar field ratio, α.

The question of whether the emission occurs on open field lines originating from the polar cap, or on closed field lines, is important to the emission physics. The polar-cap angular radius at altitude r is ${\theta }_{\mathrm{pc}}\simeq \sqrt{r/R_{\mathrm{LC}}}$ with R_LC = cP_*/2π the light-cylinder radius¹⁰. If the colatitude of the emission region is such that θ < θ_pc, then the emission originates from the open field-line region. In our posterior sample, the probability that it occurs is ≃0.2% for A and ≃15% for B. The rest of the sample is well within the closed-field-line region. For A, emission in the open region only happens for the highest altitudes ∼1000R_* while it occurs at nearly all altitudes for B thanks to correlations with spin and colatitudes. Although this would appear to disfavour models where emission occurs on open field lines, we caution that our assumption of an aligned dipole prevents us from coming to any conclusions, unless perhaps we were to assume a priori a model requiring near-alignment (Beniamini & Kumar 2025).

Both sources, interestingly, share the same beam opening angle, Ω ≃ 0.05 rad. Since the proper frame of emission is assumed to be ultra-relativistic with respect to the observer, the beam angle is the combination of the relativistic beaming angle and of the streamline divergence. Neglecting the latter, one can estimate a lower limit on the opening angle such that Ω ≳ 1/γ, which translates into a minimum plasma Lorentz factor of γ ≳ 20. This is compatible with the limit of γ ≳ 930(ϵ/0.01)^−1/2(r/R_*)^−1/2P_*^1/2 obtained in Beniamini & Kumar (2025), where ϵ quantifies the polarisation degree, which results from constraints on beam coherence within a polar-cap geometrical model that makes different assumptions from ours.

On the other hand, assuming that the value of Ω is dominated by streamline divergence, we can estimate the characteristic transverse size of the emission region (see Sect. 2); that is, the size corresponding to an opening of one Ω. Due to the Gaussian nature of the angular emission profile, an observer may receive the emission from at most a few Ω. Approximating the transverse cross-section of the emission to an ellipse, we find semi-major axes that typically do not exceed 10% of the altitude of the emission site (Table 2 and Fig. 6). This validates the assumption of locality of the emission region.

In absolute terms, the maximum characteristic transverse size is about a ≲ 2 × 10³ km for A and a ≲ 4 × 10² km for B (assuming R_* = 12 km), with a 99% probability in both cases (see also Table 2). We conclude that the transverse size of the emission regions is generally smaller than a few thousand kilometres. This is compatible with the results in Nimmo et al. (2025), which used scintillation to constrain the transverse extent of the emission region to be smaller than 3 × 10⁴ km.

4.4. Spin and effect of spin-down

The spin periods found in this work, $2.3_{-0.5}^{+0.5}$ s for A and $0.8_{-0.2}^{+0.1}$ s for B (Table 2, Fig. 8), place these sources among the fastest known magnetars. Indeed, only two neutron stars with magnetar-like behaviours have a spin period below 2 s (Olausen & Kaspi 2014)¹¹, namely PSR J1846−0258 (Livingstone et al. 2011) and Swift J1818.0−1607 (Esposito et al. 2020) at ∼0.3 s and ∼1.4 s, respectively. A large majority of magnetars spin with periods between 5 and 12 s (Olausen & Kaspi 2014). PSR J1846−0258 is peculiar as it is a very young neutron star, ∼800 yr, which mostly behaves like a pulsar with occasional magnetar-like outbursts, and may be seen as a transition object. It is remarkable that, for B, the ∼1 s spin period matches the one inferred from the empirical scaling law relating sub-pulse quasi-periods to spin periods exhibited in Kramer et al. (2023).

The spin period decay timescale, $\tau =P_{\ast }/2{\dot{P}}_{\ast }$, caused by the electromagnetic braking of the dipolar magnetic field, is given by τ ≃ 46 yr P_*²(B_*/3 × 10¹⁰ T)⁻², where B_* is the surface value of the dipole (e.g. Lyne & Graham-Smith 2012). This leads to τ ≃ 550 yr for A and τ ≃ 66 yr for B, respectively, using median values of Table 2. Here we used a fiducial value of B_* about the median of known magnetars (Olausen & Kaspi 2014). Since burst morphology and in particular the downward-drifting rate depends on the spin period, the model predicts that morphology will evolve on the timescale given by τ. Such effects might therefore be observable in the near future. In addition, burst morphology may also evolve due to magnetic reconfigurations.

On the other hand, if the progenitors are somewhat low-field magnetars (see Sect. 4.2), their characteristic ages might be two orders of magnitude larger. In this case, we do not expect the effect of spin-down on burst morphology to be significant in the coming decades. We note that the three known low-field magnetars have much longer characteristic ages than the above estimates, of the order of 10⁵ to 10⁶ years (Rea et al. 2014, 2013; Rodríguez Castillo et al. 2016). However, it has been proposed that low-field magnetars can be born as such when the magnetic field is amplified by the Tayler-Spruit dynamo, which is triggered when the star is spun up by fallback material (Barrère et al. 2022; Igoshev et al. 2025).

4.5. Radius-to-frequency mapping

Thorsett (1991) showed that the broadband evolution of the pulse profile can be fit by a continuous law of f ∝ 1/(r − r₀)^β, where r₀ is a constant and β ∼ 1, although with variations across pulsars. Thus choosing β = 1 corresponds to the low-frequency asymptote of this empirical law. However, we caution that this type of radius-to-frequency mapping cannot explain the profile of a non-negligible number of pulsars (Posselt et al. 2021) and thus is not general.

If one considers that the emission results from curvature radiation then f ∝ 1/ρ, where ρ is the curvature radius (e.g. Ruderman & Sutherland 1975, Eq. 18). In a magnetic field dominated by its toroidal component then ρ ∼ r as field lines asymptotically tend to be circular. Of course, this assumes that all particles have the same energy, which might be locally true, but is certainly too naive in general.

A larger value of β would likely lead to emission taking place very close to the stellar surface, at a few stellar radii, because of the steep decrease in the emission frequency with distance. That behaviour was seen during the inconclusive run with f = f_*(R_*/r)^p (see Sect. 2.2). In principle, a larger distance would nonetheless be possible if it were associated with extreme values of the surface emission frequency, f_*. Either way, low-altitude emission or very large f_* both appear rather unlikely from a physical point of view, as was discussed above in Sect. 4.1.

Most importantly, further work should aim to constrain the spectral part of the model, at least by fitting for the parameter β and possibly using different types of laws to connect spectral properties with geometry. This will, however, likely require the use of a broadband dataset.

4.6. Caveats

The model in this paper relies on many assumptions that will need to be tested. The assumption that emission tangents magnetic field lines is uncertain. It appears to be correct at least for those pulsars for which the rotating vector model functions. However, it does not always work in pulsars, and the emission mechanism might be different in FRB sources. In this respect, an analysis of polarisation with our model might bring useful complementary information (Voisin 2023).

It is in principle perfectly possible within the model to use a geometry not based on magnetic field lines; however, this might require going beyond geometry and solving the plasma physics. If using the magnetic geometry is correct, the particular geometry we have used in this work might nonetheless be vastly oversimplified. The dipole is unlikely to be aligned, the toroidal field is usually concentrated near the equator in known solutions, multi-polar components might be prominent close enough to the star, and plasma effects can distort the field. Implementing more realistic geometries is certainly a way to improve the model.

In this work the spectrum of an emitting element is essentially a dirac delta function at the frequency given by radius-to-frequency mapping. One expects a broader spectrum such that several altitudes would blend at a given frequency. In addition, radius-to-frequency mapping might not apply, although one can hope for it to be a reasonable approximation locally. In this respect, the model could be improved using physically motivated emission and propagation models.

Finally, the locality assumption might fail in some circumstances. This might be the case for example for FRB20191221A, for which a 217 ms periodic signal has been observed for 3 s (Andersen et al. 2022). If this were to be associated with a rotating polar-cap emission, then the locality approximation would not apply due to the quickly varying geometry around the pole (or a sequence of local models would be necessary). Another possibility is that emission remains local, but several distinct regions contribute at the same time. This is a rather straightforward evolution of the model.

5. Conclusion

In this work, we have fitted the data from the first CHIME/FRB catalogue for the two repeating sources FRB20180814A (A) and FRB20180916B (B) with the geometrical model presented in Voisin (2023). This model associates burst morphologies with the general properties of the source, such as its spin.

For both sources we have inferred a fast spin, of the order of ∼1 − 2 s, and a high toroidal magnetic field, from a few to a hundred times the magnitude of the dipolar field at the surface, decaying as ∼1/r⁴. Based on these results, one may speculate that the relatively high degree of repetition of these sources is to be associated with a relatively young star, as can be deduced from its fast spin, and with a particularly twisted magnetosphere that is very actively releasing its energy.

As a by-product, we obtained the location and size of the emission regions, typically located at ∼100R_* and with characteristic sizes of ≲1000 km. Such size is compatible with recent observations (Nimmo et al. 2025). We also obtained a conservative lower limit on the Lorentz factor of the emitting plasma of γ ≥ 20.

Assuming a typical magnetar dipolar field, we obtained a spin-down age of the order of a few centuries and possibly as short as a few decades. In the latter case, a change in burst morphology could be detectable in the near future. However, the large toroidal-to-poloidal ratio as well as the typical emission altitudes are evocative of low-field magnetars, in which case the characteristic spin-down ages could be up to two orders of magnitude larger. This exploratory work can be improved by using a more realistic model of the magnetic configuration, using a larger dataset with broad-band data in order to constrain radius-to-frequency mapping, and including polarisation as is proposed in Voisin (2023).

Data availability

Codes, Multinest analysis results, and additional residual plots are available on the Zenodo repository at https://dx.doi.org/10.5281/zenodo.16421753

Acknowledgments

The authors would like to thank Dr. Fabrice Mottez and Dr. Marilyn Cruces for helpful discussions and reviews. The authors thank the referee for their insightful review. This work was granted access to the HPC resources of MesoPSL financed by the Region Ile de France and the project Equip@Meso (reference ANR-10-EQPX-29-01) of the programme Investissements d’Avenir supervised by the Agence Nationale pour la Recherche.

References

Appendix A: Correlation plots of posterior probabilities

Fig. A.1. Correlation plot of the posterior distribution of parameters for source A.

Fig. A.2. Correlation plot of the posterior distribution of parameters for source B.

Appendix B: Modes

All Tables

All Figures

Fig. 1. Illustration of the model, from Voisin (2023). The emitter travels from left to right, emitting in a forward cone (the cone is illustrative, in practice we use a smooth angular profile in this work). Emission frequency varies by Δf as it travels radially by Δr. Dotted lines show the observer’s direction, with γ the angle with respect to the path. The two thick ticks delimit the visible segment beyond which no emission can be seen from within the cone. The characteristic emission frequency, fc, maps to the radius, rc, at the centre of the segment. The visible segment being in a rotating frame, it varies with time until it vanishes.

In the text

Fig. 2. Residuals of a particular event from B, FRB20190519B, using best-fit parameters. Left-hand side: CHIME/FRB dynamic spectrum (middle), intensity curve (top), and spectrum (left). The one-Ω envelope (or emission angle iso-contour) of the emission beam is shown in blue, with a dashed purple line showing its characteristic frequency as a function of time. The vertical dashed red lines mark the injection time for each sub-burst. In the spectrum and intensity plots the blue lines represent the data, and the orange lines the model. Right-hand side: Same but showing the residuals between the data and the model.

In the text

Fig. 3. Diagram of the two-stage fitting process. Global parameters are iterated using nested sampling. For a fixed set of global parameters, the ne envelopes (one per event) each have a set of local parameters comprising two envelope parameters, f0, T0 (see main text), and parameters for each of the n(k) sub-bursts (also denoted as nb/e in the main text), with k the envelope index. Local parameters were fitted independently to the dynamic spectrum of each event using least-square minimisation. Initial guesses were obtained thanks to pre-processing.

In the text

	Fig. 4. Posterior distribution of χ2 for source A (left) and B (right). As a reference, the orange lines represent the theoretical χ2 distribution for a Gaussian posterior with 7 degrees of freedom.
In the text

Fig. 5. Correlation plot of the posterior sample of source A (a) and B (b) marginalised over all but three parameters: the surface emission frequency, log10(f*), the surface toroidal-to-poloidal magnetic-field ratio, log10(α), and the radial power-law exponent of the toroidal field, p. A Gaussian mixture model has been fitted to this marginalised distribution. On the 2D correlation plots, grey iso-contours of the model are shown together with coloured ellipses representing the one-standard-deviation area of each mode. The same colours are used on the 1D plots to represent the parts of the histogram attributed to each mode, and the line represents the model. In particular, the low mode is blue, the medium mode is orange, and the high mode is green.

In the text

Fig. 6. Posterior distributions of transverse sizes: semi-major axis, a, and semi-minor axis, b, of the approximating ellipses, relative to the radius, r. The first row represents values for A, and the second row values for B. For better visibility the left-hand histograms are limited to values of a/r < 0.25, which excludes from the plot about 0.1% of all samples for A and 0.004% for B. The colours represent modes according the convention of Fig. 5.

In the text

	Fig. 7. Posterior distributions of radius r, colatitude θ, and toroidal-to-poloidal magnetic ratio of emission site for source A, top row, and source B, bottom row. The colours represent modes according the convention of Fig. 5.
In the text

	Fig. 8. Posterior distributions of the spin periods for FRB20180814A (A), left-hand side, and FRB20180916B (B), right-hand side. The colours represent modes according to the convention of Fig. 5.
In the text

	Fig. A.1. Correlation plot of the posterior distribution of parameters for source A.
In the text

	Fig. A.2. Correlation plot of the posterior distribution of parameters for source B.
In the text