% [y1k y2k] = block_separator(mix1k, mix2k, wave, 0, g_max, str_path_WC )
%
% Written by: Mariana S.C. Almeida , Luís B. Almeida
% e-mail: mariana.almeida@lx.it.pt , luis.almeida@lx.it.pt
%
% Paper: [1] M. S. C. Almeida and L. B. Almeida, "Wavelet-based nonlinear separation
% of two page images", Neurocomputing, 2007.
%
% BLOCK_SEPARATOR : Separates front and back image pages (mix1 and mix2).
% Separation is separatly applyed to block images (mixk2,mixk1).
% Block images must overlap 256 pixels. 
% Should be used when of the two image dimention is larger than 1024
% (calls: "separator" function)
%
% wave  - String indicating which wavelet trnasform is used:
%         'inv' - Discrete stationary wavelet transform is used
%         'complex' -  complex wavelet transform is used
%                      (http://taco.poly.edu/WaveletSoftware/)
%         Default: 'inv'
% init  - Indicates if initial linear pre-processing is prefromed:
%         1 - with pre-processing
%         1 - without pre-processing
%         Default: 1
% g_max - value of compasatory gain (see [1] for details).
%         Default: 3.
% str_path_WC - String indicating the path where is the code for complex wavelet transform
%               (availavle at: http://taco.poly.edu/WaveletSoftware/)
%               Default: runing path.


 [y1k y2k] = block_separator(mix1k, mix2k, wave, init, g_max, str_path_WC );

dim_bloco=1024;
overlap = 256;
[idim jdim] = size(mix1);

[nblocos_i overlaps_i] = no_blocos(idim,overlap, dim_bloco);
[nblocos_j overlaps_j] = no_blocos(jdim,overlap, dim_bloco);

% % indices !!
extra_i = floor(min(overlaps_i)/2)-1;
for k=1:nblocos_i
    aux = 1 + (k-1)*dim_bloco - sum(overlaps_i(1:k-1)) ;
    indis_i(k,:) = (aux:aux+dim_bloco-1);
    
    if indis_i(k,1)==1      aux1=1;    else aux1=indis_i(k,1)+extra_i ; end
    if indis_i(k,end)==idim aux2=idim; else aux2=indis_i(k,end)-extra_i ; end
    indis_i2{k} = (aux1:aux2);
end 
extra_j = floor(min(overlaps_j)/2)-1;
for k=1:nblocos_j
    aux = 1 + (k-1)*dim_bloco - sum(overlaps_j(1:k-1)) ;
    indis_j(k,:) = (aux:aux+dim_bloco-1);
    
    if indis_j(k,1)==1      aux1=1;    else aux1=indis_j(k,1)+extra_j ; end
    if indis_j(k,end)==jdim aux2=jdim; else aux2=indis_j(k,end)-extra_j ; end
    indis_j2{k} = (aux1:aux2);
end 
% % (end indices)

if init
    X_M(1,:,:)= double(mix1);  X_M(2,:,:)= double(mix2);
    [X_out Mf] = esf_sim(reshape(X_M,2,[]));
    [auxi,idim,jdim]=size(X_M);
    X_M = reshape(X_out,2,idim,jdim);
    clear X_out
else
    [auxi,idim,jdim]=size(X_M);
end

y1=[]; y2=[]; y1k=[]; y2k=[];

for i = 1 : nblocos_i
    for j = 1 : nblocos_j

        mix1k = mix1( indis_i(i,:) , indis_j(j,:) );
        mix2k = mix2( indis_i(i,:) , indis_j(j,:) );
        %         figure(5), imagesc(mix1k-mix2k), colormap gray , axis image, title('alinhamento do bloco')
        save auxi_ext y1 y2 mix1 mix2,
        clear y1 y2 mix1 mix2 y1k y2k

        fprintf('Aligment of block : ( %d , %d ) ',i,j);   
        [y1k y2k] = separator(mix1k, mix2k, wave, 0, g_max, str_path_WC );
        % %         denois_inv_g_k;
        load auxi_ext
        
        y1k_=y1k; y2k_=y2k;
        if i==1 & i ~= nblocos_i
            y1k = y1k(1:end-extra_i,:); y2k = y2k(1:end-extra_i,:);
        elseif i~=1 & i == nblocos_i
            y1k = y1k(extra_i+1:end,:); y2k = y2k(extra_i+1:end,:);
        end
        if j==1 & j ~= nblocos_j
            y1k = y1k(:,1:end-extra_j); y2k = y2k(:,1:end-extra_j);
        elseif j~=1 & j == nblocos_j
            y1k = y1k(:,extra_j+1:end); y2k = y2k(:,extra_j+1:end);
        end       
        
        y1( indis_i2{i} , indis_j2{j} ) = y1k;
        y2( indis_i2{i} , indis_j2{j} ) = y2k;

%         figure(10), imagesc(y1k), colormap gray , axis image, title('Den 1')
%         figure(11), imagesc(fliplr(y2k)), colormap gray , axis image, title('Den 2')

    end
end

% figure, imagesc(mix1_den), colormap gray , axis image, title('Den 1 - Final ')
% figure, imagesc(fliplr(mix2_den)), colormap gray , axis image, title('Den 2 - Final')
