## Abstract

We perform Monte Carlo light scattering simulations to study the angular distribution of the fluorescence emission from turbid media and compare the results to measured angular distributions from fluorescing white paper samples. The angular distribution of fluorescence emission is significantly depending on the concentration of fluorophores. The simulations show also a dependence on the angle of incidence that is however not as evident in the measurements. A detailed analysis of the factors affecting this angular distribution indicates that it is strongly correlated to the mean depth of the fluorescence process. The findings can find applications in fluorescence spectroscopy and are of particular interest when optimizing the impact of fluorescence on e.g. the appearance of paper as the measured values are angle dependent.

© 2015 Optical Society of America

## 1. Introduction

Fluorescent dyes are widely used to improve or affect the appearance of textiles and paper products. Daylight fluors are typically used for attention-grabbing applications, whereas a special type of fluorescent dyes is used to increase the whiteness of paper and textiles. These dyes are referred to as fluorescent whiteness agents (FWAs) or optical brightening agents (OBAs). They absorb UV radiation and re-emit light in the visible part of the electromagnetic spectrum. The appearance of a reflecting surface depends on its optical properties, on the observer, and on the spectral power distribution (SPD) of the light source. For fluorescing turbid media, the radiance factor itself depends on the SPD [1] in addition to the medium optical properties. In order to predict the appearance in different illuminations, bispectral radiance data relating spectral radiance to individual excitation wavelengths is required. This matrix of bispectral radiance factors is known as the Donaldson matrix [2] and can be determined with the double monochromator method [3, 4]. These measurements are generally performed in a fixed measurement geometry following one of the recommended geometries from the CIE [5]. Extension to other measuring geometries or viewing environments relies on the assumption of a Lambertian distribution of the luminescent radiance factor.

The angular distribution of the fluorescence emission from turbid media has gained attention in recent years. Holopainen et al. [6] studied the fluorescence emission from solid amorphous material and showed a clear non-Lambertian behaviour. On the other hand, Johansson et al. [7] found that the fluorescence emission of uncoated paper was much less anisotropic and Tominaga et al. [8] reported nearly Lambertian luminescent radiance factor from coloured fluorescing papers and tapes. For the case of a green fluorescing paper, the fluorescence even increased with increasing viewing angles, as opposed to the measurements reported by Holopainen et al.

General radiative transfer theory has been shown to successfully describe the angular reflectance of turbid media [9]. Gebhart et al. [10] performed Monte Carlo simulations of fluorescence in turbid media and reported that the fluorescence is Lambertian regardless of the angle of incidence or of the optical properties of the medium. However, the simulations were performed with parameters typical for tissue with low scattering coefficients, and with constant absorption coefficient and asymmetry factor. In this paper we study the anisotropic behaviour of the fluorescence from turbid media by comparing Monte Carlo simulations to the measured bidirectional luminescence distribution function (BLDF) of paper samples with different amount of FWAs. The Donaldsson matrix of the samples is known from previous measurements with a bispectrophotometer [11], and is used to estimate the optical parameters needed in the simulations. It is important to characterise the dependence of the anisotropy on the medium optical parameters and angle of incidence since this provides information on the underlying factors that affect the angular distribution of fluorescence. This is of particular importance when determining the quantum efficiency of fluorescing turbid media or the impact of fluorescence on appearance since these will be significantly affected by the measurement geometry.

This paper is structured as follows: Section 2 describes a method used to measure the relative BLDF with a single monochromator goniospectrophotometer and a UV filter, the Monte Carlo simulation model and a parameter estimation method developed to retrieve the model parameters from bi-spectrophotometer measurements. The results presented in Section 3 show the impact of the optical properties on the angular distribution of the fluorescence, and that there is a strong dependence of measured relative BLDF on fluorophore concentration. Section 4 provides the analysis of the findings and underlines the possible reasons for discrepancies regarding observed angular dependence of fluorescence reported in the literature. We conclude and propose future work in Section 5.

## 2. Method

