© M. H. Xu et al. 2022

1 Introduction

In fundamental astronomy, two cardinal improvements in the precision of astrometric measurements of positions of celestial objects have been achieved by very long baseline interferometry (VLBI) at radio wavelengths and by the European Space Agency mission Gaia¹ (Gaia Collaboration 2016) at optical wavelengths. A detailed comparison between the third realization of the International Celestial Reference Frame (ICRF) (ICRF3; Charlot et al. 2020) and Gaia Early Data Release 3 (EDR3; Gaia Collaboration 2021) shows that the median difference of radio and optical source positions for more than 2000 common sources is about 0.5 milliarcsecond (mas), and for the sources with optical G magnitudes <18.0 mag the median difference is on the order of 0.3 mas (Xu et al. 2021b). On the other hand, for the sources that have radio to optical position differences larger than 3σ, which is ~24% of the common sources, it has become more and more convincing that the position differences are parallel to the directions of the radio jets (Kovalev et al. 2017; Petrov & Kovalev 2017; Plavin et al. 2019a; Xu et al. 2021b). Our recent work has shown that these sources are more likely to have extended structure at radio wavelengths (Xu et al. 2021b). Gaia and VLBI already have the potential to detect astrophysical properties of the radio and optical emission at mas scales for these objects (Plavin et al. 2019a; Petrov et al. 2019), which are mostly active galactic nuclei (AGNs).

In the last several decades, celestial reference frame (CRF) sources have predominantly been observed in the S and X bands by geodetic VLBI to derive their positions for building the ICRFs (Ma et al. 1998; Fey et al. 2015; Charlot et al. 2020). The geodetic VLBI observations have been coordinated by the International VLBI Service for Geodesy and Astrometry (IVS²; Schuh & Behrend 2012; Nothnagel et al. 2017) since 1999. In order to achieve 1 mm station position accuracy and 0.1 mm yr^-1 velocity stability on global scales, the IVS is developing a next-generation, broadband geodetic VLBI system, known as the VLBI Global Observing System (VGOS; Niell et al. 2007; Petrachenko et al. 2009). While the primary goals are to greatly increase the precision of Earth orientation parameters and the International Terrestrial Reference Frame (ITRF; Altamimi et al. 2016) using VGOS, the CRF, as an integral part of geodetic VLBI, is also expected to be improved in the VGOS era. At the same time, because the source positions in the ICRF3 have a median formal error of 0.1 mas in right ascension and 0.2 mas in declination, they need to be improved by one order of magnitude to an accuracy of about 1/30 of the interferometer beamwidth in order to achieve the ambitious goals of VGOS. A considerable effort is thus required to mitigate the systematic errors in source positions due, for instance, to source structure. Based on 21 actual VGOS sessions, it was demonstrated that the measurement noise in VGOS group delay observables is at the level of ~2ps (Xu et al. a2021). A network of nine VGOS stations carried out biweekly, 24-hour observing sessions in 2019 and 2020 and will start to observe on a weekly basis in the near future; this network is expanding globally, as planned (Behrend et al. 2020), with another eleven stations built and nine stations in the planning stage as of January 2021. Now VGOS observations are used to build the latest realization of the ITRF³.

Currently, the VGOS system observes simultaneously in four 512 MHz wide bands centered at 3.3, 5.5, 6.6, and 10.5 GHz, which are labeled as bands A, B, C, and D, respectively, with 32 recording channels (see Niell et al. 2018 and Table C.1 for the channel frequencies used in the current VGOS observations), while the legacy VLBI system observes in the S (~2.2GHz) and X (~8.6GHz) bands. There is a substantial difference in the observing frequencies and thus in the received radio signals from AGNs between the legacy system and the VGOS system. Typically, the celestial sources in geodetic VLBI observations have extended structure with a core-jet morphology at mas scales, where the jet is the highly collimated outflow of magnetized plasma and the core is a compact, bright, and flat spectrum feature in the upstream part of the flow. These celestial sources have different structures in different bands, and the structure changes over time at sub-mas angular scales; therefore, we cannot assume that the locations of the radio emission of a source in different frequency bands are at the same position or have stable relative positions between the four bands at angular scales of ~0.1 mas to 1.0 mas. The fact that source positions change with frequency is one of the reasons why the ICRF3 has three separate catalogs in the three frequency bands (Charlot et al. 2020). The fundamental question regards the position of a source determined by VGOS observations with respect to the locations of the radio emission in the four bands. We aim to address this question here.

On the other hand, because the impact of source structure on VGOS group delay observables is about one order of magnitude larger than the random noise level, as shown in Xu et al. (2021a), the IVS is making an effort to derive source-structure corrections for geodetic solutions by collecting information of antenna system temperatures and gain curves and producing images from VGOS observations. For now, before VGOS source-structure corrections are generated, it is necessary that the images in the four bands be aligned with respect to one another. The impact of the potential errors in that alignment needs to be studied. This is equivalent to the previous question that we aim to discuss.

The purposes of this study are thus twofold: (1) to investigate the potential impacts of source structure on the source positions as determined by VGOS observations if these effects are not modeled and (2) to demonstrate the importance of image alignment over frequency when one wants to correct for these effects. The paper is structured as follows. We introduce in Sect. 2 the changes of source positions with respect to frequency based on the ICRF3 and images obtained based on actual VGOS data. In Sect. 3 we first describe how the variations in the channel phases affect the VGOS observables and then derive the formula of the position estimate from VGOS observations as a function of the locations of the radio emission in the four bands. We present a discussion in Sect. 4 and our conclusion in Sect. 5.

