EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 634 (February 2020) A&A, 634 (2020) A87 Full HTML
Open Access
Issue
A&A
Volume 634, February 2020
Article Number A87
Number of page(s) 11
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/201936564
Published online 12 February 2020

© R. Lico et al. 2020

Licence Creative Commons
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Open Access funding provided by Max Planck Society.

1. Introduction

Blazars are the most extreme objects in the family of active galactic nuclei (AGNs). With their jets closely aligned to our line of sight, they represent the best target for the study of both the physics of particle acceleration and the role played by magnetic fields in these extreme plasma environments (e.g. Blandford & Königl 1979; Marscher et al. 2008). Blazars include both flat spectrum radio quasars (FSRQs – with an optical spectrum dominated by strong emission lines) and BL Lac objects (BL Lacs – showing almost featureless spectra with occasional weak optical emission lines, Stickel et al. 1991). Blazar emission spans the entire electromagnetic spectrum from radio up to TeV γ-ray energies, and is mostly non-thermal in nature and linearly polarised, providing us with important insights into the magnetic field structure (e.g. Gabuzda et al. 1994; Gómez et al. 2008; Orienti et al. 2011; Hovatta et al. 2012; Casadio et al. 2019). Blazars tend to show erratic variability across the electromagnetic spectrum on different timescales, ranging from years (Ulrich et al. 1997) down to a few minutes (e.g. Aharonian et al. 2007; Albert et al. 2007a.) A small subset of blazar sources exhibit possible quasi-periodic variability across the radio, optical, and X-ray emission bands (e.g. Villata et al. 2004; Hovatta et al. 2008; Valtonen et al. 2008; Wiita 2011; King et al. 2013). Over the last few years, thanks to the Fermi Large Area Telescope (Fermi-LAT) continuous monitoring of the sky in the MeV-GeV energy range (Atwood et al. 2009), quasi-periodic variability has been investigated at γ-ray energies in a number of blazars (e.g. Prokhorov & Moraghan 2017). Such high-energy periodicity might be related to jet precession and/or modulation of the accretion rate onto the central engine(s). As such, the study of quasi-periodic variability at γ-ray energies can help shed light on fundamental issues such as the disc–jet connection and the nature of the jet’s magnetic fields, and could provide further insight into gravitational wave production in binary super massive black hole (SMBH) systems (Abbott et al. 2016).

In this context, the BL Lac object PG 1553+113 has exhibited complex high-energy variability and has been detected at MeV/GeV γ-ray energies by the Fermi-LAT at a significance level above 10 σ (Abdo et al. 2009). PG 1553+113 is the first blazar for which a quasi-periodic pattern with a period of ∼2.18 ± 0.08 yr (Ackermann et al. 2015) has been observed in its γ-ray light curve, providing us with a unique opportunity to investigate this peculiar behaviour.

PG 1553+113 exhibits an almost featureless optical spectrum (Miller & Green 1983; Falomo & Treves 1990) and the emission at all wavelengths, from radio to very high energies (VHE; E > 100 GeV), can be attributed to the non-thermal jet emission, without significant contributions from other substructures (e.g. the disc, the corona, or the broad line region). The host galaxy is yet undetected and the redshift determination is still uncertain. Using a method based on Bayesian statistics and the extragalactic background light absorption effects on the VHE spectrum, Abramowski et al. (2015) obtained a value of z = 0.49 ± 0.04. Additional indirect methods, such as the Lyα forest method and the absence of a break in the VHE spectrum, constrain the redshift to the range 0.3–0.6 (e.g. Qin et al. 2018; Landoni et al. 2014; Danforth et al. 2010; Prandini et al. 2010; Mazin & Goebel 2007). This source is classified as a high synchrotron peaked (HSP) blazar (Giommi et al. 1995) and was detected by the HESS and MAGIC γ-ray Cherenkov telescopes in 2005 (Aharonian et al. 2006; Albert et al. 2007b).

The quasi-periodic trend found in the γ-ray light curve of PG 1553+113, although still debated (e.g. Covino et al. 2019), is matched by similar behaviour at lower frequencies. In addition to the periodicity in the γ-ray emission, hints of periodicities have been found in the optical emission light curve spanning ∼10 yr (2004–2014) obtained within the Tuorla optical monitoring program1 and in the 0.3–10 keV light curve obtained by XRT on board the Neil Gehrels Swift satellite, although at a lower significance level due to the sparse observations (Tavani et al. 2018).

To investigate the radio properties of the innermost regions of the jet, where the γ rays are thought to be produced, high-resolution very long baseline interferometry (VLBI) observations are needed. On VLBI scales, the source shows a compact nuclear region with jet structure in the northeast direction (position angle ∼50°), that is well detected out to roughly 20 pc from the core region (Rector et al. 2003). Beyond this distance the jet becomes very faint and diffuse. Although several multi-wavelength (MWL) observing campaigns of PG 1553+113 (e.g. Raiteri et al. 2015, 2017; Hovatta et al. 2016; Itoh et al. 2016) were conducted in recent years, only sparse VLBI observations of PG 1553+113 can be found in the literature, including the 15 GHz observations within the MOJAVE database2 and the TeV Blazars very long baseline array (VLBA) monitoring program led by B. G. Piner and P. G. Edwards3. For these reasons, we conducted a systematic multi-frequency (15, 24, and 43 GHz) VLBA monitoring of PG 1553+113 covering a period between two consecutive maxima in the γ-ray activity during 2015–2017. The main purpose of this observing campaign is to investigate and characterise the parsec-scale source properties and their evolution with time as well as the possible connections with the observed flux density modulation.

