| Issue |
A&A
Volume 705, January 2026
|
|
|---|---|---|
| Article Number | A85 | |
| Number of page(s) | 15 | |
| Section | Extragalactic astronomy | |
| DOI | https://doi.org/10.1051/0004-6361/202453625 | |
| Published online | 09 January 2026 | |
Jet reorientation revealed by intermittent jet activity in radio galaxy 0954+556
1
Shanghai Astronomical Observatory, Chinese Academy of Sciences, Shanghai 200030, China
2
School of Physical Science and Technology, ShanghaiTech University, Shanghai 201210, China
3
Instituto de Astrofísica de Andalucía-CSIC, Glorieta de la Astronomía s/n, E-18008 Granada, Spain
4
School of Astronomy and Space Science, University of Chinese Academy of Sciences, Beijing 100049, China
5
Department of Space, Earth and Environment, Chalmers University of Technology, Onsala Space Observatory, SE-439 92 Onsala, Sweden
6
School of Physics and Astronomy, Yunnan University, Kunming 650500, China
7
National Astronomical Observatories, Chinese Academy of Sciences, Beijing 100101, China
8
Xinjiang Astronomical Observatory, Chinese Academy of Sciences, 150 Science-1 Street, Urumqi 830011, China
9
Chinese Academy of Sciences South America Center for Astronomy, National Astronomical Observatories, CAS, Beijing 100101, China
10
Instituto de Estudios Astrofísicos Facultad de Ingeniería y Ciencias Universidad Diego Portales, Av. Ejército 441, Santiago, Chile
★ Corresponding authors: This email address is being protected from spambots. You need JavaScript enabled to view it.
Received:
29
December
2024
Accepted:
17
October
2025
Context. Intermittent jet activity of active galactic nuclei (AGNs) is a common phenomenon, whereas significant jet reorientation during episodic jet activity in relatively young radio galaxies is rarely reported. The quasar 0954+556 at z = 0.903 is an intriguing source exhibiting an unusual radio jet structure with significantly different jet directions at kiloparsec (kpc) and parsec (pc) scales. At kpc scales, images from the Very Large Array (VLA) exhibit a bright core, a linear jet extending ∼24 kpc to the northwest, and a discrete jet component ∼16 kpc to the northeast. At pc scales, images from the Very Long Baseline Array (VLBA) show a two-component structure with a projected separation of ∼360 pc in the north–south direction.
Aims. The peculiar structure of 0954+556 might result from jet reorientation. Here, our aim was to investigate the possible mechanism via multiscale and multifrequency deep radio images.
Methods. We performed VLA and VLBA observations of 0954+556. Together with some existing data in the NRAO data archive, we made multiple VLA images at 1.4–22 GHz and VLBA images at 1.7–43 GHz for various image analyses of the jet structure.
Results. We identified the location of the radio core at pc scales, detected the faint counter-jets at both pc and kpc scales for the first time, and revealed a diffuse emission region connecting pc- and kpc-scale forward jets. Our spectral index distribution and spectral aging analysis indicate that 0954+556 might undergo at least two episodes of jet activity during the current AGN phase. Moreover, pc-scale polarization maps display a well-resolved spine-sheath polarization structure.
Conclusions. It seems that the jet direction of 0954+556 changed significantly during intermittent jet activity. This may explain the different jet orientations and spectral ages observed from kpc to pc scales. The research provides a strong case that AGN jet direction might change rapidly on timescales of one million years.
Key words: instrumentation: high angular resolution / techniques: interferometric / galaxies: active / galaxies: jets
© The Authors 2026
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.
This article is published in open access under the Subscribe to Open model. This email address is being protected from spambots. You need JavaScript enabled to view it. to support open access publication.
1. Introduction
It is widely accepted that nearly all galaxies in the Universe host supermassive black holes (SMBHs) at their centers (Richstone et al. 1998). Most galaxies are dormant due to the lack of fuel-feeding central black holes. Active galactic nuclei (AGNs, Urry & Padovani 1995) are a small fraction that produce high luminosity in a concentrated region at the center of a galaxy where the SMBH is accreting significant amounts of gas. As a result of accretion, the magnetic field anchored in the rotating disk drives winds that carry energy and angular momentum away from the disk, leading to jet formation and outflows (BP model, Blandford & Payne 1982). Alternatively, the energy and angular momentum are extracted from a spinning black hole by a purely electromagnetic mechanism (BZ model, Blandford & Znajek 1977).
It was a question of whether the active state of a galaxy is a repeated phenomenon or a one-time and long-term activity, after which the galaxy never reactivates. Hatziminaoglou et al. (2001) studied all galaxies, and found they go through different phases of nuclear activity over their lifetimes, which can be supported by the observational fact that there are more young sources than sources with extended and old structures (see Section 1 in Czerny et al. 2009). The recurrent activity can be recognized when new activity occurs and the emissions of past episodes remain for a long time. Consequently, a young radio source embedded in a relic structure is expected.
Intermittent jet activity has been reported many times. Double-double radio galaxies (DDRG, Schoenmakers et al. 2000) and X-shaped radio galaxies (XRG, Cheung 2007) serve as prime examples of such intermittent phenomena. A DDRG consists of a pair of double lobes with a common center, the outer lobes being the relic of past activity and the inner lobes being a manifestation of new jet activity. To date, four cases of triple-double radio galaxies, which contain three pairs of double lobes resulting from three separate episodes, have been identified: B0925+420 (Brocksopp et al. 2007), J1409-0302 (Hota et al. 2011), J1216+0709 (Singh et al. 2016), and J1225+4011 (Chavan et al. 2023). An XRG displays a conventional pair of active lobes and an additional pair of low surface brightness “wings” that are considered as remnants from a rapid realignment. In addition to these two types, signs of intermittent activity may be present in young radio galaxies that are in the nascent stage of radio galaxy evolution, reflecting the initial phases of their radio emission.
In addition to the observational evidence of intermittent jet activity, many theoretical studies have proposed possible explanations for AGNs switching on and off. Mergers (Merritt & Ekers 2002; Gopal-Krishna et al. 2003; Blundell 2008) may play a significant role on long timescales (106–108 years). The infalling galaxy passes through the host galaxy and turns around every ∼107 yr, causing instabilities in the host galaxy during each pass. After several encounters, the colliding galaxies merge completely, usually taking up ∼108 yr. Ionization instabilities in a narrow and unstable zone of the accretion disk can propagate throughout the disk, resulting in intermittent activity on timescales of 106–108 yr (Siemiginowska et al. 1996). Additionally, the short timescale activity is related to radiation pressure instabilities in the disk and mainly operates on timescales 104–106 yr (Janiuk et al. 2002; Czerny et al. 2009; Grzędzielski et al. 2017).
For this work we focused on 0954+556 (also named 4C +55.17, z = 0.903,1Ahumada et al. 2020), which was suggested as a case of intermittent activity due to its peculiar morphology (Rossetti et al. 2005). 0954+556 has been categorized as a flat-spectrum radio quasar (FSRQ) because of the presence of broad optical emission lines in its spectrum, high optical/UV core luminosity (absolute B-band magnitude, MB < −23), and high γ-ray luminosity (Lγ ≃ 1047 erg s−1) (Wills et al. 1995; Véron-Cetty & Véron 2006; Adelman-McCarthy et al. 2008; McConville et al. 2011). At parsec (pc) scales, two components separated by ∼50 mas were detected by the Very Long Baseline Array (VLBA) at 1.7 and 4.8 GHz, resembling a double oriented in the north–south direction. The total extent corresponds to a projected linear size of 640 pc (80 mas). In contrast, at kiloparsec (kpc) scales, 0954+556 presents a triple structure with a projected linear size up to 36 kpc, i.e., 4.5 arcsec, under the detection of the Very Large Array (VLA) at 1.4 GHz. At 5 GHz, the triple structure is resolved into a bright core, a linear jet extending northwest for approximately 3 arcsec, and a discrete component 2 arcsec northeast of the core. As the parsec-scale structure coincides with the kiloparsec-scale core, the morphology of 0954+556 has naturally raised a question about the reason for the jet changing direction between the two scales.
The classification of 0954+556 is controversial. It meets the criteria of compact symmetric objects (CSOs, Readhead et al. 1996; Kiehlmann et al. 2024), featuring a very low radio variability and a symmetric double structure with a linear size less than 1 kpc in VLBA images, similar to the classical doubles of Fanaroff-Riley Class II. Accordingly, it has been classified as a case of a CSO at parsec scales, and it was considered to be a case of a medium symmetric object (MSO, the larger version of CSOs, with a linear size ranging from 1 to 20 kpc) at kiloparsec scales, as the radio flux is dominated by the core region, i.e., the CSO part. McConville et al. (2011) further argued the preference of 0954+556 as a CSO. As a subclass of young radio galaxies, CSOs could be interpreted in several ways (O’Dea & Saikia 2021): as sources in their early stages of evolution, as sources confined by a dense interstellar medium, or as sources experiencing intermittent jet activity. Therefore, whether this source arises from intermittent activity deserves further investigation.
Previous radio studies of 0954+556 focused primarily on total intensity maps, due to a lack of data, thus leaving the above questions unanswered. For this work we further investigated the source by analyzing the total intensity, spectral index distribution, spectral aging, and polarized emissions, using data from VLA and VLBA over a wide frequency range from 1.4 to 43 GHz. From the perspective of CSOs, we explored the intermittent jet activity of 0954+556 and how the intermittent activity affects source morphology.
The paper is organized as follows. Section 2 summarizes the observations and data reduction. Section 3 presents the obtained multifrequency images and results of analysis on variability at γ-ray and radio wavelengths. In Section 4 the results of the spectral analysis are presented. We also present the results of VLBI polarimetry at 1.7, 5.0, and 8.4 GHz in Section 5. Section 6 includes the discussion and Section 7 gives a summary. Cosmological parameters are adopted as follows: H0 = 67.4 km s−1 Mpc−1, ΩM = 0.315, ΩΛ = 0.685 (Planck Collaboration VI 2020), and therefore 1 mas ∼ 8 pc. Throughout the paper, Sν ∝ ν−α, in which Sν represents the flux density, ν is frequency, and α is the spectral index.
2. Observation and data reduction
2.1. Observations
Over the last two decades, we have observed 0954+556 multiple times either as the target source or as the phase-reference calibrator, thus a great deal of VLA and VLBA data have been accumulated by now, including the following: two VLA datasets at 8.4 and 22.5 GHz from December 2000 (observation code: AH0721); one VLBA dataset at 1.7 GHz from February 2000 (BH065), two VLBA datasets at 5 and 8.4 GHz acquired in July 2002 (BH096), all designed for polarimetry; and three VLBA datasets at 5, 22, and 43 GHz in mid-2016 (BZ061, BL222). Furthermore, we obtained two additional VLA datasets at 1.4 and 5 GHz in November 2000 (AM0672), one VLA dataset at 15 GHz in October 1996 (TESTT), as well as five VLBA datasets at 5 GHz between September 2002 and August 2005 (BM170, BM208) from the NRAO archive2. Table 1 lists the information of each observation, such as frequency, observation code, array, and observation date.
Image parameters.
To make a brief analysis of the variability of 0954+556, we collected γ-ray photon flux data from the Fermi-LAT Light Curve Repository (LCR)3 (Abdollahi et al. 2023) and the total flux at radio band from the Owens Valley Radio Observatory (OVRO) 40 m Telescope Monitoring Program (Richards et al. 2014). The γ-ray data were divided into seven-day bins, covering a period between 8 August 2008 and 6 December 2024; the data represent the average integrated photon flux with energies from 0.1 to 100 GeV. In the OVRO program 0954+556 was observed at 15 GHz every four days from 19 April 2009 to 30 December 2011.
2.2. Data reduction
The data were calibrated using the NRAO Astronomical Image Processing System (AIPS) (Greisen 2003) following the standard procedures described in the AIPS cookbook4. The data from VLA were inspected and edited using LISTR, QUACK, and UVFLG, followed by flux and phase calibration. The flux calibrator was set to 3C 48 or 3C 286 in our calibration by SETJY and CALRD. Finally, the task CALIB was used to calibrate amplitudes and phases.
Following Appendix C in the AIPS cookbook, we calibrated data from VLBA in AIPS. On each track of VLBA data at 1.7, 5, and 8.4 GHz, we performed global fringe fitting twice for 0954+556 to correct the structure-induced phase errors, considering the pc-scale structure of 0954+556 at these frequencies is extended. For the first global fringe fitting, we solved the residual delay, rate, and phase by assuming a point source model. The output uv data was transformed into images with a clear structure of 0954+556 by the iteration of CLEAN and phase self-calibration in Difmap (Shepherd 1997). The second global fringe fitting was performed by using the cleaned images. Then we applied the resulting solutions containing corrections for the structure-induced phase errors to 0954+556. Once the initial calibration was completed in AIPS, self-calibration (both phase and amplitude) and imaging loops were performed using Difmap to obtain images with high dynamic ranges.
The 22 GHz VLBA data were imaged twice after global fringe-fitting. First, we imaged the data in the full uv range (2.3–636 Mλ), and only diffuse emissions were detected marginally. The ratio of peak flux density to the noise level in the cleaned image was so low (∼6) that it was infeasible to perform self-calibration. In the second trial, the maximum uv range was limited to 48.0 Mλ (see detailed description in Section 2.4). On the resultant image, the peak-to-noise ratio was raised to 20. Phase-only self-calibration was applied very carefully for several iterations to recover the flux density.
At 43 GHz, in global fringe fitting, solutions with a signal-to-noise ratio (S/N) higher than 5 (the default S/N threshold in the AIPS task FRING) only took 10%, so we adjusted the threshold of S/N to 4, and then the good solutions increased to 20%. The good solutions were applied to the data, and a tentative 43 GHz VLBI image was produced. As the peak-to-noise ratio of the resultant image was only 10, no self-calibration was performed.
The resulting images are presented in Figs. 1 and 2. The parameters of the images (e.g., frequency, array, date, uv range, size of the synthesized beam, rms noise level in a clean image, peak intensity, and total clean flux density) are listed in Table 1.
![]() |
Fig. 1. VLA flux density images of 0954+556 at kpc scale. The contours are drawn at –1, 1, 2, 4, 8, …, of the first contour level (3σrms). The synthesized beams are plotted in the bottom left corner of each image. In the top left panel, the red circles overlapping the contours present the results of modelfit. In the top middle panel, the 5 GHz image shows a weak counter-jet. To better reveal the weak counter-jet, the 5 GHz image was restored with a circular beam of 900 mas, as shown in the top right panel. |
![]() |
Fig. 2. VLBA flux density images of 0954+556 at pc scale. The contours are drawn at –1, 1, 2, 4, 8, …, of the first contour level (3σrms). The synthesized beams are also plotted in the bottom left corner of each image. The image in the top right corner is the tapered image at 1.7 GHz (UVTAPER 0.5, 10). |
2.3. Polarization calibrations
We performed polarization calibration in AIPS on VLBA observations designed for polarimetry, i.e., 1.7 GHz data obtained in February 2000 (BH065) and 5 GHz and 8.4 GHz data obtained in July 2002 (BH096). First, we used the AIPS task RLDLY to solve the R-L delays (right circular and left circular polarization). The instrumental polarization (D-terms) was solved by running the AIPS task LPCAL on a weakly polarized source, which is DA 193 in project BH065, and OQ 208 in project BH096.
The electric vector position angle (EVPA) of the linear polarization was calibrated by comparing the integrated EVPA of a calibrator obtained in our observation with the reference EVPA of the calibrator reported by the Very Large Array polarization monitoring program5 in nearby epochs. For the project BH096 observed on 20 July 2002, J2253+1608 (3C 454.3) was used as the calibrator for EVPA calibration at 5 and 8.4 GHz (Hu et al. 2024a). The nearest epoch of J2253+1608 in the Very Large Array polarization monitoring program was on 10 July 2002, and the differences between the reference EVPA and the integrated EVPA of J2253+1608 were 29.8° ±2.6° at 5 GHz and −41.9° ±2.0° at 8.4 GHz.
As explained in detail by Hu et al. (2024b), there is no information close to the observational date of the project BH065, so it is impossible to recover the absolute EVPA distribution of 0954+556 at 1.7 GHz by using an EVPA calibrator. While the absolute EVPAs at 1.7 GHz are arbitrary, the relative EVPAs across the entire total intensity structure are accurate. Our results reveal the EVPA distributions of 0954+556 at 5 and 8.4 GHz are consistent (Section 5), indicating that the rotation measure (RM) is close to zero. In addition, at kpc scales, the RM of 0954+556 was estimated to be +1 ± 1 rad m−2 in Simard-Normandin et al. (1981). Assuming that the polarization remained stable from 2000 to 2002, and given the RM value, we can recover an approximate absolute EVPA distribution at 1.7 GHz by aligning it with those at higher-frequency maps. Specifically, we rotated the EVPAs at 1.7 GHz to match the orientation at 5 GHz, using the southern lobe as a reference. The resulting EVPA distributions at three frequencies are presented in Section 5, and the image at 1.7 GHz generally represents the absolute EVPA distribution at this frequency.
2.4. Preparations for spectral analysis
The spectral index distribution map requires the alignment of images. The core of an extragalactic jet is a region with an optical depth τs ≈ 1. Due to synchrotron self-absorption, the core position varies with the observing frequency and follows the relation rcore ∝ ν−1/kr (Blandford & Königl 1979), which is a phenomenon known as the core-shift effect (Blandford & Königl 1979; Kovalev et al. 2008; Sokolovsky et al. 2011). In contrast, the jet part is typically optically thin, showing a steep spectrum. As the optically thin jets are assumed not to change with frequency, it is used to align images at different frequencies. Aligning optically thin compact components across images is the most common way to correct for the core shift (Hovatta et al. 2014). In terms of sources lacking compact features and showing diffuse structures, alignment is instead performed using the two-dimensional cross-correlation method (Croke & Gabuzda 2008), which adopts spatial correlations of optically thin extended jet structure to register images.
To avoid errors brought by different uv samplings and resolutions, the images used in spectral index maps and the spectral aging analysis were produced with data in the same uv range, in which the minimum uv distance is the smallest value of the highest frequency and the maximum is the largest value of the lowest frequency. In addition, images were set to the same map size and cell size, and restored with a beam of the lowest frequency. Pixels with flux density below 5.4σrms, alongside those with spectral index uncertainty exceeding 40%, were blanked when we produced the spectral index maps.
To produce the kpc-scale spectral index map, the 1.4 and 5 GHz VLA data were restricted to a common uv range of 5.1 − 168.7 kλ and were inverse Fourier transformed into images with a beam of 1449 × 1941 mas, −85.13° to obtain the 1.4–5 GHz spectral index map. For the 5 and 8.4 GHz VLA data in a uv range of 11.3 − 561.2 kλ, the inverse Fourier transform was used, and the resultant images were restored with a beam of 445.3 × 587.2 mas, −81.7°.
To produce pc-scale spectral index maps, 1.7, 5 (20 July 2002), and 8.4 GHz VLBA data were used. For the data in a uv range of 2.3 − 48.0 Mλ, the inverse Fourier transform was used and the images were restored with a beam of 5.4 × 7.3 mas, −15.5° to produce 1.7–5 GHz index map; for data in a uv range of 5.4 − 39.5 Mλ, the inverse Fourier transform was used and the images were restored with a beam of 4.3 × 7.0 mas, −26.0° to produce 5–8.4 GHz index map.
To produce kpc-scale images for spectral fitting in Section 4.2, VLA data in a uv range of 30.8 − 168.6 kλ were utilized, and the resulting images were restored with a beam of 1257 × 1676 mas, −84.4°. At pc scales, VLBA data at 1.7, 5, 8.4, 22 GHz incorporating the 15 GHz VLBA data from Rossetti et al. (2005) in a uv range of 7.1 − 48.0 Mλ were utilized, and the resulting images were restored with a beam of 5.0 × 6.5 mas, −18.3°. As listed in Table 1, the uv range of data at 43 GHz is 22.3 − 1105.4 Mλ. If the 43 GHz data is added into spectral analysis, the uv range used in the spectral analysis will be 22.3 − 48.0 Mλ. Many visibilities in every dataset will be flagged. Thus, data at 43 GHz was not used in the spectral analysis.
3. Results
3.1. Morphology
Figure 1 shows the kpc-scale morphology of 0954+556 at 1.4, 5, 8.4, 15, and 22 GHz. At 1.4 GHz, two extended structures are observed to the northeast and northwest of the brightest emission region (hereafter component NE, component NW, and kpc-scale core, respectively), resembling a triple structure. The clear triple structure is revealed by 5 and 8.4 GHz images with higher angular resolution. The northwestern extension is resolved into a linear jet structure (hereafter referred to as linear jet) extending at a position angle of −60° for about 3 arcsec, apparently brightened at its end. The northeastern extension is resolved into an individual knot located about 2 arcsec away from the kpc-scale core. On our 15 and 22 GHz VLA images, the kpc-scale core and the innermost 0.5 arcsec of the jet were detected. The end of the jet was only detected as diffuse emissions. The rest of the jet and the northeastern knot were resolved and were not detected due to the sensitivity limitation.
In addition, the images at 1.4 and 5 GHz are likely to reveal a counter-jet. First, we performed MODELFIT in Difmap to fit 1.4 GHz data with several circular Gaussian models, which are superimposed on the contours in Figure 1. Table 2 lists the parameters of each Gaussian, with uncertainties estimated following Lee et al. (2008). Notably, a circular Gaussian located southeast of the core may signal a counter-jet. Then, at 5 GHz, a weak emission is detected southeast of the core above 3σrms. This feature is considered reliable, as it is clearly visible in the residual map. To emphasize the low surface brightness emission, we restored the 5 GHz image with a circular beam of 900 mas, approximately twice the size of the initial beam, as presented in the top right panel of Figure 1. With a larger beam, the weak emission exceeds 6σrms, enhancing its reliability. Therefore, it probably represents the counter-jet at kpc scale.
Parameters of the fitted circular Gaussian models for the VLA image at 1.4 GHz.
At pc scales, 0954+556 mainly exhibits two components positioned in the north–south direction, as shown in Figure 2, consistent with previous studies. According to Rossetti et al. (2005), the two components were named component N and component S, and there was no indication of where the core is. In addition, Rossetti et al. (2005) detected the diffuse tail-like emission (hereafter referred to as tail) extending from component S toward the west at 1.7 GHz. In our VLBA image at 1.7 GHz, the tail is traced farther, extending more than 100 mas to the northwest, roughly aligned with the kpc-scale linear jet. We also detected some emission extending northeastward from component N for a few tens of mas, likely an indication of a counter-jet (CJ). These structures are more clearly revealed by the slightly tapered image (top right panel of Figure 2).
In the VLBA images at 5 GHz, Rossetti et al. (2005) and McConville et al. (2011) detected component N and component S, while our observation revealed a greater number of faint features. As mentioned in Section 2.1, we accumulated seven epochs of 5 GHz VLBA data over a period of 14 years from July 2002 to July 2016. In Figure 2, we present the image of the latest epoch with the highest dynamic range. Due to an on-source time of ∼3 hours, the CJ was more clearly detected, extending first to the northwest and then bending to the northeast. The bright portion of the tail was also detected. Because of the limited on-source time and the lack of antennas from other 5 GHz datasets, to reconstruct images with higher sensitivity we concatenated uv data from 20 July 2002 and 22 September 2002; from 14 June 2004 and 3 November 2004; and from 4 April 2005 and 18 August 2005. Since the source is stable over seven epochs, and the two uv datasets used for concatenation in each year are close in time, our results are reliable. Three reconstructed images (in 2002, 2004, and 2005) are shown in Figure 3, and the image parameters are listed in Table 3. In Figure 3, the image in 2002 only show the emission of the CJ extends northwest, while the images in 2004 and 2005 show the CJ extending northwest and bending northeast, consistent with the image in 2016. Thus, we consider the reliable CJ to be the counter-jet at pc scale.
![]() |
Fig. 3. VLBA flux density images of 0954+556 at 5 GHz. Left panel: Reconstructed image of the visibilities from a combination of data on 20 July 2002 and 22 September 2002. Middle panel: Reconstructed image of the visibilities from a combination of data on 14 June 2004 and 3 November 2004. Right panel: Reconstructed image of the visibilities from a combination of data on 4 April 2005 and 18 August 2005. The contours are drawn at –1, 1, 2, 4, 8, …, of the first contour level (3σrms). The synthesized beams are plotted in the bottom left corner of each image. |
Image parameters of concatenated 5 GHz VLBA datasets.
At 8.4 GHz, the brightest emitting region in component N elongates in the north–south direction and possibly turns southwestward, indicating an inner jet in this region. Most of the emissions from component S were resolved out, and only the brightened region and its southern edge were detected, implying a spectral steepening between 5 and 8.4 GHz in this component. The Radio Fundamental Catalog (RFC) provides the highly accurate coordinates of distant extragalactic radio sources derived from X–S dual-frequency astrometric analysis, helping build the International Celestial Reference Frame (ICRF). As 0954+556 is registered with ICRF, the image center of 0954+556 is aligned with the position provided by RFC. We used the latest position from RFC (rfc_2025a6, Petrov & Kovalev 2025) to roughly represent the brightest position in the 8.4 GHz VLBA image, i.e., the image center, and used the optical position in Gaia Data Release 37 (Gaia Collaboration 2023). Following Xu et al. (2021), the angular offset between the Gaia and VLBI positions is 2.58 ± 0.22 mas. This small offset supports the identification of component N as the region with the radio core, and the CJ component as the counter-jet.
The innermost core region of 0954+556 is revealed at higher frequencies. The 15 GHz VLBA image in Rossetti et al. (2005) displayed a well-defined jet in component N, initially extending southeast before bending southwest. In our work, 0954+556 is imaged for the first time at 22 and 43 GHz at pc scale. The 22 GHz map produced with data in the uv range between 2.3 and 48.0 Mλ is presented in the lower left panel of Figure 2. The structure is consistent with the bending observed at 15 GHz. The lower right panel of Figure 2 shows the tentative 43 GHz image of 0954+556. The detected emission is compact, slightly extending southward, and can be fitted with a circular Gaussian with a radius ∼78 μas (0.624 pc).
The jet morphology of 0954+556 exhibits a complex bending pattern both at pc and kpc scales. At pc scales, the jet initially extends southeast and then bends southwest, as shown in the 15 and 22 GHz VLBA images. The lower frequency images (1.7, 5, and 8.4 GHz) indicate that the jet possibly bends back toward the southeast, reaching the edge of component N, and continues southwestward, connecting to component S. The pc-scale counter-jet initially extends northwest, then turns toward the northeast. To connect the structure from pc to kpc scales, the observed tail emission might suggest that the pc-scale forward jet undergoes a single bend toward the northwest to align with the kpc-scale forward jet. However, it seems that the pc-scale counter-jet possibly bends southeast to align with the kpc-scale counter-jet, and then bends again to the north to reach component NE.
3.2. Variability
Figure 4 shows the γ-ray light curve and the 15 GHz radio light curve of 0954+556 produced with the data introduced in Section 2.1. No flare is seen in either of the light curves. To quantify the variability of 0954+556, we utilized the fractional variability amplitude (Fvar; Edelson et al. 2002; Vaughan et al. 2003), and the related equations are given here.
![]() |
Fig. 4. Light curve of 0954+556. The top panel shows a Fermi-LAT γ-ray light curve of 0954+556 from 8 August 2008 to 6 December 2024, divided into seven-day bins. All points are plotted along with their statistical errors. The bottom panel is the radio light curve obtained from OVRO at 15 GHz, from 19 April 2009 to 30 December 2011, with an interval of four days. The gray area in the top panel corresponds to the bottom panel. |
The mean square error is defined as
where N is the number of data points, and σerr, i is the error on the data point. In addition, the sample variance is
where xi is the value of the data point, and
is the mean of xi. Therefore, the fractional variability amplitude is given by
and the uncertainty of fractional variability amplitude ΔFvar is calculated as
The Fvar derived from the LCR dataset and the 15 GHz OVRO data is 6.9 ± 2.5% and 1.6 ± 0.3%, respectively, and are considerably lower than those observed in blazars, which are known for their strong variability. For example, in the LCR dataset, the blazar PKS 0805-07 exhibits an Fvar above 75% (Akbar et al. 2024), and the blazar 3C 279 shows an Fvar of 100% (Rani et al. 2018). In the radio band, the Fvar of blazar Mrk421 obtained by OVRO is 15% (Schleicher et al. 2019), while blazar 3C 273 reaches an Fvar of 29% at 15 GHz from the data from the University of Michigan Radio Astronomy Observatory (UMRAO) (Soldi et al. 2008). Compared with typical blazars, 0954+556 is not significantly variable, thus indicating stable emission zones for 0954+556. Notably, although 0954+556 is classified as a FSRQ in the Fermi catalog, Ballet et al. (2020) reported that while the blazars all have fractional variability exceeding 10%, 0954+556 shows a flux decrease of about 13% from August 2008 to August 2018, without any flare. To assess the stability of the structure, we analyzed all the total flux density images at 5 GHz, from 2002, 2004, 2005, and 2016. All of the 5 GHz images were convolved with a common beam corresponding to the observation in 2016, and an average image was constructed by taking the mean of the four images. The residual maps were produced by subtracting each individual image from the average, as shown in Figure A.1. The intensity differences across these residuals range from –15 mJy/b to 15 mJy/b, which is negligible compared to the intensity of the source. Thus, the structure of 0954+556 remained stable over the 14-year span, without significant variability.
4. Spectral analysis
4.1. Spectral index distribution
The kpc-scale spectral index maps were derived from the VLA images at 1.4 and 5 GHz and 5 and 8.4 GHz, respectively. At 1.4 and 5 GHz, the NE component was fitted with circular Gaussian models with flux densities of 61.4 ± 14 mJy and 24.6 ± 5 mJy, respectively, yielding α ∼ 0.7 ± 0.2 and suggesting an optically thin synchrotron emission between the two frequencies. Analogously, at 5 and 8.4 GHz, the NE component was modeled with circular Gaussian models of 22.5 ± 4.3 and 10.5 ± 3.1 mJy, respectively, resulting in α ∼ 1.5 ± 0.7, further supporting its optically thin nature. Thus, its positions measured at each frequency pair should be identical with each other (Kovalev et al. 2008). We therefore aligned the images by the centroid of component NE to produce kpc-scale spectral index maps.
The top panel of Figure 5 shows the spectral index distribution at kpc scales, overlaid with contours representing the 1.4 and 5 GHz VLA images. In the 1.4–5 GHz spectral index map, the spectral index distribution of the core, the component NE and NW tend to be flat, with α above 0. While the component NE and NW tend to have approximate spectral index distribution, the core displays a flatter spectral index. In the 5–8.4 GHz spectral index map, the steepest spectrum appears in component NE. The spectrum at the kpc-scale core remains flat. Going outward along the linear jet, the spectrum steepens and the α varies. At the jet terminus, the spectrum is flattened, and the α decreases. The brightening and spectral flattening in the jet terminus indicate a local acceleration of emitting particles.
![]() |
Fig. 5. Spectral index distribution maps of 0954+556. From left to right, the panels show the distribution of spectral index between 1.4 and 5 GHz at kpc scales, between 5 and 8.4 GHz at kpc scales, between 1.7 and 5 GHz at pc scales, and between 5 and 8.4 GHz at pc scales. The overlapped contour images are a VLA image at 1.4 GHz with a uv range 5.1 − 168.7 Kλ, a VLA image at 5 GHz with a uv range 11.3 − 561.2 Kλ, a VLBA image at 1.7 GHz with a uv range 2.3 − 48.0 Mλ, and a VLBA image at 5 GHz with a uv range 5.4 − 39.5 Mλ, respectively. All counters start at 3σrms and increase by a factor of 2. The color represents the value of spectral index α. Pixels with flux density below 5.4σrms, alongside those with spectral index uncertainty exceeding 40%, are blanked. |
At pc scales, we used 1.7 (7 February 2000), 5 (20 July 2002), and 8.4 (20 July 2002) GHz VLBA images to produce spectral index maps. Since no compact optically thin component was identified at the pc scales, we utilized a 2D cross-correlation method (Fromm et al. 2013; Kim & Trippe 2014) to determine the core-shift between two frequencies, and thus to align the images. Given that the 1.7 GHz and 5 GHz observations are not simultaneous, we attempted to construct a spectral index map by combining the 1.7 GHz data (7 February 2000) with three 5 GHz datasets (20 July 2002, 3 November 2004, and 4 April 2005). The three resulting spectral index maps are very similar, and thus we chose the closest epoch. In addition, the source has very weak variability according to Section 3.2, making our spectral index map reliable and convincing. The bottom panels of Figure 5 present the 1.7–5 GHz spectral index map overlapping on the 1.7 GHz VLBA contour image and the 5–8.4 GHz spectral index map overlapping on the 5 GHz VLBA contour image.
The 1.7–5 GHz spectral index map shows that the spectral index distribution of components N and S are very similar; both exhibit flat spectra across most of their regions and steeper spectra toward the edges. This similarity makes it challenging to identify the core region based solely on the spectral index map at these frequencies. However, in the 5–8.4 GHz spectral index map, component S shows a significantly steeper spectrum with α ∼ 2, while component N retains a relatively flat region, suggesting continued particle acceleration in component N rather than S.
4.2. Radiative age
An aging synchrotron spectrum from a relativistic electron population with an initial injected power-law energy distribution is described as a power-law steepening at high frequencies due to particle energy losses (Kardashev 1962; Pacholczyk 1970; Murgia et al. 1999). The break frequency where the spectrum starts steepening depends on the elapsed time since the particles were injected, and on the intensity of the magnetic field.
To perform spectral aging analysis, we utilized a Python-based package SYNCHROFIT8 (Turner et al. 2018; Turner 2018; Quici et al. 2022), which uses an adaptive maximum likelihood algorithm to fit the observed radio spectrum. SYNCHROFIT provides continuous-injection on (CI-on), continuous-injection off (CI-off), Jaffe-Perola (JP), and Kardashev-Pacholczyk (KP) models (Kardashev 1962; Pacholczyk 1970; Jaffe & Perola 1973; Murgia et al. 1999; Orienti et al. 2007) to estimate the radiative age of jet components. In the CI-on and CI-off models the electron population is injected continuously, so electrons with different ages are mixed. The injection keeps taking place in the CI-on model, whereas the continuous injection was kept for a long time and switched off in the CI-off model. In the JP and KP models the particles were injected in a single pass. The difference between the JP and KP models is whether electron pitch-angle scattering is considered or not. Among the four models, the break frequency νbr and the non-aged spectral index αinj are constrained by fitting the flux density across multiple frequencies. In addition, the CI-off model gives the remnant fraction T, i.e., the fractional time spent in an inactive phase to the total source age.
The radiative age τsyn is derived as
in which z is the redshift of the source, B is the magnetic field, BIC = 3.18 × 10−6(1 + z)2 gauss is the magnetic field equivalent to the microwave background, and v is a proportional constant. The details are in Section 2.2.1 of Turner et al. (2018).
We used the equipartition magnetic field strength (Beq) representing the condition of minimum total energy density as the magnetic field to derive τsyn. The following equation given by Miley (1980) was adopted to calculate the magnetic field strength,
Here k = 1 is the ratio of heavy particle energies to electron energies, η = 1 is the filling factor, and z is the redshift of a target. Based on assuming cylindrical symmetry, θx and θy represent the size of a component in a plane perpendicular to the line of sight, and s is the length of a component along the line of sight; ϕ is the angle between the magnetic field and the direction of the line of sight; F0 in Jy is the flux density at frequency v0; and α is the spectral index.
The model assessment is conducted by a reduced-chi square equation,
where N is the number of observation data points, P is the number of fitted parameters, N − P represents the degree of freedom, xi are the observed values, μi are the expected values based on the model, and σi is the uncertainty of the observed values.
The integral flux density of each component was extracted with AIPS task TVSTAT. Its uncertainty was estimated as
(Bruni et al. 2021), where Nregion is the number of pixels of the chosen emission region, and Nbeam is the number of pixels of the resolution beam. The estimated flux density with the uncertainty of each component is listed in Table 4.
Total flux density of each component for spectral analysis.
4.2.1. Spectral aging at kpc scale
Section 2.4 described how we prepared the images for spectral fitting, and the resultant images are displayed in Figure A.2. After being convolved with the beam at the lowest frequency, the linear jet becomes a circular Gaussian, so the overall structure can be recognized as the core, component NE, and component NW (Figure A.2). The fitted values of αinj, νbr, and T are listed in Table 5. The fitting curves overlapped by data points are presented in Figure 6.
![]() |
Fig. 6. Four different fit results for each component of 0954+556. The red, blue, black, and yellow dashed lines represent the JP, KP, CI-on, and CI-off models respectively. The VLA core and component N are fitted by a simple power law. For component N, the data at 43 GHz was not used in the spectral analysis as mentioned in Section 2.4. |
The flux density of the core component is well fit by a simple power law with an index of ∼0.44. For component NW, the fitted injection index from the CI-on model is αinj = 0.54, while the other models yield a slightly lower value of αinj = 0.49. The CI-on and the KP models both give a νbr within the observed frequency range, whereas the CI-off and JP models produce a νbr that is a bit higher than the highest-frequency point. The quality of fits is good, with the χred2 of 1.1 for the CI-on model and ∼0.6 for the other models. To estimate Beq, according to the morphology in Figure A.2, we assumed a circular geometry for component NW with a diameter of 3.5 arcsec, resulting in the corresponding magnetic field strength as 4.4 × 10−5 gauss. The resultant τsyn for component NW ranges from 0.3 to 1.4 Myr. It is notable that above 15 GHz, the emissions of NW are dominated by the hot spot at the jet terminus, and thus the derived radiative age should be much lower than the actual source age. Due to the lack of data points, fitting the spectrum of component NE was not performed.
4.2.2. Spectral aging at pc scale
The integral flux densities with uncertainties of component N and component S were extracted in the same way, as described in Sections 2.4 and 4.2. The derived values and images are listed in Table 4 and Figure A.3. The fitted values of αinj and νbr are listed in Table 5. The fitting curves overlapped by data points are presented in Figure 6.
Results of radiative age.
The emission of component N between 1.7 and 22 GHz can only be fit with a simple power law with an index of 0.38. The observational spectrum of component S is curved, and apparently our data disfavor the CI-on model as the CI-on model steepens too gradually. Among all models, the CI-on model yielded the highest χred2, more than twice that of the other models. For the remaining models, the fitting values of νbr ∼ 5 GHz are reasonable since they are consistent with what is revealed by the total intensity maps. Coupled with the brightening at the southwestern edge of component S below 8.4 GHz, a local acceleration in this region is strongly indicated. Furthermore, the spectrum steepening sharply beyond 5 GHz suggests that the acceleration might be historical and had been quenched for a long time.
The geometric size of component S was estimated using a circular region with a diameter of 40 mas (Figure A.3). The magnetic field strength in component S was estimated to be 2.6 × 10−3 gauss, assuming energy equipartition. Therefore, in the CI-off, JP and KP models the value of τsyn was derived to be ∼3.6 kyr, ∼2.7 kyr, and ∼1.8 kyr, respectively. The actual age of S should be older than these estimated radiative ages, considering the historical local acceleration in this region. The separation between the brightest part of component N and component S is ∼45 mas, i.e., 360 pc. Assuming a typical advance speed of 0.03–0.1c for CSOs, the estimated dynamic age of component S is 11.7–39.1 kyr, close to its radiative age.
5. Polarization
The 1.7 GHz VLBA observation conducted in February 2000, and the 5 and 8.4 GHz VLBA observations conducted in July 2002 were designed for polarimetry (Section 2.1). The results of the polarimetry are presented in Figure 7. The fractional polarization (m = P/I) is shown as pseudo-color, overlapping contours of the total intensity. The observed EVPAs are indicated by vectors. Following the description in Hovatta et al. (2012), we estimated the uncertainties of Stokes I, Q, U images, and considered the errors associated with antenna feed leakage. Then we estimated the uncertainties of P, m, and EVPA, and clipped pixels with I < 3σI, P < 3σP, and m < 3σm, as described in Hu et al. (2024b). Table 6 lists the statistics of the above images: total polarized flux density, peak polarized intensity, rms noise of the polarization image, and the average fractional polarization (
). From Figure 7 we found that the EVPA distributions at 5 and 8.4 GHz are nearly identical, and its RM at kpc scales was estimated to be +1 ± 1 rad m−2 (Simard-Normandin et al. 1981), allowing us to obtain the EVPA distribution at 1.7 GHz, as mentioned in Section 2.3.
Parameters of polarization images.
![]() |
Fig. 7. Polarization images at different frequencies. The restoring beam is indicated by a gray ellipse in the lower left corner of each image. The contours are drawn at –1, 1, 2, 4, 8, …, of the first contour level (3σrms). Superimposed on the contours of flux density are vectors whose lengths are proportional to the polarized intensity; the directions indicate the direction of the E-vectors, and the colors show the fractional linear polarization. In all polarization images, the 1 mas line is 1 mJy/b. We clipped pixels with I < 3σI, P < 3σP, and m < 3σm. |
0954+556 is a rare case in that the detected polarized emissions cover most of the jet structure at pc scales. In general, the jet at pc scales propagates toward the south (from component N to S) and then bends to the northwest (the tail). The polarization structure of 0954+556 reveals a clear spine-sheath configuration throughout the entire structure (Pushkarev et al. 2005; Gabuzda et al. 2014; Gabuzda 2021). In the central spine region, the EVPAs are aligned with the jet axis, whereas the EVPAs are roughly perpendicular to the jet direction at the edges, i.e., the sheath. This means the inferred B-field at the inner jet is oriented orthogonally to the jet and at the jet edges along the jet (Pacholczyk 1970), which is a strong indication of the existence of the helical magnetic field at pc scales within the jet (Gabuzda 2021). Notably, the EVPAs in the spine region basically follow the jet even at the bending, further supporting the idea that the jet has bent a few times to reach the brightest region of component S, as mentioned in Section 3.1.
The average fractional polarization
is 3.8%, 5.2%, and 7.5% at 1.7 GHz, 5 GHz, and 8.4 GHz, respectively. Interestingly, the fractional polarization m increases from the central spine toward the outer sheath. At all three frequencies, component N shows an increase in m from ∼0.05 in the central region to ∼0.3 at the jet edge, with the localized edge reaching m ∼ 0.7. Component S exhibits the same trend at 1.7 and 5 GHz as component N, and shows a high m of 0.4–0.7 at 8.4 GHz. The tail observed at 1.7 GHz shows high values of m ranging from 0.3 to 0.7. This is also an indication of the existence of a helical magnetic field (Gabuzda 2021).
6. Discussion and conclusion
6.1. Classification
As mentioned in Section 1, the classification of 0954+556 is controversial; it is either a FSRQ in the γ-ray frequency range or a CSO in the radio frequency range. A blazar is usually referred to as an AGN with a powerful relativistic jet aligning well with our line of sight, characterized by extremely high brightness temperatures (Tb) over the threshold of inverse Compton catastrophe (1012 K), apparent superluminal motions observed at pc scales, and fast and drastic variations of flux density over a wide range of frequencies (Homan et al. 2015; Liodakis et al. 2018; Homan et al. 2021). However, for 0954+556, the Tb of the brightest emission component detected by VLBA at 5 and 15 GHz was only at the magnitude of 108 K (Rossetti et al. 2005; McConville et al. 2011). Following the method in Homan et al. (2021), we estimated the Tb of the brightest component at 8.4 GHz to be ∼6 × 108 K. No superluminal motion was reported. No evidence of variability was found in γ-rays over the 19-month LAT observing period from August 2008 to March 2010, and its optical flux density was reported to change by only 30% over 7 years from 2002 to 2009 (McConville et al. 2011), nor was any variability at radio wavelengths found (Fan et al. 2007). Our analysis also supports that the variation of 0954+556 is much lower than typical blazars at γ-ray and radio bands (see Section 3.2), indicating that 0954+556 is a non-aligned AGN.
The structure of 0954+556 detected by VLBA resembles a compact double with a projected separation of 360 pc, leading to 0954+556 constituting a CSO in previous studies. However, the core position of this CSO was not identified. In our results, unlike that of component N, the spectral index distribution of component S shows a significant decline in its flux density over 5 GHz, indicating a difference between the two components. As analyzed in Section 4.2.2, the radiative age analysis reveals that there is no particle acceleration in component S, whereas component N harbors the radio core in which the acceleration of the particles is still taking place. There is evidence that the radio core is located in the brightest emission region of component N. First, the optical position of this source lies just 2.58 ± 0.22 mas from the peak of the 8.4 GHz VLBA image (Section 3.1). Given that the angular size of component N is approximately 25 mas, this close alignment strongly suggests that the optical and radio cores coincide. Second, our observations detected the counter-jet at 1.7 and 5 GHz (Figure 2). Third, a jet extending southward is observed in the central part of the structure (Sections 3 and 5). At 15 GHz, this structure appears as a well-defined jet initially extending southeast and subsequently turning southwest. The southeast segment of this jet is also visible at 22 GHz. The images at 1.7, 5, and 8.4 GHz show that the jet continues to elongate to the southeast. To connect to the brightest region of component S, the jet then extends southwest again. The EVPAs within the jet are along the jet, while they are perpendicular to the direction of jet propagation at the edges of the structure (Figure 7), indicating a spine-sheath polarization structure. Additionally, the fractional polarization increases from the central spine toward the outer sheath, from ∼0.05 to ∼0.7. Finally, the radio spectrum of component N between 1.7 and 22 GHz can only be fit by a power law with an index of 0.38 (Section 4.2.2), implying a continuous supply of relativistic electrons. Thus, at pc scales we identify that the core of this CSO is located at component N.
Although the radio core is near the position of peak brightness, it is not identified in our VLBA images. The core usually exhibits a flat spectrum, and our spectra index maps lack the resolution to distinguish a flat-spectrum core feature. At 22 and 43 GHz where the core might be dominant, the structure of 0954+556 is highly resolved on long baselines. Therefore, we disfavor identifying the emission seen at 22 and 43 GHz as the core. The core might be intrinsically faint, blended with nearby jet emission, or self-absorbed even at higher frequencies.
6.2. Intermittent activity
The main morphology of 0954+556 in this study is consistent with previous work (Reid et al. 1995; Rossetti et al. 2005; McConville et al. 2011), a triple structure at kpc scales, and a two-component structure located in the north–south direction at pc scales, which naturally raised a question of jet realignment. Regarding this concern, Rossetti et al. (2005) proposed that 0954+556 could be a case of a smothered source or of recurrent activity, involving two distinct phases responsible for the kpc-scale and the pc-scale structures, respectively, as inferred from the total intensity structure. Our spectral analysis provides clear evidence that the observed emitting components in 0954+556 emerged from independent episodes of activity, rather than the pc-scale jet is now feeding the kpc-scale structure to form a continuous jet extending from pc to kpc.
0954+556 has experienced at least two episodes of jet activity, forming the triple structure at kpc scales and a two-component structure at pc scales. The two-sided jets are clearly detected at both spatial scales. The kpc-scale jets are collimated with a position angle of ∼–60°. At pc scales, the high-resolution VLBA images revealed the north–south structure with a bending counter-jet and a bending forward jet. The estimated radiative ages of the outermost and innermost regions differ by over two orders of magnitude, approximately 0.3–1.4 Myr for component NW and 1–4 kyr for component S. Three prominent components, i.e., component NE at the kpc scales, component NW at the kpc scales, and component S at the pc scales, are notably misaligned, and the counter-jets corresponding to components NW and S are detected, so it might lead us to assume that component NE represents a third jet activity. However, due to insufficient data, we were unable to perform a spectral fitting to derive a reliable radiative age for component NE. Furthermore, from 1.4 to 5 GHz, the two components NE and NW have similar spectral index distributions. At the 5–8.4 GHz spectral index map, although the component NE is slightly steeper than component NW, the opposite kpc-scale radio jets or lobes in double-lobed galaxies and quasars often have different spectral index distributions (Liu & Pooley 1991; Dennett-Thorpe et al. 1997). Thus, we consider that components NE and NW originated from the same jet activity. The radiative age of the kpc-scale structure estimated by component NW reflects that the jet activity is on a short timescale, which allows us to detect the remaining emission from the earlier jet activity, i.e., the tail extending westward in VLBA images at 1.7 and 5 GHz.
In addition, the low RM in 0954+556 proves the intermittent activity indirectly. Young radio galaxies such as CSOs are typically embedded in the dense interstellar medium of their host galaxies, exhibiting relatively high RM between 102 and 104 rad m−2 (An et al. 2010; Rastello et al. 2016; Tremblay et al. 2016) and low fractional polarization of less than 1%9 (O’Dea & Saikia 2021). In contrast, the RM detected in 0954+556 is close to zero, and its m is relatively high, increasing from 3.8% at 1.7 GHz to 7.5% at 8.4 GHz. The negligible RM indicates that the ambient medium surrounding the source is extremely tenuous, limiting strong low-frequency depolarization, and the high m is therefore plausible. Such conditions can arise naturally if the source undergoes recurrent jet activity. In the case that the re-collapse phase of the radio structure lasts longer than the repetition timescale, the repetitive shocks may maintain the interstellar medium warm and efficiently drive the gas out of the host galaxy. The process was described in detail by Czerny et al. (2009) and Siemiginowska et al. (2010). It is more likely that 0954+556 is a CSO maintained by intermittent jet activity rather than a single nascent outburst confined by a dense interstellar medium.
6.3. Jet reorientation in intermittent activity
0954+556 is a rare case in which the jet direction changed with different restarting activity. The earliest jet activity produced the forward linear jet and the corresponding counter-jet, which is 60 degrees north by west at kpc scales, and the discrete component NE, which is northeast of the core. The counter-jet and forward jet with a roughly north–south direction at pc scales were produced by the latest active phase. There are only a few reported cases of jet reorientation, all of which are related to a radio galaxy that harbors a blazar core. PBC J2333.9-2343 is the first one with significant jet reorientation ascribed to restarting activity, transitioning from a radio galaxy to a blazar (Hernández-García et al. 2017, 2023), but no reasonable explanation has been proposed for the new jet pointing toward us by chance. Subsequently, Pajdosz-Śmierciak et al. (2018) revealed that the change in jet direction of B1646+499 could be explained by episodic jet activity with a precessing axis. In addition, a pilot study (Pajdosz-Śmierciak et al. 2022) supplemented seven blazars characterized by exceptionally extended radio galaxy morphology several hundred kiloparsecs in size, and considered that a sudden jet realignment could be caused by a galaxy merger for 5BZU J1238+5325, by jet precession for 5BZU J1345+5332, and by the realignment of the SMBH spin as an interaction with the spin of the outer accretion disk for the remaining sources.
Two broad-line systems were not found in the Sloan Digital Sky Survey (SDSS) optical spectrum of 0954+55610. Nevertheless, this does not definitively exclude the possibility of a recent merger as a binary SMBH system could remain undetected if one SMBH lacks an emission-line system, if the emission lines from the two SMBHs are blended into a single, unusually broad profile, or if obscuration hides one system (Landt et al. 2010). Meanwhile, the morphology of a non-blazar source is not affected strongly by jet precession since a small intrinsic bend is amplified by a small viewing angle along the line of sight. It was proposed that the intermittent jet activity on a timescale between 103 and 106 years is caused by a radiation pressure instability within an accretion disk (Janiuk et al. 2002; Czerny et al. 2009; Siemiginowska et al. 2010). Under this hypothesis, the instability occurs when the radiation pressure overtakes the gas pressure in the disk. The duration of the activity and the outburst separation are determined by the mass of the central BH, the accretion rate, and the viscosity parameter. The virial mass of the supermassive black hole (MBH, vir) in 0954+556 is log(MBH, vir/M⊙) = 8.97 ± 0.4, and its Eddington ratio, which is the ratio of bolometric luminosity to Eddington luminosity, has been estimated to be log(Lbol/LEdd) = − 0.6 by fitting the profile of the MgII line (Shen et al. 2011, The Sloan Digital Sky Survey Data Release 7). In such a condition, if the viscosity parameter takes a value between 0.2 and 0.02, the activity will last for 103 − 104 years and the separation time between episodes will be 104.0 − 104.5 years (Czerny et al. 2009). Considering that our analysis gives the lower limits, the implied timescales are consistent with our results in the order of magnitude.
7. Conclusions
We investigated the complex jet morphology of AGN 0954+556, uncovering evidence of intermittent jet activity. Our analysis, based on multifrequency data from the VLA and VLBA, has led to several significant findings.
1. We prefer to classify 0954+556 as a CSO rather than a FSRQ. In our results, 0954+556 shows weak variability in the radio and γ-ray bands, supported by Fermi-LAT data over 15 years, OVRO archival data over 2 years, and our multiepoch VLBA observations during the period from 2002 to 2016.
2. The potential pc-scale radio core of 0954+556 is most likely near the position of the peak feature in the VLBA maps. This is in agreement with the optical centroid reported by Gaia, the two-sided multifrequency jet structure, the spectral index distribution, and the spine-sheath polarization structure.
3. We detected diffuse counter-jet components at kpc scales and at pc scales and a pc-to-kpc-scale bridge emission region for the first time. This confirms that the forward jet bends about 120 degrees from pc to kpc scales.
4. Together with the spectral index map and spectral aging analysis, we interpreted the apparent complex jet structure from pc to kpc scales as a consequence of at least two episodes of jet activity and rapid jet reorientation on timescales of one million years, which is possibly caused by accretion disk instabilities.
Data availability
The reduced images associated with this work are available at the CDS via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/705/A85
Acknowledgments
We thank the anonymous referee for constructive comments that helped us improve the manuscript. This work supported by the National Key R&D Program of China (2018YFA0404602). This work has made use of data from VLA and VLBA, both operated by The National Radio Astronomy Observatory, which is a facility of the National Science Foundation operated under a cooperative agreement by Associated Universities, Inc., as well as the light curve data from The Fermi Large Area Telescope (LAT) Light Curve Repository (LCR) (Abdollahi et al. 2023), on behalf of NASA and the Fermi Large Area Telescope Collaboration. This research has made use of vlpy (Li et al. 2018) which is a Python package used for VLBI data analysis and NASA/IPAC Extragalactic Database NED which is operated by the JPL, California Institute of Technology, under contract with the National Aeronautics and Space Administration.
References
- Abdollahi, S., Ajello, M., Baldini, L., et al. 2023, ApJS, 265, 31 [NASA ADS] [CrossRef] [Google Scholar]
- Adelman-McCarthy, J. K., Agüeros, M. A., Allam, S. S., et al. 2008, ApJS, 175, 297 [NASA ADS] [CrossRef] [Google Scholar]
- Ahumada, R., Allende Prieto, C., Almeida, A., et al. 2020, ApJS, 249, 3 [NASA ADS] [CrossRef] [Google Scholar]
- Akbar, S., Shah, Z., Misra, R., & Iqbal, N. 2024, ApJ, 977, 111 [Google Scholar]
- An, T., Hong, X. Y., Hardcastle, M. J., et al. 2010, MNRAS, 402, 87 [NASA ADS] [CrossRef] [Google Scholar]
- Ballet, J., Burnett, T. H., Digel, S. W., & Lott, B. 2020, arXiv e-prints [arXiv:2005.11208] [Google Scholar]
- Blandford, R. D., & Königl, A. 1979, ApJ, 232, 34 [Google Scholar]
- Blandford, R. D., & Payne, D. G. 1982, MNRAS, 199, 883 [CrossRef] [Google Scholar]
- Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433 [NASA ADS] [CrossRef] [Google Scholar]
- Blundell, K. M. 2008, ASP Conf. Ser., 386, 467 [Google Scholar]
- Brocksopp, C., Kaiser, C. R., Schoenmakers, A. P., & de Bruyn, A. G. 2007, MNRAS, 382, 1019 [NASA ADS] [CrossRef] [Google Scholar]
- Bruni, G., Gómez, J. L., Vega-García, L., et al. 2021, A&A, 654, A27 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Chavan, K., Dabhade, P., & Saikia, D. J. 2023, MNRAS, 525, L87 [NASA ADS] [CrossRef] [Google Scholar]
- Cheung, C. C. 2007, AJ, 133, 2097 [NASA ADS] [CrossRef] [Google Scholar]
- Croke, S. M., & Gabuzda, D. C. 2008, MNRAS, 386, 619 [CrossRef] [Google Scholar]
- Czerny, B., Siemiginowska, A., Janiuk, A., Nikiel-Wroczyński, B., & Stawarz, Ł. 2009, ApJ, 698, 840 [NASA ADS] [CrossRef] [Google Scholar]
- Dennett-Thorpe, J., Bridle, A. H., Scheuer, P. A. G., Laing, R. A., & Leahy, J. P. 1997, MNRAS, 289, 753 [NASA ADS] [Google Scholar]
- Edelson, R., Turner, T. J., Pounds, K., et al. 2002, ApJ, 568, 610 [CrossRef] [Google Scholar]
- Fan, J. H., Liu, Y., Yuan, Y. H., et al. 2007, A&A, 462, 547 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Fromm, C. M., Ros, E., Perucho, M., et al. 2013, A&A, 557, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Gabuzda, D. C. 2021, Galaxies, 9, 58 [NASA ADS] [CrossRef] [Google Scholar]
- Gabuzda, D. C., Reichstein, A. R., & O’Neill, E. L. 2014, MNRAS, 444, 172 [Google Scholar]
- Gaia Collaboration (Vallenari, A., et al.) 2023, A&A, 674, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Gopal-Krishna, Biermann, P. L., & Wiita, P. J. 2003, ApJ, 594, L103 [NASA ADS] [CrossRef] [Google Scholar]
- Greisen, E. W. 2003, Astrophys. Space Sci. Lib., 285, 109 [NASA ADS] [Google Scholar]
- Grzędzielski, M., Janiuk, A., Czerny, B., & Wu, Q. 2017, A&A, 603, A110 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Hatziminaoglou, E., Siemiginowska, A., & Elvis, M. 2001, ApJ, 547, 90 [Google Scholar]
- Hernández-García, L., Panessa, F., Giroletti, M., et al. 2017, A&A, 603, A131 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Hernández-García, L., Panessa, F., Bruni, G., et al. 2023, MNRAS, 525, 2187 [CrossRef] [Google Scholar]
- Homan, D. C., Lister, M. L., Kovalev, Y. Y., et al. 2015, ApJ, 798, 134 [Google Scholar]
- Homan, D. C., Cohen, M. H., Hovatta, T., et al. 2021, ApJ, 923, 67 [NASA ADS] [CrossRef] [Google Scholar]
- Hota, A., Sirothia, S. K., Ohyama, Y., et al. 2011, MNRAS, 417, L36 [NASA ADS] [CrossRef] [Google Scholar]
- Hovatta, T., Lister, M. L., Aller, M. F., et al. 2012, AJ, 144, 105 [Google Scholar]
- Hovatta, T., Aller, M. F., Aller, H. D., et al. 2014, AJ, 147, 143 [Google Scholar]
- Hu, X.-Z., Hong, X., Zhao, W., et al. 2024a, ApJ, 965, 74 [Google Scholar]
- Hu, X.-Z., Hong, X., Zhao, W., et al. 2024b, ApJ, 974, 111 [Google Scholar]
- Jaffe, W. J., & Perola, G. C. 1973, A&A, 26, 423 [NASA ADS] [Google Scholar]
- Janiuk, A., Czerny, B., & Siemiginowska, A. 2002, ApJ, 576, 908 [NASA ADS] [CrossRef] [Google Scholar]
- Kardashev, N. S. 1962, Soviet Ast., 6, 317 [Google Scholar]
- Kiehlmann, S., Lister, M. L., Readhead, A. C. S., et al. 2024, ApJ, 961, 240 [NASA ADS] [CrossRef] [Google Scholar]
- Kim, J.-Y., & Trippe, S. 2014, J. Korean Astron. Soc., 47, 195 [Google Scholar]
- Kovalev, Y. Y., Lobanov, A. P., Pushkarev, A. B., & Zensus, J. A. 2008, A&A, 483, 759 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Landt, H., Cheung, C. C., & Healey, S. E. 2010, MNRAS, 408, 1103 [NASA ADS] [CrossRef] [Google Scholar]
- Lee, S.-S., Lobanov, A. P., Krichbaum, T. P., et al. 2008, AJ, 136, 159 [Google Scholar]
- Li, X., Mohan, P., An, T., et al. 2018, ApJ, 854, 17 [NASA ADS] [CrossRef] [Google Scholar]
- Liodakis, I., Hovatta, T., Huppenkothen, D., et al. 2018, ApJ, 866, 137 [Google Scholar]
- Liu, R., & Pooley, G. 1991, MNRAS, 249, 343 [NASA ADS] [Google Scholar]
- McConville, W., Ostorero, L., Moderski, R., et al. 2011, ApJ, 738, 148 [Google Scholar]
- Merritt, D., & Ekers, R. D. 2002, Science, 297, 1310 [NASA ADS] [CrossRef] [Google Scholar]
- Miley, G. 1980, ARA&A, 18, 165 [Google Scholar]
- Murgia, M., Fanti, C., Fanti, R., et al. 1999, A&A, 345, 769 [NASA ADS] [Google Scholar]
- O’Dea, C. P., & Saikia, D. J. 2021, A&ARv, 29, 3 [Google Scholar]
- Orienti, M., Dallacasa, D., & Stanghellini, C. 2007, A&A, 461, 923 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Pacholczyk, A. G. 1970, Radio Astrophysics. Nonthermal Processes in Galactic and Extragalactic Sources (Freeman) [Google Scholar]
- Pajdosz-Śmierciak, U., Jamrozy, M., Soida, M., & Stawarz, Ł. 2018, ApJ, 868, 64 [Google Scholar]
- Pajdosz-Śmierciak, U., Śmierciak, B., & Jamrozy, M. 2022, MNRAS, 514, 2122 [Google Scholar]
- Petrov, L. Y., & Kovalev, Y. Y. 2025, ApJS, 276, 38 [NASA ADS] [CrossRef] [Google Scholar]
- Planck Collaboration VI. 2020, A&A, 641, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Pushkarev, A. B., Gabuzda, D. C., Vetukhnovskaya, Y. N., & Yakimov, V. E. 2005, MNRAS, 356, 859 [Google Scholar]
- Quici, B., Turner, R. J., Seymour, N., et al. 2022, MNRAS, 514, 3466 [NASA ADS] [CrossRef] [Google Scholar]
- Rani, B., Jorstad, S. G., Marscher, A. P., et al. 2018, ApJ, 858, 80 [Google Scholar]
- Rastello, S., Dallacasa, D., & Orienti, M. 2016, Astron. Nachr., 337, 42 [Google Scholar]
- Readhead, A. C. S., Taylor, G. B., Xu, W., et al. 1996, ApJ, 460, 612 [NASA ADS] [CrossRef] [Google Scholar]
- Reid, A., Shone, D. L., Akujor, C. E., et al. 1995, A&AS, 110, 213 [NASA ADS] [Google Scholar]
- Richards, J. L., Hovatta, T., Max-Moerbeck, W., et al. 2014, MNRAS, 438, 3058 [Google Scholar]
- Richstone, D., Ajhar, E. A., Bender, R., et al. 1998, Nature, 385, A14 [Google Scholar]
- Rossetti, A., Mantovani, F., Dallacasa, D., Fanti, C., & Fanti, R. 2005, A&A, 434, 449 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Schleicher, B., Arbet-Engels, A., Baack, D., et al. 2019, Galaxies, 7, 62 [NASA ADS] [CrossRef] [Google Scholar]
- Schoenmakers, A. P., de Bruyn, A. G., Röttgering, H. J. A., van der Laan, H., & Kaiser, C. R. 2000, MNRAS, 315, 371 [Google Scholar]
- Shen, Y., Richards, G. T., Strauss, M. A., et al. 2011, ApJS, 194, 45 [Google Scholar]
- Shepherd, M. C. 1997, ASP Conf. Ser., 125, 77 [Google Scholar]
- Siemiginowska, A., Czerny, B., & Kostyunin, V. 1996, ApJ, 458, 491 [Google Scholar]
- Siemiginowska, A., Czerny, B., Janiuk, A., et al. 2010, ASP Conf. Ser., 427, 326 [Google Scholar]
- Simard-Normandin, M., Kronberg, P. P., & Button, S. 1981, ApJS, 45, 97 [NASA ADS] [CrossRef] [Google Scholar]
- Singh, V., Ishwara-Chandra, C. H., Kharb, P., Srivastava, S., & Janardhan, P. 2016, ApJ, 826, 132 [NASA ADS] [CrossRef] [Google Scholar]
- Sokolovsky, K. V., Kovalev, Y. Y., Pushkarev, A. B., & Lobanov, A. P. 2011, A&A, 532, A38 [CrossRef] [EDP Sciences] [Google Scholar]
- Soldi, S., Türler, M., Paltani, S., et al. 2008, A&A, 486, 411 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Tremblay, S. E., Taylor, G. B., Ortiz, A. A., et al. 2016, MNRAS, 459, 820 [NASA ADS] [CrossRef] [Google Scholar]
- Turner, R. J. 2018, MNRAS, 476, 2522 [CrossRef] [Google Scholar]
- Turner, R. J., Shabala, S. S., & Krause, M. G. H. 2018, MNRAS, 474, 3361 [NASA ADS] [CrossRef] [Google Scholar]
- Urry, C. M., & Padovani, P. 1995, PASP, 107, 803 [NASA ADS] [CrossRef] [Google Scholar]
- Vaughan, S., Edelson, R., Warwick, R. S., & Uttley, P. 2003, MNRAS, 345, 1271 [Google Scholar]
- Véron-Cetty, M. P., & Véron, P. 2006, A&A, 455, 773 [Google Scholar]
- Wills, B. J., Thompson, K. L., Han, M., et al. 1995, ApJ, 447, 139 [NASA ADS] [CrossRef] [Google Scholar]
- Xu, M. H., Lunz, S., Anderson, J. M., et al. 2021, A&A, 647, A189 [EDP Sciences] [Google Scholar]
Appendix A: Extra figures
![]() |
Fig. A.1. Subtraction maps to assess structural changes in 0954+556. The images in Figure 3 were convolved to a common restoring beam, corresponding to the beam from the observation on 13 June 2016. An average image was generated by taking the mean of the four epochs (2002, 2004, 2005, and 2016). Each individual image was then subtracted from the average image to highlight deviations. The color bar represents the intensity in units of mJy/b. The subtraction maps shows the intensity differences range approximately from -15 mJy/b to 15 mJy/b, indicating a stable structure of 0954+556. |
![]() |
Fig. A.2. Images used when analyzing radiative age at kpc scale. The contours are drawn at -1, 1, 2, 4, 8, ..., of the first contour level (3σrms). The synthesized beams are plotted in the bottom left corner of each image. Three components used in the spectral analysis, i.e., NE, core, and NW, are distinguished by two gray dashed lines. |
![]() |
Fig. A.3. Images used when analyzing radiative age at pc scale. The contours are drawn at -1, 1, 2, 4, 8, ..., of the first contour level (3σrms). The synthesized beams are plotted in the bottom left corner of each image. The division criteria for the north and south components are represented by a gray dashed line. |
All Tables
All Figures
![]() |
Fig. 1. VLA flux density images of 0954+556 at kpc scale. The contours are drawn at –1, 1, 2, 4, 8, …, of the first contour level (3σrms). The synthesized beams are plotted in the bottom left corner of each image. In the top left panel, the red circles overlapping the contours present the results of modelfit. In the top middle panel, the 5 GHz image shows a weak counter-jet. To better reveal the weak counter-jet, the 5 GHz image was restored with a circular beam of 900 mas, as shown in the top right panel. |
| In the text | |
![]() |
Fig. 2. VLBA flux density images of 0954+556 at pc scale. The contours are drawn at –1, 1, 2, 4, 8, …, of the first contour level (3σrms). The synthesized beams are also plotted in the bottom left corner of each image. The image in the top right corner is the tapered image at 1.7 GHz (UVTAPER 0.5, 10). |
| In the text | |
![]() |
Fig. 3. VLBA flux density images of 0954+556 at 5 GHz. Left panel: Reconstructed image of the visibilities from a combination of data on 20 July 2002 and 22 September 2002. Middle panel: Reconstructed image of the visibilities from a combination of data on 14 June 2004 and 3 November 2004. Right panel: Reconstructed image of the visibilities from a combination of data on 4 April 2005 and 18 August 2005. The contours are drawn at –1, 1, 2, 4, 8, …, of the first contour level (3σrms). The synthesized beams are plotted in the bottom left corner of each image. |
| In the text | |
![]() |
Fig. 4. Light curve of 0954+556. The top panel shows a Fermi-LAT γ-ray light curve of 0954+556 from 8 August 2008 to 6 December 2024, divided into seven-day bins. All points are plotted along with their statistical errors. The bottom panel is the radio light curve obtained from OVRO at 15 GHz, from 19 April 2009 to 30 December 2011, with an interval of four days. The gray area in the top panel corresponds to the bottom panel. |
| In the text | |
![]() |
Fig. 5. Spectral index distribution maps of 0954+556. From left to right, the panels show the distribution of spectral index between 1.4 and 5 GHz at kpc scales, between 5 and 8.4 GHz at kpc scales, between 1.7 and 5 GHz at pc scales, and between 5 and 8.4 GHz at pc scales. The overlapped contour images are a VLA image at 1.4 GHz with a uv range 5.1 − 168.7 Kλ, a VLA image at 5 GHz with a uv range 11.3 − 561.2 Kλ, a VLBA image at 1.7 GHz with a uv range 2.3 − 48.0 Mλ, and a VLBA image at 5 GHz with a uv range 5.4 − 39.5 Mλ, respectively. All counters start at 3σrms and increase by a factor of 2. The color represents the value of spectral index α. Pixels with flux density below 5.4σrms, alongside those with spectral index uncertainty exceeding 40%, are blanked. |
| In the text | |
![]() |
Fig. 6. Four different fit results for each component of 0954+556. The red, blue, black, and yellow dashed lines represent the JP, KP, CI-on, and CI-off models respectively. The VLA core and component N are fitted by a simple power law. For component N, the data at 43 GHz was not used in the spectral analysis as mentioned in Section 2.4. |
| In the text | |
![]() |
Fig. 7. Polarization images at different frequencies. The restoring beam is indicated by a gray ellipse in the lower left corner of each image. The contours are drawn at –1, 1, 2, 4, 8, …, of the first contour level (3σrms). Superimposed on the contours of flux density are vectors whose lengths are proportional to the polarized intensity; the directions indicate the direction of the E-vectors, and the colors show the fractional linear polarization. In all polarization images, the 1 mas line is 1 mJy/b. We clipped pixels with I < 3σI, P < 3σP, and m < 3σm. |
| In the text | |
![]() |
Fig. A.1. Subtraction maps to assess structural changes in 0954+556. The images in Figure 3 were convolved to a common restoring beam, corresponding to the beam from the observation on 13 June 2016. An average image was generated by taking the mean of the four epochs (2002, 2004, 2005, and 2016). Each individual image was then subtracted from the average image to highlight deviations. The color bar represents the intensity in units of mJy/b. The subtraction maps shows the intensity differences range approximately from -15 mJy/b to 15 mJy/b, indicating a stable structure of 0954+556. |
| In the text | |
![]() |
Fig. A.2. Images used when analyzing radiative age at kpc scale. The contours are drawn at -1, 1, 2, 4, 8, ..., of the first contour level (3σrms). The synthesized beams are plotted in the bottom left corner of each image. Three components used in the spectral analysis, i.e., NE, core, and NW, are distinguished by two gray dashed lines. |
| In the text | |
![]() |
Fig. A.3. Images used when analyzing radiative age at pc scale. The contours are drawn at -1, 1, 2, 4, 8, ..., of the first contour level (3σrms). The synthesized beams are plotted in the bottom left corner of each image. The division criteria for the north and south components are represented by a gray dashed line. |
| 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.











![$$ \begin{aligned} \begin{aligned} B_{\mathrm{eq}} = 5.69 \times 10^{-5}\left[\frac{(1+\mathrm{k})}{\eta }(1+\mathrm{z})^{3-\alpha }\right. \frac{1}{\theta _x \theta _y\ \mathrm{s} \sin ^{3 / 2} \phi } \\ \left.\times \frac{F_0}{v_0^\alpha } \frac{v_2^{\alpha +1 / 2}-v_1^{\alpha +1 / 2}}{\alpha +\frac{1}{2}}\right]^{2 / 7} \text{ gauss}. \end{aligned} \end{aligned} $$](/articles/aa/full_html/2026/01/aa53625-24/aa53625-24-eq8.gif)





