Barycentric Parameterizations for Isotropic BRDFs - Semantic Scholar
1
Barycentric Parameterizations for Isotropic BRDFs Michael M. Stark1 , James Arvo1,2 , and Brian Smits2
Abstract— A bidirectional reflectance distribution function (BRDF) is often expressed as a function of four real variables: two spherical coordinates in each of the the “incoming” and “outgoing” directions. However, many BRDFs reduce to functions of fewer variables. For example, isotropic reflection can be represented by a function of three variables. Some BRDF models can be reduced further. In this paper we introduce new sets of coordinates which we use to reduce the dimensionality of several well-known analytic BRDFs as well as empirically measured BRDF data. The proposed coordinate systems are barycentric with respect to a triangular support with a direct physical interpretation. One coordinate set is based on the BRDF model proposed by Lafortune. Another set, based on a model of Ward, is associated with the “halfway” vector common in analytical BRDF formulas. Through these coordinate sets we establish lower bounds on the approximation error inherent in the models on which they are based. We present a third set of coordinates, not based on any analytical model, that performs well in approximating measured data. Finally, our proposed variables suggest novel ways of constructing and visualizing BRDFs. Index Terms— Computer Graphics; Picture/Image Generation; Display algorithms; Color, shading, shadowing, and texture; Reflectance
I. I NTRODUCTION
T
HE faithful simulation of surface reflection is an essential element of accurate image synthesis [1], [2]. The most fundamental relationship governing surface reflection can be stated in terms of radiometric quantities associated with a differential surface area: the differential radiance reflected from a given incident direction to a given outgoing direction is proportional to the irradiance due to a differential solid angle about the incident direction. In symbols, dL(v) ∝ L(u) d ω (1) where L denotes radiance, ω denotes solid angle, and u and v are unit vectors of the incident and reflected directions, respectively. The constant of proportionality depends on both directions, so (1) may be written as dL(v) = ρ (u, v)L(u) d ω .
(2)
The function ρ (u, v) is known as the bidirectional reflectance distribution function, or BRDF, and carries physical units of inverse steradians. This work focuses on the BRDF at a single surface point and at a single wavelength; the concepts are easily generalized to include spatiallyvarying BRDFs and wavelength dependence. Manuscript submitted January 31, 2004 of California, Irvine 2 Pixar Animation Studios 1 University
BRDFs are typically treated as real-valued functions ρ : Ω × Ω → R, where Ω denotes the unit hemisphere. In this representation a BRDF can be expressed as a function of four real variables: ρ (θi, φi , θo , φo ), where θi and φi are polar and azimuthal angles of the incident direction, respectively; θo and φo , of the reflected direction. A BRDF function can be represented more succinctly as a function ρ (u, v) of two (unit) vectors u and v, representing the incoming and outgoing direction vectors on the unit hemisphere. Because a BRDF models the physical process of light reflection, not every function f : Ω × Ω → R represents a physically plausible BRDF. In particular, a BRDF is subject to physical constraints: the second law of thermodynamics requires a BRDF to obey the reciprocity principle
ρ (θi , φi , θo , φo ) = ρ (θo , φo , θi , φi ), ρ (u, v) = ρ (v, u),
(3) (4)
and the first law requires it to conserve energy: Z
Ω
ρ (ω , v) cos θi d ω ≤ 1
(5)
for any fixed v. Reciprocity ensures that a BRDF is symmetric with respect to its incoming and outgoing directions; energy conservation ensures that no more energy is scattered out than comes in. Outside of these constraints, BRDFs can have a wide variety forms. Numerous mathematical BRDF models have found application in a variety of fields, including radiative transfer, remote sensing, computer vision, and computer graphics [3]–[7]. In addition, empirical BRDF measurements are made using a gonioreflectometer, a device that measures directionally-dependent reflection from a physical sample [8]. In this work we are primarily concerned with isotropic BRDFs, which are invariant under rotation about the surface normal and under reflection in the incidence plane. Formally a BRDF ρ is isotropic if ρ (θi , φi , θo , φo ) = ρ (θi, φi + ψ , θo , φo + ψ ) for any ψ , and ρ (θi, 0, θo , φo ) = ρ (θi, 0, θo , −φo ). An isotropic BRDF depends only on the relative position of the azimuthal angles, and can therefore be represented by a function of only three real variables θ i , θo , and |φi − φo | ≡ ∆φ . A BRDF that is isotropic, or can be treated as isotropic, is simpler to represent numerically: the number of tabular values or measurements needed is reduced from O(N 4 ) to O(N 3 ), N being the number of values in each angle. Furthermore, a gonioreflectometer that measures only isotropic reflection is mechanically simpler. It is natural to wonder whether there is a useful class of BRDFs that are representable by fewer than three variables. Indeed, mathematical models of BRDFs frequently depend only superficially on all of their arguments, and are better thought of as functions of fewer variables. For example, simple Phong-like BRDFs
2
reduce to functions of only a single variable, the deviation from the direction of mirror reflection. As we shall see, several isotropic BRDF models which initially appear to be functions of three variables readily reduce to only two. A. Goals of this Work The problem of determining the proper arity of a function, the number of arguments that it depends on, is subtle. First, such a simplification must entail more than merely representing a function as a composition of lowerarity functions, which is theoretically always possible. This follows from an intriguing theorem of Kolmogorov [9] which states that essentially all continuous functions of more than one argument can be expressed as compositions of single-argument functions and addition. Indeed, it is a straightforward exercise to express most BRDF models as compositions of addition with single-argument elementary functions such as sin, exp, etc., although doing so confers no particular advantage. In contrast, what we pursue is a representation of isotropic BRDFs having the form
ρ (u, v) = f (ζ1 (u, v), ζ2 (u, v))
(6)
where ζ1 and ζ2 are relatively simple real-valued functions. It is in this sense that we claim certain isotropic BRDFs are intrinsically functions of two variables rather than three. When such a representation is possible, it generally leads to a simpler expression for the BRDF, and also reveals something of its limitations. For example, we will show that any two-variable representation necessarily incurs some approximation error. The idea of choosing a more suitable parameterization for a BRDF has been investigated previously by Marschner [8], Rusinkiewicz [10], and subsequently utilized by Ramamoorthi and Hanrahan [11]. In these previous works the goal was to assist approximation by selecting a parameterization that produced the smallest variation of the BRDF. Kautz and McCool [12] took a different approach, approximating BRDFs with a sum of strictly twovariable BRDF lobes. Our work differs from these in two fundamental aspects. First, our goal is in part to exactly represent the BRDF where possible, using fewer variables. Second, our interest is in how well we can approximate a general isotropic BRDF with only two variables (as verified empirically and visually). In particular, we wish to preserve the Fresnel behavior of increasing specularity near grazing reflection. We use the term “variable” in the context of a BRDF function to refer to a quantity that is a function of the direction vectors u and v, as in ζ1 and ζ2 in (6). Most BRDF models depend on a number of other surface-dependent parameters, such as a Lambertian reflectance, a measure of surface roughness, etc. However, these are constant for any fixed BRDF function in the model, and are therefore not considered in our reduction methods. B. Outline We begin in Section II by examining the BRDF models of Phong and Lafortune, and through these introduce our
TABLE I S UMMARY OF N OTATION AND VARIABLE NAMES n u v r b h
= = = = = =
unit surface normal incident (incoming) direction reflected (scattered, or outgoing) direction reflection of u through n: r = 2(u·n)n − u bisector vector: b = (u + v)/2 normalized bisector, or “halfway” vector: h = b/kbk
θu = angle between u and the surface normal: cos θu = u·n θv = angle between v and the surface normal: cos θv = v·n θh = “halfway” angle; angle between b or √ h and the surface normal: cos θh = b·n/kbk = h·n (= 1 − δ )
ηuv = angle between u and v: cos ηuv = u·v (= 1 − 2α ) ηr = angle between v and the reflection of u, and vice versa: cos ηr = 2(u·n)(v·n) − u·v (= 1 − 2β ) α β δ ψ γ ξ σ
= = = = = = =
sin2 η2uv = 12 (1 − cos ηuv ) sin2 η2r = 12 (1 − cos ηr ) sin2 θh 1−γ −δ cos θu cos θv = 1 − α − β = 1 − δ − ψ (cos θu + cos θv )2 /4 = (1 − δ )(1 − α ) 1 − α − ξ = δ (1 − α )
first set of barycentric variables. Section III contains a similar treatment of the Blinn and Ward BRDF models and presents the second variable set. In Section IV we consider the implications of approximating a general BRDF using an arbitrary two-dimensional parameterization. In Section V we turn to practical approximation of measured BRDF data, and then introduce our third parameterization that overcomes some of the inadequacies of the first two. Finally, Section VI describes several applications of our two-dimensional parameterizations, including new BRDF models. II. P HONG AND L AFORTUNE BRDF S In 1973 Bui-Toc Phong introduced a simple cosine falloff heuristic for adding highlights to computer generated images [13]. Although Phong did not phrase it as such, his model can be cast in the language of BRDFs. The corresponding “Phong BRDF” is written as c ρP (u, v) = d + cs (cos ηr )k (7) π where cd and cs are surface-dependent constants, the diffuse and specular coefficients, respectively, and ηr is the angle v makes with the reflection of u. (Negative values of cos η r are clamped to 0.) The exponent k is a surface-dependent constant that controls the sharpness of the specular peak. In vector notation, cos ηr = v·r, where r ≡ 2(u·n)n − u is the reflection vector of u through the surface normal, and therefore cos ηr = 2(u·n)(v·n) − u·v. (8)
The Phong model is thus expressed in vector notation, but in essence it still depends only on the single variable η r . Any BRDF that depends only on ηr is said to be radially symmetric, or symmetric about the mirror reflection direction. Radially symmetric BRDFs such as Phong’s model have some limitations: they are inherently isotropic, and can exhibit neither off-specular peaks nor Fresnel behavior because the polar angles are not considered.
3
α
1
ng en
tangent retro-reflection
o
tt
=
φ
γ
(∆
ce
=
rfa
90 ◦ )
su
β
e
α
th
B. Barycentric Coordinates The domain of the function ρL (x, y) in (12) is the set of all pairs (x(u, v), y(u, v)) ∈ R2 for u and v on the hemisphere. It can be shown that this domain is the triangle bounded by x, y ≤ 1 and 0 ≤ x + y (Figure 1), which suggests the use of barycentric coordinates. We propose the following: ηuv 1 α (u, v) ≡ sin2 (13) = (1 − cos ηuv ) 2 2 ηr 1 β (u, v) ≡ sin2 = (1 − cos ηr ) (14) 2 2 γ (u, v) ≡ cos θu cos θv = 1 − α − β . (15)
ta
0
A. Lafortune’s Model Let us now consider the BRDF model proposed by Lafortune et al. [14], which hinges on a small but important change to the Phong model. The Phong BRDF can be written succinctly in matrix form: c ρP (u, v) = d + cs [−uT H v]k (9) π where H ≡ I − 2nnT is the Householder reflection matrix based on the surface normal n [15]. Lafortune observed that the model can be generalized by replacing the Householder matrix by a matrix that, in effect, allows for polar angles to be taken into account. In the isotropic case, Lafortune’s BRDF has the matrix form c ρL (u, v) = d + cs [−uT A v]k (10) π where A performs a scaling in the direction of the normal vector n. The matrix A can be expressed as QT SQ, where Q transforms into the local coordinates of the BRDF, and S = diag(a, a, b). (Lafortune’s model actually uses a sum of “lobes” having the form of the second term in (10).) To see the relation to Phong’s model, note that the H matrix can be expressed as QT diag(−1, −1, 1)Q. Lafortune’s model therefore possesses one additional variable that the Phong model does not: namely a/b. In terms of u and v, the Lafortune model becomes c ρL (u, v) = d + cs [(a − b)(u·n)(v·n) − au·v]k , (11) π which shows explicitly that the model can be expressed as a function of two variables c ρL (x, y) = d + cs [(a − b)x − ay]k (12) π under the coordinate transformation x(u, v) = (u·n)(v·n) and y(u, v) = u·v.
v
=
Fig. 1. The values of cos ηuv and cos ηr are restricted to a triangle, which is naturally parameterized in barycentric coordinates (α , β , γ ).
or
v
1
γ
sθ co
β
diffuse reflection (mostly)
u 0,
+
γ
=
u
x = cos ηr
α
grazing reflection
1
γ
sθ co (−1, −1)
α (1, 1)
β = 0, v = r, mirror reflection
y = cos ηuv β
α = 0, u = v; retro-reflection normal reflection
1
β
Fig. 2. Barycentric coordinates α , β , and γ and some geometric interpretations.
For fixed u and v, each of α , β , and γ lies in the unit interval [0, 1]. Also, α + β + γ = 1 so one of the coordinates is redundant—any two are sufficient to parameterize Lafortune’s isotropic BRDF model, as well as any radiallysymmetric BRDF. In terms of the barycentric coordinates, the Phong and Lafortune BRDFs are (omitting the explicit dependence on u and v) c ρP (α , β , γ ) = d + cs (1 − 2β )k (16) π c ρL (α , β , γ ) = d + cs [a(α − β ) − bγ ]k , (17) π respectively. We refer to the parameterization of (13)–(15) as the “αβ ” parameterization. These variables are similar to a subset of a parameterization proposed by Marschner [8] (Appendix II). The mapping of (u, v) to (α , β , γ ) is not injective: in general an infinite collection of (u, v) pairs will map to the same barycentric position. The αβ parameterization is therefore not sufficient to parameterize a general isotropic BRDF. Nevertheless, the barycentric variables α , β , and γ contain some useful geometric information. For example, any (u, v) pair that maps to α = 0 satisfies u = v, and by continuity α near 0 implies u is near v. In other words, the region of the triangle where α is nearly zero corresponds to backscattering. There are a number of such relationships:
α ≈0 β ≈0 α ≈1 β ≈1 γ ≈0 γ ≈1 α ≈β
⇔ ⇔ ⇔ ⇔ ⇔ ⇔ ⇔
near retro-reflection (u ≈ v) near mirror reflection (v ≈ r) near grazing mirror reflection near tangent retro-reflection one of u, v near tangent near normal reflection ∆φ near ±90◦
Figure 2 illustrates this graphically in the αβ -plane. III. B LINN AND WARD BRDF S AND THE B ISECTOR The Phong BRDF is purely heuristic in that there is no physical basis for the cosine falloff. Blinn [16] used
4
a different approach based on the observation that u and v are nearly mirror reflections if, and only if, the bisector vector b = 12 (u + v) is nearly aligned to the surface normal n. Blinn’s BRDF uses the angle b makes with n, which is called the “halfangle” θh , to measure the distance to mirror reflection. Blinn’s BRDF is thus c ρB (u, v) = d + cs cos p θh (18) π where p is the Blinn exponent, and cos θh =
b·n u+v = h·n = ·n. kbk ku + vk
(19)
The vector h = b/kbk, which is the normalized bisector, is often referred to as the “halfway” vector. The behavior of Blinn’s model is close to that of Phong’s near normal incidence, but becomes significantly different toward grazing reflection. In particular, cos θh becomes increasingly sensitive to deviation from the plane of incidence due to the bisector normalization. If θu = θv = θ is fixed, then 1 cos θh ≈ q . (20) 1 + 14 (∆φ − π )2 tan2 θ As θ → 90◦ , tan θ grows without bound and thus the specular peak becomes increasingly narrow in the azimuthal direction, approaching a singularity at θ → 90◦ . This gives a bisector lobe a distinct shape, analogous to the radial symmetry of a Phong lobe. The halfangle θh appears in a number of BRDF models. For example Ward’s model [6] has the isotropic form cd cs − tan2 θh √ (21) ρW = + exp π 4π a2 cos θu cos θv a2 where a is a (constant) measure of surface roughness. Note that (21) is a function of only two values: cos θu cos θv and θh . It is a simple exercise to show that these two values are fundamentally different from the Lafortune variables. The value of cos ηuv , for example, is not determined by a given (ψ , δ ), just as θh is not by (α , β ). Ward’s BRDF suggests an alternate bivariate parameterization. In fact, the variables can be chosen so that their domain is the same unit triangle as that of the αβ parameterization, and thus can be expressed in an analogous set of barycentric coordinates. We use
ψ (u, v) ≡ 1 − γ − δ = cos2 θh − cos θu cos θv δ (u, v) ≡ sin2 θh = 1 − cos2 θh γ (u, v) ≡ cos θu cos θv
(22) (23) (24)
and refer to this as the “ψδ ” parameterization. As functions of u and v, the (ψ , δ , γ ) variables have physical interpretations similar to those of (α , β , γ ) shown in Figure 2. The two parameterizations even share the same complementary barycentric coordinate γ . It can be shown (see Appendix I) that (ψ , δ , γ ) are indeed barycentric coordinates. As with the αβ parameterization, the ψδ parameterization is not sufficient to represent an arbitrary isotropic BRDF, although it does parameterize
Blinn’s model and Ward’s isotropic model exactly. However, any three of α , β , ψ and δ provide a complete (unique) isotropic BRDF parameterization (Appendix I). In the ψδ parameterization the Blinn and Ward BRDF models are c ρB = d + s(1 − δ ) p/2 (25) π c cs δ −1 ρW = d + (26) √ exp 2 π 4π a γ δ a2 respectively. Note that Ward’s BRDF becomes somewhat simpler, as it is expressed without the tangent function. IV. A PPLICATION TO G ENERAL BRDF S Having seen that the well-established Lafortune and Ward analytical BRDF models are inherently bivariate functions, but of different variables, we now consider the implications. Two questions immediately arise. First, what is lost when a two-variable model is used to approximate a general BRDF model or measured data? Second, are other two-variable coordinate systems better suited to BRDF approximation than the αβ and ψδ systems? A. The Projection Concept As noted, the mappings (u, v) 7→ (α , β ) and (u, v) 7→ (ψ , δ ) are not injective: each fixed (α , β ) or (ψ , δ ) generally has an infinite collection of (u, v) vector pairs that map to it. Each mapping thus partitions Ω×Ω into a collection of different equivalence classes: two vector pairs (u 1 , v1 ) and (u2 , v2 ) belong to the same class if, and only if, they map to the same (α , β ) or (ψ , δ ) value. The equivalence classes are fixed for a given coordinate set, but different coordinate sets may have very different equivalence classes. A Phong or Lafortune BRDF is necessarily constant on each αβ equivalence class, as is a Blinn or Ward BRDF on each ψδ class, but a general BRDF function need not be constant on either. However, any two-variable coordinate set can be used to approximate a general isotropic BRDF, by treating the BRDF as constant on each equivalence class. We refer to such an approximation as a projection of a BRDF onto the coordinate set. The idea is not new; in fact, the isotropic assumption is just such a projection. A real surface is not likely to have an entirely isotropic BRDF. However if it has little variation under rotation it can be approximated with an isotropic BRDF, by assuming it is constant on all vector pairs that can be rotated or reflected to each other about the surface normal, i.e., over each isotropic equivalence class. Isotropy can thus be interpreted as a projection of a general four-variate BRDF onto a three-variable coordinate set. In fact, spatial invariance of a BRDF can even be viewed as the projection of a sixvariable function onto four coordinates. Projection onto a two-variable set such as αβ or ψδ is merely a further application of this concept. The constant value of the BRDF for each equivalence class could conceivably be chosen in a variety of ways. We generally use an averaging approximation described in Section V. Whatever the selection method however, the
5
α
ψ
φ = 90◦ φ = 180◦
φ
(θ u , ∆ φ )
(θv , 0)
θminθmax
φ = −90◦
90 ◦
θ
=
(θu , −∆φ )
=
45 ◦
max deviation
θ
β
δ
Fig. 3. A visualization of αβ (left) and ψδ (right) equivalence classes. For each fixed (α , β ) or (ψ , δ ), the equivalence class consists of an infinite collection of vector pairs, represented here as pairs of point on a flat projection of the hemisphere. Paired points have the same color. The points are plotted at uniform values of θv ∈ [θmin , θmax ]. The plots are arranged over the triangle centered at the point for which the equivalence class is illustrated. The expanded plot, for which α = 0.9 and β = 0.05, is more labeled. The dashed circles in the ψδ plots show the value of θ h for the class.
success of the approximation depends on how much the BRDF varies over each equivalence class.
B. A Lower Bound on the Projection Error The diameter of the set of values in each equivalence class is a fundamental quantification of BRDF variation over the classes. For a given equivalence class E ⊂ Ω × Ω, the diameter is formally defined as diamE (ρ ) ≡
sup ((u,v),(u0,v0 ))∈E
|ρ (u, v) − ρ (u0, v0 )|
= sup(ρ ) − inf(ρ ). E
E
1 sup diamE (ρ − ρ ∗ ) 2 all E
OF
G ENERAL BRDF S
In this section we apply the two-variable parameterizations to the practical problem of approximating measured BRDF data. Before doing so, however, we examine the αβ and ψδ equivalence classes, because their structure that ultimately controls the how well the parameterizations can represent real data. Moreover, analyzing the weaknesses leads to a new two-dimensional parameterization.
(27) A. Structure of the Equivalence Classes (28)
The constant value chosen to represent the BRDF on the equivalence class is at best midway between sup(ρ ) and inf(ρ ). Therefore, if ρ ∗ represents a projected BRDF, kρ − ρ ∗ k ∞ ≥
V. P RACTICAL P ROJECTIONS
(29)
where k f k∞ denotes the L∞ norm of f (the essential supremum). Note that (29) is a lower bound on the error incurred in projecting a BRDF, and as such it provides a theoretical limit on the approximation ability of a particular twovariable representation. For example, it implies that even an infinite sum of Lafortune lobes as in (10) cannot converge to an arbitrary isotropic BRDF, because the equivalence classes remain the same for all the lobes.
The equivalence classes, being subsets of the fourdimensional space Ω × Ω, are not easy to visualize directly. One visualization approach is to plot for a fixed (α , β ) or (γ , δ ), pairs of points on the hemisphere that represent pairs of vectors in the class. Figure 3 illustrates such plots for the αβ and ψδ parameterizations. The plots are arranged in a panorama over the triangle to help illustrate how the classes change with the variables. A BRDF projection is necessarily constant on each of its equivalence classes, so the shapes of the classes provides a rough idea of what reflection phenomena the parameterization can and cannot capture. Of particular interest is the structure near grazing reflection, which is around the triangle’s upper vertex. Notice that the αβ classes narrow in θ in this region while the narrowing of the ψδ classes is in φ . This behavior has a very definite impact on the approximation ability of the parameterizations toward grazing.
6
α β
ψ δ
α σ
α = sin2 η2uv σ = kbT k ξ = kbn k
n
h
u bn p ξ
√ α
b
ηuv 2
ηuv 2
√
bT
v
σ
Fig. 4. Lines of constant parameter value on a plane with near grazing reflection. At each point on the plane, the value of the variables are determined by the direction to the black dot and the direction to the eye. Lines of constant α and β (left) have almost identical shape that causes a round highlight. Similarly, ψ and δ lines (middle) produce a highly elongated highlight shape. Lines in the proposed ασ parameterization (right) remain fundamentally different across the plane.
B. Problems with the Parameterizations Some analytic limitations on the parameterizations can be developed. In the ψδ parameterization, γ = cos θ u cos θv and δ = sin2 θh both near zero imply near grazing reflection. However, γ does not depend on ∆φ , and consequently a ψδ projected BRDF depends on ∆φ only in δ . This means that a near-grazing azimuthal slice of a ψδ -projected BRDF will have essentially the same shape as that of a Blinn BRDF. The αβ parameterization has a similar restriction. If θ u = θv are fixed near grazing, it can be shown that to first order cos∆φ =
β −α β ≈ 1−2 . β +α α
(30)
Similarly, if we fix ∆φ = 180◦ , then
β cos ηr ≈ 1 − 2 . (31) α The specular falloff is therefore approximately the same in both the polar and azimuthal directions, and consequently an αβ -projected BRDF retains a Phong-like shape even toward grazing reflection. Figure 4 illustrates the problems at grazing reflection more concretely, by displaying lines of constant parameter value. C. A New Parameterization In summary, both the αβ and ψδ parameterizations have limited approximation ability near grazing reflection. Figure 4 suggests that choosing one variable from each of the αβ and ψδ parameterizations would allow greater flexibility near grazing. Indeed, α and δ produce a useful bivariate parameterization: α = sin2 η2uv measures the separation of u and v, while δ = sin2 θh measures the deviation from mirror reflection. However, the bisector normalization instability is still present, and the parameterization is highly nonuniform, despite the fact that the domain is a square. Yet another drawback arises in measured BRDF data. For each θu and θv , measurements are typically taken at uniformly spaced ∆φ values. In the αβ and ψδ parameterizations these points conveniently lie on straight lines of constant γ = cos θu cos θv , but they lie along curves in an αδ parameterization. To overcome these difficulties, we propose a barycentric variant of an αδ parameterization comprised of α and the
Fig. 5. Geometric interpretations of the ασ parameterization: σ and ξ are the squares of the components of the bisector vector b in the bn-plane, α is the square of the distance from b to u or v.
unnormalized bisector vector. Let ηuv α (u, v) ≡ sin2 2 σ (u, v) ≡ 1 − α − ξ 2 cos θu + cos θv ξ (u, v) ≡ . 2
(32) (33) (34)
Here α = sin η2uv as usual; the variable ξ = (b·n)2 is related to, but not equivalent to cos2 θh , because the bisector vector b need not be a unit vector. Moreover, ξ depends only on θu and θv so that variation in ∆φ occurs along straight lines of constant ξ , a property useful for measured BRDF interpolation. The variables σ and ξ have a geometric interpretation: as ξ is the square of the length of the bisector projected onto the surface normal, it can be shown that σ is that of the bisector projected onto the surface tangent plane (Figure 5). More precisely, if bn = (b·n)n, and bT = b − bn , then
α = ku − bk2 = kv − bk2 σ = kbT k2 ξ = kbn k2 .
(35) (36) (37)
That is, σ and ξ are the squares of the components of b (unnormalized) in the bn-plane, while α gives some measure of how far u and v extend out of this plane. This interpretation is useful, because the bisector vector b provides much information about the reflection geometry. For example, u and v are mirror reflections if, and only if, b lies along the normal (i.e., σ = 0). Similarly, kbn k2 = ξ = 1 implies normal reflectance, and grazing reflection occurs only when both σ and ξ are nearly zero (α = 1). Thus (α , σ , ξ ) have the same approximate geometric interpretations as (α , β , γ ) and (ψ , δ , γ ) illustrated in Figure 2. The new “ασ ” parameterization is, however, fundamentally different from the αβ and ψδ parameterizations. In particular, the reflection angle ηr cannot be determined from α , σ , and ξ , although the halfangle θh can be computed from σ , and ξ from cos2 θh =
ξ . 1−α
(38)
The ασ parameterization can thus capture, but does not impose, the bisector normalization.
7
α
Krylon Blue Data 10 diamE
5
Garnet Red Data 10 diamE
5
Equivalence class panorama for the ασ parameterization.
Figure 6 illustrates the equivalence classes of the ασ parameterization. Note that the equivalence class branches are sometimes disconnected, in such a way that allows an ασ -projected BRDF to better distinguish on and off specular behavior. As we shall see, the improved flexibility near grazing of the ασ parameterization significantly enhances its approximation ability.
D. Approximating Measured Data To test their approximation efficacy, we used our parameterizations to approximate isotropic measured BRDF data provided by Cornell University [17]. The measured data for each surface consists of spectral BRDF measurements taken at sample points in θu , θv , and ∆φ , which we interpolate to obtain BRDF function values. We use the CIE XY Z tristimulus form of the data for simplicity, although the approach immediately generalizes to purely spectral data. Our interpolation method is trilinear in θu , θv , and ∆φ with some extra derivative continuity adjustments in critical places. We interpolate 1/X, 1/Y , and 1/Z rather than X, Y , and Z directly because doing so reduces the interpolation error in undersampled regions. The measured data extend only to θ = 80◦ ; we extrapolate values beyond this by smoothly interpolating to a zero value at 90◦ . The interpolation method is imperfect and causes obvious artifacts in some of the renderings. Smoothing the artifacts might improve the aesthetic appearance of the images, but it also introduces an extra approximation step which can only complicate our evaluation of the parameterizations’ ability to represent the data. Our working assumption is that a good approximation should be able to capture all BRDF behavior, including any artifacts it may have. In Section VI we consider the problem of fitting measured data in a twovariable coordinate system directly.
10
β
10
β
diamE
5
ψ
diamE
α
δ
10
5
α
Fig. 6.
diamE
5
α
σ
10
σ
diamE
5
ψ
δ
α
σ
Fig. 7. Plots of diamE of the measured data. The z-axis is linear, and the values are from the Y component. (The surfaces have been smoothed somewhat for clarity.)
For projection, we use the sampled mean value functional
ρ ∗ (x, y) =
1 N ∑ ρ (θk , 0, θu (θk , x, y), ∆φ (θk , x, y)) N k=1
(39)
with θk uniformly spaced in [θmin , θmax ]. The values of ρ in the summation are interpolated as described above. Formulas for θmin , θmax , and θu (θk , x, y), ∆φ (θk , x, y) are given in Appendix I. Figure 7 illustrates diamE (ρ ), plotted as a height field over the triangle. The ασ parameterization clearly provides the best approximation for this data. E. Rendered Images In computer graphics, the effectiveness of a BRDF approximation ultimately lies in its visual appearance. Figures 8–10 contain ray-traced images rendered using our three bivariate approximations, with global illumination computed by path tracing. The “projection” images are rendered with the approximate BRDF of each parameterization, projected using the mean value functional of (39) with N = 16. The “measured” images are reference images rendered directly from the interpolated measured data. Our first test scene (Figure 8) consists of a single sphere on a textured Lambertian surface, after a figure of Ashikhmin and Shirley [7]. Of all our scenes, this is the most forgiving of poor BRDF approximation. In fact, the images are nearly identical. This is explained by the large light source and limited grazing reflection, although Fresnel behavior is evident around the edges of the sphere. The second scene (Figure 9) is a more demanding test of grazing reflection. Note the αβ parameterization fails to capture the apparent elongation of the highlight, as explained above. Note that only the ασ -projected approximation captures any of the glare-like artifacts. The reflections of the rightmost two sources rely on extrapolated values. Our final test scene (Figure 10) consists of a row of spheres illuminated by a single point source. This arrangement demonstrates reflection in off-specular directions. For example, the central highlights on the illuminated limbs of
8
measured
αβ projection
ψδ projection
ασ projection
Fig. 8. A test scene rendered using measured BRDF data and its projection onto our barycentric coordinate sets. The sphere surface substrate is the Krylon Blue (top) and Garnet Red (bottom). The scene is illuminated by a large spherical source. Some differences in the images can be discerned in the reflection of the white circle and near the edges of the sphere. The differences however, are quite subtle.
measured
αβ projection
ψδ projection
ασ projection
Fig. 9. A test scene illustrating grazing point source reflection off of a halfplane of the Krylon Blue and Garnet Red substrates. Note how the αβ parameterization fails to capture the elongation of the highlight. The ray-like effects emanating from the highlights might pass for glare, but they are actually caused by derivative discontinuities in linear interpolation. The extended highlight on the Garnet Red surface is artificially narrow also due to an interpolation problem. Note that the ασ parameterization is able to capture some of the ray-like artifacts.
spheres illustrate the behavior at near-tangent incidence and near-normal exitance, and vice versa. The αβ parameterization cannot distinguish between this u/v configuration and others away from the highlight along the bright limb of the sphere. Consequently, the highlight is spread out in that direction. The ψδ parameterization has a similar problem, but manifested in the reflection remaining too bright toward the edges. Note that the greatest approximation inaccuracies of both the αβ and ψδ parameterizations correspond to the region in the domain near the middle of the hypotenuse of the triangle.
VI. A PPLICATIONS The advantages of an exact BRDF parameter reduction, i.e., one free of approximation error, have already been noted. In this section we discuss applications of the approximate reductions developed in this paper. A. Visualization of Measured BRDF Data General BRDF functions are difficult to visualize. However, a two-variable function is readily displayed as a surface height field. Similarly, measured data can be displayed by mapping the u and v directions of each data point to
9 measured
measured
αβ projection
αβ projection
ψδ projection
ψδ projection
ασ projection
ασ projection
Fig. 10. Renderings of a row of spheres lit by a point source. Note how the αβ projections spread out the specular highlight over the limbs, and the ψδ projections push brightness too close to the edges. The ασ projection however remains faithful to the measured data.
the bivariate parameterization and plotting the value as a point above the plane. Figure 11 illustrates this for the Krylon Blue and Garnet Red paint data under the ασ parameterization; Figure 12 shows near-grazing closeups of the Krylon Blue data under all three parameterizations. These sample plots provide a direct illustration of the effect of the projection on the data: if the points lie on a twodimensional surface, the parameterization provides an exact representation; the approximation weakens as the samples diverge into a point cloud. B. Fitting Measured Data Practical interpolation of measured BRDF data is difficult as the artifacts in our rendered images attest. Multilinear interpolation generally introduce derivative discontinuities; higher-order interpolation can cause undulations that result in extraneous highlights; smoothing methods, such as Bspline fitting, can improperly alter the shape of the BRDF lobes. All of these are exacerbated by the extreme variation in the measured values, which spans several orders of magnitude even for a moderately specular surface. To make matters worse, the data points can be highly scattered and undersampled in critical regions due to physical limitations
10
10
1
1
0.1
0.1
0.01
0.01
σ
X Y Z Krylon Blue Data
α
α σ Garnet Red Data
Fig. 11. Plots of the Krylon Blue and Garnet Red BRDF data samples. The angular positions of each sample is mapped to ασ coordinates where each CIE XY Z component is plotted as a (logarithmic) height. Note the rapid increase of the data points toward grazing reflection (α = 1).
on the measurement apparatus. When the data can be reduced in dimensionality, the interpolation problem becomes much simpler. Isotropic gonioreflectometer data is typically given as a function of θv , θu , and ∆φ . Interpolating directly in these variables is easiest with respect to the data, but has numerous boundary constraints. An additional advantage to fitting BRDF data in our barycentric parameterizations is that most of the boundary constraints are enforced automatically by
10
1 Y
β
1 Y
α αβ parameterization
δ
ψ ψδ parameterization
1 Y
TABLE II I MPROVED MODEL PARAMETER VALUES
σ
α ασ parameterization
Fig. 12. Closeups of the Krylon Blue samples in the region near grazing reflection under each parameterization. In these plots the height is the inverse of the CIE Y (luminance) component. Under the ασ projection the sample points form a recognizable surface near σ = 0, under the ψδ projection the surface is distinct but more complicated, but under the αβ projection the plotted points vary significantly near α = 1.
the projection. For example, because β corresponds to the cosine of ηr the specular lobe will have a zero derivative in ηuv at ηuv = 0, unless the surface has a singularity or a vertical tangent in β at 0. Similar conditions hold for other variables; all are result of the fact that the variables are smoothly-varying functions of the direction cosines of u and v rather than their spherical coordinates. C. Empirical BRDF Models We leave the problem of fitting general measured BRDF data for future work. As a proof of concept, we present some two-variable BRDF models based on the shape of measured data in the ασ parameterization. We emphasize that these models are purely empirical: they are constructed by approximating the shape of existing measured data without regard to any analytical model or physical interpretation beyond the basic thermodynamic constraints. The two-dimensional BRDF projections all have the same basic characteristics: the surfaces are flattest near the middle of the triangle, drop off toward the hypotenuse, and have a sharp specular upslope near σ = 0 that increases toward grazing reflection (as α or ψ approaches 1). It is our goal to create a model with minimal relative error, particularly near grazing reflection. Accordingly, we fit to 1/ρ ∗ , which approaches zero near σ = 0, rather than to ρ directly. Moreover, we use the ασ parameterization because all the measured data becomes approximately linear near σ = 0 under this parameterization. That is,
ρ ∗ (α , σ ) =
1 + f (α , σ ) m(α )σ + c(α )
(40)
where f is a bounded function that remains small as σ → 0. The αβ and ψδ parameterizations do not exhibit the property of (40) for the Cornell BRDF measurements. D. A Phong-Like Model We propose first a simple model in the spirit of Phong’s BRDF, suggested by (40): SC , (41) C+σ where D is akin to a Lambertian component. The second term is similar to the Phong specular component: S is a specular coefficient and the term is maximal at σ = 0, the
ρ1∗ (α , σ ) = Dλ +
D S0 S1 C0 C1
Krylon Blue X Y 0.0046 0.004 2500 2600 1000 1000 13 12 0.05 0.05
Z 0.023 2100 1000 14 0.02
D S0 S1 C0 C1
Garnet Red X Y 0.0035 0.001 1620 1800 800 800 20 20 0.1 0.1
Z 0.0003 2340 800 20 0.1
point of mirror reflection. The value of C is analogous to the Phong exponent in that it controls the specular falloff. However, the falloff in this model is fundamentally different from Phong’s cosine falloff, particularly toward grazing reflection. This model is similar to the approach Schlick used to approximate the Phong cosine falloff with rational functions [18]. The simple BRDF model of (41) works well away from grazing angles (as do most BRDF models). To match the Krylon Blue data we use D = [0.0062, 0.0057, 0.023], S = 0.05 and C = 0.01; for the Garnet Red, D = [0.005, 0.002, −0.0005], S = 0.045, and C = 0.15. The Z component value of D is negative for the Garnet Red, because unlike the other components, its surface curves away from the plane and this behavior can only be matched with a negative D value. In this sense D is not exactly a Lambertian term. E. An Improved Model The simple model of (41) has some obvious limitations. Most notably, there is not enough specular sharpening near grazing. We propose an improved BRDF model of the following form: 1 ρ2∗ (α , σ ) = Dλ + , (42) Sλ (α )σ +Cλ (α ) where S, and C are smooth, slowly changing functions of α . Unlike the simple model, the C and S are spectrally dependent. We found that interpolating between S and C values fit near α = 0, and near α = 1 produces acceptable results. If these fit values are S0,λ , C0,λ and S1,λ , C1,λ respectively, our S and C functions are Sλ (α ) = (1 − t3 (α ))S0,λ + t3 (α )S1,λ
Cλ (α ) = (1 − t2 (α ))C0,λ + t2 (α )C1,λ
(43) (44)
where t2 and t3 are the blending polynomials t2 (α ) = 1 − (1 − α )2, and t3 = (3 − 2α )α 2, employed to provide better shape approximation. Figure 13 illustrates the fitting properties of the two BRDF models as surface plots on the ασ parameterization. Figures 14 and 15 contain the test images rendered using the simple and the improved model along with reference images. Table II lists the parameter values we used to fit the improved model to the Krylon Blue and Garnet Red data. VII. C ONCLUSION In this paper we have introduced a collection of variables, each a function of the incident and outgoing direction vec-
11
1/ρ1∗
ρ1∗
α
ρ2∗
α
σ σ
α σ
Fig. 13. Surface plots of the simple BRDF model ρ1∗ , and the improved model ρ2∗ with the Krylon Blue data points.
measured
measured
simple model
improved model
simple model
improved model
measured
measured
simple model
simple model
improved model
improved model
Fig. 14. Test scene rendered using measured data, the simple model ρ1∗ model, and the improved model ρ2∗ .
tors u and v, that when taken in pairs are useful in approximately parameterizing isotropic BRDFs. This approach allows us to conveniently express several isotropic BRDF models in terms of two-variables rather than three, and leads to a viable new method of approximately representing measured BRDF data. At the same time, the reduction in the number of variables places theoretical bounds on the attainable accuracy of the approximation. Moreover, we have shown that projecting measured BRDF data onto our proposed two-variable ασ domain, which is based on the coordinates of the unnormalized bisector vector, results in little loss of the information contained in the original datasets. Finally we have demonstrated that a relatively simple two-variable BRDF model can capture nearly all the reflection effects exhibited by a collection of measured data. Avenues we wish to explore further include BRDF visualization based on twodimensional projections, better fitting methods in the twodimensional domains, and the construction of new physically plausible BRDFs.
Fig. 15. Test scene rendered using measured data, the simple model ρ1∗ model, and the improved model ρ2∗ .
A PPENDIX I A DDITIONAL F ORMULAS 1) αβ Parameterization: Suppose α and β are fixed, and
θv is known. Then
cos θu =
The measured BRDF data used in this work was provided courtesy of Cornell University. The authors wish to thank Michael Ashikhmin, Steve Marschner, Peter Shirley, Ken Torrance, and Steve Westin for their insights, encouragement, and many useful suggestions. This work was funded in part by an NSF CARGO grant (DMS-0138440).
cos ∆φ =
β −α . sin θu sin θv
(A.1)
For a fixed α and β , θv (and θu ) lies in the interval [θmin , θmax ], where
θmin =
ACKNOWLEDGMENTS
γ , cos θv
|ηuv − ηr | , 2
θmax =
ηuv + ηr . 2
(A.2)
θv can take on any value in this interval, and thus parameterizes the equivalence class of (α , β ) using the above expressions for θu and ∆φ . 2) ψδ Parameterization: Given u and v, the value of δ = sin2 θh can be computed from its definition in (19), or from cos θh =
cos θu + cos θv cos θu + cos θv √ = . 2 cos η2uv 2 1−α
(A.3)
12
R EFERENCES
The expression
ψ=
∆φ ∆φ 1 cos2 2 (1 − cos2 D) + sin2 2 (1 − cos2 S) 2 cos2 ∆φ (1 + cosD) + sin2 ∆φ (1 + cosS) 2
(A.4)
2
where D = θu − θv and S = θu + θv , implies the barycentric relationship of ψ , δ , and γ , as 0 ≤ ψ ≤ 1 and may take on any value therein. For the inverse parameterization, suppose γ and δ are fixed, as is θv . Then γ = 1 − ψ − δ , and γ cos θu = (A.5) cos θv (cos θu + cos θv )2 − 2(1 + γ )(1 − δ ) cos ∆φ = . (A.6) 2(1 − δ ) sin θu sin θv
θu and θu lie in [θmin , θmax ], where 1 θmin = θh − arccos(1 − 2ψ ) (A.7) 2 1 θmax = θh + arccos(1 − 2ψ ). (A.8) 2 3) ασ Parameterization: Given α , σ , and θv , ξ = 1− α − σ and p cos θu = 2 ξ − cos θv (A.9) 1 − 2α − cos θu cos θv . (A.10) cos ∆φ = sin θu sin θv
The equivalence classes are parameterized by θ in ( |θh − η2uv | if θh + η2uv ≤ π2 , p θmin = (A.11) arccos(2 ξ ) otherwise; ηuv π θmax = min(θh + , ) (A.12) 2 2 p where cos θh = ξ /(1 − α ), and cos ηuv = 1 − 2α as usual. 4) Other Relationships: Any three nonredundant variables of the three coordinate sets we have presented can parameterize any isotropic BRDF. The angles of u and v can be recovered from p p cos θu|v = ξ ± ξ − γ2 (A.13) cos ∆φ
=
q
β −α
1 2 4 (1 + γ ) − ξ
(A.14)
with the collective equation (1 − α )(1 − 2δ ) + β + γ = 2ξ .
R ELATIONSHIP
A PPENDIX II TO E XISTING PARAMETERIZATIONS
Rusinkiewicz [10] proposed the use of four spherical coordinates (θh , φh , θd , φd ) for parameterizing BRDF functions. (θh , φh ) are the spherical coordinates of the (normalized) halfway vector h, (θd , φd ) are the spherical coordinates of u after is rotated about the z-axis by −φh then about the y-axis by −θh (this transformation brings h to the surface normal). In our notation, α = sin2 θd and δ = sin2 θh , but there is no direct correspondence with σ and ξ , because Rusinkiewicz’s parameters depend on the normalized halfway vector h rather than the unnormalized bisector vector b. Marschner [8] proposed the mapping (u, v) 7→ (x, y, z) = (sin θu sin θv cos∆φ , sin θu sin θv sin ∆φ , cos θu cos θv ) to facilitate BRDF resampling. The αβ parameterization is related to this mapping as 2α = 1 − x − z and γ = z, hence the αβ equivalence classes correspond to line segments parallel to the y-axis under Marschner’s mapping.
[1] F. E. Nicodemus, J. C. Richmond, J. J. Hisa, I. W. Ginsberg, and T. Limperis, Geometrical Considerations and Nomenclature for Reflectance, ser. Monograph number 160. Washington DC: National Bureau of Standards, 1977. [2] J. T. Kajiya, “Radiometry and photometry for computer graphics,” in Advanced Topics in Ray Tracing, SIGGRAPH 90 Course Notes, vol. 24, Aug. 1990. [3] K. E. Torrance and E. M. Sparrow, “Theory for off-specular reflection from roughened surfaces,” Journal of the Optical Society of America, vol. 57, no. 9, pp. 1105–1114, Sept. 1967. [4] R. L. Cook and K. E. Torrance, “A reflectance model for computer graphics,” ACM Transactions on Graphics, vol. 1, no. 1, pp. 7–24, Jan. 1982. [5] C. Schlick, “An inexpensive BRDF model for physically-based rendering,” Computer Graphics Forum, vol. 13, no. 3, pp. 233– 246, 1994. [6] G. J. Ward, “Measuring and modeling anisotropic reflection,” Computer Graphics, vol. 26, no. 2, pp. 265–272, July 1992. [7] M. Ashikhmin and P. Shirley, “An anisotropic Phong BRDF model,” J. Graphics Tools, vol. 5, no. 2, pp. 25–32, 2000. [8] S. Marschner, “Inverse rendering for computer graphics,” Ph.D. dissertation, Cornell University, Ithaca, NY, 1998. [9] G. G. Lorentz, Approximation of Functions. New York: Chelsea Publishing Company, 1986. [10] S. Rusinkiewicz, “A new change of variables for efficient BRDF representation,” in Proceedings of the Ninth Eurographics Workshop on Rendering, Vienna, Austria, June 1998. [11] R. Ramamoorthi and P. Hanrahan, “Frequency space environment map rendering,” ACM Transactions on Graphics, vol. 21, no. 3, pp. 517–526, 2002. [12] J. Kautz and M. D. McCool, “Approximation of glossy reflection with prefiltered environment maps,” in Graphics Interface, 2000, pp. 119–126. [13] B. T. Phong, “Illumination for computer generated images,” Ph.D. dissertation, University of Utah, Salt Lake City, UT, 1973. [14] E. P. F. Lafortune, S.-C. Foo, K. E. Torrance, and D. P. Greenberg, “Non-linear approximation of reflectance functions,” in Computer Graphics Proceedings, ser. Annual Conference Series, ACM SIGGRAPH, Aug. 1997, pp. 117–126. [15] J. Arvo, “Analytic methods for simulated light transport,” Ph.D. dissertation, Yale University, 1995. [16] J. F. Blinn, “Models of light reflection for computer synthesized pictures,” in Computer Graphics (SIGGRAPH ’77 Proceedings), J. George, Ed., vol. 11, no. 2, July 1977, pp. 192–198. [17] S. H. Westin, “Measurement data, cornell university program of computer graphics.” [Online]. Available: http://www.graphics.cornell.edu/online/measurements/ [18] C. Schlick, “A fast alternative to phong’s specular model,” in Graphics Gems IV, 1994, pp. 385–387. [19] J. Arvo, K. Torrance, and B. Smits, “A framework for the analysis of error in global illumination algorithms,” in Computer Graphics Proceedings, ser. Annual Conference Series, ACM SIGGAPH, July 1994, pp. 75–84. [20] J. Kautz and M. D. McCool, “Interactive rendering with arbitrary BRDFs using separable approximations,” in Rendering Techniques ’99. Springer, 1999. [21] R. R. Lewis, “Making shaders more physically plausible,” Computer Graphics Forum, vol. 13, no. 2, pp. 109–120, Jan. 1994.
Barycentric Parameterizations for Isotropic BRDFs Michael M. Stark1 , James Arvo1,2 , and Brian Smits2
Abstract— A bidirectional reflectance distribution function (BRDF) is often expressed as a function of four real variables: two spherical coordinates in each of the the “incoming” and “outgoing” directions. However, many BRDFs reduce to functions of fewer variables. For example, isotropic reflection can be represented by a function of three variables. Some BRDF models can be reduced further. In this paper we introduce new sets of coordinates which we use to reduce the dimensionality of several well-known analytic BRDFs as well as empirically measured BRDF data. The proposed coordinate systems are barycentric with respect to a triangular support with a direct physical interpretation. One coordinate set is based on the BRDF model proposed by Lafortune. Another set, based on a model of Ward, is associated with the “halfway” vector common in analytical BRDF formulas. Through these coordinate sets we establish lower bounds on the approximation error inherent in the models on which they are based. We present a third set of coordinates, not based on any analytical model, that performs well in approximating measured data. Finally, our proposed variables suggest novel ways of constructing and visualizing BRDFs. Index Terms— Computer Graphics; Picture/Image Generation; Display algorithms; Color, shading, shadowing, and texture; Reflectance
I. I NTRODUCTION
T
HE faithful simulation of surface reflection is an essential element of accurate image synthesis [1], [2]. The most fundamental relationship governing surface reflection can be stated in terms of radiometric quantities associated with a differential surface area: the differential radiance reflected from a given incident direction to a given outgoing direction is proportional to the irradiance due to a differential solid angle about the incident direction. In symbols, dL(v) ∝ L(u) d ω (1) where L denotes radiance, ω denotes solid angle, and u and v are unit vectors of the incident and reflected directions, respectively. The constant of proportionality depends on both directions, so (1) may be written as dL(v) = ρ (u, v)L(u) d ω .
(2)
The function ρ (u, v) is known as the bidirectional reflectance distribution function, or BRDF, and carries physical units of inverse steradians. This work focuses on the BRDF at a single surface point and at a single wavelength; the concepts are easily generalized to include spatiallyvarying BRDFs and wavelength dependence. Manuscript submitted January 31, 2004 of California, Irvine 2 Pixar Animation Studios 1 University
BRDFs are typically treated as real-valued functions ρ : Ω × Ω → R, where Ω denotes the unit hemisphere. In this representation a BRDF can be expressed as a function of four real variables: ρ (θi, φi , θo , φo ), where θi and φi are polar and azimuthal angles of the incident direction, respectively; θo and φo , of the reflected direction. A BRDF function can be represented more succinctly as a function ρ (u, v) of two (unit) vectors u and v, representing the incoming and outgoing direction vectors on the unit hemisphere. Because a BRDF models the physical process of light reflection, not every function f : Ω × Ω → R represents a physically plausible BRDF. In particular, a BRDF is subject to physical constraints: the second law of thermodynamics requires a BRDF to obey the reciprocity principle
ρ (θi , φi , θo , φo ) = ρ (θo , φo , θi , φi ), ρ (u, v) = ρ (v, u),
(3) (4)
and the first law requires it to conserve energy: Z
Ω
ρ (ω , v) cos θi d ω ≤ 1
(5)
for any fixed v. Reciprocity ensures that a BRDF is symmetric with respect to its incoming and outgoing directions; energy conservation ensures that no more energy is scattered out than comes in. Outside of these constraints, BRDFs can have a wide variety forms. Numerous mathematical BRDF models have found application in a variety of fields, including radiative transfer, remote sensing, computer vision, and computer graphics [3]–[7]. In addition, empirical BRDF measurements are made using a gonioreflectometer, a device that measures directionally-dependent reflection from a physical sample [8]. In this work we are primarily concerned with isotropic BRDFs, which are invariant under rotation about the surface normal and under reflection in the incidence plane. Formally a BRDF ρ is isotropic if ρ (θi , φi , θo , φo ) = ρ (θi, φi + ψ , θo , φo + ψ ) for any ψ , and ρ (θi, 0, θo , φo ) = ρ (θi, 0, θo , −φo ). An isotropic BRDF depends only on the relative position of the azimuthal angles, and can therefore be represented by a function of only three real variables θ i , θo , and |φi − φo | ≡ ∆φ . A BRDF that is isotropic, or can be treated as isotropic, is simpler to represent numerically: the number of tabular values or measurements needed is reduced from O(N 4 ) to O(N 3 ), N being the number of values in each angle. Furthermore, a gonioreflectometer that measures only isotropic reflection is mechanically simpler. It is natural to wonder whether there is a useful class of BRDFs that are representable by fewer than three variables. Indeed, mathematical models of BRDFs frequently depend only superficially on all of their arguments, and are better thought of as functions of fewer variables. For example, simple Phong-like BRDFs
2
reduce to functions of only a single variable, the deviation from the direction of mirror reflection. As we shall see, several isotropic BRDF models which initially appear to be functions of three variables readily reduce to only two. A. Goals of this Work The problem of determining the proper arity of a function, the number of arguments that it depends on, is subtle. First, such a simplification must entail more than merely representing a function as a composition of lowerarity functions, which is theoretically always possible. This follows from an intriguing theorem of Kolmogorov [9] which states that essentially all continuous functions of more than one argument can be expressed as compositions of single-argument functions and addition. Indeed, it is a straightforward exercise to express most BRDF models as compositions of addition with single-argument elementary functions such as sin, exp, etc., although doing so confers no particular advantage. In contrast, what we pursue is a representation of isotropic BRDFs having the form
ρ (u, v) = f (ζ1 (u, v), ζ2 (u, v))
(6)
where ζ1 and ζ2 are relatively simple real-valued functions. It is in this sense that we claim certain isotropic BRDFs are intrinsically functions of two variables rather than three. When such a representation is possible, it generally leads to a simpler expression for the BRDF, and also reveals something of its limitations. For example, we will show that any two-variable representation necessarily incurs some approximation error. The idea of choosing a more suitable parameterization for a BRDF has been investigated previously by Marschner [8], Rusinkiewicz [10], and subsequently utilized by Ramamoorthi and Hanrahan [11]. In these previous works the goal was to assist approximation by selecting a parameterization that produced the smallest variation of the BRDF. Kautz and McCool [12] took a different approach, approximating BRDFs with a sum of strictly twovariable BRDF lobes. Our work differs from these in two fundamental aspects. First, our goal is in part to exactly represent the BRDF where possible, using fewer variables. Second, our interest is in how well we can approximate a general isotropic BRDF with only two variables (as verified empirically and visually). In particular, we wish to preserve the Fresnel behavior of increasing specularity near grazing reflection. We use the term “variable” in the context of a BRDF function to refer to a quantity that is a function of the direction vectors u and v, as in ζ1 and ζ2 in (6). Most BRDF models depend on a number of other surface-dependent parameters, such as a Lambertian reflectance, a measure of surface roughness, etc. However, these are constant for any fixed BRDF function in the model, and are therefore not considered in our reduction methods. B. Outline We begin in Section II by examining the BRDF models of Phong and Lafortune, and through these introduce our
TABLE I S UMMARY OF N OTATION AND VARIABLE NAMES n u v r b h
= = = = = =
unit surface normal incident (incoming) direction reflected (scattered, or outgoing) direction reflection of u through n: r = 2(u·n)n − u bisector vector: b = (u + v)/2 normalized bisector, or “halfway” vector: h = b/kbk
θu = angle between u and the surface normal: cos θu = u·n θv = angle between v and the surface normal: cos θv = v·n θh = “halfway” angle; angle between b or √ h and the surface normal: cos θh = b·n/kbk = h·n (= 1 − δ )
ηuv = angle between u and v: cos ηuv = u·v (= 1 − 2α ) ηr = angle between v and the reflection of u, and vice versa: cos ηr = 2(u·n)(v·n) − u·v (= 1 − 2β ) α β δ ψ γ ξ σ
= = = = = = =
sin2 η2uv = 12 (1 − cos ηuv ) sin2 η2r = 12 (1 − cos ηr ) sin2 θh 1−γ −δ cos θu cos θv = 1 − α − β = 1 − δ − ψ (cos θu + cos θv )2 /4 = (1 − δ )(1 − α ) 1 − α − ξ = δ (1 − α )
first set of barycentric variables. Section III contains a similar treatment of the Blinn and Ward BRDF models and presents the second variable set. In Section IV we consider the implications of approximating a general BRDF using an arbitrary two-dimensional parameterization. In Section V we turn to practical approximation of measured BRDF data, and then introduce our third parameterization that overcomes some of the inadequacies of the first two. Finally, Section VI describes several applications of our two-dimensional parameterizations, including new BRDF models. II. P HONG AND L AFORTUNE BRDF S In 1973 Bui-Toc Phong introduced a simple cosine falloff heuristic for adding highlights to computer generated images [13]. Although Phong did not phrase it as such, his model can be cast in the language of BRDFs. The corresponding “Phong BRDF” is written as c ρP (u, v) = d + cs (cos ηr )k (7) π where cd and cs are surface-dependent constants, the diffuse and specular coefficients, respectively, and ηr is the angle v makes with the reflection of u. (Negative values of cos η r are clamped to 0.) The exponent k is a surface-dependent constant that controls the sharpness of the specular peak. In vector notation, cos ηr = v·r, where r ≡ 2(u·n)n − u is the reflection vector of u through the surface normal, and therefore cos ηr = 2(u·n)(v·n) − u·v. (8)
The Phong model is thus expressed in vector notation, but in essence it still depends only on the single variable η r . Any BRDF that depends only on ηr is said to be radially symmetric, or symmetric about the mirror reflection direction. Radially symmetric BRDFs such as Phong’s model have some limitations: they are inherently isotropic, and can exhibit neither off-specular peaks nor Fresnel behavior because the polar angles are not considered.
3
α
1
ng en
tangent retro-reflection
o
tt
=
φ
γ
(∆
ce
=
rfa
90 ◦ )
su
β
e
α
th
B. Barycentric Coordinates The domain of the function ρL (x, y) in (12) is the set of all pairs (x(u, v), y(u, v)) ∈ R2 for u and v on the hemisphere. It can be shown that this domain is the triangle bounded by x, y ≤ 1 and 0 ≤ x + y (Figure 1), which suggests the use of barycentric coordinates. We propose the following: ηuv 1 α (u, v) ≡ sin2 (13) = (1 − cos ηuv ) 2 2 ηr 1 β (u, v) ≡ sin2 = (1 − cos ηr ) (14) 2 2 γ (u, v) ≡ cos θu cos θv = 1 − α − β . (15)
ta
0
A. Lafortune’s Model Let us now consider the BRDF model proposed by Lafortune et al. [14], which hinges on a small but important change to the Phong model. The Phong BRDF can be written succinctly in matrix form: c ρP (u, v) = d + cs [−uT H v]k (9) π where H ≡ I − 2nnT is the Householder reflection matrix based on the surface normal n [15]. Lafortune observed that the model can be generalized by replacing the Householder matrix by a matrix that, in effect, allows for polar angles to be taken into account. In the isotropic case, Lafortune’s BRDF has the matrix form c ρL (u, v) = d + cs [−uT A v]k (10) π where A performs a scaling in the direction of the normal vector n. The matrix A can be expressed as QT SQ, where Q transforms into the local coordinates of the BRDF, and S = diag(a, a, b). (Lafortune’s model actually uses a sum of “lobes” having the form of the second term in (10).) To see the relation to Phong’s model, note that the H matrix can be expressed as QT diag(−1, −1, 1)Q. Lafortune’s model therefore possesses one additional variable that the Phong model does not: namely a/b. In terms of u and v, the Lafortune model becomes c ρL (u, v) = d + cs [(a − b)(u·n)(v·n) − au·v]k , (11) π which shows explicitly that the model can be expressed as a function of two variables c ρL (x, y) = d + cs [(a − b)x − ay]k (12) π under the coordinate transformation x(u, v) = (u·n)(v·n) and y(u, v) = u·v.
v
=
Fig. 1. The values of cos ηuv and cos ηr are restricted to a triangle, which is naturally parameterized in barycentric coordinates (α , β , γ ).
or
v
1
γ
sθ co
β
diffuse reflection (mostly)
u 0,
+
γ
=
u
x = cos ηr
α
grazing reflection
1
γ
sθ co (−1, −1)
α (1, 1)
β = 0, v = r, mirror reflection
y = cos ηuv β
α = 0, u = v; retro-reflection normal reflection
1
β
Fig. 2. Barycentric coordinates α , β , and γ and some geometric interpretations.
For fixed u and v, each of α , β , and γ lies in the unit interval [0, 1]. Also, α + β + γ = 1 so one of the coordinates is redundant—any two are sufficient to parameterize Lafortune’s isotropic BRDF model, as well as any radiallysymmetric BRDF. In terms of the barycentric coordinates, the Phong and Lafortune BRDFs are (omitting the explicit dependence on u and v) c ρP (α , β , γ ) = d + cs (1 − 2β )k (16) π c ρL (α , β , γ ) = d + cs [a(α − β ) − bγ ]k , (17) π respectively. We refer to the parameterization of (13)–(15) as the “αβ ” parameterization. These variables are similar to a subset of a parameterization proposed by Marschner [8] (Appendix II). The mapping of (u, v) to (α , β , γ ) is not injective: in general an infinite collection of (u, v) pairs will map to the same barycentric position. The αβ parameterization is therefore not sufficient to parameterize a general isotropic BRDF. Nevertheless, the barycentric variables α , β , and γ contain some useful geometric information. For example, any (u, v) pair that maps to α = 0 satisfies u = v, and by continuity α near 0 implies u is near v. In other words, the region of the triangle where α is nearly zero corresponds to backscattering. There are a number of such relationships:
α ≈0 β ≈0 α ≈1 β ≈1 γ ≈0 γ ≈1 α ≈β
⇔ ⇔ ⇔ ⇔ ⇔ ⇔ ⇔
near retro-reflection (u ≈ v) near mirror reflection (v ≈ r) near grazing mirror reflection near tangent retro-reflection one of u, v near tangent near normal reflection ∆φ near ±90◦
Figure 2 illustrates this graphically in the αβ -plane. III. B LINN AND WARD BRDF S AND THE B ISECTOR The Phong BRDF is purely heuristic in that there is no physical basis for the cosine falloff. Blinn [16] used
4
a different approach based on the observation that u and v are nearly mirror reflections if, and only if, the bisector vector b = 12 (u + v) is nearly aligned to the surface normal n. Blinn’s BRDF uses the angle b makes with n, which is called the “halfangle” θh , to measure the distance to mirror reflection. Blinn’s BRDF is thus c ρB (u, v) = d + cs cos p θh (18) π where p is the Blinn exponent, and cos θh =
b·n u+v = h·n = ·n. kbk ku + vk
(19)
The vector h = b/kbk, which is the normalized bisector, is often referred to as the “halfway” vector. The behavior of Blinn’s model is close to that of Phong’s near normal incidence, but becomes significantly different toward grazing reflection. In particular, cos θh becomes increasingly sensitive to deviation from the plane of incidence due to the bisector normalization. If θu = θv = θ is fixed, then 1 cos θh ≈ q . (20) 1 + 14 (∆φ − π )2 tan2 θ As θ → 90◦ , tan θ grows without bound and thus the specular peak becomes increasingly narrow in the azimuthal direction, approaching a singularity at θ → 90◦ . This gives a bisector lobe a distinct shape, analogous to the radial symmetry of a Phong lobe. The halfangle θh appears in a number of BRDF models. For example Ward’s model [6] has the isotropic form cd cs − tan2 θh √ (21) ρW = + exp π 4π a2 cos θu cos θv a2 where a is a (constant) measure of surface roughness. Note that (21) is a function of only two values: cos θu cos θv and θh . It is a simple exercise to show that these two values are fundamentally different from the Lafortune variables. The value of cos ηuv , for example, is not determined by a given (ψ , δ ), just as θh is not by (α , β ). Ward’s BRDF suggests an alternate bivariate parameterization. In fact, the variables can be chosen so that their domain is the same unit triangle as that of the αβ parameterization, and thus can be expressed in an analogous set of barycentric coordinates. We use
ψ (u, v) ≡ 1 − γ − δ = cos2 θh − cos θu cos θv δ (u, v) ≡ sin2 θh = 1 − cos2 θh γ (u, v) ≡ cos θu cos θv
(22) (23) (24)
and refer to this as the “ψδ ” parameterization. As functions of u and v, the (ψ , δ , γ ) variables have physical interpretations similar to those of (α , β , γ ) shown in Figure 2. The two parameterizations even share the same complementary barycentric coordinate γ . It can be shown (see Appendix I) that (ψ , δ , γ ) are indeed barycentric coordinates. As with the αβ parameterization, the ψδ parameterization is not sufficient to represent an arbitrary isotropic BRDF, although it does parameterize
Blinn’s model and Ward’s isotropic model exactly. However, any three of α , β , ψ and δ provide a complete (unique) isotropic BRDF parameterization (Appendix I). In the ψδ parameterization the Blinn and Ward BRDF models are c ρB = d + s(1 − δ ) p/2 (25) π c cs δ −1 ρW = d + (26) √ exp 2 π 4π a γ δ a2 respectively. Note that Ward’s BRDF becomes somewhat simpler, as it is expressed without the tangent function. IV. A PPLICATION TO G ENERAL BRDF S Having seen that the well-established Lafortune and Ward analytical BRDF models are inherently bivariate functions, but of different variables, we now consider the implications. Two questions immediately arise. First, what is lost when a two-variable model is used to approximate a general BRDF model or measured data? Second, are other two-variable coordinate systems better suited to BRDF approximation than the αβ and ψδ systems? A. The Projection Concept As noted, the mappings (u, v) 7→ (α , β ) and (u, v) 7→ (ψ , δ ) are not injective: each fixed (α , β ) or (ψ , δ ) generally has an infinite collection of (u, v) vector pairs that map to it. Each mapping thus partitions Ω×Ω into a collection of different equivalence classes: two vector pairs (u 1 , v1 ) and (u2 , v2 ) belong to the same class if, and only if, they map to the same (α , β ) or (ψ , δ ) value. The equivalence classes are fixed for a given coordinate set, but different coordinate sets may have very different equivalence classes. A Phong or Lafortune BRDF is necessarily constant on each αβ equivalence class, as is a Blinn or Ward BRDF on each ψδ class, but a general BRDF function need not be constant on either. However, any two-variable coordinate set can be used to approximate a general isotropic BRDF, by treating the BRDF as constant on each equivalence class. We refer to such an approximation as a projection of a BRDF onto the coordinate set. The idea is not new; in fact, the isotropic assumption is just such a projection. A real surface is not likely to have an entirely isotropic BRDF. However if it has little variation under rotation it can be approximated with an isotropic BRDF, by assuming it is constant on all vector pairs that can be rotated or reflected to each other about the surface normal, i.e., over each isotropic equivalence class. Isotropy can thus be interpreted as a projection of a general four-variate BRDF onto a three-variable coordinate set. In fact, spatial invariance of a BRDF can even be viewed as the projection of a sixvariable function onto four coordinates. Projection onto a two-variable set such as αβ or ψδ is merely a further application of this concept. The constant value of the BRDF for each equivalence class could conceivably be chosen in a variety of ways. We generally use an averaging approximation described in Section V. Whatever the selection method however, the
5
α
ψ
φ = 90◦ φ = 180◦
φ
(θ u , ∆ φ )
(θv , 0)
θminθmax
φ = −90◦
90 ◦
θ
=
(θu , −∆φ )
=
45 ◦
max deviation
θ
β
δ
Fig. 3. A visualization of αβ (left) and ψδ (right) equivalence classes. For each fixed (α , β ) or (ψ , δ ), the equivalence class consists of an infinite collection of vector pairs, represented here as pairs of point on a flat projection of the hemisphere. Paired points have the same color. The points are plotted at uniform values of θv ∈ [θmin , θmax ]. The plots are arranged over the triangle centered at the point for which the equivalence class is illustrated. The expanded plot, for which α = 0.9 and β = 0.05, is more labeled. The dashed circles in the ψδ plots show the value of θ h for the class.
success of the approximation depends on how much the BRDF varies over each equivalence class.
B. A Lower Bound on the Projection Error The diameter of the set of values in each equivalence class is a fundamental quantification of BRDF variation over the classes. For a given equivalence class E ⊂ Ω × Ω, the diameter is formally defined as diamE (ρ ) ≡
sup ((u,v),(u0,v0 ))∈E
|ρ (u, v) − ρ (u0, v0 )|
= sup(ρ ) − inf(ρ ). E
E
1 sup diamE (ρ − ρ ∗ ) 2 all E
OF
G ENERAL BRDF S
In this section we apply the two-variable parameterizations to the practical problem of approximating measured BRDF data. Before doing so, however, we examine the αβ and ψδ equivalence classes, because their structure that ultimately controls the how well the parameterizations can represent real data. Moreover, analyzing the weaknesses leads to a new two-dimensional parameterization.
(27) A. Structure of the Equivalence Classes (28)
The constant value chosen to represent the BRDF on the equivalence class is at best midway between sup(ρ ) and inf(ρ ). Therefore, if ρ ∗ represents a projected BRDF, kρ − ρ ∗ k ∞ ≥
V. P RACTICAL P ROJECTIONS
(29)
where k f k∞ denotes the L∞ norm of f (the essential supremum). Note that (29) is a lower bound on the error incurred in projecting a BRDF, and as such it provides a theoretical limit on the approximation ability of a particular twovariable representation. For example, it implies that even an infinite sum of Lafortune lobes as in (10) cannot converge to an arbitrary isotropic BRDF, because the equivalence classes remain the same for all the lobes.
The equivalence classes, being subsets of the fourdimensional space Ω × Ω, are not easy to visualize directly. One visualization approach is to plot for a fixed (α , β ) or (γ , δ ), pairs of points on the hemisphere that represent pairs of vectors in the class. Figure 3 illustrates such plots for the αβ and ψδ parameterizations. The plots are arranged in a panorama over the triangle to help illustrate how the classes change with the variables. A BRDF projection is necessarily constant on each of its equivalence classes, so the shapes of the classes provides a rough idea of what reflection phenomena the parameterization can and cannot capture. Of particular interest is the structure near grazing reflection, which is around the triangle’s upper vertex. Notice that the αβ classes narrow in θ in this region while the narrowing of the ψδ classes is in φ . This behavior has a very definite impact on the approximation ability of the parameterizations toward grazing.
6
α β
ψ δ
α σ
α = sin2 η2uv σ = kbT k ξ = kbn k
n
h
u bn p ξ
√ α
b
ηuv 2
ηuv 2
√
bT
v
σ
Fig. 4. Lines of constant parameter value on a plane with near grazing reflection. At each point on the plane, the value of the variables are determined by the direction to the black dot and the direction to the eye. Lines of constant α and β (left) have almost identical shape that causes a round highlight. Similarly, ψ and δ lines (middle) produce a highly elongated highlight shape. Lines in the proposed ασ parameterization (right) remain fundamentally different across the plane.
B. Problems with the Parameterizations Some analytic limitations on the parameterizations can be developed. In the ψδ parameterization, γ = cos θ u cos θv and δ = sin2 θh both near zero imply near grazing reflection. However, γ does not depend on ∆φ , and consequently a ψδ projected BRDF depends on ∆φ only in δ . This means that a near-grazing azimuthal slice of a ψδ -projected BRDF will have essentially the same shape as that of a Blinn BRDF. The αβ parameterization has a similar restriction. If θ u = θv are fixed near grazing, it can be shown that to first order cos∆φ =
β −α β ≈ 1−2 . β +α α
(30)
Similarly, if we fix ∆φ = 180◦ , then
β cos ηr ≈ 1 − 2 . (31) α The specular falloff is therefore approximately the same in both the polar and azimuthal directions, and consequently an αβ -projected BRDF retains a Phong-like shape even toward grazing reflection. Figure 4 illustrates the problems at grazing reflection more concretely, by displaying lines of constant parameter value. C. A New Parameterization In summary, both the αβ and ψδ parameterizations have limited approximation ability near grazing reflection. Figure 4 suggests that choosing one variable from each of the αβ and ψδ parameterizations would allow greater flexibility near grazing. Indeed, α and δ produce a useful bivariate parameterization: α = sin2 η2uv measures the separation of u and v, while δ = sin2 θh measures the deviation from mirror reflection. However, the bisector normalization instability is still present, and the parameterization is highly nonuniform, despite the fact that the domain is a square. Yet another drawback arises in measured BRDF data. For each θu and θv , measurements are typically taken at uniformly spaced ∆φ values. In the αβ and ψδ parameterizations these points conveniently lie on straight lines of constant γ = cos θu cos θv , but they lie along curves in an αδ parameterization. To overcome these difficulties, we propose a barycentric variant of an αδ parameterization comprised of α and the
Fig. 5. Geometric interpretations of the ασ parameterization: σ and ξ are the squares of the components of the bisector vector b in the bn-plane, α is the square of the distance from b to u or v.
unnormalized bisector vector. Let ηuv α (u, v) ≡ sin2 2 σ (u, v) ≡ 1 − α − ξ 2 cos θu + cos θv ξ (u, v) ≡ . 2
(32) (33) (34)
Here α = sin η2uv as usual; the variable ξ = (b·n)2 is related to, but not equivalent to cos2 θh , because the bisector vector b need not be a unit vector. Moreover, ξ depends only on θu and θv so that variation in ∆φ occurs along straight lines of constant ξ , a property useful for measured BRDF interpolation. The variables σ and ξ have a geometric interpretation: as ξ is the square of the length of the bisector projected onto the surface normal, it can be shown that σ is that of the bisector projected onto the surface tangent plane (Figure 5). More precisely, if bn = (b·n)n, and bT = b − bn , then
α = ku − bk2 = kv − bk2 σ = kbT k2 ξ = kbn k2 .
(35) (36) (37)
That is, σ and ξ are the squares of the components of b (unnormalized) in the bn-plane, while α gives some measure of how far u and v extend out of this plane. This interpretation is useful, because the bisector vector b provides much information about the reflection geometry. For example, u and v are mirror reflections if, and only if, b lies along the normal (i.e., σ = 0). Similarly, kbn k2 = ξ = 1 implies normal reflectance, and grazing reflection occurs only when both σ and ξ are nearly zero (α = 1). Thus (α , σ , ξ ) have the same approximate geometric interpretations as (α , β , γ ) and (ψ , δ , γ ) illustrated in Figure 2. The new “ασ ” parameterization is, however, fundamentally different from the αβ and ψδ parameterizations. In particular, the reflection angle ηr cannot be determined from α , σ , and ξ , although the halfangle θh can be computed from σ , and ξ from cos2 θh =
ξ . 1−α
(38)
The ασ parameterization can thus capture, but does not impose, the bisector normalization.
7
α
Krylon Blue Data 10 diamE
5
Garnet Red Data 10 diamE
5
Equivalence class panorama for the ασ parameterization.
Figure 6 illustrates the equivalence classes of the ασ parameterization. Note that the equivalence class branches are sometimes disconnected, in such a way that allows an ασ -projected BRDF to better distinguish on and off specular behavior. As we shall see, the improved flexibility near grazing of the ασ parameterization significantly enhances its approximation ability.
D. Approximating Measured Data To test their approximation efficacy, we used our parameterizations to approximate isotropic measured BRDF data provided by Cornell University [17]. The measured data for each surface consists of spectral BRDF measurements taken at sample points in θu , θv , and ∆φ , which we interpolate to obtain BRDF function values. We use the CIE XY Z tristimulus form of the data for simplicity, although the approach immediately generalizes to purely spectral data. Our interpolation method is trilinear in θu , θv , and ∆φ with some extra derivative continuity adjustments in critical places. We interpolate 1/X, 1/Y , and 1/Z rather than X, Y , and Z directly because doing so reduces the interpolation error in undersampled regions. The measured data extend only to θ = 80◦ ; we extrapolate values beyond this by smoothly interpolating to a zero value at 90◦ . The interpolation method is imperfect and causes obvious artifacts in some of the renderings. Smoothing the artifacts might improve the aesthetic appearance of the images, but it also introduces an extra approximation step which can only complicate our evaluation of the parameterizations’ ability to represent the data. Our working assumption is that a good approximation should be able to capture all BRDF behavior, including any artifacts it may have. In Section VI we consider the problem of fitting measured data in a twovariable coordinate system directly.
10
β
10
β
diamE
5
ψ
diamE
α
δ
10
5
α
Fig. 6.
diamE
5
α
σ
10
σ
diamE
5
ψ
δ
α
σ
Fig. 7. Plots of diamE of the measured data. The z-axis is linear, and the values are from the Y component. (The surfaces have been smoothed somewhat for clarity.)
For projection, we use the sampled mean value functional
ρ ∗ (x, y) =
1 N ∑ ρ (θk , 0, θu (θk , x, y), ∆φ (θk , x, y)) N k=1
(39)
with θk uniformly spaced in [θmin , θmax ]. The values of ρ in the summation are interpolated as described above. Formulas for θmin , θmax , and θu (θk , x, y), ∆φ (θk , x, y) are given in Appendix I. Figure 7 illustrates diamE (ρ ), plotted as a height field over the triangle. The ασ parameterization clearly provides the best approximation for this data. E. Rendered Images In computer graphics, the effectiveness of a BRDF approximation ultimately lies in its visual appearance. Figures 8–10 contain ray-traced images rendered using our three bivariate approximations, with global illumination computed by path tracing. The “projection” images are rendered with the approximate BRDF of each parameterization, projected using the mean value functional of (39) with N = 16. The “measured” images are reference images rendered directly from the interpolated measured data. Our first test scene (Figure 8) consists of a single sphere on a textured Lambertian surface, after a figure of Ashikhmin and Shirley [7]. Of all our scenes, this is the most forgiving of poor BRDF approximation. In fact, the images are nearly identical. This is explained by the large light source and limited grazing reflection, although Fresnel behavior is evident around the edges of the sphere. The second scene (Figure 9) is a more demanding test of grazing reflection. Note the αβ parameterization fails to capture the apparent elongation of the highlight, as explained above. Note that only the ασ -projected approximation captures any of the glare-like artifacts. The reflections of the rightmost two sources rely on extrapolated values. Our final test scene (Figure 10) consists of a row of spheres illuminated by a single point source. This arrangement demonstrates reflection in off-specular directions. For example, the central highlights on the illuminated limbs of
8
measured
αβ projection
ψδ projection
ασ projection
Fig. 8. A test scene rendered using measured BRDF data and its projection onto our barycentric coordinate sets. The sphere surface substrate is the Krylon Blue (top) and Garnet Red (bottom). The scene is illuminated by a large spherical source. Some differences in the images can be discerned in the reflection of the white circle and near the edges of the sphere. The differences however, are quite subtle.
measured
αβ projection
ψδ projection
ασ projection
Fig. 9. A test scene illustrating grazing point source reflection off of a halfplane of the Krylon Blue and Garnet Red substrates. Note how the αβ parameterization fails to capture the elongation of the highlight. The ray-like effects emanating from the highlights might pass for glare, but they are actually caused by derivative discontinuities in linear interpolation. The extended highlight on the Garnet Red surface is artificially narrow also due to an interpolation problem. Note that the ασ parameterization is able to capture some of the ray-like artifacts.
spheres illustrate the behavior at near-tangent incidence and near-normal exitance, and vice versa. The αβ parameterization cannot distinguish between this u/v configuration and others away from the highlight along the bright limb of the sphere. Consequently, the highlight is spread out in that direction. The ψδ parameterization has a similar problem, but manifested in the reflection remaining too bright toward the edges. Note that the greatest approximation inaccuracies of both the αβ and ψδ parameterizations correspond to the region in the domain near the middle of the hypotenuse of the triangle.
VI. A PPLICATIONS The advantages of an exact BRDF parameter reduction, i.e., one free of approximation error, have already been noted. In this section we discuss applications of the approximate reductions developed in this paper. A. Visualization of Measured BRDF Data General BRDF functions are difficult to visualize. However, a two-variable function is readily displayed as a surface height field. Similarly, measured data can be displayed by mapping the u and v directions of each data point to
9 measured
measured
αβ projection
αβ projection
ψδ projection
ψδ projection
ασ projection
ασ projection
Fig. 10. Renderings of a row of spheres lit by a point source. Note how the αβ projections spread out the specular highlight over the limbs, and the ψδ projections push brightness too close to the edges. The ασ projection however remains faithful to the measured data.
the bivariate parameterization and plotting the value as a point above the plane. Figure 11 illustrates this for the Krylon Blue and Garnet Red paint data under the ασ parameterization; Figure 12 shows near-grazing closeups of the Krylon Blue data under all three parameterizations. These sample plots provide a direct illustration of the effect of the projection on the data: if the points lie on a twodimensional surface, the parameterization provides an exact representation; the approximation weakens as the samples diverge into a point cloud. B. Fitting Measured Data Practical interpolation of measured BRDF data is difficult as the artifacts in our rendered images attest. Multilinear interpolation generally introduce derivative discontinuities; higher-order interpolation can cause undulations that result in extraneous highlights; smoothing methods, such as Bspline fitting, can improperly alter the shape of the BRDF lobes. All of these are exacerbated by the extreme variation in the measured values, which spans several orders of magnitude even for a moderately specular surface. To make matters worse, the data points can be highly scattered and undersampled in critical regions due to physical limitations
10
10
1
1
0.1
0.1
0.01
0.01
σ
X Y Z Krylon Blue Data
α
α σ Garnet Red Data
Fig. 11. Plots of the Krylon Blue and Garnet Red BRDF data samples. The angular positions of each sample is mapped to ασ coordinates where each CIE XY Z component is plotted as a (logarithmic) height. Note the rapid increase of the data points toward grazing reflection (α = 1).
on the measurement apparatus. When the data can be reduced in dimensionality, the interpolation problem becomes much simpler. Isotropic gonioreflectometer data is typically given as a function of θv , θu , and ∆φ . Interpolating directly in these variables is easiest with respect to the data, but has numerous boundary constraints. An additional advantage to fitting BRDF data in our barycentric parameterizations is that most of the boundary constraints are enforced automatically by
10
1 Y
β
1 Y
α αβ parameterization
δ
ψ ψδ parameterization
1 Y
TABLE II I MPROVED MODEL PARAMETER VALUES
σ
α ασ parameterization
Fig. 12. Closeups of the Krylon Blue samples in the region near grazing reflection under each parameterization. In these plots the height is the inverse of the CIE Y (luminance) component. Under the ασ projection the sample points form a recognizable surface near σ = 0, under the ψδ projection the surface is distinct but more complicated, but under the αβ projection the plotted points vary significantly near α = 1.
the projection. For example, because β corresponds to the cosine of ηr the specular lobe will have a zero derivative in ηuv at ηuv = 0, unless the surface has a singularity or a vertical tangent in β at 0. Similar conditions hold for other variables; all are result of the fact that the variables are smoothly-varying functions of the direction cosines of u and v rather than their spherical coordinates. C. Empirical BRDF Models We leave the problem of fitting general measured BRDF data for future work. As a proof of concept, we present some two-variable BRDF models based on the shape of measured data in the ασ parameterization. We emphasize that these models are purely empirical: they are constructed by approximating the shape of existing measured data without regard to any analytical model or physical interpretation beyond the basic thermodynamic constraints. The two-dimensional BRDF projections all have the same basic characteristics: the surfaces are flattest near the middle of the triangle, drop off toward the hypotenuse, and have a sharp specular upslope near σ = 0 that increases toward grazing reflection (as α or ψ approaches 1). It is our goal to create a model with minimal relative error, particularly near grazing reflection. Accordingly, we fit to 1/ρ ∗ , which approaches zero near σ = 0, rather than to ρ directly. Moreover, we use the ασ parameterization because all the measured data becomes approximately linear near σ = 0 under this parameterization. That is,
ρ ∗ (α , σ ) =
1 + f (α , σ ) m(α )σ + c(α )
(40)
where f is a bounded function that remains small as σ → 0. The αβ and ψδ parameterizations do not exhibit the property of (40) for the Cornell BRDF measurements. D. A Phong-Like Model We propose first a simple model in the spirit of Phong’s BRDF, suggested by (40): SC , (41) C+σ where D is akin to a Lambertian component. The second term is similar to the Phong specular component: S is a specular coefficient and the term is maximal at σ = 0, the
ρ1∗ (α , σ ) = Dλ +
D S0 S1 C0 C1
Krylon Blue X Y 0.0046 0.004 2500 2600 1000 1000 13 12 0.05 0.05
Z 0.023 2100 1000 14 0.02
D S0 S1 C0 C1
Garnet Red X Y 0.0035 0.001 1620 1800 800 800 20 20 0.1 0.1
Z 0.0003 2340 800 20 0.1
point of mirror reflection. The value of C is analogous to the Phong exponent in that it controls the specular falloff. However, the falloff in this model is fundamentally different from Phong’s cosine falloff, particularly toward grazing reflection. This model is similar to the approach Schlick used to approximate the Phong cosine falloff with rational functions [18]. The simple BRDF model of (41) works well away from grazing angles (as do most BRDF models). To match the Krylon Blue data we use D = [0.0062, 0.0057, 0.023], S = 0.05 and C = 0.01; for the Garnet Red, D = [0.005, 0.002, −0.0005], S = 0.045, and C = 0.15. The Z component value of D is negative for the Garnet Red, because unlike the other components, its surface curves away from the plane and this behavior can only be matched with a negative D value. In this sense D is not exactly a Lambertian term. E. An Improved Model The simple model of (41) has some obvious limitations. Most notably, there is not enough specular sharpening near grazing. We propose an improved BRDF model of the following form: 1 ρ2∗ (α , σ ) = Dλ + , (42) Sλ (α )σ +Cλ (α ) where S, and C are smooth, slowly changing functions of α . Unlike the simple model, the C and S are spectrally dependent. We found that interpolating between S and C values fit near α = 0, and near α = 1 produces acceptable results. If these fit values are S0,λ , C0,λ and S1,λ , C1,λ respectively, our S and C functions are Sλ (α ) = (1 − t3 (α ))S0,λ + t3 (α )S1,λ
Cλ (α ) = (1 − t2 (α ))C0,λ + t2 (α )C1,λ
(43) (44)
where t2 and t3 are the blending polynomials t2 (α ) = 1 − (1 − α )2, and t3 = (3 − 2α )α 2, employed to provide better shape approximation. Figure 13 illustrates the fitting properties of the two BRDF models as surface plots on the ασ parameterization. Figures 14 and 15 contain the test images rendered using the simple and the improved model along with reference images. Table II lists the parameter values we used to fit the improved model to the Krylon Blue and Garnet Red data. VII. C ONCLUSION In this paper we have introduced a collection of variables, each a function of the incident and outgoing direction vec-
11
1/ρ1∗
ρ1∗
α
ρ2∗
α
σ σ
α σ
Fig. 13. Surface plots of the simple BRDF model ρ1∗ , and the improved model ρ2∗ with the Krylon Blue data points.
measured
measured
simple model
improved model
simple model
improved model
measured
measured
simple model
simple model
improved model
improved model
Fig. 14. Test scene rendered using measured data, the simple model ρ1∗ model, and the improved model ρ2∗ .
tors u and v, that when taken in pairs are useful in approximately parameterizing isotropic BRDFs. This approach allows us to conveniently express several isotropic BRDF models in terms of two-variables rather than three, and leads to a viable new method of approximately representing measured BRDF data. At the same time, the reduction in the number of variables places theoretical bounds on the attainable accuracy of the approximation. Moreover, we have shown that projecting measured BRDF data onto our proposed two-variable ασ domain, which is based on the coordinates of the unnormalized bisector vector, results in little loss of the information contained in the original datasets. Finally we have demonstrated that a relatively simple two-variable BRDF model can capture nearly all the reflection effects exhibited by a collection of measured data. Avenues we wish to explore further include BRDF visualization based on twodimensional projections, better fitting methods in the twodimensional domains, and the construction of new physically plausible BRDFs.
Fig. 15. Test scene rendered using measured data, the simple model ρ1∗ model, and the improved model ρ2∗ .
A PPENDIX I A DDITIONAL F ORMULAS 1) αβ Parameterization: Suppose α and β are fixed, and
θv is known. Then
cos θu =
The measured BRDF data used in this work was provided courtesy of Cornell University. The authors wish to thank Michael Ashikhmin, Steve Marschner, Peter Shirley, Ken Torrance, and Steve Westin for their insights, encouragement, and many useful suggestions. This work was funded in part by an NSF CARGO grant (DMS-0138440).
cos ∆φ =
β −α . sin θu sin θv
(A.1)
For a fixed α and β , θv (and θu ) lies in the interval [θmin , θmax ], where
θmin =
ACKNOWLEDGMENTS
γ , cos θv
|ηuv − ηr | , 2
θmax =
ηuv + ηr . 2
(A.2)
θv can take on any value in this interval, and thus parameterizes the equivalence class of (α , β ) using the above expressions for θu and ∆φ . 2) ψδ Parameterization: Given u and v, the value of δ = sin2 θh can be computed from its definition in (19), or from cos θh =
cos θu + cos θv cos θu + cos θv √ = . 2 cos η2uv 2 1−α
(A.3)
12
R EFERENCES
The expression
ψ=
∆φ ∆φ 1 cos2 2 (1 − cos2 D) + sin2 2 (1 − cos2 S) 2 cos2 ∆φ (1 + cosD) + sin2 ∆φ (1 + cosS) 2
(A.4)
2
where D = θu − θv and S = θu + θv , implies the barycentric relationship of ψ , δ , and γ , as 0 ≤ ψ ≤ 1 and may take on any value therein. For the inverse parameterization, suppose γ and δ are fixed, as is θv . Then γ = 1 − ψ − δ , and γ cos θu = (A.5) cos θv (cos θu + cos θv )2 − 2(1 + γ )(1 − δ ) cos ∆φ = . (A.6) 2(1 − δ ) sin θu sin θv
θu and θu lie in [θmin , θmax ], where 1 θmin = θh − arccos(1 − 2ψ ) (A.7) 2 1 θmax = θh + arccos(1 − 2ψ ). (A.8) 2 3) ασ Parameterization: Given α , σ , and θv , ξ = 1− α − σ and p cos θu = 2 ξ − cos θv (A.9) 1 − 2α − cos θu cos θv . (A.10) cos ∆φ = sin θu sin θv
The equivalence classes are parameterized by θ in ( |θh − η2uv | if θh + η2uv ≤ π2 , p θmin = (A.11) arccos(2 ξ ) otherwise; ηuv π θmax = min(θh + , ) (A.12) 2 2 p where cos θh = ξ /(1 − α ), and cos ηuv = 1 − 2α as usual. 4) Other Relationships: Any three nonredundant variables of the three coordinate sets we have presented can parameterize any isotropic BRDF. The angles of u and v can be recovered from p p cos θu|v = ξ ± ξ − γ2 (A.13) cos ∆φ
=
q
β −α
1 2 4 (1 + γ ) − ξ
(A.14)
with the collective equation (1 − α )(1 − 2δ ) + β + γ = 2ξ .
R ELATIONSHIP
A PPENDIX II TO E XISTING PARAMETERIZATIONS
Rusinkiewicz [10] proposed the use of four spherical coordinates (θh , φh , θd , φd ) for parameterizing BRDF functions. (θh , φh ) are the spherical coordinates of the (normalized) halfway vector h, (θd , φd ) are the spherical coordinates of u after is rotated about the z-axis by −φh then about the y-axis by −θh (this transformation brings h to the surface normal). In our notation, α = sin2 θd and δ = sin2 θh , but there is no direct correspondence with σ and ξ , because Rusinkiewicz’s parameters depend on the normalized halfway vector h rather than the unnormalized bisector vector b. Marschner [8] proposed the mapping (u, v) 7→ (x, y, z) = (sin θu sin θv cos∆φ , sin θu sin θv sin ∆φ , cos θu cos θv ) to facilitate BRDF resampling. The αβ parameterization is related to this mapping as 2α = 1 − x − z and γ = z, hence the αβ equivalence classes correspond to line segments parallel to the y-axis under Marschner’s mapping.
[1] F. E. Nicodemus, J. C. Richmond, J. J. Hisa, I. W. Ginsberg, and T. Limperis, Geometrical Considerations and Nomenclature for Reflectance, ser. Monograph number 160. Washington DC: National Bureau of Standards, 1977. [2] J. T. Kajiya, “Radiometry and photometry for computer graphics,” in Advanced Topics in Ray Tracing, SIGGRAPH 90 Course Notes, vol. 24, Aug. 1990. [3] K. E. Torrance and E. M. Sparrow, “Theory for off-specular reflection from roughened surfaces,” Journal of the Optical Society of America, vol. 57, no. 9, pp. 1105–1114, Sept. 1967. [4] R. L. Cook and K. E. Torrance, “A reflectance model for computer graphics,” ACM Transactions on Graphics, vol. 1, no. 1, pp. 7–24, Jan. 1982. [5] C. Schlick, “An inexpensive BRDF model for physically-based rendering,” Computer Graphics Forum, vol. 13, no. 3, pp. 233– 246, 1994. [6] G. J. Ward, “Measuring and modeling anisotropic reflection,” Computer Graphics, vol. 26, no. 2, pp. 265–272, July 1992. [7] M. Ashikhmin and P. Shirley, “An anisotropic Phong BRDF model,” J. Graphics Tools, vol. 5, no. 2, pp. 25–32, 2000. [8] S. Marschner, “Inverse rendering for computer graphics,” Ph.D. dissertation, Cornell University, Ithaca, NY, 1998. [9] G. G. Lorentz, Approximation of Functions. New York: Chelsea Publishing Company, 1986. [10] S. Rusinkiewicz, “A new change of variables for efficient BRDF representation,” in Proceedings of the Ninth Eurographics Workshop on Rendering, Vienna, Austria, June 1998. [11] R. Ramamoorthi and P. Hanrahan, “Frequency space environment map rendering,” ACM Transactions on Graphics, vol. 21, no. 3, pp. 517–526, 2002. [12] J. Kautz and M. D. McCool, “Approximation of glossy reflection with prefiltered environment maps,” in Graphics Interface, 2000, pp. 119–126. [13] B. T. Phong, “Illumination for computer generated images,” Ph.D. dissertation, University of Utah, Salt Lake City, UT, 1973. [14] E. P. F. Lafortune, S.-C. Foo, K. E. Torrance, and D. P. Greenberg, “Non-linear approximation of reflectance functions,” in Computer Graphics Proceedings, ser. Annual Conference Series, ACM SIGGRAPH, Aug. 1997, pp. 117–126. [15] J. Arvo, “Analytic methods for simulated light transport,” Ph.D. dissertation, Yale University, 1995. [16] J. F. Blinn, “Models of light reflection for computer synthesized pictures,” in Computer Graphics (SIGGRAPH ’77 Proceedings), J. George, Ed., vol. 11, no. 2, July 1977, pp. 192–198. [17] S. H. Westin, “Measurement data, cornell university program of computer graphics.” [Online]. Available: http://www.graphics.cornell.edu/online/measurements/ [18] C. Schlick, “A fast alternative to phong’s specular model,” in Graphics Gems IV, 1994, pp. 385–387. [19] J. Arvo, K. Torrance, and B. Smits, “A framework for the analysis of error in global illumination algorithms,” in Computer Graphics Proceedings, ser. Annual Conference Series, ACM SIGGAPH, July 1994, pp. 75–84. [20] J. Kautz and M. D. McCool, “Interactive rendering with arbitrary BRDFs using separable approximations,” in Rendering Techniques ’99. Springer, 1999. [21] R. R. Lewis, “Making shaders more physically plausible,” Computer Graphics Forum, vol. 13, no. 2, pp. 109–120, Jan. 1994.
