Two-scale Monte Carlo ray tracing for canopy-leaf vector radiative transfer coupling

Two-scale Monte Carlo ray tracing for canopy-leaf vector radiative transfer coupling

Journal Pre-proof

Two-scale Monte Carlo ray tracing for canopy-leaf vector radiative transfer coupling Abdelaziz Kallel PII: DOI: Reference:

S0022-4073(19)30685-5 https://doi.org/10.1016/j.jqsrt.2019.106815 JQSRT 106815

To appear in:

Journal of Quantitative Spectroscopy & Radiative Transfer

Received date: Revised date: Accepted date:

23 September 2019 19 December 2019 24 December 2019

Please cite this article as: Abdelaziz Kallel, Two-scale Monte Carlo ray tracing for canopy-leaf vector radiative transfer coupling, Journal of Quantitative Spectroscopy & Radiative Transfer (2019), doi: https://doi.org/10.1016/j.jqsrt.2019.106815

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier Ltd.

Highlights • Vector radiative transfer to model leaf-canopy polarized reflectance is proposed • Decomposing incident light into a basis of polarization allows to derive the Stokes phase function at leaf level • Discretizing wavelength domain, incident and exiting angles allows to create a database of Stokes phase function • Monte Carlo weighted sampling allows to predict the ray tracing from one wavelength to another on leaf and canopy levels • Polarization and 3-D structure of leaves are not negligible

1

2

Two-scale Monte Carlo ray tracing for canopy-leaf vector radiative transfer coupling

3

Abdelaziz Kallela,∗

1

a Centre

4

de Recherche en Num´ erique de Sfax, Technopole de Sfax, PO Box 275, 3021 Sfax, Tunisia

5

Abstract

6

Canopy reflectance simulation was widely developed using Monte Carlo (MC) ray tracing. Nevertheless, the

7

proposed models assume leaf as Bi-Lambertian medium. Based on the recently developed leaf MC ray tracing

8

code, this study proposes to couple MC ray tracing simulations on leaf and canopy levels. Firstly, at the

9

first level, the bidirectional scattering distribution function (BSDF) is computed for unpolarized incoming

10

rays as well as some common states of polarization producing together a polarization decomposition basis,

11

allowing to derive the Stokes phase function. Moreover, discretizing both incident and scattering angles over

12

the sphere allows to produce a database of all possible Stokes phase functions. Secondly, at the canopy

13

level, the reflectance is simulated using ray tracing technique, when ray is intercepted, the scattering is done

14

considering the appropriate Stokes phase function in the leaf database. Simulation of multiple wavelengths

15

is accelerated based on a new MC weighted sampling technique permitting to consider the same tracing for

16

all the wavelengths together. Simulation results show the relevance of such modeling compared to traditional

17

models. Indeed, from one side the canopy bidirectional reflectance depends on the leaf BSDF thus if leaves are

18

assumed Bi-Lambertian surfaces leads to inaccurate results, and from the other side, the reflectance is sensitive

19

to polarization and neglecting it affects the results mainly when the incident light is polarized. Comparison

20

with actual polarized reflectance measurements presented in literature shows good agreement with same trends

21

and variation ranges. Quantitative validation of our canopy level model is done using the ROMC web-tool.

22

Reflectance RMSE between our simulations and the ROMC reference is lower than 0.02.

23

Keywords: Canopy polarized reflectance, Leaf polarized scattering, Vector radiative Transfer, Ray tracing,

24

Monte Carlo simulation.

25

1. Introduction

26

Nowadays, using remote sensing technology to monitor crops and forest growth, phenology and health has

27

become common practice. Semi-empirical and physical based models are used to link the observations to the ∗ Correspending author. Preprint Journal of Quantitative Spectroscopy & Radiative Emailsubmitted address: to [email protected] (Abdelaziz Kallel) Transfer

December 30, 2019

1

INTRODUCTION

2

1

canopy and then leaf features. Although the latter are the most complex, they provide generally the most

2

accurate results with easier generalization from a study case to larger or other area because there is no need

3

to calibrate empirical parameters and only a priori about the physical properties of the scene are needed.

4

At canopy scale, radiative transfer (RT) theory is widely used to model the physical phenomena such

5

as extinction and scattering. This theory has been developed the first time to model rotationally invariant

6

medium such as dust having spherical form, with an extinction coefficient independent of the propagation

7

direction and scattering coefficient depending only on the angle between incidence and scattering [1]. Then,

8

in [2], authors proposed the first adaptation of RT to vegetation case. By assuming that leaves are null-size

9

bi-Lambertian facets with random angular distribution, it was possible to redefine the extinction and the

10

bidirectional scattering coefficients that depend on the direction of propagation and both incidence and exiting

11

directions, respectively [3]. This development was not enough to mimic all the features of the bidirectional

12

reflectance distribution. In particular, with the apparition of the multi-angular sensors such as POLDER and

13

PROBA-V that showed that there is an increase in the reflectance in both backward and forward directions,

14

that cannot be explained based on the original RT adaptation [4, 5]. In backscattering, a sharp peak around

15

the sun direction is observed, it is called the hot spot effect and it happens when the sun, sensor and the

16

observed scene are aligned. In this case, all the observed components are sunlit and scene appears bright,

17

then the further apart are the sun-sensor orientations, the more shadows are observed and the darker is the

18

scene. The hot spot effect is due to the non-null size of the medium particles, it was modeled coupling RT and

19

geometrical models [6, 7]. Concerning the canopy remarkable brightness in the forward direction, it is mainly

20

due to the specular effect. In fact, it appears clearly for low sensor scanning when the scattering phase angle

21

(i.e. angle between the sensor and the sun) is too high, forming therefore a bright spot around the forward

22

direction. In these directions, green vegetation among other canopy materials can appear white [8]. To take

23

this effect into account, leaves are assumed partially scattering in bi-Lambertian fashion and partially reflecting

24

in the specular direction [8]. Moreover, it is well-known that the specular light is highly polarized [8, 4, 9],

25

therefore to be taken into account, the vector RT (VRT) theory allowing to deal with the light Stokes vector

26

has been adapted [10, 11].

27

Besides, first RT models cannot allow to mimic spatial reflectance variation in very high-resolution image

28

where it is possible to distinguish individual trees in forest or rows of crops in fields. In fact, these models could

29

not allow to take into account the 3-D heterogeneity of the vegetation medium, since (semi-)analytical solution

1

INTRODUCTION

3

1

of the RT equation with variable density of vegetation is impractical. To overcome this problem, numerical

2

solutions based on Monte Carlo ray tracing [12, 13] or iterative discrete ordinate [14] have been proposed.

3

At leaf scale, first models inspired from the Kubelka-Munk (K-M) theory, are called the plate models

4

[15, 16]. They assumed that leaves are composed of N compact layers and they derived the reflectance and

5

transmittance coefficients by linking the K-M scattering coefficients to the medium refractive index variation.

6

Then, the concentration of the different components was taken latter into account in the derivation of the leaf

7

scattering coefficients in the well-known PROSPECT model [17]. Such a model produces just one coefficient

8

of reflectance and another for transmittance which does not allow to model their bidirectional variation as

9

shown in experimental results [18]. Therefore, as for vegetation canopy, the MC ray tracing techniques have

10

been applied to produce more realistic simulations of the radiative regime [19, 20, 21]. Moreover, like for

11

vegetation, laboratory experiments [22, 9] showed that bidirectional polarization effect is too important for

12

all the scattering directions and not just a point-like in the forward direction. These results were recently

13

confirmed based on VRT theory [23]. It was particularly shown that the difference between scalar and vector

14

RT is important in the case of polarized incident light.

15

In order to understand the relationship between the remote sensing observation and the leaf properties,

16

it is important to couple radiative transfer model within leaf and canopy. Generally leaf 1-D model such as

17

PROSPECT is coupled with 1-D or 3-D models [24, 11]. As PROSPECT does not specify the scattering

18

bi-directional distribution, it is assumed to be bi-Lambertian. In order to take into account the specular effect,

19

just a peak of reflectance is added in the specular direction [10]. Such assumptions neglect the leaf scattering

20

anisotropy that was identified in actual measurements [18, 9] and proved theoretically [25, 23]. Therefore, the

21

existing coupled models are not perfectly representing the radiative regime within the canopy. This issue can

22

prevent accessing to interior leaf mechanisms such as photosynthesis, as they need precise interpretation of the

23

canopy reflectance spectrum [26].

24

In this paper, it is proposed to couple leaf [23] and canopy [10] 3-D VRT models. At the leaf level,

25

such coupling requires to estimate the leaf polarized bidirectional scattering function allowing to derive the

26

Stokes phase function, and to estimate the latter for any incident and scattering directions. For this purpose,

27

optimized techniques are proposed to discretize the bidirectional scattering over the sphere. At the canopy level,

28

the ray tracing technique proposed in [10] is adapted to take into account the leaf bidirectional scattering and

29

polarization variations. Moreover, to allow simulation of multiple wavelengths, the running time is optimized

2

1

THEORETICAL BACKGROUND

4

by developing a new MC sampling technique that simulates ray tracing on only one wavelength.

2

This paper is organized as follows. Section 2 recalls the theoretical background. Then, Section 3 outlines

3

our approach. After that, Section 4 shows our simulation results and validation. It is followed by Section 5

4

discussing the new outcomes. Finally, Section 6 summarizes our conclusions and gives perspectives.

5

2. Theoretical background

6

In this section, we introduce the light representation using Stokes vector then we summarize RTM at canopy

7

and leaf levels that will be used to develop our coupled model.

8

2.1. Light representation

9

b , the component on n b For a transverse electric (TE) plane wave E that propagates in a given direction, n

10

is null. Its state of polarization can be represented within a coordinate system by the four-component Stokes

11

vector I defined as follows [27]





 I       Q    I= ,    U      V

12

(1)

13

where I corresponds to the wave radiative intensity, Q represents the degree of parallel (i.e. Q > 0) and

14

perpendicular (i.e. Q < 0) polarizations, U models the degree of diagonal polarization, and V stands for the

15

left-hand (i.e. V < 0) and right-hand (i.e. V > 0) circular polarizations.

16

17

18

19

20

21

22

23

24

b−n b . Then, when reaching the In radiative transfer I is written in the spherical coordinate system, θb − φ

