Adaptive Sampling for Final Gather

Personal Information

Name: Razvan Ranca
E-mail address: ranca.razvan@gmail.com
IRC: I'm frequently on #yafaray as Razvan

Idea to implement

Synopsis and Goals:

In most global illumination acceleration techniques, such as irradiance caching, almost all of the work goes into caching the irradiance integrals for areas of high radiance volatility. The compromise needed to be made between speed and quality is obvious, as a method that relies on interpolating nearby samples to obtain the irradiance at a certain point will invariably cause the destruction of small scale details.

The present solution relies on the fact that, due to parallax, when one is moving objects that are farther away seem to move slower than objects nearby. Therefore for any 2 nearby point’s irradiance due to objects farther away will tend to change significantly less than irradiance due to objects nearby. As such it makes sense to apply irradiance caching(or similar techniques) for the objects that are farther from our point than some threshold x and, for the rest, to use a method specially suited to computing irradiance due to nearby objects

By using specific methods for the two different types of objects we can take advantage of each method's strength while minimizing it's cost(avoiding unnecessarily dense storing of irradiance value for volatile areas for instance). Thus we can obtain results very similar to classical irradiance caching while improving speed by an order of magnitude [1]

Detailed description:

We need to compute the radiosity(or diffusely reflected light) at a point P. This can be expressed as :

$B(P) = BSDF(P) \times \int_\Omega \! cos\theta \times L(P,\omega) \, d\omega$

where $\Omega$ is the points' incident hemisphere, $\theta$ is the angle between the points' normal and $\omega$ and $L(P,\omega)$ is the incident radiance at point P from direction $\omega$ .

Our goal is to break this integral into two parts, representing irradiance due to distant and near objects, this is easily achieved like so:

$B(P) = BSDF(P) \times \left( \int_{\Omega_D} \! cos\theta \times L_D(P,\omega) \, d\omega + \int_{\Omega_N} \! cos\theta \times L_N(P,\omega) \, d\omega \right)$

where the new notations have the obvious meanings. Additionally, we shall henceforth refer to the threshold value which defines which objects are near and which are far as $\alpha$

We will notice however that computing the distant part of the hemisphere, $\Omega_D$ , turns out to be more difficult than computing the near part, as such we will rewrite the previous equation like so:

$B(P) = BSDF(P) \times \left( \int_{\Omega} \! cos\theta \times L_D^*(P,\omega) \, d\omega + \int_{\Omega_N}\left( \! cos\theta \times \left( L_N(P,\omega) - L_D^*(P,\omega) \, d\omega\right) \right)\right)$

Here, in the first term, we compute the incoming distant radiance all over the hemisphere, ignoring occlusions due to nearby objects and then, in the second term, besides computing radiance due to the nearby objects we also correct our previous calculation by removing the surplus of distant radiance. For ease of use from now on i shall refer to the first integral in the above equation as $I_D(P)$ and to the second one as $I_N(P)$ .

As mentioned earlier, calculation of $I_D(P)$ at any given shading point relies on a technique similar to irradiance caching, only instead of caching the whole integral( $I_D(P)$ ) we shall instead cache only the incoming radiance( $L_D^*(P,\omega)$ ). This is based on a few different reasons ([1], [3]) including that it allows to more easily compute $I_N(P)$ . Moreover we will store our radiances as just 9 spherical harmonic coefficients (27 actual numbers, 9 for each RGB color), as can be shown ([2]) the average error thus introduced is less than 3%. From these coefficients we directly compute the irradiance on any point by calculating:

$c_1L_{22}(x^2-y^2) + c_3L_{20}z^2 + c_4L_{00} - c_5L_{20} + 2c_1(L_{2-2}xy + L_{21}xz + L_{2-1}yz) + 2c_2(L_{11}x + L_{1-1}y + L_{10}z)$

where x, y, z are the Cartesian coordinates of the surface normal at point P, $c_a$ 's are hard wired constants and $L_{xy}$ 's are the lighting coefficients(the values we store in the radiance map). This calculation needs to be done separately for every RGB color. Thus all we have to do in order to compute $I_D(P)$ is to interpolate the values of nearby radiance samples as to obtain the radiance value at point P and then to apply the aforementioned formula.

