# Integration of a Normal Field without Boundary Condition

Jean-Denis Durou, Frédéric Courteille

► To cite this version:

Jean-Denis Durou, Frédéric Courteille. Integration of a Normal Field without Boundary Condition. 1st International Workshop on Photometric Analysis For Computer Vision - PACV 2007, Oct 2007, Rio de Janeiro, Brazil. pp.1-8. inria-00264852

HAL Id: inria-00264852

https://inria.hal.science/inria-00264852

Submitted on 18 Mar 2008

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# Integration of a Normal Field without Boundary Condition

Jean-Denis DUROU and Frédéric COURTEILLE

IRIT, UMR CNRS 5505, Toulouse, France

{durou,courteil}@irit.fr

## Abstract

We show how to use two existing methods of integration of a normal field in the absence of boundary condition, which makes them more realistic. Moreover, we show how perspective can be taken into account, in order to render the 3D-reconstruction more accurate. Finally, the joint use of both these methods of integration allows us to obtain very satisfactory results, from the point of view of CPU time as well as that of the accuracy of the reconstructions. As an application, we use this new combined method of integration of a normal field in the framework of photometric stereo, a technique which aims at computing a normal field to the surface of a scene from several images of this scene illuminated from various directions. The performances of the proposed method are illustrated on synthetic, as well as on real images.

## 1. Introduction

Photometric stereo is a computer vision technique which uses several images of a scene illuminated from various directions. Usually, this technique consists in two stages: the normal to the surface is first computed in each visible point [15]; the normal field is then integrated, using either the calculus of variations [6], direct integration [17] or frequency-domain methods [5]. In this paper, we take an interest in the second stage. More precisely, we improve two existing methods of integration in two ways: on the one hand, we show that the knowledge of the height on the boundary, a knowledge which is usually not available, is not necessary; on the other hand, we show how to extend these methods to perspective projection. This theoretical study is propounded in Section 2. In Section 3, a new method of normals integration is designed from this theoretical study, which combines both improved methods. In order to validate this new method, we test it in the framework of photometric stereo, on synthetic, as well as on real images (Section 4). In Section 5, we summarize the main contributions of the paper.

## 2. Relation between Normal and Gradient

Due to lack of space, no state-of-the-art on the integration of a normal field is given here (see e.g. [8, 10, 7, 1]). Suppose that, in each point $Q = (x, y)$ in the image of a surface $S$, we know the unit outgoing normal $\vec{n}(x, y)$ to $S$:

$$ \vec{n}(x, y) = \begin{bmatrix} n_X(x, y) \\ n_Y(x, y) \\ n_Z(x, y) \end{bmatrix}. \qquad (1) $$

The vectorial function $\vec{n}$ is called a “normal field”. The problem of integrating a normal field consists in searching for a shape $S$ i.e., for three functions $X, Y$ and $Z$ ($Z$ is called the “height”), such that the object point $P$ conjugated with $Q$ has $X(x, y)$, $Y(x, y)$ and $Z(x, y)$ as coordinates. It can be solved only if the model of projection is known.

### 2.1. Orthographic Projection

<figure><img src="image_2.png" /><figcaption>Figure 1. Two models of projection: orthographic (points $P_1$ and $Q_1$) and perspective (points $P_2$ and $Q_2$).</figcaption></figure>

Under the assumption of orthographic projection (cf. points $P_1$ and $Q_1$ in Fig. 1), $X(x, y)$ and $Y(x, y)$ are worth:

$$ \begin{cases} X(x, y) = x d/f, \\ Y(x, y) = y d/f. \end{cases} \qquad (2) $$

Let us denote $g$ the transverse magnification $f/d$. The product of both partial derivatives $\partial_x \vec{P}$ and $\partial_y \vec{P}$ is normal to
the surface (the dependences in $(x, y)$ are omitted, for seek of simplicity):

$$ \partial_x \vec{P} \wedge \partial_y \vec{P} = \begin{bmatrix} -Z_x/g \\ -Z_y/g \\ 1/g^2 \end{bmatrix}, \qquad (3) $$

where $Z_x$ and $Z_y$ denote the partial derivatives of function $Z$. The product of this vector by $\vec{n}(x, y)$ is null, which gives us the three following equations:

