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

SCM Repository

[diderot] Annotation of /branches/vis12/test/lzc2d-brute-BUG1.diderot
ViewVC logotype

Annotation of /branches/vis12/test/lzc2d-brute-BUG1.diderot

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1953 - (view) (download)

1 : glk 1953 // lzc2d-brute
2 :     //
3 :     // Demo of laplacian zero-crossing edges via brute-force search
4 :     // and then Illinois false position
5 :     //
6 :     // Process output with:
7 :     // unu jhisto -i pos.nrrd -b 512 512 -min 0 0 -max 1 1 | unu 2op neq - 0 | unu quantize -b 8 -o lzc2d.png
8 :    
9 :     int gridSize = 300;
10 :     field#2(2)[] F = c4hexic ⊛ image("../data/ddro-64.nrrd");
11 :     int stepsMax = 30;
12 :     real epsilon = 0.00001;
13 :    
14 :     strand sample (int ui, int vi) {
15 :     output vec2 pos = [lerp(0.0, 1.0, -0.5, real(ui), real(gridSize)-0.5),
16 :     lerp(0.0, 1.0, -0.5, real(vi), real(gridSize)-0.5)];
17 :     vec2 posA=pos;
18 :     vec2 posB=pos;
19 :     int stepsNum = 0;
20 :     int side = 0;
21 :     real step = 0.004;
22 :     real fA=0.0;
23 :     real fB=0.0;
24 :     real A=0.0;
25 :     real B=1.0;
26 :     bool searching = true;
27 :     update {
28 :     // We bail if we're no longer inside or taken too many steps.
29 :     if (!inside(pos, F) || stepsNum > stepsMax) {
30 :     die;
31 :     }
32 :     if (searching) {
33 :     posA = pos;
34 :     real ll = trace(∇⊗∇F(pos));
35 :     vec2 dir = normalize(∇F(pos));
36 :     real sgn = 1.0 if ll > 0.0 else -1.0;
37 :     posB = pos + sgn*step*dir;
38 :     searching = trace(∇⊗∇F(posB))*ll > 0.0;
39 :     if (!searching) {
40 :     /* set up next phase */
41 :     fA = trace(∇⊗∇F(posA));
42 :     fB = trace(∇⊗∇F(posB));
43 :     }
44 :     } else {
45 :     pos = lerp(posA, posB, fA, 0.0, fB);
46 :     real fs = trace(∇⊗∇F(pos));
47 :     if (0.0 == fs) {
48 :     stabilize;
49 :     }
50 :     if (fs*fB > 0.0) { /* not between s and b */
51 :     posB = pos;
52 :     fB = fs;
53 :     if (1 == side) {
54 :     fA /= 2.0;
55 :     }
56 :     side = 1;
57 :     } else { /* not between a and s */
58 :     posA = pos;
59 :     fA = fs;
60 :     if (-1 == side) {
61 :     fB /= 2.0;
62 :     }
63 :     side = -1;
64 :     }
65 :     if (|posA - posB| < epsilon ) {
66 :     /* BUG: is this the long-standing well-known control-flow graph bug?
67 :     compiling generates:
68 :    
69 :     error in cvtMethod(Update, ...)
70 :     uncaught exception Fail [Fail: unimplemented]
71 :     raised at common/phase-timer.sml:76.50-76.52
72 :     raised at translate/translate.sml:384.110-384.112
73 :     raised at translate/translate.sml:324.29-324.31
74 :     raised at translate/translate.sml:324.29-324.31
75 :     raised at translate/translate.sml:324.29-324.31
76 :     raised at translate/translate.sml:297.29-297.49
77 :    
78 :     but this is a place where the choice between stabilize and die
79 :     (and one of the two will certainly be called) depends on a test.
80 :     Can be re-written without an else branch but seems cleaner this way.
81 :     */
82 :     if (|∇F(pos)| > 0.0001) {
83 :     stabilize;
84 :     } else {
85 :     die;
86 :     }
87 :     }
88 :     }
89 :     stepsNum += 1;
90 :     }
91 :     }
92 :    
93 :     initially { sample(ui, vi) | vi in 0..(gridSize-1), ui in 0..(gridSize-1) };

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