Thank you very much for your valuable and constructive comments. As you have remarked, there is a number of improvements which we will be happy to consider for a future second re-processing. In the following, we would like to comment on your remarks and evaluate them w.r.t. the present version of our data set.

We have numbered the reviewer’s comments consecutively (as RC1.1 and so on), and repeated them in italic. Our replies are given below each comment.

Responses to the specific comments

RC1.1:

The four input solutions represent different realizations of the IGb14 frame and it is natural that the authors decided to harmonize the frame differences. Some of the input solutions were globally aligned to the IGb14 frame while others were aligned using a regional subset of the IGb14 frame. A regional frame alignment is a very effective way to remove common-mode noise from the daily position time series, and conversely, a global frame alignment is effective in introducing frame noise into the daily position time series. In addition, the quality of the frame alignment could change over time differently among the solutions, especially for those based on sparse regional subsets. The frame noise in the input solutions can be minimized by applying a new regional daily frame alignment to all the solutions, so that all the aligned series have a consistent level of frame noise before their combination. However, the authors decided to remove the long-term alignment bias only, instead of removing an alignment bias on a daily basis (the frame jitter) that would account for both the frame bias and the frame noise altogether.

Of course, the reviewer is right in the sense, that by using the selected IGS sites with their IGb14 frame solution the long-term alignment bias is considered. However, we think that by using daily Helmert parameters to align the individual solutions (cf. preprint line 305) the frame jitter is reduced at least to a certain extent. As you remarked, an internal regional frame alignment would further reduce the frame jitter. This will be an issue to be followed in a subsequent re-processing. Please also see our answer further below.

RC1.2:

Related to the previous point, the authors decided to align the position series to the IGb14 frame instead of the latest IGS20 frame. Since the objective is to provide a self-consistent and complete dataset of the highest quality in terms of series noise, the choice of the frame is mostly irrelevant for the correct combination of the position series. Mostly irrelevant, but maybe not totally, because the reported position discontinuities for many reference stations in Antarctica are not the same in these two reference frames (see for instance the changes in DUM1, MAW1 or DAV1). 

Of course, you are right to use the latest reference frame solution. However, at the time of the re-processing the IGS20 solution has not been available yet. In a next reprocessing we will use the IGS20/IGb20 (or whichever version will be the most recent one).

RC1.3:

Although this may have an impact on the terrestrial frame that was realized by the combined series, it might be possible that the use of a small set of stations has even a larger impact than the choice of one of these two reference frames. In this regard, a daily frame that is obtained from the combination of the four solutions themselves, without any external reference, could improve the regional frame alignment (and reduction of the frame noise) by including the maximum number of sites, including those that may not be included in the IGb14. Even the temporal sites could possibly contribute, though this may introduce a seasonal frame noise signature in the series with a frame systematically better in summer. In addition, a self-consistent combined frame would account for systematic coordinate variations that are not considered in the secular IGb14 frame, such as the periodic variations and the impact of ice inside some of the antennas of the network, among other effects. Once the frame noise is minimized and the series are combined, they could be translated to the IGb14 or the IGS20 using a selected subset of reference stations, as the authors did.

Thank you for this consideration. It is an issue we have discussed a lot in our group, with all pros and cons. We will definitely consider your suggested two-step approach for a future re-processing. We believe that the "internal" alignment of the individual solutions could be enhanced by considering an increased number of sites or even all possible sites, as you suggested. In any case, it is a comment which we regard very much valuable for consideration in a next reprocessing.

RC1.4:

Related to the previous point, it may be useful for the non-specialist reader adding a paragraph to inform that the frame of the combined dataset may not be free of secular biases and that caution must be taken when using the combined series (and their velocities) in an absolute sense. This is particularly relevant in Antarctica because, in addition to the difficulty of realizing a frame in a continent with sparse geodetic infrastructure, the stability of the origin of the global reference frame along the polar axis is affected by relatively large uncertainties. For instance, the translation drift between IGb14 and IGS20 is 0.2 mm/yr along the polar axis, i.e., a direct bias of the vertical velocities in Antarctica. The actual uncertainty is likely larger than this bias.