interface between two media, and in order to compute the reflected and refracted radiations, it is comment to b , where b write I in the parallel-perpendicular coordinate system, b l − rb − n l and rb are in the plane of scattering b and the normal to the interface) and perpondicular to it. (it contains n

b−n b ) and the second (b b ) is η, thus Assume that the rotation angle between the first system (θb − φ l − rb − n

to pass from the Stokes vector expression written in the first (Iθ− b φ b ) to the second (Ib l−b r ), one has to use the Stokes rotation matrix [28]

Ibl−b b φ b, r = L(η)Iθ− where L is the Stokes rotation matrix with expression given in [28].

(2)

2

1

2

THEORETICAL BACKGROUND

5

Note that to pass from the second to the fist coordinate sytems, on has to use L(−η). p For a monochromatic TE plane wave, one has I = Q2 + U 2 + V 2 and the wave is called perfectly polarized.

3

Conversely, for non-monochromatic wave, one has I 2 ≥ Q2 + U 2 + V 2 . Thus, one can define the degree of

4

polarization, P , as P =

5

6

7

8

9

10

11

12

p

Q2 + U 2 + V 2 . I

(3)

2.2. Interface scattering properties Let us consider two homogeneous media (medium 1 and medium 2) with different refractive indexes (n1 b l . An incident wave, Ii , propagates within and n2 , respectively), separated by a plane interface of normal n b i and reaches the interface with medium 2 with zenith angle θi . It produces two the medium 1 in direction n new waves, the reflected one, Ir , that propagates back into the first medium in the specular direction (b nr ) and

b t . The relationship linking the the transmitted one, It , that propagates in the second medium in direction n reflected and transmitted Stokes vector to the incident one are [23],

Ir

= Fr (θi )Ii ,

It

=

13

n2 Ft (θi )Ii , n1

(4)

14

where Fr (resp. Ft ) is the reflection (resp. transmission) Stokes scattering function whose expression is given

15

in [23].

16

Fr is defined for Ii and Ir written in the coordinate system relative to the parallel/perpendicular compo-

17

nents. In radiative transfer, Ii and Ir are written in the spherical coordinate system. Therefore, the Stokes

18

scattering function is replaced by the Stokes phase function, Zr , written as

19

br, n b l ) = L(−ηr )Fr (θi )L(ηi ). Zr (b ni , n

(5)

20

where L is the Stokes rotation matrix [cf. Eq. (2)], ηi (resp. ηi ) is the angle between the spherical coordinate

21

system and the incident (resp. reflectd) wave parallel-perpendicular coordinate system.

22

For details about the geometrical representation of the different angles, reader is referred to [23].

23

The same expression can be found for the Stokes scattering function relative to transmittance Zt .

2

1

THEORETICAL BACKGROUND

6

2.3. Radiative transfer modeling

2

Polarized VRT (PolVRT) [10, 23] is a MC ray tracing model used to simulate radiative regime within leaf

3

and canopy. It consists of propagating a large number of rays on a given mock-up and then tracking them

4

interception after interception until exiting the medium or vanishing. The mock-up can be a vegetation medium

5

[10] or only a single leaf [23]. In both cases, the mock-up is decomposed into voxels to accelerate the tracing. In

6

fact, the time-consuming task is the interception finding since the code has to test every possible intersection

7

between rays and components, therefore it is better to do the search over few components within small voxels

8

rather than within all the medium.

9

2.3.1. Canopy level

10

PolVRT assumes that the vegetation layer is composed of disk-shaped leaves. Leaf positions can be gen-

11

erated randomly or exported to describe homogeneous or heterogeneous vegetation layer. Other elements can

12

be added to the medium such as soil, tree trunk and branch, etc. Sun radiation reaching the scene is simu-

13

lated as a large number of rays reaching the top of the scene. Each ray is traced voxel by voxel until finding

14

the first interception by the different medium elements then a part of the light is scattered. It is done in a

15

traditional Lambertian way if the scatterer is not a leaf. Otherwise, two phenomena are taken into account:

16

(i) bi-Lambertian diffuse reflectance and transmittance; (ii) specular reflection. The first is due to rays that

17

are scattered many times inside the leaf, whereas the second is due to its epidermis scattering. In terms of ray

18

tracing, these two phenomena compete so that after interception, one ray is propagated either by Lambertian

19

scattering or specular reflection. If the former case happens, ray is propagated randomly according to the

20

Lambertian law, whereas if the latter happens, ray is propagated exactly in the specular direction. In order to

21

compute reflectance efficiently, after each interception a secondary ray is sent in direction of the sensor. If it

22

is not intercepted again it is accounted in the radiance. Such technique is called ‘photon spreading’ [29] and

23

it can be directly used only in the bi-Lambertian case as the probability of scattering ray in any direction is

24

non-null. However, in the second case, the exiting direction is fixed by the Fresnel principle, and generally it

25

is not possible to send rays in the sensor direction.

26

PolVRT showed that in the forward direction, canopy polarizes horizontally light since leaves responsible of

27

specular effect are horizontal and the polarization is parallel to the scatterer surface. In addition, it showed that

28

in inclined forward direction, diagonal polarization can be observed. In fact, leaves in this case are inclined

2

THEORETICAL BACKGROUND

7

(b) Epidermis

(a) Leaf mockup: 228µm × 194µm × 175µm.

(c) Spongy

(d) Palisade

Figure 1: Leaf structure mock-ups [23].

1

in order to fit the sensor direction, and so the polarization state. For that neither horizontal nor diagonal

2

components are null.

3

The limit of this model is that it does not take into account leaf polarized bidirectional scattering distribu-

4

tion.

5

2.3.2. Leaf level

6

7

8

9

10

11

12

13

PolVRT was recently adapted to the leaf medium case. As proposed firstly in [19] and then revisited in [21], leaf is composed principally by 4 layers. They are presented in Figure 1: • Top and bottom epidermis: they are similar layers composed of compact cells having ellipsoidal shape, which fit one another so that there is no intercellular space; • Palisade parenchyma: it is below the top epidermis and composed of narrow vertical cells of cylindrical form; • Spongy mesophyll: it is between bottom epidermis and palisade layers. It is a heterogeneous layer filled by spherical cells of different diameters occupying randomly almost the half of the allocated space.

2

THEORETICAL BACKGROUND

8

Table 1: Geometrical parameters used to generate leaf mock-up.

Tissue Epidermis

Palisade

Spongy

Parameter Cell radius Cell oblateness Cell wall thickness Cell radius Cell height Cell wall thickness Chloroplast number per cell Chlorophyll concentration Chloroplast size (axis (x, y, z)) Vacuole radius Vacuole height Maximum Cell radius Maximum Cell wall thickness Faction of airspace Chloroplast number per cell Chlorophyll concentration Chloroplast size (axis (x, y, z)) Maximum Vacuole radius Tissues thickness

Value 12.5 µm 0.25 1 µm 7.5 µm 70 µm 1 µm 243 0.7 × 106 g.m−3 2 × 1.35 × 0.4µm3 5.4 µm 65.8 µm 12 µm 1 µm 0.45 µm 102 106 gm−3 1.8 × 1.8 × 0.5µm3 9.6 µm 70 µm

1

The cells of the different layer are covered by thin walls. The cell interior, called vacuole, filled by water

2

and it occupies the most important part of the cell volume (from 40 up to 90%). Moreover, palisade and

3

spongy cells contain the different pigments such as chlorophyll and carotenoid. They belong to chloroplasts

4

which are ellipsoidal shaped organelles. Figures 1(c) & 1(d) show that they are randomly distributed within a

5

thin structure between cell wall and vacuole called cytoplasm.

6

Table 1 presents the geometrical properties of the different layers of the leaf mock-up. The oblateness of

7

the epidermis cell corresponds to the ratio between the cell vertical and horizontal radii. In our case, we choose

8

a small value (0.25) corresponding to smooth surface approximating plate epidermis [25].

9

Two physical phenomena are taken into account when ray is traveling within a given leaf; (i) decrease

10

by absorption within the different tissues that is governed by the absorption coefficients; (ii) reflection or

11

refraction in interface between two different tissues that are governed by the Stokes phase functions [Zr and

12

Zt , cf. Eq. (5)]. The reflection and refraction modify the light polarization since the parallel and perpendicular

13

components of the wave are not scattered with the same proportion. Rays travel generally a large number

14

of times from a tissue to another before exiting the medium, leading to large running time, different MC

15

optimization techniques are proposed. In particular, the average ray effect after long trip is predicted rather

16

than computed. Also, a new MC weighted sampling technique is proposed to simulate close wavelengths using

17

the same tracing. Alike, multi-angular simulations are gathered using the same tracing. Experimental results

3

PROPOSED APPROACH

9

1

show that leaves are not bi-Lambertian surfaces with thin peak of specular reflection, but the reflectance

2

increases with the sun zenith angle with large peak of specular reflection in the forward direction. Moreover,

3

as in canopy case, leaves highly polarize light in horizontal and horizontal-diagonal directions in the forward

4

and inclined-forward directions, respectively.

5

For a given incident irradiance of Stokes vector Is of direction Ωs , PolVRT allows to simulate exitance and

6

radiance within a given sphere tessellation {∆Ωn }n=1,...,Nt . The exitance set is called {Id (∆Ωn )}n=1,...,Nt and

7

the corresponding radiance set is {L(∆Ωn )}n=1,...,Nt

8

3. Proposed approach

9

Our approach uses PolVRT at two scales separately. At leaf scale, it produces a data set of Stokes scattering

10

functions that cover all possible cases of incident light polarization and orientation. Then, at canopy scale,

11

when radiation is intercepted by a given leaf, the latter data set is used to produce the exiting radiation as

12

well as to estimate radiance in sensor direction.

13

3.1. Leaf polarization modeling

14

Using the PolVRT model at leaf scale allows to estimate the scattering distribution function. The letter

15

depend on the incident light polarization. As there are infinite possibilities of polarizations varying Q, U and

16

V , it is important to be able to deal with all of them. To this end, the Stokes scattering function, F, and the

17

Stokes scattering function relative to radiance, FL , are required. Appendix A explains their derivation.

18

F and FL are the only radiative parameters at leaf level that will be used in canopy level to simulate ray scat-

19

tering by leaves. In geometrical point of view and in terms of incidence, they depend only on the ray incidence

20

angle, θ, as leaf surface is assumed rotationally invariant. In terms of exiting, as the scattering directions are dis-

21

cretized over a tessellation of the sphere {∆Ωn }n=1,...,Nt , the MC ray tracing simulation at a given wavelength,

22

λ, and incidence angle, θ, gives rise to the two sets {F(λ, θ, ∆Ωn )}n=1,...,Nt and {FL (λ, θ, ∆Ωn )}n=1,...,Nt .

23

Simulation at vegetation level requires Stokes scattering function at different incidence angles and wave-

24

lengths. The computation at leaf level is time consuming as it needs MC ray tracing, making the global

25

simulation impracticable since it requires MC ray tracing at two levels. To solve this problem, we propose

26

to generate a database containing a fine discretization of all the possible configurations in terms of incidence

27

angles and wavelengths.

3

PROPOSED APPROACH

10

1

When ray tracing is simulated at vegetation level, the ray incidence on leaf (θ) is random and can be any

2

value from 0◦ to 90◦ , and the exiting direction could be any one within the sphere. F and FL are already

3

discretized on the scattering direction, which is not the case for the incidence. We propose therefore to sample

4

it with a step of 2◦ from 0◦ to 88◦ , producing the set {θq }q=1,...,Nθ . To be fast, MC ray tracing is done only

5

for the incidence set {0◦ , 15◦ , 30◦ , 45◦ , 60◦ , 75◦ , 88◦ }, then the results are extended to all the set using the

6

physically-based interpolation technique described in [23].

7

To mimic actual observations, simulated canopy reflectance must cover in general a large number of wave-

8

lengths, therefore it is necessary to simulate leaf bi-directional scattering for the same wavelengths. In our

9

case, as the simulation at canopy level are not known a priori, therefore we propose a fine simulation of the

10

wavelengths with step equal to 1nm from 400nm to 2400nm, producing the set {λp }p=1,...,Nλ . For running time

11

optimization sake, MC ray tracing within the leaf is done with step equals 100nm, then using MC weighted

12

sampling at leaf scale the results are interpolated to the required spectral resolution.

13

At the end, databases of exitance and radiance Stokes scattering functions are produced

14

{F(λp , θq , ∆Ωm )}p=1,...,Nλ ;q=1,...,Nθ ;m=1,...,Nt and {FL (λp , θq , ∆Ωm )}p=1,...,Nλ ;q=1,...,Nθ ;m=1,...,Nt , respectively.

15

3.2. Vegetation MC ray tracing

16

3.2.1. Stokes phase functions

17

PolVRT at canopy level propagates rays in a given wavelength, when it reaches the leaf the Stokes scattering

18

functions are applied to derive the exitance and radiance. As F and FL are not necessary defined for this

19

wavelength, the closest discrete wavelength in the database is considered. Similarly, for the incidence angle,

20

the closest discrete angle is considered. The exiting direction is chosen randomly among the discrete direction

21

relative to the sphere tessellation ({∆Ωm }m=1,...,Nt ).

22

b i of wavelength λ, Ii (λ, n b i ), is reaching a leaf of normal n b l . To compute Assume that a ray of direction n

23

the scattered Stokes vector, Id , one needs to derive the Stokes phase function. For that, one has to find

24

the corresponding Stokes scattering function (F). Therefore, F discrete wavelength λp∗ ∈ {λp }p=1,...,Nλ and

25

incidence θq∗ ∈ {θq }q=1,...,Nθ have to be determined. λp∗ that is the closest to λ is considered, similarly θq∗

26

(angle between incidence and the leaf normal) that is the closest to the actual incident angle [i.e. angle between

27

b i )] is taken. The exiting directions remain discrete the direction of propagation and the normal to the leaf (b nl , n

3

1

PROPOSED APPROACH

over the sphere tessellation, the exitance is given by

b s , ∆Ωm , n b l )Ii (λ, n b i ), m = 1, . . . , Nt . Id (λ, ∆Ωm ) = Z(λ, n

2

3

11

(6)

where Z is the exitance Stokes phase function

b s , ∆Ωm , n b l ) = L(ηm )F(λp∗ , θq∗ , ∆Ωm )L(ηi ), Z(λ, n

4

(7)

5

where L(ηi ) and L(ηm ) are the Stokes rotation matrix for incident and exiting directions, allowing to pass from

6

the spherical coordinate system to the leaf one and conversely, respectively [cf. Eq. (5)], and ηi and ηm are

7

the corresponding rotation angles, respectively.

8

9

b o , L(λ, n b o ), to the incident In the same way, it is possible to link the exiting radiance in direction n

irradiance,

b o ) = ZL (λ, n b s, n b o, n b l )Ii (λ, n b i ), Ld (λ, n

(8)

b s, n b o, n b l ) = L(ηo )FL (λp∗ , θq∗ , ∆Ωm∗ )L(ηi ), ZL (λ, n

(9)

10

11

where ZL is the radiance Stokes phase function

12

13

14

15

b o ∈ ∆Ωm∗ and ηo the rotation angle from the leaf coordinate system to the spherical m∗ is chosen such as n one.

3.2.2. Radiance derivation

16

In the following, we explain how to derive the radiance angular distribution using Monte Carlo forward ray

17

tracing. We start by estimating irradiance reaching leaves at each transition and then compute the produced

18

radiance. (i)

19

Assume that Nr ray tracings are performed, the incident ray Stokes vector set is {Is }i∈{1,...,Nr } . They

20

are assumed starting from the top of canopy at positions {b rs }i∈{1,...,Nr } , respectively. In each tracing, T (i)

21

transitions are done before exiting the medium or energy becomes negligible. For a given tracing i and transition

22

t, the Stokes vector is noted It

23

(i)

(i)

(i)

and it starts after t − 1 interceptions at position rbt−1 and continues until (i)

b t represents the direction of the ray in transition t, it is the next interception at position rbt . Henceforth, n

3

1

2

3

4

5

6

7

8

PROPOSED APPROACH

12 (i)

(i)

(i)

b t is no other than the sun ray direction therefore the unit vector linking rbt−1 to rbt . Note that, when t = 1, n (i)

bs . n

(i)

(i)

(i)

(i)

(i)

Let It (λ, rb1 , . . . , rbt ) be the radiation obtained by multiple scattering of the incident ray Is (λ, rbs ) (i)

reaching initially the canopy at position rb1

(i)

= rbs .

At the end of transition t, the ray is intercepted

b l,t , then a part of the ray intensity is absorbed and another one is scattered ranby a leaf of normal n

domly among the possible discrete directions {∆Ωn }n=1,...,Nt and according to the Stokes phase functions (i) b t−1 , ∆Ωn , n b l,t )}n=1,...,Nt . To do so, all the Stokes vectors, I˜t (b {Z(λ, n nt−1 , ∆Ωn ), corresponding to all the

discrete directions, {Ωn }n=1,...,Nt , are computed,

(i) (i) (i) (i) b t−1 , ∆Ωn , n b l,t )It−1 (λ, rb1 , . . . , rbt−1 ), n = 1, . . . , Nt . I˜t (b nt−1 , ∆Ωn ) = Z(λ, n

9

(10)

10

The choice among the directions is done randomly and proportional to the intensity. Moreover, to conserve

11

energy, the chosen radiation intensity is put as the global scattering energy. The complete derivation of

12

13

(i)

(i)

(i)

It (λ, rb1 , . . . , rbt ) is given in Appendix B.

According to the ‘photon spreading’ technique, at each transition t, a secondary ray is sent to the sensor (i)

14

direction [b no (θo , φo )], if it is not intercepted it contributes to the reflectance. The radiance Stokes vector Lt

15

at a given wavelength λ reaching the top of the canopy is

(i)

17

18

19

(i)

(i)

(i)

b o ) = 1n rt ) Lt (λ, rb1 , . . . , rbt , n b o (b

16

|b nt .b no | (i) (i) (i) b t, n b o, n b l,t )It (λ, rb1 , . . . , rbt ), ZL (λ, n cos(θo )

(11)

(i)

where ZL is the scattering phase function relative to the radiance and 1n rt ) is the gap indicator function b o (b (i)

(equal to 1 if rbt

b o and null otherwise). is observed at the top of the canopy from direction n

If all the rays are considered together, the average radiance is given by (i)

Nr X T 1 X (i) (i) (i) b o) = b o ). L(λ0 , n L (λ0 , rb1 , . . . , rbt , n Nr i=1 t=1 t

20

21

In the following and without lack of generality, T = maxi∈{1,...,Nr } {T (i) } transitions will be assumed for (i)

22

all the tracings. To be consistent, It

23

written as

24

(12)

(i)

and thus Lt

are assumed null for t > T (i) . Eq. (12) can be therefore

T Nr 1 XX (i) (i) (i) b o) = b o ). L(λ0 , n L (λ0 , rb1 , . . . , rbt , n Nr t=1 i=1 t

(13)

3

PROPOSED APPROACH

13

1

The latter expression will be useful in the optimization of multi-wavelength simulation.

2

3.2.3. Multi-wavelength simulation optimization

3

MC ray tracing at a given wavelength, λ0 , is time consuming as for each transition t of a ray i, the

4

algorithm has to find the next interception. Reaching transition t corresponds to the creation of all the ray path:

5

6

7

(i)

(i)

(i)

(i)

(i)

rb1 , . . . , rbt . Given this path, it is possible to compute the corresponding exitance, It (λ0 , rb1 , . . . , rbt ), using (i)

(i)

(i)

b o ). Eqs. (10)(B.2) iteratively. Applying afterwards Eq. (11) allows to derive the radiance Lt (λ0 , rb1 , . . . , rbt , n Exitance and radiance computations are too fast since they are composed of few mathematical operations at

8

each transition. Noteworthy that if the path is known, it is possible to compute exitance and radiance of each

9

wavelength λ, substituting Z(λ, . . .) [resp. ZL (λ, . . .)] for Z(λ0 , . . .) [resp. ZL (λ0 , . . .)] in Eq. (10) [resp. Eq.

10

(11)]. Lt (λ, . . .) derivation will be used in the following to simulate the radiance distribution at wavelength

11

λ given ray tracing at λ0 .

(i)

12

If MC ray tracing is done to simulate hyperspectral data separately, the global running time will increase

13

linearly as a function of the wavelength number. To limit this time, we propose in this study to use the same

14

tracing (path find) for all the wavelengths together, then compute exitance and radiance separately.

15

If two different wavelengths (λ0 , λ) are considered, it is clear that at each interception, they have gen-

16

erally different bidirectional scattering probabilities (called πλ0 and πλ ), different exitances and different ra-

17

diances in sensor direction. In the particular case, πλ0 = πλ , the MC path sampling done at wavelength

18

19

b o ), one has just to compute the set λ0 remains valid for λ. In this case, to compute the radiance L(λ, n (i)

(i)

(i)

(i)

(i)

b o )}i=1,...,Nr ;t=1,...,T relative to the different paths {b {Lt (λ, rb1 , . . . , rbt , n r1 , . . . , rbt }i=1,...,Nr ;t=1,...,T already

20

computed for λ0 [cf. Eq. (13)]. Now, as generally πλ0 6= πλ , the path sampling at λ0 if used to approximate

21

ray propagation at λ, produces bias in the global canopy reflectance distribution. In this work, we will try to

22

correct such an artifact. This way, it will be possible to simulate all the wavelengths together.

23

24

25

26

27

28

We assume in the following that ray tracing at canopy scale is done at wavelength λ0 , and we try to estimate the reflectance at wavelength λ 6= λ0 . As in Eq. (13), if all the rays are considered together, the average radiance at wavelength λ0 and after t transitions, Lt , is given by b o) = Lt (λ0 , n

Nr 1 X (i) (i) (i) b o ). L (λ0 , rb1 , . . . , rbt , n Nr i=1 t

(14)

b o ) can Using a weighted MC sampling technique, it is shown in Appendix C that for wavelength λ, Lt (λ, n

3

1

PROPOSED APPROACH

b o ). In this case, Lt (λ, n b o ) is written as be derived using the same ray tracing as Lt (λ0 , n b o) = Lt (λ, n

2

3

6

7

8

9

10

11

Nr 1 X (i) (i) (i) (i) (i) b o )ωλ0 ,λ (b L (λ, rb1 , . . . , rbt , n r1 , . . . , rbt ), Nr i=1 t

(15)

where (i)

(i)

ωλ0 ,λ (b r1 , . . . , rbt ) =

4

5

14

(i) (i) t Y bp ) πλ (b np−1 , n (i)

p=2

(i)

bp ) πλ0 (b np−1 , n

,

(16)

b p is the direction of propagation at the pth transition and πλ0 (b b p ) and πλ (b b p ) are the corren np−1 , n np−1 , n

sponding transition probabilities respectively at wavelengths λ0 and λ, and computed as explained by Eq. (B.1). b ) using the same tracing ωλ0 ,λ can be estimated from the tracing, therefore it is possible to derive Lt (λ, n

b ), based on the weighted sampling presented in Eq. (15). The same processing can be done for all as Lt (λ0 , n

b ). This way, it transition numbers t = 1, . . . , T . Therefore, it is possible to estimate the global radiance L(λ, n is possible to simulate all the wavelengths together.

12

The weighted sampling increases the MC noise as the different samples do not have the same weight,

13

therefore the result accuracy decreases as a function of the transition number. As in vegetation canopy near-

14

infrared is the less absorbent domain, it needs the maximal number of transitions and therefore weighted

15

sampling in this domain is the less accurate. To prevent this problem, the ray tracing is performed in this

16

domain reducing therefore the MC noise.

17

Noteworthy that in [30], it was proposed a similar technique to simulate multiple wavelengths together.

18

They are based on BRDF, ρso , decomposition as a weighted sum of the multiple scattering defined as the

19

power of the leaf albedo, w. They write ρso =

20

to the vegetation 3-D structure. In [30], it was showed that doing ray tracing on a given wavelength allows

21

to retrieve {ai }i=1,2,... and thus to derive ρso for any wavelength. Compared to our derivation of Eq. (15),

22

they neglect the ωλ0 ,λ terms. Same assumptions are done in [31] that in MC ray tracing, the energy update

23

after each interception and scattering is done without corrections. It means that both works assume that the

24

reflectance and transmittance probabilities are independent of the wavelength. According to the PROSPECT

25

simulations [24], this assumption is acceptable in the bi-Lambertian case. However, in our case, it does not

26

hold and specifically the specular effect impact on the scattering probability distribution is the highest in the

27

absorbent band creating different scattering distribution from a given wavelength to another.

P

i

ai wi , where ai , i = 1, 2, . . . are geometrical coefficients linked

4

1

2

3

EXPERIMENTAL RESULTS

15

Our model will be called the Leaf Canopy Vector Radiative Transfer model (LCVRT).

4. Experimental results In this section, simulations at leaf level will be summarized followed by simulations at canopy level. Then,

4

the model is validated using a web-based tool.

5

4.1. Leaf Level

6

7

Although results at leaf level are already presented in [23], here some simulations will be summarized allowing the interpretation of the results at canopy level.

8

Figure 2 shows the variation of the leaf spectrum with sun zenith angle varying from nadir to 75◦ . CChl

9

and CCar stand for chlorophyll and carotenoid contents, respectively. Absorption is high in the visible domain

10

due to the chlorophyll absorption in addition to carotenoid before 500nm. Peaks of absorption are observed in

11

shortwave infrared due to water. The comparison between the different curves show that reflectance increases

12

as a function of the sun zenith angle due to the increase of the specular effect. Conversely, the transmittance

13

decreases. For very low incidence (i.e. θs = 75◦ ), leaf becomes too bright as just small portion of the radiation

14

penetrates the leaf interior reducing therefore the pigment absorption effect.

15

16

Note that carotenoid content will be maintained 8µg/cm2 all along the experimental results and henceforth will be omitted. 0

1 0.8

0.2

0.6

0.4

0.4

0.6

0.2

0.8

0 400

800

1200

1600

2000

1 2400

Figure 2: Leaf spectrum for different incidence angles, θs , with values shown in legend. Solid and dashed lines correspond to the reflectance and transmittance, respectively. CChl = 80µg/cm2 and CCar = 8µg/cm2 .

17

In order to understand the effect of the specular effect, Figure 3 shows the variation of the different

18

components of the reflectance (V = 0 and it is omitted). When radiation comes from nadir, reflectance

4

EXPERIMENTAL RESULTS

16

1

approaches the Lambertian case with almost constant I and null values of Q and U . When θs increases the

2

reflectance intensity increases in the forward direction due the specular effect. This is confirmed by Q which

3

decreases to reach non-negligible negative value, corresponding to horizontal polarization. Indeed, specular

4

reflection creates light polarized parallel to the leaf surface. Moreover, Q does not present a sharp peak in the

5

exact forward direction because the epidermis surface producing the specular effect is not plane but ellipsoidal

6

shaped. This is also observed in terms of U that becomes non-null in the inclined-forward directions with

7

positive values in the left-forward and right-forward directions, respectively. This diagonal polarization is

8

created by the inclined surface of the leaf that creates inclined polarization.

9

10

Note that here, we omitted the bidirectional transmittance because it is almost Lambertian and unpolarized for all the incidence cases [23].

11

Figure 4 shows chlorophyll variation impact on the spectrum. As expected it increases absorption in the

12

visible domain. It is important to note here that absorption does not increases in the blue and red bands

13

(∼450nm & ∼650nm). Indeed, in this case, chlorophyll absorbs almost all the light reaching the interior of

14

the leaf for CChl = 50µg/cm2 , therefore the increase of the content has no effect. Nevertheless, light is not

15

completely absorbed because a part of it is reflected by the epidermis before reaching the mesophyll structure.

16

17

In order to understand the effect of the incident light polarization on the spectrum, Figure 5 shows its

18

variation for unpolarized, horizontal and vertical polarized light beams. For nadir illumination, almost all the

19

curves are superimposed. Indeed in all these cases, the electric field is parallel to the leaf surface and there

20

is no orthogonal component, therefore light is scattered in the same way [10]. Then, the curve differences

21

increase with the incidence. In particular, reflecance (resp. transmittance) increases in case of horizontal (resp.

22

vertical) polarization. Indeed, specular effect is the highest for the horizontal case increasing the reflectance

23

and decreasing therefore the transmittance [23].

24

4.2. Canopy level

25

At the beginning of this section, we present the general features and a validation of the MC weighted

26

sampling approach, then we study the different effects of leaf, vegetation, soil and light polarization on the

27

polarized reflectance distribution.

28

All along this section, vegetation layer is assumed of high 2m, composed of disk-shaped leaf of radius 5cm

29

homogeneously distributed across the layer with a given leaf inclination distribution functions (LIDF) varying

4

EXPERIMENTAL RESULTS

17

(a) I, θs = 0◦

(b) I, θs = 45◦

(c) I, θs = 60◦

(d) Q, θs = 0◦

(e) Q, θs = 45◦

(f) Q, θs = 60◦

(g) U , θs = 0◦

(h) U , θs = 45◦

(i) U , θs = 60◦

Figure 3: Leaf polarized BRDF for different incident angles, θs , shown in legend and unpolarized incident light. CChl = 80µg/cm2 . 0◦ and 180◦ correspond to the backward and forward directions, respectively.

4

EXPERIMENTAL RESULTS

18 0

1 0.8

0.2

0.6

0.4

0.4

0.6

0.2

0.8

0 400

800

1200

1600

2000

1 2400

Figure 4: Leaf spectrum for different incidence agles, θs , shown in legend and two values of chlorophyll contents: 80µg/cm2 and 50µg/cm2 plotted without and with markers, respectively

0

1 0.8

0.2

0.6

0.4

0.4

0.6

0.2

0.8

0 400

800

1200

1600

2000

1 2400

Figure 5: Leaf spectrum for different incidence angles, θs , plotted in legend and different polarizations: unpolarized, horizontal and vertical plotted without marker, with big circular and small triangular markers, respectively. CChl = 80µg/cm2 .

1

from a simulation to another. The leaf azimuthal distribution is assumed uniform whereas the zenithal one, gL ,

2

is assumed varying. As suggested in [32], it is parametrized using the coefficient set {ag , bg , cg , dg } as follows gL (θL ) =

3

4

2 [ag + bg cos(2cg θL )] + dg sin(θL ), π

where θL is the zenith angle.

5

In our study, the following particular distributions are considered:

6

• Planophile distribution: ag = 1, bg = 1, cg = 1, dg = 0,

7

• Erectophile distribution: ag = 1, bg = −1, cg = 1, dg = 0,

(17)

4

EXPERIMENTAL RESULTS

19

1

• Uniform distribution: ag = 0, bg = 0, cg = 0, dg = 1.

2

The leaf mock-up presented in Section 2.3.2 is used to simulate radiative transfer at leaf level. Chlorophyll

3

concentration is maintained 80µg × cm−2 all along this section, except in Section 4.2.6 where this content effect

4

is studied. It will be therefore omitted except in the pointed out section.

5

4.2.1. General feature

6

7

In this section, we try to study the canopy reflecance angular distribution as well as its variation as a function of the wavelength. Our model LCVRT will be used to couple 3D VRT at leaf and canopy scales.

8

Figure 6 shows an example of canopy polarized reflectance variation as a function of the sensor orientation

9

for two wavelengths in the red (650nm) and near-infrared (950nm) bands simulated using LCVRT. In both

10

wavelengths, the intensity that corresponds to the scalar reflectance shows firstly a peak in the exact backward

11

direction corresponding to the hot spot effect, then an increase in the forward direction that is due to the

12

specular effect. This is confirmed by the polarization Q that shows an important decrease in the forward

13

direction. Indeed, remember that single leaves produce specular reflected light polarized mainly perpendicular

14

to their surfaces as explained in the last section. Therefore, horizontal leaves as well as slightly inclined ones

15

are producing horizontal polarization (Q < 0) in the forward direction. The polarization U is null in the exact

16

forward direction, and positive and negative in the right-forward and left-forward directions, respectively.

17

Assuming that leaves are plane and produce specular effect as well as Lambertian diffuse scattering, this effect

18

is already simulated in [10]. It was due to the specular effect produced by inclined leaf reflect light in directions

19

of sensor in the inclined-forward directions. The polarization is parallel to the leaf surface therefore it is

20

not perfectly horizontal inducing non-null diagonal polarization. However, compared to our simulations, the

21

diagonal polarization in the latter study is less pronounced. Indeed, in addition to the leaf inclination that

22

produces diagonal polarization, in our case also single leaves produce it (cf. Sec. 4.1).

23

Recent polarization measurements on leaf and vegetation cover have been done in [9]. Distributions of I,

24

Q and U have been presented for the Murraya euchrestifolia Hayata plant assumed presenting planophile leaf

25

distribution. Comparing their measurements to our simulations shows the same order and the same trend in

26

terms of intensity (cf. Fig. 4 & 5 of [9]). In particular, the specular effect is well shown in measurements as in

27

simulations. In terms of polarization Q and U , same behaviors as simulations are also shown; i.e. an increase

28

of the absolute value of Q and U in the forward direction and the inclined-forward direction, respectively (cf.

4

EXPERIMENTAL RESULTS

20

1

Fig. 7, 9, 10 & 14 of [9]). It is noticeable that the measurements present specular effect and polarization

2

amplitudes different from our simulations. This is explained by the fact that these effects depend mainly on

3

the leaf smoothness; the more the leaf is smooth, the higher are specular effect and polarization.

4

The comparison between red and near-infrared bands shows an increase in intensity but not in polarization

5

in the latter case. In fact, chlorophyll absorption effect is responsible of the decrease of the reflectance in

6

the red band. Conversely, in near-infrared less absorption is shown, it is only due to the cell wall [23]. The

7

multiple scattering within leaf [23] and vegetation cover [33] are therefore more pronounced in the latter band

8

that increases the reflectance. At the same time, after multiple scattering rays are polarized in a random

9

direction and by averaging the polarization approaches zero. This explain why polarization does not increase

10

in near-infrared. Such stability of canopy polarization is also shown in actual measurement presented in Fig.

11

7, 9, 10, 13 & 14 of [9]. Moreover, measurements on leaves in Fig. 15 of the same work show that polarization

12

is almost wavelength-invariant.

13

It is also noted that relatively to the specular peak, the hot spot shows less intensity in the red band.

14

Indeed, in this domain the absorption caused by leaf interior is too high producing low scattering particularly

15

in the hot spot direction. However, the specular reflection is not affect as it is produced by leaf surface at the

16

top of the canopy. Conversely, in the near-infrared, the hot spot and the specular peaks are of the same order

17

since the specular effect conserves its amplitude whereas the scattering increases as the absorption is too low

18

in this band.

19

Figure 7 plots the canopy directional-hemispherical spectra in terms of reflectance and degree of polariza-

20

tion [cf. Eq. (3)]. As already claimed, absorption in the visible domain is due to chlorophyll. In the shortwave

21

infrared, the absorption is due to the water content. In the near infrared, less absorption is observed. Re-

22

garding the degree of polarization, one can see that it is inversely proportional to the reflectancce. Indeed,

23

the polarization is mainly due to the specular effect of the first order of scattering because diffuse scattering

24

provides light with random polarization that leads to unpolarized result by averaging. Thus, the higher the

25

effect of the high order scattering is, the lower the polarization is. Therefore, when the leaf scattering increases,

26

the canopy diffuse scattering increases and polarization degree decreases.

27

Besides measurements on actual vegetation cover presented in Fig 6 of [9] show the same trends as the

28

simulations. Nevertheless, the degree of polarization is the highest in measurement due to the smoothness of

29

Murraya euchrestifolia Hayata leaves.

4

EXPERIMENTAL RESULTS

21

(a) I, λ = 650nm

(b) Q, λ = 650nm

(c) U , λ = 650nm

(d) I, λ = 950nm

(e) Q, λ = 950nm

(f) U , λ = 950nm

Figure 6: Canopy polarized BRDF for unpolarized incident radiation, solar zenith angle θs = 60◦ , LAI = 3 and Planophile LIDF.

0.5

0.15

0.1

0.3

0.05 0.1 500

900

1300

1700

2100

0

Figure 7: Canopy directional-hemispherical spectrum of reflectance (solid line) and degree of polarization (dashed line) for unpolarized incident light, solar zenith angle θs = 60◦ , LAI = 3 and Planophile LIDF.

1

4.2.2. Multi-wavelength validation

2

Using weighted MC sampling in multi-wavelength simulations allows to use the same ray tracing to simulate

3

different wavelengths together. It has been already validated for ray tracing within leaf [23]. Remember that

4

the weighted MC algorithm accuracy decreases with the number of scattering, therefore the leaf case is more

5

complex than the canopy one since the high order scattering in leaves is more significant than in vegetation

4

1

EXPERIMENTAL RESULTS

22

cover.

2

Figure 8 shows the original spectrum simulated using the weighted MC with original ray tracing simulation

3

at λ0 = 1050nm as well as separate simulations from 550nm to 2350nm with step equal to 100nm. The com-

4

parison between the two kinds of simulations shows very close results confirming the validity of the developed

5

approach. It is important to note that reflectances are the closest in absorbent bands in visible and shortwave

6

domain unless they are far from the original wavelength and the scattering bidirectional distributions are dif-

7

ferent. It is due to the low impact of multiple scattering on reflectance that decreases accuracy. Conversely,

8

the degrees of polarization are the farthest in absorbent bands. Even though, the error remains proportional

9

to the degree of polarization which means that the relative accuracy remains almost the same.

0.5

0.15

0.1 0.3 0.05 0.1 500

900

1300

1700

2100

0

Figure 8: Canopy directional-hemispherical spectrum of reflectance and degree of polarization using weighted MC with λ0 = 1050nm (solid lines) and separate simulations (‘◦’ for reflectance and ‘’ for polarization) for solar zenith angle θs = 60◦ , LAI = 3 and planophile LIDF.

10

4.2.3. Light polarization effect

11

Figure 9 plots the different components of the reflectance variation as functions of the sensor orientation for

12

different incident light polarization and different leaf bidirectional scattering cases in the vector and scalar RT

13

cases, respectively. Note that the scalar RT concerns just the canopy level but the leaf level remains simulated

14

using vector RT. Nevertheless, according to [23], if the leaf incident light is unpolarized, the exiting intensity is

15

slightly dependent on the polarization effect within the leaf, thus it means that when scalar RT is assumed in

16

canopy it is equivalent to scalar RT in the coupled system (i.e. leaf and canopy). Moreover, in the scalar case,

17

three kinds of leaf scattering are considered; (i) Bidirectional scattering in incident and exiting directions as

18

in the vector RT case; (ii) Directional incidence and bi-Lambertian exiting from the two faces (reflectance and

19

transmittance); (iii) Lambertian incidence and bi-Lambertian exiting. Bi-Lambertian scattering is simulated by

4

EXPERIMENTAL RESULTS

23

1

integrating exitance over the two hemispheres and computing the equivalent Lambertian scattering: reflectance

2

if the exitance is in the same hemisphere as the incoming ray and conversely transmittance. Lambertian

3

incidence assumes that the incident ray is coming from nadir, i.e. normal to the leaf surface.

4

In the VRT case, the intensity comparison reveals that the horizontal polarization presents the highest

5

specular effect whereas the vertical case is the lowest. Indeed, as explained at leaf level (cf. Sec . 4.1), the

6

horizontal component of the specular reflection is higher than the vertical one [10, 23]. Unpolarized results are

7

intermediate, moreover they are close to the scalar RT with bidirectional leaf scattering. In fact, for both cases,

8

it was proved in [23] that the first order reflectance is the average between horizontal and vertical cases. For

9

high order of scattering, these two cases behave differently, however as the observation is done in the red band

10

(650nm), high order scattering is negligible. Besides, the comparison between the different leaf scattering cases

11

in scalar RT shows that the specular effect is removed in both directional-bi-Lambertian and Lambertian-bi-

12

Lambertian cases. Indeed, this effect is due the directional effect at the leaf level, thus when the exiting light

13

is bi-Lambertian, this effect is particularly removed. Moreover, the hot spot peak appears the largest in the

14

bidirectional scattering case. Indeed, leaves are not actually Lambertian and scattering at large zenith angle is

15

high since epidermis cells are ellipsoidal shaped [23]. Therefore, the large peak in the bidirectional case is not

16

only due to hot spot but also caused by the non-Lambertian effect. Comparing the Lambertian-bi-Lambertian

17

to the directional-bi-Lambertian case shows that the former is the lowest. Indeed, in this case, the light is

18

assumed coming from nadir that corresponds to the less diffusing case as shown in Figure 2.

19

In terms of polarization, Q is more important when the incident light is polarized, since after first order of

20

scattering it remains polarized. Moreover, the scattered light conserves its incident polarization. In particular,

21

Q is negative and positive in the forward direction for the horizontal and vertical cases, respectively. It is the

22

highest in the former case since the horizontal coefficient of scattering is larger than the vertical one. For the

23

same reason, U is the highest in the horizontal case. For horizontal polarization, U is positive (resp. negative)

24

in the right-forward (resp. left-forward) case since horizontal (resp. vertical) polarization turns in the direction

25

of the first (resp. second) diagonal. This is inverted in case of vertical polarization. Indeed, when scattering is

26

inclined in the right-forward direction, the vertical polarization turns in the direction of the second diagonal.

27

28

29

It is noted that in the scalar RT case, the polarization is not taken into account for that Q and U are null (results are omitted). In order to see the variation of the polarization effect as a function of the wavelength, Figure 10 draws the

4

EXPERIMENTAL RESULTS

24

(a) I, Unpolarized

(b) Q, Unpolarized

(c) U , Unpolarized

(d) I, Scalar RT, Bidirec

(e) I, Scalar RT, Dir-Bi-Lamb

(f) I, Scalar RT, Lamb-Bi-Lamb

(g) I, Vertical

(h) Q, Vertical

(i) U , Vertical

(j) I, Horizontal

(k) Q, Horizontal

(l) U , Horizontal

Figure 9: Variation of the canopy polarized BRDF as functions of the incident light polarization (unpolarized, vertical and horizontal) and RT modeling (vector for rows 1, 3 and 4 or scalar for row 3). In the scalar case, three cases of leaf directional scattering are considered: Bidirectional (Bidir), Directional-bi-Lambertian (Dir-Bi-Lamb) and Lambertian-bi-Lambertian (LambBi-Lamb). θs = 60◦ , LAI = 3 and λ = 650nm.

4

EXPERIMENTAL RESULTS

25

1

spectra of reflectance and degree of polarization. By comparison between the cases of polarized incident ray,

2

one can see that the reflectance intensity is the highest in the horizontal polarization case because the coefficient

3

of horizontal reflection is the highest. Unpolarized and Scalar RT with bidirectional leaf case produce close

4

results. Indeed, as we already explained, in terms of first order of reflectance, they give theoretically the same

5

results. Then, for the high orders, they behave differently as the scalar RT does not take into account the

6

polarization. Nevertheless, for diffuse light, in general there is no favored direction of polarization and it is

7

therefore almost unpolarized and therefore vector RT results approach the scalar ones.

8

Scalar RT with Directional-bi-Lambertian leaves produces close results to the scalar RT case as the difference

9

is only in the probability distribution of the exiting radiation but with the same energy. Nevertheless, in the

10

visible domain, the original scalar RT case is the lowest. Indeed, first order of scattering is mainly oriented in the

11

forward direction in this case, that increases the probability of interception that is followed by high absorption

12

in this domain. In opposite, scalar RT with Lambertian-bi-Lambertian leaves assumes that radiation is coming

13

from nadir that corresponds to the less scattering case and therefore leads to the lowest spectrum of reflectance.

14

In terms of degree of polarization, the polarized incident cases are the highest and particularly the horizontal

15

one. Indeed polarization is mainly produced by the first order of scattering, and if the incident ray is polarized

16

the scattered one remains. For the scalar RT, the polarization is null.

0.3 0.5 0.2 0.3 0.1 0.1 400

800

1200

1600

2000

0 2400

Figure 10: Canopy directional-hemispherical spectrum of reflectance and degree of polarization for different incident light polarization and taking or not the polarization effect into account. LAI = 3, θs = 60◦ and Planophile LIDF. UNP, HOR and VER correspond to unpolarized, horizontal and vertical incident polarized lights, respectively. SRT DIR-L and L-L are the simulation using scalar RT considering bidirectional, directional-bi-Lambertian, Lambertian-bi-Lambertian leaf scattering, respectively. θs = 60◦ and LAI = 3.

4

17

EXPERIMENTAL RESULTS

26

4.2.4. Canopy structural parameters

1

4.2.4.1. LAI variation

2

Figure 11 draws the variation of the different polarized BRDF components as functions of LAI. The polarization

3

U is omitted as it follows the same behavior as Q. In terms of intensity, one can see that it increases as a function

4

of LAI since the incident ray interception and therefore scattering increase when foliage density increases. Plots

5

corresponding to LAI≥ 3 are too close in the forward direction as reflectance reaches its asymptotic value as

6

a function of LAI too fast. Indeed, it depends mainly on the first order specular reflection that is independent

7

on deep leaves and thus does not increase with LAI. For the same plots, a non-negligible increase is observed

8

in the backward direction and the hot spot appears clearer than in case LAI=2. This is explained by the fact

9

the in the near-infrared domain, the reflectance depends on the multiple scattering and therefore it continues

10

to increase with deep leaves and high LAI [33]. Conversely, the polarization Q is almost constant functions

11

of LAI. Indeed, as multiple scattering produces almost unpolarized light, polarization increases only by the

12

increases of the first order of scattering, for that the polarization saturation is reached earlier than intensity.

(a) I, LAI = 2

(b) I, LAI = 3

(c) I, LAI = 5

(d) Q, LAI = 2

(e) Q, LAI = 3

(f) Q, LAI = 5

Figure 11: Variation of the canopy polarized BRDF as functions of LAI for solar zenith angle θs = 60◦ , planophile LIDF and λ = 950nm.

13

Figure 12 presents the canopy reflectance and degree of polarization, P , spectra for different values of LAI.

4

EXPERIMENTAL RESULTS

27

14

As already explained, the reflectance increases as a function of LAI. The increase is inversely proportional to

1

absorption. In fact, the higher the scattering is, the more important the effect of the diffuse scattering is.

2

Similarly, P variation is shown in non-absorbent wavelengths. Particularly, it shows a decrease as a function

3

of LAI due to the increase of diffuse scattering that corresponds to rays traveling along random paths and

4

therefore reaching the observation with random polarization leading in sum to unpolarized light.

0.15

0.5

0.1 0.3 0.05 0.1 400

800

1200

1600

2000

0 2400

Figure 12: Canopy directional-hemispherical spectrum of reflectance and degree of polarization of different LAI values and for solar zenith angle θs = 60◦ and Planophile LIDF.

5

4.2.4.2. θs variation

6

Figure 13 plots the variation of the polarimetric reflectance components I and Q as functions of the observation

7

angle, for different solar incidences. Note that U is omitted as it shows the same trends as Q. The comparison

8

between the solar angles affirms that the intensity is an increasing function of the solar zenith angle. Indeed,

9

when sun is highly inclined, the optical path increases, the probability of interception increases accordingly and

10

therefore the reflectance. Particularly, the hot spot effect appears the brightest for θs = 60◦ . However, when

11

θs increases, it becomes less and less distinguishable compared to the other observation directions. Indeed, the

12

high inclination reflectance approaches its upper limit and therefore becomes close to the hot spot amplitude.

13

In terms of polarization, the simulations reveal that its effect increases as a function of the solar incidence

14

angle and the specular peak appears clearly in the forward direction only for θs = 60◦ . This variation is

15

explained by the fact that the specular effect increases with the sun-sensor angular distance. In particular,

16

Q polarization is well noticeable respectively in the forward and inclined forward directions for the highest

17

θs value (θs = 60◦ ). Then, when θs decreases to 30◦ , the polarization in the forward direction decreases but

18

remains distinguishable. This effect almost vanishes when the sun is at the nadir because the sun-sensor angle

4

19

EXPERIMENTAL RESULTS

28

becomes less than 90◦ reducing therefore the specular effect.

(a) I, θs = 60◦

(b) I, θs = 30◦

(c) I, θs = 0◦

(d) Q, θs = 60◦

(e) Q, θs = 30◦

(f) Q, θs = 0◦

Figure 13: Variation of the canopy polarized BRDF as functions of the solar zenith angle for LAI = 3, planophile LIDF and λ = 650nm.

1

Figure 14 compares the reflectance and polarization spectra of the different solar zenith angles. The

2

reflectance increases as a function of the incidence since the optical path increases. Polarization increases, too,

3

since it is mainly caused by the specular effect in the forward direction and the latter is an increasing function

4

of θs .

5

The study of the polarization as a function of the sun zenith angle has been also presented in [9]. Same

6

trends as our simulations are observed. In particular, the increase of polarization in the forward direction

7

as a function of the zenith angle are plotted in Fig. 9 & 10 of [9]. Moreover, the decrease of the degree of

8

polarization spectrum as a function of the sun zenith angle is also confirmed by measurements (cf. Fig. 6 of

9

[9]).

10

4.2.4.3. LIDF variation

11

As for the variation of θs , Figure 15 plots the different reflectance components for different leaf distributions.

12

In terms of intensity, and far from the forward direction, the planophile distribution presents the lowest

13

variation since it approaches horizontal distribution and therefore the Lambertian case. The intensity is also

4

EXPERIMENTAL RESULTS

29

0.5

0.15

0.1

0.3

0.05 0.1 400

800

1200

1600

2000

0 2400

Figure 14: Canopy directional-hemispherical spectrum of reflectance and degree of polarization of different solar incident angle. LAI = 3 and Planophile LIDF.

14

the highest (resp. the lowest) in the planophile (resp. the erectophile) case mainly close to nadir. Indeed,

1

erectophile distribution approaches vertical distribution therefore scattering in high elevated sensor orientation

2

is the lowest. The specular effect is limited to high sensor inclination (≥ 80◦ ) for uniform and erectophile

3

leaf distributions because in these cases, leaves are highly inclined and the probability to find leaves allowing

4

specular effect in the forward direction is low. In particular, in the exact forward direction, specular effect is not

5

observed in the erectophile case. This is confirmed in terms of polarization, where the forward effect is almost

6

null in this particular direction. The specular effect increase in the inclined-forward direction is explained by

7

the fact that leaves allowing this effect are more inclined than those required in exact forward. and in the

8

erectophile case the probability of inclined leaves is much higher than of horizontal ones.

9

10

Noteworthy that the decrease of intensity and polarization in the forward direction in case of erectophile distribution has been already observed in actual measurements shown in Fig. 4, 5, 13 & 14 of [9].

11

Figure 16 compares the spectra relative to the different LIDF. It is shown that the reflectance decreases with

12

the leaf distribution average inclination. Indeed, the more horizontal the leaf are the more light they diffuse in

13

high elevated orientations (corresponding to the lowest optical depth and so the highest probability to exit to

14

canopy), and therefore the highest the reflectance is. The spectra of polarization confirms that it dramatically

15

decreases in case of uniform and erectophile leaf distributions. Besides, the decrease in the spectrum from

16

planophile to erectophile distribution is already shown in the actual measurements in Fig. 6 of [9].

4

EXPERIMENTAL RESULTS

30

(a) I, LIDF=Planophile

(b) I, LIDF=Uniform

(c) I, LIDF=Erectophile

(d) Q, LIDF=Planophile

(e) Q, LIDF=Uniform

(f) Q, LIDF=Erectophile

Figure 15: Variation of the canopy polarized BRDF as functions of LIDF for LAI = 3, planophile LIDF and λ = 650nm.

17

4.2.5. Soil effect

1

All along the last sections, the soil was assumed black in order to focus on the effect of the different

2

vegetation properties on the polarized reflectance. Here, we complement our study by adding a soil. To do it,

3

a dry (bright) soil is added to our vegetation layer mock-up. It is assumed depolarizing.

4

Figure 17 plots the variation of I and Q as functions of the observation angle for two values of LAI. U

5

is also omitted in this figure as it behaves as Q. The same results are drawn in Figure 11 but for black soil.

6

When soil is added, the intensity increases mainly for observation angles close to nadir as the optical path

7

is the lowest in this case and it is likely for incident radiation reaching the soil to exit in these directions.

8

In terms of polarization, soil does not affect the results. Indeed, the considered soil is assumed depolarizing,

9

therefore radiation being interacted with it will be unpolarized and then if interacted with leaves will have

10

random polarization. In both cases (intercepted by leaves or not), radiation scattered by the soil will produce

11

an average null polarization. This explains why polarization is independent from soil reflectance.

12

Figure 18 compares the spectra of reflectance and degree of polarization of black and bright soils for two

13

values of LAI. In terms of reflectance, it is the highest for bright soil as the canopy reflectance increases with

14

soil brightness. It is noted that canopy reflectance increase is the highest for the lowest LAI value (i.e. LAI=2)

4

EXPERIMENTAL RESULTS

31

0.5

0.15

0.1

0.3

0.05 0.1 400

800

1200

1600

2000

0 2400

Figure 16: Canopy directional-hemispherical spectrum of reflectance and degree of polarization of different LIDF. LAI = 3 and θs = 60◦ .

(a) I, LAI=2

(b) Q, LAI=2

(c) I, LAI=3

(d) Q, LAI=3

Figure 17: Variation of the canopy polarized BRDF as functions of LAI, for solar zenith angle θs = 60◦ , planophile LIDF, soil reflectance ρs = 0.3 and wavelength λ = 650nm.

15

because the optical depth is the lowest in this case and it is more likely for incident radiation to reach the

1

soil and then exit the canopy. In case of bright soil, reflectance is the highest for the lowest LAI in absorbent

2

domains (visible and shortwave infrared) and for high LAI otherwise (near infrared). This is can be explained

3

by the fact that when soil reflectance is higher than the leaf albedo, it is better to have less foliage to increase

4

soil impact on canopy reflectance. Conversely, when albedo is too high, it is better to increase the foliage

5

density.

6

Degree of polarization decreases with soil brightness since the intensity increases with soil brightness,

7

whereas polarizations Q and U remain constant for such a depolarizing soil. Moreover, the degree of polarization

8

decrease is the highest for the lowest LAI value because the soil contribution to reflectance is the highest in

9

this case.

4

EXPERIMENTAL RESULTS

32

0.15

0.5

0.1 0.3 0.05 0.1 400

800

1200

1600

2000

0 2400

Figure 18: Canopy directional-hemispherical spectrum of reflectance and degree of polarization of different LAI for black and dry soil. For dry soil, the latter ‘s’ is added in legend. ρs is the bare soil reflectance. Planophile LIDF and θs = 60◦ .

(a) I, CChl = 80µg/cm2

(b) Q, CChl = 80µg/cm2

(c) I, CChl = 50µg/cm2

(d) Q, CChl = 50µg/cm2

Figure 19: Canopy polarized BRDF for unpolarized incident radiation, solar zenith angle θs = 60◦ , LAI = 3 and Planophile LIDF, λ = 550nm.

10

4.2.6. Leaf chlorophyll effect

1

In leaf scale, we will study the chlorophyll variation effect on the polarized reflectance. Figure 19 compares

2

two chlorophyll cases 50µg/cm2 and 80µg/cm2 in terms of reflectance distribution in the green domain (550nm).

3

This wavelength is prefered to 660nm because chlorophyll variation has no effect in the latter case as drawn

4

in Figure 4. It is shown that reflectance decreases as a function of the chlorophyll content due to its high

5

absorption. Polarization Q remains almost the same as they depend on the leaf surface but not the interior.

6

4.3. Validation

7

Validation will be besed on the results of the third phase of the RAdiation transfer Model Intercomparison

8

(RAMI) [32]. RAMI proposed to benchmark the existing radiative transfer model and produced a reference.

9

The latter is derived based on the surrogate truth 3D Monte Carlo models such as DART [14], FLIGHT [12]

10

and Rayspread [34] that display generally close reflectance curves. RAMI benchmark has been made useful

11

to validate model through a web-based tool called On-line Model Checker (ROMC) exercise [35]. Current

4

12

1

EXPERIMENTAL RESULTS

33

ROMC version does not take into account the leaf scattering bidirectional distribution nor the polarization effect. Therefore, it will allow to validate only the canopy level of our model.

2

Three canopy structure are proposed in ROMC: homogeneous (HOM), heterogeneous (HET) and a combi-

3

nation between them (HETHOM). For each structure, two vegetation kinds are assumed: turbid (TUR) and

4

discrete (DIS). The leaf radius is assumed null for the turbid case and equal to 0.05m or 0.1m for the discrete

5

one. Three leaf optical properties are suggested corresponding to two wavelengths (red and near infrared) and

6

a purist corner cases. They are called RED, NIR and NR1, respectively. The corresponding leaf reflectance

7

and transmittance (ρ, τ ) are assumed equal (0.0546,0.0149), (0.4957,0.4409) and (0.5,0.5), respectively. The

8

soil optical properties are varieted in accordance with leaf. Its reflectance is assumed equal to 0.127, 0.159 and

9

1, respectively. Finally three LIDF are considered: planophile (PLA), uniform (UNI) and erectophile (ERE).

10

In this work, as our model allows to simulate 3-D heterogeneous scene both HOM and HET scenes are

11

considered. Moreover, since we are interesting to realistic canopy, only discrete leaf cases (DIS) are presented.

12

Table 2 summarizes the five tested homogeneous simulation experiments and the heterogeneous one, HET01.

13

In the latter casen, the scene is composed of fifteen spheres of radius 10m covering the soil background. The

14

mock-ups corresponding to HOM and HET scenes are drawn in [32].

15

ROMC proposes eleven measurements. One can cite the total BRDF in principal (brfpp) and cross (brfop)

16

planes as well as the single collided (brfpp_co_sgl/brfop_co_sgl: BRDF in principarl and cross planes

17

due to only one scattering by vegetation), uncollided (brfpp_uc_sgl/brfop_uc_sgl: BRDF due to only one

18

scattering by the soil), multiple-scattered (brfpp_mlt/brfop_mlt: BRDF due to multiple scattering in the

19

scene: vegetation and soil) and finally the spectral albedo of the canopy, i.e., the Directional-Hemispherical

20

reflectance and foliage absorption.

21

In our validation, we are interested in the eight BRDF terms: global, co_sgl, uc_sgl and mlt in both

22

principal (brfpp) and cross (brfop) planes. Table 3 presents the performances of our model for the considered

23

scenes HOM and HET and the different measurements in terms of percentage of agreement with the RAMI

24

reference. The results are avereaged by scene. In sum, our model, LCVRT, presents satisfactory results. The

25

co_sgl results are generally the best because the first order reflectance is produced by the first ray/medium

26

interseption that is almost modeled in the same way in all the models. In brfpp_co_sgl, our model shows

27

better results than in brfop_co_sgl because in the first case appears the hot spot that is well simulated by

28

LCVRT. Conversely, the worst results correspond generally to co_sgl. Indeed, for dense vegetation, it is

4

EXPERIMENTAL RESULTS

34

Table 2: List of the ROMC tested experiences.

Experiment

θs

rl /h

LAI

LIDF

Wavelength

HOM11

0◦ , 30◦ , 60◦

0.05/1

1

ERE, PLA

NR1

HOM12

0◦ , 30◦ , 60◦

0.05/1

2

ERE, PLA

NR1

HOM13

20◦ , 50◦

0.1/2

3

PLA, UNI

NIR,RED

HOM03

20◦ , 50◦

0.05/2

3

ERE

NIR,RED

HOM15

0◦ , 30◦ , 60◦

0.05/1

5

ERE,PLA

NR1

HET01

20◦ , 50◦

0.1/10

5

UNI

NR1,NIR,RED

Table 3: Percentage of agreement of LCVRT with the ROMC reference. Measurement

brfpp

brfpp_co_sgl

brfpp_uc_sgl

brfpp_mlt

brfop

brfop_co_sgl

brfop_uc_sgl

brfop_mlt

HOM11_DIS

98.109

98.275

64.955

97.696

63.318

95.280

64.015

97.536

HOM12_DIS

98.081

98.885

61.972

96.459

75.409

86.442

59.347

93.439

HOM13_DIS

99.328

99.418

71.652

92.834

70.844

57.611

75.444

68.770

HOM03_DIS

99.218

99.635

86.864

97.867

87.244

88.200

92.761

98.263

HOM15_DIS

98.510

98.584

64.131

62.884

77.736

89.858

63.704

51.629

HET01_DIS

99.469

99.256

95.740

99.868

96.346

92.696

95.668

99.794

29

unlikely for radiation to reach the bottom of the medium and come back to the top, therefore a very large

1

number of rays is needed to accurately simulate the soil contribution. Nevertheless, its contribution to the

2

total BRDF is too small. Indeed, brfpp and brfop results remain too high even when it is not the case for

3

co_sgl. mlt performances are lower than co_sgl because the former is created by multiple scattering that is

4

implemented in different ways in the different models, in addition the MC noise increases with the number of

5

interceptions reducing the performances.

6

To analyze in depth the ROMC results, we present in Figures 20 & 21 the different reflectance bidirectional

7

distributions in the principal plane related to a homogeneous, HOM15_DIS_ERE_NR1_60, and a heterogeneous,

8

HET01_DIS_UNI_NR1_50, scenes, respectively. Note that results for the cross plane are omitted as they show the

9

same trends as the principal one. In general, our curves are within the gray-colored envelops around the ROMC

10

reference, corresponding to error lower than 5%. The sole exception is for uc_sgl_HOM15_DIS_ERE_NR1_60.

11

Indeed, in this case LAI is too high (i.e. LAI=5), the probability that incident rays be intercepted by the

4

EXPERIMENTAL RESULTS

(a) Global, (0.0188)

(b) co sgl+, (0.0025)

35

(c) uc sgl, (0.0006)

(d) mlt, (0.0069)

Figure 20: Bidirectaional reflectance in principal plane for ROMC experiment HOM15, in the purist corner case with observation zenith angle θs = 60◦ and erectophile leaf angle distribution. Bracket value is the RMSE.

(a) Global, (0.0147)

(b) co sgl, (0.0050)

(c) uc sgl, (0.0046)

(d) mlt, (0.0180)

Figure 21: Bidirectaional reflectance in principal plane for ROMC experiment HET01, in the purist corner case with observation zenith angle θs = 50◦ and uniform leaf angle distribution. Bracket value is the RMSE.

5

DISCUSSION

36

12

soil background and reach again the canopy top is too low, this fact increases the MC noise and decreases the

1

performances. Nevertheless, the reflectance in this case is lower than 0.2% that means that this result does

2

not impact the global canopy reflectance. Note that in HOM15_DIS_ERE_NR1_60, LAI is also equal to 5, even so

3

our model performs better than the latter case. This is explained by the fact that in heterogeneous case, the

4

soil is well sunny and observed, therefore the number of incident rays intercepted by it and then reaching the

5

sensor is too high, reducing therefore the MC variance. This is confirmed by the reflectance amplitude that

6

exceeds 35%. Finally, the hot spot peak is well preserved using our model.

7

5. Discussion

8

Classical radiative transfer models on canopy scale assume that leaves are bi-Lambertian and neglect the

9

scattering bidirectional distribution and polarization effects. Nevertheless, recent research tried taking into

10

account them assuming that a part of the leaf scattering is due to the specular reflection that is polarized.

11

Comparing these model simulations [10] to actual measurements such as those presented in [9] shows significant

12

differences in leaf and then vegetation cover scales. Indeed, actual leaves are not Lambertian nor presenting

13

a sharp peak of specular reflection in the forward direction. The leaf reflectance distribution is rather smooth

14

with a continuous increase in the forward direction due to the specular effect produced by its rough surface.

15

The corresponding polarization distribution shows a wide spot in the forward direction. In addition to hori-

16

zontal polarization, diagonal polarization is shown in the inclined-forward directions. It was possible to mimic

17

such a scattering behavior using MC ray tracing on leaf level [23]. Then, at canopy level, the polarized re-

18

flectance results highly depend on the leaf scattering and polarization bidirectional distributions. Therefore, it

19

is important to simulate RT at canopy scale tacking into account the non-Lambertian scattering properties of

20

leaves.

21

The developed model in this paper couples 3-D vector RT models at leaf and canopy scales. The results are

22

well correlated to the measurements shown in [9] in terms of bidirectional distributions of I, Q and U as well

23

as spectrum of degree of polarization. This confirms the validity of our model and its capability to simulate

24

polarization at canopy scale taking into account the radiative regime at leaf level. Nevertheless, simulated

25

and measured reflectance amplitudes are not too close at canopy scales (nor at leaf level as shown in our last

26

study [23]). This is mainly due to the considered leaf mock-up that does mimic the actual considered leaf. In

27

particular, to mimic Q and U , firstly epidermis surface geometry should be well known as it is the main origin of

5

DISCUSSION

37

28

polarization angular distribution. For instance, if leaf oblateness is increased, leaf specular effect approaches the

1

ideal mirror case as presented in our study [10] and sharp polarization in the forward direction will be observed.

2

Conversely, if oblateness decreases, polarization effect decreases and becomes flat around the forward direction

3

[23]. Secondly, the polarization depends on the epidermis wall refraction index, the higher it is, the higher the

4

specular effect is. In sum, to mimic actual leaf polarization, epidermis oblateness and cell wall refraction index

5

should well known. For instance, they can be inverted using actual leaf polmarization measurements. Leaf

6

scattering intensity (reflectance and transmittance) depends on the specular effect that is already explained as

7

well as the radiation penetrating the leaf interior. In turns, the latter depends on the epidermis transmittance

8

and therefore its geometry and refractive index in addition to the mesophyll structure, e.g. palisade and spongy

9

cell spatial organization, chloroplast distribution within each cell, etc. All these properties are varying from a

10

given plant spacies to another, this means that to be accurate different mock-ups should be created to mimic

11

the different kinds of leaves. This is an actual need of the current and future remote sensing community as

12

amoung the hot research topic is to estimate internal leaf mechanisms such as photosunthetic activity based on

13

xanthophyll cycle and the sun induced fluorescence (SIF) using hyperspectral and fluorescence sensors [26, 36]

14

that measure scattering spectrum with high precision.

15

It is noted that polarization Q (resp. U ) can reach around 0.1 (resp. 0.05) in the forward (resp. inclined-

16

forward) direction. Although they are not too high, they are measurable by existing polarimetric sensor. In

17

this direction, the developed approach can be helpful to model and interpret such sensor measurements. This

18

permits better estimation of the vegetation properties.

19

Sensitivity study reveals that degree of polarization increases with leaf absorption, decreases with average

20

leaf inclination (i.e. highest values in planophile case and lowest in erectophile) and remains almost constant as

21

a function of LAI. Therefore, polarization can be used to make distinction between different cultures as they can

22

be different in terms of leaf absorption and angular distribution. It is also possible to add the leaf polarization

23

bidirectional scattering distribution as characteristic since it is varying from a spacies to another. Such a

24

classification is robust with respect to plant growth since polarization is slightly sensitive to LAI. Another

25

application of the polarization can be the derivation of the average leaf inclination distribution parameter.

26

Indeed, for a given plant species for which the leaf scattering properties are known, canopy polarization is

27

mainly function of the latter parameter.

28

As canopy polarized reflectance depends on the incident light polarization, remote sensing measurements

6

CONCLUSION

38

29

can not be interpreted regardless the atmosphere polarization effect [37]. In particular, intensity increases with

1

horizontal polarization and conversely decreases with the vertical one, it is important to couple VRT at canopy

2

and atmosphere scales to take into account the atmosphere light polarization.

3

6. Conclusion

4

This work depicts canopy polarized reflectance simulation based on leaf and canopy 3-D vector radiative

5

transfer model coupling. Leaf 3-D mock-up that is presented in my recent work [23] has been used. Acceleration

6

method on leaf radiative transfer code allowing to simulate multi-wavelengths and multi-angular scattering

7

distribution has been adopted in this scale and then extended to the canopy one. It consists of using the

8

same Monte Carlo ray tracing to simulate different wavelengths, unless they do not necessarily have the same

9

directional scattering probabilities. To do it, the MC averaging is weighted by corrective coefficients to transfer

10

knowledge on bidirectional reflectance distribution from a given wavelength to another.

11

Polarization of an incident light on a given leaf within the canopy is not known a priori and it changes

12

continuously after scattering from a leaf to the next. Thus to deal with polarization, leaf scattering function

13

is simulated over a canonical basis (unpolarized, horizontal, diagonal and circular polarizations) allowing the

14

construction of the Stokes phase function.

15

16

Simulation results validate firstly the weighted MC sampling. Variation of the polarized reflectance as a function of scene geometry, leaf and canopy properties and incident light polarization are studied.

17

Scene geometry as well as the sensor and sun angular variations are studied. As in recent literature, high

18

horizontal polarization is observed in the forward direction and diagonal polarization is observed in the inclined

19

forward-direction, they are caused by the specular effect. Moreover, it was shown that polarization impact

20

is the highest for absorbent medium since diffuse light (that is generally unpolarized) effect is negligible with

21

respect to the first order scattering which is polarized. Besides, intensity and polarization are highly varying

22

functions of the sun incidence, respectively. The intensity is the highest for high sun zenith angle since the

23

optical path lengthens. Likewise, the polarization is an increasing function since the probability of specular

24

reflection in the forward direction increases with the sun-sensor angular distance.

25

In terms of vegetation cover and leaf properties three effects were studied: LAI, LIDF and leaf chlorophyll

26

content variations. When LAI increases, the reflectance intensity increases whereas the polarization remain

27

almost constant. In fact, LAI increases diffuse scattering that increases intensity with unpolarized light.

6

CONCLUSION

39

28

Concerning the leaf inclination distribution, when it varies the reflectance and intensity vary accordingly. In

1

particular, the more horizontal the leaves are, the higher the intensity and polarization are. The interception

2

increases in this case, increasing therefore the intensity. The probability of specular reflection in the forward

3

direction and therefore the degree of polarization are the highest for horizontal leaves. Besides leaf absorption

4

in the visible domain is mainly linked to the chlorophyll content. It is pointed out that reflectance in the

5

considered domain decreases as a function of the chlorophyll content, whereas Q and U remain the same as

6

they depend only on the leaf surface geometrical and optical properties. In sum, chlorophyll increases the

7

polarization degree.

8

Soil impact is also studied. Depolarizing soil has been added to our mock-up. It was shown that reflectance

9

increases with soil brightness, polarization components remains the same as radiation interacting with soil is

10

unpolarized. Consequently, the degree of polarization decreases with the soil brightness.

11

The way the radiative transfer is modeled was deeply studied. In particular, neglecting light polarization

12

and leaf bidirectional scattering were discussed. Comparison between Vector and scalar RT was especially

13

pointed out. It was shown that they provide contrasted results if incident light is polarized and mainly in

14

the near-infrared domain where the diffuse scattering is important. This proves the importance of both using

15

vector RT and taking into account the canopy incident light polarization. The letter is due to the atmospheric

16

effect. Besides neglecting the leaf bidirectional scattering leads to inaccurate simulations of the reflectance

17

proving therefore the importance of using the 3-D RT modeling within leaf.

18

Qualitative validation was done comparing our simulations to laboratory measurements presented in [9]. In

19

general, similar trends are shown for different configurations: sun and sensor angular variation, LIDF variation

20

(planophile and erectophile) and wavelength variation (red and near-infrared). However, our polarization

21

simulations are slightly lower than measurements since measured plant leaves are too bright in the specular

22

direction due to the cuticle film covering their surfaces. Such a film can be added to our leaf mock-up to

23

simulate this kind of leaves in future work.

24

Quantitative validation of our model was performed using the web-based tool, ROMC, that compares our

25

simulations to benchmark references. The results show a good agreement between our simulations and the

26

references. Even though ROMC does not take into account light polarization, it allows to validate the scalar

27

part of our model that contains scene 3-D representation, ray tracing, MC simulation and reflectance derivation.

APPENDIX A LEAF POLARIZATION MODELING

28

1

2

3

4

5

Appendix A. Leaf polarization modeling In order to derive VRT at canopy scale it is important to estimate the Stokes scattering function, F, at leaf level that allows to link the exiting radiation to the incident one. To do it, it will be shown that one has to estimate the Stokes vector corresponding to a canonical representation. As canonical decomposition, four polarizations are considered, 

IN

6

7

8

11















 1   1   1   1                   0   0   1   0          =   , IQ =   , IU =   , IV =   ,          0   1   0   0                  1 0 0 0

(A.1)

where IN , IQ , IU and IV correspond to unpolarized, perpendicular, diagonal and circular states, respectively. This way, a ray of Stokes vector I of components I, Q, U and V can be written as

I = (I − Q − U − V ) × IN + Q × IQ + U × IU + V × IV .

9

10

40

(A.2)

For given incidence and scattering directions, the exitance Stokes vector is computed using RT model at leaf level for each state of polarization. One finds 

Id,N

12

 Id,N    Qd,N  =   Ud,N   Vd,N





  Id,Q       Qd,Q    , Id,Q =      Ud,Q     Vd,Q





  Id,U       Qd,U    , Id,U =      Ud,U     Vd,U





  Id,V       Qd,V    , Id,V =      Ud,V     Vd,V



     ,    

(A.3)

13

where Id,N , Id,Q , Id,U and Id,V are the scattered exitances when the incident irradiances are IN , IQ , IU and

14

IV , respectively.

15

16

17

Due to the polarization state linearity, one can write the total exitance, Id , corresponding to the incident irradiance I as Id = (I − Q − U − V ) × Id,N + Q × Id,Q + U × Id,U + V × Id,V .

(A.4)

(I)

APPENDIX B IT

18

2

3

41

Therefore, one can find the matrix representation 

1

DERIVATION

 Id    Qd     Ud   Vd





  Id,N       Qd,N   =     Ud,N     Vd,N |

Id,Q − Id,N

Id,U − Id,N

Qd,Q − Qd,N

Qd,U − Qd,N

Ud,Q − Ud,N

Ud,U − Ud,N

Vd,Q − Vd,N

Vd,U − Vd,N {z F

5

Appendix B. It

8

9

10

11

I    Q     U    V

(A.5)

Finally, note that if one replaces exitance by radiance [cf. Eq. (24) of [23]], it will be possible to find the Stokes scattering function relative to radiance FL

7

Id,V − Id,N      Qd,V − Qd,N     Ud,V − Ud,N     Vd,V − Vd,N }



with F the Stokes scattering function.

4

6



(i)

derivation (i)

(i)

(i)

Given the Stokes vector, It−1 (λ, rb1 , . . . , rbt−1 ), produced at transition t − 1 and at the end of transition t, (i)

(i)

(i)

one has to derive It (λ, rb1 , . . . , rbt ).

(i) All the Stokes vectors corresponding to all the discrete directions, {I˜t (b nt−1 , ∆Ωn )}n=1,...,Nt , are computed

(i) (i) in Eq. (10). The choice among the directions is done therefore according to the intensity of I˜t , noted I˜t .

The probability of transition to the solid angle ∆Ωn , pλ (b nt−1 , ∆Ωn ) is therefore defined as follows

pλ (b nt−1 , ∆Ωn ) =

(i) I˜t (b nt−1 , ∆Ωn ) , n = 1, . . . , Nt . N t X (i) I˜t (b nt−1 , ∆Ωm )

(B.1)

m=1

12

13

14

(i)

b t = ∆Ωn∗ is therefore chosen randomly according to pλ . It A direction n

polarization state is the same as

(i) (i) I˜t (b nt−1 , ∆Ωn∗ ), however they do not have the same intensity. In fact, in order to conserve energy, It takes Nt X (i) all the intensity exiting the leaf, i.e. I˜t (b nt−1 , ∆Ωm ). Therefore, one finds, m=1

15

(i)

(i)

(i)

It (λ, rb1 , . . . , rbt ) =

Nt X

(i) I˜t (b nt−1 , ∆Ωm )

m=1 (i) I˜t (b nt−1 , ∆Ωn∗ )

(i)

16

Using Eq. (Appendix B) recursively allows to find It

17

in the following to derive the radiance angular distribution.

(i) I˜t (b nt−1 , ∆Ωn∗ ).

(B.2)

(i)

as a function of the incident irradiance Is . It allows

APPENDIX C

18

1

2

DERIVATION OF LT (λ) BASED ON RAY TRACING AT λ0

Appendix C. Derivation of Lt (λ) based on ray tracing at λ0 b o ) is the average radiance at transition t over all Bases on Eq. (14), it is possible to assume that Lt (λ0 , n

the possible positions of interceptions from the first to t, {b r1 , . . . , rbt }.

b o ) = Eλ0 ,b b1 , . . . , rbt , n b o )] , Lt (λ0 , n r1 ,...,b rt [Lt (λ0 , r

3

4

5

6

7

8

locations within canopy.

For the wavelength λ0 , let fλ0 ,t (b r1 , . . . , rbt |b ns , rbs ) be defined as the distribution of the ray interceptions

b s , and initial position at the top of the canopy, rbs . It can from the first to t given the incident ray direction, n

be approximated as

fλ0 ,t (b r1 , . . . , rbt |b ns , rbs ) ≈

12

13

14

15

16

Nr Y t 1 X δ (i) (b rp ), Nr i=1 p=1 rbp

where δrb(i) (b rp ) is the Dirac function, i.e. for a given space S within the canopy, C, p

11

(C.1)

b1 , . . . , rbt over the possible where Eλ0 ,b r1 ,...,b rt is the expectation operator varying the interseption positions r

9

10

42

(i)

and only if rbp ∈ S.

(C.2) R

δ

(i) rp ∈S r bp

(b rp )db rp = 1, if

It is therefore possible to write Eq. (C.1) as

b o) = Lt (λ0 , n

Z

b1 ,...,b r rt ∈C

(i)

b o )fλ0 ,t (b Lt (λ0 , rb1 , . . . , rbt , n r1 , . . . , rbt |b ns , rbs )db r1 . . . db rt .

(C.3)

It is clear that replacing fλ0 ,t by its expression given by Eq. (C.2) in Eq. (C.3) allows to retrieve the exb o ) given by Eq. (14), proving the validity of the derivation. Eq. (C.3) will be therefore used pression of Lt (λ0 , n b o ) based on the Monte Carlo ray tracing simulation done for the wavelength in the following to derive Lt (λ, n

REFERENCES

17

λ0 .

b o) = Lt (λ, n =

= = = 1

2

3

4

5

6

43

(i)

(i)

b1 , . . . , rbt , n b o )] Eλ,b r1 ,...,b rt [Lt (λ, r Z (i) b o )fλ,t (b Lt (λ, rb1 , . . . , rbt , n r1 , . . . , rbt |b ns , rbs )db r1 . . . db rt b1 ,...,b r rt ∈C

}| { fλ,t (b r1 , . . . , rbt |b ns , rbs ) b o) Lt (λ, rb1 , . . . , rbt , n fλ0 ,t (b r1 , . . . , rbt |b ns , rbs )db r1 . . . db rt fλ0 ,t (b r1 , . . . , rbt |b ns , rbs ) b1 ,...,b r rt ∈C

b1 , . . . , rbt , n b o )ωλ0 ,λ (b Eλ0 ,b r1 , . . . , rbt )] r1 ,...,b rt [Lt (λ, r

Nr 1 X (i) (i) (i) (i) (i) b o )ωλ0 ,λ (b L (λ, rb1 , . . . , rbt , n r1 , . . . , rbt ). Nr i=1 t (i)

11

12

(i)

(i)

b ), it remains to determine ωλ0 ,λ and therefore fλ,t (same expression for in Section 3.2, thus to find Lt (λ, n

b p−1 to n b p , p = 2, . . . , t. For a given value of p, the fλ0 ,t ). The latter is the probability of transitions from n b p ). Thus, it transition probability depends only on the intercepting leaf scattering phase function, πλ (b np−1 , n

is clear that varying p, the πλ probabilities are independent. Thus, one can write

fλ,t (b r1 , . . . , rbt |b ns , rbs ) =

t Y

p=2

b p ). πλ (b np−1 , n

(C.5)

The same expression can be written for fλ0 ,t , replacing πλ by πλ0 , therefore ωλ0 ,λ can be written as follows

(i)

(i)

ωλ0 ,λ (b r1 , . . . , rbt ) =

9

10

(i)

(C.4)

b o ) can be estimated like Lt (λ0 , rb1 , . . . , rbt , n b o ) using the ray tracing as explained Lt (λ, rb1 , . . . , rbt , n

7

8

ωλ0 ,λ (b r1 ,...,b rt )

z

Z

(i) (i) t Y bp ) πλ (b np−1 , n (i)

p=2

(i)

bp ) πλ0 (b np−1 , n

.

(C.6)

b ) using the same tracing ωλ0 ,λ can be estimated from the tracing, therefore it is possible to derive Lt (λ, n

b ), based on the weighted sampling presented in Eq. (C.4). as Lt (λ0 , n References

13

[1] S. Chandrasekhar, Radiative Transfer, Dover, New-York, 1950.

14

[2] J. Ross, T. Nilson, Concerning the theory of plant cover radiation regime, Investig Atm Phys Acad Sci

15

16

17

ESSR 4 (1963) 42–64. [3] J. Ross, T. Nilson, A mathematical model of the radiation regime of plant cover, Actinometry and atmospheric optics (1968) 263–281.

REFERENCES

18

1

44

[4] F. Nadal, F.-M. Br´eon, Parameterization of surface polarized reflectance derived from POLDER spaceborne measurements, IEEE Transactions on Geoscience and Remote Sensing 37 (1999) 1709–1718.

2

[5] J.-L. Roujean, J. Leon-Tavares, B. Smets, P. Claes, F. C. D. Coca, J. Sanchez-Zapero, Surface albedo and

3

toc-r 300 m products from PROBA-V instrument in the framework of Copernicus Global Land Service,

4

Remote Sensing of Environment 215 (2018) 57–73, doi:https://doi.org/10.1016/j.rse.2018.05.015.

5

6

7

8

9

10

[6] A. Kuusk, The hot-spot effect of a uniform vegetative cover, Soviet Journal of Remote Sensing 3 (4) (1985) 645–658. [7] A. Kallel, T. Nilson, Revisiting the vegetation hot spot modeling: Case of Poisson/Binomial leaf distributions, Remote Sensing of Environment 130 (2013) 188 – 204, ISSN 0034-4257. [8] V. Vanderbilt, S. Ustin, Polarization of Light by Vegetation, chap. Monte Carlo Methods, Springer-Verlag, ISBN 9783642753916, 191–228, 1991.

11

[9] Z. Sun, Z. Peng, D. Wu, Y. Lv, Photopolarimetric properties of leaf and vegetation covers over a wide

12

range of measurement directions, Journal of Quantitative Spectroscopy and Radiative Transfer 206 (2018)

13

273–285.

14

[10] A. Kallel, J. P. Gastellu-Etchegorry, Canopy polarized BRDF simulation based on non-stationary Monte

15

Carlo 3-D vector RT modeling, Journal of Quantitative Spectroscopy and Radiative Transfer 189 (2017)

16

149–167.

17

[11] J. P. Gastellu-Etchegorry, N. Lauret, T. Yin, L. Landier, A. Kallel, Z. Malenovsk, A. A. Bitar, J. Aval,

18

S. Benhmida, J. Qi, G. Medjdoub, J. Guilleux, E. Chavanon, B. Cook, D. Morton, N. Chrysoulakis,

19

Z. Mitraka, DART: Recent Advances in Remote Sensing Data Modeling With Atmosphere, Polarization,

20

and Chlorophyll Fluorescence, IEEE Journal of Selected Topics in Applied Earth Observations and Remote

21

Sensing 10 (6) (2017) 2640–2649.

22

23

24

25

[12] P. North, Three-Dimensional Forest Light Interaction Model Using a Monte Carlo Method, IEEE Transactions on Geoscience and Remote Sensing 34 (946–956). [13] P. Lewis, Three-dimensional plant modelling for remote sensing simulation studies using the Botanical Plant Modelling System, Agronomie-Agriculture and Environment 19 (1999) 185–210.

REFERENCES

26

1

2

3

4

5

6

7

45

[14] J. Gastellu-Etchegorry, V. Demarez, V. Pinel, F. Zagolski, Modeling radiative transfer in heterogeneous 3-D vegetation canopies, Remote Sensing of Environment 58 (2) (1996) 131–156. [15] W. A. Allen, H. W. Gausman, A. J. Richardson, J. R. Thomas, Interaction of Isotropic Light with a Compact Plant Leaf∗, J. Opt. Soc. Am. 59 (10) (1969) 1376–1379. [16] W. A. Allen, H. W. Gausman, A. J. Richardson, Mean Effective Optical Constants of Cotton Leaves∗, J. Opt. Soc. Am. 60 (4) (1970) 542–547. [17] S. Jacquemoud, F. Baret, PROSPECT: A model of leaf optical properties spectra, Remote Sensing of Environment 34 (2) (1990) 75–91.

8

[18] L. Bousquet, S. Lach´erade, S. Jacquemoud, I. ¨ ol Moya, Leaf BRDF measurements and model for specular

9

and diffuse components differentiation, Remote Sensing of Environment 98 (2) (2005) 201 – 211, ISSN

10

11

12

13

14

0034-4257. [19] Y. M. Govaerts, S. Jacquemoud, M. M. Verstraete, S. L. Ustin, Three-dimensional radiation transfer modeling in a dicotyledon leaf, Appl. Opt. 35 (33) (1996) 6585–6598. [20] S. L. Ustin, S. Jacquemoud, Y. Govaerts, Simulation of photon transport in a three-dimensional leaf: implications for photosynthesis, Plant, Cell & Environment 24 (10) (2001) 1095–1103.

15

[21] Y. Xiao, D. Tholen, X.-G. Zhu, The influence of leaf anatomy on the internal light environment and photo-

16

synthetic electron transport rate: exploration with a new leaf ray tracing model, Journal of Experimental

17

Botany 67 (21) (2016) 6021–6035.

18

[22] C. L. Patty, D. A. Luo, F. Snik, F. Ariese, W. J. Buma, I. L. ten Kate, R. J. van Spanning, W. B. Sparks,

19

T. A. Germer, G. Garab, M. W. Kudenov, Imaging linear and circular polarization features in leaves with

20

complete Mueller matrix polarimetry, Biochimica et Biophysica Acta (BBA) - General Subjects 1862 (6)

21

(2018) 1350 – 1363.

22

23

24

[23] A. Kallel, Leaf Polarized BRDF Simulation based on Monte Carlo 3-D Vector RT Modeling, Journal of Quantitative Spectroscopy and Radiative Transfer 221C (2018) 202–224. [24] S. Jacquemoud, W. Verhoef, F. Baret, C. Bacour, P. J. Zarco-Tejada, G. P. Asner, C. Franois, S. L.

REFERENCES

25

1

2

3

46

Ustin, PROSPECT+SAIL models: A review of use for vegetation characterization, Remote Sensing of Environment 113 (2009) S56–S66, imaging Spectroscopy Special Issue. [25] L. Bousquet, Mesure et mod´elisation des propri´et´es optiques spectrales et directionnelles des feuilles, Ph.D. thesis, UNIVERSITE PARIS 7 - DENIS DIDEROT UFR DE PHYSIQUE, 2007.

4

[26] N. Vilfan, C. V. der Tol, P. Yang, R. Wyber, Z. Malenovsk, S. A. Robinson, W. Verhoef, Extending

5

Fluspect to simulate xanthophyll driven leaf reflectance dynamics, Remote Sensing of Environment 211

6

(2018) 345 – 356.

7

8

9

10

[27] K. Liou, An Introduction to Atmospheric Radiation, International geophysics series, Academic Press, ISBN 9780124514515, 2002. [28] M. Mishchenko, L. Travis, A. Lacis, Multiple Scattering of Light by Particles: Radiative Transfer and Coherent Backscattering, Cambridge University Press, ISBN 9780521834902, 2006.

11

[29] R. L. Thompson, N. S. Goel, Two models for rapidly calculating bidirectional reflectance of complex

12

vegetation scenes: Photon spread (PS) model and statistical photon spread (SPS) model, Remote Sensing

13

Reviews 16 (3) (1998) 157–207.

14

15

[30] M. Disney, P. Lewis, P. North, Monte Carlo ray tracing in optical canopy reflectance modelling, Remote Sensing Reviews 18 (2-4) (2000) 163–196, doi:10.1080/02757250009532389.

16

[31] F. Zhao, Y. Li, X. Dai, W. Verhoef, Y. Guo, H. Shang, X. Gu, Y. Huang, T. Yu, J. Huang, Sim-

17

ulated impact of sensor field of view and distance on field measurements of bidirectional reflectance

18

factors for row crops, Remote Sensing of Environment 156 (2015) 129 – 142, ISSN 0034-4257, doi:

19

https://doi.org/10.1016/j.rse.2014.09.011.

20

[32] J.-L. Widlowski, M. Taberner, B. Pinty, V. Bruniquel-Pinel, M. Disney, R. Fernandes, J.-P. Gastellu-

21

Etchegorry, N. Gobron, A. Kuusk, T. Lavergne, S. Leblanc, P. Lewis, E. Martin, M. M˜ ottus, P. J. R. North,

22

W. Qin, M. Robustelli, N. Rochdi, R. Ruiloba, C. Soler, R. Thompson, W. Verhoef, M. M. Verstraete,

23

D. Xie, The third RAdiation transfer Model Intercomparison (RAMI) exercise: Documenting progress in

24

canopy reflectance modelling, Journal of Geophysical Research 112, doi:10.1029/2006JD007821.

REFERENCES

47

25

[33] A. Kallel, W. Verhoef, S. L. H´egarat-Mascle, C. Ottl´e, L. Hubert-Moy, Canopy bidirectional reflectance

1

calculation based on Adding method and SAIL formalism: AddingS/AddingSD, Remote Sensing of Envi-

2

ronment 112 (9) (2008) 3639–3655.

3

[34] J.-L. Widlowski, T. Lavergne, B. Pinty, M. Verstraete, N. Gobron, Rayspread: A Virtual Laboratory

4

for Rapid BRF Simulations Over 3-D Plant Canopies, in: F. Graziani (Ed.), Computational Methods in

5

Transport, vol. 48 of Lecture Notes in Computational Science and Engineering, Springer Berlin Heidelberg,

6

ISBN 978-3-540-28122-1, 211–231, 2006.

7

[35] J.-L. Widlowski, M. Robustelli, M. Disney, J.-P. Gastellu-Etchegorry, T. Lavergne, P. Lewis, P. North,

8

B. Pinty, R. Thompson, M. Verstraete, The RAMI On-line Model Checker (ROMC): A web-based bench-

9

marking facility for canopy reflectance models, Remote Sensing of Environment 112 (3) (2008) 1144–1150,

10

doi:DOI: 10.1016/j.rse.2007.07.016.

11

[36] G. H. Mohammed, R. Colombo, E. M. Middleton, U. Rascher, C. van der Tol, L. Nedbal, Y. Goulas,

12

O. P´erez-Priego, A. Damm, M. Meroni, J. Joiner, S. Cogliati, W. Verhoef, Z. Malenovsk´ y, J.-P. Gastellu-

13

Etchegorry, J. R. Miller, L. Guanter, J. Moreno, I. Moya, J. A. Berry, C. Frankenberg, P. J. Zarco-Tejada,

14

Remote sensing of solar-induced chlorophyll fluorescence (SIF) in vegetation: 50 years of progress, Remote

15

Sensing of Environment 231 (2019) 111 – 177.

16

[37] S. Y. Kotchenova, E. F. Vermote, R. Matarrese, J. Frank J. Klemm, Validation of a vector version of

17

the 6S radiative transfer code for atmospheric correction of satellite data. Part I: Path radiance, Applied

18

Optics 45 (26) (2006) 6762–6774.

1 Author declaration 1. Conflict of Interest Potential conflict of interest exists: No conflict of interest exists. I wish to confirm that there are no known conflicts of interest associated with this publication and there has been no significant financial support for this work that could have influenced its outcome. 2. Funding

No funding was received for this work.

3. Intellectual Property I confirm that I have given due consideration to the protection of intellectual property associated with this work and that there are no impediments to publication, including the timing of publication, with respect to intellectual property. In so doing I confirm that I have followed the regulations of my institutions concerning intellectual property.

5. Authorship

I confirm that I am the sole author of the manuscript

I the undersigned agree with all of the above.

Author’s name (Fist, Last)

1. ___Abdelaziz Kallel

Signature

Date

December 19th , 2019

Declaration of interests ☐x The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. ☐The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: