DEVELOPMENT OF A TIME-DEPENDENT 3-PARAMETER HELMERT DATUM TRANSFORMATION MODEL: A CASE STUDY FOR MALAYSIA

This paper aims to develop a time-dependent 3-parameter Helmert datum transformation model for Malaysia as a proposed solution to the current non-geocentric issue of the Geocentric Datum of Malaysia 2000 (GDM2000). Methodologically, the datum transformation models is categorised into three parts; firstly, the time-dependent aspect of the datum transformation model is determined using the tectonic motion velocities computed from linear least squares regression of the long-term time series of MyRTKnet stations positions from year December 2004 to 2014; whereby the station positions are obtained from high-precision daily double-difference processing of MyRTKnet and IGS stations via Bernese 5.0. Secondly, the 3 Helmert translation-only parameters, are derived between the original GDM2000 and GDM2000@2013 – the new datum coordinates which refers to ITRF2008 at epoch 3/7/2013 – via Bernese 5.0 software. Thirdly, a distortion model is computed in order to minimise the coordinate residuals between the ‘processed’ and ‘transformed’ new datum. The datum transformation model is then validated to determine the reliability of the model. The validation results show that the datum transformation model is within centimetre-level accuracy, i.e., below 3 cm, over Malaysia for forward transformations to year 2014 and 2015. Therefore, this study anticipates that it will contribute as a feasible solution for the GDM2000 issue with consideration of the core concern: the complex tectonic motion of Malaysia.


INTRODUCTION
The Geocentric Datum of Malaysia 2000 (GDM2000) is established by the Department of Surveying Malaysia (DSMM) via the Malaysia Real-time Kinematic Network (MyRTKnet) geodetic infrastructure.Presently, the GDM2000 datum coordinates are referenced to ITRF2000 at epoch 2006.On the other hand, Malaysia has been affected by earthquakes over the years, primarily after the December 26, 2004 9.2 Mw Sumatra-Andaman earthquake; followed by the 8.6 Mw Nias-Simeulue earthquake on March 28, 2005, the 8.5 Mw Bengkulu earthquake on September 12, 2007 and the 8.6 Mw Northern Sumatra, or Indian Ocean, earthquake on April 11, 2012(USGS, 2016).
As a result, the tectonic motion of Malaysia has displaced the datum coordinates; causing GDM2000 to be non-geocentric (see Shariff et al., 2014;Gill et al., 2015).Furthermore, the current International Terrestrial Reference Frame (ITRF) is ITRF2008 (with respect to the period in which this study was carried out) which is, in fact, an enhancement to ITRF2000.Hence, these would create issues such as mismatch of base maps and geospatial database, coordinate bias due to inconsistency with satellite orbits, and decreased accuracy of reference stations coordinates which in turn affects the legal traceability of the coordinates.
Therefore, GDM2000 needs to be updated in order for the issues to be rectified.However, due to tectonic motion, the geodetic datum would need to be updated over time.This then would result in an array of epochs for the geodetic datum causing confusion at the user level.Hence, a time-dependent datum transformation model would be most appropriate -as users would update the geodetic datum to a desired epoch with reference to the official datum epoch; thus, the aim of this study: to develop a time-dependent datum transformation model for Malaysia.
In this study, the well-known Helmert datum transformation model is chosen with only its 3 translation parameters selected.The translation parameters are exclusively chosen as the transformation is between satellite datums, i.e., GDM2000 and GDM2000@2013 -both referenced to the ITRF; hence, theoretically, there should not be significant rotation and scale values, whereby the latter is predominantly dependent on the changes in ellipsoid.Moreover, a 3-parameter transformation is simple to realise and update, especially in terms of its timedependent aspect, of which is represented by the tectonic motion velocities of Malaysia in this study.It should be noted that the measured tectonic motion is with regard to the crustal motion, i.e., within the Sunda block.This paper will firstly describe the methodological facet in developing the time-dependent 3-parameter Helmert transformation model for Malaysia.Next, the tectonic motion analysis and datum transformation model's validation results are discussed.This paper will attempt to further improve the tectonic motion analysis by Gill et al. (2015) with the addition of 2014 MyRTKnet data; however, less emphasis will be given as it slightly deviates from the aim of this paper.Lastly, the conclusion and final remarks are outlined.All in all, this study anticipates that it will contribute as a solution for the GDM2000 issue with consideration of the core concern: the complex tectonic motion of Malaysia.

