FR2823344A1 - Procede d'obtention d'une representation numerique, simulee de l'etat radiatif d'une scene tridimensionnelle heterogene - Google Patents

Abstract

L'invention concerne un procédé d'obtention d'une représentation numérique simulée de l'état énergétique radiatif d'une scène (5) tridimensionnelle hétérogène découpé en une pluralité de cellules (71, 72, 73) élémentaires. On définit une pluralité de secteurs d'incidence en nombre inférieur aux directions incidentes. On calcule, cellule après cellule, une énergie interceptée de secteur égale à la somme des énergies interceptées des rayons incidents dans le secteur d'incidence. On utilise une fonction de transfert de diffusion de secteur permettant de calculer une énergie diffusée de secteur au point origine équivalent de diffusion dans la direction de ré-émission OMEGA vj . On calcule l'énergie volumique diffusée comme la somme des énergies diffusées de secteur.

Description

<Desc/Clms Page number 1>

PROCEDE D'OBTENTION D'UNE REPRESENTATION NUMERIQUE SIMULEE DE L'ETAT RADIATIF D'UNE SCENE TRIDIMENSIONNELLE HETEROGENE
L'invention concerne un procédé d'obtention d'une représentation numérique simulée de l'état énergétique radiatif, dans le domaine visible et/ou en dehors du domaine visible-notamment dans le domaine compris de l'ultraviolet à l'infrarouge thermique-, d'une scène tridimensionnelle hétérogène, notamment une portion d'espace susceptible de s'étendre dans le champ de vision d'un capteur satellitaire au-dessus du sol d'une planète telle que la Terre.

Une telle représentation numérique simulée permet soit de déterminer l'état énergétique radiatif en tout point de la scène, soit de réaliser des images numériques de la scène vue chacune d'un point extérieur de la scène (par exemple un capteur satellitaire ou aéroporté) situé selon une direction d'observation (différentes images simulées vues selon différentes directions d'observation sont en général produites simultanément).

L'obtention de telles représentations numériques simulées capables de fournir avec une grande précision radiométrique l'état radiatif d'une scène tridimensionnelle (et non pas seulement l'aspect visuel général comme dans les techniques de synthèse d'images virtuelles) est importante dans certaines applications, par exemple pour le calcul précis de variables biophysiques (biomasse, étude de la végétation...) ou biochimique à partir de mesures satellitaires, ou pour l'étalonnage au sol de capteurs satellitaires à haute résolution spatiale (pouvant maintenant atteindre de l'ordre de lm ou moins).

La publication"Modelling Radiative Transfer in Heterogeneous 3-D Végétation Canopies"J-P. Gastellu-Etchegorry et al Remote Sens. Environ. 58 : 131-156 (1996) décrit un procédé d'obtention de telles images numériques. Il repose sur une approche itérative : à l'itération n, il y a diffusion de l'énergie Wout (O) selon les directions (Q) de l'espace, compte tenu de la distribution spatiale de l'énergie interceptée Wint, et non absorbée, à l'itération précédente n-1. L'itération n=l comprend l'éclairement par une source externe

<Desc/Clms Page number 2>

mono-directionnelle (telle que le soleil), alors que l'itération 2 comprend la diffusion de l'énergie interceptée à l'itération 1, et l'émission thermique, s'il y a lieu. Ce procédé souffre de deux inconvénients :
1) Comme l'indique cette publication, les temps de calcul restent extrêmement longs et les capacités mémoires nécessaires sont extrêmement importantes par rapport aux dimensions de la scène. Ainsi, toute simulation appliquée à un paysage forestier de 50mx50x20m avec des cellules de 20cm, avec 150 directions et avec une précision radiométrique compatible avec la précision usuelle des capteurs (satellitaires et autres) requiert un temps de calcul de l'ordre de la semaine, selon le type de station de travail utilisée. Ce temps calcul est essentiellement du à la méthode des harmoniques sphériques initialement employées pour gérer les diffusions multiples. Cette méthode est en effet très longue et est appliquée à toute cellule et à tout rayon, à chaque itération.

2) Les diffusions multiples sont mal prises en compte, affectant grandement la précision du modèle. La principale cause d'imprécision est la méconnaissance de l'anisotropie du rayonnement incident qui est diffusé.

