On the performance of position-domain sidereal filter for 30-s kinematic GPS to mitigate multipath errors

The noise level of kinematic Global Positioning System (GPS) coordinates is much higher than static daily coordinates. Therefore, it needs to be improved to capture details of small sub-daily tectonic deformation. Multipath is one of the dominant error sources of kinematic GPS, which the sidereal filter can mitigate. With increasing interest in applying kinematic GPS to early postseismic deformation studies, we investigate the characteristics of multipath errors and the performance of the position-domain sidereal filter using 30-s kinematic coordinates with a length of nearly 5 days. Experiments using three very short baselines mostly free from atmospheric disturbances show that multipath signature in position-domain has better repeatability at longer periods, and sidereal filtering without low-pass filtering yields a lift of power spectral density (PSD) at periods shorter than 200 s. These results recommend an empirical practice of low-pass filtering to a sidereal filter. However, a moderate cut-off period maximizes the performance of the sidereal filter because of the smaller multipath signature at longer periods. The amplitude of post-sidereal-filtered fluctuation is less than 6 mm in standard deviation, which demonstrates the nearly lowest noise level of kinematic GPS used for postseismic and other tectonic deformation studies. Our sidereal filter is proven to mitigate several peaks of power spectral density at periods up to 100,000 s, but the period dependency of PSD is not fully alleviated by sidereal filtering, which needs future investigation.


