Table des matières Python

Introduction à l'algorithme du lancer de rayons

1. Introduction

L'algorithme du lancer de rayons ([1][2]) permet d'effectuer la synthèse d'images de scènes tridimensionnelles. Il consiste à lancer des rayons depuis le point d'observation vers la scène, afin de déterminer la couleur à affecter à chaque point de l'image. Il est utilisé lorsqu'un grand réalisme est recherché (animation, architecture, etc).

Un autre algorithme utilisé pour obtenir des images de scènes est l'algorithme de rendu par objet, qui consiste à projeter les objets sur le plan image et à utiliser un tampon de profondeur (Z-buffer). Celui-ci est beaucoup plus rapide, ce qui explique son utilisation prépondérante dans le rendu en temps réel, par exemple pour les jeux vidéo. Il a cependant l'inconvénient de rendre difficile le rendu des éclairages complexes (multiples sources, diffusion de la lumière, ombres, réflexion, etc). L'algorithme du lancer de rayons est au contraire tout à fait adapté au rendu de l'éclairage.

Ce document explique le principe de l'algorithme de lancer de rayons, sous sa forme élémentaire, et présente un exemple en python.

2. Lancer de rayons

2.a. Principe optique

La formation d'une image par un système optique (œil, caméra, etc) obéit aux lois de la perspective, expliquées dans Projection perspective. Le système est constitué d'un centre optique C et d'un plan de projection (ou plan image) Π.

La lumière provenant d'un point M d'un objet suit le rayon MC pour parvenir au plan image.

figureA.svgFigure pleine page

La lumière qui part de M pour parvenir au plan image provient des sources de lumière. Elle peut en provenir directement, ou par l'intermédiaire de la diffusion et de la réflexion effectuée par les autres objets.

Le lancer de rayons consiste à considérer le trajet de la lumière dans le sens inverse du sens réel, en partant du centre optique C vers le point M. L'éclairement (et la couleur) du point image P (intersection de ce rayon avec le plan image) est égal au flux lumineux le long du rayon MP. Ce flux lumineux provient de

  • la lumière provenant directement des sources,
  • la lumière diffusée par d'autres objets.

On fait une distinction entre la réflexion spéculaire, qui suit la loi de Descartes de la réflexion, et la réflexion diffuse. Cette dernière est de loin la plus difficile à simuler, en raison du grand nombre de rayons qu'elle implique (en principe une infinité). Dans cette introduction, nous nous limitons à la première cause d'éclairement. Pour obtenir un rendu réaliste, il faut tenir compte de la lumière diffusée, qui est un phénomène très important dans la réalité. On suppose de plus que les sources de lumière sont ponctuelles.

Un autre phénomène à prendre en compte est la réflexion et la transmission éventuelles de la lumière au point M. Nous nous limiterons à la réflexion, mais la transmission se traite de manière analogue.

Le point de départ de l'algorithme est le lancer d'un rayon dans la direction CP, où P est le point de l'image dont on cherche l'éclairement. On commence par déterminer le point M d'intersection de ce rayon avec le premier objet rencontré. C'est ce point qui apparaît sur l'image en P.

Pour calculer le flux lumineux provenant directement des sources, on lance des rayons du point M vers les sources, qui sont ici ponctuelles. Ces rayons sont appelés rayons d'éclairement (représentés en bleu sur la figure). La figure montre deux rayons d'éclairement : le premier, lancé vers la source 1, l'atteint effectivement. Ce rayon contribue donc à l'éclairement de la surface au point M. Le second, lancé vers la source 2, ne l'atteint pas à cause de la présence d'un objet sur son chemin. La surface de l'objet au point M peut être plus ou moins réfléchissante. Pour prendre cela en compte, on lance à partir du point M un rayon réfléchi. Si ce rayon réfléchi rencontre un objet en un point M', la lumière provenant de M' parvient au point P par la réfléxion en M : l'éclairement du point M' contribue donc à l'éclairement du point P sur le plan image.

On voit bien la récursivité de l'algorithme : le flux lumineux qui quitte M' en direction de M se calcule de la même manière, en lançant des rayons vers les sources, et en lançant un rayon réfléchi si la surface est réfléchissante au point M'.