En effet, toute diffusion dépend de la distribution angulaire du rayonnement incident. Celle-ci n'est bien connue qu'à l'ordre 1, c'est à dire à l'itération 1, car la direction solaire est unique. Pour les itérations n > l, le rayonnement incident peut provenir de n'importe quelle direction et est a priori anisotrope. La mémorisation de la distribution angulaire du rayonnement incident à chaque itération, pour chaque cellule, requiert de très importants volumes mémoire (de l'ordre de la centaine de giga-octets pour un paysage de 50mx50mx20m). En théorie, la décomposition du rayonnement diffusé sous forme d'harmoniques sphériques est supposée résoudre ce problème. En pratique, outre les problèmes de temps calcul et de volume mémoire de stockage sus-indiqués, cette solution est imprécise, car il ne peut être stocké qu'un nombre fini de coefficients de décomposition. En outre, le fait de multiplier une dimension de la scène d'un facteur m conduit à

6 multiplier la quantité d'informations à traiter d'un facteur ni.

Les inconvénients indiqués ci-dessus empêchent donc toute application industrielle pratique.

<Desc/Clms Page number 3>

L'invention vise donc de façon générale à pallier ces inconvénients en proposant un procédé qui permet d'obtenir des représentations numériques simulées de scènes de grandes dimensions, ayant un grand nombre de cellules, avec une précision radiométrique suffisante, et avec des temps de calcul et des capacités de mémorisation compatibles avec une exploitation industrielle pratique.

L'invention a plus particulièrement pour objet un tel procédé permettant de simuler des images de télédétection planétaire-notamment terrestre-depuis l'espace, pour toute configuration expérimentale (paysage hétérogène ou non, directions d'éclairement et d'observation quelconques, etc.).

L'invention a aussi plus particulièrement pour objet un tel procédé permettant d'obtenir des images radiométriques permettant d'étalonner les capteurs satellitaires modernes au sol avant leur lancement.

L'invention vise en particulier à permettre de mieux comprendre les mesures satellitaires et par suite d'améliorer leur conception, leur exploitation et finalement contribuer à une meilleure évaluation des possibilités offertes par les capteurs satellitaires.

L'invention a aussi pour objectif de fournir une représentation tri-dimensionnelle simulée du bilan radiatif du paysage étudié, information essentielle à de nombreuses applications environnementales telles que l'étude du fonctionnement de la végétation (activité photosynthétique, flux gazeux, etc.).

Pour ce faire, l'invention concerne un procédé d'obtention d'une représentation numérique simulée de l'état énergétique radiatif d'une scène tridimensionnelle hétérogène, dans lequel : - on découpe la scène en une pluralité de cellules élémentaires comprenant des cellules, dites cellules turbides, ayant des propriétés radiatives connues supposées homogènes dans l'ensemble de chaque cellule, - on mémorise pour chaque cellule turbide des fonctions de transfert permettant de déterminer, dans chaque direction d'une pluralité discrète de Ny directions de l'espace, dites directions de ré-émission Oy, prédéterminées dans un repère de l'espace indépendant des cellules, une énergie ré-émise par la

<Desc/Clms Page number 4>

cellule turbide en fonction d'un rayonnement intercepté dû à un rayonnement incident qu'elle reçoit formé d'au moins un rayon incident défini par une direction incidente #sk appartenant à une pluralité Ns de directions de l'espace, dites directions incidentes Os, et une énergie incidente Win (#sk), dont une fraction, dite énergie interceptée Wint (ask) est interceptée par la cellule turbide, - pour au moins une portion, dite portion de sol, de la

scène, on calcule, cellule après cellule, une énergie ré-émise Wout (Qyj) par chaque cellule turbide dans chacune des directions de ré-émission Ovj, cette énergie ré-émise Wout (Ovj) étant considérée comme issue d'une énergie volumique diffusée Wdiff (qui) dans la direction de ré-émisssion Qyj par i l'ensemble de la cellule turbide, cette énergie volumique diffusée W diff (Ovj) étant calculée en considérant qu'elle est émise totalement à partir d'un point unique de

la cellule turbide, dit point origine équivalent de diffusion Ms, k (Ovj), e r2 j caractérisé en ce que : - on définit dans un repère de l'espace indépendant des

cellules une pluralité de M secteurs, dits secteurs d'incidence AQM, i, en nombre M inférieur aux directions incidentes Qs, et contenant toutes les Ns directions incidentes Os , - pour calculer, au moins pour la portion de sol, l'énergie volumique diffusée Wdiff (nvj) dans chaque direction de ré-émission Ovj à chaque cellule turbide : . on calcule, pour chaque secteur d'incidence sum, i, une énergie, dite énergie interceptée, de secteur Wmt(##M,i), égale à la somme des énergies interceptées Wint (Qsk, i) des rayons incidents sur la cellule turbide dans des directions incidentes Qsk appartenant au secteur d'incidence ##M,i, . on utilise une fonction, dite fonction de transfert de diffusion de secteur T(##M,i,#vj) permettant de calculer pour chaque secteur d'incidence Ami et chaque direction de ré-émission Qvj, une énergie, dite

énergie diffusée de secteur Wdis (AQM, ivj), diffusée au point origine équivalent

<Desc/Clms Page number 5>

de diffusion Ms, (Qvj) dans la direction de ré-émission Qvj à partir de l'énergie interceptée de secteur Wint (AQM, i) du secteur d'incidence AQM, i, . pour chaque direction de ré-émission Ovj, on calcule l'énergie volumique diffusée W din (Qvj) par la cellule turbide comme la somme des énergies diffusées de secteur Wdi QM i, Qvi) diffusées par chaque secteur d'incidence ilOM, i dans cette direction de ré-émission Qvj.

Avantageusement et selon l'invention, la fonction de transfert de diffusion de secteur T (AQM,Qvj) est une moyenne pondérée de fonctions de transfert de diffusion T (Qsk, Ovj) mémorisées pour chaque direction incidente Qsk et chaque direction de ré-émission Qvj initialement définies. Avantageusement et selon l'invention, la fonction de transfert de diffusion de secteur T (ans, i, Qvj) est définie selon la formule 1 suivante :

et où K (i) est le nombre de directions incidentes #sk dans le secteur ##M,i. i.

Les fonctions sont pré-calculées afin d'éviter le calcul répétitif de celles-ci.

Avantageusement et selon l'invention, le nombre M de secteurs d'incidence ##M,i est supérieur à 2 et inférieur à 20, notamment compris entre 4 et 14, par exemple et de préférence égal à 6, et très inférieur au nombre Ns de directions incidentes #s qui est supérieur à 100, notamment compris entre 100 et 1000. En outre, avantageusement et selon l'invention, les directions incidentes Qs possibles des rayons incidents sont déterminées et choisies dans la même pluralité discrète de directions de l'espace que les directions de ré-émission #v. Cette pluralité de directions est déterminée de façon absolue dans un repère de l'espace, indépendamment des cellules.

Avantageusement et selon l'invention, pour chaque rayon Qgk incident sur la cellule turbide, on détermine un point origine équivalent

<Desc/Clms Page number 6>

Ms, k (Ovj) de la diffusion du rayon incident par la cellule turbide dans la direction de ré-émission Ovj. L'indice"k"indique que la position du point Ms, k (Ovj) dépend de la direction incidente Qsk. Les cellules turbides étant formées de polyèdres pour chaque rayon incident selon une direction Osk et qui entre en un point Pin de l'une de ses faces, dite face d'entrée, le point équivalent origine de la diffusion Ms, k (nj) selon la direction de ré-émission Qyj est déterminé comme le point situé dans la cellule turbide, sur la direction incidente Qsb et à une distance Llrk (Osk) de la face d'entrée par la formule II suivante :

u où k est l'angle formé par la direction incidente ask avec une direction fixe 0 prédéterminée de l'espace, dite direction de référence, O, j est l'angle formé par la direction de ré-émission Qyj avec la direction de référence, u. sk=cos6sk, vj=cosevj, et Alk (Qsk) est le trajet du rayon incident à travers la cellule, gdçlf) G (ask) -. jQgk. fj. df est un coefficient, dit facteur de projection 2n perpendiculaire à la direction Qsk. Cette projection s'applique à la matière comprise dans la cellule turbide sachant que cette matière est orientée selon la direction Qr, avec Or compris dans le secteur 2n, et que la probabilité d'avoir 1,.. r\ giQr) 0 JgiQr) dr\ 1 l'orientation Qf est . On a : p. dQf-1.

27T J 271 : 2n

<Desc/Clms Page number 7>

G(Qyj)=. jQvj. Qfj. dnf est un facteur de projection perpendiculaire à la 2n 2n direction nj.

De la sorte, le point de diffusion Ms, k (nj) dépend des directions de diffusion Qyj et incidentes Qgk, ce qui améliore considérablement la précision des résultats.

Avantageusement et selon l'invention, on calcule et on utilise en tant que point origine équivalent de diffusion, pour chaque rayon

incident, au moins un point origine moyen de diffusion Msk indépendant de la direction de ré-émission et on calcule l'énergie diffusée de secteur Wdiff (, MQvj) à partir de chaque rayon incident de direction incidente Qsk comme diffusée à ce point origine moyen de diffusion Msk. De la sorte, on diminue encore les temps calcul et l'espace mémoire nécessaire.

Avantageusement et selon l'invention, on calcule, et on utilise en tant que point origine équivalent de diffusion, pour chaque rayon incident, deux points origine moyens de diffusion : un point Msk (Qt) associé à une diffusion de rayonnement selon un premier hémisphère prédéterminé de l'espace, dit hémisphère supérieur, et un point Msk(##) associé à la diffusion de rayonnement selon un deuxième hémisphère prédéterminé de l'espace, dit hémisphère inférieur, les deux hémisphères supérieur et inférieur étant choisis de façon que la direction de référence soit orthogonale au plan diamétral séparant ces deux hémisphères. On réalise ainsi un bon compromis entre d'une part une bonne précision de résultats, et d'autre part des temps de calcul et un espace mémoire nécessaire raisonnables. Avantageusement et selon l'invention, la direction de référence est la direction orthogonale au plan diamétral séparant les hémisphères supérieur et inférieur.

Dans le cas où la représentation numérique simulée comprend au moins une image simulée vue de l'extérieur de la scène selon une direction d'observation orientée vers un observateur extérieur, fixe et prédéterminée par rapport à la scène, avantageusement et selon l'invention, la direction d'observation est orientée dans l'hémisphère supérieur.

<Desc/Clms Page number 8>

Avantageusement et selon l'invention, la direction de référence orthogonale au plan diamétral séparant les deux hémisphères supérieur et inférieur est la verticale d'une zone de sol d'une planète, la scène s'étendant à partir de ladite zone de sol de la planète. Avantageusement et selon l'invention, le procédé permet d'obtenir simultanément Nv images vues selon Nv directions d'observations orientées vers l'hémisphère supérieur.

Avantageusement et selon l'invention, on choisit

Msk (nt) =Ms, k (vj (6vj =20 )) et Msk (Q) =Ms, k (vj (Ovj =160 )). Ces choix sont en général de très bons compromis, surtout pour la simulation d'images de télédétection.

Dans la plupart des applications pratiques, avantageusement et selon l'invention, on calcule deux points origine moyens Ms (Qt) et Mg (Q) respectivement associés à la diffusion selon les directions montantes appartenant audit hémisphère supérieur et respectivement les directions descendantes appartenant audit hémisphère inférieur, ces deux points origine moyens Mg (out) et Ms(##) étant utilisés en tant que points origines équivalents de diffusion pour l'ensemble des rayons incidents sur la cellule turbide. De la sorte, le calcul de l'énergie ré-émise Wout(#vj) reste simple même lorsque les cellules reçoivent plusieurs rayons incidents, comme cela est généralement le cas. Pour chaque cellule turbide, les points Ms(##) et Ms(##) sont respectivement calculés en tant que barycentres énergétiques des points origine moyens de diffusion Msk (Ot) et Msk (O), calculés sur les différents rayons incidents sur la cellule turbide, avec les relations itératives :

où ts (out) indique les coordonnées x, y et z du point Ms (out), tsk (Ot) indique les coordonnées x, y et z du point Msl"), et

<Desc/Clms Page number 9>

W, nt est l'énergie interceptée sans compter l'énergie Wmt (Osk) qui vient d'être interceptée. Après cette dernière interception, l'énergie

interceptée est donc : Wint + Wit (Msk)
Ce calcul est effectué de manière séquentielle lors de toute interaction avec une cellule.

Avantageusement et selon l'invention, le procédé est en outre caractérisé en ce qu'avant de calculer les énergies ré-émises Wout (ove) par chaque cellule, chaque point origine moyen Ms(##) et Ms(##) est remplacé par un point d'émission équivalent qui est le point Pe le plus proche choisi parmi une pluralité de points d'émission Pe possibles régulièrement répartis dans chaque cellule, pour lesquels toutes les trajectoires possibles des rayons issus de tous les points Pe selon toutes les directions de ré-émission Qvj possibles ont été préalablement calculées et enregistrées.

Ainsi, avantageusement et selon l'invention, le calcul préalable de toutes les trajectoires possibles des rayons dans la scène est effectué pour une grille de plusieurs points à l'intérieur de la cellule origine, et non pas uniquement pour le centre de la cellule. Ainsi, tout trajet de rayon issu du point Pe d'une cellule quelconque de coordonnées (X, Y, Z) dans la scène, déterminées à partir d'une cellule origine (0,0, 0), est obtenu par translation (X, Y, Z) du trajet d'un rayon parallèle issu du point Pe correspondant dans la cellule origine (0,0, 0).

Lors de toute diffusion, il est choisi les deux points Pe les plus proches des points origine moyens Ms ('n'1') et Ms (Q).

Avantageusement et selon l'invention, le procédé est caractérisé en ce qu'on calcule par itérations successives un rayonnement ré-émis par chaque cellule comprenant au moins un rayon ré-émis dans au moins une direction de ré-émission Qvj, en calculant et mémorisant à chaque itération : - le rayonnement incident sur la cellule issu des rayons ré- émis par les autres cellules adjacentes comme déterminé à l'itération précédente,

- le rayonnement ré-émis par la cellule issu notamment de la diffusion du rayonnement intercepté dû au rayonnement incident sur la cellule.

<Desc/Clms Page number 10>

En outre, avantageusement et selon l'invention, le procédé est caractérisé en ce que l'on définit préalablement dans la scène, indépendamment des cellules, notamment par triangularisation, des éléments géométriques qui représentent différents éléments du paysage, et on calcule l'énergie ré-émise Wout (Qvj) en prenant en compte ces éléments géométriques en tant qu'éléments principalement réfléchissants et/ou absorbants.

L'invention permet ainsi pour la première fois de réaliser des images simulées prenant en compte simultanément et avec une grande précision radiométrique, d'une part, les milieux translucides (par exemple la végétation, les herbes, les feuilles des arbres...) grâce au découpage en cellules turbides, mais aussi les milieux opaques (murs, toits, sols, étendues d'eau, reliefs naturels...) grâce aux éléments géométriques. Les mécanismes de diffusion, absorption et émission sont systématiquement pris en compte à l'aide de modélisations physiques très précises (fonctions de transfert pour les milieux translucides et modèles paramétriques pour les milieux opaques).

Par ailleurs, dans certaines applications avantageuses de l'invention, la scène comprend une portion de sol d'une planète, notamment la Terre, et une portion d'espace s'étendant au-dessus de cette portion de sol.

L'invention permet alors de prendre en compte avec précision la nature tridimensionnelle de l'atmosphère et des surfaces planétaires pour ce qui est de la propagation du rayonnement. Cet aspect est essentiel pour obtenir des images avec une bonne précision radiométrique. En effet, les mécanismes d'interaction "atmosphère-terre"sont extrêmement complexes, en particulier du fait de la nature tridimensionnelle de l'atmosphère et des surfaces terrestres. Ainsi, avantageusement et selon l'invention, la scène comprenant une portion de sol d'une planète, on définit en partie supérieure de la scène, au-dessus de la portion de sol, une colonne atmosphérique découpée en cellules, dites cellules atmosphériques ayant, comme les cellules turbides, des propriétés radiatives connues supposées homogènes dans l'ensemble de chaque cellule atmosphérique. Toutes les cellules atmosphériques de même altitude, c'est-à-dire situées dans un même plan perpendiculaire à la direction de référence, ont les mêmes dimensions et propriétés physiques (densité, coefficients d'extinction, etc.).

<Desc/Clms Page number 11>

Avantageusement et selon l'invention, les dimensions des cellules atmosphériques croissent avec l'altitude. En effet, les gaz et aérosols étant moins denses en altitude, la dimension des cellules peut augmenter sans nuire à la précision radiométrique.

En outre, avantageusement et selon l'invention, le procédé est caractérisé en ce que la colonne atmosphérique comprend une zone inférieure, dite basse atmosphère, au contact de la portion de sol (définissant une zone de paysage jusqu'au sol planétaire), et une zone supérieure, dite atmosphère haute et intermédiaire. De la sorte, on constate que l'on arrive à prendre en compte le comportement radiatif (diffusion, absorption, émission thermique, transmittance) de l'atmosphère avec une bonne précision, et ce malgré la complexité de ses mécanismes tridimensionnels. Ainsi, à partir d'un rayonnement solaire prédéterminé défini par un rayon incident sur la face supérieure de chaque cellule de la couche supérieure de la colonne atmosphérique, on réalise des exécutions successives suivantes :
1) on détermine le rayonnement incident sur la face supérieure des cellules de la basse atmosphère issu du rayonnement solaire après propagation dans l'atmosphère haute et intermédiaire (dans le domaine des grandes longueurs d'onde, le rayonnement incident comprend aussi le rayonnement thermique atmosphérique).

2) on détermine le rayonnement (rayonnement thermique et rayonnement réfléchi dû au rayonnement incident déterminé à l'exécution 1)) montant issu de la basse atmosphère,
3) on détermine le rayonnement rétro diffusé par l'atmosphère haute et intermédiaire vers la basse atmosphère à partir du rayonnement réfléchi déterminé à l'exécution 2),
4) on détermine le rayonnement réfléchi par la basse atmosphère vers l'atmosphère haute et intermédiaire à partir du rayonnement rétrodiffusé déterminé à l'exécution 3),
5) on détermine le rayonnement final montant dans la colonne atmosphérique-notamment au niveau de sa face supérieure-issu de la

