Event Horizon Telescope observations of the jet launching and collimation in Centaurus A

Very-long-baseline interferometry (VLBI) observations of active galactic nuclei at millimetre wavelengths have the power to reveal the launching and initial collimation region of extragalactic radio jets, down to 10–100 gravitational radii (rg ≡ GM/c2) scales in nearby sources1. Centaurus A is the closest radio-loud source to Earth2. It bridges the gap in mass and accretion rate between the supermassive black holes (SMBHs) in Messier 87 and our Galactic Centre. A large southern declination of −43° has, however, prevented VLBI imaging of Centaurus A below a wavelength of 1 cm thus far. Here we show the millimetre VLBI image of the source, which we obtained with the Event Horizon Telescope at 228 GHz. Compared with previous observations3, we image the jet of Centaurus A at a tenfold higher frequency and sixteen times sharper resolution and thereby probe sub-lightday structures. We reveal a highly collimated, asymmetrically edge-brightened jet as well as the fainter counterjet. We find that the source structure of Centaurus A resembles the jet in Messier 87 on ~500 rg scales remarkably well. Furthermore, we identify the location of Centaurus A’s SMBH with respect to its resolved jet core at a wavelength of 1.3 mm and conclude that the source’s event horizon shadow4 should be visible at terahertz frequencies. This location further supports the universal scale invariance of black holes over a wide range of masses5,6. The millimetre image of the Centaurus A nucleus by the Event Horizon Telescope reveals a highly collimated, asymmetrically edge-brightened jet. The source’s event horizon shadow should be visible at terahertz frequencies, consistent with the universal scale invariance of black holes.

