Digital representation of real 3D surfaces has a crucial importance in a variety of cutting-edge applications, such as autonomous navigation (Huang J. et al., 2021), UAV fleets (Ji et al., 2021), extended reality streaming, or telesurgery (Huang T. et al., 2021). Point clouds represent 3D surfaces by means of a set of 3D locations of points on the surface. In general, those points can be acquired by active or passive techniques (Chen et al., 2021; Rist et al., 2021), in presence of random errors, and they may be associated with color and texture information as well. Point cloud denoising can in general be applied as an enhancement stage at the decoder side of an end-to-end communication system, involving volumetric data, for e.g., for extended reality or mixed reality services. Although lossless compression of point clouds is feasible (Ramalho et al., 2021), color point cloud lossy coding based on 2D point cloud projection (Xiong et al., 2021) is increasingly relevant both in sensor networks (de Hoog et al., 2021) and autonomous systems (Sun et al., 2020). Nonlocal estimation solutions (Zhu et al., 2022) or color-based (Irfan and Magli 2021b) solutions as well as solutions for point cloud sequences (Hu et al., 2021a) have been proposed.

The extraction of visually relevant features on point cloud is needed for tasks as pattern recognition, registration, compression and quality evaluation (Yang et al.s, 2020; Diniz et al., 2021), and semantic segmentation. Point cloud (PC) processing has been widely investigated, and many of the proposed processing methods are based on the geometric properties (Hu et al., 2021b; Erçelik et al., 2021). Still, feature extraction on PC is mostly focused on information related to the point cloud shape. In addition, new PC acquisition systems for surveillance (Dai et al., 2021) or extended reality (XR) (Yu et al., 2021) require processing tools operating both on geometry and texture. Shape and texture processing needs development of new tools because of the non-Euclidean nature of the real surfaces modeled by a point cloud. In this direction, few studies in the literature simultaneously leverage both geometry and texture information. Furthermore, in the context of classical image processing, several effective tools have been inspired by the human visual systems (HVSs), which process both texture and shape, and it is sensitive to angular patterns such as edges, forks, and corners (Beghdadi et al., 2013). In particular, two-dimensional circular harmonic functions (CHFs) have been investigated for visually driven image processing and specifically for angular pattern detection. CHFs have been successfully applied to interpolation (Colonnese et al., 2013), deconvolution (Colonnese et al., 2004), and texture synthesis (Campisi and Scarano 2002). On the contrary, point cloud processing lacks HVS-inspired processing tools, which can, in principle, provide alternative perspectives.

In this article, we leverage a class of point cloud multiscale anisotropic harmonic filters (MAHFs) inspired by HVS. MAHFs were recently introduced in our conference article (Conti et al. (2021)). First, we recall the MAHF definition and describe their local anisotropic behavior that highlights directional components of the point cloud texture or geometry. In addition, we show their applicability to both geometric and textured PC data. Second, we illustrate how MAHF can be applied to visually driven PC denoising problems. Denoising is a crucial preprocessing step for many further point cloud processing techniques. In real acquisition scenarios, the perturbations on the PC vertices or on the associated signal severely affect the PC usability. MAHF is used to drive an iterative denoising algorithm so as to adapt the restoration to the local information. The proposed method differs from other competitors (Zhu et al., 2022 and the references) in linking the denoising with visually relevant features, as estimated by suitable anisotropic filtering in the vertex domain. We test the performance of the visually driven point cloud denoising algorithm (VIPDA) on synthetic and real data from the public database (Turk and Levoy 1994; d’Eon et al., 2017). Specifically, considering different signal-to-noise-ratios by adding Gaussian noise to the original data, we verify that our method outperforms state-of-the-art alternatives in denoising data.

The structure of the article is as follows. In Section 2, we review a particular class of HVS-inspired image filters, namely the circular harmonic filters, which are needed to introduce our point cloud filtering approach. In Section 3, we present a class of multiscale anisotropic filters, formerly introduced in Conti et al. (2021), and we illustrate their relation with visually driven image filters. In Section 5, we present the visually driven point cloud denoising algorithm (VIPDA) based on the proposed manifold filters. In Section 6, we show by numerical simulations that the VIPDA outperforms state-of-the-art competitors. Section 7 concludes the article.

