VLBI and GPS inter- and intra-technique combinations on the observation level for evaluation of TRF and EOP

We study the effects of combination on the observation level (COL) of different space-geodetic techniques and of networks of the same technique and present the corresponding improvement for the determination of station positions and earth orientation parameters. Data from the continuous geodetic very long baseline interferometry (VLBI) campaign CONT17 are used in a batch least-squares (LSQ) estimator. This campaign includes 15 days of observations with two legacy S/X networks, namely Legacy-1 (L1) and Legacy-2 (L2). For this study the VLBI L1 network is used as the base and reference solution. Data from the L1 network are combined first with data from co-located Global Positioning System (GPS) stations by estimating common tropospheric parameters. The derived station positions repeatabilities of the VLBI and GPS networks are evaluated with respect to single-technique solutions. In terms of precision, we find a 25% improvement for the vertical repeatability of the L1 network, and a 10% improvement for the horizontal one. The GPS network also benefits by 20% and 10% in the horizontal and vertical components, respectively. Furthermore, a combined solution using data of the L1 and L2 network is performed by estimating common earth orientation parameters. The combined L1&GPS and L1&L2 solutions are compared to the reference solution by investigating UT1 and polar motion estimates. UT1 is evaluated in terms of mean bias and formal errors with respect to the International Earth Rotation Service (IERS) C04 products which were used as a priori values. The L1&GPS solution has the lowest formal error and mean bias for UT1 with a 30% improvement. The weighted root mean square (WRMS) and weighted mean offset (WMO) differences between the obtained polar motion estimates and the ones derived by the International GNSS Service (IGS) are also compared. We find that the L1&GPS solution gives the lowest WRMS and WMO, exhibiting an average 40% improvement with respect to the reference solution. The presented results highlight the potential of COL for ongoing transition to multi-space geodetic analysis, e.g., Global Navigation Satellite Systems (GNSS) with the next-generation VLBI system.


Introduction
Two-week continuous geodetic VLBI campaigns have been organized by the International VLBI Service for Geodesy and Astrometry (IVS, Nothnagel et al. 2017) every third year since 2002. These continuous (CONT) VLBI campaigns provide a great opportunity to not only investigate the current performance of the VLBI systems (Nilsson et al. 2014), Earth rotation (Artz et al. 2010), as well as ionosphere and troposphere effects (Teke et al. 2011), but also test the feasibility of combining geodetic VLBI with other space-geodetic techniques such as GNSS or SLR (Thaller et al. 2007). The objective of such a combination is to provide a robust estimation and physically consistent treatment of common parameters of interest such as Earth orientation parameters (EOP) or zenith wet delays (ZWD), in order to ensure the highest quality and homogeneity of the final products (Rothacher et al. 2011). In addition, space-geodetic techniques have different strengths and weaknesses for recovering global geodetic parameters. Their combination is thus beneficial to fully exploit the strengths of each of them and overcome the technique-specific weaknesses (Artz et al. 2012), assuming that the combination is properly constructed, suitable weights are applied and good quality ties at co-location sites are available. One such example could be the geocenter coordinates that can be transferred from SLR to the GNSS network, in the case of a combined SLR-GNSS solution (Sośnica et al. 2019). The inclusion of geodetic VLBI into the solution on the other hand provides the full EOP information (Sovers et al. 1998), which can augment the GNSS or SLR data analysis. Utilization of many space-geodetic techniques also improves the global geometrical coverage, the quantity of observations and their density, and helps in the reduction of correlations between the estimated parameters.
The combination of observed data originating from different space-geodetic techniques is usually carried out using three different strategies in order to estimate common parameters of interest: (1) CRL (combination on the result level) refers to a case where target parameters are estimated firstly from single-technique solutions and the results are combined in the subsequent stage to obtain one set of estimates (Seitz 2015).
(2) CNL (combination on the normal equation level) implies that parameters of interest such as EOP, station coordinates or common troposphere models are estimated using technique-specific datum-free normal equation (NEQ) systems while including necessary local tie information (Thaller et al. 2007;Rothacher et al. 2011). (3) COL (combination on the observation level) implies that common inter-technique as well as techniquespecific parameters are estimated simultaneously, preserving all correlations between the estimated parameters (Coulot et al. 2007).
COL has already been performed ), but still has not been yet fully exploited (Coulot et al. 2007). The advent of the next-generation VLBI system, known as VGOS (VLBI Global Observing System), as well as the full deployment of few Global Navigation Satellite Systems (GNSS) like Galileo and BeiDou gives rise to the exciting new prospects which come from (a) the densification of observations and (b) the densification of networks. In this study, we employ COL and use the CONT17 data set to examine those two densification modes with respect to primary parameters of interest, namely the terrestrial reference frame (TRF) and EOP. Regarding the densification of observations, data from GPS stations which are co-located to the VLBI stations are used. As a consequence, the significantly improved spatio-temporal observation resolution allows to derive troposphere-related parameters with improved quality, affecting positively also the derived station positions and EOP. Regarding the densification of networks, the two S/X VLBI networks that participated in CONT17 are combined via the common access to EOP. The resulting extended network with superior global geometrical coverage improves the sensitivity to EOP determination.
We provide the description of the input data used in this study in "Input data" section. "Methods" section discusses various aspects of the combination of VLBI and GPS as well as multiple VLBI networks on the observation level. "Results" section presents the obtained results, where parameters of interest are derived from single-and multi-technique solutions. In "Discussion" section we discuss the challenges related to the combination of both space-geodetic techniques and comment on the discrepancies between the parameters of interest derived from single-technique solutions. Finally, "Conclusions" section provides the reader with the summary and conclusions, and outlines future work concerning this topic.

