function FDI(caso,dim,k,simulacoes,salvar)

%
% Função Distribuição da Informação
%
% Esta rotina estima a função de distribuição da informação nos pesos sinápticos 
% dos modelos de memória associativa auto e hetero-associativas.
%
% Comando: FDI(caso,dim,k,simulacoes,salvar);
%
% Entradas: 
%       Caso - Caso auto-associativo (caso=1) ou caso heteroassociativo (caso=2);
%              Padrão: Caso = 1;
%        dim - Dimensão dos padrões de entrada e saída ;
%              Padrão: dim = [25,20];
%          k - Número de pares de padrões armazenados nas memórias associativas;
%              Padrão: k = 3;
% simulacoes - Número de simulações (Padrão: simulacoes = 100)
%     salvar - Salva os resultados no arquivo fdi_auto.mat 
%                ou fdi_hetero.mat quando salvar=1 (padrão: salvar=1);
%

if nargin<5
   salvar=1;
end

if nargin<4
   simulacoes=100;
end

if nargin<3
   k=3;
end

if nargin<2
   dim = [25,20];
end

if nargin<1
   caso=1;
end

n=dim(1);
if length(dim)<2
   m=dim(1);
else
   m=dim(2);
end

for s=1:simulacoes
   disp('');
    tic
    % Gerar as matrizes contendo as memórias fundamentais;
    clear X
    X=sign(2*rand(n,k)-1); % Padrões bipolares
    clear Y
    Y=sign(2*rand(m,k)-1); % Padrões bipolares
    
    % Memória Associativa de Hopfield;
    if caso==1
       % Matriz dos pesos sinápticos:
       W=corr0(X);
       
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(W);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(W);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Todas as conexões sinápticas, recordação com uma única iteração!
       D_hopfield(s,1)=(norm(X-sign(W*X))<1.e-4);
       % Excluindo as conexões sinápticas:
       
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             W(ind_i(ex(ind)),ind_j(ex(ind)))=0;
          end
          D_hopfield(s,i+1)=(norm(X-sign(W*X))<1.e-4);
       end
    end
    
    % BAM;
    if caso==2
       % Matriz dos pesos sinápticos:
       W=corr(X,Y);
              
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(W);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(W);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Todas as conexões sinápticas, recordação com uma única iteração!
       D_bam(s,1)=(norm(Y-sign(W*X))<1.e-4);
       % Excluindo as conexões sinápticas:
       
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             W(ind_i(ex(ind)),ind_j(ex(ind)))=0;
          end
          D_bam(s,i+1)=(norm(Y-sign(W*X))<1.e-4);
       end
    end
    
    % Memória Associativa de Personnaz
    if caso==1
       % Matriz dos pesos sinápticos:
       W=proj(X,X);
       
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(W);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(W);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Todas as conexões sinápticas,
       D_personnaz(s,1)=(norm(X-sign(W*X))<1.e-4);
       % Excluindo as conexões sinápticas:
       
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             W(ind_i(ex(ind)),ind_j(ex(ind)))=0;
          end
          D_personnaz(s,i+1)=(norm(X-sign(W*X))<1.e-4);
       end
    end
    
    % Memória Associativa de Kanter-Sompolinsky
    if caso==1
       % Matriz dos pesos sinápticos:
       W=proj0(X);
       
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(W);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(W);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Todas as conexões sinápticas, uma única iteração!
       D_kanter(s,1)=(norm(X-sign(W*X))<1.e-4);
       % Excluindo as conexões sinápticas:
       
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             W(ind_i(ex(ind)),ind_j(ex(ind)))=0;
          end
          D_kanter(s,i+1)=(norm(X-sign(W*X))<1.e-4);
       end
    end
    
    % ABAM
    if caso==2
       % Matriz dos pesos sinápticos:
       W=proj(X,Y);
              
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(W);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(W);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Todas as conexões sinápticas, uma única iteração!
       D_abam(s,1)=(norm(Y-sign(W*X))<1.e-4);
       % Excluindo as conexões sinápticas:
       
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             W(ind_i(ex(ind)),ind_j(ex(ind)))=0;
          end
          D_abam(s,i+1)=(norm(Y-sign(W*X))<1.e-4);
       end
    end
    
    % ECAM    
    if caso==1
       % As matrizes dos pesos sinápticos são:
       %		a matriz contendo as chaves-fundamentais e sua transposta!
       % As duas matrizes contendo as conexões sinápticas são convertidas 
       %    em um vetor-linha: vec.
       clear vec
    	 vec=[reshape(X,1,n*k),reshape(X,1,n*k)];
       
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(vec);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(vec);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Transformando o vetor em matrizes das conexões sinápticas,
       clear W1; clear W2;
       W1=reshape(vec(1:k*n),n,k); W2=reshape(vec(k*n+1:(2*n)*k),n,k)';
       % Todas as conexões sinápticas, uma única iteração!
       D_ecam(s,1)=(norm(X-sign(W1*exp(W2*X)))<1.e-4);
       % Excluindo as conexões sinápticas:
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             vec(ind_i(ex(ind)),ind_j(ex(ind)))=0;
          end
          % Transformando o vetor em matrizes das conexões sinápticas,
          clear W1; clear W2;
          W1=reshape(vec(1:k*n),n,k); W2=reshape(vec(k*n+1:(2*n)*k),n,k)';
          % Calculando a distribuição da informação;
          D_ecam(s,i+1)=(norm(X-sign(W1*exp(W2*X)))<1.e-4);
       end
    end
    
    % BECAM    
    if caso==2
       % As matrizes dos pesos sinápticos são:
       %		a matriz contendo as chaves-fundamentais 
       %    e a transposta da matriz das recordações-fundamentais!
       % As duas matrizes contendo as conexões sinápticas são convertidas 
       %    em um vetor-linha: vec.
       clear vec
    	 vec=[reshape(Y,1,m*k),reshape(X,1,n*k)];
       
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(vec);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(vec);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Transformando o vetor em matrizes das conexões sinápticas,
       W1=reshape(vec(1:k*m),m,k); W2=reshape(vec(k*m+1:(m+n)*k),n,k)';
       % Todas as conexões sinápticas, uma única iteração!
       D_becam(s,1)=(norm(Y-sign(W1*exp(W2*X)))<1.e-4);
       % Excluindo as conexões sinápticas:
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             vec(ind_i(ex(ind)),ind_j(ex(ind)))=0;
          end
          % Transformando o vetor em matrizes das conexões sinápticas,
          W1=reshape(vec(1:k*m),m,k); W2=reshape(vec(k*m+1:(m+n)*k),n,k)';
          % Calculando a distribuição da informação;
          D_becam(s,i+1)=(norm(Y-sign(W1*exp(W2*X)))<1.e-4);
       end
    end
    
    X=(X+1)/2; % Padrões binários (modelos morfológicos)
    Y=(Y+1)/2; % Padrões binários (modelos morfológicos)
    
    % Memória Auto-associativa Morfológica Wxx
    if caso==1
       % Matriz dos pesos sinápticos:
       clear W
       W=boxmin(X,-X');
       
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(W>-10);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(W>-10);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Todas as conexões sinápticas, uma única iteração!
       D_Wxx(s,1)=(norm(X-boxmax(W,X))<1.e-4);
       % Excluindo as conexões sinápticas:
       
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             W(ind_i(ex(ind)),ind_j(ex(ind)))=-100;
          end
          D_Wxx(s,i+1)=(norm(X-boxmax(W,X))<1.e-4);
       end
    end
    
    % Memória Associativa Morfológica Wxy
    if caso==2
       % Matriz dos pesos sinápticos:
       clear W
       W=boxmin(Y,-X');
       
       % Número de conexões sinápticas não-nulas;
    	 clear nz;
    	 nz=nnz(W>-10);
       clear ind_i; clear ind_j
       % Localizando os índices das conexões sinápticas não nulas; 
    	 [ind_i,ind_j]=find(W>-10);
       clear ex;
       % Gerando a sequência das conexões sinápticas que serão excluídas, 
    	 ex=randperm(nz);
    	 clear passo;
       passo=fix(nz/100);
        
       % Calculando a distribuição da informação:
       % Todas as conexões sinápticas, uma única iteração!
       D_Wxy(s,1)=(norm(Y-boxmax(W,X))<1.e-4);
       % Excluindo as conexões sinápticas:
       
       ind=0;
       for i=1:100
          for j=1:passo
             ind=ind+1;
             W(ind_i(ex(ind)),ind_j(ex(ind)))=-100;
          end
          D_Wxy(s,i+1)=(norm(Y-boxmax(W,X))<1.e-4);
       end
    end
    
    % Memória Associativa Morfológica de Duas Camadas
    % - Caso Auto-associativo:
    if caso==1
       % Matrizes dos Pesos;
    	 Z=eye(k);
    	 Wzy=boxmin(X,-Z');
       Mxxz=boxmin(boxmax(Z,-X'),boxmax(X,-X'));
       % Vetor contendo as conexões sinápticas;
    	 vec=[reshape(Wzy,1,n*k),reshape(Mxxz,1,n*k)];
    
    	 % Número de conexões sinápticas;
    	 clear nnz;
    	 nz=nnz(-10<vec<10);
    	 clear ind_i; clear ind_j
    	 [ind_i,ind_j]=find(-10<vec<10);
    	 clear ex;
    	 ex=randperm(nz);
    	 clear passo;
    	 passo=fix(nz/100);
    
    	 % Calculando a Distribuição da Informação;
    	 ind=0;
    	 Wzy=reshape(vec(1:n*k),n,k);
    	 Mxxz=reshape(vec(n*k+1:(2*n)*k),k,n);
    	 D_mam2c1(s,1)=(norm(X-boxmax(Wzy,ceil(boxmin(Mxxz,X)/100)))<1.e-4);
       for i=1:100
       	for j=1:passo
         	ind=ind+1;
            vec(ind_i(ex(ind)),ind_j(ex(ind)))=-1*sign(k*n-ex(ind))*100;
        	end
     	   Wzy=reshape(vec(1:n*k),n,k);
         Mxxz=reshape(vec(n*k+1:(2*n)*k),k,n);
         D_mam2c1(s,i+1)=(norm(X-boxmax(Wzy,ceil(boxmin(Mxxz,X)/100)))<1.e-4);
       end
    end
    
    % Memória Associativa Morfológica de Duas Camadas    
    % - Caso Hetero-associativo:
    if caso==2
       % Matrizes dos Pesos;
    	 Z=eye(k);
    	 Wzy=boxmin(Y,-Z');
       Mxxz=boxmin(boxmax(Z,-X'),boxmax(X,-X'));
       % Vetor contendo as conexões sinápticas;
    	 vec=[reshape(Wzy,1,m*k),reshape(Mxxz,1,n*k)];
    
    	 % Número de conexões sinápticas;
    	 clear nnz;
    	 nz=nnz(-10<vec<10);
    	 clear ind_i; clear ind_j
    	 [ind_i,ind_j]=find(-10<vec<10);
    	 clear ex;
    	 ex=randperm(nz);
    	 clear passo;
    	 passo=fix(nz/100);
    
    	 % Calculando a Distribuição da Informação;
    	 ind=0;
    	 Wzy=reshape(vec(1:m*k),m,k);
    	 Mxxz=reshape(vec(m*k+1:(m+n)*k),k,n);
    	 D_mam2c2(s,1)=(norm(Y-boxmax(Wzy,ceil(boxmin(Mxxz,X)/100)))<1.e-4);
       for i=1:100
       	for j=1:passo
         	ind=ind+1;
            vec(ind_i(ex(ind)),ind_j(ex(ind)))=-1*sign(k*m-ex(ind))*100;
        	end
     	 	Wzy=reshape(vec(1:m*k),m,k);
       	Mxxz=reshape(vec(m*k+1:(m+n)*k),k,n);
         D_mam2c2(s,i+1)=(norm(Y-boxmax(Wzy,ceil(boxmin(Mxxz,X)/100)))<1.e-4);
       end
    end
    
    disp(sprintf('Faltam %d simulações. O tempo estimado para completar a execução é %d min.',simulacoes-s,ceil((simulacoes-s)*toc/60)));
end
disp('')

figure(1)
clf
% Construção dos gráficos!
if caso==1
   clear Dm_hopfield; Dm_hopfield=mean(D_hopfield);
   clear Dm_personnaz; Dm_personnaz=mean(D_personnaz);
   clear Dm_kanter; Dm_kanter=mean(D_kanter);
   clear Dm_ecam; Dm_ecam=mean(D_ecam);
   clear Dm_Wxx; Dm_Wxx=mean(D_Wxx);
   clear Dm_mam2c1; Dm_mam2c1=mean(D_mam2c1);
   x=[0:100];
   plot(x,Dm_hopfield,x,Dm_personnaz,x,Dm_kanter,x,Dm_ecam,x,Dm_Wxx,x,Dm_mam2c1);
   legend('MA Hopfield','MA Personnaz','MA Kanter-Somp.','ECAM','MAM Wxx','MAM 2 Camadas');
   xlabel('% Conexões Sinapticas Excluídas','FontSize',12);
   ylabel('FDI','FontSize',12);
   
   if salvar==1
      save fdi_auto Dm_hopfield Dm_personnaz Dm_kanter Dm_ecam Dm_Wxx Dm_mam2c1 dim k simulacoes;
      disp('Confira os resultados no arquivo fdi_auto.mat')
   end   
end   

if caso==2
   clear Dm_bam; Dm_bam=mean(D_bam);
   clear Dm_abam; Dm_abam=mean(D_abam);
   clear Dm_becam; Dm_becam=mean(D_becam);
   clear Dm_Wxy; Dm_Wxy=mean(D_Wxy);
   clear Dm_mam2c2; Dm_mam2c2=mean(D_mam2c2);
   x=[0:100];
   plot(x,Dm_bam,x,Dm_abam,x,Dm_becam,x,Dm_Wxy,x,Dm_mam2c2);
   legend('BAM','ABAM','BECAM','MAM Wxy','MAM 2 Camadas');
   xlabel('% Conexões Sinapticas Excluídas','FontSize',12);
   ylabel('FDI','FontSize',12);
   
   if salvar==1
      save fdi_hetero Dm_bam Dm_abam Dm_becam Dm_Wxy Dm_mam2c2 dim k simulacoes;
      disp('Confira os resultados no arquivo fdi_hetero.mat')
   end
end