Accuracy Analysis of Point Velocities Determined by Different Software Packages and GNSS Measurement Processing Methods

: The paper analyzes the geodynamic network of the City of Zagreb for periodic campaigns carried out from 2006 to 2009 which were processed by different software packages. The first computations and processing results were obtained by using the scientific software Gamit/Globk and indicate the ongoing tectonic activity of the area. In this paper, all calculations were performed by using the scientific software Bernese. Processing strategies and applied error models in the Gamit and Bernese solution are analyzed. Th e results of the previous analyses show the need to perform GNSS measurements at intervals of up to one year, which is necessary f or understanding the mechanism of the structural frame of the wider Zagreb area. The research and analysis performed in this paper indicate certain uncertainties in determining the velocities from periodic one -year GNSS measurements. When periodic GNSS observations are analyzed at time intervals shorter than 2.5 years, annual signals can cause significant errors in determining p oint velocities. The accuracy of determining velocities between annual time intervals depends on a number of factors: noise level in GNSS measurements, measurement sessions quantity and applied computation strategies. Previously , the t ime series analysis of observations was a key procedure in the context of geodynamic and geokinematic research and the FODITS algorithm was used for the analysis. A noise analysis on the daily time series of coordinates was performed for the purpose of understanding all influences on geodynamic points. Moreover, a correlation between the time series of observations was determined in order to estimate the final velocity uncertainty error. The purpose of this paper is to demonstrate the applicability of the methods and procedures used to determine the coordinates and velocities of points that can be reliably used for ge odynamic and geokinematic analys es and consequently, timely responses to various geophysical phenomena due to earthquakes


INTRODUCTION
Previous calculations for establishing the reference frame and investigating geodynamics have made it possible to determine the character and parameters of the displacement of the Earth's crust in the territory of the Republic of Croatia [1][2][3][4][5]. However, measurement techniques such as the GNSS method are also a powerful tool for monitoring displacements and deformations in the local areas of interest, respectively for geodynamic and geokinematic analyses, and increasingly in the study of seismic activities. Because permanent stations are not densely distributed in the areas of interest, non-permanent, periodical observations on GNSS points are performed. An example of such a network that is analyzed in this paper is the geodynamic network of the City of Zagreb. The results of data processing indicate a constant presence of tectonic activity in the wider area of the City of Zagreb. The analysis of the processing results of the geodynamic network of the City of Zagreb [2] over a period of 20 years shows that, in the case of multi-year time intervals between GNSS campaigns, there is an averaging of the shifts or losses of certain information on the annual shifts, which are extremely valuable data needed for a reliable analysis. The results of previous analyses show the need to perform measurements at intervals of up to one year, which is necessary for analyzing the quality of the obtained results from the GNSS campaign and understanding the mechanism of the structural frame of the wider Zagreb area. The investigations carried out in this paper relate to examining the impact of the chosen error modeling method and the GNSS measurement processing strategy, while analyzing the time series of observations and the influence of time intervals between the observations on the accuracy and the precision of the points' coordinates and velocities. The data processing of geodynamic campaigns is influenced by specific deviations due to the use of different processing parameters (antenna and receiver parameters, software changes), but also due to random errors (receiver noise and unmodified atmospheric effects). Coordinate and velocity variations may occur considering the network adjustment strategy used to detect local geodynamic impacts. Such solutions are ambiguous and dependent on the strategy employed [6]. Based on the collected observations of the Geodynamic Network of the City of Zagreb (GNCZ), the application of various error modeling methods and processing strategies for GNSS measurements has been analyzed, with the possibility of eliminating their systematic or random influences, which will ultimately enable their reliable use for geodynamic and geokinematic analyses.
The paper analyzes the time series of coordinates and velocities for discontinuities and jumps, which is today a key procedure for geodynamic and geokinematic analyses. A large number of discontinuities is caused by model updates (e.g. antenna calibration models) and changes in the GNSS data processing strategies or new ITRF implementations. Additionally, there are other sources that cause them, such as earthquakes that can cause discontinuities and velocity changes in the points near the epicenter. Discontinuities may also be present for some unknown reasons (construction of new objects around the antenna, dropping accumulated snow on the antenna, etc.), and they usually remain for longer periods (weeks) and, if identified, can be determined with accuracy below mm, and they can be applied to the time series of coordinates. The jumps are present in the time series of coordinates mainly due to systematic errors, e.g. due to bad atmospheric conditions (snow storms), anthropogenic influences, etc. Finding jumps depends on the assumptions determined by the statistical test and the level found to be significant for the analysis of the datasets [7]. Jumps occur in a shorter period of time and the values are usually greater than those of discontinuities. By analyzing the daily solutions of the coordinates' time series, significant changes in the data in the predefined epochs are checked, i.e. discontinuities and jumps are determined. Multi-day coordinate solutions show fewer RMS values (Root Mean Square) than the daily ones, but they also use reduced data, hence it is important to use the daily solutions to determine discontinuities and jumps. A data analysis without averaging the time series of coordinates provides insight into the displacements and deformations of the location of interest, and the best insight can be obtained from continuous data analyses at permanent stations. The detection of geo-kinematic phenomena requires a data analysis of multiple daily GNSS solutions since determining displacement parameters requires a previous statistical analysis of the daily solutions that provide control for jumps. Daily solutions that have jumps, and remain in the processed dataset, can deform the established parameters, e.g. seismic deformations affect the time series of coordinates, accumulate over time and can significantly affect coordinates, not only from large earthquakes, but also from the accumulation of many small earthquakes. Seismic deformation modeling helps detect discontinuities in the time series of GNSS coordinates, which is one of the main sources of error in determining the ITRF frame today [8]. Numerous discontinuities induced by earthquakes are too small for visual detection due to seasonal variations and the GNSS noise that disable their identification. However, if not taken into account, they have a great influence on the determination of point velocities, considering the precision requirements for geodynamic and geokinematic analyses.

GNSS MEASUREMENT NOISE ANALYSIS
The noise analysis of GNSS measurements mainly assumes the presence of white noise in the time series of coordinates. Several previous studies [9] have shown that if correlated noise is neglected and if only the white noise model is used, then the uncertainty rate is significantly underestimated. Considering the combination of noise affecting GNSS measurements, it is difficult to accurately determine the source of the noise and, consequently, the reliability of the estimated time series of coordinates. Understanding noise content is very important in order to identify realistic uncertainties for the parameters to be determined. When applying the least squares adjustment, daily position errors are considered to be independent. However, using this method to calculate point velocities from time series coordinates gives some uncertainty. Studies have shown that the assumption that measurement errors are random and uncorrelated from an epoch-to-epoch (white noise) is unrealistic for GNSS measurements. Sources of time-correlated noise (colored noise) in the GNSS data include: orbits, atmospheric effects, stabilization shifts, etc. [10,11]. For computing realistic errors, the approach of [11] is selected, who developed an empirical model to determine the GNSS error (σr) for individual velocity components (N, E, H) using the time series of coordinates with the presence of white and colored noise. The final velocity uncertainty error is represented by the formula (1):   2  2  2  2  3  2   12   f  w  rw  r  b   a  a  gT  gT  g T   σ  σ  σ  σ ≅   +  +  (1) where: g is the number of measurements per year, T is the total time range of observations (T = period (in years)), and a and b are empirical constants: a = 1.78; b = 0.22. White noise is a random noise (Sw "white"). There are two main types of colored noise: pink (Sf "flicker") inversely proportional to frequency and red (Srw "random walk") inversely proportional to the square of frequency. In Fig. 1, the pink line represents the mean of white and pink noise, the blue dashed line is the maximum amplitude of white, pink and red noise [10], whereas the red dashed line is the minimum amplitude of the white and pink noise calculated from the time series of coordinates [12]. The magnitudes for white and pink noise are expressed in mm and for red noise in mm/yr.

Figure 1
Analysis of measurement noise using the example of the GNSS network in the western Alps [12] By investigating and analyzing the time series of coordinates according to [11], it was concluded that white and pink noise dominate the GNSS noise spectrum for the time series coordinates, while red noise is relatively small. Pink noise is spatially correlated and has a clear dependence on geodetic latitude. Although the pink noise amplitude decreases with time, it is still the most dominant. The vertical component is three times larger than the horizontal component. The red noise component is very difficult to detect since it requires a much longer time span of observations. Stabilization noise is characterized as red noise [10] and occurs with respect to the cycles of rain, freezing and atmospheric effects on rocks and soil. If the red noise is very small, the velocity error depends strongly on the time range of the observations (T) and less on the number of the observations (g). However, several studies have demonstrated the importance of red noise recognition in the GNSS data. According to [9] cumulative changes in the earth and atmosphere, stabilization shifts with respect to the deeper layers of the crust. The most common types of stabilization in geodynamic networks are: metal pillars and deep concrete pillars. According to [10], there is no clear distinction between the deep concrete pillars and ordinary metal pillars related to the stabilization noise. Moreover, according to [10], larger noises do not occur from stabilization instability. Another study showed that deep concrete pillar stabilizations have a lower colored noise, but also slightly smaller velocity residuals [13] (Fig. 2). to the type of stabilization [13] According to [11], a strong linear correlation between the weighted RMS values of the GNSS coordinate time series and the determined noise amplitudes for white and pink noise was observed with the linear equations shown in Tab. 1. Table 1 The linear correlation of RMS and measurement noise [11] Parameter  Theoretical velocity error from annual sinusoidal signals versus the time period [14] Although previous analyses have shown that white and pink noise are dominant in the noise model for the coordinate time series, there are also various other sources of error, such as periodic (seasonal) signals that can be imposed as dominant, such as localized deformations with respect to changes in the groundwater (unknown colored noise). Therefore, the parameter of the annual and semi-annual sinusoidal signal must be introduced into the process of determining point velocities, which must be determined simultaneously with velocities.
Seasonal signals contribute to the pole and ocean tide, which are usually corrected in GNSS processing by using the available models (the nutation model and sub-day pole model, ocean tide model and coefficients). Additionally, the earth's potential and solid earth tide are modeled. Errors caused by seasonal signals, which are not usually corrected in GNSS processing, and which must be taken into account are: atmospheric pressure, snow and soil moisture and the non-tidal ocean mass. It should also be emphasized that the accuracy and precision of coordinate and velocity determinations are also influenced by different computation or modeling strategies using different software packages for GNSS processing. Only with a complete understanding of the full spectrum of GNSS errors [15] (Tab. 2) could we obtain accurate and precise results that can be used in geodynamics and geokinematics. The Western Alps are one of the most researched areas related to monitoring the mechanism of the structural frame and detailed geophysical surveys [12], [16]. The GPS ALPES group, made up of French, Swiss and Italian researchers, monitors the GPS network in the Western Alps. In the period from 1993 to 1998, more than 60 stations were observed with 3-day sessions at each station, and the data were processed by different software packages. The calculated solutions with the GAMIT/GLOBK software had a horizontal repeatability (N-E) of 4-7 mm in the year 1993 and 2-3 mm in the year 1998. A comparison of coordinates from the years 1993 and 1998 shows that the residual of GPS points velocities is less than 2 mm/year. Solutions were compared with the BERNESE software and the differences for the 1998 dataset are quite small, reaching an average of about 2 mm for the horizontal components and 5 to 10 mm for heights. The differences are slightly larger for the 1993 dataset with 4-7 mm for the horizontal components and 10-20 mm for heights. At several points observed in the year 1993, slightly higher values (greater than 10 mm) were obtained on the eastern component and up to several tens of mm in height. This can be attributed to differences in the ambiguity resolution strategy of the two software solutions and it indicates the locations where shifts need to be carefully analyzed. However, for most sites, the differences generally remain within repeatability [16]. According to [17] GNSS, processing software packages play an important role in geodetic and geophysical studies. Error estimates or standard deviations do not have to be realistic values given their daily RMS repeatability. They present standard deviations as formal errors (FEs) and scale them according to RMS errors. The scale factor (SF) is defined as the ratio between the RMS error and the formal error.

