Generate the diagnostic and antidiagnostic faces, the same process used in the experimental part that located the bubbles at '1' positions, at each SF band.

bubbles.m

load allResults.mat


ySize = 256;
xSize = 256;

imSubSample=0;
nBands = 6; % number of bandpass (cannot be changed at the moment)
fGauss=mkGaussian(11);
fGauss1=fGauss(:,1);
fGauss1=sqrt(2)*fGauss1/sum(fGauss1);
clear fGauss;

%% Gaussian windows
nPeriod = 3;
nZero = 2.18;
stdev = nPeriod * 2^(5); %stdev is equal to constant * half a period of the octave-spaced spatial frequencies
maxHalfSize5 = round(stdev * nZero); %cuts the gaussian after 2.18 * stdev
gauss5 = zeros(2*maxHalfSize5,2*maxHalfSize5);
[y,x] = meshgrid(-maxHalfSize5:maxHalfSize5,-maxHalfSize5:maxHalfSize5);
gauss5 = exp(-(x.^2/stdev^2)-(y.^2/stdev^2));
clear x, y;

% method = 'bicubic';
method = 'nearest';

stdev = nPeriod * 2^(4); %stdev is equal to half a period of the octave-spaced spatial frequencies
maxHalfSize4 = round(stdev * nZero); %cuts the gaussian after 2.18 * stdev
gauss4 = double(imresize(gauss5,[(2*maxHalfSize4+1) (2*maxHalfSize4+1)], method));

stdev = nPeriod * 2^(3); %stdev is equal to half a period of the octave-spaced spatial frequencies
maxHalfSize3 = round(stdev * nZero); %cuts the gaussian after 2.18 * stdev
gauss3 = double(imresize(gauss5,[(2*maxHalfSize3+1) (2*maxHalfSize3+1)], method));

stdev = nPeriod * 2^(2); %stdev is equal to half a period of the octave-spaced spatial frequencies
maxHalfSize2 = round(stdev * nZero);  %cuts the gaussian after 2.18 * stdev
gauss2 = double(imresize(gauss5,[(2*maxHalfSize2+1) (2*maxHalfSize2+1)], method));

stdev = nPeriod * 2^(1); %stdev is equal to half a period of the octave-spaced spatial frequencies
maxHalfSize1 = round(stdev * nZero);  %cuts the gaussian after 2.18 * stdev
gauss1 = double(imresize(gauss5,[(2*maxHalfSize1+1) (2*maxHalfSize1+1)], method));

%%%% make_bubs_std;      %%or you can calls this function to define the basics bubble sizes

cacaPlane = double(zeros(ySize,xSize));
winPlane = double(zeros(ySize,xSize));
winImage = double(zeros(xSize,ySize));
tempPlane = double(zeros(ySize+2*maxHalfSize5,xSize+2*maxHalfSize5));
[pyr,pind] = buildLpyr(double(W1),6,fGauss1); % W1 is one of the faces used in the expriment
for i = 1:5
    tempPlane = zeros(ySize+2*maxHalfSize5,xSize+2*maxHalfSize5);
    nameGauss = sprintf('gauss%d',i); % size of Guassian approriate to  SF band. 
	nameMax = sprintf('maxHalfSize%d',i);
    nameH = allDiagnostic(:,:,i); %% compute the diagnostic mask  
    %nameH = allAntiDiagnostic(:,:,i); %% uncomment this line to compute the anti_diagnostic mask

     for y = 1:256
         for x = 1:256
             if (nameH(y,x) == 1)
                 xMax = x;
                 yMax = y;
                 tempPlane(yMax:yMax+2*eval(nameMax),xMax:xMax+2*eval(nameMax)) = tempPlane(yMax:yMax+2*eval(nameMax),xMax:xMax+2*eval(nameMax)) + eval(nameGauss);
             end
         end
     end
     cacaPlane = double(tempPlane(eval(nameMax):ySize+eval(nameMax)-1,eval(nameMax):xSize+eval(nameMax)-1));% cuts the crap
     winPlane = double(sign(round(cacaPlane / 2)));
     winPlane = double(cacaPlane .* (1- winPlane) + winPlane);
        
     winImage = double(double(winImage) + double(reconLpyr(pyr,pind,[i],fGauss1) .* winPlane));

end

winImage = double(double(winImage) + double(reconLpyr(pyr,pind,[6],fGauss1)));
winImage = uint8(winImage);
imshow(winImage);