# How to make an inkblot?

How to effectively create a polygon that looks like a realistic inkblot? So far, I could come up with this (borrowing from Ed Pegg Jr.'s Rorschach demonstration):

RandomBlot[num_, opts___] := Module[{pts},
pts = RandomReal[{0, 1}, {num, 2}];
pts = pts[[Append[Last@FindShortestTour[pts, Method -> "TwoOpt"], 1]]];
Graphics[{Polygon@
Table[BSplineFunction[pts, SplineKnots -> "Clamped"][a], {a, 0, 1, 0.001}]}, opts]
];

blot = RandomBlot[40, ImageSize -> 200, AspectRatio -> 1]


There are two problems with this:

• it is a bit slow due to FindShortestTour
• the blot contains corners being too sharp and has 'bays' reaching too far inward.

Compare it to a desired outcome:

Update:

I should mention that it is not necessary to actually create a Rorschach-like image, so mirroring is not a requirement.

-
Is it important to have concave parts of in the shape? Without concave parts it's easy ... – Szabolcs Mar 22 '12 at 15:11
@Szabolcs, I noticed. ComputationalGeometryConvexHull does it quite well. – rcollyer Mar 22 '12 at 15:16
@Szabolcs I would say: let it be a parameter of the model, that specifies the amount of concavity/convexity. If the valu of this parameter is zero, the output could be an ellipse. – István Zachar Mar 22 '12 at 15:19
Can anyone post photographs of actual ink blots on paper, so we can judge what is realistic? – Colonel Panic Mar 23 '12 at 0:25
There's also the ink blot demonstration by Andy Ross link – jrhodin Mar 23 '12 at 21:16

A bit of image processing:

Table[
Blur[
Dilation[
Graphics[
Table[
Rotate[
Disk[RandomReal[{-10, 10}, {2}], {RandomReal[{1, 5}],RandomReal[{1, 5}]}],
RandomReal[{0, 3.14}]
],
{40}
]
],
DiskMatrix[20]
], 20
]// Binarize,
{3}, {3}
] // Grid


Lots of parameters to play with...

Now these are bitmaps and if vector graphics are required (the question seems to imply that) we can adapt a bit of Vitaly's code from here:

img = Thinning@EdgeDetect@p;
points = N@Position[ImageData[img], 1];
pts = FindCurvePath[points] /. c_Integer :> points[[c]];
Graphics[{EdgeForm[Directive[Dashed, Thick, Red]],FilledCurve@({Line@#} & /@ pts)}]


with p our blob bitmap. (The contour is dashed to better show that we're dealing with vector graphics here).

-
Nice, these inkblots look more realistic than the others – R. M. Mar 22 '12 at 21:12
Ah, I see fish bones ... What does that mean? (+1) – rcollyer Mar 22 '12 at 21:25
I love parameters! – István Zachar Mar 22 '12 at 21:58

While I enjoyed fiddling with each of the beautiful solutions you gave, I chose two that are the closest to what I needed in form, splatter-distribution, parameterization and speed. As a token of my appreciation I've reworked them into a dynamic demo, showcasing Szabolcs's and Sjoerd's solutions. This does not mean that the other solutions could not be included: I think all of them could be easily extended to comply with the specified parameters. I simply don't have more time. But if anyone feels like doing it, please go ahead, and edit this post!

Both methods are wrapped in a smoothing function (Blur & Binarize), and then in a "fractalization" function that detects the edge and applies some noise to it in the form of black and white disks (idea coming from Sjoerd's solution). This can be done recursively, with disks of smaller and smaller sizes, adding more subtle details.

Options[RandomBlot] =
Join[Options@Graphics, {RandomSeed -> Automatic, Elevation -> 2,
EdgeRecursion -> 2, EdgeResolution -> 300, EdgeSmoothing -> 7,
Method -> Automatic}];

RandomBlot[bulk_: .1, pat_: 0, smo_: .2, opts : OptionsPattern[]] :=
Module[{ratio, seed, ele, rec, res, rad, method, edgeNoise, range},
ratio = OptionValue@AspectRatio /. Automatic -> 1;
ele = OptionValue@Elevation /. Automatic -> 2;
rec = OptionValue@EdgeRecursion /. Automatic -> 2;
res = OptionValue@EdgeResolution /. Automatic -> 300;
rad = OptionValue@EdgeSmoothing /. Automatic -> 5;
method = OptionValue@Method /. Automatic -> "Szabolcs";
seed =
OptionValue@RandomSeed /.
Automatic -> RandomInteger@{0, 99999999999999};

edgeNoise[img_, lev_, num_: 300] :=
Module[{pt =
N@Position[
Reverse[Transpose@ImageData@Thinning@EdgeDetect@img, {2}],
1], new},
pt = Take[RandomSample@pt, Min[num, Length@pt]];
new =
Show[img,
Graphics[({RandomChoice@{Black, White},
Rotate[Disk[#, RandomReal[{0, 15/lev}, {2}]],
RandomReal@{0, \[Pi]}]} & /@ pt)]];
Blur[new, rad] // Binarize];

BlockRandom[SeedRandom@seed;
Show[
Fold[edgeNoise[#1, #2, res] &,
Blur[Switch[method,
"Szabolcs", range = 10;
With[{fun = Exp[-Round[pat*100] #.#] &,
pts = Transpose@{RandomReal[{-range, range}, {Round[100*bulk]}],
RandomReal[{-range, range}*ratio, {Round[100*bulk]}]}},
With[{fc = Compile[{xl, yl},
Total[fun[# - {xl, yl}] & /@ (pts*.9)] > 1/ele]},
RegionPlot[
fc[x, y], {x, -range, range}, {y, -range*ratio,
range*ratio}, PlotStyle -> Black, BoundaryStyle -> Black,
Frame -> False]]],

"Sjoerd",
Dilation[Graphics[{
Black,
Table[Rotate[
Disk[{RandomReal@{-10, 10},
RandomReal@({-10, 10}*ratio)},
RandomReal[{.1, 4}, {2}]],
RandomReal@{0, \[Pi]}], {Round[100*bulk]}],
White,
Table[Rotate[
Disk[{RandomReal@{-10, 10},
RandomReal@({-10, 10}*ratio)},
RandomReal[{.1, 2}, {2}]],
RandomReal@{0, \[Pi]}], {Round[100*bulk*pat]}]
}], DiskMatrix@ele]
], 100*smo] // Binarize,
Range@rec],
FilterRules[{opts}, Options@Graphics]
]]];

Manipulate[
RandomBlot[bulk, pat, smo, AspectRatio -> ratio, RandomSeed -> seed,
ImageSize -> size, Elevation -> ele, EdgeRecursion -> rec,
EdgeResolution -> res, EdgeSmoothing -> rad, Method -> method],
{{seed, 0},
Button["randomize", seed = RandomInteger@{0, 99999999999999}] &},
{{method, "Szabolcs",
"method"}, {"Szabolcs" -> "RegionPlot (Szabolcs)",
"Sjoerd" -> "Disks (Sjoerd)"}},
Delimiter,
{{bulk, .1, "bulkiness"}, 0, 3, Appearance -> "Labeled"},
{{pat, .05, "patchiness"}, 0, 1, Appearance -> "Labeled"},
{{smo, .6, "smoothness"}, 0, 1, Appearance -> "Labeled"},
{{ele, 2, "elevation"}, 0, 10, Appearance -> "Labeled"},
Delimiter,
{{rec, 0, "edge recursion"}, 0, 3, 1, Appearance -> "Labeled"},
{{res, 300, "edge resolution"}, 0, 600, 10, Appearance -> "Labeled"},
{{rad, 7, "edge smoothness"}, 0, 30, Appearance -> "Labeled"},
Delimiter,
{{ratio, 1, "ratio"}, 0, 1, Appearance -> "Labeled"},
{{size, 300, "size"}, 100, 1000, 1, Appearance -> "Labeled"}
]


A collection of blots:

-

Hopefully, this entry is not too late. As with Heike, I also work with Perlin noise, but here I restrict myself to the one-dimensional version, and consider what happens when I treat it as a polar function:

fBm = With[{permutations =
Apply[Join, ConstantArray[RandomSample[Range[0, 255]], 2]]},
Compile[{{x, _Real}},
Module[{xf = Floor[x], xi, xa, u, i, j},
xi = Mod[xf, 16] + 1;
xa = x - xf; u = xa*xa*xa*(10.0 + xa*(xa*6.0 - 15.0));
i = permutations[[permutations[[xi]] + 1]];
j = permutations[[permutations[[xi + 1]] + 1]];
(2 Boole[OddQ[i]] - 1)*xa*(1.0 - u) +
(2 Boole[OddQ[j]] - 1)*(xa - 1)*u], "CompilationTarget" -> "WVM",
RuntimeAttributes -> {Listable}]];

With[{a = 1, b = 3/4, c = 30, d = 15},
Graphics[Cases[
ParametricPlot[(a + b fBm[c + d t/(2 Pi)]) {Cos[t], Sin[t]}, {t,
0, 2 Pi}], _Line, Infinity] /. Line -> Polygon]]


Here are a few more I got:

Play around with various values of a,b,c,d (and maybe tweak SeedRandom[] while you're at it) and see what you get.

Here's how to generate an entire zoo of them:

Graphics[Table[
With[{h = RandomReal[100], k = RandomReal[100],
a = RandomInteger[{1, 6}], b = RandomInteger[{1, 6}],
c = RandomInteger[255], d = RandomInteger[{3, 30}],
n = RandomInteger[{3, 30}]},
Cases[ParametricPlot[
{h, k} + (a + b fBm[c + d t/(2 Pi)]) {Cos[t], Sin[t]},
{t, 0, 2 Pi}], _Line, Infinity] /. Line -> Polygon], {30}]]


If one wants inkblots with inherent bilateral symmetry, it is a simple matter to modify the code for generating them:

With[{a = 4, b = 3, c = 4, d = 6, n = 1},
Graphics[Cases[
ParametricPlot[(a + b fBm[c + d Sin[n t]]) {Cos[t], Sin[t]}, {t,
0, 2 Pi}], _Line, Infinity] /. Line -> Polygon]]


Here are a few more:

Again, tweak the parameters to taste.

Here, we do something similar to cormullion's ink splotches:

Graphics[Table[
With[{h = RandomReal[100], k = RandomReal[100],
a = RandomInteger[{4, 6}], b = RandomInteger[{1, 3}],
c = RandomInteger[255], d = RandomInteger[{3, 30}],
n = RandomInteger[{3, 30}]},
Cases[ParametricPlot[
{h, k} + (a + b fBm[c + d Sin[n t]]) {Cos[t], Sin[t]},
{t, 0, 2 Pi}], _Line, Infinity] /. Line -> Polygon], {30}]]


-
Never too late, and vector-graphics are always welcome as they are usually much faster than image-processed solutions. – István Zachar May 2 '12 at 7:30
You already got my vote, what more do you want? :) – rcollyer May 15 '12 at 14:44
@rcollyer: nothing else, really. :) But I do have a few more inkblot-related ideas on the pipeline... which I will add to this answer later when I've polished them up. – J. M. May 15 '12 at 14:48
Your last figure looks like plankton. – Szabolcs Jun 2 '12 at 10:10
@Szabolcs Or viruses. – Mechanical snail Dec 2 '12 at 1:08

This approach is based on a random walk of a shrinking disk. Several of these are combined and a Gaussian filter is used to smooth it out. Optionally the smoothed image can be multiplied by the original to restore the tiny "droplets" that are wiped out by the smoothing. There is a streakiness parameter which biases the random walk in a particular direction.

randomstep := RandomReal[{0,1}] Through[{Cos,Sin}[RandomReal[{0,2Pi}]]];

rndwalk[numpts_, streakiness_, numruns_] := Module[{streak}, Table[
streak = streakiness randomstep;
RandomChoice[{Identity, Reverse}]@
NestList[# + streak + 0.1 randomstep &, randomstep, numpts]
, {numruns}]];

Graphics[
Range[(Length@# - 1), 0, -1]/(10. (Length@# - 1))]] & /@
points], 50, 1];

imageprocess[pic_, filterwidth_, threshold_, droplets_, reflect_] :=
Module[{smoothed, combined},
smoothed = Binarize[GaussianFilter[pic, filterwidth], threshold];
combined = If[droplets, ImageMultiply[smoothed, pic], smoothed];
If[reflect, ImageMultiply[combined, ImageReflect[combined, Left]],
combined]];

Manipulate[
SeedRandom[seed];
imageprocess[spatter[rndwalk[numpts, streakiness, numspatters]],
filterwidth, threshold, droplets, reflect],
{{seed, 0}, 0, 10^6, 1},
{{numpts, 100}, 10, 300, 1},
{{streakiness, 0}, 0, 0.05},
{{numspatters, 10}, 1, 20, 1},
{{filterwidth, 10}, 1, 20},
{{threshold, 0.6}, 0, 1},
{{droplets, True}, {True, False}},
{{reflect, True}, {True, False}}]


-
Very realistic and fast, I love it! – István Zachar Apr 14 '12 at 9:34
+1, those are awesome! You get a lot of alien looking pictures, though. – rcollyer May 2 '12 at 3:49
This is amazing! – Eli Lansey May 2 '12 at 13:37
@rcollyer, I'm glad I'm not the only one seeing aliens! I wasn't sure if it was some deep-seated psychological issue... – Simon Woods May 2 '12 at 20:47
Here's two I made that I particularly liked: i.stack.imgur.com/aQjWa.png and i.stack.imgur.com/2ljZO.png – rcollyer May 2 '12 at 21:48

This solution uses Perlin noise to generate the blobs.

To generate the noise we use the following function. Here, range is the domain on which we will generate the noise, res is the number of points in x and y direction, and seed is the seed for the random number generator.

perlinNoise[range_: {{0, 1}, {0, 1}}, res_: {30, 30}, seed_: 1] :=
With[{
grid = (SeedRandom[seed]; Table[With[{t = RandomReal[2 Pi]},
{Cos[t], Sin[t]}], {res[[1]]}, {res[[2]]}]),
intf = 3 #^2 - 2 #^3 &},
Function[{x, y},
Module[{ind, xr, yr, pmat},
(* ind=={{xmin,xmax},{ymin,
ymax}} *)
{xr, yr} =
MapThread[Rescale[#1, #2, {1, #3}] &, {{x, y}, range, res}];
ind = Floor[{xr, yr}];
ind = MapThread[Min[#1, #2 - 1] &, {ind, res}];
{xr, yr} -= ind;

pmat = {{{xr, yr}.grid[[ind[[1]], ind[[2]]]],
{xr, yr - 1}.grid[[ind[[1]], ind[[2]] + 1]]},
{{xr - 1, yr}.grid[[ind[[1]] + 1, ind[[2]]]],
{xr - 1, yr - 1}.grid[[ind[[1]] + 1, ind[[2]] + 1]]}};

{1 - intf[xr], intf[xr]}.pmat.{1 - intf[yr], intf[yr]}]]]


Next we're going to superpose some noise on top of a two-dimensional Gaussian function and we use ListContourPlot to plot the region where the resulting function is larger than some level (note that we could use RegionPlot as well, but ListContourPlot is faster). There are a lot of parameters to play with here, such as the resolution of the noise, the ratio between the noise and Gaussian surface, and the level of the contour. For example for

res = 30; seed = 3; level = .6; ratio = .15;


we get

f = perlinNoise[{{-1, 1}, {-1, 1}}, {res, res}, seed];
tab1 = Table[Exp[-(x^2 + y^2)] + ratio*f[x, y], {x, -1, 1, 1/(2 res)}, {y, -1, 1,
1/(2 res)}];
pl = ListContourPlot[ArrayPad[tab1, {{1, 1}, {1, 1}}],
InterpolationOrder -> 2, Contours -> {level}, ContourShading -> None,
Frame -> False];
pl/. {___, a__Line} :> FilledCurve[Thread[{{a}}]]


By increasing the ratio get more splattering:

res = 30; seed = 3; level = .6; ratio = .8;


And by lowering the resolution you get smoother blobs:

res = 8; seed = 3; level = .6; ratio = .8;


-
For the last one, I see a butterflied chicken, or possibly a face. Guess that indicates that I'm crazy. – rcollyer Mar 22 '12 at 23:53
@rcollyer First fish bones, now a butterflied chicken. Maybe you're just hungry. – Heike Mar 22 '12 at 23:57
Not anymore, just had pie, or if you prefer Pi/4. – rcollyer Mar 22 '12 at 23:58
There's a Perlin noise example on the doc page of Compile, if you're interested – Szabolcs Mar 23 '12 at 9:26

This solution is based on a generation of blots that are deformations of a circle:

drop[
tan_ /; NumberQ[N[tan]] && NonNegative[tan],
n_Integer /; Positive[n], num_Integer /; Positive[num]] :=
tangentialAmplitudes, tmp},
radialDirection[phi_] := N[{Cos[phi], Sin[phi]}];
tangentialDirection[phi_] := N[{-Sin[phi], Cos[phi]}];
radialAmplitudes = Table[RandomReal[{-1, 1}], {n}];
tangentialAmplitudes = Table[RandomReal[{-1, 1}], {n}];
tmp = Table[radialDirection[phi] + Apply[Plus,
radialDirection[# phi]/#^2] &, Range[3, n + 2]]] + Apply[Plus,
Map[N[tan*tangentialAmplitudes[[# - 2]]*
tangentialDirection[# phi]/#^2] &, Range[3, n + 2]]],
{phi, 0, N[2 Pi - Pi/num], N[2 Pi/num]}]; Append[tmp, First[tmp]]
];


Then ink splatters can be generated in this way:

inksplatter = Image[Graphics[
Table[With[{loc = {RandomReal[{-5, 5}], RandomReal[{-5, 5}]}},
Scale[Polygon[Map[(# + loc) &,
drop[RandomReal[{.5, 1.5}], RandomReal[{.5, 1.5}], 10, 100]]],
RandomReal[{.1, 1.2}]]
], {20}], PlotRange -> {{-8, 8}, {-8, 8}}, AspectRatio -> 1]]


that produces this:

and a nice inkblot:

ImageAdjust[ImageMultiply[inksplatter, ImageReflect[inksplatter, Left]], {0, .2}]


-
The last one's definitely a panda with mutton chops! – R. M. Mar 23 '12 at 13:37
@R.M possibly a pig. – rcollyer Mar 23 '12 at 18:37

An inkblot used to look like this, in the days when I used fountain pens and indian ink, rather than Mathematica:

blot = Image[BubbleChart[RandomReal[1, {20, 3}] , Axes -> None,
Frame -> None, ColorFunction -> Function[Black],
BubbleSizes -> {.001, .3}, Background -> LightGray,
ChartElementFunction -> "NoiseBubble", ImageSize -> 400]]


Edit: The paper is then to be folded in half to stimulate those revealing subconscious thoughts:

ImageAdjust[ImageMultiply[blot, ImageReflect[blot, Left]], {0, .2}]


-
Now you just need to factor in the porosity of the blotting paper, throw in some surface tension and capillary action and you'll be all set =) – R. M. Mar 22 '12 at 19:48
Blotting paper was not yet invented in those days – belisarius has settled Mar 22 '12 at 21:13
But since benches in those times were not flat but slanted any ink splattered on the paper immediatly ran off. – István Zachar Mar 22 '12 at 22:07
It would be nice to extract/reproduce the actual method "NoiseBubble" uses, as this could be a wrapper applied to the other solutions to give more realistic blot-edges. I guess it is just a radial positive noise over points of the circle. – István Zachar Mar 23 '12 at 9:38
To continue the food theme, the last one looks like a crab. – rcollyer Mar 24 '12 at 2:46

Here's a slow and concave version:

blot[smoothness_: 20, points_Integer: 10] :=
With[
{fun = Exp[-smoothness #.#] &, pts = RandomReal[1, {points, 2}]},
RegionPlot[
Total[fun[# - {x, y}] & /@ pts] > .5, {x, -.5, 1.5}, {y, -.5, 1.5},
Frame -> False, PlotStyle -> Black, BoundaryStyle -> Black]
]

Grid@Table[blot[], {3}, {3}]


Per Leonid's suggestion, here's a considerably faster version using "just in time" compiling:

blotc[smoothness_: 20, points_Integer: 10] :=
With[{fun = Exp[-smoothness #.#] &, pts = RandomReal[1, {points, 2}]},
With[{fc = Compile[{xl, yl}, Total[fun[# - {xl, yl}] & /@ pts] > .5]},
RegionPlot[fc[x, y], {x, -.5, 1.5}, {y, -.5, 1.5},
Frame -> False, PlotStyle -> Black, BoundaryStyle -> Black]
]
]


Thanks to the speed of the Mathematica compiler, this will speeds it up about 5 times on my computer.

Here's a fast but always convex version:

<< ComputationalGeometry
pts = With[{points = RandomReal[1, {20, 2}]}, points[[ConvexHull[points]]]]
Graphics@FilledCurve[BSplineCurve[pts, SplineClosed -> True]]


-
I'm not sure if this is desired but sometimes the slow code will generate multiple disconnected regions: i.stack.imgur.com/NAu8N.png Upvote regardless! Haha you edited your answer while I was writing my comment. – s0rce Mar 22 '12 at 15:26
@s0rce Yep, but real inkblots can be disconnected too :) It's a good very good question how to get rid of disconnected parts. One possibility is MorphologicalComponents. If we're doing image processing, we might as well try implementing this in terms of image processing functions: Binarize@ImageAdjust@ ImageCorrelate[ ImagePad[ Image@SparseArray@Thread[RandomInteger[{1, 100}, {10, 2}] -> 1], 20], GaussianMatrix[40]] (It's a start, but it still has problems) – Szabolcs Mar 22 '12 at 15:33
Next question could be leopard spots :D Maybe someone will write a full reaction-diffusion system simulation – Szabolcs Mar 22 '12 at 15:35
@Szabolcs, I've already done that but to make a 2D pattern out of an ODE system is not trivial :) The answer is really nice, I think that disconnected spots/holes make it more realistic! – István Zachar Mar 22 '12 at 15:38
Nice! +1. You can speed up your first function 8-10x by using With[{fc = Compile[{xl, yl}, Total[fun[# - {xl, yl}] & /@ pts] > .5]}, RegionPlot[fc[x, y],...]] in place of your RegionPlot call (compiled to MVM target, not C). I think this is a good example of JIT compilation in action for mma. – Leonid Shifrin Mar 22 '12 at 16:00