We will include a respective short paragraph to better explain this problem of the secular bias and stability of our reference frame solution. 

The uncertainty caused by a (secular) drift of the z-axis realization would add to the overall uncertainty of inferred coordinate change rates, esp. the local vertical component.

RC1.5:

On a similar point related to the noise in the position series, the NEWC group is the only one that applied non-tidal atmospheric loading corrections (probably at the observation level?).

The impact of NTAL in Antarctica is not negligible with detrended and deseasoned variability (RMS) in the range of 2-4 mm in the vertical component, mostly depending on the distance to the coast. It may be appropriate to apply this correction to the series from the other groups in order to have a more consistent dataset before the combination, especially since the weight of each solution was obtained from level of noise in each series.

Yes, of course, you are completely right. The NTAL correction was only applied in the NEWC processing, namely on the observational level (based on VMF1 atmospheric pressure grids), hence it might lead to less “noisy” daily coordinate solutions. 

We decided to address this issue already in the current data set in order to minimise the inconsistencies as far as possible, even if a re-computation was not explicitly required. Therefore, the NTAL effect will be added back (restored) to the NEWC solution before inferring the combined time series product. This will increase the consistency between the individual solutions. In the paper revision, we will add respective statements in Section 3.1 to clarify the procedure. Consequently, we will adjust respective figures (namely Figs. 3-5, 7, 9, 10) and some numbers (e. g. Tab. 1). However, we do not expect a qualitative change in the statements in the manuscript. 

Also, we decided to additionally correct the final data set of combined coordinate time series for NTAL. This decision was made since the NTAL correction will further decrease the noise in the combined time series (which better facilitates a subsequent inference of coordinate change rates), and the resulting coordinate time series will better suit the needs and expectations of the users who will use our product for subsequent analyses in solid Earth science (especially GIA).

Both adjustments will also require to update the final coordinate time series at PANGAEA, which can be done without problems (as informed by the responsible colleague of PANGAEA). The format of the dataset will not change.

RC1.6:

The TUD group is the only one that applied corrections for the orientation of the ground antenna patterns. Some of the IGS stations in Antarctica report antennas not oriented to north (DAV1, OHI3, maybe others). However, the antenna orientations were not considered to compute the IGb14 reference frame. Due to the small number of reference stations used for the alignment, I wonder if the inconsistency of the TUD antenna corrections could partly explain the different frame bias of this solution shown in Fig. 3. In other words, the level of agreement with the reference frame needs to be considered for each selected station used in the alignment.

Thank you for this comment which underlines the need to clarify this issue also in the text. In order to be consistent with the conventions used to generate the IGb14, the antenna orientations of the used IGS reference sites were set to zero. Thus, we do not think that a frame bias originates from this issue. We will add some sentences to discuss the details of the TUD processing.

RC1.7:

The authors use the Lomb-Scargle periodogram to validate the overall quality of the combined series. However, it is not clear if the standard normalized LS periodogram was used. The standard LS periodogram is normalized by the series variance, and thus it provides information about the relative distribution of power at different periods, which is very convenient for stacking series of different quality. However, for the same reason, there is not much one can say about the absolute noise levels or absolute peak amplitudes when comparing normalized periodograms from different datasets, unless they have the same overall variance. The later might be the case for the individual solutions if their noise levels are similar, but this may not be the case according to the previous comments (frame noise and NTAL). Note that not only the noise levels between globally and regionally aligned series can differ, but also the shape of their power spectra. The same applies between series corrected for NTAL or not. In any case, since the authors claim that the combined series have smaller variance, which would be expected, the normalized periodogram of the combined solution is not fully comparable to the others. As an example, for systematic periodic variations that have not been affected by the combination, like the fortnightly aliasing and the draconitic harmonics, a smaller variance in the series would make them look apparently stronger.