METHODOLOGY FOR THE DEVELOPMENT OF A TIME-DEPENDENT 3-PARAMETER HELMERT DATUM TRANSFORMATION MODEL
Several steps were implemented for the development of the model.Firstly, the time-dependent aspect of the datum transformation model was determined using the tectonic motion velocities computed from linear least squares regression of the long-term time series of MyRTKnet stations positions from year December 2004 to 2014; whereby the station positions were obtained from high-precision daily double-difference GPS processing of MyRTKnet and IGS stations via Bernese 5.0.The time series was plotted and linear least squares regression was computed via GPS Interactive Time Series Analysis software (GITSA) (Goudarzi et al., 2013).A velocity vector map was then plotted in order to visualise the tectonic motion of Malaysia via Generic Mapping Tools (GMT) (Wessel and Smith, 1998).Secondly, the 3 translational parameters, were derived between the original GDM2000 and GDM2000@2013 -the new datum coordinates which refers to ITRF2008 at epoch 3/7/2013 -via Bernese 5.0 software (Dach et al., 2007).
Thirdly, a distortion model was computed in order to minimise the coordinate residuals between the original and new datum.The datum transformation model was then validated with coordinates at epoch 9/7/14, 1/4/15, and 1/7/15 in order to determine the reliability of the model.

High-precision GPS processing
To estimate daily solutions of the MyRTKnet stations, Bernese version 5.0 was used by employing its double difference quasiionosphere free (QIF) strategy.65 MyRTKnet stations and 24 IGS stations were chosen with GNSS data spanning from December 2004 to December 2014.Only 15 out of the 24 IGS stations were selected (see figure 1) as fiducial stations for datum definition as they represented stable motions throughout the data time span.The processing strategy and parameters adopted are given in table 1.

Velocity vector estimation and mapping
After the daily solutions were estimated, a time series of daily solutions for the selected MyRTKnet stations were plotted using GITSA, a software developed by Goudarzi et al. (2013) for time series analysis using MATLAB.Once each station's time series was plotted, the outliers were removed at 99% confidence level within GITSA as it may affect the linear regression line later for estimating the velocity vectors.The determination of velocity vectors from linear least squares regression must fulfil two criteria: (1) minimum of 2.5 years solution in order to reduce annual and semi-annual effects in geodetic time series, of which will cause biased estimated velocities (Blewitt and Lavallée, 2002), and (2) time series with long data gaps, i.e., few months, are not chosen to estimate the velocity vectors.Hence, steadystate periods of at least 2.5 years were chosen to compute the velocity vectors.These steady-state periods are: the 2008 to 2011 period which represents the post-seismic tectonic motion after the pivotal 2004 Sumatra-Andaman earthquake, and the October 2012 to 2014 period which represents the post-seismic motion after the 2012 Northern Sumatra earthquake.Refer to figure 2 for visualisation and explanation.Once the velocity vectors for each station is computed, they were plotted on map using GMT.It is noted that the velocities are in centimetre per year (cm/yr) and in Cartesian system, but converted to metre per year (m/yr) for the datum transformation model.
There are a few reasons why the tectonic velocity vectors are required for the datum transformation model; primarily, it is to compute the time-dependent aspect, other reasons include: (1) determination of size of the area of transformation (see Mitsakaki, 2004), and (2) common points' selection based on conformality (see Collier, 2002).These will be discussed later.

Time-dependent aspect estimation for the 3-parameter datum transformation model
The concept of utilising tectonic motion for the time-dependent aspect was adapted from Stanaway and Roberts (2009) who used Euler pole rotations for their 3-parameter datum transformation model for Australia.The time-dependent aspect is commonly known as the parameters' rates of change.The rates of change were determined by computing the average velocities of the MyRTKnet stations that were selected as common points for the transformation.

