© ESO, 2016

1. Introduction

Radio-loud active galactic nuclei (AGN) and blazars in particular – constituting their beamed population – are among the most copious and variable emitters of radiation in the Universe. They intrinsically feature double, opposite-directed fast plasma outflows. These are most likely driven by the conversion of gravitational into electromagnetic and kinetic energy, in the immediate vicinity of spinning supermassive black holes (SMBHs; Blandford & Znajek 1977) or from the accretion disk surrounding them (Blandford & Payne 1982; for reviews cf. Meier et al. 2001; Blandford 2001; Meier 2003) The linear extent of these jets greatly exceeds the dimensions of their hosts and undisrupted can reach up to a few million parsecs (pc).

Their extreme phenomenology includes fast superluminal motion at parsec and sub-parsec scales, high degree of radio and optical linear polarization, along with rapid broadband flux density and polarization variability at time scales down to minutes. Their double-humped spectral energy distributions (SEDs) first peak between radio and soft X-rays because of synchrotron emission, and then at keV to GeV energies, likely owing to inverse Compton (IC) scattering (see e.g. Ulrich et al. 1997; Urry 1999). However, the detailed processes giving rise to blazar characteristics are still under intense debate. Outstanding questions include (i) the mechanism(s) that drive their rapid variability across the whole electromagnetic spectrum; (ii) the location in the jet where the high-energy emission originates; and (iii) the target photon field for IC up-scattering. Possibilities for the latter include the dusty torus and the ambient infrared (IR) photon field, the broad-line region (BLR), and/or accretion disk photons. Through their observational predictions, different scenarios can be tested. Theoretical considerations indicate that the bulk of high-energy γ-ray emission probably emanates from regions of low ambient density, otherwise the pair production process γ + γ → e⁻ + e⁺, would prohibit γ-ray photons from escaping the MeV/GeV production region, introducing a high γ-ray opacity environment (Dermer et al. 2012; Tavecchio & Ghisellini 2012). Nonetheless, observational evidence is still contradictory, indicating that GeV-emission can be produced close to the central engine, within the broad-line region (BLR), or even in the vicinity of the accretion disk some 100 Schwarzschild radii from the SMBH (e.g. Blandford & Levinson 1995). Other findings lend support to the “larger-distance scenario”, whereby MeV/GeV photons originate from shock-shock interaction or turbulent cells within the plasma flow (see e.g. Valtaoja & Teräsranta 1995; Marscher 2014).

Correlated variability and flare timing analyses between radio and γ-ray bands place the production site of γ rays a few pc upstream of the mm-band unit-opacity surface (Fuhrmann et al. 2014), while other findings point towards travelling or standing shocked regions downstream of the jet (Valtaoja & Teräsranta 1995). With the bulk of blazar activity taking place in regions close to the central engine, very-long-baseline interferometry (VLBI) at mm wavelengths is essential for addressing these questions. This technique is able to deliver the highest angular resolution attainable and to penetrate the opacity barrier that renders these regions inaccessible to lower frequency observations. As an example, work based on 7 mm VLBI and flux monitoring of OJ 287 places the γ-active region more than 14 pc away from its central engine (Agudo et al. 2011).

In 2008 August, the Fermi-GST (Gamma-ray Space Telescope, hereafter Fermi, Ritz 2007; Atwood et al. 2009) detected the blazar PKS 1502+106, which showed a rapid and strong γ-ray outburst followed by high and variable flux over the next months (Ciprini 2008; Abdo et al. 2010). This activity at high energies triggered an intensive multi-wavelength campaign, covering the radio, optical, UV, and X-ray bands. The outburst was accompanied by a simultaneous flare at optical/UV/X-rays (Pian et al. 2011), with a significantly delayed counterpart at radio bands (Fuhrmann et al. 2014).

The powerful blazar PKS 1502+106 (OR 103, S3 1502+10) at a redshift of z=1.8385, D_L=14176.8 Mpc (Adelman-McCarthy et al. 2008) is classified as a flat-spectrum radio quasar (FSRQ) with a mass of its central engine of ~10⁹M_⊙ (Abdo et al. 2010, and references therein). In X-rays, it is known as a significantly variable source (George et al. 1994). Radio interferometric observations with the VLA at 1.4 GHz reveal the large-scale morphology of PKS 1502+106, with a straight jet at a position angle (PA) of about −160° (Cooper et al. 2007). Previous VLBI findings suggest a curved, asymmetric jet and a multi-component, core-dominated source (Fey et al. 1996) with very fast apparent superluminal motion of up to (37 ± 9) c (An et al. 2004).

In the following, we present a comprehensive VLBI study of PKS 1502+106 at three frequencies: 15, 43, and 86 GHz (wavelengths of 20, 7, and 3 mm, respectively). Our millimeter-VLBI monitoring of PKS 1502+106 was triggered by the flaring activity that the source underwent and started in early 2008. Data at 43 and 86 GHz span the time period between 2009 May and 2012 May, comprising six observing epochs. VLBI monitoring at mm wavelengths is complemented by nineteen epochs of VLBA observations at 15 GHz from the MOJAVE program (Lister et al. 2009).

Throughout this paper we adopt S ∝ ν^{+ α}, where S is the radio flux density, ν the observing frequency, and α the spectral index, along with the following cosmological parameters: H₀ = 71 km s^-1 Mpc^-1, Ω_m = 0.27, and Ω_Λ = 0.73.

At the redshift of PKS 1502+106, an angular separation of 1 milliarcsecond (mas) translates into a linear distance of 8.53 pc and a proper motion of 1 mas yr^-1 corresponds to an apparent superluminal speed of 79 c.

The paper is structured as follows. In Sect. 2 we introduce the technical aspects of the observations and the data reduction techniques. In Sect. 3 we present our findings concerning the phenomenology of PKS 1502+106, while in Sect. 4 we show the deduced physical parameters. Section 5 discusses our findings in a broader context and in conjunction with previous work by other authors. In Sect. 6 we present our summary and conclusions.

2. Observations and data reduction

2.1. High-frequency 43 GHz and 86 GHz GMVA data

The mm-VLBI observations with the Global Millimeter VLBI Array (GMVA¹) comprise six epochs spanning the period between 2009 May 7 and 2012 May 18 with observations obtained twice per year, except for 2010 October. The experiments included the following European and US stations: Effelsberg (100 m, Ef), Onsala (20 m, On), Pico Veleta (30 m, Pv), Plateau de Bure interferometer (6×15 m phased array, Pb), Metsähovi (14 m, Mh), and Yebes (40 m, Ys) in Europe, and the Very Long Baseline Array (VLBA) stations that are equipped with 86 GHz receiving systems (8×25 m) in the US. For a complete presentation of the observational setup we refer to Tables 1 and 2. The GMVA experiments were carried out in dual-polarization mode at the stations where this setup is available (Ys and On only record left circular polarisation) at a bit rate of 512 Mbits/s. The observing strategy at 86 GHz comprised seven-minute scans of PKS 1502+106, bracketed by scans of calibrators for the full duration of the session. Calibrator sources included 3C 273, 3C 279, 3C 454.3, and PKS 1510-089. The 43 GHz scans, using only the VLBA, were obtained in gaps when European stations were performing pointing calibration. Data were correlated at the MK4 correlator of the MPIfR in Bonn, Germany.

