How to create a new "person curve"?

Question

Wolfram|Alpha has a whole collection¹ of parametric curves that create images of famous people. To see them, enter WolframAlpha["person curve"] into a Mathematica notebook, or person curve into Wolfram|Alpha. You get a mix of scientist, politicians and media personalities, such as Albert Einstein, Abraham Lincoln and PSY:

The W|A parametric people curves are constructed from a combination of trigonometric and step functions. This suggests that the images might have been created by parametrising a sequence of contours... which is backed up by some curves being based of famous photos, e.g., the W|A curve for PAM Dirac:

is clearly based on the Dirac portrait used in Wikipedia:

Here's a animation showing each closed contour of Abraham Lincoln's portrait as the plot parameter $t$ increases by $2\pi$ units:

Since the functions are so complicated, I can't believe that they were manually constructed. For example, the function to make Abe's bow tie is (for $8\pi < t < 10\pi$)
The full parametric curve for Abe has 56 such curves tied together with step functions and takes many pages to display.

So my question is:

How can I use Mathematica to take an image and produce a good looking "people curve"?

Answers can start from line art and just automatically parametrise the lines or they can start from a picture/portrait and identify a set of contours that are then parametrised. Or any other (semi)automated approach that you can think of.

¹ _{^{At the time of posting this question, it has 37 such curves.}}

Answer 1

This now has been discussed in Wolfram blog posts by Michael Trott:

Part 1: Making Formulas… for Everything — From Pi to the Pink Panther to Sir Isaac Newton
Part 2: Using Formulas... for Everything — From Complex Analysis Class to Political Cartoons to Music Album Covers
Part 3: Even More Formulas… for Everything—From Filled Algebraic Curves to the Twitter Bird, the American Flag, Chocolate Easter Bunnies, and the Superman Solid

Here is one of the example apps from blog - go read it in full - fun! Don't miss the link to download the notebook with complete code and apps at the end of the blog.

Answer 2

This shows a way to parametrise a line using the method suggested by Rahul Narain in a comment, i.e. using Fourier to approximate the data with a set of sinusoids. I use Rationalize to convert all the reals back to rationals, this isn't necessary but it makes the expression look more like those used in Wolfram Alpha.

param[x_, m_, t_] := Module[{f, n = Length[x], nf},
  f = Chop[Fourier[x]][[;; Ceiling[Length[x]/2]]];
  nf = Length[f];
  Total[Rationalize[
     2 Abs[f]/Sqrt[n] Sin[Pi/2 - Arg[f] + 2. Pi Range[0, nf - 1] t], .01][[;; Min[m, nf]]]]]

tocurve[Line[data_], m_, t_] := param[#, m, t] & /@ Transpose[data]

tocurve will take a Line, a number of modes m and a symbolic parameter t and return a parametrisation of the line data. Because of the implied periodicity of the data in Fourier it will only work properly on closed lines.

The hard part is getting a good set of lines from the image of a person. Here's a much simpler example using ListContourPlot to extract the outline of a silhouette.

First load an image and do a bit of preprocessing to ensure a nice contour:

img = Import[
   "http://catclipart.org/wp-content/uploads/2012/11/elephant-silhouette-clip-art.gif"];

img = Binarize[img~ColorConvert~"Grayscale"~ImageResize~500~Blur~3]~Blur~3;

Now extract contours and plot the parametrised curve with 500 modes:

lines = Cases[Normal@ListContourPlot[Reverse@ImageData[img], Contours -> {0.5}], _Line, -1];

ParametricPlot[Evaluate[tocurve[#, 500, t] & /@ lines], {t, 0, 1}, Frame -> True, Axes -> False]

With fewer modes the detail starts to smooth out. Here's the 30 mode curve:

The parametrisation consists of sinusoids:

curves // Short

Answer 3

This was supposed to be a comment to Simon's answer, but it's gotten too long. Still, I wanted to share a somewhat cleaned-up version of Simon's Fourier-fitting function param[] (which I have renamed to FourierCurve[]):

FourierCurve[x_, m_, t_, tol_: 0.01] := Module[{rat = Rationalize[#, tol] &, fc},
  fc = Take[Chop[Fourier[x, FourierParameters -> {-1, 1}]], Min[m, Ceiling[Length[x]/2]]];
  2 rat[Abs[fc]].Cos[Pi (2 Range[0, Length[fc] - 1] t - rat[Arg[fc]/Pi])]]

This has the virtue of returning a function that genuinely closes up; more precisely, if f[t_] = FourierCurve[pts, modes, t], then f[0] == f[1]. (The indiscriminate use of Rationalize[] in the earlier version prevented a nice closure of the resulting curve.)

As Rahul alludes to in his comment, this is more or less the "epicycle" approach of Ptolemy for determining the paths of planetary orbits.

---

Of course, Fourier fitting can also be applied to space curves as well as plane curves. Here's an example:

{f[t_], g[t_], h[t_]} = FourierCurve[#, 20, t] & /@
                        KnotData["FigureEight", "SpaceCurve"]["ValuesOnGrid"];

ParametricPlot3D[{f[t], g[t], h[t]}, {t, 0, 1}, Axes -> None, Boxed -> False,
                 Method -> {"TubePoints" -> 20}, PlotStyle -> Blue, ViewPoint -> Top] /. 
Line[pts_, rest___] :> Tube[pts, 1/8, rest]

Since most of the knots given in KnotData[] have their space curves given as InterpolatingFunction[] objects, you can use this approach if you prefer to have explicit parametric expressions for those knots.

Answer 4

Simon Woods' cool answer is based on two independent discrete Fourier transforms: one for $x$, one for $y$. This answer differs in that it uses only one, but on complex values. It's it not a very big difference that is a required step for this kind of things.

img = Import["https://i.sstatic.net/wtJoA.png"];
img = Binarize[img~ColorConvert~"Grayscale"~ImageResize~500~Blur~3];
pts = DeleteDuplicates@Cases[Normal@ListContourPlot[Reverse@ImageData[img], 
       Contours -> {0.5}], _Line, -1][[1, 1]];

z = pts[[All, 1]] + I*pts[[All, 2]];
m = 10;
n = Length@z;
cn = 1/n*Table[Sum[z[[k]]*Exp[-I*i*k*2 Pi/n], {k, 1, n}], {i, -m, m}];
{f[t_], g[t_]} = {Re@#, Im@#} &@
    Sum[cn[[i + m + 1]]*Exp[I*i*t], {i, -m, m}] // ComplexExpand;
ParametricPlot[{f[t], g[t]}, {t, 0, 2 Pi}]

• $m=10$:

• $m=20$:

To make the difference clear: Simon Woods' code with m = 3 gives

{{4341/8 + 1593/8 Sin[37/11 + (44 t)/7] + 58 Sin[1/3 + (88 t)/7], 
  3225/11 - 867/8 Sin[6/5 - (44 t)/7] + 142/3 Sin[45/23 + (88 t)/7]}}

while the above code gives:

{23.2813 - 38.9784 Cos[t] + 16.5986 Cos[2 t] - 17.5251 Cos[3 t] - 
  195.161 Sin[t] + 55.6546 Sin[2 t] - 15.5714 Sin[3 t], -101.483 - 
  101.745 Cos[t] + 44.3024 Cos[2 t] - 34.0713 Cos[3 t] + 
  37.4785 Sin[t] - 16.4218 Sin[2 t] + 5.10139 Sin[3 t]}

The base functions are the same for x and y here.

Note If someone manages to compute my cn using Fourier I'd be happy!