function [p,proj] = CS_EX1_projections(KB)
% CS_EX1_projections - distance from H_n center for projections:
%     (1) onto axes; (2) onto pairise 45 degrees; (3) all axes 45 degrees
% On input:
%     KB (KB struct vector): KB
% On output:
%     p (2n+2(n choose k)+2 x1 vector): distance of proj pt to H_n center
% Call:
%     p = CS_EX1_projections(KB3);
% Author:
%     T. Henderson
%     UU
%     Fall 2022
%

h = CS_KB2hp(KB);
[mp2n,np1] = size(h);
n = np1 - 1;
num_proj = 2*n + 2*nchoosek(n,2) + 2;
p = zeros(num_proj,1);
proj = zeros(num_proj,n);
xmin = zeros(n,1);
xmax = ones(n,1);
options = optimoptions('linprog','Algorithm','interior-point',...
    'Display','none');
pc = 0.5*ones(n,1);

% project onto coordinate axes
for d = 1:n
    index = 2*(d-1) + 1;
    e = zeros(n,1);
    e(d) = 1;
    e = e - pc;
    x = linprog(e,-h(:,1:end-1),h(:,end),[],[],xmin,xmax,options);
    if ~isempty(x)
        p(index) = norm(pc-x);
        proj(index,:) = x';
    end
    e(d) = -1;
    x = linprog(e,-h(:,1:end-1),h(:,end),[],[],xmin,xmax,options);
    if ~isempty(x)
        p(index) = norm(pc-x);
        proj(index,:) = x';
    end
end

% project onto pairwise 45 degree axes
pairs = nchoosek([1:n],2);
[num_pairs,~] = size(pairs);
for k = 1:num_pairs
    index = 2*(k-1) + 1 + 2*n;
    e1 = pairs(k,1);
    e2 = pairs(k,2);
    e = zeros(n,1);
    e(e1) = 1;
    e(e2) = 1;
    e = e - pc;  % new idea for direction
    e = e/norm(e);
    x = linprog(e,-h(:,1:end-1),h(:,end),[],[],xmin,xmax,options);
    if ~isempty(x)
        p(index) = norm(pc-x);
        proj(index,:) = x';
    end
    e = -e;
    x = linprog(e,-h(:,1:end-1),h(:,end),[],[],xmin,xmax,options);
    if ~isempty(x)
        p(index+1) = norm(pc-x);
        proj(index+1,:) = x';
    end
end

% project onto all-axis 45 degree axis
e = ones(n,1);
e = e/norm(e);
x = linprog(e,-h(:,1:end-1),h(:,end),[],[],xmin,xmax,options);
p(end-1) = norm(pc-x);
proj(end-1,:) = x';
e = -e;
x = linprog(e,-h(:,1:end-1),h(:,end),[],[],xmin,xmax,options);
p(end) = norm(pc-x);
proj(end,:) = x';

tch = 0;
