ANALYSIS AND MANAGEMENT OF GNSS REFERENCE STATION DATA

Original scientific paper Global Navigation and Satellite System (GNSS) is used in a wide range of applications from low accuracy navigation to high accuracy point positioning. GNSS reference stations are required for using the GNSS applications effectively. GNSS reference stations at the Yildiz Technical University in IstanbulTurkey (YLDZ) and at the Afyon Kocatepe University in Afyon-Turkey (AFKU) were established within the scope of the research projects supported by their universities in 2013. In this study, quality control procedure, data process, analysis and management of YLDZ reference station data are explained. Besides, test results about the Real Time Kinematic (RTK) positioning method of the reference station are given. The reference station YLDZ meets the International GNSS Service (IGS) standards based on the processing results. According to the RTK results, coordinates can be determined with 1,5 cm and better accuracy on the field by using YLDZ station.


Introduction
The first GPS reference stations were set up to transmit DGPS corrections to improve the ship navigation accuracy on North America.Depending on the developments on the accuracy of GNSS measurement techniques, reference stations have been established all over the world.They are used to monitor the Earth's crust, to provide geodetic control, to support surveying, engineering, GIS data collection, machine control and precise positioning.
Centimeter-level accuracy positioning in real-time based on GNSS was developed in the mid-1990s and is currently named as RTK positioning.It involves a reference receiver transmitting its raw measurements or observation corrections to a rover receiver by a communication link.The data processing at the rover receiver includes ambiguity resolution of the differenced carrier phase data and coordinate estimation of the rover position.In this single base RTK method, the maximum distance between reference and rover receiver must not exceed 10 to 20 kilometers for determining the point coordinates reliably.This limitation is caused by distancedependent biases such as orbit error, ionospheric and tropospheric signal refraction.Continuous GNSS reference stations were developed in order to model and, accordingly, to reduce the distant-related errors that significantly affect the accuracy of RTK positioning.
Typical continuously operating GNSS reference station consists of a GNSS receiver and antenna mounted on a stable mast.The receiver operates continuously to store the raw GNSS observations, and stream the RTK and DGPS corrections to any of the RTK, GIS and GNSS navigation applications.The receiver is usually controlled by the remote computer system, which checks its health and functionality.The observed GNSS data is stored on central FTP server where the users can access their accounts.Currently, numerous GNSS reference stations networks are operating all over the world such as International GNSS Service (IGS Network), GNSS Earth Observation Network System (GEONET) in Japan, European Reference Frame (EUREF).On the other hand, The International GNSS Service (IGS)-USA, Centro di Geodesia Spaziale G. Colombo (ASI), Matera -Italy, Bundesamt für Kartographie und Geodäsie (BKG) -Germany, Institut Géographique National (IGN) -France, Slovak University of Technology (SUT) Slovakia, GNSS Analysis Center at Shanghai Astronomical Observatory (SHAO) -China, Geospatial Information Authority of Japan (GSI)-Japan are the examples of the analysis centers that analyze the GNSS reference station data.
Reference stations at Yildiz Technical University in Istanbul (YLDZ) and Afyon Kocatepe University in Afyon (AFKU) were established as a part of the two research projects supported by their universities.These two stations are planned as the first nodes of new reference stations of GNSS Geo-Sensor Network (GGN-TR) in South West of Turkey.The main purpose of these two active stations is streaming the GNSS corrections for RTK and DGPS applications such as monitoring the crustal movements and being GNSS data source for engineering application in the regions.In the scope of the research project, GNSS data collected by YLDZ station is processed and analyzed continuously and published on a web site for the users automatically.The ultimate aim of the project is to operate on an international GNSS data analysis centers.In this study, data processing and management scheme of GNSS Reference Station YLDZ is explained.

