Computational Visual Media Computational Adaptive slices for acquisition of anisotropic BRDF Adaptive slices for acquisition of anisotropic BRDF Adaptive slices for acquisition of anisotropic BRDF

BRDF continues to be used as a fundamental tool for representing material appearance in computer graphics. In this paper we present a practical adaptive method for acquisition of anisotropic BRDF, based on sparse adaptive measurement of the complete four-dimensional BRDF space by means of one-dimensional slices, which form a sparse four-dimensional structure in the BRDF space, and can be measured by continuous movements of a light source and sensor. Such a sampling approach is advantageous especially for gonioreﬂectometer-based measurement devices where the mechanical travel of a light source and a sensor imposes a signiﬁcant time constraint. In order to evaluate our method, we have performed adaptive measurements of three materials and we simulated adaptive measurements of thirteen others. This method has one quarter the reconstruction error of that resulting from regular non-adaptive BRDF measurements using the same number of measured samples. Our method is almost twice as good as a previous adaptive method, and it requires from two to ﬁve times fewer samples to achieve the same results as alternative approaches.


Introduction
Realistic appearance of spatially homogeneous materials is usually represented by means of a bidirectional reflectance distribution function (BRDF) as introduced in Ref. [1]. Precise measurement of the BRDF is time consuming 1 Institute of Information Theory and Automation of the CAS, Pod Vodárenskou věží4, 185 08, Praha 8, Czech Republic. E-mail: filipj@utia.cas.cz ( ). Manuscript received: 2017-09-07; accepted: 2017-  due to the very high number of samples of the function that have to be taken. In this paper, we introduce a method for adaptive measurement of the BRDF that provides precise results even when using a limited number of samples. The method does not rely on a database of already measured materials and scales well, so it is applicable where precision is important. We build upon the method in Ref. [2] which measures the BRDF by one-dimensional slices, but substantially extend its contributions.
We note that the BRDF is a four-dimensional vector-valued function f r (θ i , θ v , ϕ i , ϕ v ) of the illumination direction ω i = [θ i , ϕ i ] and the viewing direction ω v = [θ v , ϕ v ]; it defines how light is reflected from a material (see Fig. 1). The isotropic BRDF is a three-dimensional simplification of the BRDF which can represent only a subclass of spatially homogeneous materials such as plastics or paints. In contrast, the full-dimensional BRDF is anisotropic. Anisotropic materials have variable reflectance when rotated around a surface normal, as commonly found for many real-world materials that contain directional elements such as fabric threads or wood grain.
The proposed method is based on our findings concerning the typical behavior of anisotropic BRDFs. We assume that each two-dimensional BRDF subspace with fixed θ i , θ v (see Fig. 2) and its important features can be captured by several  Visualization of a BRDF with highlighted axial (red), diagonal (blue), horizontal (green), and vertical (cyan) slices. Note that the BRDF is a four-dimensional continuous function and that its axes are orthogonal.
diagonal and anti-diagonal cross-sections (so-called slices), i.e., one-dimensional signals. Sampling the BRDF subspace only along the slices greatly reduces the number of samples needed (see Ref. [2]). Another saving is made by choosing only certain subspaces (e.g., by limiting values of θ i , θ v to multiples of 15 • : see Fig. 2(left)). Subsequently, elevation-dependent behavior is captured by another type of slice: see Fig. 2(right). The slices form a sparse 4D structure in BRDF space. Their ability to capture important features of the BRDF depends on their density. The more dense the sampling, the more precise reconstruction that can be achieved, but at the cost of more samples. In order to balance measurement time and reconstruction quality, we thus provide optimal configurations of such a sampling structure with respect to a chosen number of samples.
Although the structure itself substantially reduces the number of BRDF samples in comparison to dense regular sampling, we must reduce further the number of measured values to make the method practical. This is done by using sparse, adaptive sampling along individual slices.
Although there are approaches to efficient material appearance measurement using specific lighting systems, using, e.g., basis-illumination [3] or linear light sources [4], we focus on a sequential measurement process as typically used by goniometric devices having a single point-light source and a single accurate sensing device.
The main contribution of this paper is, to the best of our knowledge, the first adaptive method for precise measurement of the full 4D anisotropic BRDF which does not rely on a database of already measured materials.
Section 2 outlines our research in the context of related work, while Section 3 explains the principle of the BRDF slices. Section 4 proposes a heuristic algorithm for adaptive measurement of values along them. Missing data reconstruction from the slices is explained in Section 5 and optimal placement of the slices in BRDF space is analyzed in Section 6. Section 7 compares the proposed method with two interpolation methods that work with nearly-uniformly distributed samples, using simulated-and real-measurement scenarios. We also compare our method with a previous adaptive method [2]. Finally, Section 8 discusses advantages and limitations of the method and Section 9 concludes the paper.