In this paper we describe the observing campaign and the data analysis in Sect. 2, while the main results are presented in Sect. 3. A general discussion and summary of our findings are presented in Sects. 4 and 5, respectively. Throughout the paper we use a ΛCDM cosmology with h = 0.71, Ωm = 0.27, and ΩΛ = 0.73 (Komatsu et al. 2011). The spectral index α is defined as Sν ∝ ν−α, and all angles are measured from north to east. At a redshift of z = 0.49 (Abramowski et al. 2015) 1 mas corresponds to ∼5.2 pc.

2. Observations and data analysis

We monitored PG 1553+113 with the VLBA at 15, 24, and 43 GHz in full polarisation from February 2015 to June 2017. The observations are based on two related observing projects: two six-hour observing sessions during the period of its maximum γ-ray activity in 2015 (BL214) and a full two-year monitoring program with roughly bimonthly runs (BL215, 54-h in total). Table 1 reports the log of the observations. In some epochs one or more VLBA stations were missing or flagged due to technical problems (as reported in the last column in Table 1); when a station is missing only for a specific observing band, it is indicated in brackets. We were able to obtain absolute electric vector position angle (EVPA) orientations for the four observing epochs (MJD 57055, 57078, 57650, 57755) for which we had quasi-simultaneous Karl G. Jansky very large array (JVLA) observations at similar frequencies. We used the source J1310+3220 as the instrumental polarisation calibrator. We also made use of the 15 GHz observations provided by the Owens Valley Radio Observatory (OVRO) blazar monitoring program (Richards et al. 2011). PG 1553+113 has been monitored with OVRO since 2008 with a cadence of about two observations per week. In this work we use the observations over the MJD range 54696-581534.

Table 1.

Details about the different VLBA observing sessions.

We used the software package Astronomical Image Processing System (AIPS; Greisen 2003) for the data calibration, the fringe-fitting, and the detection of cross-polarised fringes. We determined the instrumental polarisation leakages (D-terms) with the AIPS task LPCAL. For the final images we used the CLEAN and self-calibration procedures in the DIFMAP software package (Shepherd 1997). The core flux densities were obtained by fitting elliptical Gaussian components to the calibrated visibilities for each epoch at each frequency in DIFMAP. In this work the core is identified as a bright and stationary feature at the upper end of the jet, with a roughly flat spectrum and high brightness temperature (of the order of 1010 K, see Sect. 3.5).

For producing the polarisation images presented in this work we made use of IDL routines developed by J. L. Gómez and for the ridge-line analysis the HEADACHE5 python package developed by J. Liu was used.

3. Results

3.1. Images

In Fig. 1 we show the natural weighted images for PG 1553+113 obtained during our 2015–2017 VLBA observing campaign at 15 GHz (left panel), 24 GHz (middle panel), and 43 GHz (right panel). The source shows a compact core region with a northeast jet structure detected up to ∼1.5 mas (∼7.8 pc in linear size) at 15 GHz and ∼1 mas (∼5.2 pc) at 43 GHz. This morphology is in agreement with 15 GHz MOJAVE images, as well as with the 22 and 43 GHz results by Piner & Edwards (2018). The jet position angle (PA) has an average value of ∼50°. Further details and constrains about the PA across the different observing epochs are provided in Sect. 3.3. Colours represent the linear polarised emission, which is mainly detected in the core region, while the black contours represent the total intensity emission.

thumbnail Fig. 1.

Natural weighted 15 GHz (left panel), 24 GHz (middle panel), and 43 GHz (right panel) VLBA images of PG 1553+113 from the first observing epoch (top image in each panel) to the fourteenth epoch (bottom image in each panel). The vertical separation between images is not proportional to the time difference between epochs. The images at each frequency have been convolved with a common beam with FWHM 0.6 mas × 1.2 mas, 0.4 mas × 0.8 mas and 0.3 mas × 0.7 mas at 15, 24 and 43 GHz, respectively, as shown in the bottom-left corner in each panel. The overlaid lowest total intensity contour is at 0.4%, 1.3%, and 1.8% of the peak at 15, 24, and 43 GHz (see Table 2), respectively, with the following contours a factor of two higher. The colour scale represents the linearly polarised intensity.
Open with DEXTER

3.2. Light curves

In Fig. 4 we show the core region total (upper image) and polarised (lower image) intensity light curves. In each image results at 15, 24, and 43 GHz are reported in single frames from top to the bottom, respectively.

We monitored PG 1553+113 between February 2015 to June 2017 with a cadence of about 2 months, and with an approximately six-month gap between the second and third observing runs. Flux densities are reported in Table 2. We find two maxima in the total intensity emission (upper image in Fig. 4) around MJD 57312 (Oct 2015) and MJD 57916 (June 2017), that is, at the beginning and the end of the observing period, and a minimum around 57558 (Jun 2016). The total intensity light curves are found to vary in the range 120–200 mJy at 15 GHz, 100–170 mJy at 24 GHz, and 70–170 mJy at 43 GHz. The variability is higher at 43 GHz, where the fractional variability amplitude (Fvar, Vaughan et al. 2003) is 0.26 ± 0.03, while at 24 and 15 GHz Fvar is 0.16 ± 0.03 and 0.13 ± 0.03, respectively. Red triangles in the first frame in the top image of Fig. 4 represent the 15 GHz OVRO light curve.

Table 2.

