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

SCM Repository

[diderot] View of /examples/iso2d-demo/iso2d.diderot
ViewVC logotype

View of /examples/iso2d-demo/iso2d.diderot

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2079 - (download) (annotate)
Mon Nov 5 12:50:51 2012 UTC (6 years, 7 months ago) by jhr
File size: 1583 byte(s)
  Working on demos
// iso2d
//
// Demo of finding isocontours via Newton-Raphson method.
// Initializes positions on a grid, and each update applies one
// step of Newton-Raphson.
//
// Process output with:
// unu jhisto -i iso2d.txt -b 512 512 -min 0 0 -max 1 1 | unu 2op neq - 0 | unu quantize -b 8  -o iso2d.png

int gridSize = 300;
image(2)[] diderot = image("data/ddro-80.nrrd");
field#1(2)[] F = ctmr ⊛ diderot;
input int stepsMax = 50;
input real stepScale = 1.0;
real epsilon = 0.0001;

strand sample (int ui, int vi) {
  // world is 1x1 centered at (0.5, 0.5)
    output vec2 pos = [lerp(0.0, 1.0, -0.5, real(ui), real(gridSize)-0.5),
                       lerp(0.0, 1.0, -0.5, real(vi), real(gridSize)-0.5)];
  // set the isvalue to 50, 30, or 10, depending on whichever we're closest to
    real isoval = 50.0 if F(pos) >= 40.0 
                       else 30.0 if F(pos) >= 20.0
                                 else 10.0;
    int steps = 0;
    update {
      // We bail if we're no longer inside or taken too many steps.
        if (!inside(pos, F) || steps > stepsMax) {
	    die;
        }
        vec2 grad = ∇F(pos);
        if (|grad| == 0.0) {    // can't compute step if |∇F|, so have to bail
	    die;
        }
        vec2 norm = normalize(grad);
        vec2 delta = -((F(pos) - isoval)/|grad|)*norm;  // Newton-Raphson step
        if (|delta| < epsilon) {    // we've converged if step is small enough
	    stabilize;
        }
        pos += stepScale * delta;
        steps += 1;
    }
}

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