Input data
The VLBI data set used in this study covers 15 days of CONT17 from a fourteen-station S/X network, referred to as the Legacy-1 (L1) network and a thirteen-station S/X network referred to as the Legacy-2 (L2) network . We focus on the first network as the reference for this study since all involved sites have both geodetic VLBI and GPS instrumentation, which is not the case for the second one. For this reason, the L2 network serves only as a complementary network to achieve improved global coverage and its contribution is measured via the effect on EOP estimation. Thus, station positions for the L2 network are neither presented nor evaluated. The considered time period spans from November 28 at 00:00 UTC to December 12 at 23:59 UTC. We use the data from the co-located GPS stations, which are part of the International GNSS Service (IGS). Table 1 shows the participating sites of the Legacy-1 network and their VLBI and GPS stations, while Fig. 1 depicts their geographical distribution.
The CONT17 campaign lasted from November 28th to December 13th, 2017. Its concurrency with other types of observing sessions, like, e.g., the intensive sessions between Wettzell and Kokee, and station-specific technical problems, like, e.g., the Fortaleza station being out-of-order for approximately 3 days, meant that certain data-gaps occurred which are presented in Fig. 2.
The geodetic VLBI data are stored in the vgosDb format (Gipson 2014), where each 24-h set of observations corresponds to a single vgosDb database. We used the analysis-ready version of the databases, commonly known as version-4 databases, where the group delay ambiguities are already resolved, ionosphere delays estimated and prospective clock breaks detected. The geodetic VLBI analysis was carried out using the S/X source catalogue of the 3 rd realization of the International Celestial Reference Frame (ICRF3, Charlot et al. , 2020). In order to analyze the GPS observations in the static Precise-Point Positioning (PPP) mode, we used observation data in the Receiver INdependent EXchange (RINEX) format from the same period and IGS final products comprising satellite ephemerides and both satellite and IGS station clocks (International GNSS Service 2020).

Methods
Combination of space-geodetic techniques promises to provide a consistent set of nuisance parameters such as station clocks, ZWD and related tropospheric gradients in north (GRN) and east (GRE) directions. The utilization of a common troposphere model at co-located stations and the inclusion of inter-technique weighting allow to identify the presence of technique-specific errors and deficiencies in the modeling as well as outliers and artifacts that might be present in the data. As outlined before, the combination on the observation level (COL) is the most promising approach for this purpose. However, before the combination takes place, it is necessary to first assess the single-technique solutions and study the behavior and quality of the common parameters of interest.

Common troposphere parameters
With the co-located GPS and VLBI stations it is possible to use one troposphere model that can be applied to both techniques, provided that there is a small spatial  distance between the co-located GPS and VLBI reference points, i.e., they share the same atmospheric conditions. The common troposphere model requires consistent modeling of the hydrostatic delays. This is achieved by applying the GPT3 model (Lagler et al. 2013) and using model-based pressure values that are calculated for a given height of the reference points of the co-located stations. The height differences also imply slightly different ZWD for the co-located stations. However, the effect was found to be rather small for the 14 sites used in the study, mainly below the standard deviation of the estimated ZWD, and thus is neglected so far in the study. The values can be calculated according to Rothacher et al. (2011) and are listed in Table 2 along with the RMS differences between the VLBI-derived and GPS-derived tropospheric parameters (ZWD, GRN, GRE) at the co-located stations.

Common clock parameters
In principle, it is possible to estimate a common clock correction for both considered techniques (with some inter-technique clock offsets) as usually reference signals from Hydrogen masers (H-masers) provide the same frequency standard to both GPS receivers and VLBI backends at the core sites .
To assess the feasibility of using a common clock model for our purposes, an analysis of the clock behavior of all participating GPS receivers was carried out using the GIPSY/OASIS II software package (Webb and Zumberge 1997).
The Modified Allan Standard Deviation (MDEV) calculated for three badly performing and a well performing GPS stations is shown in Fig. 3. The unusual behavior of some of the GPS receivers includes the case of ZECK, where a loss of lock to the reference signal was detected for the GPS receiver on day 7 of CONT17. In  Table 2 The expected ZWD differences ( δZWD ) due to the height difference between the VLBI and GPS reference points (Rothacher et al. 2011), shown together with the WRMS differences and mean biases between the VLBI-derived and GPSderived ZWD addition, clock breaks were detected for stations BADG and HRAO. Moreover, not all of the GPS receivers at the CONT17 sites were connected to the H-maser. Therefore, the combination of clocks was not pursuit in this study and clocks for the co-located VLBI and GPS stations were modeled as separate parameters. Information on the VLBI and GPS clock parameterization used in this study is given in "Parameter estimation" section.

Combination of different VLBI networks
The L1 and L2 networks do not have co-located stations for any site, apart from the WE, and thus common troposphere resolution cannot be obtained. The strength of this combination approach comes from the fact that two different geometries are combined into one and hence EOP estimates of an improved quality are possible. However, a special consideration must be given when formulating the problem with respect to how minimal constraints should be applied in order to remove the rank deficiency present in that case.

Application of minimal constraints
Geodetic observations usually prevent from a reliable definition of all components of the coordinate system (origin, orientation, and scale) with respect to which the station positions are estimated. This issue, in the geodetic network adjustment, leads to a rank-deficient system of NEQ. The most common approach for addressing this issue is to apply the so-called minimal constraints (MC) realized in the form of No-Net-Translation (NNT), No-Net-Rotation (NNR) or No-Net-Scale (NNS) conditions (Altamimi et al. 2002). The philosophy behind MC is that the number of the constraints is exactly equal to the rank defect of the normal matrix. The normal matrix which is constructed using observations from a single VLBI network has a rank defect of 6, in the case when both EOP and station positions are estimated. Thus MC are fulfilled by applying NNT and NNR, i.e., specifying the origin and orientation of the realized frame to be that of the a priori frame. When combining VLBI networks and since there is no direct connection with common observations between stations of different networks, there will be still two distinct frame realizations with different origins, thus applying NNT constraints separately on the two networks. For consistency in EOP determination, the orientation of the two frames must be common and thus a single NNR constraint is applied on the combined L1&L2 network. Furthermore, MC are sensitive to noisy observations and inaccurate a priori information on station coordinates, so-called data noise and datum noise effect, respectively (Kotsakis 2018). Thus a correct choice of stations to apply NNT/NNR is important for the quality of the solution. In our case, the NNT/NNR constraints were imposed on twelve out of fourteen VLBI stations on the L1 network and on the complete L2 network. In the case of L1, the KASHIM11 station was excluded from the datum constraint as it experienced severe radio frequency interference at X-band, which led to the exclusion of four X-band channels at the correlation stage. The ZELENCHK station was also excluded as it was found that the a priori position was suboptimal due to significant discrepancy on the (fixed) station velocity information with respect to the co-located ZECK GPS station.

Parameter estimation
The CONT17 data were analyzed as 15 daily solutions independent with respect to each other with no global parameters valid for the whole 15-day period. The analysis was carried out using the multi-technique spacegeodetic analysis software c5++ (Hobiger et al. 2010). An overview of the design of the software is presented in Fig. 4. The inclusion of local ties, while feasible, was not pursued in this study.
Four different analyses were carried out, namely an VLBI L1-network only (sol-V), a GPS-only (sol-G) a combined VLBI L1&GPS (sol-VG) and a combined VLBI L1&L2 (sol-VV). All available VLBI data were used (scheduled with an elevation cutoff angle of 5 • ). The GPS data however, were decimated from the original 30-s sampling to the 5-min sampling and an elevation cut-off angle of 10 • was applied to all GPS stations. A two-stage parameter estimation process for all solutions was employed. First, the data were used to estimate a second-order polynomial for station clocks. For L1 and L2 networks in particular, one station clock per network was fixed and used as reference. Then using the polynom as a priori information and thus removing the trend from the clock parameters a second-stage estimation was performed. Station clocks, zenith wet delays and tropospheric gradients were modeled as continuous piecewise linear (PWL) offsets, while EOP were estimated as constant daily offsets. The time resolution of the PWL parameters varied according to the temporal resolution of the dataset and is presented in detail in Table 3.
In the processing of GPS data, no integer-fixing of phase ambiguities was performed, and instead they were resolved as real numbers. In addition, an elevationdependent weighting was used where formal errors of each observation are multiplied with the wet mapping function evaluated for a given observation elevation angle. For geodetic VLBI data, an elevation-dependent weighting was also applied, with the difference that it is baseline-based. This means that, besides the formal uncertainty derived at the post-correlation analysis stage, the observation error in our analysis includes quantities evaluated for two stations forming a baseline, i.e., wet mapping functions computed for the given elevation angles. No correlations between the individual observations were considered in our case. In the case of VLBI, baseline-based variance component estimation (VCE) of the Helmert type is employed. In the case of GPS, the VCE is computed per station and per type of observable (code-or carrier-phase). A general description of the VCE estimator is presented in Eq. 1, where v i is the residual vector, P i the inverse of the variance matrix and n i the number of data-points for the i th group of observations, while n is the total number of observations, and u is the total number of unknowns. For more details, see, e.g., Bähr et al. (2007).
Inter-technique weighting was applied as an extension of the single-technique VCE where the total number of observations was considered to be the sum of both VLBI and GPS and the same applied for the number of (1) The basic concept of c5++ (figure from Hobiger and Otsubo (2014)). The c5++ analysis software allows for COL using common geophysical models for all the supported space-geodetic techniques, namely SLR, GNSS and VLBI unknowns. As a result of these adjustments, the chisquared statistic for both single-technique and combined solutions was approximately equal to one. The a priori positions at each observation epoch were calculated using station positions and velocities expressed in ITRF2014, including the post-seismic displacement model (Altamimi et al. 2016) for sol-VG and sol-G, while the VTRF2019 was used for sol-V and sol-VV. In terms of the periodic displacements triggered by geodynamical phenomena, the modeling approach applied here followed the IERS conventions (Petit et al. 2010). The tidal S1-S2 atmospheric loading (Ray and Ponte 2003) was also used and applied to the positions of both VLBI and GPS reference points. EOP were estimated as corrections to the a priori IERS-14-C04 (Bizouard et al. 2019) values. The whole list of the a priori used is presented in Table 3.

Results
The effect of combination on the L1 network by employing COL of VLBI with co-located GPS stations is examined in the following context. First, the derived station positions are evaluated w.r.t. the results of single-technique solutions; during the performed COL the common troposphere model (ZWD, GRN, GRE) was derived based on both GPS and VLBI observations. The EOP determination is then compared w.r.t. the single-technique solution and to a solution where network densification is achieved via the combination of both interferometric networks L1&L2 (sol-VV).

