DART : a 3D model for remote sensing images and radiative budget of earth surfaces

Modeling the radiative behavior and the energy budget of land surfaces is relevant for many scientific domains such as the study of vegetation functioning with remotely acquired information. DART model (Discrete Anisotropic Radiative Transfer) is developed since 1992. It is one of the most complete 3D models in this domain. It simulates radiative transfer (R.T.) in the optical domain: 3D radiative budget and remote sensing images (i.e., radiance, reflectance, brightness temperature) of vegetation and urban Earth surfaces, for any atmosphere, wavelength, sun/view direction, altitude and spatial resolution. It uses an innovative multispectral approach (flux tracing, exact kernel, discrete ordinate techniques) over the whole optical domain. Here, its potential is illustrated with the case of urban and tropical forest canopies. Moreover, three recent improvements in terms of functionality and operability are presented: account of Earth/Atmosphere curvature for oblique remote sensing measurements, importation of 3D objects simulated as the juxtaposition of triangles with the possibility to transform them into 3D turbid objects, and R.T. simulation in landscapes that have a continuous topography and landscapes that are non repetitive. Finally, preliminary results concerning two application domains are presented. 1) 2D distribution of the reflectance, brightness temperature and radiance measured by a geostationary satellite over a whole continent. 2) 3D radiative budget of natural and urban surfaces with a DART energy budget (EB) component that is being developed. This new model, called DARTEB, is intended to simulate the energy budget of land surfaces.