Before the introduction of the MAHFs, a step back is necessary in order to contextualize the research problems by investigating other filter methods in the Euclidean domain.

In several important applications in the field of image processing, circular harmonic functions (CHFs) have been used (Panci et al., 2003; Colonnese et al., 2010). As mentioned previously, CHFs have been widely applied in image processing applications because they are able to detect relevant image features, such as edges, lines, and crosses, i.e., they perform the analysis in an analogous way to the behavior of the HVS during the pre-attentive step. It is important to remark here that the results of the order-1 CHFs are complex images, in which the module corresponds to the edge magnitude while the phase describes the orientation. Taken together, this filtering procedure returns precious information about the structures of the output image; in fact, it underlines the edges by simultaneously measuring their intensity and direction. The interest in CHFs also stems from the fact that they can be integrated within an invertible filter bank, thereby being exploited for suitable processing, for e.g., image enhancement, in the CHF-transformed domain (Panci et al., 2003).

CHFs’ properties relate to the specific way in which they characterize the information belonging to two points. In fact, they encode the distance as well as the geometric direction that joins them.

These aspects are evident in the mathematical formulation of CHFs. Let us consider the 2D domain of the continuous CHF described by the polar coordinates (r, ϑ) that, respectively, represent the distance from the origin and the angle with the reference x axis. The CHF of order k is the complex filter defined as: h k r , ϑ = g k r e j k ϑ , (1) where the influence of the radial (r) and the angular (ϑ) contributions are separated by the two factors. With the aim of preserving the isomorphism with the frequency space, the functions g _k (r) in Eq. 1 are usually isotropic Gaussian kernels. The variable k defines the angular structure of the model. For k = 0, the zero-order CHF returns output as a real image, represented by the low-pass version of the original one. As a general consideration, when the order k increases, CHFs are able to identify more and more complex structures on the images, such as edges (for k = 1), lines (for k = 2), forks (for k = 3), and crosses (for k = 4). We can see the effect of increasing the k order of the heat kernel on a sphere in Conti et al. (2021).

Based on the definition of CHF and the introduction of a scale parameter α, circular harmonic wavelet (CHW) (Jacovitti and Neri 2000) of order k can be introduced and can be typically applied in the context of multi-resolution problems.

Finally, it is worth observing that the CHF output has been shown to be related to the Fisher information of the input w.r.t rotation and translation parameters. The Fisher information of an image w.r.t. shift/rotation estimation is associated with the power of the image first derivative w.r.t. the parameter under concern (Friedlander 1984). The CHF in the 2D account for a local derivative of the signal has been shown to be related to the Fisher information w.r.t. localization parameters (Neri and Jacovitti 2004).

In order to show a visual example of the application of CHF on real data, we consider the “cameraman image,” and we apply first-order CHF as an example to represent the effect of CHF on a real image. We report the results in Figure 1, in which we can see the module in panel A and the phase in panel B.

Stemming on these studies on directional harmonic analysis of 2D signals, we introduce in the following sections the multiscale anisotropic harmonic filters to be adopted in the non-Euclidean manifold domain.

Although the HVS is very complex in nature, its low-level behavior, as determined by the primary visual cortex, is well characterized by being bandpass and orientation selective (Wu et al., 2017). Therefore, the CHF mimics these features on 2D images, and in this section, we show how to extend this behavior to manifold filters. Specifically, in this section, we describe a new class of visually driven filters operating on a manifold in the 3D space; the preliminary results on such filters appear in Conti et al.(2021). We extend the presentation in Conti et al. (2021) by an in-depth analysis of their relation with the CHF and by providing new results about their applications to point cloud filtering.

Our general idea consists in the extension of CHFs to 2D manifolds embedded in 3D domains. In this direction, the two key points to adapt to this new scenario are as follows: 1) we need to define a smoothing kernel that corresponds to the isotropic Gaussian smoothing in the 2D case; and 2) we have to identify an angular measurement on the surface of the manifold in the 3D space.

