The polygon midpoint puzzle

Suppose I assigned the following task:

Write a program to create a polygon from a random set of point (any sort of craziness allowed with crossing lines, concave shapes, etc.!)

Show what happens to the polygon as you iteratively replace update the polygon by replacing each point with its midpoint.

Step 1: Write code to plot a polygon.

In the first step of our solution, we'll write code to plot a polygon.

In [6]:
using Plots # this is a plotting package for Julia
n = 11 # generate 11 random points
x = randn(n)
y = randn(n)

scatter(x,y) # this plots each point
for i=1:n-1
    plot!([x[i],x[i+1]], [y[i],y[i+1]], color=:lightblue)
end
plot!()
Out[6]:

That's a mess! Let's add some commands to clean up the plot and actually plot that last line.

In [8]:
plot!([x[n],x[1]],[y[n],y[1]],color=:lightblue) # plot the last line
plot!(legend=false)
Out[8]:

Now we are going to wrap this into a function. This will make it easier to reuse.

In [10]:
"""
`plot_polygon` A function to plot a set of points as a polygon.
"""
function plot_polygon(x,y)
    if maxabs(x) > 0.5
        x = 0.5*(x/maxabs(x)); # scale so x max is 0.5;
    end
    if maxabs(y) > 0.5
        y = 0.5*(y/maxabs(y)); # scale so x max is 0.5;
    end
    
    scatter(x,y); # draw the points

    # draw the lines
    npts = length(x);
    for i=1:npts-1
        plot!([x[i], x[i+1]], [y[i], y[i+1]], color=:lightblue);
    end
    # draw the last line
    plot!([x[length(x)], x[1]], [y[length(y)], y[1]], color=:lightblue);

    # set the limits on the plot, and adjust properties
    plot!(xlim=(-0.5,0.5), ylim=(-0.5,0.5), legend=false)
end

n = 11;           # use a 11 point polygon.
x = randn(n,1)
y = randn(n,1)
plot_polygon(x,y)
Out[10]:

Step 2 write code to take the midpoints

In the next step, we need to see what happens when we take the midpoints of each line of the polygon. This will generate a new set of points from the original set of points.

In [51]:
xnew = zeros(n,1)
ynew = zeros(n,1)
for i=1:n-1
    xnew[i] = (x[i] + x[i+1])/2.
    ynew[i] = (y[i] + y[i+1])/2.
end
# handle the stinker at the end
xnew[n] = (x[n]+x[1])/2
ynew[n] = (y[n]+y[1])/2

# Show the midpoints without the final lines, this is a
# really simple version of plot_polygon.
scatter(x,y)
plot!(x,y)
scatter!(xnew,ynew)
Out[51]:
In [ ]:
function poly_midpoints(x,y)
    n = length(x)
    xnew = zeros(n,1)
    ynew = zeros(n,1)
    for i=1:n-1
        xnew[i] = (x[i] + x[i+1])/2.
        ynew[i] = (y[i] + y[i+1])/2.
    end
    # handle the stinker at the end
    xnew[n] = (x[n]+x[1])/2
    ynew[n] = (y[n]+y[1])/2
    return xnew,ynew
end
In [7]:
## With a for loop and normalized to keep it from shrinking
n = 11
x = randn(n,1)
y = randn(n,1)
x = x - sum(x)/n # center on the origin
y = y - sum(y)/n
x = x/norm(x) # rescale
y = y/norm(y)
niter = 30
anim = @animate for iter=1:niter
    
    xnew,ynew = poly_midpoints(x,y)
    xnew = xnew/norm(xnew)
    ynew = ynew/norm(ynew)

    plot_polygon(xnew,ynew)

    x = xnew;
    y = ynew;

end

gif(anim, "tmp.gif", fps = 10)

## It turns out you can prove this ALWAYS
# looks like an ellipse at 45 degs!
# http://www.cs.cornell.edu/cv/ResearchPDF/EllipsePoly.pdf
INFO: Saved animation to /home/juser/CS314-2016/tmp.gif
Out[7]:
In [ ]: