The origin of levelling in Algeria dates back to the French colonial period. It was confided to the Geographical Service of the army (GSA), then taken in charge by the National Geographical Institute – France. Provided with a similar material to that used by the Service of the General Levelling of France (NGF), its operators elaborated from 1887 a geometric levelling network of limited precision to first and second order. Between 1953 and 1954, the levelling works were led in parallel with the re-observation of the North parallel and the 1st order complementary of the coastline. The lines of spirit levelling executed since 1976 by the National Institute of Cartography and Remote Sensing (INCT) followed in parallel the works of the geodesy by the densification of the North parallel. Nowadays, the Algerian levelling network includes, 33,757 km levelled. This coverage remains insufficient for some areas, particularly in South of Algeria and where new measurements are needed to accomplish a homogeneous coverage.

3.1 Presentation of the gravimetric geoid model

3.1.1 Gravity anomalies

Gravity anomalies used in this work were supplied by the Bureau gravimétrique international (BGI). These measurements with an initial precision of 1 m Gal were expressed in the Geodetic Reference System GRS67 and were transformed into the GRS80 system. Moreover, we applied an atmospheric correction recommended by the International Association of Geodesy (IAG, 1971) in order to eliminate the influence of atmospheric masses. All data have been checked carefully to remove duplicate points, and validated for gross errors by applying the Least Square Collocation method. The error ratio detected using this technique represents 1.21%, which remains small with regard to the number of measures used, and proves that the gravity data measurements provided by BGI are of good quality and subsequently, can be used for the Algerian geoid computation. Fig. 1 shows the geographical distribution of the available gravity data. From this figure, it becomes clear that the coverage with gravity observations is not sufficient for some land areas and new measurements are needed to accomplish a homogeneous coverage.

3.1.4 Geoid computation

The computations were made using the “GRAVSOFT” package, developed during a number of years at the National Survey and Cadastre of the Copenhagen University. For the actual solution, the Fast collocation technique and the Remove–Restore procedure were used to compute the geoid heights. The procedure is summarised as follows:

• • Available free air gravity anomalies ΔgOBS are reduced for the reference geopotential model ΔgGM and terrain effects ΔgRTM, yielding smooth residual anomalies Δgr.

  Δgr=ΔgObs−ΔgGM−ΔgRTM

  (3)

• • Gridding of Δgr to produce a 5’ × 5’ grid of residual gravity anomalies and then converted it into residual geoid undulations by Fast Collocation method.
• • Restore the effects removed in the first step (Reference geopotential model and terrain effects) on the prediction grid to get the final geoid model.

N=Nr+NGM+NRTM

(4)

The geopotential model OSU91A developed to degree and order 360 has been used as a reference to remove the long wavelength component of the gravity data. However, in order to smooth the anomalous gravity field of the test area, i.e. to remove the irregular variations of the topography, the Residual Terrain Model (RTM) reduction is applied, which takes into account only the short-wavelengths of the topography. The method employs a smooth reference surface which may be derived either from a spherical harmonic expansion of the global topography or by averaging the local elevation data. The computation of the effects of the topography according to the RTM reduction modelling method is based on two elevation grids: 1 km × 1 km digital terrain model for the North of Algeria; and a global topographic model ETOPO5 of 5’ × 5’ for the rest of the area which was used up to a distance of 200 km. The reference surface of 15’ × 15’ needed for the RTM reduction is formed from the 5’ × 5’ blocks of ETOPO5. This grid is smoothed further by taking moving averages over the 3 × 3 adjacent blocks. The statistics, in mGals, of the reduced gravity data are regrouped in Table 1.

Table 1

Statistiques, dans la zone de calcul, des données de gravité réduites.

Anomalies	Minimum	Maximum	Mean	Standard deviation
Δ g obs	−82.59	136.20	4.89	26.84
Δ g obs − Δ g OSU 91 A	−67.75	123.46	0.12	12.31
Δ g r = Δ g obs − Δ g OSU 91 A − Δ g RTM	−34.17	110.38	3.09	13.00

