% Third Implementation (Search Subprogram)
% May 29-31, 2007
% (c) 2007 Christopher Mitchell
% All rights reserved

% First we need to load a bunch of images from file.
% Load them from the preexisting training file
load eigenCells eigenCells
eigenSize = size(eigenCells);
eigenN = eigenSize(1);
eigenX = input('Width of training images (px): ');
eigenY = input('Height of training images (px): ');
threshold = input('Threshold: ');
if eigenX*eigenY ~= eigenSize(2)
    disp('Dimensional mismatch in training set.');
    return  %This would be a dim mismatch
end

% Set up matrices beforehand for speed
s = input('Dimensions match.  Now we need the test image: ','s');
testImage = imread(s);
[IdimY,IdimX] = size(testImage);
matches = zeros(IdimY-eigenY+1,IdimX-eigenX+1);
subimage = zeros(eigenY,eigenX);
subimagevector = zeros(1,eigenX*eigenY);
cellMatch = zeros(eigenN);

% Now do the tests for each possible subimage.  This should be a LONG
% process since we're operating on a 250,000-element matrix. :P
tic;
for X = 1:size(matches,2)
    for Y = 1:size(matches,1)
        subimage = testImage(Y:Y+eigenY-1,X:X+eigenX-1);
%        if max(max(subimage)) ~= min(min(subimage))
%            subimage = subimage./(max(max(subimage))-min(min(subimage)));
%        end
%        subimage = subimage - min(min(subimage));
        for i = 1:eigenY
            for j = 1:eigenX
                subimagevector(1,((i-1)*eigenX)+j) = subimage(i,j);
            end
        end
        %now do the facespace transformations
        
        %% Commented-out percentage printing code (needs buffer flushing) %%
        %if floor((X*size(matches,1)+Y)/1000) == (X*size(matches,1)+Y)/1000
        %    s = sprintf('%d%% complete',floor(100*((X*size(matches,1)+Y)/(size(matches,1)*size(matches,2)))));
        %    disp(s);
        %end
        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        
        %To face space
        for i=1:eigenN
            cellMatch(i) = inner(subimagevector,eigenCells(i,:));
        end
        
        %Back to image space
        subimagevector2 = zeros(1,eigenX*eigenY);
        for i=1:eigenN
            subimagevector2 = subimagevector2 + cellMatch(i)*eigenCells(i,:);
        end
        
        %Find the MSE
        if max(max(subimagevector2))~=min(min(subimagevector2))
            matches(Y,X) = sum(((subimagevector2-min(min(subimagevector2)))/(max(max(subimagevector2))-min(min(subimagevector2)))-subimagevector).^2);
        else
            matches(Y,X) = threshold.^2;
        end
    end
end
toc;

matchesBox = ones(IdimY-eigenY+1,IdimX-eigenX+1);
maxval = max(max(matches));

RGBimg = zeros(size(testImage,1),size(testImage,2),3);
RGBimg(:,:,1) = testImage;
RGBimg(:,:,2) = testImage;
RGBimg(:,:,3) = testImage;

while min(min(matches)) < threshold
    [Y,I] = min(matches);
    [X,J] = min(Y);
    coordX = J;
    coordY = I(J);
    s = sprintf('(%d,%d) @ MSE threshold of %d',coordX,coordY,X);
    disp(s);
    for X = coordX:coordX+eigenX-1
        for Y = coordY:coordY+eigenY-1
            if (X>0 && Y>0 && X<=size(RGBimg,2) && Y<=size(RGBimg,1)) && (X == coordX || X == coordX+eigenX-1 || Y == coordY || Y == coordY+eigenY-1)
               RGBimg(Y,X,2) = 255;
               RGBimg(Y,X,1) = 0;
               RGBimg(Y,X,3) = 0;
            end
        end
    end
    for X = coordX-floor(eigenX/2):coordX+floor(eigenX/2)
        for Y = coordY-floor(eigenY/2):coordY+floor(eigenY/2)
            if X>0 && Y>0 && X<=size(matches,2) && Y<=size(matches,1)
               matches(Y,X) = maxval;
            end
        end
    end
end

matches = matches-min(min(matches));
matches = matches./(max(max(matches)) - min(min(matches)));
imshow(RGBimg./255);

%All done here!