Computing $I_N(P)$ is even easier. We shall estimate its value at each shading point by utilizing an aggressive heuristic, namely we consider all back facing triangles and all triangles which are entirely beneath P's hemisphere as not visible and all the rest as visible. It turns out that this surprisingly simple heuristic works quite well for real life scenes, and doesn't perform too badly in contrived scenes either [1]. Thus our equation is:

$I_N(P) = \sum_{t \in T} \left(FF(t,P) \left(\frac {B(t)}{\pi} - L_D^*(P,\omega_t)\right)\right)$

where T is the set of viable triangles, that is triangles that pass our heuristic test and whose center is closer to P than $\alpha$ (the threshold), FF(t,P) is triangle's t form factor with respect to point P, B(t) is it's course radiosity which is computed in a preprocess step and $L_D^*(P,\omega_t)$ is the incoming distant radiance from the direction of triangle t's center. Since all triangles being evaluated are already considered visible the form factor computation becomes trivial:

$FF(t,P) = \frac n {2\pi} \sum_{e\inE} \frac {e \times v_e} {|e \times v_e|}$

where n is the surface normal at point P, E is the set of vectors running along triangle t's edges and $v_e$ is the vector going from point P to the start vertex of vector e. The incoming radiance, $L_D^*(P,\omega_t)$ is once more computed easily from the cached radiance map by first rewriting it in regards to the spherical coordinates on the unit sphere, as $L_D^*(\theta,\phi)$ and then applying the formula

$L_D^*(\theta,\phi) = \sum_{l,m}\left(L_{lm}Y_{lm}(\theta,\phi)\right)$

where (for our approximation) 0<=l<=2 and -l<=m<=l. Here $L_{lm}$ are the stored spherical harmonic coordinated for point P and $Y_{lm}(\theta,\phi)$ is the spherical harmonic function, which, in our case, since we only care about it's first 9 values becomes a polynomial of degree <=2.

All that remains is the preprocess step in which we split triangles which are too large, compute a course radiosity for every triangle and create the actual radiance map. The first two tasks are fairly trivial. Our threshold, $\alpha$ is also the maximum dimension for any triangle, triangles larger than that shall be split. This is because a number of suppositions we make in this approach rely on the assumptions that triangles are small enough to be treated more or less like points (they have a single course radiosity value and so on). The course radiosity is simply computed for the center of the triangle from a course global illumination solution that must have been previously run (e.g.: photon map). The last item is slightly trickier.

