University of Huddersfield Repository Curvature based sampling of curves and surfaces Computer Aided Geometric Design Curvature based sampling of curves and surfaces ✩

Eﬃcient sampling methods enable the reconstruction of a generic surface with a limited amount of points. The reconstructed surface can therefore be used for inspection purpose. In this paper a sampling method that enables the reconstruction of a curve or surface is proposed. The input of the proposed algorithm is the number of required samples. The method takes into account two factors: the regularity of the sampling and the complexity of the object. A higher density of samples is assigned where there are some signiﬁcant features, described by the curvature. The analysed curves and surfaces are described through the B-splines spaces. The sampling of surfaces generated by two or more curves is also discussed.


Introduction
Acquiring points from a freeform surface is usually a long process. Fast and reliable inspection of the product has become an important factor due to the market requirements, such as low volume and high customised products (Mackelprang et al., 2015). An efficient sampling method enables the acquisition of small amount of point to reconstruct the surface and verify if the product is within tolerance.
Although triangular and quadrilateral meshes can represent a general freeform surface, it is not possible to represent the class of surfaces generated by the Cartesian product of two or more curves. In this paper non-uniform rational B-splines (NURBS, Piegl and Tiller, 1997) curves and surfaces are used to describe the shape of an object. NURBS are the standard surface representation in Computer Aided Geometric Design (CAGD). If only the point cloud is available there are some techniques to reconstruct the surface using interpolation (Piegl and Tiller, 1997) or approximation methods (Eck and Hoppe, 1996). When the nominal shape is known, using its NURBS representation has several advantages also in finite element method simulation (Cottrell et al., 2009). Sampling of parametric curve or surface is usually carried out by adding sampled points where critical points are found. Hernández-Mederos and Estrada-Sarlabous (2003) proposed a sampling method using the re-parametrisation of the curve based on a mixture between the arc length, or uniform, parametrisation and the value of the bending energy of the curve section. The authors showed that with the parametrisation it is possible to add samples points where there are abrupt changes of the curve.
A seminal work on sampling of freeform surfaces can be found in the paper of Elkott et al. (2002). The authors have developed an automatic sampling plan for freeform surfaces. The developed algorithm selects the sampling points according ✩ This paper has been recommended for acceptance by Konrad Polthier. to the complexity of the NURBS surface, which is represented by the curvature changes and the size of the patches. Starting from the size of each patch of the NURBS surface and the mean of the Gaussian curvature, the patches are ranked according to a linear combination of the two computed values. The number of sampled points are distributed among the surface according to the rank.
Another sampling method can be found in Obeidat and Raman (2009), where it is proposed to sample the points with higher curvature for inspection purpose. The algorithm take as input a NURBS surface, it is divided into patches based on the knots vectors and ranked based on their size. The points are sampled according to the higher values of the surface curvature.
A sampling method that is not based on the curvature of the surface can be found in Yu et al. (2013). The authors developed an adaptive sampling method based on the distance between reconstructed geometry and the CAGD model.
A review of intelligent sampling techniques can be found in Wang et al. (2012). Among the reviewed techniques, two adaptive sampling methods have been analysed: triangle and rectangle patch adaptive subdivision. Both the algorithm are based on the iterative subdivision of the of the surface domain in triangle (rectangle) until the reconstruction error is smaller than a defined threshold.
In the previous papers sampling methods for complex curves and freeform surfaces have been analysed. Most of them are based on the incremental reconstruction of the surface, i.e. samples are incrementally added until the reconstruction error is lower than a pre-defined threshold. If the number of sample points it is set in advance, instead of the threshold, it is possible to use a non greedy algorithm to select the samples.
A class of freeform surfaces that have not been investigated in the aforementioned papers regards the surfaces generated from two or more curves such as ruled surface, skinned surfaces, etc. (Piegl and Tiller, 1997). A sampling method based on this curves, called generatrices, is proposed and the reconstruction performance are compared to the other methods.
This paper is construct as follows: in Section 2 a sampling method for NURBS curves is presented, in Section 3 a sampling method for freeform and generated surfaces is developed and Section 4 summarise the paper.