Station position repeatabilities
The 3D position repeatabilities are shown both (a) for the co-located L1 and GPS networks (Fig. 5) and (b) for station-specific cases (Figs. 6,7,8). In the first case, the improvement that COL induces is visible between the sol-V, sol-G and sol-VG solutions by the different color opacity. More specifically, the analysis indicates that for the L1 network t, on average, the repeatability in the up component is reduced by 25% In the case of the horizontal components, we note an improvement of 16% and 6% in the north and east components, respectively. The GPS network on the other hand experiences the biggest improvement in the north component by 40%, followed by the up component by 10%, while the east shows a slight degradation of 5%. Looking at station-specific results, where V stands for VLBI and G for GPS stations, direct inspection of Fig. 8 shows the relative changes. In 80% of the cases a significant improvement is visible, varying between 10 and 60%. For 20% of the cases, a small degradation is visible mainly in the GPS network with the exception of HA and YE sites, where the degradation is visible on the co-located VLBI station as well. Mostly, it does not exceed 18%, except for the HA site where the GPS station shows a significant increase in the up component. Whether these relative changes are significant or not can be examined looking at the absolute level of the repeatabilities shown in Fig. 7. The L1 network is also examined in terms of the baseline length repeatabilities. The combined solution improves this metric considerably as can be seen in  . Evaluated at 6000 km, COL reduces the baseline repeatability by 35% from 6.5 to 4.2 mm. Further elaboration will follow in "Discussion" section.

Zenith wet delays
As shown in Table 4, the mean WRMS difference between the ZWD derived from sol-VG and sol-G is approximately 2 mm. The average bias approaches 1 mm. As sol-VG follows closely sol-G, it is expected that the differences of the troposphere derived from those two solutions w.r.t. sol-V will be similar. This effect is illustrated in Fig. 9, where the dominant role of sol-G to the sol-VG derived ZWD is evident, due to the characteristics of GPS observations, i.e., simultaneous observations of many satellites distributed on the local skies and a good temporal resolution of observations in relation to the VLBI technique.

Frame-defining parameters
A 7-parameter Helmert transformation was performed to examine the effect of sol-VG on the frame-defining parameters of the L1-network and the results are presented in Table 5. It is observed that the induced frame bias, as indicated by the parameters themselves, is all but eliminated for the determination of the origin and the scale of the frame and is approximately one order of magnitude less for its rotation. The sensitivity, as indicated by the uncertainty levels, becomes higher with the benefit being more pronounced on the translation and scale parameters where the uncertainty levels drop by an order of magnitude. It should be noted though that  they are already sufficiently small on the sub-millimeter level. On the other hand, the gain on the uncertainty level of the rotation parameters, ranging from 40 to 75%, while smaller may have positive impact on the EOP determination.

EOP repeatabilities
EOP were evaluated augmenting the information available from the L1 network in two different ways: (a) densifying the interferometric network itself by incorporating L2 network into the solution and (b) better resolving the troposphere by the combined L1&GPS solution. One should keep in mind that EOP products were estimated via the L1 and/or L2 networks. The co-located GPS network was only used to better resolve troposphere and thus its effect is only indirect on the EOP estimation. While EOP products can be accessed through global  solutions of GPS networks, those solutions routinely involve a dense network of stations, something that cannot be achieved with the geometry available from the CONT17 campaign and the combination schemes pursued in this study. The process of evaluating EOP products starts with looking at UT1-UTC estimates. These estimates which were obtained from the 15 daily solutions should in principle constitute a zero-mean set (Thaller 2008). Figures 11, 12 with the corresponding weighted mean offset (WMO) and formal errors show that WMO w.r.t. the a priori values are about −19 µ s for both sol-V and sol-VV, while sol-VG is causing it to drop to −12.8 µ s. These results show that both combination schemes do not induce any additional bias in the solution and thus can be used to evaluate polar motion estimates. Further elaboration follows in "Discussion" section.
Polar motion (PM) estimates (x p , y p ) are evaluated with the aid of IGS PM products. These products have been obtained through global solutions on a dense GPS network and thus provide an excellent external reference in order to quantify the effect of the combination. The differences of the PM estimates to the IGS PM products were calculated and presented in the form of WRMS and WMO. Figures 13, 14 show the effect that the combination schemes have on these two metrics. More specifically, augmentation of the L1 with the L2 network shows the biggest effect on the mean bias of x p which lowers by 60%, while an improvement of 10% is also visible in the WRMS metric. The x p component shows a 19% improvement to the WRMS metric. COL with the co-located GPS stations also affects most the WMO of x p which lowers by 85%. In this case, y p shows a bigger improvement than before, with the WMO and WRMS both decreased by 50%. The celestial pole offsets showed non-significant changes between the different combinations and thus are not discussed in this study.

