% Code for rendering three dimensional results generated by the gpca% algorithm.  Users should call plot_subspacesfunction pl = plot_subspaces3(X, labels, groupBases, basisDimensions)GRID_LINE_COUNT = 15;% Generate a list, existList, of labels that have had points assigned to them.highestLabel = max(labels);existList = zeros(highestLabel, 1);  % existList will have at most highestLabel elements.groupCount = 0;for labelIndex = 1:highestLabel,    if (sum(labels == labelIndex) > 0) % If the group has at least one element in it,        groupCount = groupCount + 1;        existList(groupCount) = labelIndex;    end;end;existList = existList(1:groupCount);markers = ['.ox+*svp^h<>']; markerCount = length(markers); markerIndex = 1;colors = ['rbgmcyk']; colorCount = length(colors); colorIndex = 1;% Iterate through every group that contains at least one element.for groupIndex = 1:groupCount,            basis = groupBases{existList(groupIndex)}(:,1:basisDimensions(existList(groupIndex)));     groupData = X(:, labels == existList(groupIndex));    if (size(basis,2) == 1)        % Plot the one dimensional subspace (line)        projected_data = groupData' * basis;        upper_bound = max(projected_data) * basis;        lower_bound = min(projected_data) * basis;        l = line([upper_bound(1) lower_bound(1)], [upper_bound(2) lower_bound(2)], [upper_bound(3) lower_bound(3)]);        set(l, 'Color', colors(colorIndex));        set(l, 'LineWidth', 3);            elseif (size(basis,2) == 2)        % Plot the plane as a surface object                % Compute a set of basis vectors that align with the plot axes....        % (hack)        row_d(1) = norm(basis(1,:));        row_d(2) = norm(basis(2,:));        row_d(3) = norm(basis(3,:));        [ignore sort_index] = sort(row_d);        sort_index = fliplr(sort_index);        basis(sort_index,:) = rref(basis(sort_index,:)')';        basis = normalize(basis);                % Compute Bounds for the data.            projected_data = (X') * basis;                        % The data bounds in the new basis.            max_pd = max(projected_data, [], 1);            min_pd = min(projected_data, [], 1);%             %             % The bounds as vectors in 3-space.%             ub1 = max_pd(1) * basis(:,1);%             ub2 = max_pd(2) * basis(:,2);%             lb1 = min_pd(1) * basis(:,1);%             lb2 = min_pd(2) * basis(:,2);                % Compute an array of values in three space in an order that can be        % used for generating a surface plot.        ub1 = max_pd(1);        ub2 = max_pd(2);        lb1 = min_pd(1);        lb2 = min_pd(2);        tempArray = zeros(GRID_LINE_COUNT, GRID_LINE_COUNT, 3);  % Pages of this array are X, Y, and Z values.        AIndex = 1;        for A = linspace(lb1, ub1, GRID_LINE_COUNT),            BIndex = 1;            for B = linspace(lb2, ub2, GRID_LINE_COUNT),                tempArray(AIndex, BIndex, :) = A * basis(:,1) + B * basis(:,2);                BIndex = BIndex + 1;            end            AIndex = AIndex + 1;        end                % Define vectors for the various group colors.        % colors = 'rbkmgcy'        groupColors = [1 0 0; 0 1 0; 0 0 0; 1 0 1; 0 1 0; 0 1 1; 1 1 0];        %CData = groupColors(colorIndex,:);        %        h = surf(tempArray(:,:,1), tempArray(:,:,2), tempArray(:,:,3));        h = surf(tempArray(:,:,1), tempArray(:,:,2), tempArray(:,:,3));        set(h,'FaceColor',colors(colorIndex),'EdgeColor','black')        set(h, 'EdgeAlpha', .2);        alpha .2;        %         for i = 0:0.1:1%             %% Draw the horizontal lines%             l1 = line([ub1(1)+i*ub2(1)+(1-i)*lb2(1) lb1(1)+i*ub2(1)+(1-i)*lb2(1)], ... %                 [ub1(2)+i*ub2(2)+(1-i)*lb2(2) lb1(2)+i*ub2(2)+(1-i)*lb2(2)], ...%                 [ub1(3)+i*ub2(3)+(1-i)*lb2(3) lb1(3)+i*ub2(3)+(1-i)*lb2(3)]);%             set(l1, 'Color', colors(colorIndex));%             %set(l1, 'LineWidth', 2.5);    %             %% Draw the vertical lines%             l2 = line([i*ub1(1)+(1-i)*lb1(1)+ub2(1) i*ub1(1)+(1-i)*lb1(1)+lb2(1)], ...%                 [i*ub1(2)+(1-i)*lb1(2)+ub2(2) i*ub1(2)+(1-i)*lb1(2)+lb2(2)], ...%                 [i*ub1(3)+(1-i)*lb1(3)+ub2(3) i*ub1(3)+(1-i)*lb1(3)+lb2(3)]);%             set(l2, 'Color', colors(colorIndex));%             %set(l2, 'LineWidth', 2.5);    %         end    end        % Advance to the next color and marker.    markerIndex = markerIndex + 1; if (markerIndex > markerCount) markerIndex = 1; end;    colorIndex = colorIndex+1; if (colorIndex > colorCount) colorIndex = 1; end;end;xlabel('axis 1')ylabel('axis 2')zlabel('axis 3')pl = [];