PROCESSING AND ANALYSIS OF THE GEODYNAMIC NETWORK OF THE CITY OF ZAGREB
The aim of this study was to examine the application of appropriate error models and processing strategies for GNSS measurements by testing different software and processing methods that will allow reliably determined point coordinates and velocities to be used in geodynamic and geokinematic analyses. The impact of using a different reference frame in the adjustment (local vs. global network) was also examined. The GNSS measurements on the Geodynamic Network of the City of Zagreb (GNCZ) have been carried out for more than 20 years and they confirm the constant presence of tectonic activity in the area around Medvednica. The stabilization of geodynamic points was performed at most points in the form of deeply stabilized concrete pillars. All GNSS measurements campaigns of GNCZ [2] were processed in the same way using the scientific GNSS software GAMIT/GLOBK ver. 10.34 i 10.6 [19], which uses a Kalman filter to determine velocities from time-separated coordinates. When processing GNSS measurements of GNCZ by using the GAMIT / GLOBK software, all parameters were set for the local network, and the two most stable points in the study area were selected as the reference: the city's permanent CAOP (metal pillar on the building) GPS station and ZZF (metal pillar on the building) GPS station. The offsets of all other points in the network refer to the vector between the mentioned two points. The points of the baseline also had the possibility of a relative shift, but according to [18], they proved to be very stable. In this paper, the Bernese scientific software ver. 5.2 [20] was used for data processing and analysis of GNCZ for periodic campaigns from 2006 to 2009, using the least squares adjustment method. The results of the GNSS processing with the Bernese software were compared with the results of the processing obtained by the GAMIT / GLOBK software used in the realized campaigns. Using the prescribed standards for processing and analyzing GNSS data (IGS, EUREF), an investigation of the impact of various parameters on final processing results was performed. When processing and analyzing GNCZ with the Bernese software, IGS stations were selected as reference points: GRAZ (active since 1993), MATE (active since 1992) and MEDI (active since 1995). All three points are included in the IGS08 reference frame used in the adjustment. It was important to choose permanent stations that have a continuum of observations given the long period of the observation of the GNCZ. These stations are close to the area of interest (Eurasian plate), but far enough from the regional / local structures and they are not located in tectonically very active areas (to reduce the possibility of reference points moving). Additionally, these stations are under the constant control of IGS (International GNSS service). Connecting to the ITRF frame can be done as a global or local/regional connection. The use of global IGS stations as a reference frame is considered to be the most accurate method, while the local/regional connection (e.g. with CROPOS stations) is achieved by an indirect connection with the ITRF frame. Each software package (Gamit and Bernese) was used with its recommended settings, i.e. with the implemented independent GNSS data processing strategies. At the same time, important processing parameters and common standard procedures prescribed by the IGS and EUREF were used, such as: applying the PPP (Precision Point Positioning) and DD (Double Difference) calculation methods, solving residuals with the outlier ejecting, using global tropospheric and ionospheric models, determination of weights in elevation-dependent observations, ocean correction model usage (FES2004), while the atmospheric correction model was not available. The IGS final orbits and corresponding Earth orientation parameters were used. Other parameters taken into account relate to the limitations and specificities of each software and the applied processing strategies (ambiguity solution method, adjustment method, application of different reference points), which can influence the final solution and must be taken into account (Tab. 3). One significant difference in the adjustment method is that the GLOBK uses a Kalman filter (equivalent to sequential least squares if there are no stochastic parameters) that calculates with covariance matrices rather than normal equations (with which Bernese calculates), and it therefore requires an a priori limit for each estimated parameter. GAMIT produces a solution, which then uses GLOBK to calculate velocities with "weak" constraints on all parameters [19]. The Gamit solution aims to apply "weak" constraints to the entire reference frame: EOP (Earth Orientation Parameters), SV (Space Vehicle) and coordinate parameters, while the realization of the reference frame is implemented through seven Helmert transformation parameters. Bernese, on the other hand, uses the normal equations expanded with the point coordinate parameters. In the Bernese solution, the equation with minimum constraints is performed. The parameters are fixed: EOP and SV to IGS standard products and apply "weak" constraints to coordinates only. The reference frame is implemented through four parameters of the Helmert transformation (without translation conditions).   The differences in the solutions obtained for the two-year period from 2006 to 2008 were also compared, and it can be seen that the differences between the solutions of the two software programs are significantly smaller, and an absolute decrease of the displacement amplitudes for both solutions can be observed (see Fig. 5).
It can be seen from the graph in An additional analysis was performed by using a network of CROPOS stations (CROatian POsitioning System). As CROPOS stations were not available for the entire observation period of the periodic campaigns from 2006 to 2009, as the system started with the operation in December 2008, the use of the CROPOS reference points for comparisons and analyses was performed only for the year 2009. Although geophysical impacts are also present in global solutions and usually have more noise in observations than the local or regional network solutions [21], the regional and local impacts such as atmospheric or hydrological effects are usually not a priori modeled in global analyzes, but can be modeled for a regional or local station network if meteorological and hydrological data are available. Annual amplitudes of atmospheric influences can reach 4 mm for the radial component and are usually lesser than 0.5 mm for the horizontal component [15]. Global seasonal cycles in temperature and humidity depend on latitude, and the amplitudes of the north-south gradients are 40% larger than

