Skip to content

Latest commit

 

History

History

halftone

halftone.diderot: Half-toning via interacting repulsive particles

This example is based on the sphere.diderot example, which in turn is based on the circle.diderot example; both of these should be read and understood first. While in those cases a surface (circle or sphere) is sampled uniformly, in this case the flat surface of an image is sampled non-uniformly, according to the image intensity. This creates a result similar to the Electrostatic Halftoning method of Schmaltz et al. The very regular hexagonal grid patterns that this program generates by energy minimization in uniform areas may be sub-optimal for artistic or signal-processing purposes; this program lacks the jittering that better haltoning methods add to avoid this. The regular sampling generated by this program could help subsequent meshing.

The standard test for this kind of program is a linear ramp, which is available via:

../fs2d/fs2d-scl -which x -width 2 -size0 401 -size1 401 |
  unu crop -min 0 100 -max M M-100 |
  unu affine -1 - 1 0 1 |
  unu pad -min -2 -2 -max M+2 M+2 -o img.nrrd
rm -f out.nrrd

The domain of this img.nrrd is [-1,1] along X, and [-0.5,0.5] along Y, and which varies from 0 at X=-1 to 1 at X=1. Two extra samples are added at the boundaries to ensure that even at the boundary, the inside() test passes. With this proxy image in place, we can compile:

diderotc --snapshot --exec halftone.diderot

We use snapshots to monitor the progress of computation. To make NN initial positions with random number seen RNG that fit within the ramp image domain:

NN=300
RNG=5
echo 0 0 | unu pad -min 0 0 -max M $((NN-1)) |
  unu 1op rand -s $RNG | unu affine 0 - 1 -1 1 -o vec2.nrrd
echo 1 0.5 | unu 2op x vec2.nrrd - -o vec2.nrrd

Then to run with snapshots saved every iteration (-s 1), but limiting the program to 800 iterations (-l 800), as well as cleaning results from previous run:

rm -f pos-????.{png,nrrd} pos.nrrd
./halftone -s 2 -l 800 -radmm 0.04 1 -eps 0.0001 -pcp 2

Running with this large value (0.04) of minimum radius (considering the image domain) is good for giving a visual impression of how the particle system populates the domain. Next, some unu hacking makes images of the evolving system, which are then turned into an animated ramp.gif (compare to `ramp-ref.gif).

SZ=200
OV=2
export NRRD_STATE_VERBOSE_IO=0
for PIIN in pos-????.nrrd; do
IIN=${PIIN#*-}; II=${IIN%.*}
   echo "post-processing $PIIN to pos-$II.png ... "
   unu jhisto -i $PIIN -min -1 -0.5 -max 1 0.5 -b $((OV*SZ*2)) $((SZ*OV)) |
     unu resample -s /$OV /$OV -k bspln5 -t float |
     unu quantize -b 8 -min 0 -max $(echo "0.15 / ($OV * $OV)" | bc -l) -o pos-$II.png
done
convert -delay 6 pos-*.png ramp.gif

On the other hand, to quantitatively check that the particle density is as it should be, we run with many more particles (which can take a few minutes to finish). The following uses -radmm 0.004 1 (versus -radmm 0.04 1 above) for the particle computation, so that inter-particle distance may be a tenth of what it was previously (resulting in up to 100 times greater particle density). A histogram of the X positions, plotted above the rasterized image of the particle positions, provides visual confirmation that horizontal (X position) particle density approaches the linear ramp expected from the linear ramp of the underlying image.

rm -f {hp,pos}-????.{png,nrrd} pos.nrrd
./halftone -s 10 -l 800 -radmm 0.004 1 -eps 0.00004 -pcp 2
SZ=200
OV=2
export NRRD_STATE_VERBOSE_IO=0
for PIIN in pos-????.nrrd; do
   IIN=${PIIN#*-}; II=${IIN%.*}
   echo "post-processing $PIIN to pos-$II.png ... "
   unu jhisto -i $PIIN -min -1 -0.5 -max 1 0.5 -b $((OV*SZ*2)) $((SZ*OV)) |
     unu resample -s /$OV /$OV -k bspln3 -t float |
unu quantize -b 8 -min 0 -max $(echo "2 / ($OV * $OV)" | bc -l) -o pos-$II.png
   unu slice -i $PIIN -a 0 -p 0 |
     unu histo -min -1 -max 1 -b $((SZ/3)) |
     unu dhisto -h $((SZ/3)) -nolog |
     unu resample -s $((SZ*2)) = -k box |
     unu join -i - pos-$II.png -a 1 -o hp-$II.png # combine histogram and position image
done
convert -delay 6 hp-*.png hp.gif

The top part of the resulting hp.gif (compare to hp-ref.gif) shows the expected convergence to a linear variation in particle density, even though there are too many particles to individually distinguish in the rasterized image.

Finally, to have some fun with a picture of Diderot himself, we invert the image intensity and darken a bit (tending to have more space between particles), recompile with the new proxy image (the array axis sizes changed), and then run with new initial positions (without snapshots this time):

unu 2op - 1 ../data/ddro.nrrd | unu gamma -g 0.75 -o img.nrrd
diderotc --exec halftone.diderot
NN=1000
RNG=5
echo 0 0 | unu pad -min 0 0 -max M $((NN-1)) |
  unu 1op rand -s $RNG | unu affine 0 - 1 -1 1 -o vec2.nrrd
./halftone -l 800 -radmm 0.01 0.2 -eps 0.000009 -pcp 1
echo 1 -1 |
unu 2op x pos.nrrd - |
unu save -f text |
sed -e s/$/\ 0.005\ 0\ 360\ arc\ closepath\ fill/ |
cat head.eps - tail.eps > ddro.eps
epstopdf ddro.eps

The ddro.eps is generated with the help of head.eps and tail.eps, which assume that the points live inside the domain [-1,1]x[-1,1]. The final command for generating ddro.pdf is the epstopdf that comes with some LaTeX installations, but there are many other ways to convert from eps to pdf. The ddro.pdf result (compare to drro-ref.pdf) should look something like the ddro.png produced by

unu quantize -b 8 -i ../data/ddro.nrrd -o ddro.png

Note that the minimum particle radius specified to -radmm was 0.01, hence at their tightest packing particles were 0.01 away from each other (they sat in each other's potential wells, at radius 0.01). Drawing each point with a circle of radius 0.005 (via sed above) creates, where the field is highest (where ddro.nrrd is darkest), a dense packing of mutually tangent circles, which is visible in the output.