3-parameter estimation of the Helmert datum transformation model:
As for the 3-parameter estimation, there are a few components that must be taken into account before estimating them.Firstly, the new datum coordinates must be realized.The method to realize a datum is slightly complex.
There is no specific length of time for realizing a new datum, it can be from a few years (see Mateo and Mackern, 2012;Habrich, 2007) to a week or two (see Gianniou, 2010).In this study, 1 week data, i.e., GPS week 1747, was processed and combined to realise the new datum at epoch 2013.5047 or 3/7/2013 -termed GDM2000@2013.This week was selected as it was a stable week, i.e., coordinates repeatability below 1 cm and 3 cm for the horizontal and vertical components, respectively.The datum is realized with the same procedure as section 2.1, i.e., minimum constrained to the IGS coordinates.The IGS coordinates employed was from the published IGS weekly datum coordinates which are in SINEX format and available at ftp://cddis.gsfc.nasa.gov/gps/products.Hence, the new datum refers to the IGb08 frame which is basically a subset of ITRF2008.Refer Rebischung (2012) for more information on IGb08 and its weekly published coordinates.Reliability of the new datum coordinates was achieved by obtaining close-to-zero values for the 7 datum transformation parameters between IGS processed/combined coordinates and published weekly coordinates (Habrich, 2007) (shown later).
Once GDM2000@2013 was realised, the next component was the common points' selection.The selected common points must be conformal; particularly, with respect to this study, they must have almost similar velocity vectors in order to minimise the abortion of errors into the parameters.With the common points' selected, the 3 datum transformation parameters was then derived between GDM2000 and GDM2000@2013 via Bernese's Helmert program.

Distortion model determination
A distortion model that accommodates the coordinate residuals between the processed/combined coordinates and transformed coordinates of GDM2000@2013 was then estimated.The transformed coordinates are basically coordinates obtained after applying the 3 transformation parameters to the original GDM2000 coordinates.In other words, the distortion model functions to shift the transformed coordinates into the processed/combined datum coordinates so as to avoid any large discrepancies in the results.The distortion model is constant, i.e., it can be applied after the transformation at any epoch.

The datum transformation model validation
Three validation epochs were chosen to realize a new set of coordinates following the same datum realization procedure as section 2.4.These epochs are 2014.5912(9/7/14), 2015.2479(1/4/15), and 2015.4973(1/7/15).The model from section 2.6 is employed to transform the original GDM2000 coordinate to the chosen validation epochs and compared with the processed/combined coordinates.The average Helmert residuals in local (N, E, U) coordinate system and the 3 Helmert parameters will be used as the results for the validation.It is noted that only epoch 2014.5912uses 1 week data, i.e., GPS week1800, for it coordinates determination.The 2015 MyRTKnet data was not available during this study; hence, only 1 day data was utilised.The 2 days MyRTKnet data for the 2015 validation was downloaded from the E-biz JUPEM Geoportal.

RESULTS AND DISCUSSION
The results and discussion are categorised to three parts: (1) the tectonic motion analysis, (2) the time-dependent 3-parameter Helmert datum transformation model estimation, and (3) the validation of the time-dependent 3-parameter Helmert datum transformation model.

The tectonic motion analysis of Malaysia
From the Bernese processing of daily solutions between years December 2004 and 2014, most daily solutions provided excellent ambiguity resolution, on average, above 75%.Positional RMS error for MyRTKnet stations were mostly below 1 mm for the horizontal component and 1.5 ±0.5 mm for the height component.
Sarawak consists stations SEMA, UMAS, SARA, KAPI, BIN1, NIAH, and MRDI, while Sabah consists stations BEAU, KENI, TMBN, UMSS, RANA, BELU, MRDU, KUDA, MTAW, DATU, and SEMP.The division into regions is in fact due to the resulting complex motion after the 2004 Sumatra-Andaman earthquake.The average horizontal velocities is tabulated in table 2. The considerable difference between the North-west region and the rest of Peninsular Malaysia may be due to a so-called 'pull' towards the epicentre of the 2004 earthquake, causing the velocities to decelerate.In general, from 2008 to 2011, Peninsular and East Malaysia moves south-east at an average horizontal velocity of 1.91 cm/yr and 2.72 cm/yr, respectively.As a result of the 2004 earthquake, Peninsular Malaysia decelerated and experienced a velocity vector gradient.Next, the average horizontal velocities after the 2012 Northern Sumatra earthquake, i.e., the Oct-2012 to 2014 period, is given in table 3. The directions of the Peninsular Malaysia regions are almost uniform with East Malaysia.The North-west region velocity vectors also seems to be more consistent with the rest of Malaysia.In terms of velocities, the velocity vectors gradient along Peninsular Malaysia is much less visible; though there are slight differences -~3 mm on average -between the regions (see table 3) -evidently, there still exists a minor deceleration in the North-west region.This suggests that the post-seismic effect of 2004 earthquake is still present.As for East Malaysia, the tectonic motion continues a uniform motion alike the previous period, but with small differences in direction.Only 2 years and 3 months data were used which would not generate concrete results as the time series contains seasonal effects -as can be seen from the increase in the velocities standard error.In general, from Oct-2012 to 2014, Peninsular and East Malaysia moves south-east at an average horizontal velocity of 2.72 cm/yr and 2.93 cm/yr, respectively.Nonetheless, this length is sufficient for the study on the datum transformation model.It is stressed here that this finding does not assume that the postseismic deformation of Oct-2012 to 2014 is minimum and can be neglected; the post-seismic motion is still on going and it is adopted as part of the time-dependent transformation model.

The time-dependent 3-parameter Helmert datum transformation model
Firstly, the new datum is realized with reference to IGb08 at epoch 2013.5027(3/7/2013).Table 4 provides the agreement between the processed/combined IGS coordinates and the weekly IGb08 frame coordinates at epoch 2013.5027.Therefore, due to the close-to-zero values for the transformation parameters, it is evident that the MyRTKnet stations are successfully tied to the IGb08 frame at epoch 2013.5027 as well.Note that all 7 parameters are used to validate the new datum as the realization of international reference frames takes into account all these parameters.

Parameters Peninsular Malaysia East Malaysia
Translation in X (mm) -0.9 ±1.7 1.9 ±2.0 Translation in Y (mm) -1.1 ±1.5 -1.9 ±.7 Translation in Z (mm) 0.5 ±1.9 0.4 ±2.2 Rotation around X-axis (seconds) 0.00004 ±0.00008As for the size of transformation area, it is evident that Peninsular and East Malaysia were moving differently before the 2012 Northern Sumatra earthquake based on section 3.1.Therefore, two separate datum transformation models are required.As for the selection of common points, Peninsular Malaysia has an nonuniform pattern, which results in tedious selection of common points.Initially, stations from the central region of Peninsular Malaysia was selected as common points due to the fact it has the most stations available that are conformal; then, by trialand-error, the common points with the least Helmert residuals and transformation parameters standard errors were selected.East Malaysia followed the same trial-and-error process for common points' selection.Note that the common points' selection was based on the transformation between GDM2000 and GDM2000@2013.Next, the 3 transformation parameters are derived.
The distortion model is then computed, i.e., the coordinates discrepancy between the processed/combined coordinates and the transformed coordinates of GDM2000@2013 for each MyRTKnet station.The distortion model is absolutely necessary due to the large coordinate discrepancy of the North-west stations which had discrepancies up to 5 cm for station ARAU and UUMK for the X component.The distortion model contains the coordinate differences in Cartesian system for each MyRTKnet station; it is not an average value due to the large coordinate differences between the some regions of Peninsular Malaysia, namely the North-west region.
Lastly, the 3-parameters' rates of change are determined.The rates of change for Peninsular Malaysia is the average of the selected common points' velocities.

