Omni-directional light in raytracing program gives wrong render c++

Question

I am trying to implement an omni-directional light source (a.k.a., point light source) in my raytracing program in C++. I am not getting the expected results, but I can't figure out the problem. Maybe someone can see what I am doing wrong. I have included the two functions that are responsible for raytracing and the light. The ClosestIntersection function finds the closest intersection and a triangle. That is used later in the DirectLight function. I would really appreciate any help.

#include <iostream>
#include <glm/glm.hpp>
#include <SDL.h>
#include "SDLauxiliary.h"
#include "TestModel.h"
#include "math.h"

using namespace std;
using glm::vec3;
using glm::mat3;

// ----------------------------------------------------------------------------
// GLOBAL VARIABLES

const int SCREEN_WIDTH = 500;
const int SCREEN_HEIGHT = 500;
SDL_Surface* screen;
int t;
vector<Triangle> triangles;
float focalLength = 900;
vec3 cameraPos(0, 0, -4.5);

vec3 lightPos(0.5, 0.5, 0);
vec3 lightColor = 14.f * vec3(1,1,1);

// Translate camera
float translation = 0.1;        // use this to set translation increment

// Rotate camera
float yaw;
vec3 trueCameraPos;

const float PI = 3.1415927;

// ----------------------------------------------------------------------------
// CLASSES

class Intersection;

// ----------------------------------------------------------------------------
// FUNCTIONS
void Update();
void Draw();
bool ClosestIntersection(vec3 start, vec3 dir, const vector<Triangle>& triangles,
    Intersection& closestIntersection);
vec3 DirectLight(const Intersection& i);
// ----------------------------------------------------------------------------
// STRUCTURES
struct Intersection
{
    vec3 position;
    float distance;
    int triangleIndex;
};

float m = std::numeric_limits<float>::max();

int main(int argc, char* argv[])
{
    LoadTestModel(triangles);

    screen = InitializeSDL(SCREEN_WIDTH, SCREEN_HEIGHT);
    t = SDL_GetTicks(); // Set start value for timer.

    while (NoQuitMessageSDL())
    {
        Update();
        Draw();
    }

    SDL_SaveBMP(screen, "screenshot.bmp");
    return 0;
}

void Update()
{
    // Compute frame time:
    int t2 = SDL_GetTicks();
    float dt = float(t2 - t);
    t = t2;
    cout << "Render time: " << dt << " ms." << endl;
    }
}

void Draw()
{
    if (SDL_MUSTLOCK(screen))
        SDL_LockSurface(screen);

    for (int y = 0; y<SCREEN_HEIGHT; ++y)
    {
for (int x = 0; x < SCREEN_WIDTH; ++x)
{
    vec3 start = cameraPos;
    vec3 dir(x - SCREEN_WIDTH / 2, y - SCREEN_HEIGHT / 2, focalLength);
    Intersection intersection;
    if (ClosestIntersection(start, dir, triangles, intersection))
    {
        //vec3 theColor = triangles[intersection.triangleIndex].color;
        vec3 theColor = DirectLight(intersection);
        PutPixelSDL(screen, x, y, theColor);
    }
    else
    {
        vec3 color(0, 0, 0);
        PutPixelSDL(screen, x, y, color);
    }
}
    }

    if (SDL_MUSTLOCK(screen))
        SDL_UnlockSurface(screen);

    SDL_UpdateRect(screen, 0, 0, 0, 0);
}

bool ClosestIntersection(vec3 s, vec3 d,
    const vector<Triangle>& triangles, Intersection& closestIntersection)
{
    closestIntersection.distance = m;
    for (size_t i = 0; i < triangles.size(); i++)
    {
        vec3 v0 = triangles[i].v0;
        vec3 v1 = triangles[i].v1;
        vec3 v2 = triangles[i].v2;
        vec3 u = v1 - v0;
        vec3 v = v2 - v0;
        vec3 b = s - v0;
        vec3 x;

        // Determinant of A = [-d u v]
        float det = -d.x * ((u.y * v.z) - (v.y * u.z)) -
            u.x * ((-d.y * v.z) - (v.y * -d.z)) +
            v.x * ((-d.y * u.z) - (u.y * -d.z));

        // Cramer'r Rule for t = x.x
        x.x = (b.x * ((u.y * v.z) - (v.y * u.z)) -
            u.x * ((b.y * v.z) - (v.y * b.z)) +
            v.x * ((b.y * u.z) - (u.y * b.z))) / det;

        if (x.x >= 0)
        {
            // Cramer'r Rule for u = x.y
            x.y = (-d.x * ((b.y * v.z) - (v.y * b.z)) -
                b.x * ((-d.y * v.z) - (v.y * -d.z)) +
                v.x * ((-d.y * b.z) - (b.y * -d.z))) / det;

            // Cramer'r Rule for v = x.z
            x.z = (-d.x * ((u.y * b.z) - (b.y * u.z)) -
                u.x * ((-d.y * b.z) - (b.y * -d.z)) +
                b.x * ((-d.y * u.z) - (u.y * -d.z))) / det;

            if (x.y >= 0 && x.z >= 0 && x.y + x.z <= 1 && x.x < closestIntersection.distance)
            {
                closestIntersection.position = x;
                closestIntersection.distance = x.x;
                closestIntersection.triangleIndex = i;
            }
        }

    }
    //end of for loop

    if (closestIntersection.distance != m)
    {
        return true;
    }
    else
    {
        return false;
    }

}