Related work
The main objective of this paper is an approach for adaptive acquisition of the BRDF that uses a limited set of samples. Our methods of adaptive measurement and data interpolation are closely related, as our task is to reconstruct an unknown function well by optimal placement of novel samples based on previously measured values. This optimal placement depends on the chosen interpolation method. Good candidates are global interpolation methods [5] based on radial basis functions (RBF) or kriging. Although these methods might seem suited to solving the problem of adaptive sampling for measurement of material appearance, they have high computational demands which become intractable when the number of samples exceeds several thousands. Adaptive measurement of material appearance is investigated in Ref. [6], but only in two dimensions; the approach samples the reflectance field. An adaptive approach for image-based BRDF measurement is proposed in Ref. [7], with planning of viewing and illumination directions based on minimization of uncertainty in parameters of an analytical model. Nauyoks et al. [8] fit six isotropic BRDF models to the measured data. They iteratively add new samples by including illumination and viewing directions where the models predict different outputs, i.e., the models disagree. Some goniometers [9] use a two pass adaptive method, first making a full pass over the outgoing hemisphere, and then making a fine pass around the detected specular peak in the form of a spiral pattern.
Matusik et al. [10] represent isotropic BRDFs using a wavelet basis or a linear combination of BRDFs.
Similarly, Nöll et al. [11] represent the same BRDFs using basis functions; however, deviations of their reconstruction from the reference are approximated by a basis of correction functions. Nielsen et al. [12] present an approach to reconstructing isotropic BRDFs from basis functions using extremely sparse measurements, while Vavra and Filip [13] have extended their approach to anisotropic BRDFs. These approaches optimize a set of appropriate directional samples based on information in a database; however, they are not adaptive (except for Ref. [13]) according to the measured materials. As a result, they cannot precisely measure materials with features not present in the database.
None of the methods above provide a reliable way to effectively acquire dense BRDF data without knowledge of previously recorded BRDFs. On the other hand, several interpolation methods applicable to measured samples exist, based on, for instance, barycentric interpolation [14], thin plate splines [15], or a Gaussian pyramid in the so-called pullpush algorithm [16]. Recently, Ward et al. [9] reconstruct uniformly sampled BRDF measurements by interpolated radial basis functions using a masstransport solution.
Another option is displacement interpolation [17] based on solving the generalized masstransport optimization problem. It performs well, but its computational demands increase dramatically when the number of samples is in the order of thousands, and therefore factorization is required. An alternative approach to BRDF data interpolation is to use data factorization by Tucker decomposition [18]. This can represent both isotropic and anisotropic BRDFs with quality comparable to analytic BRDF models in reasonable reconstruction time.
As these interpolation methods differ in quality and speed, we have tested several of them and selected two as references for evaluation of our proposed approach based on adaptive measurement.
Our paper builds on Ref. [2] which uses slices to capture the main features of 2D subspaces of the BRDF. Values on the slices are measured adaptively to minimize the number required. Although individual subspaces can be captured very precisely, the whole BRDF is obtained by simple interpolation and is therefore rather inaccurate. In this paper, we introduce two new types of slices that enable us to capture all important features and to achieve better BRDF reconstruction quality. The adaptive sampling algorithm is enhanced to enable very efficient sampling of the slices based on a chosen count of samples. This allows control over the time needed by the measurement process, which is very important in industrial practice. We also consider the optimal placement of the slices, which is very important for achieving good results. This topic is completely overlooked in Ref. [2].

