Yuan HUANG,Feipeng DA

1School of Automation,Southeast University,Nanjing 210096,China

2Shenzhen Research Institute of Southeast University,Shenzhen 518000,China

Abstract: When obtaining three-dimensional (3D) face point cloud data based on structured light, factors related to the environment,occlusion,and illumination intensity lead to holes in the collected data,which affect subsequent recognition.In this study, we propose a hole-filling method based on stereo-matching technology combined with a B-spline. The algorithm uses phase information acquired during raster projection to locate holes in the point cloud, simultaneously extracting boundary point cloud sets. By registering the face point cloud data using the stereo-matching algorithm and the data collected using the raster projection method, some supplementary information points can be obtained at the holes. The shape of the B-spline curve can then be roughly described by a few key points, and the control points are put into the hole area as key points for iterative calculation of surface reconstruction.Simulations using smooth ceramic cups and human face models showed that our model can accurately reproduce details and accurately restore complex shapes on the test surfaces. Simulation results indicated the robustness of the method, which is able to fill holes on complex areas such as the inner side of the nose without a prior model.This approach also effectively supplements the hole information, and the patched point cloud is closer to the original data. This method could be used across a wide range of applications requiring accurate facial recognition.

Key words:Three-dimensional(3D)point cloud;Hole filling;Stereo matching;B-spline

1 Introduction

Computer vision technology is one of the key core technologies of artificial intelligence (Russell et al., 2010). Meanwhile, visual knowledge (Pan,2019) represents knowledge in a new form, which can be applied in the fields of security, finance, hardware, marketing, automobile design, and medicine.Within the field of computer vision, facial recognition algorithms based on fundamental issues of visual recognition(Pan,2021)have received extensive attention in recent years, and many two-dimensional (2D)and three-dimensional (3D) facial recognition algorithms have been produced over a short period(O’Toole et al., 2012; Black et al., 2016; Gilani and Mian, 2016). At present, 2D facial recognition technology is present in a large number of commercial products. However, compared with 2D recognition,3D facial recognition technology,which involves more extracted information and higher accuracy, is still in the research and development stage.Collection of data and vulnerability to ambient light are two important factors that constrain further development. Threedimensional face data can usually be acquired by laser 3D measurement or raster projection measurement(Schaffer et al., 2011; Liu K et al., 2014; Stone and Skubic,2015).However,if a recognition system uses raster projection measurement for data collection,holes in the cloud point data may appear on the inner side of the nose and at the corners of the eyes in the early stages of data storage.This is due to the need to project auxiliary light during measurement, as well as the degree of cooperation from the user. Holes seriously affect automatic data storage and create issues related to subsequent identification.Therefore,filling point cloud holes is a crucial step in the preprocessing of data.

There are two major difficulties related to filling holes in point cloud data. The first is rapid location and extraction of hole boundaries. Its accuracy directly affects effective filling.Existing hole-boundary extraction methods are based mainly on 3D data.These are based on a point cloud grid(Orriols and Binefa,2003;Jun, 2005) or on the geometric features of boundary points(Floater and Reimers,2001;Kurlin,2014;Nguyen et al., 2015). The former method usually extracts smoother point cloud boundaries, but requires a large amount of computation.The latter is less computationally intensive,but is susceptible to noise and outliers.

The second difficulty,and also the more critical,concerns producing fill-in data that is close to the original data. It is well known that filling a hole in a flat-panel point cloud data set is different from filling a point cloud data set describing a table that contains a bump (the hole is located at the bump). Existing hole-filling algorithms directly employ information related to the neighborhood around the hole, directly using methods such as a spline curve (Chui and Lai,2000; Bendels et al., 2006) or a neural network (Carr et al., 2001) for interpolation. This type of algorithm is simple, easy to use, and does not require any additional information. However, for a hole area with large curvature variation or complicated shape, the filling may not be ideal because of excessive information loss. Pernot et al. (2007) proposed a method of filling holes that involved inserting a topological mesh, which was then deformed.In this approach,the shape of the inserted mesh is the result of minimizing the quadratic function based on the linear mechanical model. This model is used to approximate the curvature variation between the inner mesh and surrounding mesh. It is also possible to specify other geometric constraints to further shape the inserted mesh. Another method has been proposed based on moving least squares to fill point cloud holes (Wang and Oliveira, 2007); it can be used to adaptively interpolate geometric and shadow information.

