/*
 $DDRO_HOME/branches/vis12/bin/diderotc --exec bug026.diderot \
 && ./bug026 \
 && unu jhisto -i x.nrrd -b 800 800 -min 0 0 -max 1 1 \
     | unu 2op neq - 0 \
     | unu quantize -b 8 -o bug026.png
*/

input real isoval = 0.5;
field#1(2)[] F = ctmr ⊛ image("../data/ddro.nrrd") - isoval;

int grid = 150;
int stepsMax = 10;
real epsilon = 0.000001;

strand RootFind(vec2 x0) {
   output vec2 x = x0;
   int steps = 0;
   update {
      // Stop if we're no longer inside or taken too many steps.
      if (!inside(x, F) || steps > stepsMax)
         die;
      /* BUG: this code generates
         error converting OffsetField
         uncaught exception Fail [Fail: bogus operator OffsetField]
           raised at common/phase-timer.sml:76.50-76.52
           raised at common/phase-timer.sml:76.50-76.52
           raised at high-to-mid/high-to-mid.sml:297.85-297.87
           raised at high-to-mid/high-to-mid.sml:294.32-294.78
        The code can be fixed by replacing all "∇F(x)" with uses of a
        new variable "vec2 grad = ∇F(x)", but why should that be needed? */
      if (|∇F(x)| == 0.0)
         die;
      // the Newton-Raphson step
      vec2 delta = normalize(∇F(x)) * F(x)/|∇F(x)|;
      // we've converged if the change is small enough
      if (|delta| < epsilon)
         stabilize;
      x -= delta;
      steps += 1;
   }
}

initially { RootFind([lerp(0, 1, -0.5, ui, grid-0.5),
                      lerp(0, 1, -0.5, vi, grid-0.5)])
            | vi in 0..(grid-1), ui in 0..(grid-1) };

/*
 different in 3D:
 vec2 -> vec3
 field#1(2)[] F -> field#1(3)[] F
 initially over 3D grid
*/

<script type="text/javascript">//<![CDATA[
		function toggle_ffErrors() {
			var errorsblock = document.getElementById("ffErrorsBlock");
			var errorsgroup = document.getElementById("ffErrors");
			if (errorsblock.style.display == "none") {
				errorsblock.style.display = "block";
				errorsgroup.style.right = "10px";
			} else {
				errorsblock.style.display = "none";
				errorsgroup.style.right = "300px";
			}
		}
//]]></script>
<div id="ffErrors">
<a href="javascript:toggle_ffErrors();">Click to toggle</a>
<div id="ffErrorsBlock">
	<div class="error">does not end with &lt;/html&gt; tag</div>
	<div class="error">does not end with &lt;/body&gt; tag</div>
	<div class="info">The output has ended thus:  };

/*
 different in 3D:
 vec2 -&gt; vec3
 field#1(2)[] F -&gt; field#1(3)[] F
 initially over 3D grid
*/</div>
</div></div>
</body></html>