Thank you for discussing this issue. We understand that it is necessary to clarify what was done in this respect, since the ways of scaling the PSD or LS periodograms (LSP), respectively, also depends on the used software package. In order to infer a realistic noise level we would not prefer a normalized LSP which gives values between 0 and 1. We used the un-normalized LSP, where the PSD takes values similar to the FFT analysis, but differing by a factor of 2. We will include a short sentence for clarification.

RC1.8:

It would be interesting to add error bars to the values of the parameter rates shown in Fig. 3a and maybe discuss about the significance of the differences.

We estimated the formal errors from the time series of Helmert parameters. The formal errors (standard deviations) are very small and would barely be visible in the Figure 3a. Instead, we will give additional information on the derived uncertainties in the captions.

RC1.9:

Eq. 3: when computing the scaled variances, the length and period considered for each series must be similar among the solutions. This is because the series are not stationary (as shown in Fig. 7) with a noise variance that increases with the length considered. The noise level can also change over time, especially for series starting in the 90s.

Of course, you are correct. We want to emphasize the empirical character of this approach. The time series provided by the individual ACs will differ in temporal coverage and the overall time base, this in turn affects the residual RMS. We will edit the section and further emphasize on the limitations of our approach.

RC1.10:

Line 358-359: this sentence may not be clear for the reader. The scaling factor of the variances is only relevant for the weighted combination. The variances of the estimated combined values, in a station-by-station basis, should not be affected by optimistic or pessimistic variance factors.

Thank you for the comment. You are correct and we will revise the text accordingly.

We want to thank the reviewer very much for the constructive remarks and helpful comments. In the following, we will respond to the detailed comments, and how to consider them for a revision of the paper. Also, your comments will be considered for a future second reprocessing, which is already in the planning stage.

We have numbered the reviewer’s comments consecutively (as RC2.1 and so on), and repeated them in italic. Our replies are given below each comment.

Responses to the specific comments

RC2.1:

The aim of the paper was to process the GPS data consistently during the considered time span and obtain consistent combined position time series. However, the differences in processing strategies could be noticed between the individual AC solutions, which could potentially decrease their mutual consistency. For example, NEWC AC applied non-tidal atmospheric loading displacements for stations; these corrections were not applied by the other three ACs. TUD AC was the only one to apply atmospheric tidal loading displacements and correct the for the antenna misalignments towards true north (correcting for these misalignments is also not consistent with the IGS14/IGb14 framework used in this study). Regarding an elevation mask, TUD used 5 degrees while UTAS and NEWC used 10 degrees (for OSU the value is not provided). Also, two ACs used orbit products from the IGS repro3 project together with the igs14.atx antenna model and IGb14 terrestrial reference frame, which are not consistent with the antenna model and reference frame used for deriving the IGS repro3 products. Maybe the authors could add a comment on these discrepancies and their potential impact on the results. Nevertheless, for future reprocessing it would be important to harmonize the processing strategies among the participating ACs. It could be also reasonable if all the analysis were carried out in the latest and most precise IGS framework (IGS20/igs20.atx). It could also allow for a consistent extension of the product by next (operational) solutions.

You are absolutely right, there are differences in the processing choices of the individual ACs. In future reprocessing efforts, we will focus on a more harmonized strategy to further improve the consistency between the solutions and, in turn, the quality of the combined data product.

We regard the reduction of S1/S2 atmospheric tidal loading effects (as in the TUD solution) as small, since the amplitudes are well below 0.3 mm for the polar regions (according to the model of Ray and Ponte (2003) we applied). Likewise, the correction would be omitted in a future reprocessing to enhance the consistency between the individual solutions.

We will clarify the section regarding the correction of the antenna azimuth in the TUD solution. You and RC1 noted correctly, that this correction is not consistent with IGb14 conventions. Therefore, in the TUD solution the antenna azimuth correction was applied only for those sites, which have not been used to define the reference frame in the solution.

W.r.t. the NTAL inconsistency, we refer to the reply to RC1 (RC1.5) who also raised this very valid point.

The inconsistency of the NEWC solution with regard to the use of repro3 and igs14.atx/ the IGb14 will probably result in a small scale effect.

