LOFAR Detection of 110–188 MHz Emission and Frequency-dependent Activity from FRB 20180916B

On 128 occasions between 2019 June 6 and 2020 August 26, FRB 20180916B was observed with LOFAR (van Haarlem et al. 2013) using its beam-formed modes (Stappers et al. 2011) for 112 hr. The Cobalt correlator and beam former (Broekema et al. 2018) coherently combined signals from the HBAs of the LOFAR Core stations to create a tied-array beam pointing to the best-known position of FRB 20180916B. The best CHIME/FRB localization of α_J2000 = 01^h57^m43^s, ${\delta }_{{\rm{J}}2000}=+65^\circ 42^{\prime} 00^{\prime\prime}$ (estimated uncertainty of $\sim {2}^{{\prime} };$ derived from baseband data) was used until 2019 August, while the much more precise position from the EVN localization, α_J2000 = 01^h58^m007502, ${\delta }_{{\rm{J}}2000}=+65^\circ 43^{\prime} 00\buildrel{\prime\prime}\over{.} 3152$ (2.3 mas uncertainty), from Marcote et al. (2020) was used afterward. The angular separation between the two positions is $2\buildrel{\,\prime}\over{.} 1$. As all observations used either the innermost 22 or 20 LOFAR Core stations, the ∼3' FWHM of the tied-array beam was slightly offset from the actual celestial position of FRB 20180916B in the observations prior to 2019 August, which led to a factor of ∼2 lower sensitivity for those early observations. Observations were obtained with the source at altitudes ranging from 37° up to culmination at 77°, with 86% of the observations being obtained at altitudes of 60° or higher.

The earliest LOFAR observations (up to approximately source activity Cycle 26; see Figure 1) were scheduled in response to CHIME/FRB detections. After the identification of a 16.3 day activity period by CHIME/FRB, observations were scheduled close to the peak in the CHIME/FRB-derived activity window (Cycle 26 and later). The observations in Cycle 43 were intended to broadly cover the full range of activity phases.

For all observations, 400 subbands of 195.3125 kHz each were recorded, covering observing frequencies of 110–188 MHz for a total bandwidth of 78.125 MHz. For 56 hr of the observations, the Cobalt correlator and beam former generated total intensity Stokes I filter-bank data, with a frequency and time resolution of 3.05 kHz and 983.04 μs, respectively. Nyquist-sampled, dual-polarization complex voltages ("CV" data, with 195.3125 kHz frequency resolution and 5.12 μs time resolution) were recorded for the other 56 hr of observations. For the dispersion measure (DM) of FRB 20180916B, DM = 348.772 ± 0.006 pc cm⁻³ (Nimmo et al. 2021), the dispersion delay from infinite frequency to the top of the observed LOFAR band (188 MHz) is 40.9 s, and the delay within the band, from 188 MHz down to 110 MHz, is 78.7 s.

The CV data allow for coherent dedispersion, a technique that we are employing to search for FRBs and millisecond pulsars using LOFAR and our DRAGNET GPU cluster (Bassa et al. 2017a, 2017b). In this paper, we apply coherent dedispersion for burst characterization (see Section 3), but, given the large burst widths and to allow for a homogeneous burst search of both the Stokes I and the CV data sets, we chose not to coherently dedisperse the CV data for the purposes of searching for signals. Instead, the CV data were channelized and time-averaged and the polarizations were squared and summed offline to form Stokes I filter banks with the same time and frequency resolution as the Stokes I filter banks (3.05 kHz and 983.04 μs) generated by Cobalt in real time. To reduce their size, these filter banks were first averaged in frequency by a factor of 16 with digifil (van Straten & Bailes 2011) using incoherent dedispersion to DM = 350 pc cm⁻³. With this setup, dispersive smearing due to incoherent dedispersion at the DM of FRB 20180916B varies from 1.3 to 6.7 ms over the LOFAR band. This temporal smearing is negligible compared to the burst widths. Radio frequency interference (RFI) was identified using the rfifind tool from the Presto software (Ransom 2001; Ransom et al. 2002) and replaced with random noise of the appropriate mean and standard deviation. Next, dedispersed time series were generated between DMs of 300 and 400 pc cm⁻³ with steps of 0.1 pc cm⁻³ using the GPU-accelerated Dedisp dedispersion library (Barsdell et al. 2012).

To maximize sensitivity toward possibly narrowband radio bursts (motivated by Gourdji et al. 2019; Kumar et al. 2021), these time series were created for the full observing band (110–188 MHz), as well as three overlapping halves (110–149, 130–169, and 149–188) and seven overlapping quarters of the band (110–129 MHz, 120–139 MHz, etc.). These time series were cross-correlated with top-hat functions with widths up to 150 ms using a GPU-accelerated version of single_pulse_search.py from Presto to search for bursts. All candidate burst events with S/N > 7 were visually inspected to distinguish bursts from residual RFI. For the 12 data sets with known bursts, we verified that the burst selection is complete by assessing the burst candidates using the FETCH deep learning–based classifier (Agarwal et al. 2020). All of the single pulses identified were grouped using Single-pulse Searcher (Michilli & Hessels 2018), and redundant burst candidates were eliminated before putting them through the FETCH classifier.

