Morgane Gerardin
email: gerardin.morgane@gmail.com

Pauline Martinetto

Nicolas Holzschuch

Dry powders reflectance model based on enhanced back-scattering: case of hematite α -Fe 2 O 3

By performing Bidirectionnal Reflectance Distribution Function (BRDF) measurements, we have identified back-scattering as the main phenomena involved in the appearance of dry nano-crystallized powders. We introduce an analytical and physically-based BRDF model that relies on the enhanced back-scattering theory to accurately reproduce BRDF measurements.

These experimental data were performed on optically thick layers of dry powders with various grains' morphologies. Our results are significantly better than the ones obtained with previous models. Our model has been validated against the BRDF measurements of multiple synthesized nano-crystallized and mono-disperse α -Fe 2 O 3 hematite powders. Finally, we discuss the ability of our model to be extended to other materials or more complex powder morphologies.

Introduction

Accurate reproduction of powdered material appearance is an important issue in many domains, such as environmental sciences, teledetection, digital prototyping, entertainment and art preservation. Powder materials are ubiquitous in real life. They can be used with a binder like pigments in most painting techniques and in some cosmetics e.g. lipstick. But they are also used without any binder. Dry powders are found as spices, dry soils, dust, and dry pigments for example.

The latter are particularly important as they are involved in parietal paintings from paleolithic to neolithic era. In some cases, the pigments were used pure and dry on the wall. Other application techniques required the pigments to be dispersed in an organic binder [START_REF] Vignaud | Le groupe des «bisons adossés» de lascaux. étude de la technique de l'artiste par analyse des pigments[END_REF][START_REF] Casoli | Research on the organic binders in archaeological wall paintings[END_REF]. However, only few traces of these binders are left nowadays, and they are not the main component responsible for the painting color. The parietal paintings that are visible nowadays can therefore be described as a dry pigment layer on top of a rock substrate. These cultural heritage artefacts have been degraded over the years by both environmental and anthropogenic causes, so that some of them are nowadays almost illegible [START_REF] Vignaud | Le groupe des «bisons adossés» de lascaux. étude de la technique de l'artiste par analyse des pigments[END_REF]. Images of such heavily deteriorated artefacts are recorded using image enhancement techniques [START_REF] Brady | Documenting and analyzing rock, paintings from torres strait, ne australia, with digital photography and computer image enhancement[END_REF][START_REF] Domingo | Evaluating conventional and advanced visible image enhancement solutions to produce digital tracings at el carche rock art shelter[END_REF]. Time, weathering and human activity are responsible for the peeling off, the erasure and possible color changing of the original paintings.

Computer Graphics is a real asset for the preservation of parietal paintings. Photo-realistic image renderings using accurate appearance models would allow the reproduction of such artefacts and the reconstruction of the most erased ones. It will also promote their broadcasting to a large audience.

The colour of a powder material depends on both its chemical composition and on its morphometry, which refers to the size and the shape of its grains. This hypothesis is well known and particularly well justified in the case of pure α -Fe 2 O 3 hematite powders that can show various shades from light red to purple [START_REF] Chalmin | Les blasons de lascaux[END_REF][START_REF] Gerardin | Influence of micro-structural features on the colour of nanocrystallised powders of hematite and visible-nir reflectance spectra simulations[END_REF]. The chemical composition of the material is involved in the absorption phenomena, while the morphology of the grains influences the way light is scattered by the material. Such light/matter interactions at the nanometric scale result in the macroscopic material appearance an observer perceives.

At the macroscopic scale, light scattering by a material is modeled using a Bidirectional Scattering Distribution Function (BSDF) model. When there is no light transmission and only light scattered in reflection is at stake, the BSDF is reduced to Bidirectional Reflectance Distribution Function (BRDF). Most models assume that the light is mainly scattered at the surface of an object, around the specular direction ie. in the direction opposite to the incident direction. Depending on the roughness of the surface, a diffuse look is obtained. According to these surface models, the color of the rendered materials only depends on their chemical composition through the optical constants. These models are not suitable in the case of dry powders for two reasons. First, powders are a medium composed of a collection of randomly oriented scatterers, so light is scattered multiple times in the volume of the medium. To account for volume scattering, one may use the Mie scattering theory [START_REF] Bohren | Absorption and scattering of light by small particles[END_REF][START_REF] Frisvad | Computing the scattering properties of participating media using lorenz-mie theory[END_REF]. This is relevant for powders diluted into a binder, and takes into account volume scattering depending on the size and shape of the grains. However, in the absence of a binder, the grains composing the powder are too close to each other, making the Mie scattering model unsuitable as shown by [START_REF] Gerardin | Influence of micro-structural features on the colour of nanocrystallised powders of hematite and visible-nir reflectance spectra simulations[END_REF]. The second issue is that only few of surface scattering models account for scattering in the backward direction [START_REF] Belcour | Bidirectional reflectance distribution function measurements and analysis of retroreflective materials[END_REF].

In this paper, we show that contrary to the prediction of most BRDF models, the major contribution to a dry powder BRDF is not the specular peak, but the back-scattering peak. We propose an analytical and physically-based BRDF model involving six parameters able to simulate the appearance of optically thick layers of dried pigments. We separate the model into two independent contributions: volume scattering and surface scattering. The former decomposes into single and multiple scattering that are responsible for the sharp back-scattering lobe, while the latter accounts for the specular scattering. This model has been developed and experimentally validated for seven hematite α -Fe 2 O 3 powders of different morphologies. Although our model still tends to slightly underestimate the back-scattering contribution to the BRDF, it gives a good intuition on the color variations of the studied samples. It is ready-to-use to reproduce the color of mutliple studied hematite powders through renderings.

In Section 2 we will review previous works on BRDF models suitable for the description of dry powder materials, and more particularly on models that take into account the back-scattering contribution. The theoretical background necessary for the modeling of light scattering in a particle medium is introduced in Section 3. The α -Fe 2 O 3 samples and their measurements are described in Section 4. The description of our appearance model and its optimization procedure is done on Section 5 whereas in Section 6 we present its results, showing that it fits well the experimental BRDF of the studied samples. We conclude and give some future work directions in Section 7.

Previous work

Dry powders and pigments are complex materials composed of a dense collection of particles with nanometric to micrometric grains size. In the absence of a binder, the dry grains are too close to each other to be modeled using Mie scattering theory [START_REF] Gerardin | Influence of micro-structural features on the colour of nanocrystallised powders of hematite and visible-nir reflectance spectra simulations[END_REF]. Their appearance is the result of absorption and scattering phenomena, so its description suggests the full understanding of the interactions of light with the material, at both its surface and in its volume.

We have identified two models that could be used for dry powders' BRDF simulations: the spectralon model and the RTLSR model. Firstly, they both take into account the back-scattering, a phenomena which has already been observed on dry powders such as meteorite soils [START_REF] Beck | Photometry of meteorites[END_REF].

Then, although these models have not been used for reproducing the appearance of dry powders or pigments, they have been introduced for describing particulate media without binder. These two reasons make them worthy to study in the case of dry powders.

An empirical model has been developed by [START_REF] Lévesque | Measurement and modeling of the spectralon spectro-polarimetric bidirectional reflectance distribution function (brdf)[END_REF] to reproduce the BRDF of the Spectralon®, a material that resembles a powder since it is made out of a collection of compressed polymer particles. This model proposes a decomposition into four independent scattering phenomena.

These contributions are by order of importance: Deep scattering D, which accounts for volume scattering at small incidence and observation angles, Back-scattering B and Specular scattering S, to account for the backward and forward surface scattering, and Subsurface scattering SS, which accounts for volume scattering at more grazing incidence/observation geometries. Each of these contributions is modeled using very simple analytical expression.

BRDF(θ i , θ o ) = D(θ i , θ o ) + B(θ i , θ o ) + S(θ i , θ o ) + SS(θ i , θ o ) (1)
One should note that this model has been developed in the incident plane. To model the complete BRDF of a material some extension is necessary to account for the different azimuth angles. We propose such an adaptation in Supplementary work.

Volume scattering contributions: For small incident angles θ i , there is a higher probability for the light to penetrate deeply into the medium before it exits with a direction close to the surface normal, with a small observation angle θ o . In this case, the entry and the exit points of the light are quite close. Such a phenomenon is called Deep scattering D and is identified as the main contribution to the BRDF. Its modeling in the incident plane is inspired from an empirical expression from [START_REF] Walthall | Simple equation to approximate the bidirectional reflectance from vegetative canopies and bare soil surfaces[END_REF] and involves five material dependent fitting parameters d j :

D(θ i , θ o ) = d 0 + d 1 (θ 2 i + θ 2 o ) + d 3 (θ 4 i + θ 4 o ) + d 3 θ 2 i θ 2 o + d 4 θ 4 i θ 4 o (2)
Deep scattering results in a large scattering lobe that describes the general shape of the measured BRDF.

For larger incident angles, light penetrates much less deeply and propagates right underneath the medium surface to an exit further away from the entry point which is modeled by the Subsurface scattering SS. Two models were proposed by [START_REF] Lévesque | Measurement and modeling of the spectralon spectro-polarimetric bidirectional reflectance distribution function (brdf)[END_REF] involving an exponential term so that the more grazing the angle, the more important this contribution is. For simplicity, we only give here the expression of the simplest, which involves two empirical parameters ss 0 and ss 1 :

SS(θ i , θ o ) = ss 0 exp(ss 1 |θ i + θ o |) (3)
Surface scattering contributions: The surface scattering phenomena are responsible for the finer details of the measured BRDF. The Back scattering is modeled using a Lorentzian function with amplitude b 0 and width b 1 , whereas the Specular scattering is described as a Gaussian function with width s 1 and whose amplitude depends on two other parameters s 0 and s 2 .

As the Spectralon model involves a large number of empirical but material-dependent parameters optimized in the specific case of the Spectralon material, it can be adapted to reproduce with a high fidelity a large range of materials. Nevertheless, these parameters have no explicit relationship with the chemical composition nor the physical parameters of the considered material so they cannot be predicted, but only numerically optimized from measurement fittings.

RTLSR

The RTLSR (Ross-Thick Li-Sparse Reflectance) model [START_REF] Roujean | A bidirectional reflectance model of the earth's surface for the correction of remote sensing data[END_REF][START_REF] Lucht | An algorithm for the retrieval of albedo from space using semiempirical brdf models[END_REF], is a semi-empirical model used in teledetection for simulating the BRDF of forest canopies. It decomposes the BRDF into two contributions: volume scattering contribution R vol deduced from an approximation of the radiative transfer theory (RTT) [START_REF] Ross | The radiation regime and architecture of plant stands. the hague: Dr. w[END_REF] and accounting for both single and multiple scattering between the leaves, and surface scattering R geo which is inferred according to geometrical considerations [START_REF] Li | Geometric-optical bidirectional reflectance modeling of a conifer forest canopy[END_REF] and accounts for shadowing and occlusions of leaves on each others. An optimization factor α allows to control the proportion of each of the contributions.

BRDF(θ i , θ o , ϕ i , ϕ o ) = αR vol (θ i , θ o , ϕ i , ϕ o ) + (1 -α)R geo (θ i , θ o , ϕ i , ϕ o ) (4)
Initially developed for thick and sparse canopies, extensions of this model have been proposed

for thin and dense media [START_REF] Wanner | On the derivation of kernels for kernel-driven models of bidirectional reflectance[END_REF]. [START_REF] Jiao | Development of a snow kernel to better model the anisotropic reflectance of pure snow in a kernel-driven brdf model framework[END_REF] brought modifications to the model for the study of snow, and multiple works [START_REF] Li | Geometric-optical bidirectional reflectance modeling of a conifer forest canopy[END_REF][START_REF] Nilson | A reflectance model for the homogeneous plant canopy and its inversion[END_REF][START_REF] Ni | An analytical hybrid gort model for bidirectional reflectance over discontinuous plant canopies[END_REF] agree that it can be used for bare soils although they do not explicitly introduce the changes that should be done in such a case. We propose an adaptation of this model for dry powders in Supplementary work. Its unique optimization factor makes this model too rigid and remains unable to reproduce fine appearance variations.

Background

As a dry powder layer involves a large number of scatterers with relatively high density, it is complex to describe the scattering of light in such a medium. When there is no interference phenomena at stake, it can be done using the RTT [START_REF] Chandrasekhar | Radiative transfer[END_REF][START_REF] Ishimaru | Wave propagation and scattering in random media[END_REF]. As it relies on the assumption that the scattered fields are uncorrelated, RTT by itself is unable to predict the appearance of materials which exhibits a strong back-scattering contribution.