BERNESE vH
GAMIT vH the amplitudes of the east-west gradients. Longitudedependent variations resulting from the local and regional effects show smaller mean amplitudes of east-west gradients compared to north-south gradients. According to [23], the amplitudes of annual and semi-annual signals in tropospheric gradients are generally smaller for the east-west component (mean absolute value of 0.17 mm) versus the north-south gradients (mean absolute value of 0.30 mm). However, the impact on heights is much greater than that on horizontal coordinates and the height is up to several tens of mm [23]. The analysis showed that the impact of the choice of reference stations (global/regional/local) on the calculated point coordinates was significant, mostly for the height component, where the differences of 2 cm were obtained between the global and local system (Fig. 6), while for the horizontal components (E and N), the differences were up to 1.5 cm. The differences between the IGS and CROPOS reference stations in the height component are about 5 mm. Therefore, it can be concluded that the selection of reference stations (global vs local) significantly affects the height component, while positionally, this influence is smaller. Figure 6 Global vs. regional vs. local reference stations' network differences

Noise Analysis of the Points of the Geodynamic Network of the City of Zagreb
A noise analysis was carried out for the campaigns from 2006 to 2009, and according to the formulas of [11] respectively, the values of σ w and σ f for all components were   campaigns. Moreover, this point shows higher amplitudes of displacement in that period. Therefore, the amplitudes of displacements at the points for periodic campaigns should be interpreted with attention, because higher measurement noise also affects the higher final error of velocity uncertainty.