Introduction
Global Navigation Satellite System (GNSS) observation is a powerful tool to study the surface deformation of the Earth. While studies on tectonics-related phenomena often employ daily static solutions of GNSS, raw observations, namely, carrier-phase of microwave transmitted from satellites, are usually recorded at a much shorter interval, for instance, 30, 15, or 1 s (e.g., Bock et al. 2000;Genrich and Bock 1992;Larson et al. 2003), or sometimes even shorter (e.g., Galetzka et al. 2015;Genrich and Bock 2006). Contrary to daily static solutions, kinematic analysis determines coordinates of antenna position at every observation epoch. Kinematic GNSS has been exploited to capture rapid ground motion during the passage of seismic waves (e.g., Bock et al. 2004;Jiang et al. 2021;Larson et al. 2003;Miyazaki et al. 2004) and postseismic deformation following great earthquakes (e.g., Jiang et al. 2021;Kato et al. 2016;Milliner et al. 2020;Miyazaki and Larson 2008;Morikami and Mitsui 2020;Munekane 2012;Twardzik et al. 2019). Daily static coordinates miss postseismic deformation less than 24 h after the mainshock as well as fast deformation rate change soon after the mainshock ("early" postseismic deformation), but kinematic GNSS coordinates can capture such early postseismic deformation in detail. Because postseismic deformation is a relaxation process of coseismic stress change, involving aseismic fault creep (afterslip) and viscoelastic relaxation of the upper mantle (e.g., Wang et al. 2012), this initial postseismic deformation contains critical information of the frictional characteristics of megathrust fault and the rheological characteristics of the upper mantle. Improvement of kinematic coordinates, hence, is crucial to gain more insights into postseismic deformation.
The nominal error of kinematic GNSS coordinates is usually on the order of centimeters or more (e.g., Bock et al. 2000;Jiang et al. 2021;Twardzik et al. 2019), considerably larger than daily static coordinates (~ a few mm; e.g., Jiang et al. 2021;Milliner et al. 2020). This larger error with kinematic GNSS coordinates needs mitigating to capture postseismic displacements, which are as small as a few millimeters at sites far from the mainshock. Many previous studies have been devoted to mitigating error sources in the GNSS time series, including atmospheric and ionospheric disturbances, multipath, satellite orbit, and Earth tides (e.g., Bock and Melgar 2016 and references therein). Among these error sources, multipath is the delay of carrierphase of microwave arrival caused by reflection at any objects around an antenna. Therefore, multipath is inherent and has been usually mitigated by datadriven approaches. Previous studies have demonstrated that the sidereal filter is useful to mitigate multipath errors in position coordinates by taking advantage of the repeatability of satellite constellation (e.g., Genrich and Bock 1992). As long as the multipath environment (i.e., arrangement of objects around each antenna) is kept unchanged over time, the geometrical relation responsible for multipath errors is also kept unchanged. Therefore, we can largely mitigate multipath errors by taking a difference of two periods apart by the repeat time of satellite constellation. A key parameter for sidereal filtering is the repeat time of satellite constellation, which has been of primary interest in previous studies. Different system of GNSS has a different repeat time. For example, the repeat time of Global Positioning System (GPS) satellites has been known to be nearly a sidereal day (i.e., 86,164 s = 1 day minus 236 s) (e.g., Bock et al. 2004;Genrich & Bock 1992;Nikolaidis et al. 2001). However, later studies have revealed that the repeat time is different among satellites and variable with time (e.g., Larson et al. 2007;Ragheb et al. 2007).
Various methods have been established for sidereal filtering, but they are primarily classified into two approaches, an observation-domain and a positiondomain filtering. The observation-domain filtering eliminates the multipath signal from raw carrier-phase observations before the positioning analysis (e.g., Atkins and Ziebart 2016; Iwabuchi et al. 2004;Ragheb et al. 2007). Two methods to construct the filter have been known in this approach. In one method, a sidereal filter is constructed from carrier-phase residuals of another day with a priori position from the static analysis and the constructed filter is subtracted from data during periods of interest with the appropriate time shift. This method allows for consideration of different orbit periods of each GNSS satellite (e.g., Atkins and Ziebart 2016; Ragheb et al. 2007;Wang et al. 2018), so it is applicable to mitigate multipath in both single-GNSS and multi-GNSS kinematic positioning analysis (e.g., Geng et al. 2017Geng et al. , 2018 even though different GNSS has largely different repeat periods of their satellites. The other method of sidereal filter construction in the observation-domain is the hemispherical mapping (e.g., Dong et al. 2016;Fuhrmann et al. 2015Fuhrmann et al. , 2021Iwabuchi et al. 2004;Moore et al. 2014;Zheng et al. 2019). Here, the sidereal filter is constructed by mapping the carrier-phase residuals on a topocentric hemisphere above each site as a function of the elevation angle and the azimuth of satellites. As this approach is free from various orbit repeat period of each satellite, it can be easily implemented to real time GNSS positioning (Dong et al. 2016;Fuhrmann et al. 2015) combined in the antenna PCV file (Moore et al. 2014). In either method, the features raised so far are essential advantages over the position-domain filtering approach in which we are allowed to assume only one representative repeat time Itoh and Aoki Earth, Planets and Space (2022) 74:23 of multipath signature (e.g., Bock et al. 2000;Dai et al. 2014;Genrich and Bock 1992;Ragheb et al. 2007). This is because, in the position-domain filtering approach, the filter is constructed from the position time series. Yet, previous studies show that the use of one representative repeat time is a reasonable approximation in their GPS data analysis Ragheb et al. 2007). The position-domain approach cannot handle multiple GNSS with largely different repeat time appropriately, so sidereal filtering for multi-GNSS positioning has been limited to the observation-domain approach (e.g., Geng et al. 2018;Zheng et al. 2019). In both observation-and position-domain approaches, the representative orbit or constellation repeat period of GPS is ~ 10 s shorter than a sidereal day, with its deviation typically within 10 s. This sub-sidereal-day repeat time yields better performance of sidereal filtering than using the repeating time of exactly one sidereal day (e.g., Agnew and Larson 2007;Choi et al. 2004;Larson et al. 2007;Ragheb et al. 2007).
Another practice in improving the performance of sidereal filtering is the removal of short-period fluctuations from a sidereal filter regardless of the observation-or the position-domain filtering. Here, previous studies empirically suggest that sidereal filtering does not work well to get rid of shorter period fluctuations and hence they are removed prior to the sidereal filter application (e.g., Atkins and Ziebart 2016; Dai et al. 2014;Geng et al. 2017Geng et al. , 2018Larson et al. 2007). Dong et al. (2016) pointed out that the hemispherical mapping method, a class of the observation-domain approaches, functionally removed the shorter period fluctuations because of the discretization of hemisphere. Most previous experiments, however, are more or less biased by other error sources other than multipath, such as atmospheric disturbances, because they have employed Precise Point Positioning (PPP) or double difference analysis with long baselines, in which situation they do not cancel (e.g., Atkins and Ziebart 2016; Bock et al. 2000). Only a few previous studies have conducted experiments with baselines short enough to cancel the atmospheric disturbances out (e.g., Bock et al. 2000;Dai et al. 2014;Dong et al. 2016;Geng et al. 2017), but they have not carried out systematic tests on the role of removal of the shorter period fluctuations in mitigating multipath errors. Wang et al. (2018) carried out a systematic experiment on the role of shorter period removal in sidereal filtering using a very short baseline (~ 12.5 m), but they employed only one baseline prepared for this purpose. Hence, there is still a need to investigate the behavior of multipath errors and the performance of sidereal filtering at different periods using multiple baselines in an environment free from atmospheric disturbances and with other primary error sources eliminated as much as we can.
Most previous studies on multipath and sidereal filtering have used 1-s kinematic coordinates by explicitly or implicitly assuming an application to seismic wave measurements (e.g., Atkins and Ziebart 2016; Geng et al. 2017Geng et al. , 2018Larson et al. 2007;Wang et al. 2018), and a few studies have used 30-s kinematic coordinates (e.g., Bock et al. 2000;Fuhrmann et al. 2015). As kinematic GNSS positioning can capture sub-daily postseismic deformation in detail, it is vital to investigate multipath errors and the optimum sidereal filtering for 30-s kinematic coordinates with a length of several days, which are yet to be fully understood. Investigation of multipath errors at longer periods than achieved by previous studies would also provide useful information on the use of kinematic GNSS positioning to monitor anomalous crustal activities such as earthquakes (e.g., Kawamoto et al. 2017;Melgar et al. 2020).
In this study, we investigate multipath errors and the performance of position-domain sidereal filtering to mitigate them from 30-s kinematic GPS coordinates. The new feature of this study arises from combination of (1) investigation of the multipath errors at up to the time scale of early postseismic deformation; (2) systematic experiments on removal of the relatively short-period fluctuations from the sidereal filter series, which have been empirically done, and (3) use of the environment mostly free from a substantial portion of atmospheric and ionospheric disturbances as well as orbit errors by constructing three very short baselines ranging from 3 to 65 m. The second and third features are more or less overlapping with Wang et al. (2018)'s focus, but our study covers the wider period range using the longer data length of nearly 5 days. Our study can be interpreted as a reexamination of their findings using the different baselines, which provides us an opportunity to further improve our understanding of noise characteristics of kinematic GPS.

