% Creates a Gaussian pyramid from the source image in by 
% repeated applying a Gaussian low-pass kernel.  nlvl specifies
% the desired number of pyramid levels.  The most blurry level
% is last entry in pyr, and the original image is the first entry.

function [pyr] = make_gpyr(in, nlvl)

% Create output cell array
pyr = cell(nlvl,1);

% Copy source image as first level
pyr{1} = double(in);

% Build each successive level by blurring the previous one
% and reducing its size.
for i=2:nlvl
  %pyr{i} = imresize(gp_blur(pyr{i - 1}), 0.5);
  % Don't downsample
  %pyr{i} = gp_blur(pyr{i - 1});
  pyr{i} = new_blur(in,i);
end

% Blurs an image with the same kernel used in when creating
% a Gaussian pyramid, but does not reduce the image.  Returns
% blurred image as a double type, not uint8.

function [out] = gp_blur(in)

Gk = [.05 .25 .4 .25 .05];

out = imfilter(double(in), Gk'*Gk, 'replicate');

function [out] = new_blur(in, lvl)
boxlen = 2^(lvl-1)+1;
out = conv2(in,ones(boxlen),'same');
border = (boxlen-1)/2;
[y,x] = size(in);
out(1:border,:) = NaN;
out(y-border+1:y,:) = NaN;
out(:,1:border) = NaN;
out(:,x-border+1:x) = NaN;