In recent years, Liu YJ et al. (2016) proposed a classification and repairing strategy based on a fuzzy inference system. This approach involves an algorithm that divides holes into five types,each with corresponding filling strategies. It also considers the recovery of sharp features of the model, based on different hole types. At the same time, researchers have suggested that data can be processed in combination with additional information and surrounding structures. Quinsat and Lartigue (2015) proposed a digital hole-filling method based on mesh deformation.This method relies on a numerical model as a priori knowledge of the nominal mesh.It uses an a priori computer aided design(CAD) model and the continuity of the neighborhood around the hole as a constraint for recovering data.Panchetti et al. (2010) used 2D information related to the corresponding picture and the curvature of the point cloud mesh to repair holes. The employed algorithm relates the position of the mesh vertices to the light intensity of the image pixels. The objective function based on a mechanical model of a bar network is used to simulate the evolution of curvature over the mesh. These methods can satisfactorily fill complex areas.However,no matter whether these establish an a priori CAD model or provide registration between 2D pixels and 3D point clouds, the overall complexity of these algorithms is high. Moreover, the conditions for obtaining the prior model are demanding and not universal.

Based on the above reasons and methodological extension, and considering the measurement principle of a raster projection measurement system, in this study we propose a 3D face point cloud hole-filling algorithm, using binocular stereo matching and a Bspline. The cross-source point cloud data is obtained using the structure from motion (SFM) algorithm of stereo matching technology (Jeong et al., 2011; Chen et al., 2016). The supplementary information points in the hole area are extracted to assist surface iterative computation of the B-spline curve. Pure SFM data cannot be compared with raster data in accuracy, but when the missing information is too much, certain known information is definitely better than the repair of completely unknown information. Combining spline iteration by extracting observation points, our final repair data can be close to the data of pure raster projection. Experimental findings indicate high reliability and potential for practical applications of our approach.

2 Method

The setup of the proposed facial recognition system is shown in Fig.1.

Fig.1 Hardware setup for the facial recognition system

While acquiring point cloud data using a monocular measurement system consisting of a digital light projector (DLP) and charged-coupled device (CCD),two wide dynamic 2D cameras are added for ranging.Simultaneous collection and backup of the 2D picture information can also be achieved. Since the person who is being photographed is required to turn his or her head slightly to the left and right during acquisition,sufficient 2D picture information can be obtained.

A flow chart describing the facial recognition algorithm proposed in this study is shown in Fig.2.

The algorithm described in Fig. 2 may appear rather complicated because it converts the point cloud data back and forth between 2D and 3D representations. However, for the repair of data related to a complex facet without a priori models, this step cannot be avoided. The 3D face data is quickly reconstructed using a standard SFM algorithm and registered with the measured data.The 2D phase information is then used to extract the boundary point sets on the hole areas of the point cloud. The required supplementary information points are extracted according to the boundary information of the registration results. Finally, the hole is repaired using a B-spline on the boundary point data sets with added supplementary points.

Fig. 2 Flow diagram describing the proposed algorithm for facial recognition

The SFM data sets are obtained using the Visual-SFM (Zheng and Wu, 2015) toolkit. There are two main advantages in employing this toolkit. First, it uses the patch-based multi-view stereo (PMVS)(Furukawa and Ponce, 2010) algorithm for dense reconstruction. A point cloud that has been densely reconstructed by PMVS has a relatively high reconstruction accuracy and integrity on the test site (Goesele et al., 2006). Some results are even comparable to laser-scanned data(Shi et al.,2011).In previous experimental findings, the reconstruction integrity of the method was above 90%.The second advantage of using this toolkit is that the reconstruction efficiency is high,and the time complexity of the main algorithm scales asO(n).

3 Data preprocessing

3.1 Acquisition of phase information

We use a raster projection measurement method based on light mechanism. It processes the collected sinusoidal phase shift fringe image sequence to calculate the wrapped phase.We assume that the gray values of the four phase shift images are as follows:

whereIi(u,v) is the gray value of theithimage,I'is the background value of the fringe light intensity,I'is the modulation intensity,andθis the phase field to be obtained.The wrap phase value is