CROPOS STATION COORDINATE TIME SERIES ANALYSIS
The GNSS permanent station data have been widely used for positioning and navigation, as well as for studying geodynamics. In the last few years, the GNSS data have been used extensively for the exploration and analysis of geodynamic and geokinematic shifts. In determining the velocities of points from GNSS measurements, the existence of continuous and long-lasting GNSS observations is of great importance for quality analyses and the evaluation of the accuracy of calculated velocities. The analysis of the observations of these points enables, in addition to the detection of the shift trend, the detection of periodic phenomena as well as the phenomena that have nonlinear structures in the time series of coordinates. The analysis of changes in the coordinate time series is being increasingly used as an accurate and reliable source of information on local geodynamics. An example is the permanent IGS and EPN networks used for the global, regional and local geodynamics research. The GNSS networks of various institutions have been established for the research of local geodynamics from the analysis of the time series of coordinates: Geospatial Information Authority of Japan -GSI, Southern California Integrated GPS Network -SCIGN, Jet Propulsion Laboratory -JPL, Massachusetts Institute of Technology -MIT, Scripps Orbit and Permanent Array Center -SOPAC etc. [22]. With the establishment of the CROPOS station network at the end of 2008, the state network of the permanent stations of the Republic of Croatia was realized; while for the area of the City of Zagreb, such a network of permanent stations was not feasible for economic reasons. Therefore, periodic measurements (24 h observations with two to three sessions) on the points of the City of Zagreb network were used for geodynamic research. The closest permanent station in the surveyed wider area of the City of Zagreb was the CAOP station, which is no longer in operation as a permanent station, and since the establishment of CROPOS, the closest stations have been ZAGR at the Faculty of Geodesy and ZABO at the cadastral office building located on the north side of Medvednica. It should be said that the CROPOS stations are stabilized on the roofs of buildings on fixed metal pillars. For the ZAGR and ZABO stations, daily sessions were calculated on a monthly basis for 2018 and for a 10-year period (2009-2018), where an observation period was selected in the campaigns of the geodynamic network of the City of Zagreb (June-July) to eliminate the seasonal impact. The same strategy and reference network (IGS stations) were used for the processing and adjustment as for the processing of the Geodynamic Network of the City of Zagreb. The adjustment statistics for the CROPOS stations ZAGR and ZABO in the processing of annual sessions show that the RMS value is between 1.12 -1.13 mm, while the "Chi-square" test for all years is between 1.39 -1.79.
Previous research and analyses on the permanent GNSS stations that are established worldwide show annual variations with maximum amplitudes of ±1 cm for the horizontal components and ±2 cm for the vertical components [24]. The analyses performed in the daily sessions over a 10-year period (2009-2018) show that the two CROPOS stations ZAGR and ZABO are very stable. An analysis along the horizontal axis N shows a linear trend of increasing in the north direction by 2 cm for ZABO and 4 cm for ZAGR in a 10-year period, and there are no major jumps at points, hence we can conclude that the points are stable in terms of jumps and discontinuities (Fig. 8). An analysis along the horizontal axis E shows a linear trend of increasing in the east direction by 1 cm, which means that the points move in the northeast direction without significant jumps and discontinuities at the ZAGR and ZABO locations in a 10-year period ( Fig. 9). At ZAGR, a slight deviation was observed in 2012, which corresponds to the fact that in that year, September was chosen for observations, which is slightly later than for other sessions (atmospheric conditions are different), which were observed in June or July, which could be the reason for the differences of 4 mm.
The height axis -H shows displacements of up to 2 cm in the 10-year period, where the linear trend is a slight increase in height for the ZAGR and ZABO points (Fig. 10). For the ZAGR point, changes in the direction (jumps) were seen in the years 2011, 2012 and 2016, while variations for the ZABO point were less expressed in the same time period.
If we look at the one-year period in which the monthly data for the two stations ZAGR and ZABO for 2018 were processed, and which may contain seasonal variations during the year, we can see a slight increase along the N axis during one year and about 5 mm variations on both stations. A slight increase is also visible along the E axis with a variation of 7 mm and a jump (5 mm) in May with a similar trend on both stations. The H axis shows a linear trend of subsidence for both the ZAGR and ZABO stations by 1 cm, where a jump in April (7 mm) is visible for the ZAGR point, while monthly coordinate variations (changes of direction by 3 mm) are visible for the ZABO point. It should be emphasized here that the ZABO point is located on a high metal pillar which can cause these small variations, as well as seasonal conditions (atmosphere). The spatial distribution of the parameters of periodic signals and their correlation with changes in the atmospheric factors and other geophysical factors is necessary to consider the whole spectrum of errors. Atmospheric data were not determined at the ZAGR and ZABO stations, nor were any other geophysical influences recorded, however, the analysis of the time series of coordinates was performed, and the FODITS algorithm of the Bernese software was used for the analysis. The analysis of time series of coordinates enables the detection of periodic signals (amplitude and phase of the signal), as well as other elements that can affect the coordinates of points. For the ZAGR and ZABO stations, the amplitude and phase of the periodic function were determined, as well as the basic components that were added to the functional model, and which correspond to the coordinates of the points and their linear velocities. Using the FODITS algorithm, new elements were determined in order to improve the functional model representing the time series of coordinates. For the ZAGR station, new elements of discontinuity were recorded for the observations in 2011, 2012, 2014 and 2016, which are visible primarily at the height component (Fig. 10). New elements of discontinuity for 2011, 2014 and 2016 have also been determined for the ZABO station, while for 2014, a new element for linear velocity was also determined. No significant deviations were recorded for any component, and no shifts were recorded due to earthquakes, equipment changes or other events. With the obtained new elements of the functional model, the coordinates and velocities in the adjustment process are recalculated. Periodic variations of the coordinates of up to 1 cm (in a 10-year period) were recorded for the height component, and according to [24], such variations are standard values for stable permanent networks in Europe. No correlation was noticed between the two stations ZAGR and ZABO in the time series of coordinates and related to the spatial component (distance of about 10 km), except for the already known fact that the points move in the northeast direction. In the height component, a descent was recorded for the two stations ZAGR and ZABO for the year 2018, however, for a ten-year period, a slight ascent was recorded for those two stations.

CROPOS Station Velocity Determination Analysis
An analysis of the displacements on the geodynamic network of the City of Zagreb from different epochs showed that it is necessary to carry out GNSS campaigns at the intervals of one year, since in the case of long-term intervals, an averaging of displacements occurs. However, the accuracy of determining the velocities from one-year periods depends on a number of factors: the intensity of the measurement noise, measurement sessions quantity and the applied computation strategies. Point velocities for periodic GNSS measurements correspond to the estimated linear rate of the daily solutions of two or more measurements. Because periodic measurements rely on two or more measurements with several days of data, periodic signals and jumps cannot be determined with certainty. These parameters represent the uncertainty of determining velocities from periodic measurements. Most campaigns are observed over the same time period to minimize the impact of the annual and semiannual signals (seasonal impact), as it was the case with the City of Zagreb Geodynamic Network. However, according to [14], when processing GNSS observations over a period of less than 2.5 years, annual signals can cause significant errors in the determination of point velocities. When a period of more than 4.5 years between the observations is available, the velocity determination error drops significantly to negligible values [14]. Therefore, it is recommended that a minimum period of 2.5 years be used as the standard for determining the velocity solution, so that the velocities of the points are reliably determined. Moreover, the velocity data for geophysical interpretations obtained in a shorter period of time should be taken with caution. According to [20], this recommendation was adopted for processing with the Bernese v5.2 GNSS software, and it is used in the FODITS module.
The uncertainty of velocities decreases with the expression: 1 / √ the number of samples, but is not a function of time (duration of measurement) that depends on the time span. The uncertainty of 1 mm / yr can be achieved at best in a one-year period, but also at worst in a 4.5-year period.
From one-year periods between observations, the velocities of point cannot be reliably determined on the basis of a linear velocity trend. By analyzing the calculated annual and multi-year velocities for the years from 2009 to 2018, results were obtained for the CROPOS stations ZAGR and ZABO that correspond to previous studies [14]. On the example of the point ZAGR, determining the velocity in the height component for a one-year period versus multi-year solutions shows amplitudes of up to 1.4 cm (Fig. 11), while in the horizontal components (E and N), they are up to 7 mm. These differences can be attributed to periodic and stochastic influences on the determined velocity. For this purpose, it is necessary to determine the annual parameters of velocity corrections from the previous analysis of the time series of coordinates for annual periods. In this way, kinematic parameters containing periodic and stochastic influences can be determined, which enable the definition of real values of velocities for annual periods. According to [14], the errors in determining velocity considering the time component were   It can therefore be concluded that the research conducted in this paper shows similar results as in [14]. It can be seen from Fig. 11 that for periods longer than 2.5 years, the differences in velocities begin to equalize, while after 4.5 years, they are negligible. The variation of the annual component mainly comes from hydrological and atmospheric influences. Therefore, if annual signals are not taken into account, they can significantly reduce the accuracy of determining the point velocities used for highly accurate purposes such as plate tectonics or reference frames, and especially for geodynamic and geokinematic research.