Discussion
COL was employed to augment VLBI observations performed with the use of the VLBI CONT17 Legacy-1 network (L1) in two distinctive manners. First, data from the co-located GPS stations were used and common tropospheric parameters were estimated. The reasoning behind this approach was that GPS data with their superior temporal and spatial resolution and the simultaneous tracking of multiple satellites with a multitude of elevation and azimuth angles would allow for a better resolution of the tropospheric elongation. Correlation between this parameter and station positions would therefore result in a lower uncertainty for the latter, in principle for both co-located stations. In general, an improved ZWD estimation affects total tropospheric delay considerably and subsequently the local height of a station. This effect is visible in the significant improvement that L1 network experiences in the up component. All 14 participating stations show a reduction of the repeatability between 10 and 50%. The components of the horizontal plane in the local topocentric frame were also evaluated. Improved estimation of tropospheric gradients should allow for lower repeatabilities on the north and east components. This is confirmed for most of the cases of the L1 network with the exception of FO (north), HA (north and east), MA (east) and YE (north). Out of those, the degradation for HA is less than 5% while for YE and FO the repeatability of the horizontal position components is already significantly lower than the vertical one, more than a factor of two. The east component of the MA site also follows the behavior of the horizontal position components of YE and FO and its degradation is counter-balanced by the more significant improvement in its north position component. Overall, a clear improvement can thus be observed for the L1 network station positions as an effect of the COL.
The results for the co-located GPS stations could also be evaluated. The improvement of the repeatabilities due to the presence of quasar observations may not be so visible since tropospheric estimates of the combined solution tend to be dominated by the presence of GPS. However, since the COL is employed, the correction to the a priori EOP, which was performed on the VLBI module of the software, is shared across all techniques meaning that the GPS component has also access and should benefit by it. One can see a clear benefit on the up component in 10 out of the 14 stations where an improvement between 10 and 50% is observed while the total improvement on the vertical component is approximately 15%. The total repeatability of the east position component for the GPS stations, co-located with the L1 network, shows a slight degradation and this is visible in the station-specific evaluation as well, where a small increase in the repeatabilities is observed for 9 out of the 14 stations. The north component, on the other hand, fares a lot better with 12 of the 14 stations showing an improvement while the total improvement on this component is approximately 50%, the biggest among the three. This can be attributed to two reasons. First, during sol-VG, the GPS module inherits the EOP corrections from the VLBI module, and thus using values no longer fixed to a priori, one would expect to see the biggest benefit on the horizontal component of the repeatabilities. Second, besides their general uniform distribution, natural radio sources used in VLBI are also observable in the northern (or southern) parts of the local skies where GPS satellites, due to the orbit inclinations, are not visible. Thus augmentation with quasar observations for the GPS stations allows for a total improvement in the resolution of the tropospheric parameters but there might be cases where station-specific biases of the co-located VLBI station propagate into the troposphere estimates and may cause a slight perturbation on the expected improvement of station position repeatabilities.
As already seen in Fig. 5, sol-V for the L1 network shows lower sensitivity on the up component with the WRMS slightly above 8 mm, a factor of four more than the horizontal components. This translates to baseline length repeatabilities having a monotonic quasi-linear behavior due to the major contribution of the local topocentric up components on the baseline vector, as shown in Fig. 10. The effect of COL here is visible in two manners, namely the total level of the WRMS and the scatter around the fitted polynom (spread). The reduced total level w.r.t. the single-technique solution can be deduced from what was previously discussed, namely the effect of better resolved tropospheric parameters in all, but especially, the up component of the station positions. The scatter shows the presence of system-specific biases, like, e.g., the gravi-thermal deformation, that may affect a station position estimate. More specifically, a systematic error in a station position results in clusters of baselines to that station showing a systematically higher WRMS than other of comparable length. For example, the VLBI station at ZE is the common factor on the cluster of baselines between 2000 and 4000 km that are visibly separated from the rest in sol-V. A similar situation occurs in the case of the VLBI station at KS, resulting in separated clusters on the 8000 to 12,000 km baseline range. For this station there were known issues with the presence of radio frequency interference in 4 X-band channels and subsequent exclusion from the fringe-fitting stage. COL is able, thus to mitigate such effects resulting in a smaller scatter around the fitted polynom and visible declustering.
Evaluating EOP estimation and especially UT1-UTC can be challenging as there is no other technique apart from VLBI which provides a direct estimate. Corrections to the C04 a priori values, obtained from the 15 daily solutions should have a WMO close to zero. The existence of unmodeled system-specific biases like, e.g., the gravi-thermal deformation of the antennas, the sparse resolution of troposphere or clock errors can propagate into UT1-UTC estimates and as is the case in this analysis where the WMO is close to −18 µ s. The feasibility of an attempted COL can be evaluated investigating the effect that it has on the WMO of the UT1-UTC estimates. A meaningful combination scheme should result in a WMO of a lower or equivalent level (Thaller 2008). Figures 11, 12 show that this condition is fulfilled with sol-VV resulting in a similar level and sol-VG lowering the WMO by 30%, namely to −12 µs.
PM estimates, on the other hand, have an external reference to compare to in the form IGS PM products. Sol-V shows a WRMS between 60 and 80 µ as which is comparable to the 40-70 µ as range shown in, e.g., Nilsson et al. (2019). One would hope that both combination schemes would augment the sensitivity of sol-V for both polar motion components. Sol-VV adds more baselines in both east-west and north-south directions, while sol-VG augments the information through correlation to atmospheric parameters. In addition, the solutions containing interferometric networks only, exhibit the so-called data-noise effect in the geodetic datum determination (Kotsakis 2018). More specifically, station-specific datagaps and observation noise cause errors that propagate into the frame-defining parameters via the selection of sites used for NNR/NNT constraints. The COL using GPS, since it contains co-located stations, intervenes on the station-specific level to stabilize the solution. This stabilization is visible in Table 5. The sol-VG, in fact, not only removes systematic bias from the frame realization but also makes the solution more sensitive to EOP determination, as shown by the lower uncertainty level of the rotation parameters. Thus, a second correlation gets augmented in sol-VG, that between frame-defining parameters and EOP. Based on the inspection of the effect on both WRMS and WMO, it is visible that sol-VV is more beneficial on the x p estimation. This is due to the fact that all but one stations of the L2 network are located in the northern hemisphere and hence the combination adds a lot of baselines in the east-west direction, improving the sensitivity to x p . On the other hand, sol-VG shows an improvement on both x p and y p , where a consistent decrease of WMO values on the level of 50-70% is visible for both components as well as the WMO for x p which is almost neglegible in our combination approach. As mentioned before, COL of two distinct VLBI networks, as pursued here, is carried out using two networks that do not share one or multiple stations. The connection between them is established by estimating common ERP and imposing the NNR constraint based on a set of stations chosen from both networks. The presence of a common set of stations would additionally allow for adding a single NNT constraint, strengthening the combination further and allowing for enhanced mitigation of networkspecific biases.

