Irregular Grid Interpolation using Radial Basis Function for Large Cylindrical Volume

Irregular grid interpolation is one of the numerical functions that often used to approximate value on an arbitrary location in the area closed by non-regular grid pivot points. In this paper, we propose a method for achieving efﬁcient computation time of radial basis function-based non-regular grid interpolation on a cylindrical coordinate. Our method consists of two stages. The ﬁrst stage is the computation of weights from solving linear RBF systems constructed by known pivot points. We divide the volume into many subvolumes. In the second stage, interpolation on an arbitrary point could be done using weights calculated on the ﬁrst stage. At ﬁrst, we ﬁnd the nearest point with the query point by structuring pivot points in a K-D tree structure. After that, using the closest pivot point, we could compute the interpolated value with RBF functions. We present the performance of our method based on computation time on two stages and its precision by calculating the mean square error between the interpolated values and analytic functions. Based on the performance evaluation, our method is acceptable.


Introduction
Interpolation is a numerical method of estimating value on an in-between data point among known data points. Methods of interpolation can be categorized based on how the known data points arrangement. If the arrangement of known data points is scattered, then we need the interpolation method on an irreg-18 Jurnal Ilmu Komputer dan Informasi (Journal of Computer Science and Information), volume 13, issue 1, February 2020 ular grid. Distortion on known points arranged in a regular grid at a cylindrical volume could cause the movement of pivot points, and the arrangement of pivot points become irregular. There are several methods for multivariate irregular interpolation: nearest-neighbor interpolation, inverse distance weighing [1], kriging [2], and Radial Basis Function (RBF) based interpolation [3]. The most widely used of irregular grid interpolation is RBF-based interpolation. Depth on RBF definition and application can be found in [3]. RBF-based interpolation is a viable choice for approximating interpolation on irregular geometry because it has high order accuracy and computational stability for large numbers of scattered data points even in high dimension [4]. Previous studies show many implementation of RBF-based interpolation can be found in various fields: reconstruction from scattered data in computer vision and graphics [5], [6] and [7], mesh deformation on electro-magnetic problem [8], and signal localization approximation [9].
Cylindrical-coordinate usually used in numerical/analytic approaches for many physics problems. One of the problems is approximation charge distribution as a function at a hollow cylindrical detector such as Time Projection Chamber (TPC) at A Large Ion Collider Experiment (ALICE) CERN [10]. Assume we have regular pivot points at the cylindrical volume when it is distorted due to space charge [11], then the pivot points become irregular grid. Hence, we need irregular interpolation methods for function approximation at the cylindrical volume. Interpolation in cylindrical volume methods and application can also be found in thermal and mechanical response analysis [12], near-field milimeter-wave cylindrical scanning [13], and other problems, such as [14], [15], [16]. Most of them provide interpolation methods for regular grid interpolation. While [15] presents a paper on the evaluation of seven spectral methods for interpolation, differ with [15], our main focus here is how to adapt RBF-based interpolation on irregular grid pivot points.
In this paper, we present a method for implementing RBF-based interpolation on largely known data points set specifically in 3D cylindrical volume. This problem arises when large regular data points are deformed due to certain of the case, or the known points are scattered in the volume. At the initialization phase, we first compute RBF weights for all pivot points and its neighbor pivot points. We also constructed a KD-tree structure from all pivot points [17]. The KD-tree is used to find the nearest pivot point of the point we want to interpolate. Once the nearest point is found, we could apply interpolation to approximate the function. This paper is organized as follows: in Section 2 we present the definition of cylindrical coordinates, KD-tree, and RBF-based interpolation, in Section 3 we present our proposed RBF-based interpolation on an irregular grid in a cylindrical algorithm. Implementation and discussion on the performance of our proposed algorithm are given in Section 4. The paper is concluded with a conclusion given in Section 5. Using KD-tree structure and RBF weights from the initialization phase, we could approximate the interpolation function in an arbitrary point P using Equation 2. First, we find the nearest pivot point P from the query point P. After that, invoke w from w (calculated RBF-weights from Algorithm 2). Last, we calculate P.

Definitions
In this section, we present the definition of cylindrical coordinates, KD-tree structure, and RBFinterpolation.