In the following sections, we elaborate on the filters description first in the case of a 2D manifold defined in a continuous 3D domain and then in the case of its discretized version, as represented by a point cloud. The main notation is reported in Table 1.

Let us consider a real-valued D-dimensional signal on the point cloud vertices s(p _i ) ∈ R ^D , i = 0, …, N − 1. Applying k-th-order MAHF for the vertex domain signal on graph filtering obtains the output point cloud signal r(p _i ) as: r p i = ∑ j = 0 N − 1 φ k p i , p j s p j . (7)

The filtering realized by the MAHFs performs an anisotropic harmonic angular filtering of the signal defined on the point cloud. The period of the harmonic analysis decreases as the filter order increases, and MAHF of different orders k are matched to different angular patterns of the input signal s(p _i ) i = 0, …, N − 1.

MAHF applies to both texture and shape signals, depending on the choice of the input signal s(p _i ) i = 0, …, N − 1. For video point clouds, the signal s(p _i ) can represent the luminance and the chrominances observed at the point p _i . If this is the case, the MAHF output highlights texture patterns on the surface. On the other hand, the signal s(p _i ) i = 0, …, N − 1 can be selected so as to represent geometric information. A relevant case is when the signal represents the normal to the point cloud surface at each vertex p _i as follows: s p i = n x p i n y p i n z p i = n p i , i = 0 , … , N − 1 . (8)

The application of MAHF to the signal defined as in Eq. 8 will be exploited in the following derivation of VIPDA.

To sum up, the MAHFs can be applied to different kinds of data defined on point clouds, and they provide a way to extract several point-wise shape and texture point cloud features for different values of the order k. A schematic representation of MAHF application to texture (luminance) and shape (point cloud normals) information is illustrated in Figure 3.

In this section, we illustrate the VIPDA approach, based on application of the aforementioned MAHF to the problem of PC denoising.

Let us consider the case in which the point cloud vertices are observed in presence of an additive noise. Thereby, the observed coordinates are written as q i = p i + w i , i = 0 , … , N − 1 , (9) where w _i is an i.i.d. random noise. The denoising provides an estimate p ̂ i of the original locations p _i . Let us remark that this problem is different from recovery of a signal defined at the vertices, which will be addressed in the future work.

Point cloud denoising algorithms typically leverage 1) data fidelity (Irfan and Magli 2021a), 2) manifold smoothness (low-rankedness) (Dinesh et al., 2020), and 3) local or cooperative averaging (Chen et al., 2019) objectives.

Here, the local manifold smoothness is accounted by adapting the estimator p ̂ i to the local shape variability as estimated at the first algorithm stage. Specifically, the MAHF is applied to the point cloud estimated normals n ⁽⁰⁾(q _j ) as: ν 0 q i = ∑ j = 0 N − 1 φ k q i , q j n 0 q j . (10)

Thereby, each point q _i is assigned a weight related to the normal variations in its neighborhood. The key idea is that when fast variations of the normal are observed around q _i , the surface smoothness is reduced, and the set of neighboring points to be exploited to compute the estimate p ̂ i should be reduced accordingly. Therefore, the size K _i of the neighborhood of the point q _i to be used in the estimation stage is selected based on the MAHF-filtered signal ν (q _i ). Specifically, K _i is selected as a function of the mean square value of the MAHF output, namely K i = K ‖ ν q i ‖ 2 , where K is an integer function defined on R. This is exemplified in Figure 4 (left), which illustrates a point cloud (namely a low-resolution version of the Stanford bunny in Turk and Levoy (1994)) and two neighborhoods of different sizes K _i and K _j around points characterized by different values of the mean square value of the MAHF output, namely ‖ ν (q _i )‖², plotted in Figure 4 (right).

After the size of the estimation window at each point is given, the denoising algorithm iteratively alternates 1) the computation of a candidate estimate of the point location based on spatially adaptive averaging over the K _i -size neighborhood of the i-th vertex and 2) the update of the current estimate, along the direction of the normal to the surface.