The validation of the time-dependent 3-parameter Helmert datum transformation model
The validation is carried out at three epochs: 2014.5912(9/7/14), 2015.2479(1/4/15), and 2015.4973(1/7/15).Each epoch is aligned to IGb08, i.e., processed/combined coordinates.Then, compared to the transformed coordinates via the model in equation ( 1).The comparison is carried out by analysing the average Helmert residuals which gives the average coordinates differences after 'alignment' between the processed/combined and transformed coordinates.Table 7. Averaged velocities in Cartesian system for Peninsular and East Malaysia.These averaged velocities will then be directly used in equation ( 1), the time-dependent 3-parameter datum transformation model, in section 2.6

Station
The 'alignment' value is the 3 translation parameters, whereby the smaller the better -good agreement between the 2 coordinate sets.This 'alignment' uses the same common points from table 5 and 6 for Peninsular and East Malaysia, respectively.Generally, if the average Helmert residuals is large, it indicates that some stations does not 'align' well between both coordinate sets.Thus, though the alignment may have close-to-zero translation values, it does not mean all stations, primarily the non-common points, would agree well.

Figure
Figure 1.Fiducial and non-fiducial IGS stations selected for Bernese processing

Figure 2 .
Figure 2. Time series of station UPMS (University Putra Malaysia) depicting the 2008-2011 and Oct-2012-2014 steady-state periods, and the post-seismic and co-seismic effects due to earthquakes.UPMS was chosen as it is affected by all the aforementioned earthquakes.The North and East components are in topocentric (local) system, while the Up component is the ellipsoidal heights.Though the time series is self-explanatory, for a more complete elaboration, refer toGill et al. (2015) The equation for employment of the 3 datum transformation parameters, their rates of change and the distortion model is as follows: Whereby, the T and S subscripts depicts the Target and Source 3D Cartesian coordinates, and ( , , ) are the translation components, ( , , ) are the averaged station velocities (rates of change), ( , , ) are the distortion model for each MyRTKnet station, and and denotes the target and reference, i.e., 2013.5027,epoch in decimal years, respectively.
post-seismic motion due to the 2012 N.S. earthquake

Figure 3 .
Figure 3. MyRTKnet velocity vector map for the 2008 to 2011 period.Note that stations are displayed as the station names and arrows will overlap one another.As for East Malaysia, most stations began operation in April 2009; hence, they are not selected for velocity vector determination.Station AMAN has a different direction and velocity as it experiences severe subsidence parameters between processed/combined IGS coordinates and the weekly IGb08 frame coordinates at epoch 2013.5027Secondly, since the datum transformation model covers the transformation from epoch 2006 to 2013, the 2008-2011 period has a significant role in determining the size of transformation area and selection of common points.

Table 3 .
Averaged horizontal velocities of Malaysia by region and their directions for the Oct-2012-2014 period.
Based on table 3 and figure4, the tectonic motion of Peninsular Malaysia is almost uniform compared to the 2008-2011 period.

Table 6 .
The time-dependent 3-parameter datum transformation model values for East Malaysia with its selected common points Table8and 9 provides the average (avg.)and maximum (max.)Helmert residual in local system (sys.)and the translation values for each of the validation epochs for Peninsular and East Malaysia, respectively.The maximum Helmert residuals are given based on the station with the highest value for each local system component.Note that HR is the acronym for Helmert Residuals in table 8. Based on table 8, Peninsular Malaysia stations have significantly small average Helmert residuals which suggests that most stations agree well with datum transformation model coordinate results.However, it seems station TOKA has local issues which results in large Helmert residuals.If TOKA is considered an outlier, the maximum residual would be 1.2 cm due to station CENE at epoch 2015.4973.As for the translation parameters, only epoch 2015.2479 has a large TY value.This stems from the processing quality, whereby the ambiguity resolution was below 70%, i.e., lower accuracy results.Note that ambiguity resolution should be above 75% for good accuracy results.All in all, the model's accuracy over Peninsular Malaysia is within centimetre-level, i.e., taking into account the maximum Helmert residuals.