A signal-processing framework for inverse rendering.
Realism in computer-generated images requires accurate input models for lighting, textures and BRDFs. One of the best ways of obtaining high-quality data is through measurements of scene attributes from real photographs by inverse rendering. However, inverse rendering methods have been largely limited to settings with highly controlled lighting. One of the reasons for this is the lack of a coherent mathematical framework for inverse rendering under general illumination conditions. Our main contribution is the introduction of a signal-processing framework which describes the reflected light field as a convolution of the lighting and BRDF, and expresses it mathematically as a product of spherical harmonic coefficients of the BRDF and the lighting. Inverse rendering can then be viewed as deconvolution. We apply this theory to a variety of problems in inverse rendering, explaining a number of previous empirical results. We will show why certain problems are ill-posed or numerically ill-conditioned, and why other problems are more amenable to solution. The theory developed here also leads to new practical representations and algorithms. For instance, we present a method to factor the lighting and BRDF from a small number of views, i.e. to estimate both simultaneously when neither is known.
Introduction
To create a realistic computer-generated image, we need both
an accurate, physically-based rendering algorithm and a detailed
model of the scene including light sources and objects specified
by their geometry and material properties-texture and reflectance
(BRDF). There has been substantial progress in the development
of rendering algorithms, and nowadays, realism is often limited by
the quality of input models. As a result, image-based rendering is
becoming widespread. In its simplest form, image-based rendering
uses view interpolation to construct new images from acquired
images without constructing a conventional scene model.
The quality of view interpolation may be significantly improved
if it is coupled with inverse rendering. Inverse rendering measures
rendering attributes-lighting, textures, and BRDF-from
photographs. Whether traditional or image-based rendering algorithms
are used, rendered images use measurements from real ob-
jects, and therefore appear very similar to real scenes. Measuring
scene attributes also introduces structure into the raw imagery,
making it easier to manipulate the scene. For example, an artist can
change independently the material properties or the lighting.
Inverse rendering methods such as those of Debevec et al. [6],
Marschner et al. [21], and Sato et al. [32], have produced high
Real Photograph Rendered Image

Figure