2 Changes in source positions over frequency

Based on globally absolute astrometric observations by VLBI in the S/X, K, and X/Ka bands, the ICRF3 was established independently at these three frequencies (Charlot et al. 2020). The two catalogs in the K and X/Ka bands were aligned to the S/X catalog by applying no-net-rotation constraints, leading to insignificant spin rotations between them. However, the glide deformation between the S/X and K catalogs is 10-30 microarc-seconds (Mas), and that between the S/X and X/Ka catalogs is as large as ~300 Mas. These deformations were believed to be due mainly to the different observing networks of the three bands – for example, the X/Ka observations were made from a network of four stations – instead of the intrinsic source position differences. After applying transformations, including both the first-degree and the second-degree deformation parameters fitted from the position differences between these three catalogs (16 parameters in total), the median angular separation of the K-band positions relative to the S/X-band positions is ~0.2 mas, the same level as that for the X/Ka-band positions; about 6% of the common sources in the K and S/X catalogs have position differences significant at the 3σ confidence level, and 11% for the common sources in the X/Ka and S/X catalogs (Charlot et al. 2020).

Figure 1 shows the radio positions in the three catalogs of the ICRF3 and the Gaia EDR3 position for source 1157-215 overlaid on its radio images at 8.7 GHz from the Astrogeo image database⁴ and at 15.3 GHz from the Monitoring Of Jets in Active galactic nuclei with VLBA Experiments (MOJAVE) program (Lister et al. 2018). As shown in the images, the ratio of the flux densities between the southeastern component and the northwestern component largely increases, from 15.3 to 8.7 GHz; this shifts the radio positions toward the southeast direction when frequency decreases. The northwestern component has an inverted spectrum and is more likely to be the core than the brightest component in these images. For source 1157-215, there are angular separations of the positions in the K, X/Ka, and Gaia catalogs relative to that in the S/X catalog that are all significant at the 3σ confidence level (see Tables 14-16 in Charlot et al. 2020). Figure 1 demonstrates that source structure and source positions change significantly with frequency and that for this source the radio positions move in the direction of the optical position when the radio frequency increases. As the jet components are optically thin (i.e., dimmer at higher frequencies), the identity of the peak component can even change between the core and the jet components as the frequency band changes, which may happen to some radio sources in VGOS observations (e.g., 3C418, 0642+449, 1030+415, and 2000+472). This leads to large position offsets between different frequencies, as the source position referred to in VLBI observations is typically dominated by the position of the peak.

In general, the radio sources targeted by VGOS are chosen to be as compact as possible. These sources, however, are not without structure. They are often marginally resolved and exhibit “invisible” structure. This structure has an extent on the order of that of the interferometer beam or smaller. While such a structure can be almost invisible in the images, its group delay and phase observables correspond closely to that of a point source located at the position of the brightness centroid, which can change with frequency. Moreover, the interferometer beam of geodetic VLBI that has antennas on Earth typically has a diameter on the order of 0.5–1.5 mas in the VGOS frequency bands, and the source position accuracies required by VGOS are about 30 times smaller than this (i.e., 0.02 mas). Therefore, the invisible structure is as important as the visible structure in VGOS observations. To investigate the general cases of compact sources, we used source models constructed from actual VGOS observations by: (1) deriving closure images using the method in Xu et al. (2021d); (2) calibrating VGOS observations based on the closure images; and (3) performing model fitting using difmap to model the structure with a small number of two-dimensional Gaussian components.

Figure 2 plots such reconstructed models for source 0016+731 in the four frequency bands based on VGOS observations VO0051, and the Gaussian parameters are listed in Table 1. At the highest three frequencies, the structure of source 0016+731 can be consistently modeled with three components: (1) a compact core (C), which is the northwestern component, (2) the first jet component located at the center (J1), and (3) a weaker jet component (J2), which is the southeastern component. However, only two components are needed to represent the structure in the lowest frequency band due to its significantly lower angular resolution. By comparing the Gaussian components in the four bands, it is very likely that components C and J1 are combined into a single component at 3.3 GHz. Even though the structure at 3.3 GHz is different from that at 5.5 GHz – due to the different spectral indices of the core and the jets changing the relative flux densities and the core shift changing the relative positions – we can assume that its structure at 3.3 and 5.5 GHz is the same when illustrating the potential impact. In this case, since at 5.5 GHz the angular separation between C and J1 is 0.42 mas with a flux density ratio of 0.89, the position of the modeled core at 3.3 GHz would be shifted toward the east by about 0.20 mas – that is, 0.42 mas x (0.89/(1.0+0.89)) – due to the contribution of J1, which is only visible in the highest three frequency bands. As shown by this example, jet components may significantly affect the position of the modeled core at the lower frequencies due to lower resolution and cause misalignments between the bands despite source structure being taken into account. Furthermore, if the effects of source structure are not modeled, the source position shifts due to the contribution of the jet components will happen in all the frequency bands, therefore further increasing the image misalignment between the bands.