<Desc/Clms Page number 12>

somme des rayonnements réfléchis déterminés aux exécutions 2) et 4) après propagation à travers l'atmosphère haute et intermédiaire.

De préférence, l'atmosphère haute et intermédiaire comprend une zone supérieure, dite atmosphère haute, qui dans l'exécution 3) ne diffuse du rayonnement que vers le haut, tout le rayonnement arrivant dans cette atmosphère haute étant diffusé vers le haut. Dans l'exécution 3), seule l'atmosphère intermédiaire (zone de l'atmosphère haute et intermédiaire située sous l'atmosphère haute) est donc susceptible de diffuser du rayonnement vers la basse atmosphère.

L'invention s'étend à un produit logiciel apte à être chargé dans la mémoire vive d'un système informatique pour la mise en oeuvre d'un procédé selon l'invention de façon à réaliser un dispositif informatique selon l'invention.

Ce produit logiciel est un programme d'ordinateur adapté pour être chargé en mémoire vive interne d'un ordinateur. Il peut être proposé sous la forme d'un support d'enregistrement ou mis à disposition pour téléchargement sur un réseau de transmission d'informations (internet).

L'invention s'étend aussi à un support d'enregistrement adapté pour pouvoir être lu par un lecteur relié à un système informatique, caractérisé en ce qu'il comprend un programme enregistré adapté pour pouvoir être chargé en mémoire vive du système informatique et le programmer pour mettre en oeuvre un procédé selon l'invention. Un support d'enregistrement selon l'invention comprend donc un produit logiciel selon l'invention et permet, lorsqu'il est lu sur un lecteur d'un ordinateur, le produit logiciel étant chargé en mémoire vive, l'obtention de représentations numériques simulées-notamment d'images simulées-de l'état énergétique radiatif d'une scène tridimensionnelle hétérogène.

L'invention s'étend à un procédé, un produit logiciel et un support d'enregistrement caractérisés en combinaison par tout ou partie des caractéristiques apparaissant ci-dessus ou ci-après.

D'autres buts, caractéristiques et avantages de l'invention apparaissent de la description suivante qui se réfère aux figures annexées dans lesquelles :

<Desc/Clms Page number 13>

- la figure 1 est un organigramme fonctionnel général d'un procédé selon l'invention, - la figure 2 est un schéma illustrant un exemple de la scène tridimensionnelle d'un procédé selon l'invention,

- la figure 3 est un schéma illustrant un exemple de portion de sol et, partiellement, une portion atmosphérique d'une scène tridimensionnelle hétérogène terrestre découpée en cellules parallélépipédiques dans un procédé selon l'invention, - les figures 4a à 4c illustrent des exemples de secteurs d'incidence définis dans un procédé selon l'invention, - la figure 5 est un schéma illustrant les différentes étapes de prise en compte du couplage atmosphère-sol dans un procédé selon l'invention, - la figure 6 est un schéma de principe illustrant la détermination d'un point d'émission, - les figures 7a et 7b sont des schémas illustrant la détermination des points origine moyens de diffusion d'une cellule turbide dans un procédé selon l'invention, - la figure 8 est un organigramme chronologique général des étapes successives d'un procédé selon l'invention,

- la figure 9 est un organigramme chronologique du calcul de points d'émission des cellules turbides lors de la première itération d'un procédé selon l'invention, - la figure 10 est un organigramme chronologique du calcul d'un point d'émission d'une cellule opaque dans un procédé selon l'invention, - la figure 11 est un organigramme chronologique du calcul de points d'émission des cellules turbides aux itérations subséquentes (autres que la première) du procédé selon l'invention, - la figure 12 est un organigramme chronologique des étapes mises en oeuvre par le module atmosphère d'un procédé selon l'invention.

La figure 1 représente un schéma des différents modules fonctionnels d'un procédé selon l'invention. Sur cette figure, on a représenté en

<Desc/Clms Page number 14>

pointillés les modules de saisie 1, 2 et 3 qui permettent à l'utilisateur de saisir les différents paramètres de configuration selon l'application envisagée.

Le module de saisie 1 permet de définir les caractéristiques géométriques de la scène tri-dimensionnelle que l'on souhaite représenter. Dans l'exemple représenté, la scène comprend une portion de sol 6 d'une planète, notamment la Terre, et une portion d'atmosphère 7 s'étendant à partir des éléments de la portion de sol jusqu'à une altitude prédéterminée.

Le module de saisie 2 permet de définir les caractéristiques optiques (réflectance foliaire pru transmittance foliaire Ta, etc. ) et thermiques (Ti) de la scène tri-dimensionnelle que l'on souhaite représenter.

A partir des paramètres saisis dans les modules de saisie 1 et 2, un module, dit module maquette 4, crée une maquette informatique 5 représentative de la scène, représentée figures 2 et 3. Dans cette maquette 5, la scène est découpée en une pluralité de cellules élémentaires qui dans l'exemple représenté et de préférence, sont cubiques. En variante, ces cellules peuvent être parallélépipédiques.

