Torrance-sparrow off-specular reflection and linear subspaces for object recognition

Abstract

The Torrance-Sparrow model of off-specular reflection is recast in a significantly simpler and more transparent form in order to render a spherical-harmonic decomposition more feasible. By assuming that a physical surface consists of small, reflecting facets whose surface normals satisfy a normal distribution, the model captures the off-specular enhancement of the reflected intensity distribution often observed at large angles of incidence and reflection, features beyond the reach of the phenomenological broadening models usually employed. In passing we remove a physical inconsistency in the original treatment, restoring reciprocity and correcting the dependence of reflectance on angle near grazing incidence. It is noted that the results predicted by the model are relatively insensitive to values of its one parameter, the width of the distribution of surface normals.

Description

BACKGROUND OF THE INVENTION

[0001] 1. Field of Invention

[0002] The present invention relates generally to machine vision, and more particularly, to a method of modeling the off-specular reflection of light to allow a machine to recognize the image of an object illuminated by varying unknown lighting conditions.

[0003] 2. Description of Related Art

[0004] It is clear that changes in lighting can have a large effect on how an object, such as a person, looks. Understanding this is essential to building effective object recognition systems.

[0005] While the Lambertian (or uniform) model of the scattering of light from real surfaces is central to any study of object illumination, more sophisticated modeling is essential. In general, a variety of conditions combine to introduce anomalies into the simple model, especially those related to specular reflection (strong scattering over a narrow range of angles) and geometrically induced shadows. Recently, consideration has been given to the examination of the harmonic components of broadened, specular reflection, and the harmonics of cast shadows.

[0006] U.S. patent application Ser. No. 09/967,206 filed Sep. 28, 2001, entitled “Broadened-Specular Reflection and Linear Subspaces for Object Recognition”, the disclosure of which is hereby incorporated by reference, considered harmonic components of broadened, specular reflection in the Phong model.

[0007] While the Phong model is mathematically contrived, the Torrance-Sparrow model upon which the present invention is based, is a physically motivated model.

[0008] In addition, while the Phong model could be put into a reasonably tractable form by means of standard harmonic analysis, this is an unreasonable goal for the Torrance-Sparrow model. In accordance with the teachings of the present invention, a modified harmonic analysis was used which still enabled an approximation much more readily in the state space of spherical harmonics than in pixel space as graphics people usually do. In the previous invention, the higher harmonics were discarded. In the present invention the higher harmonics are kept in a modified (approximate) form, more approximate the less important the harmonics are.

[0009] Most importantly, the present invention makes use of simpler functions where arguments capture more of the angular dependencies than is possible with the usual expansions.

[0010] U.S. patent application Ser. No. 09/705,507, filed Nov. 3, 2000, entitled “Lambertian Reflectance and Linear Subspaces”, the disclosure of which is hereby incorporated by reference, considered the relationship between the function that describes the lighting intensity, and the reflectance function that describes how much light an object reflects as a function of its surface normal, under a given lighting condition. Representing these functions as spherical harmonics, it was shown that for Lambertian reflectance the mapping from lighting to reflectance is a convolution with a nearly perfect low-pass filter. High-frequency components of the lighting hardly affect the reflectance function. Therefore, nearly all Lambertian reflectance functions could be modeled as some linear combination of nine spherical harmonic components.

[0011] These low-dimensional approximations have two key advantages. First, they are described with few parameters, thus optimization is simpler. Moreover, the images produced under low-frequency lighting are inherently immune to small changes in orientation. A small rotation in the surface normal is equivalent to a small rotation in the light. However, low-frequency light does not appreciably change with a small rotation. The analysis shows this, in that an image produced by an nth order harmonic light component is computed by taking an nth degree polynomial of the components of the surface normal.

[0012] While promising, the invention described in U.S. patent application Ser. No. 09/705,507 considered only Lambertian reflectance. Specular reflectance, described in U.S. patent application Ser. No. 09/967,206, can significantly affect how an object appears. For example, a human face may have some oiliness, particularly on the forehead, which reflects light specularly. This is seen as highlights on the object.

[0013] While the invention described in U.S. Ser. No. 09/705,507 involved only Lambertian reflection, the Phong model is a simple mathematical rendering of diffuse specular reflectance. The Torrance-Sparrow model underlying the present invention is a fundamental improvement to make the reflectance physically accurate, a necessity for face recognition.

[0014] Specular reflectance cannot be approximated well by only considering low-dimensional linear subspaces which capture only low-frequency lighting. Intuitively, the set of possible reflectance functions for a perfect mirror is the same as the set of possible lighting functions. The question remains how to cope with specularity in object recognition using models that will inevitably be imperfect.

[0015] Specularities are often thought of as mainly present in spatially localized highlights, which occur when a bright light source is reflected by a shiny object. Such highlights can potentially be identified in the image and processed separately. For example, the algorithms described by Coleman, et al., in “Obtaining 3-Dimensional Shape of Textured and Specular Surfaces Using Four-Source Photometry”, Physics-Based Vision Principles and Practice, Shape Recovery, pp. 180-199 (1992), and by Brelstaff, et al., in “Detecting Specular Reflection Using Lambertian Constraints”, International Conference on Computer Vision, pp. 297-302 (1988), treat specular pixels as statistical outliers that should be discarded. However, specularities due to low-frequency lighting can be spread out throughout the entire image. Since these specularities affect every pixel in the image, their effect must be better understood; they cannot simply be discarded.

[0016] Some previous works have considered specular effects also as convolutions. Basri, et al. in NEC Research Institute Technical Report No. 2000-112, entitled “Lambertian Reflectance and Linear Subspaces” (August 2000) as revised in Technical Report No. 2000-172R (December 2000) considered a version of the Phong model in which the specular lobe is a rotationally symmetric shape centered about a peak that is a 180-degree rotation of the direction of the incoming light about the surface normal. In such a model, specular reflection doubles the frequency of the incoming light, and then acts as a convolution. However, this model is heuristic, and not physically based. It does not extend to realistic models.

[0017] Ramamoorthi, et al., in “On the Relationship between Radiance and Irradiance: Determining the Illumination From Images of a Convex Lambertian Object”, Journal of the Optical Society of America A, Vol. 18, No. 10, pp. 2448-2459 (2001), considered the effects of arbitrary bi-directional reflectance distribution functions (BRDFs) on the light field. They were able to develop algorithms for recovering the BRDF and/or the lighting given a measured light field and an object with a known 3D structure. Though interesting, these results are not generally applicable to machine vision because, without special sensing devices, the light field is not available. An image of an object is produced by the reflectance function that describes how the surface normals reflect light in one fixed direction. Therefore, an image corresponds to a 2D slice through the light field. For an arbitrary BRDF, the 2D reflectance function for a fixed viewpoint no longer results from a convolution of the lighting.

BRIEF SUMMARY OF THE INVENTION

[0018] The goal of the present invention is to allow a machine, such as a computer, to match an input image of an object, such as a human face, to an object model stored in a database associated with the machine, and to perform this matching with a high degree of accuracy.

[0019] The steps of the recognition method of the present invention include providing a database of three-dimensional object models, providing an input image, positioning each three dimensional model relative to the input image, determining an optimal rendered image of the model that is most similar to the input image, computing a measure of similarity between the input image and the optimal rendered image of each model, and selecting the object whose optimal rendered image is most similar to the input image. The step of determining an optimal rendered image includes deriving a reflectance function for each object that includes the effects of both broadened-specular and diffuse reflection, and optimizing that reflectance function to most closely match the input image.

[0020] The described method of the present invention includes rendering the images produced by an object model when illuminated only by the fundamental along with each harmonic component of light to create the effect that there are no negative regions, or some other combinations which still cover all the harmonics of interest. The summation of these components defines the reflectance function for the object model. The derived reflectance function takes into account the intensity of light components incident on the object model, the intensity of light diffusely reflected from the model, and the intensity of off-light specularly reflected from the model, the last component being dependent on the angle between the direction of the observer and the direction of maximum specular reflection.

[0021] The advance in the art provided by the present invention resides in the fact that unlike the Lambertian case in which the reflectance depends only on the surface normal, and only harmonics of the surface normal are considered, the reflectance expressed by any model including off-specular reflection involves two coupled dependencies, that on the surface normal and that on the direction of the reflected light. The number of terms in such an expansion grows only as n2 for the Lambertian, but as n4 for diffuse, specular reflection, and in addition n is larger for the latter case. The term “n” is the harmonic which is sometimes referred to as “l” in the text below. This complexity is avoided by using knowledge of physically correct, specular reflection to preserve both dependencies in a single harmonic expansion by using harmonics of a modified direction dependent on both the surface normal and the reflected direction of interest, returning the complexity to n2. This result is somewhat analogous to having a spatial Fourier series in a plane (2D). For arbitrary surfaces one has products of Fourier components in the x and the y directions. However, if one knows the structure is roughly circular, or colinear, or some other symmetry, a 1D modified expansion may suffice, not simply x or y, but in a combination appropriate to the problem.

[0022] A principal object of the present invention is therefore, the provision of a method of modeling the off-specular reflection of light from an object to allow a machine to recognize the image of an object illuminated by varying lighting conditions.

[0023] Another object of the present invention is the provision of a method of modeling the off-specular reflection of light from an object using a modified Torrance-Sparrow model of the object to allow a machine to recognize the image of an object illuminated by varying lighting conditions.

[0024] A further object of the present invention is the provision of a method of modeling the off-specular reflection of light from an object using a modified Torrance-Sparrow model of the object to allow a machine to recognize the image of an object illuminated by varying lighting conditions for use in face recognition.

[0025] Further and still other objects of the present invention will become more clearly apparent when reading the following description of the invention

DETAILED DESCRIPTION OF THE INVENTION

[0026] Although detailed and/or microscopic models of reflection exist, these apply principally to homogeneous, rigid, material surfaces. For reflection from inhomogeneous, flexible surfaces such as faces, where most of the features of interest are configurational, a more robust model that captures the most important aspects of the scattered light is more appropriate. Thus the Lambertian model, which simply scatters incident radiation equally into the hemisphere above each surface element, is usually chosen.