2.1.1. A priori calibration

After correlation, the raw data were loaded in and inspected with the National Radio Astronomy Observatory’s (NRAO) Astronomical Image processing System (AIPS) software (Greisen 1990). Phase calibration, the most crucial step because of atmospheric and instrumental instabilities, was made using the task FRING within AIPS in the following manner. First, an initial manual phase calibration was applied incrementally to the data (Martí-Vidal et al. 2012). That is, calculating single subband phases and delays, first for the VLBA antennas and then for European stations, and applying them to the data. These phases and delays were estimated with respect to a reference antenna using a set of visibility data from a bright source from among those listed in Sect. 2.1. Subsequently, global fringe fitting was performed, and phases, delays, and rates were deduced and applied to the multi-band data. We performed the initial a priori amplitude calibration within AIPS and applied an atmospheric opacity correction to the data. The visibility amplitudes were calibrated using the measured system temperatures and telescope gain curves. Finally, fully calibrated visibility data were exported to DIFMAP for further analysis.

2.1.2. Imaging and model fitting

For the subsequent imaging and analysis we used the DIFMAP software (Shepherd 1997; Shepherd et al. 1994). After importing the visibility data, continuous deconvolution cycles were performed with the application of the CLEAN algorithm (Högbom 1974), followed by phase self-calibration. Amplitude self-calibration, with decreasing integration time, was performed after every sequence of CLEAN and phase self-calibration and only when a representative model for the source had been obtained. Extreme caution was exercised to avoid artificially and/or irreversibly altering the morphology of the source by introducing or eliminating spurious/fake or real source features. This was accomplished by comparing the model with the visibilities and by preventing the standard deviation of visibility amplitudes from varying by more than 10%.

To parameterize each image and quantitatively handle the relevant quantities, the MODELFIT algorithm within DIFMAP was employed. To simplify the analysis and reduce the number of unconstrained parameters, we only fit circular Gaussian components. The fit parameters are the component flux density, its size expressed as the full width at half maximum (FWHM) of the two-dimensional Gaussian and its position angle (PA) expressed in degrees (between 0° and 180° in the clockwise and 0° to −180° in the counterclockwise direction). We modeled the data separately at each frequency and each observing epoch. The number of components was initially unrestricted. The limit to the number of components was set by testing whether the addition of a new one caused the χ² value for the fit to decrease significantly.

We adopted the following uncertainties of the model fitting. For the component flux density an uncertainty of 10% is a reasonable and widely adopted estimate (Lister et al. 2009, 2013). The uncertainty in the position is set to onefifth of the beam size if the component is unresolved, otherwise we take onefifth of the component FWHM as an estimate of the uncertainty. This is increased to half a beam for unresolved knots at 86 GHz. The beam size for each observing epoch is calculated through (Lobanov 2005). The uncertainty in the PA is the direct translation of the positional uncertainty (ΔX) from mas to degrees using the simple trigonometric formula ΔPA = arctan(ΔX/r), with r the radial separation in mas.

To reliably define the jet components we adopted the following criterion. A physical region of enhanced emission within the jet should repeatedly appear in three consecutive observing epochs. We relaxed this criterion for 86 GHz data and only when images in specific intermediate epochs suffered from a low dynamic range that prevented us from discerning a component whose presence was indicated by data in adjacent epochs. In the process of obtaining the final kinematical model at each frequency, the component separation, size, and flux density were used as parameters combined with the requirement that they vary smoothly, in a way compatible with physical processes.

2.2. 15 GHz MOJAVE survey data

We also made use of the publicly available data from the Monitoring Of Jets in Active galactic nuclei with VLBA Experiments (MOJAVE²) program (Lister et al. 2009). Data on PKS 1502+106 comprise nineteen observing epochs that span the period between 2002 August 12 and 2011 August 15. Except for a few occasions, all ten VLBA stations participated in the observations (see Table 3). Data are provided as fully self-calibrated visibilities and as such, imaging and model fitting were performed, following the steps described in Sect. 2.1.2.

2.3. Single-dish 15, 43, and 86 GHz F-GAMMA data

Single-dish multi-frequency monitoring data were obtained within the framework of the long-term Fermi-GST AGN Multi-frequency Monitoring Alliance (F-GAMMA³) program (Fuhrmann et al. 2007; Angelakis et al. 2010). F-GAMMA is the coordinated effort for simultaneous observations at 11 bands in the range between 2.64 GHz and 345 GHz. Regular monthly observations including the Effelsberg 100-m and IRAM 30-m telescopes are closely coordinated to ensure maximum coherency. The observing strategy uses cross-scans, that is, slewing over the source position both in azimuth and elevation. After initial flagging, pointing offset, opacity, sensitivity, and gain-curve corrections are applied to the data (for details see Fuhrmann et al. 2014; Angelakis et al. 2015; Nestoras et al., in prep.).

From the broad F-GAMMA frequency coverage we only used the bands that matched our mm-VLBI frequencies, which are those at 14.60, 43.00 and 86.24 GHz. A detailed single-dish multi-wavelength study using the full spectral coverage of F-GAMMA as well as optical/X-ray data will appear in a forthcoming article.

2.4. Fermi/LAT γ-ray data

Fermi is the latest-generation space-bound telescope built for and dedicated to detailed observations of the high-energy sky. Since its launch in 2008, its two instruments, the Large Area Telescope (LAT, Atwood et al. 2009) and the Gamma-ray Burst Monitor (GBM), scan the full sky about once every three hours in the energy range from ~10 keV to > 300 GeV, providing unprecedented temporal and spectral coverage.

Here, we employed the monthly binned LAT γ-ray light curve of PKS 1502+106 presented in Fuhrmann et al. (2014), which covers the time period between MJD 54 707 (2008.66, 2008 Aug. 29) and MJD 55 911 (2011.96, 2011 Dec. 16). Details on the γ-ray flux (E> 100 MeV) light curve and analysis are given in the aforementioned reference (see also Abdo et al. 2010).

3. PKS 1502+106 jet phenomenology

3.1. Parsec-scale jet morphology

MODELFIT VLBI maps at 15, 43, and 86 GHz are presented in Figs. 1 through 3, respectively. PKS 1502+106 exhibits a pronounced core-dominated morphology with a continuous one-sided jet at parsec scales.

The entire series of VLBI images – at all frequencies – features a distinct region of enhanced emission responsible for the bulk of the observed flux density (see Sect. 3.3). This brightest feature is referred to as the core. For the purposes of our analysis we assumed it to be stationary for the full length of our observations and across all frequencies. This is consistent with the fact that a putative opacity shift is smaller than our positional accuracy (see Sect. 5.2).

