sphere and disk area lights mostly ready (I think :S).

Makefile

@@ -1,6 +1,9 @@
 CXX = g++
 TARGET = ray
-OBJECTS = main.o sampling.o camera.o environment.o disk.o plane.o sphere.o phong_brdf.o hsa_brdf.o directional_light.o point_light.o spot_light.o scene.o tracer.o path_tracer.o whitted_tracer.o
+OBJECTS = main.o sampling.o camera.o environment.o disk.o plane.o sphere.o \
+          phong_brdf.o hsa_brdf.o directional_light.o point_light.o \
+          spot_light.o sphere_area_light.o disk_area_light.o scene.o tracer.o \
+          path_tracer.o whitted_tracer.o
 DEPENDS = $(OBJECTS:.o=.d)
 CXXFLAGS = -ansi -pedantic -Wall -DGLM_FORCE_RADIANS -fopenmp
 LDLIBS = -lfreeimage -ljson_spirit

TODO.org

@@ -12,8 +12,8 @@
  - [X] Directional lights
  - [X] Point lights
  - [X] Spotlights
- - [ ] Area lights
-   - [ ] Sphere lights
+ - [-] Area lights
+   - [X] Sphere lights
    - [ ] Box/planar lights
  - [X] Phong shading
  - [X] Specular reflections

directional_light.cpp

@@ -12,11 +12,11 @@ using glm::dot;
 using glm::pow;
 using glm::max;
 
