**1. Introduction**

Scattered data interpolation and approximation are still active research topics in computer-aided design (CAD) and geometric modeling [1–9]. This is because engineers and scientists often face the problem of how to produce smooth curves and surfaces for the raw data obtained from experiments or observations. This is where scattered data interpolation can be used to assist them. To construct smooth curves and surfaces, some mathematical formulations are required. This can be achieved using functions which are well-established, such as the Bézier, B-spline, and radial basis functions (RBFs). All these methods are guaranteed to produce smooth curves and surfaces.

The formulation problem in scattered data interpolation can be described as follows:

Given functional data

$$(x\_{i\prime}y\_{i\prime}z\_i)\_{\prime\prime} \text{ } i = \mathbf{1}, \mathbf{2}, \dots, N$$

construct a smooth *C*<sup>1</sup> surface *z* = *F*(*x*, *y*) such that

$$z\_i = F(\mathbf{x}\_{i\prime} y\_i)\_{\prime} \quad i = 1, 2, \dots, N$$

To solve the above problem, there are many methods that can be used, such as meshless methods (e.g., radial basis functions (RBFs) and many types of Shepard's families). However, some meshless schemes are global. Fasshauer [10] gave details on many meshless methods to solve the problems arising in scattered data interpolation and approximation, as well as partial differential equations. Beyond that, another approach that can be used to solve the problem is the triangulation of the given data points. Then, the Bézier or spline triangular can be used to construct a piecewise smooth surface with some degree of smoothness, such as *C*<sup>1</sup> or *C*2. The Shepard triangular can also be used to produce a continuous surface from irregular scattered data. For instance, Cavoretto et al. [6], Dell'Accio and Di Tommaso [11], and Dell'Accio et al. [12,13] have discussed the application of the Shepard triangular for surface reconstruction. However, their schemes require more computation time in order to produce the interpolated surfaces.

Crivellaro et al. [14] applied RBFs to reconstruct 3D scattered data via new algorithms, which involves an adaptive multi-level interpolation approach based on implicit surface representation. The least squares approximation is used to remove the noise that appears in the scattered data. Chen and Cao [15] employed neural network operators of a logistic function through translations and dilation. Meanwhile, Bracco et al. [2] considered scattered data fitting using hierarchical splines where the local solutions are represented in variable degrees of the polynomial spline. Zhou and Li [16] studied scattered noise data by extending the weighted least squares method via triangulating the data points. Zhou and Li [17] discussed the scattered data interpolation for noisy data by using bivariate splines defined on triangulation. Qian et al. [18] also considered scattered data interpolation by using a new recursive algorithm based on the non-tensor product of bivariate splines. Liu [19] constructed local multilevel scattered data interpolation by proposing a new idea (i.e., nested scattered data sets), and scaled the compactly supported RBFs. Borne and Wende [3] also considered the meshless scheme based on definite RBFs for scattered data interpolation. In their study, they applied the domain decomposition methods to produce a symmetric-saddle point system. Joldes et al. [20] modified the moving least squares (MLS) methods by integrating the polynomial bases to solve the scattered data interpolation problem. Brodlie et al. [5] discussed the constrained surface interpolation by using the Shepard interpolant. The solution to the problem is obtained by solving some optimization. Lai and Meile [21] discussed scattered data interpolation by using nonnegative bivariate triangular splines to preserve the shape of the scattered data. Schumaker and Speleers [22] considered the nonnegativity preservation of scattered data by using macro-element spline spaces including Clough–Tocher macro-elements. Furthermore, they also give general results for range-restricted interpolation. Karim et al. [23] discussed the spatial interpolation for rainfall data by employing cubic Bézier triangular patches to interpolate the scattered data. Karim et al. [24] have constructed a new type of cubic Bézier-like triangular patches for scattered data interpolation. Karim and Saaban [25] constructed the terrain data using cubic Ball triangular patches [23]. In this study, they show that the scattered data interpolation scheme by Said and Rahmat [26] is not *C*<sup>1</sup> everywhere. Thus, a new condition for *C*<sup>1</sup> continuity is derived. The final surface is *C*<sup>1</sup> and provides a smooth surface. Feng and Zhang [27] proposed piecewise bivariate Hermite interpolations based on triangulation. They applied the scheme for large scattered data sets to produce high-accuracy surface reconstruction. Sun et al. [28] constructed bivariate rational interpolation defined on a triangular domain for scattered data lying on a parallel line. They only considered a few data sets, and it was not tested for large data sets. By using a rational spline, the computation time increases. Bozzini et al. [4] proposed a polyharmonic spline to approximate the noisy scattered data.

The main motivation of the present study is described in the following paragraphs. In triangulationbased approaches to scattered data interpolation, cubic Bézier triangular or quintic Bézier triangular patches are the common methods. The quartic Bézier triangular has received less attention due to the need to solve optimization problems in order to calculate the Bézier ordinates. This approach increases the computation time. There are four steps in constructing a surface using a triangulation method: (a) triangulate the domain by using Delaunay triangulation; (b) specify the first partial derivative at the data points (sites); (c) assign the control points or ordinates for each triangular patch; and finally (d) the surface is constructed via a convex combination scheme. Goodman and Said [29] constructed a suitable *C*<sup>1</sup> triangular interpolant for scattered data interpolation using a convex combination scheme between three local schemes. Their work is different from that of Foley and Opitz [30]. However, both studies developed a *C*<sup>1</sup> cubic triangular convex combination scheme. Said and Rahmat [26] constructed a scattered data surface using cubic Ball triangular patches [31,32] with the same approach as in Goodman and Said [29]. Based on the numerical results, their scheme gave the same results as cubic Bézier triangular patches. The main advantages of cubic Ball triangular patches are that the required computation is 7% less when compared with the work of Goodman and Said [29]. This is what has been claimed by References [26,29]. However, in the work of Karim and Saaban [25], it was proved that Said and Rahmat [26] is not *C*<sup>1</sup> continuous everywhere, and Karim and Saaban [25] found that the [26] scheme produced the same surface for scattered data interpolation when the inner coefficient was calculated by using Reference [29]. Hussain and Hussain [33] proposed the rational cubic Bézier triangular for positivity-preserving scattered data interpolation. They claimed that their proposed scheme is *C*<sup>1</sup> positive everywhere. However, from their results, it is possible that their scheme may not be positive everywhere. Chan and Ong [7] considered range-restricted interpolation using a cubic Bézier triangular comprising three local schemes. All the schemes were implemented by estimating the partial derivatives at the respective knots using the method proposed by Goodman et al. [34].

Other than the use of cubic Ball and cubic Bézier triangular patches for scattered data interpolation, there are some studies that have utilized quartic Bézier triangular and rational quartic Bézier triangular patches for scattered data interpolation. For instance, Saaban et al. [35] constructed *C*<sup>1</sup> (or *G*1) scattered data interpolation based on the quartic Bézier triangular. Piah et al. [36] considered *C*<sup>1</sup> range-restricted positivity-preserving scattered data interpolation by using the quartic Bézier triangular. They employed an optimization method (i.e., the minimized sum of squares) to calculate the inner Bézier points proposed in Saaban et al. [35]. Hussain et al. [37] extended this idea to construct convexity-preserving scattered data interpolation. Hussain et al. [38] constructed a new scattered data interpolation scheme by using the rational quartic Bézier triangular. They applied it to positivity-preserving interpolation. However, to achieve *C*<sup>1</sup> continuity, we still need to solve some optimization problems. This is the main weakness of quartic Bézier triangular patches when applied to scattered data interpolation. Some good surveys on scattered data interpolation can be found in [39–43].

The present study aims to answer the following research questions:

a. Can we construct a scattered data interpolation scheme by using quartic triangular patches but without an optimization method?

b. How can we produce a *C*<sup>1</sup> surface (everywhere)?

c. Is the proposed scheme better than some existing schemes in terms of CPU time, coefficient of determination (R2), and maximum error?

To answer these research questions, we will use the quartic triangular basis initiated by Zhu and Han [44]. The main advantage of using this quartic basis is that it only requires ten control points to construct one triangular patch. This is the same as the number of control points in the cubic Bézier triangular patch. Thus, in order to construct *C*<sup>1</sup> scattered data interpolation using the quartic spline triangular, we can employ the Foley and Opitz [30] cubic precision scheme to calculate the inner ordinates. With this, the optimization problem required in a quartic triangular basis will be avoided. Hence, this will show that the proposed scheme is local. Furthermore, the proposed scheme is different from the works of Lai and Meile [21] and Schumaker and Spellers [22], even though all schemes required triangulation of the given data in the first step.

Some contributions from the present study are described below:

1. The proposed scattered data interpolation scheme produces a *C*<sup>1</sup> surface without any optimization method like Piah et al. [36], Saaban et al. [35] and Hussain et al. [37,38].

2. The proposed scheme is local; meanwhile, the schemes presented in Piah et al. [36], Saaban et al. [35] and Hussain et al. [37,38] are global.

3. Based on the CPU time needed to construct the surface, the proposed scheme is faster than quartic Bézier triangular patches. Thus, the reconstruction of scattered surfaces from large data sets can be performed in less time.

4. Furthermore, the proposed positivity-preserving scattered data interpolation is capable of producing a better interpolated surface than quartic Bézier triangular patches. This lies in contrast to scattered data schemes by Ali et al. [1], Draman et al. [9] and Karim et al. [24], which are not positivity-preserving interpolations.

This paper is organized as follows: In Section 2 we give a review of the triangular basis initiated by Zhu and Han [44], and the derivation of the quartic triangular basis with ten control points. Some graphical results are presented, as well as the construction of a local scheme comprising convex combination between three local schemes. The numerical results and the discussion are given in Section 3 with various numerical and graphical results, including a comparison with some existing schemes. Error analysis is also investigated in this section. The construction of the positive scattered data interpolant is discussed in Section 4. Meanwhile, numerical results for positivity-preserving scattered data interpolation are shown in Section 5. Conclusions and future recommendations are given in the final section.