In addition to the core, the 15 GHz jet can be decomposed into three to four distinct MODELFIT components at each observing epoch. At 43 GHz the jet is better represented by three components, while at 86 GHz we used between two to four MODELFIT components because of the varying dynamic range of the images and the intrinsic structural and flux density changes of the source. Nevertheless, we are able to cross-identify two components in all three frequencies and a third one only between 3 and 86 GHz. In our analysis we use the following nomenclature: at 15 GHz, five jet features can be identified (Fig. 1). The three earliest visible are labelled Ca, Cb, and Cc. This selection is due to their estimated ejection times, which are consistent with a separation from the 15 GHz core, before or close in time to the first observing epoch we employ (see Table 4 and Sect. 3.2 on obtaining the kinematical parameters). We use arithmetic labelling (C1, C2, C3) for components that can be positively cross-identified between at least two frequencies. More specifically, C1 and C2 are present at 15, 43, and 86 GHz, while component C3 is only visible at the latter two frequencies due to blending effects at 15 GHz.

Overall, the maximum angular extension of the jet at 15 GHz reaches about 4 mas towards the south-east with the inner jet, up to a distance of ~2 mas, lying at a smaller position angle (PA) and more towards the east (see Fig. 1). A misalignment between the inner and outer jet is evident from the maps at 43 and 86 GHz as seen in Figs. 2 and 3. This is further supported by the relative RA and Dec of fitted components shown in Fig. 4 at the two high frequencies, where regions closer to the core are resolved. Consequently, observations at 43 and 86 GHz reveal the inner jet with highest detail.

3.2. VLBI kinematics at 15, 43, and 86 GHz

To obtain each knot’s proper motion, μ, through the narrow curved jet of PKS 1502+106, we perform a weighted linear least-squares fit to their epoch-to-epoch radial separation from the respective core at each frequency. Component radial separations from the core are shown in Fig. 5 for 15, 43, and 86 GHz, respectively, where the linear fits are also visible. The inferred values of apparent speed in mas/yr and in units of speed of light (β_app = υ_app/c), along with ejection dates for all robust jet features are reported in Tables 4 through 6.

Compared to the kinematical results of the MOJAVE program presented in Lister et al. (2013), our component identification is slightly different. It is also based on the cross-identification with the high-frequency data, but the overall picture does not differ significantly and knot speeds along with their ejection dates agree very well.

The apparent velocity, β_app, of each component is given by (1)with

• μ

  the proper motion in rad s^-1,

• D _L

  the luminosity distance in m,

• c

  the speed of light in m s^-1, and

• z

  the redshift of the source.

As evident from Fig. 5 and Tables 4 through 6, PKS 1502+106 is characterized by extremely fast superluminal motion of features within the jet. Measured apparent velocities are in the range of ~10 to 22 c at 15 GHz, while at 43 GHz the range is 2–11 c and from 7 to 20 c for 86 GHz. The maximum β_app per frequency is , , , all characterizing the motion of the superluminal feature C1.

The independently deduced kinematical models per frequency agree in general, with only a few differences that need to be discussed. These are mainly encountered at our mid-frequency of 43 GHz for components C1 and C2. Here, C1 appears to be travelling downstream the jet with only about half the speed observed at 15 and 86 GHz, where the measurements agree within the observational uncertainties. This is attributed to the large positional uncertainties, but also to the small number of high-frequency VLBI observations (especially at 43 GHz), which are also limited in uv-coverage. Consequently, the absence of even a single data point can drastically affect the quality of the fit and the deduced kinematical parameters.

Inferred apparent component speeds allow for a study of the β_app profile along the radial direction of the jet. This is presented in Fig. 6. Up to a distance of about 1.5 mas from the core the apparent speeds observed follow an almost monotonically increasing trend and peak at about 1 mas at all frequencies. Beyond 1.5 mas downstream, the trend seems to break with only one but reliable measurement at 15 GHz for the historical component Ca, which only has a β_app = (9.7 ± 1.2) at a radial separation of ~2.8 mas, indicating a change of physical conditions within the flow, such as a deceleration and/or change of the viewing angle. In the immediate vicinity of the core, knot C3 moves outwards with an apparent superluminal speed of (5.0 ± 3.7) c and (7.1 ± 0.7) c at 43 and 86 GHz, respectively.

The XY (RA–Dec) position and the radial separation plots presented in Figs. 4 and 5 reveal the erratic path of the superluminal feature C2 at 43 and 86 GHz. The emerging pattern is consistent with and might even support the scenario that C2 is following a helical trajectory. This scenario cannot be corroborated because of the limited number of observations and the consequently low sampling rate of the underlying motion.

Other physical effects may also contribute to small observed differences. For example, the small angles under which blazars are viewed have the potential of greatly magnifying relativistic aberration effects. Additionally, the core can constitute a non-stationary feature with its position subject to erratic motion. However, such an effect is not clearly seen by the motion of other superluminal components. An additional effect with increasing resolution is that more compact regions are picked up by the interferometer, thus shifting the centroid of a fitted component within a larger emission region (as seen at lower frequencies).

For those jet components that can be identified across observing frequencies, ejection dates again agree well, with the exception of C2, whose wobbly motion does not allow a good estimation of the ejection date at 43 GHz. Components C1, C2, and C3 have been separated from the core in the time period covered by the VLBI monitoring. C1 and C2 were ejected between 2005 and 2006, and C3 in ~2008, at a time close to the onset of the multi-frequency outburst of PKS 1502+106.

3.3. Radio flux density decomposition

The flare in PKS 1502+106 is dominated by a single, clearly visible event from γ-ray energies down to radio frequencies, between MJD 54 689 and MJD 55 295. The upper panel of Fig. 7 features the Fermi/LAT γ-ray light curve, and the following panels depict the total intensity radio data, both from the filled-aperture F-GAMMA observations and from cm/mm-VLBI monitoring. The total single-dish flux density is decomposed into light curves of the individual VLBI core and jet components. We show those of the positively cross-identified components C1, C2, and C3, omitting the historical components at 15 GHz and not reliably identified features.

The core dominates the total flux density at VLBI scales at all frequencies. Radio flux density decomposed into core and MODELFIT components accounts for the largest portion of the total single-dish flux density and nicely follows the flare evolution. The sum of the flux density of all MODELFIT components deviates from the total single-dish flux density only in a few cases. This is attributed to the low cadence of mm-VLBI monitoring and the source’s extreme variability behaviour. The substructure in the well-sampled F-GAMMA and γ-ray light curves provides evidence for this fast variability.

Our mm-VLBI monitoring reveals that the radio flare, for its whole duration, is dominated by emission from the core, which at 15 GHz has varied in flux density by a factor of ~3 with a time scale of approximately two years (see Fig. 7). The same behaviour persists also at 43 and 86 GHz, where its flux density varies by a factor of ~6–7 and ~4–5, respectively, within our GMVA monitoring period.

At 43 and 86 GHz, where resolution allows, it is clear that the flare is not an attribute of the core alone. Here, the most recently ejected component C3 appears to be in a decaying phase during the first mm-VLBI observing epoch, but still shows significantly higher flux densities compared to all later epochs (see middle panels of Figs. 7 and A.1). Despite the limited number of data for C3, it can be clearly seen that except for the core, it is the only moving jet feature that is flaring, with its flux density dropping by a factor of ~20 at both frequencies. For completeness, all individual component light curves are shown in Appendix A.

