function out=miaOrdinal(N,dispresults)
% Function which computes the gamma index
%
%  Required input arguments:
%
%       N    :    Contingency table.
%                 Matrix or Table.
%                 Matrix or table which contains the input contingency
%                 table (say of size I-by-J).
%
%  Optional input arguments:
%
%  dispresults : Boolean. If dispresults is true the routine shows the
%               number of concordances, the number of discordances and the
%               value of the gamma index. The default value of dispresults
%               is false
%
%   Output   
%
% out : scalar. The value of the gamma index.
%
% Examples:

%{
    %%  miaOrdinal with just one input argument.
    % Rows of N indicate the results of a written test with levels:
    % 'Sufficient' 'Good' Very good'
    % Columns of N indicate the results of an oral test with levels:
    % 'Sufficient' 'Good' Very good'
          N=[20    40    20;
            10    45    45;
             0     5    15];
    out=miaOrdinal(N);
%}

%{
    %%  miaOrdinal with just two input arguments.
    % Rows of N indicate the results of a written test with levels:
    % 'Sufficient' 'Good' Very good'
    % Columns of N indicate the results of an oral test with levels:
    % 'Sufficient' 'Good' Very good'
          N=[20    40    20;
            10    45    45;
             0     5    15];
    out=miaOrdinal(N,true);
%}


%% Beginning of code
% Se il numero di argomenti di input è minore di 2 allora l'output non
% viene mostrato
if nargin<2
    dispresults=false;
end

% Se N è una table allora viene trasformata in array di doubles
if istable(N)
    N=table2array(N);
end



% Estrai il numero di righe e di colonne delle tabella di input N
[I,J] = size(N);

% Trovo le coppie discordanti e concordanti 
    C=0;
    D=0;
    for i=1:I
        for j=1:J
            if i<I && j<J
                xsel=N(i+1:I,j+1:J);
                C=C+N(i,j)*sum(xsel(:));
            end
            if j>1 && i<I
                xsel=N(i+1:I,1:j-1);
                D=D+N(i,j)*sum(xsel(:));
            end
        end
    end

% Gamma index
out = (C - D)/(C + D);

if dispresults ==true
    disp(['Il numero totale di concordanze è: ' num2str(C)])
    disp(['Il numero totale di discordanze è: ' num2str(D)])
    disp(['Il valore dell''indice gamma è: ' num2str(out)])
end

%% Il codice per calcolare gli altri indici e metterli tutti dentro una struct è riportato di seguito

% n = sum(sum(N)); %sample size
% % totpairs= total number of pairs
% totpairs=0.5*n*(n-1);
% % nidot = row sums
% nidot=sum(N,2);
% % ndotj = columns sums
% ndotj=sum(N,1);
% 
% % Tau-a index
% taua = (C - D)/totpairs;
% % Taub-b index
% % Tx = number of pairs tied on X
% Tx=sum(nidot.*(nidot-1)/2);
% % Ty = number of pairs tied on Y
% Ty=sum(ndotj.*(ndotj-1)/2);
% taub=(C - D)/sqrt( (totpairs-Tx)*(totpairs-Ty) );
% % Tau-c index
% m=min(I,J);
% tauc=m*2*(C-D)/( n^2 *(m-1));
% % Somers index
% wr=n^2-sum(nidot.^2);
% som = 2*(C-D)/wr;
% 
% out=struct;
% out.gam=gam;
% out.taua=taua;
% out.taub=taub;
% out.tauc=tauc;
% out.som=som;

end