Changes in source positions with frequency are mainly due to two factors: extended source structure and the frequency dependence of the core position, the so-called core shift (Blandford & Königl 1979). This is supported by the ICRF3 multi-frequency catalogs and the optical positions from Gaia. It has been shown that the differences between the radio positions and the optical positions are typically parallel to the jet directions (see, e.g., Lambert et al. 2021; Liu et al. 2021). These two effects in general are larger in the lower frequency bands (S/X and the VGOS frequency bands) than in the higher frequency bands (K and X/Ka bands). It is expected, however, that the impact of the effects in the S band on the previous ICRF realizations is significantly reduced. According to the linear combination of the delay observables in the S and X bands, which is done in geodetic data analysis after the fringe fitting process to remove the ionospheric effects, the source position determined from VLBI observations in these two bands, denoted by k_S/X, is given approximately by $$k_{S/X}=1.07k_X-0.07k_S,$$(1)

where k_X and k_S are the group-delay source positions in the two bands (see, e.g., Porcas 2009). When a source is point-like in the S and X bands, these two positions coincide with the jet base (Porcas 2009); in general cases (i.e., extended sources), they are, however, pointing in different directions. As shown in Eq. (1), the contribution of the S-band positions, including potential variations, to the S/X band positions are reduced by a factor of 14. We should note that the effects of source structure in the S band are scaled down by the same factor. In the legacy VLBI system, the source positions in the higher band thus dominate over that in the lower band. In the VGOS system, however, observations are made in four frequency bands instead of two bands, and the ionospheric effects are fitted in the fringe fitting process, through which group delay observables are determined. As a consequence, the contribution of the four individual bands to the VGOS source positions needs to be studied thoroughly, which is the purpose of the next section.

Fig. 1 Radio positions from the ICRF3 in the S/X, K, and X/Ka bands and the Gaia position for source 1157-215, on its radio images at 8.7 GHz from Astrogeo (left) and at 15.3 GHz from MOJAVE (right). Overlay contours are shown at levels of peak percentage, specified at the bottom of the plot. Since these radio images do not have information regarding the absolute positions, the peaks are formally assumed to be located in the S/X-band positions. Error bars are the uncertainties from the four catalogs. The K- and X/Ka-band positions are the ones after the full deformation transformations (Charlot et al. 2020). The angular separations of the K- and X/Ka-band positions and the Gaia position relative to the S/X-band position for this source are 2.29 ± 0.31 mas, 3.54 ± 0.32mas, and 3.97 ± 0.13 mas, respectively.

Fig. 2 Images of source 0016+731 from VGOS session VO0051 (February 20, 2020) at 3.3 GHz (upper left), 5.5 GHz (upper right), 6.6 GHz (bottom left), and 10.5 GHz (bottom right). The ellipses indicate the Gaussian components, the parameters of which are reported in Table 1. Overlay positive contours are shown at nine levels of peak percentage, starting from 0.2% and increasing by a factor of 2, while the negative contour, indicated as dash lines, has a peak percentage of – 0.2%. These images are constructed from the VGOS observations, which are calibrated based on the images derived from closure phases and closure amplitudes (Xu et al. 2021d). Therefore, the pixel flux densities have arbitrary units. The beam size is shown as a filled ellipse in the bottom-left corner of each plot.

3 Simulation for VGOS source positions

3.1 Changes in broadband observables due to phase variations

We recall the model of using 32 channel phases to determine the broadband observables as done in the VGOS post processing using fourfit⁵. Following Cappallo (2016), the visibility phase ϕ_νi at the channel with frequency v_i can be expressed as $${\varphi }_{v_i}={\tau }^g\left(v_i-v_0\right)+{\varphi }_0-\frac{k\delta \text{TEC}}{v_i},$$(2)

where τ^g, ϕ₀, and δTEC are the broadband group delay, the broadband phase at the reference frequency ν0, which is 6.0 GHz in the current data processing, and the differential total electron content (TEC) along the ray paths of the radio emission from a source to the two stations of a baseline, respectively; these three quantities are the broadband observables simultaneously fitted in fourfit. The constant k is equal to 1.3445 when phases are in units of turns of a cycle, delays in units of nanoseconds (ns), frequencies in units of GHz, and δTEC in units of TECU⁶. We can derive the phase delay via τ^p = (ϕ₀ + N)/v₀, where N is the integer number of phase turns and τ^p is in units of ns.

We denote the observational equation as follows, $$l=Ax+\sigma ,$$(3)

where l is the vector of the visibility phases observed in each of the 32 channels, A is the design matrix, x is the vector of the three unknowns (i.e., τ^g, ϕ₀, and δTEC), and σ is the noise vector. In the present study, the vector l consists of the simulated channel phases based on delay offsets or source position offsets. The design matrix A, with a dimension of 32 × 3, consists of the partial derivatives of the visibility phase with respect to the three unknowns, which can be derived based on Eq. (2). By assuming that the visibility amplitudes over 32 channels are flat (i.e., the phases ϕ_νi are given equal weights), the normal matrix can be derived as A^TA. Via least squares fitting (LSF), the changes in the broadband observables due to changes in the channel phases can be determined.

The results shown in Table 2 were obtained from LSF for seven possible combinations of bands with a 1 ps delay offset applied to the channel phases in these bands based on the channel frequencies listed in Table C.1. The four scenarios labeled as “Case 1” correspond to the cases where only one of the four bands has a delay offset that causes variations in the eight channel phases. The normal matrix in the LSF process of estimating the broadband observables is independent of channel phases, and the estimates of the broadband observables are linearly dependent on them. Therefore, there are two features in this broadband fitting process. First, the coefficients of propagating the delay offsets in individual bands to the broadband observables are independent of the magnitudes of those delay offsets: they are invariable for a given set of channel frequencies. Second, the results of all other possible combinations of bands with delay offsets can also be derived from the four basic scenarios in Case 1 through a linear combination. For instance, the first scenario in Case 2 is the summation of the results of the first two scenarios in Case 1, and Case 3 is the summation of that of the four scenarios in Case 1.

For general scenarios, where the delay offsets in the four bands are Δτ_Α, Δτ_B, Δτ_C, and Δτ_D, respectively, the change in the broadband group delay, denoted by Δτ^g can be obtained from the results in Table 2 as $$\text{Δ}{\tau }^{\text{g}}=+0.505\text{Δ}{\tau }_{\text{A}}-1.448\text{Δ}{\tau }_{\text{B}}-1.458\text{Δ}{\tau }_{\text{C}}+3.401\text{Δ}{\tau }_{\text{D}},$$(4)

and the change in the phase delay, denoted by Δτ^P as $$\text{Δ}{\tau }^{\text{p}}=-0.883\text{Δ}{\tau }_{\text{A}}+1.729\text{Δ}{\tau }_{\text{B}}+2.044\text{Δ}{\tau }_{\text{C}}-1.889\text{Δ}{\tau }_{\text{D}},$$(5)

This development is similar to what was done for the investigation of the impact of constant instrumental delays between different bands on VGOS broadband delays by Corey & Himwich (2018). Equation (4) can be used to recover the results of the 15 possible scenarios in Fig. 1 and Table 2 of their VGOS memo⁷.

3.2 Source positions determined by VGOS observations

The delay offset in an individual band can change with time and be baseline- and source-dependent due to some physical causes, for instance a source position offset. We denote the (phase-delay) source position in band A by k_A = k₀ + Δk_A, where k₀ is a reference direction and Δk_A is a position offset with respect to the reference direction, and k_B = k₀ + Δk_B, k_C = k₀ + Δk_C, and k_D = k₀ + Δk_D for the other three bands, respectively. These source positions are defined at the center frequencies of the individual bands because source structure and core shift evolve with frequency even within a band.

The (phase) delay offset in band A due to the source position offset Δk_A is computed as $$\text{Δ}{\tau }_{\text{A}}=-\frac{B\cdot \text{Δ}{k}_{\text{A}}}{c},$$(6)

where B is the baseline vector and c is the speed of light. In this calculation, there is no need for spherical trigonometry when the position offset is small. This is always possible because we would expect the differences in the locations of the radio emission of the CRF sources between the four bands to be at the mas level or even smaller. By applying this common delay-position relation to the four delay terms in Eq. (4), the VGOS group delay position with respect to the reference direction, denoted by Δk^g, is jointly determined by the locations of the radio emission in the four bands as follows: $$\text{Δ}{k}^{\text{g}}=+0.505\text{Δ}{k}_{\text{A}}-1.448\text{Δ}{k}_{\text{B}}-1.458\text{Δ}{k}_{\text{C}}+3.401\text{Δ}{k}_{\text{D}}.$$(7)

Similarly, the VGOS phase-delay position with respect to the reference direction, denoted by Δk^p, is given by $$\text{Δ}{k}^{\text{p}}=-0.883\text{Δ}{k}_{\text{A}}+1.729\text{Δ}{k}_{\text{B}}+2.044\text{Δ}{k}_{\text{C}}+1.889\text{Δ}{k}_{\text{D}}.$$(8)

The summation of the four coefficients on the right-hand side of Eq. (7) is equal to unity, which is also the case for Eq. (8) (the round-off error of the displayed coefficients notwithstanding). Taking it one step further, the VGOS group delay absolute position k^g = k₀ + Δk^g and phase-delay absolute position k^p = k₀ + Δk^p can be obtained from $$\begin{array}{l}{k}^{\text{g}}=k_0+0.505\text{Δ}{k}_{\text{A}}-1.448\text{Δ}{k}_{\text{B}}-1.458\text{Δ}{k}_{\text{C}}+3.401\text{Δ}{k}_{\text{D}}\hfill \\ =+0.505k_{\text{A}}-1.448k_{\text{B}}-1.458k_{\text{C}}+3.401k_{\text{D}}\hfill \end{array}$$(9) $$\begin{array}{l}{k}^{\text{p}}=k_0-0.883\text{Δ}{k}_{\text{A}}+1.729\text{Δ}{k}_{\text{B}}+2.044\text{Δ}{k}_{\text{C}}-1.889\text{Δ}{k}_{\text{D}}\hfill \\ =-0.883k_{\text{A}}+1.729k_{\text{B}}+2.044k_{\text{C}}-1.889k_{\text{D}},\hfill \end{array}$$(10)

respectively. As we can see from these equations, they are independent of the reference direction k₀.

A simulation based on the VGOS session VO0051 (February 20, 2020) was run to demonstrate the results for group delays. Figure 3 shows two cases, where we assume that the position of source 0016+731 in band B is offset by 0.2 mas in declination with respect to that in the other three bands (in the same position and formally selected as k₀) and in band D by 0.1 mas in declination. By referring to k₀, the position-offset-induced phases at the 32 channels of each individual observation were calculated, from which the broadband observables were fitted using Eq. (2). The calculation was done for all the observations of source 0016+731 in the session one by one. The source positions determined by these simulated broadband group delays for the 1308 observations of the source in the session are (0.000, −0.289) mas and (0.000, 0.340) mas for the two cases, respectively. These agree with the predictions from Eq. (7): in the first case the shift is about 1.45 times that for band B and in the opposite direction, while in the second case it is about 3.4 times the shift for band D in the same direction. This can also be seen from the relative magnitude and sign of the blue curves in Fig. 3. In the simulation, the position offset of 0.2 mas causes delays on a 9000 km baseline with a magnitude of ~30 ps, which is smaller than the phase-delay ambiguity spacings at the four bands and thus does not cause an issue of 2π ambiguities in channel phases when doing LSF. This issue could arise in the actual observations if position offsets are larger.