Summary of the total and polarised intensity parameters shown in Fig. 4 and in the lower frame of Fig. 6.

As shown in the lower image in Fig. 4 we detect polarised emission from the source core region with an average value of 2.6, 2.7, and 2.5 mJy at 15, 24, and 43 GHz, respectively. The polarised emission is also found to be variable, with Fvar ≃ 0.40 ± 0.05 at all three observing frequencies. During MJD 57276 a peak in the polarised emission of 5.6 ± 0.7, 5.5 ± 0.6, and 4.8 ± 0.6 mJy is detected at 15, 24, and 43 GHz, respectively. The fractional polarisation varies in the range 0.7−3.1% at 15 GHz, 0.6−3.5% at 24 GHz, and 1.1−3.7% at 43 GHz. The average fractional polarisation is 1.7%, 2.0%, and 2.2% at 15, 24, and 43 GHz, respectively.

To asses the significance of the variability observed both in the total and polarised intensity light curves, we performed a χ2 analysis by testing the null hypothesis of non-variability (with the flux density being constant around the average value). In all cases the null hypothesis can be rejected with a confidence level above 3σ.

3.3. Jet total intensity ridge line

It is apparent upon visual inspection of the source images (Fig. 1) that the jet direction varies across epochs. Given that the extended emission of the source on VLBI scales is not prominent and well-defined, it is not straightforward to represent the jet brightness distribution by means of Gaussian model-fit components. In order to constrain the jet direction and to further quantify its variation with time, we calculate the total intensity ridge line in each epoch at 15 GHz after convolving the images with a circular beam with a 0.6 mas radius. We use the method proposed by Pushkarev et al. (2017) by taking azimuthal slices along the jet in a polar coordinate system centred on the core. We look for the half-point that divides the intensity integrated along the arc into two equal sections. The slices along the jet direction are taken in the image plane in steps of 0.03 mas and only pixels with a signal-to-noise ratio of S/N > 10 are taken into account. This approach is particularly suitable in the case of a jet showing a limb-brightened structure in total intensity, that is a distinctive feature of HSP blazars (e.g. Giroletti et al. 2008; Piner et al. 2010; Lico et al. 2014). As an example, in the left panel of Fig. 2 we present the 15 GHz total intensity image obtained during MJD 57276, convolved with a 0.6 mas × 0.6 mas circular beam, with the ridge line (yellow line) overlaid. The limb brightened structure is clearly visible. The same image is shown in polar coordinates in the right panel of Fig. 2, where the total intensity limb-brightened structure is also clearly visible (blue colour).

thumbnail Fig. 2.

Evolution with time of the total intensity (upper image) and polarised flux density (lower image) in the core region of PG 1553+113 during the MJD 57055–57916 period. In each image we report the 15 GHz (upper frame), 24 GHz (middle frame), and 43 GHz (lower frame) results. The overlaid red triangles in the upper image (top panel) represent the 15 GHz OVRO light curve in the MJD range 56974–57990.
Open with DEXTER

In the upper panel of Fig. 3 we present the ridge lines obtained for all of the 14 observing epochs in a single plot with different colours, where the jet wandering can clearly be seen. Overall, the ridge line extends from 0.9 to ∼1 mas from the core region. The jet direction in each epoch is determined by averaging all the angles of each ridge-line point between 0.6 mas from the core region (i.e. beyond the beam radius) and the outermost ridge-line point (≲1 mas). All angles thus obtained are shown in the lower panel of Fig. 3 and reported in Table 3. In principle we could also extract information about the ridge-line curvature. However, given that the jet of PG 1553+113 between 15 and 43 GHz is not particularly extended and is not well defined (core-dominated source), any estimates regarding the curvature would be largely uncertain.

thumbnail Fig. 3.

Left image: 15 GHz total intensity (contours and colour scale) image of PG 1553+113 during the third observing epoch (MJD 57276), convolved with a 0:6 mas 0:6 mas circular beam. The overlaid lowest total intensity contour is at 0.6% of the peak, with the following contours a factor of two higher. The ridge line is overlaid as a yellow line. Right image: same 15 GHz image in polar coordinates. The colour scale represents the total intensity emission and the white line represents the ridge line.
Open with DEXTER

Table 3.

Jet PA values with standard errors (σPA), as reported in the bottom panel in Fig. 3, and core region spectral indexes between 15 and 43 GHz (SI15−43) with the uncertainties σSI15−43.

To estimate the angle of the funnel in which the jet wobbles, we first align the 15 GHz images (convolved with a common beam with FWHM 0.6 mas × 1.2 mas) according to the position of the fitted core component, and after subtracting the core component emission we stack all of the residual images. Only those pixels with a S/N > 8 are taken into account. In the residual stacked image (Fig. 5), the jet-emitting region extends up to ∼1 mas (∼5.2 pc) in a funnel with an angle of ϕ ∼ 100 deg. A total intensity limb-brightened structure is clearly visible, with the southeast limb being brighter (peak flux density ∼3.7 mJy beam−1) than the northern one (peak flux density ∼2.5 mJy beam−1).

3.4. Intrinsic polarisation angle and rotation measure

Due to the lack of polarisation calibrators on VLBI scales, in order to obtain the absolute orientation of the EVPAs, a comparison with quasi-simultaneous single-dish or JVLA observations is required. The final EVPA values for the four epochs for which close in time JVLA observations were available (see Sect. 2) are shown in the bottom panel of Fig. 6 for 15 (black circles), 24 (red triangles), and 43 (blue stars) GHz. The EVPAs are found to be variable across the different observing epochs, and change from being roughly aligned (∼50°) to approximately transverse (∼150°) to the jet axis.

