close all;
clear i;

echo on; clc;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% cs515 Demo: Smoothness of f ~ decay of f^
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Useful data for first two examples:

n = [-128:127];              % range of fourier coefficients
t = n / 128 * pi;            % define functions on [-pi,pi)
ttt = [t-2*pi,t,t+2*pi];     % time axis for periodic extension
N = length(n);    

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Ex 1: A piecewise constant function related to H(t)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

f = (t<0) - (t>=0);                      % evaluate f(t)
% F = (2*pi)*fftshift(fft(f))/N;         % FFT evaluation: F(n) = f^(n)
F = (2i*pi) * sin(n*pi/2) .* sinc(n/2);  % Use analytic form of f^(n)

echo off;
figure;
subplot(2,1,1);
plot(t,f,t,zeros(1,N),':');
axis([-pi,pi,-1.5,1.5]);
title('Piecewise constant f(t) := H((t/pi+1)/2)');
set(gca,'YTick',[-1,0,1]);
set(gca,'XTick',[-pi,0,pi]);
set(gca,'XTickLabel','-pi|0|pi');
icons = ['text(.5,.5,''Next>>'',''Horiz'',''center'')'];
callbacks = ['uiresume'];
btngroup('GroupID', 'TestGroup', 'ButtonID', ['Next>>'], 'Callbacks', callbacks,...
    'Position', [0.45 .01 .1 .05], 'IconFunctions', icons, 'PressType', 'flash');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