[0027] The method of the present invention involves a much more complicated, physical model of reflectance developed by Torrance and Sparrow. The Torrance-Sparrow model of off-specular reflection is recast in a significantly simpler and more transparent form in order to render a spherical-harmonic decomposition more feasible. By assuming that a physical surface consists of small, reflecting facets whose surface normals satisfy a normal distribution, the model captures the off-specular enhancement of the reflected intensity distribution often observed at large angles of incidence and reflection, features beyond the reach of the phenomenological broadening models usually employed. The results predicted by the new model are relatively insensitive to values of its one parameter, the width of the distribution of surface normals. This is accomplished using both physical intuition and mathematical simplification. The result is an extraction of the principal dependences in the harmonic representation with a minimum of dependence on the details of the model.

[0028] Torrance-Sparrow Reflection Model (TSRM)

[0029] Torrance and Sparrow presented a physical model which captured several seemingly anomalous features of broadened, specular reflection from surfaces of metals and insulators. These features included off-specular maxima in the reflected intensity distribution and comparatively larger relative reflectance at larger angles of incidence. Neither of these features could be included in phenomenological broadening models such as the simpler Phong model which lacked both.

[0030] The model may be described as follows: the surface is taken to comprise a distribution of small, randomly-oriented, mirror-like facets, and the effects of shadowing and masking of facets by adjacent facets are included. Such roughened surfaces are the usual result of fabricating non-single-crystal samples. The off-specular maxima and enhanced reflectance result from the larger number of facets available for scattering at large angles of incidence. For example, in the limit of grazing incidence, the facet normal can vary by π radians in an arc orthogonal to both the plane of incidence and the nominal surface leaving the incident and reflected directions invariant.

[0031] In the following simplified derivation, there is a direct finding of Blinn's compact expressions for the shadowing/masking factor G and a slight altering of the facet-orientation distribution function. A detailed description of Blinn's work is found in the article J. F. Blinn, “Models of Light Reflection for Computer Synthesized Pictures”, Computer Graphics 11 (2), pp. 192-198 (1977). The masking factor G is described in detailed below. The present model removes a physical inconsistency in the original treatment, one which has been carried forward verbatim in nearly all subsequent discussions. This is significant since in its present form the TSRM cannot account even qualitatively for the experimentally-observed, off-specular reflection at high angles of incidence. A primary concern in the present invention is to represent the reflection in terms of spherical harmonics of the incoming illumination, analogous to prior treatments.

[0032] The Shadowing-Masking Factor G

[0033] As one would expect of scattering from a rough, faceted surface, a good deal of local masking and shadowing would occur between the facets. While some of this contributes to the diffuse (Lambertian) scattering observed, most is light scattered only once. The TSRM treats this by assuming that each facet of normal nγ is paired with a complementary facet of normal n′γ(n′γ=2(n·nγ)n−nγ) forming a V-shaped, trough-like local structure. In general, only a fraction (G) of the light falling on a trough can reach a reflecting surface, be reflected, and exit the trough without being shadowed on entering or masked on exiting. This purely geometric factor is G:0≦G≦1, and is independent of the depth (or size) of the trough.

[0034] Torrance and Sparrow provide a complicated, trigonometric derivation whose result was greatly compacted by Blinn, and is given in below. We provide a simple, direct derivation which includes other useful forms for G as well.

[0035] Consider a ray of light ne reflected off the surface of a facet nγ at a distance of x (along the facet) from the bottom of the trough, and exiting the trough by just grazing the upper edge of the adjoining facet of normal n′γ (at a distance y (along the facet) from the bottom of the trough). Since any light reflected in direction ne forms an interior angle of β with facet nγ and α with facet n′γ,sinβ=ne·nγ and sinα=−ne·n′γ. But by the law of sines x/y=sinα/sinβ: hence x/y=−ne·n′γ/ne·nγ, and therefore

G b =(ne·nγ+ne·n′γ)/ne·nγ

[0036] But n′γ=2(n·nγ)n−nγ(n′γ is just the reflection of nγ about n, the nominal surface normal). Hence, Gb=2(n·nγ)(n·ne)/nγ·ne. Gb=2(n·nγ)(n·n′e)/nγ·n′e is found similarly, or just by reciprocity, replacing the outgoing ne with the incoming n′e. We note in passing that (nγ+n′γ)=(2+2nγ·n′γ)1/2n, rendering

G b =(2+2nγ·n′γ)1/2n·ne/nγ·ne

[0037] which is yet another form for Gb. Strictly speaking, ne includes a component parallel to the lateral (nγ×n′γ) axis of the trough, and although neither nγ nor n′γ has a corresponding component, this component does enter the normalization of ne. However, as only the ratio of the two dot products enters, use of the full ne and n′e is justified.

[0038] Broadened, Specular Reflection in TSRM

[0039] Given a distribution of facets f(nγ), where nγ is the facet normal, the specular reflection rTS of incident lumination L(ni) can be written as

r TS (ne,n)=u(n·ne)∫dniL(ni)∫dnγf(nγ)G(n,nγ,ne)δ(ni−n′e)

[0040] where ni is the (negative) direction of the incident light, ne is the direction of interest of the reflected light, n is the normal of the nominal surface in the region of interest, and n′e=(2nγ(nγ·ne)−ne) is the incident spectral direction corresponding to ne for a facet whose normal is nγ. The term u(x) is the unit step function: u(x)=0 for x<0, u(x)=1 for x>0. The term u(n·ne) ensures that only regions of the surface capable of reflecting light into ne are considered. (The corresponding criterion for incident light, u(n·ni), is included in the shadowing/masking factor G.) The Fresnel factor included by Torrance and Sparrow has been omitted. All vectors nk are normalized.