We note that source structure and core shift evolve continuously with frequency, thus leading to additional phase variations within individual bands. These intraband phase variations have much narrower frequency ranges compared to the phase variations from band to band; they cause a second-order impact on broadband VGOS observables and are not accounted for in this study.

Two other kinds of simulations were performed, the results of which are shown in Appendices A and B. It is possible in the future that channel frequencies in VGOS observations will be changed to a wider bandwidth and a broader frequency range than the current settings. The results for two other possible sets of channel frequencies are given in Appendix A. The main conclusion from this simulation is that, with a broader frequency range, the coefficients in Eq. (7) will be reduced significantly. This also means that if the VGOS observing frequencies are changed – for instance if the upper limit of frequency is extended to 14 GHz, if frequency setups are adjusted in realtime to avoid radio frequency interference at individual stations, or if one band at a station is missing because of radio frequency interference or har dware problems – the source position estimates may change substantially. This may have a significant impact on the creation and use of the VGOS CRF. The results from the simulation with different reference frequency values for the broadband phase, ν0 in Eq. (2), are discussed in Appendix B. The magnitudes of the four coefficients can be significantly reduced for VGOS phase-delay positions by increasing the reference frequency for the phase observables in the VGOS data processing. On the other hand, the coefficients for the other two observables do not change.

Fig. 3 Simulation of broadband group delays of the baseline GGAO12M–ISHIOKA by assuming the position of source 0016+731 in band B to be offset by +0.2mas in declination (left) and the position in band D to be offset by +0.1 mas in declination (right). In these two cases, the reference position is selected to be the positions in the other three bands. The red dots, corresponding to the 47 observations for the source on the baseline in session VO0051, show the delay offsets induced by the applied position shifts, in band B for the left panel and in band D for the right panel. The blue dots are the broadband group delays based on the simulation for the 47 observations. The models based on these position estimates are shown as blue curves.

4 Discussion

4.1 Variations in VGOS Source Positions

The source positions determined by the VGOS observables are linearly dependent on the locations of the four-band radio emission, as shown in Eq. (7) for group delays and in Eq. (8) for phase delays. The summations of the four linear coefficients are unity for both cases of group delay and phase-delay source positions; however, three of them have absolute values larger than 1, and two have negative values. A major consequence is that if there are position variations due to structure evolution in the three highest bands, they would cause variations in the estimates of the VGOS source positions that are larger than the actual position variations in individual bands. Figure 4 shows a one-dimensional diagram of the VGOS group delay positions in the four scenarios of position offsets between various bands. In three of these four simplest scenarios, the VGOS group delay positions are located far away from the area of the actual radio emission – k₀ to k_A|B|C|D. This will add complexities to the understanding of the radio and optical position differences in the future.

Based on the routine geodetic data analysis of session VO0051 using vSolve⁸,25 of the 63 sources with more than three observables available have position adjustments with respect to their S/X positions with magnitudes larger than three times the respective uncertainties estimated from the solution. The median angular separation of the positions from VGOS and S/X observations for these 25 sources is 0.387 mas with an uncertainty of 0.056 mas. Since a VGOS CRF is not yet available, it is a common strategy when analyzing VGOS observations to estimate source positions, which is usually unnecessary for the analysis of routine legacy S/X observations (Sergei Bolotin, priv. comm.).

As shown by Fomalont et al. (2011) and based on multifrequency observations for four close CRF sources, the ICRF position of a source can be dominated by a jet component displaced from the radio core by ~0.5 mas and moving with a velocity of 0.2 mas yr⁻¹. This offset and the change will be significantly amplified in VGOS position estimates, leading to much larger position variations in a VGOS CRF than in the previous ICRF realizations. However, it should be noted that this study is for four ICRF sources only. Larger statistics would be necessary to further assess this point.

Fig. 4 One-dimensional diagram of the relation between the VGOS group delay positions (blue dots) and the locations of the radio emission in individual bands (red rhombuses). Four scenarios are shown: the location of the radio emission in band A, B, C, or D as marked by kA|B|C|D is offset by 0.1 mas with respect to the locations in the other three bands, which are located at the origin marked by k0. The VGOS group delay positions are shown as blue dots: for example, $k_D^{\text{VGOS}}$ is where the VGOS group delay position is located when the location of the radio emission in band D is offset by 0.1 mas with respect to the other three bands.

4.2 Aligning the images in the four bands

