% 
% Aliasing Reduction in Clipped Signals
% Fabian Esqueda, Stefan Bilao and Vesa Välimäki
% Aalto University, Dept. of Signal Processing and Acoustics
% 
% MATLAB code used to implement Fig. 10
% 
% August 2015
% 
%%
close all; clear all; clc;

% INITIAL SETUP
step = 1/1000;
d = 0:step:1-step;
dl = length(d);
z = zeros(1,4*dl);

% SETUP FIGURE
figure(1)
set(1, 'Color', 'w');
% set(1,'units','normalized','outerposition',[0 0 0.38 0.6])
fontsize = 16;
fontname = 'Times';
markersize = 16;
linewidth = 2;

% polySINC
polySINC_L_L = d.^3/6; 
polySINC_L = - d.^3/2 + d.^2/2 + d/2 + 1/6;
polySINC_R = d.^3/2 - d.^2 + 2/3;
polySINC_R_R = - d.^3/6 + d.^2/2 - d/2 + 1/6;
polySINC = [polySINC_L_L polySINC_L polySINC_R polySINC_R_R];

% PLOT polySINC
subplot(2,2,1), plot(polySINC, 'k','LineWidth',2);
set(gca,'fontsize',fontsize,'fontname',fontname);
xlabel({'Time (s)','(a)'}');
set(gca,'XTick', [1, dl 2*dl 3*dl 4*dl]); 
set(gca,'XTicklabel',{'-2T', '-T', '0', 'T', '2T'});
set(gca,'YTick', [0, 0.5 1]); 
set(gca,'YTicklabel',{'0', '0.5', '1'});
axis([1 length(polySINC) -0.2 1.2])
xL = xlim;
line(xL, [0 0], 'color','k');  %x-axis

% PLOYBLEP
polyBLEP_L_L = d.^4/24;
polyBLEP_L = - d.^4/8 + d.^3/6 + d.^2/4 + d/6 + 1/24;
polyBLEP_R = d.^4/8 - d.^3/3 + (2*d)/3 + 1/2;
polyBLEP_R_R = - d.^4/24 + d.^3/6 - d.^2/4 + d/6 + 23/24;
polyBLEP = [polyBLEP_L_L polyBLEP_L polyBLEP_R polyBLEP_R_R];

% PLOT POLYBLEP
subplot(2,2,2), plot(polyBLEP, 'k','LineWidth',2);
set(gca,'fontsize',fontsize,'fontname',fontname);
xlabel({'Time (s)','(b)'}');
set(gca,'XTick', [1, dl 2*dl 3*dl 4*dl]); 
set(gca,'XTicklabel',{'-2T', '-T', '0', 'T', '2T'});
set(gca,'YTick', [0, 0.5 1]); 
set(gca,'YTicklabel',{'0', '0.5', '1'});
axis([1 length(polyBLEP) -0.2 1.2])
xL = xlim;
line(xL, [0 0], 'color','k');  %x-axis

% PLOYBLAMP
polyBLAMP_L_L = d.^5/120;
polyBLAMP_L = -d.^5/40 + d.^4/24 + d.^3/12 + d.^2/12 + d./24 + 1/120;
polyBLAMP_R = d.^5/40 - d.^4/12 + d.^2/3 + d./2 + 7/30;
polyBLAMP_R_R = - d.^5/120 + d.^4/24 - d.^3/12 + d.^2/12 + (23*d)/24 + 121/120;
polyBLAMP = [polyBLAMP_L_L polyBLAMP_L polyBLAMP_R polyBLAMP_R_R];

% PLOT POLYBLAMP & TRIVIAL RAMP
subplot(2,2,3), plot(polyBLAMP, 'k','LineWidth',2), hold on;
set(gca,'fontsize',fontsize,'fontname',fontname);
set(gca,'XTick', [1, dl 2*dl 3*dl 4*dl]);
set(gca,'XTicklabel',{'-2T', '-T', '0', 'T', '2T'});
set(gca,'YTick', [0, 1 2]); 
set(gca,'YTicklabel',{'0', '1', '2'});
xlabel({'Time (s)','(c)'}');
axis([1 length(polySINC) -0.2 2.2]);
xL = xlim;
line(xL, [0 0], 'color','k');  %x-axis

% TRIVIAL RAMP
line([2000 3200],[0 1.2],'color','k','LineWidth',2,'LineStyle','--');
line([1 2000],[0 0],'color','k','LineWidth',2,'LineStyle','--');


% PLOYBLAMP
polyBLAMP_L_L = d.^5/120;
polyBLAMP_L = -d.^5/40 + d.^4/24 + d.^3/12 + d.^2/12 + d./24 + 1/120;
polyBLAMP_R = d.^5/40 - d.^4/12 + d.^2/3 - d./2 + 7/30;
polyBLAMP_R_R = -d.^5/120 + d.^4/24 - d.^3/12 + d.^2/12 - d./24 + 1/120;
polyBLAMP = [polyBLAMP_L_L polyBLAMP_L polyBLAMP_R polyBLAMP_R_R];

% PLOT POLYBLAMP RESIDUAL
subplot(2,2,4), plot(polyBLAMP, 'k','LineWidth',2);
set(gca,'fontsize',fontsize,'fontname',fontname);
set(gca,'XTick', [1 dl 2*dl 3*dl 4*dl]); 
set(gca,'XTicklabel',{'-2T', '-T', '0', 'T', '2T'});
set(gca,'YTick', [0, 7/30]); 
set(gca,'YTicklabel',{'0', '7/30'});
xlabel({'Time (s)','(d)'}');
axis([1 length(polySINC) -0.02 0.27]);
xL = xlim;
line(xL, [0 0], 'color','k');  %x-axis

% print -deps2 Fig-polyBLAMP-steps-N3