30-s kinematic analysis
We employ differential positioning rather than the precise point positioning (PPP) to focus on the multipath effect as error sources. We construct three baselines using preexisting GNSS sites at various places in the world (Table 1) to explore common characteristics of the multipath and sidereal filter. The approximate lengths of the baselines used, POTM-POTS in Germany, GODE-GODN in the USA, and TSK2-TSKB in Japan, are 3 m, 65 m, and 36 m, respectively, which are short enough to cancel most atmospheric and ionospheric disturbances, the primary sources of positioning uncertainties, as well as satellite orbit errors (e.g., Bock et al. 2000). Throughout this study, codes of two sites in each baseline are connected with "-" in the order of the rover (kinematic) and reference sites (Table 1). These sites are located in governmental research institutes, so enigmatic changes of multipath environment in a short time, especially change of object arrangement around the antennas, are less anticipated than sites in public places such as parks and schools. Antennas and receivers of the two sites for each baseline are not the same in most cases (Table 1); this factor could be an additional positioning uncertainty source (e.g., Park et al. 2004). As different GNSS sites are commonly equipped with different antennas and receivers, we do not attempt to prepare new sites with the same equipment for this study. Yet, we attempt to mitigate effects inherent to different equipment types between two sites by implementing the phase center variation correction and the receiver correlation type information.
We employ the TRACK package (version 1.41) of the GAMIT/GLOBK program (Herring et al. 2015(Herring et al. , 2018a to carry out the kinematic analysis of the three baselines. We retrieve daily Receiver Independent Exchange Format (RINEX) files between 01 January 2019 and 31 December 2019. RINEX files on the day of year (DOY) 169 and 176 are not available at TSK2 and TSKB and the RINEX file at TSK2 is not available on DOY 339. The sampling interval of the RINEX files and hence obtained kinematic GNSS coordinates is 30 s, except for the 1-Hz analysis shown in sections "1-Hz kinematic analysis" and Noise characteristics at shorter periods than 200-500 s. We process a one-day session for each TRACK run. As our baselines are considered short enough to cancel the ionospheric delay out, we do not take the linear combination of L1 and L2 carrier-phase observations, but independently use the two observation types to fix the integer ambiguities (e.g., Schaffrin and Bock 1988) by employing the "short" mode of TRACK. We only use observations from GPS satellites above their elevation angle of 15 degrees. We always exclude PRN4 because the assigned space vehicle changed a few times during the analysis period. We do not use observables from other satellite systems (e.g., GLONASS) because they have quite different recurrence periods (e.g., Geng et al. 2018). We use the precise orbit product of the International GNSS Service (i.e., IGS final). As all the baselines are short, we do not estimate the atmospheric delay nor model tidal response. We implement differential code bias information and phase center models of satellite and stations (i.e., the ANTEX information) aligned to ITRF2014 (Altamimi et al. 2016). The radome information at TSKB and GODE happens to be unavailable in the ANTEX table, so only the antenna information is implemented for these two sites. A priori coordinates of these sites are constrained to daily static coordinates in ITRF2014 processed by Nevada Geodetic Laboratory, University of Nevada, Reno (Blewitt et al. 2018). When the daily coordinates are missing, we interpolate daily coordinates on their neighbor days. The uncertainty of the a priori coordinates is set to 10 cm. The process noise of kinematic coordinate estimates is set to 5mm/ √ 30s.

1-Hz kinematic analysis
We carry out 1-Hz kinematic analysis for GODE-GODN to examine the effect of sampling interval on the performance of sidereal filtering, as will be discussed in section "Noise characteristics at shorter periods than 200-500 s". The analysis setting is mostly same as the 30-s analysis except for the following points. First, the data length processed by one TRACK run is 6 h (c.f., 24 h in the 30-s analysis). Second, the number of 6-h-long data sets we process is only six to save the computation time (Additional file 2: Table S1) (c.f., all the days in 2019 in the 30-s analysis). Finally, the process noise of kinematic coordinate estimates is set to 0.91mm/ √ 1s , equivalent to that for the 30-s analysis ( 5mm/ √ 30s).

Postprocessing to remove outliers
In this study, we analyze the change in the relative position of the rover site from the reference site. Given the very short baseline lengths, it is reasonable to assume no changes of baseline lengths throughout the analysis period. Hence, any variation of obtained coordinates should predominantly represent effects of systematic errors inherent in the positioning technique such as the multipath. The obtained kinematic coordinates include epochs largely deviating from the average position and not repeating over  Fig. 1 and Additional file 1: Fig. S1); such epochs should be removed as outliers because they are not the phenomena of our interest. We design the following outlier removal procedure consisting of five steps: (1) Remove epochs estimated from a small number of double differences (DDs) and with large root mean square (RMS) of post-fit DD residuals. The minimum number of DDs for each epoch to retain is set to 8 for all the baselines, while the maximum value of allowable RMS of post-fit DD residuals is variable with baselines (Table 2; Additional file 1: Fig. S2).
(2) Remove epochs deviated from the average of coordinates of all the remained epochs by 4 times their standard deviation (Table 2). Each epoch is retained at this and next step only when coordinates of all the three (north, east, and up) components are below the threshold. (3) Remove epochs deviated from the average of coordinates within a specified length block by 4 or 5 times their standard deviation ( Table 2). The data length for each block is 430,700 s. This block length is a common least multiple of the sampling interval (30 s) and the representative repeat period of multipath (86,154 s; Ragheb et al. 2007). We explain the reason of this choice later in this section. This step removes epochs with local, relatively smaller deviations, which cannot be removed at Step 2. (4) Linear interpolation of the omitted epochs at Steps 1 to 3. (5) Split the year-long time series obtained at Step 4 into 430,700-s-long blocks.
Quality of observations and resultant positions are different among the baselines, so there is no reasonable common threshold at Stage 1 and 3; performance of a threshold designed for the noisier site would be obviously worse at the less noisy site. Hence, we determine the thresholds used at these steps by trial-and-error for each site not to remove too many epochs (Table 2). Eventually, 1.0, 1.4, and 4.2% of the total epochs are generated by interpolation for baselines POTM-POTS, GODE-GODN, and TSK2-TSKB, respectively. Although we do not include these interpolated epochs to assess the performance of the sidereal filter, we still need to interpolate removed epochs because low-pass filtering and power spectral density computation requires equally sampled time series. We present details of these analyses in section "Evaluation of performance of position-domain sidereal filter".
We obtain 73 430,770-s blocks and use them as a minimum unit to carry out all the examinations with 30-s kinematic coordinates. The assumed repeat time, 86,154 s (Ragheb et al. 2007), is slightly different in other studies and varies with time and each GPS satellite, but the deviation from 86,154 s is known much smaller than the sampling rate of 30 s (e.g., Agnew and Larson 2007; Choi et al. 2004;Larson et al. 2007;Ragheb et al. 2007). A repeat period of 86,160 s (Bock et al. 2000), the closest multiple of the sampling interval (i.e., 30 s) to 86,154 or 86,164 s, could also be the block length as a reasonable approximation of a repeat period. Yet, we need to prepare time series longer than the sidereal day to investigate multipath characteristics at the time scale of early postseismic deformation (e.g., Tsang et al. 2019;Twardzik et al. 2019). Then, the repeat time of 430,700 s is the best option because it naturally comprises the sampling interval, the representative repeat period, and the time scale of early postseismic deformation (i.e., ~ 5 days). This constellation repeat time is shorter than a sidereal day (~ 86,164 s) and our filter length is about 5 times of the standard sidereal filter, but we keep using the word "sidereal filter" to express our filter to remove multipath errors for simplicity, in contrast to previous studies which defined new terms (e.g., Choi et al. 2004).
Prior to Step 1, we concatenate the daily session of kinematic GPS coordinates obtained by the TRACK runs to construct the year-long coordinate series. Analyzing another 24-h session from 12:00:00 on 1 January 2019 to 11:59:30 on 2 January 2019 using a concatenated RINEX files (Additional file 1: Fig. S3) of POTM and POTS confirms that the concatenation in the positiondomain does not introduce artificial gaps or trends; the difference between the observation-domain and the position-domain concatenation is less than a few mm within ~ 30 min from the day boundary.

Repeatability of coordinate fluctuations
We first assess the repeatability of the coordinate fluctuations, which are fundamental in applying the sidereal filter, quantitatively by taking block-by-block correlation coefficients without the epochs filled by the interpolation (section Postprocessing to remove outliers). We define the correlation coefficient (CC) as follows: where d (i) k and d (i) are a coordinate at the kth epoch and an average of the coordinates of the ith block, respectively. N ij is the total number of epochs used for the computation of CC ij , the correlation coefficient between the , Itoh and Aoki Earth, Planets and Space (2022) 74:23 ith and jth blocks. As stated above, coordinates filled by the interpolation are excluded for those terms. Then, we also compute the correlation coefficients after low-pass filtering to examine if longer-period fluctuations have better temporal repeatability than shorter-period fluctuations. We first demean and detrend the time series and then apply the low-pass filtering by the lp command of Seismic Analysis Code (SAC) (Goldstein and Snoke 2005;Helffrich et al. 2013). The low-pass filtering is designed as a second-order Butterworth filter with various cut-off periods. We test the cut-off periods of 100, 200, 500, 1000, 2000, 5000, 10,000, and 20,000 s.

Sidereal filtering
Next, we apply our 430,770-s length sidereal filtering to evaluate its performance. We make pairs of two neighboring blocks and take a difference between the two blocks. We do not use pairs in which two blocks are not neighbors because it is practically common to generate  Fig. S2; see section "Postprocessing to remove outliers"). Blue dots indicate coordinates with 30-s intervals after the interpolation of neighboring epochs for outliers removed from the green dots with the whole-year and the local standard deviation criteria (see section "Postprocessing to remove outliers"). The same plots but for the east and up components are provided in Additional file 1: Fig. S1 a sidereal filter from coordinates of a period as close to the period of interest as possible. We compare the performance of unfiltered and low-pass filtered sidereal filters, which is evaluated by variance reduction (VR) defined as: where σ B,i is a standard deviation of coordinates of the block i before sidereal filtering, and σ A,i,j is that after sidereal filtering using the block j as a sidereal filter. We exclude the epochs filled by the interpolation (section "Postprocessing to remove outliers") when computing σ A,i,j and σ B,i . We obtain 144 VR values from 72 pairs of neighboring blocks because VR i,j is asymmetric with respect to i and j.