Introduction
The BRDF is a four-dimensional vector-valued function of the elevation angles θ i , θ v and the azimuthal angles ϕ i , ϕ v . To enable efficient sampling of the entire BRDF, we propose to use four types of one-dimensional slices. Two types, axial and diagonal slices, were introduced in Ref. [19], where only one slice of each type per 2D BRDF subspace is used. We propose to use up to dozens of slices per subspace to capture subtle details of the BRDF. Moreover, we suggest extending the concept to additional horizontal and vertical slices.

Axial and diagonal slices
Axial and diagonal slices are made in the 2D BRDF subspaces defined by fixed θ i and θ v . These slices are designed to optimize capture of specular reflections and anisotropic reflections. Axial slices s A (red in all figures) are perpendicular to anisotropic reflections, while diagonal slices s D (blue in all figures) are perpendicular to specular reflections (see Fig. 3): where α or β determines the position of the slice in the 2D subspace chosen by elevation angles θ i , θ v . Each subspace is typically measured using

Fig. 3
A schema of axial (red, perpendicular to anisotropic highlights) and diagonal (blue, perpendicular to specular highlights) slices placed in the 2D subspace (fixed both θi = 60 • and θv = 60 • ) and their function values. Note that the subspace is periodic. several axial and several diagonal slices to accurately capture all reflections and their shapes (see Fig. 3).

Horizontal and vertical slices
Horizontal and vertical slices are 1D subspaces of the BRDF. They are designed to capture changes in reflectance values when the viewing elevation angle θ v (horizontal slice, green) or the illumination elevation angle θ i (vertical slice, cyan) is changed while other parameters are fixed (see Fig. 2(right)): , so only, e.g., the horizontal slices need to be measured. All four types of slices together represent a sparse 4D structure in the four-dimensional BRDF space that effectively captures the main visual features of the measured BRDF (see Fig. 2(right)).

Adaptive sampling along the slices
Each slice can be interpreted as an unknown onedimensional signal that we need to measure and reconstruct. For axial and diagonal slices, the signal is periodic with a period of 360 • . To adaptively place samples, we use an enhanced version of the heuristic which was first introduced in Ref. [2].
Any slice can be sampled uniformly with a defined step (e.g., 1 • ), or adaptively, the latter decreasing the number of samples on one hand and increasing reconstruction accuracy in areas with high variance of the signal on the other. As the behavior of the signal is unknown, the adaptive algorithm must rely on already measured samples, adding new samples in areas where doing so improves the accuracy of the reconstructed signal. When samples are taken in appropriate directions and their count is sufficient, values for the entire slice can be interpolated sufficiently well using, e.g., piecewise cubic splines.
Our proposed heuristic in Algorithm 1 enables very efficient adaptive sampling of the slices using a given count of samples n. The algorithm is based on a simple assumption. If the value of a sample can be predicted well by neighboring samples, Algorithm 1 Algorithm for adaptive sampling along the slices Input: n, k, p1, p2. Output: T , a table of measured directions and their values.

1:
Measure values of all n0 samples at all intersections of the BRDF slices to create T .

2:
p1-times perform: the neighborhood of the sample can probably be predicted well too. Therefore, there is no need to place new samples there. On the other hand, if the value of any sample cannot be predicted well by the neighboring samples, it is possible that the neighborhood of the sample cannot be predicted well even by the sample itself together with its neighbors. Therefore, we should place new samples there.
The algorithm comprises three steps as shown in Fig. 4. In the first step, the algorithm samples the signal only in directions corresponding to intersections of slices, to collect initial information about the signals. Then, several iterations of adaptive measurement are performed relying on the ASI (adaptive sampling iteration) function, which adaptively selects and measures n i sampling candidates. See Algorithm 2.
The number of samples n i in each iteration depends on parameters k, p 1 , p 2 . The parameter k defines the ratio of the count of samples measured in the second step of Algorithm 1 to the count of samples taken in the third step. We have found experimentally that a reasonable choice is k = 0.9. The parameters p 1 and p 2 define how many iterations are performed in each step; our default choice is p 1 = p 2 = 5. Using these parameters, the algorithm forces division even of intervals on the slices where the error is lower in the second step, as there could be some hidden variation of the signal. In the third step, only intervals still having high cross-validation error are divided. Although the algorithm is not very

Fig. 4
Individual steps of the adaptive sampling algorithm, and the cross-validation procedure (right).

Algorithm 2 Adaptive sampling iteration function
Evaluate the leave-one-out cross-validation error for each direction in T (for each slice independently): interpolate the value of a sample in T based on neighboring samples and evaluate the error by comparison with the true value of the sample.

2:
Make a list of potential directions (see crosses in Fig. 4) ordered by descending weight. A potential direction lies at the center of each pair of neighboring already measured directions in T . Its weight is the maximal error at the neighboring directions.

3:
Measure values of the first n i directions in the list and append them to T .
sensitive to the parameters, one should not decrease the total number of iterations p 1 + p 2 too much, as intervals on the slices would then be insufficiently divided. However, using too many iterations (up to addition of only one new sample during an iteration) is not recommended either, as the algorithm would become too focused on the specular highlights and might overlook variations of the signal in unexplored intervals.
The ASI function requires the table T of already measured directions and their values and the count n i of samples to be identified and measured by the function. The output of the function is the appended table T . The function itself has three steps. First, the leave-one-out cross-validation error is evaluated for each already measured sample. The error is evaluated independently for each slice, so the sample at the intersection of the slices has several error values. Evaluation for each sample is done by excluding the evaluated sample from the dataset, linear interpolation of a value at the location of the evaluated sample from those of the neighboring samples, and computation of the difference between the interpolated value and the actual value of the evaluated sample (see Fig. 4(right)).
As a distance measure, we use the maximum difference over all color channels. Then, in the second step of the ASI function, a list of all directions where new measurements might be performed is prepared. Each direction is assigned a weight equal to the maximum of the cross-validation errors of its already measured neighbors (see Fig. 4(right)). The list is then sorted in descending manner, and values for the first n i directions from the list are measured. The table T is expanded by appending the newly measured directions and their values.

Reconstruction of BRDF
Once the measurement is done, we want to reconstruct values of the four-dimensional BRDF. In this section, we show how to rewrite equations introduced by Filip et al. [2] which serve for reconstruction of a two-dimensional BRDF subspace as a step in computing the four-dimensional BRDF space. Then, equations for interpolation of any value in the BRDF space from the sparse 4D structure (in the form of the four types of BRDF slices) can be easily derived.
We start with the same data as in Ref. [2]: the four values at the corners of the square formed by slice intersections cxȳ, wherex,ȳ ∈ {0, 1}, the two values p x0 , p x1 on the axial slices and the two values q 0y , q 1y on the diagonal slices. First, we apply bilinear interpolation of the corner values: (1) Values of the slices are interpolated linearly and differences are computed as ∆p xy = p xy − c xy ∆q xy = q xy − c xy The final value of the reconstructed function is r xy = max(c xy + ∆p xy + ∆q xy , min xy ) min xy = min(p x0 , p x1 , q 0y , q 1y ) Extension to four dimensions is straightforward and is shown in Section S4 in the Electronic Supplementary Material (ESM). We apply multilinear interpolation to all sixteen corners of the fourdimensional hyper-cube cxȳzw, wherex,ȳ,z,w ∈ {0, 1}. Then, values of the slices are interpolated in the three remaining dimensions and the final value of the reconstructed function is computed.

Basics
Although values along individual slices are measured adaptively, the positions of the slices in BRDF space must be determined in advance. Here we investigate how to place the slices optimally for a given number of samples. Due to a lack of available anisotropic BRDF measurements with high angular density and high accuracy, we performed a study on ten materials represented by the state-of-the-art analytical BRDF model of Kurt et al. [20]. We selected ten materials, four of them measured in Ref. [21] and fitted in Ref. [20], and the remaining materials came from our own measurements. All of these materials are anisotropic and include fabrics, brushed aluminum, and untreated wood. Their fitted BRDFs are shown in Fig. 5 and the model parameters are included in Section S2 in the ESM.

Placement theory
The axial and diagonal slices are placed into 2D subspaces selected by θ i and θ v . The most straightforward approach is to deploy the slices regularly into the subspaces using: α, β ∈ {ka s ; k ∈ {0, . . . , a n − 1} , a n = 2π/a s } where a n is the number of slices of one type in one subspace and a s is the azimuthal step between two consecutive slices. The single parameter a s has to be chosen. While the positions of the axial and diagonal slices might be chosen arbitrarily, for simplicity, we limit ourselves to regular distribution of both types of slices using a common step parameter. Similarly, 2D subspaces are selected by a single parameter, the elevation step e s : where e n = θ max /e s + 1 is the number of elevations used and θ max is the maximal elevation. Measured subspaces are those resulting from all combinations of θ i , θ v . Finally, horizontal and vertical slices are placed through every intersection of the axial and diagonal slices (see Fig. 2).

Analysis of optimal placement
We performed an experiment to determine the optimal positions of the slices based on a userchosen number of samples n. We used all ten BRDFs of the ten materials (Fig. 5)  To evaluate the quality of the reconstructed BRDFs, we computed the mean relative error (MRE) between the reference BRDF f r (θ i , θ v , ϕ i , ϕ v ) and its reconstruction f r (θ i , θ v , ϕ i , ϕ v ): is the number of data points and λ ∈ {R, G, B} is a color channel. Individual color channels are treated separately and the results summed. We express MRE as a percentage.
We computed graphs of MRE as a function of the number of adaptive samples n. This resulted in |a s | × |e s | × |m| = 5 × 8 × 10 = 400 graphs of err m,a s ,e s (n), where m stands for one of 10 materials. Each graph captures the reconstruction error for a large range of samples n. Figure 6 plots some of those graphs averaged across all the materials. The first group of graphs in Fig. 6(a) shows the error on fixing  the parameter e s = 14 • and varying the parameter a s . The second group of graphs in Fig. 6(b) shows the error when fixing the parameter a s = 20 • and varying the parameter e s . The fast convergence of the individual graphs confirms the efficiency of the adaptive sampling algorithm. From the graphs, it is possible to conclusively select the best combination of the parameters a s and e s for a selected number of samples n, for a given material m. Unfortunately, these values differ according to material, as each exhibits individual behavior with respect to changes in the azimuthal and elevation angles. Therefore, we select optimal values a s , e s of these parameters for a given number of samples n so that the sum of errors across all materials relative to the achievable error is minimized using The resulting optimal values for a s , e s are summarized in Table 1; when a new material is measured, one should select appropriate values of the parameters for the chosen number of samples.

Results
This section presents results of using the proposed method. Firstly, due to lack of reliable data, we use synthetic data generated by the BRDF model (see Fig. 5). Thus, we can easily and quickly obtain a BRDF value in any direction and results of the experiment are not influenced by errors caused by a measurement process. To evaluate the method on real data, we used a 3D scene (see Fig. 11) and measured all data needed to visualize the scene using a gonioreflectometer. Finally, we compare our method with the previous adaptive method in Ref. [2].

Simulated measurement experiment
We evaluated accuracy of the proposed method in comparison to use of uniformly distributed samples, taken in directions according to one of the thirty sampling schemes we designed (see Fig. 7 and Section S1 of the ESM). These schemes produced in total from n = 435 to n = 354,061 reciprocal samples. Values of the samples were interpolated using barycentric [14] or RBF [5] interpolation. Note that the second method is global while the first is local, and is therefore suitable for fast rendering on a GPU. Both methods compute results for each color channel separately. We interpolated the BRDF to a four-dimensional array using a uniform step of 2 • and a maximal elevation of 80 • , so the dimensions of the array were |θ i |×|θ v |×|ϕ i |×|ϕ v | = 41×41×180×180. Note that the error values can be compared only for the given number of samples as the barycentric and RBF interpolation methods operated on the 30 predefined sampling schemes. The average MRE (over all materials and schemes) for barycentric interpolation was almost 7.5 times as large as for the BRDF slices. RBF interpolation was better, but its average MRE was still almost 3.9 times as large as for the BRDF slices.  Figure 8 shows reconstruction error as a function of the number of samples for all tested methods, for two materials (for other materials see Section S3 of the ESM). Convergence of the proposed algorithm to low MRE values as the number of samples increases is very fast (see the green line in Fig. 8). While the RBF method performs well for lower number of samples, the proposed method has better accuracy for higher number of samples. With more than 5000 samples, our method achieves high quality results that are significantly better than those achieved using uniform interpolation methods.
To evaluate quality of the reconstructed BRDF, we computed the MRE between the reference BRDF and its reconstruction at N = 3 × 41 × 41 × 180 × 180 =163,393,200 data points. Improvement in the MRE due to usage of the proposed method instead of barycentric or RBF interpolation for sampling schemes 14 (8911 samples) and 19 (18,721 samples) is shown in Fig. 9. In Fig. 10, we compare all three methods with a reference rendering for three materials and 8911 reciprocal samples. We used the grace environment represented by means of 256 lights. (Results for additional materials are shown in Section S3 in the ESM.) The results achieved show that our method provides the best reproduction of the specular reflections and the anisotropic highlights.

Practical measurement experiment
The previous experiment was performed using reference BRDFs represented by an analytical model. However, such an approach sacrifices some visual features of the original reflectance behavior that cannot be reliably represented by the model. Therefore, we performed practical BRDF measurement experiments, collecting 8911 samples. Firstly, we recorded these samples uniformly (using sampling scheme 14, see Fig. 7 and Section S1 of the ESM) and interpolated them using barycentric and RBF methods. Then, we recorded the same number of samples adaptively using our proposed method. Altogether, we measured two datasets for each of three materials (fabric112, fabric135, fabric136) publicly available at http://btf.utia.cas.cz.    All measurements were performed using the UTIA gonioreflectometer [2] that allows for the placement of an almost arbitrary combination of illumination and viewing directions with high angular accuracy. Acquisition time for the adaptive method was about 10% longer than for uniform sampling of the same number of samples due to the data processing overhead required by the adaptive algorithm. It took approximately 25 hours to acquire 8911 samples.
We compared the data in an applied situation using a 3D scene comprised of four spheres illuminated by a single point-light source as shown in Fig. 11. Therefore, many combinations of illumination and view directions were covered providing a comprehensive visualization of the properties of these materials. The rendered images were divided into a sparse raster with only 6195 occupied pixels representing the directions reachable by the gonioreflectometer; therefore, only this number of BRDF values was measured for the three materials. The entire scene was then rendered using those pixels, which we call control samples; we used them as our reference. Note that pixels representing directions unreachable by the gonioreflectometer due to occlusion of view of the camera by the light source are not included into the control samples (resulting in a white spot in the difference images in Fig. 12, on the leftmost sphere). Their values were interpolated from regular measurements for purposes of visualization.
The insufficient precision of analytical models, and also of non-adaptive measurement, is shown by our experiment in Fig. 12. Renderings using tabulated data (rows 2, 3, and 4) outperform analytical models (rows 5 and 6) in terms of precision. Therefore, tabulated functions should be used whenever high accuracy is important.
We also rendered the same scene using values obtained by the three methods to be compared and evaluated their results at all control samples. Figure 12 shows the reference scene side-by-side with its reconstruction using all tested methods as well as results of two tested analytical anisotropic models due to Kurt et al. [20] and Ward [22]. From the accompanying error values, we conclude that our method has by far the best accuracy for real BRDF data.
In addition, we show how our method compares  12 All reachable directions in the virtual scene (see Fig. 11) measured and compared with barycentric and RBF interpolations of uniformly measured data, fits of two analytical models ( [20,22]), and our proposed adaptive measurement. 8911 samples were used. Difference values are given as CIE ∆E / RMSE / PSNR (dB). The last row shows results of dense BRDF measurement [23].
to an approach to dense BRDF measurements in Ref. [23]. It relies on spatial homogeneity of measured material and slight variations of illumination and viewing angles within each captured image. A total of 8505 images capture nearly three million BRDF samples. After interpolation of missing samples, the data are uniformly resampled with a sampling step of 2 • in all four dimensions; these data currently represent one of the best anisotropic BRDF datasets available. The last row of Fig. 12 uses such data performance for our test scene, giving difference images and CIE ∆E, RMSE, and PSNR values.
When compared to the proposed method (using 8911 samples) in the previous row, one can observe the comparable accuracy of both methods, but while the sequential measurement time for over 8000 images is comparable, the proposed method does not assume largely homogeneous material with the associated challenges related to sample preparation.