Because of the effect of Faraday rotation, which occurs when a polarised wave propagates through a magnetised plasma, the observed polarisation angle (χobs), at a given observing wavelength λ is rotated with respect to the intrinsic angle (χint). Taking this effect into account is essential for obtaining the intrinsic polarisation angle as well as the intrinsic magnetic field orientation. For a Faraday rotating medium external to the emitting region χobs = χint + RM ×  λ2, with RM being the rotation measure defined as:

(1)

where ne is the electron density (cm−3), B∥ is the component of the magnetic field parallel to our line of sight (mG), and dl the path length (pc). It is evident, within this theoretical framework, that a linear dependence exists between χobs and λ2. We therefore compute linear fits to our observed EVPAs at 15, 24, and 43 GHz to obtain estimates of χint and RM.

The values of RM and χint obtained for these four observing epochs are presented in Table 4 and in Fig. 6 (top and middle panels, respectively). The average RM value is ∼−1.0 ± 0.4 × 104 rad m−2, varying between −1.3 and −0.8 ×  104 rad m−2. The intrinsic polarisation angle is also variable across the different observing epochs, varying between 135° ± 7° (i.e. roughly transverse to the jet axis) and 210° ± 7° (i.e. roughly parallel to the jet axis).

Table 4.

Rotation measure and intrinsic polarisation angle values for the core region obtained from the linear fits of EVPAs vs. λ2.

3.5. Brightness temperature measurements

By fitting the core brightness distribution in the uv-plane with an elliptical Gaussian component, via the modelfit procedure in DIFMAP, we can determine the observed rest-frame core brightness temperature (e.g. Tingay et al. 1998):

(2)

where Score corresponds to the fitted core flux density (Jy) for a given observing frequency ν (GHz), θmaj and θmin are the FWHM (mas) of the major and minor axes of the elliptical Gaussian core component, and z is the redshift. The resulting values are reported in Table 5, where , and represent the maximum, the minimum, and average values at each frequency across the different observing epochs.

Table 5.

Core region brightness temperature values in units of 1010 K. For more details see Sect. 3.5.

We also estimate the rest-frame core variability brightness temperature using a variability timescale of τ ∼ 300 days, which is the average time range between the maximum and minimum values observed in the core total intensity light curve:

(3)

with ΔSmax being the difference between the maximum and the minimum core flux density values and dL (m) the luminosity distance (e.g. Liodakis et al. 2017). The values at the three observing frequencies are reported in Table 5 (fourth line).

Since blazar jets are closely aligned to our line of sight, the effects of Doppler boosting must be taken into account when investigating the intrinsic source properties (Urry & Padovani 1995). and are related to the intrinsic brightness temperature via the following equations:

(4)

(5)

where δ = [γ(1−β cos θ)]−1 is the Doppler factor (with θ being the viewing angle between the jet and our line of sight, β the jet speed in units of the speed of light, and γ = (1−β2)−1/2 the bulk Lorentz factor).

4. Discussion

In recent years, several multi-frequency observing campaigns have been devoted to the study of the quasi-periodic variations (on timescales of ∼2 yr) detected in the γ-ray light curves of the PG 1553+113 (Ackermann et al. 2015). While optical periodic variations on different timescales have been extensively investigated in blazars (e.g. Valtonen et al. 2006; Li et al. 2009; Britzen et al. 2018), at γ-ray energies this has been possible only with the advent of the Fermi-LAT in 2008 thanks to the continuous and systematic sky monitoring in the MeV-GeV energy range. A similar approximately two-year quasi-periodicity (with significance above >3σ) has also been claimed in a few other blazars such as PKS 0537-441 (Sandrinelli et al. 2016), BL Lacertae (Sandrinelli et al. 2017), PKS 2155-304 (Zhang et al. 2017a), PKS 0301-243 (Zhang et al. 2017b), and J1043+2408 (Bhatta 2018).

However, while some dedicated studies confirm the approximately two-year periodicity in the γ-ray light curve of PG 1553+113 (e.g. Ackermann et al. 2015; Prokhorov & Moraghan 2017; Ait Benkhali et al. 2020; Yan et al. 2018), there are other studies which demonstrate that there is no solid evidence supporting such hints of year-long periodicities in blazars (e.g. Covino et al. 2019). We emphasise that caution is required when claiming quasi-periodicities on timescales of a few years in blazar γ-ray emission. This is because of the limited γ-ray light-curve duration (based on ∼10 yr of LAT operational activity), and the red-noise contamination that can easily mimic such short-period variability patterns (e.g. Covino et al. 2019).

One approach to further substantiate the existence of γ-ray periodicity is the use of comparisons to other wavelengths. In the case of PG 1553+113, the optical light curve shows a variability trend in good agreement with the γ rays, confirming the approximately two-year periodic pattern. Moreover, hints of possible periodic behaviour are also found at X-rays, although this is not statistically significant due to the very sparse temporal sampling (Ackermann et al. 2015; Tavani et al. 2018). In this work we extend the investigation to the radio emission from the innermost jet region by means of a dedicated multi-frequency VLBA observing campaign. It is apparent from the 15, 24, and 43 GHz VLBA light curves shown in Fig. 4 that two periods of enhanced radio emission occur around October 2015 (MJD 57312) and June 2017 (MJD 57916), with a minimum appearing around June 2016 (MJD 57558). We note that the spectral index flattens during the periods of enhanced activity, with α15−43 ∼ 0.2 ± 0.1, while it is steeper during the low-activity state, with α15−43 ∼ 0.5 ± 0.1 (see Table 3). No clear repeating patterns are discernible in the light curve during the period 2015–2017. Furthermore, the detection of multiple cycles on a longer observing period would be necessary to claim quasi-periodicity in the parsec-scale radio emission.

thumbnail Fig. 4.

Upper image: ridge lines of PG 1553+113 for all 14 observing epochs at 15 GHz. Each epoch is represented by a different colour and is labelled with a progressive number. The two dashed arcs represent 0.6 mas and 1 mas from the core region. Lower image: jet direction across the observing epochs represented by the average ridge-line angle (see Sect. 3.3 for more details).
Open with DEXTER

We note that the 15 GHz VLBA flux density trend over the entire observing period is in good agreement with the light curve obtained from the OVRO monitoring (red triangles overlaid in the top panel of Fig. 4). Indeed, we find Fvar ∼ 0.13 ± 0.03 for both light curves and the Pearson correlation coefficient between our 15 GHz VLBA and the closest OVRO observations (in time) is r = 0.90, implying a correlation above a 3σ level. These findings allow us to argue that: (1) the overall variability in the radio emission is driven by changes occurring within the VLBI core region, and (2) that the OVRO light curve is representative of the core region flux density modulation. Based on these assumptions, it would appear that despite the presence of rapid variability in the OVRO light curve since 2008 (lower frame in Fig. 7) there is no clear and well-defined periodicity in the radio emission (in contrast to the well-defined periodicity in the γ-ray emission).

This is confirmed by applying the weighted wavelet z-transform (WWZ) method (Foster 1996), which is widely used for testing and quantifying quasi-periodic oscillations in blazar light curves (e.g. Hovatta et al. 2008; Ackermann et al. 2015; Bhatta 2017; Tavani et al. 2018; Ait Benkhali et al. 2020). We adopt a WWZ period search window in the frequency range 5.5−20.0 × 10−4 day−1, with a frequency step of 3 × 10−5 day−1, resulting in approximately 50 test frequencies. We choose the minimum frequency in such a way that at least two complete cycles could be detected throuth the full observing time range that covers ∼10 yr, while the maximum frequency corresponds to a timescale of 500 days (i.e. about 80 times the mean time separation). The wavelet window width is defined by using a decay constant c = 0.0125. In Fig. 7 we show the OVRO light curve WWZ power (colour scale in the top-right corner) as a function of time (x-axis) and period (y-axis) over the MJD range 54696-58153. From this analysis, no variability pattern is seen in the OVRO light curve over the entire time range, in contrast to that seen in the γ-ray light curve (Ackermann et al. 2015; Tavani et al. 2018). A possible quasi-periodic modulation is found over a limited MJD range (57300-58153) with a period of ∼930 days (∼2.5 yr).

We assess the statistical significance of the WWZ variability timescale by means of Monte Carlo simulations, following the approach proposed by (and references therein Emmanoulopoulos et al. 2013), which takes into account the so-called ‘red noise leakage’ and aliasing effects due to the sampling properties of the real data sets (i.e. finite length and finite time resolution). We first generate 2000 light curves with the same power spectrum density (PSD) and power density function (PDF) as the actual observations. We then calculate the WWZ periodograms for the simulated light curves using the same parameters as for the OVRO 15 GHz observations. The averaged WWZ is then compared with the one obtained from the observed light curve. The significance intervals are over-plotted as white contours in Fig. 7, in which the dotted, dashed, and solid lines represent 1σ, 2σ, and 3σ, respectively.

These results confirm the lack of a clear periodic variability or pattern that is stable over time at radio frequencies.

4.1. A wobbling jet

In general, quasi-periodic MWL emission modulation in blazars has one of two possible origins: (i) The first is a geometrical origin with periodic variations in the Doppler beaming factor produced by periodic changes of the angle between the jet axis and our line of sight. These changes in orientation may be related to the jet precessing or rotational motion, and/or helical structure within relativistic jets (e.g. Camenzind & Krockenberger 1992; Abraham 2000; Stirling et al. 2003; Nakamura & Meier 2004; Rieger 2004; Raiteri et al. 2015). (ii) The second is an intrinsic origin, with quasi-periodic plasma injection into the relativistic jet that is due to quasi-periodic instabilities in the accretion flow, producing the observed variability patterns in the emitted radiation (e.g. Honma et al. 1992; Tchekhovskoy et al. 2011; Pihajoki et al. 2013). In the case of a geometrical effect, the quasi-periodic variability is only observed in the frame of the observer, while in the case of pulsational instabilities in the accretion flow, such quasi-periodic emission variations are also present in the jet comoving frame. In both scenarios the presence of a binary SMBH system is often invoked (e.g. Cavaliere et al. 2017; Sobacchi et al. 2017; Tavani et al. 2018; Caproni et al. 2017; Britzen et al. 2018). As suggested by Begelman et al. (1980), binary SMBHs may produce wiggles, helical motion, and periodic variability in the jets of AGNs as a kinematical consequence of their orbital motion and of jet precession.