The simulation so far does not explicitly model the effects of relative source structure – the apparent two-dimensional distribution of emission on the sky as a function of frequency, as opposed to an ideal point source per frequency. Relative source structure affects the phases of individual frequency channels in a frequency-, baseline-length-, and baseline-orientation-dependent manner, and these phases in turn affect the broadband group delay and phase. Given a model of the relative source structure, such as images of a source in the four VGOS bands, the frequency channel phases can be corrected to mimic the phases of a point source. However, as images derived using self-calibration (e.g., Wilkinson et al. 1977; Cornwell & Wilkinson 1981; Pearson & Readhead 1984; Thompson et al. 2017) and closure-based imaging (Chael et al. 2016, 2018) techniques lose absolute source position information, it is necessary to properly align the images in the four bands in order to generate coherent source-structure phase corrections across all four bands. Otherwise, a misalignment of the images in the four bands will introduce a change in the position estimate from the VGOS observables with these corrections applied. This has been demonstrated by Xu et al. (2021d). The exact formulas of the impact of the misalignment are given above as Eq. (7) for group delays and Eq. (8) for phase delays. These results are important for deriving VGOS source-structure corrections.

Because the core of source 0016+731 at 3.3 GHz as modeled by difmap is shifted toward the east by about 0.20 mas, the corresponding impact on the VGOS source position can thus be calculated according to Eqs. (7) and (8); it is about 0.10 mas toward the east for group delays and about 0.18 mas toward the west for phase delays. Such effects are five to ten times larger than the accuracy of 0.02 mas required to achieve the VGOS goal of 1 mm station position accuracy, and we must emphasize that the discussion here assumes that the effects of source structure are already corrected based on the aligned images through the modeled cores in the four bands. The impact will be even more complex and larger if they are not corrected.

The third type of VGOS broadband observables, δTEC, may help determine the differences in source positions in the four bands or align the images over frequency because these observables are affected by source position offsets as well. Its change due to position offsets among the four bands, denoted by ΔTEC, is given by $$\text{ΔTEC\hspace{0.17em}}=-B\cdot \left(-118.7\text{Δ}{k}_{\text{A}}+168.3\text{Δ}{k}_{\text{B}}+200.7\text{Δ}{k}_{\text{C}}-250.3\text{Δ}{k}_{\text{D}}\right)/c,$$(11)

where ΔTEC is in units of TECU given that the baseline vector is in units of km, the speed of light is in units of km s⁻¹, and the position offsets are in units of mas. This equation is derived in the same way as those for broadband group and phase delays in Sect. 3, but with a scaling factor applied to the coefficients for ΔTEC in Table 2 to be in accordance with the units stated here. The summation of the four coefficients on the right-hand side is zero; therefore, ΔTEC is independent of the reference source position selected, which is canceled out if introduced explicitly in Eq. (11) by replacing Δk_A by k_A – k₀, Δk_B by k_B – k₀, and so forth. In the case where the position offsets are a function of the frequency to the power of −1, as discussed in Sect. 4.3, ΔTEC is equal to zero. With external information about the ionospheric effects for VGOS observations available, this quantity may help validate the image alignment in the presence of core shift.

If there is a position offset in only one of the four bands with a magnitude of 0.1 mas, ΔTEC for the observations of a 6000 km baseline will have a pattern of a sinusoidal wave over 24 hours of Greenwich mean sidereal time with a magnitude in the range 0.2–0.5 TECU based on Eq. (11), depending on which band the position offset occurs in. Through the comparison of δTECs from VGOS observations and global TEC maps⁹ derived from global navigation satellite system observations, the precision of δTEC from global TEC maps is probably on the order of 1–2 TECU (Zubko et al., in prep.). Therefore, position offsets among the four bands larger than 0.4 mas may be detectable by comparing δTEC estimates from VGOS and global TEC maps.

4.3 Differences in VGOS group delay and phase-delay source positions

Here we consider the scenario where the position change of a source over frequency ν is a power-law function ΔX/ν, where ΔX is a constant position shift, following Marcaide & Shapiro (1984) and Lobanov (1998). This change corresponds to the core shift effect with a frequency power law of −1. In the band with a center frequency v_center, which can be calculated from Table C.1, the source position offset then is ΔX/ν_center. Applying the position offsets in this scenario to Eqs. (7) and (8) gives Δk^g = 0 and Δk^p = ΔX/6.0. The results show that if core shift is a function of the frequency to the power of −1, group delays refer to the position of the AGN jet base, which is the position at the infinite frequency, but the phase delays refer to the position of the actual radio emission at the reference frequency. This result was originally derived by Porcas (2009). We can see that Eqs. (7) and (8) describe the differences in the estimated positions between group delays and phase delays for general scenarios. This generalization is necessary because astronomical results such as those from Fromm et al. (2013) and Plavin et al. (2019b) demonstrate that the position dependence on frequency is often not to the power of −1, and in fact the power index can be variable in time. Instead of considering core shift as a function of frequency, we can describe it directly as position offsets in individual bands.

As shown in Eqs. (7) and (8), the position offset in one band leads to position changes in the opposite directions for the broadband group delays and phase delays. These two kinds of broadband observables thus actually refer to different directions, with possibly a large separation on the sky. If the phase observables from VGOS can be used in geodetic solutions, as was recently done for very short baselines (see, e.g., Varenius et al. 2021; Niell et al. 2021; Xu et al. 2021c), the source position differences need to be addressed. Based on Eqs. (7) and (8), the difference in source positions between VGOS group delays and phase delays, denoted by Δk^g-p, is given by $$\text{Δ}{k}^{\text{g-p}}=+1.388\text{Δ}{k}_{\text{A}}-3.177\text{Δ}{k}_{\text{B}}-3.501\text{Δ}{k}_{\text{C}}+5.290\text{Δ}{k}_{\text{D}}.$$(12)

Since the sum of the four coefficients is zero, Δk^g-p is a quantity independent of the reference position. The difference in source positions between VGOS group delays and phase delays measures an absolute-scalar product as a combination of the four-band position offsets, including core shift. The measurement noise levels of both types of VGOS observables are 1−2 ps or even smaller, which allows this position product to be detected at the few μas level (Xu et al., in prep.).

5 Conclusion

We have derived the formulas of the source position estimates from VGOS broadband group delays and phase delays as a function of the locations of the radio emission in the VGOS frequency bands. The resolution across the source is very different in different VGOS bands (a factor of three within the current frequency range and nearly a factor of five between 3 and 14 GHz for the future range), and parts of the jet components may contribute to the flux density and position of the modeled core at the lower frequencies, in which case this will move the “core” position down the jet at these frequencies, especially if such jet components have steeper spectra. Source position offsets between various VGOS bands are expected to be common. The variations in the VGOS source positions will be significantly larger than the actual changes in the emission locations due to structure evolution if they happen in one of the three highest frequency bands. We should expect VGOS sources to wander around parallel to the jet directions as a function of time with much larger magnitudes than the S/X positions because of variations in source structure and core shift being amplified for the broadband observables, as demonstrated by Eqs. (1) and (9). The only way to mitigate these frequency-dependent impacts on VGOS source positions is to measure source structure and core shift (see, e.g., Kovalev et al. 2008; Fomalont et al. 2011; Sokolovsky et al. 2011) and then correct for these effects. This is critical in order to make full use of the high quality data from VGOS, which show a random noise level of ~2ps.

If we want to derive source-structure corrections for VGOS observations, aligning the images in the four bands is essential: a misalignment can translate to a larger offset in the position estimate from the VGOS observations with these corrections applied than the magnitude of the misalignment itself. The difference between source position estimates from VGOS group delays and phase delays is also discussed; it depends on the effect of core shift and, if the corrections for source structure are applied to the observables, the image alignment over frequency.

Acknowledgements

We are grateful to the IVS VGOS stations at GGAO (MIT Haystack Observatory and NASA GSFC, USA), Ishioka (Geospatial Information Authority of Japan), Kokee Park (U.S. Naval Observatory and NASA GSFC, USA), McDonald (McDonald Geodetic Observatory and NASA GSFC, USA), Onsala (Onsala Space Observatory, Chalmers University of Technology, Sweden), Westford (MIT Haystack Observatory), Wettzell (Bundesamt für Kartographie und Geodäsie and Technische Universität München, Germany), and Yebes (Instituto Geográfico Nacional, Spain), to the staff at the MPIfR/BKG correlator center and the MIT Haystack Observatory correlator for performing the correlations and the fringe fitting of the data, to the NASA GSFC VLBI group for doing the geodetic solutions, and to the IVS Data Centers at BKG (Leipzig, Germany), Observatoire de Paris (France), and NASA CDDIS (Greenbelt, MD, USA) for the central data holds. The authors are in particular grateful to the reviewer who helped to improve the manuscript. We would like to thank Bill Petrachenko for his helpful comments on this work. We used the Astrogeo VLBI FITS image database for our work, and specifically we thank Bo Zhang for providing the image of source 1157-215 at 8.7 GHz. This research has made use of data from the MOJAVE database that is maintained by the MOJAVE team (Lister et al. 2018). The research was supported by the Academy of Finland project No. 315721, National Natural Science Foundation of China (No. 11973023), and the German Research Foundation grants HE5937/2-2 and SCHU1103/7-2.

Appendix A Simulations with Different Settings of Channel Frequencies

The IVS proposes observing with 1 GHz wide bands and eventually extending to higher frequencies, as originally planned (Petrachenko et al. 2009).

We first performed the simulation for the frequency range 3.0-11.2 GHz but with 992 MHz wide bands, as listed in Table C.2 (Bill Petrachenko, private communication). The equations equivalent to Eqs. 7 and 8 for this wider bandwidth are given by $$\text{Δ}{k}_{\text{992,11}}^{\text{g}}=+0.399\text{Δ}{k}_{\text{A}}-1.382\text{Δ}{k}_{\text{B}}-1.359\text{Δ}{k}_{\text{C}}+3.342\text{Δ}{k}_{\text{D}}.$$(A.1)

and $$\text{Δ}{k}_{\text{992,11}}^{\text{p}}=-0.796\text{Δ}{k}_{\text{A}}+1.177\text{Δ}{k}_{\text{B}}+2.018\text{Δ}{k}_{\text{C}}-1.955\text{Δ}{k}_{\text{D}}.$$(A.2)

We then performed the simulation for the frequency range 3.0−14.0 GHz, as listed in Table C.3 (Bill Petrachenko, private communication). The results are given as follows: $$\text{Δ}{k}_{\text{992,14}}^{\text{g}}=+0.248\text{Δ}{k}_{\text{A}}-1.008\text{Δ}{k}_{\text{B}}-0.995\text{Δ}{k}_{\text{C}}+2.755\text{Δ}{k}_{\text{D}}.$$(A.3)

and $$\Delta k_{992,14}^{\text{p}}=-0.633\Delta k_{\text{A}}+1.663\Delta k_{\text{B}}+2.078\Delta k_{\text{C}}-2.108\Delta k_{\text{D}}.$$(A.4)

