Text archives Help

---

Author: arobison
Date: Tue Feb 12 01:02:17 2008
New Revision: 2051

Modified:
   trunk/Model/Primitives/KenslerShirleyTriangle.cc
   trunk/Model/Primitives/KenslerShirleyTriangle.h
   trunk/Model/Primitives/MeshTriangle.cc
   trunk/Model/Primitives/MeshTriangle.h
Log:
Adding a MeshTriangle method for computing smooth surface derivatives.
Currently only KenslerShirley triangles use it.  There's also a small
bug somewhere flipping some of my barycentric coordinates around, that
will be fixed shortly.



Modified: trunk/Model/Primitives/KenslerShirleyTriangle.cc
==============================================================================
--- trunk/Model/Primitives/KenslerShirleyTriangle.cc    (original)
+++ trunk/Model/Primitives/KenslerShirleyTriangle.cc    Tue Feb 12 01:02:17 
2008
@@ -46,6 +46,22 @@
   context.done();
 }
 
+void KenslerShirleyTriangle::computeSurfaceDerivatives(const RenderContext& 
context, RayPacket& rays) const
+{
+  for(int i = rays.begin(); i != rays.end(); ++i) {
+    float a = rays.getScratchpad<float>(SCRATCH_U)[i];
+    float b = rays.getScratchpad<float>(SCRATCH_V)[i];
+
+    Vector dPdu, dPdv;
+    computeSplineSurfaceDerivatives(a,b, dPdu,dPdv);
+
+    rays.setSurfaceDerivativeU(i, dPdu);
+    rays.setSurfaceDerivativeV(i, dPdv);
+  }
+
+  rays.setFlag(RayPacket::HaveSurfaceDerivatives);
+}
+
 void KenslerShirleyTriangle::computeTexCoords2(const RenderContext&, 
RayPacket& rays) const
 {
   const int which = myID*3;

Modified: trunk/Model/Primitives/KenslerShirleyTriangle.h
==============================================================================
--- trunk/Model/Primitives/KenslerShirleyTriangle.h     (original)
+++ trunk/Model/Primitives/KenslerShirleyTriangle.h     Tue Feb 12 01:02:17 
2008
@@ -38,6 +38,8 @@
       computeTexCoords2(context, rays);
     }
 
+    virtual void computeSurfaceDerivatives(const RenderContext&, RayPacket& 
rays) const;
+
     void intersect(const RenderContext& context, RayPacket& rays) const;
 
     void computeNormal(const RenderContext& context, RayPacket &rays) const;

Modified: trunk/Model/Primitives/MeshTriangle.cc
==============================================================================
--- trunk/Model/Primitives/MeshTriangle.cc      (original)
+++ trunk/Model/Primitives/MeshTriangle.cc      Tue Feb 12 01:02:17 2008
@@ -58,3 +58,165 @@
   }
   rays.setFlag( RayPacket::HaveSurfaceDerivatives );
 }
+
+void MeshTriangle::computeSplineSurfaceDerivatives(float alpha, float beta,
+                                                   Vector& dPdu, Vector& 
dPdv) const
+{ 
+  // This code based on the ATI PN Triangle paper
+  // http://portal.acm.org/citation.cfm?id=364387
+
+  // The basic idea is to fit a three-sided cubic bezier patch to the
+  // triangle and then pull the surface derivatives out of the
+  // analytic surface at the hit point as defined by the given
+  // barycentric coordinates.  Once the surface derivatives have been
+  // found with respect to the barycentric coordinates, we use the
+  // inverse function theorem to find the derivatives with respect to
+  // the UV parameterization of the triangle as defined by the texture
+  // coordinates.
+
+  // TODO(arobison): a *lot* of this code can be precomputed and stored
+  // in the triangle, this probably suggests a new MeshTriangle type
+  // that supports spline fitting, or just a modification of the existing
+  // types to store the spline control points and some of the constant 
+  // derivatives.
+
+  // Grab triangle components
+
+  const int which = myID*3;
+
+  const unsigned int uv0_idx = mesh->texture_indices[which+0];
+  const unsigned int uv1_idx = mesh->texture_indices[which+1];
+  const unsigned int uv2_idx = mesh->texture_indices[which+2];
+
+  const Vector uv0 = mesh->texCoords[uv0_idx];
+  const Vector uv1 = mesh->texCoords[uv1_idx];
+  const Vector uv2 = mesh->texCoords[uv2_idx];
+
+  const unsigned int vec0_idx = mesh->vertex_indices[which+0];
+  const unsigned int vec1_idx = mesh->vertex_indices[which+1];
+  const unsigned int vec2_idx = mesh->vertex_indices[which+2];
+
+  const Vector vec1 = mesh->vertices[vec0_idx];
+  const Vector vec2 = mesh->vertices[vec1_idx];
+  const Vector vec3 = mesh->vertices[vec2_idx];
+
+  const unsigned int n0_idx = mesh->normal_indices[which+0];
+  const unsigned int n1_idx = mesh->normal_indices[which+1];
+  const unsigned int n2_idx = mesh->normal_indices[which+2];
+
+  const Vector n1 = mesh->vertexNormals[n0_idx];
+  const Vector n2 = mesh->vertexNormals[n1_idx];
+  const Vector n3 = mesh->vertexNormals[n2_idx];
+
+  // Compute spline control vertices
+
+  const Vector b300 = vec1;
+  const Vector b030 = vec2;
+  const Vector b003 = vec3;
+
+  const Real w12 = Dot((vec2 - vec1), n1);
+  const Real w21 = Dot((vec1 - vec2), n2);
+  
+  const Real w23 = Dot((vec3 - vec2), n2);
+  const Real w32 = Dot((vec2 - vec3), n3);
+  
+  const Real w13 = Dot((vec3 - vec1), n1);
+  const Real w31 = Dot((vec1 - vec3), n3);
+
+  const Vector b210 = (2.*vec1 + vec2 - w12*n1)/3.;
+  const Vector b120 = (2.*vec2 + vec1 - w21*n2)/3.;
+  const Vector b021 = (2.*vec2 + vec3 - w23*n2)/3.;
+  const Vector b012 = (2.*vec3 + vec2 - w32*n3)/3.;
+  const Vector b102 = (2.*vec3 + vec1 - w31*n3)/3.;
+  const Vector b201 = (2.*vec1 + vec3 - w13*n1)/3.;
+
+  const Vector E = (b210 + b120 + b021 + b012 + b102 + b201)/6.;
+  const Vector V = (vec1 + vec2 + vec3)/3.;
+
+  const Vector b111 = E + (E - V)/2.;
+
+  // The spline surface is defined by the following equation:
+  //
+  // b(alpha, beta) = sum {i+j+k=3} (b_ijk * 3!/(i!j!k!) alpha^i * beta^j * 
gamma^k
+  // where gamma = 1 - alpha - beta
+  // or:
+  //
+  // b(u,v,w) = b300*u^3 + b030*v^3 + b003*w^3 +
+  //            3*b210*u^2*v + 3*b120*u*v^2 +
+  //            3*b021*v^2*w + 3*b012*v*w^2 +
+  //            3*b201*u^2*w + 3*b102*u*w^2 +
+  //            6*b111*u*v*w
+  //
+  // so the partials with respect to the barycentric coordinates are:
+  //
+
+  const Real gamma = 1 - alpha - beta;
+
+#if 0
+  // from the paper's incorrect expansion of the sum
+  const Vector dPdalpha = 3*b030*alpha*alpha + b210*3*gamma*gamma + 
2*b120*3*gamma*alpha +
+    2*b021*3*alpha*beta + b012*3*beta*beta + b111*6*gamma*beta;
+  const Vector dPdbeta = 3*b003*beta*beta + b201*3*gamma*gamma + 
2*b021*3*alpha*alpha +
+    2*b102*3*gamma*beta + 2*b012*3*alpha*beta + b111*6*gamma*alpha;
+  const Vector dPdgamma = 3*b300*gamma*gamma + 2*b210*3*gamma*alpha + 
b120*3*alpha*alpha +
+    2*b201*3*gamma*beta + b102*3*beta*beta + b111*6*alpha*beta;
+#else
+  // correct expansion
+  const Vector dPdalpha = 3*b300*alpha*alpha + 
+                          6*b210*alpha*beta  + 3*b120*beta*beta +
+                          6*b201*alpha*gamma + 3*b102*gamma*gamma +
+                          6*b111*beta*gamma;
+
+  const Vector dPdbeta = 3*b030*beta*beta +
+                         3*b210*alpha*alpha  + 6*b120*alpha*beta +
+                         6*b021*beta*gamma   + 3*b012*gamma*gamma +
+                         6*b111*alpha*gamma;
+
+  const Vector dPdgamma = 3*b003*gamma*gamma +
+                          3*b021*beta*beta   + 6*b012*beta*gamma +
+                          3*b201*alpha*alpha + 6*b102*alpha*gamma +
+                          6*b111*alpha*beta;
+#endif
+
+  // We know that:
+  // dPdu = dPdalpha*dalphadu + dPdbeta*dbetadu + dPdgamma*dgammadu
+  // dPdv = dPdalpha*dalphadv + dPdbeta*dbetadv + dPdgamma*dgammadv
+
+  // So we need to compute d{alpha,beta,gamma}d{u,v}. We use the jacobian 
and invert it
+  // to get the local inverse derivatives
+
+  // Since u(a,b) = u0*a + u1*b + u2*(1-a-b) (similarly for v(a,b))
+  // we know the derivatives are:
+  const Real dudalpha = uv0.x()-uv2.x();
+  const Real dudbeta = uv1.x()-uv2.x();
+  const Real dvdalpha = uv0.y()-uv2.y();
+  const Real dvdbeta = uv1.y()-uv2.y();
+
+  // Now invert the matrix
+  const Real determinant = dudalpha * dvdbeta - dudbeta * dvdalpha;
+
+  if(Abs(determinant < 1e-8)) {
+    dPdu = Vector(1,0,0);
+    dPdv = Vector(0,1,0);
+    return;
+  } 
+
+  const Real inv_det = 1./determinant;
+
+  // And pull out the inverted derivatives
+  const Real dalphadu = inv_det * dvdbeta;
+  const Real dalphadv = inv_det * -dudbeta;
+  const Real dbetadu = inv_det * -dvdalpha;
+  const Real dbetadv = inv_det * dudalpha;
+
+  // These follow from: gamma = (1 - alpha - beta)
+  const Real dgammadu = -dalphadu - dbetadu;
+  const Real dgammadv = -dalphadv - dbetadv;
+
+  // Now that we have the correct local derivative information, change from 
+  // barycentric coordinates to UV space
+  
+  dPdu = (dPdalpha * dalphadu) + (dPdbeta * dbetadu) + (dPdgamma * dgammadu);
+  dPdv = (dPdalpha * dalphadv) + (dPdbeta * dbetadv) + (dPdgamma * dgammadv);
+
+}

Modified: trunk/Model/Primitives/MeshTriangle.h
==============================================================================
--- trunk/Model/Primitives/MeshTriangle.h       (original)
+++ trunk/Model/Primitives/MeshTriangle.h       Tue Feb 12 01:02:17 2008
@@ -36,6 +36,9 @@
     virtual void computeSurfaceDerivatives(const RenderContext& context,
                                            RayPacket& rays) const;
 
+    void computeSplineSurfaceDerivatives(float alpha, float beta,
+                                         Vector& dPdu, Vector& dPdv) const;
+
     //List of the triangle classes that implement MeshTriangle.
     enum TriangleType {
       WALD_TRI,

---

• [Manta] r2051 - trunk/Model/Primitives, Austin Robison, 02/12/2008

---

Archive powered by MHonArc 2.6.16.