Dans la maquette 5, il est également défini grâce au module de saisie 1, indépendamment des cellules, des éléments géométriques qui représentent différents éléments de paysage principalement réfléchissants et/ou absorbants-notamment le relief du sol (y compris les surfaces de sol nu, les surfaces aquatiques et les routes, mais à l'exclusion des surfaces végétales), des murs et des toits de bâtiments, et des troncs d'arbres-. Ces éléments géométriques sont définis par triangularisation du paysage, c'est à dire par un maillage représentant le paysage en tant que juxtaposition de polygones (triangles et parallélogrammes) selon les coordonnées des points extrêmes des éléments géométriques saisis et selon la finesse du maillage spatial choisie dans le module de saisie 1.

Le module de saisie 3 permet de définir dans un repère de l'espace indépendant des cellules le mode de discrétisation de différentes directions de l'espace selon lesquelles le rayonnement peut se propager. On définit au moins une pluralité de directions et leur répartition spatiale. Pour chaque pluralité de directions, on définit le nombre N de directions Q, et

<Desc/Clms Page number 15>

éventuellement des densités de distribution spécifiques de ces directions dans certains secteurs angulaires. On peut ainsi définir Ns directions incidentes Os, et Nv directions de ré-émission Qv.

De préférence, les directions incidentes Os possibles des rayons incidents sont déterminées dans la même pluralité discrète de directions de l'espace que les directions de ré-émission Qv (Ng=Nv=N). De surcroît, on définit une direction particulière supplémentaire correspondant à la direction, dite

direction incidente initiale 0', d'une source de lumière principale telle que le s soleil. Il y a donc N+1 directions définies.

En conséquence, toute maquette 5 contient : - des cellules turbides, chaque cellule turbide ayant des propriétés radiatives connues supposées homogènes dans toute la cellule turbide ; en particulier, des cellules turbides, dites cellules foliaires 71, représentent la végétation (herbes, cultures, prairies, couronnes d'arbres,...) et des cellules turbides, dites cellules atmosphériques 72, représentent l'atmosphère ; -des cellules dites cellules opaques 73 contenant au moins une partie d'un élément géométrique tel que décrit ci-dessus, et éventuellement de la matière turbide (végétation ou atmosphère).

Comme on le voit figure 2, la portion d'atmosphère comprend avantageusement plusieurs niveaux d'atmosphère discrets : - un niveau d'atmosphère basse 11, dont les cellules 72 ont la même dimension que les cellules de la portion de sol 6 et qui s'étendent audessus de cellules foliaires 71 ou opaques 73 de la portion de sol 6 ; - un niveau d'atmosphère intermédiaire 12 s'étendant audessus de l'atmosphère basse 11, et dont les cellules 72 ont une dimension plus grande ; - un niveau d'atmosphère haute 13 dont les cellules 73 ont une dimension encore plus grande-notamment occupent chacune toute une couche horizontale de la scène-.

Les cellules turbides 71,72 sont caractérisées par des propriétés optiques et géométriques spécifiques définies grâce aux modules de

<Desc/Clms Page number 16>

saisie 1,2. Ainsi, pour chaque type de végétation, il est défini une orientation angulaire foliaire représentée par un coefficient LAD et une biomasse foliaire représentée par un coefficient LAI. Pour l'atmosphère, il est défini pour les gaz (exposant m) et les aérosols (particules solides ou liquides ; exposant p) des profils verticaux (variant selon l'altitude z) d'albédo de diffusion simple corn (z) et

(OP (z), de coefficient d'extinction (Xm p de coefficient d'extinction par e (Z) et (Xe (Z) l absorption a (z) et c (z), de fonction de transfert de diffusion Tm (Qsk, Ov) et Tp(#sk,#vj), ainsi qu'un profil de température T (z).

Les coefficients d'extinction am, a, a et a sont a OED déterminés à toute altitude z à partir de la connaissance des épaisseurs optiques -C a atmosphériques totales associées à l'absorption gazeuse m, à là diffusion gazeuse Tn', à l'absorption des aérosols 1 et à là diffusion des aérosols Tp. Pour cela, les d d profils de densité volumique des gaz et des aérosols sont supposés décroître exponentiellement avec l'altitude. L'ozone est le seul gaz qui fait exception. Par suite :

où le facteur d'échelle H spécifie le type de profil exponentiel.

Il est utilisé un facteur d'échelle Hm pour les gaz uniformément mélangés dans l'atmosphère, c'est à dire tous les gaz autres que l'ozone et la vapeur d'eau. Il est utilisé un facteur d'échelle HHo pour la vapeur d'eau. Il est utilisé une loi analytique pour représenter le profil de l'ozone, car celui-ci est surtout centré autour de 20-30km d'altitude. Les gaz uniformément mélangés sont essentiellement N2, Oz, CO2, CH4 et N2O. D'autre part, il est utilisé un facteur d'échelle Hp pour représenter le profil de densité volumique des aérosols.

Les épaisseurs optiques Td (X) et la (X) sont calculées à d a toute longueur d'onde X à partir d'une table d'épaisseurs optiques T (i) et T (i) à a des longueurs d'onde Ài pré-définies. Les valeurs Td (i) et Tm (ki) sont calculées d a avec des modèles atmosphériques comme 6-S (Vermote E. F., Tanré D., Deuzé J. L., Herman M. and Morcrette J.-J.,"Second simulation of the satellite signal in the solar spectrum, 6S : An overview", IEEE Trans. Geosci. Remote Sensing, vol.

<Desc/Clms Page number 17>

35, pp. 675-686, May 1997) pour des conditions atmosphériques standards telles qu'une atmosphère tropicale. L'utilisateur peut donner des valeurs particulières, si nécessaire. La forte variabilité spatiale et temporelle de la vapeur d'eau atmosphérique explique que la table stocke l'épaisseur optique Tm (i, dk) de H2O a pour K valeurs d'épaisseur d'eau pouvant être précipitée dk. Par suite, le module maquette 4 calcule Tm (?, d) par interpolation sur la quantité d'eau atmosphérique a pouvant être précipitée d et sur la longueur d'onde 1.

La distribution angulaire de la diffusion des gaz est modélisée avec la fonction de phase de Rayleigh (Lenoble J., 1993, Atmospheric radiative transfer, ISBN 0-937194-21-2, Deepak Publishing, Hampton, Virginia, USA) :

où 8 est le facteur de dépolarisation. Il dépend de l'anisotropie des molécules et est toujours petit. Pour une atmosphère terrestre moyenne ##0.0279. l'angle #sv est l'angle entre les directions Os et n,. Il vérifie

la relation : cossv= [cosOv. coses+sinev. sinOs. cos (0v-0s)] L'épaisseur optique des aérosols est calculée à toute longueur d'onde À par la loi empirique d'Angstrom : () To)- ['] . Le e e coefficient d'Angstrom ss dépend du type d'aérosols présents. Il est en général compris entre 0. 5 et 2. D'autre part, l'on a : T) = (Dp. T () et T) = (l-û) p). T), d e a e où cop est l'albédo de diffusion simple. Le calcul du taux de diffusion et d'absorption par les aérosols atmosphériques à toute altitude z nécessite donc la connaissance de 4 paramètres : l'épaisseur optique T (o) à une longueur d'onde e o pré-définie, le coefficient ss, l'albédo de diffusion simple cop et le facteur d'échelle Hp.

<Desc/Clms Page number 18>

La distribution angulaire de la diffusion des aérosols est modélisée avec la somme de deux fonctions PP (IF) de Henyey-Greenstein (Lenoble, 1993) :

où le facteur d'asymétrie g représente le surplus de diffusion avant par rapport à la diffusion arrière. Il est donc nul si PP (IF) est symétrique. Il dépend du type d'aérosols et de la longueur d'onde.

Un module géométrique 9 calcule et mémorise les directions incidentes Os et les directions de ré-émission Ov en fonction des valeurs saisies dans le module de saisie 3. En outre, le module géométrique 9 calcule et mémorise, dans le repère de l'espace indépendant des cellules, une pluralité de M secteurs volumiques, dits secteurs d'incidence AM, en nombre M inférieur aux Ns directions incidentes Os mais contenant toutes les Ns directions incidentes Os, complémentaires les uns des autres dans l'angle solide plein et formant une partition de l'angle solide plein.

De préférence, le nombre Ns de directions incidentes (le nombre N de directions dans le cas où Ng=Nv=N) est supérieur à 100, notamment compris entre 100 et 1000.

De préférence, le nombre M de secteurs d'incidence, saisi grâce au module de saisie 3, est supérieur à 2 et inférieur à 20. La figure 4a représente un exemple avec M=6 ; la figure 4b représente un exemple avec M=10 ; la figure 4c représente un exemple avec M= 14. La valeur M=6 représente un bon compromis pour les scènes terrestres.

Un module, dit module optique 8, calcule et mémorise pour chaque type de cellule turbide quatre fonctions de transfert à partir des propriétés optiques et géométriques des éléments contenus dans la cellule saisies par les modules de saisie 1,2 : (1) une fonction de transfert de transmittance T (d2), par unité de déplacement et par unité de densité volumique de surface de la matière contenue dans la cellule. Cette fonction permet de calculer l'énergie transmise

<Desc/Clms Page number 19>

Wtr (Q) à travers la cellule turbide, et donc aussi l'énergie interceptée Wjnt (Q) par la cellule turbide, selon chacune des Ns+ 1 directions incidentes Os pré-définies, (cf. publication"Modelling Radiative Transfer in Heterogeneous 3-D Végétation Canopies"J-P. Gastellu-Etchegorry et al Remote Sens. Environ. 58 : 131-156 (1996)) ; (2) une fonction de transfert de diffusion initiale

