Home My Page Projects Code Snippets Project Openings diderot
Summary Activity Tracker Tasks SCM

SCM Repository

[diderot] View of /tests/vis15-bugs/dvresys.diderot
ViewVC logotype

View of /tests/vis15-bugs/dvresys.diderot

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4820 - (download) (annotate)
Tue Nov 1 21:29:40 2016 UTC (2 years, 9 months ago) by cchiw
File size: 7440 byte(s)
added varcount
/*

try compiling with (vis15):

diderotc --exec dvresys.diderot

*/

input image(3)[] vol ("volume dataset to render") = image("dvresys-vol.nrrd");
input vec3 camEye ("camera look-from point") = [6, 9, 2];
input vec3 camAt ("camera look-at point") = [0, 0, 0];
input vec3 camUp ("camera pseudo-up vector") = [0, 0, 1];
input real camNear ("relative to look-at point, distance to near clipping plane (where rays start from)") = -3;
input real camFar ("relative to look-at point, distance to far clipping plane") = 3;
input real camFOV ("field-of-view angle (in degrees) subtended vertically by view plane") = 15;
input bool camOrtho ("whether to use orthographic, instead of perspective, projection") = false;
input int iresU ("# samples on horizontal axis of image") = 640;
input int iresV ("# samples on vertical axis of image") = 480;
input real rayStep ("inter-sample distance along view direction") = 0.1;

input real refStep ("reference (unit) step length, for normalizing opacities") = 0.1;
input real transp0 ("transparency close enough to 0 to terminate ray") = 0.005;
input real isoval ("isovalue at which to render soft isosurface") = 0;
input real thick ("approximate thickness (in world-space) of soft isosurface") = 0.1;
input real maxAlpha ("maximum opacity on rendered surface") = 1;
input real phongKa ("Blinn-Phong ambient component") = 0.1;
input real phongKd ("Blinn-Phong diffuse component") = 0.7;
input real phongKs ("Blinn-Phong specular component") = 0.3;
input real phongSp ("Blinn-Phong specularity exponent") = 60;
input vec3 litdir ("direction (non-normalized) towards light source, in (U,V,N) view-space") = [-1, -2, -1];
input vec3 mcnear ("material color at near clipping plane (for depth cuing)") = [1,1,1];
input vec3 mcfar ("material color at far clipping plane") = [1,1,1];
/* this fog has no practical value except as a subtle indicator
   of the extent and shape of the volume data domain */
input vec4 fog ("fog RGBA inside volume data domain") = [1,1,1,0];

/* Convolve volume data with one of various possible kernels;
   see ../mip/mip.diderot for more info */
//field#1(3)[] V = c1tent ⊛ vol;
//field#1(3)[] V = ctmr ⊛ vol;
field#2(3)[] V = bspln3 ⊛ clamp(vol);
//field#4(3)[] V = c4hexic ⊛ vol;

/* create a field to render from the original volume data field */
field#2(3)[] F = V - isoval;

// (boilerplate) computing ray parameters and view-space basis
vec3 camN = normalize(camAt - camEye);  // N: away from eye
vec3 camU = normalize(camN × camUp);    // U: right
vec3 camV = camN × camU;                // V: down (right-handed frame)
real camDist = |camAt - camEye|;
real camVmax = tan(camFOV*π/360)*camDist;
real camUmax = camVmax*iresU/iresV;
real camNearVsp = camNear + camDist; // near clip, view space
real camFarVsp = camFar + camDist;   // far clip, view space
// convert light directions from view-space to world-space
vec3 litwsp = transpose([camU,camV,camN])•normalize(litdir);

/* 2-D opacity function, after Levoy. The 1.4 is a trick to ensure
  that there is a segment of positions receiving maximum opacity. */
function real alpha(real v, real g) = maxAlpha*clamp(0, 1, 1.4 - |v|/(g*thick));

function real alphx(real x) = maxAlpha*clamp(0, 1, 1.4 - |x|/thick);

/*
function tensor[3,3] inverse(tensor[3,3] A) {
  real dd = (  A[0,0]*A[1,1]*A[2,2]
             + A[1,0]*A[2,1]*A[0,2]
             + A[2,0]*A[0,1]*A[1,2]
             - A[2,0]*A[1,1]*A[0,2]
             - A[1,0]*A[0,1]*A[2,2]
             - A[0,0]*A[2,1]*A[1,2]);
  tensor[3,3] ret = [[ (A[1,1]*A[2,2] - A[1,2]*A[2,1]),
                      -(A[0,1]*A[2,2] - A[0,2]*A[2,1]),
                       (A[0,1]*A[1,2] - A[0,2]*A[1,1])],
                     [-(A[1,0]*A[2,2] - A[1,2]*A[2,0]),
                       (A[0,0]*A[2,2] - A[0,2]*A[2,0]),
                      -(A[0,0]*A[1,2] - A[0,2]*A[1,0])],
                     [ (A[1,0]*A[2,1] - A[1,1]*A[2,0]),
                      -(A[0,0]*A[2,1] - A[0,1]*A[2,0]),
                       (A[0,0]*A[1,1] - A[0,1]*A[1,0])]];
  ret /= dd;
  return ret;
}
*/

