A smart stochastic approach for manifolds smoothing
Eurographics Symposium on Geometry Processing 2008 Pierre Alliez and Szymon Rusinkiewicz (Guest Editors)
Volume 27 (2008), Number 5
A smart stochastic approach for manifolds smoothing A. F. El Ouafdi†1 , D. Ziou1 and H. Krim2 1
Départment d’informatique, Université de Sherbrooke, Québec, Canada Engineering Department, North Carolina State University, USA
2 Computer
Abstract In this paper, we present a probabilistic approach for 3D object’s smoothing. The core idea behind the proposed method is to relate the problem of smoothing objects to that of tracking the transition probability density functions of an underlying random process. We show that such an approach allows for additional insight and sufficient flexibility compared with existing standard smoothing techniques. In particular, we are able to propose a newer, faster, and simpler smoothing approach that retains and enhances important manifold features. Furthermore, it is demonstrated to improve performance over existing smoothing techniques.
1. Introduction Recent developments in acquisition technology make it possible to obtain objects with highly detailed features. During all stages of the object’s acquisition process, noise is introduced by the inaccuracy of a sensor’s calibration, imprecision in registration procedure, limited sampling resolution etc. The noise is usually described by a probabilistic model and the degradation often yields the resulting observation object, most commonly presented in additive form U = U0 + η
(1)
where U is the noisy object, U0 is the corresponding noisefree object, and η is noise mostly assumed to be Gaussian. The challenge facing the smoothing process consists of recovering the underlying object U0 from U while preserving the object’s crucial features such as edges and corners. 1.1. Related work In recent years, various approaches have been proposed to tackle the problem of object smoothing while preserving crucial characteristics. Most smoothing methods can be classified into two main categories, depending on whether the smoothing process is performed under a probabilistic or deterministic framework. Both categories are tightly connected † Corresponding authors. [email protected], [email protected]
E-mail:
c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation ° Published by Blackwell Publishing, 9600 Garsington Road, Oxford OX4 2DQ, UK and 350 Main Street, Malden, MA 02148, USA.
(a)
(b)
Figure 1: Smoothing a scanned model: On the left is the input Chinese dragon model and on the right is the result of the proposed smoothing method.
since they are based directly or indirectly on the heat diffusion principle. In what follows, we briefly review some examples of smoothing methods in both categories. For a recent review, we refer the interested reader to [BPK∗ 07]. Deterministic methods: In this category, one can distinguish two families of smoothing approaches that either discretize a PDE expressing the heat diffusion principle or perform the smoothing process along the normal vector. The key step in the PDE based smoothing method is the PDE discretization step usually performed by the finite element FE , finite volume FV, or finite difference FD methods. Such discretization procedures replace the PDE by a linear system
1358
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
of algebraic equations. An implicit or semi-implicit time discretization is generally adopted to achieve a stable numerical schema. In this setting, [CDR04] discretize the PDE of the heat diffusion by a FE method and estimate the diffusion tensor by fitting the local surface quadratically. Such a smoothing method enhances the triangular surface in the principal directions of curvature. Bajaj and Xu [BX03] generalize the work of Claranz et al. to a higher dimension and combine a limit function of loop subdivision with a PDE to smooth higher order functions on the surface while fairing it. Using a level set method, Tasdizen et al. [TWBO03] propose a twostep smoothing method by solving simultaneously two PDEs on level set surfaces. Very recently, [EZ08] decompose the heat diffusion principle into basic laws and derive the numerical schema which allows one to perform the smoothing process by discretizing the basic laws using computation algebraic topological tools. In the particular case when the diffusion tensor is a scalar and locally constant function, the PDE that expresses the heat diffusion principle is reduced to the mean curvature evolution equation along the normal vector. Along this vein, the mean curvature flow smoothing method proposed by [MDSB03] adjusts the vertex position along the surface normal with speed proportional to the mean curvature flow. The important features are preserved by a weighting parameter which depends on the principal curvature values. [HP04] propose a prescribed mean curvature flow method by preconditioning the anisotropic mean curvature vector. The diffusion tensor is estimated by a normal cycle approach. The recent work of [YBS06] extended the non-local image filter to meshes by computing a local radial basis function approximation to define the similarity measure. Probabilistic methods: Due to the nature of the noise, the smoothing process may also be interpreted from a statistical point of view. The noise is considered as a stochastic process with known parameters. In this setting, [PSZ01] suggest an adaptative Wiener filtering smoothing technique for semi-regular meshes. A pre-processing step is needed to convert the meshes to multi-resolution representations before smoothing. Based on robust statistics and local statistical predictors of surfaces, [JDD03] propose a non-iterative smoothing method. The surface features are preserved by means of a robust influence weight function together with a predictor for vertex positions based on the surface’s tangent plane. A similar statistical predictor also appears in the bilateral mesh smoothing approach proposed by [FDCO03]. Recently, [BK04] provide a stochastic interpretation of the standard [PM90] diffusion method for image filtering. Also they propose a two-sided gradient edges detector to filter a noisy input signal while preserving important features. Intuitively, the authors model the nonlinear diffusion by a differential stochastic equation SDE. Next they construct a probability density of the SDE to resolve the nonlinear or constrained diffusion problem.
1.2. Contribution The anisotropic diffusion is a powerful smoothing approach that performs the directional smoothing and enhancing of the manifold features while preserving adequately curvatures, edges and corners. To the best of our knowledge, the existing anisotropic smoothing methods discretize the PDE heat diffusion by a FE or FD methods [Wei94, BX03, CDR04]. The derivation of a stable numerical schema from the discretization of the PDE by the such as discretization methods requires the resolution of a linear system at each vertex at each iteration, which induces numerical error, time consuming and computationally expensive. Closely related to the intuitive idea of the image smoothing method proposed by [BK04]. We construct in this work an anisotropic Gaussian kernel to smooth an arbitrary manifold initially corrupted with natural or artificial noise. Such kernel is the transition probability density of an constrained diffusion process. The proposed method is simpler than solving the heat diffusion numerically. Additionally, the proposed approach is fast, unconditionally stable and able to retain and even enhance important manifold features such as edges and corners. 1.3. Paper Organization In the next section we begin with a brief description of anisotropic diffusion process on manifold, and then, we deduce a short time solution of the anisotropic diffusion problem. In Section 3, we discretise the diffusion process on triangular meshes and we propose a new algorithm for 3D objects smoothing whose evaluation through numerous substantiating examples is provided in Section 4. We finally give some concluding remarks in Section 5. 2. Diffusion on a Riemannian manifold In this section, we begin with a brief description of the constrained diffusion process from the stochastic point of view. Then, we deduce a short time solution of anisotropic diffusion on a manifold. 2.1. Constrained diffusion process on a manifold Assume that the object to be processed can be represented by an n-dimensional smooth and compact manifold M embedded in Rn+1 . The diffusion processes constitute a powerful and a fundamental tool for smoothing manifolds. They often arise as physical processes in the real world such as the heat diffusion. The diffusion process, known in physics as the transfer of energy from the more energetic to the less energetic particles, can be assimilated to the random motions of a virtual particles. For instance, the particle trajectory X(t) = y at the point y at time t is a Markov process governed by the following stochastic differential equation SDE dX(t) = µ(X(t)) dt + σ(X(t)) dB(t)
(2)
c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
1359
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
X(0) = x ∈ Ω ⊆ M t ≥ 0, where Ω is a smooth compact domain of M, B(t) is mdimensional Brownian motion that generally models the noise term in (1), µ is the drift velocity vector field and Σ = 12 σσT is the diffusion tensor assumed to be symmetric and positive definite [IW81]. For an unconstrained diffusion process, the drift velocity vector field and the diffusion tensor are proportional to the first and second moments of random changes in coordinates, respectively. To generate a diffusion path that obeys the desired diffusion behaviour, a constrained or anisotropic diffusion process is achieved by modifying the data-driven control parameter µ and Σ [Mor04]. In computer vision, several smoothing methods may be interpreted as constrained diffusion processes. For instance, the standard gradient based smoothing method proposed by [PM90] can be interpreted as a constrained diffusion process with vanishing drift velocity coefficient, and a diffusion coefficient equal to a conductance function that varies with the local gradient amplitude, which allows to preserve important sharp features such as edges in images [BK04]. More generally, most smoothing methods that perform the directional smoothing of the manifold constrain implicitly the diffusion process through a diffusion tensor based on the principal curvatures ki and the principal direction of curvatures ei , for i = 1, · · · n. Such a diffusion tensor is generally given by Σ˜ = UVU T
(3)
where U = (e0 | · · · |en |N), V is a diagonal matrix with diagonal term (G(k0 ) · · · G(kn ), 1), N is the unit normal vector, ki and ei are respectively the eigenvalues and eigenvectors of the shape operator matrix, and G is a conductance function. When constrained, the diffusion process endows implicitly the manifold with a metric equal to the inverse of the diffusion tensor (Σ−1 ) in order to constrain the particles dynamics [Kam86]. In our case, as the diffusion tensor (3) is symmetric and positive definite for a bounded eigenvalues. Its inverse is also a symmetric and definite positive. Thus, the tensor Σ˜ −1 is a natural choice to define a new Riemannian metric on M to achieve a directional smoothing of the manifold. One may notice that the inverse tensor metric on manifold induces differential operators like the gradient and Laplace-Beltrami operators which are generally encountered in the manifold smoothing methods.
In many practical problems, the determination of the transition probability function of the diffusion process X(t) is of interest [Gar96]. The transition probability function on Ω, noted pΩ (t, x, y), describes the probability for a particle starting initially at the location X(0) = x to reach X(t) = y ∈ Ω at time t. For a diffusion process on a Riemannian manifold, it can be shown (see for example [Mol75, Hsu02]), that the transition probability function pΩ (t, x, y) is a minimal positive solution of the following backward Fokker-PlanckKolmogorov FPK problem in Ω ∂u(x,t) = (∆˜ M + µ)u(x,t), (x,t) ∈ Ω × (0, ∞), ∂t u(x,t) = 0 (x,t) ∈ ∂Ω × (0, ∞)
(4)
u(x, 0) = u0 (x) x ∈ Ω, with u0 (x) = δx , where δx is the Dirac function. Here ∆˜ M denote the Laplace-Beltrami operator associated with the metric tensor Σ˜ −1 defined on the local region Ω. For a short time period and for vanishing drift term µ = 0, the transition probability function, minimal positive solution of (4), has the following asymptotic expansion [Dod83] ˜
2
pΩ (t, x, y) ≈ Φ(x, y)e−d(x,y)
/4t
(5)
and fulfills asymptotically the normalization condition Z
Ω
pΩ (t, x, y)dy = 1,
(6)
where Φ(x, y) is a continuous function on Ω × Ω such that ˜ y) = (x − y)T Σ˜ −1 (x − y) Φ(x, y) ≈ 1 for x close to y, d(x, is the Riemannian distance between x and y relative to the metric tensor Σ˜ −1 . In a general case, when the initial condition u0 (x) is a continuous function on Ω, the FPK problem (4) admits an unique solution of the form Z
u(x,t) =
Ω
pΩ (t, x, y)u0 (y)dy.
(7)
The above equation may be interpreted as an anisotropic convolution between the Gaussian kernel pΩ and the continuous function u0 . Note that the function u0 may also be interpreted as a quantity represented on the region Ω at time t = 0 (e.g., initial temperature, gray levels, colors, vertices), in which case, u(x,t) represents the anisotropic diffusion of u0 at time t; this is a behaviour similar to what is required in anisotropic smoothing. 3. Diffusion on triangular meshes
2.2. A short time solution of the anisotropic diffusion problem on manifold In general, it is difficult to find an explicit solution X(t) of SDE (2). Such a solution is only possible when (2) can be transformed into a linear equation of low dimension. In addition, an accurate numerical solution requires the generation of a large number of Brownian motion paths with a small time step, which is computationally expensive [KP92]. c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
Discretization is the basis of many numerical procedures, this also holds continuous for the diffusion process. The object of discretization is the Riemannian manifold M and the continuous time axis. For our purposes, let us review this discretization on a 2-manifold M2 embedded in R3 . To solve the diffusion process problem numerically, the 2dimensional spatial domain M2 on which the diffusion process takes place is discretized by a polygonal or triangular
1360
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
mesh M2 = (V, E, S), where V = {v1 , · · · , vm } is the set of vertices, E = {ei j } is the set of edges, and S = {S1 , · · · , Sn } is the set of triangles. Two vertices vi , v j ∈ V are adjacent( written vi ∼ v j ) if they are connected by an edge ei j ∈ E. The neighborhood of a vertex vi is the set Vi∗ = {v j ∈ V : vi ∼ v j } ∪ {vi }. We denote by Si∗ the set of all triangles incident to the vertex vi . To each vertex vi , we associate a triangular local surface Ωi = ∪S j ∈Si∗ S j formed by the union of all triangles in Si∗ . The time discretization is performed with a constant time step. We suppose that we have nt equal time intervals [tk ,tk+1 ] of length τ = tk+1 − tk . In what follows, we consider that tk represents the initial time of the diffusion. After the space and time discretization, the second step of the discretization procedure is to build from (5) a normalized transition probability function on the polygonal local surface Ωi . Thus, let us consider the truncated normalized transition probability function on Ωi given by pΩi (τ, x, y) . Ωi pΩi (τ, x, z)dz
(8)
p˜Ωi (τ, x, y) = R
The denominator term is introduced to correct the truncation error and to enforce the normalization condition (6) of p˜Ωi on Ωi . Such denominator term may be expressed at x = vi as Z
Z
Ωi
pΩi (τ, vi , z)dz =
∑
S j ∈Si∗ S j
e
˜ i ,z)2 /4τ −d(v
dz.
(9)
Using the standard integrals approximation techniques, we can approximate the integration on each triangle S j formed by the vertices v j0 , v j1 and v j2 as Z
˜
Sj
2
e−d(vi ,y)
/4τ
dy ≈
| S j | 2 −d(v ˜ i ,v j )2 /4τ r e . ∑ 3 r=0
By substituting the above expression in (9), we reformulate the normalized transition probability function p˜Ωi as ˜
p˜Ωi (τ, vi , v j ) =
2
e−d(vi ,v j ) ∑vl ∈Vi∗ αl
/4τ
˜ i ,vl )2 /4τ e−d(v
,
(10)
T
|S |
where αl = ∑ p 3p , S p ∈ Sl∗ Si∗ where Sl∗ is the set of ˜ i, v j ) = all incident triangles to the vertex vl ∈ Vi∗ and d(v T ˜ −1 (vi − v j ) Σ (vi − v j ). Having the normalized transition probability (10) at hand, we can now express the local solution u(vi ,tk+1 ) (7) at vertex vi at time tk+1 with initial condition u(y,tk ) on Ωi as Z
u(vi ,tk+1 ) =
∑
S j ∈Si∗ S j
p˜Ωi (τ, vi , y)u(y,tk )dy.
(11)
As before, using the standard integrals approximation techniques, we can rewrite u(vi ,tk+1 ) in the above expression as a normalized weighted sum of u(v j ,tk ), v j ∈ Vi∗ as u(vi ,tk+1 ) =
∑
v j ∈Vi∗
ωi j u(v j ,tk ),
(12)
with ˜
ωi j =
2
α j e−d(vi ,v j ) ∑vl ∈Vi∗ αl
/4τ
˜ i ,vl )2 /4τ . e−d(v
(13)
It is easy to verify that ∑v j ∈Vi∗ ωi j = 1. To smooth anisotropically the mesh M2 , we consider that u(v,t) in (12) represents the evolution of the vertex v at time t, that is u(v,t) = v(t). We can thus filter anisotropically each vertex vi of M2 according to the update rule: vi (tk+1 ) ←
∑
v j ∈Vi∗
ωi j v j (tk ),
(14)
where τ is the time of the diffusion. The smoothing process is summarized in the algorithm 1. Basically, the smooth-
Anisotropic Smoothing(M2 ,n) for k ← 0 to n do for i ← 0 to m do Estimate the diffusion tensor Σ˜ (3) at vi . Update the vertex vi according to (14). end end Algorithm 1: The probabilistic smoothing algorithm. The parameters are: n the number of smoothing iterations and m the number of vertices in M2 . ing method described by the algorithm 1 adjusts a vertex vi successively in short time steps. At each smoothing iteration, the new position of vi is determined by a weighted average of the old positions of vertices in the set Vi∗ . The contribution of a vertex v j is controlled by the Riemannian distance d˜ between v j and vi . In particular, the larger the distance, the less the contribution. Therefore, if vi lies on a sharp edge, the neighboring vertices on the edge has much greater influence on the new location of vi . Intuitively, this is how sharp features will be preserved. On the other hand, noise may also produce sharp features, but they tend be appear equally likely around vi that their effects cancel out. This intuitively explains why smoothing can be achieved. The eigenvalues (k1 , k2 ) and eigenvectors (e1 , e2 ) of the shape operator that appear in the diffusion tensor (3) are estimated using the approach proposed by [EZ08], and the conductance function G is the Cauchy weight function given by ( 1 if | s |≤ λ, G(s) = (15) ( 1+(1s )2 otherwise, λ
where λ plays the role of the threshold parameter, such that the vertices on the edges are characterized by k1 > λ or k2 > λ, and the vertices on corners are characterized by k1 > λ and k2 > λ. c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
(a)
(b)
(a) Noisy object
1361
(b) PSM
Figure 2: Original Objects used for experimentations: (a) fandisk, (b) Venus head.
4. Results and discussion In this section, we present comparative results of the proposed probabilistic smoothing method PSM with the anisotropic geometric diffusion method AGD [CDR04], the mean curvature flow method MCF [MDSB03], the prescribed mean curvature flow method PMC [HP04] and the non-iterative smoothing method NIM [JDD03]. To test the robustness of the smoothing methods to reduce the noise, we generate a zero-mean Gaussian noise where the standard deviation σ is proportional to the average mesh edge length e, ¯ i.e. σ = α e. ¯ The noise is then added to each vertex component along the normal and tangential directions. We use two objects in our comparison: the fandisk and the Venus head objects (Figure 2). For these objects, Table 1 presents timing results and parameter settings used for our method and the comparative smoothing methods. To be objective in our comparison’s results with the other smoothing approaches, we used the public implementations of the smoothing approaches [MDSB03], [HP04] and [JDD03]. Also, we tried to choose parameter settings producing the best results. In what follows, we present a qualitative and quantitative comparison of the five smoothing methods.
(c) AGD
(d) MCF
(e) PMC
(f) NIM
Figure 3: Comparison of the smoothing methods on the fandisk object.
Qualitative evaluation Figure 3 presents the smoothing results for the synthetic fandisk object characterized by curved sharp edges and a vanishing ridge. To this object, we added a zero-mean Gaussian noise with a standard deviation σ = 0.1 e. ¯ The MCF, PMC smoothing methods reduce the noise in the non characteristic regions but fail to smooth adequately some edges and corners, this is probably due to the fact that the smoothing process is performed only along the normal vector at each vertex, which allows to reduce the noise along the normal vector while the edges and corners are not adequately smoothed. The AGD method was able to reduce a noise but also damaged some corners and edges. This is essentially due to the estimation of the curvature tensor that determines the obc 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
ject features by fitting quadratically a local surface, since the quadratic local surface interpolation does not ensure a good representation of the object features. In addition, the higher order derivatives of the quadratic local surface represent a factor of noise enhancing. The NIM preserves adequately the characteristics during the smoothing process. However, the object is under-smoothed. This is closely related to the determination of the optimal value of the two smoothing parameters σs and σg . For the proposed method PSM, it can be clearly seen that due to its anisotropic directional smoothing of the object features, it is able to faithfully preserve crucial curved sharp edges and corners. Figure 4 shows smoothing results for Venus head object, which is a real scanned object.
1362
(a) Noisy object
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
(a)
(b)
(c)
(d)
(b) PSM
Figure 5: Enlarged view of the fandisk and Venus head objects after smoothing iteration by the PSM method. (a) - (c) Original objects, (b) - (d) smoothed by the PSM method.
Quantitative evaluation (c) AGD
(d) MCF
(e) PMC
(f) NIM
Figure 4: Comparison of the smoothing methods on the Venus head object.
To which, we added a zero-mean Gaussian noise, where the standard deviation was also set to σ = 0.1e. ¯ It can be seen that the PSM method more effectively removes noise while preserving critical features. In Figure 5, we depict an enlarged view of the fandisk and Venus head objects in order to show the better performance of the proposed method. In particular, the sharp features of the fandisk object and the mouth and the engravings on Venus head are very well preserved by our method.
Let Ω be a reference 3D object (assumed to be sufficiently dense) and let Ωk be the object obtained from a reference object Ω by adding noise and applying k iterations of a smoothing process. To quantify the better performance of the proposed approach in comparison with the AGD, MCF and PMC smoothing methods, we compute the vertex-position error metric Ev (Ω, Ωk ) proposed by [BO03]. We plot and analyze the graph of the error metric for the smoothing results applied to the fandisk and Venus head objects, setting the maximum number of smoothing iterations at 40. According to the error measure graphs in Figure 6, we see that the results converge to optimally smoothed objects with minimal error value between 4 and 10 iterations. Beyond these optimal results, the errors are growing and the objects are over smoothed. We see that the proposed method has a better performance than the other smoothing methods, and so, it gives the best results indicating the consistency with the qualitative comparisons. We have applied our smoothing technique to filter a large real-world scanned object which is naturally corrupted with a noise. Figures 1 and 7 shows the smoothing result for the Chinese dragon and Gargoyle part objects. The timing and the parameters setting are reported in Table 4. Note that sharp features are well preserved while the noise is reduced. In all the experiments, we observe that the proposed technique is able to reduce the noise while preserving important geometric structure of the objects, in a very fast and efficient way. Discussion The proposed method works better than previous approaches because an exact short time solution of the anisotropic heat c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
1363
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
Objects
fandisk object
Vertices error measure Ev
0.06
fandisk
PSM MCF AGD PMC
0.05
0.04
Venus head
0.03
0.02
0.01
0
0
5
10
15
20
25
30
35
40
Iterations (a)
Venus head object
Vertices error measure Ev
0.4
PSM MCF AGD PMC
0.35
0.3
Methods AGD PMCF MCF NIM PSM AGD PMCF MCF NIM PSM
0.2
6 7 10 10 6
λ1 0.005 0.01 0.01 2.7 0.01 0.1 0.01 0.01 2.5 0.5
λ2 0.85 0.055 0.4 0.85 10 0.45 1 1.25
Time 10s 6s 4s 11s 2s 50s 18s 10s 56s 5s
Table 1: Parameter setting and timing results. Here n represents the number of iterations. For the AGD and our method, λ1 and λ2 denote the time step and the threshold value for the principal curvature, respectively. For the PMC, λ1 and λ2 denotes the step width (explicit scheme) and the feature’s detection parameter, respectively. For the MCF, λ1 represents the step width (implicit scheme). For the NIM, λ1 and λ2 denote the step widths (σs ) and the influence weight Gaussian (σg ), respectively. Object Chinese dragon Gargoyle part
0.25
n 6 10 10
n 6 6
λ1 0.2 0.1
λ2 0.85 0.85
Time 239s 21s
Table 2: Parameter setting and timing results for smoothing the real scanned models.
0.15
0.1
0.05
0
5
10
15
20
25
30
35
40
Iterations (b)
Figure 6: Evolution of the vertices error measure in terms of the smoothing iterations.
(a)
(b)
Figure 7: Smoothing a scanned model: On the left is the input Gargoyle part model and on the right is the result of the proposed smoothing method
c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
diffusion problem (4) is discretized on triangular mesh to smooth noisy objects. In addition, the method smoothes the edges tangentially ("along themselves"), which allows to reduce the noise along the edges and preserves the corners adequately. The parameter λ in (15) determines also the amount of smoothing along the sharp features. For instance, smoothing a geometrical object with large values of λ tends to over smooth the object because the sharp features are not detected, while smoothing with small values of λ will tends to under smooth the object, since the artificial edges introduced by the noise are preserved. In this work, a short time step size assumption is essentially motivated by a theoretical and practical requirements : 1) Theoretically, a short time step allows to drive a simplified approximation of the transition probability function (5) fundamental solution of the FPK problem (4). 2) Practically, the short time step permits to reduce the effect of the area and the volume shrinkage problem during the smoothing process. By taking a short time step size to smooth noisy objects in our experimentations, the objects error measure in Figure 6 that reflects implicitly the deformation of the smoothed objects is minimal and stable during the smoothing iterations by the proposed smoothing approach. If a large amount of smoothing is desired, large time step values should be used, in such as case, the algorithm remains unconditionally stable
1364
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
because at each smoothing iteration, each vertex is updated by means of normalized weight sum of the vertices in its neighborhood. As a linear combination of the bounded real values by a normalized weight is a bounded real value. 5. Conclusion In this paper, we presented a novel probabilistic method for 3D object smoothing. The core idea behind the proposed method is to relate the problem of smoothing objects to that of tracking the transition probability density functions of an underlying stochastic process. We showed that such an approach allows additional insight and sufficient flexibility compared to the standard smoothing techniques. The experimental results clearly showed a much improved performance of the proposed method in comparison with the current approaches used in 3D objects smoothing. For future work, we propose to incorporate such probabilistic framework to smooth a 3D object with texture. Acknowledgements The completion of this research was made possible thanks to Natural Sciences and Engineering Research Council of Canada (NSERC), Bell Canada’s support through its Bell University Laboratories RδD program. The authors wish to thank anonymous reviweres for many helpful comments. Chinese dragon and Gargoyle part objects are courtesy of AIM@SHAPE repository. References [BK04] BAO Y., K RIM H.: Smart nonlinear diffusion: A probabilistic approach. IEEE Trans. Pattern Anal. Mach. Intell. 26, 1 (2004), 63–72. [BO03] B ELYAEV A., O HTAKE Y.: A comparison of mesh smoothing methods. In Israel-Korea BiNational Conference on Geometric Modeling and Computer Graphics (2003), pp. 83–87. [BPK∗ 07] B OTSCH M., PAULY M., KOBBELT L., A L LIEZ P., L ÉVY B., B ISCHOFF S., R ÖSSL C.: Geometric modeling based on polygonal meshes. In SIGGRAPH Course Notes (2007), ACM. revised course notes. [BX03] BAJAJ C. L., X U G.: Anisotropic diffusion of surfaces and functions on surfaces. ACM Trans. Graph 22, 1 (2003), 4–32. [CDR04] C LARENZ U., D IEWALD U., RUMPF M.: Processing textured surfaces via anisotropic geometric diffusion. IEEE Trans. Image Process. 13, 2 (2004), 248–261. [Dod83] D ODZIUK J.: Maximum principle for parabolic inequalities and the heat flow on open manifolds. Indiana ˝ Univ. Math. J. 32, 5 (1983), 703U716.
[FDCO03] F LEISHMAN S., D RORI I., C OHEN -O R D.: Bilateral mesh denoising. ACM Trans. Graph. 22, 3 (2003), 950–953. [Gar96] G ARDINER C. W.: Handbook of Stochastic Methods: For Physics, Chemistry and the Natural Sciences. Springer, November 1996. [HP04] H ILDEBRANDT K., P OLTHIER K.: Anisotropic filtering of non-linear surface features. Comput. Graph. Forum 23, 3 (2004), 391–400. [Hsu02] H SU E. P.: Stochastic Analysis on Manifolds. American Mathematical Society, 2002. [IW81] I KEDA N., WATANABE S.: Stochastic differential equations and diffusion processes. North-Holland Publishing Co., Amsterdam, 1981. [JDD03] J ONES T., D URAND F., D ESBRUN M.: Noniterative, feature-preserving mesh smoothing. In ACM SIGGRAPH ’03 (2003), pp. 943–949. [Kam86] K AMPEN N. G. V.: Brownian motion on a manifold. Journal of Statistical Physics 44, 1/2 (1986), 1–24. [KP92] K LOEDEN P. E., P LATEN E.: Numerical Solution of Stochastic Differential Equations. Springer, November 1992. [MDSB03] M EYER M., D ESBRUN M., S CHRÖDER P., BARR A. H.: Discrete differential-geometry operators for triangulated 2-manifolds. In Visualization and Mathematics III. 2003, pp. 35–57. [Mol75] M OLCHANOV S. A.: Diffusion processes and riemannian geometry. RUSS. MATH. SURV. 30, 1 (1975), 1–63. [Mor04] M ORSE D.: Theory of constrained brownian motion. Advances in Chemical Physics 128 (2004), 65–189. [PM90] P ERONA P., M ALIK J.: Scale-space and edge detection using anisotropic diffusion. IEEE Trans. Pattern Anal. Mach. Intell 12, 7 (1990), 629–639. [PSZ01] P ENG J., S TRELA V., Z ORIN D.: A simple algorithm for surface denoising. In IEEE Vis. ’01, (2001), pp. 107–112. [TWBO03] TASDIZEN T., W HITAKER R., B URCHARD P., O SHER S.: Geometric surface processing via normal maps. ACM Trans. Graph. 22, 4 (2003), 1012–1033. [Wei94] W EICKERT J.: The thoretical foundations of anisotropic diffusion in image processing. In Theoretical Foundations of Computer Vision (1994), pp. 221–236. [YBS06] YOSHIZAWA S., B ELYAEV A. G., S EIDEL H.P.: Smoothing by example: Mesh denoising by averaging with similarity-based weights. In IEEE Shape modelling international (2006).
[EZ08] E LOUAFDI A. F., Z IOU D.: A global physical method for manifold smoothing. In IEEE Shape modelling international (2008). c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
Volume 27 (2008), Number 5
A smart stochastic approach for manifolds smoothing A. F. El Ouafdi†1 , D. Ziou1 and H. Krim2 1
Départment d’informatique, Université de Sherbrooke, Québec, Canada Engineering Department, North Carolina State University, USA
2 Computer
Abstract In this paper, we present a probabilistic approach for 3D object’s smoothing. The core idea behind the proposed method is to relate the problem of smoothing objects to that of tracking the transition probability density functions of an underlying random process. We show that such an approach allows for additional insight and sufficient flexibility compared with existing standard smoothing techniques. In particular, we are able to propose a newer, faster, and simpler smoothing approach that retains and enhances important manifold features. Furthermore, it is demonstrated to improve performance over existing smoothing techniques.
1. Introduction Recent developments in acquisition technology make it possible to obtain objects with highly detailed features. During all stages of the object’s acquisition process, noise is introduced by the inaccuracy of a sensor’s calibration, imprecision in registration procedure, limited sampling resolution etc. The noise is usually described by a probabilistic model and the degradation often yields the resulting observation object, most commonly presented in additive form U = U0 + η
(1)
where U is the noisy object, U0 is the corresponding noisefree object, and η is noise mostly assumed to be Gaussian. The challenge facing the smoothing process consists of recovering the underlying object U0 from U while preserving the object’s crucial features such as edges and corners. 1.1. Related work In recent years, various approaches have been proposed to tackle the problem of object smoothing while preserving crucial characteristics. Most smoothing methods can be classified into two main categories, depending on whether the smoothing process is performed under a probabilistic or deterministic framework. Both categories are tightly connected † Corresponding authors. [email protected], [email protected]
E-mail:
c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation ° Published by Blackwell Publishing, 9600 Garsington Road, Oxford OX4 2DQ, UK and 350 Main Street, Malden, MA 02148, USA.
(a)
(b)
Figure 1: Smoothing a scanned model: On the left is the input Chinese dragon model and on the right is the result of the proposed smoothing method.
since they are based directly or indirectly on the heat diffusion principle. In what follows, we briefly review some examples of smoothing methods in both categories. For a recent review, we refer the interested reader to [BPK∗ 07]. Deterministic methods: In this category, one can distinguish two families of smoothing approaches that either discretize a PDE expressing the heat diffusion principle or perform the smoothing process along the normal vector. The key step in the PDE based smoothing method is the PDE discretization step usually performed by the finite element FE , finite volume FV, or finite difference FD methods. Such discretization procedures replace the PDE by a linear system
1358
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
of algebraic equations. An implicit or semi-implicit time discretization is generally adopted to achieve a stable numerical schema. In this setting, [CDR04] discretize the PDE of the heat diffusion by a FE method and estimate the diffusion tensor by fitting the local surface quadratically. Such a smoothing method enhances the triangular surface in the principal directions of curvature. Bajaj and Xu [BX03] generalize the work of Claranz et al. to a higher dimension and combine a limit function of loop subdivision with a PDE to smooth higher order functions on the surface while fairing it. Using a level set method, Tasdizen et al. [TWBO03] propose a twostep smoothing method by solving simultaneously two PDEs on level set surfaces. Very recently, [EZ08] decompose the heat diffusion principle into basic laws and derive the numerical schema which allows one to perform the smoothing process by discretizing the basic laws using computation algebraic topological tools. In the particular case when the diffusion tensor is a scalar and locally constant function, the PDE that expresses the heat diffusion principle is reduced to the mean curvature evolution equation along the normal vector. Along this vein, the mean curvature flow smoothing method proposed by [MDSB03] adjusts the vertex position along the surface normal with speed proportional to the mean curvature flow. The important features are preserved by a weighting parameter which depends on the principal curvature values. [HP04] propose a prescribed mean curvature flow method by preconditioning the anisotropic mean curvature vector. The diffusion tensor is estimated by a normal cycle approach. The recent work of [YBS06] extended the non-local image filter to meshes by computing a local radial basis function approximation to define the similarity measure. Probabilistic methods: Due to the nature of the noise, the smoothing process may also be interpreted from a statistical point of view. The noise is considered as a stochastic process with known parameters. In this setting, [PSZ01] suggest an adaptative Wiener filtering smoothing technique for semi-regular meshes. A pre-processing step is needed to convert the meshes to multi-resolution representations before smoothing. Based on robust statistics and local statistical predictors of surfaces, [JDD03] propose a non-iterative smoothing method. The surface features are preserved by means of a robust influence weight function together with a predictor for vertex positions based on the surface’s tangent plane. A similar statistical predictor also appears in the bilateral mesh smoothing approach proposed by [FDCO03]. Recently, [BK04] provide a stochastic interpretation of the standard [PM90] diffusion method for image filtering. Also they propose a two-sided gradient edges detector to filter a noisy input signal while preserving important features. Intuitively, the authors model the nonlinear diffusion by a differential stochastic equation SDE. Next they construct a probability density of the SDE to resolve the nonlinear or constrained diffusion problem.
1.2. Contribution The anisotropic diffusion is a powerful smoothing approach that performs the directional smoothing and enhancing of the manifold features while preserving adequately curvatures, edges and corners. To the best of our knowledge, the existing anisotropic smoothing methods discretize the PDE heat diffusion by a FE or FD methods [Wei94, BX03, CDR04]. The derivation of a stable numerical schema from the discretization of the PDE by the such as discretization methods requires the resolution of a linear system at each vertex at each iteration, which induces numerical error, time consuming and computationally expensive. Closely related to the intuitive idea of the image smoothing method proposed by [BK04]. We construct in this work an anisotropic Gaussian kernel to smooth an arbitrary manifold initially corrupted with natural or artificial noise. Such kernel is the transition probability density of an constrained diffusion process. The proposed method is simpler than solving the heat diffusion numerically. Additionally, the proposed approach is fast, unconditionally stable and able to retain and even enhance important manifold features such as edges and corners. 1.3. Paper Organization In the next section we begin with a brief description of anisotropic diffusion process on manifold, and then, we deduce a short time solution of the anisotropic diffusion problem. In Section 3, we discretise the diffusion process on triangular meshes and we propose a new algorithm for 3D objects smoothing whose evaluation through numerous substantiating examples is provided in Section 4. We finally give some concluding remarks in Section 5. 2. Diffusion on a Riemannian manifold In this section, we begin with a brief description of the constrained diffusion process from the stochastic point of view. Then, we deduce a short time solution of anisotropic diffusion on a manifold. 2.1. Constrained diffusion process on a manifold Assume that the object to be processed can be represented by an n-dimensional smooth and compact manifold M embedded in Rn+1 . The diffusion processes constitute a powerful and a fundamental tool for smoothing manifolds. They often arise as physical processes in the real world such as the heat diffusion. The diffusion process, known in physics as the transfer of energy from the more energetic to the less energetic particles, can be assimilated to the random motions of a virtual particles. For instance, the particle trajectory X(t) = y at the point y at time t is a Markov process governed by the following stochastic differential equation SDE dX(t) = µ(X(t)) dt + σ(X(t)) dB(t)
(2)
c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
1359
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
X(0) = x ∈ Ω ⊆ M t ≥ 0, where Ω is a smooth compact domain of M, B(t) is mdimensional Brownian motion that generally models the noise term in (1), µ is the drift velocity vector field and Σ = 12 σσT is the diffusion tensor assumed to be symmetric and positive definite [IW81]. For an unconstrained diffusion process, the drift velocity vector field and the diffusion tensor are proportional to the first and second moments of random changes in coordinates, respectively. To generate a diffusion path that obeys the desired diffusion behaviour, a constrained or anisotropic diffusion process is achieved by modifying the data-driven control parameter µ and Σ [Mor04]. In computer vision, several smoothing methods may be interpreted as constrained diffusion processes. For instance, the standard gradient based smoothing method proposed by [PM90] can be interpreted as a constrained diffusion process with vanishing drift velocity coefficient, and a diffusion coefficient equal to a conductance function that varies with the local gradient amplitude, which allows to preserve important sharp features such as edges in images [BK04]. More generally, most smoothing methods that perform the directional smoothing of the manifold constrain implicitly the diffusion process through a diffusion tensor based on the principal curvatures ki and the principal direction of curvatures ei , for i = 1, · · · n. Such a diffusion tensor is generally given by Σ˜ = UVU T
(3)
where U = (e0 | · · · |en |N), V is a diagonal matrix with diagonal term (G(k0 ) · · · G(kn ), 1), N is the unit normal vector, ki and ei are respectively the eigenvalues and eigenvectors of the shape operator matrix, and G is a conductance function. When constrained, the diffusion process endows implicitly the manifold with a metric equal to the inverse of the diffusion tensor (Σ−1 ) in order to constrain the particles dynamics [Kam86]. In our case, as the diffusion tensor (3) is symmetric and positive definite for a bounded eigenvalues. Its inverse is also a symmetric and definite positive. Thus, the tensor Σ˜ −1 is a natural choice to define a new Riemannian metric on M to achieve a directional smoothing of the manifold. One may notice that the inverse tensor metric on manifold induces differential operators like the gradient and Laplace-Beltrami operators which are generally encountered in the manifold smoothing methods.
In many practical problems, the determination of the transition probability function of the diffusion process X(t) is of interest [Gar96]. The transition probability function on Ω, noted pΩ (t, x, y), describes the probability for a particle starting initially at the location X(0) = x to reach X(t) = y ∈ Ω at time t. For a diffusion process on a Riemannian manifold, it can be shown (see for example [Mol75, Hsu02]), that the transition probability function pΩ (t, x, y) is a minimal positive solution of the following backward Fokker-PlanckKolmogorov FPK problem in Ω ∂u(x,t) = (∆˜ M + µ)u(x,t), (x,t) ∈ Ω × (0, ∞), ∂t u(x,t) = 0 (x,t) ∈ ∂Ω × (0, ∞)
(4)
u(x, 0) = u0 (x) x ∈ Ω, with u0 (x) = δx , where δx is the Dirac function. Here ∆˜ M denote the Laplace-Beltrami operator associated with the metric tensor Σ˜ −1 defined on the local region Ω. For a short time period and for vanishing drift term µ = 0, the transition probability function, minimal positive solution of (4), has the following asymptotic expansion [Dod83] ˜
2
pΩ (t, x, y) ≈ Φ(x, y)e−d(x,y)
/4t
(5)
and fulfills asymptotically the normalization condition Z
Ω
pΩ (t, x, y)dy = 1,
(6)
where Φ(x, y) is a continuous function on Ω × Ω such that ˜ y) = (x − y)T Σ˜ −1 (x − y) Φ(x, y) ≈ 1 for x close to y, d(x, is the Riemannian distance between x and y relative to the metric tensor Σ˜ −1 . In a general case, when the initial condition u0 (x) is a continuous function on Ω, the FPK problem (4) admits an unique solution of the form Z
u(x,t) =
Ω
pΩ (t, x, y)u0 (y)dy.
(7)
The above equation may be interpreted as an anisotropic convolution between the Gaussian kernel pΩ and the continuous function u0 . Note that the function u0 may also be interpreted as a quantity represented on the region Ω at time t = 0 (e.g., initial temperature, gray levels, colors, vertices), in which case, u(x,t) represents the anisotropic diffusion of u0 at time t; this is a behaviour similar to what is required in anisotropic smoothing. 3. Diffusion on triangular meshes
2.2. A short time solution of the anisotropic diffusion problem on manifold In general, it is difficult to find an explicit solution X(t) of SDE (2). Such a solution is only possible when (2) can be transformed into a linear equation of low dimension. In addition, an accurate numerical solution requires the generation of a large number of Brownian motion paths with a small time step, which is computationally expensive [KP92]. c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
Discretization is the basis of many numerical procedures, this also holds continuous for the diffusion process. The object of discretization is the Riemannian manifold M and the continuous time axis. For our purposes, let us review this discretization on a 2-manifold M2 embedded in R3 . To solve the diffusion process problem numerically, the 2dimensional spatial domain M2 on which the diffusion process takes place is discretized by a polygonal or triangular
1360
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
mesh M2 = (V, E, S), where V = {v1 , · · · , vm } is the set of vertices, E = {ei j } is the set of edges, and S = {S1 , · · · , Sn } is the set of triangles. Two vertices vi , v j ∈ V are adjacent( written vi ∼ v j ) if they are connected by an edge ei j ∈ E. The neighborhood of a vertex vi is the set Vi∗ = {v j ∈ V : vi ∼ v j } ∪ {vi }. We denote by Si∗ the set of all triangles incident to the vertex vi . To each vertex vi , we associate a triangular local surface Ωi = ∪S j ∈Si∗ S j formed by the union of all triangles in Si∗ . The time discretization is performed with a constant time step. We suppose that we have nt equal time intervals [tk ,tk+1 ] of length τ = tk+1 − tk . In what follows, we consider that tk represents the initial time of the diffusion. After the space and time discretization, the second step of the discretization procedure is to build from (5) a normalized transition probability function on the polygonal local surface Ωi . Thus, let us consider the truncated normalized transition probability function on Ωi given by pΩi (τ, x, y) . Ωi pΩi (τ, x, z)dz
(8)
p˜Ωi (τ, x, y) = R
The denominator term is introduced to correct the truncation error and to enforce the normalization condition (6) of p˜Ωi on Ωi . Such denominator term may be expressed at x = vi as Z
Z
Ωi
pΩi (τ, vi , z)dz =
∑
S j ∈Si∗ S j
e
˜ i ,z)2 /4τ −d(v
dz.
(9)
Using the standard integrals approximation techniques, we can approximate the integration on each triangle S j formed by the vertices v j0 , v j1 and v j2 as Z
˜
Sj
2
e−d(vi ,y)
/4τ
dy ≈
| S j | 2 −d(v ˜ i ,v j )2 /4τ r e . ∑ 3 r=0
By substituting the above expression in (9), we reformulate the normalized transition probability function p˜Ωi as ˜
p˜Ωi (τ, vi , v j ) =
2
e−d(vi ,v j ) ∑vl ∈Vi∗ αl
/4τ
˜ i ,vl )2 /4τ e−d(v
,
(10)
T
|S |
where αl = ∑ p 3p , S p ∈ Sl∗ Si∗ where Sl∗ is the set of ˜ i, v j ) = all incident triangles to the vertex vl ∈ Vi∗ and d(v T ˜ −1 (vi − v j ) Σ (vi − v j ). Having the normalized transition probability (10) at hand, we can now express the local solution u(vi ,tk+1 ) (7) at vertex vi at time tk+1 with initial condition u(y,tk ) on Ωi as Z
u(vi ,tk+1 ) =
∑
S j ∈Si∗ S j
p˜Ωi (τ, vi , y)u(y,tk )dy.
(11)
As before, using the standard integrals approximation techniques, we can rewrite u(vi ,tk+1 ) in the above expression as a normalized weighted sum of u(v j ,tk ), v j ∈ Vi∗ as u(vi ,tk+1 ) =
∑
v j ∈Vi∗
ωi j u(v j ,tk ),
(12)
with ˜
ωi j =
2
α j e−d(vi ,v j ) ∑vl ∈Vi∗ αl
/4τ
˜ i ,vl )2 /4τ . e−d(v
(13)
It is easy to verify that ∑v j ∈Vi∗ ωi j = 1. To smooth anisotropically the mesh M2 , we consider that u(v,t) in (12) represents the evolution of the vertex v at time t, that is u(v,t) = v(t). We can thus filter anisotropically each vertex vi of M2 according to the update rule: vi (tk+1 ) ←
∑
v j ∈Vi∗
ωi j v j (tk ),
(14)
where τ is the time of the diffusion. The smoothing process is summarized in the algorithm 1. Basically, the smooth-
Anisotropic Smoothing(M2 ,n) for k ← 0 to n do for i ← 0 to m do Estimate the diffusion tensor Σ˜ (3) at vi . Update the vertex vi according to (14). end end Algorithm 1: The probabilistic smoothing algorithm. The parameters are: n the number of smoothing iterations and m the number of vertices in M2 . ing method described by the algorithm 1 adjusts a vertex vi successively in short time steps. At each smoothing iteration, the new position of vi is determined by a weighted average of the old positions of vertices in the set Vi∗ . The contribution of a vertex v j is controlled by the Riemannian distance d˜ between v j and vi . In particular, the larger the distance, the less the contribution. Therefore, if vi lies on a sharp edge, the neighboring vertices on the edge has much greater influence on the new location of vi . Intuitively, this is how sharp features will be preserved. On the other hand, noise may also produce sharp features, but they tend be appear equally likely around vi that their effects cancel out. This intuitively explains why smoothing can be achieved. The eigenvalues (k1 , k2 ) and eigenvectors (e1 , e2 ) of the shape operator that appear in the diffusion tensor (3) are estimated using the approach proposed by [EZ08], and the conductance function G is the Cauchy weight function given by ( 1 if | s |≤ λ, G(s) = (15) ( 1+(1s )2 otherwise, λ
where λ plays the role of the threshold parameter, such that the vertices on the edges are characterized by k1 > λ or k2 > λ, and the vertices on corners are characterized by k1 > λ and k2 > λ. c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
(a)
(b)
(a) Noisy object
1361
(b) PSM
Figure 2: Original Objects used for experimentations: (a) fandisk, (b) Venus head.
4. Results and discussion In this section, we present comparative results of the proposed probabilistic smoothing method PSM with the anisotropic geometric diffusion method AGD [CDR04], the mean curvature flow method MCF [MDSB03], the prescribed mean curvature flow method PMC [HP04] and the non-iterative smoothing method NIM [JDD03]. To test the robustness of the smoothing methods to reduce the noise, we generate a zero-mean Gaussian noise where the standard deviation σ is proportional to the average mesh edge length e, ¯ i.e. σ = α e. ¯ The noise is then added to each vertex component along the normal and tangential directions. We use two objects in our comparison: the fandisk and the Venus head objects (Figure 2). For these objects, Table 1 presents timing results and parameter settings used for our method and the comparative smoothing methods. To be objective in our comparison’s results with the other smoothing approaches, we used the public implementations of the smoothing approaches [MDSB03], [HP04] and [JDD03]. Also, we tried to choose parameter settings producing the best results. In what follows, we present a qualitative and quantitative comparison of the five smoothing methods.
(c) AGD
(d) MCF
(e) PMC
(f) NIM
Figure 3: Comparison of the smoothing methods on the fandisk object.
Qualitative evaluation Figure 3 presents the smoothing results for the synthetic fandisk object characterized by curved sharp edges and a vanishing ridge. To this object, we added a zero-mean Gaussian noise with a standard deviation σ = 0.1 e. ¯ The MCF, PMC smoothing methods reduce the noise in the non characteristic regions but fail to smooth adequately some edges and corners, this is probably due to the fact that the smoothing process is performed only along the normal vector at each vertex, which allows to reduce the noise along the normal vector while the edges and corners are not adequately smoothed. The AGD method was able to reduce a noise but also damaged some corners and edges. This is essentially due to the estimation of the curvature tensor that determines the obc 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
ject features by fitting quadratically a local surface, since the quadratic local surface interpolation does not ensure a good representation of the object features. In addition, the higher order derivatives of the quadratic local surface represent a factor of noise enhancing. The NIM preserves adequately the characteristics during the smoothing process. However, the object is under-smoothed. This is closely related to the determination of the optimal value of the two smoothing parameters σs and σg . For the proposed method PSM, it can be clearly seen that due to its anisotropic directional smoothing of the object features, it is able to faithfully preserve crucial curved sharp edges and corners. Figure 4 shows smoothing results for Venus head object, which is a real scanned object.
1362
(a) Noisy object
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
(a)
(b)
(c)
(d)
(b) PSM
Figure 5: Enlarged view of the fandisk and Venus head objects after smoothing iteration by the PSM method. (a) - (c) Original objects, (b) - (d) smoothed by the PSM method.
Quantitative evaluation (c) AGD
(d) MCF
(e) PMC
(f) NIM
Figure 4: Comparison of the smoothing methods on the Venus head object.
To which, we added a zero-mean Gaussian noise, where the standard deviation was also set to σ = 0.1e. ¯ It can be seen that the PSM method more effectively removes noise while preserving critical features. In Figure 5, we depict an enlarged view of the fandisk and Venus head objects in order to show the better performance of the proposed method. In particular, the sharp features of the fandisk object and the mouth and the engravings on Venus head are very well preserved by our method.
Let Ω be a reference 3D object (assumed to be sufficiently dense) and let Ωk be the object obtained from a reference object Ω by adding noise and applying k iterations of a smoothing process. To quantify the better performance of the proposed approach in comparison with the AGD, MCF and PMC smoothing methods, we compute the vertex-position error metric Ev (Ω, Ωk ) proposed by [BO03]. We plot and analyze the graph of the error metric for the smoothing results applied to the fandisk and Venus head objects, setting the maximum number of smoothing iterations at 40. According to the error measure graphs in Figure 6, we see that the results converge to optimally smoothed objects with minimal error value between 4 and 10 iterations. Beyond these optimal results, the errors are growing and the objects are over smoothed. We see that the proposed method has a better performance than the other smoothing methods, and so, it gives the best results indicating the consistency with the qualitative comparisons. We have applied our smoothing technique to filter a large real-world scanned object which is naturally corrupted with a noise. Figures 1 and 7 shows the smoothing result for the Chinese dragon and Gargoyle part objects. The timing and the parameters setting are reported in Table 4. Note that sharp features are well preserved while the noise is reduced. In all the experiments, we observe that the proposed technique is able to reduce the noise while preserving important geometric structure of the objects, in a very fast and efficient way. Discussion The proposed method works better than previous approaches because an exact short time solution of the anisotropic heat c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
1363
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
Objects
fandisk object
Vertices error measure Ev
0.06
fandisk
PSM MCF AGD PMC
0.05
0.04
Venus head
0.03
0.02
0.01
0
0
5
10
15
20
25
30
35
40
Iterations (a)
Venus head object
Vertices error measure Ev
0.4
PSM MCF AGD PMC
0.35
0.3
Methods AGD PMCF MCF NIM PSM AGD PMCF MCF NIM PSM
0.2
6 7 10 10 6
λ1 0.005 0.01 0.01 2.7 0.01 0.1 0.01 0.01 2.5 0.5
λ2 0.85 0.055 0.4 0.85 10 0.45 1 1.25
Time 10s 6s 4s 11s 2s 50s 18s 10s 56s 5s
Table 1: Parameter setting and timing results. Here n represents the number of iterations. For the AGD and our method, λ1 and λ2 denote the time step and the threshold value for the principal curvature, respectively. For the PMC, λ1 and λ2 denotes the step width (explicit scheme) and the feature’s detection parameter, respectively. For the MCF, λ1 represents the step width (implicit scheme). For the NIM, λ1 and λ2 denote the step widths (σs ) and the influence weight Gaussian (σg ), respectively. Object Chinese dragon Gargoyle part
0.25
n 6 10 10
n 6 6
λ1 0.2 0.1
λ2 0.85 0.85
Time 239s 21s
Table 2: Parameter setting and timing results for smoothing the real scanned models.
0.15
0.1
0.05
0
5
10
15
20
25
30
35
40
Iterations (b)
Figure 6: Evolution of the vertices error measure in terms of the smoothing iterations.
(a)
(b)
Figure 7: Smoothing a scanned model: On the left is the input Gargoyle part model and on the right is the result of the proposed smoothing method
c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
diffusion problem (4) is discretized on triangular mesh to smooth noisy objects. In addition, the method smoothes the edges tangentially ("along themselves"), which allows to reduce the noise along the edges and preserves the corners adequately. The parameter λ in (15) determines also the amount of smoothing along the sharp features. For instance, smoothing a geometrical object with large values of λ tends to over smooth the object because the sharp features are not detected, while smoothing with small values of λ will tends to under smooth the object, since the artificial edges introduced by the noise are preserved. In this work, a short time step size assumption is essentially motivated by a theoretical and practical requirements : 1) Theoretically, a short time step allows to drive a simplified approximation of the transition probability function (5) fundamental solution of the FPK problem (4). 2) Practically, the short time step permits to reduce the effect of the area and the volume shrinkage problem during the smoothing process. By taking a short time step size to smooth noisy objects in our experimentations, the objects error measure in Figure 6 that reflects implicitly the deformation of the smoothed objects is minimal and stable during the smoothing iterations by the proposed smoothing approach. If a large amount of smoothing is desired, large time step values should be used, in such as case, the algorithm remains unconditionally stable
1364
A. F. El Ouafdi, D. Ziou & H. Krim / A smart stochastic approach for manifolds smoothing
because at each smoothing iteration, each vertex is updated by means of normalized weight sum of the vertices in its neighborhood. As a linear combination of the bounded real values by a normalized weight is a bounded real value. 5. Conclusion In this paper, we presented a novel probabilistic method for 3D object smoothing. The core idea behind the proposed method is to relate the problem of smoothing objects to that of tracking the transition probability density functions of an underlying stochastic process. We showed that such an approach allows additional insight and sufficient flexibility compared to the standard smoothing techniques. The experimental results clearly showed a much improved performance of the proposed method in comparison with the current approaches used in 3D objects smoothing. For future work, we propose to incorporate such probabilistic framework to smooth a 3D object with texture. Acknowledgements The completion of this research was made possible thanks to Natural Sciences and Engineering Research Council of Canada (NSERC), Bell Canada’s support through its Bell University Laboratories RδD program. The authors wish to thank anonymous reviweres for many helpful comments. Chinese dragon and Gargoyle part objects are courtesy of AIM@SHAPE repository. References [BK04] BAO Y., K RIM H.: Smart nonlinear diffusion: A probabilistic approach. IEEE Trans. Pattern Anal. Mach. Intell. 26, 1 (2004), 63–72. [BO03] B ELYAEV A., O HTAKE Y.: A comparison of mesh smoothing methods. In Israel-Korea BiNational Conference on Geometric Modeling and Computer Graphics (2003), pp. 83–87. [BPK∗ 07] B OTSCH M., PAULY M., KOBBELT L., A L LIEZ P., L ÉVY B., B ISCHOFF S., R ÖSSL C.: Geometric modeling based on polygonal meshes. In SIGGRAPH Course Notes (2007), ACM. revised course notes. [BX03] BAJAJ C. L., X U G.: Anisotropic diffusion of surfaces and functions on surfaces. ACM Trans. Graph 22, 1 (2003), 4–32. [CDR04] C LARENZ U., D IEWALD U., RUMPF M.: Processing textured surfaces via anisotropic geometric diffusion. IEEE Trans. Image Process. 13, 2 (2004), 248–261. [Dod83] D ODZIUK J.: Maximum principle for parabolic inequalities and the heat flow on open manifolds. Indiana ˝ Univ. Math. J. 32, 5 (1983), 703U716.
[FDCO03] F LEISHMAN S., D RORI I., C OHEN -O R D.: Bilateral mesh denoising. ACM Trans. Graph. 22, 3 (2003), 950–953. [Gar96] G ARDINER C. W.: Handbook of Stochastic Methods: For Physics, Chemistry and the Natural Sciences. Springer, November 1996. [HP04] H ILDEBRANDT K., P OLTHIER K.: Anisotropic filtering of non-linear surface features. Comput. Graph. Forum 23, 3 (2004), 391–400. [Hsu02] H SU E. P.: Stochastic Analysis on Manifolds. American Mathematical Society, 2002. [IW81] I KEDA N., WATANABE S.: Stochastic differential equations and diffusion processes. North-Holland Publishing Co., Amsterdam, 1981. [JDD03] J ONES T., D URAND F., D ESBRUN M.: Noniterative, feature-preserving mesh smoothing. In ACM SIGGRAPH ’03 (2003), pp. 943–949. [Kam86] K AMPEN N. G. V.: Brownian motion on a manifold. Journal of Statistical Physics 44, 1/2 (1986), 1–24. [KP92] K LOEDEN P. E., P LATEN E.: Numerical Solution of Stochastic Differential Equations. Springer, November 1992. [MDSB03] M EYER M., D ESBRUN M., S CHRÖDER P., BARR A. H.: Discrete differential-geometry operators for triangulated 2-manifolds. In Visualization and Mathematics III. 2003, pp. 35–57. [Mol75] M OLCHANOV S. A.: Diffusion processes and riemannian geometry. RUSS. MATH. SURV. 30, 1 (1975), 1–63. [Mor04] M ORSE D.: Theory of constrained brownian motion. Advances in Chemical Physics 128 (2004), 65–189. [PM90] P ERONA P., M ALIK J.: Scale-space and edge detection using anisotropic diffusion. IEEE Trans. Pattern Anal. Mach. Intell 12, 7 (1990), 629–639. [PSZ01] P ENG J., S TRELA V., Z ORIN D.: A simple algorithm for surface denoising. In IEEE Vis. ’01, (2001), pp. 107–112. [TWBO03] TASDIZEN T., W HITAKER R., B URCHARD P., O SHER S.: Geometric surface processing via normal maps. ACM Trans. Graph. 22, 4 (2003), 1012–1033. [Wei94] W EICKERT J.: The thoretical foundations of anisotropic diffusion in image processing. In Theoretical Foundations of Computer Vision (1994), pp. 221–236. [YBS06] YOSHIZAWA S., B ELYAEV A. G., S EIDEL H.P.: Smoothing by example: Mesh denoising by averaging with similarity-based weights. In IEEE Shape modelling international (2006).
[EZ08] E LOUAFDI A. F., Z IOU D.: A global physical method for manifold smoothing. In IEEE Shape modelling international (2008). c 2008 The Author(s) ° c 2008 The Eurographics Association and Blackwell Publishing Ltd. Journal compilation °
