Radial basis function interpolation of non-matching grids surfaces for volume calculation.

Prada, Marcela ; Teusdea, Alin Cristian ; Fetea, Ioana 等

1. INTRODUCTION

Volume calculation starting from two horizontal walls with grid (i.e. planar projections) matching is the classical well known method. The non-matching grids are an issue for the surveying software, mostly when the scanned point clouds are huge. In these cases, the coverage surfaces imply a huge number of manually trimming of undesired vertex (points) link that passes through the volume.

In this paper, is presented the volume calculation starting from two vertical diabaz open mine walls with non-matching grids (i.e. the vertical planar projections). The vertical walls have irregular shapes and are 3D scanned in different moments. From the 3D scanned point clouds are extracted the 3D contours of the vertical walls within are done the radial basis function (RBF) 3D interpolations. The results are two 3D surfaces with higher resolutions but with non-matching grids (connections) or surrounding contours.

The next step is to build the coverage surfaces (the upper, the lower, the right and the left part) in order to calculate the total volume between the initial walls. The coverage surfaces are built considering the conjugate (connecting) parts of the walls contours and by RBF 3D interpolation. The end of this stage is a 3D hollow object with a surrounding surface, built up with points, that has same resolution over the height.

In the final stage, the 3D object is sliced over the discrete values domain of the height generating a set of sliced surfaces. The volume is calculated by multiplying the sliced area and the height resolution over the entire set of slices.

2. METHODS AND SAMPLES

2.1 Samples

As a case study, for this paper the diabaz open mine from Mures valley, Romania, was included. First, there were scanned two certain surface walls of the same mine site (S1 at initial time (figure 1) and S2 at final time) using the smart station GPS R3 Total Station Trimble 5600 DR 200+, in direct reflex mode. The mean error for point positioning is 39.82mm (dx=11mm, dy=13mm, dz=36mm).

[FIGURE 1 OMITTED]

The point clouds post processing of the two walls consists in trimming the undesired point connections and then in capturing the contours (figure 1--the big dotted lines for left and right parts, boxed line for the upper part and diamonded line for the lower part).

In order to calculate the volume of the 3D hollow object bordered by these two vertical surfaces, one have to interpolate them over a fixed vertical resolution. This method enables the 3D hollow object to be regularly sliced and to provide a continuously sliced contour in order to get the area of the slice.

2.2 Radial basis function interpolation

Radial basis function (RBF) interpolation consists in finding the coefficients, [lambda] = ([[lambda].sub.1], ..., [[lambda].sub.N]), for a base of radial functions and the coefficients, c = ([c.sub.1], ..., [c.sub.1]), for a set of fitting polynomial, p = {[p.sub.1], ..., [p.sub.1]}, so that this interpolation function s(x) defined below (Carr et al., 2003; Boer et al., 2007)

