Actually, the algorithm is straightforward, it is perfectly well explained in the video you linked, and I'm not entirely sure why you haven't simply tried it yourself. So what I will provide is a version for grayscale images and I'll leave it to you to extend this for colour images.
The algorithm itself is highly modular and you can work on each subproblem almost separately. In the video, it was explained that you put nails around the image and then use one long dark string that goes from nail to nail to resemble the image.
We start by using the pixel matrix of an image which is a matrix. Then we extract the nail positions which are all pixels around the image. Extracting the positions for all four sides can be done with the dimension of the matrix
getNailPositions[{nx_, ny_}] := Join[
Table[{x, 1}, {x, nx - 1}],
Table[{nx, y}, {y, ny - 1}],
Table[{x, ny}, {x, nx, 2, -1}],
Table[{1, y}, {y, ny, 2, -1}]
]
As explained in the video, we start at one nail and we need to make a line to all other nails and check, which of these lines is the darkest one. One prerequisite is an algorithm that takes two points and gives us all integer positions on the line between the points. I have written here about the Bresenham algorithm which is stunningly simple. For the answer here, I'm going to re-implement my solution using Compile
so that we can check all lines in parallel
bresenhamC = Compile[{{p0, _Integer, 1}, {p1, _Integer, 1}},
Module[{dx = 0, dy = 0, sx = 0, sy = 0, err = 0, newp = p0,
result = Internal`Bag[p0]},
{dx, dy} = Abs[p1 - p0];
{sx, sy} = Sign[p1 - p0];
err = dx - dy;
While[newp =!= p1,
newp = With[{e2 = 2 err}, {
If[e2 > -dy, err -= dy; Compile`GetElement[newp, 1] + sx,
Compile`GetElement[newp, 1]],
If[e2 < dx, err += dx; Compile`GetElement[newp, 2] + sy,
Compile`GetElement[newp, 2]]
}];
Internal`StuffBag[result, Compile`GetElement[newp, 1]];
Internal`StuffBag[result, Compile`GetElement[newp, 2]];
];
Partition[Internal`BagPart[result, All], 2]
],
CompilationTarget -> "C",
Parallelization -> True,
RuntimeAttributes -> {Listable}
];
But that's not all, because it is not sufficient to get all the positions. What we need is the mean darkness along the lines, as this is our selector that indicates which string we draw next. So having two nails we want to check, we need to get all points along the line in between and then we have to calculate the mean pixel values in the image on this line. There are other measures that could be used when you want to implement an adaption for colour images.
With[{bresenhamC = bresenhamC},
pixelMeasure =
Compile[{{img, _Real, 2}, {p0, _Integer, 1}, {p1, _Integer, 1}},
Module[{points = bresenhamC[p0, p1], v = p1 - p0},
Mean[Table[Compile`GetElement[img, p[[1]], p[[2]]], {p, points}]]
],
CompilationOptions -> "InlineCompiledFunctions" -> True,
Parallelization -> True,
RuntimeAttributes -> {Listable},
CompilationTarget -> "C"
]
]
If we have found the best next nail (p1
) to our starting nail p0
, we need to brighten the image on this line. Otherwise, the algorithm would end up using the same connections. I will make this in place which means that we have one starting image that gets continually brighter. To make this work, I'm using HoldFirst
to get the symbol that contains the image matrix instead of getting the matrix itself
SetAttributes[weakenLine, {HoldFirst}];
weakenLine[img_, {p0_, p1_}] := With[{pts = bresenhamC[p0, p1]},
Do[
Part[img, Sequence @@ p] += 0.03
, {p, pts}
]
]
One step of the algorithm works as follows
- Take the image matrix, a starting nail point, and a list of all possible nail points
- Select the next best nail point that has the lowest pixel measure alone the line
- Make this line weaker in the image matrix
- Return the new next nail position
That gives us
SetAttributes[step, {HoldFirst}];
step[img_, p0_, nails_] := Module[{points, sel},
sel = First@Ordering[pixelMeasure[img, p0, nails]];
weakenLine[img, {p0, nails[[sel]]}];
nails[[sel]]
]
We are ready to put everything together to calculate the series of nail positions we need to visit. I left out one little detail: It is probably a good idea to pad the image matrix with one pixel-layer of white. This will prevent the algorithm from going along the image border
calculateLines[img_, n_] := Module[
{
data = ArrayPad[ImageData[img, "Real"], 1, 1.0],
positions,
lines
},
positions = With[{pos = getNailPositions[Dimensions[data]]},
Table[pos[[i]], {i, 1, Length[pos]}]
];
Print[Dynamic[Image[data]]];
NestList[step[data, #, positions] &, {1, 1}, n]
]
Now it is time for a first test. The dynamic output I put inside the above function is to show you how the original image slowly degrades. Three parameters influence the degradation:
- The number of lines
n
you create
- How much white gets added in the
weakenLine
function
- How many nailpoints you have choosen inside the
calculateLines
function (currently all pixel positions around the image)
Here is how the code can be used:
img = ColorConvert[Import["https://i.sstatic.net/YRCQD.png"], "Grayscale"];
lines = calculateLines[img, 2000];
This is the degraded image after 2000 lines. As you can see nothing much happend and with the current settings, we aim for a very detailed image that needs way more lines.
To show this in a Graphics
(which we won't do later), we need to massage the coordinates of the lines a bit as they are in pixel coordinates which use a different coordinate-system:
Graphics[{Opacity[0.1, Black], AbsoluteThickness[1],
Line[{1, -1}*Reverse[#] & /@ lines]}]
Since we have already have implemented the Bresenham, a better idea is to use it and draw the lines onto a completely white image. That however, is again very simple. We need a more general function of our weakenLine
algorithm that applies a function to all pixels on a line between points p0
and p1
SetAttributes[applyLine, {HoldFirst}];
applyLine[img_, {p0_, p1_}, f_Function] :=
With[{pts = bresenhamC[p0, p1]},
Do[
Part[img, Sequence @@ p] = f@Part[img, Sequence @@ p]
,{p, pts}
]
]
Now, we calculate 10000 lines and after that, we draw the lines onto the white canvas
lines = calculateLines[img, 10000];
Time for tea.. and when this is done
Dynamic[Image[out]]
out = ConstantArray[1.0, Reverse[ImageDimensions[img]] + 2];
Do[
applyLine[out, line, Function[val, val - .028]],
{line, Partition[lines, 2, 1]}
]