Curve sampling
The methods proposed in this paper is based on NURBS curves and surfaces, which can represent both the nominal geometry of the CAGD model and the reconstruction of a general freeform surface. The construction of the B-splines and NURBS curves is now introduced. Given a non decreasing sequence of knots A d-th degree B-splines curve is defined as where p i is the set of control points. A d-th degree NURBS curve as where w i is the weight of the i-th B-splines basis function. The parametrisation proposed in this paper is similar to the one analysed in Hernández-Mederos and Estrada-Sarlabous (2003). The authors have proposed to re-parametrise the curve according the arc length and the square of the curvature. The arc length parametrisation of a generic parametric curve can be computed as (do Carmo, 1976) where • is the l 2 norm. This is also called uniform parametrisation because if a curve is re-parametrised according to l(t) its speed is unitary. An example of a B-splines curve and a uniform arc length sampling is shown in Fig. 1a. The parameter value as a function of the computed parametrisation is reported in Fig. 1b, this function is used to perform the points sampling. It is possible to observe that the greater difference between the uniform and the initial parametrisation are located near the boundaries, i.e. the distances between the sampled points is smaller with arc length sampling.
The curvature is a geometric quantity that measures the rate of change of the unit tangent vector of the curve, it is usually computed to describe the complexity of a curve (Hernández-Mederos and Estrada-Sarlabous, 2003). The curvature of a parametric curve can be computed as (do Carmo, 1976) Considering that it is an intrinsic value of the curve, a curvature parametrisation can be computed as where the integral is computed according to the infinitesimal arc length dl(u). Compared to the parametrisation proposed in Hernández-Mederos and Estrada-Sarlabous (2003), in this paper the square of the curvature is not used and the integral is performed on the curve, i.e. weighting the integrand function with r (u) . It should be noted that otherwise the value of k p (t) depends on the curve parameter t and not only on the geometry of the curve.
The cumulative integral in Equation (7) is shown in Fig. 2b. The curvature varies where there are high changes of the tangent vector, but where the curve is flat it is almost constant, such as in the end of the graph. This leads to an high density of sampled points in some parts of the curves and where the curve is flat there are few sampled points (see Fig. 2a).
To overcome the lack of points in the flat zones, it is possible to compute a mixed parametrisation, similar to the one in Hernández-Mederos and Estrada-Sarlabous (2003), as In this paper it is assigned an equal weight both to the arc length and the curvature parametrisation. In Fig. 3b the parametrisation computed with Equation (8) is shown. Since it is a mean of the two quantities in the flat zones the slope of the cumulative integral is higher compared to the curvature one; where there are big changes of the tangent vector the slope is higher than the arc length one. This parametrisation is a compromised between the uniform sampling and a sampling based only on the complexity of the curve. A 50 points uniform sampling is shown in Fig. 3a, the density of the sampled points is higher where the curvature is high but there are points also in the flat zones.

Reconstruction results
In this Section the reconstruction performances based on different sampling methods are evaluated. Three curves are analysed: the curve used in the previous paragraph, a curve representing a simulation of a milled profile and a profile of a real surface. The root mean square error (RMSE) is used as a performance indicator where n is the number of the points, d is the dimension of the B-splines curve or surface, r ij (s) is the true value and r ij (s) is estimated point with the reconstruction model. The reconstruction algorithms used in this paper are the piecewise cubic splines interpolation implemented in MATLAB ® (2014) through the cubicinterp function and the Akima splines interpolation (Akima, 1970;Henning, 2012). The piecewise cubic spline interpolation is a globally C 2 cubic function, the second derivative is usually set to zero, this lead to the so called natural cubic spline. The Akima interpolation is third degree function with a globally C 1 continuity. Since there are less constraints in the Akima spline, it can easily describe curves with abrupt changes. The reconstruction is then performed with the method that allows to achieve the better performance, i.e. the Akima interpolation when the curve has some jumps while the piecewise interpolation when it is "smoother". In this analysis it is assumed that the true value of the parameter of the curve is known, i.e. the point cloud parametrisation is not taken into account.
Five sampling methods are investigated: • uniform arc length sampling, arc_length; • uniform mixed arc length and curvature sampling, al_curv; • uniform parametric sampling, uniform; • random latin hypercube sampling (LHS), lhs (Santner et al., 2003); • the method implemented in Hernández-Mederos and Estrada-Sarlabous (2003), al_curv2. The first analysed curve is the curve described in the previous Section; 30 points sampled with the mixed parametrisation are shown in Fig. 4a. The RMSEs of the analysed models with the piecewise cubic reconstruction are shown in Fig. 4b. It is possible to observe that the reconstruction with the proposed mixed parametrisation has a smaller reconstruction error.
In The values of the RMSEs are shown in Figs. 6b and 7b. The reconstruction, based on the Akima interpolation, with the proposed methods produces better reconstruction results. It should be noted that, due to the reconstruction algorithm, the parametrisation of the initial point cloud is close to the arc length parametrisation, consequently the arc_length and the uniform have similar performances.
This Section has shown that the reconstruction based on the proposed sampling method outperforms the other analysed models. The mean of the arc length and curvature parametrisations is a good compromise between the uniformity and the complexity of the curve. In the next Section the curve sampling is extended to the surface scenario.

Surface sampling
A simple sampling method based on the vector product of two sampling in u and v direction is presented. After a brief introduction of the NURBS surfaces, the proposed sample method is presented and the reconstruction is compared with some common sampling techniques. An investigation on the sampling of some surfaces generated from curves is finally performed; the surface sampling methods are compared with generatrices curves-based samplings.
A NURBS surface can be computed through the tensor product of two NURBS basis function in two orthogonal direction u and v. A NURBS surface of degree p in u direction and degree q in v direction is a piecewise rational function described by (Piegl and Tiller, 1997) where p ij are the control points, w ij is the weights of the product of the B-splines basis B i (u) and B j (v), u and v are the knots vectors of the B-splines basis in u and v direction. An example of a B-splines surface is shown in Fig. 8.
In order to present the proposed sampling method, the first step is to compute the area of a parametric surface (do Carmo, 1976) then it is possible to calculate the marginal cumulative area along the u direction as The values A • (•) are used to compute the uniform sampling along the u and v directions. The sampling on the surface can then be computed as the cross product of the two marginal samples. This is not the exact equal area sampling, but it uses the same principle of the NURBS surface construction. An example of the B-splines surface with the two marginal cumulative areas is shown in Fig. 9b, while a 10 × 10 uniform sampling along the marginal areas is reported in Fig. 9a.
As for the curve sampling a value that measures the complexity of the surface has to be chosen. An important shape property is the mean curvature, that is computed as (do Carmo, 1976) k (u, v) where k 1 (u, v) and k 2 (u, v) are the principal curvatures of the surface. The mean curvature is used instead of the Gaussian one because if one of the principal curvatures is null, the Gaussian curvature vanish to 0. When a principal curvature is 0 and if it is not possible to find the right parametrisation, a sampling model based on the Gaussian curvature cannot "see" that the curvature in one direction is changing. Similar to the curve case, the cumulative integral of the mean curvature on the surface is computed. The marginal cumulative integral along the u direction is and along v direction is These quantities can be used to derive other two sampling strategies: one is based on the curvature and the other is based on a mixed measure between the area and the curvature. The mixed marginal parametrisations are where the area and the curvature have the same importance. The mixed marginal parametrisations of the example surface are shown in Fig. 10b, while a sample of 10 × 10 points is drawn in Fig. 10a.

Reconstruction results
In this Section three freeform surfaces are analysed. The first and the second are surfaces with some, respectively, slow and rapid changes of the shape; the third surface is a structured surface. The following sampling strategies are analysed: • uniform based on the cumulative marginal areas, area; • uniform based on the mixed cumulative marginal areas and mean curvatures, area_curv; • uniform on the parameters space u-v, uniform; • latin hypercube sampling, lhs (Santner et al., 2003); • Hammersley sampling, hammersley (Wong et al., 1997); • Halton sampling, halton (Wong et al., 1997); • triangular patch sampling implemented in Shih et al. (2008), tri_patch.
The performances of the parametrisation based only on the mean curvature are not evaluated because, as in the curve scenario, it does not add points on the flat parts of the surface. The analysed adaptive sampling is the triangular patch sampling since it has achieved good reconstruction performances in Wang et al. (2012) and Shih et al. (2008). The reconstruction is performed through the piecewise cubic interpolation implemented in the fit function of MATLAB ® (2014).
As for the curve test cases the u and v parameters during the reconstruction are assumed known, i.e. the parametrisation problem in not considered. The RMSEs of the surface in Fig. 8 are shown in Fig. 11. The uniform, area and area_curv sampling have comparable reconstruction errors, which are smaller than the errors of the other sampling strategies.   12 shows the sampled points and the colour map of the reconstruction error. Reconstruction with both uniform and area have the same disadvantages of, respectively, the uniform and the arc_lenght parametrisations. With the area_curv parametrisation the errors appear almost constant along the surface. The three random sample techniques analysed, i.e. hammersley, halton and lhs, generate some cluster of points in the surface domain, so the reconstruction is poor in regions where the points' density is not high. The reconstruction with the samples of the tri_patch method has the worst performance in some parts of the surface with small curvature, bad not flats.
Since the previous surface has a simple behaviour, the performances of a more complex surface are evaluated (X3D, 2015). The plot of the surface with a sampling of 10 × 10 points with the area_curv method and the RMSEs values are shown in Fig. 13. When the surface has a complex shape both the area and the mixed parametrisation can ensure better reconstruction results. The third surface analysed is shown in Fig. 14a. It represents a structured surface. The values of the RMSE are reported in Fig. 14b. With the proposed sampling method it is possible to add more samples on the critical points of the surface and this leads to better reconstruction results.
In this Section it has been shown that a reconstruction with samples that take into account the complexity of the surface through the mean curvature leads to a smaller RMSE when the complexity of the surface increases.

Generated surfaces
After the analysis of some general freeform surfaces, the reconstruction of surfaces generated from two or more curves is analysed. These surfaces can be classified into: prism surface, ruled surface, surface of revolution, swung surface, skinned surface, swept surface and coons surface; for a better explanation regarding their construction see Piegl and Tiller (1997).
During the design phase the generatrices are known, it is therefore possible to use these information to extract the sampled points. The curves based sampling is build independently along the generatrices curve; the samples on the surface are, in most of the cases, the cross product of the two curves based samples. In this paper this sampling method is called curves.
Each type of generated surface is briefly recalled the construction and a sampling method is proposed. The reconstruction results are finally compared with the surface based sampling.

Prism surface
A prism surface is a generalisation of the cylinder, a curve r 1 (u) on the x-y plane is swept to a distance d. The surface can be represented as Since a prism surface is generated extruding a curve along an axis, the proposed sampling method along the curves is: • a sampling based on the mixed curve parametrisation along the swept curve; • a uniform sampling along the extrusion direction.
The surface sampling is then computed as the cross product of the two curves' samplings.
An example of an extruded surface is show in Fig. 15a with the curves based sampling of 30 × 3 points.

Ruled surface
Let r 1 (u) and r 2 (u) be two B-splines curves, the surface representing the linear interpolation between the two curves is called ruled surface.
Let n u and n v be the samples along u and v directions, the sampling is computed as: • sampling along the generatrices curves n u samples with the mixed sampling method; • for each pair of the samples along the u direction compute n v samples through linear interpolation.
A ruled surface with the two generating curves is show in Fig. 15b. In this Figure there is the sampling of 30 × 3 points with the curves sampling method.

Surface of revolution
A surface of revolution, see Fig. 15c, is a surface that is build through the rotation of a curve, called generatrix, along an axis. The proposed sampling is: • compute n u , uniformly spaced, samples along a circle; • compute n v sample along the generatrix with the mixed sampling method; • compute the samples on the surface as the cross product of the samples along the curves.

Swung surface
A swung surface is a generalisation of the surface of revolution. Let r 1 (u) be a curve defined on the x-z plane, called profile curve and r 2 (v) a curve on the x-y plane, called trajectory The swung surface is a surface defined as where and α is a scaling parameter. The sampling strategy can be defined similarly to the surface of revolution case: • compute n u samples along the profile curve with the mixed sampling method; • compute n v samples along the trajectory curve with the mixed sampling method; • compute the samples on the surface as the cross product of the samples along the curves.

Skinned surface
Let {r k (u)} be a set of curves, called section curves. These curves represent the cross section of the surface along the u direction, the curves along the v direction are defined through an interpolation with a desired continuity. The sampling strategy can be defined similarly to the ruled surface one: • sampling along the section curves n u samples with the mixed sampling method; • for each tuple of the samples along the u direction perform an interpolation with the same continuity of the skinned surface in v direction; • for each interpolated curve compute n v sampling through the arc length based sampling method.
In the last step the sampling method based on the arc length is preferred because it leads to a better reconstruction result.
An example of a skinned surface build from four profile curves is shown in Fig. 15e, the surface is of class C 2 along the v direction. The sample size used for the reconstruction is equal to n × 5.