The value range of the wrapped phaseφ(u,v) is[-π, π], as shown in Fig. 3. The true phase value should be

Fig.3 The wrapped phase of a line in an image

In Eq. (3),k(u,v) is an integer, representing the number of fringe cycles where the image point (u,v)is located.We determine the value ofkand expand the wrapped phase into a continuously distributed absolute phase,called the unwrapped phase.The gray code pattern is a set of binary fringe pattern sequences with the number of cycles multiplied, and the maximum number of cycles is consistent with the sinusoidal fringe pattern.We arrange the acquired gray code fringe images in time series,decode the pixel points according to the gray value to obtain a binary number, and uniquely determine the cycle numberk.The gray code encoding principle is shown in Fig.4.

Fig.4 Gray code pattern coding principle

The above calculation process of obtaining the absolute phase is independent of the processing of different pixels, and only the phase value of the key point can be calculated, reducing the time consumed by the algorithm.We take a set of horizontally coded fringe images of Venus’ head, and obtain the phase field results as Fig.5.

Fig.5 Horizontal absolute phase of the Venus surface

3.2 Hole boundary extraction of the point cloud

In order for point cloud data to be repaired, the 2D phase information obtained is used to extract the hole boundary.Our method uses mainly the feature that hole boundary points are easy to identify and extract in a 2D phase map.

For non-closed point cloud data, the boundary points are composed of outer contour points and inner hole area boundary points.In a 2D phase diagram,the point cloud boundary points appear as phase jumps.To eliminate the influence of noise points and remove outliers, the boundary points of the point cloud are defined as phase jumps.The phase of the point itself is not 0,and the phases of the points in the neighborhood are also not equal to 0. In the experiment, the number of neighborhood points is set to eight.

Let the size of the 2D phase map beM×N, andp[m][n] is a point on the phase map. The set of all boundary points of the point cloud isB. The set of outer contour points isBE.The set of boundary points in the inner hole area isHp. Whenp[m][n] satisfies Eq.(4),p[m][n]∈B:

where ?is the threshold for a phase jump, set to 5 in the experiment.

After extracting all the boundary points, further processing of the data is needed.On a 2D phase image,the outer contour points are characterized by the initial and final non-zero points of each row in the point cloud boundary point setB. Whenp[m][n] satisfies Eq.(5),p[m][n]∈BE:

wherem=1, 2, …,M, and hole boundary pointsHp={p[m][n]|p[m][n]∈B-BE}, || is the logical judgement for or. It can be seen from Shi et al. (2011) that the 3D coordinates [XW,YW,ZW]Tof the point cloud data measured using raster projection have the following relationship with corresponding pointp[m][n]=θon the phase map:

whereθis the phase value,a1-a8are system parameters,ρis the scale factor,andAcis the internal parameter matrix of the 3×3 camera.a1-a8andAccan be obtained through the calibration of the system.

Therefore, the 2D phase coordinates of the boundary points of the hole are substituted. The corresponding 3D coordinate pointH={hi|i=1, 2, …,k}of the hole boundary can then be calculated using Eqs.(6)and(7).

Fig. 6 is a schematic describing the extraction of hole boundary points in a face point cloud. The green points are the extracted boundary points. Fig. 6a shows the original face point cloud map, Fig. 6b is the schematic of the boundary points, and Fig. 6c is the schematic after the removal of outer boundary points.

Fig. 6 Hole boundary extraction in a face point cloud:(a) initial point cloud; (b) boundary points (green); (c) point cloud without the outer boundary points

3.3 Point cloud registration

To obtain the supplementary points for the point cloud data set acquired using SFM in the hole area, it is necessary to register point cloud sets acquired by raster projection and SFM.In this paper,the initial registration is achieved using random sample consensus(RANSAC) (Buchin et al., 2009) based on a coplanar four-point set. The improved iterative closest point(ICP) provided in our previous work (Huang Y et al.,2015; Huang Y and Da, 2019) is then employed to accurately register the point cloud. For two crosssource point clouds (Huang XS et al., 2017) of the same object,the data can be converted using