Comparison with a previous adaptive method
We have compared our method with the method of Filip et al. [2] using synthetic data and our own reimplementation of that method. Its results were much worse than ours when we use 12 axial and 12 diagonal slices per subspace as suggested in Ref. [2]. Moreover, its parameter t must be tuned carefully to achieve the chosen number of samples, which would not be possible in a real measurement scenario. Thus, we have also enhanced the method with our approach to adaptive sampling (see Section 4) which enables us to choose the number of samples easily. Also, we have taken our conclusions on the optimal placement of samples (see Section 6) into account, which helps greatly. This allows us to evaluate the additional contribution due to the horizontal and vertical slices over and above that of the axial and diagonal slices. Results of the comparison are summarized in Fig. 13. Our method on average has only about half the error of even the enhanced version of the method of Filip et al. [2].

Validation on specular anisotropic materials
Finally, we validated our method on several highly specular anisotropic materials. As such ground truth measurements are not available with sufficient angular density, we have to rely on an anisotropic BRDF model again. We selected three highly reflective self-adhesive aluminum foils with directional finishing. We captured their BRDFs uniformly using 81 samples per hemisphere [24] and fitted this data using an anisotropic BRDF model [20]. Photographs of the foils and the captured BRDFs are shown in Fig. 14.
The fitted BRDF model was used to generate dense ground truth data for testing of the compared methods. Figure 15, last column, shows that the proposed method clearly performs better than uniform sampling with barycentric or RBF interpolation. One can observe the contribution of  the proposed method in comparison to its variant without horizontal and vertical slices [2].
Similar results are obtained for the two sampling schemes when the MREs are computed between the ground truth and the values reconstructed by the tested methods, as shown in Table 2. Also see additional results in Section S6 of the ESM. of the proposed method, and gives timings for the interpolation methods used. Furthermore, we discuss reliability of the proposed approach that finds the optimal placement of the slices.

Advantages
In contrast to competitive methods, the main advantage of our proposed adaptive sampling approach is its faster decrease in reconstruction errors with increasing number of samples (especially for lower counts, see Fig. 6). On average, we achieve almost four times lower relative errors given the same number of samples, and between two to five times fewer samples are necessary to achieve the same reconstruction errors. Additionally, the proposed adaptive sampling method operates along the onedimensional slices in the BRDF space, allowing data acquisition by continuous movement of a light source and a sensor. It is especially beneficial for accurate gonioreflectometer-based measurement devices.

Limitations
Common to all tested methods is their decreasing improvement in accuracy as increasing number of samples, approaching an asymptotic value. One reason for this behavior is lack of samples at locations representing specular and anisotropic highlights. They are not sampled properly by the uniform sampling used in the RBF and barycentric methods, nor by the uniform positioning of the slices in BRDF subspaces by the proposed method. We believe that the development of a better parameterization in combination with an adaptive placement of the slices would further improve performance. As the proposed method is particularly designed to represent features perpendicular to the slices, it can suffer from improper representation of curved highlights or visual features occurring between the slices. Such behavior is typical for subspaces with a large difference between viewing and illumination elevation angles, as shown in Fig. 16. These errors can be suppressed by an elevation-angle dependent interpolation method [25] that interpolates data along the anisotropic highlights. Unfortunately, this method is not yet elaborate enough to reconstruct the whole 4D BRDF space. a uniform step of 2 • requires about 45 seconds using barycentric interpolation, whereas RBF interpolation takes about 22 minutes; neither depends on the number of samples. Reconstruction of the whole BRDF array by the proposed method using our MATLAB implementation typically takes 4 minutes. All timings were obtained using a single core of an Intel Xeon E5-2643 at 3.3 GHz. Implementation of the method for interpolation of an arbitrary value on a GPU is described in Section S5 of the ESM.