Cylindrical Coordinates
Our proposed method is to work in a cylindrical coordinate. Cylindrical-coordinate is a generalization of two-dimensional polar coordinates to three dimensions by superposing a height (z) axis. It provides a natural extension of polar coordinates to three-dimension. In the cylindrical coordinate system, a point in space is represented by an ordered triple (r, θ, z) where (r, θ) represents the polar coordinates of the point's projection in the xy-plane and z represents the point's projection onto the z-axis. Suppose there are two points in a cylindrical volume p 0 = (r 0 , θ 0 , z 0 ) and p 1 = (r 1 , θ 1 , z 1 ). Distance between p 0 and p 1 is define as follows [18]: We use distance definition in a cylindrical coordinate given in Equation 1 extensively in RBF-based interpolation.

KD-Tree Structure for Finding Nearest Point
KD-tree is a data structure used for organizing points in K-dimensional space to perform a nearest neighbor search. KD-Tree partitions each point (called a node) in the dataset into axisaligned in a hierarchical manner. Once KD-tree constructed, finding the nearest pivot point has complexity O(logn). More on how to construct KD-tree and how to use KD-tree for finding the nearest pivot point based on a query point can be found in [19]. Several parallel solution have been done in [20] and [21] to reduce the time complexity.
Pivot points (list of known points) in threedimensional space can be arrange in regular grids or in irregular grids as given in Figure 1 and Figure  2. Finding the nearest point in the regular grid is straight forward, however in an irregular grid is not so straight forward. One of the methods is to construct a KD-tree structure as a look-up tree.

Radial Basis Function-based Interpolation
One of methods for multivariate interpolation on scattered data is RBF-based interpolation. Radial basis function ϕ(|p 0 −p|) is a real function that based only by distance from a pivot point. Some of RBF functions shown in Table 1 [3].
Assume we have a set of pivot points {p 0 , · · · , p N } where p = (r, θ, z) is a coordinate in cylindrical volume. RBF-based interpolation use the pivot points, distance function defined in Equation 1 and a selected RBF function. The interpolation function is defined as following equation: where w is computed as such the interpolation function f (p) interpolate all pivot points. Hence, we have to solve linear system in Equation 3.
There are several ways for solving Equation 3 such as using Singular Value Decomposition technique such as in [22].

Irregular Grid Interpolation for Large Data Points on Cylindrical Coordinate
In this section we present our proposed method for implementation RBF-based interpolation in large cylindrical volume. We divide our method into two phases namely intitialization and interpolation. In the first phase we construct a KD-tree structure for looking-up nearest pivot point and calculated weights of RBF-interpolation system as in Equation 3 for all pivot points (where size of neighborhood we choose as a parameter). In the following phase, we can approximate any arbitrary point in the cylindrical volume by using Equation 2 after a KD-tree is constructed and weights is computed.

Initialization Phase
In the initialization phase, we compute the RBFbased interpolation weights in Equation 3 by traversing all pivot points in the cylindrical volume. At first the neighborhood size of pivot points is chosen (s r , s θ , s z ). Then, r 0 which is used in the RBF function at Table 1 is calculated by GetRadius0RBF() function. r 0 is the average distance of all pivot 20 Jurnal Ilmu Komputer dan Informasi (Journal of Computer Science and Information), volume 13, issue 1, February 2020 . v is a vector of the interpolation function of pivot points. Hence, we have a linear system A w = v with w is unknown.

Interpolation Phase
Using KD-tree structure and RBF weights from initialization phase, we could approximate the interpolation function in an arbitrary point P = (r, θ, z) using Equation 2. First, we find nearest pivot point p i,j,k from the query point P . After that, we invoke w i,j,k from w (calculated RBF-weights from Algorithm 2). Last, we calculate f (P ) = srs θ sz n=0 w i,j,k [n]ϕ( p n − P ).

Experiment and Performance Analysis
We implemented our method in C++ and compiled with GNU GCC. The experiment conducted at a CPU based machine whose specification as follows RAM 16 GB and a core of Intel-Xeon E5-2695 2.10 GHz. Purposes of the experiments are to evaluate its error performance and its computation time for each phase (initialization and interpolation).

