A Sparse Non-parametric BRDF Model

The bidirectional reflectance distribution function [Nicodemus et al. 1977] describes how light scatters at the surfaces of a scene, depending on their material characteristics. The Bidirectional Reflectance Distribution Function (BRDF) is a 4D function parameterized by the incident and exitant scattering angles and can be described using either parametric models [Blinn 1977; Cook and Torrance 1982; Ashikhmin and Shirley 2000; Walter et al. 2007; Löw et al. 2012] or data-driven models [Lawrence et al. 2004; Jakob et al. 2014; Bilgili et al. 2011; Bagher et al. 2016; Sun et al. 2018; Tongbuasirilai et al. 2019]. Parametric models present great artistic freedom and the possibility to interactively tweak parameters to achieve the desired look and feel. However, most analytical models are not designed for efficient and accurate representation of the scattering properties of measured real-world materials. Data-driven models, however, enable the use of measured BRDFs and real-world materials directly in the rendering pipeline, and are commonly used in computer vision applications [Romeiro et al. 2008]. Data-driven models can represent BRDFs in many different ways. Iterative-factored representations approximate BRDFs with multiple low-rank components [Lawrence et al. 2004; Bilgili et al. 2011; Tongbuasirilai et al. 2019], while hybrid analytical data-driven models [Bagher et al. 2016; Sun et al. 2018] rely on non-parametric components or basis functions computed using specific weighting and optimization schemes.

The efficiency, or performance, of a non-parametric model is typically measured in terms of the number of variables/coefficients required to represent a BRDF at a given quality. The number of coefficients, in turn, depends on the efficacy of the underlying basis functions, often called dictionaries, in modeling the BRDF. Most, if not all, existing methods either sacrifice the model accuracy to achieve fast reconstruction for real-time applications, or aim for high image fidelity leading to increased storage and computational requirements. Another important aspect is the complexity of the basis functions used in the representation. At one end of the spectrum, we have analytical basis functions such as spherical harmonics and wavelets [Ramamoorthi and Hanrahan 2001; Claustres et al. 2003], which do not require the storage of the basis functions, at the cost of low approximation accuracy, unless a large number of coefficients is used. On the other end of the spectrum, we have decomposition-based methods [Bilgili et al. 2011] that model the BRDF as a multiplication of a set of coefficients and a basis matrix/tensor computed from data. Unfortunately, these approaches require a computationally expensive decomposition, e.g., Principal Component Analysis (PCA) or Singular Value Decomposition (SVD), for each BRDF individually and suffer from a high storage cost for the basis itself. Another problem is that the expressiveness of existing bases/decomposition methods is limited. Except for a few, they are in most cases not designed for BRDF data, hence requiring high numbers of coefficients for accurate BRDF representation.

The goal in this article is to develop a new data-driven BRDF model using a set of trained basis functions that enable highly accurate representations with a minimal number of coefficients to represent any BRDF. To solve this challenge, we derive a model that in essence relies on decomposing BRDFs into a coefficient–basis pair, but uses machine learning to adapt the basis to the space of BRDFs to provide maximally sparse coefficients. Sparse BRDF modeling is achieved using a novel BRDF dictionary ensemble and a novel model selection algorithm to efficiently represent a wide range of real-world materials. The learned BRDF dictionary ensemble consists of a set of basis functions trained such that they guarantee a very sparse BRDF representation and near optimal signal reconstruction. Moreover, our model takes into account the multidimensional structure of measured BRDFs (e.g., 3D or 4D depending on the parameterization) and can exploit the information redundancy in the entire BRDF space to reduce the number of coefficients. While in this article, we focus on isotropic materials, our method can be readily applied to anisotropic materials, given a data set of such materials.

The learned ensemble is versatile and can be trained only once to be reused for representing a wide range of previously unseen materials. Additionally, the dictionary ensemble is not restricted to a single BRDF transformation as previous models. Instead, multiple BRDF transformations can be included in the ensemble training such that for each individual BRDF the best representation can be automatically selected and used. This is achieved by a novel model selection method to pick a dictionary in the ensemble that leads to the sparsest solution, the smallest reconstruction error, and the most suitable transformation with respect to rendering quality. For the experiments and evaluations presented here, we use the MERL [Matusik et al. 2003] and RGL-EPFL [Dupuy and Jakob 2018] databases, which are divided into a training set and a test set used for evaluation.