-inline vec3 DirectionalLight::direction(vec3 point) {
+vec3 DirectionalLight::direction(vec3 point) const {
   return m_position;
 }
 
-inline float DirectionalLight::distance(vec3 point) {
+float DirectionalLight::distance(vec3 point) const {
   return numeric_limits<float>::max();
 }
 

directional_light.hpp

@@ -14,8 +14,8 @@ public:
 
   DirectionalLight(vec3 _p, vec3 _d, vec3 _s): InfinitesimalLight(normalize(_p), _d, _s) { }
 
-  virtual vec3 direction(vec3 point);
-  virtual float distance(vec3 point);
+  virtual vec3 direction(vec3 point) const;
+  virtual float distance(vec3 point) const;
   virtual vec3 diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
   virtual vec3 specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
 };

disk.cpp

@@ -1,6 +1,12 @@
-#include "disk.hpp"
+#include <glm/gtc/constants.hpp>
 
+#include "disk.hpp"
+#include "sampling.hpp"
+
+using glm::cos;
+using glm::sin;
 using glm::dot;
+using glm::pi;
 
 bool Disk::intersect(Ray & r, float & t) const {
   float _t;
@@ -15,3 +21,18 @@ bool Disk::intersect(Ray & r, float & t) const {
 
    return false;
 }
+
+vec3 Disk::sample_at_surface() const {
+  float theta = random01() * (2.0f * pi<float>());
+  float r = glm::sqrt(random01() * m_radius);
+  float x = r * cos(theta);
+  float y = r * sin(theta);
+  float z = 0.0f;
+  vec3 sample = vec3(x, y, z);
+  rotate_sample(sample, m_normal);
+  return sample;
+}
+
+void Disk::calculate_inv_area() {
+  m_inv_area = 1.0f / pi<float>() * (m_radius * m_radius);
+}

disk.hpp

@@ -13,15 +13,25 @@ class Disk : public Plane {
 public:
   float m_radius;
   
-  Disk(Material * mat = NULL): Plane(mat), m_radius(1.0f) { }
+  Disk(Material * mat = NULL): Plane(mat), m_radius(1.0f) {
+    calculate_inv_area();
+  }
 
-  Disk(float x, float y, float z, float nx, float ny, float nz, float _r, Material * mat = NULL): Plane(x, y, z, nx, ny, nz, mat), m_radius(_r) { }
+  Disk(float x, float y, float z, float nx, float ny, float nz, float _r, Material * mat = NULL): Plane(x, y, z, nx, ny, nz, mat), m_radius(_r) {
+    calculate_inv_area();
+  }
 
-  Disk(vec3 _p, vec3 _n, float _r, Material * mat = NULL): Plane(_p, _n, mat), m_radius(_r) { }
+  Disk(vec3 _p, vec3 _n, float _r, Material * mat = NULL): Plane(_p, _n, mat), m_radius(_r) {
+    calculate_inv_area();
+  }
 
   virtual ~Disk() { }
 
   virtual bool intersect(Ray & r, float & t) const;
+  virtual vec3 sample_at_surface() const;
+
+protected:
+  virtual void calculate_inv_area();
 };
 
 

disk_area_light.cpp

@@ -0,0 +1,37 @@
+#include <iostream>
+#include <cassert>
+
+#include "disk_area_light.hpp"
+#include "ray.hpp"
+
+using glm::normalize;
+
+const float BIAS = 0.000001f;
+
+vec3 DiskAreaLight::diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const {
+  float d, att;
+  vec3 l_dir, ref;
+
+  l_dir = normalize(direction(i_pos));
+  d = distance(i_pos);
+  att = 1.0f / (m_const_att + (m_lin_att * d) + (m_quad_att * (d * d)));
+
+  return (att * m.m_brdf->diffuse(l_dir, normal, r, i_pos, m_diffuse)) / m_figure->pdf();
+}
+
+vec3 DiskAreaLight::specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const {
+float d, att;
+  vec3 l_dir, ref;
+
+  l_dir = normalize(direction(i_pos));
+  d = distance(i_pos);
+  att = 1.0f / (m_const_att + (m_lin_att * d) + (m_quad_att * (d * d)));
+
+  return (att * m.m_brdf->specular(l_dir, normal, r, i_pos, m_specular, m.m_shininess)) / m_figure->pdf();
+}
+
+void DiskAreaLight::sample_at_surface(vec3 point) {
+  Disk * d = static_cast<Disk *>(m_figure);
+  m_last_sample = m_figure->sample_at_surface();
+  m_n_at_last_sample = d->m_normal;
+}

disk_area_light.hpp

@@ -0,0 +1,26 @@
+#pragma once
+#ifndef DISK_AREA_LIGHT_HPP
+#define DISK_AREA_LIGHT_HPP
+
+#include "light.hpp"
+#include "disk.hpp"
+
+class DiskAreaLight: public AreaLight {
+public:
+  float m_const_att;
+  float m_lin_att;
+  float m_quad_att;
+
+  DiskAreaLight(Disk * _s, float _c = 1.0, float _l = 0.0, float _q = 0.0):
+    AreaLight(static_cast<Figure *>(_s)),
+    m_const_att(_c),
+    m_lin_att(_l),
+    m_quad_att(_q)
+  { }
+
+  virtual vec3 diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
+  virtual vec3 specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
+  virtual void sample_at_surface(vec3 point);
+};
+
+#endif

figure.hpp

@@ -13,7 +13,7 @@ class Figure {
 public:
   Material * m_mat;
 
-  Figure(Material * mat = NULL) {
+  Figure(Material * mat = NULL): m_inv_area(0.0f) {
     m_mat = mat == NULL ? new Material() : mat;
   }
   
@@ -21,8 +21,18 @@ public:
     delete m_mat;
   }
 
+  virtual float pdf() const {
+    return m_inv_area;
+  }
+
   virtual bool intersect(Ray & r, float & t) const = 0;
   virtual vec3 normal_at_int(Ray & r, float & t) const = 0;
+  virtual vec3 sample_at_surface() const = 0;
+
+protected:
+  float m_inv_area;
+
+  virtual void calculate_inv_area() = 0;
 };
 
 #endif

light.hpp

@@ -2,12 +2,14 @@
 #ifndef LIGHT_HPP
 #define LIGHT_HPP
 
-#include <glm/vec3.hpp>
+#include <glm/glm.hpp>
 
+#include "figure.hpp"
 #include "material.hpp"
 #include "ray.hpp"
 
 using glm::vec3;
+using glm::length;
 
 class Light {
 public:
@@ -19,53 +21,60 @@ public:
 
   Light(): m_position(vec3(0.0f)), m_diffuse(vec3(1.0f)), m_specular(vec3(1.0f)) { }
 
-  Light(vec3 _p, vec3 _d, vec3 _s): m_position(_p), m_diffuse(_d), m_specular(_s) { }
+  Light(ltype_t _t): m_position(vec3(0.0f)), m_diffuse(vec3(1.0f)), m_specular(vec3(1.0f)), m_type(_t) { }
+
+  Light(vec3 _p, vec3 _d, vec3 _s, ltype_t _t): m_position(_p), m_diffuse(_d), m_specular(_s), m_type(_t) { }
   
   virtual ~Light() { }
 
   virtual ltype_t light_type() {
-    return type;
+    return m_type;
   }
   
-  virtual vec3 direction(vec3 point) = 0;
-  virtual float distance(vec3 point) = 0;
+  virtual vec3 direction(vec3 point) const = 0;
+  virtual float distance(vec3 point) const = 0;
   virtual vec3 diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const = 0;
   virtual vec3 specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const = 0;
 
 protected:
-  ltype_t type;
+  ltype_t m_type;
 };
 
 class InfinitesimalLight: public Light {
 public:
-  InfinitesimalLight(): Light() {
-    type = INFINITESIMAL;
-  }
+  InfinitesimalLight(): Light(INFINITESIMAL) { }
   
-  InfinitesimalLight(vec3 _p, vec3 _d, vec3 _s): Light(_p, _d, _s) {
-    type = INFINITESIMAL;
-  }
+  InfinitesimalLight(vec3 _p, vec3 _d, vec3 _s): Light(_p, _d, _s, INFINITESIMAL) { }
   
-  virtual vec3 direction(vec3 point) = 0;
-  virtual float distance(vec3 point) = 0;
+  virtual vec3 direction(vec3 point) const = 0;
+  virtual float distance(vec3 point) const = 0;
   virtual vec3 diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const = 0;
   virtual vec3 specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const = 0;
 };
 
 class AreaLight: public Light {
 public:
-  AreaLight(): Light() {
-    type = AREA;
+  Figure * m_figure;
+
+  AreaLight(): Light(AREA), m_figure(NULL), m_last_sample(vec3()), m_n_at_last_sample(vec3()) { }
+
+  AreaLight(Figure * _f): Light(AREA), m_figure(_f), m_last_sample(vec3()), m_n_at_last_sample(vec3()) { }
+
+  virtual vec3 direction(vec3 point) const {
+    return m_last_sample - point;
   }
 
-  AreaLight(vec3 _p, vec3 _d, vec3 _s): Light(_p, _d, _s) {
-    type = AREA;
+  virtual float distance(vec3 point) const {
+    return length(m_last_sample - point);
   }
 
-  virtual vec3 direction(vec3 point) = 0;
-  virtual float distance(vec3 point) = 0;
   virtual vec3 diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const = 0;
   virtual vec3 specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const = 0;
+  virtual void sample_at_surface(vec3 point) = 0;
+
+protected:
+  vec3 m_last_sample;
+  vec3 m_n_at_last_sample;
 };
 
 #endif

path_tracer.cpp

@@ -15,8 +15,9 @@ vec3 PathTracer::trace_ray(Ray & r, Scene * s, unsigned int rec_level) const {
   Figure * _f;
   vec3 n, color, i_pos, ref, sample, dir_diff_color, dir_spec_color, ind_color, amb_color;
   Ray mv_r, sr, rr;
-  bool vis;
+  bool vis, is_area_light = false;
   float kr, r1, r2;
+  AreaLight * al;
 
   t = numeric_limits<float>::max();
   _f = NULL;
@@ -35,8 +36,19 @@ vec3 PathTracer::trace_ray(Ray & r, Scene * s, unsigned int rec_level) const {
     i_pos = r.m_origin + (t * r.m_direction);
     n = _f->normal_at_int(r, t);
 
+    is_area_light = false;
+    // Check if the object is an area light;
+    for (vector<Light *>::iterator it = s->m_lights.begin(); it != s->m_lights.end(); it++) {
+      if ((*it)->light_type() == Light::AREA && static_cast<AreaLight *>(*it)->m_figure == _f)
+	is_area_light = true;
+    }
+
+    // If the object is an area light, return it's emission value.
+    if (is_area_light) {
+      return _f->m_mat->m_emission;
+
     // Check if the material is not reflective/refractive.
-    if (!_f->m_mat->m_refract) {
+    } else if (!_f->m_mat->m_refract) {
       // Calculate the direct lighting.
       for (size_t l = 0; l < s->m_lights.size(); l++) {
 	// For every light source
@@ -44,21 +56,39 @@ vec3 PathTracer::trace_ray(Ray & r, Scene * s, unsigned int rec_level) const {
 
 	if (s->m_lights[l]->light_type() == Light::INFINITESIMAL) {
 	  // Cast a shadow ray to determine visibility.
-	  sr = Ray(s->m_lights[l]->direction(i_pos), i_pos + n * BIAS);
+	  sr = Ray(s->m_lights[l]->direction(i_pos), i_pos + (n * BIAS));
+
 	  for (size_t f = 0; f < s->m_figures.size(); f++) {
 	    if (s->m_figures[f]->intersect(sr, _t) && _t < s->m_lights[l]->distance(i_pos)) {
 	      vis = false;
 	      break;
 	    }
 	  }
-	} else if (s->m_lights[l]->light_type() == Light::AREA) {
-	  // Area lights not supported with Whitted ray tracing.
-	  vis = false;
-	}
 
 	// Evaluate the shading model accounting for visibility.
 	dir_diff_color += vis ? s->m_lights[l]->diffuse(n, r, i_pos, *_f->m_mat) : vec3(0.0f);
 	dir_spec_color += vis ? s->m_lights[l]->specular(n, r, i_pos, *_f->m_mat) : vec3(0.0f);
+
+	} else if (s->m_lights[l]->light_type() == Light::AREA) {
+	  // Cast a shadow ray towards a sample point on the surface of the light source.
+	  al = static_cast<AreaLight *>(s->m_lights[l]);
+	  al->sample_at_surface(i_pos);
+	  sr = Ray(al->direction(i_pos), i_pos + (n * BIAS));
+
+	  for (size_t f = 0; f < s->m_figures.size(); f++) {
+	    // Avoid self-intersection with the light source.
+	    if (al->m_figure != s->m_figures[f]) {
+	      if (s->m_figures[f]->intersect(sr, _t) && _t < al->distance(i_pos)) {
+		vis = false;
+		break;
+	      }
+	    }
+	  }
+
+	  // Evaluate the shading model accounting for visibility.
+	  dir_diff_color += vis ? s->m_lights[l]->diffuse(n, r, i_pos, *_f->m_mat) : vec3(0.0f);
+	  dir_spec_color += vis ? s->m_lights[l]->specular(n, r, i_pos, *_f->m_mat) : vec3(0.0f);
+	}
       }
 
       // Calculate indirect lighting contribution.

plane.cpp

@@ -19,3 +19,11 @@ bool Plane::intersect(Ray & r, float & t) const {
 vec3 Plane::normal_at_int(Ray & r, float & t) const {
   return vec3(m_normal);
 }
+
+vec3 Plane::sample_at_surface() const {
+  return vec3(0.0f);
+}
+
+void Plane::calculate_inv_area() {
+  m_inv_area = 0.0f;
+}

plane.hpp

@@ -14,17 +14,28 @@ public:
   vec3 m_point;
   vec3 m_normal;
   
-  Plane(Material * mat = NULL): Figure(mat), m_point(vec3(0.0f)), m_normal(vec3(0.0f, 0.0f, 1.0f)) { }
+  Plane(Material * mat = NULL): Figure(mat), m_point(vec3(0.0f)), m_normal(vec3(0.0f, 0.0f, 1.0f)) {
+    calculate_inv_area();
+  }
 
-  Plane(float x, float y, float z, float nx, float ny, float nz, Material * mat = NULL): Figure(mat), m_point(vec3(x, y, z)), m_normal(normalize(vec3(nx, ny, nz))) { }
+  Plane(float x, float y, float z, float nx, float ny, float nz, Material * mat = NULL):
+    Figure(mat),
+    m_point(vec3(x, y, z)),
+    m_normal(normalize(vec3(nx, ny, nz)))
+  {
+    calculate_inv_area();
+  }
 
   Plane(vec3 _p, vec3 _n, Material * mat = NULL): Figure(mat), m_point(_p), m_normal(normalize(_n)) { }
 
   virtual ~Plane() { }
 
   virtual bool intersect(Ray & r, float & t) const;
-
   virtual vec3 normal_at_int(Ray & r, float & t) const;
+  virtual vec3 sample_at_surface() const;
+
+protected:
+  virtual void calculate_inv_area();
 };
 
 

point_light.cpp

@@ -11,11 +11,11 @@ using glm::dot;
 using glm::pow;
 using glm::max;
 
-inline vec3 PointLight::direction(vec3 point) {
+vec3 PointLight::direction(vec3 point) const {
   return normalize(m_position - point);
 }
 
-inline float PointLight::distance(vec3 point) {
+float PointLight::distance(vec3 point) const {
   return length(m_position - point);
 }
 
@@ -23,9 +23,8 @@ vec3 PointLight::diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const {
   float d, att;
   vec3 l_dir, ref;
 
-  l_dir = m_position - i_pos;
-  d = length(l_dir);
-  l_dir = normalize(l_dir);
+  l_dir = normalize(direction(i_pos));
+  d = distance(i_pos);
   att = 1.0f / (m_const_att + (m_lin_att * d) + (m_quad_att * (d * d)));
 
   return att * m.m_brdf->diffuse(l_dir, normal, r, i_pos, m_diffuse);
@@ -35,9 +34,8 @@ vec3 PointLight::specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const
   float d, att;
   vec3 l_dir, ref;
 
-  l_dir = m_position - i_pos;
-  d = length(l_dir);
-  l_dir = normalize(l_dir);
+  l_dir = normalize(direction(i_pos));
+  d = distance(i_pos);
   att = 1.0f / (m_const_att + (m_lin_att * d) + (m_quad_att * (d * d)));
 
   return att * m.m_brdf->specular(l_dir, normal, r, i_pos, m_specular, m.m_shininess);

point_light.hpp

@@ -14,8 +14,8 @@ public:
 
   PointLight(vec3 _p, vec3 _d, vec3 _s, float _c, float _l, float _q): InfinitesimalLight(_p, _d, _s), m_const_att(_c), m_lin_att(_l), m_quad_att(_q) { }
 
-  virtual vec3 direction(vec3 point);
-  virtual float distance(vec3 point);
+  virtual vec3 direction(vec3 point) const;
+  virtual float distance(vec3 point) const;
   virtual vec3 diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
   virtual vec3 specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
 };

sampling.cpp

@@ -14,8 +14,8 @@ using glm::pi;
 
 const float PDF = (1.0f / (2.0f * pi<float>()));
 
-float random01() {
-  return static_cast<float>(rand()) / static_cast<float>(RAND_MAX);
+inline float random01() {
+  return static_cast<float>(rand() % 1024) / 1025.0f;
 }
 
 vec2 sample_pixel(int i, int j, float w, float h, float a_ratio, float fov) {
@@ -49,7 +49,7 @@ vec3 sample_hemisphere(const float r1, float r2) {
   return vec3(x, r1, z);
 }
 
-void rotate_sample(vec3 & sample, vec3 & n) {
+void rotate_sample(vec3 & sample, const vec3 & n) {
   vec3 nt, nb;
   mat3 rot_m;
 

sampling.hpp

@@ -14,6 +14,6 @@ extern float random01();
 extern vec2 sample_pixel(int i, int j, float w, float h, float a_ratio, float fov);
 extern void create_coords_system(const vec3 &n, vec3 &nt, vec3 &nb);
 extern vec3 sample_hemisphere(const float r1, float r2);
-extern void rotate_sample(vec3 & sample, vec3 & n);
+extern void rotate_sample(vec3 & sample, const vec3 & n);
 
 #endif

scene.cpp

@@ -18,6 +18,8 @@
 #include "directional_light.hpp"
 #include "point_light.hpp"
 #include "spot_light.hpp"
+#include "sphere_area_light.hpp"
+#include "disk_area_light.hpp"
 
 using std::cerr;
 using std::endl;
@@ -44,6 +46,7 @@ static const string DLT_KEY = "directional_light";
 static const string PLT_KEY = "point_light";
 static const string SLT_KEY = "spot_light";
 static const string SAL_KEY = "sphere_area_light";
+static const string DAL_KEY = "disk_area_light";
 static const string PAL_KEY = "planar_area_light";
 
 static const string ENV_TEX_KEY = "texture";
@@ -137,6 +140,12 @@ Scene::Scene(const char * file_name, int h, int w, float fov) {
 	else if ((*it).name_ == SLT_KEY)
 	  m_lights.push_back(read_light((*it).value_, SPOT));
 
+	else if((*it).name_ == SAL_KEY)
+	  m_lights.push_back(read_area_light((*it).value_, SPHERE));
+
+	else if((*it).name_ == DAL_KEY)
+	  m_lights.push_back(read_area_light((*it).value_, DISK));
+
 	else
 	  cerr << "Unrecognized key \"" << (*it).name_ << "\" in the input file." << endl;
       }
@@ -475,5 +484,31 @@ Light * Scene::read_light(Value & v, light_type_t t) {
     return static_cast<Light *>(new SpotLight(position, diffuse, specular, const_att, lin_att, quad_att, spot_cutoff, spot_exp, spot_dir));
   else
     throw SceneError("Unknown light type.");
+
   return NULL;
 }
+
+Light * Scene::read_area_light(Value & v, light_type_t t) {
+  Disk * d;
+  Sphere * s;
+  DiskAreaLight * dal;
+  SphereAreaLight * sal;
+  Light * l;
+
+  if (t == SPHERE) {
+    s = static_cast<Sphere *>(read_sphere(v));
+    sal = new SphereAreaLight(s);
+    m_figures.push_back(static_cast<Figure *>(s));
+    l = static_cast<Light *>(sal);
+
+  } else if (t == DISK) {
+    d = static_cast<Disk *>(read_disk(v));
+    dal = new DiskAreaLight(d);
+    m_figures.push_back(static_cast<Figure *>(d));
+    l = static_cast<Light *>(dal);
+
+  } else
+    throw SceneError("Unknown area light type");
+
+  return l;
+}

scene.hpp

@@ -19,8 +19,6 @@ using std::vector;
 using std::runtime_error;
 using json_spirit::Value;
 
-typedef enum LIGHT_TYPE { DIRECTIONAL, POINT, SPOT } light_type_t;
-
 class SceneError: public runtime_error {
 public:
   explicit SceneError(const string & what_arg): runtime_error(what_arg) { }
@@ -28,6 +26,8 @@ public:
 
 class Scene {
 public:
+  typedef enum LIGHT_TYPE { DIRECTIONAL, POINT, SPOT, SPHERE, DISK, PLANE } light_type_t;
+
   vector<Figure *> m_figures;
   vector<Light *> m_lights;
   Environment * m_env;
@@ -45,6 +45,7 @@ private:
   Figure * read_plane(Value & v);
   Figure * read_disk(Value & v);
   Light * read_light(Value & v, light_type_t t);
+  Light * read_area_light(Value & v, light_type_t t);
 };
 
 #endif

scenes/scene3.json

@@ -30,7 +30,7 @@
 	}
     },
 
-    "sphere":{
+    "sphere_area_light":{
 	"position": [1.0, 0.0, -3.25],
 	"radius": 1.5,
 	"material": {

scenes/scene5.json

@@ -7,15 +7,17 @@
 	"position": [0.0, 0.0, -2.0],
 	"radius": 0.5,
 	"material": {
-	    "diffuse": [0.5, 0.5, 0.5]
+	    "diffuse": [1.0, 1.0, 0.0],
+	    "transmissive": true,
+	    "ref_index": 1.33
 	}
     },
 
-    "sphere": {
-	"position": [0.0, 2.0, -2.0],
-	"radius": 0.25,
+    "sphere_area_light": {
+	"position": [0.0, 1.0, -2.0],
+	"radius": 0.15,
 	"material": {
-	    "emission": [10.0, 10.0, 10.0]
+	    "emission": [1.0, 1.0, 1.0]
 	}
     },
 

scenes/scene6.json

@@ -0,0 +1,90 @@
+{
+    "sphere_area_light": {
+	"position": [0.0, 0.75, -1.0],
+	"radius": 0.15,
+	"material": {
+	    "emission": [1.0, 1.0, 1.0]
+	}
+    },
+    
+    "sphere": {
+	"position": [0.2, 0.0, -0.75],
+	"radius": 0.25,
+	"material": {
+	    "diffuse": [1.0, 1.0, 1.0]
+	}
+    },
+
+    "sphere": {
+	"position": [-0.5, -0.5, -1.5],
+	"radius": 0.5,
+	"material": {
+	    "diffuse": [0.0, 0.0, 0.0],
+	    "rho": 1.0
+	}
+    },
+
+    "sphere": {
+	"position": [-0.5, -0.5, 0.6],
+	"radius": 0.5,
+	"material": {
+	    "diffuse": [1.0, 1.0, 0.0],
+	    "transmissive": true,
+	    "ref_index": 1.33
+	}
+    },
+
+    "plane": {
+	"position": [0.0, -1.0, 0.0],
+	"normal": [0.0, 1.0, 0.0],
+	"material": {
+	    "diffuse": [0.0, 1.0, 0.0],
+	    "specular": [0.0, 0.0, 0.0]
+	}
+    },
+
+    "plane": {
+	"position": [-2.0, 0.0, 0.0],
+	"normal": [1.0, 0.0, 0.0],
+	"material": {
+	    "diffuse": [1.0, 0.0, 0.0],
+	    "specular": [0.0, 0.0, 0.0]
+	}
+    },
+
+    "plane": {
+	"position": [2.0, 0.0, 0.0],
+	"normal": [-1.0, 0.0, 0.0],
+	"material": {
+	    "diffuse": [0.0, 0.0, 1.0],
+	    "specular": [0.0, 0.0, 0.0]
+	}
+    },
+
+    "plane": {
+	"position": [0.0, 1.0, 0.0],
+	"normal": [0.0, -1.0, 0.0],
+	"material": {
+	    "diffuse": [0.0, 1.0, 1.0],
+	    "specular": [0.0, 0.0, 0.0]
+	}
+    },
+
+    "plane": {
+	"position": [0.0, 0.0, -2.0],
+	"normal": [0.0, 0.0, 1.0],
+	"material": {
+	    "diffuse": [1.0, 0.0, 1.0],
+	    "specular": [0.0, 0.0, 0.0]
+	}
+    },
+
+    "plane": {
+	"position": [0.0, 0.0, 1.1],
+	"normal": [0.0, 0.0, -1.0],
+	"material": {
+	    "diffuse": [1.0, 1.0, 0.0],
+	    "specular": [0.0, 0.0, 0.0]
+	}
+    }
+}

sphere.cpp

@@ -1,4 +1,8 @@
+#include <iostream>
+#include <glm/gtc/constants.hpp>
+
 #include "sphere.hpp"
+#include "sampling.hpp"
 
 using namespace glm;
 
@@ -38,3 +42,23 @@ vec3 Sphere::normal_at_int(Ray & r, float & t) const {
   vec3 i = vec3(r.m_origin + (t * r.m_direction));
   return normalize(vec3((i - m_center) / m_radius));
 }
+
+vec3 Sphere::sample_at_surface() const {
+  float theta;
+  float u, sqrt1muu, x, y, z;
+
+  // Sampling formula from Wolfram Mathworld:
+  // http://mathworld.wolfram.com/SpherePointPicking.html
+  theta = random01()* (2.0f * pi<float>());
+  u = (random01() * 2.0f) - 1.0f;
+  sqrt1muu = glm::sqrt(1.0f - (u * u));
+  x = m_radius * sqrt1muu * cos(theta);
+  y = m_radius * sqrt1muu * sin(theta);
+  z = m_radius * u;
+
+  return vec3(x, y, z) + m_center;
+}
+
+void Sphere::calculate_inv_area() {
+  m_inv_area = 1.0f / (4.0 * pi<float>() * (m_radius * m_radius));
+}

sphere.hpp

@@ -13,17 +13,26 @@ public:
   vec3 m_center;
   float m_radius;
   
-  Sphere(Material * mat = NULL): Figure(mat), m_center(vec3(0.0f)), m_radius(0.5f) { }
+  Sphere(Material * mat = NULL): Figure(mat), m_center(vec3(0.0f)), m_radius(0.5f) {
+    calculate_inv_area();
+  }
 
-  Sphere(float x, float y, float z, float r, Material * mat = NULL): Figure(mat), m_center(vec3(x, y, z)), m_radius(r) { }
+  Sphere(float x, float y, float z, float r, Material * mat = NULL): Figure(mat), m_center(vec3(x, y, z)), m_radius(r) {
+    calculate_inv_area();
+  }
 
-  Sphere(vec3 _c, float r, Material * mat = NULL): Figure(mat), m_center(_c), m_radius(r) { }
+  Sphere(vec3 _c, float r, Material * mat = NULL): Figure(mat), m_center(_c), m_radius(r) {
+    calculate_inv_area();
+  }
 
   virtual ~Sphere() { }
 
   virtual bool intersect(Ray & r, float & t) const;
-
   virtual vec3 normal_at_int(Ray & r, float & t) const;
+  virtual vec3 sample_at_surface() const;
+
+private:
+  virtual void calculate_inv_area();
 };
 
 #endif

sphere_area_light.cpp

@@ -0,0 +1,45 @@
+#include "sphere_area_light.hpp"
+#include "ray.hpp"
+
+using glm::max;
+using glm::dot;
+using glm::normalize;
+
+vec3 SphereAreaLight::diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const {
+  float d, att, ln_dot_d, d2, g;
+  vec3 l_dir, ref;
+
+  l_dir = normalize(direction(i_pos));
+  ln_dot_d = dot(-m_n_at_last_sample, l_dir);
+  if (ln_dot_d > 0.0f) {
+    d2 = glm::distance(m_last_sample, i_pos);
+    d2 *= d2;
+    g = ln_dot_d / d2;
+    d = distance(i_pos);
+    att = 1.0f / (m_const_att + (m_lin_att * d) + (m_quad_att * (d * d)));
+    return (att * m.m_brdf->diffuse(l_dir, normal, r, i_pos, m_diffuse) * g) / m_figure->pdf();
+
+  } else
+    return vec3(0.0f);
+}
+
+vec3 SphereAreaLight::specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const {
+  float d, att, ln_dot_d;
+  vec3 l_dir, ref;
+
+  l_dir = normalize(direction(i_pos));
+  ln_dot_d = dot(-m_n_at_last_sample, l_dir);
+  if (ln_dot_d > 0.0f) {
+    d = distance(i_pos);
+    att = 1.0f / (m_const_att + (m_lin_att * d) + (m_quad_att * (d * d)));
+    return (att * m.m_brdf->specular(l_dir, normal, r, i_pos, m_specular, m.m_shininess)) / m_figure->pdf();
+
+  } else
+    return vec3(0.0f);
+}
+
+void SphereAreaLight::sample_at_surface(vec3 point) {
+  Sphere * s = static_cast<Sphere *>(m_figure);
+  m_last_sample = m_figure->sample_at_surface();
+  m_n_at_last_sample = normalize(vec3((m_last_sample - s->m_center) / s->m_radius));
+}

sphere_area_light.hpp

@@ -0,0 +1,26 @@
+#pragma once
+#ifndef SPHERE_AREA_LIGHT_HPP
+#define SPHERE_AREA_LIGHT_HPP
+
+#include "light.hpp"
+#include "sphere.hpp"
+
+class SphereAreaLight: public AreaLight {
+public:
+  float m_const_att;
+  float m_lin_att;
+  float m_quad_att;
+
+  SphereAreaLight(Sphere * _s, float _c = 1.0, float _l = 0.0, float _q = 0.0):
+    AreaLight(static_cast<Figure *>(_s)),
+    m_const_att(_c),
+    m_lin_att(_l),
+    m_quad_att(_q)
+  { }
+
+  virtual vec3 diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
+  virtual vec3 specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
+  virtual void sample_at_surface(vec3 point);
+};
+
+#endif

spot_light.cpp

@@ -17,9 +17,8 @@ vec3 SpotLight::diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const {
   float d, att, spot_effect;
   vec3 l_dir, ref;
 
-  l_dir = m_position - i_pos;
-  d = length(l_dir);
-  l_dir = normalize(l_dir);
+  l_dir = normalize(direction(i_pos));
+  d = distance(i_pos);
   spot_effect = dot(m_spot_dir, -l_dir);
   
   if (acos(spot_effect) < radians(m_spot_cutoff)) {
@@ -35,9 +34,8 @@ vec3 SpotLight::specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const {
   float d, att, spot_effect;
   vec3 l_dir, ref;
 
-  l_dir = m_position - i_pos;
-  d = length(l_dir);
-  l_dir = normalize(l_dir);
+  l_dir = normalize(direction(i_pos));
+  d = distance(i_pos);
   spot_effect = dot(m_spot_dir, -l_dir);
 
   if (acos(spot_effect) < radians(m_spot_cutoff)) {

spot_light.hpp

@@ -17,10 +17,20 @@ public:
 
   SpotLight(): PointLight(), m_spot_cutoff(45.0f), m_spot_exponent(0.0f), m_spot_dir(vec3(0.0f, -1.0f, 0.0f)) { }
 
-  SpotLight(vec3 _p, vec3 _d, vec3 _s, float _c, float _l, float _q, float _sco, float _se, vec3 _sd): PointLight(_p, _d, _s, _c, _l, _q),
-												       m_spot_cutoff(_sco),
-												       m_spot_exponent(_se),
-												       m_spot_dir(normalize(_sd)) { }
+  SpotLight(vec3 _p,
+	    vec3 _d,
+	    vec3 _s,
+	    float _c,
+	    float _l,
+	    float _q,
+	    float _sco,
+	    float _se,
+	    vec3 _sd):
+    PointLight(_p, _d, _s, _c, _l, _q),
+    m_spot_cutoff(_sco),
+    m_spot_exponent(_se),
+    m_spot_dir(normalize(_sd))
+  { }
 
   virtual vec3 diffuse(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;
   virtual vec3 specular(vec3 normal, Ray & r, vec3 i_pos, Material & m) const;

whitted_tracer.cpp

@@ -1,3 +1,4 @@
+#include <iostream>
 #include <limits>
 
 #include <glm/gtc/constants.hpp>
@@ -14,8 +15,9 @@ vec3 WhittedTracer::trace_ray(Ray & r, Scene * s, unsigned int rec_level) const
   Figure * _f;
   vec3 n, color, i_pos, ref, dir_diff_color, dir_spec_color;
   Ray mv_r, sr, rr;
-  bool vis;
+  bool vis, is_area_light;
   float kr;
+  AreaLight * al;
 
   t = numeric_limits<float>::max();
   _f = NULL;
@@ -34,8 +36,19 @@ vec3 WhittedTracer::trace_ray(Ray & r, Scene * s, unsigned int rec_level) const
     i_pos = r.m_origin + (t * r.m_direction);
     n = _f->normal_at_int(r, t);
     
+    is_area_light = false;
+    // Check if the object is an area light;
+    for (size_t l = 0; l < s->m_lights.size(); l++) {
+      if (s->m_lights[l]->light_type() == Light::AREA && static_cast<AreaLight *>(s->m_lights[l])->m_figure == _f)
+	is_area_light = true;
+    }
+
+    // If the object is an area light, return it's emission value.
+    if (is_area_light) {
+      return clamp(_f->m_mat->m_emission, 0.0f, 1.0f);
+
     // Check if the material is not reflective/refractive.
-    if (!_f->m_mat->m_refract) {
+    } else if (!_f->m_mat->m_refract) {
 
       // Calculate the direct lighting.
       for (size_t l = 0; l < s->m_lights.size(); l++) {
@@ -51,9 +64,22 @@ vec3 WhittedTracer::trace_ray(Ray & r, Scene * s, unsigned int rec_level) const
 	      break;
 	    }
 	  }
+
 	} else if (s->m_lights[l]->light_type() == Light::AREA) {
-	  // Area lights not supported with Whitted ray tracing.
-	  vis = false;
+	  // Cast a shadow ray towards a sample point on the surface of the light source.
+	  al = static_cast<AreaLight *>(s->m_lights[l]);
+	  al->sample_at_surface(i_pos);
+	  sr = Ray(al->direction(i_pos), i_pos + (n * BIAS));
+
+	  for (size_t f = 0; f < s->m_figures.size(); f++) {
+	    // Avoid self-intersection with the light source.
+	    if (al->m_figure != s->m_figures[f]) {
+	      if (s->m_figures[f]->intersect(sr, _t) && _t < al->distance(i_pos)) {
+		vis = false;
+		break;
+	      }
+	    }
+	  }
 	}
 
 	// Evaluate the shading model accounting for visibility.