The uGMRT (Gupta et al. 2017) observations of FRB 20180916B were carried out on three different days: 2020 May 29, 2020 June 15, and 2020 July 1. The observations were intentionally scheduled close to the peak of the CHIME/FRB-derived activity window, and the telescope was phased up toward the EVN position of FRB 20180916B during all observations. On 2020 May 29, we used both Band 2 (180–280 MHz) and Band 3 (250–350 MHz) simultaneously with two phased-array beams employing two subarrays: one using 10 of the available 29 central-square antennas at Band 2 (180–280 MHz) and another at Band 3 (250–500 MHz) using the remaining 19 antennas. We used the uGMRT Wideband Backend (Reddy et al. 2017) to record coherently dedispersed total intensity filter-bank data with the passband split into 2048 channels with 327.68 μs sampling time in both phased-array beams. The observations were divided into eight 20 minute sessions with a provision to rephase the subarrays in between the sessions to account for the temporal instrumental gain and ionospheric changes. We used the same setup on June 15 with five 20 minute sessions, but the Band 2 data were not usable due to strong RFI. On 2020 July 1, we observed only in Band 3 set to 250–450 MHz using 24 antennas in four 40 minute sessions. The coherently dedispersed filter-bank data were recorded with 81.92 μs time resolution and 2048 frequency channels across the 200 MHz band. On all 3 days, we also recorded "ON" source and "OFF" source data on 3C48 for calibration purposes.

To prepare the data for searching, we identified and mitigated narrowband RFI and broadband time-domain RFI using gptool (Susobhanan et al. 2020). The single-pulse search was carried out using a machine-learning technique based on Zhang et al. (2018). We have trained our convolutional neural network model using archival uGMRT data and simulated CHIME/FRB-like FRBs with various burst morphologies. The details of our implementation will be published elsewhere. This particular search is tuned for 200 pc cm⁻³ ≤ DM ≤ 500 pc cm⁻³ and is sensitive to pulse widths up to 256 ms. The dynamic spectra of candidates were visually examined to distinguish astrophysical signals from spurious RFI.

CHIME/FRB searches intensity data from 1024 stationary synthesized beams for dispersed transients over an ∼200 deg² field of view in the 400–800 MHz octave in real time (CHIME/FRB Collaboration 2018). The intensity data have a 0.98304 ms time resolution and 16,384 frequency channels.

The median daily exposure of the experiment to the sky position of FRB 20180916B is 746 s (i.e., the time spent by the source within the FWHM of the synthesized beams at 600 MHz). The transit time of the source through the primary beam of the instrument, however, is much longer (∼40 minutes), albeit with significant variation in sensitivity that is still in the process of being quantified outside of the FWHM of the synthesized beams.

A total of 18 bursts were detected using LOFAR: 14 were found in the full-bandwidth data, and four fainter bursts were identified in the time series generated from the half- and quarter-bandwidth data segments. The dynamic spectra of these bursts are presented in Figure 2. Four of the bursts occurred during CV observations and thus have available polarimetric information. These are labeled CV1–CV4, while the bursts obtained in total intensity Stokes I data are labeled I01–I14.

Zoom In Zoom Out Reset image size

Bursts I01–I08 correspond to the bursts L02–L09 presented by Pastor-Marazuela et al. (2020). Our analysis did not recover their L01 burst due to periodic baseline variations affecting 2.6 hr of the LOFAR observations. These variations were the result of a rare temporary clock skew at one of the LOFAR stations, where some network packets of this station arrived out of sync at the COBALT beam former (Broekema et al. 2018) and were discarded. Hence, the number of stations added coherently in the tied-array beam varied at the ∼1 s beam former block size of COBALT. The resulting baseline variations increased the noise in the dedispersed time series of these observations, placing the S/N of this burst below our detection threshold.

FETCH positively identified the same set of 18 bursts found through visual inspection, including the four faint bursts found only in the subband search, using the full-bandwidth data. The 14 bursts identified by visual inspection in the full-bandwidth data were classified as astrophysical pulse candidates with probabilities greater than 83% by all 11 available models (labeled a–k) of FETCH. The three fainter bursts, I01, I02 and I11, were positively identified by at least nine of the 11 models with probabilities greater than 64%. The faintest burst, CV2, could only be identified by FETCH models e and h, each with a probability of more than 93%. No additional bursts were found in these data sets.

The LOFAR burst properties were determined by fitting Gaussian profiles to the time and frequency averages of the dynamic spectra to obtain their arrival times, temporal and spectral widths (FHWM), and emission frequencies (Table 1). All observed bursts are band-limited, with spectral widths (defined as FWHM) ranging from 20 to 50 MHz and the majority of the bursts, 15 out of 18, peaking in brightness above 160 MHz. The temporal width of the bursts varies between 40 ms for bursts peaking in the top of the LOFAR band (∼180 MHz) and 160 ms near the bottom of the LOFAR band (∼120 MHz).

To calibrate the LOFAR detections, we have subtracted the mean of an off-burst region in each frequency channel and divided by the standard deviation of the off-burst region. We then convert each channel i from S/N to flux units using the radiometer equation,

$$\begin{eqnarray}&&{\rm{\Delta }}{S}_{i}=\displaystyle \frac{{T}_{\mathrm{sys},i}}{{G}_{i}\sqrt{{n}_{{\rm{p}}}{\rm{\Delta }}{\nu }_{i}{t}_{{\rm{s}}}}},\end{eqnarray} \tag{ 1 }$$

where T_sys,i is the system temperature (receiver and sky), G_i is the gain, n_p = 2 is the number of summed polarizations, Δν_i is the channel bandwidth, and t_s is the sampling time. For the calculation of T_sys and G, we take into account the number of core stations used and correct for the zenith-angle dependence, as described by Kondratiev et al. (2016). We calculate these values at the times of the burst detections. In the center of the band (∼149 MHz), T_sys ≈ 1090 K (note that FRB 20180916B is at a Galactic latitude of b = 37 with T_sky ≈ 700 K at 150 MHz) and G ≈ 4.3 K Jy⁻¹, on average. We repeat this measurement using six different independent but adjacent ∼100 ms off-pulse regions in the same observation as the burst and quote the average (band-averaged) fluence and peak flux of the bursts, as well as the standard deviation across the six samples, in Table 1. Note that the systematic uncertainty in the flux measurements could be larger than reported in Table 1 due to unaccounted-for contributions from the ionosphere and the Galactic plane or other bright sources in the primary beam of the HBAs. As the systematic errors likely underestimate system noise, the fluence measurements are likely also underestimated.