8
Knowledge of the radiative behavior and the energy budget of land surfaces is essential 9 for studying the functioning of natural and urban surfaces with remotely acquired 10 information. Account of their 3D nature is often essential because in most cases these 11 surfaces are not isotropic. For example, it has long been known that the albedo of a 12 canopy with an anisotropic Bi-directional Reflectance Factors (BRF) may be 13 underestimated by as much as 45% if it is computed with nadir reflectance only (Kimes 14 and Sellers, 1985). Radiative transfer (R.T.) models have the potential for correcting this 15 type of error provided they account for the three dimensional (3D) nature of Earth 16 surfaces. Neglect of the 3D structure of canopies can lead to large errors on the 3D 17 radiation budget and remote sensing measurements. For example, for vegetation BRF and 18 directional brightness temperature (DTDF) distribution functions, errors can be as large as 19 50%, depending on instrumental (e.g., view and sun directions) and experimental (e.g., 20 vegetation heterogeneity) conditions (Gastellu-Etchegorry et al., 1999). The problem is 21 similar for urban canopies due to their strong spatial heterogeneity. The application of 22 R.T. modeling to urban surfaces is important in the context of the advent of satellite 23 sensors with spatial and spectral resolutions that are more and more adapted to urban 24 characteristics such as building dimensions and temperature spatial variability. It explains 25 the numerous works conducted in the field of remote sensing of urban surfaces (Soux et 26 al., 2004;Voogt and Oke, 1998). The use of descriptions with qualitatively based land use 27 data instead of more fundamental surface descriptors is a source of inaccuracy for 28 modeling BRFs and DTDFs (Voogt and Oke, 2003). 29 R.T. models are essential tools for assessing accurately radiative quantities such as the 30 exitance, the irradiance and remote sensing measurements in the optical and thermal 31 domains. However, in order to meet this objective, models must account for the three 32 dimensional (3D) nature of Earth surfaces. Here, we consider vegetation canopies and urban 33 canopies. This consideration of the 3D architecture of Earth surfaces is possible with the so-34 called 3D models. Generally speaking, the latter ones are intended to be accurate, robust and 35 more comprehensive than other models. Ideally, they should be used in place of other 36 models. However, they are often more difficult to manage, both in terms of computation 37 time and landscape description. Moreover, when dealing with specific situations, one needs 1 that the model be accurate and robust, but one does not necessarily need that the model be 2 comprehensive. This explains that in many cases, the objective of 3D models is to calibrate 3 models that are simpler to manage in terms of landscape description, computation time, etc. 4 Once calibrated, these models can meet the required accuracy levels. 5 These remarks stress the usefulness of 3D R.T. models. A number of 3D models is being 6 developed by the scientific community (Widlowski et al., 2007). Usually, they are designed 7 for a given type of landscape (e.g., natural or urban), with or without topography and 8 atmosphere, for a specific type of application (e.g., remote sensing or radiation budget), and 9 for a given spatial resolution of analysis (e.g., simulation of a tree crown with or without 10 branches). Moreover, remote sensing models are usually designed for a given spectral 11 domain (e.g., sun reflected spectral domain or thermal infrared spectral domain). 12 A common problem when designing a R.T. model is to assess to which level of detail 13 landscapes must be taken into account. This is very especially important in term of 14 landscape simulation. However, it affects also the mathematical formulation of R.T. 15 Generally speaking, one should take into account all landscape elements that have a 16 significant influence on the application we are dealing with (i.e., landscape radiative budget 17 or remote sensing measurement). In practice, the answer can be complex. For example, 18 when simulating remote sensing measurements of heterogeneous Earth surfaces in the 19 visible spectral domain, is it necessary to simulate the atmosphere within a unique "Earth -20 Atmosphere" system in order to simulate with a good accuracy its complex interaction with 21 earth surfaces? 22 Some of these aspects are discussed here with the brief presentation of DART (Discrete 23 Anisotropic Radiative Transfer) model. This is one of the most complete 3D models 24 designed for simulating the radiative budget and the satellite observations of the land 25 surfaces in the visible, near infrared and thermal infrared of land surfaces. It was originally 26 developed (Gastellu-Etchegorry et al., 1996) for simulating remote sensing images of 3D 27 vegetation canopies in the visible / near infrared (NIR) spectral domain. Afterwards, it was 28 extended to the thermal infrared domain and to the simulation of any landscape: urban or 29 natural, with atmosphere and topography. As a result, the present DART model simulates 30 the radiation budget and remote sensing images of vegetation and urban canopies, for any 31 experimental (sun direction, canopy heterogeneity, topography, atmosphere, etc.) and 32 instrumental (view direction, spatial resolution, etc.) configuration. 33 After a brief presentation of DART, two types of applications are discussed: urban and 34 forest canopies. This is followed by the presentation of three recent improvements: reflectance models (e.g., Flight (North, 1996), Sprint (Thompson and Goel, 1998) DART simulates R.T. in heterogeneous 3-D landscapes with the exact kernel and discrete 20 ordinate methods. It uses an iterative approach: radiation intercepted in iteration "i" is 21 scattered in iteration "i+1". Any landscape is simulated as a rectangular matrix of 22 parallelepipedic cells. Figure 1  For the 3 modes, the atmosphere can be treated as a propagating medium or as an interface. 25 In any case, landscape irradiance has 2 components: direct sun W( s ,x,y) and atmospheric 26 W a ( n ,x,y) source vectors. W( s ) propagates along direction ( s ). W( s ) and W a ( n ) are 27 simulated from a fictitious cell layer on top of the scene (Figure 1), with values equal to: 28 W( s ) = E s ( s ).| s |.x.y and W a ( n ) = L a ( n ).| n |.x.y. n 29 where x.y is the area of the cell face,  s =cos s ,  n =cos n , E s ( s ) is the solar constant at the 30 top of the scene, and  s denotes the solar incident direction. L a ( n ) is the atmospheric 31 radiance along direction ( n ), with n [1 N'], where N' is the number of downward discrete 32 directions. It is null at the top of the atmosphere. 33 DART landscape modeling is as independent as possible from the RT modeling in order to 34 allow DART to simulate RT on landscapes simulations that are generated by any other 35 model. It can combine them with its own simulated landscapes. Imported landscapes and 36 landscape elements can be edited to some extent. Geometric transformations (i.e., 3D 37 translation, 3D homothety, 3D rotation) can be applied and optical properties can be 38 assigned.

39
DART uses 2 complementary approaches for simulating landscapes: 40 -Juxtaposition of cells that contain one or several turbid medium (i.e., cloud of infinitely 41 small planar elements). This is useful for simulating volumes of foliar elements such as 42 grass and tree crown. A turbid medium is characterized by its volume density, an 43 angular distribution and optical properties (i.e., abaxial reflectance, adaxial reflectance, 44 transmittance). 45 -Juxtaposition of translucent triangles. This is useful for simulating the ground, the 1 branches, the urban surfaces (i.e., walls and roofs) and also foliar elements. A single cell 2 can contain several turbid medium and several triangles or part of them. 3 In addition to the atmosphere and to the ground and its topography, DART simulates 4 4 types of landscape elements: 5 -Trees with exact or random locations and specific optical properties. Each tree is made 6 of a trunk, possibly with branches, simulated as triangles, and a tree crown simulated as 7 the juxtaposition of turbid cells. Tree crown can have a number of predefined shapes 8 (e.g., ellipsoid, cone, trapezoid, etc.), with specific vertical and horizontal distributions 9 of leaf volume density. Trees with different geometric and optical properties can be 10 mixed.

11
-Grass or crops. They are simulated as a volume of turbid medium. This volume can 12 located anywhere in space (x,y,z). 13 -Urban elements (i.e., houses, roads,...). The basic element is a house with walls 14 characterized by the location of their 4 upper corners and roofs characterized by their 4 15 upper corners. 16 -Water elements (i.e., river, lake). They are simulated as surfaces with any optical 17 property (e.g., anisotropic reflectance possibly with specular component). 18 Generally speaking, two types of radiation interaction take place. (1)  simulations where DART input variables A 1 , A 2 , A 3 ,.. take N 1 , N 2 , N 3 ,…, respectively. 5 Results are stored in a LUT (Look Up Table) for further display and / or processing. 6 -Manipulation of DEM (Digital Elevation Model): this is used for importing, creating 7 and resampling DEMs. 8 -Simulation of foliar spectra with the Prospect model (Jacquemoud and Baret, 1990). 9 -Simulation of scene spectra (reflectance, brightness temperature, radiance). It can be 10 computed using a single DART simulation that is conducted with N spectral bands, or 11 with the help of the sequencer for running N DART simulations with 1 spectral band 12 each. 13 -Simulation of broadband reflectance, brightness temperature, radiance, irradiance,… It 14 is the sum of a few DART simulated narrow spectral bands, possibly weighted by 15 sensor spectral sensitivity. 16 -Importation of land cover maps for a direct simulation of DART scenes. 17 DART models and tools are managed with a user friendly Graphic User Interface (GUI: 18 Figure 2) to input all necessary parameters (e.g., view and illumination conditions) and to 19 specify the required products. They can be also managed with command lines such as 20 scripts written in Python programming language. DART computation code is written in 21 C++ language (more than 300 000 lines of code). The GUI is written in Java language.  Measurements are essentially images for the "Flux tracking" mode and waveforms for 31 the "Lidar" mode. Images are simulated in the focal plane of the satellite sensor. The 32 cross section of the emitters and scatterers at the origin of the signal are taken into 33 account for improving the image quality, especially for scenes with marked 3D 1 architectures (urban elements, topography). A bi-linear interpolation method is used for 2 projecting the horizontal upper grid of the scene onto an over sampled grid in the 3 sensor plane, at any altitude (BOA to TOA). 4 -Radiative budget: 3D, 2D and 1D distribution of the radiation that is intercepted, 5 scattered, emitted and absorbed. 6 7 DART potential for simulating remote sensing images of urban and forest canopies is 8 illustrated here with 2 examples. 9

Urban canopy 10
The example shown here is derived from the CAPITOUL project of Meteo France ( that urban elements in the data base are not individual houses or buildings but unrelated 20 walls and roofs was a difficulty. district. They are simulated for a sensor at the bottom of the atmosphere (i.e., BOA image) 28 and for a sensor at the top the atmosphere (TOA). The bluish tone of the TOA image, 29 compared to the BOA image, is due to the fact that atmosphere scatters more in the blue 30 than in the red spectral domain. DART images realism is very useful for verifying that the 31 land surface is correctly simulated. Moreover, it helps also for testing the coherence of the 32 RT. For example, at the 1 st scattering order, the reflectance of shadows is null. 33

Forests 34
One considers here DART simulations of the reflectance of a tropical forest (Sumatra, 35 Indonesia) in the visible (VIS), near infrared (NIR) and short wave infrared (SWIR) spectral 36 domains. Optical and geometric characteristics are given in Gastellu-Etchegorry et al., 1998. 37 Sun direction is characterized by a 35° off-nadir angle ( s ) and a 200° azimuth angle ( s ). 38 Figure 5 shows DART VIS images for three viewing directions: nadir viewing direction 39 ( v =0°), sun backscattering viewing direction ( v =35°,  v =200°), also called hot spot 40 configuration, and specular configuration (   Figure 6 shows DART NDVI (Normalized Vegetation Index) and VIS / NIR / SWIR 23 reflectance values in the principal and perpendicular solar planes, for several sun off-nadir 24 directions. Reflectance values share some similar characteristics: (1) a well marked bowl 25 shape in the principal solar plane with a spectrally dependent minimum for a direction 26 between the specular and nadir directions, (2) a strong maximum in the hot spot direction, 27 more marked in the VIS than in the NIR and SWIR, (3) a systematic increase of reflectance if 28 view zenith angles exceed a threshold value (e.g. 50° in the VIS if  s =0°), which depends on 29 the spectral domain and the sun direction, (4) an azimuth symmetry relative to the principal 30 solar plane, and (5) a relatively small variability in the perpendicular solar plane. 31 Reflectance is maximal at hot spot ( vis =0.073,  nir =0.615,  swir =0.324) and minimal for 32 directions between nadir and specular configurations ( vis =0.025,  nir =0.354,  swir =0.12), 33 whereas NDVI is minimal (0.79) in the hot spot direction and maximal (0.87) in the 34 specular direction. 35 In addition to DART simulations, Figure 6 shows reflectance simulations carried out with 36 the well known SAIL model (Verhoef, 1984). For this model, vegetation is an homogeneous 37 turbid medium. Figure 6 allows one to stress the role of canopy architecture on forest 38 reflectance. Its neglect by the SAIL model explains that DART reflectance is much smaller. 39 Reflectance differences depend a lot on view direction. Smaller differences occur in the hot 40 spot directions because no shadows occur for these viewing directions. Outside the hot spot 41 configuration, for  s =0°, mean relative reflectance difference is around 60-70% in the VIS, 42 50% in the SWIR and 25% in the NIR. Differences tend to increase with sun off-nadir angle. 43 The bowl shape of VIS and SWIR reflectance is more marked with SAIL than with DART. It 44 corresponds to the fact that in the forward principal plane, for large off-nadir viewing 45 angles  v , an increase of  v implies that DART reflectance increases less than SAIL 46 reflectance, especially for large off-nadir sun angles. Indeed, the canopy 3-D heterogeneous 1 structure ensures that an increase of  v leads to an increase of the proportion of shaded tree 2 crowns that are viewed in the principal plane. This effect tends to be less marked in the NIR 3 than in the VIS because the role of shadows is less marked in the NIR, due to the increased 4 occurrence of multiple scattering processes. Moreover, the presence of shadows explains 5 that with oblique sun illumination, minimal values of DART reflectance are more shifted 6 towards the specular direction than those of SAIL. This is also true for NDVI. Image simulation is very useful for understanding forest reflectance behavior with 7 experimental and instrumental parameters. Here, this is shown for sky radiation. Figure 7  8 shows DART nadir NIR images of part of the tropical forest, for 2 extreme atmosphere 9 conditions: SKYL equal to 0 and 1, with  Atmosphere irradiance SKYL Total irradiance . We note that some 10 tree crowns are invisible with SKYL=0 (bottom of Figure 7), because they are shaded, and 11 12 13  For an off-nadir angle , an atmosphere path length between altitudes z a and z b is smaller 39 than in an horizontal atmosphere (i.e., z/|µ|, with µ=cos and z = |z a -z b |). Moreover, 40 gas and aerosol vertical distributions are not constant. Here, atmosphere is assumed to be 41 spherical ( Figure 27). {µ < 0, z < 0} for the 2 nd one. 10 The optical depth of path AB is: 11 . , with t the altitude relative to z a 12 and l the path length from A, we have: = − . μ + . μ + t. (t + 2 )l if µ > 0 and : The computation of (z a ,z,µ) could be solved with an integration by parts: 15 This approach requires the derivative '(z) of the extinction coefficient. However, '(z) is 1 known only for an ideal case such as an exponential atmosphere. In order to work with any 2 atmosphere (e.g., non exponential vertical profile of O 3 ), (z a ,z,µ) is computed with: 3 Δ ( , Δ , μ) = ( ) . = ( ) . | Δ |, with ( ) the mean extinction coefficient within [z a 4 z a +z] 5 Thus, before simulating the atmosphere RT, DART computes µ sph (z i ,µ j ) for all I atmosphere 6 layers and all J discrete directions, including the sun direction, with i  [1 I] and j  [1 J]. 7 The use of µ sph (z i ,µ j ) allows one to treat the atmosphere as an horizontal plane for simulating 9 the atmosphere RT: any path r(z i ,µ j ) that is computed for an horizontal atmosphere layer 10 We have:

Transformation of "3D triangle objects" into "3D turbid objects" 13
As already mentioned DART can import 3D objects ( Figure 9) that are simulated as groups 14 of triangles. In many cases, especially for vegetation, these objects are simulated with 15 tremendous numbers of triangles (e.g., 10 6 triangles). This is very costly is terms of 16 computation time and volume memory, especially if we work with forests… In order to 17 solve this problem, we designed a module that transforms 3D objects simulated as triangles 18 (i.e., 3D triangle objects) into 3D objects simulated as turbid medium (i.e., 3D turbid objects). 19 Figure 10 shows the schematic approach. It is reminded that a cell filled with turbid medium 20 is characterized by its LAI (Leaf Area Index), its leaf angle distribution (LAD) and the leaf 21 optical properties (e.g., transmittance and adaxial and abaxial reflectance, possibly with 22 specular parameters). The "transformation" module computes the LAI and LAD of all leaf 23 24 25 26 Fig. 9. Examples of 3D trees that are imported by DART. 27 elements within each cell of a 3D cell matrix. This is done for all or part of the groups of 1 triangles of the "3D triangle" object. Optical properties of the turbid medium are those of the 2 triangles. On the other hand, the LAD is either the one that is computed or a predefined 3 LAD. The possibility to use a predefined LAD is well adapted to the case of 3D triangle 4 objects with low volume densities of triangles. 5 6 7 Fig. 10. Schematic representation of the transformation of a "triangle cell" into a "turbid 8 cell". 9 Here, the transformation of a 3D triangle scene into a turbid 3D turbid scene is illustrated by 10 DART color composite images of the citrus tree "3D triangle scene" and "3D turbid scene" 11 ( Figure 11). The associated 2D reflectance polar plots are shown also. For these 12 2 cases, RT was simulated in the blue, green and red spectral bands. Results are very 13 encouraging: reflectance values are very close. Actually, reflectance values of the 14 3D triangle and turbid scenes are much closer than the 3D triangle objects contain a lot of 15 triangles. It can be noted that DART images in Figure 11 b and e are duplicated 16 2 x 2 times. This mode of representation of images is often useful for better interpreting 17 simulated images where objects (e.g., trees) and their shadows cross the scene boundaries. 18 Compared to the usual simulation of trees with classical tree crown shapes such ellipsoids 19 or cones, the transformation of 3D triangle tree crowns into 3D turbid tree crowns is very 20 interesting for keeping the 3D architecture of trees. reflectance polar plots. 10

Finite landscapes and continuous infinite landscapes topography 11
DART was designed to operate with infinite scenes that are made of a DART simulated 12 pattern that is periodic. When a ray exits the scene, it re-enters the scene by the scene 13 opposite side. This is the so-called "repetitive topography" method. This approach is used in 14 most 3D models. It works very well with landscapes without topography or if the 15 topography is identical on the landscape opposite sides. Thus, it is erroneous in presence of 16 any topography. It is also erroneous if one wants to simulate a finite landscape without 17 interaction with their neighborhood. We solved these 2 problems by introducing 2 new 18 landscape modeling methods, called "Continuous topography" and "Isolated landscape", 19 respectively. This implied to adapt the 3 DART RT modeling modes (i.e., flux tracing, Monte 20 Carlo and lidar). In short, landscapes can be simulated with 3 methods (Figure 12): 21 -Repetitive topography: the landscape and the topography are periodic. It that case, a 1 landscape with a simple slope is actually simulated as a series of periodic slopes, which 2 implies undesirable illumination (Figure 12.a) and view (Figure 12.b)   continuous, due to the discontinuity of the slope) Isolated scene (no adjacency effects). 7

Radiative budget 8
Different modules and tools were added to DART for better simulating and managing the 9 3D radiative budget. DART computes the different terms of the radiative budget: absorbed 10 energy per cell, intercepted energy per cell, scattered energy per cell and downward and 11 upward energy on top cell faces. Here, this is illustrated by a schematic tree landscape 12 ( Figure 14). DART simulation was conducted in the near infrared with a SKYL equal to 0.2 13 and a sun direction ( s = 160°,  s = 90°). The possibility to obtain 2D and 3D displays is very 14 helpful for understanding the radiation interaction within the tree canopy. For example, 15 here, larger downward radiation occurs in the air cells that are directly illuminated by the 16 sun. This quantity can be larger than 1 because tree crowns scatter radiation through these 17 air cells. Larger interception occurs on the illuminated tree crowns and the ground that is 18 directly illuminated by the sun. Larger scattering occurs where interception is maximal and 19 where local reflectance / albedo is maximal. This explains that scattering by illuminated tree 20 crowns is larger than scattering by the illuminated ground. Indeed, in the near infrared, leaf 21 albedo is larger than ground reflectance. It is interesting to note that the "Broadband" 22 module of DART computes broadband radiation budget with the DART simulated narrow 23 band radiation budgets. Carlo method that is adapted for taking into account the usually strong anisotropy of the 8 phase function of landscape elements. In short, the occurrence probabilities of scattering 9 events that have the same order of magnitude are grouped for obtaining groups that have 10 cumulated probabilities with the same order of magnitude (Gastellu-Etchegorry et al., 2010). 11 The Monte Carlo mode in the DART model was initially developed for assessing the 12 accuracy of DART flux tracing method, using the same simulations of landscapes. Indeed, 13 flux tracing RT modelling requires some simplifying hypotheses for representing multiple 14 scattering. The associated inaccuracy depends on the trade-off between the expected 15 accuracy and computational time of simulations. The advantage of the Monte Carlo 16 approach is to simulate multiple scattering processes as a succession of exactly modelled 1 single scattering processes. 2 In order to simulate Lidar waveforms, the DART Lidar module works with a Gaussian 3 spatial distribution for illumination, and a Gaussian (time) laser pulse. It is suited for small 4 and large footprints. It is very flexible because it inherits major features of the DART model: 5 urban and natural scenes, with topography, atmosphere, etc. Examples of waveforms are 6 presented here for urban and natural scenes. 7 Figure 15 illustrates lidar modelling of a urban area that is the St Sernin district (Figure 3). It 8 shows the image that is simulated with the flux tracing mode and the waveform that is 9 simulated with the lidar module. As expected, the horizontal axis of the waveform gives the 10 altitude of scene elements (buildings, vegetation, ground). The 2 images are very similar. Actually, the degree of similarity depends on the number of 22 photons that are used. As expected, the Monte Carlo mode is usually much more expensive 23 in terms of computation time. An interesting point is that DART simulates the images of the 24 illuminated and view footprints ( Figure 20 d and e). This is very useful for interpreting the 25 simulated waveform (Figure 20 f). As expected, the latter one shows a peak that corresponds 26 to the ground and a peak that corresponds to the tree crowns. 27 The citrus waveform (Figure 16.f) being related to the LAI (Leaf Area Index) vertical profile, 28 Ueberschlag (2010) assessed the potential of DART for retrieving the LAI of forests with an 29 inversion procedure. Results were encouraging. As expected, similarly to the lidar response, 30 the LAI retrieval depends a lot on the location of trees within the footprint, except if the 31 lidar signal is uniform. For a Gaussian lidar signal, Figure 17 shows that in the case of a tree 32 cover the LAI of which is 0.5, the retrieved LAI can vary from 0.33 up to 0.70, depending on 33 the tree location within the footprint. 34 1 Fig. 16. Lidar simulation of citrus trees (a). Images simulated with the Monte Carlo mode 2 (b) and the flux tracing mode (c). Images of Lidar ground (d) and view area (e) footprint. 3 Waveform (f). The lidar has a 3ns pulse duration, a 0.5ns acquisition rate, a 4m footprint 4 radius and a 0.368 Gaussian illumination parameter. The potential of altitude mapping was assessed. Figure 18.a shows an altitude map where 14 the altitude of each pixel is the local higher altitude. This is called the Lidar first return 15 altitude map. Figure 18.b is an interpolation of Figure 18.a for obtaining a 3D display. 16

Radiance and radiative budget at continental scale 4
A method was developed to create automatically for every spectral band, any date and any 5 land area, maps of radiometric products (i.e., radiance, reflectance, brightness temperature) 6 at a continental scale. For that, it realizes a spatial interpolation on a set of georeferenced 7 DART products that are created by running the DART "sequencer" module with time, 8 wavelength and land surface location used as variable parameters. 9 Results shown here are for the African continent, with a geostationary satellite (35800km 10 altitude, 0° N, 17° E) for 2 spectral bands (550nm: Figure 19 and Figure  Winter. 8 In a 1 st step, the tool "Sequencer" created a grid of 20x20 DART simulations at 5 wavelengths 1 (0.4μm, 0.55μm, 0.67μm, 0.9μm, 1.65μm, 11μm), with an automatic computation of the sun 2 and view sensor angles in the local reference system, as a function of the date, time and 3 coordinates. The grid covers almost all of the African continent: "40° S -44° N" in latitude 4 (step of 4.2°) and "20° W -54° E" in longitude (step of 3.7°). In a 2 nd step, the grid of DART 5 products is interpolated (Python script and Numpy module). Results are displayed with a 6 Python script and the module Basemap. 7 Logically, TOA radiances differ a lot from BOA radiances, because of the atmosphere. 8 Moreover, one can note: 9 -BOA radiances are maximal in the East in the morning and West in the evening. During 10 summer, they are larger in the Northern Hemisphere and lower in the Southern 11 hemisphere. The situation is reversed in winter. 12 -TOA radiances tend to be maximal in the East in the morning and in the West in the 13 evening. This effect is less clear than for BOA radiances. During summer, TOA 14 radiances are larger in the Northern Hemisphere and lower in the southern hemisphere. 15 The situation is reversed in winter. 16 The impact of atmosphere on TOA radiances depends on its optical thickness and phase 17 function and on the sun and view directions. The analysis of DART images allows 18 one to verify that this impact is not symmetric from the point of view of the satellite 19 sensor. 20 The Figure 11 and Figure 12 show the maps of luminance to 900nm of the African 21 continent to 4 dates (spring, summer, autumn, winter) and 3:00 (UTC 8h, 12h UTC 16h 22 UTC) defined above. As Figure 9 and Figure 10, this is 100 * 100 maps obtained by 23 interpolation of DART simulations performed for 400 geographical sites, with a 24 Lambertian soil ( = 900 0351: "Brown to dark brown gravelly loam" of USDA Soil 25 Conservation service) and an atmosphere "US Standard" with aerosols 23km visibility. 26 Logically, the TOA radiances and reflectances are much less affected by the atmosphere in 27 the near infrared than in the field of green. The spatial variability of radiances is mainly 28 due to the spatial variability of the illumination of the land surface. The geometric 29 configuration "Sun -Earth" plays a much more important role than the geometric 30 configuration "Sensor -Sun-Atmosphere". This is particularly the case of the phase 31 function of gases and aerosols. This explains the much larger symmetry of the near-32 infrared maps of luminance. 33 Figure 23 shows the annual evolution of BOA (a) and TOA (b) radiances of 6 cities: Algiers 34 (36. Luanda (8.50 ° S, 13.14 ° E, UTC +1). The evolution of the radiance L ref (t) at the nadir of the 37 satellite sensor is used as a reference. Any radiance L(t) larger than L ref (t) indicates that the 38 sensor receives light above this reference. 39 Figure 24 shows the 11µm TOA radiance, and its associated Brightness temperature, of the 40 African continent. It was obtained with a ground surface characterized by an emissivity 41 equal to 1 and a 300K thermodynamic temperature, a US standard atmosphere and an 42 atmosphere water vapor thickness equal to 1.4cm. The simplicity of this configuration 43 explains the spatial symmetry of brightness and temperature. Indeed, one should take into 44 account the actual local emissivity and thermodynamic temperature.

DARTEB energy budget simulation 11
Energy budget modeling is essential for many application domains such as the functioning 12 of land surfaces. It is also essential for obtaining realistic simulations of satellite 13 measurements in the thermal infrared. Indeed, it allows one to obtain the 3D distribution of 14 thermodynamic temperature. This explains present efforts for developing a model, called 15 DARTEB (DART Energy Budget), that simulates the 3D radiative budget of urban and 16 vegetation surfaces, possibly with topography and atmosphere. DARTEB uses the 3-D 17 DART radiative budget and it models all physical phenomena, other than radiation, that 18 contribute to the energy budget: heat conduction, turbulent momentum and heat fluxes, 19 water reservoir evolution, vegetation photosynthesis, evapotranspiration. The example 20 shown here is for a urban canopy. In that case, non-radiative mechanisms that contribute to 21 the energy budget are simulated with the equations of the TEB urban surface scheme 22 (Masson, 2000). This scheme works with a canyon geometry. For example, turbulent fluxes 23 and conduction are computed with classical boundary-layer laws using equations of TEB. 24 Brightness temperature (K) TOA radiance (W/m 2 /sr/µm Conversely to the TEB scheme, DARTEB uses a 3-D cell discretization, in addition to the 1 layer discretization of roofs, walls and roads: modeling is conducted on a DART cell per cell 2 basis. As a result, fluxes are computed for each point of the 3-D scene. The transfer 3 coefficients for turbulent heat and moisture fluxes are identical; they differ from the transfer 4 coefficients for momentum fluxes. For DARTEB, the urban canopy is simulated as the 5 juxtaposition of urban street canyons. Here, we worked with a single urban canyon, for 6 remaining in the validity domain of TEB equations (Masson, 2000). Each surface type (wall, 7 soil, roof) is discretized into several layers for simulating the conduction fluxes to or from 8 the ground and building interiors. The number of layers for road, wall and roof can differ. A 9 minimal number of three layers is advised because temperature gradients can be large and 10 because of the multi-layer structure of the walls and roofs. 11 DARTEB uses a prognostic approach for assessing the 3D radiative budget distribution, and 12 consequently the 3D temperature distribution. Temperature values at time "k -1" are used 13 for computing the 3D TIR and energy budgets at time "k", which allows one to compute the 14 3D temperature distribution at time "k", using the 3D visible and NIR radiation budget at 15 time "k" (Figure 25). DART  West orientation) in Toulouse. 30 The simulated and measured road temperature curves are very similar (Figure 26.a) proportions of the 2 components of the canyon illumination: sun and sky illumination. 19 Here, these components are driven by the atmosphere optical depth and sun zenith angle. 20 However, in the absence of measurements, the atmosphere optical depth is assumed to be 21 constant. 22 The wall (Figure 26. It is also more and more used for radiation budget in urban and natural environments, with 45 the objective to couple it with functioning models. Moreover, DART is also more and more 1 used in the field of education for teaching the physical bases of remote sensing and radiative 2 budget. Generally speaking, the possibility to create very easily urban and natural 3 landscapes is very interesting for scientists who want assess remote sensing measurements 4 in any experimental and instrument conditions. 5 Today, work continues for improving DART, for both the physics, functionalities and 6 computer science aspects. Present ongoing improvements concern the storing and 7 manipulation of LUTs as actual databases, the inversion of satellite images, the Lidar and 8 the design of a new GUI. Planned improvements concern the modeling of the Earth-9 Atmosphere system in presence of clouds. DART will not simulate RT within clouds. These 10 will be considered as interfaces with specific optical properties and geometric dimensions. 11 Other major planned improvements will be the simulation of RT adapted to water bodies 12 and to fire. 13

14
This work was supported by the Centre National d'Etudes Spatiales (CNES) in the frame of 15 the project "High spatial resolution of the land surfaces with a geostationary satellite". The 16 authors are also thankful to team of scientists the Space Center for Biosphere Studies 17 (CESBIO) who contribute to the development of DART. 18 19 For an Earth target T and a satellite sensor S, DART computes the local incidence angle , 20 the sensor off-nadir direction  and the sun azimuth angle in the Earth surface reference 21 system as a function of the locations (latitude, longitude, altitude) of T and S, with the 22 assumption that the Earth is an ellipsoid. 23

Annex: Computation of view incidence angle
Let us note O the Earth center and Oz the "South -North" axis perpendicular to the 24 equatorial plane. The axis Ox is perpendicular to Oz and is in the plane {OT , Oz}. The plane 25 that contains the target T and the vertical Oz is an ellipse defined by Let R target and R sat be the Earth radius at points t and s that are the vertical projection on the 33 Earth ellipsoid of target T and satellite S respectively. Points t and s are defined by: 34 -t: Lat t , Long t . Target T altitude is H target above the Earth surface 35 -s: Lat sat , Long sat . Satellite altitude is H sat above the Earth surface 36 Zenith angles of vectors OS and OT , with s and t the vertical projections of points S and T, 37 are: 38  sat = -Lat sat and  t = -Lat t respectively (in radians) In the pre-defined Earth system, point t is defined by OT {R t ( t ).sin( t ), 0, R t ( t ).cos( t )}. The 13 unit vector W ,∥ of the tangent to the ellipse along the axis Ox at point t is defined by: 14  This allows one to compute the incidence angle , the sensor off-nadir angle  and the 9 satellite azimuth  sat in the local system: