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

SCM Repository

[diderot] Annotation of /tests/glk-tests/back/fsamp.diderot
ViewVC logotype

Annotation of /tests/glk-tests/back/fsamp.diderot

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4640 - (view) (download)

1 : glk 4640
2 :     input int sz0 = 30;
3 :     input int sz1 = 30;
4 :     input int sz2 = 30;
5 :     input real width = 4;
6 :     input int which = 1;
7 :     vec3 spc = [width/(sz0-1), width/(sz1-1), width/(sz2-1)];
8 :    
9 :     input vec3 pp = [0.2,0.2,0.2];
10 :     vec4 qq = normalize([pp[0],pp[1],pp[2],1.0]);
11 :     tensor[3,3] rot = [[qq[0]*qq[0] + qq[1]*qq[1] - qq[2]*qq[2] - qq[3]*qq[3],
12 :     2*(qq[1]*qq[2] - qq[0]*qq[3]),
13 :     2*(qq[1]*qq[3] + qq[0]*qq[2])],
14 :     [2*(qq[1]*qq[2] + qq[0]*qq[3]),
15 :     qq[0]*qq[0] - qq[1]*qq[1] + qq[2]*qq[2] - qq[3]*qq[3],
16 :     2*(qq[2]*qq[3] - qq[0]*qq[1])],
17 :     [2*(qq[1]*qq[3] - qq[0]*qq[2]),
18 :     2*(qq[2]*qq[3] + qq[0]*qq[1]),
19 :     qq[0]*qq[0] - qq[1]*qq[1] - qq[2]*qq[2] + qq[3]*qq[3]]];
20 :     vec3 xdir = rot•[spc[0], 0.0, 0.0];
21 :     vec3 ydir = rot•[0.0, spc[1], 0.0];
22 :     vec3 zdir = rot•[0.0, 0.0, spc[2]];
23 :     vec3 orig = -(xdir*(sz0-1) + ydir*(sz1-1) + zdir*(sz2-1))/2;
24 :    
25 :     function real func(vec3 pos) {
26 :     if (1 == which) {
27 :     return pos•pos/2;
28 :     } else if (2 == which) {
29 :     // a C1 union of a downward parabala (constant hessian)
30 :     // and a gaussian (hessian goes to zer0)
31 :     real sig = 0.7;
32 :     real aa = 0.1209853623/(sig^3);
33 :     real bb = 3*0.1209853623/sig;
34 :     real rr = |pos|;
35 :     real cap = -aa*rr*rr + bb;
36 :     real gss = exp(-rr*rr/(2*sig*sig))/(sqrt(2*π)*sig);
37 :     // the rr < sig conditional won't work, so use lerp instead
38 :     real ww = clamp(0, 1, lerp(0, 1, -1, 100*(rr-sig), 1));
39 :     real ret = (1-ww)*cap + ww*gss;
40 :     return ret;
41 :     } else if (3 == which) {
42 :     real x = pos[0];
43 :     real y = pos[1];
44 :     return atan2(y,x);
45 :     } else {
46 :     return 2*pos[0] + erf(10*(0.8-|pos|));
47 :     }
48 :     }
49 :    
50 :     strand sample(int xi, int yi, int zi) {
51 :     output real vv = 0;
52 :     update {
53 :     if (0 == xi && 0 == yi && 0 == zi) {
54 :     print("NRRD0004\n");
55 :     print("# Complete NRRD file format specification at:\n");
56 :     print("# http://teem.sourceforge.net/nrrd/format.html\n");
57 :     print("type: float\n");
58 :     print("dimension: 3\n");
59 :     print("sizes: ", sz0, " ", sz1, " ", sz2, "\n");
60 :     print("kinds: space space space\n");
61 :     print("endian: little\n");
62 :     print("encoding: raw\n");
63 :     print("space dimension: 3\n");
64 :     print("space directions: (", xdir[0], ",", xdir[1], ",", xdir[2],
65 :     ") (", ydir[0], ",", ydir[1], ",", ydir[2],
66 :     ") (", zdir[0], ",", zdir[1], ",", zdir[2], ")\n");
67 :     print("space origin: (", orig[0], ",", orig[1], ",", orig[2], ")\n");
68 :     print("data file: vv.nrrd\n");
69 :     print("byte skip: -1\n");
70 :     }
71 :     vec3 xyz = orig + xi*xdir + yi*ydir + zi*zdir;
72 :     vv = func(xyz);
73 :     stabilize;
74 :     }
75 :     }
76 :     initially [ sample(xi, yi, zi) | zi in 0..(sz2-1), yi in 0..(sz1-1), xi in 0..(sz0-1) ];

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