# Fitting a curve using a spline, Fourier transform, etc.

1. Mar 5, 2013

### Jamin2112

1. The problem statement, all variables and given/known data

Just wondering if my output seems wrong. The interpolating polynomial looks like it's way off, though I've looked over my code many times and it seems right (?).

clc
clear all
format long

x1=[1:1/10:4];
y1=zeros(1,length(x1));
y1(find(x1<=3))=floor(x1((find(x1<=3))).^2);
y1(find(x1>3))=9*cos(3*x1(find(x1>3))-9).^3;
xx=[1:1/1000:4];
yy=zeros(1,length(xx));
yy(find(xx<=3))=floor((xx(find(xx<=3))).^2);
yy(find(xx>3))=9*cos(3*xx(find(xx>3))-9).^3;
p=polyfit(x1,y1,length(x1)+1);
y30=polyval(p,xx);
ffty1=fft(y1);
fsum=0;
for k=0:15
fsum=fsum+(ffty1(k+1))*(k+1)*exp(2*pi*i*k*(xx-1)/3);
end
for k=16:30
fsum=fsum+(ffty1(k+1))*(k+1)*exp(2*pi*i*(k-31)*(xx-1)/3);
end
f31=fsum/31;
ynakstruct=spline(xx,yy);
ynakcoeffs=ynakstruct.coefs;
ynak=zeros(1,length(xx));
ynak(1)=yy(1);
for k=2:length(ynak)
ynak(k)=ynakcoeffs(k-1,4)+ynakcoeffs(k-1,3)*(xx(k)-xx(k-1))+ynakcoeffs(k-1,2)*(xx(k)-xx(k-1))^2+ynakcoeffs(k-1,1)*(xx(k)-xx(k-1))^3;
end
yclstruct=spline(xx,[0 yy 0]);
yclcoeffs=yclstruct.coefs;
ycl=zeros(1,length(xx));
ycl(1)=yy(1);
for i=2:length(ycl)
ycl(k)=yclcoeffs(k-1,4)+yclcoeffs(k-1,3)*(xx(k)-xx(k-1))+yclcoeffs(k-1,2)*(xx(k)-xx(k-1))^2+yclcoeffs(k-1,1)*(xx(k)-xx(k-1))^3;
end
figure;
subplot(2,2,1)
hold on
plot(xx,yy,xx,y30,'-')
legend('yy', 'y30')
e1=(1/3001)*norm(y30-yy,1);
title(sprintf('Error: %d', e1))
subplot(2,2,2)
hold on
plot(xx,yy,xx,f31,'-')
legend('yy', 'f31')
e2=(1/3001)*norm(f31-yy,1);
title(sprintf('Error: %d', e2))
subplot(2,2,3)
hold on
plot(xx,yy,xx,ynak,'-')
legend('yy', 'ynak')
e3=(1/3001)*norm(ynak-yy,1);
title(sprintf('Error: %d', e3))
subplot(2,2,4)
hold on
plot(xx,yy,xx,ycl,'-')
legend('yy', 'ycl')
e4=(1/3001)*norm(ycl-yy,1);
title(sprintf('Error: %d', e4))