Stability of placement of slices
To test stability of the proposed placement optimization, we performed the following experiment. We computed the optimal values of the parameters m times using Eq. (3), each time with one material left out of the computation. This provides m different results, which we graphically compare in Fig. 17 to the results in Table 1. The first line in Fig. 17 Fig. 17 Top row: The color coded optimal placement of slices in terms of es/as for various numbers of samples. Other rows: Stability test of optimal placement of slices: in each row, one material is left out from computation of the optimal placement.
As one can see, boundaries of individual regions (representing combinations of parameters e s /a s ) vary slightly but not dramatically. When we leave out Red velvet, there is a green region (14/20) which does not appear on other lines substituting purple (10/36) and blue (8/36) regions. When we leave Brushed alum, fabric041, or wood01 out, the blue region (8/36) is missing. Note that usage of the parameters for neighboring regions leads to only slightly worse reconstruction errors.
We evaluated the increase in MRE due to usage of the optimal placement of the BRDF slices according to Table 1, instead of usage of the best placement for the tested material (which is unknown for a newly measured material).
Evaluation was performed across a broad range from 55 samples (resulting from the number of slice intersections in the sparsest configuration e s = 28 • , a s = 180 • ) to 354, 061 samples (resulting from the densest tested uniform distribution of samples across a hemisphere, see Fig. 7). Results are summarized in Table 3. While the maximal increase in the error is almost 1.5 times, the average MRE across all materials is only 1.2 times worse than if we use the best possible placement of the slices. We conclude that usage of the optimal placement according to Eq. (3) works well in practice for measurement of an unknown material.

Conclusions
Our paper deals with efficient sampling and reconstruction of anisotropic BRDFs using a chosen number of samples. Firstly, we sample a fourdimensional dataspace by means of one-dimensional slices, with adaptive density of samples along the slices according to material properties. To find the optimal placement of the slices in the 4D BRDF space, we perform an extensive study on ten materials, and we provide parameters for optimal placement based on the required number of measured samples.
We have compared the accuracy of our method with barycentric and RBF interpolation approaches as applied to 30 uniform sampling schemes with varying numbers of samples. In simulated measurements that use an analytical anisotropic BRDF model as a source of reference data, the proposed method achieves on average almost fourtimes lower reconstruction errors than uniform sampling combined with either of two interpolation methods. Alternatively, for a given reconstruction error, our method requires between two to five times fewer samples than the competing approaches.
In real measurement of the anisotropic BRDFs, our method achieves half the reconstruction error of other approaches, and our further experiments suggest that the captured BRDF data are amongst the best publicly available anisotropic BRDFs. We have also verified that the proposed method performs almost twice as well as a previous adaptive method [2], mainly due to its use of two additional types of BRDF slices. Moreover, the proposed method readily allows reconstruction of a BRDF on a GPU.
In our future work, we plan to enhance the interpolation method proposed in Ref. [25] to eliminate the main disadvantage of our method, reconstruction artifacts. Also, we plan to extend our method to a more time-and resource-efficient measurement of spatially-varying representations such as SVBRDFs or BTFs. schemes, fitted parameters of BRDF models, and addition results for all tested materials are available in the online version of this article at Open Access The articles published in this journal are distributed under the terms of the Creative Commons Attribution 4.0 International License (http:// creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.
Other papers from this open access journal are available free of charge from http://www.springer.com/journal/41095. To submit a manuscript, please go to https://www. editorialmanager.com/cvmj.