From the results of Table 1, it is obvious that the OSU91A reference field fits well the gravity in the area under consideration, and the smoothing of the gravity data is considerable after the removal of the topographic effect if we take into account only the mean and standard deviation values. The standard deviation of data reduces from 26.84 mGals to 13.00 mGals (a reduction of 50%).

The gridding of Δgr to produce a 5’ × 5’ grid of residual gravity anomalies between the limits 18° ≤ φ ≤ 38° and −9° ≤ λ ≤ 11° was done with the GEOGRID program using a fast quadrant-search collocation algorithm. The prediction of the grid gravity value was computed from the five closest points in each quadrant around the prediction point, using a second-order Markov model covariance function fitted to local data (correlation length x_1/2 = 25 km), and taking into account individual variance of the gravity data.

However, the first step in the use of the fast collocation for residual geoid heights estimates from the grid of residual gravity anomalies resides in the determination of the empirical covariance function and subsequently the selection of its corresponding analytic representation capable to describe the local behaviour of the gravity field.

The empirical covariance function of the BGI residual gravity anomalies was computed by EMPCOV program (Tscherning et al., 1992) using the following formula:

CδΔgi,δΔgj=C(ψ)=1M∑δΔgi.δΔgj

(5)

where M is the number of combinations, ψ is the spherical distance between Q_i and Q_j. δΔg is the BGI free-air residual gravity anomalies. The summation was made for all the combinations of the data points Q_i and Q_j where the interstation distances in the range between (ψ−Δψ/2) and (ψ+Δψ/2), and Δψ denotes the sampling interval size.

These empirical values were used in order to estimate the parameter values of a local covariance function model. As a model of the covariance function, we have adopted, in our case, a classical analytic covariance function for the disturbing potential proposed by Tscherning and Rapp (1974), which is expressed by a sum of series of Legendre polynomials (Moritz, 1980):

K(ψ)=KG(ψ)+KL(ψ)

(6)

with:

KG(ψ)=αGMRE2.∑l=2NεlRE2rr′l+1Pl(cosψ)KL(ψ)=∑i=N+1+∞A(i−1)(i−2)(i+4)RB2rr′i+1Pi(cosψ)

where ɛ_l: error degree variances of the reference field, N: the maximum degree of the reference field coefficients (here 360), α: scale factor of errors associated with the reference field, R_B: radius of a Bjerhammar sphere, r and r’: radial distances of the points, R_E: mean earth radius, A: scale factor of the degree variances.

The actual use of this model as a local covariance function requires the estimation of three parameters: the radius of the sphere of Bjerhammar (R_B) and two scale factors α and A. These parameters are estimated by fitting this model to a number of empirical covariance values using the least square adjustment. The results of the adjustment of the empirical covariance function on the Tscherning & Rapp model were obtained by the program COVFIT8 (Knudsen, 1987). Fig. 2 shows the empirical and analytic covariance functions of the residual anomalies. The agreement is optimum between them up to a spherical distance of around 0.30° (∼ 30 km) which is much greater that the correlation distance..

The gridded residual gravity anomalies have subsequently been converted to geoid heights by using the Fast Collocation method (Bottoni and Barzaghi, 1993). The Collocation was carried out on a grid of 241 × 241 points. The final geoid was obtained by adding the geopotential model and the residual terrain effect on the 5’ × 5’ residual geoid grid. The major contributions to the final geoid are coming from the OSU91A geopotential model. Fig. 3 shows a map of the geoid solution in Algeria contoured with 1 m intervals.

3.2 GPS/Levelling data

For this investigation, 247 precise GPS levelled points have been used, of which ten are benchmarks of the first order levelling network, and the others belong to the second order levelling network. All of these points are located in the North of Algeria. The geographical distribution of the available GPS/levelling data is shown in Fig. 4.

For our study and in order to assess the precision of the gravimetric geoid estimate, this model has been compared with 247 GPS/Levelling derived geoid undulations. All of these points are located in the North of Algeria and distributed as in Fig. 4 where most are close to the station of Arzew. The statistics of the absolute differences in benchmarks between the gravimetric geoid model and the GPS/Levelling derived geoid height are summarised in Table 2 and they indicate that the overall agreement between the gravimetric and GPS/Levelling derived geoid undulation is about 68 cm in terms of standard deviation and that there exist systematic biases between these two kinds of geoid representations with a mean value of about −1.481 m. These systematic biases are due to the datum inconsistencies inherent among the height types, each of which refers to a different reference surface and long wavelength systematic errors in N, distortions in the orthometric height datum due to an over-constrained adjustment of the levelling network, effect of various geodynamic effects, improper or non-existent terrain/density modelling in the geoid modelling and orthometric heights and negligence of the sea surface topography at the tide gauges, poorly modeled GPS errors (i.e. tropospheric effects), and assumptions and theoretical approximations made in processing observed data. The major part of these discrepancies is usually attributed to long wavelengths errors introduced by the use of a global model in gravimetric geoid computation that is not necessarily optimal in the area of interest, to the systematic errors and to datum inconsistencies inherent among the height types which can be reduced by the use a corrector surface model (Fotopoulous, 2003).

In practice, to minimize these deviations, we have tested a four and a third polynomial transformation models in order to find the best model able to represent the correction that should be applied to the gravimetric geoid.

The four-parameters model is the most commonly used in such adjustments and is given by the following formula (Heiskanen and Moritz, 1967):

The third order polynomial model may be adequate to absorb non-periodic large scale systematic errors, such as systematic levelling errors:

Among the 247 GPS levelled stations, only 16 well-distributed GPS points (see Fig. 2) are used as benchmarks in a least-squares adjustment to solve for the unknown parameters of the parametric model transformation, and the remaining points were excluded from the combined adjustment in order to estimate the real accuracy given by cross-validation. The fitting of the gravimetric geoid model to a set of levelled GPS reference points has been carried out by a program, called ADJ_GLG (Benahmed Daho, 2004) developed at Geodetic Laboratory of the National Centre of Space Techniques for combined adjustment between the gravimetric geoid model and GPS/Levelling data using the Hybrid approach (Kotsakis and Sideris, 1999). The statistics of the differences after fitting out the systematic biases and tilts between the gravimetric geoid and the GPS/Levelling data using the previous parameter models transformation are summarised in Table 3. These statistics in benchmarks show that a good fit between the Algerian geoid and GPS/Levelling has been reached when the third order polynomial was used as corrector surface. The standard deviation of these discrepancies between these two representations at benchmarks is about ± 1.4 cm after fitting. This discrepancy depends upon the quality and the density of data integrated in the combined adjustment, and represents the internal accuracy of the adjustment. In addition, and before any special tests are performed for checking particular features of the general adjustment model, some global statistical test must be used first. Levelled GPS points with global residual module larger that four times RMS were rejected. No outliers can be identified. This does not necessarily mean that outliers do not exist, but more exhaustive statistical means must be used to identify the possible blunders on the gravimetric geoid, levelling and GPS networks by using their individual component residuals. Since this is not the purpose of this work, it was assumed that all GPS/Levelling data points were valid. Furthermore, the Chi² (χ²) test with significance level α = 1% was employed. This test is positive and indicates that the model is not problematic. The next step is to test if a parameter or set of parameters in the model are significant. The application of the Student test on the parameters values shows that the estimated parameters are statistically significant.

However, the previous RMS value is not a real accuracy of the determined geoid; this provided proof that the combined adjustment can optimally fit the gravity geoid to the GPS levelling point in the least square sense. So, and in order to pronounce on the conclusion drawn on the parametric model choice and to assess the real accuracy of the adjustment, we have applied the cross-validation procedure. Its principal consists to compare the GPS/Levelling derived geoid undulations at control points with adjusted ones computed using the estimated parameters for the each model transformation. The statistics of these differences for tested parameter models are summarized in Table 4. We can see that the combination of the gravimetric geoid model with GPS/Levelling gives the best results in experimental area when the third order polynomial model was used as a corrector surface. The standard deviation of the differences between these two representations at control points is about ± 5.5 cm after fitting. In most points (89%), the differences are smaller than 5 cm. There are very few stations where really large differences occur. However, the obtained accuracy cannot be generalized for whole country, considering the number of GPS\Levelling points and their distribution.

One of the main goals of a combined adjustment between levelling, GPS and the geoid is to substitute expensive levelling through cheaper GPS measurements at least in regions where the highest accuracy is not necessary. In order to test the capability of GPS for height determination, a special GPS campaign was carried out in 2004 for the altimetric auscultation of an oil storage reservoir located at 40 km from Oran (in the West of Algeria). A total of about 132 new GPS stations, for which fifty points are considered as benchmarks, were established. The GPS observations were performed using ASHTECH Z-12 dual frequency receivers with a observation period between 15 mn and 45 mn and were processed with the Bernese GPS software version 4.2 using the precise ephemerides supplied by IGS and the TYRGEONET station (Arzew) as reference point. Code pseudo ranges were used to estimate receiver a priori co-ordinates and clock synchronisation with respect to GPS time. All GPS stations have been connected to the national height system through spirit levelling, which consists of orthometric heights. The comparisons between GPS and levelling were performed by transforming the GPS ellipsoidal heights into orthometric ones using the adjusted gravimetric geoid heights. These differences on the benchmarks and their statistics, after fitting out the systematic biases and tilts using the optimal four-parameter transformation model, are summarised in Tables 5 and 6, respectively. Their analysis permits us to conclude that the orthometric heights on benchmarks can be deducted from the GPS observations with a precision better that 2 cm. It is about 3 cm in the remaining points. We note that this precision remains not sufficient for GPS/Levelling purpose but it can be used in low accuracy scientific applications.

A detailed study of the possibility of substituting classic spirit levelling by GPS/geoid derived orthometric heights was presented. The main goal of this work is to assess the precision of the transformation of the ellipsoidal heights from GPS into orthometric heights in the northern part of Algeria. The field experiment introduced above clearly proves that the combination of GPS/Levelling and gravimetric geoid model computed from the gravity data supplied by BGI is capable of producing orthometric heights with an acceptable accuracy for the low order levelling network densification with regard to the actual national levelling network status. The accuracy obtained is thought to be quite sufficient for many applications compared to the alternative method of spirit levelling, which is time consuming and impractical from the financial point of view. Unfortunately, the non-availability of GPS levelling data in the whole of the country, with a homogeneous distribution and sufficient density, does not allow a more reliable assessment at the national scale of the quality of the used gravimetric geoid model.

The analysis of the results shows that the third order polynomial model transformation used for modelling the deviations between the gravimetric geoid and the available GPS/Levelling points gives the best results in the experimental area. It also shows that the four-parameter transformation model usually employed in studies dealing with the combined adjustment of GPS, levelling and geoid data is not always sufficient for removing the systematic errors introduced by the datum differences of the three height data sets. Nevertheless, we note that these results remain dependent on the quality and the density of data used, particularly the levelled GPS data of which the present configuration is incomplete from the point of view of density and quality.

Finally, the obtained results were satisfactory, so in the near future, new solutions will be proposed and the additional comparisons should be made for a complete error assessment of the our local geoid model. These will include new gravity data, topographic information, and the new data of GPS/Levelling in order to reach an acceptable accuracy on all the country.

Nowadays, it does not seem to be possible to replace precise levelling by GPS/Levelling, at least not on short distances, but with refined observation and data analysis techniques, it will certainly be appropriate in the future.