s(x)=p(x)+[N.summation over (i=1)][[lambda].sub.i] x [phi]([absolute value of x-[x.sub.i]), x [member of] [R.sup.N] (1)

has to pass through the values of definition (Carr et al., 2001)

s([x.sub.i]) = [y.sub.i], i = [bar.1N] and [N.summation over (j=1)] [[lambda].sub.j] x p([x.sub.j]) = 0, (2)

where ([x.sub.i]; [y.sub.i]) are the coordinates of N known points.

The thin plate radial function, [phi](r) = [r.sup.2] x ln(r), was chosen for the studied case. These conditions, under the matrix form, can be written the following form (Carr et al., 2003)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (3)

where we have: [R.sub.i,j] = [phi]([absolute value of [x.sub.i] - [x.sub.j]]), [P.sub.t,l] = [p.sub.l]([x.sub.i]), [Y.sub.i] = [y.sub.i], i, j = [bar.1, n], l = [bar.1,m]. The generated equations system has the solution given by

c = [[([P.sup.T] x [R.sup.-1] x P).sup.-1]] x ([P.sup.T] x [R.sup.-1] x Y), [lambda] = ([R.sup.-1] x Y) - ([R.sup.-1] x P) x [[([P.sup.T] x [R.sup.-1] x P).sup.-1] x ([P.sup.T] x [R.sup.-1] x Y)]. (4)

3. RESULTS

One of the initial conditions--small roughness dynamic range for the walls--of the data, involve neglecting the fitting polynomial. The 3D RBF interpolation within a certain contour for the S1 wall is presented in figure 2--and the S2 wall 3D RBF interpolation within its contour is done in the same way.

In order to generate the four coverage surfaces one must connect the upper S1 with the upper S2 contour parts and provide the upper coverage surface contour and forth for the lower, right and left ones. The next step is to make the 3D RBF interpolation within the contours to generate the four coverage surfaces (figure 3) at 0.1m resolution.

The final result of the six RBF 3D interpolated surfaces is a 3D hollow object (figure 4) that encapsulates the extracted diabaz volume between the considered moments, in local 3D coordinates, with 0.1m resolution for all the coordinates.

[FIGURE 2 OMITTED]

[FIGURE 3 OMITTED]

[FIGURE 4 OMITTED]

Slicing the 3D object at a fixed height provides a set of 2D points which consists of a polygonal surface. The polygonal area equation is given by (Cha Zhang & Tsuhan Chen, 2001)

S = 1/2 [[summation over (i,j)][([x.sub.i] x [y.sub.i+1] - [y.sub.i] x [x.sub.i+1])], i, j = [bar.1,M]. (5)

where M is the number of the polygonal surface points.

In table 1 are presented the results of diabaz mine extracted volume. There are showed three ways to calculate the same volume with the slice method (i.e. the Simpson rule). The numerical results are very close and they generate a very low standard deviation and relative error.

4. CONCLUSION

This work presents the numerical algorithm that provides the volume calculation starting with two vertical 3D scanned walls with non-matching vertical planar projections (grids). The basic methods used are the 3D RBF interpolation within the surface contour and the slice method of volume calculation (i.e. the Simpson rule). All these are applied on data taken from a diabaz open mine from Mures valley, Romania. The numerical results provided with the presented algorithm are accurate regarding that the standard deviation of the calculated volume is smaller than the theoretical one ([+ or -] 0.173205)([m.sup.3]).

The presented issue can not be solved with the actual surveying software in real time and in an economical efficient way, fact that qualify the presented method as novelty.

5. REFERENCES

Carr, J.C.; Beatson, R.K.; McCallum, B.C.; Fright, W.R.; McLennan, T. J. & Mitchell, T. J. (2003). Smooth surface reconstruction from noisy range data, ACM GRAPHITE 2003, pp. 119-126, ISBN 1-58113-578-5, Melbourne, Australia, February 2003, ACM, NY USA

Cha Zhang & Tsuhan Chen, (2001). Efficient feature extraction for 2D/3D objects in mesh representation, IEEE Conference on Image Processing ICIP 2001, Vol. 3, pp. 935-938, ISBN 0-7803-6725-1, Thessaloniki, Greece, October 2001, IEEE

Boer, A. de; Schoot, M.S. van der & Bijl, H. (2007). Mesh deformation based on Radial Basis Function Interpolation, Computers & Structure. Fourth MIT Conference on Computational Fluid and Solid Mechanics, Vol. 85, Issues 11-14, June-July 2007, pp.784-795, ISSN 0045-7949

Sirakow, N.M.; Iwanowski, M.; Hack, D.R. & Feves, M.L. (2003). Morphological approach to volume calculation of complex 3D geological objects, APCOM 2003: 31st International Symposium on Application of Computers and Operations Research in the Mineral Industries, pp.329-334, ISBN 1-919783-46-6, Cape Town, SA, May 2003, South African Institute of Mining and Metallurgy, Johannesburg

Tab. 1. Volume calculation results

Calculus formulas (Sirakov et al., 2003), k = [bar.1,N]

V=[[summation over (k)][S.sub.k]]x[DELTA]z             2262.503

V=[[summation over (k)][1/2([S.sub.k]+[S.sub.k+1])
x([z.sub.k+1] - [z.sub.k])]+[([S.sub.N]+[S.sub.N-1])
x([z.sub.N]-[z.sub.N-1])]                              2262.006

V=[[summation over (k)][([S.sub.k]+[S.sub.k+1])x 1/2
([z.sub.k+1]-[z.sub.k])]]                              2262.132

V=(2622.213667 [+ or -] 0.149169)([m.sup.3])

[[epsilon].sub.V] = 0.0066%