Dans certains cas, le tracé récursif des rayons réfléchis peut conduire à une boucle sans fin, ou du moins limitée par le nombre de récursions autorisée. Il faut donc limiter volontairement le nombre de récursion, ce qui revient à limiter le nombre de réflexions multiples.

2.b. Modèle de surface

Voyons comment une source de lumière ponctuelle contribue au flux lumineux. La figure suivante montre le rayon incident, le rayon envoyé vers la source et la normale à la surface au point M.

figureB.svgFigure pleine page

Il faut bien noter que le sens réel de la lumière est inverse de celui représenté par les flèches, qui correspond aux sens des rayons calculés. On cherche le flux lumineux envoyé à partir du point M dans la direction du rayon incident. En général, ce flux dépend des angles entre les deux rayons et de la normale n à la surface au point M.

Le modèle le plus simple pour représenter les propriétés optiques d'une surface est le modèle de Lambert. Il consiste à attribuer le facteur cos(θ) pour l'éclairement du point M, et à supposer que le flux lumineux envoyé dans la direction du rayon incident est indépendant de l'angle qu'il fait avec la normale. On introduit aussi un coefficient de diffusion pour la surface, et le flux s'écrit donc :

I est l'intensité de la source.

Les grandeurs Φ, I et kd dépendent en fait de la longueur d'onde λ. La couleur de l'objet dépend de la fonction kd(λ). La source elle-même peut déliver une intensité I(λ) fonction de la longueur d'onde. Une manière simple de traiter l'influence de la longueur d'onde est de considérer trois composantes RVB (rouge,vert, bleu) pour ces grandeurs. Par exemple, la surface d'un objet rouge aura une composante rouge de son coefficient de diffusion plus grande que les composantes bleue et verte. Ce mode de calcul permet d'obtenir directement les composantes RVB du point P sur l'image finale.

Dans la suite de cette introduction, on se contente de calculer des flux monochromes, sans tenir compte de l'influence de la longueur d'onde.

En ce qui concerne le rayon réfléchi, on doit définir un coefficient de réflexion. Physiquement, ce coefficient de réflexion dépend de l'angle entre le rayon incident et la normale. Pour simplifier, on considèrera un coefficient de réflexion indépendant de l'angle.

2.c. Algorithme

L'algorithme est principalement constitué d'une fonction lancer_rayon(rayon) qui cherche le point de la scène rencontré par un rayon et calcule le flux lumineux le long de ce rayon. La fonction s'appelle elle-même récursivement dans le cas où le point rencontré est sur une surface réfléchissante.

Voici la séquence des opérations effectuées par cette fonction :

  • Calcul des intersections avec les tous les objets de la scène. Lorsqu'il y a une intersection, on stocke la distance t (analogue à un temps de parcours) entre le point de départ et le point d'intersection.
  • Si aucun objet n'est rencontré, renvoi de l'éclairement du fond.
  • L'intersection la plus proche est retenue, celle dont la valeur de t est la plus petite.
  • Le flux est initialisé à zéro.
  • Pour chaque source ponctuelle, envoi d'un rayon vers la source. Les intersections de ce rayon avec les objets de la scène sont recherchées. Si aucune intersection n'est trouvée, la contribution de cette source au flux lumineux est ajoutée, en utilisant un modèle de surface.
  • Si la surface est réfléchissante, calcul du rayon réfléchi en utilisant la loi de Descarte de la refléxion spéculaire. La fonction lancer_rayon(rayon) est appelée (appel récursif) et le résultat est ajouté au flux, avec un facteur de multiplication (coefficient de réflexion).
  • Renvoi du flux.

Pour effectuer la synthèse d'une image complète, il faut échantillonner le plan image, par exemple avec un quadrillage régulier, puis appeler la fonction lancer_rayon(rayon) pour chaque rayon CP et attribuer l'éclairement obtenu au point P du plan image. Une manière simple de le faire est de lancer un rayon pour chaque pixel de l'image finale.

2.d. Exemple de calcul d'intersection

Le calcul de l'intersection d'un rayon avec un objet de la scène dépend évidemment de la forme de l'objet. Nous allons considérer le cas simple d'un objet sphérique.