Even after dedispersion to the best-fit DM from Nimmo et al. (2021), the brightest LOFAR bursts show residual time delays toward lower observing frequencies, with one burst (CV4) broadening with decreasing frequency. Due to the absence of visible burst substructure, it is unclear if these effects are due to DM underestimation, multipath scattering, or the frequency drifts of unresolved subbursts (i.e., the "sad trombone" effect; Hessels et al. 2019). Deviations from the canonical Δt ∝ ν⁻² dispersion relation are expected at low radio frequencies (see Hassall et al. 2012, and references therein), and it is also possible that the residual time delays in CV4 result from such deviations.

To measure a scattering timescale, we use a least-squares fitting routine, previously developed for CHIME/FRB bursts (e.g., CHIME/FRB Collaboration et al. 2019a; Josephy et al. 2019), to model the dynamic spectra of the LOFAR bursts. We analyze the three brightest bursts that can be coherently dedispersed: bursts CV1, CV3, and CV4, where dspsr (van Straten & Bailes 2011) was used to generate coherently dedispersed single-pulse profiles. We fixed DM = 348.772 pc cm⁻³ (Nimmo et al. 2021) and assumed that the bursts consist of only one component. We referenced the scattering timescale to 150 MHz and fixed the scattering index to −4. The results are provided in Table 2 and Figure 3. The measured scattering timescale for the two bursts with the highest S/N (CV1 and CV4) is ∼50 ms at 150 MHz. The other measured timescale, for burst CV3, is more uncertain due to the low-detection S/N and narrow bandwidth of the burst.

Zoom In Zoom Out Reset image size

We recorded polarization, using LOFAR's orthogonal linear feeds, for bursts CV1–4; three of these bursts are bright enough to perform polarimetric analyses. We corrected for azimuth and elevation-dependent gain using the LOFAR beam model, as implemented within the dreamBeam package. ²⁶ No additional polarization calibration, beyond that already implemented to form the individual station and Core tied-array beams (Stappers et al. 2011), was performed; thus, a slight degradation of the polarization fraction (≲5%) could arise, e.g., from the thermal expansion of the cables. The absence of an absolute polarization calibration also prevents us from comparing the PA between different bursts. In this analysis, the frequency resolution was increased to ∼12 kHz (16 channels synthesized within each 195 kHz subband) to have an intrachannel depolarization smaller than 1% across the whole band at the value of the rotation measure (RM) reported by CHIME/FRB Collaboration et al. (2019a). We measured the RM value of each burst by using RM synthesis ²⁷ (Burn 1966; Brentjens & de Bruyn 2005) and a deconvolution algorithm (Heald 2009); the Stokes parameters of each burst have been corrected for the resulting value.

We report the polarization profiles for the three bursts in Figure 4 using the coherently dedispersed pulse profiles from Section 3.1. The linear polarization fractions of CV1, CV3, and CV4 are roughly 70%, 60%, and 30%, respectively, much lower than the ∼100% reported at higher frequencies of 300–1700 MHz (CHIME/FRB Collaboration et al. 2019a; Chawla et al. 2020; Nimmo et al. 2021). Noutsos et al. (2015) performed a long observation of PSR B2217+47 to study the depolarization fraction as a function of the hour angle (HA) for LOFAR. They found a depolarization fraction of <10% for ∣HA∣ < 6 hr and zenith angle <50°. Our observations of FRB 20180916B are well within this range, with a maximum HA of 2.1 hr and a maximum zenith angle of 203. Given the absence of any visible artifacts in the Faraday dispersion function that could indicate the presence of signal leakage or uncorrected delay between the polarization channels (such as a peak at RM = 0 rad m⁻² or symmetric peaks around RM = 0 rad m⁻²), we conclude that the observed depolarization toward lower radio frequencies is predominantly astrophysical in origin.

Zoom In Zoom Out Reset image size

Part of the depolarization at lower frequencies could therefore be intrinsic to the source emission or related to extrinsic propagation effects. In particular, CV4 is observed at the lowest frequencies and has the lowest linear polarization fraction, which is compatible with depolarization due to scattering, as observed in some pulsars (e.g., Noutsos et al. 2015; Xue et al. 2019). As at higher frequencies, the circular polarization fraction here is consistent with 0%. The PA curve is nearly flat at higher frequencies, though subtle variations are seen on short timescales (Nimmo et al. 2021). While the PA curve of CV4 is consistent with being flat (${\chi }_{\mathrm{red}}^{2}=1.1$), the PA curve of the brightest burst, CV1, shows a hint of an increase at later times, with ${\chi }_{\mathrm{red}}^{2}=2.2$ with respect to a straight line. Day et al. (2020) similarly found evidence for a time-varying PA in the wide and potentially scatter-broadened FRB 190608.

The lack of refined polarization calibration is not expected to affect the RM values measured by LOFAR (e.g., Sobey et al. 2019), which are reported in Table 2. We compare the new LOFAR RM values obtained here with other measurements presented by CHIME/FRB Collaboration et al. (2019a) and Chawla et al. (2020). The ionospheric contribution to the observed burst RMs was determined using ionFR ²⁸ (Sotomayor-Beltran et al. 2013), which utilizes data from the International Geomagnetic Reference Field and IONosphere map EXchange (IONEX) global ionospheric maps. The resulting values are shown in Figure 5 as a function of date and source activity phase, as measured by PR3 and refined in Section 3.5. There is a hint of RM values changing systematically as a function of activity phase, as measured by different instruments in different activity cycles. However, the RM variations may simply be stochastic; more detections are needed to investigate this further.

Zoom In Zoom Out Reset image size

We detected five bursts in the uGMRT data: two in the 250–350 MHz band on 2020 June 15 and three in the 250–450 MHz band on 2020 July 1. We did not detect any bursts on 2020 May 29. The burst dynamic spectra are presented in Figure 6. To calibrate the bursts, we have used the counts per jansky estimated for every clean frequency channel from the "ON" and "OFF" data of 3C 48. The conversion factors are multiplied with the filter-bank counts to get the calibrated data. The measured peak flux densities and fluences are presented in Table 1.