In order to study the impact of fluorophore concentration on the angular distribution of the fluorescence emission, we measure the angle-resolved radiance in the emission wavelength band of paper samples of different FWA concentrations when illuminated with monochromatic UV excitation at two different angles of incidence. The goniospectrophotometer we use does not have a monochromator on the detector side. We therefore first show that using a UV-blocking filter in front of the detector enables determining the relative bi-directional luminescence distribution function (BLDF) in the emission band, provided that the angular distribution of the fluoresced light is independent of the emission wavelength. We make use of a Monte Carlo light scattering simulation tool to model the scattering and fluorescence process with radiative transfer theory and compute the BLDF from 360 nm to 440 nm. In order to determine the scattering and absorption coefficients at the excitation wavelength and the quantum efficiency, we use measurements from a bi-spectrophotometer with fixed geometry and a parameter estimation method suitable for opaque samples. We then use the Monte Carlo model to study the impact of the optical parameters on the BLDF and verify the assumption made in the BLDF measurements. More details are given in the following subsections.

#### 2.1. Bi-spectrophotometer measurements

We use four 80 g/m^{2} unfilled sulphate bleached chemical pulp papers with FWA concentrations 3, 9, 18 and 36 kg/T of dried pulp. The Donaldson matrix of the samples was measured with a custom-made bispectrophotometer at the University of Eastern Finland in Joensuu [12]. The samples are a subset of the samples characterized by Coppel et al. [11] who describe the detailed measurement setup and sample manufacturing. The bispectrophotometer measurements were performed in a 0/45 geometry in the range 280 to 700 nm. In the present study, only the bispectral radiance factors *D*(360, 360) and *D*(360, 440) are used. This combines maximum light absorption by the FWA at 360 nm and maximum fluorescence at 440 nm. The measurements were performed on a single 80 g/m^{2} sheet over a black background (*D*_{0}) and on an opaque pad of similar sheets (*D*_{∞}), except for the 36 kg/T sample for which only *D*_{∞} had been measured. Also, at this FWA concentration, the luminescent radiance factor of one single sheet is very close to that of an opaque sheet.

#### 2.2. Goniospectrophotometric measurements

The samples were measured with a Lambda 1050 spectrophotometer (Perkin Elmer) equipped with an ARTA goniophotometer accessory (OMT Solutions BV). The sample is illuminated with collimated monochromatic light, and the reflected light is detected with a broadband photomultiplier tube (PMT). In order to filter out the fluorescence from the total radiation, an UV-blocking filter with a cut-off wavelength of 415 nm (UV-dichroic blocking filter LP415-35.5, Midwest Optical Systems Inc.) is mounted in front of the detector. The sample is placed on a motorized rotation stage, and the angle of the incident light is varied by rotating the sample. The detector is mounted on a second motorized rotation stage, allowing the detector to be positioned at any angle relative the normal of the sample surface, except for ±30° in proximity to the light source. A detailed description of the instrument is given in [7].

Measurements are carried out with an excitation wavelength of 360 nm and with angle of incidence 0° and 45°. We limit our analysis to measurements made in the plane of incidence and to viewing polar angles ranging from 30 to 75°. The size of the illuminated area is 5 x 8 mm at normal incidence, and the solid angle subtended by the detector when viewed from the sample plane is 0.0049 sr (detection half angle 0.8°). The goniophotometer is a double-beam instrument allowing simultaneous determination of the incident Φ* _{i}* and reflected Φ

*power.*

_{r}The bidirectional reflectance distribution function (BRDF) *f _{r}* is defined as [13]

*L*(

_{r}*θ*,

_{i}*ϕ*,

_{i}*θ*,

_{r}*ϕ*) is the differential radiance reflected in direction (

_{r}*θ*,

_{r}*ϕ*) and differential solid angle

_{r}*dω*, originating from the differential irradiance d

_{r}*E*(

_{i}*θ*,

_{i}*ϕ*) impinging on an area element d

_{i}*A*from direction (

*θ*,

_{i}*ϕ*), as illustrated in Fig. 1.

_{i}As in any practical measurement situation, the differential quantities have to be approximated by their finite counterparts, i.e *f _{r}* ≈

*L*. The measurements are conducted with a setup where the entire illuminated spot

_{r}/E_{i}*A*is contained within the detected area

_{i}*A*, as illustrated in Fig. 2.

_{f}The reflected radiance is, according to the radiance invariance principle, equal to the detected radiance, which in terms of measured power is given by

where*A*is the area of the detector aperture and Ω

_{d}*is the projected solid angle defined by the detector’s field of view. Since the measured power is contained within the solid angle*

_{f}*ω*=

*A*cos

_{i}*θ*

_{r}/r^{2}, where

*r*is the distance between the detector and the sample surface, Ω

*can be replaced by*

_{f}*ω*, which gives us

*A*is the size of the illuminated area at normal incidence. The final expression for the BRDF, in terms of measured quantities, is then

_{s}*λ*

_{1}is the excitation wavelength and

*λ*

_{2}is the emission wavelength. However, since there is no monochromator at the detection side, the output voltage

*V*(

*θ*) of the detector is in this case related to the detected power according to

_{r}*λ*) is the sensitivity of the detector and Ψ(

*λ*) is the transmittance of the UV-blocking filter. Upon excitation at 360 nm the measured output is thus the total contribution of the emitted light at all wavelength above the cut-off wavelength and the BLDF at one specific emission wavelength

*λ*

_{2}cannot be retrieved. However, assuming that the angular distribution of the emitted light is independent of the wavelength one can write where

*f*

_{1}(

*θ*) is the angular distribution of the emitted power and

_{r}*f*

_{2}(

*λ*

_{2}) its spectral distribution. By substituting Eq. (8) into Eq. (7), the output voltage is given by where

*C*is independent of the viewing angle

*θ*. Thus the measured voltage is then proportional to the emitted power at a specific emission wavelength

_{r}*λ*

_{2}. The BLDF calculated according to Eq. (6), although not equal to the correct BLDF in absolute values, can therefore after normalization be used for comparative measurements.

#### 2.3. Simulation model

Radiative transfer (RT) theory describes the propagation of light in turbid media such as paper characterised by their scattering (*s*) and absorption (*k*) coefficients and a phase function describing the angular distribution of the single scattering events. The single scattering phase function is modelled with the Henyey-Greenstein function [15] controlled by the asymmetry factor *g*. When fluorescence occurs, the radiation absorbed is re-emitted isotropically at a longer wavelength according to a fluorescence matrix *Q* determining the probability for emission at different wavelength upon absorption at a lower wavelength. Within the turbid medium, the probability of travelling a distance *t* before scattering or absorption is given by [16]

*a*) and the mean free path (

*l*

_{e}) expressed as where 1/

*l*is the extinction coefficient. Since the samples are characterised by their basis weight (80 g/m

_{e}^{2}) and their thicknesses are unknown, the scattering and absorption coefficients are reported in m

^{2}/g and

*l*in g/m

_{e}^{2}. The Open PaperOpt simulation tool [17] is used to solve the RT equation accounting for fluorescence using a Monte Carlo method. Wave packets are sent onto the top surface. Their initial position and direction follow the light source distributions. The path of each wave packet is recorded until it eventually leaves the simulation volume or is absorbed. The path length probability between events is governed by the extinction coefficient. The probability of scattering or absorption is then given by the single scattering albedo, and the wave packet may finally change wavelength upon absorption (i.e. fluorescence) according to the quantum efficiency

*Q*. When determining the scattering and absorption at the emission wavelength and the albedo at the excitation wavelength, fluorescence is not of concern and a faster discrete-ordinate solver, DORT2002 [18], is then used.

#### 2.4. Parameter estimation

We assume that the asymmetry factor *g* does not depend on wavelength and it is therefore the same at the excitation and emission wavelengths. A value of 0.8 has been previously shown to represent measured BRDF data of uncoated paper reasonably well for all wavelengths [9] but values of 0.6 and 0.4 are also used here to study the impact on the angular dependency of the fluorescence.

Once the asymmetry factor *g* and basis weight are defined, the scattering and absorption coefficients can be estimated by minimizing the square difference between simulated and measured reflectance of one single sheet and of one opaque pad of identical sheets. This method is used to determine *s* and *k* at 440 nm. However, due to the large absorption by the FWA in the UV part of the electromagnetic spectrum this method cannot be used to determine *s* and *k* at 360 nm since the reflectance of one sheet does not significantly differ from the reflectance of an opaque pad. A modification of the method proposed by Coppel et al. [11] is used instead. This method makes use of the fact that although the single sheet is opaque in the UV band of the electromagnetic spectrum, the luminescent radiance factor is larger for an opaque pad than for a single sheet. The method we use here is as follows.

- The ratio
*k/s*(or the single scattering albedo*a*) is first estimated from the reflectance factor of the opaque pad only, thus optimizing*k/s*to minimize the difference between simulated and measured*D*_{∞}(360, 360). *s*and*Q*are then optimized to minimize the difference between measured and simulated bispectral reflectance factor of a single sheet [*D*_{0}(360, 440)] and of an opaque pad of sheets [*D*_{∞}(360, 440)]. This optimization is performed here in two steps. First, different*s*values are generated and the corresponding*k*values are computed according to the estimated*k/s*ratio. For each*s*(360) and*k*(360) pair,*Q*is optimized to minimize the difference between simulated and measured bispectral reflectance factor. This optimization can be formulated as$$\underset{Q(360,440)}{\text{min}}{\left[{{D}_{0}^{\prime}(360,440)|}_{s(360),k(360)}-{D}_{0}^{\u2033}(360,440)\right]}^{2},$$where*D′*_{0}and*D″*_{0}are the simulated and measured bispectral reflectance factors, respectively. This gives three optimized parameters,*s*,*k*, and*Q*, for each trial*g*. In the second step,*D*_{∞}(360, 440) is simulated for each set of parameters and the optimal set can be determined graphically by fitting simulated and measured*D*_{∞}(360, 440). Note that varying*s*with constant*k/s*ratio is equivalent to varying the mean free path*l*_{e}while keeping the single scattering albedo*a*constant.

*D*

_{0}is not available, the scattering coefficients at excitation and emission wavelengths as well as the quantum efficiency are estimated by extrapolation of the values at lower FWA amount. The absorption coefficients are then determined to fit the measured

*D*

_{∞}.

## 3. Results

#### 3.1. Parameter estimation

The scattering and absorption coefficients and the quantum efficiency of the paper samples estimated from the measurements with the 0/45 bi-spectrophotometer are given in Table 1 for three different asymmetry factors. Figure 3 shows the difference between simulated and measured bispectral reflectance of an opaque pad [Δ*D*_{∞}(360, 440)] of the 9 kg/T FWA paper sample with *g* = 0.8 and for different mean free paths *l*_{e}. The single scattering albedo at the excitation wavelength is kept constant while *Q* is fitted in order to match simulated and measured bispectral reflectance for the case of a single 80 g/m^{2} paper sheet over a black background. Since each *l*_{e}(360) and *Q*(360, 440) pair impacts differently on the bispectral reflectance of the opaque pad, they can be estimated directly from Fig. 3, and since the albedo is held constant we obtain *s* and *k* at the excitation wavelength. Note that the estimated scattering coefficients depend strongly on the assumed asymmetry factor *g*, since more forward scattering (larger *g*) requires a larger scattering coefficient to result in the same reflectance [19].

#### 3.2. Effect of mean free path at excitation wavelength and of asymmetry factor on BLDF

Figure 4 shows the simulated scaled BLDF at 440 nm of an opaque pad of 18 kg/T FWA samples at 0° angle of incidence for different mean free paths *l*_{e}(360) and different asymmetry factor *g*. As in the parameter estimation, the single scattering albedo *a*(360) is held constant for each *g* to fit the measured *D*_{∞}(360, 360), which is equal to *D*_{0}(360, 360) since the single sheet is opaque at the excitation wavelength. Note that *a*(360) is on the other hand dependent on *g* in order to keep *D*_{∞}(360, 360) constant. The estimated *l*_{e}(360) that fit both *D*_{0}(360, 440) and *D*_{∞}(360, 440) correspond to the ○ markers. Similar curves were obtained for the 3 and 9 kg/T FWA samples. The shape of the scaled BLDF depends significantly on the assumed asymmetry factor and on *l*_{e}(360). Smaller *l*_{e} at excitation wavelength and larger *g* lead to a less Lambertian BLDF. The corresponding simulated BRDFs at the 360 nm excitation wavelength are shown in Fig. 5. The simulated BRDFs all cross at 45° viewing angle since the single scattering albedo *a*(360) is optimised to fit the measured *D*_{∞}(360, 360) in a 0/45 geometry. We observe that, as opposed to their impact on the BLDF, the mean free path *l*_{e}(360) has a negligible impact on the BRDF at excitation wavelength and that the asymmetry factor *g* only affects the BRDF at larger viewing angles.

#### 3.3. BLDF versus emission wavelength

Figure 6 shows the scaled BLDF at different emission wavelengths for the 3 kg/T FWA sample illuminated at 45° and 360 nm. The angular distribution of the fluoresced light is independent of the wavelength for all emission wavelengths above 420 nm. On the other hand, the fluoresced light at 400 nm is much more Lambertian than at the other wavelengths. However, as shown in Fig. 7, the light absorption is large and the quantum efficiency small between 400 and 420 nm. This means that very little light is fluoresced below 420 nm (the absolute BLDF is much larger at wavelengths above 420 nm). Hence, the assumption of constant angular distribution made in Eq. (8) holds as the contribution from wavelengths below 420 nm is negligible. When the absorption at the emission wavelength is large, only light fluoresced close to the top surface has a probability large enough to reach the surface and the fluoresced light becomes more Lambertian.

#### 3.4. Measured and simulated fluorescence

Figure 8 shows the measured scaled BLDFs and the simulated ones for 0° and 45° angles of incidence, and with estimated scattering and absorption coefficients and quantum efficiencies for *g* = 0.8 according to Table 1. Both the measured and simulated BLDFs are strongly non-Lambertian at low FWA amount and decrease significantly with increasing viewing angle. The BLDFs become more Lambertian with increasing amount of FWA and at 36 kg/T concentration of FWA, the measured BLDF is close to Lambertian for 45° angle of incidence. According to the simulation, the anisotropy depends also on the angle of incidence and is larger for 0° than for 45°. This behavior is not evident in the measurements, as it will be discussed in the next section.

## 4. Discussion

Both measurements and simulations show that light fluoresced in turbid media is non-Lambertian and that the BLDF decreases with increasing viewing angles. The shape of the BLDF is significantly dependent on the FWA concentration and thus on the absorption coefficient in the excitation band of the FWA. At high FWA concentrations, the BLDF is more Lambertian than at lower concentrations. The Monte Carlo simulations and parameter estimation show that the scattering coefficient and the asymmetry factor in the excitation band also play a significant role. The BLDF is less Lambertian for larger mean free paths at constant single scattering albedo in the excitation band, and larger asymmetry factor increases the non-Lambertianity. Although the dependence on the angle of incidence is only clear at the largest FWA concentration for the measurement, the simulated BLDFs are more Lambertian at 45° than at 0° angle of incidence.

All of these observations can be related to the mean depth of the fluorescence process. Smaller mean free path at constant single scattering albedo means larger scattering and absorption coefficients, making the fluorescence process occur closer to the top surface. Increasing the FWA concentration (larger *k*) or the angle of incidence also leads to the light fluorescing closer to the medium surface. Although the fluorescence emission is isotropic, fluoresced light emitted at oblique angles has a larger probability to be scattered before it exits the medium. If the fluorescence process takes place closer to the medium surface, more light exits the medium at oblique angles. This is in analogy with a similar effect found in the anisotropic reflectance of turbid media with relatively large absorption [20, 21].

Our results are in accordance with the measurement reported by Holopainen et al. [6] who stated that the fluorescent emission from solid amorphous material is clearly non-Lambertian. The different results reported by Johansson et al. [7], who found a nearly Lambertian fluorescent emission from paper samples, can be explained by different concentration of fluorophores (FWAs) and angles of incidence, as Holopainen et al. used 0° and Johansson et al. 45° angle of incidence.

We notice that as the FWA concentration increases, the absorption coefficient naturally increases but the scattering coefficient decreases, especially for smaller values of *g*. This interdependence of the scattering and absorption coefficients has been the subject of debate in the paper optics community. Our results are in agreement with the the findings of Neuman et al. [22], who showed that most of the decrease of *s* with increasing *k* can be attributed to anisotropic reflectance and single scattering albedo.

Although the general trend is the same for the simulated and measured BLDFs in Fig. 8, the simulations with *g* = 0.8 show a larger decrease of the fluorescence with viewing angle than the measurements. Using smaller values of the asymmetry factor *g* leads to better agreement between simulations and measurements, and therefore a more accurate estimation of *g* is required to predict the angular distribution of fluorescence. The measurements are more accurate at smaller viewing angles than at larger viewing angles because of e.g. the effect of finite-size detection [23]. Another aspect neglected in this study is the effect of surface reflection that can potentially change the angular distribution of the fluoresced light. We believe however that this will have a minor impact on non coated papers, with no well defined surfaces, as compared to e.g. coated papers, which is the subject for future work.

## 5. Conclusion

We have measured the angular distribution of the fluorescence from paper samples with different amount of fluorescing whitening agents (FWA). We showed that the fluorescent emission goes from strongly anisotropic to almost Lambertian depending in the concentration of FWA. Using Monte Carlo light scattering simulations, we characterised the dependence of the fluorescence anisotropy on the optical parameters and angle of incidence. This analysis showed that the fluorescence angular distribution depends on the single scattering albedo and mean free path at excitation angle, and on asymmetry factor and angle of incidence, which explains the different results reported in the literature. Furthermore, the results further show that the angular distribution is strongly correlated to the mean depth of the fluorescence process. This means that the angular distribution of fluorescence must be taken into account when determining the quantum efficiency of fluorescing turbid media. Also, when the optical parameters are known, the depth of the fluorophores may be deduced from the angular distribution. A thorough understanding of these issues are important in the application of fluorescent dyes in paper, textile, and other appearance-enhancing or attention-grabbing usage of fluorescence. Future work includes accurate determination of the asymmetry factor and its wavelength dependence and modelling the effect of surface reflection on the BLDF.

## Acknowledgments

This work was supported by the Marie Curie Initial Training Networks (ITN) CP7.0 N-290154 funding and the Kempe Foundations SMK-1252, who are gratefully acknowledged. We acknowledge the valuable feedback from anonymous reviewers.

## References and links

**1. **L. G. Coppel, M. Andersson, O. Norberg, and S. Lindberg, “Impact of illumination spectral power distribution on radiance factor of fluorescing materials,” in *Proceedings of IEEE Colour and Visual Computing Symposium* (IEEE, 2013), pp. 1–4.

**2. **R. Donaldson, “Spectrophotometry of fluorescent pigments,” Br. J. Appl. Phys. **5**, 210–214 (1954). [CrossRef]

**3. **J. C. Zwinkels, D. S. Gignac, M. Nevins, I. Powell, and A. Bewsher, “Design and testing of a two-monochromator reference spectrofluorimeter for high accuracy total radiance factor measurements,” Appl. Opt. **36**, 892–902 (1997). [CrossRef] [PubMed]

**4. ** ASTM-E2153-01, “Standard practice for obtaining bispectral photometric data for evaluation of fluorescent color,” (2011).

**5. ** CIE, “CIE 15: Technical report: Colorimetry, 3rd edition,” Tech. Rep. (2004).

**6. **S. Holopainen, F. Manoocheri, and E. Ikonen, “Non-Lambertian behaviour of fluorescence emission from solid amorphous material,” Metrologia **46**, S197 (2009). [CrossRef]

**7. **N. Johansson and M. Andersson, “Angular variations of reflectance and fluorescence from paper - the influence of fluorescent whitening agents and fillers,” in *Proceedings of 20th Color and Imaging Conference* (IS&T, 2012), pp. 236–241.

**8. **S. Tominaga, K. Hirai, and T. Horiuchi, “Measurement and modeling of bidirectional characteristics of fluorescence objects,” in *Lecture Notes in Computer Science*, A. Elmoataz, O. Lezoray, F. Nouboud, and D. Mammass, eds. (Springer International Publishing, 2014), pp. 35–42. [CrossRef]

**9. **M. Neuman, L. G. Coppel, and P. Edström, “Angle resolved color of bulk scattering media,” Appl. Opt. **50**, 6555–6563 (2011). [CrossRef] [PubMed]

**10. **S. C. Gebhart, A. Mahadevan-Jansen, and W.-C. Lin, “Experimental and simulated angular profiles of fluorescence and diffuse reflectance emission from turbid media,” Appl. Opt. **44**, 4884–4901 (2005). [CrossRef] [PubMed]

**11. **L. G. Coppel, M. Andersson, and P. Edström, “Determination of quantum efficiency in fluorescing turbid media,” Appl. Opt. **50**, 2784–2792 (2011). [CrossRef] [PubMed]

**12. **P. Turunen, J. Kinnunen, and J. Mutanen, “Modeling of fluorescent color mixing by regression analysis,” in *Proceedings of Fourth European Conference on Colour in Graphics, Imaging, and Vision* (IS&T, 2010), pp. 94–100.

**13. **F. E. Nicodemus, J. C. Richmond, J. J. Hsia, I. Ginsberg, and T. Limperis, *Geometrical Considerations and Nomenclature for Reflectance* (National Bureau of Standards, 1977).

**14. **B. Bernad, A. Ferrero, A. Pons, M. L. Hernanz, and J. Campos, “Upgrade of goniospectrophotometer GEFE for near-field scattering and fluorescence radiance measurements,” Proc. SPIE **9398**, 93980E (2005).

**15. **L. Henyey and J. Greenstein, “Diffuse radiation in the galaxy,” Astrophys. J. **93**, 70–83 (1941). [CrossRef]

**16. **L. Wang, S. L. Jacques, and L. Zheng, “MCML–Monte Carlo modeling of light transport in multi-layered tissues,” Comput. Meth. Prog. Bio. **47**, 131–146 (1995). [CrossRef]

**17. **L. G. Coppel, P. Edström, and M. Lindquister, “Open source Monte Carlo simulation platform for particle level simulation of light scattering from generated paper structures,” in *Proceedings of paper making research symposium*, E. Madetoja, H. Niskanen, and J. Hämäläinen, eds. (Kuopio University, 2009).

**18. **P. Edström, “A fast and stable solution method for the radiative transfer problem,” SIAM Rev. **47**, 447–468 (2005). [CrossRef]

**19. **M. Neuman, L. G. Coppel, and P. Edström, “Point spreading in turbid media with anisotropic single scattering,” Opt. Express **19**, 1915–1920 (2011). [CrossRef] [PubMed]

**20. **M. Neuman and P. Edström, “Anisotropic reflectance from turbid media. I. Theory,” J. Opt. Soc. Am. A **27**, 1032–1039 (2010). [CrossRef]

**21. **M. Neuman and P. Edström, “Anisotropic reflectance from turbid media. II. Measurements,” J. Opt. Soc. Am. A **27**, 1040–1045 (2010). [CrossRef]

**22. **M. Neuman, L. G. Coppel, and P. Edström, “A partial explanation of the dependence between light scattering and light absorption in the Kubelka-Munk model,” Nord. Pulp Pap. Res. J. **27**, 426–430 (2012). [CrossRef]

**23. **N. Johansson, M. Neuman, M. Andersson, and P. Edström, “Influence of finite-sized detection solid angle on bidirectional reflectance distribution function measurements,” Appl. Opt. **53**, 1212–1220 (2014). [CrossRef] [PubMed]