Pour chaque rayon lancé, il faut rechercher l'intersection avec les N objets de la scène. Pour les scènes très complexes, cela fait beaucoup de calculs, d'autant plus qu'il y a au moins autant de rayons que de pixels dans l'image. On est alors amené à utiliser une méthode d'accélération, qui consiste à éliminer le plus rapidement possible les objets qui ne sont pas coupés par le rayon. Nous n'aborderons pas ce sujet dans cette introduction. On note seulement qu'une méthode d'accélération (parmis d'autres) consiste à définir pour chaque objet une sphère englobante. On commence par rechercher l'intersection du rayon avec cette sphère. Si cette intersection existe, on peut aborder le calcul de l'intersection avec l'objet, qui peut être beaucoup plus complexe.

Un rayon est défini par un point initial P et un vecteur unitaire U. Une sphère est définie par son centre C et son rayon r. Voyons comment déterminer l'intersection éventuelle du rayon avec la sphère. Une première idée est de considérer l'équation paramétrique du rayon et de rechercher le paramètre t pour l'intersection avec la sphère : cela donne une équation du second degré pour t. L'inconvénient de cette méthode est qu'elle nécessite beaucoup de calculs (en particulier le calcul d'une racine carrée), même lorsque le rayon ne coupe pas la sphère, ce qui est finalement le cas le plus fréquent.

Une méthode géométrique ([1]) permet d'éliminer rapidement le cas du rayon ne rencontrant pas la sphère, et de calculer le paramètre t seulement lorsque cela est nécessaire.

figureC.svgFigure pleine page

On commence par comparer la distance PC avec le rayon r pour savoir si le point de départ du rayon se trouve dans ou hors la sphère. Soit H le projeté orthogonal de C sur le rayon. La distance algébrique de P à H est

Si cette distance est négative et si P est en dehors de la sphère, alors le rayon ne coupe pas la sphère.

Le carré de la distance entre le centre de la sphère et le point H est :

Le carré de la distance MH est :

Si ce carré est négatif, c'est que le point H est hors de la sphère et donc que le rayon ne coupe pas la sphère. On peut finalement calculer la distance t=PM. Si le point P est hors de la sphère, on a :

Si le point P est dans la sphère, on a :

Si la valeur de t est nulle, c'est que le point initial du rayon se trouve déjà sur la surface. Dans ce cas, on considère qu'il n'y a pas d'intersection. Ce cas doit en principe se produire lorsqu'on lance un nouveau rayon à partir de M (rayon d'éclairement ou rayon réfléchi). Cependant, des erreurs d'arrondis dans les calculs peuvent placer le point M légèrement à l'intérieur de la sphère, ce qui conduira à l'obtention d'une intersection avec la sphère, qui n'a physiquement pas lieu puisque le point se trouve en fait sur la surface. Pour éviter ce problème, on modifie légèrement la valeur de t afin que le point M soit de manière certaine hors de la sphère, par exemple en multipliant t par 0.999.

3. Implémentation Python

