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

SCM Repository

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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4820 - (view) (download)

1 : glk 4788 /*
2 :    
3 :     try compiling with (vis15):
4 :    
5 :     diderotc --exec dvresys.diderot
6 :    
7 :     */
8 :    
9 :     input image(3)[] vol ("volume dataset to render") = image("dvresys-vol.nrrd");
10 :     input vec3 camEye ("camera look-from point") = [6, 9, 2];
11 :     input vec3 camAt ("camera look-at point") = [0, 0, 0];
12 :     input vec3 camUp ("camera pseudo-up vector") = [0, 0, 1];
13 :     input real camNear ("relative to look-at point, distance to near clipping plane (where rays start from)") = -3;
14 :     input real camFar ("relative to look-at point, distance to far clipping plane") = 3;
15 :     input real camFOV ("field-of-view angle (in degrees) subtended vertically by view plane") = 15;
16 :     input bool camOrtho ("whether to use orthographic, instead of perspective, projection") = false;
17 :     input int iresU ("# samples on horizontal axis of image") = 640;
18 :     input int iresV ("# samples on vertical axis of image") = 480;
19 :     input real rayStep ("inter-sample distance along view direction") = 0.1;
20 :    
21 :     input real refStep ("reference (unit) step length, for normalizing opacities") = 0.1;
22 :     input real transp0 ("transparency close enough to 0 to terminate ray") = 0.005;
23 :     input real isoval ("isovalue at which to render soft isosurface") = 0;
24 :     input real thick ("approximate thickness (in world-space) of soft isosurface") = 0.1;
25 :     input real maxAlpha ("maximum opacity on rendered surface") = 1;
26 :     input real phongKa ("Blinn-Phong ambient component") = 0.1;
27 :     input real phongKd ("Blinn-Phong diffuse component") = 0.7;
28 :     input real phongKs ("Blinn-Phong specular component") = 0.3;
29 :     input real phongSp ("Blinn-Phong specularity exponent") = 60;
30 :     input vec3 litdir ("direction (non-normalized) towards light source, in (U,V,N) view-space") = [-1, -2, -1];
31 :     input vec3 mcnear ("material color at near clipping plane (for depth cuing)") = [1,1,1];
32 :     input vec3 mcfar ("material color at far clipping plane") = [1,1,1];
33 :     /* this fog has no practical value except as a subtle indicator
34 :     of the extent and shape of the volume data domain */
35 :     input vec4 fog ("fog RGBA inside volume data domain") = [1,1,1,0];
36 :    
37 :     /* Convolve volume data with one of various possible kernels;
38 :     see ../mip/mip.diderot for more info */
39 :     //field#1(3)[] V = c1tent ⊛ vol;
40 :     //field#1(3)[] V = ctmr ⊛ vol;
41 :     field#2(3)[] V = bspln3 ⊛ clamp(vol);
42 :     //field#4(3)[] V = c4hexic ⊛ vol;
43 :    
44 :     /* create a field to render from the original volume data field */
45 :     field#2(3)[] F = V - isoval;
46 :    
47 :     // (boilerplate) computing ray parameters and view-space basis
48 :     vec3 camN = normalize(camAt - camEye); // N: away from eye
49 :     vec3 camU = normalize(camN × camUp); // U: right
50 :     vec3 camV = camN × camU; // V: down (right-handed frame)
51 :     real camDist = |camAt - camEye|;
52 :     real camVmax = tan(camFOV*π/360)*camDist;
53 :     real camUmax = camVmax*iresU/iresV;
54 :     real camNearVsp = camNear + camDist; // near clip, view space
55 :     real camFarVsp = camFar + camDist; // far clip, view space
56 :     // convert light directions from view-space to world-space
57 :     vec3 litwsp = transpose([camU,camV,camN])•normalize(litdir);
58 :    
59 :     /* 2-D opacity function, after Levoy. The 1.4 is a trick to ensure
60 :     that there is a segment of positions receiving maximum opacity. */
61 :     function real alpha(real v, real g) = maxAlpha*clamp(0, 1, 1.4 - |v|/(g*thick));
62 :    
63 :     function real alphx(real x) = maxAlpha*clamp(0, 1, 1.4 - |x|/thick);
64 :    
65 :     /*
66 :     function tensor[3,3] inverse(tensor[3,3] A) {
67 :     real dd = ( A[0,0]*A[1,1]*A[2,2]
68 :     + A[1,0]*A[2,1]*A[0,2]
69 :     + A[2,0]*A[0,1]*A[1,2]
70 :     - A[2,0]*A[1,1]*A[0,2]
71 :     - A[1,0]*A[0,1]*A[2,2]
72 :     - A[0,0]*A[2,1]*A[1,2]);
73 :     tensor[3,3] ret = [[ (A[1,1]*A[2,2] - A[1,2]*A[2,1]),
74 :     -(A[0,1]*A[2,2] - A[0,2]*A[2,1]),
75 :     (A[0,1]*A[1,2] - A[0,2]*A[1,1])],
76 :     [-(A[1,0]*A[2,2] - A[1,2]*A[2,0]),
77 :     (A[0,0]*A[2,2] - A[0,2]*A[2,0]),
78 :     -(A[0,0]*A[1,2] - A[0,2]*A[1,0])],
79 :     [ (A[1,0]*A[2,1] - A[1,1]*A[2,0]),
80 :     -(A[0,0]*A[2,1] - A[0,1]*A[2,0]),
81 :     (A[0,0]*A[1,1] - A[0,1]*A[1,0])]];
82 :     ret /= dd;
83 :     return ret;
84 :     }
85 :     */
86 :    
87 :     // how to render ray through (rayU,rayV) on view plane
88 :     strand raycast(int ui, int vi /* real rayU, real rayV */) {
89 :     // cell-centered sampling of view plane (intersects look-at)
90 :     real rayU = lerp(-camUmax, camUmax, -0.5, ui, iresU-0.5);
91 :     real rayV = lerp(-camVmax, camVmax, -0.5, vi, iresV-0.5);
92 :     // creation of per-strand ray state based on ../mip/mip.diderot
93 :     real rayN = camNearVsp;
94 :     vec3 UV = rayU*camU + rayV*camV;
95 :     vec3 rayOrig = camEye + (UV if camOrtho else [0,0,0]);
96 :     vec3 rayVec = camN + ([0,0,0] if camOrtho else UV/camDist);
97 :     /* alphaFix is used for opacity correction, so that attenuation
98 :     happens in world-space rather than ray sample index space.
99 :     The actual distance between samples on this ray is |rayVec|*rayStep */
100 :     real alphaFix = |rayVec|*rayStep/refStep;
101 :     vec3 eyeDir = -normalize(rayVec);
102 :     // output for this ray
103 :     output vec4 rgba = [0,0,0,0];
104 :     // state for this ray is current color ...
105 :     vec3 rgb = [0,0,0];
106 :     // ... and current tranparency
107 :     real transp = 1;
108 :     /* example of turning on debuging for one strand (pixel);
109 :     can can then have "if (verb) { print(...); }"
110 :     bool verb = 369==ui && 242==vi; */
111 :    
112 :     update {
113 :     rayN += rayStep; // increment ray position
114 :     if (rayN > camFarVsp) { // ray passed the far clipping plane
115 :     stabilize;
116 :     }
117 :     vec3 pos = rayOrig + rayN*rayVec; // pos is ray sample position
118 :     if (!inside(pos,F)) { // If not inside field domain,
119 :     continue; // then move onto next iteration
120 :     }
121 :     // compute fog contribution
122 :     real aa = fog[3];
123 :     if (aa > 0) {
124 :     aa = 1 - (1 - aa)^alphaFix;
125 :     rgb += transp*aa*[fog[0],fog[1],fog[2]];
126 :     transp *= 1 - aa;
127 :     }
128 :     // compute data contribution
129 : cchiw 4820 field#0(3)[3,3] hhF = ∇⊗∇F;
130 :     field#0(3)[3,3] hiF = inv(hhF);
131 :     field#0(3)[3] deltaF = hiF•∇F;
132 :     vec3 delta = deltaF(pos);
133 :     real{3} eval = evals(hhF(pos));
134 : glk 4788 real sthr = 0.001;
135 :     if (eval{0} < 0 && eval{1} < 0 && eval{2} < 0 && eval{0} < -sthr) {
136 :     aa = alphx(|delta| - 0.3);
137 :     }
138 :     if (aa == 0) {
139 :     continue;
140 :     }
141 :     vec3 grad = -delta;
142 :     aa = 1 - (1 - aa)^alphaFix;
143 :     // Note: not standard Phong diffuse (no dark hemisphere)
144 :     real dcomp = lerp(0, 1, -1, normalize(grad)•litwsp, 1)^2;
145 :     // If phongKs=0, this is conditional expression speeds things slightly
146 :     real scomp = max(0,normalize(grad)•normalize(eyeDir+litwsp))^phongSp
147 :     if phongKs != 0 else 0.0;
148 :     // simple depth-cueing
149 :     vec3 dcol = lerp(mcnear, mcfar, camNearVsp, rayN, camFarVsp);
150 :     // lots of things could be basis for material color; this
151 :     // is contrived to show off some univariate mapping
152 :     vec3 mcol = [1,1,1];
153 :     // light color is currently [1,1,1]
154 :     rgb += transp*aa*((phongKa + phongKd*dcomp)*modulate(dcol,mcol)
155 :     + phongKs*scomp*[1,1,1]);
156 :     transp *= 1 - aa;
157 :     if (transp < transp0) { // early ray termination
158 :     transp = 0;
159 :     stabilize;
160 :     }
161 :     }
162 :     stabilize {
163 :     if (transp < 1) { // un-pre-multiply opacities
164 :     real aa = 1-transp;
165 :     rgba = [rgb[0]/aa, rgb[1]/aa, rgb[2]/aa, aa];
166 :     }
167 :     }
168 :     }
169 :    
170 :     initially [ raycast(ui, vi)
171 :     | vi in 0..iresV-1, // slower
172 :     ui in 0..iresU-1 ]; // faster

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