Power spectral density
A comparison of power spectral density (PSD) of kinematic coordinates before and after sidereal filtering is useful for assessing characteristic periods of multipath errors and whether the sidereal filter can mitigate them. First, we compute the Fourier spectrum of each block's coordinates by Fast Fourier Transform (FFT) using the fft command of SAC and then convert them to PSD by taking the squares of the amplitude of the Fourier spectrum. Prior to applying FFT to the coordinates, we demean, detrend, and taper the time series. We then apply a Hanning window to coordinates within 5% from each end to taper the time series. We obtain 72 PSDs from the coordinates before the sidereal filtering. The number of PSDs obtained from coordinates after the sidereal filtering depends on whether the low-pass filtering is applied to the sidereal filter; we derive 72 and 144 PSDs from the time series after the sidereal filtering without and with the low-pass filtering, respectively. We take an average of these PSDs when we discuss their characteristics; however, for averaging, we exclude PSDs derived from blocks containing 600 or more consecutive interpolated epochs because the linearly interpolated epochs would slightly lift the average PSD at long periods.

Repeatability of coordinate fluctuations
Kinematic GPS coordinates after the outlier removal show a notable periodic fluctuation pattern over time, but its period is shorter than the solar day seen as slant stripe patterns when the fluctuation amplitude is drawn by color (Fig. 2). When the plot is rearranged by the 430,770-s blocks, vertical stripes of fluctuation amplitude emerge (Fig. 3), strongly suggesting that these fluctuations are predominantly caused by multipath. It also suggests that the accumulation of the deviation of the "true" repeat time from the assumed representative value in this study (86,154 s) is not significant, at least, with our sampling interval (30 s). On the other hand, the stripe pattern is not uniform over time, exhibiting that the amplitude of multipath evolves over time. It can be interpreted as the temporal change of multipath environment. It is not always easy to identify their origins because it is highly case dependent, but they are possibly due to, for example, change in the number of visible satellites (e.g., Larson et al. 2007) and the recurrence position of satellites due to the maneuver (e.g., Choi et al. 2004) and/or the arrangement of reflective bodies around the antenna and their surface reflectivity (e.g., Elósegui et al. 1995). Typical nominal errors of the coordinates obtained in this study are 3-5 mm and 4-6 mm for the horizontal and vertical components, respectively. The standard deviation of fluctuation of each block is ranging 1.6-3.7 mm, 1.0-4.7 mm, and 2.4-7.6 mm for north, east, and up components, respectively. Correlation coefficients of 72 pairs of neighboring blocks are mostly between 0.5-1, 0.3-0.8, and 0.5-1, for North, East, and Up components, respectively, showing good repeatability in the 430,770-s window (Fig. 4). As a rule of thumb, the repeatability of GODE-GODN is better than POTM-POTS, which is still better than TSK2-TSKB. By extracting longer-period components with low-pass filtering, the correlation coefficient of these pairs is improved (Fig. 4), consistent with the idea of empirically applying low-pass filtering in constructing sidereal filter (e.g., Atkins and Ziebart 2016; Geng et al. 2017Geng et al. , 2018Larson et al. 2007). However, the correlation between two blocks becomes worse with cut-off periods longer than 20,000 s (Fig. 4). Low-pass filtering with a cut-off period of 20,000 s does not yield the maximum correlation coefficient for most pairs (Additional file 1: Fig. S4). Low-pass filtering with a cut-off period of 20,000 s yields a much smaller amplitude of fluctuation than that with a shorter cut-off period (Additional file 1: Figs. S5, S6, and S7). Hence, it is reasonable to  Left (a, d, and g), center (b, e, and h), and right (c, f, and i) columns indicate north, east, and up components, respectively. Neither low-pass filtering (LPF) nor sidereal filtering (SRF) is applied. Green color indicates epochs filled by interpolation of neighboring epochs during the data cleaning process Itoh and Aoki Earth, Planets and Space (2022) 74:23 conclude that random fluctuation has bigger effects than the periodic multipath signature at longer periods, which degrades the correlation coefficient values. This period dependency holds when we compute correlation coefficients of all the possible pairs of blocks (Additional file 1: Fig. S8). Up to a certain cut-off period of low-pass filtering, the correlation coefficient improves even with pairs of temporarily apart blocks. This suggests that static multipath environment (e.g., antenna height and nearby buildings) causes a relatively longer-period component of multipath signature. The correlation of temporarily farther blocks is worse at the same period band, and it is worse at relatively shorter periods (Additional file 1: Fig. S8). As already discussed above, the accumulation of deviation of the "true" repeat time from the assumed value should emerge in all the baselines because the orbit repeat time changes over time among all satellites, which are not consistent with our correlation coefficient diagrams (Additional file 1: Fig. S8). Temporal change of local multipath environment would  Left (a, d, and g), center (b, e, and h), and right (c, f, and i) columns indicate north, east, and up components, respectively. Neither low-pass filtering (LPF) nor sidereal filtering (SRF) is applied. Green color represents epochs filled by interpolation of neighboring epochs during the data cleaning process. The same plots but those after LPF are provided in Additional file 1: Figs. S5, S6, and S7 Itoh and Aoki Earth, Planets and Space (2022) 74:23 have larger impacts on multipath signature in the position-domain because the correlation coefficients of three baselines decrease with time differently (Additional file 1: Fig. S8).

Performance of position-domain sidereal filtering
Our sidereal filtering reduces the coordinate fluctuation efficiently when low-pass filtering with a cut-off period of 500 s is applied to the sidereal filter ( Fig. 5 and Additional file 2: Fig. S9). The standard deviation of fluctuation after the filtering is in a range of 1.0-2.6 mm, 0.8-4.6 mm (but only four examples exceed 3.5 mm), 1.5-5.8 mm for north, east, and up components, respectively. Even without low-pass filtering to the sidereal filter, the post-sidereal-filtered coordinates look similar (time series C and D in Fig. 6, Additional file 2: Figs. S10 and S11). The difference between these two cases is not clearly discernable, apparently suggesting that low-pass filtering to a sidereal filter is unnecessary. However, the difference between them is clearly found as the difference of variance reduction (Fig. 7). The most frequent VR range of the horizontal components improved by up to 20% when low-pass filtering with a cut-off period of 500 s is applied to the sidereal filter (Fig. 7a, b, d, e, g, h). The range of the most frequent VR of the vertical component did not improve, but it is clear that more pairs have larger VR values with the low-pass filtered sidereal filter (Fig. 7c, f, i).
On the other hand, applying a low-pass filtered sidereal filter with a cut-off period of 5000 or 20,000 s clearly leaves periodic fluctuations in post-sidereal-filtered coordinates (Fig. 6, Additional file 2: Figs. S10, S12, S13, S14, and S15). The reason for this observation is obvious; with a longer cut-off period, the amplitude of  Fig. S10), so we cannot mitigate moderately correlated multipath errors due to the removal of them. Accordingly, VR decreases to 30% or further worse in most cases when fluctuation at a period longer than 20,000 s is only used for a sidereal filter (Fig. 7). This competition finds that a cut-off period of 500 s yields the maximum VR in most cases for POTM-POTS and GODE-GODN (except for the east component) after testing the eight cut-off periods (Fig. 8a-d and f ). The maximum cut-off period with the maximum VR for TSK2-TSKB and east component of GODE-GODN has a broader distribution between 500 to 2000s ( Fig. 8e and g-i). These results imply that there is no magic number of the cut-off period of lowpass filtering for improving the performance of sidereal filtering, although 200 or 500 s might indicate the lowest limit of candidate cut-off period for the 30-s data.

Fig. 5
Coordinate fluctuation of three baselines as labeled in each block after sidereal filtering (SRF). The sidereal filter is made of each succeeding block. Low-pass filtering (LPF) with a cut-off period of 500 s is applied to coordinates used as sidereal filter. Left (a, d, and g), center (b, e, and h), and right (c, f, and i) columns indicate north, east, and up components, respectively. Green color represents epochs filled by interpolation of neighboring epochs during the data cleaning process. The same plots but using a different sidereal filter are provided in Additional file 2: Figs. S9, S11, S12, S13, S14, and S15 This conclusion is consistent with a simulation showing that the period of multipath signature depends on the site configuration and the satellite elevation angle (e.g., Larson et al. 2007). Trial-and-errors for the cut-off period with checking correlation coefficient and variance reduction would be important to improve the performance of the position-domain sidereal filter.

Power spectral density
PSD of coordinate fluctuations before sidereal filtering contains many peaks at periods of up to 100,000 s ( Fig. 9). Previous studies have not clearly identified some of these peaks of PSD at the longest period, especially longer than > 50,000 s (e.g., Bock et al. 2000;Geng et al. 2017Geng et al. , 2018Larson et al. 2007), likely because of their block length, quality, or both. Sidereal filtering mitigates most of these peaks (Fig. 9). On the other hand, interestingly, the PSD at the shortest periods (< ~ 200 s) somewhat increases by sidereal filtering without low-pass filtering in all the baselines tested (Fig. 9). The lift of PSD at the shortest period can be easily avoided by low-pass-filtered sidereal filter because it contains little fluctuation at these periods (Fig. 9). Not only the variance reduction, but also our PSD analysis also recommends the low-pass filtering practice to improve the position-domain sidereal filter. The qualitatively same conclusion has been drawn in Wang et al. (2018), which applied the sidereal filtering to 1-Hz GPS data in the observation-domain. We further discuss this phenomena at the shortest period in section "Noise characteristics at shorter periods than 200-500 s".

