function inverted = NE_inversion_experiment1(max_count)
% NE_inversion_experiment1 - tile PD with regular pentagon
%    put pentagon in center of PD and invert through sides recursively
% On input:
%     max_count (int): maximum number of pentagons to create (non-uniquely)
% On output:
%     inverted (nx10 array): linearized set of unique pentagons created
% Call:
%     inv = NE_inversion_experiment1(40);
% Author:
%     T. Henderson
%     UU
%     Summer 2023
%

ZERO_THRESH = 0.00001;

poly = NE_gen_pentagon(0);
corners_PH = poly(1).corners2;
[num_corners,~] = size(corners_PH);
corners_PD = zeros(num_corners,2);
for k = 1:num_corners
    [c1,c2] = NE_H2D(corners_PH(k,1),corners_PH(k,2));
    corners_PD(k,1) = c1;
    corners_PD(k,2) = c2;
end
% boundaries of square
sq_pts1 = NE_PD_seg_pts([-sqrt(2)/2,sqrt(2)/2],[sqrt(2)/2,sqrt(2)/2]);
sq_pts2 = NE_PD_seg_pts([sqrt(2)/2,sqrt(2)/2],[sqrt(2)/2,-sqrt(2)/2]);
sq_pts3 = NE_PD_seg_pts([sqrt(2)/2,-sqrt(2)/2],[-sqrt(2)/2,-sqrt(2)/2]);
sq_pts4 = sq_pts2;
sq_pts4(:,1) = -sq_pts4(:,1);
t = [0:0.1:2*pi,2*pi];
x = cos(t);
y = sin(t);

clf
plot(x,y,'k');
hold on
plot(sq_pts1(:,1),sq_pts1(:,2),'b.');
plot(sq_pts2(:,1),sq_pts2(:,2),'b.');
plot(sq_pts3(:,1),sq_pts3(:,2),'b.');
plot(sq_pts4(:,1),sq_pts4(:,2),'b.');
axis equal
xx = 0.4143*cos(t);
yy = 0.4143*sin(t);
plot(xx,yy,'k');
plot(xx,yy,'k.');

inverted = [];
queue(1).corners = corners_PD;

done = 0;
count = 0;
while done==0
    count = count + 1
    pent = queue(1).corners;
    queue(1) = [];
    len_queue = length(queue);
    new_pent = 1;
    if NE_already_inverted_pent(pent,inverted)
        new_pent = 0;
    else
        inverted = [inverted;pent(:)'];
    end
    
    if new_pent==1
        % draw current pentagon
        p_pts1 = NE_PD_seg_pts(pent(1,:),pent(2,:));
        p_pts2 = NE_PD_seg_pts(pent(2,:),pent(3,:));
        p_pts3 = NE_PD_seg_pts(pent(3,:),pent(4,:));
        p_pts4 = NE_PD_seg_pts(pent(4,:),pent(5,:));
        p_pts5 = NE_PD_seg_pts(pent(5,:),pent(1,:));
        plot(p_pts1(:,1),p_pts1(:,2),'r');
        plot(p_pts2(:,1),p_pts2(:,2),'r');
        plot(p_pts3(:,1),p_pts3(:,2),'r');
        plot(p_pts4(:,1),p_pts4(:,2),'r');
        plot(p_pts5(:,1),p_pts5(:,2),'r');
        drawnow;
        
        % inversion through face 1
        circle1 = NE_orthog_circle(pent(1,:),pent(2,:));
        p1i = NE_inversion_all(pent(1,:),pent(2,:),circle1,pent);
        len_queue = len_queue + 1;
        queue(len_queue).corners = p1i;
        
        % inversion through face 2
        circle2 = NE_orthog_circle(pent(2,:),pent(3,:));
        p2i = NE_inversion_all(pent(2,:),pent(3,:),circle2,pent);
        len_queue = len_queue + 1;
        queue(len_queue).corners = p2i;
        
        % inversion through face 3
        circle3 = NE_orthog_circle(pent(3,:),pent(4,:));
        p3i = NE_inversion_all(pent(3,:),pent(4,:),circle3,pent);
        len_queue = len_queue + 1;
        queue(len_queue).corners = p3i;
        
        % inversion through face 4
        circle4 = NE_orthog_circle(pent(4,:),pent(5,:));
        p4i = NE_inversion_all(pent(4,:),pent(5,:),circle4,pent);
        len_queue = len_queue + 1;
        queue(len_queue).corners = p4i;
        
        % inversion through face 5
        circle5 = NE_orthog_circle(pent(5,:),pent(1,:));
        p5i = NE_inversion_all(pent(5,:),pent(1,:),circle5,pent);
        len_queue = len_queue + 1;
        queue(len_queue).corners = p5i;
        
        if count>=max_count
            done = 1;
        end
        tch = 0;
    end
end

tch = 0;
