clear;
%import traced data
data_quad1 = readtable('quadrant1.csv');
data_quad2 = readtable('quadrant2.csv');

num_transects_1 = size(data_quad1,2)/2;
num_transects_2 = size(data_quad2,2)/2;

%scientific notation correction factor
factor = 1000;

%store xy transects in cells
transects_1 = cell(num_transects_1, 1);
transects_2 = cell(num_transects_2,1);
max_length = 0;

%read from .csv 
for i = 1:num_transects_1
    col1 = i*2-1;
    col2 = i*2;

    data1 = table2cell(data_quad1(:,col1));
    data2 = table2cell(data_quad1(:,col2));
    
%check if data is separated by space or by column, read data for xy values from each case 
    if isnan(data2{4})
        sep = true;
    else
        sep = false;
    end

    if (sep)
        data1 = data1(3:end);
        x_1 = [];
        y_1 = [];
        for j = 1:length(data1)
            c_data_1 = data1{j};
            c_data_1 = convertCharsToStrings(c_data_1);
            if (strcmp(c_data_1,""))
                break;
            end
            str_split = strsplit(c_data_1,' ');
            c_x_1 = str2double(str_split(1));
            c_y_1 = str2double(str_split(2));

            x_1 = [x_1; c_x_1];
            y_1 = [y_1; c_y_1];
        end
        x_1 = x_1*factor;
        y_1 = y_1*factor;
    else
        data1 = data1(3:end);
        data2 = data2(3:end);
        data1 = cell2mat(data1)*factor;
        data2 = cell2mat(data2)*factor;

        x_1 = data1;
        y_1 = data2;
        sel = ~isnan(x_1);
        x_1 = x_1(sel);
        y_1 = y_1(sel);
    end
%store xy values in cell array 
    transects_1{i} = [x_1,y_1];
end

%repeat for second quadrant data
for i = 1:num_transects_2
    col1 = i*2-1;
    col2 = i*2;

    data3 = table2cell(data_quad2(:,col1));
    data4 = table2cell(data_quad2(:,col2));

    if isnan(data4{4})
        sep = true;
    else
        sep = false;
    end

    if (sep)
        data3 = data3(3:end);
        x_2 = [];
        y_2 = [];
        for j = 1:length(data3)
            c_data_2 = data1{j};
            c_data_2 = convertCharsToStrings(c_data_2);
            if (strcmp(c_data_2,""))
                break;
            end
            str_split = strsplit(c_data_2,' ');
            c_x_2 = str2double(str_split(1));
            c_y_2 = str2double(str_split(2));

            x_2 = [x_2; c_x_2];
            y_2 = [y_2; c_y_2];
        end
        x_2 = x_2*factor;
        y_2 = y_2*factor;
    else
        data3 = data3(3:end);
        data4 = data4(3:end);
        data3 = cell2mat(data3)*factor;
        data4 = cell2mat(data4)*factor;

        x_2 = data3;
        y_2 = data4;
        sel = ~isnan(x_2);
        x_2 = x_2(sel);
        y_2 = y_2(sel);
    end

    transects_2{i} = [x_2,y_2];
end

% extract all x and y values into separate NxM zero-padded matrices of x and y-coordinates of M circle segments 
M = size(transects_1,1);
maxN= size(transects_1{1}(:,1),1);
for i = 1:num_transects_1
    if size(transects_1{i}(:,1),1)>maxN
    maxN=size(transects_1{i}(:,1),1);
    elseif size(transects_2{i}(:,1),1)>maxN
    maxN=size(transects_2{i}(:,1),1);
    end
end
N=maxN;

all_x_1=nan(N,M);
all_x_2=nan(N,M);
all_y_1=nan(N,M);
all_y_2=nan(N,M);

for i=1:num_transects_1
    all_x_1(1:size(transects_1{i}(:,1),1),i)=transects_1{i}(:,1);
end

for i=1:num_transects_2
    all_x_2(1:size(transects_2{i}(:,1),1),i)=transects_2{i}(:,1);
end

%zero-padded NxM matrix of all x values 
all_x = [all_x_1;all_x_2];


for i=1:num_transects_1
    all_y_1(1:size(transects_1{i}(:,2),1),i)=transects_1{i}(:,2);
end

for i=1:num_transects_2
    all_y_2(1:size(transects_2{i}(:,2),1),i)=transects_2{i}(:,2);
end

%zero-padded NxM matrix of all x values 
all_y=[all_y_1;all_y_2];

centered_x=all_x-4.3445e+03;
centered_y=all_y-4.0627e+03;

for i = 1:size(centered_x,2)
    c_x = centered_x(:,i);
    c_y = centered_y(:,i);

    r = sqrt(c_x.^2 + c_y.^2);
    theta = atan2(c_y, c_x);

    [theta, idx] = sort(theta);
    r = r(idx);

    c_x = r.*cos(theta);
    c_y = r.*sin(theta);
    centered_x(:,i) = c_x;
    centered_y(:,i) = c_y;
end



% save('processedXY.mat','all_x','all_y')
save('CenteredProcessedXY.mat','centered_x','centered_y')

figure(1)
clf;

areas = NaN(size(centered_x,2),1);
ravg=[];
for i = 1:size(centered_x,2)
    x = centered_x(:,i);
    y = centered_y(:,i);

    sel = ~isnan(x);
    x = x(sel);
    y = y(sel);
    
    x = [x; x(1)];
    y = [y; y(1)];

    r = sqrt(x.^2+y.^2);
    ravg(i) = mean(r);


    hold on;
    p = plot(x,y,'k','LineWidth',1);
    hold off;

    c_a = polyarea(x,y);
    areas(i) = c_a;
   

end

ravg=ravg';
axis equal
axis ij
grid on