[0041] Assume a facet distribution of normals of the form

f (nγ)=exp(γn·nγ)/N

N= 2π(eγ−1)/γ

[0042] The normalization N is such that the integral of f(nγ) over the hemisphere where n·nγ>0 is unity. This differs only slightly in form from the Gaussian distribution in the angle α (n·nγ=cosα≈1−α2/2) used in the TSRM. This change, though minor numerically, offers significant algebraic simplification. (The heart of our distribution is therefore exp(−γα2/2) while that of TSRM is exp (−c2α2). The TSRM c=0.05/degree works out to c2=8.207/radian2, corresponding to a γ of 16.414. In the present invention γ=5, 10 and 20 were chosen for the numerical work.)

[0043] As discussed above, the shadowing/masking factor G is given by 1$\begin{array}{c}G\left(n,n_{\gamma },n_e\right)=\max \left(0,\min \left(1,G_b,G_c\right)\right)\\ G_b=2n·n_{\gamma }\left(\frac{n·n_e}{n_{\gamma }·n_e}\right),G_c=2n·n_{\gamma }\left(\frac{n·n_e^{\prime }}{n_y·n_e^{\prime }}\right)\end{array}$

[0044] Clearly 0≦G≦1: G=1 when the reflection is such that there is neither masking nor shadowing of the incident or reflected light; G=Gb when some of the reflected light is masked, which can only occur when n·ne<n·n′e; G=Gc when some of the incident light is shadowed, which can only happen when n·n′e<n·{circumflex over ( )}ne:Gb<Gc and Gc<Gb, respectively. (Since nγ is the normal of the specularly reflecting facet, nγ·ne=nγ·n′e.) (It should be noted that whereas most treatments of TSRM consider the reflected distribution in ne for fixed ni, our primary interest is the reflected light in direction ne for a given distribution in ni. This means that the most strongly reflected light will come from ni where n·ni<n·ne. By reciprocity, this concurs with the usual case where light reflected from a fixed ni to some ne is maximized for n·ni>n·ne.)

[0045] The expression for rTS can be readily cast into a spherical-harmonic form: First recall 2$\delta \left(n_i-n_e^{\prime }\right)={\Sigma }_{lm}{Y}_{lm}\left(n_e^{\prime }\right)Y_{lm}^{*}\left(n_i\right)$

[0046] It is also advantageous to change the variable of the remaining integration from nγ to ne. As in Torrance and Sparrow, this is done to express the integrand in terms of the solid angle of the reflected (TS) or the incident (here) light. Torrence and Sparrow refer to previous work of W. A. Rense in the Journal of Optical Society of America, 40, pp 55ff (1950) for this conversion. However, it is readily possible to derive the result by using n′e=2nγ(nγ·ne)−ne as follows: we can change nγ to nγ+dnγ in only two dimensions, either perpendicular to nγ and ne, dnγ⊥, or in the plane of nγ and ne (the plane of reflection) but perpendicular to nγ,dnγ∥. (As nγ is of unit length, being a normal vector, nγ·dnγ=0.) Thus dn′e⊥=2dnγ⊥(nγ·ne)+2nγ(dnγ⊥·ne)=2(nγ·ne)dnγ⊥ as dnγ⊥·ne=0. Similarly dn′e∥=2dnγ∥. This is more readily seen by noting that in the plane of reflection the change in the angle between incident and reflected rays is double the change between incident and normal. Thus dnγ=dn′e/4(nγ·ne), and rTS takes the form

r TS (ne,n)=u(ne·n)ΣlmLlmRlm(ne,n)

R lm (ne,n)=∫dn′eYlm(n′e)G(n,nγ,ne)eγn nγ4N(nγ·ne)

L lm =∫dn i L (ni)Y*lm(ni)

[0047] Here of course, nγ satisfies n′e=(2nγ(nγ·ne)−ne).

[0048] Missing in the above expression is the (nγ·n′e)/(n·ne) factor present in Torrance and Sparrow, and in Blinn. The factor arises because in passing from intensities to fluxes and then back to intensities, Torrance and Sparrow perform the first conversion as if the reflection were due to individual facets, which is correct, but the second conversion is performed as if the reflection were due to the average surface, which is incorrect. The correct factor is, therefore, (nγ·n′e)/(nγ·ne) which, of course, is unity. The present invention avoids this unnecessary complication and the result satisfies reciprocity while TSRM's does not.

[0049] Spherical-Harmonic Factor Rlm (ne,n)

[0050] For any reasonable γ (e.g. γ>10), the [G exp(γn·nγ)/(nγ·ne)]=F factor reaches a maximum for n·nγ≈1, and the maximum is comparatively sharp, at least for the relatively low values of lm of interest. The change of variables from nγ to n′e ensures that the domain of n′e, dn′e, is well-contained even for relatively small (nγ·ne): recall dn′e=4(nγ·ne)dnγ. Also, it is the diminishing value of the (nγ·ne) factor which leads to larger reflectance for larger angles of incidence: for smaller (nγ·ne), a larger fraction of the facets contribute to the reflection.

[0051] Now, first calculate the value nγc of nγ which maximizes F. The result is that there are two regions in (n·ne) to be considered: n·ne>n·neg where the maximum occurs for G=1; n·ne<n·neg where the maximum occurs for Gc=1 and nγc=nγg. The critical value n·neg depends only on γ; the larger γ, the larger n·neg. From nγc we can obtain n′ec=2nγc(nγc·ne)−ne and take Ylm(n′ec) out of the integral in Rlm(ne,n). Then the remaining integral can be approximated as a two-dimension Gaussian integral and evaluated.

[0052] First calculate the facet normal nγc which maximizes F for G=1. One method is to differentiate with respect to n·nγ, the component of nγ parallel to n (from which nγ can be determined, since by symmetry nγc lies in the (n,ne) plane). There is an easier way, which, moreover, will be of value when evaluating the integral over n′e. Suppose m is a vector in the (n,ne) plane. Let η be the angle between m and nγ. Then m·nγ=cos η. If now one changes nγ to nγ+dnγ∥, this corresponds to rotating nγ in the (n,ne) plane. Since d(m·nγ)/dη=−sin η=−(1−cos2η)1/2=−(1−(m·nγ)2)1/2, and since nγ has unit length (so that dη=dnγ∥), it follows that d(m·nγ)/dnγ∥=−(1−(m·nγ)2)1/2. Hence, ∂(n·nγ)/∂nγ∥=−(1−(n·nγ)2)1/2, ∂(nγ·ne)/∂nγ∥=−(1−(nγ·ne)2)1/2.

[0053] Varying F with respect to nγ∥ shows that nγc satisfies

γ(1−(n·nγ)2)1/2(nγ·ne)=(1−(nγ·ne)2)1/2

[0054] and

γ(1−n·nγ2)1/2(n·nγn·ne−(1−n·nγ2)1/2(1−n·ne2)1/2)=n·nγ(1−n·ne2)1/2+n·ne(1−n·nγ2)1/2

[0055] (The notation has been simplified so that n·nγ2=(n·nγ)2 and n·ne2=(n·ne)2.) Solving for (1−n·nγ2)1/2 the result is: 3${\left(1-n·n_{\gamma c}^2\right)}^{1/2}=\frac{n·n_e\left(n·n_{\gamma }-{\gamma }^{-1}\right)-{\left(n·{n_e^2\left(n·n_{\gamma }-{\gamma }^{-1}\right)}^2-4\left(1-n·n_e^2\right)n·n_{\gamma }/\gamma \right)}^{1/2}}{2{\left(1-n·n_e^2\right)}^{1/2}}$

[0056] By taking n·nγ=1 on the right-hand side of the equation, and knowing n·ne and γ, a very accurate value for n·nγc is obtained on the left-hand side of the equation, which can be improved by iteration if necessary. From n·nγ we have nγ·ne=n·nγn·ne−(1−n·nγ2)1/2(1−n·ne2)1/2, and hence F for nγ=nγc and G=1. (When solving the quadratic equation, the minus sign is chosen in order to obtain the correct limit as n·ne→1.)

[0057] Now it is necessary to ascertain for which value n·neg of n·ne Gc=1. For n·ne<n·neg, G<1 and F will decrease more rapidly than any corresponding increase in the other factors. Hence, Gc=1 fixes n·nγc for n·ne<n·neg. Thus setting Gc=2n·nγn·n′e/nγ·n′e=1 results in

n γ ·n e =2n·nγ(2n·nγnγ·ne−n·ne),

[0058] from which it follows that

(1−n·nγg2)1/2=n·nγn·ne(1−2(4n·nγ2−1)−1)(1−n·ne2)−1/2

[0059] which as above can be solved for n·nγg in a few iterations. Here nγg is nγ for Gc=1 given ne. Equating the expression for nγc and nγg, results in an expression for neg

n·n eg =(1+γ(2n·nγ(4n·nγ2−1)−1−γ−1)(1−2(4n·nγ2−1)−1))−1/2

[0060] where n·nγ=n·nγg for n·ne=n·neg. Since n·nγg is nearly unity,

n·n eg ≈(1+(2γ/3−1)/3)−1/2

[0061] These expressions enable us to evaluate Fmax for all (n,ne) pairs. Values for Fm/F0 are given in Table 1, where Fm=Fmax and F0 is F at the angle for specular reflection (n·nγ=1) . In Table 1, θe is the angle of reflection of interest relative to the nominal-surface normal (n·ne=cosθe), θi is the incident angle which maximizes the reflectance into θe(n·n′e(nγe)=cosθi), and Fm/F0 is the ratio of the reflectance from light incident at θi to that incident at θe, i.e., nominal specular reflection. The smaller γ, the greater the fluctuation permitted in the facet normals, and hence the greater angle as well as magnitude of off-specular scattering which results. The horizontal lines in the table distinguish regions where the maximum reflection occur for Gc>1 (smaller angles) and Gc=1 (larger angles). For the latter, the location of the maxima is independent of γ, the key parameter in Gaussian distributions. Unlike TSRM, due to masking and shadowing the reflection goes rapidly to zero as the incident angle exceeds θi. Also, the magnitude of the enhancement is significantly reduced. The former is expected from physics, and is absent in the TSRM. In the TSRM it results from incorrectly calculating the scattered flux. The latter presumable results from not including the appropriate Frensnel factor.

	TABLE 1
	γ = 5		γ = 10		γ = 20