Conclusions
The presented results indicate that the combination on the observation level tends to be beneficial for the analysis of space-geodetic data. Such an approach leads to improved results for parameters of paramount interest in space-geodesy, i.e., TRF and EOP. Both inter-and intratechnique combinations of space-geodetic observations can be performed with this approach.
The inter-technique aspect was addressed by combining data from different techniques at co-located sites, i.e., VLBI and GPS, by estimating common tropospheric parameters. VLBI and GPS data from the CONT17 campaign were used for the presented approach. This combination strategy led to improvements concerning determination of VLBI station positions and related baseline length repeatabilities by up to 25%.
The intra-technique aspect was addressed by combining simultaneously observed data of the same spacegeodetic technique, e.g., data from two different VLBI networks operating at the same time. Two VLBI networks observing during the CONT17 campaign were combined by estimating a common set of EOP. This combination strategy led to an improvement in precision of the derived polar motion and UT1 values by 20% to 30%, illustrated by a better agreement of our polar motion estimates with independent GPS-based estimates.
As an outlook for the future, the plan is to extend the described COL approach within the c5++ analysis software and to combine a multitude of space-geodetic techniques to address the topic in a comprehensive manner: • Including more GPS stations than used in the presented study, not necessarily all co-located with VLBI station but allowing for an adequate density and global distribution of the GPS network. This may allow for determination of common EOP through direct inference from both space-geodetic techniques, i.e., VLBI and GPS, with the possibility of extending the data analysis with nutation components derived with the use of VLBI observations. • Including additionally all the other existing GNSS, i.e., GLONASS, Galileo and Beidou, in the combinations is a logical subsequent step and is of course expected to benefit the tropospheric parameters at co-location stations, as well as the TRF and EOP determination due to an increased amount of the available satellites at local skies and different orbital characteristics of those navigation systems. It will