% First Implementation Attempting to Define Nucleus Space
% May 24, 2007
% (c) 2007 Christopher Mitchell
% All rights reserved

disp('Stage 1...');
% First we need to load a bunch of images from file.
% Note: this program calculates dynamic array sizes, but it assumes
%       that all images have the same (square) dimensions.  Make it so!
RGB1 = rgb2gray(imread('./v1/cell1.png'));
RGB2 = rgb2gray(imread('./v1/cell2.png'));
RGB3 = rgb2gray(imread('./v1/cell3.png'));
RGB4 = rgb2gray(imread('./v1/cell4.png'));
RGB5 = rgb2gray(imread('./v1/cell5.png'));
RGB6 = rgb2gray(imread('./v1/cell6.png'));
RGB7 = rgb2gray(imread('./v1/cell7.png'));
RGB8 = rgb2gray(imread('./v1/cell8.png'));

disp('Stage 2...');
% Next we need to create a big matrix containing all of the data we just
% loaded as 8 matrices
[dimY,dimX] = size(RGB1);
RawData = zeros(dimY*dimX,dimY*dimX);
for i = 1:dimY
    for j = 1:dimX
        RawData(1,((i-1)*dimX)+j) = RGB1(i,j);
        RawData(2,((i-1)*dimX)+j) = RGB2(i,j);
        RawData(3,((i-1)*dimX)+j) = RGB3(i,j);
        RawData(4,((i-1)*dimX)+j) = RGB4(i,j);
        RawData(5,((i-1)*dimX)+j) = RGB5(i,j);
        RawData(6,((i-1)*dimX)+j) = RGB6(i,j);
        RawData(7,((i-1)*dimX)+j) = RGB7(i,j);
        RawData(8,((i-1)*dimX)+j) = RGB8(i,j);
    end
end

disp('Stage 3...');
% Step 3: Compute and sort eigenvectors and eigenvalues
[eigenVectors,eigenValues] = eig(RawData);
eigenValues = diag(eigenValues);
eigenVectors(:,dimX*dimY+1) = eigenValues;
eigenVectors = sortrows(eigenVectors,dimX*dimY+1);
eigenVectors = eigenVectors(:,1:dimX*dimY);
disp('Stage 4...');
disp('Completed.');

% eigenCell:
for i = 1:dimY
    for j = 1:dimX
        eigenCell1(i,j) = eigenVectors(1,j+dimX*(i-1));
        eigenCell2(i,j) = eigenVectors(2,j+dimX*(i-1));
        eigenCell3(i,j) = eigenVectors(3,j+dimX*(i-1));
        eigenCell4(i,j) = eigenVectors(4,j+dimX*(i-1));
        eigenCell5(i,j) = eigenVectors(5,j+dimX*(i-1));
        eigenCell6(i,j) = eigenVectors(6,j+dimX*(i-1));
        eigenCell7(i,j) = eigenVectors(7,j+dimX*(i-1));
        eigenCell8(i,j) = eigenVectors(8,j+dimX*(i-1));
    end
end
imwrite(0.5+0.5*eigenCell1,'./v1/eig1.png','png');
imwrite(0.5+0.5*eigenCell2,'./v1/eig2.png','png');
imwrite(0.5+0.5*eigenCell3,'./v1/eig3.png','png');
imwrite(0.5+0.5*eigenCell4,'./v1/eig4.png','png');
imwrite(0.5+0.5*eigenCell5,'./v1/eig5.png','png');
imwrite(0.5+0.5*eigenCell6,'./v1/eig6.png','png');
imwrite(0.5+0.5*eigenCell7,'./v1/eig7.png','png');
imwrite(0.5+0.5*eigenCell8,'./v1/eig8.png','png');