θe	θi	Fm/F0	θi	Fm/F0	θi	Fm/F0
00.00	00.00	1.000	00.00	1.000	00.00	1.000
18.20	27.98	1.014	22.45	1.006	21.14	1.002
25.84	41.01	1.032	32.20	1.014	28.79	1.006
31.79	52.73	1.057	40.10	1.023	35.62	1.010
36.87	65.84	1.093	47.21	1.034	41.54	1.015
41.41	73.80	1.150	54.04	1.049	46.98	1.022
45.57	75.19	1.200	60.98	1.070	52.14	1.030
49.46	76.49	1.254	68.67	1.099	57.18	1.039
53.13	77.71	1.286	77.71	1.147	62.18	1.054
56.63	78.88	1.322	78.88	1.204	67.45	1.073
60.00	80.00	1.355	80.00	1.256	73.30	1.103
63.26	81.09	1.382	81.09	1.303	81.09	1.154
66.42	82.14	1.408	82.14	1.344	82.14	1.224
69.51	83.17	1.430	83.17	1.380	83.17	1.286
72.54	84.17	1.448	84.17	1.413	84.17	1.341
75.52	85.17	1.464	85.17	1.439	85.17	1.389
78.46	86.15	1.477	86.15	1.461	86.15	1.428
81.37	87.13	1.487	87.13	1.478	87.13	1.457
84.26	88.09	1.494	88.09	1.487	88.09	1.477
87.13	89.05	1.499	89.05	1.492	89.05	1.490
88.85	89.62	1.507	89.62	1.496	89.62	1.496

[0062] Now, evaluate the integral over n′e. Having factored the spherical harmonic Ylm(n′ec) out of the integral in Rlm(ne,n), the remaining integral I is

I=∫dn′ e exp(H(n·n65 ,nγ·ne)), H=γn·nγ−ln((nγ·ne)/G)

[0063] As ∂H/∂nγ∥=0 determines nγ=nγc, and ∂H/∂nγ⊥=0 for nγ=nγc by symmetry, we can write approximately that 4$I=\exp \left(H_c\right)\int d{n}_e^{\prime }\exp \left(-\frac{1}{2}\left({A_{\&LeftDoubleBracketingBar;}\left(n_{e\&LeftDoubleBracketingBar;}^{\prime }-n_{e\&LeftDoubleBracketingBar;c}^{\prime }\right)}^2+A_{\perp }{n}_{e\perp }^{\prime 2}\right)\right),$

I= 2πexp(Hc)/(A∥A⊥)1/2, Hc=H(n·nγc,nγc·ne)

[0064] Here A∥=−∂2H/∂n′e∥2, and A⊥=−∂2H/∂n′e⊥2, both evaluated for nγ=nγc. The cross term ∂2H/∂nγ∥∂nγ⊥ is zero by symmetry.

[0065] To evaluate A∥ and A⊥ we use a variety of simplifications. First, changing nγ to (nγ+dnγ⊥), it is necessary to include a factor of (1+dnγ⊥2)1/2≈1−dnn⊥/2 to preserve the normalization. Thus ∂nγ/∂nγ⊥=0 and ∂2nγ/∂2nγ⊥=−nγ/2. (Also ∂2nγ/∂nγ∥∂nγ⊥=0.) To obtain the ∂2nγ/∂n′e⊥2 recall that dn′e⊥=2(nγ·ne)dnγ⊥. Thus (for G=1) A⊥ becomes 5$\begin{array}{c}{A}_{\perp }=-\frac{1}{4}(\gamma \frac{1}{n_{\gamma }·n_e}\frac{\partial }{\partial n_{\gamma \perp }}\left(\frac{1}{n_{\gamma }·n_e}\frac{\partial n·n_{\gamma }}{\partial n_{\gamma \perp }}\right)-\\ \frac{1}{n_{\gamma }·n_e}\frac{\partial }{\partial n_{\gamma \perp }}\left(\frac{1}{n_{\gamma }·n_3}\frac{\partial \ln \left(n_{\gamma }·n_e\right)}{\partial n_{\gamma \perp }}\right))\end{array}$

[0066] where ∂(n·nγ)/∂nγ⊥=0,∂2(n·nγ)/∂nγ⊥2=−(n·nγ)/2, and similarly for (nγ·ne). Evaluating one finds A⊥=(γn·nγc−1)/8(nγc·ne)2 For A∥ recall that dn′e∥=2dnγ∥ and that ∂(n·nγ)/∂nγ∥=−(1−(n·nγ)2)1/2. A similar calculation, also for G=1, yields 6$A_{\&LeftDoubleBracketingBar;}=\frac{1}{4}\left(\gamma n·n_{\gamma c}-{\left(n_{\gamma c}·n_e\right)}^{-2}\right)$

