16
$\begingroup$

How can I find the simplest curve that gives me a scroll like this?

enter image description here

This is where I am at the moment and doing it by hand is becoming more and more cumbersome.

First, take a drawing of the scroll since the outline is already processed:

enter image description here

Then, define a spiral in a parametric plot with several parameters that move things around:

Manipulate[
 Show[ParametricPlot[
   1/d (1 + a Sin[ϕ/2 + f]/4 + a Sin[ϕ/4 + f]/4) Norm[
     2 + ϕ^1.1] {Sin[-ϕ], 0.075 + Cos[ϕ]}, {ϕ, 0, 
    5.25 π}, PlotRange -> {-1, 1}],
  Epilog -> Disk[{0, 0}, .075],
  Prolog -> {Texture[img], 
    Polygon[{Scaled[{0, 0}], Scaled[{1, 0}], Scaled[{1, 1}], 
      Scaled[{0, 1}]}, 
     VertexTextureCoordinates -> {{0 + b, c}, {1 + b, c}, {1 + b, 
        1 + c}, {0 + b, 1 + c}}]}], {a, .5, 
  2}, {b, .091}, {c, .148}, {d, 30, 40}, {f, 0, π}]

(b and c are there to center the image only, could not find a nicer way)

What I tried to do was to increase the radius of the spiral as a function of $\phi$ but I have the impression a computer would be smarter than me in finding this. This is the best I can do atm.

enter image description here

I know there is a new FindFormula but cannot figure out how to make it work for parametric functions.

I would very much like something that does not involve BSplineCurve. A piecewise function on $\phi$ is welcome. Violin shapes were made with compass and rulers so I imagine that there must be a simple function for the curvature and how it changes as the scroll turns. With simple function I mean something like this: enter image description here

In this case the radius of each segment on the scroll is given by the projection from the lower triangle. This can be taken to the continuous case as a simple function.

Using the method provided by nikie I can obtain different scrolls varying the functions I use for the fit. This is the one I most like fitFn = a + b theta + c theta^2 + d /theta + e/theta^2 + h Sin[ f (-theta)^.5 + g];

Maybe you can find out to which one it corresponds

enter image description here

$\endgroup$
5
  • 3
    $\begingroup$ This may be helpful. $\endgroup$ Apr 22, 2016 at 1:29
  • $\begingroup$ Okay I'm confused: splines are still piecewise functions after all, so you seem to have conflicting desires. Also, "with compass and rulers" can be taken to mean that the arcs could very well be made up of circle arcs and line segments suitably pasted together. $\endgroup$ Apr 22, 2016 at 1:46
  • $\begingroup$ Maybe I should rephrase it as "don't use BSplineCurve" which would be easy enough to fit with enough points. Let me try to edit the question $\endgroup$ Apr 22, 2016 at 1:50
  • $\begingroup$ With respect to "simple": the lituus seems to be a good starting point as any for a volute. $\endgroup$ Apr 22, 2016 at 2:08
  • $\begingroup$ @J.M. I tried a few Archimedes' spirals but it doesn't really work, you need to add some kind of "asymmetry" (for lack of a better word), that's what my 1 + a Sin[\[Phi]/2 + f]/4 + a Sin[\[Phi]/4 + f]/4 is doing, modulating the radius of the spiral to try to fit it the real scroll. I thought about subtracting a simple spiral to the edge and fitting the difference to a function with FindFormula but don't really know yet how to yet. $\endgroup$ Apr 22, 2016 at 2:38

2 Answers 2

20
$\begingroup$

Step 1: find the curve in the black and white image. For this, I will use a shortest path search, and to make sure it finds the path I'm looking for (instead of one of the short cuts through the labels), I will give it a few "path markers" along the way that the path should visit in order:

First, load the image, get the black pixels:

img = Import["https://i.stack.imgur.com/XBEkk.png"];    
blackPts = PixelValuePositions[Binarize[img], 0];    
nearestFn = Nearest[blackPts -> Range[Length[blackPts]]];    

Next, enter rough path markers using a LocatorPane

pathMarkers = {{373, 5}, {290, 595}, {475, 588}, {495, 371}, {382, 476}};    
LocatorPane[Dynamic[pathMarkers], 
 Dynamic@Show[img, 
   Graphics[{Red, 
     Line[blackPts[[First /@ nearestFn /@ pathMarkers]]]}]]]

enter image description here

Then, turn the black pixels into a graph, with point distances as edge weight:

vertices = Range[Length[blackPts]];
edges = Select[
   Flatten[Thread /@ Thread[vertices <-> nearestFn[blackPts, 20]]],
   #[[1]] < #[[2]] &];
edgeWeights = 
  edges /. {i_ <-> j_ :> Norm[N@blackPts[[i]] - blackPts[[j]]] - .1};
g = Graph[vertices, edges, EdgeWeight -> edgeWeights];

(Don't try to evaluate g in a notebook. It's a very large graph, displaying it crashed the kernel on my machine. But we don't have to display it, we just want to find paths in it, which is very fast.)

Then find the shortest path that touches out path markers:

pathPoints = 
  Flatten[FindShortestPath[g, ##] & @@@ 
    Partition[First /@ nearestFn /@ pathMarkers, 2, 1]];
path = blackPts[[pathPoints]];
HighlightImage[img, {Thick, Line[path]}]

enter image description here

That looks promising.

Step 2: transform to polar coordinates.

I find handling polar coordinates easier working in the complex plane:

center = {402, 465};
pathComplex = N[(# - center & /@ path).{1, I}];

The radii are obvious:

radii = Abs[pathComplex];

The branch cut takes a bit more effort:

first, multiply each point with the conjugate of the next point:

pointPairs = 
  MapThread[#2*Conjugate[#1] &, {pathComplex[[;; -2]], 
    pathComplex[[2 ;;]]}];

Now, the phase of pointPairs is the angle between each point pairs. And since the points are close to each other, this will not contain any discontinuities due to the branch cut:

anglesBetweenPointPairs = Im@Log[pointPairs];

I just have to sum these angles up to get the original phase again:

angles = Prepend[Accumulate[anglesBetweenPointPairs], 0] + 
   Im[Log[pathComplex[[1]]]];

We can check that it's the same path:

ListPolarPlot[{angles, radii}\[Transpose], Joined -> True]

enter image description here

And that there are no path cuts:

ListLinePlot[{angles, radii}\[Transpose]]

enter image description here

You can even move the center around and watch as the polar curve changes dynamically to see where the best center point might be:

center = {402, 465};
Row[{
  LocatorPane[Dynamic[center], Image[img, ImageSize -> All]],
  Dynamic@
   Module[{pathComplex, radii, pointPairs, anglesBetweenPointPairs, 
     angles},
    pathComplex = N[(# - center & /@ path).{1, I}];
    radii = Abs[pathComplex];
    pointPairs = 
     MapThread[#2*Conjugate[#1] &, {pathComplex[[;; -2]], 
       pathComplex[[2 ;;]]}];
    anglesBetweenPointPairs = Im@Log[pointPairs];
    angles = 
     Prepend[Accumulate[anglesBetweenPointPairs], 0] + 
      Im[Log[pathComplex[[1]]]];
    ListLinePlot[{angles, radii}\[Transpose], ImageSize -> 500]
    ]}]

enter image description here

(I have no idea what kind of function might have been used for the radius/angle relationship, so I'm not going to attempt the actual curve fitting. IMHO, if you don't know anything about a functional relationship, a spline usually is the best guess.)

FWIW, I get a somewhat decent approximation using these terms:

Clear[g, eq]
fitFn = a + b theta + c theta^2 + d/theta + e/theta^2;
eq[theta_] = 
  fitFn /. FindFit[{angles, radii}\[Transpose], 
    fitFn, {a, b, c, d, e, f, g, h}, theta];

Show[
 ListLinePlot[{angles, radii}\[Transpose], ImageSize -> 500],
 Plot[eq[theta], {theta, Min[angles], Max[angles]}, 
  PlotStyle -> Directive[Dashed, Red], PlotRange -> All]]

enter image description here

HighlightImage[img,
 Line[Array[center + ReIm[Exp[#*I + Log[eq[#]]]] &, 500, 
   MinMax[angles]]]]

enter image description here

$\endgroup$
6
  • $\begingroup$ Beautiful. At the moment I;m using this to fit: a + b theta + c Cos[theta + e] + f theta^2 + g Sqrt[Abs[theta]] I don't know either what function was, but I guess it should be easy to obtain with compass and ruler via different projections. Now, which kind of projections are those... that's another question. Thanks a lot! $\endgroup$ Apr 22, 2016 at 19:58
  • $\begingroup$ I've added a linear fit. No idea if this can be constructed using ruler and compass. $\endgroup$ Apr 22, 2016 at 20:17
  • 1
    $\begingroup$ This might be useful for taking care of the branch cuts in the angles. $\endgroup$ Apr 22, 2016 at 22:30
  • $\begingroup$ One of the most beautiful things of the scroll is that has a certain "movement". If you plot the difference between your last fit and the curve you obtain a very marked oscillation with a defined frequency during at least 2 cycles. I think that's what make the scrolls alive somehow. I added some figures obtained with your code. Thanks again. $\endgroup$ Apr 22, 2016 at 23:37
  • $\begingroup$ I thought the "oscillation" was because I missed the actual center point... $\endgroup$ Apr 23, 2016 at 7:06
13
$\begingroup$

I didn't get a good answer, but this is how I would do it and probably how I would evaluate the quality of an answer. Much of what I do below involves manual work, but it could be automated:

First, we need some points along the curve. I selected them manuallyenter image description here

If I hadn't manually selected them, I would have used some image processing including a thinning transform and then I would have used FindCurvePath.

We also want to select a point that seems like the "center" of the spiral.

center = {{165.92578125`, 147.453125`}};
pts = {{183.65234375`, 292.625`}, {201.4765625`, 
   287.3671875`}, {215.21484375`, 282.234375`}, {226.21484375`, 
   277.7578125`},.....};

And we want to normalize the pts so the "center" is the origin:

normalizedPts = (# - First@center) & /@ pts;

We want to work with this in polar coordinates:

ToPolarCoordinates@normalizedPts

But there's a problem. There's a branch cut. I've manually edited the results so that there isn't a branch cut. This involves subtracting or adding 2Pi to the results so that they form a continuous curve:

arc1 = {{146.25014293696273`, 
   1.4492904217728535`}, .... {82.16481360563434`, -3.113682087252578`}}

arc2 = {{79.36095240482918`, 2.9819303345670916`},... {39.305506825650895`, -2.883525354070561`}}

arc3 = {{36.32153376025393`, 3.0092472204025533`}, {19.075687437056864`, -2.8414075099042826`}}

branchCutsRemoved = Reverse/@ 
 Join[arc1, (# - 2 Pi) & /@ arc2, (# - 4 Pi) & /@ arc3];

I reversed all the points so that they're in {arg, abs} form. This means we are looking for the absolute value of the point as a function of arg. I use theta for arg or the angle.

We now have a fairly reasonable function of theta to try and approximate:

ListPlot[branchCutsRemoved, PlotRange -> All]

enter image description here

You can try to fit to this:

eq = Normal@
  NonlinearModelFit[branchCutsRemoved, 
   a + b theta + c  Cos[d theta + e] + f  Cos[d theta/2 + e], {a, b, 
    c, d, e, f, g, h}, theta]

The fit seems... okayish:

Show[ListPlot[sorted, PlotStyle -> Red], Plot[eq, {theta, -16, 5}]]

enter image description here

$\endgroup$
2
  • $\begingroup$ thanks a lot, I didn't know about the coordinate tool, quite useful. Is there a way to automatically fix the branches? Will leave the question open for a bit to see if another answer appears but this is kinda exactly what I was looking for, now I will check the coefficients for different scrolls and see how the different "schools" behave. $\endgroup$ Apr 22, 2016 at 18:56
  • $\begingroup$ Yes. You could write a function that goes through a list and when it sees a value drop by 2Pi adds 2Pi to all subsequent values. $\endgroup$
    – Searke
    Apr 22, 2016 at 19:17

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.