% we add 4 level noises namely 
e{1} = randn(129,1);
for j=2:5,
    e{j} = e{j-1}/2;
end
dx = 1/2^7;
x=0:dx:1;
[y{1},y{2},y{3},y{4}] = testgen(x);
for j=1:4,
    dev(j) = norm(y{j}(2:end,2)-y{j}(1:end-1,2));
end    
% 1) test relative L1 L2 , L_inf norm against true signals
% 2) do the same comparison for 1st and 2nd derivatives
% 3) plot the error curves
% 4) randomly sample 257 points from [0,1] repeat 1-3

% generate random points from [0,1] and compute the true values
rx = rand(257,1);
rx = sort(rx);
[ry{1},ry{2},ry{3},ry{4}] = testgen(rx);

% error test
for j = 1:5,
    for k = 1:4,
        sig = y{k}(:,2)+e{j};
        [c,D1,D0,f1,f0] = phlstfit(sig,2);
        [f{j}{k},df{j}{k},d2f{j}{k}]=optView(c,D1,D0,f1,f0,x);
        [rf{j}{k},rdf{j}{k},rd2f{j}{k}]=optView(c,D1,D0,f1,f0,rx);
        %grid points error
        errf{j}{k} = f{j}{k}-y{k}(:,2);
        errdf{j}{k} = df{j}{k}-y{k}(:,3);
        errd2f{j}{k} = d2f{j}{k}-y{k}(:,4);
        %random points error
        rerrf{j}{k} = rf{j}{k}-ry{k}(:,2);
        rerrdf{j}{k} = rdf{j}{k}-ry{k}(:,3);
        rerrd2f{j}{k} = rd2f{j}{k}-ry{k}(:,4);
        %L1,L2,Linf errors of f,df,d2f
        fL1{j}{k} = norm(errf{j}{k},1)/norm(y{k}(:,2),1);
        fL2{j}{k} = norm(errf{j}{k})/norm(y{k}(:,2));
        fLinf{j}{k} = norm(errf{j}{k},inf);
        rfL1{j}{k} = norm(rerrf{j}{k},1)/norm(ry{k}(:,2),1);
        rfL2{j}{k} = norm(rerrf{j}{k})/norm(ry{k}(:,2));
        rfLinf{j}{k} = norm(rerrf{j}{k},inf);
        dfL1{j}{k} = norm(errdf{j}{k},1)/norm(y{k}(:,3),1);
        dfL2{j}{k} = norm(errdf{j}{k})/norm(y{k}(:,3));
        dfLinf{j}{k} = norm(errdf{j}{k},inf);
        rdfL1{j}{k} = norm(rerrdf{j}{k},1)/norm(ry{k}(:,3),1);
        rdfL2{j}{k} = norm(rerrdf{j}{k})/norm(ry{k}(:,3));
        rdfLinf{j}{k} = norm(rerrdf{j}{k},inf);
        d2fL1{j}{k} = norm(errd2f{j}{k},1)/norm(y{k}(:,4),1);
        d2fL2{j}{k} = norm(errd2f{j}{k})/norm(y{k}(:,4));
        d2fLinf{j}{k} = norm(errd2f{j}{k},inf);
        rd2fL1{j}{k} = norm(rerrd2f{j}{k},1)/norm(ry{k}(:,4),1);
        rd2fL2{j}{k} = norm(rerrd2f{j}{k})/norm(ry{k}(:,4));
        rd2fLinf{j}{k} = norm(rerrd2f{j}{k},inf);
    end
end