[0067] Although the second term on the right appears to become large as n·ne decreases, it will be apparent that taking G=1 is no longer correct. As the dependence is not a particularly sensitive one compared to the spherical-harmonic dependence, approximate 7$A_{\&LeftDoubleBracketingBar;}\approx \frac{1}{4}\left(\gamma ·n_{\gamma c}-1\right)$

[0068] as an upper bound. Combining the above finally that results

R * lm (ne,n)=21/2γeγn nγcYlm(n′ec)/(γn·nγc−1)(eγ−1)

[0069] where n′ec=2nγc(nγc·ne)−ne.

[0070] Harmonic Representation

[0071] The approximation that the spherical harmonic Ylm(n′e) is relatively constant over the most significant portion of the integral over n′e is clearly not accurate for sufficiently large l. Consider the simpler case of ne=n: reflection in the normal direction. Define

I l (γ)=∫dn′eeγn nγG(n·n′e)Yl0(n′e)/4n·nγ=2π((2l+1)/4π)1/2∫2−1/21dueγuPl(2u2−1)g(u)

[0072] where Pl(x) is the Legendre function of order l:

P l (x)=dl(x2−1)l/dxl/2ll !

[0073] and where g(u)=1 for 31/2/2≦u≦1, and g(u)=2(2u2−1) for 21/2/2≦u≦31/2/2. Note that Il(γ)=NRl0(n,n). Let I*l(γ)be Il(γ) with Pl(2u2−1) set equal to its value for u=1, that is set to unity. Then the ratio

I l /I * l (γ)

[0074] gives a measure of the reduction of Rlm with increasing l.

[0075] While one can evaluate Il(γ) completely, consider the following simplification. For large γ, only the portion Pl(u) of near u=1 will matter. Now Pl(1)=1 for all l, and

P′ l (1)=lPl(1)+P′l−1(1)

[0076] according to the standard recursion relation xP′n(x)−P′n−1(x)=nPn(x) .

[0077] Thus P′l(1)=l(l+1)/2, and for l>0 we can replace Pl(x) with alx+bl, xl≦x≦1, xl=−bl/al, al=l(l+1)/2, bl=1−al, a straight line with slope al. After a little algebra one obtains

I l (γ)≈((2l+1)π)1/2(eγ/γ)(1−h−1+e−h/h)

[0078] with h=2γ/l(l+1). Comparing with

I * l (γ)=((2l+1)π)1/2(eγ−1)/γ

[0079] where Ylm(1) is pulled out of the integrand, one can replace the R*lm(ne,n) given above with

R lm (ne,n)=Il(γ)R*lm(ne,n)/I*l(γ)

[0080] to correct for the 1-dependence of Ylm(n′e) in the integral over facets. The I/I* goes as h/2(1−e−γ) for h small.

[0081] A complete evaluation of Il(γ) is performed. To do this we make use of several expressions: First, we can write any Legendre polynomial in the form 8$P_l\left(x\right)=\sum _{r=0}^{\left[1/2\right]}{\left(-1\right)}^r\left(2l-2r\right)!x^{l-2r}/2^{l}r!\left(l-r\right)!\left(l-2r\right)!$

[0082] Here [z]=the greatest integer contained in z.

[0083] Second, we will need the binomial expansion for (2u2−1) 9${\left(2u^2-1\right)}^{l-2r}=\sum _{k=0}^{l-2r}\left(l-2r\right)!{\left(2u^2\right)}^k{\left(-1\right)}^{l-k}/k!\left(l-2r-k\right)!$

[0084] Third, the indefinite integral 10$\Gamma \left(u;\gamma ,k\right)=\int du{e}^{\gamma u}{u}^{2k}=\left(e^{\gamma u}/{\gamma }^{2k+1}\right)\sum _{n=0}^{2k}{\left(-\gamma u\right)}^{2k-n}\left(2k\right)!/\left(2k-n\right)!$

[0085] Finally, defining 11$\begin{array}{c}{I}_l\left(b,a\right)={\int }_a^{b}du{e}^{\gamma u}{P}_l\left(2u^2-1\right)\\ I_{\mathrm{gl}}\left(b,a\right)={\int }_a^{b}du{e}^{\gamma u}2\left(2u^2-1\right)P_l\left(2u^2-1\right)\end{array}$

[0086] we can write Il(γ)=Il(1,31/2/2)+Igl(31/2/2,21/2/2) yielding Rlm(ne,n)≈R*lm(ne,n)Il(γ)/I*(r), with I*(γ)=I*(1,31/2/2)+I*g(31/2/2,21/2/2),I*(u)=exp(γu)/γ, and I*g(u)=exp(γu) (γ2u2−2γu+2)/γ3. After some algebra, we find the indefinite integrals Il(u) and Igl(u) in Il(b,a), respectively to be 12$\begin{array}{c}{I}_l\left(u\right)=\sum _{k=0}^l{C}_0\left(l,k\right)\Gamma \left(u;\gamma ,k\right)\\ I_{\mathrm{gl}}\left(u\right)=-2C_0\left(l,0\right)e^{\gamma u}/\gamma +2\sum _{k=1}^{l+1}{C}_1\left(l,k\right)\Gamma \left(u;\gamma ,k\right)\\ C_0\left(l,k\right)=\sum _{r=0}^{\left[\left(l-k\right)/2\right]}\frac{{\left(-1\right)}^{l+r-k}\left(2l-2r\right)!}{2^{l-k}r!\left(l-r\right)!k!\left(l-2r-k\right)!}\\ C_1\left(l,k\right)=\sum _{r=0}^{\left[\left(l-k+1\right)/2\right]}\frac{{\left(-1\right)}^{l+r-k+1}\left(2l-2r\right)!\left(l-2r+1\right)}{2^{l-k}r!\left(l-r\right)!k!\left(l-2r-k+1\right)!}\end{array}$

