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

SCM Repository

[diderot] Diff of /branches/pure-cfg/test/vr-curv-quad.diderot
ViewVC logotype

Diff of /branches/pure-cfg/test/vr-curv-quad.diderot

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 866, Tue Apr 19 20:26:44 2011 UTC revision 911, Thu Apr 21 18:37:43 2011 UTC
# Line 1  Line 1 
1  // vr-curv  // vr-curv
2  //  //
3  // curvature-based volume rendering  // Demonstration of curvature-based transfer functions in volume rendering
4    //
5    // process output with:
6    // unu reshape -i mip.txt -s 4 480 345 | overrgb -i - -b 0.1 0.15 0.2 -g 1.2 -o - | unu quantize -b 8 -min 0 -max 255  -o vr-curv.png
7    
8    image(3)[] img = load("../data/quad-patches-pad.nrrd");
9    field#2(3)[] F = img ⊛ bspln3;
10    image(2)[3] imgRGB = load("../data/hue.nhdr");
11    field#0(2)[3] RGB = imgRGB ⊛ tent;
12    
 string dataFile = "../data/quad-patches-pad.nrrd";  
13  vec3 camEye = [183.021, -411.857, 686.458];  vec3 camEye = [183.021, -411.857, 686.458];
14  vec3 camAt = [55.4037, 40.2407, 4.32732];  vec3 camAt = [55.4037, 40.2407, 4.32732];
15  vec3 camUp = [0.0, 0.0, 1.0];  vec3 camUp = [0.0, 0.0, 1.0];
16  real camNear = -42.0;  real camNear = -42.0;
17  real camFar = 42.0;  real camFar = 42.0;
18  real camFOV = 7.15;  real camFOV = 7.15;
 //int imgResU = 640;  
 //int imgResV = 480;  
19  int imgResU = 480;  int imgResU = 480;
20  int imgResV = 345;  int imgResV = 345;
21  //int imgResU = 320;  real rayStep = 0.2;
22  //int imgResV = 240;  real valOpacMin = -0.25;
23  real rayStep = 2.0;  real valOpacMax = -0.10;
 real valOpacMin = -0.1;  
 real valOpacMax = 0.1;  
 image(2)[] imgR = load("../data/hue-0.nrrd"); field#0(2)[] RR = imgR ⊛ tent;  
 image(2)[] imgG = load("../data/hue-1.nrrd"); field#0(2)[] GG = imgG ⊛ tent;  
 image(2)[] imgB = load("../data/hue-2.nrrd"); field#0(2)[] BB = imgB ⊛ tent;  
 //image(2)[3] imgRGB = load("../data/hue.nhdr"); field#0(2)[3] RGB = imgRGB ⊛ tent;  
24    
25  real camDist = |camAt - camEye|;  real camDist = |camAt - camEye|;
26  real camVspNear = camNear + camDist;  real camVspNear = camNear + camDist;
# Line 32  Line 31 
31  real camVmax = tan(camFOV*3.1415926536/360.0)*camDist;  real camVmax = tan(camFOV*3.1415926536/360.0)*camDist;
32  real camUmax = camVmax*real(imgResU)/real(imgResV);  real camUmax = camVmax*real(imgResU)/real(imgResV);
33    
34  vec3 lightVspDir = [0.9, -1.0, -2.5];  vec3 lightVspDir = [0.9, -2.5, -2.5];
35  vec3 lightDir = normalize(lightVspDir[0]*camU + lightVspDir[1]*camV + lightVspDir[2]*camN);  vec3 lightDir = normalize(lightVspDir[0]*camU +
36                              lightVspDir[1]*camV +
37  real phongKa = 0.05;                            lightVspDir[2]*camN);
38  real phongKd = 0.80;  vec3 lightRGB = [1.0, 1.0, 1.0];
39  real phongKs = 0.20;  
40  real phongSp = 30.0;  real phongKa = 0.1;
41    real phongKd = 0.7;
42  image(3)[] img = load (dataFile);  real phongKs = 0.3;
43  field#2(3)[] F = img ⊛ bspln3;  real phongSp = 45.0;
44    
45  strand RayCast (int ui, int vi)  strand RayCast (int ui, int vi) {
 {  
46      real rayU = lerp(-camUmax, camUmax, -0.5, real(ui), real(imgResU)-0.5);      real rayU = lerp(-camUmax, camUmax, -0.5, real(ui), real(imgResU)-0.5);
47      real rayV = lerp(-camVmax, camVmax, -0.5, real(vi), real(imgResV)-0.5);      real rayV = lerp(-camVmax, camVmax, -0.5, real(vi), real(imgResV)-0.5);
48      vec3 rayVec = (camDist*camN + rayU*camU + rayV*camV)/camDist;      vec3 rayVec = (camDist*camN + rayU*camU + rayV*camV)/camDist;
49    
50      real rayN = camVspNear;      real rayN = camVspNear;
51      real transp = 1.0;      real rayTransp = 1.0;
52      real gray = 0.0;      vec3 rayRGB = [0.0, 0.0, 0.0];
53      real or = 0.0;      output vec4 outRGBA = [0.0, 0.0, 0.0, 0.0];
     real og = 0.0;  
     real ob = 0.0;  
     output vec4 rgba = [0.0, 0.0, 0.0, 0.0];  
54    
55      update      update {
     {  
56         vec3 pos = camEye + rayN*rayVec;         vec3 pos = camEye + rayN*rayVec;
57         if (inside (pos,F)) {         if (inside (pos,F)) {
58            // ########## BEGIN per-sample code            real val = F(pos);
59            real val = F@pos;            vec3 grad = -∇F(pos);
           vec3 grad = -∇F@pos;  
60            vec3 norm = normalize(grad);            vec3 norm = normalize(grad);
61            if (val > valOpacMin) {  // we have some opacity            if (val > valOpacMin) {  // we have some opacity
62                tensor[3,3] H = ∇(∇F)@pos;                // begin curvature computation
63                  tensor[3,3] H = -∇(∇F)(pos);
64                tensor[3,3] P = identity[3] - norm⊗norm;                tensor[3,3] P = identity[3] - norm⊗norm;
65                tensor[3,3] G = (P•H•P)/|grad|;                tensor[3,3] G = (P•H•P)/|grad|;
66                real disc = max(0.0, sqrt(2.0*|G|^2 - trace(G)^2));                real disc = max(0.0, sqrt(2.0*|G|^2 - trace(G)^2));
67                real k1 = max(-0.19, min(0.19, (trace(G) + disc)/2.0));                real k1 = (trace(G) + disc)/2.0;
68                real k2 = max(-0.19, min(0.19, (trace(G) - disc)/2.0));                real k2 = (trace(G) - disc)/2.0;
69                vec2 kk = [k1,k2];                // finished curvature computation; begin finding sample RGBA
70                vec3 rgb = RGB@kk if inside(kk,RGB)                k1 = max(-1.0, min(1.0, 4.2*k1));
71                      else [0.0,0.0,0.0];                k2 = max(-1.0, min(1.0, 4.2*k2));
72                real opac = min(1.0, lerp(0.0, 1.0, valOpacMin, val, valOpacMax));                vec3 matRGB = RGB([k1,k2]) if inside([k1,k2],RGB)
73                real ld = lightDir • norm;                              else [1.0,0.0,1.0];
74                real mat = (  phongKa                real alpha = min(1.0, lerp(0.0, 1.0,
75                            + phongKd*(ld if ld > 0.0 else 0.0)                                           valOpacMin, val, valOpacMax));
76                            + phongKs*(ld^phongSp if ld > 0.0 else 0.0));                real ld = max(0.0, norm • lightDir);
77                gray = gray + transp*opac*mat;                real hd = max(0.0, norm • normalize(lightDir +
78                or = or + transp*opac*mat*rgb[0];                                                    normalize(camEye - pos)));
79                og = og + transp*opac*mat*rgb[1];                // Phong shading
80                ob = ob + transp*opac*mat*rgb[2];                vec3 pntRGB = (phongKa*matRGB
81                transp = transp*(1.0 - opac);                               + phongKd*ld*modulate(matRGB, lightRGB)
82                                 + phongKs*hd^phongSp*lightRGB);
83                  // composite with existing ray color and transparency
84                  rayRGB = rayRGB + rayTransp*alpha*pntRGB;
85                  rayTransp = rayTransp*(1.0 - alpha);
86            }            }
87         }         }
88         if (transp < 0.001) {  // early ray termination         if (rayTransp < 0.01) {  // early ray termination
89  //          rgba = [gray, gray, gray, 1.0];            rayTransp = 0.0;
90            rgba = [or, og, ob, 1.0];            outRGBA = [rayRGB[0], rayRGB[1], rayRGB[2], 1.0 - rayTransp];
91            stabilize;            stabilize;
92         }         }
93         if (rayN > camVspFar) {         if (rayN > camVspFar) {
94            rgba = [gray, gray, gray, 1.0-transp];            outRGBA = [rayRGB[0], rayRGB[1], rayRGB[2], 1.0 - rayTransp];
95            stabilize;            stabilize;
96         }         }
97         rayN = rayN + rayStep;         rayN = rayN + rayStep;

Legend:
Removed from v.866  
changed lines
  Added in v.911

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