Zoom In Zoom Out Reset image size

We caution the reader not to overinterpret the spectral structure of the uGMRT bursts due to the baseline variations from residual broadband RFI in the dynamic spectra. The four subbursts in burst G5, however, are likely of astrophysical origin; the dispersion delay (quadratic and ∼16 s from 450 to 250 MHz) and dispersion smearing within the frequency channels of these candidates are as expected for the DM of the source, and we find no candidate bursts at other DMs with similar properties. Using DM_phase ²⁹ (Seymour et al. 2019), we align the subbursts of G5 and find an optimum DM = 349.5 ± 0.1 pc cm⁻³. In Figure 7, we show a comparison of the dynamic spectra of the bursts, dedispersed to the fiducial and optimum DM values.

Zoom In Zoom Out Reset image size

We present 17 new bursts detected by CHIME/FRB since PR3. ³⁰ The burst dynamic spectra are presented in Figure 8. The morphologies of the bursts are comparable to those previously detected by CHIME/FRB; they exhibit narrow (50–150 MHz) bandwidths and sometimes show downward-drifting subbursts (CHIME/FRB Collaboration et al. 2019a).

Zoom In Zoom Out Reset image size

We construct burst models and measure peak fluxes and fluences as in CHIME/FRB Collaboration et al. (2020a). In summary, we fit single- or multicomponent models of dynamic spectra to the 0.98304 ms total intensity data for each burst using a least-squares algorithm. Data are flux calibrated using transits of steady sources and scaled by the beam response using the best-known location of FRB 20180916B and a model for the synthesized beams. For bursts detected outside of the FWHM of the synthesized beams at 600 MHz, peak fluxes and fluences are lower limits.

We recalculate the activity period and burst rate of FRB 20180916B using the same methods as described in PR3. We measure the activity period to be 16.33 ± 0.12 days, with a 5.2 day window, based on the now-reported 55 CHIME/FRB detections. We use reference MJD ϕ₀ = 58,369.40 to put the average arrival time of the CHIME/FRB bursts at ϕ = 0.5. We estimate the CHIME/FRB detection rate to be 0.8 ± 0.3 bursts hr^–1 above a fluence threshold of 5.1 Jy ms for a ±2.6 day interval around each cycle of activity (31 detections in 39.1 hr of exposure). In the 1σ activity window of ±0.96 days around each cycle of activity, we estimate the detection rate to be ${1.5}_{-0.6}^{+0.8}$ bursts hr^–1 (22 detections in 14.9 hr of exposure).

Both LOFAR and uGMRT predominantly observed FRB 20180916B during the 5.4 day activity window from PR3, though in later activity cycles, LOFAR also targeted activity phases outside of the activity window (see Figure 1). Figure 9 shows the barycentered arrival times of the CHIME/FRB, uGMRT, and LOFAR bursts folded on the activity period of 16.33 days as a function of observing frequency, MJD, and the cumulative fractions and histograms of the exposure of each instrument. The bursts observed with LOFAR fall in the range of activity phases of 0.53 < ϕ < 0.79, corresponding to a range of 4.1 days within the 16.33 day activity period. This observed LOFAR activity window is nominally shorter than the 5.2 day activity window width observed from the CHIME/FRB bursts, but additional LOFAR observations may show the activity window to be wider. The average activity phase of the LOFAR bursts is ϕ ∼ 0.66, which is offset by 2.6 days from the average activity phase of FRB 20180916B determined from CHIME/FRB bursts (PR3). We note that the LOFAR burst detections as a function of activity phase are not just simply a reflection of the observing exposure. Through a two-sample Kolmogorov–Smirnov test, we can reject the null hypothesis that the cumulative distribution functions of the LOFAR exposure and bursts shown in Figure 9 are drawn from the same distribution (p ∼ 10⁻¹⁰).

Zoom In Zoom Out Reset image size

We performed numerical simulations to investigate the impact of the nonuniform exposure of the LOFAR observations as a function of activity phase (see Figures 9(c) and (d)) on the properties of the activity window at LOFAR frequencies. We define an activity window by its width w and central phase ϕ₀ and compute the effective exposure of the LOFAR observations that fall within this activity window. The top panel of Figure 10 shows the LOFAR exposure given the properties of the activity window. As 18 bursts were observed with LOFAR, the exposure will provide the burst rate r as a function of w and ϕ₀. Given the activity window properties and the burst rate, we draw burst arrival times from a uniform distribution within the activity window for the activity cycles spanning the LOFAR observations (see Figure 1). For each w and ϕ₀ combination, we run multiple simulations to obtain the fraction of simulations where the simulated bursts fall within both the LOFAR observations and the observed phase range of LOFAR bursts (0.53 < ϕ < 0.79). The bottom panel of Figure 10 shows this fraction, which we treat as the probability that all simulated bursts fall within the observed LOFAR activity phase range.

Zoom In Zoom Out Reset image size

The simulations show that the observed activity window of the LOFAR bursts is delayed with respect to the activity window observed by CHIME/FRB. The best-fit parameters for w and ϕ₀ are determined through a Markov Chain Monte Carlo analysis using the emcee software (Foreman-Mackey et al. 2013) to sample and maximize the probability that all simulated bursts coincident with the LOFAR observations fall within the observed LOFAR activity range. Flat priors were used for both parameters, though the width was limited to the observed minimum range of w > 4.1 days. The posterior distributions were obtained using 32 walkers for 20,000 steps, well beyond 100 times the largest autocorrelation of the fitted parameters. After discarding a burn-in phase of 1000 steps and thinning by 35 steps, we obtained ${\phi }_{0}={0.72}_{-0.04}^{+0.07}$ and $w={5.0}_{-0.8}^{+2.3}$ days. These values correspond to a LOFAR burst rate of $r={0.32}_{-0.04}^{+0.08}$ hr⁻¹ for a fluence limit of 26 Jy ms.