[0087] From these the l-dependent factor Il(γ)/I*l(γ) follows.

[0088] The method of deriving the reflectance function described above is used in the object recognition process. To match an input image with an object model in the database, the object model must be positioned to match the position of the input image. There are standard methods of performing this function known in the art, though a preferred method is disclosed in U.S. patent application Ser. No. 09/538,209, filed Mar. 30, 2000, entitled “Method for Matching a Two Dimensional Image to One of a Plurality of Three Dimensional Candidate Models in a Database”, the disclosure of which is hereby incorporated by reference.

[0089] Once positioned, images of the model are rendered under each of one or more spherical harmonic components of light as described above. These component rendered images can be linearly combined to arrive at an optimal rendered image of the model. An optimal rendered image means a rendered image that most closely approximates the input image. The comparison of the various rendered images to the input image preferably is made by calculating the Euclidian Distance between the images, in a manner that is known in the art. The optimization is performed by a least squares fit technique between the vector that is the reflectance function and the matrix that represents the input image, also known in the art.

[0090] Other steps follow the optimization in the recognition process. One is computing some measure of similarity between the input image and the optimal rendered image for each model. Preferably, this is accomplished be measuring the Euclidian Distance between each optimal image and the input image. The object model whose optimal rendered image is nearest the input image, assuming some threshold value of similarity must be reached, is a match to the input image.

[0091] While the reflectance function of the present invention is described in application to the field of machine vision, it is not limited to that field. There are numerous parallels between this work and advances in the field of computer graphics, and our invention has applicability there as well. The field of computer graphics is continually striving to produce more realistic rendered images based on object models that exist only in the computer. The reflectance function of the present invention advances that work as well.

[0092] In summary, in accordance with the teachings of the present invention, there have been derived several harmonic expressions for the Torrance-Sparrow model of off-specular reflection. In the process there is a significantly simplification of both the derivation and the expression of the results of their model and restored reciprocity between incidence and reflection. There is noted an insensitivity of the position and magnitude of the off-specular reflection at large angles to the distribution of facet normals. This result is an advancement of the existing state of the art, since results which are too sensitive to parametric models are unreliable. Also results significantly different from those predicted would necessitate a qualitatively different model.

[0093] While there has been described a method of modeling the off-specular reflection of light to allow a machine to recognize the image of an object illuminated by varying unknown lighting conditions, it will be apparent to those skilled in the art that variations and modifications are possible without deviating from broad teachings of the present invention which shall be limited solely by the scope of the claims appended hereto.

Claims

1. A method of deriving the reflectance function of an object model under a variety of lighting conditions, the object model comprising a plurality of surfaces, each of the surfaces defining a normal vector pointing perpendicularly outward, the method comprising the steps of: (a) rendering, for a given orientation of the object model, a plurality of images produced by the object model when illuminated by each of a plurality of spherical harmonic components of incident light, wherein said rendering comprises: (i) calculating the intensity of incident light components upon the object model, relative to the normal at each surface; (ii) calculating the intensity of light components diffusely reflected by the model toward the observer; (iii) calculating the intensity of light components off-specularly reflected toward the observer, relative to the observer's angle with respect to the angle of maximum specular reflection; and (b) defining a reflectance function for the object as a linear combination of said plurality of images.

2. The method according to claim 1, wherein said plurality of spherical harmonic components include the 0th order spherical harmonic component of incident light.

3. The method according to claim 1, wherein the object model is a face.

4. The method according to claim 1, wherein said calculating the intensity of light components broadened-specularly reflected toward the observer comprises using a modified Torrance-Sparrow model of the object.

5. The method according to claim 4, wherein the object model is a face.

6. A method of recognizing objects under various lighting conditions comprising the steps of: (a) providing a database comprising a plurality of three dimensional models, (b) providing an input image, (c) positioning each three dimensional model relative to the input image, (d) determining, for each three dimensional model, a rendered image which is most similar to the input image, said determining step comprising: (i) deriving a reflectance function that describes an approximation of the set of all possible rendered images that each three dimensional model can produce under all possible lighting conditions, said rendered images including both diffusely and off-specularly reflected light, and said reflectance function employs a model of broadened-specular reflectance that accounts for the angle between the direction of maximum specular reflectance; and (ii) optimizing the reflectance function to determine rendered image of each model that is most similar to the input image; (e) computing a measure of similarity between the input image and each optimal rendered image; and (f) selecting the three dimensional model corresponding to the optimal rendered image whose measure of similarity is most similar to the input image.

7. The method of claim 6, wherein the input image is the image of a face.

8. The method of claim 7, where the three dimensional models are models of heads.

9. The method according to claim 6 wherein said off-specular reflectance model comprises a modified Torrance-Sparrow model.

10. The method of claim 9, wherein the input image is the image of a face.