T (QQyj) mémorisée pour la direction incidente initiale Qo et chaque direction s s de ré-émission Ovj initialement définies ; cette fonction indique l'énergie que tout élément volumique de la cellule diffuse à la première itération du procédé selon chacune des Nv directions de ré-émission (Qtyj), compte tenu de l'énergie

Wint (Qsk) interceptée selon la direction incidente initiale (no) ; cette fonction de s transfert de diffusion initiale T (Q, Qvj) est préalablement calculée par le module s optique 8 ; (3) une fonction de transfert de diffusion de secteur T (##M, i,#vj) de calculer pour chaque secteur d'incidence AQM, i et chaque direction de ré-émission Qvj, une énergie, dite énergie diffusée de secteur Wdiff(##M,1,#vj), diffusée au point origine équivalent de diffusion Ms,k(#vj) dans la direction de ré-émission Qvj à partir de l'énergie interceptée de secteur Wuit (AQM, j) du secteur d'incidence AQM, i. Cette fonction est calculée selon la formule 1 ci-dessus. Les fonctions de transfert de diffusion T(#sk,#vj) sont préalablement calculées par le module optique 8 ; (4) une fonction de transfert thermique T (uf, Tf, Sf, vj) qui permet de calculer l'énergie thermique ré-émise par la cellule dans la direction de ré-émission, Qvj, compte tenu de la quantité de matière Uf présente ainsi que de la température Tf et de l'émissivité #f de cette matière (Thèse Modélisation des bilans radiatifs et énergétique des couverts végétaux , P. Guillevic, Université Paul Sabatier, 13/12/1999. Paragraphe III.4.3.2 p47 à 49).

Il est à noter que dans un procédé selon l'invention, on utilise une fonction de transfert de transmittance T (fl) qui est la même pour chaque type de cellule et chaque direction Q, et qui ne dépend pas en particulier

<Desc/Clms Page number 20>

du trajet du rayon dans la cellule ni de la quantité de matière présente dans la cellule.

Par exemple, pour des cellules de végétation foliaire, la fonction de transmittance définie par la loi de Beer-Lambert s'écrit :

où Ali est la distance du trajet du rayon dans la cellule, uf (i) est la quantité de surface foliaire par unité de volume G(j, #) est un facteur de projection des feuilles perpendiculairement à la direction Q considérée, définissant un type de végétation.

Dans le procédé de l'invention, on utilise une approximation de cette fonction, et on ne mémorise que N valeurs pour chaque type j de cellule turbide, soit T =exp [-G (j, Q)]. Pour des cellules de végétation

foliaire de densité volumique foliaire uni), on a donc : T (ill, Q) = [y] ''.

Tous ces calculs étant effectués, un module, dit module radiatif 10, simule la propagation du rayonnement dans la scène par itérations successives, et, à chaque itération, par suivi de rayon à travers les cellules de la scène, comme expliqué ci-après, et fournit une représentation numérique simulée de l'état radiatif de la scène. Cette représentation peut comprendre des images 11 de télédétection, pour toutes les directions de ré-émission #vj montantes et/ou des matrices tri-dimensionnelles 12 représentatives du bilan radiatif de la scène.

La figure 8 est un organigramme général des différentes étapes mises en oeuvre dans un procédé d'obtention selon l'invention, par le module radiatif 10.

Comme cela est illustré schématiquement figure 6 dans l'étape initiale 14, le module radiatif 10 prédéfinit dans une cellule origine de coordonnées (0,0, 0) (figure 3) une pluralité de points d'émission possibles Pe, régulièrement répartis au sein de la cellule. Dans l'exemple représenté figure 6, la cellule est vue en élévation et comprend 53=125 points d'émission possibles Pe.

En outre, le module radiatif 10 calcule au préalable et enregistre toutes les trajectoires possibles à travers la scène des rayons pouvant être émis à partir de chaque point d'émission possible Pe selon les différentes directions de ré- émission Qy par suivi de rayon.

<Desc/Clms Page number 21>

Ce calcul 14 initial fournit donc n3. N liste chaînées qui sont mémorisées, où n3 est le nombre de points Pe et N est le nombre total de directions Q autres que la direction initiale Q. Une liste chaînée correspond s donc à un trajet possible selon une direction donnée et à partir d'un point Pe donné. Tout élément de liste chaînée est une structure qui contient les coordonnées du point d'entrée et la longueur du trajet dans la cellule en cours ainsi que l'indice qui indique la cellule suivante.

On détermine 15 tout d'abord l'éclairement solaire au sommet de l'atmosphère 7. Cet éclairement solaire est défini par la direction

s Ç20) comme incidente initiale no et par une énergie incidente initiale Win (Q) comme représenté figure 2.

Une première étape 16 d'un module, dit module

atmosphère, simule la propagation du rayonnement solaire W, (Qs) à travers s l'atmosphère et détermine 17 l'éclairement solaire direct incident sur la portion de sol 6 (c'est-à-dire l'énergie résultant de la transmission du rayonnement solaire à travers l'atmosphère et incidente sur la portion de sol 6). Lors de cette étape 16 du module atmosphère, il est également déterminé 18 l'éclairement diffus du paysage, c'est-à-dire le rayonnement résultant de la diffusion 19 après interception des rayonnements incidents sur les différentes cellules atmosphériques, augmenté de l'émission thermique 20 des cellules atmosphériques calculé précédemment. L'éclairement solaire direct de la portion de sol 6 est déterminé par simulation de la partie du rayonnement solaire initial non intercepté par les cellules atmosphériques.

Lors de l'étape 22 subséquente, on détermine, cellule après cellule, l'énergie ré-émise Wout (Qvj) par diffusion par chaque cellule de la portion de sol 6 dans chacune des directions de ré-émission Qyj. Pour chaque cellule, on calcule tout d'abord au moins un point d'émission Pe. On calcule d'une part 23 deux points d'émission Pend, Pet pour les cellules turbides, et d'autre part, et simultanément, lors de l'étape 24, un point d'émission Pe pour chaque élément surfacique de toute cellule opaque 73.

<Desc/Clms Page number 22>

La figure 9 représente les différentes étapes du calcul des points Pe, Pet. A partir du rayonnement solaire direct issu de l'étape 17, on calcule 25 pour la première cellule turbide 71, 72 rencontrée par chaque rayon incident selon la direction solaire ÇÊS, une énergie interceptée de secteur Wint (Q) égale à l'énergie interceptée par la cellule turbide, et ce par application de la fonction de transfert de transmittance T : W) =W (Q ;) [l-[T ()] -'], s s s où Uf (i) est par exemple la densité volumique foliaire d'une cellule de végétation et où AI, est la longueur de la trajectoire du rayon de direction Q dans la cellule s i.

On calcule ensuite 26 deux points origine moyens de diffusion M (nit) et M (Q) associés respectivement à la diffusion du s s rayonnement selon un premier hémisphère supérieur de l'espace, et, selon un deuxième hémisphère inférieur. Pour ce faire, on applique la formule (II) cidessus et on choisit par exemple et de préférence Ms0(##)=Ms,k(#vj(#v=20 )) et

M (Q) =Ms, k (Qvj (6v= 1600)). Lors de cette première itération sk=. s v k=
La cellule turbide 71,72 pouvant recevoir plusieurs rayons incidents de direction #s0, on calcule deux points origine moyens Ms (out) et Ms(##) en tant que barycentres énergétiques des points Msk (d) et respectivement Msk(##), et ce lors de l'étape 27 subséquente. Ce calcul est effectué avec les relations itératives susmentionnées.

La figure 7a représente schématiquement la détermination du point origine moyen de diffusion Msk (out) et la figure 7b représente un exemple de calcul d'un point origine moyen Ms(##), avec trois rayons incidents d'énergie incidente Wsl), Wint (Qs2), et Wint (s3) respectivement dans le cas général de trois directions incidentes distinctes #s1, #q2, #s3. Dans la première itération, le même calcul est effectué mais toutes les directions incidentes sont

confondues avec no. De façon plus générale, et lors des itérations subséquentes, s pour chaque rayon incident, on détermine les trois points origine moyens de diffusion Ms1(##), Ms2(##) et Ms3(##). On fait le barycentre de ces trois points pour trouver le point origine moyen Ms(##),. Il est à noter que les figures 7a, 7b,

<Desc/Clms Page number 23>

sont tracées dans le plan, la cellule étant vue en élévation, mais doivent être considérées dans l'espace, les cellules étant des volumes.

On réitère ces différentes étapes 25 à 27 successivement pour toutes les cellules turbides 71,72 de la portion de sol 6 comme schématisé par la flèche 28 de la figure 5.

On calcule ensuite 29 l'énergie diffusée Wdiff(#s0,#vj) diffusée au point origine moyen Ms (out) et Ms (), selon que la direction de ré- émission considérée Qvj est montante ou descendante, à partir de l'énergie interceptée de secteur Wint (Q) calculée à l'étape 25 précédente et par application

de la fonction de transfert de diffusion initiale T (Q, Qvj) précalculée, soit : 0 0 0 Wdiff ( j) = W). T (Qvj). s s
On détermine ensuite 30 quels sont les points d'émission possibles Pet et Pe les plus proches des points origine moyens Ms (out) et Ms () comme représenté schématiquement figure 6. En effet, on a préalablement calculé, pour chaque cellule, une pluralité de points d'émission Pe possibles, régulièrement répartis dans la cellule. Egalement, le module radiatif 10 a calculé préalablement et enregistré les trajectoires possibles des différents rayons issus selon les directions incidentes Qgk possibles, et les directions de ré- émission Ovj possibles, à partir des différents points d'émission Pe possibles de la cellule origine (0,0, 0) (calcul initial 14). Dans une cellule quelconque de coordonnées (X, Y, Z), tout trajet de rayon issu d'un point d'émission Pe possible est déterminé à partir des listes chaînées mémorisées pour la cellule origine (0,0, 0) par translation (X, Y, Z) du trajet (liste chaînée) d'un rayon parallèle issu du point Pe correspondant dans la cellule origine (0,0, 0). Ainsi, on calcule lors de l'étape 31 l'énergie diffusée par la cellule turbide dans chacune des directions de ré-émission Qvj, c'est-à-dire l'énergie Wout (Ovj) en supposant que cette énergie est ré-émise à partir du point Pe qui est soit le point Pet, soit le point Pe selon que la direction Qvj est montante ou descendante.

