function y = curve()
clf;

figure(gcf);
set(gca, 'XLim', [0, 10], 'YLim', [0, 10], 'DataAspectRatio', [1,1,1]);

c = getPoints()
hold on;
drawCurveAndTarget(c, 'b');
n = size(c, 2);
select = floor(n/2);
p = getTarget()
drawCurveAndTarget(c, 'b', p, 'g')
hold off;
% now solve for new points, q, with
%
% L(q) = L(c)    [keep shape the same]
% q(select) = p  [make the selected point move to location p]
% q(0) = c(0)    [and make the curve start at the same place as before]

q = findQ(c, select, p)
clf;
hold on;
drawCurveAndTarget(c, 'b', p, 'g')
set(gca, 'XLim', [0, 10], 'YLim', [0, 10], 'DataAspectRatio', [1,1,1]);
hold on
drawCurveAndTarget(q', 'b', q(select, :)', 'o')
hold off
set(gca, 'XLim', [0, 10], 'YLim', [0, 10], 'DataAspectRatio', [1,1,1]);


function q = findQ(c, select, p)
n = size(c, 2);
row = zeros(1, n);
row([1, 2]) = [-1, 1];
L = zeros(n, n)
for i = 1:n
    L(i,:) = row;
    row = circshift(row, [0 1]);
end
row = zeros(1, n)
row(1) = 1
M = [L; 50* row]; % target is now L(c); c(:, 0)'
row = zeros(1, n)
row(select) = 1
M = [M; 50* row]; % target is now L(c); c(:, 0)'; pt'

target = [L*c'; 50* c(:, 1)'; 50* p]
q = M \ target;
M*q - target
4


function c = getPoints()
[x, y] = getpts(gcf);
c =  [x, y]'

function drawCurveAndTarget(c, col, pt, col2)
xs = c(1,:);
xs = xs([1:end, 1])
ys = c(2,:);
ys = ys([1:end, 1])
plot(xs, ys, col);
if (nargin == 4)
    plot(pt(1), pt(2), 'or')
end

function y = getTarget()
[x, y] = ginput(1)
y = [x, y]