The object FRB 20180916B already held the record for the FRB with the lowest-frequency emission detected to date. Using the GBT and SRT, respectively, Chawla et al. (2020) and Pilia et al. (2020) previously presented a total of 10 burst detections at frequencies as low as 300 MHz. Here we present 18 bursts detected in the 110–188 MHz band.

Notably, the majority of these bursts are brighter in the top half of the band, but burst CV4 clearly demonstrates that emission can be detected down to at least 110 MHz (Figure 2). Searches for FRB 20180916B emission at radio frequencies <100 MHz are thus well motivated and are underway using joint observations with LOFAR and the New Extension in Nançay Upgrading LOFAR (NenuFAR; Bondonneau et al. 2020). These observations can better quantify whether the observed burst rate is systematically reduced at lower radio frequencies. While only three of the 18 bursts we have observed here are emitting predominantly in the lowest part of the LOFAR HBA band (110–140 MHz), one must take observational biases like increased sky temperature (T_sky) and larger pulse width into account (our detection metric scales as $F/\sqrt{w}$, where F and w are burst fluence and width, respectively). The current sample is too small to make robust statements about the burst rate declining toward the bottom of the LOFAR HBA band.

Similarly, the detection of FRB 20180916B at 110–188 MHz provides renewed hope for detections in wide-field, low-frequency surveys (e.g., Sanidas et al. 2019).

The simple fact that FRB 20180916B is visible at 110 MHz sets new requirements on the emission mechanism and constraints on the effect of free–free absorption in the local medium. Based on the lack of free–free absorption at 300 MHz, Chawla et al. (2020) argued that FRB 20180916B is not associated with a hypercompact H ii region or a young (<50 yr) supernova remnant. This is consistent with the lack of local Hα emission or a persistent radio counterpart, as shown by Tendulkar et al. (2021) and Marcote et al. (2020), respectively. Here we show that the circumburst environment is optically thin to free–free absorption at 110 MHz. Considering an ionized nebula of size L_pc and DM < 70 pc cm⁻³ (Marcote et al. 2020), we use the following expression for the optical depth due to free–free absorption (Condon & Ransom 2016):

$$\begin{eqnarray}\begin{array}{rcl}{\tau }_{\mathrm{ff}} & = & 1.6\times {10}^{-3}\times {\left(\displaystyle \frac{T}{{10}^{4}{\rm{K}}}\right)}^{-1.35}\times {\left(\displaystyle \frac{\nu }{1\mathrm{GHz}}\right)}^{-2.1}\\ & & \times \displaystyle \frac{1}{{f}_{\mathrm{eff}}{L}_{\mathrm{pc}}}\times {\left(\displaystyle \frac{\mathrm{DM}}{70\mathrm{pc}{\mathrm{cm}}^{-3}}\right)}^{2}\ll 1,\end{array}\end{eqnarray} \tag{ 2 }$$

where f_eff accounts for the volume-filling factor and the electron density fluctuation in the circumburst medium. For f_eff = 1, we find L ≫ 0.16 pc (T/10⁴ K)^−1.35. The Crab Nebula, for comparison, is about 1.7 pc in radius. Piro (2016) discussed how a surrounding supernova remnant can absorb bursts at low radio frequencies. In the context of that work (see their Figure 4), our detection of a burst down to 110 MHz implies a source age of at least 100–600 yr for assumed ejecta masses between 3 and 10 M_⊙. Similarly, Marcote et al. (2020) discussed how the relatively low RM and lack of a persistent radio counterpart constrain the age of FRB 20180916B to be ≥300 yr, in the context of models that describe it as a young flaring magnetar in a dense nebula (Metzger et al. 2019).

Models that describe FRBs via synchrotron maser emission from decelerating relativistic blast waves predict that the intrinsic FRB fluence is lower at lower radio frequencies (Metzger et al. 2019). In contrast, at LOFAR frequencies of 110–188 MHz, we see >10× higher fluences compared to bursts detected by CHIME/FRB at 400–800 MHz (Table 1). Similarly, the detection of 110 MHz emission also places constraints on magnetospheric emission heights and models, e.g., the curvature radiation model (Kumar et al. 2017), where a particle density of ∼10¹⁷ cm⁻³ is necessary to produce a typical FRB luminosity (L_iso ∼ 10⁴³ erg s⁻¹), whereas the particle density also cannot be much larger than this value for the emitting region to remain transparent to ∼100 MHz emission.

The spectrotemporal properties of the 17 CHIME/FRB bursts presented here (Figure 8) are consistent with those of the 38 CHIME/FRB bursts previously presented in the literature. As such, we do not discuss their properties in detail.

The five new uGMRT bursts (detected at 200–450 MHz) complement the 15 previously detected at 550–750 MHz (Marthi et al. 2020). Of the new uGMRT bursts, G5 (Figures 6 and 7) has the most intriguing spectrotemporal properties, showing substructures that do not completely match the typical "sad trombone" downward-drift features seen from repeaters (e.g., Hessels et al. 2019). A similar morphology has been observed in some bursts from FRB 20121102A (Caleb et al. 2020; Hilmarsson et al. 2021) and FRB 20180916B (Chawla et al. 2020). However, it is unclear whether these represent drifting subbursts within a single burst envelope or a closely spaced set of separate bursts peaking at different frequencies and then each individually drifting downward in frequency. The short-timescale substructure of burst G5 also allows us to better investigate possible DM variations (as we discuss in the next subsection).