Another key contribution of this article is a novel algorithm for interpolating two or more BRDFs in the sparse coefficient space. We propose to perform this task in the sparse coefficient space. The key challenge is when the two BRDFs to be interpolated do not lie in the same subspace, i.e., they use different dictionaries of the ensemble. We propose a novel formulation for interpolation of BRDFs that lie in different coordinate systems. The proposed interpolation algorithm admits the construction of a smooth surface over the coefficient space of all BRDFs represented using our model.

The main contributions of this article are as follows:

We compare the proposed non-parametric BRDF model to the current state-of-the-art methods and demonstrate that it performs significantly better in terms of rendering Signal-to-Noise Ratio (SNR) and visual quality, as well as having a lower reconstruction error. To the authors’ knowledge this is the first BRDF model based on sparse representations and dictionary learning.

Notations. Throughout the article, we use the following notational convention. Vectors and matrices are denoted by boldface lower-case (\(\mathbf {a}\)) and bold-face upper-case letters (\(\mathbf {A}\)), respectively. Tensors are denoted by calligraphic letters, e.g., \(\mathbf {\mathcal {A}}\). A finite set of objects is indexed by superscripts, e.g., \(\lbrace \mathbf {A}^{(i)}\rbrace _{i=1}^{N}\), whereas individual elements of \(\mathbf {a}\), \(\mathbf {A}\), and \(\mathbf {\mathcal {A}}\) are denoted \(\mathbf {a}_i\), \(\mathbf {A}_{i_1,i_2}\), \(\mathbf {\mathcal {A}}_{i_1,\dots ,i_n}\), respectively. The determinant of \(\mathbf {A}\) is denoted \(\mathrm{det}(\mathbf {A})\). The \(\ell _p\) norm of a vector \(\mathbf {s}\), for \(1\le p\le \infty\), is denoted by \(\Vert \mathbf {s}\Vert _p\). Frobenius norm is denoted \(\Vert \mathbf {s}\Vert _F\). The \(\ell _0\) pseudo-norm of a vector, denoted \(\Vert \mathbf {s}\Vert _0\), defines the number of non-zero elements in the vector. The n-mode product of a tensor \(\mathbf {\mathcal {S}}\) and a matrix \(\mathbf {U}\) is denoted \(\mathbf {\mathcal {S}}\times _n\mathbf {U}\).

Measured BRDFs have proven to be an important tool in achieving photo-realism during rendering [Matusik et al. 2003; Dupuy and Jakob 2018; Dong et al. 2016]. Even highly complex surfaces such as layered materials require multiple components of measured data to construct novel complex materials [Jakob et al. 2014]. Measured materials, however, are high-dimensional signals with large memory footprint and a key challenge is that small approximation errors can lead to visual artifacts during rendering. To efficiently represent such high-dimensional measured BRDF data, one can use parametric models, or data-driven models, since densely sampled BRDF data imposes a large memory footprint, making it impractical to use in many applications.

Parametric models. By careful modeling, BRDFs can be encoded with only a few parameters. The components or factors of such models are based on either assumptions describing by the physics of light—surface interactions using, e.g., microfact theory [Cook and Torrance 1982; Walter et al. 2007; Holzschuch and Pacanowski 2017], or empirical observations of BRDF behaviors [Blinn 1977; Ward 1992; Ashikhmin and Shirley 2000; Nishino and Lombardi 2011; Löw et al. 2012]. However, in many practical cases and applications, existing parametric models cannot accurately fit measured real-world data, as mentioned in Bagher et al. [2016].