// how to render ray through (rayU,rayV) on view plane
strand raycast(int ui, int vi /* real rayU, real rayV */) {
   // cell-centered sampling of view plane (intersects look-at)
   real rayU = lerp(-camUmax, camUmax, -0.5, ui, iresU-0.5);
   real rayV = lerp(-camVmax, camVmax, -0.5, vi, iresV-0.5);
   // creation of per-strand ray state based on ../mip/mip.diderot
   real rayN = camNearVsp;
   vec3 UV = rayU*camU + rayV*camV;
   vec3 rayOrig = camEye + (UV if camOrtho else [0,0,0]);
   vec3 rayVec = camN + ([0,0,0] if camOrtho else UV/camDist);
   /* alphaFix is used for opacity correction, so that attenuation
   happens in world-space rather than ray sample index space.
   The actual distance between samples on this ray is |rayVec|*rayStep */
   real alphaFix = |rayVec|*rayStep/refStep;
   vec3 eyeDir = -normalize(rayVec);
   // output for this ray
   output vec4 rgba = [0,0,0,0];
   // state for this ray is current color ...
   vec3 rgb = [0,0,0];
   // ... and current tranparency
   real transp = 1;
   /* example of turning on debuging for one strand (pixel);
      can can then have "if (verb) { print(...); }"
   bool verb = 369==ui && 242==vi; */

   update {
      rayN += rayStep;          // increment ray position
      if (rayN > camFarVsp) {   // ray passed the far clipping plane
         stabilize;
      }
      vec3 pos = rayOrig + rayN*rayVec;  // pos is ray sample position
      if (!inside(pos,F)) {              // If not inside field domain,
         continue;                       // then move onto next iteration
      }
      // compute fog contribution
      real aa = fog[3];
      if (aa > 0) {
         aa = 1 - (1 - aa)^alphaFix;
         rgb += transp*aa*[fog[0],fog[1],fog[2]];
         transp *= 1 - aa;
      }
      // compute data contribution
      field#0(3)[3,3] hhF = ∇⊗∇F;
      field#0(3)[3,3] hiF = inv(hhF);
      field#0(3)[3] deltaF = hiF•∇F;
      vec3 delta = deltaF(pos);
      real{3} eval = evals(hhF(pos));
      real sthr = 0.001;
      if (eval{0} < 0 && eval{1} < 0 && eval{2} < 0 && eval{0} < -sthr) {
         aa = alphx(|delta| - 0.3);
      }
      if (aa == 0) {
         continue;
      }
      vec3 grad = -delta;
      aa = 1 - (1 - aa)^alphaFix;
      // Note: not standard Phong diffuse (no dark hemisphere)
      real dcomp = lerp(0, 1, -1, normalize(grad)•litwsp, 1)^2;
      // If phongKs=0, this is conditional expression speeds things slightly
      real scomp = max(0,normalize(grad)•normalize(eyeDir+litwsp))^phongSp
                   if phongKs != 0 else 0.0;
      // simple depth-cueing
      vec3 dcol = lerp(mcnear, mcfar, camNearVsp, rayN, camFarVsp);
      // lots of things could be basis for material color; this
      // is contrived to show off some univariate mapping
      vec3 mcol = [1,1,1];
      // light color is currently [1,1,1]
      rgb += transp*aa*((phongKa + phongKd*dcomp)*modulate(dcol,mcol)
                        + phongKs*scomp*[1,1,1]);
      transp *= 1 - aa;
      if (transp < transp0) { // early ray termination
         transp = 0;
         stabilize;
      }
   }
   stabilize {
      if (transp < 1) {  // un-pre-multiply opacities
         real aa = 1-transp;
         rgba = [rgb[0]/aa, rgb[1]/aa, rgb[2]/aa, aa];
      }
   }
}

initially [ raycast(ui, vi)
            | vi in 0..iresV-1,   // slower
              ui in 0..iresU-1 ]; // faster

root@smlnj-gforge.cs.uchicago.edu
ViewVC Help
Powered by ViewVC 1.0.0