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

SCM Repository

[diderot] Annotation of /branches/charisee/ertest/vispaper/iso3d/giso3d.diderot
ViewVC logotype

Annotation of /branches/charisee/ertest/vispaper/iso3d/giso3d.diderot

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2938 - (view) (download)

1 : glk 2937 /*
2 :    
3 :     ../../../bin/diderotc giso3d.diderot \
4 :     && ./giso3d \
5 :     && unu crop -i giso3d.txt -min 0 0 -max 1 M | \
6 :     unu jhisto -b 600 600 -min -1.2 -1.2 -max 1.2 1.2 -t float | \
7 :     unu 1op log1p | \
8 :     unu quantize -b 8 -o giso3d.png \
9 :     && open giso3d.png
10 :    
11 :     */
12 :    
13 :     field#4(3)[] V= c4hexic ⊛ load("parab0.nrrd");
14 :    
15 :     field#3(3)[3] N = ∇V/|∇V|;
16 :     field#3(3)[3,3] P = identity[3] - N⊗N;
17 :     field#2(3)[3,3] H = ∇⊗∇V;
18 :     field#3(3)[3,3] G = -(P•H•P)/|∇V|;
19 :    
20 : glk 2938 // compilation stalls out with either of these
21 : glk 2937 field#2(3)[] F = trace(G) - -3;
22 : glk 2938 //field#2(3)[] F = G:G - 9;
23 : glk 2937
24 :     // (to try later after this is fixed: isosurfaces of flowline curvature)
25 :    
26 :     int grid = 30;
27 :     int stepsMax = 10;
28 :     real epsilon = 0.000001;
29 :    
30 :     strand RootFind(real xi,real xj, real xk) {
31 :     output vec3 x = [xi,xj,xk];
32 :     int steps = 0;
33 :     update {
34 :     // Stop if we're no longer inside or taken too many steps.
35 :     if (!inside(x, V) || steps > stepsMax)
36 :     die;
37 :     // subsequent expressions are undefined if |∇F| is zero
38 :     if (|∇F(x)| == 0.0)
39 :     die;
40 :     // the Newton-Raphson step
41 :     vec3 delta = normalize(∇F(x)) * F(x)/|∇F(x)|;
42 :     // we've converged if the change is small enough
43 :     if (|delta| < epsilon)
44 :     stabilize;
45 :     x -= delta;
46 :     steps += 1;
47 :     }
48 :     }
49 :    
50 :     initially { RootFind(lerp(-2, 2, -0.5, ui, grid-0.5),
51 :     lerp(-2, 2, -0.5, vi, grid-0.5),
52 :     lerp(-2, 2, -0.5, wi, grid-0.5))
53 :     | wi in 0..(grid-1), vi in 0..(grid-1), ui in 0..(grid-1) };

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