Data-driven models. Due to their non-parametric property, data-driven models are superior to parametric models in that the number of degrees of freedom, or implicit model parameters, is much higher. This means that the representative power is higher and the expected approximation error is lower. Factored BRDF models use decomposition techniques to factorize BRDF into several components. Matrix and tensor decompositions have been used by Lawrence et al. [2004], Soler et al.[2018], Bilgili et al. [2011], and Tongbuasirilai et al. [2019]. Moreover, factorization-based models have been used for interactive BRDF editing [Kautz and McCool 1999; Ben-Artzi et al. 2006] and real-time rendering [Soler et al. 2015]. Analytical basis functions have been used for BRDF measurements using spherical harmonics [Ghosh et al. 2010;, 2009; Tunwattanapong et al. 2013] and BRDF modeling using rational functions [Pacanowski et al. 2012].

A problem with existing factored models is that rank-1 approximations in most cases lead to inferior results. Accurate modeling requires iterative solutions with many layered factors. Analytic-data-driven BRDF models [Bagher et al. 2016; Sun et al. 2018] employ analytical models extended to a higher number of parameters fitted with measured data to achieve higher accuracy. Cooper et al. [2021] employ the Gaussian mixture model for BRDF representation and clustering that is utilized in the inverse rendering problem. The recent advancement of machine learning algorithms, in particular deep learning, have led to new research directions on BRDF-related topics [Dong 2019]. Deep learning has been used for BRDF editing [Hu et al. 2020], BRDF acquisition [Deschaintre et al. 2019;, 2018; Li et al. 2018], compact BRDF representation for importance sampling [Sztrajman et al. 2021], and inverse rendering using invertible neural BRDFs [Chen et al. 2021]. Moreover, deep learning has also been used for Bidirectional Texture Function (BTF) compression, see, e.g., Rainer et al. [2019].

Dictionary Learning. One of the most commonly used dictionary learning methods is K-SVD [Aharon et al. 2006], and its many variants [Mukherjee et al. 2015; Mazhar and Gader 2008; Rusu and Dumitrescu 2012; Marsousi et al. 2014], where a 1D signal (i.e., a vector) is represented as a linear combination of a set of basis vectors, called atoms. A clear disadvantage of K-SVD for BRDF representation is signal dimensionality. For instance, a measured BRDF in the MERL data set, excluding the spectral information, is a \(90\ \times \ 90\ \times \ 180=\text{1,458,000}\)-dimensional vector. In practice, the number of data points needed for K-SVD dictionary training should be a multitude of the signal dimensionality to achieve a high-quality dictionary. In addition to unfeasible computational power required for training, the limited number of available measured BRDF data sets renders the utilization of K-SVD impractical.

In contrast to 1D dictionary learning methods, multidimensional dictionary learning has received only little attention in the literature [Roemer et al. 2014; Hawe et al. 2013; Ding et al. 2017]. In multidimensional dictionary learning, a data point is treated as a tensor, and a dictionary is trained along each mode. For instance, given our example above, instead of training one \(\text{1,458,000}\)-dimensional dictionary for the MERL data set, one can train three dictionaries (i.e., one for each mode), where the atom size for these dictionaries are 90, 90, and 180, corresponding to the dimensionality of each mode. To the best of our knowledge, there exists only a few multidimensional dictionary learning algorithms. Our sparse BRDF model in this article is inspired by the multidimensional dictionary ensemble training proposed in Miandji et al. [2019], which has been shown to perform well for high dimensional signals such as light fields and light field videos. We will elaborate on our training scheme for BRDFs in Section 3.2.

Our non-parametric model is based on learning a set of multidimensional dictionaries, a dictionary ensemble, spanning the space of BRDFs, i.e., the space in which each BRDF is a single multi-dimensional point. Each dictionary in the ensemble consists of a set of basis functions (i.e., orthonormal matrices), representing each dimension of the BRDF space, that admits sparse representation of any measured BRDF using only a small number of coefficients as illustrated in Figure 1. The dictionary ensemble is trained only once on a given training set of measured BRDFs and can then be reused to represent a wide range of different BRDFs. This is in contrast to previous models that use tensor or matrix decomposition techniques, where the basis and the coefficients are calculated for each BRDF individually.

A major challenge when using machine learning methods, and in particular dictionary learning, on BRDFs is the high dynamic range inherent to the data. In Section 3.1, we describe two data transformations that when applied on measured BRDFs, they improve the fitting to our non-parametric model, see Section 4. The training of the multidimensional dictionaries is described in Sections 3.2 and 3.3, followed by our model selection technique in Section 3.4, where we describe a method to select the most suitable dictionary in the ensemble for any unseen BRDF such that the coefficients are maximally sparse, the modeling error is minimal, and that the data transformation used is one that leads to a better rendering quality. Finally, in Section 3.5, we present a novel algorithm for interpolating two or more measured BRDFs, represented using our model, directly in the representation space (i.e., using the sparse coefficients).

A BRDF can be parameterized in many different ways [Rusinkiewicz 1998; Löw et al. 2012; Stark et al. 2005; Barla et al. 2015]. Our dictionary learning approach does not rely on the parameterization of given BRDFs as long as the resolution of these BRDFs is the same. For simplicity and to facilitate a fair comparison to previous works, all the data sets we use here are based on the Rusinkiewicz’s parameterization [Rusinkiewicz 1998] defined by \((\theta _h, \theta _d, \phi _d)\), where \(\theta _h, \theta _d\) are elevation angles of the half vector and the difference vector, respectively, \(\phi _d\) is the azimuthal angle of the difference vector. The sample BRDFs we use here are of resolution \(90\times 90\times 180\). The data sets are sampled with MERL encoding where a square-root function is applied to \(\theta _h\); see Matusik et al. [2003].

This section presents an evaluation of the proposed BRDF model and comparisons to the current state-of-the-art models in terms of BRDF reconstruction error and rendering quality. The rendering results were generated using PBRT [Pharr and Humphreys 2010] with the Grace Cathedral environment map. The images were rendered at a resolution of \(512\times 512\) pixels using 512 pixel samples in PBRT with the directlighting surface integrator and 256 infinite light-source samples.

The BRDF dictionary was trained using materials from the MERL database [Matusik et al. 2003] and RGL-EPFL isotropic BRDF database [Dupuy and Jakob 2018]. We split the MERL and RGL-EPFL materials into a training set and a test set. The training set contains 136 materials, where 85 materials are from the MERL dataset and 51 materials are from the RGL-EPFL dataset. The test set contains 28 materials with 15 materials from the MERL dataset, 8 materials from the DTU data set [Nielsen et al. 2015] and 5 materials from RGL-EPFL. The training and test sets cover a wide range of material classes and none of the materials in the test set appear in the training set.

Each BRDF color channel is processed independently for the training and model selection. We use the Rusinkiewicz’s parameterization [Rusinkiewicz 1998], at a resolution of \(90\times 90\times 180\), i.e., we have \(m_1 = 90\), \(m_2 = 90\), and \(m_3 = 180\). For our experiments, we trained four ensembles, each with \(K = 8\) dictionaries and with training sparsities of \(\tau _l=32\), \(\tau _l=64\), \(\tau _l=128\), and \(\tau _l=256\). We then construct one ensemble by taking the union of the dictionaries in the four ensembles that were trained, as described in Section 3.3. The training BRDFs were transformed using log-plus (\(\rho _{t1}\)) and cosine-weighted-log (\(\rho _{t2}\)) functions before starting the training, hence resulting in 272 materials. The training took about 11 hours on a machine with 40 cores clocked at 2.4 GHz.

Once the ensemble is trained, we use the model selection algorithm, described in Section 3.4, to obtain the optimal sparse coefficients of each BRDF in the test set. The model selection for each BRDF and for all three color channels takes about 61 seconds using a CPU with 16 cores clocked at 3.8 GHz. Note that the model selection can be parallelized over the BRDFs in the test set, since this task is performed independently for each BRDF. Moreover, the reconstruction of all the elements of a measured BRDF, i.e., a total of \(\text{4,374,000}\) elements, takes 17 seconds using the same hardware setup. It should be noted that during rendering, only one or a few values from the measured BRDF is needed at each intersection point. In this case, we utilize the algorithm proposed in Miandji et al. [2019] to recover a single element. In addition, since the reconstruction in Equation (7) is based on the n-mode product of a sparse tensor and dense matrices, we can utilize the GPU-based algorithm in Baravdish et al. [2019] for real-time applications. We have left this application for future work. A CPU-based implementation of our method for model selection and reconstruction is provided at https://github.com/emiandji/sparse-brdf.

To evaluate our sparse BRDF model, we use two quality metrics: SNR that is calculated on the rendered images (floating-point images) and Relative Absolute Error (RAE), which is computed on linear BRDF values. Note that since the reference BRDF is in linear BRDF domain, i.e., it is not transformed, we invert the transformation of the reconstructed BRDF (which may be \(\rho _{t1}\) or \(\rho _{t2}\) depending on the model selection result) to transform it back to the linear space prior to computing the RAE. The RAE is defined aswhere \(\rho _{ref}\) is the reference BRDF, and \(\rho _{recon}\) is the reconstructed BRDF. Even though rendering SNR (or PSNR) is typically used to evaluate BRDF models in many publications, RAE is very useful to capture the model accuracy without relying on a specific rendering setup. This is because RAE is computed on all BRDF values, while rendering SNR only takes into account a portion of the BRDF values, depending on the rendering setup. Note that for rendering, we invert Equations (1) and (2) to transform the reconstructed BRDFs back to the original linear domain, which is required by the renderer.

We compare our results to Bagher et al. [2016] (naive model), Bilgili et al. [2011] (Tucker decomposition) and Tongbuasirilai et al. [2019] (rank-1 CPD decomposition with L = 1) on 15 MERL test materials. The method of Bagher et al. [2016] stores \((90+90+180+2) = 362\) coefficients per channel, Bligili et al. [2011] uses \((128+16+16+64+2) = 226\) coefficients, and the CPD decompositions from Tongbuasirilai et al. [2019] uses \((90+90+180) = 360\) coefficients per channel. Since the Tucker and CPD methods use an iterative approach, we limit our comparisons to L = 1, i.e., a single factorization was performed so that the number of coefficients used for all models were roughly the same. The CPD method was tested using two different parameterizations: the PDV [Löw et al. 2012; Tongbuasirilai et al. 2019] and HD [Rusinkiewicz 1998] parameterizations.

To the best of our knowledge, the model of Bagher et al. [2016] is the current state-of-the-art for non-parametric BRDF models. Therefore, we set the number of coefficients \(\tau _t\) for our model such that we match the storage complexity of Bagher et al. [2016]. Since our representation is sparse, we require \(1+1+1\) bytes for storing the location of each nonzero element in a sparse coefficient tensor with size \(90\times 90\times 180\), as well as 8 bytes for each nonzero value. Simple calculations show that by using \(\tau _t=262\) coefficients for our model for each color channel, we can match the storage complexity of Bagher et al. [2016], which uses 362 coefficients to model each color channel of a BRDF. Note that the discrepancy in the number of coefficients is due to the fact that we need to store the location of nonzero values, which penalizes our method to use less coefficients.

For the rendered images, shown in Figures 2, 3, and 9, we apply gamma-corrected tone-mapping. The error images, also known as false-color, produced by normalizing the absolute image-space error in the range \([0,1]\), followed by applying a jet color map using MATLAB. Normalization is done individually for each BRDF represented using different models. All the error images are multiplied by 10 to facilitate visual comparisons. Apart from what follows, additional quantitative and qualitative results are included in the supplementary materials.

Table 3 reports SNR statistics for 15 test materials in the MERL database. The average SNR of our model is about 8dB, 5dB, and 10dB higher for log-plus, cosine-weighted-log and the proposed model selection based on Gamma-mapped-MSE, respectively, when compared to Bagher et al.; moreover, our results show a smaller standard deviation on SNR. Additionally, our proposed model selection method can achieve higher SNR on average compared to the cases where we only use log-plus or cosine-weighted-log. The lower standard deviation indicates that the proposed model can represent the MERL materials more faithfully. Table 4 shows a direct comparison of our model to that of Bagher et al. [2016] for each BRDF in the MERL test set using rendering SNR and BRDF-space RAE. Here, we use our Gamma-mapped-MSE metric to choose between the transformations. Compared to the model of Bagher et al., our approach achieves significantly higher visual quality on 13 of 15 materials. In Figure 2, we present example renderings of four BRDFs in the MERL test set modeled using our method and Bagher et al. [2016]. Our results are obtained using the proposed Gamma-mapped-MSE for model selection. The HDR environment map used here is Grace Cathedral.

The cosine-weighted-log transformation suppresses the grazing angle BRDF values. For diffuse materials, this leads to better visual results and significantly higher rendering SNR. It is evident from Figure 2 that the log-plus transformation is better for glossy materials as the cosine-weighted-log transformation leads to color artifacts for some materials, e.g., gold-metallic-paint2, red-metallic-paint, and violet-acrylic. More results for further analysis is available in the supplementary material.

Table 5 below shows rendering SNR and BRDF-space RAE values for the RGL-EPFL test set using both \(\rho _{t1}\) and \(\rho _{t2}\). Our BRDF model and selection method can efficiently represent the RGL-EPFL data set with an average SNR of more than 38dB. Our model selection method on the RGL-EPFL test set missed on 2 of 5 materials, which are cc-amber-citrine-rgb and vch-dragon-eye-red-rgb. The SNR values demonstrate that our data-driven model can accurately represent and faithfully reconstruct the BRDFs in the test set. See the supplementary materials for rendered images obtained using our model applied on the RGL-EPFL data set. Although we use a fixed number of coefficients for all materials here for the sake of comparisons, i.e., \(\tau _t=262\), ideally we would like to tune this parameter for each BRDF individually. For instance, for challenging materials such as cc-amber-citrine-rgb, the model accuracy significantly increases with higher number of coefficients. We analyze the effect of \(\tau _t\) on rendering SNR and BRDF RAE in Figure 4. Rendering results for two example BRDFs when we set \(\tau _t\) to 16, 64, 256, 1,024, and 4,096 are included in Figure 5.

To demonstrate the robustness of our sparse non-parametric model for representing BRDFs in a test set, we also evaluate it using eight test samples in the DTU data set [Nielsen et al. 2015]. Note that we use the same dictionary described above and that none of the materials from the DTU data set were used in the training set. The results are summarized in Table 6, where we report rendering SNR and BRDF-space RAE for \(\rho _{t1}\), \(\rho _{t2}\), and our model selection based on Gamma-mapped-MSE. Our BRDF model and selection method can reproduce the DTU data set with average SNR of more than 48dB. Our model selection algorithm on the DTU test set missed on 2 of 8 materials, which are blue-book and green-cloth. Visual quality examples of the rendered images are presented in Figure 3. The difference between \(\rho _{t1}\) and \(\rho _{t2}\) is evident in this figure. We can see that \(\rho _{t1}\) is favored by glossy materials, while \(\rho _{t2}\) is more effective in modeling low-frequency or diffuse-like materials. For instance, we observe that the artifacts seen on the cardboard and green-cloth renderings for log-plus (\(\rho _{t1}\)) do not appear in cosine-weighted-log (\(\rho _{t2}\)) renderings. For this figure, we used the Grace Cathedral HDR environment map.

Our results confirm the discrepancy between BRDF-space error metrics (such as RAE) and rendering quality using SNR. For example, blue-metallic-paint in Table 4, cardboard in Table 6, and vch-dragon-eye-red-rgb in Table 5 demonstrate how RAE contradicts the rendering SNR. The lower the BRDF-space RAE is, the more accurate the model represents a BRDF. However, a rendered image is dependent on a variety of additional factors such as geometry of objects, lighting environment, and viewing position. As a result, the BRDF-space RAE and rendering SNR have to be considered together for the evaluation of a BRDF model.

Figure 4(a) demonstrates the effectiveness of our BRDF model on all data sets with respect to the number of coefficients, \(\tau _t\). In terms of the decline of error with respect to the number of coefficients, all data sets show a similar behavior. Moreover, all the data sets have poor RAE when \(\tau _t \lt 32\). This behavior is expected, since the training sparsity is set to \(\tau _l = 32\), \(\tau _l = 64\), \(\tau _l = 128\), and \(\tau _l = 256\). And when we have \(\tau _t \lt \tau _l\), we expect a decline in performance. Figure 4(b) indicates that our model generally produces incrementally better rendering results when the number of coefficients is increased for all the data sets. Our model performs best for the MERL data set, because most of the training data are from MERL, while EPFL has fewer number of materials in the training. Although we do not include any materials from the DTU data set in the training, the rendering SNR and BRDF RAE are competitive with MERL reconstructions.

We also evaluated our BRDF model with the PCA-based method presented in Nielsen et al. [2015]; see Figure 6. Our RAE performs similar to PCA-based method with 32 principal components and when we have \(\tau _t=32\). We see that for a lower number of coefficients the PCA-based method outperforms our method, since our minimum training sparsity is 32. It should be noted that the PCA dictionary exhibits a very high storage cost compared to our BRDF dictionary ensemble. The PCA dictionary contains \(\hbox{1,458,000}\times p\) elements, where p is the number of principal components, while our combined dictionary ensemble consists of \((90*90+90*90+180*180)*32=\text{1,555,200}\) elements. This shows that the storage complexity of the PCA dictionary is proportional to the number of principal components, whereas for our dictionary ensemble the storage complexity is constant. In Figures 7(a), 7(b), and 7(c), we present the rendering SNR and in Figures 7(d), 7(e), and 7(f), we present the logarithmic BRDF-space RAE for all three data sets in comparison to Nielsen et al. [2015], where we compensate for the size of the PCA dictionary by increasing the number of coefficients of our method, \(\tau _t\). When considering the size of the dictionaries, the PCA method with 2,3 and 4 principal components corresponds to our model with \(\tau _t= \hbox{453,601}\), \(\tau _t= \hbox{939,601}\), and \(\tau _t= \hbox{1,425,601}\), respectively. The bar plots indicate that, for all data sets, our model performs better than the PCA technique by a large margin on both quality metrics.

To demonstrate our BRDF interpolation method on the sparse coefficient space (the sparseland), we provide a number of challenging examples in Figure 8 for different material classes using Equation (19) introduced in Section 3.5. The interpolation coefficient, \(\alpha\), was set between 0.0 and 1.0 with increments of 0.1. Figure 8(a) illustrates the interpolation between two materials, i.e., red-metallic-paint and specular-yellow-phenolic, and Figure 8(b) shows the interpolation of two diffuse materials, namely, dark-red-paint and pink-fabric2. Figure 8(c) gives an example of a layered material, cc-amber-citrine-rgb, and a glossy material, silver-metallic-paint2. Finally, Figure 8(d) is an example of vinyl film, vch-dragon-eye-red-rgb, and a diffuse material, red-fabric. Note that each BRDF pair we use here share the same transformation (i.e., \(\rho _{t1}\) or \(\rho _{t2}\) as discussed in Section 3.1). As shown in the rendered images, our BRDF interpolation on the sparse coefficient space provides a smooth transition between the material pairs. For these experiments, we used \(\tau _t=512\) for the model selection algorithm. Moreover, as mentioned in Section 3.5, we run the projection step of the model selection algorithm on the interpolated BRDF to obtain a desired sparsity, in this case \(\tau _t=512\). Note that this parameter is user-defined and can be changed to obtain a trade-off between quality and storage complexity for the interpolated BRDF.

In Figure 9, we evaluate our BRDF model using the Princeton scene with the following materials: blue-metallic-paint, gold-metallic-paint2, pink-fabric2, silver-metallic-paint2, and specular-yellow-phenolic. We rendered the scene with path tracing in PBRT [Pharr and Humphreys 2010] using the uffizi environment map and with \(2^{17}\) samples-per-pixel. Figures 9(a) and 9(c) present rendered images from our model and Bagher et al., respectively. Our model achieves an 8.03dB advantage in SNR over the model of Bagher et al.