La figure 10 représente le calcul 24 du point d'émission Pe dans une cellule opaque à partir de l'éclairement solaire direct calculé à l'étape 17, dans la première itération. On détermine 33 si le rayon incident est ou non

<Desc/Clms Page number 24>

intercepté par un élément géométrique de la cellule. Si tel est le cas, on détermine à partir des propriétés optiques connues de l'élément géométrique, l'énergie interceptée par cet élément géométrique et les coordonnées du point d'intersection correspondant lors de l'étape 34. On calcule ensuite 35 le barycentre énergétique des différents points d'intersection correspondant aux différents rayons incidents sur l'élément de surface dans la cellule, ainsi que l'énergie totale interceptée correspondante. On réitère ces étapes 33 à 35 pour chaque élément de surface de chaque cellule opaque et sur toutes les cellules opaques, comme représenté schématiquement par la flèche 36. Dans le cas où aucun élément géométrique n'intercepte le rayon incident, comme déterminé lors de l'étape 33, on passe immédiatement à l'élément géométrique subséquent par itération comme représenté par la boucle 37 figure 10.

Une fois ces itérations effectuées, pour chaque cellule qui contient au moins un élément géométrique (surface), on calcule 38 pour chaque élément géométrique dont l'énergie d'interception est non nulle le point d'émission Pe possible qui est le plus proche du barycentre des points d'interception et qui est extérieur à tout volume circonscrit par des éléments géométriques.

Le module de saisie 2 a permis de définir des lois paramétriques permettant de déterminer les diffusions issues des éléments géométriques à partir des énergies interceptées par ces éléments géométriques dans les différentes directions. Il suffit donc d'appliquer ces lois paramétriques pour déterminer l'énergie diffusée Wdiff (Qv) au point d'émission Pe par chaque élément géométrique selon chaque direction de ré-émission Ovj. On calcule ensuite 31 l'énergie ré-émise Wout (Qvj) correspondant à cette énergie diffusée par

application de la fonction de transfert de transmittance Ty), si la cellule contient de la matière turbide.

Par ailleurs, à partir de l'émission thermique atmosphérique 20 permettant d'obtenir l'éclairement diffus du paysage 18 et à partir de l'émission thermique de la portion de sol qui peut aussi être calculée 39, on obtient à l'issue de la première itération Il du procédé, le rayonnement ré-émis par chacune des cellules de la maquette 5 à partir du rayonnement solaire initial.

<Desc/Clms Page number 25>

Dans cette première itération Il, on ne prend donc pas en compte les phénomènes de diffusion multiples. En effet, contrairement aux diffusions multiples, la diffusion d'ordre 1 n'est influencée que par un rayonnement Wm (vos) monodirectionnel selon la direction (). De plus, s s l'importance des diffusions multiples est beaucoup plus faible que celle des diffusions d'ordre 1.

