title: What can cause Phong specular shading to produce gamut overflows?
body: <p>I am currently implementing a basic raytracer in c++. Works pretty well so far, matte materials (with ambient and diffuse brdf) work as expected so far.</p>
  
  <p>Adding specular highlights would result in the full <a href="http://en.wikipedia.org/wiki/Phong_shading" rel="nofollow noreferrer">Phong Model</a> and that's exactly what i tried to do.</p>
  
  <p>Unfortunately, i encounter gamut overflow, with all kinds of values for the specular reflection constant ks and the exponent. Here are some examples.</p>
  
  <pre><code>// Sphere material definition:
  ka    = 0.9;
  kd    = 1.0;
  ks    = 0.3;
  exp   = 1.0;
  color = rgb(1.0, 1.0, 0.98);
  </code></pre>
  
  <p>image: <a href="http://dl.dropbox.com/u/614366/cornell_1.png" rel="nofollow noreferrer"><a href="http://dl.dropbox.com/u/614366/cornell_1.png" rel="nofollow noreferrer">http://dl.dropbox.com/u/614366/cornell_1.png</a></a></p>
  
  <pre><code>// Sphere material definition:
  ka    = 0.9;
  kd    = 1.0;
  ks    = 0.3;
  exp   = 20.0;                  // only changed exp
  color = rgb(1.0, 1.0, 0.98);
  </code></pre>
  
  <p>image: <a href="http://dl.dropbox.com/u/614366/cornell_2.png" rel="nofollow noreferrer"><a href="http://dl.dropbox.com/u/614366/cornell_2.png" rel="nofollow noreferrer">http://dl.dropbox.com/u/614366/cornell_2.png</a></a></p>
  
  <pre><code>// Sphere material definition:
  ka    = 0.9;
  kd    = 1.0;
  ks    = 0.1;                  // only changes here
  exp   = 0.1;                  // and here
  color = rgb(1.0, 1.0, 0.98);
  </code></pre>
  
  <p>image: <a href="http://dl.dropbox.com/u/614366/cornell_3.png" rel="nofollow noreferrer"><a href="http://dl.dropbox.com/u/614366/cornell_3.png" rel="nofollow noreferrer">http://dl.dropbox.com/u/614366/cornell_3.png</a></a></p>
  
  <p>and here are some relevant excerpts of the code:</p>
  
  <p>in raycast.cpp</p>
  
  <pre><code>namespace {
      const float floatmax = std::numeric_limits&lt;float&gt;::max();
  }
  
  rgb
  RayCast::trace ( const Ray&amp; ray ) const
  {
      HitRecord rec(scene_ptr_);
      float tmax = floatmax;
      float tmin = 0.0;
  
      if ( scene_ptr_-&gt;shapes.hit(ray,tmin,tmax,rec) )
      {
          rec.ray = ray;
          return rec.material_ptr-&gt;shade(rec);
      }
  
      return scene_ptr_-&gt;bgcolor;
  }
  </code></pre>
  
  <p>in phong.cpp</p>
  
  <pre><code>rgb
  Phong::shade ( HitRecord&amp; hitrec ) const
  {
      Vector wo = - hitrec.ray.dir();
  
      rgb L = ambient_brdf_ptr_-&gt;rho(hitrec,wo) *
          hitrec.scene_ptr-&gt;ambient_ptr-&gt;L(hitrec);
  
      int num_lights = hitrec.scene_ptr-&gt;lights.size();
  
      for (int i = 0; i &lt; num_lights; ++i)
      {
          Vector wi     = hitrec.scene_ptr-&gt;lights[i]-&gt;get_direction(hitrec);
          float  ndotwi = dot(hitrec.normal, wi);
  
          if ( ndotwi &gt; 0.0 )
          {
              L += ( diffuse_brdf_ptr_-&gt;f (hitrec, wo, wi)  +
                     specular_brdf_ptr_-&gt;f(hitrec, wo, wi)
                   ) * hitrec.scene_ptr-&gt;lights[i]-&gt;L(hitrec) * ndotwi;
          }
      }
  
      return L;
  }
  </code></pre>
  
  <p>in specular.cpp</p>
  
  <pre><code>namespace {
      const rgb black(0.0,0.0,0.0);
  }
  
  rgb
  Specular::f ( const HitRecord&amp; hitrec, const Vector&amp; wo, const Vector&amp; wi ) const
  {
      rgb L(0,0,0);
      float ndotwi = dot(hitrec.normal, wi);
  
      Vector r = -wi + 2.0 * hitrec.normal * ndotwi;
      float rdotwo = dot(r, wo);
  
      // reflection detected
      if ( rdotwo &gt; 0.0 )
          L = ks_ * pow(rdotwo, exp_);
  
      return L;
  }
  
  rgb
  Specular::rho ( const HitRecord&amp; hitrec, const Vector&amp; wo ) const
  {
      return black;
  }
  </code></pre>
  
  <p>in sphere.cpp</p>
  
  <pre><code>bool
  Sphere::hit ( const Ray&amp; ray, interval_t tmin, interval_t tmax, HitRecord&amp; hitrec ) const
  {
      Vector org = ray.origin() - center_;
      Vector dir = ray.dir();
      float a = dot(dir, dir);
      float b = dot(dir, org) * 2;
      float c = dot(org, org) - pow(radius_, 2);
      float discriminant = pow(b,2) - 4*a*c;
  
      if ( discriminant &gt; 0 )
      {
          discriminant = sqrt(discriminant);
          double t = ( -b - discriminant ) / ( 2*a );
  
          if ( t &lt; tmin )
              t = ( -b + discriminant ) / ( 2*a );
  
          if ( t &gt; tmin and t &lt; tmax )
          {
              // hit detected
              hitrec.t            = t;
              hitrec.hit          = true;
              hitrec.normal       = unify( t*ray.dir() + org );
              hitrec.material_ptr = material_ptr_;
              hitrec.hitpoint     = ray.origin() + t * ray.dir();
              hitrec.ray          = ray;
  
              return true;
          }
      }
  
      return false;
  }
  </code></pre>
  
  <p>do you have ideas where the error could be caused? what are possible factors that lead to such results?</p>
  
  <p>thanks in advance,
  patrick.</p>