In formulas, at the l-th iteration, the candidate estimate of the i-th point cloud vertex is computed as p ̃ i l = α 0 p ̂ i l − 1 + ∑ j ∈ η i ; K i α j p ̂ j l − 1 , where η(i; K _i ) denotes the set of K _i nearest neighbors of the i-th point cloud vertex. Then, the estimate is updated as p ̂ i l = p ̂ i l − 1 + ρ 0 l n l − 1 p ̂ i l − 1 n l − 1 p ̂ i l − 1 T p ̃ i l − p ̂ i l − 1 , where ρ ₀ ∈ (0, 1] is a parameter controlling the update rate throughout the iterations. Finally, the normals n ( l ) ( p ̂ i ( l ) ) are recomputed on the estimated point cloud p ̂ i ( l ) , i = 0 , … , N − 1 .

To sum up, the MAHF is applied once for all at the beginning of the iterations. For each vertex, the set of neighboring points is identified. Then, a candidate new point is computed as a weighted average of the neighbor and of the point itself. Then, the point estimated at the previous iteration is updated only by projection of the correction on the direction of the normal to the surface, as illustrated in Figure 5.

The normals to the surface are recomputed. The algorithm is terminated after few iterations (1-3 in the presented simulation results).

The algorithm is summarized as follows:

Input: Noisy point clouds coordinates q _i , i = 0, …, N − 1; heat kernel spread t; update rate control parameter ρ ₀.

Output: Denoised point clouds coordinates p ̂ i , i = 0 , … , N − 1 .

Inizialization:

Computation of the heat kernel. K t ( G ) ( p i , p j ) , j = 0 , … , N − 1 , i = 0 , … , N − 1

Computation of the normals n ⁽⁰⁾ and of their filtered version. ν ( p i ) = ∑ j = 0 N − 1 φ ( 1 ) ( q i , q j ) n ( 0 ) ( q j )

Computation of K i = K ‖ ν ( q i ) ‖ 2 and of the K _i nearest neighborhood η(i; K _i )

Iteration: for l = 1, …, L − 1

Computation of the candidate estimate p ̃ i ( l ) = α 0 p ̂ i ( l − 1 ) + ∑ j ∈ η ( i ; K i ) α j p ̂ j ( l − 1 ) , with α _j = (α ₀)/K _i

Update of the current estimate. p ̂ i ( l ) = p ̂ i ( l − 1 ) + ρ l [ n ( l − 1 ) ( p ̂ i ( l − 1 ) ) ] [ n ( l − 1 ) ( p ̂ i ( l − 1 ) ) ] T p ̃ i ( l ) − p ̂ i ( l − 1 )

Update of the normals n ( l ) ( p ̂ i ( l ) ) .

An overview of the VIPDA algorithm stages appears in Figure 6.

In this section, we present simulation results associated with the application of MAHF on synthetic and real PC, and we measure the performances of VIPDA. In particular, in subsec.6.1, we illustrate the MAHF behavior, also in comparison with CHF, and in subsec.6.2 we assess the performance of VIPDA, also in comparison with state-of-the-art denoising algorithms.

This work has presented a novel point cloud denoising approach, the visually driven point cloud denoising algorithm (VIPDA). The proposed method differs from other competitors in linking the denoising with visually relevant features, as estimated by suitable anisotropic angular filters in the vertex domain.

VIPDA leverages properties inspired by those of the human visual system (HVS), and it is viable for application on the texture and geometry data defined over a point cloud. The VIPDA approach leads to smooth denoised surfaces since it iteratively corrects the noise component normal to the manifold underlying the point cloud by projecting the observed noisy vertex toward the plane tangent to the underlying manifold surface. The performance of VIPDA has been numerically assessed on real open access data and compared with state-of-the-art alternatives.

The proposed algorithm is effective in denoising real point cloud data, and thereby it makes point cloud modeling more suitable for real applications. Our findings pave the way to HVS-inspired point cloud processing, both for enhancement and restoration purposes, by suitable anisotropic angular filtering.