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

SCM Repository

[diderot] Annotation of /trunk/test/lic-turb2d.diderot
ViewVC logotype

Annotation of /trunk/test/lic-turb2d.diderot

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2636 - (view) (download)

1 : glk 1391 // lic-turb.diderot
2 :     //
3 :     // demo of line integral convolution (LIC) on 2D turbulent flow
4 :     //
5 :     // First, need to create padded vector dataset v.nrrd with:
6 :     // unu pad -i ../data/turb2d.nrrd -min 0 -5 -5 -max M M+5 M+5 -b pad -v 0.0001 -o v.nrrd
7 :     // We pad with value 0.0001; we padded with zero than normalize() would fail
8 :     // But if you try take too many steps (increasing stepNum) then eventually
9 :     // you will leave the vector field domain, which will crash the program,
10 :     // so really the padding should be done with a vector field that
11 :     // pointwards inwards, towards the center of the domain, so that anything
12 :     // leaving the domain will get pushed back in. OR, really, the code
13 :     // below should be smart and use inside() or some other test to stop
14 :     // the integration when we get too close to the domain boundary.
15 :     //
16 :     // Then, need to create noise texture with:
17 :     // unu slice -i ../data/turb2d.nrrd -a 0 -p 0 | unu resample -s 1020 561 | unu 1op nrand -o R.nrrd
18 :     // where "1020 561" is copied from imgSizeX and imgSizeY below; we start
19 :     // with unu slice in order to inherit the image orientation info, to
20 :     // create a noise texture at the same world-space position as the data
21 :     //
22 :     // process output with:
23 :     // unu reshape -i mip.txt -s 1020 561 | unu quantize -b 8 -o lic.png
24 :    
25 :    
26 :     int imgSizeX = 1020;
27 :     int imgSizeY = 561;
28 :     real h = 0.005; // step size of integration
29 :     int stepNum = 10; // take this many steps both upstream and downstream
30 :    
31 : jhr 2636 field#1(2)[2] V = image("v.nrrd") ⊛ ctmr;
32 :     field#0(2)[] R = image("R.nrrd") ⊛ tent;
33 : glk 1391
34 :     strand LIC (int xi, int yi) {
35 :     real xx = lerp(0.0, 6.78, 0.0, real(xi), real(imgSizeX)-1.0);
36 :     real yy = lerp(0.0, 3.74, 0.0, real(yi), real(imgSizeY)-1.0);
37 :     vec2 pos0 = [xx,yy];
38 :     vec2 forw = pos0;
39 :     vec2 back = pos0;
40 :     output real sum = R(pos0);
41 :     int step = 0;
42 :    
43 :     update {
44 :     // Euler integration step
45 :     // forw = forw + h*V(forw);
46 :     // back = back - h*V(back);
47 :     // Midpoint method step
48 :     forw += h*normalize(V(forw + 0.5*h*normalize(V(forw))));
49 :     back -= h*normalize(V(back - 0.5*h*normalize(V(back))));
50 :     sum += R(forw) + R(back);
51 :     step += 1;
52 :     if (step == stepNum) {
53 :     // the pow() is a way to adjust how much velocity modulates
54 :     // contrast; lower values increase contrast at small velocity
55 :     sum = pow(|V(pos0)|,0.6)*sum/real(1 + 2*stepNum);
56 :     stabilize;
57 :     }
58 :     }
59 :     }
60 :    
61 :     initially [ LIC(xi, yi) | yi in 0..(imgSizeY-1), xi in 0..(imgSizeX-1) ];

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