## mip.diderot: basic maximum intensity projection (MIP) volume rendering
Maximum intensity projection is the minimalist volume visualization
tool. This implementation is about as short as possible; much of the code is
spent on specifying and setting up the ray geometry. It can also be adapted
(by changing a line or two) to do other kinds of projections. The code for
camera geometry and ray sampling will be re-used verbatim in other volume
rendering programs. This may later become a target for evolving Diderot to
support libraries that contain common functionality.
Note that the explicitly pedagogical nature of this example means that some
unusual usages are required (like recompiling the program multiple times for
different datasets). In addition, the interaction loop of changing the
program parameters via command-line option, generating (on the command-line)
output images, and viewing the images, is very clumsy. Diderot supports
compiling programs to libraries, around which one can generate nicer GUIs,
but these examples are intended to be as stand-alone and minimal as possible.
Just like [`iso2d`](../iso2d) depends on first creating a dataset with
[`fs2d-scl`](../fs2d), we need to create a volume dataset with [`fs3d-scl`](../fs3d)
in order to compile this program `mip.diderot`:
../fs3d/fs3d-scl -which 13 -width 3 -sz0 73 -sz1 73 -sz2 73 | unu save -f nrrd -o cube.nrrd
In this case the volume size is chosen to ensure that the local maxima of
this -which 13 synthetic function, a cube frame with maxima at
(x,y,)=(+-1,+1,+1), are actually hit by grid sample points, which helps
reason about subsequent debugging. We can inspect the volume by tiled slicing
in `unu` to make [cube-tile.png](cube-tile.png):
unu pad -i cube.nrrd -min 0 0 0 -max M M 79 |
unu tile -a 2 0 1 -s 10 8 |
unu quantize -b 8 -o cube-tile.png
These examples are mainly to demonstrate Diderot, but some learning some
`unu` is a useful side-effect. Type `unu` to get a list of all `unu` commands,
and type, for example, `unu pad` to get the usage info for `unu pad`.
We copy `cube.nrrd` to `vol.nrrd`, so that `mip.diderot` can see it in its
default location (in the current directory), we then we can compile
`mip.diderot`:
cp cube.nrrd vol.nrrd
../../vis12/bin/diderotc --exec mip.diderot
And then make some MIP renderings:
./mip -out0 0 -camFOV 20 -rayStep 0.03 -iresU 300 -iresV 300
unu quantize -b 8 -i out.nrrd -o cube-persp.png
./mip -out0 0 -camFOV 20 -rayStep 0.03 -iresU 300 -iresV 300 -camOrtho true
unu quantize -b 8 -i out.nrrd -o cube-ortho.png
Make sure you can run these steps to get the same
[cube-persp.png](cube-persp.png) and [cube-ortho.png](cube-ortho.png). The
only difference in the second image is the use of orthographic instead of
perspective projection, which is apparent in the result. The `-out0 0`
initialization of the MIP accumulator is safe because `unu minmax cube.nrrd`
reports that the minimum value is 0; this avoids an extra `unu 2op exists`
step.
The usual parameters in volume rendering can be set by command-line options.
This sampling of the various options is not exhaustive; you should try more
yourself. For the purposes of this demo we set a (shell) variable to store
repeatedly used options, which can be over-ridden later in the command-line.
OPTS="-out0 0 -camFOV 20 -rayStep 0.03 -iresU 300 -iresV 300"
./mip $OPTS -camAt 1 1 1; unu quantize -b 8 -i out.nrrd -o cube-111.png
./mip $OPTS -camFOV 30; unu quantize -b 8 -i out.nrrd -o cube-wide.png
./mip $OPTS -rayStep 0.3; unu quantize -b 8 -i out.nrrd -o cube-sparse.png
./mip $OPTS -camFar 0; unu quantize -b 8 -i out.nrrd -o cube-farclip.png
./mip $OPTS -camEye 7 7 7; unu quantize -b 8 -i out.nrrd -o cube-eyefar.png
Make sure that these commands generate similar [cube-111.png](cube-111.png),
[cube-wide.png](cube-wide.png), [cube-sparse.png](cube-sparse.png),
[cube-farclip.png](cube-farclip.png), and [cube-eyefar.png](cube-eyefar.png)
for you, and also make sure that you understand why the results look the way they do.
Aside from camera and ray parameters, we now have a basis for testing that
field reconstruction can be invariant with respect to the details of the
sampling grid. We first sample a quadratic function in a few different ways
on a very low resolution grid (which is translating and rotating across these
four datasets):
../fs3d/fs3d-scl -which 4 -sz0 10 -sz1 10 -sz2 10 | unu save -f nrrd -o parab0.nrrd
../fs3d/fs3d-scl -which 4 -sz0 10 -sz1 10 -sz2 10 -angle 30 -off 0.25 0.25 0.25 | unu save -f nrrd -o parab1.nrrd
../fs3d/fs3d-scl -which 4 -sz0 10 -sz1 10 -sz2 10 -angle 60 -off 0.05 0.50 0.50 | unu save -f nrrd -o parab2.nrrd
../fs3d/fs3d-scl -which 4 -sz0 10 -sz1 10 -sz2 10 -angle 90 -off 0.75 0.75 0.75 | unu save -f nrrd -o parab3.nrrd
This demo (of volume rendering with sampling grid invariance) involves a very
unusual invocation of the Diderot compiler and execution of the resulting
program, so the following steps are unfortunately rather clumsy: each
different grid orientation requires a recompilation. First, we use the
program as is, with the
field#0(3)[] F = vol ⊛ tent;
line **uncommented out**, and then run some sh/bash commands:
OPTS="-out0 0 -inSphere true -camOrtho true -camEye 8 0 0 -camFOV 15 -rayStep 0.01 -iresU 300 -iresV 300"
for I in 0 1 2 3; do
echo === $I ====
cp parab${I}.nrrd vol.nrrd
../../vis12/bin/diderotc --exec mip.diderot
./mip $OPTS
unu quantize -b 8 -min 0 -max 1 -i out.nrrd -o parab${I}-tent.png
done
This makes orthographic MIPs of the four low-res parab?.nrrd datasets, producing
images that reveal the underlying sampling grid:
[parab0-tent.png](parab0-tent.png),
[parab1-tent.png](parab1-tent.png),
[parab2-tent.png](parab2-tent.png),
[parab3-tent.png](parab3-tent.png).
The fact that these images differ means that the quadratic function sampled
by `fs3d-scl -which 4` is not exactly reconstructed by the `tent`
reconstruction kernel. If the reconstruction had been exact, then the
rendering would have been entirely determined by the underlying quadratic
function, rather than the particulars of the sampling grid. Note that the
`inSphere` parameter is important here: it ensures that the region rendered
(the sphere) has the same rotational symmetry as the quadratic function
itself.
Now **we change one line of code below**, so that only the
field#1(3)[] F = vol ⊛ ctmr;
field definition is uncommented (commenting out `field#0(3)[] F = vol ⊛ tent;`).
This changes the reconstruction kernel the Catmull-Rom cubic interpolating spline,
which can accurately reconstruct quadratic functions.
Then we re-run the commands above, but saving the results to differently-named images.
OPTS="-out0 0 -inSphere true -camOrtho true -camEye 8 0 0 -camFOV 15 -rayStep 0.01 -iresU 300 -iresV 300"
for I in 0 1 2 3; do
echo === $I ====
cp parab${I}.nrrd vol.nrrd
../../vis12/bin/diderotc --exec mip.diderot
./mip $OPTS
unu quantize -b 8 -min 0 -max 1 -i out.nrrd -o parab${I}-ctmr.png
done
This generates four new images:
[parab0-ctmr.png](parab0-ctmr.png),
[parab1-ctmr.png](parab1-ctmr.png),
[parab2-ctmr.png](parab2-ctmr.png),
[parab3-ctmr.png](parab3-ctmr.png).
As should be visually clear, these images are all the same (or very nearly so), demonstrating the
accuracy of the Catmull-Rom spline for quadratic functions.
Things to try by further modifications of this program:
* Change the rendering from being a maximum intensity, to a summation intensity
(e.g. `out += F(pos)`). Make sure `out0` is zero.
* Change rendering to be mean intensity, which requires adding a counter to
count the number of samples of current ray inside the volume.
* Change the rendering to summating of gradient magnitude
(e.g. `out += |∇F(pos)|`).
* Change the rendering to generate a `vec3` output summing gradient vectors
(e.g. `out += ∇F(pos)`). Type of `out` must be `vec3`, initialized to `[0,0,0]`.