On propose une implémentation simplifiée, avec des calculs géométriques à 2 dimensions (l'extension à 3 dimensions est immédiate). La couleur n'est pas pris en compte : le flux lumineux est une grandeur scalaire. La représentation graphique des sphères et des rayons est faite avec le module Pygame.

lancerRayons2d.py
import math
import pygame
              

Un rayon est représenté par un objet de la classe Rayon. Le constructeur prend en argument le point initial et le vecteur directeur. Celui-ci est normalisé par le constructeur. La couleur est utilisée pour la représentation graphique, laquelle est effectuée par la fonction dessiner.

class Rayon:
    def __init__(self,P,U):
        self.P = P # point initial
        norm = math.sqrt(U[0]*U[0]+U[1]*U[1])
        self.U = [U[0]/norm,U[1]/norm] # vecteur unitaire
        self.t = 0.0
        self.couleur = [255,0,0]
    def dessiner(self,screen,echelle):
        Q = [self.P[0]+self.U[0]*self.t,self.P[1]+self.U[1]*self.t]
        pygame.draw.line(screen,self.couleur,[self.P[0]*echelle,self.P[1]*echelle],[Q[0]*echelle,Q[1]*echelle])

              

La classe Intersection permet de mémoriser les propriétés d'un point d'intersection d'un rayon avec une surface : rayon incident, distance entre le point de départ du rayon et le point d'intersection, point d'intersection, normale, coefficients de réflexion et de diffusion. La fonction rayon_reflechi calcule le rayon réfléchi avec la relation :

La fonction cosinus renvoit le cosinus entre un rayon et la normale.

class Intersection:
    def __init__(self,rayon_incident,t,P,N,ref,dif):
        self.rayon_incident = rayon_incident
        self.t = t # temps
        self.P = P # point
        self.N = N # normale
        self.ref = ref # coef reflexion
        self.dif = dif # coef diffusion
    def rayon_reflechi(self):
        U = self.rayon_incident.U
        ps = self.N[0]*U[0]+self.N[1]*U[1]
        R0 = U[0]-2*ps*self.N[0]
        R1 = U[1]-2*ps*self.N[1]
        return Rayon(self.P,[R0,R1])
    def cosinus(self,rayon):
        U = rayon.U
        return self.N[0]*U[0]+self.N[1]*U[1]

               

La classe Objet est la classe générique des objets de la scène. La fonction intersection_rayon calcule l'intersection avec un rayon. S'il n'y pas d'intersection, la valeur False est renvoyée. La fonction dessiner effectue la représentation graphique de l'objet.

class Objet:
    def __init__(self):
        pass
    def intersection_rayon(self,rayon):
        # renvoit l'intersection si elle existe
        return False
    def dessiner(self,screen,scale):
        pass

              

La classe Sphere représente une sphère dont on donne le centre, le rayon, le coefficient de réflexion et le coefficient de diffusion.

class Sphere(Objet):
    def __init__(self,C,r,ref,dif):
        Objet.__init__(self)
        self.C = C
        self.r = r
        self.r2 = r*r
        self.ref = ref
        self.dif = dif
    def intersection_rayon(self,rayon):
        P = rayon.P
        U = rayon.U
        dPCx = self.C[0]-P[0]
        dPCy = self.C[1]-P[1]
        dPC2 = dPCx*dPCx+dPCy*dPCy
        if dPC2 > self.r2:
            exterieur = True
        else:
            exterieur = False
        tph = dPCx*U[0]+dPCy*U[1]
        if tph<0 and exterieur:
            return False
        D2 = dPC2-tph*tph
        tmh2 = self.r2-D2 
        if tmh2<0:
            return False
        tmh = math.sqrt(tmh2)
        if exterieur:
            t = tph-tmh
        else:
            t = abs(tph)+tmh
        if t==0:
            return False
        t = t*0.999
        I = [P[0]+U[0]*t,P[1]+U[1]*t]
        N = [I[0]-self.C[0],I[1]-self.C[1]]
        norm = math.sqrt(N[0]*N[0]+N[1]*N[1])
        N = [N[0]/norm,N[1]/norm]
        return Intersection(rayon,t,I,N,self.ref,self.dif)
    def dessiner(self,screen,echelle):
        r = self.r*echelle
        pygame.draw.ellipse(screen,[0,0,0],[self.C[0]*echelle-r,self.C[1]*echelle-r,2*r,2*r],2)
  
              

La classe Source représente une source de lumière ponctuelle, donnée par sa position et son intensité. La fonction rayon_eclairage calcule le rayon d'éclairage à lancer vers cette source depuis un point P donné.

class Source:
    def __init__(self,S,i):
        self.S = S
        self.i = i
    def rayon_eclairage(self,P):
        U = [self.S[0]-P[0],self.S[1]-P[1]]
        norm = math.sqrt(U[0]*U[0]+U[1]*U[1])
        U = [U[0]/norm,U[1]/norm]
        rayon = Rayon(P,U)
        rayon.t = norm
        rayon.couleur = [0,0,255]
        return rayon
    def dessiner(self,screen,echelle):
        r = 0.1*echelle
        pygame.draw.ellipse(screen,[0,0,255],[self.S[0]*echelle-r,self.S[1]*echelle-r,2*r,2*r],2)

              

La classe Scene permet de stocker les différents éléments de la scène : sources, objets et rayons. Elle contient la fonction lancer_rayon qui effectue le calcul du flux lumineux pour un rayon. L'argument niveau sert à contrôler le niveau de récursion. Sa valeur initiale doit être nulle.

class Scene:
    def __init__(self,fond):
        self.fond = fond
        self.liste_objets = []
        self.nombre_objets = 0
        self.liste_rayons = []
        self.liste_sources = []
        self.max_recursion = 20
    def ajouter_objet(self,objet):
        self.liste_objets.append(objet)
        self.nombre_objets += 1
    def ajouter_source(self,source):
        self.liste_sources.append(source)
    def dessiner_objets(self,screen,echelle):
        for objet in self.liste_objets:
            objet.dessiner(screen,echelle)
    def dessiner_rayons(self,screen,echelle):
        for rayons in self.liste_rayons:
            rayons.dessiner(screen,echelle)
    def dessiner_sources(self,screen,echelle):
        for source in self.liste_sources:
            source.dessiner(screen,echelle)
    def dessiner(self,screen,echelle):
        self.dessiner_objets(screen,echelle)
        self.dessiner_rayons(screen,echelle)
        self.dessiner_sources(screen,echelle)
    def lancer_rayon(self,rayon,niveau):
        if niveau > self.max_recursion:
            return self.fond
        liste_intersections = []
        for objet in self.liste_objets:
            intersection = objet.intersection_rayon(rayon)
            if intersection:
                liste_intersections.append(intersection)
        if len(liste_intersections)==0:
            return self.fond
        premiere_intersection = liste_intersections[0]
        for intersection in liste_intersections:
            if intersection.t < premiere_intersection.t:
                premiere_intersection = intersection
        intersection = premiere_intersection
        rayon.t = intersection.t
        self.liste_rayons.append(rayon)
        flux = 0.0
        for source in self.liste_sources:
            rayon_eclairage = source.rayon_eclairage(intersection.P)
            eclairage = True
            for objet in self.liste_objets:
                if objet.intersection_rayon(rayon_eclairage):
                    eclairage = False
            if eclairage:
                flux += intersection.dif * source.i * intersection.cosinus(rayon_eclairage)
                self.liste_rayons.append(rayon_eclairage)
        if intersection.ref != 0:
            rayon_reflechi = intersection.rayon_reflechi()
            flux += intersection.ref*self.lancer_rayon(rayon_reflechi,niveau+1) # appel recursif
        return flux

              

4. Exemple

Voici un exemple avec 4 sphères réfléchissantes et trois sources. Les rayons sont lançés à partir d'un point situé à l'infini, ce qui revient à faire une projection orthoscopique. Les rayons sont donc parallèles.

import math
import pygame
import numpy
from matplotlib.pyplot import *
from lancerRayons2d import *
            
scene = Scene(0.0)
scene.ajouter_objet(Sphere([5.0,5.0],1.0,0.5,1.0))
scene.ajouter_objet(Sphere([7.0,2.0],1.0,0.5,1.0))
scene.ajouter_objet(Sphere([2.0,5.0],0.7,0.5,1.0))
scene.ajouter_objet(Sphere([8.0,7.0],0.8,0.5,1.0))
scene.ajouter_source(Source([1,8],1.0))
scene.ajouter_source(Source([9,0.5],1.0))
scene.ajouter_source(Source([0.5,0.5],1.0))
N=100
x = numpy.arange(N)*10.0/N
flux = numpy.zeros(N)
for k in range(N):
    flux[k] = scene.lancer_rayon(Rayon([x[k],0],[0,1]),0)
pygame.init()
screen = pygame.display.set_mode([500,500])
pygame.display.set_caption("Rayons")
clock = pygame.time.Clock()
done = False
echelle = 50.0
while not done:
    clock.tick(5)
    for event in pygame.event.get():
        if event.type == pygame.QUIT:
            done = True
    screen.fill([255,255,255])
    scene.dessiner(screen,echelle)
    pygame.display.flip()
pygame.image.save(screen,"../../../../figures/graphie/rayons/introrayons/rayons1.png")
pygame.quit()
figure(figsize=(7,5))
plot(x,flux,'o')
xlabel("x")
ylabel("flux")
            
figAfigA.pdf

Voici le tracé des rayons : en rouge les rayons issus du centre optique (rayons parallèles) et les rayons réfléchis. Seuls les rayons rencontrant une sphère sont représentés. En bleu sont représentés les rayons d'éclairement, ceux qui parviennent effectivement à la source.

rayons
Références
[1]  R. Malgouyres,  Algorithmes pour la synthèse d'images et l'animation 3D,  (Dunod, 2055)
[2]  P. Shirley, S. Marschner,  Fundamentals of computer graphics,  (A K Peters, 2009)
Creative Commons LicenseTextes et figures sont mis à disposition sous contrat Creative Commons.