Le procédé se poursuit ensuite par itérations successives en calculant et mémorisant à chaque itération In : - lors d'une étape d'initialisation 40 de l'itération, pour chaque cellule, le rayonnement intercepté Wint (Qsk) issu de l'interception du rayonnement incident Win (Qsk), sachant que ce dernier est le rayonnement ré- émis Wout(#vj) par les autres cellules adjacentes comme déterminé à l'itération précédente,

- le rayonnement ré-émis Wout (vj) par chaque cellule issu notamment de la diffusion du rayonnement intercepté sur la cellule.

L'énergie interceptée Wint (ami) par chaque secteur est initialisée à zéro à chaque étape 40 d'initialisation d'itération.

Là encore, on calcule pour chaque cellule au moins un point d'émission du rayonnement ré-émis, lors de chaque itération In. De même que pour la première itération, on distingue un calcul 41 pour les cellules turbides et un calcul 42 pour les cellules opaques. Le calcul 41 des points d'émission Pet et Peul pour les cellules turbides est expliqué plus en détail figure 11. On choisit tout d'abord lors de l'étape 43, un rayon incident sur la cellule turbide selon la direction incidente Qsk et avec une énergie incidente Win(#sk). Lors de l'étape 44 subséquente, on détermine à quel secteur Anm, i appartient la direction incidente Qsk. On calcule ensuite 45, l'énergie interceptée Wint (Qsk) du rayon incident par la

cellule turbide par application de la fonction de transfert de transmittance T (flk) soit : Wint (Qsk) =Win (Qsk) OÙ Udi) est par exemple la densité volumique foliaire d'une cellule de végétation et où AI, est la longueur de la trajectoire du rayon de direction Qsk dans la cellule i.

<Desc/Clms Page number 26>

On incrémente ensuite 46 l'énergie interceptée Wint (AQM, i) pour chaque secteur AM. i en y ajoutant l'énergie interceptée Wffit (nk) par chaque rayon incident.

On calcule ensuite 47 deux points origine moyens de diffusion Msk (out) et Msk(##), associés respectivement à la diffusion du rayonnement selon un premier hémisphère supérieur de l'espace, et, selon un deuxième hémisphère inférieur. Pour ce faire, on applique la formule (II) ci-

dessus et on choisit Msk (qat) =Ms, k (Qvj (6v=20 )) et Msk () =Ms, k (Qvj (Oy=160 )).

Lorsque la cellule turbide reçoit plusieurs rayons incidents, on calcule deux points origine moyens Ms (nt) et Ms Ms(##), en tant que barycentres énergétiques des points Msk (Ot) et respectivement Msk (O), et ce lors de l'étape 48 subséquente. Ce calcul est effectué avec les relations itératives susmentionnées.

On fait ensuite une itération 49 sur les différents rayons incidents Msk sur la cellule turbide.

On calcule ensuite, lors de l'étape 50, les énergies diffusées de secteur Wdifî (AQM, i,vj) diffusées par chacun des M secteurs AQM, i dans les directions de ré-émission Qvj par les points origine moyens Ms (Qt) et Ms Ms(##), déterminés à l'étape 48, et ce à l'aide de la fonction de transfert de diffusion de secteur calculée selon la formule (I) ci-dessus. En appliquant cette fonction de transfert de diffusion de secteur à l'énergie interceptée de secteur

Wint (AM, i), on obtient, pour chaque direction de ré-émission Qvj, l'énergie diffusée de secteur Wdiff (AM, i, ij soit :

On effectue ensuite la somme de ces énergies diffusées de secteur pour déterminer l'énergie volumique diffusée W (Anj) par la cellule turbide dans chaque direction de ré-émission Ovj, c'est-à-dire qu'on effectue la somme des différentes énergies de secteur diffusées par les différents secteurs d'incidences AQM, i.

Lors de l'étape 51, on détermine les deux points d'émission

Pet et Pe les plus proches des points origine moyens Mg (out) et Mg (Q). En

<Desc/Clms Page number 27>

considérant que l'énergie volumique diffusée Wdiff (Anj) est émise à partir du point d'émission Pet ou Pe prédéterminé, et selon l'énergie interceptée par la matière turbide entre ce point d'émission et le point Pout (Qyj) de sortie du rayonnement correspondant de la cellule turbide, on détermine avec la fonction de transfert de transmittance T l'énergie ré-émise Wout (Ovj) On réitère le calcul 50 des énergies diffusées et ré-émises pour toutes les directions de ré- émission Qyj.

On détermine ainsi lors de l'étape 52 l'énergie ré-émise Wout (Qvj) par la cellule turbide dans la direction de ré-émission Qyj à partir des deux points d'émission Pet et Pe selon que la direction de ré-émission Ovj est montante ou descendante. On réitère ce calcul sur toutes les cellules turbides de la portion de sol 6 de la maquette 5.

Pour une cellule opaque 73, le calcul 42 du point d'émission Pe est similaire au calcul 24 décrit en relation à la figure 20 pour la première itération Il. Après ce calcul, l'énergie ré-émise Wout (Ovj) par la cellule opaque 73 est également calculée lors de l'étape 52 avec une expression paramétrique comme décrit précédemment dans l'étape 31 de la première itération Il. La différence est que cette fois-ci, l'expression paramétrique est adaptée au cas d'un rayonnement incident plus ou moins isotrope et non plus monodirectionnel selon la direction . L'ensemble du calcul est réitéré pour toutes les cellules de la portion de sol 6 et cette itération est représentée par la boucle 53 figure 4. Pour passer d'une itération à la suivante, on considère dans l'étape d'initialisation 40 que l'énergie ré-émise lors d'une itération In correspond à l'énergie incidente de l'itération subséquente Ion+1. Après un nombre n prédéterminé d'itérations, dépendant de la précision souhaitée, on connaît l'état radiatif des cellules de la portion de sol 6.

On détermine 54 ensuite si une rétrodiffusion atmosphérique doit ou non être prise en compte. Dans l'affirmative, le module atmosphère exécute une étape 55 représentée plus en détail figure 12. En effet, ce calcul de rétrodiffusion atmosphérique 55 peut ne pas être effectué à toutes les exécutions du calcul du module radiatif 10. Comme représenté schématiquement figure 5, la diffusion atmosphérique est prise en compte lors d'une première

<Desc/Clms Page number 28>

exécution El. Lors d'une deuxième exécution E2 du calcul radiatif, le module radiatif 10 calcule le rayonnement réfléchi par la portion de sol 6 vers la portion d'atmosphère 7 (correspondant à l'atmosphère haute 13 et intermédiaire 12) à partir du rayonnement incident sur la face supérieure des cellules de l'atmosphère basse 11 comme déterminé lors de la première exécution El. Lors de cette deuxième exécution E2, on n'exécute donc pas l'étape 55 de calcul de rétrodiffusion atmosphérique. Par contre, lors d'une troisième exécution E3 subséquente, on détermine le rayonnement rétrodiffusé par l'atmosphère haute 13 et intermédiaire 12 (étape 55) vers la basse atmosphère 11 (et vers la portion de sol 6) à partir du rayonnement réfléchi déterminé lors de la deuxième exécution E2. Lors d'une quatrième exécution E4, le module radiatif 10 détermine le rayonnement réfléchi par la portion de sol 6 et la basse atmosphère Il vers l'atmosphère haute 13 et intermédiaire 12, à partir du rayonnement rétrodiffusé tel que déterminé lors de la troisième exécution E3. Lors de la quatrième exécution E4, on n'exécute donc pas le calcul 55 de rétrodiffusion atmosphérique.

Lors d'une cinquième exécution E5 subséquente, on détermine le rayonnement final montant dans la colonne atmosphérique - notamment au niveau de la face supérieure de la scène, c'est-à-dire de la portion atmosphère 7-après propagation à travers l'atmosphère. Cette séquence d'exécutions du procédé selon l'invention est représentée schématiquement figure 5.

Lors de chaque calcul 55 de rétrodiffusion, comme représenté figure 12, on détermine tout d'abord 56 la réflexion du rayonnement par la portion de sol 6 et la basse atmosphère 11 selon les directions de ré- émission montantes, c'est-à-dire appartenant à l'hémisphère supérieur.

On calcule ensuite 57 la rétrodiffusion atmosphérique en tant que convolution de l'énergie montante issue de la basse atmosphère 11 avec une fonction de transfert notée FTBA-BA. Cette dernière indique l'énergie rétrodiffusée WBOA(i',j',#s#) par l'atmosphère sur chaque pixel (i', j') du haut de la basse atmosphère 11 et pour toute direction descendante Qy, sachant l'on a une

énergie unité montante WBOA (i,vt) selon une direction montante Ovt, à partir

<Desc/Clms Page number 29>

d'un pixel (ij). Par suite, la fonction de transfert est notée FTBA-BAALvls où Ai=i'-i et Aj=j-j'. L'expression de la convolution est :

La fonction de transfert est précalculée par le module optique 8 en simulant l'énergie WaoAC'j's) rétrodiffusée par l'atmosphère, en tout pixel (i', j') et pour toute direction descendante Ost sachant que la seule source d'énergie est un pixel (ij) du haut de la basse atmosphère 11 qui envoie

une énergie WBOA (ij.. vt) selon une direction montante Qyi. Ce calcul préliminaire est effectué de manière itérative sur toutes les directions montantes n, t. L'emploi de la fonction FTBA-BA diminue beaucoup les temps calcul, car il évite de simuler la rétrodiffusion atmosphérique pour chaque énergie WBOA(i,j,#v#) issue de tout pixel du haut de la basse atmosphère 11, selon chacune des directions montantes out. Le temps calcul est donc sensiblement diminué d'un facteur égal au nombre de pixels de la basse atmosphère 11, soit un terme généralement supérieur à 500, voire beaucoup plus.

On calcule ensuite 58 la réflexion par la portion de sol 6 et la basse atmosphère 11 de la rétrodiffusion atmosphérique précédemment calculée lors de l'étape 57. Lors du calcul 58 de réflexion, seules les directions de ré-émission montantes sont prises en compte. On calcule enfin lors de l'étape 59 le rayonnement transmis à tous les niveaux et notamment au sommet de l'atmosphère, selon les directions de ré-émission montantes, c'est-à-dire en direction du capteur 60 (observateur), qui peut être, par exemple, un capteur aéroporté ou satellitaire. L'énergie transmise depuis le sommet de la basse atmosphère 11 jusqu'au capteur 60 est calculée à l'aide d'une fonction de transfert

notée FTBA-capteur (Ai, Aj, Qvi-, Qst). Cette fonction indique l'énergie montante WBoA'J'st) en tout pixel (i', j') à l'altitude du capteur 60, selon toute direction montante Ost, sachant que la seule source d'énergie est due à un pixel (ij) de la basse atmosphère 11 qui envoie une énergie WBOA(i,j,#v#) selon la direction montante Qvf. Par suite, l'énergie reçue au niveau du capteur 60 est

<Desc/Clms Page number 30>

où WBoA (ij, Qvt) est l'énergie totale montante issue de tout pixel (ij) de la basse atmosphère 11. Cette énergie est donc la somme de deux termes qui représentent les énergies que la portion de sol 6 diffuse vers le haut lors des étapes 56 et 58.

L'invention peut faire l'objet de nombreuses variantes par rapport au mode de réalisation décrit ci-dessus. Il est possible de réaliser non seulement des images de télédétection, mais également un bilan énergétique tridimensionnel radiatif de la scène. Les images de télédétection peuvent être de type satellitaire ou autres. L'invention permet aussi d'obtenir des représentations numériques simulées d'autres scènes que des scènes terrestres.

Les différentes fonctions mentionnées ci-dessus peuvent être aisément mises en oeuvre par programmation informatique, c'est-à-dire par programmes d'ordinateurs. En particulier, le procédé selon l'invention a été mis en oeuvre et programmé en langage C sur station de travail IBM Risc 6000 (R) sous environnement Linux. Il a été possible de simuler des images représentant des paysages naturels et urbains avec un temps de calcul de l'ordre de 4h pour des paysages complexes. Ces temps de calcul peuvent être inférieurs à 15min, et même beaucoup moins, pour des paysages relativement simples.

Claims (20)

REVENDICATIONS

1/-Procédé d'obtention d'une représentation numérique simulée de l'état énergétique radiatif d'une scène (5) tridimensionnelle hétérogène, dans lequel : - on découpe la scène (5) en une pluralité de cellules (71, 72,73) élémentaires comprenant des cellules, dites cellules turbides (71,72), ayant des propriétés radiatives connues supposées homogènes dans l'ensemble de chaque cellule, - on mémorise pour chaque cellule turbide (71, 72) des fonctions de transfert permettant de déterminer, dans chaque direction d'une pluralité discrète de Nv directions de l'espace, dites directions de ré-émission Qy, prédéterminées dans un repère de l'espace indépendant des cellules, une énergie

1 1 ré-émise par la cellule turbide en fonction d'un rayonnement intercepté dû à un rayonnement incident qu'elle reçoit formé d'au moins un rayon incident défini par une direction incidente #sk appartenant à une pluralité discrète de Ns directions de l'espace, dites directions incidentes Qs, et une énergie incidente Wm (Qsk), dont une fraction, dite énergie interceptée Wint(#sk), est interceptée par la cellule turbide, - pour au moins une portion, dite portion de sol (6), de la

scène, on calcule, cellule après cellule, une énergie ré-émise Wout (Ovj) par chaque cellule turbide (71, 72) dans chacune des directions de ré-émission Qvj, cette énergie ré-émise Wout(#vj) étant considérée comme issue d'une énergie volumique diffusée Wdiff(Qvj) dans la direction de ré-émisssion Qvj par l'ensemble de la cellule turbide, cette énergie volumique diffusée Wdiff (vj) étant calculée en considérant qu'elle est émise totalement à partir d'un point unique de la cellule turbide, dit point origine équivalent de diffusion Ms, k (vj), caractérisé en ce que : - on définit dans un repère de l'espace indépendant des cellules (71,72) une pluralité de M secteurs, dits secteurs d'incidence ##M,i, en

<Desc/Clms Page number 32>

énergie diffusée de secteur Wff (AQM, i,vj), diffusée au point origine équivalent de diffusion Ms, k (Ovj) dans la direction de ré-émission Qvj à partir de l'énergie interceptée de secteur Wint (arm, i) du secteur d'incidence QM, i, . pour chaque direction de ré-émission Qvj, on calcule l'énergie volumique diffusée WdQyj) par la cellule turbide comme la somme des énergies diffusées de secteur Wdiff (AQMj, Qyj) diffusées par chaque secteur d'incidence AQM, i dans cette direction de ré-émission Qvj.

. on calcule, pour chaque secteur dtincidence SQM i, une énergie, dite énergie interceptée de secteur Wint (AQM, i), égale à la somme des énergies interceptées Wmt (Osk, D des rayons incidents sur la cellule turbide dans des directions incidentes Osk appartenant au secteur d'incidence ##M,i, . on utilise une fonction, dite fonction de transfert de diffusion de secteur T (AQM, i, Qyj) permettant de calculer pour chaque secteur d'incidence ##M,i et chaque direction de ré-émission #vj, une énergie, dite

nombre inférieur aux directions incidentes Qs et contenant toutes les Ns directions incidentes Qs, - pour calculer, au moins pour la portion de sol (6), l'énergie volumique diffusée Wdiff (Qvj) dans chaque direction de ré-émission Qvj par chaque cellule turbide :

2/-Procédé selon la revendication 1, caractérisé en ce

que la fonction de transfert de diffusion de secteur T (ans, i,yj) est une moyenne pondérée de fonctions de transfert de diffusion T (Qsk, Qvj) mémorisées pour chaque direction incidente Qsk, et chaque direction de ré-émission Qvj.

3/-Procédé selon l'une des revendications 1 ou 2, caractérisé en ce que le nombre M de secteurs d'incidence est supérieur à 2 et inférieur à 20.

<Desc/Clms Page number 33>

4/-Procédé selon l'une des revendications 1 à 3, caractérisé en ce que le nombre Ns de directions incidentes est supérieur à 100, notamment compris entre 100 et 1000.

5/-Procédé selon l'une des revendications 1 à 4, caractérisé en ce que les directions incidentes Qs possibles des rayons incidents sont déterminées dans la même pluralité discrète de directions de l'espace que les directions de ré-émission Qv.

6/-Procédé selon l'une des revendications 1 à 5, caractérisé en ce que les cellules turbides (71,72) étant formées de polyèdres, pour chaque rayon incident selon une direction Qsk et qui entre en un point Pin de l'une de ses faces, dite face d'entrée, on détermine le point origine équivalent Ms, k (Qvj) de la diffusion du rayon incident par la cellule turbide (71,72) dans la direction de ré-émission Ovj, comme le point situé dans la cellule turbide, sur la direction incidente Qgk. et à une distance Ark (ask) de la face d'entrée déterminée par la formule II suivante :

) D. sk=cos6sk, fvj=cosOvj, et Alk (Qsk) est le trajet du rayon incident à travers la cellule, G (ask) est un facteur de projection de la matière perpendiculairement à la direction Qsk prédéterminé pour chaque type de cellule, G (Qyj) est un facteur de projection de la matière perpendiculairement à la direction Qyj prédéterminé pour chaque type de cellule.

OÙ Ogk est l'angle formé par la direction incidente Qsk avec une direction fixe prédéterminée de l'espace, dite direction de référence, Ovj est l'angle formé par la direction de ré-émission Ovj avec la direction de référence,

<Desc/Clms Page number 34>

7/-Procédé selon la revendication 6, caractérisé en ce qu'on calcule et on utilise en tant que point origine équivalent de diffusion, pour chaque rayon incident, au moins un point origine moyen de diffusion Msk indépendant de la direction de ré-émission Qvj, et on calcule l'énergie diffusée de

secteur Wdi e (AQM, ivj) à partir de chaque rayon incident de direction incidente Qsk comme diffusée à ce point origine moyen de diffusion Msk.

8/-Procédé selon la revendication 7, caractérisé en ce qu'on calcule et on utilise en tant que point origine équivalent de diffusion, pour chaque rayon incident, deux points origine moyens de diffusion : un point Msk (d) associé à une diffusion de rayonnement selon un premier hémisphère prédéterminé de l'espace, dit hémisphère supérieur, et un point Msk (Q) associé à la diffusion de rayonnement selon un deuxième hémisphère prédéterminé de l'espace, dit hémisphère inférieur, les deux hémisphères supérieur et inférieur étant choisis de façon que la direction de référence soit orthogonale au plan diamétral séparant ces deux hémisphères.

9/-Procédé selon la revendication 8, dans lequel la représentation numérique simulée comprend au moins une image simulée vue de l'extérieur de la scène (5) selon une direction d'observation, fixe et prédéterminée par rapport à la scène, caractérisé en ce que la direction d'observation est orientée vers l'hémisphère supérieur.

10/-Procédé selon l'une des revendications 8 ou 9, caractérisé en ce que la direction de référence orthogonale au plan diamétral séparant les deux hémisphères supérieur et inférieur est la verticale d'une zone de sol d'une planète, la scène (5) s'étendant à partir de ladite zone de sol de la planète.

11/-Procédé selon la revendication 6 et l'une des revendications 8 à 10, caractérisé en ce qu'on choisit

M) =Ms, k (Qvj (Ovj=20 )) et M) =Mvj (el60 )).

12/-Procédé selon l'une des revendications 8 à 11, caractérisé en ce qu'on calcule deux points origine moyens Mg (Qt) et Ms (Q) respectivement associés à la diffusion selon les directions montantes appartenant

<Desc/Clms Page number 35>

audit hémisphère supérieur et respectivement les directions descendantes appartenant audit hémisphère inférieur, ces deux points origine moyens Ms (qat) et Ms (Q) étant communs et utilisés en tant que points origine équivalents de diffusion pour l'ensemble des rayons incidents sur la cellule turbide (71,72).

13/-Procédé selon la revendication 12, caractérisé en ce qu'on calcule les points origine moyens Ms (Qt) et Ms () en tant que barycentres énergétiques des points origine moyens de diffusion Msk (out) et respectivement Msk (), calculés sur les différents rayons incidents sur la cellule turbide.

14/-Procédé selon la revendication 13, caractérisé en ce qu'avant de calculer les énergies ré-émises Wout (Ovj) par chaque cellule (71,72, 73), chaque point origine moyen Ms(##) et Ms(##) est remplacé par un point d'émission équivalent qui est le point Pe le plus proche choisi parmi une pluralité de points d'émission Pe possibles régulièrement répartis dans chaque cellule, pour lesquels toutes les trajectoires possibles des rayons issus de tous les points Pe selon toutes les directions de ré-émission Qvj possibles ont été préalablement calculées et enregistrées.

15/-Procédé selon la revendication 14, caractérisé en ce que pour déterminer les points d'émission Pe possibles, on divise une cellule (71, 72,73), en une pluralité de zones semblables et on choisit les milieux de ces zones en tant que points d'émission Pe possibles.

16/-Procédé selon l'une des revendications 1 à 15, caractérisé en ce qu'on calcule par itérations successives un rayonnement ré-émis par chaque cellule (71,72, 73) comprenant au moins un rayon ré-émis dans au moins une direction de ré-émission Ovj, en calculant et mémorisant à chaque itération : - le rayonnement incident sur la cellule (71,72, 73) issu des rayons ré-émis par les autres cellules adjacentes (71,72, 73), comme déterminé à l'itération précédente, - le rayonnement ré-émis par la cellule issu notamment de la diffusion du rayonnement intercepté dû au rayonnement incident sur la cellule.

<Desc/Clms Page number 36>

17/-Procédé selon l'une des revendications 1 à 16, caractérisé en ce que l'on définit préalablement dans la scène (5) indépendamment des cellules (71,72, 73), des éléments géométriques qui représentent différents éléments de paysage, et en ce que l'on calcule l'énergie ré- émise Wout (Qvj) en prenant en compte ces éléments géométriques en tant qu'éléments principalement réfléchissants et/ou absorbants.

18/-Procédé selon l'une des revendications 1 à 17, caractérisé en ce que la scène (5) comprenant une portion de sol (6) d'une planète, et une portion d'atmosphère (7) s'étendant au-dessus de cette portion de sol jusqu'à une altitude prédéterminée, on définit au-dessus de la scène une colonne atmosphérique (7) découpée en cellules, dites cellules atmosphériques (72), ayant des propriétés radiatives connues supposées homogènes dans l'ensemble de chaque cellule atmosphérique (72), toutes les cellules atmosphériques (72) de même altitude ayant les mêmes dimensions et propriétés radiatives.

19/-Procédé selon la revendication 18, caractérisé en ce que les dimensions des cellules atmosphériques (72) croissent avec l'altitude.

20/-Procédé selon l'une des revendications 18 ou 19, caractérisé en ce que la colonne atmosphérique (7) comprend une zone inférieure, dite basse atmosphère, au contact de la scène et une zone supérieure, dite atmosphère haute et intermédiaire, et en ce qu'à partir d'un rayonnement solaire prédéterminé défini par un rayon incident sur la face supérieure de chaque cellule (72) de la couche supérieure de la colonne atmosphérique (7), on réalise des exécutions successives suivantes :

1) on détermine le rayonnement incident sur la face supérieure des cellules de la basse atmosphère issu du rayonnement solaire après propagation dans l'atmosphère haute et intermédiaire,

2) on détermine le rayonnement réfléchi par la zone de paysage et la basse atmosphère dans l'atmosphère haute et intermédiaire à partir du rayonnement incident déterminé à l'exécution 1),

<Desc/Clms Page number 37>

5) on détermine le rayonnement final montant dans la colonne atmosphérique-notamment au niveau de sa face supérieure-issu de la somme des rayonnements réfléchis déterminé aux exécutions 2) et 4) après propagation à travers l'atmosphère haute et intermédiaire

4) on détermine le rayonnement réfléchi par la zone de paysage et la basse atmosphère vers l'atmosphère haute et intermédiaire à partir du rayonnement rétrodiffusé déterminé à l'exécution 3),

3) on détermine le rayonnement rétrodiffusé par l'atmosphère haute et intermédiaire vers la basse atmosphère à partir du rayonnement réfléchi déterminé à l'exécution 2),