CONCLUSION
Using the GNSS observation data from the geodynamic network of the City of Zagreb for the period from 2006 to 2009, different error models and GNSS measurement processing strategies were analyzed with the possibility of eliminating systematic and random influences. It all made it possible to determine the coordinate and velocity reliability and it showed the possibilities of using it for geodynamic and geokinematic analyses. This research has shown that when applying different parameters such as a local reference solution versus a global one, the identity of the reference frame must be ensured, because otherwise, solutions can be interpreted ambiguously in a geodynamic analysis. Furthermore, certain differences in the coordinate and velocity solutions are visible when using two software programs (Bernese and GAMIT/GLOBK), which results from different adjustment and velocity calculation strategies. The research showed that when applying different software solutions for processing high-precision measurements that use different processing strategies, the differences in the velocity of up to 1 cm per the horizontal components E and N and of up to 2 cm per the height component are possible. However, the mean values are at the mm level. Uncertainties in determining velocities are visible in the annual periods since they contain a periodic and stochastic signal. It can be seen that for linear velocities in the 2.5 periods, these differences are below mm, which directly indicates the connection between a higher measurement noise and the processing strategies implemented in software solutions. Higher values were obtained at several points, and in the case of significant shifts, it is necessary to perform the analysis by using the time series analysis, because in this way, we obtain reliable mean values that are within repeatability. Therefore, in order to monitor the geodynamics of a particular area of interest, it is necessary to use a previously noise analysis on the daily time series of coordinates in the area of interest with the aim of determining the true character of the geodynamic change. The research and analyses performed in this paper indicate certain uncertainties in determining the velocities from periodic measurements. They are affected by periodic signals and can cause significant errors in determining point velocities. Therefore, several studies recommend a minimum period of 2.5 years to be used as the standard for determining the velocity solution in order to reliably determine the velocities of geodynamic points. Velocities which are used as kinematic parameters for geophysical interpretations and that are obtained from the one-year period between the observations on geodynamic points must be interpreted with a detailed analysis. On the other hand, permanent GNSS stations brings a new quality to geodynamic research and the future of global and local geodynamic research is under constant monitoring. The analysis of changes in the time series of coordinates on permanent stations was also used in this paper as an accurate and reliable source of information on local geodynamics. At shorter time intervals (such as periodic observations on geodynamic networks), these values are influenced by the errors of different periodic and stochastic components. Moreover, one of the assumptions in the GNSS observation processing is that measurement errors are random and timeless (white noise), which leads to unrealistically estimated uncertainties in determining velocities, and affects the height component (~ 1 cm) the most. Previous research shows that white and pink noise dominate the GNSS noise spectrum for the time series of coordinates, hence it is necessary to calculate the true value of noise based on the analysis of the time series of the coordinates and the application of one of the models of the spectral analysis or empirical formulas. To determine the geokinematic parameters, the data of several daily GNSS solutions of the CROPOS stations were analyzed, since the determination of the displacement parameters requires a previous statistical analysis of the daily solutions, which provides control for the jumps and discontinuities. A permanent network of the CROPOS was available from 2009, therefore, the analysis of the time series of coordinates for the two stations ZAGR and ZABO was performed by using the FODITS algorithm and new elements of the functional model were determined. Such parameters can be used to determine the true character of the displacement/velocity of the geodynamic point (GMGZ) by using one of the methods of spatial interpolation. Further determination of parameters from the analysis of the time series of the CROPOS station coordinates will enable the use of the Geodynamic Network of the City of Zagreb for geodynamic and geokinematic purposes (the monitoring of the displacement parameters due to earthquakes or from other sources). The research conducted in this paper shows that in the geodynamic analyses for the interpretation of displacements, in addition to geological or seismological data, the influence of processing parameters for estimating reliability should be determined. The GNSS processing strategies analysis, as well as the determination of local geodynamic parameters from the analysis of the time series of coordinates, should be considered. The determination of the kinematic parameters at the points of the geodynamic network of the City of Zagreb will enable a further analysis of the observation epochs and determination of the kinematic model for the monitoring of the geodynamic and geokinematic phenomena. Additionally, an atmospheric impact analysis should be considered for an additional interpretation of the geodynamic and geokinematic parameters. Furthermore, for detailed analyses, it is necessary to establish denser permanent stations (~ 5 km or more) in the area of interest (the City of Zagreb and Medvednica).