The 18 LOFAR bursts constitute our first robust view of FRB emission below 300 MHz. These bursts have large widths (40–160 ms; Table 1) compared to other FRBs in general (e.g., Figure 13 in Petroff et al. 2019) and FRB 20180916B at high frequency (Nimmo et al. 2021). This broadening is not due to dispersion smearing (see Section 2.1 and the coherently dedispersed waterfalls in Figure 3). For comparison, the typical burst width at 1.7 GHz is ∼2–3 ms, and burst substructure has been detected on a wide range of timescales, down to ∼3–4 μs (Nimmo et al. 2021). It is not completely clear whether the larger burst widths at lower frequencies are simply due to scatter broadening (extrinsic multipath propagation through the intervening material) or whether they reflect an intrinsic aspect of the emission process, like the "sad trombone" effect seen from FRB 20180916B and other repeaters, where the subburst drift to later times increases toward lower frequencies (Hessels et al. 2019).

Previously, Marcote et al. (2020) estimated a scattering time of 2.7 μs at 1.7 GHz, and Chawla et al. (2020) constrained the scattering time to be <1.7 ms at 350 MHz (this is compatible with the higher-frequency measurement). For a scattering time that scales with frequency as τ_scatt ∝ ν⁻⁴, the effect is 30× larger at 150 MHz compared to 350 MHz and 100× larger at 110 MHz. Based on the measured scattering time of Marcote et al. (2020), we expect a scattering time of ∼40 ms at 150 MHz for a ν⁻⁴ scaling. This matches well with the modeled scattering times of the high-S/N CV data bursts CV1 and CV4 (Table 2). Regardless of modeling, bursts I14 and CV4 are the brightest detected bursts in the top and bottom halves of the LOFAR HBA band, respectively, and both show clear asymmetric tails (as do other high-S/N bursts). Their measured widths of 54 ± 2 ms at 161 MHz (I14) and 158 ± 7 ms at 116 MHz (CV4) are also consistent with what one would predict by extrapolating from the scattering measurement of Marcote et al. (2020). However, the LOFAR bursts may be additionally broadened by a poorly resolved "sad trombone" effect as well. The drift toward later times increases toward lower frequencies (Hessels et al. 2019; Josephy et al. 2019) and, based on FRB 20180916B observations in the CHIME/FRB band (PR3), could be ∼10 ms per ∼50 MHz at LOFAR frequencies. We also note that many bursts from FRB 20121102A show asymmetric burst profiles regardless of scattering (Hessels et al. 2019). Scatter broadening can blend multiple intrinsic burst components (see also Section 4.3 of Day et al. 2020). Additional LOFAR observations spanning a wide range of epochs can constrain whether the scattering time varies, as is seen from the Crab pulsar (Driessen et al. 2019). Our burst detections spanning ∼10 months show no obvious evidence for this, however.

The detection of FRB 20180916B, with DM ∼ 350 pc cm⁻³, in the LOFAR HBA band contrasts with results from the LOFAR HBA (110–188 MHz) census of slow and millisecond pulsars (Bilous et al. 2016; Kondratiev et al. 2016), as well as the new and known pulsars detected in the LOFAR Tied-Array All-Sky survey (LOTAAS; Sanidas et al. 2019). None of these observations detect Galactic pulsars with DM ≳ 220 pc cm⁻³. The absence of LOFAR pulsar detections above this limit is consistent (Sanidas et al. 2019) with predictions from the Galactic scattering relations derived by Bhat et al. (2004) and Geyer et al. (2017). Indeed, extragalactic FRBs are often significantly less scattered than Galactic pulsars of comparable DM (see, e.g., Figure 16 of Cordes & Chatterjee 2019), though note that FRB 20180916B has a low Galactic latitude, b = 37, toward Galactic longitude l = 1297.

In previous simultaneous LOFAR HBA, GBT, and CHIME/FRB observations, bursts were detected at >300 MHz but not at 110–188 MHz (Chawla et al. 2020). Of the CHIME/FRB, uGMRT, and LOFAR observations presented here, the arrival times of five CHIME/FRB bursts from CHIME/FRB Collaboration et al. (2020a) overlap with LOFAR observations, of which three have been presented in Chawla et al. (2020). None of these bursts have counterparts in the 110–188 MHz band of LOFAR. Similarly, Pearlman et al. (2020) used simultaneous CHIME/FRB and DSN 70 m dish observations to demonstrate a burst detection in the CHIME/FRB band but none at 2.3 or 8.4 GHz. Clearly, FRB 20180916B bursts have a low instantaneous bandwidth, as has been seen for FRB 20121102A (Gourdji et al. 2019; Majid et al. 2020) and beautifully demonstrated for FRB 20190711A (Kumar et al. 2021). The LOFAR HBA bursts we present here show bandwidths of 20–50 MHz and are consistent with an emerging picture in which repeating FRB bursts have typical fractional bandwidths (BW/ν_obs) of ∼20% (e.g., CHIME/FRB Collaboration et al. 2019a; Gourdji et al. 2019). The large sample of FRB 20180916B bursts that are now available from 110–1700 MHz could also be stacked to determine an average spectral index, though that requires careful consideration of CHIME/FRB beam effects and other instrumental biases/selection effects and is beyond the scope of this work.

Lastly, for three LOFAR HBA bursts with full polarimetric data available, we have measured the linear and circular polarization fractions (Figure 4). We find a broad similarity to the polarimetric properties measured for FRB 20180916B bursts at 350 MHz (Chawla et al. 2020) and 1.7 GHz (Nimmo et al. 2021); i.e., the high linear polarization fraction, with flat PA during the burst, and negligible circular polarization fraction persist from 130 to 1700 MHz, almost four octaves in radio frequency. For the low-frequency burst CV4, the lower linear polarization fraction is likely due to scattering, though investigations are ongoing to determine whether this could be the result of Faraday conversion (Gruzinov & Levin 2019; Vedantham & Ravi 2019). The remarkably similar polarimetric burst profiles across a wide range of radio frequencies are dissimilar to what is seen in some pulsars (Noutsos et al. 2015) and provide a novel constraint on FRB emission theory.

