Potatro, RayMarching et DistanceFields : une histoire de SphereTracing

Le « Raymarching + DistanceFields » était la technique la plus utilisée dans les intros 4k présentées en 2009 (et continuera probablement d'avoir une forte influence en 2010). Chez Frequency, xt95 a introduit cette technique grâce au merveilleux travail de iq et de SystemKit. Lorsque wullon a demandé à Frequency de produire une petite intro pour annoncer la construction du site internet demoscene.fr, xt95 a construit un brouillon prometteur de Potatro à partir duquel j'ai pu travailler pour améliorer le design global, les formes « isos » et les scènes. Même si je n'avais pas pris le temps de comprendre le principe du « RayMarching + DistanceFields », cela ne me semblait pas si abscons.

De retour de la #main 2009, en regardant le BBS de pouet, je suis passé sur la discussion : « So, what do remote field equations look like? And how do we solve them? »(1). Question triviale, mais les premières réponses n'ont pas éclairci ma basique compréhension que j'avais auparavant, mais le but principal restait toujours obscur pour moi. Mais Ferris posta alors un lien sur l'article de John C. Hart, « Sphere tracing: a geometric method for the antialised ray tracing of implicit surfaces » de 1996 et iq se référa à deux articles, un premier article plus vieux de John C. Hart, « Ray Tracing Deterministic 3-D Fractals » de 1989 et une reprise du principe dans le « Chapitre 8. Per-Pixel Displacement Mapping With Distance Functions » du livre GPU Gems 2 par NVIDIA datant de 2004. En parcourant l'article original de John C. Hart tout devint clair ! Cet article est un essai introduisant quelques principes de base du « RayMarching + DistanceFields », plus explicitement appelés « SphereTracing ».

Je ne vais pas détailler le principe du « lancer de rayon », mais je vais juste le mettre en relation avec le « RayMarching ». Pour effectuer le rendu, le « lancer de rayon » doit calculer les intersections du rayon provenant de la caméra avec un objet de la scène (pour extraire les informations du matériel de la surface : la couleur diffuse, réfraction, diffraction, etc.) et donc, calculer les formules d'intersections pour chaque objet. Dans le cas où les objets sont modélisés à partir d'un ensemble de triangles, il est nécessaire de calculer les intersections du rayon avec chaque triangle dans la scène pour déterminer quel est l'objet concerné par le rayon. Le calcul des intersections est très facile à implémenter, les formules sont établies, donc il n'y a aucun problème, sauf du côté des performances, car nous devons calculer les intersections du rayon avec un nombre substantiel de surfaces (même s'il y a des techniques pour accélérer, comme la division de l'espace avec, par exemple, les K-d Trees).

Par contre, lorsque l'on veut afficher des surfaces implicites continues ou discontinues (aussi appelées « courbes de niveau » ou « isosurfaces »), cela devient difficile de calculer l'intersection, car cela requiert le calcul de l'intersection entre l'équation de la surface (Fig. 1) et le trajet du rayon. Ce calcul serait spécifique à chaque équation. Donc, si le principe du rendu par lancer de rayon est basé sur l'intersection d'un rayon avec les objets de la scène, comment pouvons-nous calculer efficacement l'intersection du rayon avec une surface implicite ?

Fig. 1 Définition d'une surface implicite (voir Wikipédia)

Avant de parler de « RayMarching », faisons un petit rappel pour faciliter notre discussion. Dans l'exemple donné d'une seule sphère de rayon r (et d'origine 0) : la formule mathématique est :

Sous la forme d'une distance, l'équation de la sphère peut être écrite différemment :

Pour calculer l'intersection avec sp1, une sphère de rayon r1 en utilisant la technique de « RayMarching » avec un rayon partant de p1(x1,y1,z1) et ayant une direction v1(VX1,VY1,VZ1), nous allons à partir de l'origine du rayon (x1,y1,z1), calculer la distance menant à la surface de la sphère, en remplaçant les éléments suivants :

en utilisant les propriétés suivantes pour d1 :

• si d1 = 0, alors le point (x1, y1, z1) est placé exactement sur la surface de la sphère ;
• si d1 < 0, alors le point (x1, y1, z1) est à une distance | d1 | dans la surface ;
• si d1 > 0, alors le point (x1, y1, z1) est à une distance | d1 | en dehors de la surface.

Si d1 > 0, le point p1 est en dehors de la sphère, il suffit de calculer l’intersection i1, nous devons nous déplacer à partir de p1 d'une distance de d1 dans la direction v1 : i1 = p1 + d1v1.

Dans le cas d'une seule sphère, le calcul de l'intersection en utilisant la technique de « RayMarching » est trivial (et n'est pas plus intéressant qu'une résolution analytique), mais si nous avons des formes plus complexes, nous ne pourrons pas déterminer l'intersection d'un seul coup. Prenons l'exemple d'une seconde sphère qui serait ajoutée à la scène précédente :

Nous souhaitons déterminer quelle sphère (p1, v1) est touchée par le rayon. Avec une seule sphère, nous avons déterminé l'intersection en déplaçant simplement p1 à d1.v1. Mais avec une autre sphère sp2, comment pouvons-nous avancer sans toucher sp2 ? Si nous nous déplaçons à sp1 sans vérification et que p2 est sur le chemin, nous allons sauter sp2 et louper l'intersection avec sp2.

Par contre, nous connaissons la distance entre p1 et sp2, qui est :

L'astuce est de se déplacer dans la direction v1 de la distance minimale de tous les objets de la scène. Donc nous nous déplaçons d'un pas vers l'intersection, un pas qui est égal à minD :

Tant que minD > 0, cela signifie que nous n'avons encore touché aucun objet. Nous pouvons généraliser l'avancée du rayon par :

Pour expliquer un peu mieux la progression du rayon et les interactions complexes avec les objets, vous pouvez télécharger SphereTracing2DApp.zip, une application 2D développée en C# qui permet de voir la progression du rayon et le nombre d'étapes qu'il doit effectuer avant de trouver l'intersection :

Le point blanc dans le coin bas gauche est l'origine du rayon (la caméra). La ligne rouge est la direction du rayon. Le premier cercle blanc est la distance minimale entre l'origine et tous les objets de la scène (ici, le rectangle bleu). Le rayon peut donc se déplacer de cette distance sans rencontrer un seul objet. L'opération est ainsi répétée jusqu'au point de rencontre avec un objet (où minD < epsilon). La ligne verte nous indique quel objet a été sélectionné comme « distance field » le plus proche.

Nous comprenons maintenant que le calcul d'une intersection du rayon de la caméra avec une surface iso est en réalité une approximation incrémentale d'une sphère de collision se déplaçant suivant le rayon d'intersection.

Comme nous l'avons vu, l'algorithme repose sur le déplacement d'une sphère d'un rayon équivalent à la plus petite distance entre le rayon d'origine et la surface de tous les objets de la scène. Cette technique a été nommée « SphereTracing » par John.C.Hart. Tandis que dans la demoscene, iq a introduit le terme de « DistanceFields », mais a ainsi éclipsé le terme original qui me semble plus explicite. Donc lorsque nous parlons de « RayMarching-DistanceFields » nous appliquons la technique de « SphereTracing » !

L'algorithme général de « SphereTracing » est le suivant :

ro // r0 : position initiale de la caméra
v // v : Direction de la caméra

d = 0 // d : distance vers l'objet  
p = ro // p : prochaine position de la caméra
do (  
    d = iso (r) // Distance vers objet  
    p = p + d * v;  
    step++;  
) while (d < epsilon && step < MAX_STEP)

La fonction iso retourne la distance minimale r à la surface de tous les objets de la scène. Dans le cas des deux sphères, la fonction iso aurait été :

Si aucune intersection n'est rencontrée, l'algorithme de « SphereTracing » est arrêté grâce à une distance maximale ou en limitant le nombre d'étapes (MAX_STEP).

La sphère n'est bien sûr pas le seul élément que nous pouvons utiliser pour modéliser une surface iso. Comme nous l'avons vu, n'importe quelle fonction f(x, y, z) devrait fonctionner (pourvu que la fonction ne retourne pas de valeurs introduisant trop de discontinuités). Donc, il est assez facile d'utiliser d'autres primitives (boîte, cône, tore…) comme décrit dans l'article de John C. Hart.

Pour les intros 4k, c'est impressionnant de voir la complexité des scènes qui peuvent être produites avec ce système. Jusqu'à présent je n'ai parlé que d'une seule sphère. Comment pouvons-nous alors construire des scènes complexes ? En utilisant le principe Constructive Solid Geometry (CSG) comme le montre John.C.Hart, une des plus vieilles techniques de modélisation 3D : c'est suffisant pour effectuer des « opérations logiques » sur les formules de surfaces implicites : intersection, union, complément…

Donc, étant donnés deux objets, chacun ayant sa fonction iso d1 et d2 respectivement, peuvent être calculés :

• l'union des deux objets (d1 ⋃ d2), il suffit de prendre le minimum de d1 et d2 : iso(p) = min(d1(p),d2(p)) ;
• l'intersection des deux objets (d1 ⋂ d2), il suffit de prendre le max de d1 et d2 : iso(p) = max (d1(p),d2(p)) ;
• l'inverse d'un objet (-d1), il suffit d'inverser le signe : iso(p) = -d1(p) ;
• la soustraction des deux objets (d1 - d2), il suffit d'effectuer une intersection entre d1 et l'inverse de d2 : iso(p) = max(d1(p),-d2(p)).

La translation peut être effectuée directement dans la formule iso. Pour une sphère de centre c (cx, cy, cz) et de rayon r, la distance géométrique au point p(x, y, z) est :

tout comme pour le redimensionnement, la rotation, la torsion, etc. Voici un exemple de torsion sur l'axe des y :

Dans sa présentation à la NVScene8, iq a fourni une liste complète de fonctions : modulo, mélange, bruit, etc. De multiples techniques utilisent aussi texture3D pour créer des transformations plus riches sur les surfaces implicites (voir Clouds Dream ou Rudebox #main 2009).

Une fois que vous sortez de la boucle du « SphereTracing », si d < epsilon, p est sur la surface de l'objet observé. Pour calculer la normale de la surface, nous utilisons un gradient dans p, nous calculons une approximation de la dérivée de la surface. Prenons n(p), la normale p et eps assez petit :

À partir de là, il est assez simple d'implémenter les réflexions, les ombres, les lumières, un « Object-Space-Ambiant-Occlusion », le bump mapping… cela n'est limité que par les ressources GPU !

Tout comme nous l'avons vu, le nombre d'étapes requises pour le « SphereTracing » peut être élevé, notamment lorsqu'un rayon est proche d'un objet sans le toucher (voir l'application C# pour comprendre le phénomène). Dans ce cas, le rayon avance par petits pas et le nombre de calculs devient énorme (plusieurs fois le calcul de la formule iso d'un simple point). Il y a différentes méthodes pour l'éviter, comme utiliser un pas minimum.

Pour cela, voici quelques astuces à suivre :

• minimiser les if/then/else/switch dans le code, les sous-fonctions, etc. Au mieux, éliminez-les ;
• par conséquent, évitez la gestion de la scène dans le shader (if scene == 0 alors cette fonction iso, sinon une autre), mais gérez-la en dehors du shader (par le biais d'une compilation de différents shaders en GLSL, ou l'utilisation de différentes passes/techniques en HLSL) ;
• construisez des shaders spécifiques par scène (un qui supporte la réflexion, un autre pour l'ambiant occlusion + réflexion) et mettez en inline l'initialisation et le format iso. Pour cela, aidez-vous de techniques comme la « templatization » telle que la concaténation/substitution de code dans l'initialisation du shader (ce qui a été fait dans la démonstration Rudebox de Alcatraz avec des directives préprocesseur #define, ou dans Muon baryon avec la concaténation/substitution en C) ;
• effectuez une grossière étape de précalcul de la scène dans une texture 3D afin d'accélérer la passe du « SphereTracing » puis appliquez le « SphereTracing » sur chaque pixel de l'écran lorsque nécessaire.

Nous avons vu le principe du « SphereTracing » ou du « RayMarching-DistanceFields », une technique largement à la mode dans les intros 4k récentes. Cette technique a l'avantage d'être particulièrement simple à implémenter et incroyablement efficace dans son rendu. C'est un pari sûr de dire que cette technique sera toujours présente en 2010 et que les scènes deviendront de plus en plus complexes et riches (avec, par exemple, la création d'outils spécifiques pour aider la modélisation de scènes CSG).

Merci à Alexandre Mutel de nous avoir permis de traduire son article.

Merci à Winjerome pour sa relecture et _Max_ et ClaudeLELOUP pour leur correction orthographique.