% clear all;
% close all;

% filtersSelection_codebook is for outputing the matching pursuit version of template, Frame version of template, and the deformed templates for each
% training images


function [template, templateMP, sketchTemplate] = filtersSelection_codebook(I, GaborScaleList, DoGScaleList, numOrient, numSketch, locationShiftLimit, orientShiftLimit, sx, sy, LocalNormOrNot,scaleFactor,thresholdFactor,savingFolder)

%% set parameters 

numScale = length(GaborScaleList); %number of scales for Gabor wavelets
numDoG = length(DoGScaleList);
numFilter = numScale*numOrient+numDoG;
epsilon = .1; 
Upperbound = 6.; 
SHUTUP = 0.; 

numImg = length(I);

%% read training images
% imgFold = 'Image';
% imgName = dir([imgFold '/*.jpg']);
% numImg = size(imgName, 1);
% clear I;
% for (i = 1 : numImg)
%     img = imread([imgFold '/' imgName(i).name]);
%     img = imresize(img, [sx,sy],'nearest');
%     img = im2double(rgb2gray(img));
%     img=img-mean(img(:));
%     img=img/std(img(:));
%     I{i} = img;
% end


%% filter images
tic;
[allFilterR, allFilterI, allSymbol] = MakeFilterBank(GaborScaleList, DoGScaleList, numOrient);
disp(['finished generating filter bank: '  num2str(toc) ' seconds']);

tic;
[Crr, Cri, Cir, Cii] = CorrFilter(allFilterR, allFilterI); %calculate correlation among filters
disp(['finished computing filters correlation: ' num2str(toc) ' seconds']);

tic;
[allFilteredR, allFilteredI, allFiltered] = ApplyFilter(I, allFilterR, allFilterI);
disp(['finished images filtering: ' num2str(toc) ' seconds']);

half = zeros(1, numOrient*numScale+numDoG); 
for (i=1:numOrient*numScale+numDoG)
    half(i) = (size(allFilterR{1, i}, 1)-1)/2;
end

% local normalization
% % if LocalNormOrNot
% %     % 
% %     
% %     for iI = 1 : numImg
% %       for iScale=1:numScale  
% %         h = half(1+(iScale-1)*numOrient); 
% %         
% %         [allFilteredR(iI, 1+(iScale-1)*numOrient : numOrient+(iScale-1)*numOrient)]...
% %             = LocalNormalizeDouble(allFilteredR(iI, 1+(iScale-1)*numOrient : numOrient+(iScale-1)*numOrient),[],h,round(scaleFactor*h),round(scaleFactor*h),thresholdFactor);    
% %         
% %         [allFilteredI(iI, 1+(iScale-1)*numOrient : numOrient+(iScale-1)*numOrient)]...
% %             = LocalNormalizeDouble(allFilteredI(iI, 1+(iScale-1)*numOrient : numOrient+(iScale-1)*numOrient),[],h,round(scaleFactor*h),round(scaleFactor*h),thresholdFactor);
% %       end
% %     
% %       % for DoG
% %       for iScale=1:numDoG  
% %         ind  = numOrient*numScale+iScale;
% %         h = half(ind);       
% %         [allFilteredR(iI, ind)]=LocalNormalizeDouble(allFilteredR(iI, ind),[],h,round(scaleFactor*h),round(scaleFactor*h),thresholdFactor);                
% %       end
% %     end
% %     
% %     for iI = 1 : numImg
% %         for iF = 1 : numFilter
% %              allFiltered{iI,iF}=(allFilteredR{iI,iF}.^2 + allFilteredI{iI,iF}.^2)/2;
% %         end
% %     end
% %  
% % end


%% prepare for storing experimental results
% multiple-scaled Gobar template
syma = zeros(sx, sy); 
syma_top = zeros(sx, sy); 
if numDoG>0
    numDoGTemplate=1;
else
    numDoGTemplate=0;
end
  
sym = cell(1, numScale + numDoGTemplate); 
for (i = 1 : numScale + numDoGTemplate)
   sym{i} = zeros(sx, sy); 
end
    
% recontructed images
Asym = cell(1, numImg); 
for (i = 1 : numImg)
    Asym{i} = zeros(sx, sy); 
end

% deformed Gabor templates
Asyma = cell(1, numImg); 
for (i = 1 : numImg)
    Asyma{i} = zeros(sx, sy); 
end

% deformed Gabor templates visualized by only a few of top features
Asyma_top = cell(1, numImg); 
for (i = 1 : numImg)
    Asyma_top{i} = zeros(sx, sy); 
end

% deformed DoG templates
Asym0 = cell(1, numImg*(numScale+1));
for (i = 1 : numImg*(numScale+1))
    Asym0{i} = zeros(sx, sy);  
end

filterSelected = zeros(1, numSketch);
xSelected = zeros(1, numSketch);
ySelected = zeros(1, numSketch);
numSketchSelected = 0;


%% do matching pursuit
tic;
MatchingPursuit(numImg, numOrient*numScale+numDoG, numOrient, numDoG, allFilteredR,...
                allFilteredI, allFiltered, Crr, Cri, Cir, Cii, half, sx, sy,...
                allFilterR(1, :), allFilterI(1, :), allSymbol(1, :), syma, sym,...
                Asym, Asym0, Asyma, locationShiftLimit, orientShiftLimit, numSketch,...
                epsilon, Upperbound, SHUTUP, filterSelected, xSelected, ySelected, numSketchSelected);
 
disp(['finished Matching Pursuit: ' num2str(toc) ' seconds']);
  
 
 %% deformed template for every training image          
 filterSelected = filterSelected(1:numSketchSelected); 
 xSelected = xSelected(1:numSketchSelected);
 ySelected = ySelected(1:numSketchSelected);
 
 
 %% store the matching pursuit version of template
 templateMP=struct('selectedFilter',single(filterSelected),'selectedRow',single(xSelected),'selectedCol',single(ySelected), 'numSelected',single(numSketchSelected),'sym',-syma);
 sketchTemplate=-syma;
 %%
   

 %[indEliminated, numSelectedFeatures]= createMask(filterSelected,xSelected,ySelected, sx, sy, numOrient, numScale, numDoG);
 [filterSelected, xSelected,ySelected, numSelectedFeatures]= indexTtransformation(filterSelected, xSelected, ySelected, sx, sy, numOrient, numScale, numDoG);

 %% store the Frame version of template
 template=struct('selectedFilter',single(filterSelected),'selectedRow',single(xSelected),'selectedCol',single(ySelected), 'numSelected',single(numSelectedFeatures));

  
 
%% show and save experimental results
% tic;
% ShowAndSaveResults;
% disp(['saving the filter selection results: ' num2str(toc) ' seconds']);
% 
% %% generate html to show filter selection experimental results
% GenerateHtmlForFilterSelection;