$$ \left\{ \begin{array}{l} g n_Z Z_x = -n_X, \\ g n_Z Z_y = -n_Y, \\ n_Y Z_x - n_X Z_y = 0. \end{array} \right. \qquad (4) $$

System (4), which is linear in $Z_x$ and $Z_y$, is ill-constrained if the rank of its matrix is less than 2 i.e., if $n_Z^2 = 0$ and $n_X n_Z = 0$ and $n_Y n_Z = 0$, which occurs only if $n_Z = 0$. This particular case must be considered carefully, since $\vec{n}(x, y)$ is parallel to the image plane $\Pi_i$ in this case. Under the assumption of orthographic projection, the image points which lie on an “occluding contour” satisfy this property. Thus, even if the normal is easy to determine for the image points which lie on an occluding contour (the normal is parallel to $\Pi_i$ and is also normal to the occluding contour), system (4) is ill-constrained and $Z$ cannot be computed for such points.

Now let us consider the image points which do not lie on an occluding contour. System (4) is well-constrained and its resolution gives us the gradient of $Z$:

$$ \nabla Z = 1/g \begin{bmatrix} p \\ q \end{bmatrix}, \qquad (5) $$

where the usual notations $p = -n_X/n_Z$ and $q = -n_Y/n_Z$ have been used. Thus, the problem of integrating a normal field has been reformulated as a problem of integrating the gradient of a function of two variables, whose solution is straightforward:

$$ Z(x, y) = Z(x_0, y_0) + \int_{(x_0, y_0)}^{(x, y)} \frac{p(u, v)du + q(u, v)dv}{g}. \qquad (6) $$

It results from this expression that $Z$ is computed up to an additive constant. Thus, the Root Mean Square Error in $Z$, denoted $|\Delta Z|_2$, depends on this constant. In the tests, as we will see, no boundary condition is imposed, thus this constant is computed so that $|\Delta Z|_2$ be minimal.

## 2.2. Perspective Projection

Under the assumption of perspective projection (cf. points $P_2$ and $Q_2$ in Fig. 1), $X(x, y)$ and $Y(x, y)$ are worth:

$$ \left\{ \begin{array}{l} X(x, y) = x Z(x, y)/f, \\ Y(x, y) = y Z(x, y)/f. \end{array} \right. \qquad (7) $$

The product of $\partial_x \vec{P}$ and $\partial_y \vec{P}$ is a little more complicated than in the orthographic case:

$$ \partial_x \vec{P} \wedge \partial_y \vec{P} = Z/f \begin{bmatrix} -f Z_x \\ -f Z_y \\ Z + x Z_x + y Z_y \end{bmatrix}. \qquad (8) $$

Knowing that this vector is parallel to $\vec{n}(x, y)$, this gives us the three following equations:

$$ \left\{ \begin{array}{l} f n_Z Z_x + n_X [Z + x Z_x + y Z_y] = 0, \\ f n_Z Z_y + n_Y [Z + x Z_x + y Z_y] = 0, \\ n_Y Z_x - n_X Z_y = 0. \end{array} \right. \qquad (9) $$

System (9) is not linear in $Z_x$ and $Z_y$, but homogeneous in $Z$. Thus, it is useful to introduce a change in the unknown:

$$ T = \ln |Z|. \qquad (10) $$

We obtain a new system which is linear in $T_x$ and $T_y$:

$$ \left\{ \begin{array}{l} [f n_Z + x n_X] T_x + y n_X T_y = -n_X, \\ x n_Y T_x + [f n_Z + y n_Y] T_y = -n_Y, \\ n_Y T_x - n_X T_y = 0. \end{array} \right. \qquad (11) $$

This system is ill-constrained if the rank of its matrix is less than 2 i.e., if:

$$ \left\{ \begin{array}{l} f n_Z [x n_X + y n_Y + f n_Z] = 0, \\ -n_X [x n_X + y n_Y + f n_Z] = 0, \\ -n_Y [x n_X + y n_Y + f n_Z] = 0. \end{array} \right. \qquad (12) $$

As $n_X$, $n_Y$ and $n_Z$ cannot simultaneously be equal to 0, since $\vec{n}$ is a unit vector, system (12) holds only if $x n_X + y n_Y + f n_Z = 0$. Moreover, $x n_X + y n_Y + f n_Z = \vec{C} \vec{Q} \cdot \vec{n}$, then system (12) holds only if the image point $Q$ lies on an occluding contour, as in the orthographic case.

Now let us consider the image points which do not lie on an occluding contour. System (11) is well-constrained and its resolution gives us the gradient of $T$:

$$ \nabla T = \begin{bmatrix} r \\ s \end{bmatrix}, \qquad (13) $$

where the following notations are used:

$$ \left\{ \begin{array}{l} r = - \frac{n_X}{x n_X + y n_Y + f n_Z}, \\ s = - \frac{n_Y}{x n_X + y n_Y + f n_Z}. \end{array} \right. \qquad (14) $$

Once again, the problem of integrating a normal field has been reformulated as a problem of integrating the gradient of a function of two variables, whose solution is straightforward:

$$ T(x, y) = T(x_0, y_0) + \int_{(x_0, y_0)}^{(x, y)} [r(u, v)du + s(u, v)dv]. \qquad (15) $$
From (15) and (10), we deduce:

$$Z(x, y) = Z(x_0, y_0) \exp \left\{ \int_{(x_0, y_0)}^{(x, y)} [r(u, v)du + s(u, v)dv] \right\}. \quad (16)$$

It results from this expression that Z is computed up to a multiplicative constant. It is also noticeable that (16) requires the knowledge of the focal distance f, as well as the location of the principal point O (since the coordinates x and y depend on it).

## 2.3. Integrability of a Normal Field

In order to ensure that the normal field is integrable i.e., that Eqs. (6) or (16) are independent of the integration path, it is necessary and sufficient that p and q (in the orthographic case) or r and s (in the perspective case) satisfy the Schwartz equations $\partial p/\partial y = \partial q/\partial x$ or $\partial r/\partial y = \partial s/\partial x$. In practice, a normal field is never rigorously integrable. There are two answers to this problem. The first one consists in using several integration paths between $(x_0, y_0)$ and $(x, y)$, and to mean the integrals, as Wu and Li do (see [17] and Section 3.1). The second answer considers that Eqs. (5) or (13) are least square problems, as Horn and Brooks do (see [6] and Section 3.2). Now let us give the continuous formulation of this second answer, which uses the calculus of variations.

Under the assumption of orthographic projection, the resolution of Eq. (5) in the least square sense amounts minimizing the following functional:

$$\mathcal{F}(Z) = \iint_{(x,y) \in \Omega} F(x, y, Z, Z_x, Z_y) \, dx \, dy, \quad (17)$$

where $\Omega$ denotes the “domain of reconstruction”, and:

$$F(x, y, Z, Z_x, Z_y) = [Z_x - p/g]^2 + [Z_y - q/g]^2, \quad (18)$$

which measures the departure of the normal field from integrability. The calculus of variations tells us that searching for a function Z that minimizes $\mathcal{F}(Z)$ is equivalent to the resolution of the associated Euler equation [6]:

$$\frac{\partial F}{\partial Z} - \frac{\partial}{\partial x} \left[ \frac{\partial F}{\partial Z_x} \right] - \frac{\partial}{\partial y} \left[ \frac{\partial F}{\partial Z_y} \right] = 0. \quad (19)$$

From (18) and (19), we deduce:

$$\nabla^2 Z = 1/g [p_x + q_y]. \quad (20)$$

This is a Poisson equation, which is not particularly difficult to solve, even analytically [11]. More precisely, solving Eq. (20) is really equivalent to the search for an extremum of $\mathcal{F}(Z)$ only if Z is constrained on the boundary $\partial\Omega$ of $\Omega$.

Otherwise, Eq. (20) has to be complemented with the “natural boundary equation” (see [6] and Section 3.2.2).

Under the assumption of perspective projection, the resolution of Eq. (13) in the least square sense amounts minimizing the following functional:

$$\mathcal{G}(T) = \iint_{(x,y) \in \Omega} G(x, y, T, T_x, T_y) \, dx \, dy, \quad (21)$$

where:

$$G(x, y, T, T_x, T_y) = [T_x - r]^2 + [T_y - s]^2. \quad (22)$$

The associated Euler equation is:

$$\nabla^2 T = r_x + s_y, \quad (23)$$

which is a Poisson equation again. This equation must also be complemented with the natural boundary equation if T is unconstrained on $\partial\Omega$.

## 3. A New Method of Normals Integration

In this section, which makes up the algorithmic contribution of the paper, we extend and combine two existing methods of normals integration. On the one hand, we show that the knowledge of the height on $\partial\Omega$ is not necessary. Since such a knowledge is usually not available, the need for a boundary knowledge could lead us to use erroneous data. On the other hand, we show how to extend both these methods to perspective projection. Finally, we combine them, in order to design a new method of normals integration which is fast and accurate at the same time.

### 3.1. Improvement of Wu and Li's Method

Propagation methods of integration consist in computing the height from starting points, where the height is known. Such points are usually located on $\partial\Omega$. These methods are fast but work well only if the normal field is integrable. Let us focus on Wu and Li's method [17].

#### 3.1.1 Wu and Li's Method

Wu and Li make a discrete approximation of the integral in Eq. (6). They compute the means of the discrete integrals that are computed along the two paths that are displayed in Fig. 2 and they show that the accuracy is optimal if $(x_0, y_0)$ and $(x, y)$ lie on the same diagonal. Suppose for seek of simplicity that the pixels are squared. The image points $(x_0, y_0)$ and $(x, y)$ can thus coincide with pixels $(i_0, j_0)$ and $(i, j)$. From now on, we will use the subscripted notation, e.g. $u_{i,j}$ for $u(i, j)$. Between two neighbouring pixels $(i, j)$ and $(i+1, j+1)$, a discrete approximation of (6) is:

$$Z_{i+1,j+1} \approx Z_{i,j} + \delta/(2g) [I_{1+2} + I_{3+4}], \quad (24)$$
<figure><img src="image_3.png" /><figcaption>Figure 2. The two paths used by Wu and Li: 1 followed by 2; 3 followed by 4.</figcaption></figure>

where $\delta$ denotes the pixel size ($\delta = 0.05$ in the tests), and:

$$ \begin{cases} I_{1+2} = p_{i,j} + p_{i+1,j} + q_{i+1,j} + q_{i+1,j+1}, \\ I_{3+4} = q_{i,j} + q_{i,j+1} + p_{i,j+1} + p_{i+1,j+1}. \end{cases} \quad (25) $$

Obviously, the knowledge of Z on $\partial\Omega$ is not required for the implementation of (24).

### 3.1.2 First Improvement

Since the shape can be computed only up to an additive constant, we arbitrarily fix the height at a starting point $Q_0$. Then, in order not to privilege a direction compared to another, we treat the pixels in spiral (cf. Fig. 3-a). The shape computed using the true normals of a synthetic vase is shown in Fig. 3-b: error $|\Delta Z|_2$ is equal to 0.03, but would be equal to 0.22 if we had imposed the true values of the height on $\partial\Omega$ and computed Z from the boundary to the center of the image.

<figure><img src="image_4.png" /><figcaption>Figure 3. First improvement of Wu and Li's method: (a) treatment of the pixels in spiral; (b) reconstructed shape of a vase from the true normals, without boundary condition.</figcaption></figure>

### 3.1.3 Second Improvement

Using the theoretical study of Section 2, we know that perspective integration is very similar to orthographic projection.

From this statement, we deduce the following approximation:

$$ T_{i+1,j+1} \approx T_{i,j} + \delta/2 [J_{1+2} + J_{3+4}], \quad (26) $$

with:

$$ \begin{cases} J_{1+2} = r_{i,j} + r_{i+1,j} + s_{i+1,j} + s_{i+1,j+1}, \\ J_{3+4} = s_{i,j} + s_{i,j+1} + r_{i,j+1} + r_{i+1,j+1}. \end{cases} \quad (27) $$

The main difference between (24) and (26) is that the focal length $f$ and the location of the principal point $O$ must be known to implement (26).

In Fig. 4, the graph of the true height $Z_t$ of the vase is shown i.e., the surface of equation $z = Z_t(x, y)$, which is not the surface of the vase. Contrary to the orthographic case, $Z_t$ is not uniform on $\partial\Omega$ (cf. Fig. 4-b). The graph of the height Z computed by the improved method of Wu and Li is shown in Fig. 5. It is globally satisfactory, even if $|\Delta Z|_2 = 0.11$, which is a little more than in the orthographic case. It is noticeable that the height on $\partial\Omega$ is qualitatively right (cf. Fig. 5-b), whereas no boundary condition is used.

<figure><img src="image_5.png" /><figcaption>Figure 4. Graph of the true height $Z_t$ of the vase: (a) perspective view; (b) side view.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 5. Graph of the height Z computed from the true normals by the improved method of Wu and Li: (a) perspective view; (b) side view.</figcaption></figure>

## 3.2. Improvement of Horn and Brooks' Method

The global methods of integration [6, 5] deeply differ from the propagation methods. Apart from their slowness,
they have two main advantages: they are more robust to noise; in the case where the Schwartz equation is not satisfied, they provide however an acceptable shape. Let us now focus on Horn and Brooks' method.

### 3.2.1 Horn and Brooks' Method

Horn and Brooks propose [6] a resolution of the Poisson equation (20) that comes from the following approximation of the expression (17) of $\mathcal{F}(Z)$:

$$ \mathcal{E}(\bar{Z}) = \sum_{(i,j) \in \Omega_1} \left[ \frac{Z_{i+1,j} - Z_{i,j}}{\delta} - \frac{p_{i+1,j} + p_{i,j}}{2g} \right]^2 + \sum_{(i,j) \in \Omega_2} \left[ \frac{Z_{i,j+1} - Z_{i,j}}{\delta} - \frac{q_{i,j+1} + q_{i,j}}{2g} \right]^2. \quad (28) $$

In this expression: $\Omega_1$ denotes the set of pixels $(i, j) \in \Omega$ such that $(i + 1, j) \in \Omega$; $\Omega_2$ denotes the set of pixels $(i, j) \in \Omega$ such that $(i, j + 1) \in \Omega$; $\bar{Z} = (Z_{i,j})_{(i,j) \in \bar{\Omega}}$; $\bar{\Omega}$ denotes the set of pixels $(i, j) \in \Omega$ whose four nearest neighbours are in $\Omega$. The values $Z_{i,j}$ of the pixels lying on $\partial\Omega$ are not considered as unknowns, since Horn and Brooks use a boundary condition of the Dirichlet type. For a pixel $(i, j) \in \bar{\Omega}$, we get from (28):

$$ \frac{\partial \mathcal{E}}{\partial Z_{i,j}} = \frac{2}{\delta} \left[ - \left( \frac{Z_{i+1,j} - Z_{i,j}}{\delta} - \frac{p_{i+1,j} + p_{i,j}}{2g} \right) - \left( \frac{Z_{i,j+1} - Z_{i,j}}{\delta} - \frac{q_{i,j+1} + q_{i,j}}{2g} \right) + \left( \frac{Z_{i,j} - Z_{i-1,j}}{\delta} - \frac{p_{i,j} + p_{i-1,j}}{2g} \right) + \left( \frac{Z_{i,j} - Z_{i,j-1}}{\delta} - \frac{q_{i,j} + q_{i,j-1}}{2g} \right) \right]. \quad (29) $$

The characterization $\nabla \mathcal{E} = 0$ of an extremum of $\mathcal{E}$ can thus be written [6], for $(i, j) \in \bar{\Omega}$:

$$ 4 Z_{i,j} - (Z_{i+1,j} + Z_{i,j+1} + Z_{i-1,j} + Z_{i,j-1}) + \frac{\delta}{2g} (p_{i+1,j} - p_{i-1,j} + q_{i,j+1} - q_{i,j-1}) = 0. \quad (30) $$

This equation is indeed a discrete approximation of the Poisson equation (20), if centered finite difference approximations of $p_x$ and $q_y$ are used and if the Laplacian operator of $Z$ is approximated by:

$$ \nabla^2 Z(i, j) \approx \frac{Z_{i+1,j} + Z_{i,j+1} + Z_{i-1,j} + Z_{i,j-1} - 4 Z_{i,j}}{\delta^2}. \quad (31) $$

### 3.2.2 First Improvement

In order to avoid the need for $Z$ on the boundary, it suffices to consider that all the values $Z_{i,j}$, for $(i, j) \in \Omega$, are unknowns. This implies that the equations $\partial \mathcal{E}/\partial Z_{i,j} = 0$, for

$(i, j) \in \Omega$, are not all written under the form (30). For example, if two neighbours $(i + 1, j)$ and $(i, j + 1)$ of $(i, j)$ only are in $\Omega$, then $\partial \mathcal{E}/\partial Z_{i,j} = 0$ is written:

$$ 2 Z_{i,j} - (Z_{i+1,j} + Z_{i,j+1}) = \frac{-\delta}{2g} (p_{i+1,j} + p_{i,j} + q_{i,j+1} + q_{i,j}). \quad (32) $$

The pixels which require a particular treatment all lie on $\partial\Omega$. In fact, equations such as (32) are nothing else than the discrete version of the “natural boundary equation” that has already been cited in Section 2. The results of Fig. 6 do not use any boundary condition.

<figure><img src="image_7.png" /><figcaption>Figure 6. Improvement of Horn and Brooks' method: reconstructed shapes from the true normals, without boundary condition, at iterations (a) 500 and (b) 8000.</figcaption></figure>

### 3.2.3 Second Improvement

In the perspective case, we saw in Section 2.3 that $\mathcal{F}(Z)$ must be replaced by $\mathcal{G}(T)$. The generalization of Horn and Brooks' method for perspective projection is thus straightforward and provide the results of Fig. 7.

<figure><img src="image_8.png" /><figcaption>Figure 7. Graphs of the heights $Z$ computed from the true normals by the improved method of Horn and Brooks, at iterations (a) 500 and (b) 8000.</figcaption></figure>

## 3.3. A New Method of Normals Integration

In order to take advantage of both them, we combine Wu and Li's and Horn and Brooks' methods. Wu and Li's method is fast and provides shapes of good quality, which consequently allows us to reduce the number of iterations
<figure><img src="image_9.png" /><figcaption>Figure 8. True shape of the vase: (a) perspective view; (b) side view.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 9. Reconstructed shape from the true normals by our new method of integration: (a) perspective view; (b) side view.</figcaption></figure>

of Horn and Brooks' method from 8000 to 500. Thus, the reconstruction is of good quality, even if the normal field is not integrable.

Contrary to orthographic integration, perspective integration does not allow us to directly compute the scene surface, since the graph $z = Z(x, y)$ of the computed height $Z$ is not the scene surface. But, knowing that an image point $Q = (x, y)$ is conjugated with the object point $P = (x Z(x, y)/f, y Z(x, y)/f, Z(x, y))$, it is easy to compute the scene surface from $Z$. In Fig. 8, the surface of the vase after resampling of the object points corresponding to the true height $Z_t$ is shown. Obviously, the bottom and the top of the vase are contained in two planes parallel to Cyz. The computed shape from the true normals of the vase using our new method of normals integration, after resampling, is shown in Fig. 9. This shape is visually similar to the true shape of the vase.

# 4. Application to Photometric Stereo

Among the computer vision techniques for 3D-reconstruction, the photometric methods use the relation between the greylevel and the normal to the surface. Thus, they are appropriate to test the new method of integration of a normal field described in Section 3.

## 4.1. Photometric Methods for 3D-reconstruction

### 4.1.1 Image Irradiance Equation

The image irradiance equation is the basic equation of the photometric methods for 3D-reconstruction, namely shape-from-shading and photometric stereo. This equation, which expresses the conservation of light energy, is particularly simple for a Lambertian surface illuminated by a homogeneous parallel beam [15]:

$$- \rho(P) \vec{L} \cdot \vec{n}(P) = I(Q). \quad (33)$$

In (33), an object point $P$ is characterized by its albedo $\rho(P)$ and its unit outgoing normal $\vec{n}(P)$, whereas its conjugated image point $Q$ is characterized by its greylevel $I(Q)$. Note that, for a colour image, $I(Q)$ as well as $\rho(P)$ must be considered as vectors in $\mathbb{R}^3$, and no longer as scalars. The light beam is characterized, in direction and density, by a vector $\vec{L}$. Note that both members of Eq. (33) are not homogeneous: in fact, rather than an equality, (33) is a relation of proportionality. Finally, it is noticeable that the resolution of (33) is possible only if the relation between $P$ and $Q$ is known i.e., if the model of projection is known.

### 4.1.2 Shape-from-shading and Photometric Stereo

The knowledge of $\rho(P)$, $\vec{L}$ and $I(Q)$ does not allow us to determine $\vec{n}(P)$ from Eq. (33). There is indeed an infinity of vectors which form the same angle with vector $\vec{L}$. This ambiguity, which is well-known, constitutes the main difficulty of shape-from-shading. We avoid shape-from-shading to illustrate our new method of integration, since it would be difficult to decide on its accuracy, knowing that a normal field computed through using shape-from-shading would not be accurate enough, especially in the case of real images [18, 4]. We prefer to make the tests in the framework of photometric stereo, a problem which is well-posed in most cases, contrary to shape-from-shading, since several images of the same scene, illuminated under various directions, are available. The photometric stereo technique, which is due to Woodham [15], has enjoyed some renewal in the last years [14, 9, 7, 13, 12, 16, 3, 2].

If several images of a scene are available, then in each image point $Q$, the equations (33) obtained for each image form a system, whose unknowns are the coordinates of the vector $\vec{m}(P) = \rho(P) \vec{n}(P)$. Three images at least, illuminated under directions that are non-coplanar, are thus necessary to determine $\vec{m}(P)$. The main difference between shape-from-shading and photometric stereo is that the albedo of the surface must be known, at least partially, so that shape-from-shading may become a sufficiently constrained problem, whereas this is not necessary for photometric stereo. If the albedo is known, then photometric
stereo can work with two images only, but if three images
are available, then it is not necessary to know the albedo.

## 4.2. Tests

### 4.2.1 Synthetic Images

We simulate five images of the same vase as that of Sec-
tion 3, under perspective projection (cf. Fig. 10). Note that
a coloured texture is warped on the vase. In each image
point, the albedo and the normal are computed by the res-
olution of a system of five equations and three unknowns.
The reconstructed shape using our new method of integra-
tion is displayed in Fig. 11. It is visually satisfactory, but
one may wonder whether this still holds with real images.

### 4.2.2 Real Images

We test our method of integration on three photographs of Beethoven's bust (cf. Fig. 12) that are available on the web<fnref n="1" />. The material of the bust seems to be approximately Lambertian. Moreover, estimates of the directions and densities of the light beams are given. The shape computed using our new method of integration is visually satisfactory (cf. Fig. 13).

# 5. Conclusion and Perspectives

In this paper, after a rigorous theoretical study, we show
how to improve two existing methods of normals integra-
tion. On the one hand, we show how to use them without
boundary condition, in order to make them more realistic.
On the other hand, we extend them to perspective projection
and show that this improves the accuracy of the reconstruc-
tions. Finally, we show that the joint use of both improved
methods provides very satisfactory results.

As an application, we use this new method of integration
in the framework of photometric stereo. The performances
of the method are illustrated on synthetic, as well as on real
images: they are very encouraging. Note that this applica-
tion of our work amounts taking into account perspective
in photometric stereo, a problem which has recently been
pioneered in [13].

As a perspective, we now aim at reconstructing scene
surfaces from images taken under badly-controlled condi-
tions i.e., for non-Lambertian surfaces illuminated by non-
homogeneous or non-parallel beams, a problem which has
already been addressed by Basri *et al.* [2].

## References

[1] A. Agrawal, R. Chellappa, and R. Raskar. An Algebraic Approach to Surface Reconstruction from Gradient Fields. In

Proc. 10th IEEE Int. Conf. Comp. Vis. (volume I), pages 174–181, Beijing, China, Oct. 2005. 1

[2] R. Basri, D. W. Jacobs, and I. Kemelmacher. Photometric Stereo with General, Unknown Lighting. *Int. J. Comp. Vis.*, 72(3):239–257, May 2007. 6, 7

[3] C.-P. Chen and C.-S. Chen. The 4-Source Photometric Stereo Under General Unknown Lighting. In *Proc. 9th Eur. Conf. Comp. Vis. (volume III)*, pages 72–83, Graz, Austria, May 2006. 6

[4] J.-D. Durou, M. Falcone, and M. Sagona. Numerical Methods for Shape-from-shading: A New Survey with Benchmarks. *Comp. Vis. and Im. Underst.*, 2007. (to appear). 6

[5] R. T. Frankot and R. Chellappa. A Method for Enforcing Integrability in Shape from Shading Algorithms. *IEEE Trans. Patt. Anal. Mach. Intell.*, 10(4):439–451, July 1988. 1, 4

[6] B. K. P. Horn and M. J. Brooks. The Variational Approach to Shape From Shading. *Comp. Vis., Grap., and Im. Proc.*, 33(2):174–208, Feb. 1986. 1, 3, 4, 5

[7] I. Horowitz and N. Kiryati. Depth from Gradient Fields and Control Points: Bias Correction in Photometric Stereo. *Im. and Vis. Comp.*, 22(9):681–694, Aug. 2004. 1, 6

[8] R. Klette and K. Schlüns. Height data from gradient fields. In *Proc. Machine Vision Applications, Architectures, and Systems Integration*, pages 204–215, Boston, Massachusetts, USA, Nov. 1996. 1

[9] L. Noakes and R. Kozera. Nonlinearities and Noise Reduction in 3-Source Photometric Stereo. *J. Math. Imag. and Vis.*, 18(2):119–127, Mar. 2003. 6

[10] A. Robles-Kelly and E. R. Hancock. A Graph-Spectral Method for Surface Height Recovery from Needle-maps. In *Proc. IEEE Conf. Comp. Vis. and Patt. Recog. (volume I)*, pages 141–148, Kauai, Hawaii, USA, Dec. 2001. 1

[11] T. Simchony, R. Chellappa, and M. Shao. Direct Analytical Methods for Solving Poisson Equations in Computer Vision Problems. *IEEE Trans. Patt. Anal. Mach. Intell.*, 12(5):435–446, May 1990. 3

[12] P. Tan, S. Lin, and L. Quan. Resolution-Enhanced Photometric Stereo. In *Proc. 9th Eur. Conf. Comp. Vis. (volume III)*, pages 58–71, Graz, Austria, May 2006. 6

[13] A. Tankus and N. Kiryati. Photometric Stereo under Perspective Projection. In *Proc. 10th IEEE Int. Conf. Comp. Vis. (volume I)*, pages 611–616, Beijing, China, Oct. 2005. 6, 7

[14] J. R. A. Torreão. Geometric-photometric approach to monocular shape estimation. *Im. and Vis. Comp.*, 21(12):1045–1061, Nov. 2003. 6

[15] R. J. Woodham. Photometric Method for Determining Surface Orientation from Multiple Images. *Opt. Engin.*, 19(1):139–144, 1980. 1, 6

[16] T.-P. Wu and C.-K. Tang. Dense Photometric Stereo by Expectation Maximization. In *Proc. 9th Eur. Conf. Comp. Vis. (volume IV)*, pages 159–172, Graz, Austria, May 2006. 6

[17] Z. Wu and L. Li. A Line-Integration Based Method for Depth Recovery from Surface Normals. *Comp. Vis., Grap., and Im. Proc.*, 43(1):53–66, July 1988. 1, 3

[18] R. Zhang, P.-S. Tsai, J. E. Cryer, and M. Shah. Shape from Shading: A Survey. *IEEE Trans. Patt. Anal. Mach. Intell.*, 21(8):690–706, Aug. 1999. 6

<footnote n="1">¹http://www.ece.ncsu.edu/imaging/Archives/Image DataBase/Industrial/</footnote>
<figure><img src="image_11.png" /><figcaption>Figure 10. Five images of the vase illuminated under various directions. Angles $\phi$ (azimuth) and $\theta$ (complementary of the elevation) are in degrees: (a) ($\phi_1$; $\theta_1$) = (0; 0); (b) ($\phi_2$; $\theta_2$) = (45; 30); (c) ($\phi_3$; $\theta_3$) = (135; 30); (d) ($\phi_4$; $\theta_4$) = (225; 30); (e) ($\phi_5$; $\theta_5$) = (315; 30).</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Figure 11. Reconstructed shape from the five images of Fig. 10: side view (left) and perspective view (right).</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Figure 12. Three photographs of Beethoven's bust available on the web: (a) ($\phi_6$; $\theta_6$) = (-72.48; 15.12); (b) ($\phi_7$; $\theta_7$) = (-4.85; 11.46); (c) ($\phi_8$; $\theta_8$) = (69.31; 16.11).</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Figure 13. Reconstructed shape from the three photographs of Fig. 12: side view (left) and perspective view (right).</figcaption></figure>