The experiments conducted in this article are limited to isotropic BRDFs. This is partially due to the lack in availability of anisotropic BRDF data. As part of the future work, we aim to extend our sparse models to efficiently represent also anisotropic materials. Moreover, we acknowledge the fact that the discrepancy between BRDF-space errors and the rendering quality is still an open problem. Although we showed significant improvements using our Gamma-mapped-MSE, we believe that a more sophisticated metric that takes into account the support of the BRDF function can improve our results. Our model is relatively robust to noise. However, we believe that an application of a denoising pass that is tailored to measured BRDFs, prior to training and model selection, can greatly improve our results. This is expected, since it is well-known that even a moderate amount of noise in measured BRDFs translates to lower rendering quality; and that noise reduces the sparsity of the representation, hence increasing the model complexity. An alternative to applying a denoiser is to modify the training and model selection methods to be noise-aware.

The proposed interpolation algorithm enables a variety of applications where the BRDFs processing can be performed directly in the coefficient space rather than the substantially higher-dimensional BRDF space. For instance, rendering of smooth surfaces where each point on the surface has a distinct measured BRDF (i.e., spatially varying BRDFs) can be efficiently performed using our model. Moreover, editing of measured BRDFs can be performed directly in the coefficient space by utilizing the shallow network proposed in Hu et al. [2020] that maps our sparse representation of materials to known BRDF properties. We have left these applications of our model for future work. Another interesting venue for future work is the utilization of our model for sampling BRDFs, using, e.g., a gonioreflectometer, where only a few samples are taken to reconstruct a full BRDF, similar to Nielsen et al. [2015].

Our model is GPU-friendly, because it consists of parallel operations where atoms and coefficients can be efficiently distributed and operated among GPU cores. There exists a variety of GPGPU algorithms for performing key tasks in our method; for instance, a GPU-based implementation of the n-mode product is described in Baravdish et al. [2015]. As the main focus in this article is to derive and evaluate a novel non-parametric BRDF model, we have left the GPU implementation of our method to future work. However, given the small memory footprint of our model, and fast BRDF reconstruction using Equation (7), existing real-time ray tracing systems can be modified to accommodate our model for measured materials. Moreover, since the quality of our model is solely controlled using the \(\tau _t\) parameter, one can linearly control the rendering performance using this parameter.

Quality metrics for effective BRDF fitting is an open problem [Lavoué et al. 2021]. In this article, we showed the discrepancy between the results of RAE and rendering SNR. We then propose a model selection technique to reduce the gap between the visual quality and the numerical quality metrics. We believe that our model selection algorithm can benefit from more accurate error metrics to further increase the visual quality. Due to the robustness of our model selection algorithm, any new mathematical error metric that is suitable for measured BRDFs can be used in Equation (6).

One of the limitations of the proposed model is the assumption of same resolution for the BRDFs in the training and testing sets. This is a limiting factor for BRDF data sets that are not densely sampled. One possible solution for these scenarios is that we train the model based on densely sampled BRDFs. For a given sparsely sampled BRDF to be modeled, we first create a dense BRDF by copying known values and setting unknown values to zero. We then iteratively perform the model selection algorithm on the obtained BRDF. At the first iteration, we obtain a crude approximation for the unknown values. This crude approximation is then improved iteratively by the consecutive runs of the model selection algorithm. We have left a full exploration of such method for future work.

This paper presented a novel non-parametric sparse BRDF model in which a measured BRDF is represented using a trained multidimensional dictionary ensemble and a set of sparse coefficients. We showed that with careful model selection over the space of multidimensional dictionaries and various BRDF transformations, we achieve significantly higher rendering quality and model accuracy compared to current state-of-the-art. Finally, we introduced a novel algorithm for interpolation of two or more BRDFs, modeled using our approach, directly in the coefficient space, even when the given BRDF data points used for the interpolation utilize distinct dictionaries. We evaluated the performance of our model and algorithm using three different data sets, MERL, RGL-EPFL, and one provided by Nielsen et al. [2015]. For the vast majority of the BRDFs used in the test set, we achieve a significant advantage over previous models.