Since the whole point of this algorithm is that through the dividing of objects into near and distant we can significantly reduce the number of samples needed to be stored, our samples will obviously be stored sparsely. Also, since we are only interested in distant radiance here, we will disregard any rays originating at our sample point whose intersections are closer to the point than $\alpha$ . However, a problem occurs in that if the samples happen to be poorly placed in a certain region such that many of the rays are disregarded the sample becomes poor, with high error margins. The solution is, instead of shooting a large number of random rays from a single point to determine a radiance sample we shall shoot one ray each from a large number of points and then store their average computed radiance as the sample (this sample will be stored in the average geometric position of the points that contributed to it's creation). This ensures that the sample is representative of the area it is supposed to cover, and minimizes errors due to bad placement. Additionally each such set will attempt to be maximal, so that we don't store needlessly many samples. The specific mechanism being that if any of the rays in the set of points has an intersection within the bounding sphere of the set and the bounding sphere is larger than our threshold $\alpha$ then the set is too large and we divide it into two and repeat the operation recursively. This means that the distance between our samples will be larger than $\alpha$ . Once we have our set, the incoming radiance in the sampled direction for each point in the set is obtained from the course global illumination solution.

The actual mechanics of storing the radiance value as spherical harmonic coordinates merit mentioning as well. The general formula for calculating our coordinates is:

$L_{lm} = \int_{\theta=0}^{\pi} \int_{\phi=0}^{2\pi}\!L(\theta,\phi)Y_{lm}(\theta,\phi)sin\theta \, d\theta d\phi$

where $L(\theta,\phi)$ is the incoming radiance from those coordinates in the unit sphere and $Y_{lm}(\theta,\phi)$ , as was explained earlier, is the spherical harmonic function. In our approximation this integral becomes:

$L_{lm} = \sum_{p \in S} \!L(p)Y_{lm}(\theta,\phi)sin\theta$

where S is the set of points previously mentioned, L(p) is the incoming radiance measured at that point and $\theta$ and $\phi$ are the spherical coordinates corresponding to the sampled direction for point p.

[1] http://www.okanarikan.com/papers/s2005/radiance.php

[2] http://graphics.stanford.edu/papers/envmap/

[3] http://portal.acm.org/citation.cfm?id=1401132.1401228

Workflow and Usability:
Users will be able to specify the threshold value to be used which determines the scarcity of the samples in the radiance cache and the number of rays to trace per sample which influences the quality of the samples in the radiance cache. Then they shall sit back and be amazed by yafaray's new astonishing speed :p

Model Implementation:

-- this is just for clarification of the above detailed description in a more succint manner, and also to give a general idea of the layout, what every part of the algorithm needs and what it returns, that sort of thing.

function: computeDiffuseRad(BSDF bsdf, Point p, Normal n, Double threshold, Triangle[] triangles, RadianceMap radianceMap)

return bsdf.getValue(p)*(distant(radianceMap,n,p) + near(radianceMap, triangles, threshold, n, p))

function: preProcess(threshold, &triangles, sample1, sample2, diffuseMap, &radianceMap, &courseRadiosity, Point[] shadingPoints) // sample 1 and 2 are 2D random Samples

triangles.subdiv(threshold) // subdivide triangles larger than threshold

for(triangle t in triangles)

t.setCourseRadiosity(diffuseMap.computeRadiosity(t.center)) // estimate triangle radiosity as the radiosity at it's center point

Point[][] sets = computeSets(sample1, sample2, threshold, shadingPoints) // compute sets which will create radiance samples, method presented in detailed description

for(Point[] set in sets)

radianceMap.add(computeRadiance(set),computePosition(set)) // again details of implementation in detailed description, above

function: distant(radianceMap,n,p)

return radianceMap.computeRadiosity(p, n)

function: near(radianceMap, triangles, threshold, n, p)

Triangle[] viableTriangles = checkViability(triangles, threshold, n, p) // apply described heuristics, also eliminate triangles which aren't near

Irradiance totalIrad = new Irradience()

for(triangle t in viableTriangles)

totalIrad += formFactor(n,t)*(t.getCourseRadiosity()/pi - radienceMap.computeIncomingRadiosity(p,t.center))

return totalIrad

Temptative schedule:

25.03 - 3.04	Read Physically Based Rendering
04.04 - 07.04	Understood Project idea and wrote proposal
08.04 - 19.05	Figure out details of yafaray integration, pester mentors with questions, read some more on global illumination(+ take uni exams :p)
20.05 - 23.05	Make a detailed implementation plan, taking into account both the algorithm presented above and the yafaray integration
24.05 - 05.06	Implement the distant part of the equation directly using the photon map instead of radiance cache(to test yafaray integration mostly)
06.06 - 11.06	Debug above, I'm sure i'll misunderstand some import detail of the yafaray photon map, should figure it out at this stage
12.06 - 24.06	Implement my preprocessing step
25.06 - 28.06	Modify distant term to use radiance cache, now with a very small threshold value it should act like an (inneficient) irradiance cache
28.06 - 04.07	Debug, this time it'll probably involve pestering some more mentors to tell me if my generated images look like they should(with a bad irrad cache)
05.07 - 22.07	Implement near term
14.07 - 16.07	Submit midterm evaluation
22.07 - 29.07	Final debugging
29.07 - 02.07	Writing documentation
16.08 - 18.08	Final evaluation

-- I do plan to work full time on this project over the summer.

-- Possible vacation time has been accounted for in my schedule(2 weeks)

-- Number of working hours per week is more difficult to estimate, I'm inclined to say as much as are needed. As a maximum i could get to 40 maybe even 45 a week if need be.

Bio

Personal
I'm Razvan, currently an undergraduate student in the informatics department of the University of Edinburgh(1). In my free time I'm usually either solving interesting informatics/math problems(project Euler anyone?(2)), attempting to learn new programming languages, reading webcomics, or doing long distance running. When i have more time i usually enjoy attempting nearly masochistic riddles/puzzles (such as notpron(3)) (with limited success i might add :P).
Programming Background
I'm fairly fluent in C++, Java and Haskell, i also know some C#. About C++, i did 4 years of that(was in an intensive informatics class), so i am fairly familiar with both the language and various algorithms + some STL. Also in my last year of high school i won the regional stage of the informatics Olympiad and qualified to the national one.
Graphics

Last semester a Haskell programming competition was held in which the goal was to make a program that can produce "interesting" images. Looking over submissions from the previous years i noticed some 3d pictures that i found pretty impressive(though to you they surely wouldn't be :lol: ) and discovered ray tracing. During the next 6 days i went from not having heard of the technique in my life to implementing a simple program in haskell that can do the basic stuff + texture mapping, environment cube and procedural textures(full specs(4), pictures (5),(6)). This was all tremendous fun. Then this semester another programming competition was held, this time to make a simple game. The time frame was about 3 weeks and i had no idea where to even start. In the first 10 or so days i managed to learn both Java3D and some game design theory from a couple of (large) books, then in the last 10 days i coded what is more a compilation of 3D game features than an actual game(full specs(7), picture(8)). This was not quite so much fun. I realized that games (even pseudo games like mine) require a lot of repetitive and just annoying tasks to be done besides the interesting ones. Actually the things i enjoyed most turned out to be the terrain generation , the missile trajectory modeling and the particle effects. Maybe not surprisingly these are the tasks that have a somewhat interesting mathematical foundation beneath them.

Why choose me?

I realize that I have pretty poor graphics experience per se, certainly beneath most of the other people that applied for ray tracing ideas, however hopefully i've managed to hint at two things i know. One, when a project interests me i can churn out quite a bit of work in a short period of time(reading and mostly understanding Physically based rendering in less than 10 days for instance), and two ray tracing is exactly the kind of project that interests me, clever programming based on clever math. Additionally the fact that I'm applying for this project(and only this project) in spite of my lack of experience hopefully shows that I'm not so much interested in just getting into gsoc, but more into getting in on a project I'll actually enjoy and thus be able to perform well on.

(1) University of Edinburgh: http://www.ed.ac.uk/home School of informatics: http://www.inf.ed.ac.uk/

(2) Project Euler: http://projecteuler.net/ Don't be fooled by the early problems, they get hard fast :)

(3) Notpron: http://notpron.org/notpron/

(4) My aaaamazing raytracer's specs :P :

Conic perspective
Sphere and Plane detection
Shadow rays
Lambert brightness
good ol' Perfect Reflections
Specular reflections (Phong)
Exponential tone mapping function (RGB clamping looked ugly after the introduction of specular reflections
Texture mapping(interesting math here as well)
Cube environment mapping(very badly implemented as the deadline was looming by the time i got here :) scales horribly)
Procedural textures (noise - type) - found this very interesting going from random numbers to decent looking textures - cool :D

(5) http://i6.photobucket.com/albums/y218/Haggis_McMutton/moreOrbs.jpg -- reflection coefficients set too high in this one, some anti aliasing wouldn't hurt either

(6) http://i6.photobucket.com/albums/y218/Haggis_McMutton/orbs.jpg

(7) My even more aaaamazing game's features :

Environment(on the fly) generated with diamond square algorithm, stored as quads
Mipmap texturing of environment based on quad height
Picking for the environment collision detection divides environment into optimum number of squares(based on environment size) as intermediate step, before checking individual quads, for performance gain
3rd person camera
Laser tube rotate towards selected target
Lasers check if selected target can be hit from characters position (if not intermediate obstacle is hit)
Hermite curve interpolation for the missile trajectory
Very simple enemy AI
Point particle effects
Animated explosion gifs
Minimap(generated on the fly) traces player and objectives location
Player lasers and missiles can overheat(with resulting spectacular!! particle effects :p )
Background music + sounds

(8) http://i6.photobucket.com/albums/y218/Haggis_McMutton/Screenshot-1.png

Attachment	Size
Ray Tracer - Haskell.zip	14.91 KB
3D Game - Java.zip	21.8 KB

Search

Main Menu