Very-long-baseline interferometry (VLBI) observations of active galactic nuclei at millimeter wavelengths have the power to reveal the launching and initial collimation region of extragalactic radio jets, down to 10 -100 gravitational radii (r g ≡ GM/c 2 ) scales in nearby sources 1 .Centaurus A is the closest radio-loud source to Earth 2 .It bridges the gap in mass and accretion rate between the supermassive black holes (SMBHs) in Messier 87 and our galactic center.A large southern declination of −43°has however prevented VLBI imaging of Centaurus A below λ 1 cm thus far.Here, we show the millimeter VLBI image of the source, which we obtained with the Event Horizon Telescope at 228 GHz.Compared to previous observations 3 , we image Centaurus A's jet at a tenfold higher frequency and sixteen times sharper resolution and thereby probe sub-lightday structures.We reveal a highly-collimated, asymmetrically edge-brightened jet as well as the fainter counterjet.We find that Centaurus A's source structure resembles the jet in Messier 87 on ∼ 500 r g scales remarkably well.Furthermore, we identify the location of Centaurus A's SMBH with respect to its resolved jet core at λ 1.3 mm and conclude that the source's event horizon shadow 4 should be visible at THz frequencies.This location further supports the universal scale invariance of black holes over a wide range of masses 5,6 .
In this work, we present the first image of Centaurus A (Cen A) obtained by the Event Horizon Telescope (EHT) with a nominal resolution of 25 microarcseconds (µas) at λ 1.3 mm.For a black hole mass of (5.5 ± 3) × 10 7 M 7 , we are probing jet structures down to scales of ∼ 200 r g ≈ 0.6 light days.It has recently become possible to model these scales with sophisticated general relativistic magnetohydrodynamics (GRMHD) simulations 8 , where jet ejection and their symbiotic relationship with accretion flows are simulated from first principles.We have observed Cen A in a six-hour-long track on April 10, 2017.The EHT, as a novel and heterogeneous high-frequency VLBI array, poses unique calibration challenges.To obtain robust results, independent of assumptions made during the data calibration, we base our scientific analysis on two data sets, which we obtained from two independent calibration pipelines: rPICARD 9 and EHT-HOPS 10 (Methods 1.1.2).
Figure 1 presents our reconstruction of the jet image structure derived from the EHT data using a regularized maximum-likelihood method, next to the large-scale source morphology and the similarly edge-brightened morphology of the M 87 jet on comparable gravitational scales.These images are convolved with Gaussian beams set by their respective nominal instrumental resolutions, as per standard practice in radio-interferometric imaging, to suppress possibly spurious fine-scale structures in the image model.The brightness temperatures T [K] shown are related to flux densities S in jansky (Jy) through the observing wavelength λ, Boltzmann constant k B , and angular resolution element Ω as T = λ 2 (2k B Ω) −1 S. The λ 1.3 mm Cen A jet has a narrow, collimated profile and exhibits one-sidedness, pronounced edge-brightening, and a northwest-southeast (NW-SE) brightness asymmetry.The approaching jet extends towards the north-east and the faint counterjet is directed south-westwards.The total compact flux density in our image is ∼ 2 Jy.The identification of the jet apex and black hole position (Methods 1.4.1) is shown in the unconvolved image model of Figure 2. We can use interferometer data with a high signal-to-noise ratio to super-resolve image features beyond the nominal resolution of the instrument.We therefore base our analysis on the robust features of the unconvolved image model.We have verified the robustness of the counterjet feature with synthetic data studies (Supplementary Fig. 4).The estimated jet position angle (PA) on the sky of 48°± 5°agrees with cm-wave VLBI observations 3 .The cm-band data also constrain the inclination angle of the jet axis with respect to our line of sight to θ ∼ 12°-45°, assuming the jet does not bend along the line of sight.
The Cen A λ 1.3 mm jet exhibits three types of brightness asymmetries (R): between the jet and counterjet, the sheath and spine, and the NW vs. SE ridgelines (Methods 1.4.2).We take the two bright radiating streams of the approaching-and counter-jet as jet 'arms' and denote the maximum intensity region along each arm as 'ridgeline'.The jet-to-counterjet intensity ratio R j/cj can naturally be explained for a relativistic outflow with θ = 90°, where jet emission will be Doppler boosted and counterjet emission de-boosted.We find R j/cj 5, which is in agreement with cm-wave VLBI observations 3 and suggests that the initial acceleration of the jet occurs within the inner collimation region imaged in this study.
There is no jet spine emission in our image.With synthetic data studies, we found that spine emission exceeding ∼ 20 % of the sheath radiation intensity would be detectable, i.e., R sh/sp > 5 (Methods 1.3).The intensities of the brightest, central SE and NW jet components in the unconvolved image are (32 ± 8) × 10 9 K and (20 ± 4) × 10 9 K, respectively.The brightness ratio between these components follows as R s/n = 1.6 ± 0.5.
The collimation profile of the jet width W follows a narrow expansion profile with distance to the apex z as W ∝ z k with k = 0.33 ± 0.05| stat ± 0.06| sys (Figure 3).Resolution and potentially optical depth effects prevent us from pinning down the jet opening angle ψ jet at small z, where the jet converges towards the apex.We denote the boundary between the inner convergence region and the outer jet with a clearly defined collimation and easily traceable jet ridgelines as z col .For the brighter and straighter SE arm, we have W (z col ≈ 32 µas) ≈ 25 µas, i.e., the brightest jet component marks the boundary between the convergence and strongly collimated regions here (Figure 2).If we assume the two jet ridgelines to meet at the apex, we find ψ jet 40°as a conservative estimate.Factoring in the range of possible θ values yields ψ int 10°-30°for the intrinsic, deprojected opening angle (Methods 1.4.3).
The M 87 11 (NGC 4486, 3C 274, Virgo A), Markarian 501 12 , and restarted 3C 84 jets 13 also show strong edgebrightening and large initial opening angles on comparable scales seen at similar inclination angles of ∼ 18°.The expansion profile of Cen A lies in between the parabolic profile of M 87 (k = 0.5) and the almost cylindrical profile of 3C 84 (k = 0.2), which implies a strong confinement of the 3C 84 jet by a shallow pressure gradient from the ambient medium.For the inner Cen A jet, this suggests strong magnetic collimation or the presence of external pressure and density gradients of P ext ∝ z −4k = z −1.3 and ρ ext ∝ z 1−4k = z −0.3 (Methods 1.4.4).Radiatively inefficient accretion flows alone, which are expected to operate in the M 87, 3C 84, and Cen A sub-Eddington lowluminosity AGN (LLAGN) sources, have comparatively steeper pressure and density gradients 14 .This may indicate the presence of winds, which are likely to be launched by this type of accretion flow.The noticeable similarity and prominence of edge-brightened jet emission in M 87, 3C 84, and Cen A suggests the dominance of jet sheath emission to be an emerging feature in LLAGN.In GRMHD simulations, the sheath manifests itself as interaction region between an accretion-powered outflow 15 and the fast jet spine, which is potentially powered by the black hole spin 16 .The mass-loaded sheath has a higher intrinsic emissivity compared to the evacuated spine.The same type of LLAGN-applicable GRMHD simulations also self-consistently develop a collimating helical magnetic field structure in the jet, which is confirmed observationally in many active galactic nuclei (AGN) 17 .The dominating sheath emissivity and helical magnetic field structure provides a natural intrinsic explanation for the prevailing edge-brightening in LLAGN and can also explain the NW-SE brightness asymmetry.This model and alternative geometric explanations for the brightness asymmetries are discussed in the Methods 1.4.2section.
The basic radiative properties of these jets can be analytically understood with a simple model 18 , where particle and magnetic energy density equipartition is assumed, while the particle density decays with z −2 .Under these conditions, an optically thick and self-absorbed compact feature is expected (the core), whose position z core along the jet is frequency dependent with z core ∝ ν − 1 19,20 .This radio core corresponds to the photosphere, where the optical depth τ (ν) to photons at the observing frequency ν is unity.The jet is optically thick upstream and optically thin downstream.The photosphere moves closer to the jet apex at higher frequencies, until the point where either the launching point is reached near the horizon, or particle acceleration has not yet begun 21 .The scale of a jet 'nozzle' emission cannot be smaller than the √ 27 r g photon capture radius (Methods 1.4.5).The combination of all emission regions along the jet gives rise to a flat to inverted radio spectrum, peaking at a maximum frequency ν, determined by the black hole mass M and accretion rate Ṁ and scaling as ν ∝ 17 6,19,22 .Here, D is the distance of the black hole to the observer and F r the observed radio flux density.These scaling relations follow from the assumption that the jet's internal gas and magnetic pressures are linearly coupled to the accretion rate and maintain a fixed ratio along the jets.The proportionality constant between Ṁ and F r generally depends on the jet's velocity, electron and magnetic energy densities, particle distribution spectrum, and inclination angle.Therefore, we are only able to make a first order estimates.It should further be noted that X-ray binary observations 23 have revealed a more complex relationship between ν and Ṁ , where the innermost particle acceleration zone in the jet may not remain stationary and source-specific accretion disk parameters come into play.The same effects are expected to also influence ν in AGN, which substantiates the fact that only order of magnitude estimates can be provided for ν.We assume the brightest features in our image to correspond to the radio cores at 230 GHz, which is discussed in Methods 1.5.Our assumption is affirmed by three consistent and independent measurements of ν, but future spectral information is needed for a definitive confirmation.We show that ν lies in the THz regime for Cen A based on the core shift that we can determine from our image, scaling relations with the M 87 jet, and the spectral energy distribution of Cen A.
We take the distance from the brightest pixel in the image to the estimated position of the jet apex and obtain a core shift of z core = 32 ± 11 µas.Based on this distance and the uncertain inclination angle, we estimate that an observing frequency of ν CenA ∼ 10 -60 THz (Methods 1.4.5) will reach the base of the jet at the black hole innermost stable circular orbit (ISCO).A caveat is that we do not take the effect of the uncertain ambient medium into account in this simple picture.
Independently, we can use the above scaling relations to estimate the order of magnitude of ν CenA by comparing the Cen A jet with M 87, which has ν M87 = 228 GHz 24,25 .For the cm jet radio core, a flux density of ∼ 1 Jy is measured for both sources 3,26 , which yields ṀCenA ∼ 0.1 ṀM87 for the accretion onto the black hole and therefore ν CenA ∼ 26 ν M87 ≈ 6 THz (Methods 1.4.5), in agreement with our observations and the assumed position of the black hole at the jet apex within an order of magnitude.Based on comparable jet velocities (∼ 0.3 c -0.5 c) and inclination angles (∼ 20°), we have assumed the amount of Doppler boosting to be similar in both jets.The relation of accretion rates would constrain ṀCenA to be 9 × 10 −5 M yr −1 27 or 7 × 10 −5 ṀEdd in terms of the Eddington accretion rate for an assumed radiative efficiency of 10 %.
The core spectral energy distribution of Cen A peaks at ∼ 10 13 Hz 28 , which may be the equivalent of the submillimeter bump seen in Sgr A* 29,30 , and would further support our hypothesis.
Observed correlations between the masses of accreting black holes and their X-ray and radio luminosities form the basis of a unified fundamental plane of scale-invariant black hole accretion.This scale-invariance has been derived based on stellar-mass black holes, which have a break frequency ν XRB in the near-infrared, and supermassive 10 8 M -10 10 M AGN, where ν SMBH lies in the radio to submillimeter regime 5,6 .With our observation, we demonstrate that the simple fundamental relations for the black hole jet activity still holds for a source with a mass of 5.5 × 10 7 M and ν CenA in the THz regime, in between those two types of black holes.Our method used to determine the optimal frequency to observe black hole shadows based on core shift, jet power, and source spectrum is in principle applicable to any LLAGN.
Our findings suggest that the black hole shadow 4 of Cen A would be visible in a bright, optically thin accretion flow at an observing frequency of a few THz.At this high frequency, a VLBI experiment above the Earth's troposphere would be able to resolve the 1.4±0.8µas shadow diameter with a minimal baseline length of ∼ 8000 km.
resource at PTI and TACC was allocated through AST170028.The simulations were performed in part on the SuperMUC cluster at the LRZ in Garching, on the LOEWE cluster in CSC in Frankfurt, and on the HazelHen cluster at the HLRS in Stuttgart.This research was enabled in part by support provided by Compute Ontario (http://computeontario.ca), Calcul Quebec (http://www.calculquebec.ca) and Compute Canada (http://www.computecanada.ca).We thank the staff at the participating observatories, correlation centers, and institutions for their enthusiastic support.
This paper makes use of the following ALMA data: ADS/JAO.ALMA#2016.1.01198.V. ALMA is a partnership of the European Southern Observatory (ESO; Europe, representing its member states), NSF, and National Institutes of Natural Sciences of Japan, together with National Research Council (Canada), Ministry of Science and Technology (MOST; Taiwan), Academia Sinica Institute of Astronomy and Astrophysics (ASIAA; Taiwan), and Korea Astronomy and Space Science Institute (KASI; Republic of Korea), in cooperation with the Republic of Chile.The Joint ALMA Observatory is operated by ESO, Associated Universities, Inc. (AUI)/NRAO, and the National Astronomical Observatory of Japan (NAOJ).The NRAO is a facility of the NSF operated under cooperative agreement by AUI.APEX is a collaboration between the Max-Planck-Institut für Radioastronomie (Germany), ESO, and the Onsala Space Observatory (Sweden).The SMA is a joint project between the SAO and ASIAA and is funded by the Smithsonian Institution and the Academia Sinica.The JCMT is operated by the East Asian Observatory on behalf of the NAOJ, ASIAA, and KASI, as well as the Ministry of Finance of China, Chinese Academy of Sciences, and the National Key R&D Program (No. 2017YFA0402700) of China.Additional funding support for the JCMT is provided by the Science and Technologies Facility Council (UK) and participating universities in the UK and Canada.The LMT is a project operated by the Instituto Nacional de Astrófisica, Óptica, y Electrónica (Mexico) and the University of Massachusetts at Amherst (USA), with financial support from the Consejo Nacional de Ciencia y Tecnología and the National Science Foundation.The IRAM 30-m telescope on Pico Veleta, Spain is operated by IRAM and supported by CNRS (Centre National de la Recherche Scientifique, France), MPG (Max-Planck-Gesellschaft, Germany) and IGN (Instituto Geográfico Nacional, Spain).The SMT is operated by the Arizona Radio Observatory, a part of the Steward Observatory of the University of Arizona, with financial support of operations from the State of Arizona and financial support for instrumentation development from the NSF.The SPT is supported by the National Science Foundation through grant PLR-1248097.Partial support is also provided by the NSF Physics Frontier Center grant PHY-1125897 to the Kavli Institute of Cosmological Physics at the University of Chicago, the Kavli Foundation and the Gordon and Betty Moore Foundation grant GBMF 947.The SPT hydrogen maser was provided on loan from the GLT, courtesy of ASIAA.The EHTC has received generous donations of FPGA chips from Xilinx Inc., under the Xilinx University Program.The EHTC has benefited from technology shared under open-source license by the Collaboration for Astronomy Signal Processing and Electronics Research (CASPER).The EHT project is grateful to T4Science and Microsemi for their assistance with Hydrogen Masers.This research has made use of NASA's Astrophysics Data System.We gratefully acknowledge the support provided by the extended staff of the ALMA, both from the inception of the ALMA Phasing Project through the observational campaigns of 2017 and 2018.We would like to thank A. Deller and W. Brisken for EHT-specific support with the use of DiFX.We acknowledge the significance that Maunakea, where the SMA and JCMT EHT stations are located, has for the indigenous Hawaiian people.
The grants listed above collectively fund the Event Horizon Telescope project.    1 Methods

Processing of observational data
This section describes the 2017 EHT observations of Cen A, the model-independent calibration 9 with two separate pipelines, the flux density calibration, and known measurement issues and systematics with corresponding mitigation strategies.The final data sets coming out from the two pipelines are both used for the scientific analysis as crossverification (Supplementary Fig. 5).

Data acquisition
Cen A (PKS 1322−428, hosted in the NGC 5128 elliptical galaxy, α J2000 = 13 h 25 m 27.62 s , δ J2000 = −43°1 8.81 ) was observed by the EHT in a six hour long track on April 10 2017, with a total on-source integration time of 105 minutes (Supplementary Fig. 6).The observations were carried out by the Atacama Large Millimeter/submillimeter Array (ALMA), Atacama Pathfinder Experiment (APEX), James Clerk Maxwell Telescope (JCMT), Large Millimeter Telescope Alfonso Serrano (LMT), South Pole Telescope (SPT) 32 , Submillimeter Array (SMA), and Submillimeter Telescope (SMT) 33 .For ALMA, 37 of the 12 m dishes were phased-up 34 .The EHT IRAM 30 m Telescope is not able to see Cen A jointly with rest of the array due to the low declination of the source.The data were recorded on Mark 6 VLBI recorders 35 with 2-bit sampling in two 2 GHz wide bands, 'low' and 'high', centred around 227.1 GHz and 229.1 GHz, respectively.Unless stated otherwise, results are derived using the combined low+high band data.Quarter-wave plates at each site except ALMA were used to observe circularly polarized light.The data were correlated with the DiFX software 36,37 .The PolConvert 38 software was used to convert the phased ALMA data 34 from a linear polarization basis to a circular basis after correlation, based on solutions from the calibration of the connected-element ALMA data 39 .

Data reduction pipelines
The autocorrelation normalization, feed angle rotation, fringe fitting, bandpass calibration, and a priori correction of atmospheric phase turbulence 40 were performed independently by the CASA-based rPICARD 9,41 and EHT-HOPS 10,42 pipelines.DiFX produces FITS-IDI and Mark4 data.rPICARD uses the FITS-IDI product and converts it into the Measurement Set format.EHT-HOPS uses the Mark4 data.Both software packages convert the calibrated data into the UVFITS format for further processing.
rPICARD performs an upstream correction for the feed rotation angle and uses station-based global fringe fitting based on an unpolarized point source model to correct for phases, delays, and rates consistently for the right circular polarization and left circular polarization signal paths 43 .Atmospheric phase and residual delay variations are corrected within the expected coherence time by fringe fitting segmented data of each VLBI scan.The segmentation length is set by the signal-to-noise ratio (S/N) of each baseline.
For EHT-HOPS, the feed rotation angle is corrected after the fringe-fitting together with an additional polarization calibration step, where complex R-L gain offsets are solved for.Delays and rates are found in a baseline-based fringe search and referenced to individual stations with a least-squares optimization 44 .Atmospheric phases are corrected by fitting a polynomial phase model to the data on baselines to the most sensitive reference station in each scan.A round-robin approach is used to avoid fitting to thermal noise and the degree of the polynomial is set by the S/N of the data.

Gain amplitude calibration
The flux density calibration is done based on determined station sensitivities in a common framework for the rPICARD and EHT-HOPS data 45 .The sensitivity of a station i is given by its system equivalent flux density (SEFD i ) in Jy, which takes into account the gain and total noise power along a telescope's signal chain as a function of time t and frequency ν.On a baseline i-j, correlation coefficients ξ (i,j) in units of thermal noise are calibrated to a physical radiation intensity scale of correlated flux density S (i,j) through where η Q is the quantization efficiency.For data recorded with 2-bit sampling, we have η Q ≈ 0.88.The gains of co-located stations were solved based on a contemporaneous measurement of the total flux density S 0 = 5.62 Jy of the source with the ALMA interferometer 10,39,45 .The correlated flux S measured between two co-located sites p and q should be equal to S 0 and for a third station o, we should have S po = S qo .It follows, that we can solve for the station-based amplitude gains A of p and q with a self-calibration approach.Here, the model is given by the constant flux density S 0 seen by baselines between co-located sites.No gain corrections for non-co-located ('isolated') stations are solved for.
Ad hoc correction factors are used to correct signal losses at APEX due to an injected instrumental signal and at SMA due to temporary losses of bandwidth 45 .Additionally, LMT and SPT suffered from pointing problems, which result in significant amplitude variations between and within VLBI scans.These losses cannot be estimated a priori and must be corrected with self-calibration gain solutions derived within short ∼ 10 s segments from high S/N data.The SMT station was able to track the source down to an elevation of a few degrees.Large self-calibration gain factors are therefore needed towards the end of the experiment.Besides these known data issues, gain corrections factors are well constrained within a determined a priori error budget ranging between 10 -20 % for the individual stations 45 .

Imaging
In this section, we describe how we obtained our image model from the observational data.In a first step, we have established a blind consensus between different imaging methods.Then, we have fine-tuned the parameters of one method, eht-imaging 46,47 , for the rPICARD and EHT-HOPS data to obtain final images for the analysis of the Cen A jet structure.
The highest resolution images of this southern source prior to this work were obtained within the TANAMI program 31 at 8 and 22 GHz with a maximum resolution of 400 µas, showing an extremely collimated structure with multiple distinct radio knot emission regions 3 .In a previous single-baseline non-imaging study of Cen A, a bright compact core was detected at 215 GHz 48 .

Blind challenge
Similarly to the method used when the shadow of M 87* was resolved by the EHT 25,49 , we have carried out a blind imaging challenge before proceeding to the scientific analysis of the data.In this challenge, a number of individuals have reconstructed an image of the source independently of each other.Early (not fully verified) low band data from the EHT-HOPS pipeline was used, which had slightly larger amplitude gain errors from outdated a priori calibration parameters.Out of twelve total images, six had acceptable reduced χ 2 cp < 2 for the closure phases.These images were obtained with the eht-imaging and SMILI 50,51 regularized maximum likelihood methods and the Difmap 52,53 and CASA 9,54 CLEAN methods 55,56 .The images that did not make the χ 2 cut often showed spurious emission features and strong sidelobe structures.

Final imaging method
With the imaging challenge, we have established that different methods converge towards the same robust source structure (Supplementary Fig. 5), independent of shared human bias.Further imaging analysis of the rPICARD & EHT-HOPS science release data was pursued with the final M 87* eht-imaging script 49 , which is based on application of a regularized maximum likelihood method that includes a maximum-entropy term.Using a secondmoment-based pre-calibration, LMT gains were stabilized with respect to the better constrained SMT amplitudes 57 .As Cen A is sufficiently compact within the EHT beam, the short LMT -SMT baseline measures a Gaussian-like source structure.We have performed an initial self-calibration to a Gaussian with size Θ maj × Θ min at a position angle Θ PA and with a total flux of S G .Any gains that were erroneously introduced in this process can later be reconciled in image-based self-calibration steps.To solve for the image brightness distribution I with a regularized maximum likelihood method (employed by eht-imaging), we are minimizing Here, D represents the collection of data terms, which are derived from the measured visibilities and have approximately normal noise statistics 58 : amplitudes, closure phases, and log closure amplitudes.Corresponding to each data term, we have a goodness-of-fit function χ 2 D = {χ 2 amp , χ 2 cp , χ 2 lca } and relative weighting α D = {α amp , α cp , α lca }, We have performed four incremental imaging runs with subsequent self-calibration, over which we have increased the weight of each data term: α D ∀D.Regularizer terms Λ R are included with weights β R to impose additional assumptions on the image.We have imposed two regularization parameters: One for a maximum entropy method (MEM) 59 term with weight β MEM and another one for the amount of compact flux Z 0 in the image with weight β z .The MEM term minimizes the entropy of I with respect to a prior image Φ, which results in a similarity between the two images for each pixel i.Here, we used Λ MEM = −S −1 For the MEM prior image Φ, we have chosen a Gaussian model oriented along the direction of the large-scale jet, which we also used as initialization for our imaging.It is expected that Z 0 < S 0 , as a significant portion of the flux measured by ALMA may come from different emission mechanisms and larger scales outside of the EHT field of view.In fact, the ∼ 150 m JCMT-SMA baseline sees a flux density of about 5 Jy and at 2 km, ALMA-APEX recovers only ∼ 4 Jy.For M 87*, the EHT measured Z 0 ≈ S 0 /2 49 .
The numerical values of the final imaging parameters are given in the supplementary material (Table 1).Optimal parameters were chosen based on an empirical minimization of χ 2 D , median station gains A (sc) from self-calibration, and patches of spurious flux in the image.Additionally, we took the similarity of image reconstructions from the rPICARD and EHT-HOPS data for a given set of parameters into account to avoid over-fitting to data peculiarities that result from assumptions made during the data calibration.A variety of images that can be reconstructed with various combinations of the free imaging parameters can be shared upon reasonable requests.We have chosen for an eht-imaging reconstruction of the rPICARD data for our final image, as this imaging method and data set have been studied most extensively.
Our images are shown in units of brightness temperature T [K], which is related to a flux density S in Jy through the observing wavelength λ, Boltzmann constant k B , and angular resolution element Ω as T = λ 2 (2k B Ω) −1 S. Fundamental data properties and fits of the final image model to the data are shown in Supplementary Fig. 7.In Supplementary Fig. 8, we show the measured amplitudes projected along and perpendicular to the jet PA.Along the jet axis, amplitudes fall off quickly at long projected baseline lengths, indicating the absence of substructures along the jet.Perpendicular to the jet, 'bouncing' amplitudes out to large projected baseline lengths occur, due to the strong intensity gradients across the transverse jet profile.

Synthetic data imaging tests
We have used the SYMBA software 60 to perform imaging studies based on simulated observations.Given an input source model M, SYMBA follows the entire EHT signal path to predict which source structure I would eventually be reconstructed.Thereby, we can assess how close our image reconstruction comes to the ground truth structure M of a fabricated observed source.SYMBA simulates the parameters of the Earth's atmosphere with the ATM module 61 to add sky noise, signal attenuation, and phase turbulence.Next, gain, leakage, pointing, and focus errors plus thermal noise are introduced for each telescope in the array based on known telescope properties 62 .Afterwards, the simulated corrupted data are calibrated by rPICARD in the same way as observational data.The u, v coverage and SEFD sensitivities are taken from the 2017 Cen A EHT observation track.The simulated calibrated data are then imaged with the same final eht-imaging script used to image the observational data in this work.
To assess the robustness of secondary features in our image reconstruction, we have performed three synthetic data tests (Supplementary Fig. 4).Firstly, a control study to demonstrate that the output reconstruction from SYMBA correctly matches the input model M final .Then, we have removed the counterjet and emission features at large distance to the apex z from M final to verify that these do not spuriously appear in our simulated observation.
Furthermore, we have explored the upper limit on the brightness T sp of potential emission from the jet spine by adding a weak emission component in the central jet region to M final .The goal was to find the smallest T sp , which would still be registered as an emission region in the reconstructed image I (M final + T sp ).

Jet structure analysis
This section describes how we extract fundamental jet parameters from our image based on geometric arguments.

The position of the jet apex
We can empirically determine the approximate position of the jet apex, where the jet and counterjet are being launched, from the high resolution image model shown in Figure 2. A zoomed-in version of this plot is shown in Supplementary Fig. 9, which is overlaid with visual aids for the determination of the apex location.The first consideration is that the apex should be located in the region where the streamlines of the approaching jet converge.In our image, this convergence region lies upstream of the optically thick radio cores (assumed to correspond to the brightest regions of the jet) for both the NW and SE jet arms.Here, we are limited by the resolution of our instrument, but a tentative merge of the two arms can be seen.The upper arm (region I in the figure) exhibits a strong bend, while the lower arm (IV) remains mostly straight.We note that a similar structure, where one jet arm appears to be straighter than the other one is also present in the M 87 jet 26 .The second consideration is the symmetry between the approaching jet and the counterjet.We note that there is no clear correspondence between individual features in the jet and counterjet.The counterjet appears straight with two components in the upper region (II) and one component in the lower region (III).As the apex must be upstream of the counterjet, the closest component of the receding jet to the approaching jet constrains how far upstream of the approaching jet the apex position can be.In fact, the position we assume for the apex based on the first consideration, where the streamlines of the approaching jet converge, lies halfway between the radio core in region I and the closest counterjet component in region II.It should be noted that a simple extrapolation of only the edge-brightened approaching jet would place the apex well inside the faint counterjet region.
Based on the robustness of our image reconstructions with different data sets, software packages, and imaging parameters, we assume a positional uncertainty of 5 µas for the robust features of the image model, which is in agreement with the width of the jet ridgelines.Taking all constraints on the apex location into account, we estimate an uncertainty of 10 µas on the position.
For the determination of z core , the pixel and jet apex position uncertainties are added in quadrature.Based on possible jet apex positions within the estimated uncertainty, we fit the W ∝ z k jet profile multiple times and derive a systematic error of ±0.06 on k.When we used image model convolved with the nominal resolving beam, we obtain k = 0.35 with a statistical error of ±0.2.

Brightness asymmetries
The jet-counterjet asymmetry is most likely caused by relativistic boosting.We can calculate the R j/cj brightness ratio by taking the average image flux density within 50 × 100 µas rectangular regions on opposite sides of the apex.This ratio has to be interpreted with care, since the two regions may be at different distances to the jet apex.Moreover, counterjet radiation may be absorbed by the accretion flow and intrinsic jet-counterjet differences may arise from asymmetries in the jet launching process and the ambient medium 63,64 .
If we assume the intrinsic emissivity to be the same in the jet sheath and spine, beaming effects can be invoked to explain observed differences in brightness across the jet.We note that the intrinsic emissivity of the jet sheath is likely larger than that of the spine, as mentioned in the main body of this work.The simplifying assumption of identical intrinsic emissivities can nonetheless be used to derive straightforward estimates for jet velocity components and the inclination angle θ, since Doppler boosting is expected to have a considerable contribution to the observed source structure.If the inclination angle θ is not too small, a significant portion of the spine emission may be beamed away from the line of sight.If the sheath and spine velocities are cβ sh and cβ sp , respectively, the ratio of I sh sheath and I sp spine intensities in a continuous jet follows as with α (I ∝ ν α ) as the spectral index of the optically thin jet components.Assuming a typical spectral index of α = −0.7 and identical intrinsic emissivities, we can constrain the sheath and spine velocities with Equation 3and (1−β 2 ) −0.5 (1−β cos (θ)) has a minimum of 1 − cos 2 (θ) at β = cos (θ).It follows that the sheath-spine asymmetry can only be explained via beaming for β sp > cos (θ), independent of the assumed value for the spectral index.For Cen A, the jet spine emission may be beamed away from the line of sight, when its velocity exceeds 0.7 c -0.9 c, while the sheath moves with a slower velocity.In fact, the emitting plasma of the large-scale jet was observed to move with 0.24 c -0.37 c 3 .
For a full three-dimensional picture of a jet, where we assume the sheath to be symmetric in the φ direction around the spine in a cylindrical coordinate system, different spine and sheath emissivities, due to beaming or intrinsic effects, cannot, on their own, explain edge-brightening.The reason is that the sheath emission will contribute to any sightline towards the jet.A more detailed description, where also the optical depth is taken into account, is given in the next paragraph.In the remainder of this section, we go through the different scenarios that could cause the observed edge-brightening.First, we discuss a common interpretation related to pathlength differences.As this only works in optically thin regions, we put the presence of helical magnetic fields forward as the most likely, intrinsic explanation for edge-brightening in LLAGN.We then discuss more exotic scenarios, of a rotating or asymmetric jet, which might be tested through future observations.
In the optically thin jet regions, the integrated column density along sightlines through the jet at different distances from its axis (center vs. edges) can be used to explain edge-brightening.These are sightlines that, across the transverse extent of the jet, enter the jet at different locations.The sightlines first pass through the near side of the jet and exit again at the other side of the jet, the far side.If we assume the absence of intrinsic spine emissivity (due to weak mass-loading or beaming of radiation into a narrow cone away from the line of sight), the observed radiation will be produced by a sheath of thickness ∆R.For a line sight that goes exactly through the center of the jet, we pass twice through the sheath, which would amount to a pathlength of 2∆R/ sin θ, when the pathlength is short enough to locally approximate the jet as a cylinder.For a local jet radius R j , the column density along a sightline through the edge of the jet will be larger by a factor of ∼ R j /∆R 26 .Here, we have neglected changes in emissivity as sightlines pass through material at different distances to the jet apex.This simple model is capable of explaining edge-brightening in optically thin jet regions, where radiation along longer pathlengths accumulates.For Cen A, this would imply a thin radiating sheath with ∆R < 0.04R j .
However, the edge-brightening in Cen A extends to the presumably optically thick radio core, suggesting that different physics are at play in this jet.The likely presence of a helical magnetic field 17,65,66,67 combined with a rotating sheath and the inclination angle θ, can lead to favorable/unfavorable pitch angles that maximize/minimize the synchrotron emissivity along the edges/center of the jet.For a power-law distribution of electrons, where in the rest-frame of the jet, the electron density n follows their energy E as dn ∝ E −p dE, the synchrotron emission coefficient in the rest frame scales as j ν ∝ |B sin χ| (p+1)/2 ν −(p−1)/2 68,69 .Here, B is the magnetic field strength, χ is the angle between the magnetic field and line of sight, and ν is the radiation frequency.The corresponding absorption coefficients scale as α ν ∝ |B sin χ| (p+2)/2 ν −(p+4)/2 68,69 .It can be seen that no asymmetries in χ would arise across the transverse jet profile for a purely poloidal (B z ) magnetic field.The edge-brightening is maximized for perpendicular angles χ between the line of sight and magnetic field at the jet edges, while the magnetic field is oriented parallel to the line of sight in the center of the jet.In future work, we will study the polarimetric properties of the jet with the EHT to test this hypothesis as explanation for the edge-brightening.To get a handle on χ, it will be necessary to narrow down the inclination angle θ and jet velocity with monitoring observations to take relativistic aberration into account.
For optically thick jet regions upstream of the radio core, the relativistic boosting is sensitive to the shape of the emitting region and less sensitive to the Doppler factor 18 .In the presence of a fast helical jet flow and θ > 0, part of the jet will rotate towards the observer and the other part will rotate in the opposite direction on the sky.Beyond the initial jet launching region, the jet is strongly collimated and the viewing angle to the jet edges will be very close to θ.For a flow with toroidal and poloidal components, we denote the angle of the helical velocity component β h with respect to the poloidal direction along the line of sight with φ h .For two identically shaped, optically thick radio core components of intensity I s at the SE jet edge and I n at the NW edge, we thus have 26 For a counter-clockwise jet rotation and R s/n ∼ 1.6, we get the weak constraint of When the bulk velocities of the NW and SE jet sheaths are known, φ h and subsequently β h can be determined 26 .We note that the the small linear scales resolved by the EHT in Cen A uniquely allow us to track relativistic dynamics across days in this source with future observations.In an alternative scenario, this tentative NW-SE brightness asymmetry seen in Cen A could be explained with two distinct jet components having different velocities or different inclinations angles with respect to the line of sight.
In this work, we have interpreted the edge-brightening in terms of a naturally emerging spine-sheath jet structure in LLAGN, based on results from GRMHD simulations that are applicable to those type of sources.However, the same phenomenon is also observed in more powerful AGN; for example Cygnus A 70,71 , where an accretion flow operating at ∼ 1 % of the Eddington limit is unlikely to be radiatively inefficient 72 .

Collimation profile
Following the NW and SE jet ridgelines, we bin distance values to the jet apex into intervals of 10 µas in size.Within each bin, we select the brightest pixel to obtain the central location along the ridge.We impose a statistical uncertainty of 5 µas on distances z in accordance with the width of the jet ridgelines in our image model.The width W of the jet is taken as the distance between the two jet arms.The profile of our image is shown in Figure 3 together with the corresponding average opening angle computed from the jet width as a function of distance to the apex.
Resolution limitations prevent us from tracing down the exact value of the initial jet opening angle ψ jet near the apex, where the analysis of binned distance values becomes uncertain.Nonetheless, we can derive an upper limit on ψ jet with a simple geometric argument: The jet has a clearly defined collimation region beyond some distance from the apex, at z > z col .To estimate z col , we have used the SE jet arm, as it is brighter, straighter, and has a more clearly identifiable compact brightness core.If we now assume that the jet stream converges monotonically towards the apex for z < z col and that the apex itself does not correspond to an extended region, we have If the inclination angle θ is known, the intrinsic opening angle ψ int can be computed as 73 The jet remains collimated out to kpc scales and contains multiple particle acceleration sites in a knotted structure 3,74,75,76 .The source is a well-suited laboratory for models of AGN feedback 77,78 and the creation of ultra-high-energy cosmic rays 79,80 .

Confinement by the ambient medium
Analytic theory for axisymmetric, relativistic, poynting-dominated outflows can be used to derive exact asymptotic solutions for the influence an ambient medium on the collimation of a jet.One can show that in the presence of external pressure gradient P ext (z) = P 0 z −κ , the jet expansion profile W as a function of distance along the jet axis z follows 81,82 in a simplified form, with C 1 a numerical constant.At large z and for a shallow external pressure gradient with κ < 2, we obtain 81,82 for S(z) = C 4 z 1−κ/2 −C 5 and C 2 , C 3 , C 4 , and C 5 numerical constants.Equation 10shows that the ambient pressure will confine the jet into a W ∝ z k profile with k = κ/4.Additionally, oscillations along the jet boundary can occur in a non-equilibrium state for C 2 = (2 − κ)/π, C 3 = 0 82 .

The location of the black hole
Given a measurement of the core shift z core with respect to the black hole, we can gauge the observing frequency ν, which corresponds to a small self-absorbed nozzle region at the footprint of the jet 83 .This region corresponds to a peak or break from a jet-dominated flat radio spectrum as it is the smallest region where particle acceleration can occur.The minimum scale where a jet can be launched by a black hole is given by the ISCO.The size of the emission region of this nozzle would be given by the photon capture radius.Thus, we can estimate ν as In this expression, M is the mass of the black hole and D the distance from the black hole to the observer.With the derived scaling relation of ν , we can relate the break frequencies of two sources if their accretion rates or jet properties are known.Here, Ṁ is the black hole accretion rate and F r is the observed flat-spectrum radio flux density.In particular, if we assume for two sources to share the same basic intrinsic jet properties and orientation with respect to Earth, we have While these expressions are strictly speaking only true for a filled conical jet, they appear to describe the emission from the jet sheath and its basic scaling properties reasonably well 84,85,86 and allows one to make a first order estimate of the characteristic radio frequency of near-horizon emission.
We have used the above equations to estimate the accretion rate of Cen A to the one of M 87 based on the assumption of a similar coupling between SMBH inflows and jet power.External Faraday Rotation effects and a generally variable rotation measure further complicates the assumed relation of accretion rates 87 , which should thus be taken as only an order-of-magnitude estimate.It is however worth pointing out that the black hole growth rate measured over cosmic timescales from X-ray cavity fluxes from the jet radio lobes is ∼ 10 −3 M yr −1 for both Cen A and M 87 88 .

Alternative interpretations for the brightest jet features
In this work, we have interpreted the brightest jet features as radio cores, which mark the transition region between upstream synchrotron self-absorbed jet regions and downstream optically thin areas.In our image, we are able to resolve the self-absorbed region between the putative radio core and jet apex, which coincides with the location of the SMBH and its accretion disk.With current telescopes, the radio core and upstream region remains unresolved for most AGN (see Table 2 in 1 for example).
The radio core interpretation of the brightest jet features seems most plausible given our data.Based on simple analytical jet theory, a bright radio core is expected to be present in VLBI images.Radio cores are typically seen in sources similar to Cen A and the core-shift typically follows the standard ν −1 relation in most sources 89 .In fact, special circumstances have to be invoked to explain the absence of radio cores in VLBI images.For example, obscuration by an optically thick region in the foreground.We do think that this is a likely scenario for our observation given the small scales probed in our image, the high observing frequency, and the proximity of the source.Moreover, the core shift we have computed in Cen A agrees with the core SED of the source and fundamental plane equations that relate the jet power of Cen A to the M 87 jet power.Our image is dominated by the brightest, compact jet features, which would be weakly polarized and have flat spectrum as radio cores 1 .It should be noted that wide-bandwidth ALMA interferometer data, that were taken simultaneously with the VLBI observations, show a flat spectrum between 212 GHz and 230 GHz and place a 3 σ upper limit of 0.15 % on the linear polarization fraction 87 .The ALMA measurements are however at a larger arcsecond resolution and we resolve out 64 % of the flux measured by ALMA with the EHT.We therefore need future polarimetric and spectral VLBI results for confirmation.
However, with the current observations, we cannot conclusively rule out the possibility that there is insufficient particle acceleration in the jet, such that no radio core is formed at λ 1.3 mm, while a core is present at longer wavelengths 3,76 .In this scenario, the bright jet regions would most likely correspond to a shock within the jet flow.The strongest counterargument here is that the radio spectrum of the core turns over at THz frequencies.This emission is most likely produced by the jet, ergo, particle acceleration should occur up to the energies that produce THz synchrotron emission.

Data availability
The ALMA raw visibility data can be retrieved from the ALMA data portal under the project code 2016.1.01198.V.The calibrated Stokes I VLBI visibility data of Centaurus A can be obtained from a DOI listed under https: //eventhorizontelescope.org/for-astronomers/data with the code 2021-D03-01.Image FITS files and scripts to reproduce the plots are available from the corresponding author upon reasonable request.
This work is based on the 'ER6' data production scripts, for which the 30e6ca14fb50275013c668285a3b476f9bc85436 91da63236db34f3a31b5309b18ac159128f28a35 image was used.
The docker image used here is tagged as dec65699ccc0acdc6e6ba8f218d6724537fc613a and can be found on https://hub.docker.com/r/mjanssen2308/symba.

Figure 1 :
Figure 1: The jet structure of Cen A in comparison with M 87.The left panel (a) shows the large scale jet of Cen A from an 8 GHz (λ 3.7 cm) TANAMI 31 observation in November 2011 3 on a logarithmic color scale.The middle panel (b) depicts our final EHT image from April 2017, blurred to the nominal resolution for a uniform weighting of the visibilities (the beams are shown in the bottom right corners).The reconstruction is based on the rPICARD data and is shown on a square-root scale, where values below a brightness temperature of 3 × 10 8 K are clipped, due to a lower dynamic range compared to the longer-wavelength observations.An unclipped and unconvolved version of this image is shown in Figure 2. The right panel (c) shows the M 87 jet at 43 GHz (λ 7 mm) from a Very Long Baseline Array observation in June 2013 9,26 on a logarithmic scale.North is up and east is to the left.The physical, linear scales of the full field of views shown in the three images are 2 parsec (pc) for TANAMI (left), 0.007 pc for the EHT (middle), and 0.6 pc for the VLBA (right).

Figure 2 :
Figure 2: Image analysis of the final model.The model image shown in the left panel (a) corresponds to the image from Figure 1 with a pixel size of 2 µas.Here, the tentative position of the jet apex is indicated with a circle.The size of the circle indicates the uncertainty in the apex location.The right panel (b) shows the central brightness temperatures along the jet ridgelines from the model in four quadrants as a function of distance to the jet apex.The quadrants I, II, III, and IV correspond to the similarly marked regions of the jet in the left panel.Negative values for the distance to the jet apex are assigned for the counterjet region.Brightness temperatures of the fainter NW (orange line) and brighter SE (blue line) arms are shown in the upper and lower panels, respectively.

Figure 3 :
Figure 3: Cen A collimation profile derived from the final model.The top panel shows the jet width W as a function of apex distance z, overplotted with a W (z) = Az k least-squares fit, with k = 0.33 ± 0.05| stat ± 0.06| sys (Methods 1.4.3).The systematic error has been derived based on the uncertainty in the apex position.The error bars are derived from a 5 µas one sigma statistical uncertainty on each pixel in the model.The bottom panel shows the corresponding jet opening angle between the apex and the two jet arms, which becomes unreliable at z z col = 32 µas due to instrumental resolution limitations (see text).The transition region is indicated with a vertical blue dashed line.A base-10 logarithmic scale is used for the horizontal axis.

Supplementary Figure 4 :
SYMBA synthetic data imaging tests.The left panel shows the input models M, which are based on our observational image reconstruction.The corresponding reconstructions I from the SYMBA pipeline are displayed in the right panel.(a-b) a control study, where the final model M final underlying the image reconstruction I (obs) final from the observational rPICARD data is passed through SYMBA.(c-d and e-f) same as (a-b) but the counterjet and extended jet emission features have been removed.The SYMBA reconstructions are convolved with a restoring beam, which matches the nominal resolution of the observation (shown in the bottom right corner).