close all; clear; clc;
seed = 28; rng(seed);

%% Fit window_lossfun many times from random starting points
N = 1024; n = 0:(N-1); % Window axis
Nf = 100; % Number of refits
pfs = NaN(Nf,3); fvals=NaN(Nf,1); % Preallocate
for j = 1:Nf
    p0 = rand(1,2); % Random starting point
    [pf,fvals(j,:)] = fminsearch(@window_lossfun,p0);
    [w,pf] = p_to_window(pf,N); % Convert to properly scaled window
    pfs(j,:) = pf;
    if ~mod(j,round(Nf/20)); fprintf('%02d%% done\n',round(j/Nf*100)), end
end
[fvals,I] = sort(fvals); pfs = pfs(I,:); % Sort results according to fval

%% Draw hypersurface
%
figure(1); clf; colormap(plasma);
p1s = linspace(0,1,301); p2s = linspace(0,1,300); % p1 and p2 vectors
L = NaN(numel(p1s),numel(p2s)); % Preallocate

for j = 1:numel(p1s)
    for k = 1:numel(p2s)
        L(j,k) = window_lossfun([p1s(j),p2s(k)]);
    end
    if ~mod(j,round(numel(p1s)/20))
        imagesc(p2s,p1s,L);
        axis square; grid on; caxis([-8 -4]);
        xlabel('p2'); ylabel('p1');
        drawnow;
        fprintf('%02d%% done\n',round(j/numel(p1s)*100));
    end
end

%% Cluster results into Nm groups
Nm = 4;
I = clusterdata(pfs,'distance','chebychev','maxclust',Nm);
I2 = zeros(size(I));
for j = 1:Nm
    Ifvals(j) = mean(fvals(I==j));
end
[~,II] = sort(Ifvals);
for j = 1:Nm
    I2(I==II(j)) = j;
end
I = I2;

%% Plot
x = 1:numel(fvals);
cm = lines;
cm = [cm(1,:); cm(4,:); cm(2,:); cm(5,:); cm(3,:)];

figure(2); clf;
s1 = subplot(1,2,1); hold on;
plot(x,fvals,'k-');
for j = 1:Nm
    plot(x(I==j),fvals(I==j),'.','color',cm(j,:),'MarkerSize',10);
end
grid on; box on; axis square;
s1.XTickLabel = '';
xlabel('sorted run');
ylabel('phi');

% Calculate mean windows
N = 1024; n = 0:(N-1);

b = bsxfun(@times,cos((0:2:4)'*pi*n/(N-1)),sign(mod(1:3,2)-1/2)');
w1 = mean(pfs(I==1,:))*b;
w2 = mean(pfs(I==2,:))*b;
w3 = mean(pfs(I==3,:))*b;
w4 = mean(pfs(I==4,:))*b;
w5 = mean(pfs(I==5,:))*b;

% Plot windows
s2 = subplot(1,2,2); hold on;
semilogy(n,abs(fft(flattopwin(N))),'k','LineWidth',1);
semilogy(n,abs(fft(w1)),'.-','color',cm(1,:),'LineWidth',1,'MarkerSize',10);
%semilogy(n,abs(fft(w2)),'.-','color',cm(2,:),'LineWidth',1,'MarkerSize',10);
semilogy(n,abs(fft(w3)),'.-','color',cm(3,:),'LineWidth',1,'MarkerSize',10);
semilogy(n,abs(fft(w4)),'.-','color',cm(4,:),'LineWidth',1,'MarkerSize',10);
s2.YScale = 'log';
xlim([0 10]);
grid on; box on;
s2.MinorGridLineStyle = 'none';
s2.YTickLabel = ''; s2.XTickLabel = '';
xlabel('f'); ylabel('fft(w)');
axis square;