function GetNaturalStat
% GetNaturalStat - Obtain statistics from natural images.
%

% CONSTANTS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
scale = .7;  % scale of Gabors, length = 17 pixels
norient = 16;  % number of orientations
Upperbound = 6.; % saturation value
binsize = .2;  % binsize for computing F

% TRAINING EXAMPLES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
rootFolder = '..';
negFolder = sprintf('%s/image/neg',rootFolder);  % image input folder
negTrainImages = dir(sprintf('%s/neg*.jpg', negFolder));
nimage = length(negTrainImages); % #images
nimage = min(2,nimage); % enough to use the first 2 images only
negTrainImages = negTrainImages(1:nimage);
Isize = zeros(nimage, 2);
I = cell(1,nimage);
for i = 1 : nimage
    tmpIm = imread([negFolder '/' negTrainImages(i).name]);
    if size(tmpIm,3) == 3
        tmpIm = rgb2gray(tmpIm);
    end
    I{i} = double(tmpIm);
    Isize(i, :) = size(I{i});
end

% FILTERING IMAGES %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
disp('start filtering');
tic
[allfilter allsymbol] = makefilter(scale, norient);  % generate Gabor filters 
h = (size(allfilter{1}, 1)-1)/2;  % half size of Gabor
allfiltered = cell(nimage,norient);
for i = 1:nimage
    allfiltered(i,:) = applyfilterfftsameLocalNorm(I(i), allfilter);  % filter training images
end
disp(['filtering time: ' num2str(toc) ' seconds']);

binnum = floor(Upperbound/binsize)+1;  % binnumbers
histog = zeros(binnum, 1);  % store F  

% MEX-C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
disp('start histogramming');
Chistog(nimage, norient, allfiltered, h, Isize(:,1), Isize(:,2), binsize, binnum, histog, Upperbound);
disp(['histogramming time: ' num2str(toc) ' seconds']);

% RESULTS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
r = (0:(binnum-1))*binsize;

M = 50; 
e = zeros(M, 1); 
Z = zeros(M, 1); 
for k=1:M
    lambda = (k-1.)/10.; 
    p = exp(lambda*r).*(histog'); 
    Z(k) = sum(p*binsize); p = p/Z(k); 
    %plot(r, p); 
    e(k) = sum(r.*p*binsize);
end
lam = (0:(M-1))/10.; 
lz = log(Z); 

save(sprintf('%s/result/natural_stat.mat',rootFolder), 'binnum','M',...
    'lam', 'e', 'lz', 'Upperbound',...
    'scale', 'norient', 'allfilter', 'allsymbol');