vec3 DirectLight(const Intersection& i)
{
    vec3 n = triangles[i.triangleIndex].normal;
    vec3 r = lightPos - i.position;
    float R2 = r.x * r.x + r.y * r.y + r.z * r.z;
    vec3 D = (lightColor * fmaxf((glm::dot(glm::normalize(r), n)), 0)) / (4 * PI * R2);
    return D;
}

Answer 1

This isn't a super-great answer, but I managed to make your code work without the strange shading discontinuities. The problem happens in ClosestIntersection and maybe Gareth's answer covers it. I need to stop looking at this now, but I wanted to show you what I have before I leave, and I need an Answer to post some code.

// This starts with some vec3 helper functions which make things
// easier to look at
float Dot(const vec3& a, const vec3& b) {
    return a.x * b.x + a.y * b.y + a.z * b.z;
}
vec3 Cross(const vec3& a, const vec3& b) {
   return vec3(a.y*b.z - a.z*b.y, a.z*b.x - a.x*b.z, a.x*b.y - a.y*b.x); 
}
float L2(const vec3& v) { return v.x*v.x + v.y*v.y + v.z*v.z; }
float Abs(const vec3& v) { return std::sqrt(L2(v)); }

// Here is the replacement version of ClosestIntersection
bool ClosestIntersection(vec3 cam, vec3 dir,
    const vector<Triangle>& triangles, Intersection& closestIntersection) 
{
    closestIntersection.distance = m;
    vec3 P0 = cam;
    vec3 P1 = cam + dir;
    for (size_t i = 0; i < triangles.size(); ++i) {
        vec3 v0 = triangles[i].v0;
        vec3 v1 = triangles[i].v1;
        vec3 v2 = triangles[i].v2;
        // Dan Sunday
        // http://geomalgorithms.com/a06-_intersect-2.html
        vec3 u = v1 - v0;
        vec3 v = v2 - v0;
        // w = P-v0, solve w = su +tv (s, t are parametric scalars)

        vec3 n = Cross(u, v);
        float ri = Dot(n, (v0 - P0)) / Dot(n, (P1 - P0));
        vec3 Pi = P0 + ri * (P1- P0);
        vec3 w = Pi - v0;

        // s = w . (n x v) / (u . (n x v))
        // t = w . (n x u) / (v . (n x u))

        float s = Dot(w, Cross(n, v)) / Dot(u, Cross(n, v));
        float t = Dot(w, Cross(n, u)) / Dot(v, Cross(n, u));

        if(s >= 0 && t >= 0 && s+t <= 1) {
            float dist = Abs(cam - Pi);
            if(dist < closestIntersection.distance) {
                closestIntersection.position      = Pi;
                closestIntersection.distance      = dist;
                closestIntersection.triangleIndex = int(i);
            }
        }
    }

    return closestIntersection.distance != m;
}

Good luck.

Answer 2

If I'm understanding the code in ClosestIntersection correctly, here's what it's doing for each triangle:

• Let u,v be the vectors from one vertex of the triangle to the other two vertices. Let d be (the reverse of) the direction of the ray we're considering.
• And let b be the vector from that vertex of the triangle to the camera.
• Find p,q,r so that b = pd+qu+rv (p,q,r are what your code calls x.x, x.y, x.z).
• Now the ray meets the triangle if p>0, q>=0, r>=0, q+r<=1 and the distance to the intersection point is p.

So, the conditions on q,r make sense; the idea is that b-qu-rv is the vector from the camera to the relevant point in the triangle and it's in direction d. Your distances aren't really distances, but along a single ray they're the same multiple of the actual distance, which means that this works fine for determining which triangle you've hit, and that's all you use them for. So far, so good.

But then you say closestIntersection.position = x; and surely that's all wrong, because this x isn't in the same coordinate system as your camera location, triangle vertices, etc. It's in this funny "how much of d, how much of u, how much of v" coordinate system which isn't even the same from one triangle to the next. (Which is why you are getting discontinuities at triangle boundaries even within a single face, I think.)

Try setting it to v0+x.y*(v1-v0)+x.z*(v2-v0) instead (I think this is right; it's meant to be the actual point where the ray crosses the triangle, in the same coordinates as all your other points) and see what it does.