A deeper analysis of the back-scattering phenomena (see Section 3.5) occurring in a scattering medium refers to Enhanced Back-Scattering (EBS), which describes it as the result of constructive interferences [START_REF] Tsang | Scattering of electromagnetic waves: advanced topics[END_REF][START_REF] Van Der | Light scattering in strongly scattering media: multiple scattering and weak localization[END_REF]. Contrary to the RTT, the scattered fields are supposed to be correlated. The EBS is therefore the sum of the two terms from RTT (single and multiple scattering) and a third one accounting for the constructive interferences in the backward direction. As the analytical results are deduced within the diffusion approximation, only isotropic particles are considered.

Nevertheless, Tsang & Kong [START_REF] Tsang | Scattering of electromagnetic waves: advanced topics[END_REF] extended the model to account for anisotropic scattering in the case of incident field at normal incidence.

This section gives an overview of the theoretical aspects and hypotheses used in the EBS theory.

Medium description

Let's consider a region 0 free of scatterers with wavenumber k 0 = 2π/λ, and a medium 1 characterized by an effective complex wavenumber K = K + iK containing scattering spherical inclusions with a density D. All scatterers have the same radius b and scattering properties (cf. Figure 1). Note that for clearer reading, we adopt bold notations for vectors representation

(ì v = v) in the following.
The medium is illuminated in the incoming direction K i with an angle θ i to the surface normal and azimuthal angle ϕ i by light coming from a point source at the position r i in the region 0.

The light scattered by the slab is observed at a point in region 0 in the direction K o , described by the angle θ o to the normal, and ϕ = |ϕ i -ϕ o | in the azimuthal plane at the position r o which is far enough from the interface (||r o || 1). To ease the description of the scattering inside the layer, a point scatterer approximation is made: the scattering events occur at the center point of the spherical scatterers.

Foldy's approximation

For a sufficiently small density D of scatterers, the effective wavenumber K can be expressed using the Foldy's approximation [START_REF] Foldy | The multiple scattering of waves. i. general theory of isotropic scattering by randomly distributed scatterers[END_REF]:

K ≈ 2π λ + jλDF (5)
where j 2 = -1, F is the complex scattering amplitude describing the scattering properties of the particles, and λ is the wavelength at which the light propagates.

θ i θ i θ o θ o ⃗ K i ⃗ K i ⃗ K o ⃗ K o r 1 r 2 r 3 r N -2 r N -1 r N 0, k 0 1, K x z y ψ(r i ) G 11 (r n , r m ) G 01 (r n , r o ) r i r o L Fig. 1
. N-th order multiple scattering in a medium composed of particles. The incident light follows a path going from the particle at r 1 to the one at r N (black path). Interferences occur with a second path that goes from r N to r 1 (gray path). The incoherent part of the resulting field (denoted as L N in Section 5.1) contributes to the ladder term whereas its coherent part (denoted as C N in Section 5.1) accounts for the cyclical term.

Light propagation

The propagation of light inside a medium containing scattering particles is illustrated by the Figure 1. The light wave coming from region 0 without being scattered is described by the incident field :

ψ(r i ) = exp(jK i • r i ) (6)
where K i is of magnitude K given by the Foldy's approximation, and oriented in the incident direction.

After it has been scattered by a particle, light travels in the region 1 towards another scatterer, as if it was in a space free of scatterers. Under this hypothesis, the propagation of light from a position r towards r is given by the Green's function G(r, r ), that satisfies the Helmholtz equation :

(∆ + k 2 )G(r, r ) = δ(r -r ) (7)
By averaging over the particles, the mean Green's function referring to the propagation of light from an inclusion to another is given as follows:

G 11 (r n , r m ) = exp(jK • |r n -r m |) 4π|r n -r m | (8)
where r m and r n are the positions of two distinct scatterers in the medium 1. In the same way,

G 01 (r n , r o ) = exp(jk 0 r o ) 4πr o exp(-jK o • r n ) (9)
describes the propagation from a scatterer at position r n in region 1 towards the camera at r o in region 0. The vector K o is of magnitude K and oriented towards the camera.

Light scattering by one particle

The scattering events caused by the encountering of light with the n-th scatterer is described using the transition operator T n that we introduced as :

T n (r n ) = 4πF(Θ n )δ(r in -r on ) = 4πF(Θ n ) ( 10
)
where Θ n = θ on -θ in is the scattering angle between the incident direction on the n-th scatterer θ in and the outgoing direction θ on .

The complex scattering amplitude F(Θ n ) is linked to the scattering phase function p(Θ n ) using the relationship:

ρp(Θ n ) = 4πN κ |F(Θ n )| 2 ( 11
)
where ρ is the single scattering albedo and κ the extinction rate.

Although any phase function can be used, a common one would be the Henyey-Greenstein phase function [START_REF] Henyey | Diffuse radiation in the galaxy[END_REF] (Eq. 12). It involves one parameter g known as the anisotropy and whose value lies between -1 and 1 describing either backward or forward scattering. Isotropic scattering would then be described with g = 0.

p(Θ) = 1 4π 1 -g 2 (1 + g 2 -2gcosΘ) 3/2 (12)

Enhanced back-scattering

The EBS theory offers a physical description of the back scattering phenomenon observed in dense media. It is an extension of RTT [START_REF] Chandrasekhar | Radiative transfer[END_REF][START_REF] Ishimaru | Wave propagation and scattering in random media[END_REF]: the so called ladder terms accounting for both single scattering BRDF ebs, L 1 and some multiple scattering BRDF ebs, L 2+ take into account the contribution of uncorrelated scattered fields already described by the RTT. This model also accounts for coherent scattered fields that induce interference effects as light undergoes multiple scattering events in the medium. It is described by a cyclical term BRDF ebs, C and results in a sharp peak in the back-scattering direction.

The BRDF is then composed of three contributions. Since the angular and spectral dependency of these three components is the same than the overall BRDF (θ i , θ o , ϕ i , ϕ o , λ), it is not reminded in the following equations for readability purposes.

BRDF ebs (θ i , θ o , ϕ i , ϕ o , λ) = BRDF ebs, L 1 (κ, g) + BRDF ebs, L 2+ (κ) + BRDF ebs, C (κ) (13)
This model has been developed for thick layers of a medium composed of spherical scatterers of radius b embedded in a host medium [START_REF] Tsang | Scattering of electromagnetic waves: advanced topics[END_REF][START_REF] Van Der | Light scattering in strongly scattering media: multiple scattering and weak localization[END_REF][START_REF] Bret | Multiple light scattering in porous gallium phosphide[END_REF][START_REF] Lönn | Utilizing enhanced backscattering for determination of scattering properties in turbid media[END_REF].