whereRis the rotation matrix and the orthogonal unit matrix is a 3×3 one,T=[TX,TY,TZ]Tis the 3D translation vector,sis the scaling factor,[XW,YW,ZW]is a point in the data setPWobtained by raster projection, [XN,YN,ZN] is a point in the data setPNobtained by SFM, andMis a 4×4 matrix. In the initial registration, the corresponding regions with distinct features are first marked inPWandPN, respectively, to increase the registration efficiency. The RANSAC method is used to randomly select four points of any non-collinear line,find the corresponding points inPN,and calculate the corresponding transformation matrixM.Afterriterations,a transformation matrixMthat minimizes the matching errors between the two point clouds is selected,andMis used to obtain the initial matching result. In the experiment,r=10. In precise registration,first we find the four points closest to the Euclidean distancePN, for each pointpinPW. Next, we find the center of gravity of the smallest neighborhood formed by these four points, thus forming pairs of points.The parametersRandTare then calculated using the least squares method.The final result is applied toPNto obtain the point cloud setPN'. After obtainingPN',the supplementary information points of the hole area required by the algorithm can be extracted to guide the repair of the point cloud hole.The method of screening supplementary information points is as follows:H={hi|i=1,2,…,k}are the hole boundary points of the point cloud. For each pointPNinPN', it must be determined whether it is in the bounding box ofH. If so,p∈PI;otherwise,pis not an element ofPI,wherePIis a set of supplementary information points. At this point, the point cloud setPs=PW∪PIto which supplementary information points are added can be obtained.

4 Point cloud hole filling

After obtaining point setPsconsisting of hole boundary points and supplementary information points,a B-spline surface is employed to initialize hole filling.

4.1 B-spline surface

The number of control verticesbij(i=0,1,…,m;j=0,1,…,n)is set to(m+1)(n+1),while the values of the two columns of surface nodesuandvsatisfyu0≤u1≤…≤um+k+1andv0≤v1≤…≤vm+l+1. Ak×lth-order Bspline surface is then defined using Eq.(9):

where

HereNik(u) andNjl(v) are B-spline basis functions of thekth-order andlth-order B-spline surfaces, defined by two columns of nodesuandv, respectively. The polyhedron composed of (m+1)(n+1) control verticesbijis usually referred to as the control mesh of the Bspline surface, or the control polyhedron. The shape of the control polyhedron generally reflects the shape of the B-spline surface.

4.2 Parameterization of the point set

In the process of B-spline surface fitting for the selected scattered point setPs, the parameterization of scattered points is indispensable. Commonly used parameterization methods are uniform parameterization,centripetal parameterization,and cumulative chord length parameterization. However, these methods are employed mainly for data points in an array of topological matrices. The base surface projection method is used mainly for the parameterization of scattered points.

A 3D point cloud contains information in three spatial directions:X,Y, andZ. However, in research and practical applications, we find thatXandYoften carry more detailed information because of the way in which subjects are used to having facial images taken. Therefore, for the point cloud data setPs, the four extreme points (Xmin,Ymin), (Xmin,Ymax), (Xmax,Ymin), and (Xmax,Ymax) of the 2D (X,Y) coordinates are separately selected to form a quadrilateral plane. The basic surface is then constructed using this as a plane.However, since a small number of data points may fall outside the plane after projection, it is necessary to adaptively enlarge the area of the quadrilateralABCDto ensure that all points of the point setPscan fall in the plane after being projected. The basic curved surfaceBsrequired for the parameterization process is constructed,as shown in Fig.7.

Fig. 7 Construction process for a basic curved surface Bs:(a) the quadrilateral plane is enlarged to ensure that all data points can fall within the plane; (b) the curved sur‐face is constructed

The parameterized coordinates (us,vs) ofPscan be obtained through iteration using the base projection method.AddingPsinto Eq.(9),we obtain

Finally, the least squares method is used to solve the over-constrained equations.The control vertices of the surface are determined to obtain surfaceBsof the preliminary fitting.

However, after initial fitting, the distance error between point setPsand surfaceBsis still relatively large. Therefore, it is necessary to further modify the control vertices,change the shape of the control polyhedron, and change the shape of the surface so that surfaceBsapproaches point setPs. An iterative approximation algorithm based on the Newton iteration method used to improve the fitting accuracy is described in Section 4.3.

4.3 Iterative approximation to fitting surfaces

Given a surfaceC(u,v)and a pointFof an arbitrary parameter, the distance vector between any point on surfaceC(u,v)and pointFcan be expressed as a function, with (u,v) as a parameter, as shown in Eq.(14):

To determine parameter (u0,v0) when the value of ||V(u,v)|| is the lowest, we may set the line where the vectorV(u,v)is located to bel.If linelintersects with surfaceC(u,v) at pointC(u0,v0), the vectorV(u0,v0) needs to be perpendicular to the tangent plane of the surfaceC(u,v) tangent to the pointC(u0,v0),as shown in Fig.8.

Fig. 8 The condition for the shortest distance between point F and arbitrary parameter surface C(u,v)

ThenV(u0,v0)must be perpendicular to the partial differentialsCu(u0,v0) andCv(u0,v0) of surfaceC(u,v)at(u0,v0)in each direction:

The point selected in the point cloud to fit the surface isPs, and the initially fitted surface isS(u,v).The distance vector of any point onPsandS(u,v) can be expressed as a function with(u,v)as a parameter:

and

HereSu(u,v) andSv(u,v) are the first-order partial derivatives of surfaceS(u,v) foruandv, respectively.The barycentric coordinates (us,vs)ofPsare calculated and used as the initial estimate. Eq. (17) is solved according to the Newton iteration method:

where

where δuand δvare the iterative steps in both directions,d=d(us,vs),fu,fv,gu, andgvrepresent the firstorder partial derivatives of the corresponding vector pairuandvat point (us,vs), andSuu,Suv,Svv, andSvurepresent the second-order partial derivatives of surfaceS(u,v)at point(us,vs)foruandv.

According to Eq.(18),we obtain

After solving Eq. (22), iterative step lengths δuand δvcan be determined, and the new parameter (us+δu,vs+δv) is substituted into Eq. (13) to obtain the surface control point, and the parameters of the surface are re-determined.The final value is set toγ＞0(according to the needs of the given system), and Eq. (22) is iteratively solved untiland surfaceS′(u,v)is determined.

To fill the holes in the final fitted B-spline surfaceS′(u,v), we consider a line on surfaceS′(u,v), and then points on the line(Fig.9).

Fig.9 Point collection on a surface based on its curvature

In a point cloud, curvature changes are usually relatively great when the point cloud density is higher.Therefore,the hole can be filled smoothly by selecting points based on such a parameter as curvature. In the process of obtaining the curvature, the original point cloud is first divided into grids,and then theknneighboring points of a particular point are obtained.Based on the relationship between each point and itsknneighboring points, the average curvatureρˉ0and the average point distancedˉof the original point cloud are then obtained. The average curvatureρˉsofPs(s=0,1,…,t)is similarly obtained,and the interval between collected points is set to Δω=ρˉ0dˉρˉs.When collecting points on the B-spline surfaceS′(u,v), the equal-parameter curve is obtained at an equal interval of Δωin the direction ofuon the surface. Points are then collected at an equal interval in the direction ofvon each equal-parameter curve, to obtain discrete points for filling the holes.

5 Simulation results

To verify the effectiveness of the proposed algorithm, simulations were conducted on two models: a smooth ceramic cup and a human face. The simulations were performed on a computer with Core i5, 4 GB memory, and the Windows 10 operating system. Matlab 2015b was used as the software environment.

First,the algorithm was tested on ordinary smooth surface objects. Fig. 10a shows the ceramic cup point cloud data obtained by raster projection measurement.It can be seen that since the reflectivity of the ceramic cup was too high, a hole was formed in the middle of the point cloud. Fig. 10b shows the point cloud data measured using VisualSFM. Since the smooth ceramic cup itself lacked texture, speckle information was projected to improve the reconstruction accuracy when acquiring point cloud data. Fig. 10c shows the point cloud map after the holes were repaired by the proposed algorithm. The data set generated by VisualSFM provided 96 supplementary information points for hole repair, and the number of added patch points was 404.During ceramic cup hole repairing, when the B-spline curve was used to compute the surface equation of the corresponding points,the supplementary points extracted from the point cloud data obtained by VisualSFM will provide this process with additional information.As a result, the surface equations more closely matched the original data, and the hole repair results were more accurate.

Fig. 10 Point cloud of a ceramic cup repaired using raster projection (a), VisualSFM (b), and the proposed algo‐rithm(c)