Despite the low frequencies of our LOFAR HBA burst detections, the large burst widths and likely presence of scattering make it difficult to precisely and accurately determine the burst DMs in order to investigate potential DM variability. For this reason, we have dedispersed all bursts in Figure 2 to a best-fit literature value of DM = 348.772 pc cm⁻³ (Nimmo et al. 2021). Previously, PR3 searched for DM variability, also as a function of FRB 20180916B's activity phase. Using structure-optimized DMs (Hessels et al. 2019) from four bursts detected using the CHIME/FRB baseband capture system, they found no DM variations with magnitude ≳0.1 pc cm⁻³.

Using three high-S/N bursts detected at 550–750 MHz with uGMRT and searching for a DM that maximizes burst substructure (Hessels et al. 2019), Marthi et al. (2020) also found no evidence for large DM variations. Their strongest constraint, DM = 348.8 ± 0.1 pc cm⁻³ (for their burst 11), is also consistent with the DM = 348.772 pc cm⁻³ (Nimmo et al. 2021) value we use throughout this paper. Optimizing the substructure of burst G5 (Figure 7), we find DM = 349.5(1) pc cm⁻³. However, as we noted in the previous subsection, there is some ambiguity in the interpretation of the burst structure of G5. Nevertheless, any DM variations are constrained to be ≲1 pc cm⁻³.

The LOFAR HBA burst detections do provide precise RM values. Comparing with previous GBT and CHIME/FRB measurements (CHIME/FRB Collaboration et al. 2019a; Chawla et al. 2020) and correcting for the variable RM contribution from the Earth's ionosphere, we find that FRB 20180916B shows significant variations at the level of ∼2–3 rad m⁻² (∼2% fractional; though the variable component of the RM may be much lower than the total line-of-sight RM, in which case the fractional variations could be much larger). With only seven measurements in hand, it is not yet clear whether the observed RM is correlated with the activity phase of the source or varies stochastically (Figure 5). Alternatively, the RM (and DM) could also depend on radio frequency, e.g., if there is a frequency dependence to the location of the emission region. Some pulsars are known to show variable apparent RM across their pulse profiles that indicates emission from multiple locations within the magnetosphere (e.g., Dai et al. 2015; Ilie et al. 2019). Similar effects might be seen in FRBs if their emission is magnetospheric in origin as well.

A revised analysis of the 38 previously published and 17 newly presented CHIME/FRB bursts—spanning 2 yr, and thus 45 activity cycles, from 2018 September 16 to 2020 September 19—provides a refined source activity period of P_act = 16.33 ± 0.12 days using the same methodology as in PR3. This is consistent with the previously published determination of P_act = 16.35 ± 0.15 days. The refined activity window is 5.2 days, as compared to the 5.4 day window that was measured using only 38 events.

Folding all bursts modulo P_act, we see that the 18 LOFAR HBA bursts are systematically delayed in activity phase by ∼3 days (0.2 cycles), compared to the 55 CHIME/FRB bursts (Figure 9). As we show in Figure 10, this effect is not simply a reflection of the observational exposure. LOFAR HBA bursts have been detected in four activity cycles, namely Cycles 20, 21, 35, and 37. Though the number of observations and bursts per cycle is low, these are consistent with the delayed activity being a time-invariant effect. Cycle 37 provides the best single-cycle observational coverage and burst sample between CHIME/FRB and LOFAR (Figure 1) and is also consistent with the overall picture one obtains by summing over all activity cycles. An ongoing campaign of NenuFAR, LOFAR HBA, and Effelsberg observations will better characterize this frequency-dependent activity and determine whether the low-frequency activity window has a larger duty cycle or not. A larger burst sample can also, in principle, determine the functional form of the activity delay with frequency, e.g., whether the delay as a function of frequency is linear or quadratic.

Higher-frequency detections of FRB 20180916B at the L band seem to preferentially arrive at the start of the CHIME/FRB-derived activity window (PR3; Aggarwal et al. 2020). In follow-up observations of repeaters with periodic activity, it thus seems wise to cover a broad window and not only the peak of activity derived at a different frequency.

We now consider what models can naturally accommodate the wealth of observational facts available for FRB 20180916B.

Recently, Tendulkar et al. (2021) used HST and Gran Telescopio Canarias observations to demonstrate that FRB 20180916B is significantly offset with respect to nearby star-forming regions in its host galaxy. Assuming it is a neutron star formed in one of these regions, as opposed to in situ, it is more likely to be an old (100 kyr–10 Myr) source. Coupled with the observed 16.3 day activity period, this led Tendulkar et al. (2021) to suggest that the system could be an HMXB, where the bursts are possibly generated through an interaction between the companion wind and neutron star magnetosphere (sometimes called a "cosmic comb" model ³¹ ; e.g., Zhang 2017; Ioka & Zhang 2020). Furthermore, the detection of a wide range (factor of ∼1000) of emission timescales and microsecond structure in some FRB 20180916B bursts (Nimmo et al. 2021) is more naturally explained in terms of magnetospheric emission models, as opposed to those that invoke emission in a relativistic shock far from the neutron star.

The HMXBs are relatively common; there are roughly 200 known in the Milky Way (Coleiro & Chaty 2013; Walter et al. 2015). Furthermore, studies find that binary interaction dominates the evolution of massive stars (Sana et al. 2012); hence, highly magnetized neutron stars in such binaries should not be particularly rare. In fact, HMXBs are arguably too abundant to explain repeating FRBs, unless one invokes particular evolutionary stages, sporadic emission episodes, and/or viewing geometry to explain the lack of radio burst detections from Galactic HMXBs, as we discuss further below.

The detection of low-frequency bursts with no simultaneous high-frequency emission (Pearlman et al. 2020) strongly challenges models in which the 16.3 day observed periodicity is the result of absorption by the companion wind (Lyutikov et al. 2020). Rather, our LOFAR HBA detections and the lack of any observed DM variations, with ΔDM ≲ 0.1 pc cm⁻³ throughout the active window (PR3), suggest that we have a relatively clean line of sight to the burst source itself. Nonetheless, the observed RM variations of FRB 20180916B are atypical and not seen from isolated Galactic pulsars; they suggest that it is in a special local environment.