3.4. Spectra of individual knots

While our high-frequency GMVA data at 43 and 86 GHz are concurrent, this is not the case for the 15 GHz VLBA data. To assemble quasi-simultaneous spectra of individual superluminal components within the jet of PKS 1502+106, we selected the 15 GHz epochs closest in time to each of the six 43/86 GHz GMVA observing runs. The quasi-simultaneous data forming each of the six spectra are labelled I to VI in Tables 1 through 3. Assembled spectra are shown in Fig. 8. The mean separation between the 15 GHz and 43/86 GHz observing epochs composing them is ~51 days.

Components C1 and C2 are positively cross-identified in all three observing frequencies (Sect. 3.1), and hence quasi-simultaneous spectra were assembled for them. Knot C1 is the only sufficiently sampled component that shows evidence of a synchrotron self-absorbed (SSA) spectrum featuring a spectral break at the turnover frequency, ν_m and flux density S_m. To obtain S_m, ν_m, and the optically thin spectral index we fitted a model SSA spectrum of the following form (e.g. Türler et al. 2000): (2)where

• ν _m

  is the turnover frequency in GHz,

• S _m

  is the turnover flux density in Jy,

• τ _m

  is the optical depth at ν_m,

• α _t

  is the spectral index of the optically thick part of the emission, and

• α

  is the optically thin spectral index.

The spectral index for the optically thick part was fixed to the canonical value α_t = 5 / 2. We fitted the observed spectra of C1 simultaneously and obtained S_m = (0.08 ± 0.01) Jy, ν_m = (33.2 ± 2.0) GHz and α = −(2.5 ± 0.4) as the average spectral index. The spectral index characterizing knot C1 is arguably very steep. This could indicate that the component at the highest frequency is resolved out as a result of the high resolution (i.e. small observing beam), leading to an underestimation of its flux density. An additional uncertainty in α of up to 15% is also to be expected, arising from the different uv-coverages (e.g. Fromm et al. 2013). To obtain the same parameters during the highest state of the knot, the fit was performed once more only at epoch 2009.9. Resulting parameters are as follows: S_m = 0.133 Jy, ν_m = 36.7 GHz and α = −1.98. These were used for the calculation of the magnetic field from SSA in Sect. 4.5.

For components C2 and C3, which show no sign of an SSA break, a power-law fit was performed separately at each observing epoch, and the mean optically thin spectral index ⟨α⟩, characterizing the time-averaged spectrum, along with its standard error are reported in Table 7.

C2 is also a well-sampled feature that during the whole length of our VLBI monitoring is characterized by a steep optically thin spectrum with a mean power-law index of ⟨α⟩ = −(1.0 ± 0.3) and a peak at ≤ 15 GHz. However, the situation is less constrained for component C3, which is only seen at 43 and 86 GHz as a result of the resolution limit of the 15 GHz observations. Here, the analysis yielded an optically thin two-point spectral index ⟨α⟩ = −(0.8 ± 0.3).

The core is the most compact stationary feature at each frequency, and it is characterized by a flat spectrum, that is, α ≥ −0.5. However, at epoch 2009.9 it appears to exhibit an SSA turnover at ν_m ~ 43 GHz. The mean spectral index of the core is found to be ⟨α⟩ = −(0.22 ± 0.07), indicative of unresolved substructure.

4. Inferred physical parameters

In this section, we deduce the physical parameters of the jet of PKS 1502+106, given the phenomenological characteristics presented in Sect. 3. We discuss the Doppler factor, viewing angle, and brightness temperature distribution along the jet and estimate magnetic field strengths.

4.1. Doppler and Lorentz factor estimates

As shown in Sect. 3.2, the highest apparent speed measured in the jet of PKS 1502+106 is β_app ~ 22 and corresponds to that of component C1 at 15 GHz. From here we can readily estimate the minimum Lorentz factor that characterizes the flow (see for example the review by Urry & Padovani 1995). It follows that γ ≥ γ_min for the jet with the minimum Lorentz factor given by (3)and for the observed maximum apparent velocity we obtain γ_min = 22.1. Consequently, if the jet is viewed at the critical angle (see Sect. 4.2), Doppler factors can be estimated as well, since δ_{βapp, C1} ~ β_{app, C1} = 22.1 at 15 GHz. For knots C2 and C3, it is δ_βapp ~ γ_{min, C2} = 9.3 and δ_βapp ~ γ_{min, C3} = 7.2, from the 86 GHz observations.

In addition to the Doppler factor estimate obtained from the minimum Lorentz factor, under the valid assumption of causality another estimate of the Doppler factor, δ_var, can be derived using the temporal variation of the flux density of individual components. For details on this method and for a full derivation see Jorstad et al. (2005) and Onuchukwu & Ubachukwu (2013), respectively. The variability Doppler factor is given by (4)where

• d _eff

  is the effective angular size of the component in rad and

• Δt_var

  is the variability time scale, [ = dt/ ln(S_max/S_min) ], with dt the time between S_max and S_min.

Components C1 and C3 can be used for this purpose. Both components are unresolved, that is, they are smaller than the beam size at the epoch of their highest respective flux densities at all three observing frequencies. In all cases (for C1 at 15 GHz and for C3 at 43 and 86 GHz) we used the minor axis of the restoring beam as the upper limit for the component angular size. This selection is based on the fact that the minor axis is positioned almost parallel to the jet axis, thus allowing the highest resolution to be achieved along it. Specifically, for the higher frequencies and component C3, the minor axis of the 86 GHz beam was selected, since it yields the most stringent upper limit to its size.

A simple inspection of the light curves indicates that component flux density variations are not adequately sampled, while the point of maximum flux density most probably precedes our first observing epoch. Nevertheless, based on the comparison with the single-dish data (Sect. 3.3), this point cannot be too displaced (in time), nor can it be much higher (in flux density) from the first data point in the mm-VLBI light curves. We therefore employed the method, since useful nominal values can be obtained. We calculated Δt_var, with dt and S_max obtained directly from the light curves of individual components.

Following the same argumentation, an estimate of the Lorentz factor can be obtained by (5)Results of the calculations are summarized in Table 8. Component C1 is characterized by the highest variability Doppler (δ_var ~ 51.4) and Lorentz (γ_var ~ 30.5) factors. All figures are consistent with the estimates of the minimum Lorentz factor, γ_min, presented previously. Knot C3 is also found to have consistent behaviour at the two high frequencies with and , while at 86 GHz the estimation yields and .

An additional Doppler factor estimate can be obtained from the relative strength of the synchrotron self-absorption and equipartition magnetic fields. This is discussed in Sect. 4.5.

4.2. Viewing angle towards PKS 1502+106

Using the observed apparent speed, we can set constraints on the aspect angle under which a jet is viewed in two ways.

First, by calculating the critical angle which results in the observed β_app with the minimum γ. The apparent speed is given by (Rees 1966) (6)where β is the intrinsic velocity in units of c and θ is the viewing angle. Using the definition of the Doppler factor, δ = [ γ(1 − βcosθ) ] ^-1, Eq. (6) can be rewritten as (7)and β_app is maximized (i.e. β_{app, max} = βγ) for sinθ_c = γ^-1 or sinθ_c = δ^-1, since for this angle θ_c, γ = δ. Using the minimum Lorentz factor, γ_min = 22.1, calculated above, the critical angle is (8)We can also estimate the viewing angle towards the source, based on causality arguments. This “variability viewing angle”, θ_var, is given by (9)This analysis was performed for two components with significant variability (see Fig. A.1) and yields viewing angles for each component separately. The results are shown in the last column of Table 8. For the fastest superluminal knot, C1, we obtain the smallest viewing angle θ_var ~ 0.8°. For C3 we obtain consistent values at 43 / 86 GHz of ~3.4° and ~2.9°, respectively.

Interestingly, components C1 and C3 are travelling with very different apparent speeds (Sect. 3.2) in two different regions of the jet – that is, in the outer (r> 1 mas) and inner jet (r< 0.5 mas), respectively. In combination with the difference of their inferred Doppler factors, this indicates a “two-region scenario”, wherein physical conditions are intrinsically different. This possibility is further discussed in the following.

4.3. Opening angle of the jet

In this section we investigate the apparent and intrinsic opening angles of the jet, φ_app and φ_int, by employing two different approaches. First, we use all MODELFIT components at all three frequencies and perform a simultaneous linear fit to their sizes to obtain the jet opening rate. In the second approach we deduce φ_app from each individual MODELFIT component based on its size and kinematical characteristics.

4.3.1. φ_int from simultaneously fitting the effective sizes

In Fig. 9 the effective deconvolved size, d_eff = 1.8 FWHM_d with , of all resolved MODELFIT components at all three frequencies is shown as a function of radial separation from the core. The factor 1.8 accounts for the fact that a Gaussian with FWHM equal to d roughly corresponds to a circle of angular diameter 1.8 d. We fitted d_eff with respect to the radial separation from the core, using a simple linear model assuming a constant jet opening angle. The least-squares fit yields a slope of , translating into an apparent half-opening angle of (40.4 ± 1.7)° and an apparent full-opening angle of φ_app = (80.8 ± 3.4)°. The de-projected opening angle is then given by (10)which for θ = θ_c = 2.6° (see Sect. 4.2) yields φ_int = (3.7 ± 0.2)°.

4.3.2. φ_int from individual knots

We performed the analysis in the uv-plane, using the component size, FWHM_d, and radial distance, r, from the core. First, we calculated the opening angle corresponding to each MODELFIT component separately. Assuming the component filling the entire jet cross section, φ_app is given by (11)Then, we averaged over effective size and distance, thus obtaining the average opening angle of a given component at a mean radial separation from the core. The intrinsic opening angle of the jet, making use of the critical viewing angle for each component, is given by Eq. (10). The results are shown in Table 9. The mean intrinsic opening angle for the jet, deduced by averaging the ⟨φ_int⟩ of all components, is (3.9 ± 1.0)°.

4.4. Radial distribution of the brightness temperature

The apparent brightness temperature of a given VLBI component in the source rest frame is given by (12)where

• S _ν

  is the component’s flux density in Jy,

• d _eff

  is the effective size of the emitting region in mas, and

• ν

  is the observing frequency in GHz.

In Fig. 10 we present the radial distribution of brightness temperatures along the jet of PKS 1502+106. The T_b is frequency dependent, with higher values at the lowest frequencies. Observed brightness temperatures approach and mostly exceed the equipartition limit (Readhead 1994) of about 5 × 10¹⁰ K for the core, at all frequencies (see Fig. 10). The core is at all times unresolved and thus the T_b values represent lower limits. The same trend of high T_b persists throughout the first mas downstream of the core, with a few components being resolved. Specifically, in the regions

• between 0 and 1 mas, high T_b is seen at all frequencies, shown by lower limits, with a handful of measured values of lower T_b, as observing frequency increases. They lie in the range from 10⁷ K at 86 GHz, to slightly above 10⁹ K at 15 GHz;

• between 1 and 2 mas, the division persists, indicative of the sampling of different portions of the jet. Here the measurements at 15 GHz are clustered in the range 10⁹–10¹⁰ K at 15 GHz, while at the two high frequencies the T_b is significantly lower, T_b ≥ 10⁷ K;

• beyond 2 mas, the jet is observable only at 15 GHz and the distribution of T_b is characterized by a decreasing trend from 10⁹ K to 10⁸ K at a distance of 4.6 mas.

The mean – over time – core brightness temperature at our three observing frequencies is  K,  K, and  K.

4.5. Magnetic field estimates

Estimates of the magnetic field strength for individual components can be obtained under the assumption that SSA is the dominant process shaping the observed spectra. Following Marscher (1983), the strength of the magnetic field can then be calculated as (13)where

• b(a)

  is a tabulated spectral index-dependent parameter (see Marscher 1983, Table 1),

• ϑ

  is the component size in mas,

• ν _m

  is the spectral turnover frequency in GHz, and

• S _m

  is the spectral turnover flux density in Jy.

Given the strong dependencies of B_SSA, mainly on ν_m and ϑ, we note that measurement uncertainties of these quantities translate into a large uncertainty of B_SSA, as given by Eq. (13). Nevertheless, we attempted the estimation for those components that the turnover frequency and flux density are somewhat constrained, namely knot C1 and the core. Even if a Doppler factor δ, correction is not applied, an apparent value along the jet can be obtained (Bach et al. 2005). Table 10 summarizes our findings.

C1 exhibits a strong magnetic field of ~5.0δ G. This high value could be an indication of a highly magnetized disturbance travelling along the jet. Its high apparent velocity of 22.1 c and high Doppler factor (22.1 ≤ δ ≤ 51.4) support this scenario.

The characteristic spectral flatness of the core is traditionally attributed to the superposition of a number of individual SSA components, each with a distinctly peaked spectrum. The spectral shape of the core at epoch 2009.9 (see Fig. 8) indicates an SSA peak close to a frequency of ~43 GHz. Taking this as the turnover frequency and with S_m = 2.3 Jy, we obtain an SSA magnetic field of B_SSA ~ 2δ mG.

Another constraint on the magnetic field can be set assuming equipartition between the energy of radiating relativistic particles (E_e ∝ B^-1.5) and that of the magnetic field (E_B ∝ B²). The equipartition magnetic field that minimizes the total energy content E_tot = (1 + k)E_e + E_B is (see e.g. Bach et al. 2005) (14)where

• k

  is energy ratio between electrons and heavy particles,

• f

  is a tabulated function of α, and upper and lower synchrotron frequency cutoffs, ν_a, ν_b in GHz,

• L

  is the synchrotron luminosity of the source given by , and

• R

  is the size of the knot in cm.

For k ~ 100 and f(−0.5,10⁷,10¹¹) = 1.6 × 10⁷, we obtain the expression (15)Knot C1 is characterized by an equipartition magnetic field of about 0.3 G and the core features the highest equipartition magnetic field, B_eq ~ 1.2 G. Results are shown in Table 10.