Analytical expressions of the multiple scattering terms are obtained by considering that for optically thick media with high albedo the resulting overall scattering statistically tends to be isotropic. So the phase function used in multiple scattering is supposed to be constant p(Θ) = 1 4π .

This is known as the diffusion approximation and leads to results in Eqs. 14.

                         BRDF ebs, L 1 (κ, g) = ρp(Θ) µ o + µ i 1 -exp -κd 1 µ o + 1 µ i BRDF ebs, L 2+ (κ) = 3ρ 2 1 -µ 2 i c 2 2Q 1 -µ o c - µ 2 i µ i + µ o BRDF ebs, C (κ) = 3ρ 2 [c + v(1 -exp(-2cτ o ))] 2µ i µ o cv[(c + v) 2 + u 2 ] (14a) (14b) (14c)
with µ i = cosθ i , µ o = cosθ o , d the thickness of the powder layer, and τ o the extinction ratio whose expression is given in Eq. 16.

Multiple parameters involved in Eqs 14 are introduced to make the equation more legible and depend on the geometry under which the medium is illuminated and observed. These parameters have been introduced by [START_REF] Van Der | Light scattering in strongly scattering media: multiple scattering and weak localization[END_REF] for isotropic scatterers such as:

                                   u = k 0 κ (µ i -µ o ) v = 1 2 1 µ i + 1 µ o α = k 0 κ f c = 3(1 -ρ) + α 2 Q = [1 + µ i c -(1 -µ i c)exp(-2cτ o (λ))] 4c (15a) (15b) (15c) (15d) ( 15e
)
where f is the norm of the vector f defined as difference between the projections of the normalized incoming and outgoing direction vectors (Figure 2) which allows the BRDF modelization in any azimuthal plane. Only the single scattering term involves the phase function in such a case, and therefore the anisotropy parameter g if the Henyey-Greenstein phase function is at use (see Eq. 12). Tsang & Kong [START_REF] Tsang | Scattering of electromagnetic waves: advanced topics[END_REF] introduced the consideration of anisotropic scatterers but only with normal incidence.

We will discuss the extension of the model for anisotropic scattering in the Section 5. The last parameter involved in the BRDF expressions (Eqs. 14) is the extinction ratio τ o that is defined as intrinsic to the scatterers' material. It accounts for internal reflections that may occur between the interface of the medium. A definition is given by [START_REF] Zhu | Internal reflection of diffusive light in random media[END_REF]:

τ o (λ) = 2 3 1 -R(λ) 1 + R(λ) (16) with R(λ) = 3C 2 + 2C 1 3C 2 -2C 1 + 2 and C i = ∫ π 2 0 R(θ, λ) cos i θ sin θdθ
In the previous equations, R(θ) is the Fresnel's angular reflection coefficient [START_REF] Born | Principles of optics: electromagnetic theory of propagation, interference and diffraction of light[END_REF]. The parameter τ o thereby depends on the optical constants of the material the particles are made of and on the wavelength.

Sample description and BRDF measurements

The macroscopic physical phenomena at stake in the appearance of a material can be inferred through BRDF measurements.

Sample description

For the study of dry powder appearance, we focused on synthetic nano-crystallized hematite α -Fe 2 O 3 powders with various colors, that we obtained using different routes allowing to control the nucleation and growth of α -Fe 2 O 3 crystallites in solution. These samples are presented in Various grain shapes and sizes ranging from 10 nm to 10 µm are obtained. The small grains are characterized by a lighter red color (samples 1 to 5) whereas a darker purple hue is observed for the pigments with larger grains (samples 6 and 7). We may also note that samples 5, 6 and 7 can be described by two grain sizes: either two different size populations (sample 5), or large grains that are agglomerates of smaller ones (samples 6 and 7).

BRDF measurements

The whole BRDF of each sample is measured by combining two distinct measurements, performed with two different setups depending on which scattering geometry is at stake. The first one does not account for the back-scattering, while the second one is specifically performed to account for this contribution. All the measurements are done on d = 1 mm thick powder layers.

Measurement out of the back-scattering direction

The main part of the BRDF is measured using the spectro-goniophotometer SHADOWS [START_REF] Potin | Shadows: a spectro-gonio radiometer for bidirectional reflectance studies of dark meteorites and terrestrial analogs: design, calibrations, and performances on challenging surfaces[END_REF],

for over 180 geometries: the scattered light is observed for

θ o ∈ [-70 • , 70 • ] with a step of 10 • , for incidences θ i ∈ 0 • , 20 • , 40 • , 60 • and azimuthal angles ϕ ∈ 0 • , 30 • , 60 • , 90 • , 135 • .
The spectral range has been limited to 360 nm -860 nm with a step of 20 nm. The reflectance measurements are calibrated with a Spectralon. Every point is measured with an integration time of 300 ms and averaged over 10 spectra. Such measurement takes about 18 hours.

When it comes to a geometry where incident and observation directions are identical, the instrument design causes the sensor to be obstructed by the light source. A specific device needs to be designed to evaluate the amount of light scattered in these directions.

Measurement in the back-scattering direction

Our design for measuring the back-scattering contribution to the BRDF is inspired from the work of [START_REF] Belcour | Bidirectional reflectance distribution function measurements and analysis of retroreflective materials[END_REF] and [START_REF] Lévesque | Measurement and modeling of the spectralon spectro-polarimetric bidirectional reflectance distribution function (brdf)[END_REF]. As these only deal with vertically held samples, they cannot be used for studying dry powders. We propose a new design to allow the measurement of any horizontally held sample.

We chose a halogen light source to illuminate at normal incidence the sample after being reflected by a beam splitter at 45 • . The size of the illumination spot is minimized by focusing the light with a 50 mm focal lens. An additional diaphragm is added to clean the focused beam.

The direction of the incident light θ i is adjusted by positioning the beam splitter at an angle of 45 • + θ i /2. As the light interacts with the powder, it is scattered in every direction and preferably back into the incident direction. It is then partially transmitted through the beam splitter and collected by a Thorlabs F950FC-A collimator. The latter is fixed on a vertical rotating plate so that we can scan different observation directions around the exact back-scattering one in the incident plane. Although the incident light may be partially polarized after reflecting on the beam splitter, for small enough incident angles the diffuse interaction with the sample has the property of depolarizing the light. For large incidence angles, the reflection on the surface becomes more important and so are the polarization effects. To neglect these polarization effects, we will only consider incident angles smaller than 40 • . We use a fiber spectrometer OceanOptics USB650 along with the SpectraSuite software in order to detect and visualize the measured signal.