x2=[1:1/100:4];
y2=zeros(1,length(x2));
y2(find(x2<=3))=floor(x2((find(x2<=3))).^2);
y2(find(x2>3))=9*cos(3*x2(find(x2>3))-9).^3-9;
yy=zeros(1,length(xx));
yy(find(xx<=3))=floor((xx(find(xx<=3))).^2);
yy(find(xx>3))=9*cos(3*xx(find(xx>3))-9).^3;
p=polyfit(x2,y2,length(x2)+1);
y301=polyval(p,xx);
ffty2=fft(y2);
fsum=0;
for k=0:150
fsum=fsum+(ffty2(k+1))*(k+1)*exp(2*pi*i*k*(xx-1)/3);
end
for k=151:300
fsum=fsum+(ffty2(k+1))*(k+1)*exp(2*pi*i*(k-301)*(xx-1)/3);
end
f301=fsum/301;
ynakstruct=spline(xx,yy);
ynakcoeffs=ynakstruct.coefs;
ynak=zeros(1,length(xx));
ynak(1)=yy(1);
for k=2:length(ynak)
ynak(k)=ynakcoeffs(k-1,4)+ynakcoeffs(k-1,3)*(xx(k)-xx(k-1))+ynakcoeffs(k-1,2)*(xx(k)-xx(k-1))^2+ynakcoeffs(k-1,1)*(xx(k)-xx(k-1))^3;
end
yclstruct=spline(xx,[0 yy 0]);
yclcoeffs=yclstruct.coefs;
ycl=zeros(1,length(xx));
ycl(1)=yy(1);
for k=2:length(ycl)
ycl(k)=yclcoeffs(k-1,4)+yclcoeffs(k-1,3)*(xx(k)-xx(k-1))+yclcoeffs(k-1,2)*(xx(k)-xx(k-1))^2+yclcoeffs(k-1,1)*(xx(k)-xx(k-1))^3;
end
figure(2)
subplot(2,2,1)
hold on
plot(xx,yy,xx,y301,'-')
legend('yy', 'y301')
e1=(1/3001)*norm(y301-yy,1);
title(sprintf('Error: %d', e1))
subplot(2,2,2)
hold on
plot(xx,yy,xx,f301,'-')
legend('yy', 'f301')
e2=(1/3001)*norm(f301-yy,1);
title(sprintf('Error: %d', e2))
subplot(2,2,3)
hold on
plot(xx,yy,xx,ynak,'-')
legend('yy', 'ynak')
e3=(1/3001)*norm(ynak-yy,1);
title(sprintf('Error: %d', e3))
subplot(2,2,4)
hold on
plot(xx,yy,xx,ycl,'-')
legend('yy', 'ycl')
e4=(1/3001)*norm(ycl-yy,1);
title(sprintf('Error: %d', e4))

2. Relevant equations

fft(x) returns the Fourier coefficients for the function x given by sample points
polyfit(x,y,n) returns n polynomial coefficients for a polynomial fitted to the points given by the vectors x and y
polyval(p,x) returns the values for a polynomial with coefficients given in the vector p at the points in the vector x

3. The attempt at a solution

So my spline things seem to be working out. It's just the other two that are looking way off, even though I've checked my code a bunch. Just wondering if maybe you can think of some numerical reason why this might be the case.

2. Mar 5, 2013

### Staff: Mentor

Can you get some debug data from the poly fit?
Did you test it with an easier function? Or with more sampling points? The "fit" has as many equations as the polynomial has degrees of freedom, so it is not a real "fit" here.

What happens to ffty1, and where do you plot it?

3. Mar 5, 2013

### Jamin2112

Isn't it asking for a polynomial of order 30? So I should be getting those coefficients from the smaller set of 31 data points y1, for a polynomial of order 30 with coefficients p={c0, c1, ..., c30}, and I can evaluate it at the larger set of data points using polyval like I did.

% Get the Fourier coefficients (complex #'s) for the 31 pts of data
ffty1=fft(y1);
% Use those to construct the function f31 at the pts xx={1, 1.01, 1.02, ... 4.0}
fsum=0;
for k=0:15
fsum=fsum+(ffty1(k+1))*(k+1)*exp(2*pi*i*k*(xx-1)/3);
end
for k=16:30
fsum=fsum+(ffty1(k+1))*(k+1)*exp(2*pi*i*(k-31)*(xx-1)/3);
end
f31=fsum/31;

4. Mar 6, 2013

### Staff: Mentor

Sure, but you use 31 data points to determine 31 coefficients. That is just a linear equation system, and it should have a unique solution, where the polynomial hits all 31 points.

fsum is a function? Oh, okay.
I don't understand why you need two different cases here.

5. Mar 6, 2013

### Jamin2112

Could you explain the second bullet point to me?

6. Mar 7, 2013

### Staff: Mentor

Why do you use different formulas for k<=15 and k>15?

7. Mar 7, 2013

### Jamin2112

Isn't that what it shows on the assignment?

8. Mar 8, 2013

### Staff: Mentor

Right, but I don't understand the reason there either. It might be related to the specific way the fourier transform is calculated.