In order to test the jet precession scenario in PG 1553+113, we constrain the jet position angle across epochs by means of the total intensity ridge line, as described in Sect. 3.3. The jet position angle, as shown in Fig. 3, is found to be variable across the different epochs with a range of values between ∼40° and ∼60°. We assessed the significance of the PA variability by means of a χ2 analysis, and the null hypothesis of a constant PA equal to the average value of 50° could be rejected with a confidence level above 3σ. A direct signature of the PA variation across epochs is the large angle (∼100°) of the funnel in which the jets wander. This was determined via stacking the images of the 15 GHz residuals for all the observing epochs (Fig. 5). The wobbling jet motion indicates that geometrical effects can play a role in the observed emission variability through Doppler boosting modulation. However, we do not find any clear connection between the jet position angle variation and the total intensity emission or with the approximately two-year γ-ray quasi-periodic variability pattern. For this reason, under the assumption that γ rays and radio emission are produced in the same region, we argue that the Doppler boosting modulation alone cannot account for the observed recurrent oscillations in the γ-ray light curve. Additional physical mechanisms are required. Within the framework of a binary SMBH system, one plausible scenario is that the secondary black hole crosses and perturbs the accretion disc of the primary black hole, inducing a temporary enhancement of the accretion rate, which in turn leads to increased jet emission (Sillanpaa et al. 1988; Valtaoja et al. 2000; Caproni et al. 2017; Britzen et al. 2018). Another mechanism that could be responsible for the jet wobbling is the possible shuttle of the core closer to the jet apex, as reported in several studies (e.g. Kovalev et al. 2008; Niinuma et al. 2015). This core-shuttle effect can result from either changes in the physical conditions at the jet base or the ejection of a new component blended within the unresolved core region (e.g. Hodgson et al. 2017; Lisakov et al. 2017).

thumbnail Fig. 5.

15 GHz residual stacked image after core subtraction (see Sect. 3.3). The colour scale represents the total intensity emission. The lowest contour is three times the image noise (that is ∼0.34 mJy beam−1), with the following contours being a factor of two higher. The ϕ symbol represents the angle of the funnel in which the jet wobbles across time.
Open with DEXTER

thumbnail Fig. 6.

Core region RM (upper panel), intrinsic polarisation angle (middle panel), and absolute EVPA orientation (lower panel) values during the epochs MJD 57055, 57078, 57650, and 57755. In the lower panel 15, 24, and 43 GHz are represented by black circles, red triangles, and blue star symbols, respectively.
Open with DEXTER

thumbnail Fig. 7.

Upper frame: two-dimensional WWZ power (colour scale on the right side) distribution as a function of time (x-axis) and period (y-axis) of the 15 GHz OVRO light curve over the MJD range 54696-58153. The significance intervals are plotted as white contours, in which the dotted, dashed, and solid lines represent 1σ, 2σ, and 3σ, respectively. Lower frame: OVRO light curve over the MJD range 54696-58153.
Open with DEXTER

One possible explanation for the presence of quasi-periodicity in γ rays but not in radio could come from distinct emission regions between the two frequencies. Within this context, an additional scenario for producing enhanced γ-ray emission is described by Bosch-Ramon et al. (2012) by means of the possible dynamical interaction of a relativistic jet with matter clumps or the atmosphere of an evolved star in the near vicinity of a SMBH. As the authors of this work point out, in general a considerably amount of dust, stars, and gas tend to cluster in the central region of galaxies and the relativistic jets in AGNs are expected to frequently interact with such ambient material.

4.2. Polarisation properties

The polarised emission from the core region (lower image in Fig. 4) is on average more variable (Fvar ∼ 0.40 ± 0.05) than the total intensity emission (Fvar ∼ 0.20 ± 0.03), and no variability patterns can be easily recognised during this 2015–2017 VLBA observing campaign. The most notable features in the polarised emission light curves are represented by a ∼5.0 ± 0.6 mJy peak reached at all of the observing frequencies during the third observing epoch (MJD 57276), and a second one of ∼3.5 ± 0.4 mJy during the fourteenth observing epoch (MJD 57916). Both peaks in the polarised emission are detected during the two periods in which the source is undergoing enhanced activity in the total intensity emission. This behaviour may indicate that the processes producing the increased total intensity emission activity are also responsible for the enhanced polarised emission.

The Faraday corrected EVPAs at 15, 24, and 43 GHz, as well as the intrinsic polarisation angle, obtained from the EVPAs versus λ2 linear fits during MJD 57055, 57078, 57650, and 57755 have been found to vary between being roughly parallel to the jet axis to being roughly transverse, with no clear connection with the total intensity emission modulation. This indicates that the magnetic field configuration is also variable with time. Such behaviour, together with the low degree of polarisation in the core region (on average ∼2.5 mJy), could be explained by the presence of multiple polarised subcomponents blended within the beam. The net observed polarisation properties integrated over the VLBI core region result from the sum of the relative contributions of each subcomponent, whose properties vary with time.

The linear dependence found between the observed EVPAs with λ2 indicate that the Faraday screen is mostly external to the emitting region (e.g. Burn 1966). The RM has been found to vary with time in a range between −1.3 and −0.8 ×  104 rad m−2. Given that the jet position angle wobbles across the different observing epochs, we can argue that the RM variations are produced by variations of the path length (dl) through the external Faraday screen and/or variations in the electron density distribution ne (Stirling et al. 2003). The persistent negative sign in the observed RM across the different epochs, within the scenario in which the Faraday screen is mostly represented by the hot accretion flow or wind outflow, could be explained by a misalignment of the jet axis with respect to the wind axis (Park et al. 2019). If compared to Mrk 421 (e.g. Lico et al. 2014) and Mrk 501 (e.g. Hovatta et al. 2012), the RM value found for PG 1553+113 is the highest measured in a blazar of HSP type to date. However, sample sizes are still not large enough for a statistically significant characterisation of the average RM properties of the HSP blazar population. In this regard, the results of recent and ongoing MWL polarisation VLBA observations will be presented in a further publication.

4.3. Brightness temperature

By investigating the brightness temperature properties we can obtain important information regarding the physical conditions within the jet, namely the energy balance within the emitting region. We find that the core brightness temperature decreases as the frequency increases, with . According to theoretical models (e.g. Marscher 1995), there is a critical frequency νc below (above) which the brightness temperature increases (decreases) with frequency. Within this scenario, and according to our findings, νc should be <15 GHz. Furthermore, by investigating the brightness temperature in a large sample of radio jets in a range of frequencies between 2 and 86 GHz, Lee (2014) found on average νc ≈ 9 GHz. This decreasing trend of with frequency could indicate that within the two emitting regions probed by the 15 and 43 GHz VLBA observations there is an acceleration of the jet flow with increasing distance from the central engine (see also Melia & Konigl 1989; Marscher 1995; Lee et al. 2016; Nair et al. 2019).

The variability brightness temperature (fourth line in Table 5) for each observing frequency is higher than , as expected because of the different dependence on δ. Given the limitations of the current data sets, we note that for estimating the variability timescale is assumed to be of the order of the average time range between the maximum and minimum values observed in the core total intensity light curve. Even though this approximation is not particularly accurate, it is suitable for the aims of the current analysis. Using Eqs. (4) and (5), with at the peak flux density, we can estimate both and δ:

We obtain an average value for δ of ∼1.4, which is a typical value for HSP objects as obtained from VLBI observations (e.g. Giroletti et al. 2004; Lico et al. 2012; Piner & Edwards 2018). Such a low value for the Doppler factor further supports our conclusion that the Doppler boosting modulation alone cannot account for the observed variability.

The resulting is of the order of 1.5 ×  1010 K, well below the ∼5 ×  1011−1012 K inverse Compton limit (Kellermann & Pauliny-Toth 1969), and in agreement with the typical found in HSP blazars (e.g. Piner et al. 2010; Piner & Edwards 2014; Lico et al. 2016). The physical mechanism generally invoked for explaining the values in HSPs is energy equipartition between particles and the magnetic field. Invoking equipartition yields an upper limit of ∼5 ×  1010 K (Readhead 1994). The estimated from our 15, 24, and 43 GHz observations is on average below the equipartition limit, indicating that we are probing an emitting region where the magnetic field energy density (uB) is higher than the non-thermal particle energy density (up). Using equation (5d) from Readhead (1994) we find log(up/uB) ∼ −4.6. A similar physical scenario is found in the innermost regions of the radio galaxy M 87, as reported in Kim et al. (2018).

5. Summary and conclusions

The HSP blazar PG 1553+113 has been observed intensely since an approximately two-year quasi-periodic variability pattern was recognised in its γ-ray light curve (Ackermann et al. 2015). In this work we explored the parsec-scale radio properties of the source by means of a 15, 24, and 43 GHz VLBA observing campaign in total and polarised intensity during the period 2015–2017.

Two periods of enhanced activity around October 2015 and June 2017 emerge from the total intensity light curve with a minimum around June 2016. However, in contrast to the γ-ray emission, no hints of quasi-periodic variability are found in the VLBI emission or in the 15 GHz OVRO light curve over a period covering about 10 yr.

We detect polarised emission in the core region (with a polarisation percentage varying in the range ∼1−4%) that is variable across epochs with no clear recurrent patterns in the light curve. The polarisation angle has been found to be variable across epochs as well, but without any clear connection to the total intensity or polarised emission. We also find a variable RM in the core region, ranging between −1.3 and −0.8 ×  104 rad m−2.

The core brightness temperature is found to decrease with increasing frequency, in agreement with Lee et al. (2016), likely suggesting that an acceleration of the jet flow is occurring within the emitting regions probed by the 15 and 43 GHz VLBA observations (Melia & Konigl 1989; Marscher 1995). Using both and we estimate a Doppler factor of ∼1.4 and K, indicating that within the emitting region the total energy is dominated by the magnetic field.

By means of a total intensity ridge-line analysis we constrain the jet position angle across the different observing epochs. We find that the jet direction varies in range between ∼40° and ∼60°, indicating that geometric effects could play a role in the observed emission variability through Doppler boosting modulation. However, there is no direct and clear connection between the jet wobbling motion and either the radio flux density or the γ-ray variability pattern, and additional physical mechanisms are invoked. One possible mechanism responsible for the jet wobbling is the core-shuttle effect described by Hodgson et al. (2017) and Lisakov et al. (2017). Another possibility is the presence of a binary SMBH system at subparsec separation, inducing a precessing motion in the jet as well as perturbations in the accretion disc with a consequent modulation in the accretion rate (e.g. Ackermann et al. 2015; Caproni et al. 2017; Tavani et al. 2018). While evidence of binary SMBHs systems at kiloparsec scales have been found in a few tens of objects through direct detection of the two centres, parsec or subparsec systems, which are expected to form quickly according to evolution models, are more elusive. The detection of parsec- or subparsec-scale binary SMBH systems would also be an important assessment of the prediction of coalescence of SMBHs and of consequent gravitational wave production (Begelman et al. 1980; Bhatta 2018; Ait Benkhali et al. 2020). This goal could be achieved in the next few years thanks to the efforts of the Event Horizon Telescope project (Event Horizon Telescope Collaboration 2019), which could potentially allow us to spatially resolve sub-parsec binary SMBH systems.