Only a marginal improvement in terms of reducing the magnitudes of the coefficients is gained by using a wider 992 MHz bandwidth for the current VGOS frequency range. However, significant improvement is achieved by moving the upper frequency band to 14 GHz for group delay positions; the magnitudes of the four coefficients decrease by 20% to 50% compared to those of the current frequency setting in Table C.1. This means that the impact of any position offset in a single band has less influence on the VGOS group delay position in this case.

For completeness, the equations for the δTEC in these two scenarios (equivalent to Eq. 11) are given by $${\text{ΔTEC}}_{992,11}=-B\cdot \left(-112.7\text{Δ}{k}_{\text{A}}+169.0\text{Δ}{k}_{\text{B}}+199.3\text{Δ}{k}_{\text{C}}-255.6\text{Δ}{k}_{\text{D}}\right)/c$$(A.5)

and $${\text{ΔTEC}}_{992,14}=-B\cdot \left(-96.1\text{Δ}{k}_{\text{A}}+150.4\text{Δ}{k}_{\text{B}}+194.7\text{Δ}{k}_{\text{C}}-250.0\text{Δ}{k}_{\text{D}}\right)/c.$$(A.6)

Appendix B Comparison with different values of v₀

A simulation was performed to compare the results when varying the reference frequency ν0 in Eq. 2 rather than keeping the nominal value of 6.0GHz. The relation for group delay positions, given by Eqs. 4 and 7, does not change when ν₀ is varied, nor does the relation for the ionospheric observables. However, the coefficients related to the VGOS phase-delay position (Eqs. 8 and 10) are found to decrease when ν₀ increases. When setting, for example, v₀=8.5 GHz, the VGOS phase-delay position is given by $$\text{Δ}{k}_{v_0=8.5\text{\hspace{0.17em}GHz}}^{\text{p}}=-0.475\text{Δ}{k}_{\text{A}}+0.794\text{Δ}{k}_{\text{B}}+1.014\text{Δ}{k}_{\text{C}}-0.333\text{Δ}{k}_{\text{D}}.$$(B.1)

By increasing v₀ from 6.0GHz to 8.5 GHz, τ^p, which is inversely proportional to v₀, simply decreases by 30%. Additionally, the magnitudes of the four coefficients are actually reduced by 16%, 24%, 30%, and 52%, respectively. However, changing the reference frequency from the central frequency 6.0 GHz will introduce the errors from group delay estimates into the phase delay estimates. The optimum reference frequency for phases will have to be a compromise between these two factors.

Appendix C Frequency Settings for VGOS Observations

The frequencies of the 32 channels in the current VGOS observations are shown in Table C.1. Two other possible sets of observing frequencies, which are being tested and may be used by VGOS in the future, are listed in Tables C.2 and C.3.

References

All Tables

All Figures

Fig. 1 Radio positions from the ICRF3 in the S/X, K, and X/Ka bands and the Gaia position for source 1157-215, on its radio images at 8.7 GHz from Astrogeo (left) and at 15.3 GHz from MOJAVE (right). Overlay contours are shown at levels of peak percentage, specified at the bottom of the plot. Since these radio images do not have information regarding the absolute positions, the peaks are formally assumed to be located in the S/X-band positions. Error bars are the uncertainties from the four catalogs. The K- and X/Ka-band positions are the ones after the full deformation transformations (Charlot et al. 2020). The angular separations of the K- and X/Ka-band positions and the Gaia position relative to the S/X-band position for this source are 2.29 ± 0.31 mas, 3.54 ± 0.32mas, and 3.97 ± 0.13 mas, respectively.

In the text

Fig. 2 Images of source 0016+731 from VGOS session VO0051 (February 20, 2020) at 3.3 GHz (upper left), 5.5 GHz (upper right), 6.6 GHz (bottom left), and 10.5 GHz (bottom right). The ellipses indicate the Gaussian components, the parameters of which are reported in Table 1. Overlay positive contours are shown at nine levels of peak percentage, starting from 0.2% and increasing by a factor of 2, while the negative contour, indicated as dash lines, has a peak percentage of – 0.2%. These images are constructed from the VGOS observations, which are calibrated based on the images derived from closure phases and closure amplitudes (Xu et al. 2021d). Therefore, the pixel flux densities have arbitrary units. The beam size is shown as a filled ellipse in the bottom-left corner of each plot.

In the text

Fig. 3 Simulation of broadband group delays of the baseline GGAO12M–ISHIOKA by assuming the position of source 0016+731 in band B to be offset by +0.2mas in declination (left) and the position in band D to be offset by +0.1 mas in declination (right). In these two cases, the reference position is selected to be the positions in the other three bands. The red dots, corresponding to the 47 observations for the source on the baseline in session VO0051, show the delay offsets induced by the applied position shifts, in band B for the left panel and in band D for the right panel. The blue dots are the broadband group delays based on the simulation for the 47 observations. The models based on these position estimates are shown as blue curves.

In the text

Fig. 4 One-dimensional diagram of the relation between the VGOS group delay positions (blue dots) and the locations of the radio emission in individual bands (red rhombuses). Four scenarios are shown: the location of the radio emission in band A, B, C, or D as marked by kA|B|C|D is offset by 0.1 mas with respect to the locations in the other three bands, which are located at the origin marked by k0. The VGOS group delay positions are shown as blue dots: for example, $k_D^{\text{VGOS}}$ is where the VGOS group delay position is located when the location of the radio emission in band D is offset by 0.1 mas with respect to the other three bands.

In the text