GNSS reference stations and analysis of data
Yildiz Technical University GNSS reference station (YLDZ) has been operational since July 2012.The station was established on the roof of the Civil Engineering Faculty building in Istanbul, (Fig. 1).There are not any physical obstacles affecting the signal quality around the station.Spectra Precision Ashtech Proflex 800 GNSS receiver and Aero Antenna AT 1675-120 SPKE Choke Ring GNSS antenna is used on YLDZ reference station.Proflex 800 receiver is a powerful positioning solution that delivers state-of-the-art RTK features in a rugged, highly integrated receiver design.The Z-Blade GNSScentric technology of the receiver uses all available GNSS signals equally (without preference to any particular constellation) to deliver fast and stable RTK solutions.AT1675-120 Choke ring antenna is known for very smooth, with minimal ripple, amplitude and phase pattern.Therefore, it provides excellent ability to reject multipath signal from any source.The choke ring antenna of GNSS receivers can provide millimeter precision measurements for surveying and geological measurement applications.The GNSS observation is processed by a computer program named PEGASUS that is developed by the project team members.PEGASUS is an in-house program with MATLAB scripts and runs almost full automatically in the high-capacity Linux servers without human intervention except for periodical performance upgrade maintenances.Two redundant Sun® Blade servers running PEGASUS program for data processing and data storing are hosted in climate controlled IT facility of Istanbul Kultur University, partnering institution in the projects.PEGASUS program consists of four modules, (Fig. 2).

PEGASUS -pre process
Observation rate 1Hz GNSS signals received throughout the day are stored in the receivers main storage where the raw observation data are converted to RINEX (v2.11) format at the end of the day 24:00 (UTC 22:00).Having converted to 1sec RINEX data, the PEGASUS Module 1 starts automatically at 06:00 (UTC 04:00) every morning since January 10 th , 2014 the first day the automatic process initiation.PEGASUS reaches the ftp addresses of the stations' receivers and retrieves the previous day's RINEX and navigation data.If the GNSS stations are responding and data are successfully retrieved, Module 1 applies a decimation routine to 1sec RINEX data in order to derive 10 and 30sec RINEX files.
Data quality check of GNSS permanent stations is a crucial topic for all the applications where their observations can be employed, both in post-processing and in RTK positioning.Therefore it is necessary to have analysis tools able to characterize data quality both from the point of view of the site where the station is installed (multipath and electromagnetic interferences) and from the point of view of the receiver precision.Presently the most used software to perform these investigations is TEQC (Translation, Editing and Quality Checking), a freeware package distributed by UNAVCO (UNiversity NAVstar COnsortium).TEQC translates from the binary receiver format to the RINEX format, to edit existing RINEX files, and to quality-check the data before post processing of GNSS data.The decimation routine, also developed during the project, also counts the epochs of the observations in RINEX data availability check currently and possible future developments of quality control and assurance.To that date, the quality control of the 24h RINEX data is also computed using TEQC software [1].

Figure 3 Data quality of YLDZ station
Multipath on L1 and L2 frequencies (MP1 and MP2), which is calculated by TEQC software, is a quality standard for GNSS Reference station data.The multipath values of YLDZ and the other adjacent IGS reference stations data are examined in Fig. 3 and multipath of YLDZ is found acceptable compared to some other IGS stations.

PEGASUS -main process
After the pre-process phase, all data was processed using version 10.5 of the GAMIT software [2].GAMIT is the software used for processing and analyzing the GNSS data on UNIX platforms.In PEGASUS main process phase, the 30-second interval Rinex data of the YLDZ GNSS station is processed on daily-basis since 2012.55GNSS day.Using a multistep approach used for PEGASUS main process, the ionosphere free LC linear combination of GPS L1 and L2 carrier phase observations were analyzed using the GAMIT software [2].In the first step, satellite orbits, earth orientation parameters, station coordinates, tropospheric zenith delays and tropospheric gradients were estimated from the double differenced GNSS phase observations for each day, applying only loose constraints to the estimated parameters.In this step we included data from 16 surrounding International GNSS Service (IGS) stations to link our study region (Anatolia) measurements to the International Terrestrial Reference Frame (ITRF00-05-08).
In the analysis, the antenna phase center was derived according to the height dependent model.The precise orbit information (IGS-F) in the SP3 format required by the GAMIT process was downloaded from the Scripps Orbit and Permanent Array Center (SOPAC) database [2].The Earth Rotation Parameters (ERPs) were obtained from USNO Bull.B (United States Naval Observatory Bulletin B).During the analysis, LC (L3) the ionosphereindependent linear combination of the L1 and L2 carrier waves was used.The FES2004 OTL grid was used to interpolate OTL components from a global grid [2,3,4,5].
The accuracy of the established evaluation strategy can be checked through the NRMS values obtained as a result of the evaluation.NRMS values are expected to be within the range of 0,15 ÷ 0,25.The fact that the NRMS value in the daily solutions was within this range indicates that the established model is accurate and that the noise levels of the measurements are at an acceptable level [2].In the process results, the RMS values of coordinates estimated from IGS stations are within 2 ÷ 3 mm level.

PEGASUS -analysis
Time series of the reference stations determined by using the coordinates obtained as the result of the daily processes.Linear, periodic and irregular movement of the reference stations can be determined by time series analysis, Eq. (1 x (t i ) time series t i (i = 1, 2, 3, …, N) of the reference stations can be written as below for the North, East and Up coordinates [6].In conducting the time series analysis, data must be regular and free of problems.As it can be seen in Fig. 4, firstly time-dependent graph of the data series is drawn and the series can be exercised visually.
Sometimes instant and non-recursive movements can be obtained in the series.These types of measurements are named as outliers and are either much lower or higher than the average of the measurement on the series.In this study, Addive Outlier (AO) method defined by [7] was used for determining the outliers.Additive outliers are usually caused by wrong data input and measurements and affect only one measurement.These outliers either much bigger or smaller than the expected value can be eliminated from the series.If some measurements of the x (t i ) series are affected by additive outliers AO(t i ) on t i time, the new series {y(t i )} can be defined, Eq. (2).
). ( ) ( ) ( Generally the parameter of each component in the series is estimated by Least Squares method (LS).If there is no outlier in the series, this method is useful and effective for parameter estimation.But if there is outlier in the series, LS method is defective for parameter estimation.
Because of the limitation of the LS method explained above, outliers were investigated by Robust Estimation method in this study.This method gives reliable results if there is outlier in the series.In the robust estimation method, function of the corrections is minimum instead of the root squares of the corrections.The series is analyzed by using re-weighted LS method [8,9,10,11,12,13,14]. Outliers in time series are usually determined by using suitable function.This function is sum of the linear and trigonometrically (periodical movement) functions [14].In the time series analyses of GNSS stations, linear movements obtained from trend component analyses are extremely important in terms of earth sciences.This is because "a k •t" in trend component is used to determine the plate movement (or velocity) which on the time unit and the direction [6].Also, it is presumed that this annual periodic effect, which can be seen in time series, is due to atmospheric and hydrologic loadings [15,16].It is additionally determined by [17,18] that here are annual, half-year and 3-month termed periodic movements present in the series, Eq. ( 3), due to climatic effects such as accumulation of ice on antennas and their shields.
Stochastic component was not added to function for outlier detection. ).
In Eq. ( 3), A is the design matrix and x T is the transpose of the unknown parameters of the linear trend and trigonometric function.Frequencies of the annual, semiannual seasonal periodical movements are respectively (unit: 1/day); f 1 =1/365, f 1 =1/(182,5) and f 1 =1/(91,25).The x unknown parameter is determined by re-weighted LS method iteratively with the equations given below, Eq. ( 4).

 
). ..., , Here, v is the vector of the N residuals (or errors) and W is the new weight matrix and w i is chosen as weight function in Eq. ( 4).w i is taken as 1 for LS method.
In Fig. 5, time series of East coordinate (raw data), function used for time series, parameters of the function which were determined by Bi-square weighted robust method and computed East coordinates are given.
Bi-square weighted robust predictor from M-Robust Estimators was used in this study to minimize the influence of outliers.Bi-square weights minimize a weighted sum of squares, where the weight given to each observation point depends on how far the point is from the fitted line.Points near the line get full weight whereas points far from the line get a reduced weight.In Fig. 6, calculated weight factors by the Bi-square weighting function for East coordinate were given.Depending on our experience from previous studies [6], data was accepted as outlier data if weight factor of data w i < 0,6.Robust fitting with Bi-square weights uses an iteratively reweighted LS algorithm, and follows the procedure outlined below: 1) Fit the y(t i ) function with weighted Least Squares (LS) 2) Compute and standardize the adjusted residuals.The adjusted residuals can be computed by Eq. ( 5), where v i are the usual LS residuals and h i are leverages that adjust the residuals by down weighting high-leverage observation points, which have a large effect on the LS fit.
K is a tuning constant equal to 4,685, and s is the robust variance given by MAD/0,6745 where MAD is the Median Absolute Deviation of the residuals.
3) Compute the robust weights as a function of u.The bi-square weights are given by 4) The procedure is completed when the fit converges.Otherwise, next iteration of the fitting procedure should be performed by returning to the first step.Outlier in times series are eliminated from the series by using Bi-square weighted robust predictor methods.New values are determined for the eliminated values in the series.Also loss of data can be determined by this method.Re-arranged time series for East coordinate is given in Fig. 7.
After applying this work routine for the entire coordinate components (North, East, Up), for determining the periodical movements of the GNSS station, coefficients of the linear function are estimated.In Fig. 8, linear function for the East coordinate of the YLDZ station is given.As explained in PEGASUS Pre-Process phase, after the quality checks, the GNSS data are stored for web site.Apache web server on Linux computer is used for publishing the GNSS data and other products.MySQL database and PHP web-programming language is used for the web portal design.In the quality checks, some statistical values about the data are obtained and these values are automatically stored on the MySQL database for web site.The structure of the database is given in Fig. 9.As it is explained in PEGASUS Main Process, GNSS data are processed with GAMIT software automatically on daily basis.Determined coordinate values are transferred to PEGASUS-Analyses phase immediately and final analysis results are then stored on the database.By means of these automatic work routines without human intervention, fully automation for GNSS data processing, analysis and publishing for YLDZ site is accomplished.
The web portal can be reached at http://www.gnss.yildiz.edu.tr/.System is open to public with a simple user registration procedure.When the users log in the system, GNSS data and analysis result can be reached for free.

Test of RTK positioning quality
One of the services that GNSS station provides is the RTK corrections for various engineering applications.GNSS reference stations can issue the corrections for RTK applications by using its known precise position and determined position.These correction values are transmitted to users via Internet or UHF radio modem.Rover GNSS stations on the site receive these correction values and determine their own position with cm-level accuracy in real-time.
A test mechanism was established to test the correction values issued by YLDZ station.A fixed rover named S1 station is used for this test (Fig. 10).North, East and Up coordinates of the S1 station were determined in 10 seconds interval for 24 hours.7800 measurements were collected during the test.Deviation values of the coordinates from the average value of the 7800 measurements were determined in every epoch and given in Fig. 11.
Also the deviation values of the coordinates from the known position of the S1 station were determined in every epoch and given in Fig. 12.  Distances between YLDZ and S1 were determined by using the instant S1 coordinates and known coordinates of YLDZ.The distance between the points is also known as 5,013 m.Deviation values of the distances from the known distance value are given in Fig. 13.Average deviation value is determined as 0,3 cm.Height difference between YLDZ and S1 station is determined by geometrical leveling method precisely.The determined height difference is 0,183 m.Deviation of values of the height differences which were determined by RTK method from the known height difference value

Figure 1 Figure 2
Figure 1 Reference station equipment 1) a k − The parameters of the trend component, b s , c s − The parameters of the periodical component, f s − Frequency α j − The parameters of the Autoregressive (AR(p)) model β jj − The parameters of the average movement (MA(q)) model v(t i ) − Random errors (the mean is zero and the variance is the σ).

Figure 4
Figure 4 Time series for East coordinates of YLDZ station

Figure 5
Figure 5 East coordinate and linear + trigonometric functions

Figure 6
Figure 6 Weight values of each observation day computed by the Bisquare weighting function

Figure 7
Figure 7 Re-arranged time series for east coordinate

Figure 8
Figure 8 Linear trend fit of East coordinates 3 Management of reference station data

Figure 9
Figure 9 Structure of the database for YLDZ station A web portal is designed for publishing the GNSS data and their statistical values.In the web portal, users can query the GNSS data by date, and can download the GNSS data.By automatically executing work routine by the means of PEGASUS-Pre Process, GNSS data and

Figure 10 Figure 11
Figure 10 RTK test mechanism

Figure 12
Figure 12 Deviations of the RTK coordinates from the known position of S1 point

Figure 13
Figure 13 Deviation values of the distances from the known distance value

Table 1
Deviations of the RTK coordinate