It is interesting to compare FRB 20180916B with known radio pulsar/OB star binary systems. One example is the Small Magellanic Cloud (SMC) radio pulsar/B1V companion binary system PSR J0045−7319. In this 51 day highly eccentric binary, at periastron and apastron, the pulsar approaches to within four and 34 B-star radii from the companion, respectively, probing very different regions of the stellar wind. Yet Kaspi et al. (1996) found an upper limit on DM variations of <0.9 pc cm⁻³, yielding a strong constraint on the stellar wind, <10⁻¹¹ M_⊙ yr⁻¹. Given the low metallicity of the SMC, ³² a B star with such a weak wind is not unexpected.

In a similar way, the lack of large observed DM variations (i.e., ΔDM ≲ 1 pc cm⁻³) from FRB 20180916B may not be problematic for a massive star binary model. Based on population synthesis arguments, Zhang & Gao (2020) argued that B-type stars would be the most likely companions to FRB binary sources. If the putative FRB 20180916B orbit is not very eccentric, the variation in pulsar/companion distance with orbital phase may not be large. This, together with a relatively weak companion wind, may make DM variations hard to observe. On the other hand, the Galactic eccentric radio pulsar/OB star binary PSR J1740−3052 does show DM variations of order ∼2 pc cm⁻³ near periastron in its 231 day orbit; these imply a mass-loss rate on the order of 10⁻⁹ M_⊙ yr⁻¹, still low for a Galactic O or early B star (Madsen et al. 2012). For a similar strength wind, given the 16.3 day orbital period and hence closer pulsar/companion distance, the CHIME/FRB-active window would have to occur near apastron, or the orbit would have to be fairly circular to avoid detectable DM changes. Similarly, the Galactic highly eccentric 1237 day radio pulsar/Be star binary PSR B1259−63 shows large DM, RM, and scattering-time variations in the ∼50 days near periastron (ΔDM ≃ 6–8 pc cm⁻³, ΔRM ≃ 6000 rad m⁻², and increased scattering such that 1.5 GHz pulse profiles become unobservable for this 48 ms pulsar), though these are likely due to the presence of a circumstellar disk through which the pulsar passes (Johnston et al. 1992, 1996, 2005). Thus, an emission-line star (e.g., a Be star) in the putative FRB 20180916B binary system seems implausible, because the FRB 20180916B system has a shorter orbital period, hence more compact orbit size and pulsar/star separation, together with the constraints on DM variations, absence of sizable RM variations, and given our new low-frequency detections which rule out large scattering times and are offset in phase.

In the context of interacting binary models (Ioka & Zhang 2020; Lyutikov et al. 2020; Popov 2020; Du et al. 2021), the FRBs may be produced by a highly magnetized neutron star whose magnetosphere is "combed" by the ionized wind of a massive companion star. Such interaction could lead to magnetic reconnection events, which have been proposed as a source of FRBs (Lyutikov & Popov 2020). The bursts may only be visible within a funnel where the neutron star's wind shields against the companion's wind, which is otherwise opaque to induced Compton or Raman scattering for FRB emission (Ioka & Zhang 2020; see their Figure 1). Windows of observable burst activity, lasting for ∼4–5 days in the case of FRB 20180916B, then correspond to when this funnel and the induced magnetic tail of the neutron star are pointed toward Earth. This special viewing geometry could also, in principle, explain why Galactic HMXBs are not known to be prolific sources of bright radio bursts.

The funnel and magnetic tail can also be swept back by orbital motion. Wang et al. (2019) and Lyutikov (2020) discussed a radius-to-frequency mapping model to explain the "sad trombone" effect seen from repeating FRBs on timescales of milliseconds (Hessels et al. 2019). If the radio emission frequency of bursts indeed scales with distance from the central source, where the magnetic field strength and plasma density are lower at larger distances, then this could plausibly also explain why the LOFAR HBA bursts we have observed are delayed in activity phase with respect to those seen at higher frequencies by CHIME/FRB. In such a toy model, the LOFAR bursts would originate at larger distances from the neutron star in a swept-back magnetic tail.

While a model in which FRB 20180916B is a highly magnetized neutron star in an interacting HMXB system can plausibly explain all of the observed phenomena to date, several authors have argued that the 16.3 day activity period of FRB 20180916B is the rotational period (Beniamini et al. 2020), or precession period (Levin et al. 2020; Sob'yanin 2020; Yang & Zou 2020; Zanazzi & Lai 2020), of an isolated magnetar. Tendulkar et al. (2021) argued against a young source, and the consistent polarization position angle between bursts (Nimmo et al. 2021) sets significant constraints on precession models. Nonetheless, in the context of these nonbinary models, the frequency dependence of the observed activity could also be interpreted as a radius-to-frequency mapping effect in an emission cone that is slowly sweeping past the line of sight and is perhaps swept backward at higher altitudes.

Lastly, we note that FRB 20180916B and FRB 20121102A, the two best-studied repeaters, share remarkably similar phenomenology (see Section 1), despite the fact that FRB 20121102A is hosted in a much less massive dwarf galaxy and is coincident with a compact, persistent radio source. They are almost certainly of a similar physical origin, though FRB 20121102A may be in the vicinity of an accreting massive black hole (Michilli et al. 2018). Such a Galactic center–like environment may explain why, unlike FRB 20180916B, FRB 20121102A has only been detected once at radio frequencies <1 GHz (Josephy et al. 2019). Their different Galactic latitudes (b = 37 for FRB 20180916B and b = −02 for FRB 20121102A) might also play some role. Perhaps, in the context of an HMXB model, if the radio emission frequency is tied to the instantaneous altitude of the emission site, it is possible that the characteristic emission height is tied to the companion wind strength.