Raytracer minimaliste en Haskell : Gestion des lumières

L'aliasing est un phénomène inhérent aux calculs de nature discrète. Il regroupe les deux phénomènes communs suivants :

• Le phénomène de crénelage, qui se caractérise par l'apparition d'un motif en escalier le long du bord des objets. Il est fortement visible sur toutes les images que nous avons générées. Vous trouverez d'autres exemples sur la page wikipedia.
• Un second phénomène très courant, que nous n'avons pas encore rencontré (faute de disposer d'assez d'objets), mais qui se produit trop souvent, par exemple dès qu'il est question d'appliquer une texture sur un plan, est l'apparition de motifs d'interférences. Cela est dû au fait que l'on superpose notre grille de pixels (que l'on peut considérer comme un réseau) sur un réseau formé par des objets / un motif qui se répète. Des exemples sont disponibles sur cette page ainsi que sur wikipedia.

Il existe plusieurs façons de réduire l'effet d'aliasing. L'une des plus simples consiste à lancer 4 rayons sur les sommets du pixel (vue comme un carré) plutôt qu'un unique rayon central, et de faire la moyenne. Cela revient à générer une image de Largeur + 1 par Hauteur +1 pixels. On est donc dans la famille du sur-échantillonnage, ou encore supersampling.

La solution que je vous propose est du même genre, mais légèrement plus raffinée (et malheureusement, coûteuse en termes de performances). Nous allons subdiviser chaque pixel en kitxmlcodeinlinelatexdvpn . nfinkitxmlcodeinlinelatexdvp sous-pixels, où n est le taux de sur-échantillonnage. Pour chacun de ces sous-pixels, on lancera un rayon dont la position est aléatoire, mais reste à l'intérieur de ce sous-pixel. De cette façon, on construit une "grille irrégulière", ce qui permet de réduire les interférences. Ensuite, pour obtenir la couleur du pixel, on moyenne tous les sous-pixels qui le composent.

Plutôt que de moyenner tous les sous-pixels de façon aveugle, on pourrait les pondérer selon une courbe gaussienne, ou encore courbe en cloche.

On commence par une fonction utilitaire, qui va nous permettre de subdiviser une liste infinie de nombres aléatoires en sous-listes de taille fixe.

--Coupe une liste en une liste de morceaux de longueur length
spliter :: Int -> [a] -> [[a]]
spliter n [] = []
spliter n list = first : (spliter second)
  where
    (first, second) = splitAt n list

On ajoute aussi la possibilité de produire une liste de valeurs aléatoires qui restent identique à chaque exécution (on ne souhaite pas que deux exécutions de notre programme avec la même scène nous donne deux valeurs différentes).

import System.Random

-- ...

--Produit des valeurs entre -offset et +offset
randomOffsets :: RealRep -> [RealRep]
randomOffsets offset = randomRs (0, offset) g
  where
    g = mkStdGen 42

On va vouloir manipuler des pixels aux coordonnées non entières. Puisqu'au final, nos pixels ne sont que des coordonnées sur une grille, pourquoi ne pas les stocker dans des RealReps ?

--Une coordonnée dans l'image
type PixelLocation = (RealRep, RealRep)

Les conversions de computeRay deviennent alors inutiles et doivent être modifiées.

--Construit un rayon à partir d'un pixel
computeRay :: PixelLocation -> Scene -> Ray
computeRay (x, y) (Scene camera objects _) = normalise (vx, vy, vz)
  where Camera _ distance (planWidth, planHeight) _ = camera
        halfPlanWidth = int2RealRep (planWidth `div` 2)
        halfPlanHeight = int2RealRep (planHeight `div` 2)
        vx = x - halfPlanWidth
        vy = y - halfPlanHeight
        vz = -distance

Nous allons modifier la fonction getPixel pour qu'elle prenne un paramètre entier supplémentaire : le taux de sur-échantillonnage. La valeur 1 correspondra à lancer un seul rayon, alors que n signifiera n*n rayons lancés. Dans cette première version, c'est un sur-échantillonnage où le rayon sera toujours lancé dans l'angle inférieur gauche (x minimal, y minimal) du sous pixel.

type SamplingRate = Int
getPixel :: SamplingRate -> Scene -> PixelLocation -> [Word8]
getPixel samplingRate scene pixelLocation@(px, py) = colorToPixel color
  where
    --Conversion en nombre à virgule
    samplingRate' = int2RealRep samplingRate
    --Produit l'intersection
    pixelIntersection pixelLocation = return pixelLocation
                                 >>= computeRay
                                 >>= computeIntersections
                                 >>= getClosestIntersection
    --Construit la couleur du pixel
    Scene _ _ lights = scene
    intersectionToPixel intersection = computeColor intersection $ scene
    maybePixelColor pixelLocation = fmap intersectionToPixel (pixelIntersection pixelLocation scene)
    pixelColor maybeColor = case maybeColor of
      Nothing -> (0, 0, 0)
      Just color -> color
    --Calcule l'écart entre deux sous-pixels
    offset = 1 / samplingRate'
    --Construit une liste de couleurs à partir des sous-pixels
    colorList = map (pixelColor.maybePixelColor) (pixels samplingRate)
    --Construit la liste des sous-pixels
    pixels 1 = [pixelLocation]
    pixels samplingRate = [(px + x, py + y) | x <- subCells, y <- subCells]
      where
        subCells = map (\x -> (int2RealRep x) * offset - 1/2) [0..samplingRate-1]
    --Fait la moyenne de toutes les sources lumineuses
    color = (offset*offset) *! (vectorSum colorList)

Vous pouvez maintenant tester cette anti-aliasing avec le main (main si ça veut le coup).

Pour rajouter un comportement aléatoire, nous allons légèrement modifier la fonction getPixel pour qu'elle puisse prendre en paramètre une liste d'offsets à appliquer à chaque sous-pixel (elle se voit donc renommée).

type SamplingRate = Int
getPixel :: SamplingRate -> Scene -> PixelLocation -> [Word8]
getPixel sr sc px = stochasticGetPixel sr sc px []

type RandList = [RealRep]
stochasticGetPixel :: SamplingRate -> Scene -> PixelLocation -> RandList -> [Word8]
stochasticGetPixel samplingRate scene pixelLocation@(px, py) randList = colorToPixel color
  where
    -- ...
    --Construit la liste des sous pixels
    pixels 1 = [pixelLocation]
    pixels samplingRate = [(px + x, py + y) | x <- subCells, y <- subCells]
      where
        subCells = map (\x -> (int2RealRep x) * offset - 1/2) [0..samplingRate-1]
        stocasticCells = if null randList then subCells else newSubCells
          where
            newSubCells = zipWith (+) subCells randList
    -- ...

On ajoute quelques fonctions utilitaires pour appeler plus simplement stochasticGetPixel :

getPixelMap :: SamplingRate -> Scene -> [PixelLocation] -> [Word8]
getPixelMap samplingRate scene locations = pixels
  where
    stochasticPixels = map
                       (stochasticGetPixel samplingRate scene)
                       locations
    pixels = concat $ zipWith ($) stochasticPixels (generateRandLists samplingRate)

generateRandLists :: SamplingRate -> [RandList]
generateRandLists samplingRate = spliter (samplingRate*samplingRate) $ randomOffsets (1 / samplingRate')
  where
    samplingRate' = int2RealRep samplingRate

On remplace maintenant l'appel dans le main :

tableData width height = listArray ((0,0,0), (height-1, width-1, 3)) $ pixels
  where
    pixels = getPixelMap 2 scene [(int2RealRep x, int2RealRep y)
                                 | y <- [0..height - 1], x <- [0..width - 1]]

Et c'est terminé! On se quitte sur une scène illustrant tout ce que nous avons ajouté au cours de cet article.

Un grand merci à Cédric Duprez pour sa relecture orthographique.