plot(ttt,[f,f,f],ttt,zeros(1,3*N),':');
axis([-3*pi,3*pi,-1.5,1.5]);
title('Periodic Extension of f(t)');
set(gca,'YTick',[-1,0,1]);
set(gca,'Xtick',[-3*pi,-pi,pi,3*pi]);
set(gca,'XtickLabel','-3pi|-pi|pi|3pi');
set(gca,'XGrid','on');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
stairs(n,abs(F));
axis tight; grid on;
title('Frequency domain');
legend('|f\^(n)|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
stairs(n,abs(F));
axis([n(1),n(end),0,.5]);
title('Frequency domain, zoomed');
legend('|f\^(n)|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
stairs(n,abs(F));
axis([n(1),n(end),0,.5]);
title('Frequency domain, zoomed');
index = find(n);
z = 4*ones(size(n));
z(index) = z(index) ./ abs(n(index));        % z = 4/abs(n)
hold on; stairs(n,z,'r');
legend('|f\^(n)|','const / |n|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;
btndown(1,'TestGroup', 'Next>>');
disp('Hit Enter for next example'); pause;

echo on; clc;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Ex 2: A piecewise linear function related to B2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

f = (t/pi+1).*(t<0) + (1-t/pi).*(t>=0);  % evaluate f(t) 
% F = (2*pi)*fftshift(fft(f))/N;         % FFT evaluation: F(n) = f^(n)
F = pi*sinc(n/2).^2;                     % Use analytic form of f^(n)

echo off;
figure;
subplot(2,1,1);
plot(t,f);
title('Piecewise Linear f(t) := B2(t/pi+1)');
axis([-pi,pi,0,1]);
set(gca,'YTick',[0,.5,1]);
set(gca,'XTick',[-pi,0,pi]);
set(gca,'XTickLabel','-pi|0|pi');
icons = ['text(.5,.5,''Next>>'',''Horiz'',''center'')'];
callbacks = ['uiresume'];
btngroup('GroupID', 'TestGroup', 'ButtonID', ['Next>>'], 'Callbacks', callbacks,...
    'Position', [0.45 .01 .1 .05], 'IconFunctions', icons, 'PressType', 'flash');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
plot(ttt,[f,f,f],ttt,zeros(1,3*N),':');
axis([-3*pi,3*pi,0,1]);
title('Periodic Extension of f(t)');
set(gca,'YTick',[0,.5,1]);
set(gca,'Xtick',[-3*pi,-pi,pi,3*pi]);
set(gca,'XtickLabel','-3pi|-pi|pi|3pi');
set(gca,'XGrid','on');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
stairs(n,abs(F));
axis tight; grid on;
title('Frequency domain');
legend('|f\^(n)|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
stairs(n,abs(F));
axis([n(1),n(end),0,.1]);
title('Frequency domain, zoomed');
legend('|f\^(n)|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
stairs(n,abs(F));
axis([n(1),n(end),0,.01]);
title('Frequency domain, zoomed');
legend('|f\^(n)|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2);
stairs(n,abs(F));
axis([n(1),n(end),0,.01]);
title('Frequency domain, zoomed');
index = find(n);
z = (4/pi)*ones(size(n));
z(index) = z(index) ./ abs(n(index)).^2;      % z = (4/pi)/abs(n)^2
hold on; stairs(n,z,'r');
legend('|f\^(n)|','const / |n|^2');
%waitforbuttonpress;
btndown(2,'TestGroup', 'Next>>');
disp('Hit Enter for next example'); pause;

echo on; clc;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Example 3: Reconstruction of a 'smooth' function %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

n = [-128:127];   % range of Fourier coefficients
t = n/128*pi;
N = length(n);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% A piecewise cubic related to B4 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

f = (2*t/pi+2).^3/6 .* (t>=-pi & t<-pi/2) ... % ... continued
    -(3*(2*t/pi+2).^3 - 12*(2*t/pi+2).^2+12*(2*t/pi+2)-4)/6 .* (t>=-pi/2 & t<0) ...
    +(3*(2*t/pi+2).^3 -24*(2*t/pi+2).^2 + 60*(2*t/pi+2)-44)/6 .*(t>=0 & t<pi/2) ...
    -(2*t/pi-2).^3/6 .* (t>=pi/2 & t<pi);
s = max(f); f = f/s;

F = (2*pi)*fftshift(fft(f))/N;         % FFT evaluation: F(n) = f^(n)
echo off

figure;
subplot(2,1,1)
plot(t,f)
title('Piecewise cubic f(t) = B_4(2t/pi+2)')
axis([-pi,pi,0,1]);
set(gca,'YTick',[0,.5,1]);
set(gca,'XTick',[-pi,0,pi]);
set(gca,'XTickLabel','-pi|0|pi');
icons = ['text(.5,.5,''Next>>'',''Horiz'',''center'')'];
callbacks = ['uiresume'];
btngroup('GroupID', 'TestGroup', 'ButtonID', ['Next>>'], 'Callbacks', callbacks,...
    'Position', [0.45 0.01 .1 .05], 'IconFunctions', icons, 'PressType', 'flash');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2)
stairs(n, abs(F))
axis tight; grid on;
title('Frequency domain');
legend('|f\^(n)|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2)
stairs(n, abs(F))
axis([n(1),n(end),0,2*10^(-3)]);
title('Frequency domain, zoomed');
legend('|f\^(n)|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2)
stairs(n, abs(F))
axis([n(1),n(end),0,2*10^(-5)]);
title('Frequency domain, zoomed');
legend('|f\^(n)|');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(2,1,2)
stairs(n, abs(F))
axis([n(1),n(end),0,2*10^(-5)]);
title('Frequency domain, zoomed');
index = find(n);
z = (2*pi)*ones(size(n));
z(index) = z(index) ./ abs(n(index)).^4;      % z = (4/pi)/abs(n)^4
hold on; stairs(n,z,'r');
legend('|f\^(n)|','const / |n|^4');
%waitforbuttonpress;
btndown(3,'TestGroup', 'Next>>');
disp('Hit Enter to continue'); pause;

figure;

F_3 = F.*(abs(n)<2);
ff_3 = ifft(ifftshift(F_3))*N/2/pi;
plot(t,f,t,real(ff_3))
title('Original f and Reconstructed from 3 Fourier coefficients')
axis tight
legend('f', 'reconstructed f');
icons = ['text(.5,.5,''Next>>'',''Horiz'',''center'')'];
callbacks = ['uiresume'];
btngroup('GroupID', 'TestGroup', 'ButtonID', ['Next>>'], 'Callbacks', callbacks,...
    'Position', [0.45 0.01 .1 .05], 'IconFunctions', icons, 'PressType', 'flash');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(1,1,1)
F_5 = F.*(abs(n)<3);
ff_5 = ifft(ifftshift(F_5))*N/2/pi;
plot(t,f,t,real(ff_5))
title('Original f and Reconstructed from 5 Fourier coefficients')
axis tight
legend('f', 'reconstructed f');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(1,1,1)
F_21 = F.*(abs(n)<10);
ff_21 = ifft(ifftshift(F_21))*N/2/pi;
plot(t,f,t,real(ff_21))
title('Original f and Reconstructed from 21 Fourier coefficients')
axis tight
legend('f', 'reconstructed f');

%waitforbuttonpress;
btndown(4,'TestGroup', 'Next>>');
disp('Hit Enter for next example'); pause;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Example 4: Comparison of two 'smooth' functions
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

N = 256;
a = N*pi;   % freq interval := [-a,a]

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% A piecewise cubic related to B4 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
t = 4 * linspace(0,N,N*(a/pi));  
f = (t.^3) / 6 .*  (t>=0 & t<1) ...     %  ... continued
    - (3*t.^3 - 12*t.^2 + 12*t -  4) / 6  .*  (t>=1 & t<2) ...
    + (3*t.^3 - 24*t.^2 + 60*t - 44) / 6  .*  (t>=2 & t<3) ...
    - (  t.^3 - 12*t.^2 + 48*t - 64) / 6  .*  (t>=3 & t<4); 
s = max(f); f = f / s;
F = fftshift(fft(f))/N;

% w = linspace(-a,a,N*(a/pi));   
% G = 1/(4*s)*(sin(w/8)./(w/8)).^4;      % analytic form

echo on;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% A C_inf function with compact support
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
t = linspace(0,N,N*(a/pi));
g = zeros(size(t));
index = find(t>0&t<1);
g(index) = exp(1./(t(index)-1)-1./t(index));
s = max(g); g = g / s;
G = fftshift(fft(g))/N; 

echo off;
figure;
subplot(4,1,1);
t = linspace(0,1,N);
plot(t,f(1:N));
axis tight; ax = axis;
title('Piecewise cubic f(t) := B_4(4t)');
set(gca,'XTick',[0,1]);

subplot(4,1,2);
t = linspace(0,1,N);
plot(t,g(1:N),'r');
axis(ax);
title('C-inf g(t) := 55 * exp(1/(t-1)-1/t) for t in [0,1), 0 otherwise');
set(gca,'XTick',[0,1]);
icons = ['text(.5,.5,''Next>>'',''Horiz'',''center'')'];
callbacks = ['uiresume'];
btngroup('GroupID', 'TestGroup', 'ButtonID', ['Next>>'], 'Callbacks', callbacks,...
    'Position', [0.45 0.01 .1 .05], 'IconFunctions', icons, 'PressType', 'flash');
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(4,1,3);
w = linspace(-a,a,length(F));
plot(w,abs(F));
axis tight;
set(gca,'Xtick',[-a,0,a]);
set(gca,'XtickLabel',['-',num2str(a/pi),'pi|0|',num2str(a/pi),'pi']);
legend('|f\^(w)|');
ax = axis;

subplot(4,1,4);
w = linspace(-a,a,length(G));
plot(w,abs(G),'r');
legend('|g\^(w)|');
axis(ax); 
set(gca,'Xtick',[-a,0,a]);
set(gca,'XtickLabel',['-',num2str(a/pi),'pi|0|',num2str(a/pi),'pi']);
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(4,1,3); axis([-a,a,0,1e-3]);
subplot(4,1,4); axis([-a,a,0,1e-3]);
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(4,1,3); axis([-a,a,0,1e-5]);
subplot(4,1,4); axis([-a,a,0,1e-5]);
%waitforbuttonpress;
disp('Hit Enter to continue'); pause;

subplot(4,1,3); axis([-a,a,0,2e-7]);
subplot(4,1,4); axis([-a,a,0,2e-7]);
btndown(5,'TestGroup', 'Next>>');