With the present VLBA monitoring we provide new and valuable information for MWL studies of this peculiar object, furthering our understanding of the physical mechanisms that produce the observed periodicity in the γ-ray emission. This work is intended to be part of a wider and extensive MWL observing program with regular monitoring of the source since 2015 at different frequencies, the results of which will be presented in a dedicated paper (MAGIC Collaboration et al., in prep).

Acknowledgments

We thank the anonymous referee for carefully revising the paper and for the valuable and constructive comments that improved the paper. RL gratefully acknowledges the financial support and the kind hospitality from the Instituto de Astrofísica de Andalucia (IAA-CSIC) in Granada. We that Jae-Young Kim for reading the paper and for the fruitful discussion. JL acknowledges the following grants: National Key R&D Program of China (No. 2018YFA0404602), Light in China’s Western Region program (No. 2015-XBQN-B-01), National Natural Science Foundation of China (No. 11503071). This work is based on observations obtained through the BL214 and BL215 VLBA projects, which makes use of the Swinburne University of Technology software correlator, developed as part of the Australian Major National Research Facilities Programme and operated under license (Deller et al. 2011). The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. This research has made use of data from the OVRO 40-m monitoring program (Richards et al. 2011) which is supported in part by NASA grants NNX08AW31G, NNX11A043G, and NNX14AQ89G and NSF grants AST-0808050 and AST-1109911.

References

All Tables

Table 1.

Details about the different VLBA observing sessions.

In the text
Table 2.

Summary of the total and polarised intensity parameters shown in Fig. 4 and in the lower frame of Fig. 6.

In the text
Table 3.

Jet PA values with standard errors (σPA), as reported in the bottom panel in Fig. 3, and core region spectral indexes between 15 and 43 GHz (SI15−43) with the uncertainties σSI15−43.

In the text
Table 4.

Rotation measure and intrinsic polarisation angle values for the core region obtained from the linear fits of EVPAs vs. λ2.

In the text
Table 5.

Core region brightness temperature values in units of 1010 K. For more details see Sect. 3.5.

In the text

All Figures

thumbnail Fig. 1.

Natural weighted 15 GHz (left panel), 24 GHz (middle panel), and 43 GHz (right panel) VLBA images of PG 1553+113 from the first observing epoch (top image in each panel) to the fourteenth epoch (bottom image in each panel). The vertical separation between images is not proportional to the time difference between epochs. The images at each frequency have been convolved with a common beam with FWHM 0.6 mas × 1.2 mas, 0.4 mas × 0.8 mas and 0.3 mas × 0.7 mas at 15, 24 and 43 GHz, respectively, as shown in the bottom-left corner in each panel. The overlaid lowest total intensity contour is at 0.4%, 1.3%, and 1.8% of the peak at 15, 24, and 43 GHz (see Table 2), respectively, with the following contours a factor of two higher. The colour scale represents the linearly polarised intensity.
Open with DEXTER
In the text
thumbnail Fig. 2.

Evolution with time of the total intensity (upper image) and polarised flux density (lower image) in the core region of PG 1553+113 during the MJD 57055–57916 period. In each image we report the 15 GHz (upper frame), 24 GHz (middle frame), and 43 GHz (lower frame) results. The overlaid red triangles in the upper image (top panel) represent the 15 GHz OVRO light curve in the MJD range 56974–57990.
Open with DEXTER
In the text
thumbnail Fig. 3.

Left image: 15 GHz total intensity (contours and colour scale) image of PG 1553+113 during the third observing epoch (MJD 57276), convolved with a 0:6 mas 0:6 mas circular beam. The overlaid lowest total intensity contour is at 0.6% of the peak, with the following contours a factor of two higher. The ridge line is overlaid as a yellow line. Right image: same 15 GHz image in polar coordinates. The colour scale represents the total intensity emission and the white line represents the ridge line.
Open with DEXTER
In the text
thumbnail Fig. 4.

Upper image: ridge lines of PG 1553+113 for all 14 observing epochs at 15 GHz. Each epoch is represented by a different colour and is labelled with a progressive number. The two dashed arcs represent 0.6 mas and 1 mas from the core region. Lower image: jet direction across the observing epochs represented by the average ridge-line angle (see Sect. 3.3 for more details).
Open with DEXTER
In the text
thumbnail Fig. 5.

15 GHz residual stacked image after core subtraction (see Sect. 3.3). The colour scale represents the total intensity emission. The lowest contour is three times the image noise (that is ∼0.34 mJy beam−1), with the following contours being a factor of two higher. The ϕ symbol represents the angle of the funnel in which the jet wobbles across time.
Open with DEXTER
In the text
thumbnail Fig. 6.

Core region RM (upper panel), intrinsic polarisation angle (middle panel), and absolute EVPA orientation (lower panel) values during the epochs MJD 57055, 57078, 57650, and 57755. In the lower panel 15, 24, and 43 GHz are represented by black circles, red triangles, and blue star symbols, respectively.
Open with DEXTER
In the text
thumbnail Fig. 7.

Upper frame: two-dimensional WWZ power (colour scale on the right side) distribution as a function of time (x-axis) and period (y-axis) of the 15 GHz OVRO light curve over the MJD range 54696-58153. The significance intervals are plotted as white contours, in which the dotted, dashed, and solid lines represent 1σ, 2σ, and 3σ, respectively. Lower frame: OVRO light curve over the MJD range 54696-58153.
Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.