For the more complex human face, the human face point cloud data measured by raster projection is shown in Fig. 11a. The eye and eyebrow areas contained holes due to low reflectance. There were also holes in the nose area due to its own occlusion.The point cloud data of the face measured by VisualSFM is shown in Fig. 11b, and the point cloud map after the holes were repaired by the proposed algorithm is shown in Fig. 11c. It can be seen that the repairing result of a complex area, such as the nose and eyes, accurately reflected the shape of its original form.

Fig. 11 Point cloud of a human face repaired using raster projection (a), VisualSFM (b), and the proposed algorithm(c)

To further verify the effectiveness of the method,two facial areas where there were typically data holes were selected for further analysis (Fig. 12).The hole shown in Fig. 12b is near one side of the nose.The hole in Fig.12c is between the eyebrows.

Fig. 13a shows the complete point cloud data,while Figs.13b and 13c are point cloud maps repaired by the point cloud hole repairing algorithm, based on the hole in Fig. 12b. The green part is a fixed point cloud with a total of 1612 points. A graph obtained using ICP to register added points after repairing and the original hole area points is shown in Fig. 13d. It can be seen that for a large hole with a large variation in shape, the point cloud data obtained by the algorithm fitted the original shape completely. The original facial details were accurately described.

Fig. 12 Positions of point cloud holes: (a) complete point cloud; (b) hole near one side of the nose; (c) hole between the eyebrows

Fig. 13 Simulation results of point cloud hole repairing in cheek: (a) complete point cloud; (b) left-hand view of the repaired point cloud; (c) a top view of the repaired point cloud;(d)registration of the repaired point cloud

Fig. 14a shows the complete point cloud data for the space between the eyebrows on the test face,Figs. 14b and 14c show the patched point clouds, and Fig. 14d shows the result after registration with the complete point cloud. It can be seen that the added point cloud was consistent with the original point cloud in terms of curvature variation. The repaired points clearly described the detailed shape between the eyebrows.These results indicated that the proposed algorithm can effectively repair holes in areas with complex shapes.

Fig. 14 Simulation results of point cloud hole repairing in wing of nose:(a)complete point cloud;(b)left-hand view of the repaired point cloud; (c) a top view of the repaired point cloud;(d)registration of the repaired point cloud

Table 1 is a comparison between the implemented and initial points for the holes shown in Figs. 12b and 12c. The number of ICP iterations was 10. The root mean square errors (RMSE) of the point cloud registrations in Figs. 13d and 14d were 0.019 and 0.072, respectively. The implemented and initial points were consistent in terms of both their numerical values and curvature variation,both in large areas and in areas where the curvature varied greatly. For objects with complex shapes, the proposed algorithm can effectively restore the original shape. It should be mentioned that due to the large number of steps in the overall algorithm, it took a little more time to repair a hole alone, which is a problem to be solved in the future.

Table 1 Comparison between the implemented and initial points shown

6 Conclusions

In practical 3D face recognition systems, the stored point cloud data often contains holes. In this paper,a hole-filling algorithm based on stereo matching and a B-spline is proposed, using a laboratorydeveloped identification system combined with a consideration of measurement principles. The algorithm first extracts a hole boundary using features of the 3D hole boundary points in the 2D phase. The SFM data set and the raster point cloud data set are then registered, and the supplementary information points are extracted. Finally, the B-spline curve surface is used to further repair the holes,making the repaired point cloud be more consistent with the original data. Simulations on a smooth ceramic cup and a human face show that the proposed algorithm can effectively reproduce surface details and accurately restore complex surface shapes. Since the algorithm requires a high-precision SFM data set, although the accuracy and completeness of the VisualSFMreconstructed point cloud is high for most objects,this method is not particularly effective for objects that lack texture. In this case, how to improve the reconstruction accuracy is an issue to be addressed.In addition, we will explore adaptive data repairing of complex objects without a priori models to develop the proposed algorithm for use across a wide range of facial recognition applications.

Contributors

Yuan HUANG designed the research, processed the data,and drafted the paper. Feipeng DA helped organize the paper.Yuan HUANG and Feipeng DA revised and finalized the paper.

Compliance with ethics guidelines

Yuan HUANG and Feipeng DA declare that they have no conflict of interest.