In the incident plane (ϕ = 0 • ), we scan the observation angles in the range

[θ i -10 • , θ i + 10 • ]
with a step of 2 • , for incidences 0 • , 20 • and 40 • . For each incidence angle of interest, the observation angles the furthest from the exact back scattering direction that are scanned are redundant with geometries reached by the SHADOWS instrument. This allows for verification that the two distinct measurements are consistent with one another. The reflectance measurements are calibrated with a Spectralon. The integration time is chosen to maximize the detected signal and an averaging over 5 spectra is performed to minimize the noise. Less than 20 seconds is necessary to measure one spectrum.

Powder BRDF Model

As shown by the measurements on Figure 6, powders' BRDFs mainly involve back-scattering and almost no specular scattering which makes the usual specular BRDF models ineffective to reproduce their appearance. Although our model is derived from a study on pure nano-crystallized hematite powders, we chose a physically-based approach so that it can be applied to any other pigment. It involves both chemical composition and morphological parameters to accurately account for the color variations observed on such pigments. Moreover it allows to consider either isotropic or anisotropic scatterers. Inspired by both RTLSR and Spectralon's models, our model is divided into two contributions, accounting for volume and surface scattering respectively.

BRDF(θ i , θ o , ϕ i , ϕ o , λ) = A v BRDF v (θ i , θ o , ϕ i , ϕ o , λ) + A s BRDF s (θ i , θ o , ϕ i , ϕ o , λ) (17)
These two contributions are weighted with the scalar amplitude parameters A v and A s . Because the angular and spectral dependency of the two contributions BRDF v and BRDF s is the same as the complete BRDF, it will not be reminded in the following equations for readability purposes.

The main difference with the previous models is that, based on the EBS theory the back scattering is considered as part of the volume contribution. The surface contribution is then reduced to specular scattering caused by surface roughness.

Volume scattering

As EBS theory is established within the diffusion approximation, the average scattering is assumed to be isotropic. This greatly simplifies the multiple scattering contribution by considering the scattering phase function as constant over the whole scattering hemisphere. The anisotropic scattering of each particle is accounted for by choosing the phase function parameter g involved in the single scattering contribution. The extension of the model for taking into account this anisotropy in multiple scattering has only been done for a normal incidence by Tsang et al. [START_REF] Tsang | Scattering of electromagnetic waves: advanced topics[END_REF].

In the following we expose the calculations we made in order to consider anisotropic scattering in both single and multiple scattering terms for any θ i /(θ o , ϕ o ) geometry. The final model only involves two variable parameters (phase function anisotropy parameter g, extinction rate κ) that rely on the physical and chemical properties of the material that is under study.

First order ladder term L 1 : This contribution accounts for light that is only scattered once before being observed and is similar to the one introduced in Section 3.5. The anisotropy is taken care of by the scattering phase function involved in the analytical expression of BRDF ebs, L 1 in Eq. 14.

Multiple ladder term L 2+ : The phase function is involved at every scattering event that light encounters as it travels the particulate medium. We denote P 1 and P n respectively the first and the last particles that are responsible for light scattering along a light path.

At the second scattering order, the light is first incoming on P 1 , scattered with an angle Θ 1 towards the particle P n , and finally scattered towards the observer with an angle Θ n as illustrated by the Figure 5.

For the n-th scattering order, light undergoes n scattering events. As for the second order, it undergoes a first one described by the angle Θ 1 while encountering the particle P 1 and last one with the particle P n and the angle Θ n . But it also interacts with n -2 other particles in between, resulting in n -2 scattering events respectively described by the scattering angles Θ k = θ ik -θ ok , difference between the incoming direction θ ik and outgoing direction θ ok on the k-th particle, with k = 2 ; n -1 as it is shown on the Figure 5. The first scatterer that light encounters as it travels inside the medium is denoted as P 1 and involves a scattering event described by the angle Θ 1 . The last scattering event before light goes out of the medium is caused by the scatterer P n with an angle Θ n . In between these two scatterers, light is scattered with an angle Θ k at each interaction with a particle k.

By averaging over all the possible paths that light can follow inside the medium while encountering n scatterers, we can write the n-th order ladder term as follows:

L n = (4π) n-2 (ρκ) n r 2 (18) × ∫ |p(Θ n )| 2 |G 11 (r n , r n-1 )| 2 ...|p(Θ 1 )| 2 exp κ z 1 µ o + z n µ i dr 1 d...dr n
with r the distance from the sample to the observation point.

The total ladder term accounting for multiple scattering k = 2 ; ∞ ) can be written in a recursive way, giving:

L 2+ = ∬ exp(κ( z i µ o + z o µ i )A(r i , r o )dr i dr o ( 19
)
with A(r i , r o ) defined such as:

A(r i , r o ) = ρκ 4π 2 p(Θ 1 )p(Θ n ) exp(κ|r i -r o |) |r i -r o | 2 (20) + ∫ ρκ 4π p(Θ k ) exp(κ|r i -r k |) |r i -r k | 2 A(r k , r o )dr k
Under the assumption of the diffusion approximation, this term can be simplified to average the contributions of all scatterers by extending the work from [START_REF] Tsang | Scattering of electromagnetic waves: advanced topics[END_REF]. The resulting BRDF of the ladder contribution is therefore reduced to the following analytical expression:

BRDF L 2+ (κ, g) = 3ρ 2 1 -µ 2 i c 2 2Q 1 -µ o c - µ 2 i µ i + µ o (21)
with

c = 3(1 -ρ)(1 -ρg) + α 2 (22)
and parameters α and Q defined as in Eqs 15.

One should note that the ladder term expression is similar to the one defined by the Eq. 14, except that it accounts for possible anisotropic average scattering in the medium through a varying phase function of parameter g (Eq. 22). If one wants to describe isotropic scattering, a constant phase function such as p(Θ) = p (similarly g = 0) should be used which leads back to the expression originally introduced in Section 3.5.

Multiple cyclical scattering C:

As for the ladder term, we are aiming to take into account the anisotropic scattering of particles in the cyclical term. We will use the same assumptions as illustrated by the Figure 5. We consider a beam of light that is scattered by n particles as it follows a path inside the medium. By accounting for the correlation between this path and its reversal and averaging over all the possible paths, the n-th order cyclical term can be written:

C n = (4π) n-2 (ρκ) n r 2 ∫ |p(Θ n )| 2 |G 11 (r n , r n-1 )| 2 ...|G 11 (r 2 , r 1 )| 2 |p(Θ 1 )| 2 × cos k(µ i -µ o )(z 1 -z n ) exp κ 2 1 µ o + 1 µ i (z 1 + z n ) dr 1 ...dr n ( 23
)
The complete cyclical contribution is obtained by summing over the n-th orders and can be simplified replacing G 11 by its definition (Eq. 8). By using the new expression of A(r i , r o ) given by Eq. 20:

C = ∬ cos(2k(µ i -µ o )(z o -z i ))exp κ 2 z i µ o + z o µ i A(r i , r o )dr i dr o (24)
Finally, by using the diffusion approximation, the cyclical term contributes to the complete appearance of the dry powder through the partial BRDF:

BRDF C (κ, g) = 3ρ 2 [c + v(1 -exp(-2cτ o ))] 2µ i µ o cv[(c + v) 2 + u 2 ] ( 25
)
where c is defined as in Eq. 22, and u, v as in Eqs. 15.

Once again, this expression is similar to the one obtained assuming isotropic averaged scattering in the medium (Eq. 14). However it involves the phase function parameter g to account for anisotropic scattering. The use of a constant phase function p(Θ) = p (similarly g = 0) brings back to the isotropic case already discussed in Section 3.5.

Our model's complete EBS contribution accounting for anisotropic scattering in the powder medium is the sum of the two ladder components given by first equation of Eq. 14 and Eq. 21,

and the cyclical component given by Eq. 25.

Surface scattering

In our model, surface scattering is reduced to the specular contribution. The light that is incoming on a rough surface such as the interface between air and a dry powder is scattered in a lobe around the specular direction.

Two options have been considered to include such a phenomenon in our model. The former is to model the specular peak by a gaussian function. According to the measurements performed in Section 4 the specular scattering is very weak so this simple method is reliable enough to reproduce the specular contribution of dry powders, although it is empirical. By assuming a uniform distribution of the surface normals, it requires only one parameter which is the gaussian standard deviation σ. With this modeling the surface scattering BRDF is written as:

BRDF s (σ) = 1 σ √ 2π exp( -Θ 2 2σ 2 ) (26)
The second option, which is the one we chose to use, relies on the physical approach from Harvey-Shack [START_REF] Krywonos | Predicting surface scatter using a linear systems formulation of non-paraxial scalar diffraction[END_REF] that accounts for diffraction effects caused by the surface asperities. Using this theory, the incoming light on a rough surface is then scattered around the specular direction according to the following BRDF:

BRDF s (l c ) = 4π 2 λ 4 (cos θ i + cos θ o ) 2 PSD(f) (27)
where the power spread function (PSD) models a rough surface topography. We chose a gaussian power spread function, involving a roughness parameter σ s and a correlation length l c that quantifies the mean distance for which the surface roughness of two points on the surface are completely independent.

PSD(f) = πσ 2 s l 2 c exp(-π 2 l 2 c ( f 2 x + f 2 y )) (28)

Complete model

Our complete model is the weighted sum of volume and surface contributions according to Eq. 17 respectively given by the adapted enhanced back scattering model and a specular model.

Note that the roughness parameter involved in the surface contribution has been included in the amplitude parameter A s = A s σ 2 s in order to reduce the number of optimization parameters.

BRDF(θ i , θ o , ϕ i , ϕ o , λ) =A v [BRDF L 1 (κ, g) + BRDF L 2+ (κ, g) + BRDF C (κ, g)] (29)
+ A s BRDF s (l c )
The model depends on a total of 2 empirical (A v and A s ) and 3 physical (g,κ,l c ) parameters to accurately simulate the BRDF of dry powder. As it is shown in the Section 6, our model works even better while including another physical parameter to refine which is the effective radius of the grains, defined as the mean size of the grains actually seen by the incoming light and therefore responsible for the scattering. The influence of such parameter is shown on Figure 7.

The parameters involved in our model are summarized in Table 1.

Optimization procedure

In order to reproduce the BRDF of a powder sample, the model is fitted to measurements by optimizing the variable parameters of Table 1. We used the least squares minimization of scipy along with the sequential least square programming method (SLSQP) to fit our BRDF model to the measurement data y m (Eq. 30). A first option is to optimize while using a fixed radius. The average dimensions of the grains of a powder sample is deduced from SEM images and reduced to only one dimensional parameter by computing the radius of a sphere that is equivalent to the grain. To better account for the shape anisotropy this radius is computed so that the ratio of the surface over the volume of the equivalent sphere and of the grain remains constant. The final mean radius is obtained by averaging on the radii of multiple grains that are in the field of view of a single SEM image. A second option is to optimize on the grains' effective radius, which both releases constraints on the optimization process and gets rid of the uncertainty on the experimental grain size definition.

The optimization is performed for all measurement geometries γ = (θ i , θ o , ϕ i , ϕ o ) with θ i < 60 • .

For larger incidences, we were not able to measure the back-scattering contribution so we do not account for these geometries in the optimization process.

E = min( γ,θ i ≤40 • |BRDF γ (A v , g, κ, b, A s , l c ) -y γ,m | 2 ) (30)
We only account for the measurements at the wavelength 700nm. The spectral dependency of the model will be discussed later on.

Results

Complete measurement

The complete BRDF measurements of dry hematite α -Fe 2 O 3 powders were performed in two independent steps described in Sections 4.2. By comparing the results for geometries reached by both instruments, the measurements happen to be shifted in amplitude from one another. This may be caused by the underestimation of the back scattering of the spectralon used for photometric calibration. Since the difference is smaller than 5%, a re-scaling of the back-scattering measurement is done so it matches the goniometer results.

A typical experimental BRDF is shown on Figure 6 after merging the measured data. The measurements show that the BRDF of dry hematite powders is composed of an overall diffuse contribution, and two directional scattering phenomena including mainly back-scattering and some specular scattering (Figure 6). The BRDFs of characterized samples vary in intensity and width of the specularly and back scattered lobes. We give results for the six other samples under study in Supplementary Work. The spreading of the back scattering peak with the azimuth angle is due to the colormap interpolation between the measured points.

Model validation

The ability of the model to simulate the BRDF of pigments is validated against the various hematite powders measurements by following the procedure of Section 5.4. In the following we chose to define the grains' effective radius as an optimization parameter as it leads to better results. Indeed, as it is shown in Table 2 it allows a nicer modeling of the back scattering peak amplitude, while keeping values pretty close to the radius deduced from microscopy observations.

Although the two populations are well distinguishable in samples 5 and 7 (see Figure 3), it is harder to spot on sample 6. In that case, the smaller grains actually refer to the roughness at the surface of the observed ellipsoidal grains. For illustration, the results from two optimizations involving or not the radius performed on sample 1 are shown on Figure 7.

The same optimization process is performed using both Spectralon and RTLSR models adapted to be relevant in the case of dry powders. Such adaptations to the original models are detailed in Supplementary work. To compare our results to the one obtained through these two other models allows a better evaluation of our model.

Our model is able to reproduce the measured BRDF of hematite powders with an accuracy comparable to that of the Spectralon model (< 8% of the measured reflectance factor), but with a much lower complexity: only 6 parameters related to physical properties of the samples are required against 12 empirical ones for the Spectralon's model (cf. Figures 9,10 and 11). Our model nicely fits the back-scattering peak, and is even more accurate than the other ones for reproducing the slope of the measured BRDF (cf. Figure 8). The mean squared error (MSE) from the fitting procedure is also of same magnitude for both Spectralon and our model (MSE < 0.05).

In comparison, the RTLSR gives much poorer results (MSE > 0.1). Detailed values are given for each samples in Table 2. Our model is particularly efficient in the case of samples whose morphology is easily characterized, that is to say for the samples with well-crystallized single-crystals. Both the back-scattering and the specular contributions are detected and modeled, even though the latter is sometimes too weak for the modeling to be reliable. Errors are mostly spotted around the back scattering peak because of the insufficient measurement resolution and interpolation at these geometries.

One should note that we have encountered difficulties in the modelization of samples 5 and 7's BRDFs that we attributed to the more complex morphology of these powders. Since the refined parameters permit our model to fit quite precisely with the measurements in any geometry, it appears that none of them depends on the incoming, outgoing and azimuthal directions. The model is then able to predict the appearance of the studied powder in any other geometry that hasn't been considered during its optimization. The results are convenient even for more grazing angles such as incidence θ i = 60 • as shown on Figure 7. . Optimization of our model on sample 4, at different wavelengths using the optical constants from [START_REF] Querry | Optical constants[END_REF]. For smaller wavelengths, our model is not able to accurately fit the measurements.

Discussion & Conclusion

From BRDF measurements, we highlighted that in the particular case of dry powders, the light is mainly back scattered. We proposed a new BRDF model accounting for three different For now, our model can only be used as a fitted model. The parameters describing the appearance of a specific powder are retrieved from the fitting of the model to the measured BRDF of such powder. In order to predict the BRDF of any other powder of known chemical composition and morphology that has not been measured, it would be necessary, as a future work, to establish a relationship between the optimization parameters of the model and the physico-chemical parameters of the material the powder is composed of. This requires to study much more samples with more variability in sizes.

The validity of the model wavelength dependency is also a main issue that couldn't be solved throughout this work. At smaller wavelengths, the optimization is not robust enough to lead to nice fitting results (see Figure 12). This might be caused by a poor knowledge of the optical constants of the studied material which is the only fixed parameter known as spectrally varying.

Indeed, it is mainly described by the optical constant at stake whose value appears to be quite doubtful according to the discrepancies between the results in the literature [START_REF] Querry | Optical constants[END_REF][START_REF] Triaud | Refractive indices of hematite[END_REF]. Such uncertainties are often ignored in both material science and computer graphics whereas it surely affects the results.

As a future work we would like to implement our dry powder BRDF model in a physically based renderer so that photo-realistic renderings of the powders can be made. This would first allow the validation of the model by comparing the rendered image with photographs of the actual samples. It would also permit, by letting the user choose the values of the variable parameters to either reproduce the appearance of other powder materials, or to create new ones.

Although the model has been tested for only one material, we believe it can be used for any other pure pigment. In the case of parietal paintings, ochers (mixture of one iron-based coloring pigment and colorless powders such as clays) and mixtures of different pigments were more commonly used than pure and almost pure iron oxide. In order to reproduce photo-realistic appearance of such paintings, we would like to extend our model to more complex powders.

The assumption of optically thick layer of powder is also too restrictive. An extension of the enhanced back scattering theory for thin films and isotropic scatterers has been proposed by [START_REF] Van Der | Light scattering in strongly scattering media: multiple scattering and weak localization[END_REF] and may be a good starting point for our model adaptation, leading to a better reproduction of the appearance of painting layers.

Fig. 2 .

 2 Fig. 2. Definition of the vector f.

Figure 3 .

 3 Figure 3. The samples denoted from 1 to 5 have been prepared following methods described by Schwertmann & Cornell [31]. 1, 2 and 5 have been obtained by forced hydrolysis of Fe I I I salts solutions in acidic conditions, while 3 and 4 are the results of the transformation of a ferrihydrite precipitate in aqueous suspension under basic conditions. We generated peanut-like hematite particles 6 following the routes proposed by Sugimoto et al. [32], NaH 2 PO 4 salts. Sample 7 is a commercial Iron(III) oxide powder sample purchased from Puratronic®. The color and the morphometry of each sample are respectively deduced from naked eye observations completed with optical microscopy using a Nikon Eclipse LV100ND microscope, and by scanning electron microscopy observations performed with a Zeiss Ultra+ microscope.

Fig. 3 .

 3 Fig. 3. α -Fe 2 O 3 hematite samples under study, showing various colors and grains morphology. Upper row:The photographs are calibrated using an X-rite 24 patches Classic Colorchecker. Middle row: Observations using optical microscopy. Lower row: Samples morphology observations using Scanning Electron Microscopy. The samples 1 to 4 show similar grain shape, with slightly different sizes from 10 to 50 nm in diameter, and light red hue. Sample 5 presents two size populations and shows a red-pinkish hue. Samples 6 and 7 have a more complex morphology with grains that seem to be agglomerations of smaller ones. The grains are larger and a darker purple hue is observed. According to optical microscopy, 7 is even composed of different types of grains, with different colors.

Fig. 4 .

 4 Fig. 4. Back-scattering measurement device. The beam of a white source is focused on the sample surface at an incidence θ i using a 50 mm focal lens L and a beam splitter. The reflected light is collected through a collimator C in different observation directions θ o .

Fig. 5 .

 5 Fig.5. Multiple scattering with two and n particles. The first scatterer that light encounters as it travels inside the medium is denoted as P 1 and involves a scattering event described by the angle Θ 1 . The last scattering event before light goes out of the medium is caused by the scatterer P n with an angle Θ n . In between these two scatterers, light is scattered with an angle Θ k at each interaction with a particle k.

Fig. 6 .

 6 Fig. 6. Complete BRDF measurement of the hematite powder sample 1 at the wavelength 700 nm for an incidence θ i = 40 • . The measured geometries are spotted with black dots, and the intermediate values are interpolated. Light is mostly reflected in the backward direction, and a very weak specular scattering peak may be observed as well. The spreading of the back scattering peak with the azimuth angle is due to the colormap interpolation between the measured points.

Fig. 7 .

 7 Fig. 7. Influence of the grain radius optimization on the results of sample 2, in the incident plane (ϕ = 0 • ) for incidences θ i = 0 • , 20 • and 40 • . The effective radius found by optimization (b = 17 nm) is pretty close to the theoretical result deduced from SEM observations (b = 13 ± 3.0 nm).

Fig. 8 .

 8 Fig. 8. Comparison of the optimization results of sample 2 from the three BRDF models: Our model, Spectralon model, RTLSR model. (a) Results in the incident plane (ϕ = 0) for incidence θ i = 20 • . (b) Residual error between the measurements and the fits. Our model slightly underestimates the back-scattering peak compare to the Spectralon model. However, it fits the slope of the measured BRDF much better.The residual error of these two models are similarly small. In comparison, the RTLSR model is less suitable to fit the measured data.

Fig. 9 .

 9 Fig. 9. Upper row: Fitting of the measured BRDF of sample 1 using our model, spectralon and RTLSR models adapted to dry powders, at incidence θ i = 0 • . Lower row: difference between the measurement and the fit. Our physically-based model (6 optimized parameters) gives better results than the two other ones as it is more performant than the RTLSR model (1 optimized parameter) and as performant as the Spectralon model (12 optimized parameters) with far less optimization parameters.

Fig. 10 .

 10 Fig. 10. Upper row: Fitting of the measured BRDF of sample 1 using our model, spectralon and RTLSR models adapted to dry powders, at incidence θ i = 20 • . Lower row: difference between the measurement and the fit. Our physically-based model (6 optimized parameters) gives better results than the two other ones as it is more performant than the RTLSR model (1 optimized parameter) and as performant as the Spectralon model (12 optimized parameters) with far less optimization parameters.

Fig. 11 .

 11 Fig. 11. Upper row: Fitting of the measured BRDF of sample 1 using our model, spectralon and RTLSR models adapted to dry powders, at incidence θ i = 40 • . Lower row: difference between the measurement and the fit. Our physically-based model (6 optimized parameters) gives better results than the two other ones as it is more performant than the RTLSR model (1 optimized parameter) and as performant as the Spectralon model (12 optimized parameters) with far less optimization parameters.

  scattering interactions: single scattering, incoherent multiple scattering and coherent multiple scattering responsible for the sharp peak in the backward direction. It is a semi-empirical and analytical model based on the propagation of light inside a granular medium that involves the physico-chemical properties of the powders connected to only six refinement parameters.The model has been validated against the BRDF measurements of seven α -Fe 2 O 3 hematite powders, and compared to two adapted existing models. Ours is able to reproduce the measured BRDF with an accuracy comparable to the one of the Spectralon's model, but with a much lower complexity. It can be used in its original form to reproduce the appearance of a mono disperse powder with spherical grains, or in its extended form presented in Supplementary work to account for more complex morphology involving multiple grains size populations.Our model still has some limitations. The results obtained for specular scattering are quite questionable since the specular contribution is very weak and almost unnoticeable on the BRDF measurements. A more accurate model should be established by studying powders showing shiny and sparkling reflections such as eye shadow.

Table 1 .

 1 Parameters involved in our BRDF model. The model can be optimized by refining only 5 parameters while the grain radius b is fixed, or by refining 6 parameters including the same 5 previous ones and an effective radius.

	Volume scattering
	A v	Intensity of volume scattering
	g	Scattering anisotropy
	κ	Extinction rate
	b	effective radius (optional)
	Surface scattering
	A s	Intensity of surface scattering
	l c	Correlation length

Table 2 .

 2 These powders are not monodisperse since they show two size populations and agglomerations of grains. An extension of our model is proposed in Supplementary work to deal with powders presenting two grain size populations.

	Sample Radius SEM	Effective			Ours	Spectralon	RTLSR
		(nm)	radius (nm)			MSE	MSE	MSE
	1	12 ± 1.6	16			0.035	0.043	0.17
	2	13 ± 3.0	17			0.052	0.052	0.21
	3	22 ± 3.7	16			0.044	0.046	0.25
	4	25 ± 4.1	15			0.037	0.039	0.29
	Sample Large radius	Effective	Small radius	Effective ra-
		SEM (nm)	radius (nm)	SEM (nm)	dius (nm)
	5	83 ± 34	83	47 ± 5.3	46	0.063	0.074	0.71
	6	687 ± 64.0	687	10 ± 0.4	43	0.041	0.043	0.49
	7	2222 ± 1054	2222	49 ± 13	49	0.019	0.019	0.30

Results from optimizations on radii. The radius deduced from SEM observations are compared with the effective radius obtained by the optimization of our model. For samples 5 to 7, the model extension for two size populations presented in Supplementary work is used. As it takes two different radii to optimize, two effective radii are obtained. The resulting MSE of such optimizations are also given and can be compared to MSE of optimizations over the Spectralon model and the RTLSR model.

Acknowledgments:

The authors would like to thank Alain Ibanez and Simon Giraud (Institut Néel) for the hematite samples syntheses. We also thank Corinne Félix (Institut Néel) for her help in developing the back-scattering measurement device, as well as Bernard Schmitt and Olivier Brissaud (IPAG) for their help in the reflectance measurements using SHADOWS gonio-spectrophotometer. At last, we thank Romain Pacanowski (INRIA) for the proofreading of this paper. This work was financially supported by the French National Research Agency in the framework of the Investissements d'Avenir program (ANR-15-IDEX-02, Cross Disciplinary Program Patrimalp).

Disclosures:

The authors declare no potential conflicts of interest.

Data availability: Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Supplemental document: See Supplement 1 for supporting content.