Accuracy
We define the accuracy by defining residue functions. Supposed we have an analytic function f analytic (r, θ, z) and we choose a set of irregular pivot points {p 0 , · · · , p N } in the cylindrical volume. Then we choose a set of arbitrary points x 0 , · · · , x M . We can define residue (error) as the following equation: Let we choose an analytic function as follows: where r 0 and r 1 are inner radius and outer radius of the hollow cylinder, respectively. We choose this function because it represents many conditions, such as polynomial, sinusoidal, and exponential function. Hence, the performance will be more convincing if the proposed method works in a variety of functions. Function map of Equation 6 in xy-plane and rzplane is given by Figure 3. We conduct our experiments in the six following steps.
1) Define pivot points in a regular grid with size varied from (17,17,18) to (129,129,144) where consecutively the grid size in each direction (n r , n z , n θ ) 2) Translate pivot points in random distance for each direction. The result is the pivot points would be arranged in irregular grid. 3) Set values for all pivot points using Equation 6 4) Do initialization phase for constructing KDtree structure of pivot points and look-up table for RBF weights using Algorithm 2.
5) Generate random points in the volume and approximate its value using RBF-based interpolation as in Equation 2. 6) Compute for each random point for its error as in Equation 4 and 5. [a] [b] After the initialization phase, we generate some arbitrary points (in our case 1×10 6 ). For each point, we calculate the value from test function and the value from the proposed method. The difference between both values can be regarded as a residue/error. One of the residue maps is given by Figure 4. [a] [b] Fig. 4. (a) xy-plane projection of residue between the test function and approximation by interpolation, (b) rz-plane projection of residue between the test function and approximation by interpolation Table 2 shows the results of our experiment for the proposed method. In Table 2, we compare error produced from an increasing number of grid size and present mean, max, and deviation of residues from all test points. We also include RBF-function type such given by Table 1. Here we can deduce that relative errors are smaller when the number of points is increasing. However, we should make note that our test function is smooth, and our pivot points are uniformly distributed.

Computation Time
We conduct experiments to measure the CPU time of the initialization phase by running a different number of pivot points, as shown in Table 3. Initialization has two steps: (1) constructing KD-tree structure and (2) solving RBF-weights as given in Algorithm 1 and 2. Table 3 shows that RBF-weights calculation has far more complexity than KD-tree initialization. Complexity for constructing KD-tree is N × O(N logN ) (average complexity), where N is a number of pivot points. While the complexity of computing RBF-weight for each pivot point is depending on the size of the neighborhood. We can formulate a number of FLOPS: N × O(m 3 ) where m is the size of the neighborhood of pivot point.
Computation time for the interpolation phase is the sum of two steps: (1) finding the nearest pivot point from the query point by using a KD-tree structure, (2) applying the pre-computed RBF weight to Equation 3. The complexity for finding the nearest pivot point on the KD-tree structure is O(logN ). While applying Equation 3 has constant computation time. Since step (1) and (2) always the same, the computation time for the interpolation phase is relatively short and constant of every query point. From our measurement, the interpolation phase averagely requires 6.17 × 10 −3 second, which is very small.

Conclusion
In this work, we present a method for implementing RBF-based interpolation for scattered data in cylindrical coordinates. We use a KD-tree structure as a helper for finding the nearest pivot point from a query point. To handle a large volume of data, we regard each pivot point and its neighbor as a set of scattered data, which is approximate by the RBF-based irregular interpolation defined on a pivot point and its neighbor. Our experiment shows that the average of residues of RBF-based interpolation for a cylindrical coordinate is in the acceptable range (has relative error below 1 × 10 −3 ). In the computation time aspect, our proposed method has relatively constant computation time for interpolation. However, in initialization time, it has complexity N × O(m 3 ), computation time extends in polynomial time. Our further works would be two things: (1) how to handle non-uniform distributed scattered data as pivot points and (2) how to accelerate the initialization phase. In the current form, we assume that the scattered data is distributed uniformly in 22 Jurnal Ilmu Komputer dan Informasi (Journal of Computer Science and Information), volume 13, issue 1, February 2020  a large volume. Some methods can be considered, such as RBF-based with a spectral method. There are many options to accelerate the initialization phase in our proposed method: multi-core, multi-CPU, or even using an accelerator card such as GPU.