% RES = BandH(image)
% 
% Computes the texture generation procedure of Bergen and Heeger.
%
% Copyright Paul Viola, 6/97

function matched = BandH(im)

% Initialize the filters

filts = 'sp5Filters';
[lo0filt,hi0filt,lofilt,bfilts,steermtx,harmonics] = eval(filts);
fsz = round(sqrt(size(bfilts,1))); fsz =  [fsz fsz];
nfilts = size(bfilts,2);
nrows = floor(sqrt(nfilts));
nbins = 100;

% Generate a random image and match the histograms
clf;
[N,X] = histo(im(:),nbins);
[mn, mx] = range2(im);
matched = histoMatch(rand(size(im)), N, X);
showIm(im + sqrt(-1)*matched);


%% build the intial Steerable pyramid:

nlevels = maxPyrHt(size(im),size(lofilt));
[N,X] = histo(im(:),nbins);
[ipyr,ipind] = buildSpyr(im, nlevels, filts);
[mn, mx] = range2(im);

% Iteratively match histograms

for iter = 1:5
    [rpyr,rpind] = buildSpyr(matched, nlevels, filts);
    
    target = ipyr(pyrBandIndices(ipind,1));
    current = rpyr(pyrBandIndices(rpind,1));
    [tN, tX] = histo(target, nbins/1);
    histCurrent = histoMatch(current, tN, tX);
    rpyr(pyrBandIndices(rpind,1)) = histCurrent;

    nbands = spyrNumBands(ipind);

    for level = 1:nlevels
      for band = 1:nbands
	rband = 1 + band + nbands*(level-1);
	target = ipyr(pyrBandIndices(ipind,rband));
	current = rpyr(pyrBandIndices(rpind,rband));
	[tN, tX] = histo(target, nbins/level);
	histCurrent = histoMatch(current, tN, tX);
	rpyr(pyrBandIndices(rpind,rband)) = histCurrent;
      end
    end

     res = reconSpyr(rpyr, rpind, filts); 

    matched = histoMatch(res, N, X);
end
showIm(im + sqrt(-1)*matched);