Of course, you are right to use the latest reference frame solution. However, at the time of the re-processing the IGS20 solution has not been available yet. In a next reprocessing we will use IGS20/IGb20 (or whichever version will be the most recent one).

Furthermore, remaining inconsistencies contribute to the overall "processing noise". It was the goal of this "exercise" also state where differences and inconsistencies exist. Of course, in a second reprocessing we would strive for an improved consistency. 

RC2.2:

Before a combination of AC solutions, each daily AC solution was transformed to the IGb14 to minimize the systematic differences between the solutions. The authors decided to use 6-parameter transformation, which should be sufficient to obtain consistent position time series. However, using a 7-parameter transformation could probably improve the consistency between AC solutions and also the repeatability of the position time series (especially for stations at the edges of the network). 

We did not include the scale factor as the scale parameter in the Helmert transform, since the reference to the IGb14 was done in a regional sense, and then scale factor and translation in z-direction are strongly correlated.

RC2.3:

Also, the authors did not write, if daily coordinates of IGb14 stations which were used as reference for the estimation of transformation parameters were verified against the catalogue coordinates. In case of too large differences, it would be important to exclude such stations from the estimation of transformation parameters to avoid biases in the alignment.

We verified the quality of the transformation by investigating the post-fit residuals for each site. If the residuals exceed the threshold of 15 mm, 15 mm or 30 mm for the north, east and up respectively, the site is not considered for the transformation. We will clarify our procedure in the manuscript.

RC2.4:

The authors combined the AC solutions according to the site-wise approach. The presented approach is appropriate for obtaining the combined station coordinate time series. I noticed, however, that the resulting combined coordinates were not transformed to the IGb14 (like individual solutions), so I suspect that these solutions may not be optimally aligned to the IGb14. 

This is completely correct. While the individual solutions were regionally aligned to the IGb14, the combined solution might not be perfectly aligned to the IGb14. Here we refer to a comment by RC1, who proposes a two-step approach. In a first step the individual solutions are “internally” alignment without an external reference, followed by then second step in which the solutions are then aligned to the IGb14 (or any other reference frame). We will consider this approach in a next re-processing. 

RC2.5:

In the supplement also three other variants of site-wise approach were described. The other possibility, however, could be a combination of daily AC solutions on the normal equation or covariance matrix level, together with the alignment of the resulting combined solutions to the IGS reference frame (e.g., using minimum constraints conditions). I wonder if the authors considered this approach for the creation of the combined position time series as well or see some disadvantages of it in relation to the site-wise approach.

Thank you very much for this comment. Yes, we also discussed the idea of using a normal equation stacking. One problem could be that the provided PPP solution (GIPSY) is based on a site-by-site processing and, therefore, does not provide covariance information (inter-site). Here some empirical approach would be needed. Another problem might arise in the weighting of the normal equations (NQs) from the different ACs. A big advantage of the NQs-stacking approach would be that we conserve the information on the inter-site covariances.

We will consider the NQ-stacking approach in a future re-processing and compare the results to the approach presented in this paper.

RC2.6:

Minor remarks

The TUD AC used Bernese GNSS Software version 5.2 to compute its solution. The authors wrote that the VMF3 was used for the troposphere modelling. However, as far as I know, in the official release of version 5.2 it is not possible to use VMF3; this capability is available in the next version (5.4). Please add a comment to clarify this.

Thank you for the comment. We used a technical modification to make the VMF3 product available in the Bernese Software v5.2. We will add a comment in the manuscript to make clear that this is not a standard feature of the Bernese Software v.5.2. In future, the Bernese Software v5.4 will be used where the VMF3-grids are supported by default.

In Fig. 3a the authors provided the values of rotation rates. It could be useful to add the information in the text on the used convention in Eq. 1 for the positive rotation of the coordinate system (clockwise or counterclockwise).

Will be clarified in the text as well as in the caption of Fig. 3. 

Technical corrections

In line 357 the closing bracket is missing.

Will be corrected.