function [ max_response, lowest_dist_x, lowest_dist_y, lowest_dist_scale ] = calcpatch( im1, im2 )
%calcpatch Calculates the best image patch of im2 inside im1

[rows, cols, channels] = size(im1);

im1blur(:,:,1) = filter2(fspecial('disk',10),im1(:,:,1));
im1blur(:,:,2) = filter2(fspecial('disk',10),im1(:,:,2));
im1blur(:,:,3) = filter2(fspecial('disk',10),im1(:,:,3));

%im1 = im1blur;

[im1gx,im1gy] = gradient(im1blur);

im1g = im1gx + im1gy;

lowest_dist = 0;
lowest_dist_x = -1;
lowest_dist_y = -1;
lowest_dist_scale = -1;

% edit these to change starting scale and number of iterations
starting_scale = 0.1;
num_iters = 1;

% do not change these
iter = 1;
scale = starting_scale;

while (iter <= num_iters)

%fprintf('scale = %f\n', scale);

% generate our downsampled target patch
patch = imresize(im2, scale, 'bilinear');

[patchgx,patchgy] = gradient(patch);
patchg = patchgx + patchgy;

%{
patchblur(:,:,1) = filter2(fspecial('disk',5),patch(:,:,1));
patchblur(:,:,2) = filter2(fspecial('disk',5),patch(:,:,2));
patchblur(:,:,3) = filter2(fspecial('disk',5),patch(:,:,3));

patch = patchblur;
%}

[prows, pcols, pchannels] = size(patch);

dr = normxcorr2(patch(:,:,1), im1(:,:,1));
dg = normxcorr2(patch(:,:,2), im1(:,:,2));
db = normxcorr2(patch(:,:,3), im1(:,:,3));

gdr = normxcorr2(patchg(:,:,1), im1g(:,:,1));
gdg = normxcorr2(patchg(:,:,2), im1g(:,:,2));
gdb = normxcorr2(patchg(:,:,3), im1g(:,:,3));

dsum = 1.0 * (dr + dg + db) + 1.0 * (gdr + gdg + gdb);

%{
figure
imshow(dr)
figure
imshow(dg)
figure
imshow(db)

figure
imshow(gdr)
figure
imshow(gdg)
figure
imshow(gdb)
%}

osx = round(pcols/2);
osy = round(prows/2);

dsum = dsum(osy+1:osy+rows,osx+1:osx+cols);

% Eliminate any possibilites that would result in the patch being partially
% off any of the edges
validRegion = zeros(size(dsum));
validRegion(osy:rows-osy, osx:cols-osx,:) = dsum(osy:rows-osy, osx:cols-osx,:);
dsum = validRegion;

max_response = max(dsum(:));

[y x] = find(dsum == max_response);

%max_response = max_response * sqrt(scale);

if (max_response > lowest_dist)

    lowest_dist = max_response;
    
    lowest_dist_x = x;%round(x - pcols/2);
    lowest_dist_y = y;%round(y - prows/2);
    
    lowest_dist_scale = scale;
    
    %fprintf('new highest response: %f at (%i, %i) with scale %f\n', lowest_dist, lowest_dist_x, lowest_dist_y, lowest_dist_scale);
end

%fprintf('\n');

iter = iter + 1;

scale = scale + 0.01;

end % end while (iter <




end