The two magnetic field estimates have a different dependence on the Doppler factor. On the one hand, B_SSA ∝ δ for the magnetic field strength calculated from SSA, and B_eq ∝ δ^{(2 / 7α + 1)} for the minimum equipartition magnetic field. We used these to estimate δ_eq, when possible, from the ratio . The estimated equipartition Doppler factor for C1, δ_{eq, C1} ~ 50, also agrees well with δ_var ~ 51.4. Given the assumptions, all values agree reasonably well.

5. Discussion

5.1. The γ-ray/radio flare

The Fermi/LAT light curve, adopted from Fuhrmann et al. (2014) and shown in Fig. 7 (top panel), spans the time period between MJD 54 707 (2008.66) through MJD 55 911 (2011.96). The source is already at high state since the beginning of our monitoring period. From then on, the rising trend continues until MJD 54 875 (2009.12), when the absolute maximum in the monthly averaged γ-ray photon flux (of ~1.3 × 10^-6 cm^-2 s^-1) is reached. In fact, there exists some amount of substructure in the γ-ray light curve. A second – this time local – maximum is reached on MJD 54 959 (2009.35) with comparable flux. Altogether, activity at high energies, from its start until the first and absolute minimum is reached on MJD 55 295 (2010.27), lasts for almost 650 days (shaded area in Fig. 7).

Savolainen et al. (2002) have made a general connection between radio outbursts and newly ejected jet features, with most radio flares accompanied by the ejection of one or more components. While the brightness of knots associated with γ-ray/radio flares typically decays fast, here we follow the evolution of one thanks to the enhanced resolution of mm-VLBI. As discussed in Sect. 3.3 apart from the flaring VLBI core, component C3 visible at 43 and 86 GHz only, is characterized by a decaying light curve with an initial flux density level – at epoch 2009.35 – almost ten times higher than that of the other two components (C1, C2) present in the flow and seen also at the two high frequencies. Although a certain amount of flux density blending between the core and C3 is plausible, we note the relatively high spatial separation of the two components (approximately one beam at 86 GHz in 2009.35) and the high signal-to-noise ratio of C3. These two facts add to the argument that C3 is indeed flaring/decaying, in addition to the flaring core.

Our observations are consistent with C3 having an apparent speed of 7.1 c (Sect. 3.2) and its ejection year, based on the kinematical models at 43 GHz and 86 GHz, is (2007.4 ± 2.2) and (2008.3 ± 0.3), respectively. The time of zero separation for C3 (red solid line in Fig. 7) coincides with the onset of the γ-ray flare when extrapolating the rate of γ-ray flux change, as it rises linearly, back in time. We conclude, based on the ejection date of superluminal component C3, that it is most probably associated with the γ-ray flare and its radio counterpart. Since its separation from the core, it is travelling downstream the jet with an apparent superluminal velocity (see Table 6). We note the larger uncertainty in the knot’s time of ejection at 43 GHz as compared to observation at 86 GHz. This is due to the sparser sampling at this frequency and the higher impact on the fit due to the very last data point in the radial separation plot (middle panel of Fig. 5). In any case, the two figures are consistent with each other. Assuming the speed of knot C3 to remain constant, at the times of the two γ-ray photon flux maxima, in 2009.2 and 2009.4 (at MJD 54 897 and 54 959) the knot is at a projected distance of ~0.7 and ~0.8 pc downstream the 3 mm core, respectively.

A plausible scenario is that a disturbance originating at the jet nozzle of PKS 1502+106, while still upstream of the 3 mm core – i.e. the unit-opacity surface at 86 GHz – produces the observed γ-ray emission. As it continues to move outwards, it becomes optically thin while crossing the 3 mm VLBI core, hence producing the radio flare. Alternatively, the core may represent the first conical recollimation shock of the flow (Daly & Marscher 1988; Gómez et al. 1995; Bogovalov & Tsinganos 2005; Marscher 2008), in which case the flare can be attributed to shock–shock interaction and subsequent enhancement of emission. In both cases, as the disturbance moves farther downstream, it is observed as knot C3 at its decaying flux density phase.

5.2. Distance estimates to the central engine

The debate as to where in the jet the high-energy emission originates is presently far from being resolved. In this paper, combining radio single-dish monitoring data and our findings from the VLBI analysis of PKS 1502+106, we explore and further constrain the region where the pronounced MeV/GeV activity of 2009 takes place and compare this with the SED modelling results of Abdo et al. (2010).

Under the assumption that the core takes up the entire jet cross-section and for a constant opening angle, the de-projected distance of the 86 GHz core from the vertex of the hypothesized conical jet can be estimated using the expression (16)With an average size for the core ⟨FWHM_core⟩ = 0.03 mas and a nominal opening angle φ_int = (3.8 ± 0.5)° (Sect. 4.3), the 86 GHz core is constrained to a distance of d_core ≤ 8 pc away from the jet base. This figure takes into account the uncertainty in φ_int. Here, d_core is an upper limit, as the core is at all times unresolved, rendering FWHM_core an upper limit to its angular size.

In an effort to detect and assess possible correlations between radio and γ rays, Fuhrmann et al. (2014) applied a cross-correlation function analysis between the single-dish radio light curves at 86 GHz and the Fermi γ-ray light curves of 54 blazars from the F-GAMMA sample. Radio emission is typically found to lag behind the activity in γ rays. More specifically, for PKS 1502+106, a statistically significant correlation is found between the two bands (see light curves in Fig. 7). Radio lags γ rays by 14 ± 11 days at a significance level above 99%. Using Eq. (3) in the same paper, this time lag translates into a distance of 2.1 pc upstream of the τ = 1 surface at 86 GHz for the γ-ray emission region. Assuming now that the core in our 86 GHz images coincides with the same surface of opacity transition from the optically thick to thin regime, we conclude that the γ-ray emitting region is located at ≤5.9 pc or ≤ 1.82 × 10¹⁹ cm away from the base of the conical jet. These figures also represent upper limits for the distance between the black hole and the γ-ray production region, due to the uncertain but most likely small separation between the black hole horizon and the base of the jet.

From the Mg ii line profile, the bulk BLR radius of PKS 1502+106 is estimated to be R_BLR = 2 × 10¹⁷ cm (see Abdo et al. 2010, and references therein). Assuming a stratified and extended BLR, a distance of 1.82 × 10¹⁹ cm corresponds to a region at its outer edge or farther away from it. The upper limit reported here agrees with the findings of Abdo et al. (2010), where the external radiation Compton mechanism (ERC) with the BLR photon field as a potential target for IC up-scattering significantly contributes to the high-energy part of the SED. However, we note the sharp drop of the BLR photon number density with distance (e.g. Kaspi & Netzer 1999), along with the fact that the γ-ray flare discussed in Abdo et al. (2010) precedes the one analyzed here by a few weeks, which most likely belongs to the pre-main-flare period. Consequently, the contribution of other seed photon fields (e.g. that of the dusty torus, or other) must be explored.

The precision we can obtain from the VLBI data at hand is not sufficient to estimate the magnetic field and other physical parameters with the core-shift effect (Marcaide & Shapiro 1984; Lobanov 1998). It is worth noting that the results of Pushkarev et al. (2012) place the 15 GHz core at a distance of ~8 pc from the vertex of the jet. Since the unit-opacity surface at 86 GHz ought to be closer to the jet base than the core at 15 GHz, the γ-ray emission site is very likely somewhat closer to the jet base than our upper limit of 5.9 pc suggests. Finally, we note that Eq. (16) must be used with caution in estimating the distance between the core and the jet base through VLBI; if the core is unresolved, it yields only an upper limit to d_core.

5.3. Intrinsic jet properties

Given the Doppler factor and the apparent speed of superluminal knots C1 and C3 estimated in Sects. 3.2 and 4.1, we explore the jet properties of PKS 1502+106. From the expressions of both the Doppler factor and Lorentz factor and making use of Eq. (6), we obtain the expressions for the Doppler factor as a function of apparent speed given the viewing angle, δ = f(β_app | θ), and given the Lorentz factor, δ = f(β_app | γ). The loci in Fig. 11 represent these two functions for a range of viewing angles and Lorentz factors (see caption for details). Overplotted are regions of estimated δ and β_app. These help us in constraining the intrinsic properties of the flow in terms of the viewing angle and the Lorentz factor.

Component C1 is characterized by an extreme variability Doppler factor of ~51. This figure is calculated based on causality arguments and agrees well with δ_eq (see Tables 8 and 10). The same variability analysis for knot C3 suggests comparable results from both 43 and 86 GHz, with a Doppler factor of ~12−15 (Table 8). Figure 11, where the inferred ranges of Doppler factor and apparent speed are shown, indicates the two regions present in the relativistic jet flow.

The differences in the inferred aspect angle and intrinsic speed between these two knots reflect physical differences in the jet flow of PKS 1502+106. First, the portion of the jet we observe at high frequencies, and knot C3, is possibly still in its accelerating phase (see also the apparent speed profile in Fig. 6). This could explain the aforementioned differences (in θ and δ) seen between knots C3 and C1, with the latter being some 1 mas downstream of the core. The jet, as revealed from the maps and the overall structure, appears to bend, hence the high Doppler factor farther out, traced by knot C1, could also be attributed to differential Doppler boosting. Taking only the 86 GHz data (red rectangle in Fig. 11) into account, the two regions wherein components C3 and C1 are travelling are incompatible when it comes to viewing angle and intrinsic speed.

The difference between the Doppler boosting factors in two regions at two different viewing angles θ₁ and θ₂, keeping β constant, is given by (17)Here, we test whether the change of Doppler factor of ~40 between the inner and outer jet at distances < 0.5 mas and > 1 mas can be caused only by the jet bending towards our line of sight. Components C3 and C1 are travelling in these two regions, respectively. The Doppler factor is estimated from variability for both of them and also from comparing the SSA and equipartition magnetic field for C1. C3 is characterized by δ_var ~ 12−15 while the two independent estimates agree very well for the Doppler factor of C1 to an extreme value of ~50. From Fig. 11 the intrinsic speed of C3 is constrained between β ~ 0.95 and β ~ 0.995. Taking even the highest value for β, keeping it constant, and applying Eq. (17) between θ₁ = 0.8° and θ₂ = 3.4° (the maximum range observed), we obtain that the difference in Doppler factor of about 40 between the two regions cannot be reconciled by changing the viewing angle alone. For these numbers only a δ₂ − δ₁ ~ 5 can be expected. Thus, the two regions are characterized by both a change in the viewing angle and also by acceleration – that is, by a change of Lorentz factor. This can also be seen directly from the plot in Fig. 11, since only a few loci can lead from one region to the other and only for the less constrained region for C3, coming from the kinematics at 43 GHz. Converting the aforementioned intrinsic velocities into the more convenient Lorentz factor, we obtain that even for the highest β = 0.995 (for C3) and for β = 0.999 (for C1), a change of Lorentz factor from at least γ ~ 10 to γ ~ 22 takes place in the region between roughly 0.3 to 1.1 mas. This large-scale acceleration region speaks in favour of magnetic acceleration scenarios (see e.g. Komissarov et al. 2007).

In estimating the viewing angle, it is worth noting that the critical value, θ_c, might not be a good estimate for the angle between the jet and the observer’s line of sight. This is the case for small kinematical Lorentz factors (γ_min), deduced by low observed apparent speeds. In fact, if the flow is intrinsically faster than implied by the γ_min, the same apparent velocity can be observed even at a viewing angle larger than the critical. However, the critical angle does not significantly differ from the true θ for ultra-relativistic flows. This conclusion is drawn from the behaviour of β_app with respect to the viewing angle (cf. any plot of β_app vs. θ) and the following energetics argument. An electron travelling at 0.6 c has only 1.5 times more (relativistic) kinetic energy than an electron at 0.5 c. The same electron at 0.999 c or γ ~ 22 possesses 3.4 times more energy than an e⁻ travelling at a speed of 0.99 c or γ ~ 7. Obviously, it is harder to accelerate the electron in the latter case than in the former, making the minimum Lorentz factor and consequently the critical angle good proxies for the intrinsic values when very high apparent speeds are observed.

Estimation of a variability viewing angle is also possible by combining single-dish flux density measurements and VLBI kinematics as in Hovatta et al. (2009), where the authors estimate θ_var = 4.7° for PKS 1502+106. However, this figure is based on a smaller β_app = 14.6 measured at 22 GHz and 37 GHz in the period from 1990.5 to 1995.5 (cf. Teraesranta et al. 1998).

Pushkarev et al. (2009) calculated a slightly smaller intrinsic opening angle of 3.11° for the jet. The reason is twofold. First, they obtained a slightly smaller intrinsic opening angle of 37.9°, which they subsequently de-projected using the variability viewing angle reported in Hovatta et al. (2009), thus obtaining the aforementioned value. The difference with the opening angle of (3.8 ± 0.5)° we find here is justified by the higher apparent speed of component C1 (β_app = 22.1) seen in our data.

The radial brightness temperature distribution cannot corroborate or falsify the two-region scenario, since in the first mas beyond the core most of the measurements constitute lower limits. This is because the components are unresolved at all three frequencies. The decaying trend at beyond 1 mas is clear, however (see Fig. 10). The scenario of acceleration within the first mas from the core, supported by both kinematics and variability, cannot thus be tested through the T_b. The case in which the part of the jet probed by the 86 GHz VLBI data is still at its accelerating phase (see e.g. Lee 2013), with the Doppler factor increasing downstream of the 86 GHz core, is compatible with the radial T_b profile. Additionally, we note that there is no evidence for any overall rising trend of T_b with distance as we move away from the core, but rather a steady decay seen at all frequencies.

6. Summary and conclusions