Swept surface
A swept is a surface that is generated sweeping a section curve along an arbitrary trajectory. A general swept surface is defined as where r 1 (u) is the section curve, r 2 (v) is the trajectory and M(v) ∈ R 3×3 is a rotation and scaling matrix. If M (v) is the identity matrix, the swept surface can be computed as The sampling strategy can be defined as: • compute n u samples along the section curve with the mixed sampling method; • compute n v samples along the trajectory curve with the mixed sampling method; • compute the samples on the surface as the cross product of the samples along the curves.
An example of a swept surface with the profile and the trajectory and its samplings are shown in Fig. 15f. A general swept surface can also be built through a additional rotation of the section curve in order to follow the behaviour of the trajectory curve (Piegl and Tiller, 1997). Applying these transformations to the section and trajectory above the resulting swept surface is shown in Fig. 15g.

Coons surface
A coons surface is defined through four boundary curves. Let r 1,u (u) and r 2,u (u) the two boundaries curve on the u parameter space and r 1,v (v) and r 2,v (v) the two boundaries curve on the u parameter space, the bilinearly blender coons surface is defined by: where r 1 (u, v) and r 2 (u, v) are the ruled surfaces between, respectively, the curves r 1,u (u) and r 2,u (u) and the curves r 1,v (v) and r 2,v (v) and r 3 (u, v) is a bilinear tensor product surface (plane) defined as and Since a coons surface is a sum of three surface that are linear in at least one direction, possible sampling strategy is: • compute n u samples along the two generatrices curves in u direction with the mixed sampling method; • for each couple of n u samples perform an interpolation in the u-v space; • compute n v samples along the two generatrices curves in v direction with the mixed sampling method; • for each couple of n v samples perform an interpolation in the u-v space.
An example of a coons surface is shown in Fig. 15h.

Discussion of results
In this Section reconstruction performance of the proposed sampling methods are compared. Fig. 16 shows the RMSEs of the analysed surfaces. Although some sampling have different sizes along the u and v directions, on the abscissa is reported the square root of the total number of samples to maintain the scale linear. The samples based on LHS, Halton, Hammersley and the triangular patch sampling are not shown for prism, ruled and skinned surfaces because the RMSE values are too high.
From the figures is possible to conclude that the curve based and area_curv methods achieve, compared to the analysed sampling methods, equal or better reconstruction performance in most of the analysed cases. The reconstruction of the skinned surface with the curve based method is not as good as the reconstruction with the area and mixed sampling (see Fig. 16e). The sampling based on the mixed parametrisation has some problem in the reconstruction of the revolution surface, this is due to the effect of the small radius on the flat part of the surface. The area based sampling achieves good prediction results. The uniform sampling has lower performance compared to the previous mentioned methods. Performance of LHS, Halton and Hammersley samplings, since the information of the surface shape is not included in the model, are lower compared to other sampling procedures. The method based on the triangular patches achieves oscillating performances: if the number of samples is small the reconstruction is not good as the proposed methods.

Conclusion
In this paper a curve sampling method based on a mixed parametrisation between the arc length and the curvature has been proposed. If the analysed shape is a surface a simple sampling method based on the marginal cumulative integrals along a NURBS surface has been proposed. The reconstruction with the proposed methods have been compared with other common sampling strategies. It has been shown that the proposed methods achieve comparable or better reconstruction performance if the surface to be reconstructed has a freeform shape. The proposed sampling plan based on the mixed parametrisation attains the best performances except in the reconstruction of the surface of revolution.
The sampling of freeform surfaces generated from two or more curves has been investigated. A sampling method that take into account the generatrices curves has been developed. In all the tested scenario the reconstruction with the sampling based on the curves can achieve good results. In almost all the cases it has the smallest RMSE values, this means that with the generating curves it is possible to "understand" the behaviour of the surface.
In this paper a real equal sampling on the area or on the mixed parametrisation has not been computed. It can be performed through an optimisation algorithm, but the samples of the proposed methods can be used as a starting point in order to improve the convergence of the algorithm. The critical points of the curve or surface has not been analysed, these points can be added to the resulting samples or the closet point to each critical point can be moved to achieve better reconstruction results.