HYPERBLEND: SIMULATING SPECTRAL REFLECTANCE AND TRANSMITTANCE OF LEAF TISSUE WITH BLENDER

: Remotely sensing vegetation condition and health hazards requires modeling the connection of plants’ biophysical and biochemical parameters to their spectral response. Even though many models exist already, the field suffers from lack of access to program code. In this study, we will assess the feasibility of open-source 3D-modeling and rendering software Blender in simulating hyperspectral reflectance and transmittance of leaf tissue to serve as a base for a more advanced large-scale simulator. This is the first phase of a larger HyperBlend project, which will provide a fully open-source, canopy scale leaf optical properties model for simulating remotely sensed hyperspectral images. Test results of the current HyperBlend model show good agreement with real-world measurements with root mean squared error around 1 ‰. The program code is available at https://github.com/silmae/ hyperblend .


INTRODUCTION
Simulating reflectance and transmittance spectrum of leaf tissue has been of great interest in recent decades. Understanding the connection between plant leaves' optical properties and their biochemical and biophysical properties allow monitoring of plant health and condition by means of hyperspectral remote sensing. Remote sensing allows examining plant properties of large areas and in sites that are difficult to access, such as jungle or mountains. In forestry, remote sensing has been used to identify tree species (Pölönen et al., 2018) and detect insect infestation (Näsi et al., 2015), for example.
In practice, simulators aim to produce spectral response similar to spectral measurements of leaves whose biochemical and biophysical properties are known. Inverting the problem, i.e., comparing remotely sensed vegetation images to simulation results produces predictions of the plants' biochemical and biophysical properties.
There exists a myriad of leaf optical properties models that simulate vegetation spectral response in different scales. Cell level models describe the intricate internal 3D-structure of a leaf e.g., (Ustin et al., 2001). Leaf level models consider the average concentration of photoactive molecules in leaves. Canopy scale models expand leaf level models to include statistical properties of large vegetation scenes, such as average leaf orientation and soil contribution. (Jacquemoud and Ustin, 2019) Perhaps the most used leaf-scale model to-date is the PROSPECT model, which was introduced in 1990 (Jacquemoud and Baret, 1990, Jacquemoud andUstin, 2019). The first version was capable of simulating leaf spectra using a few most important leaf pigments and water content. Many improvements and adaptations (Jacquemoud and Baret, 1990, Fourty et al., 1996, Feret et al., 2008, Féret et al., 2017, Jiang et al., 2021 have since increased the number of modeled chemical and physical constituents, as well as overall accuracy of the model. * Corresponding author A leaf scale model considers scattering and absorption of leaf tissue and provides data similar to a point spectrometer measurement. Combining it with a canopy model, such as SAIL (Verhoef, 1984), increases its usefulness in modeling whole forests (Jacquemoud, 1993, Jacquemoud et al., 2009. SAIL relies on statistical properties of trees, such as average leaf inclination angle and leaf area index. Statistical approach works well when modeling space-borne hyperspectral images where a single pixel is in the scale of dozens or hundreds of meters and can contain several tree-species, undergrowth, bare soil, and so on. Actual geometry of forests becomes more important when imaged from lower altitudes. For example, pixel size of imagers mounted on a drone may well be in centimeter scale. Canopy models that account geometry are scarce. DART (Gastellu Etchegorry et al., 1996, Gastellu-Etchegorry et al., 2004 can simulate approximate geometry as voxels (volumetric pixels) and has been used to simulate canopy reflectance of Norway spruce (Picea abies) in (Janoutová et al., 2019), for example. While free-to-use for academic research, DART's restrictive licensing terms hinder its applicability. In general, the field suffers from the lack of access to program code (Jacquemoud and Ustin, 2019). Even when the code is made available, insufficient documentation and poor user interface may render the code virtually unusable for others than the creators of the software.
In this study, we present the first phase of project HyperBlend (outlined in section 2). The final objective of the project is to provide an open-source simulator for hyperspectral images of randomly generated forests with various tree species, soil, and undergrowth. This first phase of the project assesses the feasibility of open-source 3D-modeling and rendering software Blender (blender.org) in simulating leaf spectra in the sense of reflection and transmission. In other words, we will show that it is possible to create virtual materials that can represent reflection and transmission properties of leaf tissue as if measured by a point spectrometer coupled to an integrating sphere. Results presented in section 3 show good agreement with realworld measurements. Future development of the HyperBlend simulator is outlined in section 4. Section 5 concludes the study. Thoroughly documented program code of the project is available in a public GitHub repository at https://github.com/ silmae/hyperblend.

MATERIALS AND METHODS
In 3D computer graphics, rendering refers to a process of producing a synthetic image of a virtual scene (Akenine-Möller et al., 2008). Software capable of rendering is called a render engine. Implementation of a render engine depends on intended use: in some applications, such as video games, speed is of essence and less realistic images are acceptable, whereas in animated films, photorealistic images are desired even if rendering them takes a long time.
In this study, we use open-source 3D-modeling software Blender, that ships with path tracing render engine called Cycles. Path tracing type render engine is currently a popular choice for photorealistic rendering. As a subset of ray tracing, path tracing is based on numerical approximation of Kajiya's rendering equation (Kajiya, 1986). Ray tracing methods in leaf optical models have been widely applied in recent decades (Govaerts and Verstraete, 1998, Ustin et al., 2001, Kallel, 2018, Kallel, 2020. The principle of path tracing is fairly simple. First, a ray (of light) is cast from a virtual camera pixel to the virtual scene. The ray scatters off surfaces of virtual objects until it reaches a light source when it is deleted and the illumination information it carried is saved to the pixel that cast it. The scattering event is called a light bounce. To avoid too long (or infinite) rendering time, the ray is also deleted after a certain number of bounces. The materials (surface properties) of objects are defined as BRDF/BSDF (Bidirectional Reflectance/Scattering Density Function). Multiple rays are cast from a single pixel and the final pixel illumination is the average of cast rays. Averaging reduces the noise inherent in the technique as scattering direction of each light bounce is determined by a Monte Carlo process.
Volumetric extension of path tracing takes into account scattering of rays inside an object-not just on its surface (Lafortune and Willems, 1996). In this study, we make heavy use of Cycles' volumetric scattering and absorbing capabilities, in fact, we do not use surface scattering at all. Ignoring surface scattering means that we assume leaf to be a Lambertian surface, which is a commonly used approximation (Jacquemoud and Ustin, 2019).
Cycles considers light in the fashion of computer graphics, i.e., it uses color models that mimic the way a human eye sees the light instead of using spectrum of photons of varying energy (wavelength). This prevents direct modeling of spectral phenomena. To render a spectral image, we must render each wavelength separately, as a gray-scale image as suggested in (Penttilä et al., 2021).
Regardless of this admittedly severe weakness, Blender has many positive qualities difficult to find in scientific software in general: stable, reliable and actively improved code base, complete user-interface enhanced with Python scripting, numerous online tutorials, and complete documentation. Most importantly for our final objective of producing large-scale forest scenes, it is very easy to create and manipulate 3D-objects. It is even possible to replicate objects procedurally to create non-identical, Figure 1. Schematic presentation of the rendering scene. The sunlamp is shown as a circle and its direction as a dashed line. Camera is represented by a triangle and cast rays as solid lines.
Absorbed rays end with a square and scattered rays with an arrow. The leaf volume is presented by dashed box. Top row shows the camera in reflectance imaging position and bottom row in transmittance imaging position. Left column shows leaf target imaging and right column reference imaging setup. Reflectance and transmittance white reference planes are shown as horizontal solid lines in right-hand column. more complex objects. For example, a forest can be modeled by first modeling a leaf, a branch, and a trunk, and then spawning rotated and scaled leaves on a branch that can be spawned on the trunk, which in turn can be spawned on ground. In this paper, however, we will only use simple objects.

Virtual Scene and Rendering
In this study, a leaf is modeled as a box whose volume is filled with absorbing and scattering particles. Behavior of the leaf material is controlled by four parameters: absorbing particle density ρa, scattering particle density ρs, scattering anisotropy α that is the ratio of forward scattering versus backward scattering, and mix factor β that controls the ratio of scattering and absorption events. Increasing particle density increases the probability of a ray undergoing a scattering or absorption event while traversing through the leaf. Zero particle density makes the leaf effectively transparent.
The rendering scene contains the leaf, a light source, and a camera. As a light source, we use what is called a sunlamp in Blender. It simply defines collimated light of certain intensity to every part of the 3D space, i.e., a static, homogeneous vector field. For reflection rendering, the camera is placed to the lit side of the leaf and for transmission rendering, it is placed to the unlit side. Used rendering scene and imaging setups are presented in Figure 1.
The rendered gray-scale images are saved to disk in TIFF format with 16-bit color depth. As in real-world photography, if the light is too strong, the image is overexposed. We use virtual white reference plates to manually adjust the intensity of the light to avoid overexposure while using as much light as possible for good signal-to-noise ratio. Used white references are rectangular planes of ideally diffuse (Lambertian) reflective or transmissive material. White references are also used in image normalization so that pixel values of the leaf are divided by respective pixel values of the white reference. This is done separately for reflection and transmission images.
To sum up the rendering of a reflectance-transmittance value pair, four images are rendered: leaf reflectance, reflectance white reference, leaf transmittance, and trasmittance white reference. The final reflectance value is obtained by dividing leaf reflectance image pixel-wise by reflectance white reference. The transmittance value is acquired in a similar manner, normalizing the transmittance image with the transmittance white reference.
The arrangement closely matches a spectral leaf measurement with integrating sphere. Even while the integrating sphere itself is not modeled, as long as the imaged area is smaller than the leaf object, there is no risk of light escaping the camera as the rays are cast from the camera-not from the light source. This convention is often used in ray tracing render engines as it is computationally less expensive and direction of propagation does not affect the path of light (Disney et al., 2000).

Optimization
Wavelength dependent leaf material parameters ρa, ρs, α, β are found by minimizing the difference of simulated reflectance r(ρa, ρs, α, β) to measured reflectance rm and simulated transmittance t(ρa, ρs, α, β) to measured transmittance tm. Term measured is used loosely here, as there is no difference if rm and tm originate from measurements or from some simulator that produces reflectance-transmittance pairs. The optimization loop for a single wavelength is illustrated in Figure 2. The optimizer solves the following minimization problem: min ρa,ρs,α,β (r − rm) 2 + (t − tm) 2 1/2 The sum of reflection and transmission is not allowed to exceed one to ensure conservation of energy. The lower limit of absorption particle density ρa and scattering particle density ρs excludes zero, to enforce at least some interaction between rays rm, tm r, t ρa, ρs, α, β  and the leaf material. The higher limit is not strictly needed, but fully constraining the problem allows using an optimizer that is intended for constrained problems. Scattering anisotropy α can vary between zero (maximal backscattering) and one (maximal forward scattering). Mixing factor β varies similarly from zero, when there is only absorption and no scattering, to one, when there is only scattering and no absorption.
Equation (1) defines a multiobjective optimization problem as reflectance and transmittance must be optimized simultaneously. For such a problem, there does not exist a single optimal solution. Instead, all solutions are trade-offs where increasing the optimality of one objective may decrease the optimality of the other. In this situation, a good starting guess is beneficial as not only does it converge faster to a solution, it is also more likely to converge into some certain solution. Considering the continuous nature of the problem, we would like the leaf material parameters to behave continuously as well. Now, if the starting guesses of adjacent wavelengths are close to each other, the leaf material parameters are more likely to be so too.
For finding a good starting guess, we can utilize the shape of leaves' reflection and transmission spectral curves, which tend to be similar in general. Letting approximate reflectancer and approximate transmittancet be equal, they can be expressed through absorptioñ Running the optimization forã ∈ [0, 1] and fitting a polynomial to each resulting material parameter spectra, we get estimates for leaf material parameters as shown in Figure 3. The coefficients of the polynomial fitting are used to retrieve a reasonable starting guess for any absorption value.

Experiment
For testing HyperBlend, we used reflectance and transmittance spectra of 28 silver birch (Betula pendula) leaf samples collected by (Hovi et al., 2017) in summer 2016. Each sample spectrum is an average of three spectral measurements. The data is publicly available at SPECCHIO database (Hueni et al., 2020), which is administered by Remote Sensing Laboratories, Institute of Geography, University of Zurich. Standard deviation of modeled and measured values are indicated by colored shadows and gray dashed lines, respectively. The mean of simulation aligns almost perfectly with the measurements except where transmittance is close to zero. Variation of the measured data is also captured by the model as can be seen from the alignment of the shaded area and grey dashed line.
The size of rendered images was 16 pixels (4 × 4) as that is the smallest image used Blender version (2.93.5) can produce. Wavelengths from 400 nm to 2300 nm were optimized with 5 nm spectral resolution. For 381 wavelengths and 28 samples, the optimization problem was thus solved 10668 times-each requiring dozens of render calls. Wall clock time per leaf sample was approximately 70 minutes using over 9 hours of processor time on a single 4-core (Intel i7-7700K, 4.20 GHz) central processing unit (CPU). Rendering was also done on said CPU, not on a graphical processing unit (GPU), as the overhead of CPU-GPU memory swaps outweighs GPU's faster rendering time on such tiny images with simple geometry.

RESULTS
Comparison of measured and modeled leaf reflection and transmission spectra of the 28 samples is shown in Figure 4. Modeled mean spectrum differs significantly from measured mean spectrum only in short wavelengths where the transmittance is nearly zero. Variation of the measured data is also captured by the model, as can be seen from the standard deviation shown in Figure 4. There is large variation in 500-700 nm range, which is expected, as the data contains both adaxial (upper surface) and abaxial (lower surface) leaf side measurements where great variation has been reported in visible region (Hovi et al., 2017). Our model incorporates this variation mainly in the scattering anisotropy parameter.
The total root mean square error (RMSE) over the whole spectrum and all 28 samples was 0.0010 for reflectance and 0.0008 for transmittance. Wavelength-wise errors averaged over all Figure 5. Average root mean square error between measured and simulated reflectance and transmittance over the 28 samples. The error of reflectance behaves somewhat systematically in visible spectrum and around the second water absorption peak at 1940 nm. The error of transmittance is high on short wavelengths but still less than 6‰ of the maximum transmittance of 1. On average, error of reflectance is 1‰ for reflectance and 0.8‰ for transmittance.
samples are plotted in Figure 5. The plot of reflectance error mostly shows appropriately random behavior except in the range of visible light and near the second water absorption peak at 1940 nm where systematic behavior can be seen. Error of transmittance is less systematic but several times greater than average near 400 nm, which indicates that our model cannot fully handle the situation where transmission is almost zero (see also Figure 4). One reason for this might be that, as shown in Figure 3, the starting guess does not quite catch the behavior when absorption is close to one, which may hinder the optimization process. Despite the systematic elements in the error, we consider the model accurate, as even the greatest error is less than 0.006, i.e., 6‰ of the theoretical maximum value of 1. In practice, if reflectance or transmittance reached the value one, the leaf would be either fully opaque or completely transparent.
The leaf material parameters shown in Figure 6 behave continuously and are consistent with the shape of reflectance and transmittance curves as expected. For example, the first water absorption peak near 1500 nm is mirrored as a sharp increase in absorption particle density. Mix factor shows similar shape, which is natural since greater mix factor means that absorption events are favored over scattering events. Scattering anisotropy favors forward scattering, varying around 0.8 in near infrared. Value 0.5 would mean equal amount of forward scattering and back scattering, which is the case only around 500 nm. Scattering particle density shows the least variation over the modeled spectrum, which is somewhat surprising, but it still follows the shape of the reflectance curve.

DISCUSSION
In light of the results, the proposed method can accurately model physical world in its intended use, but it is significantly slower than other existing methods. For example, generating reflectance and transmittance spectra from 400 nm to 2500 nm with 1 nm resolution using web-implementation (http: //opticleaf.ipgp.fr/) of PROSPECT takes only a few seconds compared to 70 minutes using our method. Considering that we simulated five times less wavelengths, our model is slower by four orders of magnitude.
Comparison of simulation speed between PROSPECT and current version of HyperBlend is not as meaningful as it might seem. Ray tracing is meant for tracing light paths in complex scenes with three-dimensional geometry. In this study, we have reduced ray tracing into essentially one-dimensional model for the sake of assessing whether or not it can replicate light scattering in leaf tissue accurately enough to be useful. On the other hand, PROSPECT inherently lacks any spatial information and even if extended into canopy model using, e.g., SAIL, no real geometry is included.
If we hypothesize that HyperBlend will be able to simulate hyperspectral images of forest scenes with real geometry, we can construct a more meaningful comparison with an existing model DART, which also utilizes real geometry. In (Janoutová et al., 2019) DART was used to simulate 926 wavelengths in 400 pixel (20 × 20) hyperspectral image of spruce forest scene with 10 trees. Simulation times ranged from few hours on 2 processor cores to 34 hours on 8 processor cores depending on the detail level of the scene. Based on our preliminary experiments of rendering forest scenes with almost hundred trees and leaf material like the one used in this study, render time of a single wavelength 0.5 megapixel (1028 × 512) is less than half a minute on a common desktop computer with 8-core processor and several years old consumer grade GPU. Rendering this more complex scene (in terms of number of trees) in vastly higher resolution using the same number of wavelengths would then take less than eight hours. With this, we are optimistic that HyperBlend will be at least on par with existing models in terms of simulation speed.
Regardless of the future versions of HyperBlend, current version can also be improved. Speeding up the model considerably requires omitting the slow optimization process where, in the worst case, hundreds of rendered images are needed to simulate reflectance and transmittance of a single wavelength. Instead of modeling each reflectance and transmittance value accurately, approximate values might suffice.
As the parameter spectra was shown to behave continuously when reflectance and transmittance do so too, one could generate evenly spaced reflectance-transmittance pairs and optimize respective material parameters. This would form four parameter surfaces with respect to reflectance and transmittance, one for each material parameter. By fitting a surface function to each parameter surface, the material parameter retrieval would be reduced into a function evaluation. We estimate that this approximation could make the simulation 10-100 times faster. In addition to faster parameter retrieval, this strategy would allow to render images for all wavelengths on a single Blender execution reducing processor context switches and Blender initializations from thousands to one. Accuracy of such model would depend on how well selected surface fits on each parameter surface. If the material parameters do not behave smooth enough, the surface may be significantly less accurate than the optimization strategy. In that case, training a deep neural network could be considered (Imaizumi and Fukumizu, 2019).
As a final note, binding leaf material parameters used in Blender directly to reflectance and transmittance values allows coupling any leaf model that produces reflectance-transmittance pairs with HyperBlend. It also allows mixing measured and modeled values in a single scene, i.e., one tree could have leaf material derived from PROSPECT simulation, another one could use material parameters based on a field study. By expanding HyperBlend into non-transmissive materials, one could include roofs of buildings or other man-made objects into the scene.

CONCLUSIONS
In this study, we presented the first phase of a new leaf optical properties model HyperBlend, which can simulate leaf reflectance and transmittance as if measured by a point spectrometer coupled to an integrating sphere. Moreover, we showed that open-source 3D-modeling and rendering software Blender can be used in simulating spectral data, even though it cannot inherently handle spectral phenomena. Our model used just four material parameters to depict the behavior of light interacting with leaf tissue. Simulated reflectance and transmittance showed good agreement with real-world measurements.
The performance of the model in its current state was discussed. Ideas on future development of HyperBlend and estimates of its performance were presented.

FUNDING
This study was funded by Academy of Finland (327862).
Leevi Lind for numerous read-troughs and comments on the manuscript as well as all the discussions and brainstorming sessions throughout the development of this first version of HyperBlend.