We presented a comprehensive three-frequency (15, 43, and 86 GHz) VLBI study of the γ-ray blazar PKS 1502 + 106 using GMVA and additional data from the MOJAVE monitoring program. Furthermore, we employed the densely sampled F-GAMMA single-dish light curves at matching frequencies, along with the Fermi/LAT monthly binned light curve at energies > 100 MeV. The work presented here allowed us to follow the multi-frequency flare of 2008–2010 and conclude the following:

• 1.

  PKS 1502 + 106 exhibits a compact, core-dominated morphology at all three frequencies with a one-sided, bent, parsec-scale jet. Component motion within the flow is characterized by high superluminal speeds in the range 5–22 c. The most extreme example is component C1, which travels downstream the jet at about 22 c at 15 GHz and comparable speeds at 43 / 86 GHz.

• 2.

  Doppler and Lorentz factor estimates for individual components were obtained using variability and VLBI kinematical arguments. The kinematical Doppler and Lorentz factors are δ_C1 ~ γ_{min, C1} = 22.1 at 15 GHz. For knots C2 and C3, they are δ_C2 ~ γ_{min, C2} = 9.3 and δ_C3 ~ γ_{min, C3} = 7.2, respectively. From causality, knot C1, showing significant variability at 15 GHz, is characterized by δ_var as high as 51.4 and γ_var ~ 30.5. Knot C3 is also found to have consistent behaviour at the two high frequencies with and , while at 86 GHz the estimation yields and . All figures are consistent with the estimates of the kinematical Doppler and Lorentz factors. Doppler factors are strikingly different for knot C3 at r< 0.5 mas and C1 at a distance r> 1 mas from the core.

• 3.

  An additional Doppler factor can be calculated from the different dependencies of B_SSA and B_eq on δ. The estimate for C1 agrees with those from kinematics and causality.

• 4.

  Using variability arguments, we were able to constrain the viewing angle towards the source to θ ~ 3.4° for the inner and θ ~ 0.8° for the outer portions of the jet, after about 1 mas. The calculation of the critical viewing angle for knot C1 at 15 GHz yielded θ_c = 2.6°.

• 5.

  We calculated the opening angle of the jet of PKS 1502 + 106 with two methods with very well agreeing results. The nominal de-projected opening angle is (3.8 ± 0.5)°.

• 6.

  The differences in apparent speed and viewing angle point towards a jet that bends after the first mas. This results in differential Doppler boosting – that is, in an increasing δ from ~12–15 to ~50, as traced by components C3 and C1, travelling at radial distances of < 0.5 mas and > 1 mas from the core, respectively. However, this Doppler factor gradient cannot be ascribed to jet bending alone. Acceleration must also be at play within the first mas of the jet.

• 7.

  By constraining the intrinsic speed, β, and viewing angle, θ, using the calculated Doppler factors and apparent velocities, we can clearly distinguish between the two regions of the jet in which knots C1 and C3 are travelling. Specifically, closer to the core at a distance ≤ 0.5 mas the parameters obtained for C3 at 43 GHz constrain the parameters of the flow to 0.95 ≤ β ≤ 0.995 and 1.5° ≤ θ ≤ 10°. The same component at 86 GHz sets slightly stronger constraints on the viewing angle with 3° ≤ θ ≤ 8°. On the other hand, in the region where C1 is travelling, its observed δ and β_app constrain the flow to β ≥ 0.999 and 0.8° ≤ θ ≤ 2.5°. A clear division is seen between these two jet regions. Hence, we conclude that the jet of PKS 1502 + 106 bends towards us in the region beyond ~1 mas downstream of the core while it also accelerates.

• 8.

  The radial brightness temperature profile indicates very high T_b for the core region at all three frequencies. We find a trend of decreasing T_b with higher observing frequency and with increasing distance from the core. This could indicate a magnetically dominated core region probed by the 86 GHz observations where Poynting flux has not yet been fully converted into kinetic flux.

• 9.

  Component spectra allow calculating the magnetic field strength in different regions of the jet. Component C1 exhibits a distinctive SSA peak. The core shows a flat spectrum with ⟨α⟩ ~ −0.22 but at one epoch (2009.9) signs of a turnover at ~43 GHz. By comparing the SSA and equipartition minimum magnetic fields, we conclude that for component C1, B_SSA is much stronger than B_eq, indicating a high degree of Doppler boosting. Finally, for the core it is B_SSA/δ ≥ 2 mG and B_eq = 1.2 G.

• 10.

  The radio flux density decomposition into distinct VLBI components and the comparison with the single-dish radio flux density outburst indicate that the bulk of the radio emission originates from the core at all frequencies. Except for the core at 43 / 86 GHz, knot C3 shares a significant radio flux density level that decays with time. C3 is not resolved at 15 GHz because it is blended with the core.

• 11.

  Given the radio flux density decomposition and the estimated ejection time of knot C3, which is coincident with the single-dish mm flare onset, we conclude that it is responsible for the radio flare observed in PKS 1502 + 106 during the period 2008–2010. Furthermore, the previously established correlation between the radio flare and γ rays indicates that component C3 is responsible for the γ-ray flare and its radio counterparts. Arguably, the 2008–2010 flare of PKS 1502 + 106 is an event originating in a disturbed region at the jet nozzle. While still upstream of the 3 mm core, it produces the increased γ-ray emission observed, and as it continues to move, it becomes optically thin (at 86 GHz) close to the 3 mm VLBI core, thus enhancing its brightness. It then continues on, and after the 86 GHz optically thick region is seen as knot C3 at its decaying flux density phase.

• 12.

  We conclude that the γ-ray emitting region is located at ≤ 5.9 pc or ≤ 1.82 × 10¹⁹ cm away from the base of the jet. These estimates constitute upper limits for the distance between the γ-ray production region and the SMBH itself.

We conclude that PKS 1502 + 106 is a source whose complex structural dynamics needs to be further investigated with higher cadence high-resolution imaging at mas and sub-mas scales.

---

Acknowledgments

V.K. acknowledges the help of B. Boccardi and E. Ros during discussions and for their comments that improved the manuscript. Thanks also go to the anonymous referee, whose valuable comments increased the quality of the paper. V.K., J.H., I.N., and I.M. were supported for this research through a stipend from the International Max Planck Research School (IMPRS) for Astronomy and Astrophysics at the Universities of Bonn and Cologne. This research has made use of data from the MOJAVE database that is maintained by the MOJAVE team. This research has made use of data obtained with the Global Millimeter VLBI Array (GMVA), which consists of telescopes operated by the MPIfR, IRAM, Onsala, Metsahovi, Yebes, and the VLBA. The data were correlated at the correlator of the MPIfR in Bonn, Germany. The VLBA is an instrument of the National Radio Astronomy Observatory, a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. Partly based on observations with the 100-m telescope of the MPIfR (Max-Planck-Institut für Radioastronomie) at Effelsberg and observations carried out with the IRAM 30-m Telescope. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). The single-dish millimetre observations were closely coordinated with the more general flux density monitoring conducted by IRAM (Institut de Radioastronomie Millimétrique), and data from both programs are included in this paper.

References

Appendix A: VLBI component light curves

The appendix plots feature all VLBI component light curves at all observing frequencies (Fig. A.1).

Appendix B: Additional tables

All Tables

All Figures