Discussion on the noise characteristics
Period (or frequency) dependency of PSD provides information of colored noise contaminated in the coordinate time series, which can bias extracted crustal deformation signature (e.g., Mao et al. 1999;Zhang et al. 1997). Most of our nine PSDs have a kink at 200-500 s (Fig. 9), showing different noise characteristics on both sides of it. Hence, in this section, we discuss the noise characteristics shorter and longer than this kink period (i.e., 200-500 s) separately to gain insights into their origins.
Noise characteristics at longer periods than 200-500 s Bock et al. (2000) have reported period dependency of PSDs after sidereal filtering in a 50-m baseline in California; at periods between 100 and 100,000 s, its spectral index (i.e., α in PSD ~ T α where T is period) is 0.3. The origin of the remained weak period dependency of PSD is yet to be understood, but it is much weaker than the colored noise caused by atmospheric disturbances, which have the spectral index of 1 (Bock et al. 2000) or 5/3 (Williams et al. 1998). On the longer-period side of our PSD kinks (i.e., > 200-500 s), the PSDs of POTM-POTS have similar weak period dependency (Fig. 9ac), but the PSDs of the other two baselines indicate stronger period dependency with gradually becoming weaker towards longer periods. Actually, all of them should contain little atmospheric disturbances. Thermal expansion of GNSS monuments, mainly controlled by configuration of monuments (e.g., material and structure) and local meteorological condition Fig. 6 An example of the north component of coordinates before and after sidereal filtering (SRF). a Coordinates of baseline POTM-POTS before (A) and after sidereal filtering by a sidereal filter without (C) and with low-pass filtering using a cut-off period of 500 (D), 5000 (E), and 20,000 (F) seconds. Sidereal filters used are drawn with the corresponding color for each low-pass filter setting at row (B). b Same as a but for GODE-GODN. c Same as a but for TSK2-TSKB. The same plots but for east and up components are provided in Additional file 2: Fig. S10 Itoh and Aoki Earth, Planets and Space (2022) 74:23 (e.g., temperature and insolation pattern), is known to impact the positioning analysis (e.g., Hatanaka et al. 2005;Munekane 2012;Yan et al. 2009). The local meteorological condition of two sites of each baseline should be the same in our study. The antennas of both POTM and POTS are situated on the brick pillars of the same architecture (Helmholtz Centre Potsdam GFZ German Research Centre for Geosciences 2020, while the other baselines are not (International GNSS Service 2021a, 2021b, 2021c, 2021d, so one might think that the configuration of monuments is a possible reason to make the remaining period dependency of PSD in GODE-GODN and TSK2-TSKB baselines. However, the pillar height is not high at GODE (0.5 m), GODN (1.5 m) and TSKB (1.5 m), so the un-differenced thermal expansion effect at each site is likely small because the thermal expansion is proportional to the pillar height (Yan et al. 2009). The pillar of TSK2 is covered by the concentric cover, which is known to mitigate the thermal expansion effect (Munekane 2012). Furthermore, Munekane (2012) has reported diurnal change in apparent position due to the thermal expansion, so it would only produce several peaks in corresponding PSDs. Hence, we conclude that the thermal expansion of monuments is not significantly responsible for the period dependency of PSDs found in this study. Another well-known GNSS error source is the common mode noise appearing uniformly at all the sites (e.g., Wdowinski et al. 1997). One of the origins is the instability of the reference site, but such effect is inseparable in this study because each experiment in this study uses only one reference and one kinematic site. However, if more baselines are simultaneously processed by the kinematic analysis with a common reference site, the monument effect of a reference site could be estimated by common mode filter (e.g., Wdowinski et al. 1997)  Histograms of variance reduction by sidereal filtering. a-c Variance reduction using sidereal filter without (green) and with (others) low-pass filtering (LPF) for north (a), east (b), and up (c) components of a baseline POTM-POTS. A cut-off period for LPF is 500 (red), 5000 (blue), and 20,000 (brown) seconds. d-f Same as a-c but for a baseline GODE-GODN. g-i Same as a-c but for a baseline TSK2-TSKB Itoh and Aoki Earth, Planets and Space (2022) 74:23 statistical decomposition such as principal component analysis or independent component analysis (e.g., Zheng et al. 2021). In this way, fluctuation of kinematic coordinates would be further alleviated (e.g., Larson et al. 2007). We should note that we cannot rule out the possibility that some local outliers remaining after the postprocessing ( Fig. 1 and Additional file 1: Fig. S1) as well as dynamic changes of multipath environment including ground reflectivity change impact the PSD shape. Future investigations would be necessary to pursue the origin of the period dependency of PSD at the long-period band.

Noise characteristics at shorter periods than 200-500 s
The shorter-period side of the PSD kink, where the PSD lift occurred, is close to the sampling interval (30 s). Hence, there is a need to examine the possibility that the sampling interval impacts the repeatability of the multipath signature at the shorter periods. We carry out a small experiment using 1-Hz observations with the GODE-GODN baseline. We obtain six 1-Hz kinematic GPS time series with a length of 21,600 s (i.e., 6 h; section "1-Hz kinematic analysis") and subsequently apply the sidereal filtering to the three pairs without low-pass filtering. for consistency with the 30-s data analysis (Additional file 2: Table S1). The sidereal filter again succeeds in mitigating coordinate fluctuations of all the components (Fig. 10a-f ), but the PSDs at a period shorter than 100 or 200 s lift after sidereal filtering (Fig. 10g-i). Hence, the sampling interval has negligible impact on the repeatability of multipath signature in the position-domain. Similarly, Wang et al. (2018), which employed the 1-Hz GPS observations for the similar experiment using the observation-domain sidereal-filtering approach, reported the similar PSD lift at the period shorter than 50 s. Hence, the multipath noise at the shortest period with a given sampling interval might be overwhelmed by other types of noise. As mentioned in Introduction, it is well known that each GPS satellite has a different orbital period with temporal variations (e.g., Agnew and Larson 2007;Choi et al. 2004;Larson et al. 2007;Ragheb et al. 2007). Larson et al. (2007) have already pointed out that frequent adjustment of the repeat period (i.e., the amount of temporal shift of sidereal filter to a target segment) improves the performance of sidereal filtering in terms of post-sidereal-filtering position residuals. Another approach to mitigate the noise at the shortest period is the observation-domain sidereal filtering which allows us to account for the various repeat periods of each satellite and its temporal change. Successful examples are Geng et al. (2018) and Wang et al. (2018), which mitigated the noise level down to a period of 16 s and 50 s, respectively. Hence, more precise dealing of the repeat time would provide a certain improvement of the short-period multipath removal. These examples have used 1-Hz GPS observations, so their idea would not be directly applicable to improving of the sidereal filter from kinematic GPS coordinates at a sampling rate of 30 s. The spectral index of PSD at the shortest period (i.e., 60-200 or 500 s) is ~ 2 in most of our PSDs derived from 30-s and 1-Hz data, but the period dependency continuously weakens to zero near 2 s (Figs. 9 and 10g-i). The similar change in the spectral index has also been reported in Genrich and Bock (2006) and Wang et al. (2018). This shape of PSD indicates a combination of random-walk and random noise, which might overwhelm multipath noise at this short-period band. Here, we employ a sidereal filter made by stacking multiple blocks to the GODE-GODN time series to examine whether the PSD lift is alleviated because the stacking mitigates random noise contribution. A sidereal filter cannot be made of observations after the mainshock in a practical application to postseismic deformation. Therefore, we stack two or five blocks preceding a block of interest without lowpass filtering and then subtract it from the time series in the block of interest. Average PSDs after sidereal filtering show that the PSD lift is more alleviated when more blocks are stacked (Fig. 11), demonstrating mitigation of random noise by stacking. PSD of longer periods than the kink slightly decreases as well (Fig. 11a and c). However, , and up (f) coordinate fluctuations (color) after sidereal filtering without low-pass filtering to the sidereal filter. g-i Average of power spectral density (PSD) of each block of 1-Hz kinematic coordinates. North (g), east (h), and up (i) PSD before (black) and after sidereal filtering without low-pass filtering to the sidereal filter (green). Vertical thick line indicates a period of 200 s compared to the one-block or two-block cases, the PSD happens to become worse at an intermediate band, 200 -1000 s, when five blocks are stacked ( Fig. 11d-f ). The random walk noise is not mitigated by the stacking approach, so it might contribute to the PSD lift. The random walk noise in geodetic data is often interpreted as the monument instability, but its effect is not dominant at this shortest period (e.g., Langbein and Johnson 1997). Hence, the origin of random walk noise causing this PSD lift is still elusive. A practical prescription for postseismic deformation studies using kinematic GPS coordinates with a sampling of 30 s would be application of low-pass filtering or smoothing to both sidereal filter and data of interest. It is because the shortest periods are usually not a focus of early postseismic deformation studies (e.g., Milliner et al. 2020;Tsang et al. 2019).

Discussion on applicability to the practical postseismic deformation study
Practically, very short baselines like those employed in this study are not able to capture postseismic deformation signature because the motion due to the postseismic deformation at two sites constructing each baseline are uniform. To capture the postseismic deformation, coordinates derived from the differential analysis with a much longer baseline or PPP analysis are necessary. Such setting is not comparable to our experiments using the very short baselines because the atmospheric delay more or less affects the positioning result, even though its effect is accounted for by implementing the mapping function (e.g., Boehm et al. 2006). If the atmospheric delay exhibits the sidereal periodicity, the sidereal filtering technique can mitigate the atmospheric disturbances together with the multipath, but that is not the case in most situations. Nevertheless, with careful inspection of the temporal repeatability of fluctuation in the data, the sidereal filtering operation with low-pass filtering tested in this study can still be effective to reduce the multipath effect in the practical application in postseismic deformation studies. From this viewpoint, our experiments presented in this paper demonstrate the currently achievable lowest limit of the noise level of 30-s kinematic GPS coordinates used for the postseismic deformation as well as other tectonic deformation studies. Our experiments would be useful when we design early postseismic deformation studies in the future.