1: Left: Real Photograph Right: Rendered image. The BRDF used for
the rendered image was estimated under complex unknown illumination from 3 photographs
of a cat sculpture with known geometry. Our algorithm also recovered the
lighting distribution, which consisted of two directional sources and an area source.
The images above show a new view not used in BRDF recovery; the lighting is also
new, being composed of a single directional source (with known direction) not used in
BRDF estimation. These images show that the recovered BRDF accurately predicts
appearance even under novel viewing and lighting conditions.
quality measurements. However, most previous work has been conducted
in highly controlled lighting conditions, usually by careful
active positioning of a single point source. Even methods that work
in outdoor conditions, such as those of Yu and Malik [39], Sato and
Ikeuchi [31] and Love [17], are designed specifically for natural
illumination, and assume a simple parametric model for skylight.
Previous methods have also usually tried to recover only one of the
unknowns-texture, BRDF or lighting. The usefulness of inverse
rendering would be greatly enhanced if it could be applied under
general uncontrolled lighting, and if we could simultaneously estimate
more than one unknown. For instance, if we could recover
both the lighting and BRDF, we could determine BRDFs under
unknown illumination. One reason there has been relatively little
work in these areas is the lack of a common theoretical framework
for determining under what conditions inverse problems can and
cannot be solved, and for making principled approximations.
Our main contribution is a theoretical framework for analyzing
the reflected light field from a curved convex homogeneous surface
under distant illumination. We believe this framework provides a
solid mathematical foundation for many areas of graphics. With
respect to inverse rendering, we obtain the following results:
Reflection as Convolution: It has been observed qualitatively by
Miller and Hoffman [23], Cabral et al. [3], Bastos et al. [2] and
others that the reflection operator behaves like a convolution in the
angular domain. We formalize these notions mathematically. The
reflected light field can therefore be thought of in a precise quantitative
way as obtained by convolving the lighting and BRDF, i.e.
by filtering the illumination using the BRDF. We believe this is a
useful way of analyzing many computer graphics problems. In par-
ticular, inverse rendering can be viewed as deconvolution.
Well-posedness and Conditioning of Inverse Problems: Inverse
problems can be ill-posed-there may be no solutions or several
solutions. They are also often numerically ill-conditioned, i.e.
extremely sensitive to noisy input data. From our theory, we are
able to analyze the well-posedness and conditioning of a number
of inverse problems, explaining many previous empirical observa-
tions. This analysis can serve as a guideline for future research.
New Practical Representations and Algorithms: Insights from
the theory lead to the derivation of a simple practical representation,
which can be used to estimate BRDFs under complex lighting. The
theory also leads to novel frequency space and hybrid angular and
frequency space methods for inverse problems, including two new
algorithms for estimating the lighting, and an algorithm for simultaneously
determining the lighting and BRDF. Therefore, we can
recover BRDFs under general, unknown lighting conditions.
Previous Work
To describe previous work, we will introduce a taxonomy based on
how many of the three quantities-lighting, BRDF and texture-
are unknown. To motivate the taxonomy, we first write a simplified
version of the reflection equation, omitting visibility.
Z#
Here, B is the reflected light field, expressed as a function of the
surface position x and outgoing direction # . The normal vector is
#n. For simplicity, we assume that a single texture T modulates the
BRDF. In practice, we would use separate textures for the diffuse
and specular components of the BRDF.
The integrand is a product of terms-the texture T (x), the
and the lighting L(x, # i ). Inverse rendering, assuming
known geometry, involves inverting the integral to recover
one or more of #, L, or T . If two or more quantities are unknown,
inverse rendering involves factoring the reflected light field.
One Unknown
1. Unknown Texture: Previous methods have recovered the
diffuse texture of a surface using a single point light source by
dividing by the irradiance in order to estimate the albedo at each
point. Details are given by Marschner [34] and Levoy et al. [16].
2. Unknown BRDF: The BRDF [24] is a fundamental intrinsic
surface property. Active measurement methods, known as goniore-
flectometry, involving a single point source and a single observation
at a time, have been developed. Improvements are suggested
by Ward [37] and Karner et al. [12]. More recently, image-based
BRDF measurement methods have been proposed by Lu et al. [18]
and Marschner et al. [21]. If the entire BRDF is measured, it may be
represented by tabulating its values. An alternative representation
is by low-parameter models such as those of Ward [37] or Torrance
and Sparrow [36]. The parametric BRDF will generally not be as
accurate as a full measured BRDF. However, parametric models are
often preferred in practice since they are compact, and are simpler
to estimate. Love [17] estimates parametric BRDFs under natural
assuming a low-parameter model for skylight and
sunlight. Dror et al. [7] classify the surface reflectance as one of
a small number of predetermined BRDFs, making use of assumed
statistical characteristics of natural lighting. However, the inverse
BRDF problem has not been solved for general illumination.
3. Unknown Lighting: A common solution is to use a mirrored
ball, as done by Miller and Hoffman [23]. Marschner and Greenberg
[20] find the lighting from a Lambertian surface. D'Zmura [8]
proposes, but does not demonstrate, estimating spherical harmonic
coefficients. For Lambertian objects, we [29] have shown how to
recover the first 9 spherical harmonics. Previous work has not estimated
the lighting from curved surfaces with general BRDFs.
Two Unknowns
4. Factorization-Unknown Lighting and BRDF: BRDF
estimation methods have been proposed by Ikeuchi and Sato [10]
and Tominaga and Tanaka [35] for the special case when the lighting
consists of a single source of unknown direction. However,
these methods cannot simultaneously recover a complex lighting
distribution and the object BRDF. One of the main practical contributions
of this paper is a solution to this problem for curved sur-
allowing us to estimate BRDFs under general unknown illu-
mination, while also determining the lighting. The closest previous
work is that of Sato et al. [30] who use shadows to estimate the
illumination distribution and the surface reflectance properties. We
extend this work by not requiring shadow information, and presenting
improved methods for estimating the illumination.
5. Factorization-Unknown Texture and BRDF: This corresponds
to recovering textured, or spatially-varying BRDFs. Sato
et al. [32] rotate an object on a turntable, using a single point
source, to recover BRDF parameters and texture. Yu et al. [38] recover
a texture only for the diffuse BRDF component, but account
for interreflections. Using a large number of images obtained by
moving a point source around a sphere surrounding the subject,
Debevec et al. [6] acquire the reflectance field of a human face,
and recover parameters of a microfacet BRDF model for each surface
location. Sato and Ikeuchi [31] and Yu and Malik [39] recover
BRDFs and diffuse textures under natural illumination, assuming
a simple parametric model for skylight, and using a sequence of
images acquired under different illumination conditions. Most of
these methods recover only diffuse textures; constant values, or
relatively low-resolution textures, are used for the specular param-
eters. A notable exception is the work of Dana et al. [5] who generalize
BRDFs to a 6D bi-directional texture function (BTF).
6. Factorization-Unknown Lighting and Texture: We
have shown [29] that a distant illumination field can cause only low
frequency variation in the radiosity of a convex Lambertian sur-
This implies that, for a diffuse object, high-frequency texture
can be recovered independently of lighting. These observations are
in agreement with the perception literature, such as Land's retinex
theory [15], wherein high-frequency variation is usually attributed
to texture, and low-frequency variation associated with
tion. However, note that there is a fundamental ambiguity between
low-frequency texture and lighting effects. Therefore, lighting and
texture cannot be factored without using active methods or making
further assumptions regarding their expected characteristics.
General Case: Three Unknowns
7. Factorization-Unknown Lighting, Texture, BRDF:
Ultimately, we wish to recover textured BRDFs under unknown
lighting. We cannot solve this problem without further assump-
tions, because we must first resolve the lighting-texture ambiguity.
Our approach differs from previous work in that it is derived
from a mathematical theory of inverse rendering. As such, it has
similarities to inverse methods used in areas of radiative transfer
and transport theory such as hydrologic optics [26] and neutron
scattering. See McCormick [22] for a review.
In previous theoretical work, D'Zmura [8] has analyzed reflection
as a linear operator in terms of spherical harmonics, and discussed
some resulting perceptual ambiguities between reflectance
and illumination. In computer graphics, Cabral et al. [3] first
demonstrated the use of spherical harmonics to represent BRDFs.
We extend these methods by explicitly deriving the frequency-space
reflection equation (i.e. convolution formula), and by providing
quantitative results for various special cases. We have earlier reported
on theoretical results for planar or flatland light fields [27],
and for determining the lighting from a Lambertian surface [29].
For the Lambertian case, similar results have been derived independently
by Basri and Jacobs [1] in simultaneous work on face
recognition. This paper extends these previous results to the general
3D case with arbitrary isotropic BRDFs, and applies the theory
to developing new practical inverse-rendering algorithms.
Assumptions
The input to our algorithms consists of object geometry and photographs
from a number of different locations, with known extrinsic
and intrinsic camera parameters. We assume static scenes, i.e. that
the object remains stationary and the lighting remains the same between
views. Our method is a passive-vision approach; we do not
actively disturb the environment. Our assumptions are:
Reflected radiance
Coefficients of basis-function expansion of B
Incoming radiance
lm Coefficients of spherical-harmonic expansion of L
Surface BRDF
# BRDF multiplied by cosine of incident angle
# lpq Coefficients of spherical-harmonic expansion of -
Incident elevation angle in local, global coordinates
Incident azimuthal angle in local, global coordinates
Outgoing elevation angle in local, global coordinates
Outgoing azimuthal angle in local, global coordinates
Hemisphere of integration in local,global coordinates
x Surface position
Surface normal parameterization-elevation angle
Surface normal parameterization-azimuthal angle
R # Rotation operator for surface normal (#)
D l
Matrix related to Rotation Group SO(3)
Y lm Spherical Harmonic basis function
lm
Complex Conjugate of Spherical Harmonic
# l Normalization constant,
I # -1

Figure

2: Notation
Known Geometry: We use a laser range scanner and a volumetric
merging algorithm [4] to obtain object geometry. By assuming
known geometry, we can focus on lighting and material properties.
Curved Objects: Our theoretical analysis requires curved sur-
faces, and assumes knowledge of the entire 4D reflected light field,
corresponding to the hemisphere of outgoing directions for all surface
orientations. However, our practical algorithms will require
only a small number of photographs.
Distant Illumination: The illumination field will be assumed to
be homogeneous, i.e. generated by distant sources, allowing us to
use the same lighting function regardless of surface location. We
treat the lighting as a general function of the incident angle.
Isotropic BRDFs: We will consider only surfaces having
isotropic BRDFs. The BRDF will therefore be a function of only 3
variables, instead of 4, i.e. 3D instead of 4D.
interreflection will be
ignored. Also, shadowing is not considered in our theoretical anal-
ysis, which is limited to convex surfaces. However, we will account
for shadowing in our practical algorithms, where necessary.
4 Theory of Reflection as Convolution
This section presents a signal-processing framework wherein reflection
can be viewed as convolution, and inverse rendering as de-
convolution. First, we introduce some preliminaries, defining the
notation and deriving a version of the reflection equation. We then
expand the lighting, BRDF and reflected light field in spherical harmonics
to derive a simple equation in terms of spherical harmonic
coefficients. The next section explores implications of this result.
Incoming Light (L) Outgoing Light (B)
a a

Figure

3: Schematic of reflection. On top, we show the situation with respect to
the local surface. The BRDF maps the incoming light distribution L to an outgoing
light distribution B. The bottom figure shows how the rotation # affects the situation.
Different orientations of the surface correspond to rotations of the upper hemisphere
and BRDF, with global directions corresponding to local directions (# i
4.1 Preliminaries
For the purposes of theoretical analysis, we assume curved convex
isotropic surfaces. We also assume homogeneous objects, i.e.
untextured surfaces, with the same BRDF everywhere. We parameterize
the surface by the spherical coordinates of the normal vector
(#), using the standard convention that corresponds
to the north pole or +Z axis. Notation used in this section
is listed in figure 2, and a diagram is in figure 3. We will use two
types of coordinates. Unprimed global coordinates denote angles
with respect to a global reference frame. On the other hand, primed
local coordinates denote angles with respect to the local reference
frame, defined by the local surface normal. These two coordinate
systems are related simply by a rotation, to be defined shortly.
Reflection Equation: We modify equation 1 based on our as-
sumptions, dropping the texturing term, and using the surface normal
(#) instead of the position x to parameterize B. Since
assumed to be independent of x, we write it as
Finally, ( #
can be written simply as cos # i , the
cosine of the incident angle in local coordinates.
Z#
We have mixed local (primed) and global (unprimed) coordinates.
The lighting is a global function, and is naturally expressed in a
global coordinate frame as a function of global angles. On the
other hand, the BRDF is naturally expressed as a function of the
local incident and reflected angles. When expressed in the local
coordinate frame, the BRDF is the same everywhere for a homogeneous
Similarly, when expressed in the global coordinate
frame, the lighting is the same everywhere, under the assumption
of distant illumination. The reflected radiance B can be expressed
conveniently in either local or global coordinates; we have used local
coordinates to match the BRDF. Similarly, integration can be
conveniently done over either local or global coordinates, but the
upper hemisphere is easier to express in local coordinates.
We now define a transfer 1 function -
in order to absorb
the cosine term. With this modification, equation 2 becomes
Z#
Rotations-Converting Local and Global coordinates:
Local and global coordinates are related by a rotation corresponding
to the surface normal (#). The north pole in local coor-
dinates, (0 # , is the surface normal. The corresponding global
coordinates are clearly (#). We define R# as a rotation operator
2 on column vectors that rotates
and is given by Rz is a rotation about
the Z axis and Ry a rotation about the Y axis.
We can now write the dependence on incident angle in equation 3
entirely in global coordinates, or entirely in local coordinates.
R#
R#
1 If we want the transfer function to be reciprocal, i.e. symmetric with respect to
incident and outgoing angles, we may multiply both the transfer function and the reflected
light field by cos #
2 For anisotropic surfaces, we need an initial rotation about Z to set the local tangent
frame. We would then have rotations about Z, Y and Z-the familiar Euler-
Angle parameterization. Since we are dealing with isotropic surfaces, we have ignored
this initial Z rotation, which has no physical significance. It is not difficult to derive
the theory for the more general anisotropic case.
Interpretation as Convolution: In the spatial domain, convolution
is the result generated when a filter is translated over an
input signal. However, we can generalize the notion of convolution
to other transformations Ta , where Ta is a function of a, and write
Z
When Ta is a translation by a, we obtain the standard expression
for spatial convolution. When Ta is a rotation by the angle a, the
above formula defines convolution in the angular domain.
Therefore, equations 4 and 5 represent rotational convolutions.
Equation 4 in global coordinates states that the reflected light field
at a given surface orientation corresponds to rotating the BRDF to
that orientation, and then integrating over the upper hemisphere.
The BRDF can be thought of as the filter, while the lighting is the
input signal. Symmetrically, equation 5 in local coordinates states
that the reflected light field at a given surface orientation may be
computed by rotating the lighting into the local coordinate system
of the BRDF, and then doing the hemispherical integration.
4.2 Spherical Harmonic Representation
For the translational case, the well-known frequency-space convolution
formula is given in terms of Fourier transforms. For a general
operator, an analogous formula can be obtained in terms of group
representations and the associated basis functions. For translations,
these basis functions are sines and cosines-the familiar Fourier
basis. For rotations, the corresponding basis functions are spherical
harmonics, and we now proceed to derive the frequency-space
rotational convolution formula in terms of spherical harmonics.
Inui et al. [11] is a good reference for background on spherical
harmonics and their relationship to rotations. Our use of spherical
harmonics to represent the lighting is similar in some respects
to previous methods [25] that use steerable linear basis functions.
Spherical harmonics, as well as the closely related Zernike Polyno-
mials, have been used before to represent BRDFs [3, 14, 33].
Spherical harmonics are the analog on the sphere to the Fourier
basis on the line or circle. The spherical harmonic Ylm is given by
s
(l -m)!
(l +m)!
l (cos #)e Im#
where Nlm is a normalization factor. In the above equation, the
azimuthal dependence is expanded in terms of Fourier basis func-
tions. The # dependence is expanded in terms of the associated
Legendre functions P m
l . The indices obey l # 0 and -l # m # l.
The rotation formula for spherical harmonics is
l
D l
mm #)e Im#
The important thing to note here is that the m indices are mixed-a
spherical harmonic after rotation must be expressed as a combination
of other spherical harmonics with different m indices. How-
ever, the l indices are not mixed; rotations of spherical harmonics
with order l are composed entirely of other spherical harmonics
with order l. For given order l, D l is a matrix that tells us how a
spherical harmonic transforms under rotation about the y-axis, i.e.
how to rewrite a rotated spherical harmonic as a linear combination
of all the spherical harmonics of the same order.
We begin by expanding the lighting in global coordinates.
l
m=-l
Here, the coefficients Llm can be computed in the standard way by
integrating against the complex conjugate Y #
lm
Z 2#
We now represent the transfer function -
of spherical harmonics. Note that -
# is nonzero only over the upper
hemisphere, i.e. when cos #
We are interested in isotropic BRDFs, which depend only on
This implies that the BRDF is invariant with respect
to adding a constant angle # to both incident and outgoing azimuthal
angles. It can be shown from the form of the spherical harmonics
that this condition forces all terms to vanish unless
The use of the complex conjugate for Y #
lm in the expansion above
is to make We now write
Furthermore, invariance of the BRDF with respect to negating
both incident and outgoing azimuthal angles requires that -
Finally, we use only three indices for the BRDF.
To represent the reflected light field, we define a new set of orthonormal
basis functions. The normalization and form of these
functions are derived in the appendix. In particular, the matrix
D comes from the rotation formula for spherical harmonics, equation
6. It will be convenient to first define a normalization constant.
r
r
The new basis functions can then be written
The expansion of the reflected light field is now
The translational convolution theorem expresses convolution in
frequency-space as a product of Fourier coefficients. For the rotational
case, an analogous result is derived in the appendix, using
spherical harmonics instead of complex exponentials. The
frequency-space reflection equation (or rotational convolution for-
mula) is a similar product of basis-function coefficients.
Implications
This section explores the implications of our results for problems
in inverse rendering, and works out some special cases in detail.
Our theory indicates which inverse problems are tractable, as opposed
to being ill-posed or ill-conditioned. Finally, we will use the
insights gained to develop a new practical representation.
5.1 General Observations
Inverse can be manipulated to yield
Llm
We may use any index m in inverse BRDF computation. Therefore,
BRDF recovery is well-posed unless the denominator vanishes for
all m, i.e. all terms for some order l in the spherical harmonic expansion
of the lighting vanish. In signal processing terms, if the
input signal (lighting) has no amplitude along certain modes of the
filter (BRDF), those modes cannot be estimated. BRDF recovery is
well conditioned when the lighting contains high frequencies like
directional sources, and is ill-conditioned for soft lighting.
Inverse Lighting: Equation 10 can also be manipulated to yield
Similarly as for BRDF recovery, any p, q can be used for inverse
lighting. The problem is well-posed unless the denominator -
vanishes for all p, q for some l. In signal processing terms, when
the BRDF filter truncates certain frequencies in the input lighting
signal (for instance, if it were a low-pass filter), we cannot determine
those frequencies from the output signal. Inverse lighting is
well-conditioned when the BRDF has high-frequency components
like sharp specularities, and is ill-conditioned for diffuse surfaces.
Light Field Factorization-Lighting and BRDF: We now
consider the problem of factorizing the light field, i.e simultaneously
recovering the lighting and BRDF when both are unknown.
The reflected light field is defined on a four-dimensional domain
while the lighting is a function of two dimensions and the isotropic
BRDF is defined on a three-dimensional domain. This seems to indicate
that we have more knowns (in terms of coefficients of the reflected
light field) than unknowns (lighting and BRDF coefficients).
For fixed order l, we can use known lighting coefficients Llm
to find unknown BRDF coefficients -
# lpq and vice-versa. In fact,
we need only one known nonzero lighting or BRDF coefficient to
bootstrap this process. It would appear from equation 10, however,
that there is an unrecoverable scale factor for each order l, corresponding
to the known coefficient we require. But, we can also use
reciprocity of the BRDF. To make the transfer function symmetric,
we multiply it, as well as the reflected light field B, by cos #
The new transfer function -
# is symmetric with respect to incident
and outgoing directions, and corresponding indices: -
# plq .
There is a global scale factor we cannot recover, since -
B is not
affected if we multiply the lighting and divide the BRDF by the
same amount. Therefore, we scale the lighting so the DC term
1/ (4#). Now, using equations 11, 12, and 13,
Llm
In the last line, we can use any value of m. This gives an explicit
formula for the lighting and BRDF in terms of coefficients of the
output light field. Therefore, up to global scale, the reflected light
field can be factored into the lighting and the BRDF, provided
the appropriate coefficients of the reflected light field do not vanish.
5.2 Special Cases
Mirror BRDF: The mirror BRDF corresponds to a gazing
sphere. Just as the inverse lighting problem is easily solved in angular
space in this case, we will show that it is well-posed and easily
solved in frequency space. The BRDF involves a delta function,
Note that the BRDF is nonzero only when # i #/2 and #
The coefficients for the BRDF, reflected light field, and lighting are
-0.4
Clamped
l ->
coefficient
->

Figure

4: Left: Successive approximations to the clamped cosine function by adding
more spherical harmonic terms. For we get a very good approximation. Right:
The solid line is a plot of spherical harmonic coefficients A l = # l - # l For l > 1,
odd terms vanish, and even terms decay rapidly.
The factor of (-1) q is because the azimuthal angle changes by #
upon reflection. We see that the lighting coefficients correspond in
a very direct way to the coefficients of the reflected light field. In
signal processing terminology, the inverse lighting problem is well
conditioned because the frequency spectrum of a delta function remains
constant with increasing order l, and does not decay.
Single Directional Source: For convenience, we position the
coordinate axes so that the source is located at +Z, i.e. at (0, 0).
Because the directional source is described by a delta function,
the spherical harmonic expansion coefficients are given simply by
lm (0), which vanishes for m #= 0. Thus,
In angular space, a single observation corresponds to a single
BRDF measurement. This property is used in image-based BRDF
measurement [18, 21]. We see that in frequency space, there is a
similar straightforward relation between BRDF coefficients and reflected
light field coefficients. BRDF recovery is well-conditioned
since we are estimating the BRDF filter from its impulse response.
Lambertian BRDF: For a Lambertian object, the transfer function
is a scaled clamped cosine function, since it is proportional to
the cosine of the incident angle over the upper hemisphere when
and is equal to 0 over the lower hemisphere. Plots
of spherical-harmonic fits to the clamped cosine function and the
magnitude of the coefficients are shown in figure 4. Because there
is no dependence on outgoing angle, we can drop the indices p and
q. Further, the reflected light field is now effectively the surface
radiosity function, and can be expanded 3 in spherical harmonics.
l
m=-l
BlmYlm (#)
We [29] have shown that with the definitions,
\Theta

one can derive
Blm
l
We define -
l . An analytic formula for A l may be derived
[29]. It can be shown that -
A l vanishes for odd values of l > 1,
3 The basis functions Clmpq in equation 9 become
m0 (#)e Im# if we
ignore output dependence, and set (the BRDF is azimuthally symmetric). It
can be shown that this is simply Ylm (#). Equation 15 now follows naturally
from equation 10 upon dropping indices p and q. Our previous derivation [29] was
specialized to the Lambertian case, and ignored the output dependence from the onset.
and even terms fall off very rapidly as l -5/2 . More than 99% of the
energy of the BRDF filter is captured by l # 2. Numerically,
Thus, the Lambertian BRDF acts like a low-pass filter, truncating
or severely attenuating frequencies with l > 2. Therefore,
from observations of a Lambertian surface, estimation of the illumination
is formally ill-posed, and is well-conditioned only for the
lighting coefficients with l # 2, corresponding to 9 parameters-1
and 5 for order 2 (l This explains the ill-conditioning
observed by Marschner and Greenberg [20] in trying
to solve the inverse lighting problem from a Lambertian surface.
Furthermore, for practical applications, including forward rendering
[28], the reflected light field from a Lambertian surface can
be characterized using only its first 9 spherical harmonic coef-
ficients; lighting effects cannot produce high-frequency variation
in intensity with respect to surface curvature.
Phong BRDF: The normalized Phong transfer function is
R is the reflection of the outgoing (viewing) direction about
the surface normal, #
L is the direction to the light source, and s is
the shininess, or Phong exponent. The normalization ensures the
Phong lobe has unit energy. Technically, we must also zero the
BRDF when the light is not in the upper hemisphere. However, the
Phong BRDF is not physically based, so others have often ignored
this boundary effect, and we will do the same.
We now reparameterize by the reflection vector #
R, transforming
the integral over the upper hemisphere centered on the surface normal
to an integral centered on #
R. The reflection vector takes the
place of the normal in the analysis, with (#) referring to #
R, and
. The Phong BRDF after reparameterization is mathematically
analogous to the Lambertian BRDF just discussed. In
fact, the properties of convolution can be used to show that for the
Phong BRDF, blurring the lighting and using a mirror BRDF is
equivalent to using the real lighting and real BRDF. This formalizes
the transformation often made in rendering with environment
maps [23]. Specifically, equation 15 can be written as
# l Llm
Here, L # lm is the blurred illumination and -
# l is the mirror BRDF 4 .
The BRDF coefficients depend on s, and are given by
Z #/2\Theta cos # i
This integral may be solved analytically. Formulae are in the ap-
pendix, and numerical plots are in figure 5.
-0.20.20.61l ->
coefficient
->

Figure

5: Numerical plots of the Phong coefficients # l -
l , as defined by equation 18.
The solid lines are the approximations in equation 19.
4 The formula # l -
not identical to equation 14 since we have now reparameterized
by the reflection vector. This accounts for the slightly different normalization.
For large s and l # s, a good approximation is
l 2
2s
The coefficients fall off as a gaussian with width of order # s. The
Phong BRDF behaves in the frequency domain like a gaussian fil-
ter, with the filter width controlled by the shininess. Therefore,
inverse lighting calculations will be well-conditioned only up to
order # s. As s approaches infinity, # l -
and the frequency
spectrum becomes constant, corresponding to a perfect mirror.
Microfacet BRDF: We now consider a simplified 4-parameter
Torrance-Sparrow [36] model, with parameters Kd , Ks , - and #.
This microfacet model is widely used in computer graphics.
FS
F (-,
exp
The subscript h stands for the half-way vector. F (-, # is the
Fresnel term for refractive index -; we normalize it to be 1 at normal
exitance. Actually, F depends on the angle with respect to the
half-way vector; in practice, this angle is usually very close to #
For simplicity in the analysis, we have omitted the geometric attenuation
factor G. In practice, this omission is not very significant
except for observations made at grazing angles, which are usually
assigned low confidence anyway in practical applications.
We focus on the specular component, reparameterizing by the
reflection vector, as for the Phong BRDF. It will also simplify matters
to fix the exitant direction, and focus on the frequency-space
representation of the incident-angle dependence. Precise analytic
are difficult to derive, but we can make a good approxi-
mation, as shown in the appendix. For normal exitance,
\Theta
For normal exitance, the specular part of the BRDF is a gaussian,
so equation 20 simply states that even in the spherical-harmonic
basis, the frequency spectrum of a gaussian is also gaussian, with
the frequency width related to the reciprocal of the angular width.
For non-normal exitance, microfacet BRDFs are not symmetric
about the reflection vector. Unlike for the Phong BRDF, there is a
preferred direction, determined by the exitant angle. However, the
BRDF filter is essentially symmetric about the reflected direction
for small viewing angles, as well as for low frequencies l. Hence, it
can be shown by Taylor-series expansions and verified numerically,
that the corrections to equation 20 are small under these conditions.
Finally, we approximate the effects of the Fresnel factor at non-normal
exitance by multiplying our expressions by F (-, #
With respect to the conditioning of inverse problems, equation
20 indicates that inverse lighting from a microfacet BRDF is
well-conditioned only for frequencies up to order l # -1 . Equation
also indicates that BRDF estimation is ill-conditioned under
low-frequency lighting. For low-frequency lighting, we may apply
the properties of convolution as we did for Phong BRDFs, filtering
the lighting using equations 17 and 20, while using a mirror BRDF.
Note that for frequencies l << # -1 , the effects of this filtering
are insignificant. The BRDF passes through virtually all the low-frequency
energy. Therefore, if the lighting contains only low
frequencies, the reflected light field from a microfacet BRDF
is essentially independent of the BRDF filter width # -1 ; this
makes estimation of the surface roughness # ill-conditioned.
5.3 Practical Representation
Thus far, we have presented the theoretical foundation for, and
some implications of, a frequency-space view of reflection. A signal
processing approach has been used before in some other areas
of computer graphics, notably the theory of aliasing. Just as a
frequency-space analysis of aliasing gives many insights difficult
to obtain by other means, the last two sections lead to new ways of
analyzing inverse rendering problems. However, the Fourier-space
theory of aliasing is not generally used directly for antialiasing. The
ideal Fourier-space bandpass filter in the spatial domain, the sinc
function, is usually modified for practical purposes because it has
infinite extent and leads to ringing. Similarly, representing BRDFs
purely as a linear combination of spherical harmonics leads to ring-
ing. Moreover, it is difficult to compute Fourier spectra from sparse
irregularly sampled data. Similarly, it is difficult to compute the reflected
light field coefficients Blmpq from a few photographs; we
would require a very large number of input images, densely sampling
the entire sphere of possible directions.
For these reasons, the frequency-space ideas must be put into
practice carefully. Here, we first discuss two useful practical
techniques-dual angular and frequency-space representations,
and the separation of the lighting into slow and fast-varying com-
ponents. Finally, we use these ideas, and the insights gained from
the previous subsection, to derive a simple practical model of the
reflected light field for the microfacet BRDF. This representation
will be used extensively in the practical algorithms of section 6.
Dual Angular and Frequency-Space Representations:
Quantities local in angular space have broad frequency spectra and
vice-versa. By developing a frequency-space view of reflection,
we ensure that we can use either the angular-space or frequency-space
representation, or even a combination of the two. The diffuse
BRDF component is slowly varying in angular-space, but is local
in frequency-space, while the specular BRDF component is local in
the angular domain. For representing the lighting, the frequency-space
view is appropriate for the diffuse BRDF component, while
the angular-space view is appropriate for the specular component.
Separation of slow and fast-varying lighting: For the
angular-space description of the lighting, used in computing the reflected
light field from the specular BRDF component, we separate
the lighting into a slow varying component corresponding to low
frequencies or area sources-for which we filter the lighting and
use a mirror BRDF-and a fast varying component corresponding
to high frequencies or directional sources. For the frequency-space
lighting description, used for the diffuse BRDF component, this
distinction need not be made since the formulae for the Lambertian
BRDF are the same for both slow and fast varying components.
Model for Reflected Light Field: Our model for the reflected
light field from the microfacet BRDF includes three components.
s,fast
Bd is from the diffuse component of the BRDF. B s,slow represents
specularities from the slowly-varying lighting, and B s,fast specular
highlights from the fast varying lighting component.
To write Bd , corresponding to the Lambertian BRDF compo-
nent, we use the 9 parameter frequency-space representation of the
lighting. Explicitly noting l # 2, and with E being the irradiance,
E(#) =X
+l
m=-l
LlmYlm (#)
The numerical values of # l - # l are given in equation 16.
For B s,slow , we filter the lighting, using equations 17 and 20,
and treat the BRDF as a mirror. With #
R denoting the reflected
direction, and Lslow the filtered version of the lighting, we obtain
For the fast varying portion of the lighting-corresponding to
sources of angular width #-we treat the total energy of the
source, given by an integral over the (small) solid angle subtended,
as located at its center, so the lighting is a s,fast
is given by the standard equation for the specular highlight from a
directional source. The extra factor of 4 cos # o in the denominator
as compared to equation 22 comes from the relation between
differential microfacet and global solid angles.
j,fast
The subscript j denotes a particular directional source; there could
be several. Note that L j,fast is now the total energy of the source.
For BRDF estimation, it is convenient to expand out these equa-
tions, making dependence on the BRDF parameters explicit.
6 Algorithms and Results
This section shows how the theory, and in particular the model just
derived in section 5.3, can be applied to a broad range of practical
inverse rendering problems. We present two types of methods-
algorithms that recover coefficients of a purely frequency-space description
of the lighting or BRDF by representing these quantities
as a sum of spherical harmonic terms, and algorithms that estimate
parameters corresponding to our model of section 5.3. Section 6.1
on BRDF estimation demonstrates direct recovery of spherical harmonic
BRDF coefficients, as well as estimation of parametric microfacet
BRDFs using equation 24. Similarly, section 6.2 demonstrates
direct recovery of spherical harmonic lighting coefficients,
as well as estimation of a dual angular and frequency-space lighting
description as per the model of section 5.3. Section 6.3 shows how
to combine BRDF and lighting estimation techniques to simultaneously
recover the lighting and BRDF parameters, when both are
unknown. In this case, we do not show direct recovery of spherical
harmonic coefficients, as we have thus far found this to be imprac-
tical. Finally, section 6.4 demonstrates our algorithms on geometrically
complex objects, showing how it is straightforward to extend
our model to handle textures and shadowing.
To test our methods, we first used homogeneous spheres 5 of different
materials. Spheres are naturally parameterized with spherical
coordinates, and therefore correspond directly to our theory.
Later, we also used complex objects-a white cat sculpture, and a
textured wooden doll-to show the generality of our algorithms.
Data Acquisition: We used a mechanical gantry to position an
inward-looking Toshiba IK-TU40A CCD(x3) camera on an arc of
radius 60cm. Calibration of intrinsics was done by the method of
Zhang [40]. Since the camera position was computer-controlled,
extrinsics were known. The mapping between pixel and radiance
values was also calibrated. We acquired 60 images of the target
sphere, taken at 3 degree intervals. To map from image pixels to
angular coordinates (#
used image silhouettes to
find the geometric location of the center of the sphere and its radius.
Our gantry also positioned a 150W white point source along an
arc. Since this arc radius (90 cm) was much larger than the sphere
(between 1.25 and 2cm), we treated the point source as a directional
light. A large area source, with 99% of its energy in low-frequency
modes of order l # 6, was obtained by projecting white
light on a projection screen. The lighting distribution was determined
using a gazing sphere. This information was used directly
for experiments assuming known illumination, and as a reference
solution for experiments assuming unknown illumination.
We also used the same experimental setup, but with only the
point source, to measure the BRDF of a white teflon sphere using
the image-based method of Marschner et al. [21]. This independent
measurement was used to verify the accuracy of our BRDF
estimation algorithms under complex illumination.
5 Ordered from the McMaster-Carr catalog http://www.mcmaster.com
6.1 Inverse BRDF with Known Lighting
Estimation of Spherical Harmonic BRDF coefficients:
Spherical harmonics and Zernike polynomials have been fit [14] to
measured BRDF data, but previous work has not tried to estimate
coefficients directly. Since the BRDF is linear in the coefficients
# lpq , we simply solve a linear system to determine -
# lpq .

Figure

6 compares the parametric BRDFs estimated under complex
lighting to BRDFs measured using a single point source with
the method of Marschner et al. [21]. As expected [14], the recovered
BRDFs exhibit ringing. One way to reduce ringing is to attenuate
high-frequency coefficients. According to our theory, this
is equivalent to using low frequency lighting. Therefore, as seen
in figure 6, images rendered with low-frequency lighting do not
exhibit ringing and closely match real photographs, since only the
low-frequency components of the BRDF are important. However,
images rendered using directional sources show significant ringing.
Real (Marschner) Order 12
Order 6
Images
f
slices
Real
Rendered
Low-Frequency Lighting
Rendered
Directional Source

Figure

Top: Slices of the BRDF transfer function of a teflon sphere for fixed
exitant angle of 63 # i varies linearly from 0 # to 90 # from top to bottom, and
| linearly from 0 # to 360 # from left to right. The central bright feature
is the specular highlight. Left is the BRDF slice independently measured using the
approach of Marschner et al. [21], middle is the recovered value using a maximum
order 6, and right is the recovered version for order 12. Ringing is apparent in both
recovered BRDFs. The right version is sharper, but exhibits more pronounced ringing.
Bottom: Left is an actual photograph; the lighting is low-frequency from a large area
source. Middle is a rendering using the recovered BRDF for order 6 and the same
lighting. Since the lighting is low-frequency, only low-frequency components of the
BRDF are important, and the rendering appears very similar to the photograph even
though the recovered BRDF does not include frequencies higher than order 6. Right
is a rendering with a directional source at the viewpoint, and exhibits ringing.
For practical applications, it is usually more convenient to recover
low-parameter BRDF models since these are compact, can
be estimated from fewer observations, and do not exhibit ringing.
In the rest of this section, we will derive improved inverse rendering
algorithms, assuming our parametric BRDF model.
Estimation of Parametric BRDF Model: We estimate
BRDF parameters under general known lighting distributions using
equation 24. The inputs are images that sample the reflected light
field B. We perform the estimation using nested procedures. In the
outer procedure, a simplex algorithm adjusts the nonlinear parameters
- and # to minimize error with respect to image pixels. In
the inner procedure, a linear problem is solved for Kd and Ks . For
numerical work, we use the simplex method e04ccc and linear
solvers f01qcc and f01qdc in the NAG [9] C libraries. The
main difference from previous work is that equation 24 provides a
principled way of accounting for all components of the lighting and
BRDF, allowing for the use of general illumination conditions.
We tested our algorithm on the spheres. Since the lighting includes
high and low-frequency components (a directional source
and an area source), the theory predicts that parameter estimation
is well-conditioned. To validate our algorithm, we compared parameters
recovered under complex lighting for one of the sam-
ples, a white teflon sphere, to those obtained by fitting to the full
BRDF separately measured by us using the method of Marschner et
al. [21]. Unlike most previous work on BRDF estimation, we consider
the Fresnel term. It should be noted that accurate estimates
for the refractive index - require correct noise-free measurements
at grazing angles. Since these measurements tend to be the most
error-prone, there will be small errors in the estimated values of -
for some materials. Nevertheless, we find the Fresnel term important
for reproducing accurate specular highlights at oblique angles.
Parameter Our Method Fit to Data
Reflectance 0.86 0.87
Kd/(Kd +Ks) 0.89 0.91
Ks/(Kd +Ks) 0.11 0.09
RMS 9.3% 8.5%

Figure

7: Comparison of BRDF parameters recovered by our algorithm under complex
lighting to those fit to measurements made by the method of Marschner et al. [21].
The results in figure 7 show that the estimates of BRDF parameters
from our method are quite accurate, and there is only a small
increase in the error-of-fit when using parameters recovered by our
algorithm to fit the measured BRDF. We also determined percentage
RMS errors between images rendered using recovered BRDFs
and real photographs to be between 5 and 10%. A visual comparison
is shown in the first and third rows of figure 12. All these results
indicate that, as expected theoretically, we can accurately estimate
BRDFs even under complex lighting.
6.2 Inverse Lighting with Known BRDF
Previous methods for estimating the lighting have been developed
only for the special cases of mirror BRDFs (a gazing sphere), Lambertian
BRDFs (Marschner and Greenberg [20]), and when shadows
are present (Sato et al. [30]). Previous methods [20, 30] have
also required regularization using penalty terms with user-specified
weights, and have been limited by the computational complexity
of their formulations to a coarse discretization of the sphere.
We present two new algorithms for curved surfaces with general
BRDFs. The first method directly recovers spherical harmonic
lighting coefficients Llm . The second algorithm estimates parameters
of the dual angular and frequency-space lighting model of
section 5.3. This method requires no explicit regularization, and
yields high-resolution results that are sharper than those from the
first algorithm, but is more difficult to extend to concave surfaces.
The theory tells us that inverse lighting is ill-conditioned for
high-frequencies. Therefore, we will recover only low-frequency
continuous lighting distributions, and will not explicitly account
for directional sources, i.e. we assume that B s,fast = 0. The reflected
light field is essentially independent of the surface roughness
# under these conditions, so our algorithms do not explicitly
use #. The theory predicts that the recovered illumination will be a
filtered version of the real lighting. Directional sources will appear
as continuous distributions of angular width approximately #.
Estimation of Spherical Harmonic Lighting coefficients:
We represent the lighting by coefficients Llm with l # l # , and solve
a linear least-squares system for Llm . The first term in parentheses
below corresponds to Bd , and the second to B s,slow . The cutoff l #
is used for regularization, and should be of order l # -1 .
l #
l
m=-l
Estimation of Parametric Dual Lighting Model: Another
approach is to estimate the dual angular and frequency-space lighting
model of section 5.3. Our algorithm is based on subtracting out
the diffuse component Bd of the reflected light field. After this,
we treat the object as a mirror sphere, recovering a high-resolution
angular-space version of the illumination from the specular component
alone. To determine Bd , we need only the 9 lowest frequency-space
coefficients Llm with l # 2. Our algorithm uses the following
methods to convert between angular and frequency-space:
Phase 2
Phase 1
Input
f
s,slow
dB d1
lm lm

Figure

8: Estimation of dual lighting representation. In phase 1, we use frequency-space
parameters
lm to compute diffuse component B 1
d . This is subtracted from the
input image, leaving the specular component, from which the angular-space lighting
is found. In phase 2, we compute coefficients L 2
lm which can be used to determine
d
. The consistency condition is that B 1
d
or L 1
lm
In this and all
subsequent figures, the lighting is visualized by unwrapping the sphere so # ranges in
equal increments from 0 to # from top to bottom, and # ranges in equal increments
from 0 to 2# from left to right (so the image wraps around in the horizontal direction).
1. 9 parameters to High-Resolution Lighting: The inputs to
phase 1 are the coefficients L 1
lm . These suffice to find B 1
d by
equation 21. Since we assumed that B s,fast = 0,
Lslow
We assume the BRDF parameters are known, and B is the
input to the algorithm, so the right-hand side can be evaluated.
2. High-Resolution Lighting to 9 parameters: Using the angular
space values L found from the first phase, we can easily
find the 9 frequency-space parameters of the lighting L 2
lm .
Now, assume we run phase 1 (with inputs L 1
lm ) and phase 2
(with outputs L 2
sequentially. The consistency condition is that
lm -converting from frequency to angular to frequency
space must not change the result. Equivalently, the computed diffuse
components must match, i.e. B 1
lm ). This is
illustrated in figure 8. Since everything is linear in terms of the
lighting coefficients, the consistency condition reduces to a system
of 9 simultaneous equations. After solving for Llm , we run phase
1 to determine the high-resolution lighting in angular space.

Figure

9 compares the methods to each other, and to a reference
solution from a gazing sphere. Both algorithms give reasonably
accurate results. As predicted by the theory, high-frequency
components are filtered by the roughness #. In the first method,
involving direct recovery of Llm , there will still be some residual
energy for l > l # . Since we regularize by not considering
higher frequencies-we could increase l # , but this makes the result
noisier-the recovered lighting is somewhat blurred compared to
our dual angular and frequency-space algorithm (second method).

Real (Gazing Sphere) Algorithm 1
f

Figure

9: Comparison of inverse lighting methods. From left to right, real lighting
(from a gazing sphere), recovered illumination by direct estimation of spherical harmonic
coefficients with l estimation of dual angular and
frequency-space lighting model. To make the artifacts more apparent, we have set 0
to gray. The results from the dual algorithm are sharper, but still somewhat blurred
because of filtering by #. A small amount of ringing occurs for direct coefficient re-
covery, and can be seen for l Using l makes the solution very noisy.
6.3 Factorization-Unknown Lighting and BRDF
We can combine the inverse-BRDF and inverse-lighting methods
to factor the reflected light field, simultaneously recovering the
lighting and BRDF when both are unknown. Therefore, we are
able to recover BRDFs of curved surfaces under unknown complex
illumination, something which has not previously been demon-
strated. There is an unrecoverable global scale factor, so we set
we cannot find absolute reflectance. Also, the
theory predicts that for low-frequency lighting, estimation of the
surface roughness # is ill-conditioned-blurring the lighting while
sharpening the BRDF does not significantly change the reflected
light field. However, for high-frequency lighting, this ambiguity
can be removed. We will use a single manually specified directional
source in the recovered lighting distribution to estimate #.
Algorithm: The algorithm consists of nested procedures. In
the outer loop, we effectively solve an inverse-BRDF problem-a
nonlinear simplex algorithm adjusts the BRDF parameters to
minimize error with respect to image pixels. Since Kd
and # will not be solved for till after the lighting and other BRDF
parameters have been recovered, there are only
Ks and -. In the inner procedure, a linear problem is solved to
estimate the lighting for a given set of BRDF parameters, using the
methods of the previous subsection. Pseudocode is given below.
global Binput // Input images
global Kd ,Ks ,-,# // BRDF parameters
global L // Lighting
procedure Factor
Minimize(Ks ,-,ObjFun) // Simplex Method
Figure 10, Equation 26
function ObjFun(Ks ,-)
Lighting
Field
return RMS(Binput ,Bpred ) // RMS Error
Finding # using a directional source: If a directional
source is present-and manually specified by us in the recovered
lighting-we can estimate # by equating specular components predicted
by equations 22 and 23 for the center, i.e. brightest point, of
the light source at normal exitance. An illustration is in figure 10.
total
f
tot

Figure

10: We manually specify (red box) the region corresponding to the directional
source in a visualization of the lighting. The algorithm then determines Lcen ,
the intensity at the center (brightest point), L tot the total energy integrated over the
region specified by the red box, and computes # using equation 26. The method does
not depend on the size of the red box-provided it encloses the entire (filtered) source-
nor the precise shape into which the source is filtered in the recovered lighting.
Results: We used the method of this subsection-with the dual
angular and frequency-space algorithm for inverse lighting-to factor
the light field for the spheres, simultaneously estimating the
BRDF and lighting. The same setup and lighting were used for
all the spheres so we could compare the recovered illumination.
We see from figure 11 that the BRDF estimates under unknown
lighting are accurate. Absolute errors are small, compared to parameters
recovered under known lighting. The only significant
anomalies are the slightly low values for the refractive index -
caused because we cannot know the high-frequency lighting com-
ponents, which are necessary for more accurately estimating the
Fresnel term. We are also able to estimate a filtered version of the
lighting. As shown in figure 12, the recovered lighting distributions
from all the samples are largely consistent. As predicted by
the theory, the directional source is spread out to different extents
depending on how rough the surface is, i.e. the value of #. Finally,
figure 12 shows that rendered images using the estimated lighting
and BRDF are almost indistinguishable from real photographs.
Material K d Ks - #
Known Unknown Known Unknown Known Unknown Known Unknown
Teflon
Delrin 0.87 0.88 0.13 0.12 1.44 1.35 0.10 0.11
Neoprene Rubber 0.92 0.93 0.08 0.07 1.49 1.34 0.10 0.10
Sandblasted Steel 0.20 0.14 0.80 0.86 0.20 0.19
Bronze
Painted (.62,.71,.62) (.67,.75,.64) 0.29 0.25 1.38 1.15 0.15 0.15

Figure

11: BRDFs of various spheres, recovered under known (section 6.1) and unknown (section 6.3) lighting. The reported values are normalized so K d
values are reported for colored objects. We see that Ks is much higher for the more specular metallic spheres, and that # is especially high for the rough sandblasted sphere. The
Fresnel effect is very close to 1 for metals, so we do not consider the Fresnel term for these spheres.
Lighting
Known
Lighting
Unknown
Lighting
Images
Images
Images
Bronze
Sandblasted Painted
Delrin
f
Teflon
of real lighting
Filtered version
Real
Rendered
Rendered
Recovered
Real lighting

Figure

12: Spheres rendered using BRDFs estimated under known (section 6.1) and
unknown (section 6.3) lighting. The algorithm in section 6.3 also recovers the lighting.
Since there is an unknown global scale, we scale the recovered lighting distributions
in order to compare them. The recovered illumination is largely consistent between
all samples, and is similar to a filtered version of the real lighting. As predicted by
the theory, the different roughnesses # cause the directional source to be spread out to
different extents. The filtered source is slightly elongated or asymmetric because the
microfacet BRDF is not completely symmetric about the reflection vector.
6.4 Complex Objects-Texture and Shadowing
We now demonstrate our algorithms on objects with complex ge-
ometry, and discuss extensions to handle concave surfaces and textured
objects. Although the theory is developed for homogeneous
surfaces, our algorithms can be extended to textured objects simply
by letting the BRDF parameters be functions of surface position. It
would appear that concave regions, where one part of the surface
may shadow another, are a more serious problem since our theory
is developed for convex objects and assumes no self-shadowing.
However, using our new practical model of section 5.3, we will
see that the extensions necessary mainly just involve checking for
shadowing of the reflected ray and directional sources, which are
routine operations in a raytracer.
Shadowing-Concave Surfaces: In our practical model, the
reflected light field consists of 3 parts-Bd , B s,slow , and B s,fast .
s,slow depends on Lslow ( #
R), the slowly-varying component of
the lighting evaluated at the reflection vector. Our model allows us
to approximate the effects of shadowing simply by checking if the
reflected ray is shadowed. The other components are handled in the
standard manner. To consider shadowing when computing B s,fast ,
corresponding to specularities from directional sources, we check
if these sources are shadowed. Bd depends on the irradiance E,
which should now be computed in the more conventional angular-
space way by integrating the scene lighting while considering vis-
ibility, instead of using the 9-parameter lighting approximation of
equation 21. It should be emphasized that in all cases, the corrections
for visibility depend only on object geometry, and can be
precomputed for each point on the object using a ray tracer.
For parametric BRDF estimation, we modify each component of
equation 24 to consider visibility, as discussed above. Our first inverse
lighting method, that directly recovers the coefficients Llm , is
modified similarly. In equation 25, we check if the reflected ray is
shadowed, and consider shadowing when computing the irradiance
due to each Ylm . Note that B is still a linear combination of the
lighting coefficients, so we will still solve a linear system for Llm .
However, it is difficult to extend our dual angular and frequency-space
method for inverse lighting to handle concave surfaces because
Bd no longer depends only on the 9 lighting coefficients Llm
with l # 2. For light field factorization, we simply extend both the
and inverse-lighting methods as discussed.
A white cat sculpture was used to test our algorithms on complex
geometric objects that include concavities. Geometry was acquired
using a Cyberware range scanner and aligned to the images by manually
specifying correspondences. The lighting was slightly more
complex than that for the spheres experiment; we used a second
directional source in addition to the area source.
To show that we can recover BRDFs using a small number of
images, we used only 3 input photographs. We recovered BRDFs
under both known lighting, using the method of section 6.1, and
unknown lighting-using the factorization method of section 6.3,
with the inverse lighting component being direct recovery of spherical
harmonic coefficients using l Comparisons of photographs
and renderings are in figures 1 and 13. BRDF and lighting
parameters are tabulated in figure 14. This experiment indicates
that our methods for BRDF recovery under known and unknown
lighting are consistent, and are accurate even for complex lighting
and geometry. The rendered images are very close to the original
photographs, even under viewing and lighting conditions not used
for BRDF recovery. The most prominent artifacts are because of
imprecise geometric alignment and insufficient geometric resolu-
tion. For instance, since our geometric model does not include the
eyelids of the cat, that feature is missing from the rendered images.
Textured BRDFs: Since the theory shows that factorization of
lighting and texture is ambiguous, we consider only recovery of
textured BRDFs under known lighting. It is fairly straightforward
to allow Kd (x) and Ks (x) to be described by textures that depend
on surface position x. In the inner procedure of the parametric
BRDF estimation algorithm of section 6.1, we simply solve a separate
linear problem for each point x to estimate Kd (x) and Ks(x).
As an experimental test, we used a wooden doll. We compared
the real input photographs with images rendered using the recovered
textured BRDF. We also took a photograph of the same object
under a single directional source and compared this to a rendering
using the textured BRDF recovered under complex illumination.
The results in figure 15 show that our renderings closely resemble
real photographs. The main artifact is blurring of texture because
of geometry-image misregistration.
7 Conclusions and Future Work
This paper has developed a signal-processing framework for inverse
rendering. The qualitative observation that the reflected light
field is a convolution of the lighting and BRDF has been formalized
mathematically. We have shown in frequency-space why a gazing
sphere is well-suited for recovering the lighting-the frequency
spectrum of the mirror BRDF (a delta function) is constant-and
Input
View
View
Rendered, Known Lighting Real Photograph Rendered, Unknown Lighting

Figure

13: Comparison of real photographs (middle column) to images rendered
using BRDFs recovered under known lighting (left column), and using BRDFs (and
lighting) estimated under unknown lighting (right column). The top row is one of the
3 input views. The bottom row is a new view, not used for BRDF estimation. Note
that in the top row, we have composited the left and right renderings over the same
background as the middle photograph in order to make a meaningful comparison.
Parameter Known Lighting Unknown Lighting
BRDF Parameters
K d 0.88 0.90
Ks 0.12 0.10
Lighting Coefficients (l,m)

Figure

14: BRDF and lighting parameters for the cat sculpture. We see good agreement
between BRDF parameter values recovered with known and unknown lighting,
showing our methods are consistent. Note that we normalize so Kd
We may also check the accuracy of the recovered lighting. Since there is an unknown
global scale for the recovered values, we report normalized lighting coefficient values
for the first 9 spherical harmonic coefficients (in real form), which are the most
important, because they significantly affect the diffuse component of the BRDF.
why a directional source is well-suited for recovering the BRDF-
we are estimating the BRDF filter by considering its impulse re-
sponse. The conditioning properties and well-posedness of BRDF
and lighting estimation under various conditions have been de-
rived, as well as an explicit formula for factoring the reflected light
field into the lighting and BRDF. The ill-conditioning observed by
Marschner and Greenberg [20] in estimating the lighting from a
Lambertian surface has been explained, and we have shown that
factorization of lighting effects and low-frequency texture is am-
biguous. All these results indicate that the theory provides a useful
analytical tool for studying the properties of inverse problems.
The insights gained from the theory also lead to a new practical
representation. We can numerically represent quantities in angular
or frequency space, depending on where they are more lo-
cal. This leads to new algorithms which are often expressed in a
combination of angular and frequency-space. We can determine
which BRDF and lighting parameters are important, and can handle
the various components appropriately. For BRDF estimation,
the parametric recovery algorithms of Yu and Malik [39], Sato
and Ikeuchi [31], and Love [17]-which are designed specifically
for natural lighting-can be seen as special cases of this general
Rendered Real Rendered
1 view in original input sequence Same view, novel lighting
Real

Figure

15: Recovering textured BRDFs under complex lighting. The rendered images
closely resemble the real photographs, even under novel lighting.
they treat sunlight (high-frequency) and skylight (low-
separately. We provide a general framework for arbitrary
illumination, and also determine conditions under which parameter
recovery is robust. For instance, our theory predicts that
estimation of # is ill-conditioned on a cloudy day, with only low-frequency
lighting. Our framework can also be applied to developing
new frequency-space algorithms to estimate the lighting from
objects with general BRDFs. The use of frequency-space naturally
handles continuous lighting distributions. Our dual angular
and frequency-space algorithm effectively reduces the problem for
general BRDFs to that for a gazing sphere, requires no explicit reg-
ularization, and allows much higher resolutions to be obtained than
with previous purely angular-space methods [20, 30]. Finally, we
demonstrate a method for factoring the light field to simultaneously
estimate the lighting and BRDF. This allows us to estimate BRDFs
of geometrically complex objects under unknown general lighting,
which has not previously been demonstrated.
We have only scratched the surface of possible applications. In
the future, it is likely that many more algorithms can be derived
using the basic approaches outlined here. Possible algorithmic improvements
include extending the consistency condition for inverse
lighting so we can use color-space methods [13] to help separate
diffuse and specular components for colored objects. Finally, while
we have discussed only inverse rendering applications, we believe
the convolution-based approach is of theoretical and practical importance
in many other areas of computer graphics. We have already
shown [28] how to use the 9 term irradiance approximation
for efficient forward rendering of diffuse objects with environment
maps, and we believe there are many further applications.

Acknowledgements

: We are grateful to Marc Levoy for many helpful initial
discussions regarding both the interpretation of reflection as convolution and the
practical issues in inverse rendering, as well as for reading a draft of the paper. We
thank Szymon Rusinkiewicz for many suggestions regarding this research, for reading
early versions of the paper, and for being ever willing to assist us in debugging problems
with the gantry. Steve Marschner also deserves thanks for detailed comments on
many early drafts. Jean Gleason at Bal-tec assisted in supplying us with our mirror
sphere. Finally, we thank John Parissenti at Polytec Products for much useful advice
on obtaining uniform spheres, and for getting one of our steel spheres sandblasted.
The cat is a range scan of the sculpture Serenity by Sue Dawes. The work described in
this paper was supported in part by a Hodgson-Reed Stanford graduate fellowship and
NSF ITR grant #0085864 "Interacting with the Visual World."



--R

Lambertian reflectance and linear subspaces.
Increased photorealism for interactive architectural walkthroughs.
Bidirectional reflection functions from surface bump maps.
A volumetric method for building complex models from range images.
Reflectance and texture of real-world surfaces
Acquiring the reflectance field of a human face.
Estimating surface reflectance properties from images under unknown illumination.
Computational Models of Visual Processing
Numerical Algorithms Group.
Determining reflectance properties of an object using range and brightness images.
Group theory and its applications in physics.
An image based measurement system for anisotropic reflection.
The measurement of highlights in color images.
Phenomenological description of bidirectional surface reflection.
Lightness and retinex theory.

Surface Reflection Model Estimation from Naturally Illuminated Image Sequences.
Optical properties (bidirectional reflection distribution functions) of velvet.
Spherical harmonics
Inverse lighting for photography.

Inverse radiative transfer problems: a review.
Simulated objects in simulated and real environments.
Geometric Considerations and Nomenclature for Reflectance.
Efficient re-rendering of naturally illuminated environments
Hydrologic Optics.
Analysis of planar light fields from homogeneous convex curved surfaces under distant illumination.
An efficient representation for irradiance environment maps.
On the relationship between radiance and irradiance: Determining the illumination from images of a convex lambertian object.
Illumination distribution from brightness in shadows: adaptive estimation of illumination distribution with unknown reflectance properties in shadow regions.
Reflectance analysis under solar illumination.
Object shape and reflectance modeling from observation.
A global illumination solution for general reflectance distributions.
Inverse Rendering for Computer Graphics.
Estimating reflection parameters from a single color image.
Theory for off-specular reflection from roughened surfaces
In SIGGRAPH 92
Inverse global illumination: Recovering reflectance models of real scenes from photographs.
Recovering photometric properties of architectural scenes from photographs.
A flexible new technique for camera calibration.
--TR
Bidirectional reflection functions from surface bump maps
A global illumination solution for general reflectance distributions
Determining Reflectance Properties of an Object Using Range and Brightness Images
Measuring and modeling anisotropic reflection
A volumetric method for building complex models from range images
Object shape and reflectance modeling from observation
Recovering photometric properties of architectural scenes from photographs
Increased photorealism for interactive architectural walkthroughs
Reflectance and texture of real-world surfaces
Inverse global illumination
The digital Michelangelo project
Acquiring the reflectance field of a human face
A Flexible New Technique for Camera Calibration
An efficient representation for irradiance environment maps
Estimating Reflection Parameters from a Single Color Image
Reflectance Analysis Under Solar Illumination
Inverse rendering for computer graphics

--CTR
Makoto Okabe , Yasuyuki Matsushita , Takeo Igarashi , Heung-Yeung Shum, Illumination interactive design of image-based lighting, ACM SIGGRAPH 2006 Research posters, July 30-August 03, 2006, Boston, Massachusetts
G. H. Hu , S. K. Ong , Y. P. Chen , A. Y. C. Nee, Technical Section: Reflectance modeling for a textured object under uncontrolled illumination from high dynamic range maps, Computers and Graphics, v.31 n.2, p.262-270, April, 2007
Martin Fuchs , Volker Blanz , Hendrik Lensch , Hans-Peter Seidel, Reflectance from Images: A Model-Based Approach for Human Faces, IEEE Transactions on Visualization and Computer Graphics, v.11 n.3, p.296-305, May 2005
Gustavo Patow , Xavier Pueyo , Alvar Vinacua, Technical Section: User-guided inverse reflector design, Computers and Graphics, v.31 n.3, p.501-515, June, 2007
Kuang-Chih Lee , Jeffrey Ho , David J. Kriegman, Acquiring Linear Subspaces for Face Recognition under Variable Lighting, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.27 n.5, p.684-698, May 2005
Ravi Ramamoorthi , Melissa Koudelka , Peter Belhumeur, A Fourier Theory for Cast Shadows, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.27 n.2, p.288-295, February 2005
Enhua Wu , Qimin Sun , Xuehui Liu, Recovery of material under complex illumination conditions, Proceedings of the 2nd international conference on Computer graphics and interactive techniques in Australasia and South East Asia, June 15-18, 2004, Singapore
Ning Xu , Narendra Ahuja, Shape and View Independent Reflectance Map from Multiple Views, International Journal of Computer Vision, v.73 n.2, p.123-138, June 2007
Bei Xiao , David H. Brainard, Color perception of 3D objects: constancy with respect to variation of surface gloss, Proceedings of the 3rd symposium on Applied perception in graphics and visualization, July 28-29, 2006, Boston, Massachusetts
Ravi Ramamoorthi, Analytic PCA Construction for Theoretical Analysis of Lighting Variability in Images of a Lambertian Object, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.24 n.10, p.1322-1333, October 2002
Jan Kautz , Peter-Pike Sloan , John Snyder, Fast, arbitrary BRDF shading for low-frequency lighting using spherical harmonics, Proceedings of the 13th Eurographics workshop on Rendering, June 26-28, 2002, Pisa, Italy
Raanan Fattal , Dani Lischinski , Michael Werman, Gradient domain high dynamic range compression, ACM Transactions on Graphics (TOG), v.21 n.3, July 2002
Sara Keren , Ilan Shimshoni , Ayellet Tal, Placing three-dimensional models in an uncalibrated single image of an architectural scene, Proceedings of the ACM symposium on Virtual reality software and technology, November 11-13, 2002, Hong Kong, China
Yang Wang , Dimitris Samaras, Estimation of multiple directional light sources for synthesis of augmented reality images, Graphical Models, v.65 n.4, p.185-205, July
Jefferson Y. Han , Ken Perlin, Measuring bidirectional texture reflectance with a kaleidoscope, ACM Transactions on Graphics (TOG), v.22 n.3, July
Imari Sato , Takahiro Okabe , Yoichi Sato, Appearance Sampling of Real Objects for Variable Illumination, International Journal of Computer Vision, v.75 n.1, p.29-48, October   2007
Paul Debevec , Andreas Wenger , Chris Tchou , Andrew Gardner , Jamie Waese , Tim Hawkins, A lighting reproduction approach to live-action compositing, ACM Transactions on Graphics (TOG), v.21 n.3, July 2002
Chang Ha Lee , Xuejun Hao , Amitabh Varshney, Light Collages: Lighting Design for Effective Visualization, Proceedings of the conference on Visualization '04, p.281-288, October 10-15, 2004
Jan Kautz , Peter-Pike Sloan , Jaakko Lehtinen, Precomputed radiance transfer: theory and practice, ACM SIGGRAPH 2005 Courses, July 31-August
Jaakko Lehtinen , Jan Kautz, Matrix radiance transfer, Proceedings of the symposium on Interactive 3D graphics, April 27-30, 2003, Monterey, California
Kenji Hara , Ko Nishino , Katsushi Ikeuchi, Light Source Position and Reflectance Estimation from a Single View without the Distant Illumination Assumption, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.27 n.4, p.493-505, April 2005
Ravi Ramamoorthi , Dhruv Mahajan , Peter Belhumeur, A first-order analysis of lighting, shading, and shadows, ACM Transactions on Graphics (TOG), v.26 n.1, p.2-es, January 2007
Martin Fuchs , Volker Blanz , Hendrik P.A. Lensch , Hans-Peter Seidel, Adaptive sampling of reflectance fields, ACM Transactions on Graphics (TOG), v.26 n.2, p.10-es, June 07
Srinivasa G. Narasimhan , Mohit Gupta , Craig Donner , Ravi Ramamoorthi , Shree K. Nayar , Henrik Wann Jensen, Acquiring scattering properties of participating media by dilution, ACM Transactions on Graphics (TOG), v.25 n.3, July 2006
Xuejun Hao , Thomas Baby , Amitabh Varshney, Interactive subsurface scattering for translucent meshes, Proceedings of the symposium on Interactive 3D graphics, April 27-30, 2003, Monterey, California
Ravi Ramamoorthi , Pat Hanrahan, Frequency space environment map rendering, ACM Transactions on Graphics (TOG), v.21 n.3, July 2002
Sara Keren , Ilan Shimshoni , Ayellet Tal, Placing Three-Dimensional Models in an Uncalibrated Single Image of an Architectural Scene, Presence: Teleoperators and Virtual Environments, v.13 n.6, p.692-707, December 2004
Image-based skin color and texture analysis/synthesis by extracting hemoglobin and melanin information in the skin, ACM Transactions on Graphics (TOG), v.22 n.3, July
David Akers , Frank Losasso , Jeff Klingner , Maneesh Agrawala , John Rick , Pat Hanrahan, Conveying Shape and Features with Image-Based Relighting, Proceedings of the 14th IEEE Visualization 2003 (VIS'03), p.46, October 22-24,
Imari Sato , Yoichi Sato , Katsushi Ikeuchi, Illumination from Shadows, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.25 n.3, p.290-300, March
G. Narasimhan , Shree K. Nayar, A practical analytic single scattering model for real time rendering, ACM Transactions on Graphics (TOG), v.24 n.3, July 2005
Andrew Gardner , Chris Tchou , Tim Hawkins , Paul Debevec, Linear light source reflectometry, ACM Transactions on Graphics (TOG), v.22 n.3, July
Michael D. Grossberg , Shree K. Nayar, Modeling the Space of Camera Response Functions, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.26 n.10, p.1272-1282, October 2004
Frdo Durand , Nicolas Holzschuch , Cyril Soler , Eric Chan , Franois X. Sillion, A frequency analysis of light transport, ACM Transactions on Graphics (TOG), v.24 n.3, July 2005
Wojciech Matusik , Hanspeter Pfister , Addy Ngan , Paul Beardsley , Remo Ziegler , Leonard McMillan, Image-based 3D photography using opacity hulls, ACM Transactions on Graphics (TOG), v.21 n.3, July 2002
Wei-Chao Chen , Jean-Yves Bouguet , Michael H. Chu , Radek Grzeszczuk, Light field mapping: efficient representation and hardware rendering of surface light fields, ACM Transactions on Graphics (TOG), v.21 n.3, July 2002
Ronen Basri , David W. Jacobs, Lambertian Reflectance and Linear Subspaces, IEEE Transactions on Pattern Analysis and Machine Intelligence, v.25 n.2, p.218-233, February
Ko Nishino , Shree K. Nayar, Eyes for relighting, ACM Transactions on Graphics (TOG), v.23 n.3, August 2004
Xuejun Hao , Amitabh Varshney, Real-time rendering of translucent meshes, ACM Transactions on Graphics (TOG), v.23 n.2, p.120-142, April 2004
Steven M. Seitz , Kiriakos N. Kutulakos, Plenoptic Image Editing, International Journal of Computer Vision, v.48 n.2, p.115-129, July 2002
Shree K. Nayar , Peter N. Belhumeur , Terry E. Boult, Lighting sensitive display, ACM Transactions on Graphics (TOG), v.23 n.4, p.963-979, October 2004
Hendrik P. A. Lensch , Jan Kautz , Michael Goesele , Wolfgang Heidrich , Hans-Peter Seidel, Image-based reconstruction of spatial appearance and geometric detail, ACM Transactions on Graphics (TOG), v.22 n.2, p.234-257, April
Ravi Ramamoorthi , Pat Hanrahan, A signal-processing framework for reflection, ACM Transactions on Graphics (TOG), v.23 n.4, p.1004-1042, October 2004
Rodrigo L. Carceroni , Kiriakos N. Kutulakos, Multi-View Scene Capture by Surfel Sampling: From Video Streams to Non-Rigid 3D Motion, Shape and Reflectance, International Journal of Computer Vision, v.49 n.2-3, p.175-214, September-October 2002
Hendrik P. A. Lensch , Michael Goesele , Yung-Yu Chuang , Tim Hawkins , Steve Marschner , Wojciech Matusik , Gero Mueller, Realistic materials in computer graphics, ACM SIGGRAPH 2005 Courses, July 31-August
Hendrik P. A. Lensch , Michael Goesele , Yung-Yu Chuang , Tim Hawkins , Steve Marschner , Wojciech Matusik , Gero Mueller, Realistic materials in computer graphics, ACM SIGGRAPH 2005 Courses, July 31-August
light field;irradiance;inverse rendering;illumination;BRDF;radiance;spherical harmonics;signal processing